{'fraction_non_alphanumeric': 0.03493494258141726, 'fraction_numerical': 0.01018511445522369, 'mean_word_length': 4.494404029099049, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 1, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 8, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Block Distributed Compressive Sensing Based Doubly Selective Channel Estimation and Pilot Design for Large-Scale MIMO Systems

The doubly selective (DS) channel estimation in the large-scale multiple-input multiple-output (MIMO) systems is a challenging problem due to the large number of the channel coefficients to be estimated, which requires unaffordable and prohibitive pilot overhead. In this paper, firstly we conduct the analysis about the common sparsity of the basis expansion model (BEM) coefficients among all the BEM orders and all the transmit-receive antenna pairs. Then a novel pilot pattern is proposed, which inserts the guard pilots to deal with the inter carrier interference (ICI) under the superimposed pilot pattern. Moreover, by exploiting the common sparsity of the BEM coefficients among different BEM orders and different antennas, we propose a block distributed compressive sensing (BDCS) based DS channel estimator for the large-scale MIMO systems. Its structured sparsity leads to the reduction of the pilot overhead under the premise of guaranteeing the accuracy of the estimation. Furthermore, taking consideration of the block structure, a pilot design algorithm referred to as block discrete stochastic optimization (BDSO) is proposed. It optimizes the pilot positions by reducing the coherence among different blocks of the measurement matrix. Besides, a linear smoothing method is extended to large-scale MIMO systems to improve the accuracy of the estimation. Simulation results verify the performance gains of our proposed estimator and the pilot design algorithm compared with the existing schemes.


I. INTRODUCTION
Large-scale multiple-input multiple-output (MIMO) [1], [2] attracts much academic interest and is considered as a promising technology in the incoming 5G cellular systems [3]. It enhances the data throughput and improves the link reliability of wireless communication systems by taking advantage of the high spatial multiplexing gains. In order to benefit from largescale MIMO, one must obtain accurate channel state information (CSI) which guarantees data recovery and contributes to multi-antenna array gains.
Time and frequency selective channel, which is also referred to as doubly-selective (DS) channel, is related to many wireless access links, such as high-speed trains [4] and millimeterwave communications [5], [6]. Frequency selectivity is caused by multipath propagation and time selectivity results from Doppler shift. For the DS channel estimation, inter-carrier interference (ICI) is a challenging problem, which incurs a large amount of channel coefficients to be estimated and the high complexity of the estimation schemes. In [7]- [9], basis expansion model (BEM) was proposed to simplify the estimation process. The work [7] presented the research about the optimal training for DS channel estimation. The work [8] verified several channel estimation schemes including the least squares (LS) estimator, the linear minimum mean square error (LMMSE) estimator and the best linear unbiased estimator (BLUE) combined with different BEM basis in single-input single-output (SISO) systems. And then in [9] they were extended to MIMO systems with 2 transmit antennas. However, in large-scale MIMO systems with more than a dozen or even dozens of antennas, the number of channel coefficients to be estimated increases largely. The conventional schemes are not feasible since it requires unaffordable pilot overhead and prohibitive complexity to avoid the inter-antenna interference (IAI) and process the ICI. To our best knowledge, the DS channel estimation in large-scale MIMO systems has been seldom considered in the existing works.
Compressive sensing (CS) is an important framework to decrease the pilot overhead and the complexity of channel estimation by taking advantage of the channel sparsity [10]. The channel sparsity in three domains are concerned in the existing literature, including the delay-Doppler domain [11], the beam domain [12], and the delay domain [13]. In delay-Doppler domain, a large Doppler shift incurs the leakage effects, which increases the channel sparsity and deteriorates the performance. Processing the leakage effects and enhancing the sparsity require high complexity of computation [14]. The sparsity in beam domain relies on an open and wide propagation environment with few scatters. It is usually increased in the environment with rich reflection and the compressive sensing framework is no longer applicable. The most commonly exploited sparsity presents in delay domain in the existing CS based channel estimation schemes. It is analyzed that the broadband channels present sparsity in delay domain.
Various CS based channel estimation schemes appear recently. In [15], the authors proposed compressive estimation schemes for flat fading channels in large-scale MIMO systems. It analyzed the channel sparsity among all the antennas and jointly estimated the channels of all the receive antennas in the base station. The works [16]- [20] proposed compressive CSI estimation schemes under frequency selective channels. In [16], the authors considered the channel model with the common support and the individual support. Then they proposed an algorithm which adapted to this kind of specific structure for large-scale MIMO systems. The work [17] proposed a scheme for the uplink channels in large-scale MIMO systems, which was relatively simple since it wasn't concerned with IAI. In [18], an approach with an unknown sparsity was proposed, which was more realistic in the practical systems. The works [19], [20] suggested different block-structured pursuit algorithms. In [4], the authors presented the research about DS channel estimation based on CS in SISO systems. It proposed a channel estimation scheme based on the position information in high mobility systems and optimized the corresponding pilot design. The work is extended to MIMO systems with 4 transmit antennas in [11], which is referred to as a low coherence compressed (LCC) channel estimation scheme. However, on the one hand, this scheme utilized the sparsity in delay-Doppler domain which is notably increased by a large Doppler shift, and on the other hand, it could not support much more antennas. The work [13] introduced distributed compressive sensing (DCS) to DS channel estimation in SISO systems. It utilized the channel sparsity in delay domain and formulated the estimation into a DCS framework, which guaranteed a more accurate recovery.
Mutual coherence is an important factor concerned with the accuracy of the recovery in CS framework [10]. In SISO systems, several algorithms were proposed to decrease the mutual coherence of the measurement matrix. In [13], it proposed discrete stochastic optimization (DSO) to select the suboptimal pilot positions and [11] proposed an algorithm to jointly optimize the pilot values and the pilot positions. In MIMO systems, [21] and [22] employed the optimization algorithms based on the stochastic search and the genetic algorithm. However, they were designed for orthogonal pilots, which meant the requirement of a large pilot overhead in large-scale MIMO systems. The works [12] adopted the nonorthogonal pilots with equispaced pilot positions for the consideration of reducing the correlation among different virtual channels.
In this paper, we propose a block distributed compressive sensing (BDCS) based DS channel estimation scheme for the large-scale MIMO systems and a novel pilot design algorithm corresponding to the unique structure. They reduce the pilot overhead under the premise of guaranteeing the estimation accuracy. In specific, firstly we analyze the common sparsity of the BEM coefficients among all the BEM orders and all the transmit-receive antenna pairs in delay domain. Different from the sparsity analysis in the existing literature [18]- [20], we focus on the sparsity of the BEM coefficients rather than the channel coefficients. In addition to the analysis of the common sparsity among different BEM orders for SISO systems in [13], the common sparsity among different transmit antennas is considered here as well; Then, a novel pilot pattern is proposed for DS channels in large-scale MIMO systems. It combines the property of the superimposed structure and the setting of the guard pilots. The superimposed structure reflects the superimposed pilot positions of different transmit antennas, which reduces the pilot overhead. The guard pilots are designed for the ICI of the DS channels, which avoid the contamination from the data subcarriers; Moreover, we propose a BDCS based channel estimator for DS channels in large-scale MIMO systems. The unknown BEM coefficients present block and common sparsity simultaneously. The structure benefits the localization of the nonzero elements, and then leads to the performance improvement and the pilot overhead reduction; Furthermore, taking advantage of the analyzed block sparsity, a novel pilot design algorithm, referred to as block discrete stochastic optimization (BDSO), is proposed. It also contributes to the performance gain and the spectral efficiency.
The remainder of this paper is organized as follows. Section II introduces the system model and the fundamentals of CS. In section III, we present the proposed channel estimation scheme and the proposed pilot design algorithm. In Section IV, we conduct the analysis of the complexity. In section V, simulation results verify the validity of our work. Section VI concludes this paper.
[A] S represents the selected rows of A, whose indices correspond to the set S. CN (0, σ 2 ) represents the complex Gaussian distribution with zero mean and σ 2 variance. I x means the identity matrix of order x. vec(A) denotes the column-ordered vectorization of matrix A.

II. SYSTEM MODEL AND FUNDAMENTALS
In this section, we introduce our system model, which includes the transmission model and the channel model of the large-scale MIMO-orthogonal frequency division multiplexing (OFDM) systems. Besides, the fundamental knowledge of CS and structured compressive sensing (SCS) is briefly illustrated.
A. System Model 1) Transmission Model: We consider a large-scale MIMO-OFDM system. The base station is equipped with a great many antennas. It serves a number of terminals with a single antenna. As shown in Fig. 1, the antenna array is arranged in a rectangle, which consists of N B antennas. We adopt frequency division duplex (FDD) mode in our system. In the OFDM system, for any transmit-receive antenna pair, there exist N subcarriers in a parallel transmission. A part of the subcarriers are selected as pilot subcarriers to estimate the channel coefficients and the remaining ones are responsible for data. (1) Assume that H (nB ) t ∈ C N ×N describes channel matrix in time domain, and its elements can be expressed as (2) The channel matrix in frequency domain H in which W is the discrete fourier transform (DFT) matrix and W [m, n] = N −1/2 exp (−j2πmn/N ), m, n ∈ [0, N − 1], It is found that for DS channels, we have to estimate the channel coefficients of each channel tap at each time instant. The total number of coefficients to be estimated for a transmit-receive antenna pair is N L. N is the number of the samplings in time domain, which is equal to the number of the subcarriers, and L is the number of multipaths, i.e., the length of the channel. BEM [9] is an important technique for DS channel estimation, which is always introduced to reduce the number of the coefficients to be estimated. Let h denote the channel coefficients of the l-th channel tap and the n B -th transmit antenna. Each h in which, is the BEM coefficients, and is the number of the BEM orders. Apparently the number of channel coefficients to be estimated is reduced from N L to DL for one transmit-receive antenna pair. The vector of the channel taps for the n B -th transmit antenna can be formulated as in which with represents the channel coefficients at the (n + 1)-th time instant of the n B -th antenna. θ whereΘ T as its first column [9]. Due to its circularity,Θ (nB ) d can be diagonalized as where W L denotes the submatrix that extracts the first L columns of W. Accordingly, substitute (8) into (7) and the time domain channel matrix can be denoted as (9) Substituting (9) into (3), it is not hard to find that the channel matrix in frequency domain can be expressed as in which, and ∆ (nB) is the modeling error [8].

B. CS and SCS
In this part, the basic knowledge of CS and SCS are introduced. The SCS means that the sparsity presents a certain structure, including DCS, block compressive sensing (BCS), and BDCS here.
1) CS: CS is an attractive framework, which recovers a high-dimensional sparse signal from a low dimensional observed vector. It solves the underdetermined problem in which x ∈ C Z×1 is an unknown high-dimensional vector, represents the observed low-dimensional vector, and e denotes the noise term. The theory of CS is based on two important premises: • The first one is the sparsity of the high-dimensional vector, which means that x is a sparse vector with sparsity K, i.e., x 0 = K, K ≪ Z. • The second one is that the measurement matrix A satisfies restricted isometry property (RIP) condition [10]. If these two conditions are satisfied, a high probability of the exact recovery of x can be guaranteed. But it should be noted that it is difficult to verify RIP condition due to the prohibitive complexity and the tremendous computation. In practical schemes, mutual coherence property (MCP) [23] is an important reference value of the measurement matrix, which reflects the coherence between columns.
where a i and a j denote the i-th and the j-th columns of A. Lemma 1 ( [24]): Suppose that A has MCP µ and that the sparsity of x is K with K < ( 1 / µ +1)/4. Furthermore, suppose that we obtain measurements of the form r = Ax + e. Then when the set of solutions B(r) = {z : Az − r 2 ≤ ζ}, ζ is a constant value, the solutionx obeys Lemma 1 verifies that a smaller value of the MCP will lead to a more accurate recovery of x. Basis pursuit (BP) [25] and orthogonal matching pursuit (OMP) [26] are the widely adopted recovery algorithms of CS.
2) SCS: Several models of SCS are elaborated as follows. We summarize the properties of the problems and the corresponding recovery algorithms.
• DCS: Many applications are concerned with a problem that several high-dimensional sparse vectors with the same positions of nonzero elements are compressed by a common measurement matrix, It is inefficient to recover each sparse vector separately. DCS framework is applied to jointly compress and recover the multiple correlated sparse signals. The basic form of DCS is in which, is the noise matrix. All the columns of X share the same nonzero positions. The recovery accuracy also relies on the common sparsity of X and the property of the measurement matrix A as CS above. As proved in [27], DCS provides higher accuracy of the recovery with fewer observed values than CS by utilizing the common sparsity. The reason is that multiple vectors contribute to the localization of the nonzero elements. Simultaneous-OMP (SOMP) [27] is an important algorithm for the recovery of DCS. • BCS: The BCS means that the sparsity of the unknown high-dimensional vector presents block sparsity. The basic form of the problem is in which, for better illustration of the block sparsity of x ′ , it is decomposed as and accordingly the measurement matrix A can be decomposed as We can see that the unknown high-dimensional vector x ′ is decomposed into T parts. The K sparsity indicates that K of them are nonzero blocks while the remaining T −K parts are zero blocks, K ≪ T . The nonzero block means that all the elements of the block are nonzero, and the zero block is constituted by d zero elements. We set an example for further clarification about the block sparsity with assumption of T = 6 and K = 2 as in which, x 2 , x 5 are the nonzero blocks while x 1 , x 3 , x 4 , x 6 are the zero blocks. The recovery algorithm, block orthogonal matching pursuit (BOMP), is proposed in [22] to solve the problem of BCS. The block sparsity is utilized to improve the accuracy of the nonzero elements localization and contributes to the performance gain.
• BDCS: The BDCS combines the property of the DCS and the BCS. The structure includes both the block sparsity of each unknown high-dimensional vector and the common sparsity among different unknown highdimensional vectors. The form of the BDCS is described as Similar to BCS, the measurement matrix A is decomposed as and the unknown matrix X ′ can be decomposed into T × J parts as in which, The block sparsity is reflected on the j-th column of X ′ , j ∈ [1, J]. K out of the T parts are nonzero blocks while the remaining ones are zero blocks. The common sparsity means that the position of the K nonzero blocks are the same among all the J columns of X ′ . We further explain the block and common sparsity with an example under the assumption of K = 2, T = 6, J = 3.
The recovery algorithms, structured subspace pursuit (SSP) [28] and block simultaneous orthogonal matching pursuit (BSOMP) [29], exploit both the common sparsity and the block sparsity to get the more accurate recovery performance.

III. THE PROPOSED ESTIMATOR
In this section, a novel DS channel estimator is proposed for large-scale MIMO systems by extending the DS channel estimation scheme in SISO systems [13]. In the process of extension, firstly we propose a novel pilot pattern which additionally considers the superimposed pilot structure of all the antennas. Then the estimator is formulated as a BDCS based problem which is different from the DCS based problem for SISO systems. Furthermore exploiting the block sparsity, the pilot design algorithm optimizes the coherence among different blocks of the measurement matrix in comparison to the different columns in SISO systems.

A. The common sparsity of BEM coefficients
In this part, we analyze the common sparsity of the BEM coefficients for large-scale MIMO systems. In addition to the common sparsity among different BEM orders, which is discussed in [13] for SISO systems, we consider the common sparsity among different antennas here as well.
Theorem 1: The elements of the BEM coefficients set { θ , in a large-scale MIMO system share common sparsity among different BEM orders and different transmit antennas under the condition of smax C ≤ 1 10BW , in which s max denotes the maximum distance between any two transmit antennas, C is the speed of light and BW is the signal bandwidth.
Proof: The proof is conducted by three steps. Firstly we analyze the common sparsity of the channel coefficients in delay domain among all the sampling instants. Then the common sparsity of the channel coefficients among different antennas is illustrated. Finally, exploiting the relationship between the channel coefficients and the BEM coefficients, we prove the common sparsity of the BEM coefficients among different BEM orders and different antennas. Assume that there are L channel taps for a transmitreceive antenna pair and the index set of them is denoted as for any l / ∈ {l 1 , · · · , l K }. It is not hard to find that { h and their common nonzero positions in delay domain are {l 1 , · · · , l K }.

Lemma 3 ( [15]):
In large-scale MIMO systems, all the transmit-receive antenna pairs share common sparsity in delay domain if smax C ≤ 1 10BW , in which s max denotes the maximum distance between any two transmit antennas, C is the speed of light and BW is the signal bandwidth.
As referred in [15], if smax C ≤ 1 10BW , the links of all the transmit-receive antenna pairs scatter invariantly in space. Thus the indices of their strong channel taps are the same, i.e., {l 1 , · · · , l K } ⊂ [0, L − 1]. In another word, they share common sparsity in delay domain. It is safe to assume that all the transmit-receive antenna pairs have the same nonzero positions of the channel taps in a large-scale MIMO system, which means that the elements of { h have common sparsity among different antennas. As to the condition of smax C ≤ 1 10BW , in the long term evolution (LTE) systems [5] with parameters of BW=20MHz and the center frequency of 2.6GHz, we have that the 10×10 transmit antenna array has common channel support [15].
Combined Lemma 2 and Lemma 3 together, we conclude that the elements of { h The linear relationship reflects that θ have common sparsity among different BEM orders and different antennas.

B. The Proposed Pilot Pattern
In this part, a novel pilot pattern is proposed for DS channels in large-scale MIMO systems, which combines the guard pilot design and the superimposed structure among different antennas to combat the ICI and reduce the pilot overhead.
Complex exponential basis expansion model (CE-BEM) is adopted due to its simple form. We have that the basis function of the CE-BEM v d = (1, · · · , e j 2π Fig. 2, we present the geometric expression of the CE-BEM for better illustration. A square, consisting of N rows and N columns, denotes the wireless channel between a transmit-receive antenna pair. The N rows represent the N transmitted subcarriers and the N columns correspond to the N received subcarriers. The intersection of the i-th row and the i-th column means the channel coefficient of the i-th subcarrier, i ∈ [1, N ] while the intersection of the i-th row and the j-th column means the interference to the j-th subcarrier from the i-th subcarrier, j ∈ [1, N ], j = i. Under the function of the CE-BEM with order D, the square is reduced to D diagonals as shown by the red lines in Fig.  2. To combat the ICI, which is represented by the D − 1 subdiagonals, the nonzero pilots are always accompanied by several zero guard pilots on both sides. It is derived in [7] that the optimal number of guard pilots on one side is D − 1. Besides, in [8], we can see that only the central D pilots including the nonzero one are not contaminated by the data subcarriers. To clearly elaborate the motivation of our designed pilot pattern, we introduce two existing pilot patterns as shown in Fig. 3(a) and Fig. 3(b).
• In Fig. 3(a), we present a pilot pattern proposed for DS channels in SISO system [13]. Assume that the pilots are arranged in G groups. For each group, a nonzero pilot is equipped with D − 1 zero guard pilots on each side. It selects the D central pilots of each group for channel estimation and obtains an ICI-free structure.  • Fig. 3(b) depicts a superimposed pilot pattern, which is proposed for frequency selective channels in large-scale MIMO systems [28]. The nonzero pilots occupy the same positions for each antenna to reduce the large pilot overhead brought by the increased number of antennas. The sequence of the nonzero pilots for each antenna consists of random ±1 and different antennas are distinguished by different sequences. In another word, it utilizes the differentiation of different antennas in code domain. Considering that our pilot pattern is designed for DS channels in large-scale MIMO systems, we combines the above two properties together. As shown in Fig. 3(c), the guard pilots are inserted in the superimposed structure. It gives the consideration to both the ICI avoidance and the reduction of the pilot overhead. To our best knowledge, little is concerned about the pilot pattern for DS channels in large-scale MIMO systems. In the existing literature [9], the orthogonal pattern with guard pilots, as shown in Fig. 3(d), is prepared for DS channels in MIMO systems. It cannot support more antennas since the requirement of pilot subcarriers is too large.

C. The Proposed Estimator
In this part, we propose a BDCS based estimator for DS channels in large-scale MIMO systems by properly extending the DCS based channel estimation scheme in SISO systems.
1) Review of the DCS based estimator in [13]: The work [13] proposed a DCS based estimator for DS channels in SISO systems. Its pilot pattern is depicted as Fig. 3(a). The pilots are divided into G groups and for each group, the central D pilots are selected for channel estimation to guarantee the ICIfree structure. Let S d , d ∈ [0, D − 1], represent the index set of the (d + 1)-th selected pilot subcarriers of all the groups. The index set of the nonzero pilots S D−1 2 is also denoted as S cen . Thus we have Let P SISO ∈ C G×1 denote the values of the nonzero pilots. The channel estimator, which describes the relationship between the received pilots and the CE-BEM coefficients, is derived as Lemma 4.
in which, S SISO is the transmitted signal and Υ SISO is the noise term. Combining with the formula (10), the CE-BEM decomposition of the received signal is derived as in which,S SISO =diag(S SISO ) and I <d− D−1 2 > N means that the identity matrix with order N shifts down circularly for , the received pilots can be expressed as (32) Substituting (32) into (31), it is not hard to obtain (28).
Thus, it is found that As analyzed above, the elements of {θ SISO d }, d ∈ [0, D − 1] have common sparsity among different BEM orders. Thus (33) is a typical DCS problem and can be solved by the SOMP algorithm.
2) The proposed BDCS based estimator: Firstly we introduce a key finding that the CE-BEM coefficients present both the common sparsity and the block sparsity. Then exploiting the structured sparsity, we propose a BDCS based estimator for DS channel estimation in large-scale MIMO systems.
The pilot pattern of the large-scale MIMO system is described in Fig. 3(c). For each transmit antenna, the positions of the pilot subcarriers are arranged as the same way as the SISO systems. Besides, among different antennas, we adopt the superimposed structure, in which, the pilot subcarriers of all the antennas have the same indexes.
Let P (nB ) ∈ C G×1 represent the pilot values of the n B -th antenna, n B ∈ [1, N B ]. Considering that the received signal is the addition of the signals transmitted by all the antennas on the base station, we derive the relationship between the received pilots and the CE-BEM coefficients of all the transmitreceive antenna pairs as Theorem 2.
Theorem 2: The estimation of the CE-BEM coefficients of all the transmit-receive antenna pairs {θ , is conducted by (34), in which, P (nB ) =diag(P (nB ) ).
Proof: Considering that the transmitted signals of different antennas overlap with each other in the air, we have that the received signal is in which, S (nB ) is the transmitted signal of the n B -th antenna, S (nB ) =diag(S (nB) ), n B ∈ [1, N B ], and Υ is the noise term. Then the pilot subcarriers are selected as . Thus (34) is obtained as the complete expression.
For convenience of observation, we reshape (34) in a compact form as in which, and For better illustration, Fig. 4 depicts an example of the CE-BEM coefficient matrix Λ with the channel length L = 5 and the sparsity K = 2. It is obvious that Λ has common sparsity among its columns. Besides, when we adjust the rows of Λ in the order of the channel tap, it presents the nonzero blocks and the zero blocks, referred to as the block sparsity. It can be concluded that the CE-BEM coefficient matrix Λ have both common and block sparsity since the elements of {θ have common sparsity not only among different CE-BEM orders, but also among different antennas. Now we can see that the proposed estimator (34) attributes to a BDCS problem, which can be solved by SSP or BSOMP algorithms.
The work [30] proposed an estimator for DS channels in large-scale MIMO systems as well, which can only be performed in time frequency training (TFT)-OFDM systems. It utilizes the correlation of the Pseudo-random sequence for the nonzero elements localization and is not suitable for the conventional cyclic prefix (CP)-OFDM.

D. The Proposed Pilot Design Algorithm
In this part, we propose a novel pilot design algorithm to optimize the pilot positions S cen for better performance of estimation. Considering the block structure, it decreases the coherence among different blocks of the measurement matrix, which is different from the MCP of the pilot design algorithm for SISO systems.
As analyzed above, the CE-BEM coefficient matrix Λ presents block sparsity when its rows are arranged in the order of the channel tap. For more accurate localization, the columns of the measurement matrix Z are also arranged in the order of the channel tap accordingly. Assume that A l = {l + 1, l + 1 + L, l + 1 + 2L, · · · , l + 1 + (N B − 1)L}, l ∈ [0, L − 1], and we have Apparently, Z a consists of L blocks and the coherence among them affects the localization of the nonzero blocks. Assume that and µ(Z s ) takes on the reference value of the coherence among difference blocks. In our proposed Algorithm 1, referred to as BDSO, we search for a proper nonzero pilot positions S D−1 2 in (38) with a small µ(Z s ) by the iteration process. The detailed steps are summarized as follows.

end end
We initialize the values of the nonzero pilots P (nB ) ∈ Z G×1 , n B ∈ [1, N B ] with random sequences consisting of ±1, and then generate S cen = { s 1 , s 2 , . . . , s G } randomly, which The condition guarantees the groups of pilots do not overlap with each other. In the reformulation stage,S cen is obtained by changing an element of S cen andZ and Z are reformulated accordingly as (38). Then we convert the form of the matrix Z and Z toZ s and Z s as (40) and (41) for the convenience of calculating the coherence. Thus the S cen is updated according to the comparison between µ(Z s ) and µ( Z s ). Finally, in the state transition stage, the probability of each state is calculated. Comparing the current and the last states, the S cen which corresponds to the larger probability is selected. As the iteration process converges, we obtain the optimized S cen with a steady solution.
In contrast to the pilot design algorithm [13] proposed for the SISO systems, referred to as DSO, we pay attention to the coherence among different blocks of the measurement matrix instead of the columns, which avoids a large amount of unnecessary search. The works [21] and [22] proposed the pilot design algorithms for MIMO systems as well, which utilized the genetic and the iterative search to optimize the pilot positions. However, they were designed for the orthogonal pilot patterns and assigned different pilot positions to different antennas. Thus, a large number of pilot subcarriers were required and the estimator could only support a maximum of 4 antennas.

E. Linear Smoothing
In order to relieve the performance deterioration brought by the modeling error of the CE-BEM, our previous work [29] proposed a linear smoothing method for SISO systems. In specific, we approximate the channel coefficients corresponding to the N time instants with a straight line. Here we extend it to large-scale MIMO systems by repeating the smoothing process for each transmit-receive antenna pair.
• Calculate the slope of the approximate line determined by the two central points of each h • The processed channel can be expressed aŝ Remark: As described above, we mainly focus on the discussion about the downlink DS channel estimation in the large-scale systems, which is a challenging problem involving the IAI. For the integrity of the work, the uplink DS channel estimation in the large-scale MIMO systems is briefly introduced, although it is relatively simple and not concerned with the IAI. The estimator can be expressed as (42), in which, Y ′ (nB ) represents the signal received by the n B -th antenna on the base station,P ′ =diag(P ′ ), P ′ ∈Z G×1 denotes the pilot sequence,θ It is found that the CE-BEM coefficients presents common sparsity while without block sparsity, which is different from the downlink channel estimator. The uplink channel estimator (42) attributes to a DCS problem, while the downlink channel estimator (34) is a BDCS problem. Thus the recovery algorithm SOMP is suitable for the problem, instead of SSP.

IV. COMPLEXITY
In this section, we conduct the analysis of the complexity of our proposed DS channel estimation scheme for the largescale MIMO systems. The computation complexity overhead consists of three parts, the recovery algorithm of the BDCSbased problem, the BDSO pilot design algorithm and the linear smoothing method. The complexity for the BSOMP algorithm is O(G 2 ), which has been analyzed in [29]. Then as to the proposed BDSO pilot design algorithm, the total times of multiplication for calculation of µ is GL 2 N B and the complexity is O (IterGL 2 N B ). Moreover, the complexity for the linear smoothing method is O (N KN B ). Finally, it is obtained that the complexity for our proposed channel estimation scheme is O(G 2 + IterGL 2 N B + N KN B ).

V. SIMULATION
In this section, we conduct simulations by MATLAB to verify the performance gain of our proposed BDCS based channel estimator and the proposed pilot design algorithm BDSO. The proposed channel estimator is compared with the channel estimation schemes, least squares (LS) and LCC [11]. The LS algorithm performs pseudo-inverse directly, which ignores the sparsity of the CE-BEM coefficients and the LCC scheme is CS-based, which utilizes the sparsity in the delaydoppler domain. The recovery algorithms includes SOMP, SSP and BSOMP. As to the proposed pilot design algorithm BDSO, the equidistant pilot pattern [12] and the genetic algorithm (GA) [22] are selected as the comparison schemes. The equidistant pilot pattern means that the nonzero pilots are distributed uniformly, which have a equal distant with each other and the GA algorithm exploits the genetic search   The parameters of our simulated system are given here. We consider an OFDM system with the channel model ITU-Vehicular B [31], which is also adopted in the work [30] for simulation. The number of the subcarriers N = 4096 and the number of the nonzero ones G = 192. We can see that the pilot overhead is G(2D − 1)/N = 23.4%. The channel length is L = 200 with K = 6 nonzero paths. The central frequency is 2.35GHz and the bandwidth is 20MHz. QPSK is selected  as the modulation technique.

A. The Verification of BDSO
In Fig. 5-6, we verify the effectiveness of our proposed pilot design algorithm BDSO under the number of the antennas N B =12 and the vehicular speed v=300km/h. Fig. 5 depicts the variation process of µ with the increasing of the iterations in Algorithm 1. We can see that the µ of the equidistant pilot pattern is approximately equal to 1, which is the upper bound of µ. The curve corresponding to GA declines slightly in the iteration process with 1000 times. The convergence speed is too low and the final value of µ is approximately equal to 0.6. As to our proposed BDSO, in the former 500 iterations, it presents obvious trend of declining and in the latter 500 iterations, it reaches the steady state. The final value of µ is less than 0.1.   It is known that two reference values can be utilized to measure the performance of the CS recovery algorithms. The one is the accuracy of the recovery and the other is the probability of the accurate recovery. Fig. 6 depicts the normalized mean square error (NMSE) performance of 100 repetitions for each pilot design algorithm with the recovery algorithm SSP. We can see that for the equidistant pilot pattern, the percentage of the NMSE skip points reaches almost 30% and the amplitude of the NMSE skip point reaches 5dB. Besides, the NMSE of the recovered points is just -10dB, which cannot satisfy the required accuracy of the estimation. The NMSE of the recovered points for the GA algorithm is -22dB, which is acceptable for the channel estimation. But the percentage of the NMSE skip points is 15% and their amplitude reaches -15dB. As to the proposed algorithm BDSO, the NMSE performance is -25dB without any skip points. The results of the simulation are in accordance with the conclusion that a smaller µ leads to a better recovery performance.

B. NMSE Versus SNR
In Fig. 7, we conduct the simulation to compare the performance of LS, LCC and our proposed BDCS based estimator under the number of the antennas N B = 12 and the vehicular speed v = 300km/h with the BDSO optimized pilot positions. It is found that the curves of the LS and the LCC schemes are close to 0dB and they hardly vary with the increasing of SNR. They are ineffective under the above conditions since the pilot overhead cannot support so many antennas. We adopt three recovery algorithms to solve our proposed BDCS based estimator, including SOMP, SSP and BSOMP. The SSP and the BSOMP present similar accuracy of recovery and their corresponding curves almost coincide with each other. In contrast to SOMP, the SSP and BSOMP have performance gain since they utilize both the block and the common sparsity instead of the common sparsity only. Besides, we perform the linear smoothing method combined with the three recovery algorithms and obtain the obvious performance gain. BSOMPli presents the best performance of the estimation. In Fig. 8, we verify the effectiveness of our proposed BDCS based estimator with BSOMP under N B = 16, N B = 20, and N B = 24. It is found that the performance deteriorates as the number of antennas increases since their pilot overhead is the same.

C. NMSE Versus Normalized Doppler Shift
In Fig. 9, we present the variation of the NMSE as the normalized Doppler shift increases under SN R = 30dB with the BSOMP algorithm. It is found that the curve of the LS and LCC schemes are close to 0dB since they are ineffective under the pilot overhead of 23.4%. The recovery performance deteriorates with the normalized Doppler shift increasing since the increased ICI leads to the increased modeling error. The algorithm BSOMP has performance gain compared with the SOMP and the gap is more obvious with linear smoothing method. The algorithm SSP presents the worst performance since it is sensitive to the interference.
In Fig. 10, we verify the effectiveness of our proposed BDCS based estimator with more antennas N B = 16, N B = 20, N B = 24 utilizing the algorithm BSOMP under SN R = 30dB. It is found that the performance gets worse as the number of antennas increases. All curves present the same rising tendency as the normalized Doppler shift increases. Besides, the linear smoothing method contributes to the performance gain. K   Fig. 11 depicts the variation of NMSE as the channel sparsity K increases for our proposed BDCS based estimator under SN R = 25dB exploiting the algorithm BSOMP and the linear smoothing method. We present the simulation for N B = 12, N B = 16, N B = 20 and N B = 24. Since the pilot overhead is fixed, the performance deteriorates as the sparsity increases under the same number of antennas and the performance gets worse as the number of the antennas increases under the same channel sparsity.

VI. CONCLUSION
We propose a BDCS based DS channel estimator and the corresponding pilot design algorithm BDSO in large-scale MIMO systems. It supports more antennas and guarantees the performance of the estimation with the affordable pilot overhead. In the future, we may take account of the sparsity in beam domain and improve the performance of the estimation further.