Robust Adaptive Beamforming Algorithm for Sparse Array Based on Covariance Matrix Reconstruction Technology

,


Introduction
Adaptive beamforming (ADBF) is an important branch of array signal processing and has gradually become one of the symbols of array signal processing. e essence of beamforming is to adaptively change the weighting factors of each array element according to the signal environment, thereby weighting the array elements and performing spatial filtering to achieve the purpose of enhancing the desired signal and suppressing interference. It has been widely used and developed rapidly in the fields of communications, radar, sonar, seismic prospecting, radio astronomy, and biomedical engineering [1][2][3][4]. However, theory and practice prove that beamforming is very sensitive to the problem of signal model mismatch. Insufficient number of sampling snapshots, element position errors, and different array models will cause model mismatch. As long as there is an error between the true value and the prior value, especially when the desired signal component exists in the sampling covariance matrix (SCM), the performance of the beamformer will be seriously degraded [5,6]. erefore, improving the robustness and efficiency of adaptive beamforming under different SA conditions and making it resistant to prior information is an important research direction of current beamforming. e principle of the classic Capon beamforming algorithm is to maximize the output signal to interference plus noise ratio (SINR). When the desired signal direction is known, it is the best spatial filter [7,8]. en in practical applications, we use the SCM instead of the interference plus noise covariance matrix (INCM) for signal processing [9]. In general, the number of sampling snapshots is limited and the signal component exists in the sampling covariance matrix, which results in a serious degradation of the performance of the Capon beamformer.
In order to eliminate the undesirable influence of the desired signal on the beamformer, scholars have proposed many algorithms [10][11][12][13][14][15][16][17][18][19][20]. e main idea of [10] is to add a weighted identity matrix before the inversion of the received signal SCM of the array and to improve the robustness of the beamformer by reducing the dispersion of small eigenvalues corresponding to noise. However, it is difficult to determine the optimal diagonal loading factor in different situations. In algorithm [15], the uncertainty set of the desired signal direction vector (DSDV) is set as the constraint condition, and the DSDV is estimated by maximizing the output power and constructing the optimization problem. However, performance optimization is difficult to occur in the worst case and the upper bound of the mismatch vector norm is usually unknown. In recent years, [17] proposed to replace the SCM with the reconstructed INCM, which can not only eliminate the desired signal components, but also reduce the covariance matrix mismatch caused by the small snapshot effect, so that the performance of robust beamformer can be significantly improved.
Since uniform linear array (ULA) is easy to implement and satisfies the Nyquist sampling theorem, ADBF algorithms based on ULA have been widely studied and applied. For a ULA, increasing the number of array elements to increase the array aperture and the degree of freedom (DOF) of the array will lead to hardware costs and computational complexity. As we all know, sparse arrays can obtain larger array apertures and higher DOF than traditional ULA and thus have advantages in estimating the DOF and resolution of signals. erefore, without increasing the number of array elements, we are the first choice to increase the array aperture and DOF by using a nonuniform sparse array [21][22][23][24][25][26][27]. In recent years, the research of sparse arrays has made great progress. New types of nonuniform sparse arrays, such as coprime arrays (CA), nested arrays (NA), and minimal redundant arrays (MRA), have been proposed. ese arrays have the characteristics of simple array structure and easy construction. At the same time, they are optimized. Arrays have a high degree of freedom, so array signal processing based on nonuniform sparse arrays has gradually become a hot spot. Applying the traditional RAB algorithm directly to the sparse array will result in a decrease in output performance, because the traditional RAB algorithm based on a uniform array does not consider the sparse characteristics of the sparse array elements.
Coprime array is a kind of sparse array. Its idea is to combine two sparse uniform subarrays and when the distance between the subarrays is a pair of coprime integers, that is, an integer multiple of the half-wavelength of the incident signal, an unambiguous estimation result can be obtained [28][29][30]. In [31], a compressed sensing algorithm based on the compressed sensing matrix is used to connect the virtual uniform linear array (VULA) of the coprime matrix with the physical array, avoiding the dimensional mismatch problem of calculating the weight vector directly by the VULA. In [32], in the generated VULA, the VULA spatial spectrum estimator is used to obtain the power of the signal, and the DOA estimation of the desired signal is used to directly calculate the DSDV. J. Yang et al. proposed a robust adaptive beamforming algorithm based on nested array [33] and constructed INCM by projecting the spatial smoothing matrix of virtual array into the signal subspace. In [34], a novel gridless direction finding method is proposed via the low-rank covariance matrix approximation, which is based on a reweighted nuclear norm optimization. In [35], Capon adaptive beamformer is applied to nested array by using the spatial smoothing matrix of virtual array. e above algorithms are based on a specific SA, and almost no algorithm is suitable for most SAs and has good performance. Based on this fact, this paper proposes two adaptive beamforming algorithms with strong applicability. e proposed algorithms all remove the desired signal components by reconstructing the INCM, but the reconstruction methods and conditions are different. e first method assumes that the desired signal is a non-Gaussian signal, thereby using the properties of the fourth-order cumulant to remove interference signals and noise to obtain an unbiased estimate of the DSDV. e signal power and direction vector estimated through the spatial spectrum are integrated to reconstruct the covariance matrix. e second method uses the subspace algorithm and eigenvalues to obtain the direction vector and power of the signals in turn and finally uses the estimated parameters to calculate the weight vector.
e main work and contributions are summarized as follows.
(1) is paper proposes two adaptive beamforming algorithms with strong applicability for SA. e first is that the applicable condition is that the desired signal is a non-Gaussian signal, and the second requires the SA to meet certain conditions. Simulation experiments show that both algorithms have good performance. (2) In the first proposed algorithm, we estimate DSDV through the fourth-order cumulant and reconstruct the INCM by integrating the desired signal outside the region. In the second algorithm, we derive the corresponding relationship between the eigenvalues and power under the sparse array and use the defined coefficients to obtain the one-to-one correspondence between the power and the direction vector. (3) e first algorithm still has good performance when the desired signal is correlated with an interference signal. e reason is that the algorithm directly performs INCM reconstruction, which can remove the influence of the desired signal. e second algorithm does not require complex operations such as integration, so the computational complexity of the algorithm is relatively low.
In this article, we use uppercase and lowercase bold letters to represent matrices and vectors, respectively. Given a matrix A, we use A T , A H , and A * to denote the transpose, the Hermitian transpose, and the conjugate of A, respectively. E(·) is used to express statistical expectations. Vec (.) denotes the vectorization operator which stacks the column vectors of a matrix one by one. ⊗ Stands for the Kronecker product.
We organize the rest of this article as follows: Section 2 describes the signal model. Section 3 estimates the DSDV and reconstructs the INCM. Section 4 gives the algorithm simulation experiments under different conditions; Section 5 summarizes the work of this paper.

Signal Model
Assume that M far-field narrowband signals are incident on a SA composed of N omnidirectional antennas with the number of array elements. e unit array element spacing is defined as (λ/2)(d � λ/2), where λ is the signal wavelength, and the position of the sensor can be represented by an integer set as . Taking the first element to be the reference, the output of the uniform array at time t can be expressed as where s(t) ∈ C M×1 is the vector of desired and interference signals; n(t) ∈ C N×1 is a Gaussian white noise with mean value of 0 and variance of σ 2 n , which are independent of each other. A(θ) ∈ C N×M is the array manifold matrix, a(θ i ) is the direction vector of the i th signal, and Z is the number of snapshots.
where θ i (i � 1, 2, . . . , M) is the azimuth of the i th signal. Suppose we define s 1 (t) as the desired signal and the other M − 1 as the interference signals. en the array output can be expressed as where d(t) � a(θ 1 )s 1 (t) is the desired signal vector; represents the vector of the interference signals. e sampling value of the k-th snapshot on the m-th element of the array is e output of the beamformer can be expressed as where w ∈ C N×1 is the weight vector. e output SINR of the array is used to evaluate the performance of the beamformer, which is defined as follows: where σ 2 1 is the power of the desired signal; R i+n is INCM, which is given by the following formula: where σ 2 i represents the power of the i th signal; σ 2 n is the power of the noise. I N ∈ C N×N is the identity matrix whose main diagonal is 1. Here, E(·) denotes the expectation operator.
In this paper, minimum variance distortionless response (MVDR) adaptive beamforming (MVDR) criterion is used to construct weight vector. It can be built by solving the following minimization problems: where the solution is given by and w opt is also called MVDR beamformer. It is difficult to get the INCM R i+n directly in practical application. erefore, the sampling covariance matrix R x can be used instead.
where Z denotes the number of snapshots. It should be noted that when Z is very small, the gap between R x and R i+n is relatively large, which will make the desired signal be suppressed as interference; with the improvement of signal to noise ratio (SNR), the proportion of desired signal components will become larger, which will lead to serious signal self-zeroing phenomenon, especially in the case of high SNR. erefore, it is necessary to remove the desired part when estimating the covariance matrix.

The Proposed Algorithm
In this part, we propose two adaptive beamforming algorithms under SA. From formula (9), it can be seen that the weight vector is calculated from the INCM and DSDV. erefore, the core idea of the proposed algorithm is to accurately reconstruct these two parameters, so that the performance and robustness of the algorithm are better.
In the first method, we assume that the desired signal is a non-Gaussian signal and use the fourth-order cumulant method to estimate the direction vector of the desired signal. en we use the Capon spatial spectrum and integrate outside the region where the desired signal is located to reconstruct the INCM. e second method uses the direct-MUSIC method to estimate the direction of arrival (DOA) of the signals under a specific SA and then uses the estimated DOA to reconstruct the direction vector of the signals. Subsequently, the eigenvalues are used to estimate the power of interference signals and noise and substituted into formula (7) to obtain the INCM.

Estimating DSDV through the Fourth-Order Cumulant.
e first N received data of the fourth-order cumulant of the array can be expressed as International Journal of Antennas and Propagation Since the desired signal is a non-Gaussian signal, interference and noise are both Gaussian signals, so from the nature of the fourth-order cumulant, we can get where r 4,q is the fourth-order cumulant of the desired signal and β � |a e above formula shows that C 4 is a copy of the direction vector of the desired signal, and the difference between the two is only a scalar factor of β, so C 4 can be regarded as the estimated value of the direction vector of the desired signal.

Reconstructing INCM through Integration.
From formula (7), we know that we need to estimate the power of the interference signals, noise, and the DOA of the interference signals to reconstruct the INCM. In this correspondence, we use the Capon spatial spectrum estimator to estimate the power of the signal as where u(θ) and P(θ), respectively, represent the direction vector and power of the signal with the DOA of θ. It can be seen from (7) and (15) that INCM can be reconstructed as where Θ is the angular sector where the desired signal is located, and Θ is the supplementary part of Θ. In other words, Θ ∪ Θ contains the entire spatial domain, and Θ ∩ Θ is empty. erefore, the covariance matrix R i+n1 collects all the interference and noise in Θ outside the sector. Substitute the estimated DSDV C 4 and the reconstructed INCM R i+n1 into (10) to calculate the weight vector as 3.3. Estimating DOA through Subspace Algorithm. We know that the traditional ULA requires that the distance between the array elements is not greater than half the wavelength of the incident signal. When the array element spacing is greater than half the wavelength, the array is called a SA. When traditional subspace algorithms are used for SA estimation, if the array structure is not restricted, grating lobe may appear in the result, resulting in inaccurate estimation results.
In order to solve the problem of the grating lobe, we constrained the array structure to obtain the DOA estimation of the signals using the direct-MUSIC algorithm [36]. We have proved through simulation experiments that, in a sparse array that satisfies the condition the influence of the grating lobe on the DOA estimation results can be removed, and an accurate result is finally obtained, where d i represents the unit distance of the (i + 1) th element from the reference element, and (d 1 , d 2 , . . . , d N−1 ) � 1 represents the common divisor of the other elements except the reference element is 1. e process of using the direct-MUCIC algorithm for DOA estimation is as follows.
First, perform eigenvalue decomposition on the sampling covariance matrix R x , and sort the eigenvalues from large to small, and regard the eigenvector corresponding to the largest eigenvalue equal to the number of signals M as the signal subspace and the remaining N − M eigenvalue corresponding eigenvectors as the noise subspace Among them, U S is a signal subspace, Σ S is a diagonal matrix composed of M large eigenvalues, U N is a noise subspace, and Σ N is a diagonal matrix composed of N − M small eigenvalues.
From the orthogonal relationship between the noise subspace and the signal direction vector, the spatial spectrum function of the array is obtained change the value of θ, when θ is the incident signal, a(θ) is orthogonal to the noise subspace, and Q takes the maximum value. at is, we can obtain the DOA estimation of the desired signal θ 1 and the interference signal θ i , i � 2, 3, . . . , M, through the spectral peak search. Substituting θ 1 into (3) we get the estimated DSDV We take the CA, NA, and MRA as examples to obtain the DOA of the signals. e simulation diagram is as follows.
It can be seen from Figure 1 that, under the condition of satisfying formula (19), the direct-MUSIC algorithm can remove the influence of the grating lobe and accurately estimate the DOA of the signals.

Estimating the Power through Eigenvalue Decomposition.
We can write formula (19) as It can be seen from formula (22) that the eigenvalue of the covariance matrix eigenvalue decomposition has a certain corresponding relationship with the signal power. It is available to transform formula (22): Since the direction vectors of different signals are orthogonal to each other, we have By substituting formulas (24) into (23), we can get Since the first M eigenvalues correspond to signal power, the following N − M eigenvalues correspond to noise power. erefore, the small eigenvalue of the matrix is the noise power. In practice, in order to obtain higher accuracy, we average the noise power By substituting formula (26) into (25), we can get erefore, we use the large eigenvalue to subtract the small eigenvalue corresponding to the noise power and then divide it by the number of array elements N to estimate the power of the interference signal However, there is a problem: we cannot determine the one-to-one correspondence between the direction of arrival and the power of the desired signal and interference signal. erefore, we define the correlation coefficient and we substitute the eigenvector e i into U S , and a(θ) is the direction vector composed of θ i , i � 2, 3, . . . , M, and (3). When the eigenvectors e i correspond to the direction vectors a(θ), their relations are not orthogonal to each other, so the coherence coefficient is the largest. From formula (22) we get the corresponding relationship between the eigenvector and the eigenvalue, so we get the corresponding relationship between the eigenvalue and the direction of arrival. From formula (28), we get the corresponding relationship between eigenvalue and power, and finally we get the corresponding relationship between power and direction of arrival. e estimated power and direction of arrival are substituted into (6) reconfigurable INCM as Finally, by substituting the desired signal SV a s of formula (21) and the reconstructed INCM R i+n2 of formula (30) into (9), the weight vector of the proposed algorithm can be expressed as  N 2 + 2N + 1)N). erefore, the complexity of the proposed algorithm is O((N 2 + 2N + 4M + 2)N). Table 1 shows the comparison of computational complexity of the MVDR [10], WORST-CASE [15], CMA [6], NAR [33], and proposed algorithm.

Simulation Results
In this section, we consider a CA, the element positions of the two sparse uniform arrays are located at [0 5 10 15 20 25]   In other words, the array is deployed by 10 identical omnidirectional physical array elements, and the array element spacing d is half a wavelength. Suppose there are desired signal and two uncorrelated interference signals from the assumed direction 5°and from −10°and 30°hitting the array. e interference to noise ratio (INR) in each array element is set to 30 dB. Additive noise is modeled as a complex Gaussian random process with a mean value of 0. When comparing the performance of the adaptive SNR forming algorithm in the input SNR, the number of snapshots is fixed at Z � 30. In the performance comparison between the output SINR and the number of snapshots, the SNR is set to 20 dB. For each scene, 200 Monte Carlo tests will be performed.
When comparing the output SINR performance, the proposed algorithm will be compared with the MVDR algorithm [10], the worst case performance optimization algorithm [15], the beamforming based on covariance matrix reconstruction [6], and the NA beamforming based on covariance matrix reconstruction [33]. In the MATLAB simulation diagram, these algorithms use "Proposed1", "Proposed2", "MVDR", "WORST-CASE", "CMA", and "NAR" as descriptions. Assuming that the number of signals is known as a priori condition, our simulation shows a worst case performance optimized beamformer with a parameter of 3. In the Proposed1 beamformer, we assume that the angular region of the desired signal is Θ � [0°, 10°]; therefore, Θ � [−90°, 0°t ∪ n(10°, 90°].

Common Sparse Array.
In this experiment, we considered the case where the sparse array is a CA, NA, and MRA. Simulation experiments are carried out on the influence of input SNR on output SINR under different formations. Figure 2 shows the change table curve of the output SINR with the input SNR under the CA. Figure 3 depicts the variation curve of output SINR with input SNR under NA. Figure 4 shows the change curve of output SINR with input SNR under the MRA. It can be seen from the figure that as the input SNR increases, the MVDR and WORST-CASE algorithms gradually tend to be fixed values. is is because these algorithms do not remove the desired signal in the sampling covariance matrix, which leads to excessive SNR. It is self-destructive when it is high. Although the CMA and NAR algorithms remove the desired signal by reconstructing INCM, their accuracy is lower than the proposed algorithm. Simulation experiments show that the proposed algorithm has the best performance under different formations and is close to the optimal performance curve.

Received Signal Coherence.
In this simulation experiment, we analyzed the coherence between the desired signal and the interference signal, assuming that the DOA estimation of the desired signal and interference signal has no error. Figure 5 shows the change curve of output SINR and input SNR when the desired signal 5°is coherent with the interference signal from −10°and the other signals are independent of each other. Figure 6 shows the variation curve of the number of snapshots with the output SINR in the same Step 1: estimate DSDV using fourth-order cumulants in (14) Step 2: reconstruct INCM using integration outside the desired signal area using (16) Step 3: calculate the weight vector using (17) ALGORITHM 1: e first proposed algorithm steps.
Step 1: estimate the DOA of the signals and reconstruct the direction vector by (3) Step 2: deduce the corresponding relationship between eigenvalues and power in (28) Step 3: find the corresponding relationship between power and direction vector Step 4: use (30) to reconstruct INCM and calculate the weight vector in (31) ALGORITHM 2: e second proposed algorithm steps.

Algorithms
Computational complexity e MVDR algorithm [10] O((N 2 + 2N + 1)N) e WORST-CASE algorithm [15] O((5N 2 + 2N + 1)N) e CMA algorithm [6] O(N 2 G)(G ≫ N) e NAR algorithm [33] O International Journal of Antennas and Propagation situation. It can be seen from the simulation results that when the expected signal and the interference signal are coherent, the performance of the other comparison algorithms has a large degree of decline, while the performance of the proposed algorithm has a relatively small degree of decline. Simulation results show that the proposed algorithm has a strong ability to suppress coherent signals and can always achieve near optimal performance in a wide range of input SNR.

Fixed Look Direction Mismatch.
In this simulation, we studied the situation where there is a fixed DOA error.
Assuming that there is an error of 3°, the true DOA is 5°, and our estimated DOA is 2°. Figure 7 describes the curve of the output SINR of each algorithm with the input SNR. At the same time, the curve of the output SINR with the number of snapshots is shown in Figure 8. It can be seen that as the input SNR increases, the performances of the MVDR and WORST-CASE algorithms drop significantly. e performances of the proposed algorithms and CMA algorithms all drop slightly, but the performance of the proposed algorithm is still better than other algorithms.
is proves that the proposed algorithm is robust to fixed look direction mismatch.

Random Look Direction Mismatch.
In this simulation experiment, consider the random uniform error of the desired signal steering vector. Specifically, the arrival angle    of the desired signal will have a uniform random error of [−4°, 4°], which means that the actual arrival angle of the desired signal is now uniformly distributed within [1°, 9°]. Figure 9 depicts the output SINR versus input SNR curve, and Figure 10 shows the output SINR versus the number of snapshots. It can be seen from the figure that although the performance of each algorithm has declined due to the random error of the desired signal steering vector, the performance of the proposed algorithm is still better than other comparison algorithms.

Random SV Mismatch.
In this simulation experiment, we consider the impact of an uncertain interference in the desired signal on the performance of the beamformer, that is, the performance of the algorithm when the random DSDV mismatch. In this case, the real DSDV can be expressed as where a(θ 1 ) is a hypothetical DSDV with direction of θ 1 , and δ is a random DSDV mismatch generated by the uncertainty set, as follows: where α is the norm of δ, which is generated randomly from one run to another in [0, �� � 0.5 √ ]. ϕ u (u � 1, 2, . . . , N) represents the coordinate of δ independently generated from [0, 2π] in the u th run. Figures 11 and 12 show the relationship between the input SNR, the number of snapshots, and the output SINR. Obviously, as the input SNR increases, the output SINR of the proposed algorithm is always close to the ideal value, and convergence requires less snapshots. e simulation results prove that the algorithm is effective in the case of random SV mismatch. In addition, we can see that the performance of the beamformer is significantly better than other test beamformers, and we conclude that the beamformer is more robust to random DSDV mismatches.

Conclusion
Based on the sparse array, this paper proposes two adaptive beamforming algorithms suitable for different conditions. e first algorithm uses the fourth-order cumulant to remove the interference signals and noise to obtain DSDV. Estimate the power of the signals through the spatial spectrum, and substitute it into the formula to integrate and reconstruct INCM. Considering that not all desired signal is non-Gaussian signals, we propose the second method. As long as the array meets certain conditions, we can estimate the DOA of the signal through direct-MUSIC to reconstruct the direction vector of the signals, estimate the power of the interference signals and noise through the feature value, and finally calculate the weight vector using the estimated parameters. Simulation experiments show that the algorithm has good performance and strong robustness to SA.
Data Availability e data, which are produced by simulations, used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that they have no conflicts of interest.    International Journal of Antennas and Propagation 9