{'fraction_non_alphanumeric': 0.04267361423052735, 'fraction_numerical': 0.009590333303386937, 'mean_word_length': 4.294292508917954, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 2, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 16, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Robust Beamforming with Pilot Reuse Scheduling in a Heterogeneous Cloud Radio Access Network

This paper considers a downlink ultra-dense heterogeneous cloud radio access network (H-CRAN) which guarantees seamless coverage and can provide high date rates. In order to reduce channel state information (CSI) feedback overhead, incomplete inter-cluster CSI is considered, i.e., each remote radio head (RRH) or macro base station (MBS) only measures the CSI from user equipments (UEs) in its serving cluster. To reduce pilot consumption, pilot reuse among UEs is assumed, resulting in imperfect intra-cluster CSI. A two-stage optimization problem is then formulated. In the first stage, a pilot scheduling algorithm is proposed to minimize the sum mean square error (MSE) of all channel estimates. Specifically, the minimum number of required pilots along with a feasible pilot allocation solution are first determined by applying the Dsatur algorithm, and adjustments based on the defined level of pilot contamination are then carried out for further improvement. Based on the pilot allocation result obtained in the first stage, the second stage aims to maximize the sum spectral efficiency (SE) of the network by optimizing the beam-vectors. Due to incomplete inter-cluster CSI and imperfect intra-cluster CSI, an explicit expression of each UE's achievable rate is unavailable. Hence, a lower bound on the achievable rate is derived based on Jensen's inequality, and an alternative robust transmission design (RTD) algorithm along with its distributed realization are then proposed to maximize the derived tight lower bound. Simulation results show that compared with existing algorithms, the system performance can be greatly improved by the proposed algorithms in terms of both sum MSE and sum SE.

a 1000 times capacity boost is thus expected in the fifthgeneration (5G) network compared to the current fourthgeneration network [2]. To realize this 5G vision, cloud radio access network (C-RAN) has been recognized as a promising solution [3]. In the C-RAN architecture, a baseband unit (BBU) pool with powerful computation capability acts as a cloud data center, and remote radio heads (RRHs) configured only with some radio-frequency functionalities are connected to the BBU pool through optical fiber fronthaul links. Due to the simplified functionalities, RRHs can be densely and distributedly deployed to improve network access in conventional cellular networks, especially in hot spots with a large number of user equipments (UEs) like hospitals, shopping malls, etc. In addition, different from untra-dense small cell networks that suffer from cochannel interference [4], interference mitigation can be effectively realized in a C-RAN by applying Coordinated Multi-Point Transmission/Reception (CoMP) thanks to the powerful BBU pool [5].
Since a C-RAN is mainly adopted to provide high data rates in hot spots, real-time voice service and control signalling are not efficiently supported. If UEs move fast, the switching speed of the RRH service in a C-RAN is relatively high, resulting in a high signalling exchange load. On the other hand, the fronthaul capacity of a C-RAN is usually limited, making it difficult to serve all UEs in a network. Hence, lack of high power nodes may make it hard to ensure backward compatibility with the existing cellular systems [6], [7]. In order to deliver the overall control signalling and guarantee seamless coverage, an advanced architecture, known as a heterogeneous cloud radio access network (H-CRAN), was proposed to combine the advantages of both C-RAN and heterogeneous networks [8]- [10]. In an H-CRAN, besides the BBU pool and RRHs, macro base stations (MBSs) are also included. The BBU pool and MBSs are interfaced via backhaul links for coodination. Hence, the delivery of control and broadcast signalling can be shifted from RRHs to MBSs to alleviate the capacity and time delay on the fronthaul [11], [12]. Unnecessary handover and re-association can then be avoided.
In H-CRANs, cochannel interference suppression is an important technical issue. In C-RANs, only interference among RRHs exists and it can be effectively mitigated via centralized processing at the BBU pool. While in H-CRANs, when RRHs and MBSs operate on the same time-frequency resource block (RB), a UE served by RRHs (a MBS) will not only experience intra-tier interference from RRHs (MBSs), but also suffer additional inter-tier interference from MBSs (RRHs). Hence, the signal-to-interference-plus-noise ratio (SINR) expressions of UEs contain more cochannel interference terms, making the problem more challenging. In [11]- [13], interference suppression problems in H-CRANs were studied aiming at different design metrics. Specifically, reference [11] aimed to maximize the average throughput and maintain the network stability by traffic admission control, user association and resource allocation. In [12], a contract-based interference coordination framework was proposed to mitigate the inter-tier interference between RRHs and MBSs. Reference [13] aimed to maximize the energy efficiency (EE) of an H-CRAN by resource assignment and power control. All these works, however, focused on single-input single-output (SISO) networks. By considering multi-antenna MBSs and multi-antenna RRHs, the performance of H-CRANs can be further enhanced by applying beamforming techniques.
Downlink beamforming design has been widely studied in both the conventional cellular networks [14]- [16] and C-RANs [17]- [19]. Perfect global channel state information (CSI) was assumed available for system performance analysis and optimization. However, it is in practice difficult to obtain perfect CSI of all links due to channel estimation and quantization errors. In addition, feeding back the CSI from all UEs to each RRH or MBS requires excessive overhead, which may easily overwhelm the capacity of the wireless radio interface, especially in ultra-dense C-RANs [20], [21]. Hence, transmission design based on incomplete CSI has drawn great attention recently [22]- [28]. In [22] and [23], a multi-user massive multiple-input multiple-output (MIMO) network and a full-duplex MIMO cognitive radio system are respectively considered, and beamforming vectors were designed under channel uncertainties. In [24]- [26], distributed transmit beamforming was studied with imperfect CSI under different scenarios. Specifically, references [24] and [25] considered a cognitive radio network, while [26] considered a multicell cellular system. In [27] and [28], incomplete CSI was assumed in C-RANs, i.e., each RRH only estimated the CSI from UEs in its serving cluster (named intra-cluster CSI). As for the UEs outside the serving cluster, it was assumed that the RRH only had the large-scale channel gains (named inter-cluster CSI). However, perfect intra-cluster CSI was still assumed in [27] and [28]. To this end, orthogonal training has to be adopted for channel estimation. In this case, the length of pilot overhead increases linearly with the number of UEs, which could be unaffordable for ultra-dense networks.
One promising way to reduce pilot overhead is allowing pilots to be reused among UEs. Pilot reuse design has been extensively studied in massive MIMO networks [29]- [33] as well as device-to-device (D2D) underlaid systems [34]- [36]. In particular, references [29]- [31] considered pilot reuse in multiple-cell scenarios, i.e., UEs in the same cell use orthogonal pilots, and the same set of pilots are reused in different cells. In [32] and [33], it was shown that due to the uncorrelation feature of massive MIMO antennas, the same pilot could be reused by UEs with different angular positions. In [34]- [36], pilots were allowed to be reused by D2D pairs. Since RRHs are usually located dispersively and they use low power for short-distance transmission, the pilot reuse by UEs far away from each other would cause marginal pilot contamination. This paper mainly focuses on pilot scheduling and robust beamforming design for an ultra-dense H-CRAN. To the best of the authors' knowledge, this problem has not yet been studied. The main contributions of this paper are summarized as follows: • This paper considers an ultra-dense H-CRAN, where  RRHs are mainly used to provide high data rates  • In the simulation part, the performances of the proposed algorithms are illustrated and compared in terms of both sum MSE and sum SE. Simulation results show that the sum MSE of channel estimation can be effectively suppressed by the proposed PSA. Compared with the existing schemes which assume perfect CSI, the sum SE of the network can be significantly increased by pilot reuse and the proposed RTD algorithm. The rest of this paper is organized as follows. In Section II, the signal transmission of an ultra-dense H-CRAN and the estimation of intra-cluster channels are presented. In Section III, a pilot scheduling algorithm is proposed to minimize the sum MSE of channel estimation. In Section IV, an alternative algorithm along with its distributed realization are provided to maximize the sum SE of the network by optimizing the transmit beam-vectors. Numerical results are presented in Section V before conclusions in Section VI.
This paper follows commonly used notations. R and C denote the real space and the complex space, respectively. The boldface upper (lower) case letters are used to denote matrices (vectors). I N stands for the N × N dimensional identity matrix and 0 denotes the all-zero vector or matrix. " \ " represents the set subtraction operation. Superscript (·) H denotes the conjugated-transpose operation and E{·} denotes the expectation operation. g is used to give the Euclidean norm of g.

A. Signal Transmission Model
Consider the downlink of a dense H-CRAN with an MBS, a BBU pool, K RRHs and M UEs as shown in Fig. 1, where each RRH connects with the BBU pool through an optical fiber and the MBS connects with the BBU pool through a backhaul link. The MBS and each RRH are respectively equipped with B and N antennas, and each UE has a single antenna. Denote the sets of RRHs and UEs by K and M, respectively. As discussed in [5], there are usually two types of clustering methods for RRHs to serve UEs namely disjoint clustering and user-centric clustering. In this paper, user-centric cluster method is adopted, i.e., each UE prefers to access the network via a selected subset of neighboring RRHs and different clusters for different UEs may overlap. If a UE cannot be served by any RRH, the MBS will offer network access to guarantee seamless coverage. For example, in Fig. 1, each RUE is served by RRHs inside the circle centered on this RUE. Hence, RUE 1, 2, 3, 4, 5 are served by RRHs, and BUE 1 and BUE 2 are served by the MBS. Denote the sets of RUEs and BUEs by M R and M B , respectively. Denote K i ⊆ K and M k ⊆ M R as the set of RRHs serving RUE i and the set of RUEs served by RRH k, respectively. Assume that all transmitters use the same time-frequency RB to transmit signals, leading to cochannel interference. Then, the received signal at UE m is given by where h k,m ∈ C N ×1 represents the channel vector from RRH k to UE m, w k,i ∈ C N ×1 denotes the beam-vector adopted by RRH k for transmitting signal to RUE i, and x i is the zero-mean unit-variance data symbol for RUE i. Likewise, h b,m , w b,j ∈ C B×1 and x j are similarly defined for the MBS and BUE j, and n m is the complex white Gaussian noise with variance N 0 , i.e., n m ∼ CN (0, N 0 ).

B. Channel Estimation
In most of the related literature, it is usually assumed that the global CSI over the network is available for optimization. To this end, orthogonal pilots should be adopted for channel estimation. However, in a dense H-CRAN network with large numbers of UEs, obtaining CSI of all links is almost infeasible due to limited training resources. As a result, in this paper, orthogonal pilots are adopted by BUEs and RUEs in the same RRH cluster for channel estimation, while pilot reuse among RUEs in different RRH clusters is allowed to reduce pilot overhead. For further pilot overhead reduction, the pilot sequence used by a BUE is allowed to be shared with RUEs. For example, in Fig. 1, BUE 1 and BUE 2 both access the network via the MBS. Hence, BUE 1 and BUE 2 adopt orthogonal training for channel estimation. RUE 1 and RUE 2 are both served by RRH 3 and RRH 4. Therefore, RUE 1 and RUE 2 belong to the same cluster and they should be assigned orthogonal pilots. Since the set of RRHs serving RUE 1 and the set of RRHs serving RUE 5 do not overlap, RUE 1 could reuse the pilot of RUE 5. In addition, due to the long distance between RUE 1 and BUE 2, they may be allocated the same pilot for pilot overhead reduction.
In the uplink training phase, assume that RRH k estimates the CSI from all RUEs in M k and the MBS estimates the CSI from all BUEs. The large-scale channel gains from all UEs to each RRH and to the MBS, i.e., {α k,m , ∀ k ∈ K, m ∈ M} and {α b,m , ∀ m ∈ M}, are assumed to be available at both the BBU pool and the MBS. Denote Q = {1, · · · , τ } as the available pilot set and Q = [q 1 , · · · , q τ ] ∈ C τ ×τ as the pilot matrix with orthogonal column vectors (i.e., Q H Q = I τ ). τ (|M B | ≤ τ ≤ M ) is the length of the pilots and is also the number of pilots available for channel estimation (this is the smallest amount of pilots that are required). Let P(M, Q) = {(m, π m )| m ∈ M, π m ∈ Q} denote an arbitrary pilot assignment scheme, where (m, π m ) means that UE m is allocated pilot q πm . In addition, let U π = {i| π i = π, ∀ i ∈ M R } and V π = {j| π j = π, ∀ j ∈ M B }, respectively, denote the sets of RUEs and BUEs that use pilot q π for channel estimation. Note that since BUEs apply orthogonal pilots for channel estimation, it follows that |V π | ∈ {0, 1}.
1) Channel Estimation for RUEs: Given the pilot assignment scheme P(M, Q), the N ×τ dimensional received signal matrix of pilots at RRH k can be written as where p R and p B are, respectively, the pilot transmit powers of RUEs and BUEs. N (R) k is the noise matrix which consists of independently and identically distributed (i.i.d.) Gaussian elements with zero mean and variance N 0 . Then, the minimum mean square error (MMSE) estimate of h k,i , ∀ k ∈ K, i ∈ M k is given by [38] Given the channel estimate vectorĥ k,i , the true channel vector h k,i can be expressed as h k,i =ĥ k,i +h k,i , where the error vectorh k,i represents the CSI uncertainty. Due to the property of MMSE estimation [38],h k,i is statistically independent of h k,i and it follows that CN (0, δ k,i I N ), where δ k,i is given by 2) Channel Estimation for BUEs: Similarly, the B × τ dimensional received signal matrix of pilots at the MBS can be written as where N (B) is the noise matrix which consists of i.i.d. Gaussian elements with zero mean and variance N 0 . Then, is statistically independent ofĥ b,j , where δ b,j is given by In general, sum SE is a very important metric in evaluating a wireless network's performance, and sum SE maximization has been widely studied in different kinds of networks [14], [15], [39]. Since pilot reuse is assumed in this paper, pilot contamination inevitably exists. Hence, how to effectively mitigate pilot contamination is also important. A two-stage optimization framework is thus studied in the following for network performance maximization. Specifically, the sum MSE of channel estimation is minimized by designing a pilot scheduling algorithm in Section III. Based on the pilot allocation result obtained in Section III, the sum SE of the network is then maximized by optimizing beam-vectors under imperfect CSI in Section IV.

III. STAGE I: PILOT SCHEDULING
In this stage, a pilot scheduling algorithm is designed to allocate pilots to UEs based on the metric of minimizing the sum MSE of channel estimation.

A. Problem Formulation
Since pilots are reused among UEs to shorten pilot overhead, pilot contamination inevitably exists. Considering the location dispersion of UEs, it is preferred that pilot contamination can be effectively mitigated by designing an appropriate pilot scheduling algorithm. According to (4) and (7), the sum MSE of all channel CSI is given by To distinguish the channels from different RUEs, it is assumed that any two RUEs for which the sets of serving RRHs are (partially) overlapping would need to be allocated orthogonal pilots. This can be mathematically expressed as follows Thus, the problem of minimizing the sum MSE of channel estimation can be formulated as Problem (10) is a resource allocation problem, which can be readily transformed to an equivalent mixed integer programming problem which is, however, usually difficult to solve. The optimal pilot scheduling scheme can be obtained through exhaustive search (ES). However, the complexity of ES increases exponentially with the number of UEs, becoming infeasible for a dense H-CRAN. Therefore, a low complexity pilot scheduling algorithm is proposed in the following subsection. RUE

B. Pilot Scheduling Algorithm
Constraint (9) indicates that any two RUEs served by at least one common RRH should be allocated different pilots. This constraint can be equivalently represented by an |M R |× |M R | dimensional matrix A with each element given by In matrix A, when two RUEs are served by at least one common RRH, the corresponding element is one. Otherwise, the element is zero. Obviously, to satisfy constraint (9), a minimum number of pilots t exists. In order to obtain t, an undirected graph can be constructed to describe constraint (9) based on A, where any two RUEs served by at least one common RRH are connected with each other. Then, determining t is equivalent to coloring the vertices of the undirected graph with the minimum number of different colors, which can be optimally solved by using the Dsatur algorithm proposed in [37]. Note that the pilots used by BUEs are mutually orthogonal. Therefore, the length of the pilots should satisfy max {|M B | , t} ≤ τ ≤ M . Taking the H-CRAN in Fig. 1 for example, the undirected graph consisting of all RUEs in Fig. 1 can be depicted as Fig. 2 (a), and the colored graph of Fig. 2 (a) after applying the Dsatur algorithm is shown in Fig. 2 (b). It can be seen that for the considered case, t = 3.
Since M = 7 and |M B | = 2, the pilot length should satisfy Though the Dsatur algorithm can yield the minimum number of pilots required by RUEs as well as a feasible pilot allocation solution, it only takes into account constraint (9) while ignores the objective function of problem (10). Hence, this solution may not be satisfactory. For example, in Fig. 2 (b), RUEs in the same color share the same pilot, and RUEs in different colors use orthogonal pilots for channel estimation. There may exist measurable pilot contamination between RUE 2 and RUE 4 since they are not so far away from each other. If τ = t, by exchanging the colors of RUE 2 and RUE 1 without changing t, the pilot contamination may be reduced due to the relatively longer distance between RUE 1 and RUE 4. If τ > t, the pilot contamination can be certainly decreased by replacing a reused pilot in Fig. 2 (b) with an unused one.
Since channel estimation is mainly affected by pilot contamination besides the effect of noise, it is thus of great importance to further mitigate pilot contamination. To this end, it is necessary to adjust the pilots allocated to each RUE by the Dsatur algorithm. As stated above, it is difficult to solve problem (10) in polynomial time. Hence, a heuristic lowcomplexity pilot scheduling algorithm is provided in the following. In order to implement this algorithm, it is first required to figure out how to measure the level of pilot contamination. In [31], a multi-cell massive MIMO network was considered and pilot reuse was assumed among different cells. To mitigate pilot contamination, a metric was defined to indicate the interference strength among UEs, and a graph coloring based scheme was then proposed. In [36], a continuous-valued metric was defined to evaluate the potential interference strength between any two D2D pairs when the same pilot was reused, and a pilot scheduling algorithm was then provided based on this metric. Motivated by these works, a similar metric based on large-scale channel gains is defined to measure the level of pilot contamination between an RUE and another UE when they reuse the same pilot. First, for any two unconnected RUEs in the undirected graph, define The definition of β i,i ′ is inspired by the channel estimation error in (4). Inside the ln(·) operation of (12), the second term is defined to measure the level of pilot contamination experienced by RUE i from RUE i ′ if they are assigned the same pilot, and vice versa for the third term. Obviously, a larger β i,i ′ means more severe potential pilot contamination between RUE i and RUE i ′ when they are assigned the same pilot. For any RUE i, let β i,i = 0. For any two connected RUEs, since they are assigned orthogonal pilots, there will be no potential pilot contamination between them. Hence, As for an RUE and a BUE, define Motivated by [31] and [36], a pilot scheduling algorithm (PSA) is proposed and summarized in Algorithm 1. According to Algorithm 1, it is necessary to first check whether τ is larger than or equal to the minimum required number of orthogonal pilots or not. If the number of available pilots are not enough to satisfy the scheduled constraints, τ has to be increased. Then, without loss of generality, pilot q j is allocated to BUE j, and {q 1 , · · · , q t } are randomly allocated to t clusters with RUEs in each cluster using the same pilot. In order to further mitigate pilot contamination, three steps are iteratively carried out to adjust the pilots assigned to each RUE. The basic idea is that the RUE experiencing larger pilot contamination possesses a higher priority for pilot adjustment. The main steps in each iteration can be explained as follows. First, RUE i ∈ M R \ Λ experiencing the largest interference is selected. Next, pilot π which causes the least interference to RUE i is chosen from the set of available pilots, i.e., Q \ X i . Finally, pilot q π is assigned to RUE i, and sets U πi , U π , Λ and X l , ∀ l ∈ M R are updated. The algorithm will be carried out for |M R | times until all RUEs' pilots have been adjusted.
Obtain matrix A from (11) and calculate t by using the Dsatur algorithm. If Assign pilot q j to BUE j, ∀ j ∈ M B . Divide all RUEs into t clusters by using the Dsatur algorithm, and randomly allocate pilots {q 1 , · · · , q t } to them with RUEs in each cluster using the same pilot. Let X i denote the set of different pilots allocated to the RUEs connected to RUE i, ∀ i ∈ M R . for d = 1, · · · , |M R | do Update X l , ∀ l ∈ M R . end for

C. Complexity Analysis
In this subsection, the computational complexity of Algorithm 1 is analyzed with order notation. According to [40], the Dsatur algorithm involves a complexity of O |M R | 2 . In order to adjust the pilots allocated to RUEs, |M R | iterations are carried out. In each iteration, a complexity of O |M R | 2 is required to find the RUE experiencing the largest interference, and the pilot causing the least interference to this RUE. The total complexity of the iteration process is thus O |M R | 3 . As a result, Algorithm 1 involves an overall complexity of O |M R | 3 . In contrast, to obtain the optimal pilot allocation solution, the ES scheme requires a complexity of O τ M , which increases exponentially with the number of UEs. Hence, the proposed algorithm is more efficient for practical use.

IV. STAGE II: ROBUST TRANSMISSION DESIGN
After obtaining the pilot allocation result by using Algorithm 1, the sum SE maximization problem under incomplete inter-cluster CSI and imperfect intra-cluster CSI is considered in this section. A robust transmission design algorithm and its distributed realization are proposed to solve the problem.

A. Problem Formulation
In order to simplify the expression of (1), the beamvectors from all RRHs in set K i for transmitting signal to RUE i can be merged to form a large-dimension vector represent the aggregated channel vector from all RRHs in K i to UE m. Then, (1) can be reformulated as In the following, all transmit beam-vectors are designed based on the obtained channel estimates, i.e., (3), (4), (6) and (7), and channel statistics, i.e., the large-scale channel gains. Consider the block fading model, where all channels remain unchanged over the coherence interval with length T . Then, the effective SINR and the achievable rate of RUE i are, respectively, given by , respectively, denote the aggregated channel estimation vector and the aggregated error vector from all RRHs in K i to RUE i. Since only MMSE estimates of the channel vectors and the distribution of channel estimation error are available, as in [41] and [36], the useful signal in (16) only contains ĝ H i,i w i 2 , and the terms corresponding to the channel estimation errors are regarded as Gaussian noise. The expectation operation in (16) is taken over the unknown channel estimation errorsh k,i , ∀ k ∈ K i , the inter-cluster channel vectors h k,i , ∀ k ∈ K \ K i and h b,i . Similarly, the effective SINR and the achievable rate of BUE j can be, respectively, written as where the expectation is taken over the unknown channel estimation errorh b,j and h k,j , ∀ k ∈ K.
Due to the fractional form of the SINR expressions and the log(·) operation, it is difficult to obtain explicit expressions of the achievable rate. In the following theorem, a lower bound on the achievable rate is derived.
Theorem 1: Given the SINR formulas in (16) and (17), the achievable rates R (R) i and R (B) j are, respectively, lower bounded by (18) and (19) as follows where (20) and (21) are all positive definite matrices defined in Appendix A. Proof: See Appendix A.
This paper aims to maximize the sum SE of the network by designing robust transmit beam-vectors under incomplete inter-cluster CSI and imperfect intra-cluster CSI. As discussed above, explicit expressions of achievable rates of both RUEs and BUEs are unavailable, while their lower bounds can be obtained according to Theorem 1. The tightness between achievable rates and their lower bounds is verified in the simulation part. Hence, R where w is the collection of all beam-vectors, including w k,i , ∀ k ∈ K, i ∈ M k and w b,j , ∀ j ∈ M B , and P (R) k and P (B) , respectively, denote the maximum transmit power of RRH k and the MBS.

B. Robust Transmission Design
Problem (22) has a form similar to the conventional sum SE maximization problems in downlink MISO systems, which can be effectively solved by adopting the weighted minimum mean square error (WMMSE) algorithm based on the following lemma [14], [15], [42]. For brevity, the proof of Lemma 1 is omitted.
Lemma 1: In a downlink interfering network, if MMSE receive filters are adopted for signal detection, the following relationship between the MMSE and the SINR of each link holds MMSE = 1 1 + SINR .
However, two factors make it difficult to directly exploit the WMMSE algorithm to solve problem (22). First, due to imperfect channel estimation, self-interference exists in the achievable rate expressions, i.e., terms j , which is different from the typical rate expression under perfect CSI; Second, due to the expectation operation, E b,j are all positive definite matrices rather than conjugate symmetric rank-one matrices as in popular rate expressions of a downlink MISO system. Because of these two factors, Lemma 1 cannot be directly applied. To deal with this difficulty, the following theorem is first introduced. (18) can be regarded as the rate of a mobile user in an equivalent downlink interfering MISO network. The MSE and single-tap MMSE receive equalizer of this user are, respectively, given by Similarly, r can also be regarded as the rate of a mobile user in an equivalent downlink interfering MISO network, and its MSE and single-tap MMSE receive equalizer are, respectively, given by Proof: See Appendix B. According to Lemma 1 and Theorem 2, r (R) i and r (B) j can be rewritten as and problem (22) 1 Note that in the objective function of problem (30), T −τ T is omitted and log(·) is replaced with ln(·) for the convenience of the following analysis. where Though the fractional SINR expressions have been avoided, problem (30) is still nonconvex and is generally difficult to solve. To make it tractable, the following auxiliary functions are introduced to remove the ln(·) operation in (30) where u (R) i and u (B) j are newly introduced auxiliary variables. Checking the first-order optimality condition of (31) yields and the corresponding optimal solutions Therefore, according to (32) and by substituting (31) into (30), the problem becomes where u = u . Compared with problem (30), problem (34) is much easier to solve since it is convex with respect to (w.r.t.) each of the individual variables. By alternatively optimizing w, f and u, a suboptimal solution of problem (34) can be obtained.
For fixed w and u, the optimal f can be obtained from (25) and (27). For given w and f , the optimal u can be obtained according to (33). When f and u have been determined, for notational brevity, denote β where From Theorem 1 and (36), it is known that both F (R) i and F (B) j are positive definite matrices. Hence, problem (35) is a quadratically constrained quadratic programming (QCQP), and can be optimally solved by adopting a standard convex optimization solver such as CVX, which is a toolbox developed in MATLAB for solving convex problems [43].

C. Distributed Implementation of the Proposed RTD Algorithm
Solving problem (34) using the proposed RTD algorithm requires a control center with great computation power, especially for an ultra-dense multi-cell network. Hence, it is desirable to obtain the beam-vectors in a decentralized manner with only local CSI, i.e., the BBU pool uses only to the BBU pool. Second, divide problem (35) into two subproblems with the first subproblem aiming to minimize the first term of the objection function of (35) subject to constraint (22b) and the second one aiming to minimize the second term of the objection function of (35) subject to constraint (22c). These two subproblems can be independently solved, and w (R) (d) as well as w (B) (d) can thus be respectively obtained at the BBU pool and the MBS. The BBU pool and the MBS then exchange the obtained w (R) (d) and w (B) (d) with each other, and respectively calculate (27) and (33). After checking the termination criterion, the algorithm stops if the algorithm converges. Otherwise, continue to the next iteration. The computational complexity of executing the proposed RTD algorithm is analyzed in the next subsection, from which it can be seen that, compared with the centralized way, executing the algorithm in a distributed manner can help spread out the compute task of the control center over the BBU pool and the MBS. This will help reduce the computational burden of the control center, especially for a network with multiple cells. The cost of the distributed implementation is the exchange and w (B) (d) between the BBU pool and the MBS. However, when distributed implementation is adopted, either the BBU pool or the MBS solves the corresponding problem with only local CSI, i.e., there is no need to collect the overall CSI, which improves the scalability. In addition, it is shown in Section V that the RTD algorithm converges rapidly within a few iterations. Hence, the distributed implementation is suitable for practical applications.

D. Convergence and Complexity Analysis
Since Algorithm 2 is carried out in an alternative manner, it is necessary to characterize its convergence behavior. In each iteration, the optimal w is first obtained by solving (35). Then, the optimal f and u are obtained according to (25), (27) and (33). As a result, the objective function of (34) decreases in each iteration. Due to the fact that this objective function is always lower bounded, the convergence of the proposed RTD algorithm is thus guaranteed.
Then, the computational complexity of the proposed RTD algorithm is analyzed. The complexity of this algorithm mainly lies in solving QCQP problem (35). As stated in Subsection IV-C, problem (35) can be divided into two subproblems, and according to [5], both of these subproblems can be equivalently transformed to a second-order cone programming (SOCP). The total numbers of variables in the two equivalent SOCP problems are, respectively, D 1 = i∈MR N |K i | and . When the RTD algorithm is carried out in a distributed way, the BBU poll and the MBS will, respectively, solve the two subproblems, and as the centralized way, L iterations are required for the algorithm to converge. Therefore, complexities of O LD 3.

V. SIMULATION RESULTS
In this section, representative simulation results are presented to evaluate the performance of the proposed algorithms. An isolated H-CRAN is considered with the MBS located at the center of the cell and all UEs randomly distributed. When a UE is close to the MBS, it usually prefers to access the network via the MBS. Hence, it is assumed that all RRHs are uniformly distributed in a ring area centered around the base station with the radius of the inner ring to be 200 m and the radius of the outer ring equal to the radius of the cell. Each RRH has a covering radius of D max m. If the distance between a UE and an RRH is within D max m, the UE chooses to access the network via this RRH. Assume that each RRH can simultaneously serve 3 UEs. If more than 3 UEs are associated with an RRH, the RRH will choose to serve the 3 closest UEs, and the rest UEs will be served by either the other RRHs or the MBS. If a UE cannot be served by any RRHs, the conventional cellular communication serves this UE. For the sake of brevity, equal maximum power constraint for all RRHs is assumed, i.e., P (R) k = P (R) , ∀ k ∈ K. Unless otherwise specified, the other system parameters are summarized in Table  I. All simulation results are obtained by averaging over 1000 channel realizations, and each channel realization is obtained by generating a random user distribution as well as a random set of fading coefficients.

A. Performance of the Proposed PSA alorithm
In this subsection, the performance of the proposed PSA alorithm is investigated in terms of sum MSE of the network. For comparison, the results obtained by the Dsatur algorithm and the ES scheme are taken as benchmarks. In particular, the Dsatur algorithm divides all RUEs into t clusters, and randomly allocates pilots {q 1 , · · · , q t } to RUEs with RUEs in each cluster using the same pilot. The ES scheme searches all feasible pilot allocation schemes and can always find the optimal solution with high calculation cost. Note that in the simulation process, since the UEs and RRHs are randomly generated, if the number of required pilots can not be satisfied for given τ , i.e., τ < max {|M B | , t}, set τ = max {|M B | , t}.  In Fig. 3, the sum MSE of channel estimation versus pilot length is depicted. It can be seen from this figure that the proposed PSA algorithm outperforms the Dsatur algorithm greatly in terms of sum MSE. When τ increases, as expected, the sum MSE obtained by the proposed PSA algorithm decreases with τ and approaches that obtained by the ES scheme. While for the Dsatur algorithm, the sum MSE remains unchanged since it always assigns t pilots to RUEs regardless of τ . Fig. 4 depicts the sum MSE of channel estimation versus the number of UEs. Compared with the Dsatur algorithm, the sum MSE of the network can be decreased greatly by the proposed PSA algorithm. As the number of UEs grows, more channels are required to be estimated and the probability of pilot reusing increases for given τ , leading to more pilot contamination. Hence, the sum MSE increases for all considered cases. When the number of UEs is equal to the pilot length, i.e., M = τ , orthogonal training is the optimal pilot allocation solution. In this case, there will be no pilot contamination, and channel estimation is only affected by thermal noise. Hence, the sum MSE will be very small, which can be found from both Figs. 3 and 4.

B. Performance of the Proposed RTD alorithm
In this subsection, the performance of the proposed RTD alorithm is evaluated.
First, in Fig. 5, the convergence behavior of the proposed RTD algorithm under different values of K is illustrated. It can be seen from this figure that the lower bound on sum SE of the network monotonically increases during the iterative procedure and converges rapidly after only a few iterations (within 8 iterations for all considered configurations). Fig. 5 also shows that the sum SE increases with the number of RRHs. This is because more UEs will access the network via RRHs and each RUE can be served by more RRHs as K grows. Due to short transmission distance between RUEs and RRHs, SE gains can be obtained in contrast to the conventional cellular communication.
Since explicit expressions of achievable rate of both RUEs and BUEs are unavailable, this paper aims to maximize the lower bound on sum SE of the network. It is thus necessary to verify the feasibility. In Fig. 6, the gaps between the achievable  rate of RUEs, BUEs and their corresponding lower bounds are investigated. Fig. 6 shows that the achievable sum rates of RUEs and BUEs are both close to their corresponding lower bounds, indicating that it is reasonable to solve the sum SE maximization problem based on lower bounds (18) and (19). In addition, Fig. 6 also shows that the sum SE increases with RRH and MBS antennas, which is consistent with intuition.
In Fig. 7, the effect of pilot length under different values of coherence interval is investigated. It is shown that for all considered cases, the sum SE of the network first grows and then decreases w.r.t. τ . This is because only a few orthogonal pilots are reused among UEs for a small τ . In this case, the channel estimation is significantly influenced by pilot contamination. Therefore, the sum SE can be enhanced by increasing τ . However, as τ becomes large enough, the channel estimation accuracy can be hardly improved by further enlarging τ . Counterproductively, increasing pilot length reduces the number of symbols available for data transmission, and thereby reduces the sum SE. Note that in Fig. 7  sum SE of the network can be significantly increased by pilot reuse and the proposed RTD algorithm, especially when T is small. Fig. 8 depicts the sum SE of the network versus the number of UEs under different values of K. For comparison, the results obtained by using Algorithm 2 proposed in [5] are depicted as benchmarks, which assume that perfect global CSI is available. Note that reference [5] considered a C-RAN without MBS. Hence, for the sake of fairness, the same H-CRAN network is considered when obtaining the benchmarks. From Fig. 8, several observations can be made. First, as expected, the sum SE grows with K for all considered cases. Second, as the number of UEs increases, the sum SE grows monotonically for the proposed RTD algorithm, while first increases and then slightly decreases for Algorithm 2 proposed in [5]. This is because orthogonal training is used by [5] to obtain the perfect CSI. As M grows, the number of required pilots increases, and thereby decreases the sum SE. In addition, it can also be seen that compared with Algorithm 2 proposed in [5], the sum SE of the network can be significantly improved by the proposed RTD algorithm.
In Fig. 9, the effect of covering radius of RRHs on sum SE of the network is investigated. As D max increases, more UEs will access the network via RRHs and each RUE can be served by more RRHs. Hence, for the cases with M = 13 and M = 16, the sum SE of the network grows with D max for both the proposed RTD algorithm and Algorithm 2 proposed in [5]. However, when M = 10, the sum SE first decreases and then increases for both the proposed RTD algorithm and the benchmark. This is due to the fact that as D max increases, channel conditions between RRHs and RUEs become worse. Though coverage areas of RRHs become large, the sum SE gains brought by this are limited when M is small. Hence, the sum SE of the network suffers a decrease.

VI. CONCLUSIONS
This paper has studied pilot scheduling and robust transmission design problems in an ultra-dense H-CRAN. Since pilot reuse was assumed among UEs to shorten pilot overhead, pilot contamination inevitably exists. Hence, a pilot scheduling algorithm was proposed to minimize the sum MSE of all channel estimates. Afterwards, robust transmit beam-vectors were designed to maximize the sum SE of the network. Since each RRH or the MBS only has the imperfect CSI of intra-cluster UEs and only tracks the large-scale channel gains of inter-cluster UEs, it was difficult to obtain the exact achievable rate of each link. Hence, a lower bound on each UE's achievable rate was derived and an alternative robust transmission design algorithm was proposed to maximize the lower bound on sum SE. Simulation results showed that compared with existing algorithms, the system performance can be significantly improved by the proposed algorithms in terms of both sum MSE and sum SE.

APPENDIX A PROOF OF THEOREM 1
Using the convexity of log 2 1 + 1 x (∀ x > 0) and applying the Jensen's inequality, a lower bound on R (R) i can be derived as (37) at the bottom of this page. Denote From the definition ofg i,i , it is known thatg i,i consists of i.i.d. Gaussian elements with zero mean and variance given by (4). Hence, To obtain the explicit expression of G (R) i ′ ,i , denote the set of RRHs serving RUE i ′ as K i ′ = k i ′ 1 , · · · , k i ′ |K i ′ | . Then, according to the definition of g i ′ ,i , G (R) i ′ ,i can be rewritten as , · · · , |K i ′ |} represents the oth row and zth column block matrix of G (R) i ′ ,i . From Stage I it is known that for any RUE i, only channel vectors h k,i , ∀k ∈ K i are estimated, while h k,i , ∀k ∈ K \ K i are unknown. In addition, different channel vectors are independent with each other. As a result, G otherwise. (41) Since for any RUE i, h b,i is not estimated, it follows that Substituting (39), (40) and (42) into (37), the lower bound (18) can be obtained. Similarly, for BUE j, a lower bound on R (B) j can also be derived as (43) shown at the bottom of this page. Denote Since for any BUE j, only h b,j is estimated, while channel vectors g i,j , ∀ i ∈ M R are unknown, E Substituting (45) into (43), the lower bound (19) can be obtained.
Since large-scale channel gains and the variance of channel estimation errors, i.e., (4) and (7)  In order to simplify the proof, a special case is first considered. Denote where g, w, v ∈ C C×1 , E = diag {θ 1 , · · · , θ C } ∈ R C×C + and N 0 ∈ R ++ . G ∈ C C×C is a positive definite Hermitian matrix and admits an eigen-decomposition G = U ΓU H , where U is a C × C dimensional unitary matrix with each column being the eigenvector of G and Γ is a diagonal matrix whose diagonal elements are the corresponding eigenvalues of G, i.e., Γ c,c = γ c [45]. θ c and γ c , ∀ c ∈ {1, · · · , C} are all positive real constants. Let e c denote a C dimensional vector with one in the cth position and zeros elsewhere. E and Γ can then be rewritten as E = [θ 1 e 1 , · · · , θ C e V ] and Γ = [γ 1 e 1 , · · · , γ C e V ], respectively. Denotev = U H v. It follows that (47) Obviously, (47) shows that r can be seen as the rate of a mobile user in an equivalent downlink interfering MISO network. The received signal of this user is given by where the first term denotes the desired signal, and the second to the third terms are the interference from other users. g, √ θ c e c and √ γ c e c denote the channel vectors from all transmitters to the desired receiver. w and v are the transmit beam-vectors. x 0 , x c and s c represent the i.i.d. data symbols with zero mean and unit variance. n is the complex white Gaussian noise with variance N 0 , i.e., n ∼ CN (0, N 0 ). Thus, the MSE of the considered link can be written as and the corresponding MMSE receiver which minimizes (49) is given by From (39), (40), (42) and (45), it is known that E b,j are both positive definite Hermitian matrices as G in (46). Analogously, Theorem 2 can be readily verified as above.