Sparse Array Design for DOA Estimation of Non-Gaussian Signals: From Global Postage-Stamp Problem Perspective

Key Laboratory of Dynamic Cognitive System of Electromagnetic Spectrum Space Nanjing University of Aeronautics and Astronautics, Ministry of Industry and Information Technology, Nanjing 211106, China College of Electronic and Information Engineering, Nanjing University of Aeronautics and Astronautics, Nanjing 211106, China The 28th Research Institute of China Electronics Technology Group Corporation, Nanjing 210007, China


Introduction
In array signal processing, direction of arrival (DOA) estimation has drawn considerable attention, which has been widely applied in various fields, such as communication, radar, and navigation [1,2]. To obtain improved DOA estimation performance, the sparse arrays with reduced mutual coupling and enhanced degree of freedom (DOF) have been proposed, whose adjacent sensor spacing is no longer limited to half wavelength, as compared with uniform linear array (ULA). Minimum redundancy array (MRA) [3], minimum holes array (MHA) [4], coprime array (CA) [5], and nested array (NA) [6] have been proposed, which are designed based on the second-order cumulant (SOC) for Gaussian signals. Besides, these arrays utilize second order difference coarray (2-DC) or second order difference coarray of sum coarray (2-DCSC) generated by vectorization operation to obtain virtual sensors, which can be used to construct virtual nonuniform linear array (NLA) or ULA. Once the equivalent received signal of the virtual array is constructed, the DOA estimates can be obtained by compressed sensing (CS) [7], discrete Fourier transform (DFT) algorithm [8], and other algorithms. Nevertheless, there still exist some drawbacks for these sparse arrays. Specifically, there is no closed-form expression of physical sensor positions for the MRA and MHA. The CA structures generate holes in the difference coarray which decreases the consecutive DOF, and the NA structures suffer from severe mutual coupling due to the dense part. In order to tackle these issues, some sparse arrays with enhanced DOA estimation performance based on NA or CA have been proposed, such as augmented NA (ANA) [9], the generalized CA with displaced subarrays (CADiS) [10], unfolded coprime linear array (UCLA) [11,12], and augmented CA (ACA) [13].
It is noteworthy that all aforementioned sparse arrays mentioned above are all designed based on SOC, which is available for Gaussian signals. In practice, most sources are non-Gaussian, and fourth-order cumulant (FOC) matrix is exploited instead of SOC to calculate correlation matrix [14,15]. In particular, the FOC-based methods employ second order difference coarray of sum coarray (2-DCSC) to obtain virtual sensors with vectorization operation [16,17]. As the 2-DC and second order sum coarray (2-SC) are utilized successively, it is complicated to obtain the location relationship between physical sensors and 2-DCSC. In [17], Fu et al. employed virtual nested multiple-input multipleoutput (MIMO) array [18], which utilized fourth order difference coarray of sum coarray (4-DCSC), termed as NA-MIMO-DCSC. The authors optimized the process of obtaining 4-DCSC separately, i.e., fourth order sum coarray (4-SC) and 2-DC, which further improved DOF greatly when FOC was employed in MIMO systems. Unfortunately, limited to the properties of NA, the virtual sensors still suffered from a lot of redundancy in the steps of optimizing 2-DC. More importantly, the design principle of NA-MIMO-DCSC neglected the discussion of cross-sum coarray generated in the process of obtaining 4-SC.
In this paper, we propose a sparse array design method for non-Gaussian signal from the perspective of global postage-stamp problem (GPSP) [19,20]. To further illustrate this method, the two steps of obtaining sparse array geometry based on the structure of UCLA are described. Specifically, the first step is to get the longest possible two 2-SC from given sensor number of each subarray, which can be modeled as a GPSP. The solution to GPSP can also be used to obtain initial physical sensor location sets with half-wavelength interelement spacing. Subsequently, the second step is to adjust the physical sensor location sets by multiplying the appropriate coefficients, which can be determined based on the structure of UCLA. To verify the priority of the proposed array geometry, the properties of UCLA about DOF and consecutive virtual sensors are also provided. Furthermore, the expression between the DOF of proposed geometry and physical sensors is presented, and the relationship between 2-SC and cross-sum coarray is also analyzed in Section 4.1.
In particular, the contributions of this paper are summarized as follows: (1) We investigate the array design problem of sparse array with non-Gaussian signals and tackle this problem by obtaining the longest possible virtual 2-SC and exploiting the properties of UCLA (2) We formulate the obtaining longest possible 2-SC problem as a GPSP and give a simple solution to the GPSP, which can be used to determine the location sets. Furthermore, we provide the closed-form expression of maximum DOF and consecutive virtual sensors of UCLA with 2-DC, which is used to analyze the DOF of proposed array geometry (3) We devise a sparse design principle which divides the process of obtaining sparse array location into two steps, i.e., obtaining the 2-SC from the given number of subarray and determining the location of physical sensors with multiplication operation We outline this paper as follows. In Section 2, we present the sparse array model and mutual coupling effect; besides, the properties of UCLA are also investigated. Section 3 elaborates the GPSP and proposed array geometry with enhanced DOF. Section 4 gives the performance analysis. Section 5 provides simulation results, and the conclusions are drawn in Section 6.
Notations. Vectors and matrices are represented by utilizing bold low-case and bold-case characters, respectively. ⊗ is Kronecker product, and ⊙ denotes Khatri-Rao product. ð⋅Þ T , ð⋅Þ H , ð⋅Þ −1 , and ð⋅Þ * stand for the transpose, conjugate transpose, inverse, and complex conjugation of a vector or matrix, respectively. vecð⋅Þ denotes the vectorization operation. diag ð⋅Þ is the diagonal matrix operator, and k⋅k F means the Frobenius norm.

Preliminaries
In this section, we first present the definitions of 2-DC, 2-SC, and 2-DCSC. Subsequently, considering that the 2-SC of proposed array geometry is constructed based on the UCLA, we exploit the properties of UCLA in terms of DOF and consecutive lag location. Finally, we provide the sparse array model without mutual coupling and with mutual coupling, respectively.
2.1. 2-DC, 2-SC, and 2-DCSC. For a linear array with L physical sensors, whose location set can be denoted as S = fd 1 ,⋯d i ,⋯,d L g, where d i , 1 ≤ i ≤ L represents the location of i-th sensor.
Definition 1 (2-DC). The second order difference coarray location D can be defined as [6] where D + and D − represent the positive and negative elements of D, respectively.
Definition 2 (2-SC). The second order sum coarray location S s can be defined as [17] Definition 3 (2-DCSC). The second order difference coarray of sum coarray location D dcsc can be defined as 2 Wireless Communications and Mobile Computing where S and D represent the physical sensors and 2-DC location sets, respectively.
2.2. The Properties of UCLA. Figure 1 shows the configuration of UCLA, which is composed of two sparse uniform subarrays that overlap at the origin. The sensor position of UCLA can be represented as [11,12] where M and N represent the physical sensors of subarray 1 and subarray 2, respectively. d = λ/2, and λ is the wavelength. According to Definition 1, the 2-DC of UCLA can be denoted as D UCLA , and D UCLA = D + UCLA ∪ D -UCLA . In addition, the properties of UCLA are derived in the Lemma 4. As described in Lemma 4, the DOF of UCLA can be given by 2.3. Sparse Array Model without Mutual Coupling. Assume that there are K far-field narrowband uncorrelated non-Gaussian signals impinging on a nonuniform linear array with DOAs ðθ k , k = 1, 2, ⋯, KÞ, where θ k denotes the elevation angle of the k -th target. The location set of nonuniform linear array with L physical sensors can be denoted as S = fd 1 , d 2 ,⋯,d L g; then, the received signal xðtÞ can be expressed as [21].
To obtain the equivalent received signal from the virtual array, we vectorize C 4,x as [17,22] where whose elements are constructed by a 4,x ðθ k Þ with the specific form e −jπðd i −d j Þ sin θ k , d i , d j ∈ D, and p = ½σ 2 1 , σ 2 2 ,⋯,σ 2 K T , σ 2 k ðk = 1, 2,⋯,KÞ denotes the power of k-th source. 2.4. Sparse Array Model with Mutual Coupling. In practice, each sensor will inevitably be affected by the radiation of from its adjacent sensors, and consequently, the received signal should be modified in the presence of mutual coupling. That means the mutual coupling effect must be taken into consideration. Specifically, the signal model in (6) can be refined asx where C represents the mutual coupling matrix. By utilizing the B-band model in [23], the elements in the mutual coupling matrix can be represented by where d i , d j ∈ S and c 1 = 0:3e jπ/3 , c l = c 1 e −jðl−1Þπ/8 /l, l ∈ ½2, B, B = 100 represents the maximum spacing between the  [23]. In order to measure the strength of the mutual coupling effect, the coupling leakage is introduced and is defined as Based on the mutual coupling model in (12), the equivalent signal model in (9) can be reconstructed as where C vec = ðC ⊗ C * Þ * ⊗ ðC ⊗ C * Þ.

Sparse Array Design Principle
Motivated by the advantages of the sparse arrays over the uniform arrays, several sparse array configurations as well as design methods have been proposed. Unfortunately, most existing array design methods generate the difference coarray based on the SOC while the FOC with great potentiality in improving array performance is neglected. In this part, we devise a sparse array principle based on 2-DCSC employed by vectorized FOC-based method.
In fact, because 2-DCSC performs twice 2-DC operations on physical sensors, which can also be converted into 2-SC and 2-DC operations, successively, as described in (3). The relationship between virtual sensors and physical sensors is complicated; it is difficult to obtain the closed-form expression of physical location with large DOF. We optimize the process of obtaining the location of physical sensor location into two steps, i.e., obtaining longest possible 2-SC of each subarray, constructing physical sensor location set based on the structure of UCLA. Specifically, the first step can be modeled as GPSP, whose solution can be used to determine the physical sensor set lacking interelement spacing coefficients. The second step employs the structure of UCLA to adjust the self-sum coarray set obtained by the first step, which can also be used to determine the spacing of physical sensors. As a result, the physical sensor location set can be obtained. According to (3), the 2-DCSC can be calculated from physical sensor set directly.

Global Postage-Stamp Problem (GPSP).
The GPSP is a well-known combination problem that has not yet been fully solved, which can be described as follows: for given positive integers h and k, a set with k nonnegative integers is determined by Remark 5. The elements in S k need to be summed h times to realize the consecutive numbers 0, 1, 2, ⋯, n h ðS k Þ. Besides, the value of n h ðS k Þ should be as large as possible.
The solution to the GPSP has been investigated in [19,20]. Herein, we also give a simple method to solve GPSP. Specifically, h can be considered as the maximum number of stamps on an envelope; S k is the set of stamp denominations with k non-negative value. Then, the computational complexity of obtaining all possible combinations is Oðk h Þ. An example is given to illustrate the solution method to GPSP explicitly. When S k = f1, 3g, in this case, 1 and 2 stamps are required to obtain postages 1 and 2, respectively. For postage 3, there are two combinations, i.e., one stamp with postage 3 or three stamps with postage 1. According to Remark 5, because of the existence of postage 0, the number of stamps currently required should be as small as possible. Besides, it can be concluded that the minimum k can be obtained by comparing all the postage combinations. The process of this method is summarized in Algorithm 1, where h = 2. It is noteworthy that S k is treated as prior information to obtain n h ðS k Þ. When S k is unknown, its (i + 1)-th element should satisfy where a i represents i-th element of S k , 1 ≤ i ≤ k. Then, all elements of can be obtained by enumerating all the values of a i . The result of Algorithm 1 is listed in Table 1, which is also provided in [20].

Wireless Communications and Mobile Computing
Suppose that the physical sensors are composed of multiple linear subarrays, and the number of elements in each subarray has been given. Without loss of generality, assume that the physical sensors are composed of two subarrays with M and N elements, respectively. As illustrated in Section 3, to obtain the longest possible 2-SC of each subarray, the integer sets S M and S N are solved from GPSP in (15) with k = M, N, h = 2, which can be used to determine the largest possible lags n 2 ðS M Þ and n 2 ðS N Þ. It is noteworthy that the sets S M and S N can be considered as the initial physical location sets without coefficients. Based on Remark 5, if all elements of S M are negative, the consecutive integer set obtained from M is where S M represents a set with M elements and n 2 ðS M Þ is the largest possible integer mentioned in Remark 5. Similarly, the consecutive positive integer set obtained from N is Based on the structure of UCLA, we multiply all the elements in (17) and (18) by certain coefficients d 1 and d 2 , respectively; then, the new sets can be constructed as the two cosubarrays of UCLA. Specifically, the cosubarray 1 of UCLA location set associated with (17) can be expressed as where d 1 = ðn 2 ðS N Þ + 1Þd is the interelement spacing of cosubarray 1; d = λ/2, and λ is wavelength. The cosubarray 2 location set associated with (18) can be expressed as where d 2 = ðn 2 ðS M Þ + 1Þd is the interelement spacing of cosubarray 2. Note that d 1 and d 2 should be set as coprime integers. Finally, we calculate the location set of physical sensors according to the correspondence between the physical array S and coarrays S sc . The coarrays S sc = S sM ∪ S sN determined by (19) and (20) can be considered as the 2-SC result of each subarray or the self-sum coarray of S. Considering that the 2-SC of S M and S N have been multiplied by d 1 or d 2 , then the physical sensor location S can be determined by multiplication operation, Based on the discussion described above, we give an example of the proposed geometry shown in Figure 2 to verify. For given element number in two subarrays with M = 5, N = 4, firstly, the initial consecutive integer sets S M = f−6,−5,−3,−1, 0g, S N = f0, 1, 3, 4g, and the largest possible values n 2 ðS M Þ = 12, n 2 ðS N Þ = 8 can be obtained by the solution to GPSP, which have been listed in Table 1. Subsequently, the cosubarrays spacing d 1 = 9d, d 2 = 13d can be obtained by the specifical structure of UCLA. The two 2-SC cosubarray locations can be expressed as S sM = 9d × f−12,⋯,−1, 0g, S sN = 13d × f0, 1,⋯,8g. Consequently, the physical sensor location set is represented by S = 9d × f−6,−5,−3,−1, 0g ∪ 13d × f0, 1, 3, 4g.

The Design Procedure of the Proposed Array Structure.
We summarize the detailed steps for implementation of obtaining the proposed array geometry as follows: (1) Based on Table 1

Performance Analysis
In this part, we evaluate the performance of the proposed array geometry compared to other arrays from the viewpoints of achievable DOF, mutual coupling, redundancy ratio, and Cramer-Rao Bound.

Achievable DOF.
We analyze the achievable DOF of the proposed array geometry. Considering that the S sM and S sN can be associated with the cosubarrays of UCLA S sc , then according to Lemma 4, the DOF for the 2-DC of S sc is 2ðn 2 ðS N Þ + 1Þðn 2 ðS M Þ + 1Þ − 1. Based on Definition 1,   Figure 3. Besides, the DOF of NA-MIMO-DCSC has been given in [17], which can be expressed as where DOF NA-MIMO-DCSC represents the DOF of array NA-MIMO-DCSC. As for the proposed array geometry, we can obtain the following relationship where DOF proposed represents the DOF of the proposed array geometry.
As expounded in the appendix, the 2-DC can be decomposed into self-difference coarray and cross-difference coarray; similarly, the 2-SC can also be decomposed into self-sum coarray and cross-sum coarray. For the proposed geometry, the processes of obtaining two 2-SC of the two subarrays are separated, where only self-sum coarray is utilized while the cross-sum coarray is neglected, i.e., all elements in S sc are included in the 2-SC of S. Thus, the 2-DCSC of S contains all elements in the 2-DC of S sc . Consequently, equation (23) is established. The nonnegative part of 2-DCSC of the proposed geometry is also given in the bottom part of Figure 3 to verify analysis mentioned above.
However, the DOF of proposed geometry only has lower bound without closed-form expression, because the location of cross-sum coarray of the 2-SC is not available. Interestingly, for NA-MIMO-DCSC, the 2-DCSC and the 2-DC of S sc exactly coincide.
According to (22) and (23), we have It can be concluded that the proposed array geometry can offer larger DOF than NA-MIMO-DCSC, which is also verified in Figure 4 below. Definition 6 (redundancy ratio). For a physical array with total number of virtual sensors T obtained from difference coarray and the achievable DOF is given by DOF p , we define the redundancy ratios as For a linear array with L physical sensors, in general, the number of virtual sensors is T = L 2 when the correlation matrix is calculated by SOC-based method, and when FOCbased method is exploited, T = L 4 . Redundancy ratio can be considered as an important indicator which reflects the efficiency of virtual sensors.

Cramer-Rao
Bound. The Cramer-Rao Bound (CRB) can be obtained from the inverse of the Fisher information matrix (FIM). According to [24,25], we give the FIM for the signal model in (9) via vectorization-based methods as where J denotes the number of snapshots, z = vecðC 4,x Þ and C 4,x can be obtained by (7). β = ½θ 1 , ⋯, θ K , σ 2 1 ,⋯,σ 2 K T , σ 2 1 , ⋯, σ 2 K , represent the power of non-Gaussian signals.   Consequently, the CRB of the DOA estimation with the proposed geometry via vectorization operation can be expressed as

Simulation Results
In this section, 500 Monte-Carlo simulations via CS algorithm are employed to validate the superiority of the proposed array geometry in the presence of mutual coupling. Define the root mean square error (RMSE) as where θ k denotes the true elevation of the k-th target and b θ k,i is the estimated value of θ k in the i-th ði = 1, ⋯, 500Þ Monte Carlo simulation. K is the total number of the incident signals which are uncorrelated far-field non-Gaussian. Besides, the searching interval of CS algorithm is Δ = 0:01 ∘ , and the coefficients of mutual coupling are set as shown in Section 2.4.

RMSE Performance Comparison of Different Sparse
Arrays versus SNR. In this part, we compare the RMSE performance of the different sparse arrays involved in Figure 4, including the CA [5], NA [6], ACA [13], unfolded ACA, UCLA [11], ANAI-1 [9], ANAI-2 [9], MRA [3], proposed geometry, and NA-MIMO-DCSC [17]. Besides, the location sets of physical sensors in the above arrays with 8 sensors are also given in Figure 4. Figures 5 and 6 depict the DOA estimation performance and CRB performance comparison of different arrays versus SNR, respectively, where J = 1300 and θ = ½−40 ∘ ,−30 ∘ ,−20 ∘ ,−10 ∘ , 0 ∘ , 10 ∘ , 20 ∘ , 30 ∘ , 40 ∘ . As shown in Figures 5 and 6, it is observed that the DOA estimation performance of all arrays improves with the increase of SNR, and in particular, the proposed array geometry outperforms the other sparse arrays. The ACA and CA achieve the worst RMSE performance due to the reduced DOF, and the dense part of the NA results in significantly higher mutual coupling. The UCLA enjoys the maximum DOF and the minimum coupling leakage among the arrays design based on SOC, leading to better performance. The self-sum coarray of proposed array geometry is devised to the UCLA, which takes full advantage of the UCLA. According to Figure 4, the coupling leakage of the proposed geometry is much lower than the others, and the DOF is also significantly higher than other arrays, especially those based on SOC design. Due to the optimization of DOF and mutual coupling, the proposed array geometry achieves the best RMSE performance which is superior to CA, NA, ACA, UACA, UCLA, ANAI-1, ANAI-2, MRA, and NA-MIMO-DCSC.

RMSE Performance Comparison of Different Sparse
Arrays versus Snapshots. Figures 7 and 8

Wireless Communications and Mobile Computing
It is also shown clearly in Figures 7 and 8 that the proposed geometry achieves the best DOA estimation performance among all the compared sparse arrays. The 2-DCSC of CA, NA, ACA, UACA, UCLA, ANAI-1, ANAI-2, and MRA are obtained from extending SOC to FOC directly, which perform twice 2-DC on arrays manifolds successively. As a result, lots of redundancy is generated, which means a reduction in DOF, but NA-MIMO-DCSC and the proposed geometry are all based on the solution to GPSP, which reduces redundancy ratio greatly, resulting in the superior DOA estimation and CRB performance to the other sparse arrays.

Conclusion
In this paper, a sparse array design method for non-Gaussian signals DOA estimation is presented from the perspective of GPSP, and a sparse array geometry based on UCLA is given to illustrate. For vectorized FOC-based DOA estimation algorithms, we divide the process of obtaining physical sensor location sets into two steps. Specifically, the first step of obtaining the longest possible 2-SC is modeled as a GPSP, whose solutions can also be employed to determine the initial physical location sets without spacing coefficients. Then, the second step is to multiply the physical sensor sets by appropriate coefficients based on the structure of UCLA. Besides, the properties of UCLA in terms of DOF and the location of consecutive lags are exploited. Numerical simulations corroborate the superiority of the proposed geometry in terms of DOF, mutual coupling, CRB, and DOA estimation performance. Considering the relationship between 2-DCSC and physical sensors is quite complicated, the design principle with splitting 2-DCSC into 2-SC and 2-DC will contribute to our future work effectively.

Data Availability
The data used to support the findings of this study are available from the corresponding author upon reasonable request.

Conflicts of Interest
The authors declare that they have no conflicts of interest.