Fundamental limit of sample generalized eigenvalue based detection of signals in noise using relatively few signal-bearing and noise-only samples

I Introduction

The observation vector, in many signal processing applications, can be modelled as a superposition of a finite number of signals embedded in additive noise. The model order selection problem of inferring the number of signals present is the critical first step in the subsequent signal parameter estimation problem. We consider the class of estimators that determine the model order, i.e., the number of signals, in colored noise from the sample generalized eigenvalues of the signal-plus-noise sample covariance matrix and the noise-only sample covariance matrix pair. The sample generalized eigenvalues precisely correspond to the eigenvalues of the matrix formed by “whitening” the signal-plus-noise sample covariance matrix with the noise-only sample covariance matrix (assuming that the number of noise-only samples is greater than the dimensionality of the system so that the noise-only sample covariance matrix is invertible).

Such estimators are used in settings where it is possible to find a portion of the data that contains only noise fields and does not contain any signal information. This is a realistic assumption for many practical applications such as evoked neuromagnetic experiments , geophysical experiments that employ a “thumper” or in underwater experiments with a wideband acoustic signal transducer where such a portion can be found in a data portion taken before a stimulus is applied. In applications such as radar or sonar where the signals of interest are narrowband and located in a known frequency band, snapshot vectors collected at a frequency just outside this band can be justified as having the same noise covariance characteristics assuming that we are in the stationary-process-long-observation-time (SPLOT) regime .

Our main objective in this paper is to shed new light on this age old problem of detecting signal in noise from finite samples using the sample eigenvalues alone . We bring into sharp focus a fundamental statistical limit that explains precisely when and why, in high-dimensional, sample size limited settings underestimation of the model order is unavoidable. This is in contrast to works in the literature that use simulations, as in , to highlight the chronically reported symptom of model order estimators underestimating the number of signals without providing insight into whether a fundamental limit of detection is being encountered.

In recent work , we examined this problem in the white noise scenario. The main contribution of this paper is the extension of the underlying idea to the arbitrary (or colored) noise scenario. Analogous to the definition in , we define the effective number of identifiable signals in colored noise as the number of the generalized eigenvalues of the population (true) signal-plus-noise covariance matrix and noise-only covariance matrix pair that are greater than a (deterministic) threshold that is a simple function of the number of signal-plus-noise samples, noise-only samples and the dimensionality of the system. Analogous to the white noise case, increasing the dimensionality of the system, by say adding more sensors, raises the detectability threshold so that the effective number of identifiable signals might actually decrease.

An additional contribution of this paper is the development of a simple, new, algorithm for estimating the number of signals based on the recent work of Johnstone . Numerical results are used to illustrate the performance of the estimator around the detectability threshold alluded to earlier. Specifically, we observe that if the eigen-SNR of a signal is above a critical value then reliable detection using the new algorithm is possible. Conversely, if the eigen-SNR is below the critical value then the algorithm, correctly for the reason described earlier, is unable to distinguish the signal from noise.

The paper is organized as follows. We formulate the problem in Section II and state the main result in Section III. The effective number of signals is defined in Section III-A along with a discussion on its implications for applications such as array processing, sensor networks and machine learning. A new algorithm for detecting the number of signals is presented in Section IV. Concluding remarks are offered in Section V. The mathematical proofs of the main result are provided in Section VI.

II Problem formulation

We observe $m$ samples (“snapshots”) of possibly signal bearing $n$-dimensional snapshot vectors ${\bf x}_{1},\ldots,{\bf x}_{m}$ where for each $i$, the snapshot vector has a (real or complex) multivariate normal distribution, i.e., ${\bf x}_{i}\sim\mathcal{N}_{n}(0,{\bf R})$ and the ${\bf x}_{i}$’s are mutually independent. The snapshot vectors are modelled as

where ${\bf z}_{i}\sim\mathcal{N}_{n}(0,\Sigma)$, denotes an $n$-dimensional (real or complex) Gaussian noise vector where the noise covariance $\Sigma$ may be known or unknown, ${\bf s}_{i}\sim\mathcal{N}_{k}({\bf 0},{\bf R}_{s})$ denotes a $k$-dimensional (real or complex) Gaussian signal vector with covariance ${\bf R}_{s}$, and ${\bf A}$ is a $n\times k$ unknown non-random matrix. Since the signal and noise vectors are independent of each other, the covariance matrix of ${\bf x}_{i}$ can hence be decomposed as

with ′ denoting the complex conjugate or real transpose. Assuming that the matrix ${\bf A}$ is of full column rank, i.e., the columns of ${\bf A}$ are linearly independent, and that the covariance matrix of the signals ${\bf R}_{s}$ is nonsingular, it follows that the rank of ${\bf\Psi}$ is $k$. Equivalently, the $n-k$ smallest eigenvalues of ${\bf\Psi}$ are equal to zero.

If the noise covariance matrix $\bm{\Sigma}$ were known apriori and was non-singular, a “noise whitening” transformation may be applied to the snapshot vector ${\bf x}_{i}$ to obtain the vector

which will also be normally distributed with covariance

Denote the eigenvalues of ${\bf R}_{\bm{\Sigma}}$ by $\lambda_{1}\geq\lambda_{2}\geq\ldots\geq\lambda_{n}$. Recalling the formulation of the generalized eigenvalue problem [Section 8.7], we note that the eigenvalues of ${\bf R}_{\bm{\Sigma}}$ are exactly the generalized eigenvalues of the regular matrix pair $(\widehat{{\bf R}},\widehat{\bm{\Sigma}})$. Then, assuming that the rank of $\bm{\Sigma}^{-1}\bm{\Psi}$ is also $k$, it follows that the smallest $n-k$ eigenvalues of ${\bf R}_{\bm{\Sigma}}$ or, equivalently, the generalized eigenvalues of the matrix pair $({\bf R},\bm{\Sigma})$), are all equal to $1$ so that

while the remaining $k$ eigenvalues ${\bf R}_{\bm{\Sigma}}$ of will be strictly greater than one.

Thus, if the true signal-plus-noise covariance matrix ${\bf R}$ and the noise-only covariance matrix $\bm{\Sigma}$ were known apriori, the number of signals $k$ could be trivially determined from the multiplicity of the eigenvalues of ${\bf R}_{\bm{\Sigma}}$ equalling one.

The problem in practice is that the signal-plus-noise and the noise covariance matrices ${\bf R}$ are unknown so that such a straight-forward algorithm cannot be used. Instead we have an estimate the signal-plus-covariance matrix obtained as

and an estimate of the noise-only sample covariance matrix obtained as

where ${\bf x}_{i}$ for $i=1,\ldots,m$ are (possibly) signal-bearing snapshots and ${\bf z}_{j}$ for $j=1,\ldots,N$ are independent noise-only snapshots. We assume here that the number of noise-only snapshots exceeds the dimensionality of the system, i.e., $N>n+1$, so that the noise-only sample covariance matrix $\widehat{\bm{\Sigma}}$, which has the Wishart distribution , is non-singular and hence invertible with probability 1 [12, Chapter 3, pp. 97],[13, Chapter 7.7, pp. 272-276]. Following (5), we then form the matrix

and compute its eigen-decomposition to obtain the eigenvalues of $\widehat{{\bf R}}_{\widehat{\Sigma}}$, which we denote by $\hat{\lambda_{1}}\geq\hat{\lambda}_{2}\geq\ldots\geq\hat{\lambda}_{n}$. We note, once again, that the eigenvalues of $\widehat{{\bf R}}_{\widehat{\Sigma}}$ are simply the generalized eigenvalues of the regular matrix pair $(\widehat{{\bf R}},\widehat{\bm{\Sigma}})$. Note that whenever $N<n$, the signal-plus-noise sample covariance matrix ${\bf R}$ will be singular so that the $n-N$ generalized eigenvalues will equal zero, i.e., $\hat{\lambda}_{N+1}=\hat{\lambda}_{N+2}=\ldots=\hat{\lambda}_{n}=0$. Figure 1 illustrates why the blurring of the sample eigenvalues relative to the population eigenvalues makes the problem more challenging.

In this paper, we are interested in the class of algorithms that infer the number of signals buried in arbitrary noise from the eigenvalues of $\widehat{{\bf R}}_{\widehat{\bm{\Sigma}}}$ or $\widehat{{\bf R}}_{\bm{\Sigma}}$ alone. Such algorithms are widely used in practice and arise naturally from classical multivariate statistical theory where the matrix $\widehat{{\bf R}}_{\widehat{\Sigma}}$ is referred to as the multivariate F matrix . The information theoretical approach to model order estimation, first introduced by Wax and Kailath , was extended to the colored noise setting by Zhao et al in who prove consistency of their estimator in the large sample size regime; their analysis does not yield any insight into the finite sample setting.

Consequently, research has focussed on developing sophisticated techniques for improving performance of eigenvalue based methods in the finite sample setting. Zhu et al improve the performance of their eigenvalue estimator by assuming a model for the noise covariance matrix. Stoica and Cedervall improve the performance of their estimator in two reasonable settings: one, where it is reasonable to assume that the noise covariance matrix is block diagonal or banded and two, where the temporal correlation of the noise has a shorter length than the signals. Other techniques in the literature exploit other characteristics of the signal or noise to effectively reduce the dimensionality of the signal subspace and improve model order estimation given finite samples. See for example and the references in .

Informally speaking, it is evident that performance of such model order estimation algorithms is coupled to the “quality” of the estimated signal-plus-noise and noise-only covariance matrices which in turn are dependent on the number of snapshots used to estimate them, respectively. Researchers applying these techniques have noted the absence of a mathematically rigorous, general purpose formula in the literature for predicting the minimum number of samples needed to obtain “good enough” detection accuracy (see, for example [pp. 846]. A larger, more fundamental question that has remained unanswered, till now, is whether there is a statistical limit being encountered.

We tackle this problem head on in this paper by employing sophisticated techniques from random matrix theory in . We show that in an asymptotic sense, to be made precise later, that only the “signal” eigenvalues of ${\bf R}_{\Sigma}$ that are above a deterministic threshold can be reliably distinguished from the “noise” eigenvalues. The threshold is a simple, deterministic function of the the dimensionality of the system, the number of noise-only and signal-plus-noise snapshots, and the noise and signal-plus noise covariance, and described explicitly next. Note the applicability of the results to the situation when the signal-plus-noise covariance matrix is singular.

III Main result

For a Hermitian matrix ${\bf A}$ with $n$ real eigenvalues (counted with multiplicity), the empirical distribution function (e.d.f.) is defined as

Of particular interest is the convergence of the e.d.f. of $\widehat{{\bf R}}_{\widehat{\Sigma}}$ in the signal-free case, which is described next.

Let $\widehat{{\bf R}}_{\widehat{\Sigma}}$ denote the matrix in (9) formed from $m$ (complex Gaussian) noise-only snapshots and $N$ independent noise-only (complex Gaussian) snapshots. Then the e.d.f. $F^{\widehat{{\bf R}}_{\widehat{\Sigma}}}(x)\to F^{R_{\Sigma}}(x)$ almost surely for every $x$, as $m,n(m)\to\infty$, $m,N(m)\to\infty$ and $c_{m}=n/m\to c>0$ and $c^{1}_{N}=n/N\to c_{1}<1$ where

This result was proved in . When $c_{1}\to 0$ we recover the famous Marčenko-Pastur density . ∎

The following result exposes when the “signal” eigenvalues are asymptotically distinguishable from the “noise” eigenvalues.

Let $\widehat{{\bf R}}_{\widehat{\Sigma}}$ denote the matrix in (9) formed from $m$ (real or complex Gaussian) signal-plus-noise snapshots and $N$ independent (real or complex Gaussian) noise-only snapshots. Denote the eigenvalues of ${{\bf R}}_{{\Sigma}}$ by $\lambda_{1}\geq\lambda_{2}>\ldots\geq\lambda_{k}>\lambda_{k+1}=\ldots\lambda_{n}=1$. Let $l_{j}$ denote the $j$-th largest eigenvalue of $\widehat{{\bf R}}_{\widehat{\Sigma}}$. Then as $n,m(n)\to\infty$, $n,N(n)\to\infty$ and $c_{m}=n/m\to c>0$ and $c^{1}_{N}=n/N\to c_{1}<1$ we have

for $j=1,\ldots,k$ and the convergence is almost surely and the threshold ${\rm T}(c,c_{1})$ is given by

The result follows from Theorem VI.5. The threshold ${\rm T}(c,c_{1})$ is obtained by solving the inequality

where for $j=1,\ldots,k$, $t^{\prime}$, from , is given by

Note that when $c_{1}\to 0$, ${\rm T}(c,c_{1})\to(1+\sqrt{c})$ so that we recover the results of Baik and Silverstein . ∎

Theorem III.2 brings into sharp focus the reason why, in the large-system-relatively-large-sample-size limit, model order underestimation is sometimes unavoidable. This motivates our heuristic definition of the effective number of identifiable signals below:

If we denote the eigenvalues of ${\bf R}_{\bm{\Sigma}}\equiv\bm{\Sigma}^{-1}{\bf R}$ by $\lambda_{1}\geq\lambda_{2}>\ldots\geq\lambda_{k}>\lambda_{k+1}=\ldots\lambda_{n}=1$ then we define the eigen-SNR of the $j$-th signal as $\lambda_{j}-1$ then (15) essentially states that signals with eigen-SNR’s smaller than ${\rm T}(n/m,n/N)$ will be asymptotically undetectable.

Figure 2 shows the eigen-SNR threshold ${\rm T}(c,c_{1})-1$ needed for reliable detection for different values as a function of $c$ for different values of $1/c_{1}$. Such an analytical prediction was not possible before the results presented in this paper. Note the fundamental limit of detection in the situation when the noise-only covariance matrix is known apriori (solid line) and increase in the threshold eigen-SNR needed as the number of snapshots available to estimate the noise-only covariance matrix decreases.

III-B Implications for array processing

Suppose there are two uncorrelated (hence, independent) signals so that ${\bf R}_{s}=\textrm{diag}(\sigma_{{\rm S}1}^{2},\sigma_{{\rm S}2}^{2})$. In (1) let ${\bf A}=[{\bf v}_{1}{\bf v}_{2}]$. In a sensor array processing application, we think of ${\bf v}_{1}\equiv{\bf v}(\theta_{1})$ and ${\bf v}_{2}\equiv{\bf v}_{2}(\theta_{2})$ as encoding the array manifold vectors for a source and an interferer with powers $\sigma_{{\rm S}1}^{2}$ and $\sigma_{{\rm S}2}^{2}$, located at $\theta_{1}$ and $\theta_{2}$, respectively. The signal-plus-noise covariance matrix is given by

where $\bm{\Sigma}$ is the noise-only covariance matrix. The matrix ${\bf R}_{\Sigma}$ defined in (5) can be decomposed as

so we that we can readily note that ${\bf R}_{\Sigma}$ has the $n-2$ smallest eigenvalues $\lambda_{3}=\ldots=\lambda_{n}=1$ and the two largest eigenvalues

respectively, where ${\bf u}_{1}:=\bm{\Sigma}^{-1/2}{\bf v}_{1}$ and ${\bf u}_{2}:=\bm{\Sigma}^{-1/2}{\bf v}_{2}$ . Applying the result in Theorem III.2 allows us to express the effective number of signals as

Equation (18) captures the tradeoff between the identifiability of two closely spaced signals, the dimensionality of the system, the number of available snapshots and the cosine of the angle between the vectors ${\bf v}_{1}$ and ${\bf v}_{2}$. Note that since the effective number of signals depends on the structure of the theoretical signal and noise covariance matrices (via the eigenvalues of ${\bf R}_{\Sigma}$), different assumed noise covariance structures (AR(1) versus white noise, for example) will impact the signal level SNR needed for reliable detection in different ways.

III-C Other applications

There is interest in detecting abrupt change in a system based on stochastic observations of the system using a network of sensors. When the observations made at various sensors can be modeled as Gauss-Markov random field (GMRF), as in , then the conditional independence property of GMRF’s is a useful assumption. The assumption states that conditioned on a particular hypothesis, the observations at sensors are independent. This assumption results in the precision matrix, i.e., the inverse of the covariance matrix, having a sparse structure with many entries identically equal to zero.

Our results might be used to provide insight into the types of systemic changes, reflected in the structure of the signal-plus-noise covariance matrix, that are undetectable using sample generalized eigenvalue based estimators. Specifically, the fact that the inverse of the noise-only covariance matrix will have a sparse structure means that one can experiment with different (assumed) conditional independence structures and determine how “abrupt” the system change would have to be in order to be reliably detected using finite samples.

Spectral methods are popular in machine learning applications such as unsupervised learning, image segmentation, and information retrieval . Generalized eigenvalue based techniques for clustering have been investigated in . Our results might provide insight when spectral clustering algorithms are likely to fail. In particular, we note that the results of Theorem III.2 hold even in situations where the data is not Gaussian (see Theorem VI.5) as is commonly assumed in machine learning applications.

IV An algorithm for reliable detection of signals in noise

In , Johnstone proves that in the signal-free case, the distribution of the largest eigenvalue of $\widehat{{\bf R}}_{\widehat{\bm{\Sigma}}}$, on appropriate centering and scaling, can be approximated to order $O(n^{-2/3})$ by the Tracy-Widom law . In the setting where there are signals present, we expect, after appropriate centering and scaling, the distribution of the signal eigenvalues of $\widehat{{\bf R}}_{\widehat{\bm{\Sigma}}}$ above the detectability threshold will obey a Gaussian law whereas those below the detectability threshold will obey the Tracy-Widom law as in the signal-free case. An analogous results for the signal bearing eigenvalues of $\widehat{{\bf R}}_{\bm{\Sigma}}$ was proved by Baik et al and El Karoui . Numerical investigations for (see Figure 3) corroborate the accuracy of our asymptotic predictions and form the basis of Algorithm 1 presented below for estimating the number of signals at (asymptotic) significance level $\alpha$. Theoretical support for this observation remains incomplete.

Figure 4 illustrates the accuracy of the predicted statistical limit and the ability of the proposed algorithm to reliably detect the presence of the signal at this limit.

V Conclusion

Figure 4 captures the fundamental statistical limit encountered when attempting to discriminate signal from noise using finite samples. Simply put, a signal whose eigen-SNR is below the detectability threshold cannot be reliably detected while a signal above the threshold can be. In settings such as wireless communications and biomedical signal processing where the signal power is controllable, our results provide a prescription for how strong it needs to be so that it can be detected. If the signal level is barely above the threshold, simply adding more sensors might actually degrade the performance because of the increased dimensionality of the system. If, however, either due to clever signal design or physics based modeling, we are able to reduce (or identify) the dimensionality of the subspace spanned by signal, then according to Figure 4 the detectability threshold will also be lowered. With VLSI advances making sensors easier and cheaper to deploy, our results demonstrate exactly why the resulting gains in systemic performance will more than offset the effort we will have to invest in developing increasingly more sophisticated dimensionality reduction techniques. Understanding the fundamental statistical limits of techniques for signal detection in the setting where the noise-only sample covariance matrix is singular remains an important open problem.

Acknowledgements

Raj Rao was supported by an Office of Naval Research Post-Doctoral Fellowship Award under grant N00014-07-1-0269. Jack Silverstein was supported by the U.S. Army Research Office under Grant W911NF-05-1-0244. R. R. thanks Arthur Baggeroer for encouragement and invaluable feedback. This material was based upon work supported by the National Science Foundation under Agreement No. DMS-0112069. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the National Science Foundation.

We thank Alan Edelman for his encouragement and support, Interactive Supercomputing, Inc. for providing access to the Star-P parallel computing software and Sudarshan Raghunathan of Interactive Supercomputing, Inc. for his patience and support in answering our multitude of Star-P programming queries. We remain grateful to Al Davis and Chris Hill of MIT for granting us access to the Darwin Project computing cluster. Thanks to their involvement we were able to program, debug and complete the computation needed to produce Figure 4 in 4 days! Without their gracious help, the computation would have taken 3 months on the latest single processor laptop. We thank Folkmar Bornemann for providing the MATLAB code for computing the percentiles in Table II.

VI Appendix

Let for $i,j=1,2,\ldots$, $X_{ij}$ be a collection of complex valued i.i.d. random variables with $\text{\pe E}X_{1\,1}=0$ and $\text{\pe E}|X_{1\,1}|^{2}=1$. For positive integers $n$ and $m$ let ${\bf X}_{n}=(X_{ij})$, $i=1,2,\ldots,n$, $j=1,2,\ldots,m$. Assume for each $n$ ${\bf T}_{n}$ is an $n\times n$ Hermitian nonnegative definite matrix. The matrix

where ${\bf T}^{1/2}_{n}$ is any Hermitian square root of ${\bf T}_{n}$, can be viewed as a sample covariance matrix, formed from $m$ samples of the random vector ${\bf T}^{1/2}_{n}{\bf X}_{\cdot 1}$ with ${\bf X}_{\cdot 1}$ denoting the first column of ${\bf X}_{n}$, which has ${\bf T}_{n}$ for its population covariance matrix. When $n$ and $m$ are both large and on the same order of magnitude, ${\bf B}_{n}$ will not be near ${\bf T}_{n}$, due to an insufficient number of samples required for such a large dimensional random vector. However, there exist results on the eigenvalues of ${\bf B}_{n}$. They are limit theorems as $n\to\infty$ with $m=m(n)$ and $c_{n}\equiv n/m\to c$, which provide information on the eigenvalues of ${\bf T}_{n}$. One result is on the empirical distribution function (e.d.f.), $F^{B_{n}}$, of the eigenvalues of ${\bf B}_{n}$, which throughout the paper, is defined for any Hermitian $n\times n$ matrix ${\bf A}$ as

The limit theorem is expressed in terms of the Stieltjes transform of the limiting e.d.f. of the $F^{B_{n}}$’s, where for any distribution function (d.f.) $G$ its Stieltjes transform, $m_{G}$, is defined to be

There exists a one-to-one correspondence between the distribution functions (d.f.’s) and their Stieltjes transforms, due to the inversion formula

It is more convenient to work with the eigenvalues of the $m\times m$ matrix $(1/m){\bf X}_{n}^{{}^{\prime}}{\bf T}_{n}{\bf X}_{n}$, whose eigenvalues differ from those of ${\bf B}_{n}$ by $|n-m|$ zero eigenvalues. Indeed, with $I_{A}$ denoting the indicator function on the set $A$ we have the exact relationship

In simple terms $S^{\prime}_{F^{c,H}}$ is comprised of the range of values where $x_{c,H}$ is increasing.

Another result which will be needed later is the following.

Suppose each $m$ contained in the interval $[m_{1},m_{2}]$ satisfies (1) and (2) of Lemma VI.1, and $\frac{d}{dm}x_{c,H}(m_{i})\geq 0$ for $i=1,2$. Then $\frac{d}{dm}x_{c,H}(m)>0$ for all $m\in(m_{1},m_{2})$.

Limiting eigenvalue mass at zero is also derived in . It is shown that

VI-B Support of eigenvalues

Since the convergence in distribution of $F^{B_{n}}$ only addresses how proportions of eigenvalues behave, understanding the possible appearance or non-appearance of eigenvalues in $S^{\prime}_{F^{c,H}}$ requires further work.

The question of the behavior of the largest and smallest eigenvalues when ${\bf T}_{n}={\bf I}$ has been answered by Yin, Bai, and Krishnaiah in , and Bai and Yin in , respectively, under the additional assumption $\text{\pe E}|{\bf X}_{1\,1}|^{4}<\infty$: the largest eigenvalue and $\min(n,m)^{\text{th}}$ largest eigenvalue of $(1/m){\bf X}_{n}{\bf X}_{n}^{*}$ converge a.s. to $(1+\sqrt{c})^{2}$ and $(1-\sqrt{c})^{2}$ respectively, matching the support, $[(1-\sqrt{c})^{2},(1+\sqrt{c})^{2}]$ of $F$ on $(0,\infty)$. More on $F$ when ${\bf T}_{n}={\bf I}$ will be given later.

For general ${\bf T}_{n}$, restricted to being bounded in spectral norm, the non-appearance of eigenvalues in $S^{\prime}_{F^{c,H}}$ has been proven by Bai and Silverstein in . Moreover, the separation of eigenvalues across intervals in $S^{\prime}_{F^{c,H}}$, mirrors exactly the separation of eigenvalues over corresponding intervals in $S^{\prime}_{H}$ . The results are summarized below.

Assume additionally $\text{\pe E}|{\bf X}_{1\,1}|^{4}<\infty$ and the ${\bf T}_{n}$ are nonrandom and are bounded in spectral norm for all $n$.

Let $F^{c_{n},H_{n}}$ denote the “limiting” e.d.f. associated with $(1/m){\bf X}_{n}^{*}{\bf T}_{n}{\bf X}_{n}$, in other words, $F^{c_{n},H_{n}}$ is the d.f. having Stieltjes transform with inverse (21), where $c,H$ are replace by $c_{n},H_{n}$.

(*) Interval $[a,b]$ with $a>0$ lies in an open interval outside the support of $F^{c_{n},H_{n}}$ for all large $n$.

Then \text{\pe P}(\text{no eigenvalue of{\bf B}_{n}$appears in$[a,b]$for all large$n})=1.

For $n\times n$ Hermitian non-negative definite matrix ${\bf A}$, let $\lambda_{k}^{A}$ denote the $k^{\text{th}}$ largest eigenvalue of $A$. For notational convenience, define $\lambda_{0}^{A}=\infty$ and $\lambda_{n+1}^{A}=0$.

(i) If $c(1-H(0))>1$, then $x_{0}$, the smallest value in the support of $F^{c,H}$, is positive, and with probability 1, $\lambda_{m}^{B_{n}}\to x^{0}$ as $n\to\infty$.

(ii) If $c(1-H(0))\leq 1$, or $c(1-H(0))>1$ but $[a,b]$ is not contained in $[0,x_{0}]$, then $m_{F^{c,H}}(b)<0$, and for all $n$ large there is an index $i_{n}$ for which

Then \text{\pe P}(\lambda_{i_{n}}^{B_{n}}>b\text{ and\lambda_{i_{n}+1}^{B_{n}}$for all large$n})=1.

The behavior of the extreme eigenvalues of $(1/m){\bf X}_{n}{\bf X}_{n}^{*}$ leads to the following corollary of Theorem VI.1.

If $\lambda_{1}^{T_{n}}$ converges to the largest number in the support of $H$, then $\lambda_{1}^{B_{n}}$ converges a.s to the largest number in the support of $F$. If $\lambda_{n}^{T_{n}}$ converges to the smallest number in the support of $H$, then $c\leq 1$ ($c>1$) implies $\lambda_{n}^{B_{n}}$ ($\lambda_{n}^{(1/m)X_{n}^{*}T_{n}X_{n}}$) converges a.s. to the smallest number in the support of $F$ ($F^{c,H}$).

In Theorem VI.1, Case (i) applies when $n>m$, whereby the rank of ${\bf B}_{n}$ would be at most $m$, the conclusion asserting, that with probability 1, for all $n$ large, the rank is equal to $m$. From Lemma VI.1, Case (ii) of Theorem VI.1 covers all intervals in $S^{\prime}_{F^{c,H}}$ on $(0,\infty)$ resulting from intervals on $(-\infty,0)$ where $x_{c,H}$ is increasing. For all $n$ large $x_{c_{n},H_{n}}$ is increasing on $[m_{F^{c_{n},H_{n}}}(a),m_{F^{c_{n},H_{n}}}(b)]$, which, from inspecting the vertical asymptotes of $x_{c_{n},H_{n}}$ and Lemma VI.1, must be due to the existence of $\lambda_{i_{n}}^{T_{n}}$, $\lambda_{i_{n}+1}^{T_{n}}$ satisfying (23).

Theorem VI.1 easily extends to random ${\bf T}_{n}$, independent of $\{{\bf X}_{ij}:i,j\geq 1\}$ with the aid of Tonelli’s Theorem [42, pp. 234], provided the condition (*) on $[a,b]$ is strengthened to:

(**) With probability 1 for all $n$ large $[a,b]$ (nonrandom) lies in an open interval outside the support of $F^{c_{n},H_{n}}$.

Indeed, let $T$ denote the probability space generating $\{T_{n}\}$, $X$ the probability space generating $\{X_{ij}:i,j\geq 1\}$. Let their respective measures be denoted by $\text{\pe P}_{T}$,$\text{\pe P}_{X}$, the product measure on $T\times X$ by $\text{\pe P}_{T\times X}$. Consider, for example in case (ii), we define

Let $t\in T$ be an element of the event defined in (**). Then by Theorem VI.1 $I_{A}((t,x))=1$ for all $x$ contained in a subset of $X$ having probability 1. Therefore, by Tonelli’s theorem

Consider now case (ii) of Theorem VI.1 in terms of the corresponding interval outside the support of $H$ and the $H_{n}$’s. By Lemma VI.1 and condition (*), we have the existence of an $\epsilon>0$ such that $0\notin[m_{F^{c,H}}(a)-\epsilon,m_{F^{c,H}}(b)+\epsilon]$, and for all $n$ large

Let $t_{a}=-1/m_{F^{c,H}}(a)$, $t_{b}=-1/m_{F^{c,H}}(b)$. Then by Lemma VI.1 we have the existence of an $\epsilon^{\prime}>0$ for which $t_{a}-\epsilon^{\prime}>0$ and $[t_{a}-\epsilon^{\prime},t_{b}+\epsilon^{\prime}]\subset S^{\prime}_{H_{n}}$ for all $n$ large. Moreover, by (24) we have for all $n$ large

Necessarily, $\lambda_{i_{n}}^{T_{n}}>t_{b}+\epsilon^{\prime}$ and $\lambda_{i_{n}+1}^{T_{n}}<t_{a}-\epsilon^{\prime}$.

Notice the steps can be completely reversed, that is, beginning with an interval $[t_{a},t_{b}]$, with $t_{a}>0$, lying in an open interval in $S^{\prime}_{H_{n}}$ for all $n$ large and satisfying (25) for some $\epsilon^{\prime}>0$, will yield $[a,b]$, with $a=x_{c,H}(-1/t_{a})$, $b=x_{c,H}(-1/t_{b})$, satisfying condition (*). Case (ii) applies, since $[a,b]$ is within the range of $x_{c,H}(m)$ for $m<0$. If $c(1-H(0))>1$, then we would have $a>x_{0}$.

VI-C Behavior of spiked eigenvalues

Suppose now the ${\bf T}_{n}$’s are altered, where a finite number of eigenvalues are interspersed between the previously adjacent eigenvalues $\lambda^{T_{n}}_{i_{n}+1}$ and $\lambda^{T_{n}}_{i_{n}}$. It is clear that the limiting $F$ will remain unchanged. However, the graph of $x_{c_{n},H_{n}}$ on $(-1/\lambda^{T_{n}}_{i_{n}+1},-1/\lambda^{T_{n}}_{i_{n}})$ will now contain vertical asymptotes. If the graph remains increasing on two intervals for all $n$ large, each one between successive asymptotes, then because of Theorem VI.1, with probability one, eigenvalues of the new ${\bf B}_{n}$ will appear in $S^{\prime}_{F^{c,H}}$ for all $n$ large.

Theorem VI.3 below shows this will happen when a “sprinkled”, or “spiked” eigenvalue lies in $(t_{a},t_{b})$. Theorem VI.4 provides a converse, in the sense that any isolated eigenvalue of $B_{n}$ must be due to a spiked eigenvalue, the absence of which corresponds to case (ii) of Theorem VI.1.

Theorem VI.3, below, allows the number of spiked eigenvalues to grow with $n$, provided it remains $o(n)$.

Assume in additon to the assumptions in Theorem VI.1 on the ${\bf X}_{ij}$ and ${\bf T}_{n}$:

for $t=t_{a},t_{b}$. (c) $t^{\prime}\in(t_{a},t_{b})$.

For $m\in[-1/t_{a},-1/t_{b}]\cap\{-1/t^{\prime}\}^{c}$, we have

By considering continuity points of $H$ in $(\alpha,\beta)$ we see that $H$ is constant on this interval, and consequently, this interval is also contained in $S^{\prime}_{H}$.

Because of (b) we have $\frac{d}{dm}x_{c,H}(m)\geq 0$ for $m=-1/t_{a},-1/t_{b}$ (recall (24),(25)).

By Lemma VI.2 we therefore have $\frac{d}{dm}x_{c,H}(m)>0$ for all $m\in(-1/t_{a},-1/t_{b})$. Thus we can find $[\underline{t}_{a},\underline{t}_{b}]\subset[t_{a},t_{b}]$ and $\delta>0$, such that $t^{\prime}\in(\underline{t}_{a},\underline{t}_{b})$ and for all $n$ large $\frac{d}{dm}x_{c_{n},\hat{H}_{n}}(m)\geq\delta$ for all $m\in[-1/\underline{t}_{a},-1/\underline{t}_{b}]$.

It follows that for any positive $\epsilon$ sufficiently small, there exist positive $\delta^{\prime}$ with $\delta^{\prime}\leq\epsilon$, such that, for all $n$ large, both $[-1/t^{\prime}-\epsilon-\delta^{\prime},-1/t^{\prime}-\epsilon]$, and $[-1/t^{\prime}+\epsilon,-1/t^{\prime}+\epsilon+\delta^{\prime}]$:

1) are contained in $[-1/\underline{t}_{a},-/\underline{t}_{b}]$, and

2) $\frac{d}{dm}x_{c_{n},H_{n}}(m)>0$ for all $m$ contained in these two intervals.

Therefore, by Lemma VI.1, for all $n$ large, $[x_{c_{n},H_{n}}(-1/t^{\prime}-\epsilon-\delta^{\prime}),x_{c_{n},H_{n}}(-1/t^{\prime}-\epsilon)]$ and $[x_{c_{n},H_{n}}(-1/t^{\prime}+\epsilon),x_{c_{n},H_{n}}(-1/t^{\prime}+\epsilon+\delta^{\prime})]$ lie outside the support of $F^{c_{n},H_{n}}$. Let $a_{L}=x_{c,H}(-1/t^{\prime}-\epsilon-\tfrac{2}{3}\delta^{\prime})$, $b_{L}=x_{c,H}(-1/t^{\prime}-\epsilon-\tfrac{1}{3}\delta^{\prime})$, $a_{R}=x_{c,H}(-1/t^{\prime}+\epsilon+\tfrac{1}{3}\delta^{\prime})$, and $b_{R}=x_{c,H}(-1/t^{\prime}+\epsilon+\tfrac{2}{3}\delta^{\prime})$. Then for all $n$ large

It follows then that $[a_{L},b_{L}]$, $[a_{R},b_{R}]$ each lie in an open interval in $S^{\prime}_{F^{c_{n},H_{n}}}$ for all $n$ large. Moreover $m_{F^{c,H}}(b_{R})<0$. Therefore, case (ii) of Theorem VI.1 applies and we have

Therefore, considering a countable collection of $\epsilon$’s converging to zero, we conclude that, with probability 1

By Lemma VI.1, $[t_{a},t_{b}]\in S^{\prime}_{H}$, and for a suitable positive $\epsilon$, $x_{c,H}$ is increasing on $[m_{c,H}(a)-\epsilon,m_{c,H}(b)+\epsilon]$, which does not contain 0.

Therefore $t_{a}<t^{\prime}<t_{b}$. If $c(1-H(0))>1$, that is, case (i) of Theorem VI.1 holds, then $a>x_{0}$, since $x_{0}$ is the almost sure limit of $\lambda_{m}^{B_{n}}$ so $\lambda^{\prime}$ cannot be smaller than it, and necessarily $x_{0}\in S_{F}$. Therefore $m_{c,H}(b)<0$, so that $0<t_{a}$.

It must be the case that only $o(n)$ eigenvalues of $t_{n}$ lie in $[t_{a},t_{b}]$, since otherwise $[t_{a},t_{b}]$ would not be outside the support of $H$. We have then $\hat{H}_{n}\rightarrow_{{\f D}}H$ as $n\to\infty$, so from the dominated convergence theorem we have $\frac{d}{dm}x_{c_{n},\hat{H}_{n}}(m)\to\frac{d}{dm}x_{c,H}(m)$ for all $m\in[m_{c,H}(a)-\epsilon,m_{c,H}(b)+\epsilon]$, implying for all $n$ large $\frac{d}{dm}x_{c_{n},\hat{H}_{n}}(m)>0$ for all $m\in[m_{c,H}(a),m_{c,H}(b)]$. Therefore (b) is true

VI-D Behavior of extreme eigenvalues

Consider now $t^{\prime}$ lying on either side of the support of $\hat{H}$. Let $\hat{\lambda}_{n}^{\min}$ and $\hat{\lambda}_{n}^{\max}$ denote, respectively, the smallest and largest numbers in the support of $\hat{H}_{n}$ Notice that $g_{n}(t)\equiv c_{n}\int\frac{\lambda^{2}}{(\lambda-t)^{2}}d\hat{H}_{n}(t)$ is decreasing for $t>\hat{\lambda}_{n}^{\max}$, and if $\hat{\lambda}_{n}^{\min}>0$, $g_{n}$ is increasing on $(0,\hat{\lambda}_{n}^{\min})$.

Therefore, if for all $n$ large, $t^{\prime}>\hat{\lambda}_{n}^{\max}$, it is necessary and sufficient to find a $t_{a}\in(\hat{\lambda}_{n}^{\max},t^{\prime})$ for which $g(t_{a})\leq 1$ in order for (26) to hold. Similarly, if for all $n$ large $t^{\prime}\in(0,\hat{\lambda}_{n}^{\min})$, then it is necessary and sufficient to find a $t_{b}\in(0,t^{\prime})$ for which $g_{n}(t_{b})\leq 1$ in order for (26) to hold. Notice if $c(1-H(0))>1$ then $g_{n}(t)>1$ for all $t\leq\hat{\lambda}_{n}^{\min}$ and all $n$ large.

Let for d.f. $G$ with bounded support, $\lambda_{G}^{\max}$ denote the largest number in $S_{G}$. If there is a $\tau>\lambda_{H}^{\max}$ for which $g(\tau)=c\int\frac{\lambda^{2}}{(\lambda-t)^{2}}dH(t)=1$, and if $\limsup_{n}\hat{\lambda}_{n}^{\max}<\tau$, then $\tau$ can be used as a threshold for $t^{\prime}\in(\limsup_{n}\hat{\lambda}_{n}^{\max},\infty)$. Indeed, by the dominated convergence theorem, $\lim_{n\to\infty}g_{n}(t^{\prime})=g(t^{\prime})$. Therefore, if $t^{\prime}>\tau$, conditions (b) and (c) of Theorem VI.3 hold, with $t_{a}=\tau$, and $t_{b}$ any arbitrarily large number.

Similar results can be obtained for the interval to the left of $S_{F}$

As in Theorem VI.1 Tonelli’s Theorem can easily be applied to establish equivalent results when ${\bf T}_{n}$’s are random and independent of ${\bf X}$.

VI-E The eigenvalues of the multivariate F matrix

Let ${\bf Y}_{ij}$ be another collection of i.i.d. random variables (not necessarily having the same distribution as the ${\bf X}_{ij}$’s), with $\text{\pe E}{\bf Y}_{1\,1}=0$, $\text{\pe E}|{\bf Y}_{1\,1}|=1$, $\text{\pe E}|{\bf Y}_{1\,1}|^{4}<\infty$, and independent of the ${\bf X}_{ij}$’s.

We form the $n\times N$ matrix ${\bf Y}_{n}=(Y_{ij})$, $i=1,2,\ldots,n$, $j=1,2,\ldots,N$ with $N=N(n)$, $n<N$, and $c^{1}_{n}\equiv n/N\to c_{1}\in(0,1)$ as $n\to\infty$.

Let now ${\bf T}_{n}=((1/N){\bf Y}_{n}{\bf Y}_{n}^{*})^{-1}$, whenever the inverse exists.

From Bai and Yin’s work we know that with probability 1, for all $n$ large, ${\bf T}_{n}$ exists with $\lambda_{1}^{T_{n}}\to(1-\sqrt{c_{1}})^{-2}$. Whenever $\lambda_{n}^{(1/N){\bf Y}_{n}{\bf Y}_{n}^{*}}=0$ define ${\bf T}_{n}$ to be ${\bf I}$.

The matrix ${\bf T}_{n}(1/N){\bf X}_{n}{\bf X}_{n}^{*}$, typically called a multivariate $F$ matrix, has the same eigenvalues as ${\bf B}_{n}$. Its limiting e.d.f. has density on $(0,\infty)$ given by

When $c\in(0,1]$, there is no mass at , whereas for $c>1$ $F$ has mass $(1-(1/c))$ at .

We are interested in the effect on spikes on the right side of the support of the $H_{n}$.

Because of the corollary to Theorem VI.1, we know $\lambda_{1}^{B_{n}}\to b_{2}$ a.s. as $n\to\infty$. We proceed in computing the function

We will see that it is unnecessary to compute the limiting e.d.f. of ${\bf T}_{n}$. It suffices to know the limiting Stieltjes transform of $F^{(1/N)Y_{n}Y_{n}^{*}}$.

Let $H_{1}$ denote the limiting e.d.f. of $F^{(1/N)Y_{n}Y_{n}^{*}}$. We have

We use (21) to find $m_{F^{c_{1},I_{[1,\infty)}}}$:

(the sign depending on with branch of the square root is taken).

As mentioned earlier the support of $H_{1}$ is $[(1-\sqrt{c_{1}})^{2},(1+\sqrt{c_{1}})^{2}]$. We need $g(t)$ for $t>(1-\sqrt{c_{1}})^{-2}$, so we need $m_{H_{1}}(x)$ for $x\in(0,(1-\sqrt{c_{1}})^{2})$.

Since $0\in S^{\prime}_{H_{1}}$, $m_{H_{1}}(0)$ exists and is real, which dictates what sign is taken on $(0,(1-\sqrt{c_{1}})^{2})$. We find that, on this interval

and using the fact that the discriminant equals $x^{2}-2x(1+c_{1})+(1-c_{1})^{2}$,

We therefore find that for $t>(1-\sqrt{c_{1}})^{-2}$

We see that the equation $g(t)=1$ leads to the following quadratic equation in $t$:

The positive sign in front of the square root being correct due to

Reducing further we find the threshold, $\tau$, to be

We now compute the right hand side of (26). We have for $t^{\prime}\geq\tau$

A straightforward (but tedious) calculation will yield $\lambda(\tau)=b_{2}$.

Using the results from the previous section, we have proved the following:

Assume in addition to the assumptions in Theorem VI.1 on the ${\bf X}_{ij}$

(a) the ${\bf T}_{n}$, possibly random, are independent of the ${\bf X}_{ij}$, with $F^{T_{n}}\rightarrow_{{\f D}}H$, a.s. as $n\to\infty$, $H$ being the limiting e.d.f. of $F^{((1/N)Y_{n}Y_{n}^{*})^{-1}}$, defined above.

(c) With $\hat{\lambda}_{n}^{\max}$ defined to be the largest number in the support of $\hat{H}_{n}$, with probability one, $\limsup_{n}\hat{\lambda}_{n}^{\max}<\tau$ the threshold defined in (30).

References