Covariance Matrix Reconstruction via a Subspace Method and Spatial Spectral Estimation for Robust Adaptive Beamforming

Several robust adaptive beamforming (RAB) algorithms based on interference plus noise covariance matrix (INCM) reconstruction have been recently proposed, which can enhance the robustness of beamforming algorithms when certain mismatches occur in the model. However, some approaches ignore the resolution of the Capon spectral estimator (CSE), leading to reconstruction errors. This paper proposes a novel RAB algorithm formulated using the subspace projection method and spatial spectral estimation (SSE), which is named INCM-SSE. First, without using the CSE, the subspace projection matrix (SPM) is obtained through the integral of the angular sector where the signal of interest (SOI) is located. Subsequently, after estimating the direction of arrival (DOA) of incident signals using the multiple signal classi ﬁ cation (MUSIC) algorithm, we project the sample covariance matrix (SCM) onto the SPM to eliminate the SOI in ﬂ uence. Then, the estimation method of interference powers is derived. Moreover, the INCM is reconstructed based on the estimated powers and steering vector (SV) of interferences. The SV of the SOI is optimized by solving a quadratic convex optimization problem. The INCM-SSE algorithm not only employs SSE to improve the angular resolution but also reduces the in ﬂ uence of the SOI component by using SPM. Simulation results indicate that the proposed method is robust against various types of mismatches, thus achieving superior overall performance.


Introduction
Steering array antennas forming beam patterns according to certain optimal criteria are called smart antennas; they are also regarded as adaptive array antennas [1][2][3]. This description means that a computer can control the performance of the antennas, greatly improving the performance of the array system [4]. Adaptive array processing technology can adjust the weighting vector of array antennas in real time according to the signal environment. It employs adaptive beamforming algorithms to achieve a certain gain in the signal of interest (SOI) direction and suppress interferences in their directions [5][6][7]. The technology has been widely applied to the fields of radar, sonar, microphone array speech processing, medical imaging, internet of things, wireless communication, astronomy, and seismology [8][9][10][11].
However, adaptive beamforming algorithms are extremely sensitive to signal source mismatches and array geometry errors (e.g., signal direction error, sensor position perturba-tion, amplitude, and phase error). In particular, when components of the SOI exist in the snapshots, the beamformer performance is severely degraded [12]. Therefore, improving the robustness of beamforming is necessary for adaptive array processing. Over the decades, numerous robust adaptive beamforming (RAB) algorithms (e.g., diagonal loading method, eigen-subspace approach, uncertainty set constraint technology, and covariance matrix taper algorithm [13][14][15]) have been proposed. The foregoing methods are effective for resolving signal direction error and sensor position perturbation. Moreover, they are capable of improving the signal-tointerference-plus-noise ratio (SINR) under certain mismatch conditions. However, the SOI consistently exists in the sample covariance matrix (SCM) of these algorithms, limiting their performance at a high signal-to-noise ratio (SNR).
Over the past few years, a number of RAB algorithms based on interference plus noise covariance matrix (INCM) reconstruction (in which the SOI component is effectively removed) have been proposed. The INCM reconstruction, first reported in [16], was achieved by combining the Capon spectral estimator (CSE) and integral of the region separated from the angular sector of the SOI. Subsequently, the steering vector (SV) of the SOI is modified by solving a quadratically constrained quadratic programming (QCQP) problem. The work in [17] describes an efficient RAB method based on the spatial power spectrum sampling (SPSS) approach. The computational complexity of the INCM-SPSS method is lower than that reported in [16]. However, this technique requires numerous sensors to achieve a similar performance as that in [16]. Reference [18] proposes a subspace-based RAB algorithm via residual noise elimination and interference power estimation. Although its performance is superior to other algorithms, the proposed algorithm has numerous integrals, matrix multiplications, and matrix inversion operations that increase the running time. The INCM reconstruction achieved through the interference SV lying at the intersection of two subspaces is presented in [19]. However, the foregoing requires multiple eigen-decomposition operations, increasing the complexity of the algorithm. Lowcomplexity spatial sampling processing and virtually received SV-based algorithm are proposed in [20]. The foregoing is realized by a projection matrix orthogonal to the signal subspace and retains interference plus noise in a higher dimension. This method can greatly improve the efficiency of the INCM reconstruction and beamforming robustness. The work in [21] describes a RAB method based on the reduction of interference matrix and the SV of SOI optimization; however, its performance degrades as the number of interferences increases.
Nonetheless, some of the approaches mentioned here (such as those reported in [18,19,22]) ignore the resolution of the CSE. When the direction of arrival (DOA) of the interference is close, the CSE is unable to effectively identify it, leading to reconstruction errors. Although adjacent interferences can be solved by null broadening technology [23], this algorithm type may lead to an increase in the sidelobe level or broadening of the mainlobe of the beam pattern. Hence, its robustness (particularly when various signal model mismatches occur) remains to be verified. In this paper, a novel RAB method via the subspace projection method and spatial spectral estimation (SSE) is proposed. Different from previous RAB methods, the proposed algorithm obtains the subspace projection matrix (SPM) through the integral of the angular sector where the SOI is located and employs a high-resolution SSE algorithm to estimate the DOA of signals. In addition, the SV of the SOI is estimated by solving a quadratic convex optimization problem. The remaining portions of this paper are organized as follows. In Section 2, the signal model and problem background are briefly described. Section 3 introduces the details of the proposed method. The simulation and analysis results are presented in Section 4. Finally, the conclusions of the study are summarized in the last section.

Signal Model and Problem Background
Consider a uniform linear array (ULA) composed of M omnidirectional sensors, as shown in Figure 1. Suppose that L + 1 narrowband uncorrelated sources from a far field impinge on the ULA with directions θ 0 , θ 1 , ⋯, θ L . The M × 1 array observation vector at instant k has the following general form: where x s ðkÞ = a 0 ðθ 0 Þs 0 ðkÞ, x i ðkÞ = ∑ L l=1 a l ðθ l Þs l ðkÞ; nðkÞ denote the components of the SOI, interferences, and noise, respectively; s l ðkÞ, l = 0, 1, ⋯, L, denotes the signal waveform; nðkÞ represents the additive white Gaussian noise with zero mean and variance σ 2 n ; and a 0 ðθ 0 Þ and a l ðθ l Þ are the corresponding SVs of the SOI and interferences, respectively. These can be modeled as a l θ l ð Þ = 1, e j2πd sin θ l /λ ,⋯,e j2π M−1 where d = λ/2 represents the interelement spacing, λ denotes the wavelength of incident signals, and ð⋅Þ T is the transposed operator.
The output of the beamformer at time k can be written as where w is the excitation weight vector of M sensors and ð⋅Þ H denotes the Hermitian transpose. The complex vector, w, can be solved by maximizing the output SINR of the ULA: where σ 2 0 = Efs 0 ðkÞs 0 * ðkÞg denotes the power of the SOI, ð⋅Þ * is the conjugate operator, and Eð⋅Þ represents the expectation operator. The M × M theoretical INCM is R i+n , which can be expressed as ...

Wireless Communications and Mobile Computing
where σ 2 l and I denote the power of the lth interference and M × M identity matrix, respectively. The maximization problem of (4) can be formulated as which is known as the Capon beamformer. The solution to (6) is given by The optimal SINR can be obtained by substituting (7) into (4): The beam pattern of the Capon beamformer is given as follows: where θ represents the angular sector, aðθÞ is the SV of θ, and lg ð⋅Þ denotes the logarithmic operator at base 10. However, a 0 ðθ 0 Þ may be inaccurate, and the exact covariance matrix, R i+n , is unavailable in practical applications. The two are usually replaced by the nominal SV (â 0 ðθ 0 Þ) and SCM, given byR where K is the number of snapshots received. Correspondingly, the complex weight vector can be obtained by the sample matrix inversion (SMI) beamformer expressed as The beam pattern and output power of the SMI beamformer are given by respectively.

Proposed Algorithm
3.1. SPM Construction. References [18,24] present the estimate of residual noise power using the integral of the region where noise alone is present. This noise power can be removed from the Capon spectrum to derive the true powers of sources. The accurate covariance matrix of the SOI is reconstructed by integrating the region where the SOI is located: where b σ 2 n and Θ s are the estimated noise power and SOI region, respectively;â s ðθ s Þ represents the SV of θ s ; and S is the corresponding number of sample points. However, the above approach severely complicates the computations of the algorithm. By simply employing (15),R s can be obtained: Subsequently, the SPM is constructed through the eigendecomposition ofR s :R where ε m arranged in descending order denote the eigenvalues ofR s , v m represent the eigenvectors corresponding to ε m , U is the SPM, and V is composed of eigenvectors corresponding to B largest eigenvalues ofR s . Evidently,â s ðθ s Þ is orthogonal to the subspace of U, i.e., Z norm = kUâ s ðθ s Þk 2 2 = 0 , where k⋅k 2 represents the Euclidean norm. The comparability of the SPM obtained by the integral using and without using the CSE is shown in the following simulation such that a series of complex operations caused by the CSE can be avoided.
As shown in Figures 2 and 3, the values of Z norm obtained by both methods are similar, and when B = 3, the results of Θ s are virtually zero. Therefore, matrix U obtained through (15) can also be adopted as the SPM to eliminate the SOI in the SCM.

SSE.
Most of the previous algorithms employ the CSE to estimate the power or DOA of sources for reconstructing the INCM [25]. However, the angular resolution of the CSE is unsatisfactory; consequently, it cannot distinguish the DOA when the incident angles of two signals are adjacent. Therefore, different from previous approaches, a high-resolution SSE algorithm is used here. In recent years, numerous SSE algorithms have been proposed, including the multiple signal classification (MUSIC) algorithm, estimation of signal parameters via rotational invariant technique (ESPRIT) algorithm and their variants, principle of maximum entropy power spectrum, sparse Bayesian learning, support vector machine, and discrete Fourier transform-based methods [26][27][28][29][30]. This subsection presents the MUSIC algorithm, which is classical and facile to implement. It is utilized to estimate the spatial spectrum, which can be formulated as follows: where G N represents the noise subspace composed of eigenvectors corresponding to several smaller eigenvalues ofR x . To demonstrate that the resolution of the MUSIC algorithm is higher than that of the CSE, a simulation is presented in this section. Assume that a desired 25 dB signal arrives at the ULA from −5°; three interference signals arrive at −30°, 24°, and 28°with 25 dB, K = 50, and the number of sources is supposed to be known. The spatial spectra estimated by the MUSIC algorithm and CSE are shown in Figure 4; the curves are the average of 100 Monte Carlo [31] experiments of the two algorithms. Figure 4 indicates that when the incident angles of two signals are close, the CSE cannot distinguish them. This means that the DOA estimated by the CSE inevitably leads to a mismatch in the INCM reconstruction and a reduction in the output SINR. However, the resolution of the MUSIC algorithm is higher than that of the CSE; hence, the algorithm can identify two adjacent signals to reconstruct the INCM more accurately.

INCM Reconstruction and SV Estimation
. MatrixR x can also be expressed aŝ By combining the foregoing with kU H a 0 ðθ 0 Þk 2 2 = 0, the SOI component can be eliminated through the SPM, U. Then, the following formula is derived: whereR i represents the matrix that preserves interferences and b σ 2 n is equal to the minimum eigenvalue ofR x , which can be obtained through the MUSIC algorithm without introducing further computational complexity. Functions ðU H Þ −1 and U −1 can reduce the influence of SPM on the interferences after the SOI is eliminated. By ignoring the influence of SPM on (20),R i can be approximately expressed in the following form: where A is the M × L matrix composed of a l ðθ l Þ, l = 1, 2, ⋯ , L; each diagonal element in Λ denotes the interference power. Accordingly, Equation (22)

Wireless Communications and Mobile Computing
where b Λ is a diagonal matrix that represents the estimated powers andÂ denotes the matrix composed of the estimated SV of interferences:Â In (23),â l represents the lth interference SV; it can be optimized by where δ is the bound value; e θ l , l = 1, 2, ⋯, L, represents the multiple spectral peaks obtained by searching in space using (18); andã l ð e θ l Þ denotes the SV corresponding to e θ l . Equation (24) can be solved using the Lagrange multiplier method:â where μ satisfies the following constraint: Different from the technique presented in [18], we used (22) to calculate the interference powers. The use of this equation avoids the matrix inversion operation to maintain low computational complexity. In the case of adjacent interferences, the estimated power is observed to be greater than that of two incident angles far apart, as demonstrated by the simulations reported in [22]. For RAB algorithms, the overestimated power facilitates the formation of a deeper null, which can better suppress the interference and improve the output SINR of the array. In this case, the overestimated power does not require correction. Accordingly, the INCM can be reconstructed based on (22) and (23): Next, the SV of the SOI is optimized. The convergence of nominal SVâ 0 ðθ 0 Þ to the SV of interferences can be prohibited byâ H 0 ðθ 0 ÞUU Hâ 0 ðθ 0 Þ = 0. However, the above equality constraint is nonconvex; thus, the following constraint is set: whereâ 0 ðθ 0 Þ + e ⊥ denotes the actual SV and e ⊥ (orthogonal toâ 0 ðθ 0 Þ) represents the error between the actual SV and a 0 ðθ 0 Þ. The right-hand side of (28) is a value tending to zero, ensuring that the termâ 0 ðθ 0 Þ + e ⊥ is restricted to the SOI region. The solution of e ⊥ can be derived as follows: The foregoing is a QCQP problem and can be solved by the CVX toolbox [32]. Then, the SV of the SOI and weight vector can be formulated as respectively.

3.4.
Steps of Proposed Algorithm. In our approach, the computational complexity of the SPM construction is Oðmax

Simulation Results
To evaluate the capability of the proposed algorithm, a ULA with M = 10 is considered. The SOI is presumed to impinge on the ULA from −5°, the angular sector of the SOI is set to Θ s = ½−9°,−1°, and three interferences arrive at the array from −30°, 24°, and 28°with 25 dB. When comparing the performance of robust beamformers in terms of the output SINR versus input SNR, the number of snapshots is set as K = 50. In comparing the SINR and snapshots, the SNR is fixed to 10 dB. All angular sectors are sampled with a fixed  [14], INCM-LINE method [16], INCM-SPSS approach [17], INCM-RNE beamformer [18], INCM-IMR beamformer [21], and INCM-ISV method [22]. For the EIG algorithm, the number of sources is assumed to be known. To achieve superior performance, accurate estimates of SV and INCM are necessary. In this regard, the proposed method applies a high-resolution SSE algorithm to estimate the DOA of interferences, and the INCM is reconstructed through the theoretical form. Consequently, its null depths in the interference directions are deeper than those of the other methods, and its performance is similar to that of the optimal beamformer. Although nulls can be formed by the INCM-IMR beamformer in all interference directions, their depths are shallow. Moreover, because the INCM-LINE beamformer integrates the region where the interferences are located, the reconstructed INCM contains all the information of the three interferences, ensuring that nulls are formed in all the interference directions. However, their null depths are shallower than that created by the proposed algorithm.

Example 2 (Mismatch due to Look Direction Error).
In the second example, the impact of random look direction mismatch on the performance of beamformers is considered. For each simulation run, the DOA mismatches of incident signals are uniformly distributed in ½−4°, 4°. Figure 7 illustrates the output SINR versus input SNR of all beamfor-mers. This figure shows that the proposed algorithm and INCM-IMR achieve the highest output SINR at a low SNR. Further, the performance of the EIG and INCM-IMR methods is severely degraded at a high SNR. The output SINR of the proposed method and INCM-LINE are higher than those of the other algorithms because the reconstructed Step 1. Calculate the SCM applying (10), and estimate the spatial spectrum by (18).
Step 5. Calculate the complex vector of the proposed algorithm through (31).
Algorithm 1: Steps of the proposed algorithm.

Example 3 (Mismatch due to Array Geometry Error).
In the third example, a scenario with a mismatch due to array geometry error is examined. Generally, the array geometry error is modeled as the sensor position perturbation, which is assumed to be uniformly distributed in ½−0:05λ, 0:05λ (i.e., the sensor position, d, used in each run is changed in interval ½d − 0:05λ, d + 0:05λ). Figure 9 depicts the output SINR versus input SNR of all beamformers. The figure shows that the proposed method performs best among the tested algorithms; its performance is similar to that of the INCM-IMR method at a low SNR. Moreover, the foregoing mismatch type has different extents of impact on the other algorithms (e.g., it severely degrades the performance of INCM-RNE). Figure 10 displays the curves of output SINR versus the number of snapshots. It also indicates that the proposed method has stable performance provided that the number of snapshots exceeds 30. However, the output SINR of INCM-LINE considerably fluctuates with the change in the number of snapshots.

Example 4 (Mismatch due to Gain and Phase Error).
In this example, the influence of arbitrary gain and phase perturbation on the RAB algorithms is considered. The m term in (2) can be modeled as where κ m denotes the gain error derived from the Gaussian distribution Nð1, 0:05 2 Þ and β m is the phase perturbation distributed in Nð0, ð5°Þ 2 Þ. Figures 11 and 12     4.5. Example 5 (Mismatch due to Incoherent Local Scattering). In the last example, the incoherent local scattering of the SOI is considered. Assume that the SOI has a time-varying signature. Accordingly, the model can be expressed as where θ η represents the scattering angle drawn from the Gaussian generator Nðθ 0 , 4°Þ, a η ðθ η Þ is the SV corresponding to θ η , and s η ðkÞ, η = 0, 1, 2, 3, 4, is independently drawn from Nð0, 1Þ. BecauseR s is no longer a rank-one matrix, the SINR in this scenario can be written as follows:     Wireless Communications and Mobile Computing which can be maximized by where Φf⋅g represents the eigenvector of the matrix corresponding to the largest eigenvalue. Figures 13 and 14 show the output SINR versus input SNR and output SINR versus the number of snapshots, respectively. The two figures indicate that the proposed approach outperforms the other algorithms at a high SNR and considerably approximates the optimal SINR. A slight perturbation is observed in the proposed algorithm when the number of snapshots is less than 30. However, as the number of snapshots increases, the output SINR quickly exceeds that of the INCM-LINE beamformer.

Conclusions
This paper proposes a novel RAB algorithm formulated using the subspace projection method and SSE. Without using the CSE, the proposed method first obtains the SPM through the integral of the region where the SOI is located. Based on this, the SCM is projected onto the SPM to remove the SOI component. Then, the DOA of the sources is estimated by the MUSIC algorithm, and the estimation method of interference powers is derived. Subsequently, the INCM is reconstructed by the estimated powers, and the SVs of the interferences and SOI are optimized by solving a QCQP problem. Finally, the weight vector of the proposed algorithm is obtained. Simulation results indicate that the proposed RAB algorithm is superior to the other algorithms in terms of the overall performance in the following cases: comparison of beam patterns, look direction mismatch, sensor position perturbation, gain and phase error, and incoherent local scattering. In future work, the INCM-SSE algorithm can be applied to the RAB problem for other geometries of antenna arrays such as L-shape, planar, and circular array.

Data Availability
The data used to support the findings of this study are included within the article.

Conflicts of Interest
The authors declare that there is no conflict of interests regarding the publication of this paper.