Performance of Statistical Tests for Single Source Detection using Random Matrix Theory

I Introduction

The detection of a source by a sensor array is at the heart of many wireless applications. It is of particular interest in the realm of cognitive radio where a multi-sensor cognitive device (or a collaborative network The collaborative network corresponds to multiple base stations connected, in a wireless or wired manner, to form a virtual antenna system.) needs to discover or sense by itself the surrounding environment. This allows the cognitive device to make relevant choices in terms of information to feed back, bandwidth to occupy or transmission power to use. When the cognitive device is switched on, its prior knowledge (on the noise variance for example) is very limited and can rarely be estimated prior to the reception of data. This unfortunately rules out classical techniques based on energy detection and requires new sophisticated techniques exploiting the space or spectrum dimension.

In our setting, the aim of the multi-sensor cognitive detection phase is to construct and analyze tests associated with the following hypothesis testing problem:

The standard case where the propagation channel and the noise variance are known has been thoroughly studied in the literature in the Single Input Single Output case and Multi-Input Multi-Ouput case. In this simple context, the most natural approach to detect the presence of source $s(n)$ is the well-known Neyman-Pearson (NP) procedure which consists in rejecting the null hypothesis when the observed likelihood ratio lies above a certain threshold . Traditionally, the value of the threshold is set in such a way that the Probability of False Alarm (PFA) is no larger than a predefined level $\alpha\in(0,1)$. Recall that the PFA (resp. the miss probability) of a test is defined as the probability that the receiver decides hypothesis $H_{1}$ (resp. $H_{0}$) when the true hypothesis is $H_{0}$ (resp. $H_{1}$). The NP test is known to be uniformly most powerful i.e., for any level $\alpha\in(0,1)$, the NP test has the minimum achievable miss probability (or equivalently the maximum achievable power) among all tests of level $\alpha$. In this paper, we assume on the opposite that:

the noise variance $\sigma^{2}$ is unknown,

In this context, probability density functions of the observations ${\boldsymbol{y}}(n)$ under both $H_{0}$ and $H_{1}$ are unknown, and the classical NP approach can no longer be employed. As a consequence, the construction of relevant tests for (1) together with the analysis fo their perfomances is a crucial issue. The classical approach followed in this paper consists in replacing the unknown parameters by their maximum likelihood estimates. This leads to the so-called Generalized Likelihood Ratio (GLR). The Generalized Likelihood Ratio Test (GLRT), which rejects the null hypothesis for large values of the GLR, easily reduces to the statistics given by the ratio of the largest eigenvalue of the sampled covariance matrix with its normalized trace, cf. . Nearby statistics , with good practical properties, have also been developed, but would not yield a different (asymptotic) error exponent analysis.

In this paper, we analyze the performance of the GLRT in the asymptotic regime where the number $K$ of sensors and the number $N$ of observations per sensor are large but have the same order of magnitude. This assumption is relevant in many applications, among which cognitive radio for instance, and casts the problem into a large random matrix framework.

Large random matrix theory has already been applied to signal detection (see also ), and recently to hypothesis testing . In this article, the focus is mainly devoted to the study of the largest eigenvalue of the sampled covariance matrix, whose behaviour changes under $H_{0}$ or $H_{1}$. The fluctuations of the largest eigenvalue under $H_{0}$ have been described by Johnstone by means of the celebrated Tracy-Widom distribution, and are used to study the threshold and the $p$-value of the GLRT.

In order to characterize the performance of the test, a natural approach would have been to evaluate the Receiver Operating Characteristic (ROC) curve of the GLRT, that is to plot the power of the test versus a given level of confidence. Unfortunately, the ROC curve does not admit any simple closed-form expression for a finite number of sensors and snapshots. As the miss probability of the GLRT goes exponentially fast to zero, the performance of the GLRT is analyzed via the computation of its error exponent, which caracterizes the speed of decrease to zero. Its computation relies on the study of the large deviations of the largest eigenvalue of ’spiked’ sampled covariance matrix. By ’spiked’ we refer to the case where the eigenvalue converges outside the bulk of the limiting spectral distribution, which precisely happens under hypothesis $H_{1}$. We build upon to establish the large deviation principle, and provide a closed-form expression for the rate function.

We also introduce the error exponent curve, and plot the error exponent of the power of the test versus the error exponent for a given level of confidence. The error exponent curve can be interpreted as an asymptotic version of the ROC curve in a $\log$-$\log$ scale and enables us to establish that the GLRT outperforms another test based on the condition number, and proposed by in the context of cognitive radio.

Notice that the results provided here (determination of the threshold of the GLRT test and the computation of the error exponents) would still hold within the setting of real Gaussian random variables instead of complex ones, with minor modifications Details are provided in Remarks 4 and 9..

Section II introduces the GLRT. The value of the threshold, which completes the definition of the GLRT, is established in Section II-B. As the latter threshold has no simple closed-form expression and as its practical evaluation is difficult, we introduce in Section II-C an asymptotic framework where it is assumed that both the number of sensors $K$ and the number $N$ of available snapshots go to infinity at the same rate. This assumption is valid for instance in cognitive radio contexts and yields a very simple evaluation of the threshold, which is important in real-time applications.

In Section III, we recall several results of large random matrix theory, among which the asymptotic fluctuations of the largest eigenvalue of a sample covariance matrix, and the limit of the largest eigenvalue of a spiked model.

These results are used in Section IV where an approximate threshold value is derived, which leads to the same PFA as the optimal one in the asymptotic regime. This analysis yields a relevant practical method to approximate the $p$-values associated with the GLRT.

Section V is devoted to the performance analysis of the GLRT. We compute the error exponent of the GLRT, derive its expression in closed-form by establishing a Large Deviation Principle for the test statistic $T_{N}$ Note that in recent papers , the fluctuations of the test statistics under $H_{1}$, based on large random matrix techniques, have also been used to approximate the power of the test. We believe that the performance analysis based on the error exponent approach, although more involved, has a wider range of validity., and describe the error exponent curve.

Section VI introduces the test based on the condition number, that is the statistics given by the ratio between the largest eigenvalue and the smallest eigenvalue of the sampled covariance matrix. We provide the error exponent curve associated with this test and prove that the latter is outperformed by the GLRT.

Section VII provides further numerical illustrations and conclusions are drawn in Section VIII.

Mathematical details are provided in the Appendix. In particular, a full rigorous proof of a large deviation principle is provided in Appendix A, while a more informal proof of a nearby large deviation principle, maybe more accessible to the non-specialist, is provided in Appendix B.

II Generalized Likelihood Ratio Test

In this section, we derive the Generalized Likelihood Ratio Test (section II-A) and compute the associated threshold and $p$-value (section II-B). This exact computation raises some computational issues, which are circumvented by the introduction of a relevant asymptotic framework, well-suited for mathematical analysis (Section II-C).

Denote by $N$ the number of observed samples and recall that:

and respectively, by $p_{0}({\bf Y};\sigma^{2})$ and $p_{1}({\bf Y};{\boldsymbol{h}},\sigma^{2})$ the likelihood functions of the observation matrix ${\bf Y}$ indexed by the unknown parameters ${\boldsymbol{h}}$ and $\sigma^{2}$ under hypotheses $H_{0}$ and $H_{1}$.

As ${\bf Y}$ is a $K\times N$ matrix whose columns are i.i.d. Gaussian vectors with covariance matrix ${\bf\Sigma}$ defined by:

In the case where parameters $\boldsymbol{h}$ and $\sigma^{2}$ are available, the celebrated Neyman-Pearson procedure yields a uniformly most powerful test, given by the likelihood ratio statistics $\frac{p_{1}({\bf Y};{\boldsymbol{h}},\sigma^{2})}{p_{0}({\bf Y};\sigma^{2})}$.

However, in the case where $\boldsymbol{h}$ and $\sigma^{2}$ are unknown, which is the problem addressed here, no simple procedure garantees a uniformly most powerful test, and a classical approach consists in computing the GLR:

In the following proposition, which follows after straightforward computations from and , we derive the closed form expression of the GLR $L_{N}$. Denote by $\lambda_{1}>\lambda_{2}>\dots>\lambda_{K}\geq 0$ the ordered eigenvalues of $\hat{\bf R}$ (all distincts with probability one).

where $C=\left(1-\frac{1}{K}\right)^{(1-K)N}$.

By Proposition 1, $L_{N}=\phi_{N,K}(T_{N})$ where $\phi_{N,K}:x\mapsto Cx^{-N}\left(1-\frac{x}{K}\right)^{N(1-K)}$. The GLRT rejects the null hypothesis when inequality $L_{N}>\xi_{N}$ holds. As $T_{N}\in(1,K)$ with probability one and as $\phi_{N,K}$ is increasing on this interval, the latter inequality is equivalent to $T_{N}>\phi_{N,K}^{-1}(\xi_{N})$. Otherwise stated, the GLRT reduces to the test which rejects the null hypothesis for large values of $T_{N}$:

where $\gamma_{N}=\phi_{N,K}^{-1}(\xi_{N})$ is a certain threshold which is such that the PFA does not exceed a given level $\alpha$. In the sequel, we will therefore focus on the test statistics $T_{N}$.

There exist several variants of the above statistics , which merely consist in replacing the normalized trace with a more involved estimate of the noise variance. Although very important from a practical point of view, these variants have no impact on the (asymptotic) error exponent analysis. Therefore, we restrict our analysis to the traditional GLRT for the sake of simplicity.

II-B Exact threshold and pp-values

where $p_{N}(t)$ represents the complementary c.d.f. of the statistics $T_{N}$ under the null hypothesis:

Note that $p_{N}(t)$ is continuous and decreasing from 1 to 0 on $t\in[0,\infty)$, so that the threshold $p_{N}^{-1}(\alpha)$ in (8) is always well defined. When the threshold is fixed to $\gamma_{N}=p_{N}^{-1}(\alpha)$, the GLRT rejects the null hypothesis when $T_{N}>p_{N}^{-1}(\alpha)$ or equivalently, when $p_{N}(T_{N})<\alpha$. It is usually convenient to rewrite the GLRT under the following form:

The statistics $p_{N}(T_{N})$ represents the significance probability or $p$-value of the test. The null hypothesis is rejected when the $p$-value $p_{N}(T_{N})$ is below the level $\alpha$. In practice, the computation of the $p$-value associated with one experiment is of prime importance. Indeed, the $p$-value not only allows to accept/reject an hypothesis by (10), but it furthermore reflects how strongly the data contradicts the null hypothesis .

In order to evaluate $p$-values, we derive in the sequel the exact expression of the complementary c.d.f. $p_{N}$. The crucial point is that $T_{N}$ is a function of the eigenvalues $\lambda_{1},\dots,\lambda_{K}$ of the sampled covariance matrix $\hat{\bf R}$. We have

where for each $t$, the domain of integration $\Delta_{t}$ is defined by:

and $p_{K,N}^{0}$ is the joint probability density function (p.d.f.) of the ordered eigenvalues of ${\bf\hat{R}}$ under $H_{0}$ given by:

where ${\boldsymbol{1}}_{(x_{1}\geq\dots\geq x_{K}\geq 0)}$ stands for the indicator function of the set $\{(x_{1}\dots x_{K})\>:\>x_{1}\geq\dots\geq x_{K}\geq 0\}$ and where $Z_{K,N}^{0}$ is the normalization constant (see for instance , [29, Chapter 4]).

For each $t$, the computation of $p_{N}(t)$ requires the numerical evaluation of a non-trivial integral. Despite the fact that powerful numerical methods, based on representations of such integrals with hypergeometric functions , are available (see for instance , ), an on line computation, requested in a number of real-time applications, may be out of reach.

Instead, tables of the function $p_{N}$ should be computed off line i.e., prior to the experiment. As both the dimensions $K$ and $N$ may be subject to frequent changes In cognitive radio applications for instance, the number of users $K$ which are connected to the network is frequently varying., all possible tables of the function $p_{N}$ should be available at the detector’s side, for all possible values of the couple $(N,K)$. This both requires substantial computations and considerable memory space. In what follows, we propose a way to overcome this issue.

In the sequel, we study the asymptotic behaviour of the complementary c.d.f. $p_{N}$ when both the number of sensors $K$ and the number of snapshots $N$ go to infinity at the same rate. This analysis leads to simpler testing procedure.

II-C Asymptotic framework

We propose to analyze the asymptotic behaviour of the complementary c.d.f. $p_{N}$ as the number of observations goes to infinity. More precisely, we consider the case where both the number $K$ of sensors and the number $N$ of snapshots go to infinity at the same speed, as assumed below

This asymptotic regime is relevant in cases where the sensing system must be able to perform source detection in a moderate amount of time i.e., the number $K$ of sensors and the number $N$ of samples being of the same order. This is in particular the case in cognitive radio applications (see for instance ). Very often, the number of sensors is lower than the number of snapshots, hence the ratio $c$ lower than 1.

In the sequel, we will simply denote $N,K\to\infty$ to refer to the asymptotic regime (13).

The results related to the GLRT presented in Sections IV and V remain true for $c\geq 1$; in the case of the test based on the condition number and presented in Section VI, extra-work is needed to handle the fact that the lowest eigenvalue converges to zero, which happens if $c\geq 1$.

III Large random matrices - Largest eigenvalue - Behaviour of the GLR statistics

In this section, we recall a few facts on large random matrices as the dimensions $N,K$ go to infinity. We focus on the behaviour of the eigenvalues of $\bf\hat{R}$ which differs whether hypothesis $H_{0}$ holds (Section III-A) or $H_{1}$ holds (Section III-B).

As the column vectors of ${\bf Y}$ are i.i.d. complex Gaussian with covariance matrix ${\bf\Sigma}$ given by (2), the probability density of ${\bf\hat{R}}$ is given by:

where $Z(N,K,{\bf\Sigma})$ is a normalizing constant.

As the behaviour of $T_{N}$ does not depend on $\sigma^{2}$, we assume that $\sigma^{2}=1$; in particular, ${\bf\Sigma}={\bf I}_{K}.$ Under $H_{0}$, matrix $\hat{\bf R}$ is a complex Wishart matrix and it is well-known (see for instance ) that the Jacobian of the transformation between the entries of the matrix and the eigenvalues/angles is given by the Vandermonde determinant $\prod_{1\leq i<j\leq K}(x_{j}-x_{i})^{2}.$ This yields the joint p.d.f. of the ordered eigenvalues (12) where the normalizing constant $Z(N,K,{\bf I}_{K})$ is denoted by $Z_{K,N}^{0}$ for simplicity.

then $\Lambda_{1}$ converges in distribution toward a standard Tracy-Widom random variable with c.d.f. $F_{TW}$ defined by:

where $q$ solves the Painlevé II differential equation:

and where Ai$(x)$ denotes the Airy function. In particular, $F_{TW}$ is continuous. The Tracy-Widom distribution was first introduced in as the asymptotic distribution of the centered and rescaled largest eigenvalue of a matrix from the Gaussian Unitary Ensemble.

Tables of the Tracy-Widom law are available for instance in , while a practical algorithm allowing to efficiently evaluate equation (18) can be found in .

In the case where the entries of matrix ${\bf Y}$ are real Gaussian random variables, the fluctuations of the largest eigenvalue are still described by a Tracy-Widom distribution whose definition slightly differs from the one given in the complex case (for details, see ).

III-B Behaviour under hypothesis H1H_{1}

In this case, the covariance matrix writes ${\bf\Sigma}=\sigma^{2}{\bf I}_{K}+{\bf h}{\bf h}^{*}$ and matrix ${\bf\hat{R}}$ follows a single spiked model. Since the behaviour of $T_{N}$ is not affected if the entries of ${\bf Y}$ are multiplied by a given constant, we find it convenient to consider the model where ${\bf\Sigma}={\bf I}_{K}+\frac{{\bf h}{\bf h}^{*}}{\sigma^{2}}$. Denote by

the signal-to-noise ratio (SNR), then matrix ${\bf\Sigma}$ admits the decomposition ${\bf\Sigma}={\bf U}{\bf D}{\bf U}^{*}$ where ${\bf U}$ is a unitary matrix and ${\bf D}={\rm diag}\left(\rho_{K},1,\ldots,1\right).$ With the same change of variables from the entries of the matrix to the eigenvalues/angles with Jacobian $\prod_{1\leq i<j\leq K}(x_{j}-x_{i})^{2},$ the p.d.f. of the ordered eigenvalues writes:

where the normalizing constant $Z(N,K,{\bf I}_{K}+{\bf hh^{*}})$ is denoted by $Z_{K,N}^{1}$ for simplicity, ${\bf X}_{K}$ is the diagonal matrix with eigenvalues $(x_{1},\dots,x_{K}),$ ${\bf B}_{K}$ is the $K\times K$ diagonal matrix with eigenvalues $(\frac{\rho_{K}}{1+\rho_{K}},0,\dots,0)$, and for any real diagonal matrices ${\bf C}_{K},{\bf D}_{K},$ the spherical integral $I_{K}({\bf C}_{K},{\bf D}_{K})$ is defined as

with $m_{K}$ the Haar measure on the unitary group of size $K$ (see [30, Chapter 3] for details).

We refer to $\rho$ as the limiting SNR. We also introduce

Under hypothesis $H_{1}$, the largest eigenvalue has the following asymptotic behaviour as $N,K$ go to infinity:

III-C Limiting behaviour of TNT_{N} under H0H_{0} and H1H_{1}

Gathering the results recalled in Sections III-A and III-B, we obtain the following:

Let Assumption 1 hold true and assume that $\rho>\sqrt{c}$, then:

IV Asymptotic threshold and pp-values

In Theorem 1 below, we take advantage of the convergence results of the largest eigenvalue of $\hat{\bf R}$ under $H_{0}$ in the asymptotic regime $N,K\to\infty$ to express the threshold and the $p$-value of interest in terms of Tracy-Widom quantiles. Recall that $\bar{F}_{TW}=1-F_{TW}$, that $c_{N}=\frac{K}{N}$, and that $b_{N}$ is given by (17).

Consider a fixed level $\alpha\in(0,1)$ and let $\gamma_{N}$ be the threshold for which the power of test (7) is maximum, i.e. $p_{N}(\gamma_{N})=\alpha$ where $p_{N}$ is defined by (11). Then:

The $p$-value $p_{N}(T_{N})$ associated with the GLRT can be approximated by:

Theorem 1 provides a simple approach to compute both the threshold and the $p$-values of the GLRT as the dimension $K$ of the observed time series and the number $N$ of snapshots are large: The threshold $\gamma_{N}$ associated with the level $\alpha$ can be approximated by the righthand side of (24). Similarly, equation (25) provides a convenient approximation for the $p$-value associated with one experiment. These approaches do not require the tedious computation of the exact complementary c.d.f. (11) and, instead, only rely on tables of the c.d.f. $F_{TW}$, which can be found for instance in along with more details on the computational aspects (note that function $F_{TW}$ does not depend on any of the problem’s characteristic, and in particular not on $c$). This is of importance in real-time applications, such as cognitive radio for instance, where the users connected to the network must quickly decide for the presence/absence of a source.

We are now in position to prove the theorem.

The mere definition of $\zeta_{N}$ implies that $\alpha=p_{N}(\gamma_{N})={\bar{F}}_{N}(\zeta_{N})$. Due to (27), $\bar{F}_{TW}(\zeta_{N})\to\alpha$. As $F_{TW}$ has a continuous inverse, the first point of the theorem is proved.

V Asymptotic analysis of the power of the test

In this section, we provide an asymptotic analysis of the power of the GLRT as $N,K\to\infty$. As the power of the test goes exponentially to zero, its error exponent is computed with the help of the large deviations associated to the largest eigenvalue of matrix $\hat{\bf R}$. The error exponent and error exponent curve are computed in Theorem 2, Section V-A; the large deviations of interest are stated in Section V-B. Finally Theorem 2 is proved in Section V-C.

The most natural approach to characterize the performance of a test is to evaluate its power or equivalently its miss probability i.e., the probability under $H_{1}$ that the receiver decides hypothesis $H_{0}$. For a given level $\alpha\in(0,1)$, the miss probability writes:

where $I_{\rho}^{+}$ is the so-called rate function associated to $T_{N}$. This observation naturally yields the following definition of the error exponent ${\mathcal{E}}_{T}$:

the existence of which is established in Theorem 2 below (as $N,K\to\infty$). Also proved is the fact that ${\mathcal{E}}_{T}$ does not depend on $\alpha$.

The error exponent ${\cal E}_{T}$ gives crucial information on the performance of the test $T_{N}$, provided that the level $\alpha$ is kept fixed when $N,K$ go to infinity. Its existence strongly relies on the study of the large deviations associated to the statistics $T_{N}$.

In practice however, one may as well take benefit from the increasing number of data not only to decrease the miss probability, but to decrease the PFA as well. As a consequence, it is of practical interest to analyze the detection performance when both the miss probability and the PFA go to zero at exponential speed. A couple $(a,b)\in(0,\infty)\times(0,\infty)$ is said to be an achievable pair of error exponents for the test $T_{N}$ if there exists a sequence of levels $\alpha_{N}$ such that, in the asymptotic regime (13),

We denote by ${\cal S}_{T}$ the set of achievable pairs of error exponents for test $T_{N}$ as $N,K\to\infty$. We refer to ${\cal S}_{T}$ as the error exponent curve of $T_{N}$.

The following notations are needed in order to describe the error exponent ${\mathcal{E}}_{T}$ and error exponent curve ${\cal S}_{T}$.

Denote by $\Delta(\,\cdot\mid A)$ the convex indicator function i.e. the function equal to zero for $x\in A$ and to infinity otherwise. For $\rho>\sqrt{c}$, define the function:

We are now in position to state the main theorem of the section:

For any fixed level $\alpha\in(0,1)$, the limit ${\cal E}_{T}$ in (30) exists as $N,K\to\infty$ and satisfies:

if $\rho>\sqrt{c}$ and ${\mathcal{E}}_{T}=0$ otherwise.

The error exponent curve of test $T_{N}$ is given by:

if $\rho>\sqrt{c}$ and ${\cal S}_{T}=\emptyset$ otherwise.

The proof of Theorem 2 heavily relies on the large deviations of $T_{N}$ and is postponed to Section V-C. Before providing the proof, it is worth making the following remarks.

At high SNR, this yields the following convenient approximation of the miss probability:

where $\psi(c)=e^{-(1+\sqrt{c})}(1+\sqrt{c})^{c-1}{c}^{-\frac{c}{2}}$.

V-B Large Deviations associated to TNT_{N}

For instance, if $A$ is a set such that $\inf_{\textrm{int}(A)}I=\inf_{\textrm{cl}(A)}I(=\inf_{A}I)$, (where $\textrm{int}(A)$ and $\textrm{cl}(A)$ respectively denote the interior and the closure of $A$), then (38) and (39) yield

As already mentioned above, all the probabilities of interest are rare events as $N,K$ go to infinity related to large deviations for $T_{N}.$ More precisely, Theorem 2 is merely a consequence of the following Lemma.

Let Assumption 1 hold true and let $N,K\to\infty$, then:

Under $H_{0},$ $T_{N}$ satisfies the LDP in the scale $N$ with good rate function $I_{0}^{+}$, which is increasing from 0 to $\infty$ on interval $[\lambda^{+},\infty)$.

For any bounded sequence $(\eta_{N})_{N\geq 0}$,

Let $x\in(\lambda^{+},\infty)$ and let $(x_{N})_{N\geq 0}$ be any real sequence which converges to $x$. If $\rho\leq\sqrt{c}$, then:

The proof of Lemma 1 is provided in Appendix A.

In Appendix A, we rather focus on the large deviations of $\lambda_{1}$ under $H_{1}$ and skip the proof of Lemma 1-(1), which is simpler and available (to some extent) in [29, Theorem 2.6.6] see also the errata sheet for the sign error in the rate function on the authors webpage.. Indeed, the proof of the LDP relies on the joint density of the eigenvalues. Under $H_{1}$, this joint density has an extra-term, the spherical integral, and is thus harder to analyze.

Lemma 1-(3) is not a mere consequence of Lemma 1-(2) as it describes the deviations of $T_{N}$ at the vicinity of a point of discontinuity of the rate function. The direct application of the LDP would provide a trivial lower bound ($-\infty$) in this case.

In the case where the entries of matrix ${\bf Y}$ are real Gaussian random variables, the results stated in Lemma 1 will still hold true with minor modifications: The rate functions will be slightly different. Indeed, the computation of the rate functions relies on the joint density of the eigenvalues, which differs whether the entries of ${\bf Y}$ are real or complex.

V-C Proof of Theorem 2

where $c_{N}=\frac{K}{N}$ converges to $c$ and where $\eta_{N}$ is a deterministic sequence such that

Taking the limit in both terms yields $I_{0}^{+}(x)\geq I_{0}^{+}(x+\epsilon)$ by Lemma 1, which contradicts the fact that $I_{0}^{+}$ is an increasing function. Now assume that $\gamma<x$. Similarly,

for a certain $\epsilon$ and for $N$ large enough. Taking the limit of both terms, we obtain $I_{0}^{+}(x)\leq I_{0}^{+}(x-\epsilon)$ which leads to the same contradiction. This proves that $\lim_{N}\gamma_{N}=x$. Recall that by definition (31),

VI Comparison with the test based on the condition number

This section is devoted to the study of the asymptotic performances of the test $U_{N}=\frac{\lambda_{1}}{\lambda_{K}}$, which is popular in cognitive radio . The main result of the section is Theorem 3, where it is proved that the test based on $T_{N}$ asymptotically outperforms the one based on $U_{N}$ in terms of error exponent curves.

A different approach which has been introduced in several papers devoted to cognitive radio contexts consists in rejecting the null hypothesis for large values of the statistics $U_{N}$ defined by:

under both hypotheses $H_{0}$ and $H_{1}$. Therefore, the statistics $U_{N}$ admits the following limits:

The test is based on the observation that the limit of $U_{N}$ under the alternative $H_{1}$ is strictly larger than the ratio $\lambda^{+}/\lambda^{-}$, at least when the SNR $\rho$ is large enough.

VI-B A few remarks related to the determination of the threshold for the test UNU_{N}

The determination of the threshold for the test $U_{N}$ relies on the asymptotic independence of $\lambda_{1}$ and $\lambda_{K}$ under $H_{0}$. As we shall prove below that test $U_{N}$ is asymptotically outperformed by test $T_{N}$, such a study, rather involved, seems beyond the scope of this article. For the sake of completeness however, we describe unformally how to set the threshold for $U_{N}$. Recall the definition of $\Lambda_{1}$ in (16) and let $\Lambda_{K}$ be defined as:

Then both $\Lambda_{1}$ and $\Lambda_{K}$ converge toward Tracy-Widom random variables. Moreover,

where $X$ and $Y$ are independent random variables, both distributed according to $F_{TW}$ Such an asymptotic independence is not formally proved yet for $\bf{\hat{R}}$ under $H_{0}$, but is likely to be true as a similar result has been established in the case of the Gaussian Unitary Ensemble ,..

As a corollary of the previous convergence, a direct application of the Delta method [27, Chapter 3] yields the following convergence in distribution:

In particular, $\xi_{N}$ is bounded as $N,K\rightarrow\infty$.

VI-C Performance analysis and comparison with the GLRT

As for $\mathbf{F}^{+}$, function $\mathbf{F}^{-}$ also admits a closed-form expression based on $\mathbf{f}$, the Stieltjes transform of Mar$\check{\textrm{c}}$enko-Pastur distribution (see Appendix C for details).

If $\lambda_{1}$ and $\lambda_{K}$ were independent random variables, the contraction principle (see e.g. ) would imply that the following functions

defined for each $t\geq 0$, are the rate functions associated with the LDP governing $\lambda_{1}/\lambda_{K}$ under hypotheses $H_{1}$ and $H_{0}$ respectively. Of course, $\lambda_{1}$ and $\lambda_{K}$ are not independent, and the contraction principle does not apply. However, a careful study of the p.d.f. $p_{K,N}^{0}$ and $p_{K,N}^{1}$ shows that $\lambda_{1}$ and $\lambda_{K}$ behave as if they were asymptotically independent, from a large deviation perspective:

Let Assumption 1 hold true and let $N,K\to\infty$, then:

Under $H_{0},$ $U_{N}$ satisfies the LDP in the scale $N$ with good rate function $\Gamma_{0}$.

Under $H_{1}$ and if $\rho>\sqrt{c}$, $U_{N}$ satisfies the LDP in the scale $N$ with good rate function $\Gamma_{\rho}.$

For any bounded sequence $(\eta_{N})_{N\geq 0}$,

Moreover, $\Gamma_{\rho}(\lambda^{+})=I_{\rho}^{+}(\lambda^{+})$.

Let $x\in(\lambda^{+},\infty)$ and let $(x_{N})_{N\geq 0}$ be any real sequence which converges to $x$. If $\rho\leq\sqrt{c}$, then:

In the context of Lemma 1, both quantities $\lambda_{1}$ and $\lambda_{K}$ deviate at the same speed, to the contrary of statistics $T_{N}$ where the denominator concentrated much faster than the largest eigenvalue $\lambda_{1}$. Nevertheless, proof of Lemma 2 is a slight extension of the proof of Lemma 1, based on the study of the joint deviations $(\lambda_{1},\lambda_{K})$, the proof of which can be performed similarly to the proof of the deviations of $\lambda_{1}$. Once the large deviations established for the couple $(\lambda_{1},\lambda_{K})$, it is a matter of routine to get the large deviations for the ratio $\lambda_{1}/\lambda_{K}$. A proof is outlined in Appendix B.

We now provide the main result of the section.

For any fixed level $\alpha\in(0,1)$ and for each $\rho$, the error exponent ${\cal E}_{U}$ exists and coincides with ${\cal E}_{T}$.

The error exponent curve of test $U_{N}$ is given by:

if $\rho>\sqrt{c}$ and ${\cal S}_{U}=\emptyset$ otherwise.

The error exponent curve ${\cal S}_{T}$ of test $T_{N}$ uniformly dominates ${\cal S}_{U}$ in the sense that for each $(a,b)\in{\cal S}_{U}$ there exits $b^{\prime}>b$ such that $(a,b^{\prime})\in{\cal S}_{T}$.

The proof of items (1) and (2) is merely bookkeeping from the proof of Theorem 2 with Lemma 2 at hand.

Let us prove item (3). The key observation lies in the following two facts:

where $(a)$ follows from the fact that $I^{-}(\lambda^{-})=0$ and by taking $u=x,v=\lambda^{-}$. Assume that inequality $(a)$ is strict. Due to the fact that $I_{\rho}^{+}$ is decreasing, the only way to decrease the value of $I^{+}_{\rho}(u)+I^{-}(v)$ under the considered constraint $\frac{u}{v}=\frac{x}{\lambda^{-}}$ is to find a couple $(u,v)$ with $u>x$, but this cannot happen because this would enforce $v>\lambda^{-}$ so that the constraint $\frac{u}{v}=\frac{x}{\lambda^{-}}$ remains fulfilled, and this would end up with $I^{-}(v)=\infty$. Necessarily, $(a)$ is an equality and (57) holds true.

Let us now give a sketch of proof for (58). Notice first that $\frac{dI^{+}_{0}}{du}\mid_{u=x}>0$ (which easily follows from the fact that $I^{+}_{0}$ is increasing and differentiable) while $\frac{dI^{-}}{dv}\mid_{v\nearrow\lambda^{-}}=0$. This equality follows from the direct computation:

where the last equality follows from the fact that $\frac{d{\bf F}^{-}}{dx}=-{\bf f}$ together with the closed-form expression for ${\bf f}$ as given in Appendix C. As previously, write:

Consider now a small perturbation $u=x-\delta$ and the related perturbation $v=\lambda^{-}-\delta^{\prime}$ so that the constraint $\frac{u}{v}=\frac{x}{\lambda^{-}}$ remains fulfilled. Due to the values of the derivatives of $I^{+}_{0}$ and $I^{-}$ at respective points $x$ and $\lambda^{-}$, the decrease of $I^{+}_{0}(x-\delta)$ will be larger than the increase of $I^{-}(\lambda^{-}-\delta^{\prime})$, and this will result in the fact that

which is the desired result, which in turn yields (58).

Theorem 3-(1) indicates that when the number of data increases, the powers of tests $T_{N}$ and $U_{N}$ both converge to one at the same exponential speed ${\cal E}_{U}={\cal E}_{T}$, provided that the level $\alpha$ is kept fixed. However, when the level goes to zero exponentially fast as a function of the number of snapshots, then the test based on $T_{N}$ outperforms $U_{N}$ in terms of error exponents: The power of $T_{N}$ converges to one faster than the power of $U_{N}$. Simulation results for $N,K$ fixed sustain this claim (cf. Figure 4). This proves that in the context of interest ($N,K\to\infty$), the GLRT approach should be prefered to the test $U_{N}$.

VII Numerical Results

In the following section, we analyze the performance of the proposed tests in various scenarios.

Figure 2 compares the error exponent of test $T_{N}$ with the optimal NP test (assuming that all the parameters are known) for various values of $c$ and $\rho$. The error exponent of the NP test can be easily obtained using Stein’s Lemma (see for instance ).

In Figure 3, we compare the Error Exponent curves of both tests $T_{N}$ and $U_{N}$. The analytic expressions provided in 2 and 3 for the Error Exponent curves have been used to plot the curves. The asymptotic comparison clearly underlines the gain of using test $T_{N}$.

VIII Conclusion

In this contribution, we have analyzed in detail the GLRT in the case where the noise variance and the channel are unknown. Unlike similar contributions, we have focused our efforts on the analysis of the error exponent by means of large random matrix theory and large deviation techniques. Closed-form expressions were obtained and enabled us to establish that the GLRT asymptotically outperforms the test based on the condition number, a fact that is supported by finite-dimension simulations. We also believe that the large deviations techniques introduced here will be of interest for the engineering community, beyond the problem addressed in this paper.

Acknowlegment

We thank Olivier Cappé for many fruitful discussions related to the GLRT.

Appendix A Proof of Lemma 1: Large deviations for TNT_{N}

The large deviations of the largest eigenvalue of large random matrices have already been investigated in various contexts, Gaussian Orthogonal Ensemble and deformed Gaussian ensembles . As mentionned in [21, Remark 1.2], the proofs of the latter can be extended to complex Wishart matrix models, that is random matrices $\bf{\hat{R}}$ under $H_{0}$ or $H_{1}$.

In both cases, the large deviations of $\lambda_{1}$ rely on a close study of the density of the eigenvalues, either given by (12) (under $H_{0}$) or by (19) for the spiked model (under $H_{1}$). The study of the spiked model, as it involves the study of the asymptotics of the spherical integral (see Lemma 3 below), is more difficult. We therefore focus on the proof of the LDP under $H_{1}$ (Lemma 1-(2)) and omit the proof of Lemma 1-(1). Once Lemma 1-(2) is proved, proving Lemma 1-(1) is a matter of bookkeeping, with the spherical integral removed at each step.

Recall that $\lambda_{1}\geq\cdots\geq\lambda_{K}$ are the ordered eigenvalues of $\hat{\bf R}$ and that $T_{N}$ is the statistics defined in (6).

For sake of simplicity and with no loss of generality as the law of $T_{N}$ does not depend on $\sigma,$ we assume all along this appendix that $\sigma^{2}=1.$ We first recall important asymptotic results for spherical integrals.

Consider a $K$-tuple $(x_{1},\cdots,x_{K})$ and denote by $\hat{\pi}_{K,{\bf x}}=\frac{1}{K-1}\sum_{i=2}^{N}\delta_{x_{2}}$ the empirical distribution associated to $(x_{2},\cdots,x_{K})$; let $d$ be a metric compatible with the topology of weak convergence of measures (for example the Dudley distance - see for instance ). A strong version of the convergence of the spherical integral in the exponential scale with speed $N$, established in can be summarized in the following Lemma:

Recall that the spherical integral $I_{K}$, defined in (20), appears in the joint density (19) of the eigenvalues under $H_{1}$. Lemma 3 provides a simple asymptotic equivalent $cJ_{\rho}(x)$ of the normalized integral $N^{-1}\log I_{K}$. Roughly speaking, this will enable us to replace $I_{K}$ by the quantity $e^{-N\times cJ_{\rho}(x)}$ when establishing the large deviations of $\lambda_{1}$, which rely on a careful study of density (19).

A-B Proof of Lemma 1-(2)

Condition (60) is technical (see for instance [44, Lemma 1.2.18]): Instead of proving the large deviation upper bound for every closed set, the exponential tightness (60), if established, enables one to restrict to the compact sets.

(Upper bound) For any $x$, for any $M$ such that $0<x<M,$

Due to the exponential tightness, it is sufficient to establish the upper bound for compact sets. As each compact can be covered by a finite number of balls, it is therefore sufficient to establish upper estimate (61) in order to establish the LD upper bound.

The fact that (62) implies the LD lower bound (39) is standard in LD and can be found in [44, Chapter 1] for instance.

As the arguments are very similar to the ones developed in , we only prove in detail the upper bound (61). Proofs of (60) and (62) are left to the reader.

Recall the notations introduced in (12) and (19) and let $x>\lambda^{+}$, $\delta>0$. Consider the following domain:

To proceed, one has to study the asymptotic behaviour of the normalizing constant:

and the rate function $I_{\rho}^{+}$ by the function $G_{\rho}$ defined by:

where, for any compactly supported probability measure $\mu$ and any real number $y$ greater than the right edge of the support of $\mu,$

The second term in (63) is easily obtained considering the fact that all the eigenvalues are less than $M$ so that for $1\leq j\leq K,$ $|x_{1}-x_{j}|\leq 2M,$ $x_{j}^{N-K}\leq M^{N-K}$ and $(UX_{K}U^{*})_{11}\leq M.$ Now, standard concentration results under $H_{0}$ yield that:

As $c_{N}\to c$ for $N,K\to\infty$, $c\mapsto\Phi(y,c,\mu)$ is continuous and $\mu\mapsto\Phi(y,c,\mu)$ is lower semi-continuous, we obtain:

By continuity in $u$ of the two involved functions, we finally get:

This concludes the proof of the upper bound in Lemma 1-(2). The proof of Lemma 1-(1) is very similar and left to the reader.

A-C Proof of Lemma 1-(3)

which enable us to properly separate $\lambda_{1}$ from the support of $\hat{\pi}_{K,\boldsymbol{\lambda}}$. Now, with the localisation indicated above, we have for $N$ large enough,

As previously, we consider the variables $y_{j}=\frac{N}{N-1}x_{j}$ for $2\leq j\leq K$ and obtain, with the help of Lemma 3:

We have already used the fact that the first term goes to zero when $N$ grows to infinity. Recall that the fluctuations of $\frac{1}{K-1}\sum_{j=2}^{K}\lambda_{j}$ are of order $\frac{1}{N}$, therefore the second term also goes to zero as we consider deviations of order $N^{-2/3}$. Now, $N^{2/3}(\lambda_{2}-(1+\sqrt{c_{N}})^{2})$ converges in distribution to the Tracy-Widom law, therefore the last term converges to $F_{\rm TW}\left(\eta-2+r(1+\sqrt{c})^{2}\right)<1.$ This concludes the proof.

Appendix B Sketch of proof for Lemma 2: Large deviations for UNU_{N}

As stated in Remark 10, we shall first study the LDP for the joint quantity $(\lambda_{1},\lambda_{K})$. The purpose here is to outline the following convergence:

which is an illustrative way, although informal All the statements, computations and approximations below can be made precise as in the proof of Lemma 1., to state the LDP for $(\lambda_{1},\lambda_{K})$ (see (40)).

We shall now perform the following approximations:

As $x_{1}\geq\alpha_{1}\geq\lambda^{+}$ and $x_{K}\leq\beta_{K}\leq\lambda^{-}$, the last integral goes to one as $K,N\to\infty$ and:

The term $\frac{2\log(x_{1}-x_{K})}{N}$ within the exponential in the integral accounts for the interraction between $\lambda_{1}$ and $\lambda_{K}$ and its contribution vanishes at the desired rate. In order to evaluate the two remaining integrals, one has to rely on Laplace’s method (see for instance ) to express the leading term of the integrals (replacing $KN^{-1}$ by $c$ below):

We have proved (informally) that the LDP holds true for $(\lambda_{1},\Lambda_{K})$ with rate function $I_{0/\rho}^{+}(x)+I^{-}(y)$. The contraction principle [44, Chap. 4] immediatly yields the LDP for the ratio $\frac{\lambda_{1}}{\lambda_{K}}$ with rate function:

which is the desired result. We provide here intuitive arguments to understand this fact.

Appendix C Closed-form expressions for functions 𝐟\mathbf{f}, 𝐅+\mathbf{F}^{+} and 𝐅−\mathbf{F}^{-}

Consider the Stieltjes transform $\mathbf{f}$ of Mar$\check{\textrm{c}}$enko-Pastur distribution:

We gather without proofs a few facts related to $\mathbf{f}$, which are part of the folklore.

where $\sqrt{z}$ stands for the principal branch of the square-root.

As a consequence, the following hold true:

Recall the definition (32) and (52) of function $\mathbf{F}^{+}$ and $\mathbf{F}^{-}$. In the following lemma, we provide closed-form formulas of interest.

Consider the case where $x\geq\lambda^{+}$. First write

It remains to plug this identity into (71) to conclude. The representation of $\mathbf{F}^{-}$ can be established similarly.

References