The Performance Analysis Based on SAR Sample Covariance Matrix

Multi-channel systems appear in several fields of application in science. In the Synthetic Aperture Radar (SAR) context, multi-channel systems may refer to different domains, as multi-polarization, multi-interferometric or multi-temporal data, or even a combination of them. Due to the inherent speckle phenomenon present in SAR images, the statistical description of the data is almost mandatory for its utilization. The complex images acquired over natural media present in general zero-mean circular Gaussian characteristics. In this case, second order statistics as the multi-channel covariance matrix fully describe the data. For practical situations however, the covariance matrix has to be estimated using a limited number of samples, and this sample covariance matrix follow the complex Wishart distribution. In this context, the eigendecomposition of the multi-channel covariance matrix has been shown in different areas of high relevance regarding the physical properties of the imaged scene. Specifically, the maximum eigenvalue of the covariance matrix has been frequently used in different applications as target or change detection, estimation of the dominant scattering mechanism in polarimetric data, moving target indication, etc. In this paper, the statistical behavior of the maximum eigenvalue derived from the eigendecomposition of the sample multi-channel covariance matrix in terms of multi-channel SAR images is simplified for SAR community. Validation is performed against simulated data and examples of estimation and detection problems using the analytical expressions are as well given.


Introduction
Multi-channel systems with random nature appear in a wide range of fields in the literature. For several of them, the central limit theorem applies and their random behavior can be modeled by Gaussian statistics, being thus statistically fully described by the first and second order moments. In the case of multi-channel SAR systems, the assumption of zero-mean multivariate complex Gaussian distribution is frequently valid for geophysical media, being thus fully described by the complex Hermitian covariance matrix. This is the statistical case treated in this work.
Since a SAR signal corresponds to the superposition of the scattered fields from all scatters inside a resolution cell, physical parameter estimation may be performed by exploiting the multi-channel covariance matrix. For that, eigendecomposition theorems have been widely used in the literature, decomposing the whole covariance matrix into elementary quantities, in order to provide a better physical interpretation of the data.
The eigendecomposition of the multi-channel SAR covariance matrix leads to useful information in a large range of SAR applications, including Ground Moving Target Indication (GMTI) [1], polarimetric SAR (PolSAR) [2], interferometric SAR (InSAR) [3], polarimetric interferometric SAR (PolInSAR) [4], change detection [5], target detection [6] and filtering [7,8]. Particularly, the covariance matrix maximum eigenvalue has been proved to be a key parameter in various areas. In [5] and [9], the change detection problem was elaborated using the maximum eigenvalue of the interferometric covariance matrix. The eigendecomposition of the SAR covariance matrix with application to SAR GMTI is used in and to calculate the probability of moving target detection in homogeneous and heterogeneous terrain. Regarding polarimetry, the eigendecomposition of the coherency (which is also a covariance) matrix is the most common tool among incoherent target decomposition theorems for the interpretation of the earth's surface. The maximum eigenvalue of the coherency matrix is also often used in polarimetry, due to the close relationship to the dominant scattering mechanism of the scattering process [2,10,11].
Although there is a wide area of applications regarding the maximum eigenvalue of the multi-channel covariance matrix, a certain lack of information in the literature concerning its statistical characterization in terms of multi-channel SAR may be verified. The statistical description can be useful to understand bias effects and to analyze estimation as well as detection problems. The estimated (sample) covariance matrix follows a Wishart distribution, which has been of interest along the years since its first derivation. In areas as communication systems, several works using such results have been published [12][13][14][15][16]. The application of the Wishart distribution in remote sensing was however considerably late introduced [17].
The cornerstone study considering the statistical description of the covariance matrix eigendecomposition in polarimetry has been carried out in [18]. The majority of the analysis in [18] was performed on the basis of numerical methods. In this paper the results of [18] is supported by addressing analytical solutions. Additionally, it is derived an exact closed form expressions for the Moment Generating Function (MGF) of the sample covariance matrix maximum eigenvalue.
The validation of the derived expressions is performed using simulated data, demonstrating its agreement with the theory. The effect of the number of samples and underlying correlation scenario of the sample covariance matrix on the bias in the estimation of the maximum eigenvalue is also investigated. Finally, examples of applications are as well given, in the fields of estimation and detection theory.
The next section reminds the basics of the statistical description of multi-channel SAR systems and of the covariance matrix eigendecomposition. It also includes the derived theorems presenting the statistical founds. Section 3 includes numerical examples validating the theorems as well as an analysis of their behavior. Sections 4 and 5 show their usage for the analysis of the estimation bias and for the elaboration of detection problems, respectively. Section 6 concludes the paper with discussions and directions for future work. The detailed statistical derivations can be found in the Appendix.

Preliminaries
The individual elements of a m-dimensional multi-channel system may be organized in a m-dimensional vector k. As the elements are assumed to follow zero mean complex Gaussian distributions, the vector k is said to follow a m-dimensional multivariate normal distribution, with zero mean and true covariance matrix Σ among the vector elements, and is represented by [19,20]. For zero mean Gaussian statistics the covariance matrix fully describes the data, playing a key role in several application fields. In practical situations however, the true covariance matrix Σ is unknown and has to be estimated by its maximum likelihood estimator (MLE), the sample covariance matrix Z = (1/n) ∑ n j=1 k j k j † , where n is the number of estimation samples and † is the transpose conjugate operator. In the SAR context, the number of independent samples is also called looks. The elements of Z follow a m dimensional complex Wishart distribution with n degrees of freedom and true covariance matrix Σ, represented by Z ∼ W C m (n, Σ) and defined as [20]: where Γ(·) is the gamma function and etr(·) is the exponential trace of a matrix. The spectral theorem from linear algebra allows the decomposition of the full rank m-dimensional covariance matrix in a set of m one-rank covariance matrices using its eigenvalues and eigenvectors. Accordingly, the decompositions of the true covariance matrix with their real non-negative eigenvalues l i and λ i , and respective complex eigenvectors e i , e ′ i , for i = {1, ...m}.

Sample Maximum Eigenvalue Statistical Description
The following Theorem III concerns the derived Moment Generating Function (MGF) of the maximum eigenvalue of the sample covariance matrix Z of a multi-channel statistical system, which is a critical step in removing the bias of the largest eigenvalue. The Theorem III, which is derived due to the previous works (Theorems I and II), is the main result of the paper. It is detailed in the appendix and will be used for the elaboration of illustrative estimation and detection problems in Sections 4 and 5.
Throughout the paper, | · | is the matrix determinant, ⟨X⟩ n = (1/n) ∑ n i=1 X i denotes the estimator of the random matrix X formed from a sample of size n.
Theorem I: Let k ∼ N C m (0, Σ) be a m-dimensional complex vector whose elements follow zero mean Gaussian distributions with associated m × m covariance matrix Σ. Let Σ have l m ≤ .... ≤ l 1 eigenvalues. Then the Cumulative Density Function (CDF) of the maximum eigenvalue λ max of the sample covariance matrix ⟨kk † ⟩ n , with the assumption m ≤ n, is given by , and γ is the incomplete Gamma function (Equation (2.42) in [21]).
Proof: [14] gave the closed expression of CDF of the largest eigenvalue of the complex Wishart matrices. Here, it is written in the form of SAR covariance matrix after small continuation, see also Appendix A.
Theorem II: Let k ∼ N C m (0, Σ) be a m-dimensional complex vector whose elements follow zero mean Gaussian distributions with associated m × m covariance matrix Σ. Let Σ have l m ≤ .... ≤ l 1 eigenvalues. Then the Probability Density Function (PDF) of the maximum eigenvalue λ max of the sample covariance matrix ⟨kk † ⟩ n , with the assumption of m ≤ n, is given by where x n−j , and Ψ(x) and S are defined in Equation (3).
Here it can be noted that when the true covariance matrix Σ is diagonal, the PDF of the largest eigenvalue reduces to Ermolaev and Rodyushkin result [23].
Theorem III: Let k ∼ N C m (0, Σ) be a m-dimensional complex vector whose elements follow zero mean Gaussian distributions with associated m × m covariance matrix Σ. Let Σ have l m ≤ .... ≤ l 1 eigenvalues. Then for any positive integers s, the sth moment of the maximum eigenvalue λ max of the sample covariance matrix ⟨kk † ⟩ n , with the assumption of m ≤ n, is given by (6) where S indicates the constant term as in Theorem I and F A (· · · ) is the hypergeometric function of several variables with a k = n l k , Here, the sum is computed over (m − 1)! permutations of π s . S m denotes the set of all m! permutations of the set S = {1, 2, ..., m}, and sgn(π s ) denotes the sign of the permutation π s : +1 if π s is an even permutation and −1 if it is odd.
Proof: See Appendix section B. In the appendix the moment of sample maximum eigenvalue is also given in more friendlier form for programming obtained by splitting the hypergeometric function of several variables into a sum of its variables.

Dependence of the Covariance Matrix Eigenvalues
Before starting the following analysis, it would be nice to illustrate the dependencies of the eigenvalues on the covariance matrix parameters, which is necessary for the further understanding of the addressed topics. For that, a two-dimensional covariance matrix is used, originated from the expected values of the outer product of the vector k = [k 1 k 2 ] T , i.e., where σ 2 1 , σ 2 2 are the variance or power of k 1 and k 2 , respectively, ρ is the absolute and ϕ the phase of their complex correlation coefficient, and T means transpose.
The eigenvalues of the true covariance matrix can be shown to be given by [5]: The behavior of l 1 and l 2 is presented in Figure 1, as a function of a normalized power ratio b and a correlation ρ. Note their mutual dependence, making clear that b and ρ directly determine the behavior of both eigenvalues.

Validation and Analysis of the Theoretical Expressions
This section aims to validate the theorems mentioned in the previous section using simulated data. The simulated data have been generated using different multi-dimensional configurations, where the correlation between channels have been generated using the well known Mahalanobis transformation [24]. The distribution of the maximum sample eigenvalue as a function of the number of samples, for a three-dimensional system. In both case, the powers are given by σ k 1 = σ k 2 = σ k 3 = 1 and correlations by ρ k 1 k 2 = 0, ρ k 1 k 3 = 0.8 and ρ k 2 k 3 = 0. When n → ∞, λ max = l 1 = 1.8.  (4) with simulations. The theoretical PDF curves clearly agree with the histograms obtained from simulated data. As expected, the PDFs become narrower with increasing n, indicating less variance around the true value of the maximum eigenvalue l 1 = 1.8. This behavior can be better seen in Figure 2(b) where the distribution of λ max as a function of the number of estimation samples n is presented. Note also that the expected value of the distributions seems to change for different n. Figure 3(a) shows the variation of the histogram mean as a function of n for a two-dimensional case with fix correlation ρ = 0.2. In Figure 3(a), the theoretical expected value of λ max has been also over plotted for comparison. The curves match well, which validates Theorem III. The same has been carried out for the second order moment, i.e., the variance of λ max , and is presented in Figure 3(b). The agreement between the theoretical variance of λ max and the variance of the simulated data can again be confirmed. Also, as already observed in Figure 2(a), the lower the number of samples n, the higher the variance of λ max . The impact of the correlation and the number of samples on the behavior of the third (skewness) and fourth (kurtosis) order moment of the maximum sample eigenvalue λ max is presented in Figure 4. Skewness is a measure of how symmetrical the distribution is with respect to its mean. Figure 4(a,b) indicates that the skewness of λ max converges to zero for increasing n and increasing ρ, expressing the tendency to a symmetrical distribution in that cases.
Kurtosis is a measure of the peakedness of the distribution. Figure 4(c,d) shows that for increasing n and ρ, the kurtosis of the λ max distribution tends to three, which corresponds to the kurtosis of the normal Gaussian distribution. Figure 5 shows the expected value (first order moment) of the maximum sample eigenvalue λ max as a function of the true eigenvalues l max for a fixed number of samples n = 3 and k = 1. The variation of the underlying correlation of the true covariance matrix ρ, which changes the values of l max , has been as well indicated in a color code. Observe that the higher the correlation is, the higher l max becomes. (d) Figure 5. Effects of the correlation between channels on the expected value of the maximum sample eigenvalue keeping n = 3 and k = 1.

Estimation Bias
The bias of the estimator of a certain parameter is defined as the difference between the expected value of the estimator and the true value of the parameter. In the present case, the bias is hence given by Figure 5 that the bias of λ max becomes very strong for low values of the correlation l max (or equivalently, low values of ρ).
The effect is emphasized in Figure 6(a,b), where the bias of λ max is presented as a function of l max and l 2 , for n = 3 and n = 16, respectively. The values of the underlying correlation ρ have been also indicated. The bias is small either when the correlation is high or the number of samples is sufficiently large. When the correlation between channels is low, the bias becomes very significant and a large number of samples is necessary to decrease the estimation bias.

Analysis of Detection Problems
Detection theory is a means to quantify the ability of a procedure to detect a parameter (or signal) immersed in a noise environment. A decision has to be taken, in order to say yes, there is a signal, or no, there is no signal. Such decision is usually taken under the application of a decision threshold. Noise contributes off course negatively inducing wrong decisions. For the quantification of detection accuracy, some quantities are usually defined as the probability of detection (P D) and probability of false alarm (P F A), allowing a detection problem analysis.
The maximum eigenvalue of covariance matrices has been frequently used in the literature for the elaboration of certain detection problems. In the Ground Moving Target Indication (GMTI) area, for instance, the signal to clutter plus noise ratio has been studied under the distribution of the sample eigenvalues, which allowed the implementation of a constant false alarm rate detector [9]. The eigenvalues of the covariance matrix of a SAR image pair have been also used in order to formulate a problem in the change detection area [5]. Another example is the determination of the existence of just a single dominant scattering mechanism in a SAR polarimetric acquisition, which also has been evaluated using a threshold in the maximum sample eigenvalue of the polarimetric covariance matrix [11]. Target detection and polarimetric filtering represent other fields of application in which the maximum eigenvalue can be used.
Having the closed form expressions of the PDF (Equation (4)) and/or CDF (Equation (3)) of the maximum sample eigenvalue of the covariance matrix, the P D and P F A when applying a threshold in λ max can be analytically computed, allowing a complete detection problem analysis.

Detection of a Dominant Maximum Eigenvalue
Since several detection problems rely in fact on the choice of a dominant or non-dominant l max by applying a threshold in λ max , the following problem with two hypotheses is elaborated In one application area H 0 may mean that just a single scattering mechanism is present inside a resolution cell of polarimetric SAR data, in other application H 0 can mean that a change happened or a moving target is present in the scene.
Since l max is not achievable, a threshold T has to be used in λ max , originating the following P D and P F A as a function of the decision threshold T .
Hence, for a two-dimensional system configuration, the distribution p λmax (x; H 0 ) is the distribution of λ max when l 1 = 2 and l 2 = 0. For a three-dimensional system, p λmax (x; H 0 ) is determined evaluating the distribution of λ max when l 1 = 3 and l 2 = l 3 = 0. Higher dimensional system configurations follow the same rule.
On the other hand, the distribution p λmax (x; H 1 ) is given by l 1 = l 2 = 1 and l 1 = l 2 = l 3 = 1 for a two-and three-dimensional system, respectively.
For each value of T , there exists a pair (P F A, P D). The curves of P D versus P F A are called Receiver Operating Characteristic (ROC) curves and express the detection performance. The more a ROC curve bends toward the upper left, the better is the detection performance since a higher P D and lower P F A is achieved. Figure 7(a) shows the detection performance as a function of the dimension of the multichannel system, i.e., m = 2, e.g., interferometric, m = 3, e.g., polarimetric and m = 6, e.g., polarimetric-interferometric system. For all cases, with fixed number of samples n = 6, one can realize that the performance of detection is significantly improved as the number of SAR images increases.   In order to state how correlation effects the detection performance, the detection problem is reformulated as follows In this way, every time that correlation is greater than zero the eigenvalues have different values, and one is larger than the others. For the two-dimensional case, for instance, p λmax (x; H 0 ) is the distribution of λ max when l 1 > l 2 , which changes for different values of the correlation ρ. Figure 7(c) presents the ROC curves for this case. For small correlation between channels, the detector suffers from a significant false-alarm rate. As expected, when correlation increases the ROC curve bends toward the upper left, indicating better detection performance. The limit is reached when ρ = 1, meaning that l 1 = 2 and l 2 = 0, which corresponds to the cases of the previous detection problem.

Target Detection Using Polarimetric Matched Filter
In this section the application of the expressions for the detection of specific targets using the Polarimetric Matched Filter (PMF) concept is illustrated. For that, a short review on PMF is required [6,7].
In a polarimetric acquisition, the dimension of the system corresponds to the number of channels, which are in general four, but three for the backscattering reciprocal case. The measured vector k is hence three-or four-dimensional. The elements of the vector k may be weighted in a given way, in order to satisfy a certain condition. The weighting can be accounted for by making y = h † k, where h is a complex vector with same dimension as k. A condition usually aimed to be satisfied in the literature is the maximization of the quadratic detector |y| 2 (Equation (59) in [6]). The optimal weighting vector h in order to detect a distributed target with the target vector k is dependent on the target characteristics. Accordingly, the expected value of |y| 2 is given by where Σ t = E{kk † } is the target covariance matrix. Regarding Rayleigh quotient (Theorem 15.91 in [25]), for any m-dimensional complex vector x and a given m × m Hermitian matrix A, x † Ax ≤ ∥x∥ 2 a max , where a max is the maximum eigenvalue of A. The equality is valid if x is along the direction of the a max eigenvector U max (∥U max ∥ = 1). Hence, for a given distributed target with covariance matrix Σ t , the optimum weighting vector h is given by the eigenvector of the maximum eigenvalue of Σ t . The target vector k is a random sample and has alone no physical meaning. Therefore, usually multi-look processing is performed making, for n samples (or looks) The target is assumed to be immersed in polarization independent clutter. In this way, the detection procedure is thus evaluated by choosing h as the eigenvector corresponding to the maximum eigenvalue of the covariance matrix Σ t of the target to be detected, and applying a detection threshold y 2 > T .
In the presence of target y 2 is maximum, given by y 2 = ∥h∥ 2 λ max , where h is the eigenvector of λ max . Hence, assuming without loss of generality that h has unitary length (∥h∥ 2 = 1), detection performances can be made using the derivations given in this work, as done in the previous section, where a threshold has been used in order to detect a dominant maximum sample eigenvalue (λ max = y 2 > T ).
Note that the probabilities of detection and false alarm defined in last section can be more simply evaluated by A three-dimensional polarimetric target detection problem was formulated making k = [S HH S HV S V V ] T , where S i is the polarimetric scattering matrix element in channel i, and using the target covariance matrix structure where ϵ = E{|S HV | 2 } E{|S HH | 2 } , γ = E{|S V V | 2 } E{|S HH | 2 } and ρ is the complex correlation coefficient between S HH and S V V [6,26]. Figure 8 shows the ROC curves for different number of samples for three different kinds of targets. Figure 8(a) corresponds to an azimuthal symmetric target having covariance matrix with parameters ϵ = 1, γ = 0.5 and ρ = 0. Figure 8(c) corresponds to a reflection symmetric target having covariance matrix with parameters ϵ = 1, γ = 0.8 and ρ = 0.8, while Figure 8(b) to a reflection symmetric target with parameters ϵ = 0, γ = 0.8 and ρ = 0.8 in its covariance matrix. For all three cases, S HH = 1 and the P F A was evaluated using the polarization independent clutter having ϵ = 1, γ = 1 and ρ = 0. Note that the curves vary not just with the number of used samples n but are also different for different targets having different Σ t . This means that some types of targets are easier to detect than others, when using the quadratic detector described here and when the clutter is polarization independent. When the target covariance matrix is similar to the one of the clutter, the detection performance weakens (Figure 8(a)). On the other hand, when the covariance matrix of the target is significantly different from the clutter one, the detection performance improves (Figure 8(c)).

Conclusions and Discussion
In this paper a depth statistical analysis of the maximum eigenvalue of the eigendecomposition of the sample covariance matrix in terms of SAR applications is presented. The proposed analysis are supported by simulation results via several examples. The results are based on a exact closed-form expressions of PDF, CDF and MGF. In this study, existing density functions of the sample maximum eigenvalue were extended and/or implemented into multi-channel SAR system in order to obtain a simple expression of the sample eigenvalues giving a way to fruitful applications. From these closed-form expressions, it has been possible to develop new algorithms for unbiased calculations of parameters extracted from multi-channel SAR covariance matrix. In addition to these implementations, closed-form expressions were developed for the MGF of the sample maximum eigenvalue, which can be critical in the area of bias removal and detection performance analysis. This new closed-form expressions of the MGF can be also interesting for other application areas like MIMO systems (Multiple-Input Multiple-Output). Apart from estimation theory analysis including the MGF, the detection problem of the sample maximum eigenvalue has been also discussed.
This Appendix starts with the reminder of some basic results from linear algebra and statistical theory to prove the Theorem III proposed in this manuscript.
Definition 1 (Theorem A3 in [19]): The determinant of a square m × m matrix A is where ∑ π denotes the summation over all m! permutations, and sgn(π) denotes the sign of the permutations determined by (−1) π , where π is the permutation symbol. Notice that the determinant is a combination of the determinants of its sub-matrices, and permuting the columns or rows of a matrix changes the sign of the determinant related to the permutation.

Proof of Theorem I
In the case m ≤ n, the sample covariance matrix has full rank and follows a central Wishart distribution. It can be noted that [28] is an interesting work concerning the case n ≤ m. In Wishart history, the joint PDF of the real ordered sample eigenvalues ∞ ≥ λ 1 ≥ . . . ≥ λ m ≥ 0 of Z is given by Equation (36) in [18] where Ξ = diag{λ 1 , . . . , λ m }, p Ξ (Ξ) is the joint PDF of the eigenvalues of the sample covariance matrix The CDF of the maximum eigenvalue λ max = λ 1 is obtained using where p(λ 1 , λ 2 , ..., λ m ) is the joint PDF of the eigenvalues as defined in Equation (A3), and Before substituting Equation (A3) in Equation (A4) it is desirable to write a friendlier expression for the joint PDF Equation (A3). If the constant is denoted by S as in Equation (3), then Equation (A3) can be written as Due to the definition of the Vandermonde determinant [25] the multiplication of Then the joint PDF of the sample eigenvalues can be written in the form of By applying the result of Definition 2 and Definition 1 into (A4) with Ψ i (λ j ) = λ n−j i and Φ i (λ j ) = exp , it yields Finally, solving the remaining integral using Equation (351) in [25] yields the desired closed-form expression of the sample maximum eigenvalue CDF of multi-channel covariance matrix denotes the lower incomplete gamma function.

Proof of Moments of the Sample Maximum Eigenvalue
For any positive integer s, the sth moment of the maximum eigenvalue is obtained using where p λmax (x) is the PDF of λ max . It has been shown in Equation (4) that where Ω(x) and Ψ(x) are m × m matrices with ijth element Ω(x) i,j = exp x n−j and , and S denotes the constant terms given in Equation (3). It is desirable to write a friendlier expression for Equation (A11) to be useful for deriving the moments of the sample maximum eigenvalue. Due to the definition of the inverse of the m × m matrix A −1 = C T ij /|A|, Equation (A11) can be written in the following form where C T (x) is the transpose of the cofactor matrix Ψ(x). Since the trace of a square matrix is the sum of the diagonal elements of the matrix, Equation (A11) can be written as the sum of the diagonal elements of D(x) = (C T (x)Ω(x)).
Now the problem in Equation (A11) is the evaluation of the diagonal elements of the D(x). The matrix D(x) is the product of the two square matrices and the trace of the product of two matrices is C ji Ω ji = C 11 Ω 11 + · · · + C m1 Ω m1 + · · · + C 1m Ω 1m + · · · + C mm Ω mm shown as where C ji is the determinant of the sub-matrix M ji obtained from Ψ(x) by removing its jth row and ith column. Then, the j, i th cofactor of Ψ ji , denoted by C ji , is C ji = (−1) i+j |M ji |. Now, Equation (A14) is: where C ji can be extracted from Ψ(x) by Definition I as follows, Here, the sum is computed over (m − 1)! permutations of π s . S m denotes the set of all m! permutations of the set S = {1, 2, ..., m}, and sgn(π sub ) denotes the sign of the permutation π s : +1 if π s is an even permutation and −1 if it is odd. After interchanging the order of the summations and integral in Equation (A16), the integration reduces to the following form It can be seen from Equation (A18) that the integration includes the multiplication of (m − 1) times incomplete gamma functions with exponential and powers. After applying identities (ex. a solution of integration in Equation (A18) may be obtained by using the following formula where M = u 1 +· · ·+u m and A = (1/2)(a 1 +· · ·+a m ) are valid for ℜ(M +v) > 0, and ℜ(b± 1 2 a 1 ±· · ·± 1 2 a n ) > 0. Here M k,m (x) is a Whittaker function of the first kind Equation (7.622) in [25] and F A (α : β 1 , . . . , β m ; γ 1 , . . . , γ m : z 1 , . . . , z m ) is a hypergeometric function of several variables Equation (9.19) in [25]. After replacing the integration part of Equation (A18) by Equation (A21), it can easily be shown that for m dimensional SAR system with n ≥ m, the moments of the sample maximum eigenvalue is ∑ m−1 k=1 a k . Notice that k and π sub (k) takes values from (k, π sub (k) = {k, π sub (k) ∈ {1, 2, . . . , m} , k ̸ = j ∧ π sub (j) ̸ = i}) related to given i, j.
Although Equation (A22) is an exact closed-form expression for the MGF of the sample maximum eigenvalue, the hypergeometric function of several variables is in general very difficult to evaluate, and it is better to express Equation (A22) in a more tractable form for future analysis. Regarding Equation (2.53) in [21], the incomplete gamma function can be written in the closed form γ(n + 1 − π sub (k), x n l k ) = (n − π sub (k))!
for n ∈ Z + , which is a constraint fulfilled by our problem. For the sake of simplicity, f (k, π sub(k) ) = (A24) and the integration can be easily implemented into Equation (A18) as follows (n − π sub (k))!