VIEWS: 13 PAGES: 5 CATEGORY: Lifestyle POSTED ON: 8/14/2010 Public Domain
LARGE DIMENSIONAL RANDOM MATRIX THEORY FOR SIGNAL DETECTION AND ESTIMATION IN ARRAY PROCESSING∗ J. W. Silverstein† and P. L. Combettes‡ † Department of Mathematics, North Carolina State University, Raleigh, NC 27695, USA. ‡ Department of Electrical Engineering, City College and Graduate School, City University of New York, New York, NY 10031, USA. ABSTRACT sensors are modeled as observations of the random vec- tor X(t) = AS(t) + N (t), t ∈ [0, +∞[, where A is a In this paper, we bring into play elements of the spec- p × q complex matrix depending on the geometry of tral theory of large dimensional random matrices and the array and the parameters of the signals, and is demonstrate their relevance to source detection and assumed to have rank q. bearing estimation in problems with sizable arrays. The detection problem is to estimate q from the ob- These results are applied to the sample spatial covari- servation of n snapshots (X(ti ))1≤i≤n of the data pro- ance matrix, R, of the sensed data. It is seen that de- cess. Under the above assumptions, the random vec- tection can be achieved with a sample size considerably tors (X(t))t∈[0,+∞[ are i.d. with spatial covariance ma- less than that required by conventional approaches. As trix R = EX(0)X(0)∗ = ARS A∗ + σ2 Ip , where Ip de- regards to determining the directions of arrivals, it is notes the p × p identity matrix. Moreover, the p − q argued that more accurate estimates can be obtained smallest eigenvalues of R are equal to σ2 . These eigen- by constraining R to be consistent with various a pri- values are referred to as the noise eigenvalues and the ori constraints, including those arising from large di- remainder of the spectrum is referred to as the signal mensional random matrix theory. A set theoretic for- eigenvalues. Since R is not known its spectrum must malism is used to formulate this feasibility problem. be inferred from observing that of the sample covari- ance matrix R = (1/n) i=1 X(ti )X(ti )∗ . Loosely n Unsolved issues are discussed. speaking, one must then decide where the observed PROBLEM STATEMENT spectrum splits into noise and signal eigenvalues. The estimation problem is to determine the direction We consider the problem of detecting the number q of arrivals (θi )1≤i≤p of the sources. Under standard of sources impinging on an array of p (q < p) sen- hypotheses, this problem can be solved via the MUSIC sors as well as their directions of arrival, when p is method [11], which requires only the knowledge of R. large. The model for the data formation mechanism is In practice, however, only R is available, which leads the following. At each time t, the j-th signal present to poor estimates if n is not suﬃciently large. in the scene, the additive noise at the i-th sensor, and the received data at the i-th sensor are respec- In this paper, we bring into play elements of the spec- tively represented by the square-integrable complex- tral theory of large dimensional random matrices and valued r.v.’s Sj (t), Ni (t), and Xi (t). The random vec- demonstrate their relevance to source detection and tors (S(t) = [S1 (t) . . . Sq (t)]T )t∈[0,+∞[ are i.d. with estimation. ES(0) = 0 and nonsingular spatial covariance matrix RS = ES(0)S(0)∗ . Moreover, it is assumed that the SPECTRAL THEORY OF LARGE r.v.’s (Ni (t) | 1 ≤ i ≤ p , t ∈ [0, +∞[) are i.i.d. with DIMENSIONAL RANDOM MATRICES EN1 (0) = 0 and E|N1 (0)|2 = σ2 , where σ2 is unknown, and independent from the r.v.’s (Sj (t) | 1 ≤ j ≤ q , t ∈ Let M be an m × m random matrix with real-valued [0, +∞[). Let N (t) = σW (t) = σ[W1 (t) . . . Wp (t)]T eigenvalues (Λi )1≤i≤m . Then, the empirical distribu- (so that the Wi (t)’s are standardized) and X(t) = tion function (d.f.) of (Λi )1≤i≤m is the stochastic pro- [X1 (t) . . . Xp (t)]T . The data collected by the array of cess m ∗ 1 The work of the ﬁrst author was supported by NSF (∀x ∈ IR) F M (x) = 1]−∞,x] (Λi ). (1) grant DMS-8903072. m i=1 We now review the main result, a limit theorem found eigenvalues of the sample spatial covariance matrix R in [15]. to each other. This would require an extraordinar- Theorem 1. [15] Let (Yij )i,j≥1 be i.i.d. real-valued ily large sample size (sometimes unattainable) when r.v.’s with E|Y11 − EY11 |2 = 1. For each m in IN∗ , let the number of sources is large in order to obtain a Ym = [Yij ]m×n , where n = n(m) and m/n → y > 0 as good estimate. Under additional assumptions on the m → +∞, and let Tm be an m×m symmetric nonneg- signals (including the independence of the snapshots), ative deﬁnite random matrix independent of the Yij ’s Theorem 1 is used in [14] to show that, for p and n for which there exists a sequence of positive numbers suﬃciently large, with high probability, the empirical (µk )k≥1 such that for each k in IN∗ d.f. F R is close to the d.f. F of Theorem 1 for m = p, ∗ 2 +∞ 1 y = p/n, and H = F ARS A +σ Ip . k a.s. xk dF Tm (x) = trTm −→ µk as m → +∞ (2) Further analysis when H is of this form shows that a 0 m value y can be computed for which y ∈]0, y[ if and only ˜ ˜ and where the µk ’s satisfy Carleman’s suﬃciency con- if the support of F splits into at least two intervals, −1/2k dition, k≥1 µ2k = +∞, for the existence and the with the leftmost interval having mass (p − q)/p. For uniqueness of the d.f. H having moments (µk )k≥1 . Let instance, in the simulations performed in [14], p = 50 T Mm = (1/n)Ym Ym Tm . Then, a.s., (F Mm )m≥1 con- ˜ and y is found to be 1.058, which can allow a rel- verges weakly to a nonrandom d.f. F having moments atively small sample size. However, the simulations k suggest something much stronger than the splitting of k! (∀k ∈ IN∗ ) νk = yk−w µm1 · · ·µmw , the support of F is occurring, namely exact splitting w=1 m1 !· · ·mw !w! 1 w of the eigenvalues. Thus, upon mathematical veriﬁca- (3) tion of this phenomena the following appears valid: R where the inner sum extends over all w-tuples of non- with a sample size on the same order of magnitude as w negative integers (m1 , . . . , mw ) such that i=1 mi = the number of sensors will have eigenvalues noticeably w k −w +1 and i=1 imi = k. Moreover, these moments split into two groups, the group lying to the left corre- uniquely determine F .1 3 sponding to the multiplicity of the smallest eigenvalue The following theorem pertains to the case when Tm of the true covariance matrix. Detection can thus be is a multiple of the identity matrix. achieved with a sample size considerably less than that required by previous approaches. Theorem 2. When Tm = σ2 Im , F is known, having an algebraic density on the positive reals with support √ √ APPLICATION TO BEARING [σ2 (1 − y)2 , σ2 (1 + y)2 ] [7], [8], [10]. The largest eigenvalue of Mm converges almost surely [respectively ESTIMATION √ in probability] to σ2 (1+ y)2 as m → +∞ if and only if EY11 = 0 and E|Y11| < +∞ [respectively x4 P{|Y11| ≥ 4 Under our basic assumptions, the directions of arrival x} → 0 as x → +∞] [1], [6], [13], [16]. Moreover, if can be calculated from the spatial covariance matrix Y11 is standardized Gaussian, the smallest eigenvalue R via the MUSIC algorithm [11]. In practice, short √ of knowing R, one must carry out the computation of Mm converges almost surely to σ2 (1 − y)2 when y < 1 [12]. 3 2 based on an observation of the sample covariance ma- trix R. Because R is often a poor approximation of Several additional results on F are mentioned in [14], R, the method can be improved by replacing R by among them the convergence of F to H as y → 0, and a matrix that satisﬁes all the a priori constraints on a way to compute the support of F from y and H. the problem before applying MUSIC. By invoking the formalism of set theoretic estimation [3] and denoting APPLICATION TO SIGNAL DETECTION by (Ψi )1≤i≤M the constraints, this feasibility problem Existing approaches, such as those based on informa- ˜ can be formulated as that of ﬁnding a matrix R in the tion theoretic criteria, rely on the closeness of the noise matrix subset 1 The proof in [15] can easily be modiﬁed to allow M complex-valued entries in Ym and Tm , giving the same S= Si where Si = {Q | Q satisﬁes Ψi }. (4) result, provided Tm is Hermitian and we take Mm = i=1 ∗ (1/n)Ym Ym Tm . 2 Results on the extreme eigenvalues have been veriﬁed In general, ﬁnding a point in S directly is impossi- for Y11 real-valued but, again, the proofs can be extended ble. Let Πi be the projection map onto Si , i.e., Πi (Q) to the complex case. is the set of matrices in Si which lie closest to Q in terms of a distance (for computational tractability, we second absolute moment σ 2 . A direct application of shall take the Frobenius distance). Then, under cer- the analysis of [5] to their sample mean leads to a set tain conditions on the sets and the initial point Q0 , of the type the sequence (Qn )n≥0 where Qn+1 ∈ Πin (Qn ) with 2np ˜ in = n (modulo M ) + 1 will converge to a point Q in S ˜ Si = {H | | ˜ yk (H)| ≤ δi }, (5) [4]. In this scheme, the sets (Si )1≤i≤M are activated k=1 in a cyclic manner and the update is any projection of the current estimate onto the next set.3 ˜ where y(H ) is the vector obtained by stacking the real First of all we can pose the problem in the R-space ˜ and imaginary parts of the entries of Y (H ). Sets based and construct sets of estimates of the true covariance ˜ on other statistics of the entries of Y (H) can be ob- matrix. Under the above assumptions, an obvious a tained in a like-manner. This H-space framework also priori constraint about R is that its rank is q. One makes the use of constraints arising from the prop- can therefore consider the (closed, nonconvex) set S1 erties of large dimensional random matrices possible. of matrices whose rank is at most q. Other constraints Indeed, according to Theorem 2, a bound can be ob- may arise from the geometry of the array. Thus, if the tained on the largest (in the Gaussian case also on the array is linear with equispaced sensors, R will have a smallest) singular value of Y (H). In the case of the Toeplitz structure and one can take S2 to be the sub- largest singular value, one obtains the set space of Toeplitz matrices. The use of a projection √ onto the set S2 has been reported in several array pro- Si = {H | ||Y (H)||2 ≤ n(σ2 (1 + y)2 + i )}, ˜ ˜ S (6) cessing applications, e.g., [9] and is often referred to as Toeplitzation. The joint use of S1 and S2 via alternat- where · S denotes the spectral norm and i reﬂects ing projections was proposed in [2]. It should be noted a conﬁdence limit. Of course, all the constraints on that in such a process, loss of positive deﬁniteness may ARS A∗ mentioned previously (rank, structure (e.g., occur. Therefore, a third set should be added, namely Toeplitz), etc.) can still be used in the H-space. For that of positive deﬁnite matrices. In simulations, no- example, the set associated with the Toeplitz con- ticeable improvements have been reported by using straint reads ˜ the constrained covariance matrix R instead of its raw Si = {H | H H ∗ is Toeplitz }. ˜ ˜ ˜ (7) counterpart R, especially if the number of samples n and the signal-to-noise ratio are low. OPEN PROBLEMS In the above approach, one seeks to estimate R di- rectly, which limits the exploitation of information that There are several mathematical questions that need may be available about the noise. An alternative is to to be addressed concerning the behavior of large di- estimate the noiseless p×n data matrix H = n−1/2 AS mensional random matrices in connection with appli- ˜ in the model X = AS + N .4 An estimate H of H can cations to the above array processing problems. The be synthesized by projection onto various constraint three most relevant are outlined below. sets. One can then form the constrained estimate of ARS A∗ , i.e., R = HH ∗ , and apply MUSIC to it. Let ˜ Extending Theorem 1. The application of Theo- us now consider the constraints that can be imposed rem 1 requires two assumptions on the formation of on H. To this end, deﬁne the residual matrix for a the signal vector S(t). The ﬁrst is that S(t) = CV (t), ˜ ˜ given estimate H of H to be Y (H) = X − n1/2 H. ˜ where C is a ﬁxed q × q nonsingular matrix, and V (t) Note that we have Y (H) = N . Therefore, all the sta- is made up of i.i.d. random variables with same dis- tistical information pertaining to N can be imposed tribution as the noise components. Since signal and ˜ on H and several sets can be constructed according to noise components are typically assumed to be Gaus- this principle [5]. For instance, with our assumptions sian, this does not appear to be a major problem. The ˜ on the noise, the entries of Y (H) should look like i.i.d. second assumption is the independence of the signal samples from a zero mean complex population with vector across snapshots. This is more serious, even though independent samples are assumed in several 3 In particular, if all the sets are closed and convex in a mathematical treatments (for example, the work on Hilbert space, convergence to a feasible point is guaranteed computing q from information theoretic criteria), and for any Q0 . in most simulations found in the literature. The possi- 4 In this model, the ith (1 ≤ i ≤ n) column of the matri- bility of extending Theorem 1 to Ym having stationary ces X, S, and N are X(ti), S(ti ) and N (ti ), respectively. columns needs to be investigated. Eigenvalue Splitting. In the detection problem, the [2] J. A. Cadzow, “Signal Enhancement - A Composite observance of exact eigenvalue splitting in the simula- Property Mapping Algorithm,” IEEE Transactions tions is striking. A limit property stronger than weak on Acoustics, Speech, and Signal Processing, vol. 36, convergence of d.f.’s must be in eﬀect, and a proof no. 1, pp. 49-62, January 1988. of this phenomenon should be pursued. The result is [3] P. L. Combettes and M. R. Civanlar, “The Founda- basically an extension of Theorem 2 on the extreme tions of Set Theoretic Estimation,” ICASSP Proceed- T ings, pp. 2921-2924. Toronto, Canada, May 14-17, eigenvalues of (1/n)Ym Ym . 1991. Rate of Convergence. On top of these problems is the general question of how fast the quantities are [4] P. L. Combettes and H. J. Trussell, “Method of Suc- cessive Projections for Finding a Common Point of approaching their limiting values. The simulations Sets in Metric Spaces,” Journal of Optimization The- performed in [14] show that for p = 50 separation of ory and Applications, vol. 67, no. 3, pp. 487-507, De- the noise and signal eigenvalues of R are in agreement cember 1990. with the support separation in F . Preliminary anal- T [5] P. L. Combettes and H. J. Trussell, “The Use of Noise ysis on F (1/n)YmYm suggests a rate of convergence of Properties in Set Theoretic Estimation,” IEEE Trans- 1/m, supporting the view that limiting behavior can actions on Signal Processing, vol. 39, no. 7, pp. 1630- be evidenced for m not exceedingly high. 1641, July 1991. Two additional problems are mentioned here. [6] S. Geman, “A Limit Theorem for the Norm of Ran- Computation of the Projections. A shortcoming dom Matrices,” The Annals of Probability, vol. 8, no. 2, pp. 252-261, April 1980. of the proposed set theoretic approach to the deter- mination of the directions of arrivals is the numerical [7] U. Grenander and J. W. Silverstein, “Spectral Anal- tedium sometimes involved in computing the projec- ysis of Networks with Random Topologies,” SIAM tions at each iteration. In general, the sets are given Journal on Applied Mathematics, vol. 32, no. 2, pp. 499-519, March 1977. in the form Si = {Q | gi (Q) ≤ δi }, where gi is a given functional. A projection of a matrix Q onto Si is ob- [8] D. Jonsson, “Some Limit Theorems for the Eigenval- tained by solving the minimization problem5 ues of a Sample Covariance Matrix,” Journal of Multi- variate Analysis, vol. 12, no. 1, pp. 1-38, March 1982. min Q − Q subject to gi (Q) = δi , (8) [9] J. P. Lecadre and P. Lopez, “Estimation d’une Ma- Q e trice Interspectrale de Structure Impos´e,” Traite- which can be approached via the method of Lagrange e ment du Signal, vol. 1, pp. 4-17, D´cembre 1984. multipliers. However, in the case when Si is not con- c [10] V. A. Marˇenko and L. A. Pastur, “Distribution vex, local minima may occur. In such cases, eﬃcient of Eigenvalues for Some Sets of Random Matrices,” global minimization methods should be developed to Mathematics of the USSR - Sbornik, vol. 1, no. 4, pp. solve (8). 457-483, 1967. [11] R. O. Schmidt, “Multiple Emitter Location and Signal Convergence to a Feasible Point. Due to the pres- Parameter Estimation,” IEEE Transactions on An- ence of nonconvex sets, the convergence of the succes- tennas and Propagation, vol. AP-34, no. 3, pp. 276- sive projection algorithm to a feasible point cannot be 280, March 1986. guaranteed for any initial estimate [4]. Although start- [12] J. W. Silverstein, “The Smallest Eigenvalue of a Large ing the iterations at the point provided by the data Dimensional Wishart Matrix,” The Annals of Proba- (i.e., R0 = R in the R-space approach or H0 = n−1/2 X bility, vol. 13, no. 4, pp. 1364-1368, November 1985. in the H-space approach) constitutes a sensible choice, it does not always ensure convergence. Therefore, the [13] J. W. Silverstein, “On the Weak Limit of the Largest Eigenvalue of a Large Dimensional Sample Covariance convergence issue deserves further investigation. Matrix,” Journal of Multivariate Analysis, vol. 30, no. 2, pp. 307-311, August 1989. REFERENCES [14] J. W. Silverstein and P. L. Combettes, “Signal Detec- tion via Spectral Theory of Large Dimensional Ran- [1] Z. D. Bai, J. W. Silverstein, and Y. Q. Yin, “A Note dom Matrices,” IEEE Transactions on Signal Process- on the Largest Eigenvalue of a Large Dimensional ing, vol. 40, no. 8, August 1992. Sample Covariance Matrix,” Journal of Multivariate Analysis, vol. 26, no. 2, pp. 166-168, August 1988. [15] Y. Q. Yin, “Limiting Spectral Distribution for a Class of Random Matrices,” Journal of Multivariate Anal- 5 As mentioned before, computational tractability seems ysis, vol. 20, no. 1, pp. 50-68, October 1986. to favor the use of the Frobenius norm. [16] Y. Q. Yin, Z. D. Bai, and P. R. Krishnaiah, “On the Limit of the Largest Eigenvalue of the Large Dimen- sional Sample Covariance Matrix,” Probability Theory and Related Fields, vol. 78, no. 4, pp. 509-521, Au- gust 1988.