Central limit theorems for functionals of large sample covariance matrix and
mean vector in matrix-variate location mixture of normal distributions

arXiv:1602.05522v2 [] 13 Jul 2017

Taras Bodnara,1 , Stepan Mazurb , Nestor Parolyac

a
b

Department of Mathematics, Stockholm University, SE-10691 Stockholm, Sweden

Department of Statistics, School of Business, Örebro University, SE-70182 Örebro, Sweden
c

Institute of Statistics, Leibniz University of Hannover, D-30167 Hannover, Germany

Abstract
In this paper we consider the asymptotic distributions of functionals of the sample covariance matrix and the sample mean vector obtained under the assumption that the matrix
of observations has a matrix-variate location mixture of normal distributions. The central
limit theorem is derived for the product of the sample covariance matrix and the sample
mean vector. Moreover, we consider the product of the inverse sample covariance matrix
and the mean vector for which the central limit theorem is established as well. All results
are obtained under the large-dimensional asymptotic regime where the dimension p and the
sample size n approach to infinity such that p/n → c ∈ [0, +∞) when the sample covariance
matrix does not need to be invertible and p/n → c ∈ [0, 1) otherwise.

ASM Classification: 62H10, 62E15, 62E20, 60F05, 60B20
Keywords: Normal mixtures, skew normal distribution, large dimensional asymptotics, stochastic representation, random matrix theory.

1
Corresponding author. E-mail address: taras.bodnar@math.su.se. The second author appreciates the financial support of the Swedish Research Council Grant Dnr: 2013-5180 and Riksbankens Jubileumsfond Grant Dnr:
P13-1024:1

1

1

Introduction

The functions of the sample covariance matrix and the sample mean vector appear in various
statistical applications. The classical improvement techniques for the mean estimation have already been discussed by Stein (1956) and Jorion (1986). In particular, Efron (2006) constructed
confidence regions of smaller volume than the standard spheres for the mean vector of a multivariate normal distribution. Fan et al. (2008), Bai and Shi (2011), Bodnar and Gupta (2011),
Cai and Zhou (2012), Cai and Yuan (2012), Fan et al. (2013), Bodnar et al. (2014a), Wang et al.
(2015), Bodnar et al. (2016a) among others suggested improved techniques for the estimation
of covariance matrix and precision matrix (the inverse of covariance matrix).
In our work we introduce the family of matrix-variate location mixture of normal distributions (MVLMN) which is a generalization of the models considered by Azzalini and Dalla-Valle
(1996), Azzalini and Capitanio (1999), Azzalini (2005), Liseo and Loperfido (2003, 2006), Bartoletti and Loperfido (2010), Loperfido (2010), Christiansen and Loperfido (2014), Adcock et al.
(2015), De Luca and Loperfido (2015) among others. Under the assumption of MVLMN we
consider the expressions for the sample mean vector x and the sample covariance matrix S.
In particulary, we deal with two products l> Sx and l> S−1 x where l is a non-zero vector of
constants. It is noted that this kind of expressions has not been intensively considered in the
literature, although they are present in numerous important applications. The first application
of the products arises in the portfolio theory, where the vector of optimal portfolio weights is
proportional to S−1 x. The second application is in the discriminant analysis where the coefficients of the discriminant function are expressed as a product of the inverse sample covariance
matrix and the difference of the sample mean vectors.
Bodnar and Okhrin (2011) derived the exact distribution of the product of the inverse sample
covariance matrix and the sample mean vector under the assumption of normality, while Kotsiuba and Mazur (2015) obtained its asymptotic distribution as well as its approximate density
based on the Gaussian integral and the third order Taylor series expansion. Moreover, Bodnar
et al. (2013, 2014b) analyzed the product of the sample (singular) covariance matrix and the
sample mean vector. In the present paper, we contribute to the existing literature by deriving
the central limit theorems (CLTs) under the introduced class of matrix-variate distributions in
the case of the high-dimensional observation matrix. Under the considered family of distributions, the columns of the observation matrix are not independent anymore and, thus, the CLTs
cover a more general class of random matrices.
Nowadays, modern scientific data include large number of sample points which is often
comparable to the number of features (dimension) and so the sample covariance matrix and
the sample mean vector are not the efficient estimators anymore. For example, stock markets
include a large number of companies which is often close to the number of available time points.
In order to understand better the statistical properties of the traditional estimators and tests
based on high-dimensional settings, it is of interest to study the asymptotic distribution of the
above mentioned bilinear forms involving the sample covariance matrix and the sample mean
vector.
The appropriate central limit theorems, which do not suffer from the “curse of dimension-

2

ality” and do not reduce the number of dimensions, are of great interest for high-dimensional
statistics because more efficient estimators and tests may be constructed and applied in practice. The classical multivariate procedures are based on the central limit theorems assuming that
the dimension p is fixed and the sample size n increases. However, numerous authors provide
quite reasonable proofs that this assumption does not lead to precise distributional approximations for commonly used statistics, and that under increasing dimension asymptotics the better
approximations can be obtained [see, e.g., Bai and Silverstein (2004) and references therein].
Technically speaking, under the high-dimensional asymptotics we understand the case when
the sample size n and the dimension p tend to infinity, such that their ratio p/n converges to
some positive constant c. Under this condition the well-known Marčenko-Pastur and Silverstein
equations were derived [see, Marčenko and Pastur (1967), Silverstein (1995)].
The rest of the paper is structured as follows. In Section 2 we introduce a semi-parametric
matrix-variate location mixture of normal distributions. Main results are given in Section 3,
where we derive the central limit theorems under high-dimensional asymptotic regime of the
(inverse) sample covariance matrix and the sample mean vector under the MVLMN. Section 4
presents a short numerical study in order to verify the obtained analytic results.

2

Semi-parametric family of matrix-variate location mixture of
normal distributions

In this section we introduce the family of MVLMN which generalizes the existent families of
skew normal distributions.
Let


x11 . . . x1n
 .
.. 
..
..
X=
.
. 

 = (x1 , ..., xn ) ,
xp1 . . . xpn
be the p × n observation matrix where xj is the j th observation vector. In the following, we
assume that the random matrix X possesses a stochastic representation given by
d

X = Y + Bν1>
n,

(1)

where Y ∼ Np,n (µ1>
n , Σ ⊗ In ) (p × n-dimensional matrix-variate normal distribution with mean
>
matrix µ1n and covariance matrix Σ ⊗ In ), ν is a q-dimensional random vector with continuous
density function fν (·), B is a p × q matrix of constants. Further, it is assumed that Y and ν
are independently distributed. If random matrix X follows model (1), then we say that X is
MVLMN distributed with parameters µ, Σ, B, and fν (·). The first three parameters are finite
dimensional, while the fourth parameter is infinite dimensional. This makes model (1) to be of
a semi-parametric type. The assertion we denote by X ∼ LMN p,n;q (µ, Σ, B; fν ). If fν can be
parametrized by finite dimensional parameter θ, then (1) reduces to a parametrical model which
is denoted by X ∼ LMN p,n;q (µ, Σ, B; θ). If n = 1, then we use the notation LMN p;q (·, ·, ·; ·)

3

instead of LMN p,1;q (·, ·, ·; ·).
From (1) the density function of X is expressed as
Z
fX (Z) =
Rq

∗
∗
fNp,n (µ1>
(Z − Bν ∗ 1>
n )fν (ν )dν .
n ,Σ⊗In )

(2)

In a special case when ν = |ψ| is the vector formed by the absolute values of every element in
ψ where ψ ∼ Nq (0, Ω), i.e. ν has a q-variate truncated normal distribution, we get
Proposition 1. Assume model (1). Let ν = |ψ| with ψ ∼ Nq (0, Ω). Then the density function
of X is given by




e −1 Φq 0; −DEvec(Z − µ1> ), D φpn vec(Z − µ1> ); 0, F
fX (Z) = C
n
n

(3)

> −1
−1 − E> DE)−1 , and
where D = (nB> Σ−1 B + Ω−1 )−1 , E = 1>
n ⊗ B Σ , F = (In ⊗ Σ
1/2

1/2

e −1 = C −1 |F| |D|
C
with C = Φq (0; 0, Ω).
|Ω|1/2 |Σ|n/2
The proof of Proposition 1 is given in the Appendix.
It is remarkable that model (1) includes several skew-normal distributions considered by
Azzalini and Dalla-Valle (1996), Azzalini and Capitanio (1999), Azzalini (2005). For example,
in case of n = 1, q = 1, µ = 0, B = ∆1p , and Σ = (Ip − ∆2 )1/2 Ψ(Ip − ∆2 )1/2 we get
d

X = (Ip − ∆2 )1/2 v0 + ∆1p |v1 |,

(4)

where v0 ∼ Np (0, Ψ) and v1 ∼ N (0, 1) are independently distributed; Ψ is a correlation matrix
and ∆ = diag(δ1 , ..., δp ) with δj ∈ (−1, 1). Model (4) was previously introduced by Azzalini
(2005).
Moreover, model (1) also extends the classical random coefficient growth-curve model (see,
Potthoff and Roy (1964), Amemiya (1994) among others), i.e., the columns of X can be rewritten
in the following way
d

xi = µ + Bν + εi , i = 1, . . . , n ,

(5)

where ε1 , . . . , εn are i.i.d. Np (0, Σ), ν ∼ fν and ε1 , . . . , εn , ν are independent. In the random
coefficient growth-curve model, it is typically assumed that εi ∼ Np (0, σ 2 I) and ν ∼ Nq (0, Ω)
(see, e.g., Rao (1965)). As a result, the suggested MVLMN model may be useful for studying
the robustness of the random coefficients against the non-normality.

4

3

CLTs for expressions involving the sample covariance matrix
and the sample mean vector

The sample estimators for the mean vector and the covariance matrix are given by
n

1X
1
x=
xi = X1n
n
n

n

and

i=1

1 X
S=
(xi − x)(xi − x)> = XVX> ,
n−1
i=1

>
2
where V = In − n1 1n 1>
n is a symmetric idempotent matrix, i.e., V = V and V = V.
The following proposition shows that x and S are independently distributed and presents
their marginal distributions under model (1). Moreover, its results lead to the conclusion that
the independence of x and S could not be used as a characterization property of a multivariate
normal distribution if the observation vectors in data matrix are dependent.

Proposition 2. Let X ∼ LMN p,n;q (µ, Σ, B; fν ). Then
(a) (n−1)S ∼ Wp (n−1, Σ) (p-dimensional Wishart distribution for p ≤ n−1 and p-dimensional
singular Wishart distribution for p > n − 1 with (n − 1) degrees of freedom and covariance
matrix Σ),

(b) x ∼ LMN p;q µ, n1 Σ, B; fν ,
(c) S and x are independently distributed.
Proof. The statements of the proposition follow immediately from the fact that
n

x̄ =

1X
xi = ȳ + Bν and S = XVX> = YVY> .
n

(6)

i=1

Indeed, from (6), the multivariate normality of Y and the independence of Y and ν, we get that
x and S are independent; S is (singular) Wishart distributed with (n − 1) degrees of freedom
and covariance matrix Σ; x has a location mixture of normal distributions with parameters µ,
1
n Σ, B and fν .
For the validity of the asymptotic results presented in Sections 3.1 and 3.2 we need the
following two conditions
(A1) Let (λi , ui ) denote the set of eigenvalues and eigenvectors of Σ. We assume that there
exist m1 and M1 such that
0 < m1 ≤ λ1 ≤ λ2 ≤ ... ≤ λp ≤ M1 < ∞
uniformly in p.
(A2) There exists M2 < ∞ such that
>
|u>
i µ| ≤ M2 and |ui bj | ≤ M2 for all i = 1, ..., p and j = 1, ..., q

uniformly in p where bj , j = 1, ..., q, are the columns of B.
5

Generally, we say that an arbitrary p-dimensional vector l satisfies the condition (A2) if |u>
i l| ≤
M2 < ∞ for all i = 1, . . . , p.
Assumption (A1) is a classical condition in random matrix theory (see, Bai and Silverstein
(2004)), which bounds the spectrum of Σ from below as well as from above. Assumption (A2) is
a technical one. In combination with (A1) this condition ensures that p−1 µ> Σµ, p−1 µ> Σ−1 µ,
p−1 µ> Σ3 µ, p−1 µ> Σ−3 µ, as well as that all the diagonal elements of B> ΣB, B> Σ3 B, and
B> Σ−1 B are uniformly bounded. All these quadratic forms are used in the statements and the
proofs of our results. Note that the constants appearing in the inequalities will be denoted by
M2 and may vary from one expression to another. We further note, that assumption (A2) is
automatically fulfilled if µ and B are sparse (not Σ itself). More precisely, a stronger condition,
which verifies (A2) is ||µ|| < ∞ and ||bj || < ∞ for j = 1, . . . , q uniformly in p. Indeed, using
2
the Cauchy-Schwarz inequality to (u>
i µ) we get
2
2
2
2
(u>
i µ) ≤ ||ui || ||µ|| = ||µ|| < ∞ .

(7)

It is remarkable that ||µ|| < ∞ is fulfilled if µ is indeed a sparse vector or if all its elements are of

order O p−1/2 . Thus, a sufficient condition for (A2) to hold would be either sparsity of µ and B
or their elements are reasonably small but not exactly equal to zero. Note that the assumption of
sparsity for B is quite natural in the context of high-dimensional random coefficients regression
models. It implies that if B is large dimensional, then it may have many zeros and there exists
a small set of highly significant random coefficients which drive the random effects in the model
(see, e.g., Bühlmann and van de Geer (2011)). Moreover, it is hard to estimate µ in a reasonable
way when ||µ|| → ∞ as p → ∞ (see, e.g., Bodnar et al. (2016c)). In general, however, ||µ||
and ||bj || do not need to be bounded, in such case the sparsity of eigenvectors of matrix Σ may
guarantee the validity of (A2). Consequently, depending on the properties of the given data set
the proposed model may cover many practical problems.

3.1

CLT for the product of sample covariance matrix and sample mean vector

In this section we present the central limit theorem for the product of the sample covariance
matrix and the sample mean vector.
Theorem 1. Assume X ∼ LMN p,n;q (µ, Σ, B; fν ) with Σ positive definite and let p/n =
c + o(n−1/2 ), c ∈ [0, +∞) as n → ∞. Let l be a p-dimensional vector of constants that satisfies
condition (A2). Then, under (A1) and (A2) it holds that
√

where



D
−1
nσν
l> Sx − l> Σµν −→ N (0, 1) for p/n → c ∈ [0, +∞) as n → ∞ ,

µν
2
σν

= µ + Bν,


tr(Σ2 ) >
>
= µν Σµν + c
l Σl + (l> Σµν )2 + l> Σ3 l .
p

(8)

(9)
(10)

Proof. First, we consider the case of p ≤ n−1, i.e., S has a Wishart distribution. Let L = (l, x)>
6

e = LSL> = {S
e ij }i,j=1,2 with Se11 = l> Sl, Se12 = l> Sx, Se21 = x> Sl, and Se22 = x> Sx.
and define S
e = LΣL> = {Σ
e ij }i,j=1,2 with Σ
e 11 = l> Σl, Σ
e 12 = l> Σx, Σ
e 21 = x> Σl, and
Similarly, let Σ
e 22 = x> Σx.
Σ


1
Since S and x are independently distributed, S ∼ Wp n − 1,
Σ and rank L =
n−1
e
2 ≤ p with probability
one, we get from Theorem 3.2.5 of Muirhead (1982) that S|x
∼

1 e
W2 n − 1,
Σ . As a result, the application of Theorem 3.2.10 of Muirhead (1982) leads
n−1
to


1 e
−1 e
e
e
e
e
e
Σ11·2 S22 ,
S12 |S22 , x ∼ N Σ12 Σ22 S22 ,
n−1
e 11·2 = Σ
e 11 − Σ
e 2 /Σ
e 22 is the Schur complement.
where Σ
12
e 22 , then
Let ξ = (n − 1)Se22 /Σ


>

l Sx|ξ, x ∼ N


ξ >
ξ
>
>
>
2
[x Σxl Σl − (x Σl) ] .
l Σx,
n−1
(n − 1)2

From Theorem 3.2.8 of Muirhead (1982) it follows that ξ and x are independently distributed
and ξ ∼ χ2n−1 . Hence, the stochastic representation of l> Sx is given by
ξ >
l Sx =
l Σx +
n−1
>

d

r

ξ
z0
(x> Σxl> Σl − (l> Σx)2 )1/2 √
,
n−1
n−1

(11)


where ξ ∼ χ2n−1 , z0 ∼ N (0, 1), x ∼ LMN p;q µ, n1 Σ, B; fν ; ξ, z0 and x are mutually independ

dent. The symbol ”=” stands for the equality in distribution.
It is remarkable that the stochastic representation (11) for l> Sx remains valid also in the
case of p > n − 1 following the proof of Theorem 4 in Bodnar et al. (2014b). Hence, we will
make no difference between these two cases in the remaining part of the proof.
From the properties of χ2 -distribution and using the fact that n/(n − 1) → 1 as n → ∞, we
immediately receive
√


n

ξ
−1
n



D

−→ N (0, 2) as n → ∞ .

(12)

√
√
We further get that n(z0 / n) ∼ N (0, 1) for all n and, consequently, it is also its asymptotic
distribution.
Next, we show that l> Σx and x> Σx are jointly asymptotically normally distributed given
ν. For any a1 and a2 , we consider
>

>

a1 x Σx + 2a2 l Σx = a1
where x̃ = x̄ +

a2
a1 1



a2
x+ l
a1

>


a2
a2
a2
Σ x + l − 2 l> Σl = a1 x̃> Σx̃ − 2 l> Σl ,
a1
a1
a1



a2
and x̃|ν ∼ Np µa,ν , n1 Σ with µa,ν = µ + Bν + l. By Provost and
a1

7

Rudiuk (1996) the stochastic representation of x̃> Σx̃ is given by
p

1X 2
x̃ Σx̃ =
λ i ξi ,
n
>

d

i=1

√ −1/2 >
where ξ1 , ..., ξp given ν are independent with ξi |ν ∼ χ21 (δi2 ), δi = nλi
ui µa,ν . Here, the
2
2
symbol χd (δi ) denotes a chi-squared distribution with d degrees of freedom and non-centrality
parameter δi2 .
Now, we apply the Lindeberg CLT to the conditionally independent random variables Vi =
2
λi ξi /n. For that reason, we need first to verify the Lindeberg’s condition. Denoting σn2 =
P
V( pi=1 Vi |ν) we get
σn2 =

p
X

V



 X
p

1
λ4i
λ2i
2(1 + 2δi2 ) = 2 2tr(Σ4 ) + 4nµ0a,ν Σ3 µa,ν
ξi ν =
2
n
n
n
i=1

i=1

We need to check if for any small ε > 0 it holds that
p
i
1 X h
2
E
(V
−
E
(V
))
1
ν
−→ 0 .
i
i
{|Vi −E(Vi )|>εσn }
2
n→∞ σn

lim

i=1

First, we get
p
X

h

E (Vi − E(Vi |ν))2 1{|Vi −E(Vi |ν )|>εσn } ν

i

i=1
Cauchy−Schwarz

≤

Chebyshev

≤

p
X
i=1
p
X
i=1

h

i

E1/2 (Vi − E(Vi |ν))4 ν P1/2 { |Vi − E(Vi |ν)| > εσn | ν}
λ4i
n2

q
σi
12(1 + 2δi2 )2 + 48(1 + 4δi2 )
εσn

8

(13)

with σi2 = V(Vi |ν) and, thus,

≤

=
≤
≤
≤
≤

p
i
1 X h
2
E
(V
−
E
(V
|ν))
1
ν
i
i
{|Vi −E(Vi |ν )|>εσn }
σn2
i=1
q
Pp
4 12(1 + 2δ 2 )2 + 48(1 + 4δ 2 ) σi
λ
i
i σn
1 i=1 i
4
3
>
ε
2tr(Σ ) + 4nµa,ν Σ µa,ν
q
P
√
p
4 (5 + 2δ 2 )2 − 20 σi
λ
i
3 i=1 i
σn
4
3
>
ε
tr(Σ ) + 2nµa,ν Σ µa,ν
Pp
√
4
2 σi
3
i=1 λi (5 + 2δi ) σn
3
ε tr(Σ4 ) + 2nµ>
a,ν Σ µa,ν
√
3
3 5tr(Σ4 ) + 2nµ>
a,ν Σ µa,ν σmax
3
ε tr(Σ4 ) + 2nµ>
σn
a,ν Σ µa,ν
!
√
4
3
σmax
+1
3
4
>
ε
σn
1 + 2nµa,ν Σ µa,ν /tr(Σ )
√
5 3 σmax
.
ε σn

Finally, Assumptions (A1) and (A2) yield
2
2
supi λ4i + 2nλ3i (u>
σmax
supi σi2
supi λ4i (1 + 2δi2 )
i µa,ν )
=
−→ 0 ,
=
=
3
3
σn2
σn2
tr(Σ4 ) + 2nµ>
tr(Σ4 ) + 2nµ>
a,ν Σ µa,ν
a,ν Σ µa,ν

(14)

which verifies the Lindeberg condition since
2
(u>
i µa,ν )

=
=
(A2)

≤
=


2
>
>
> a2
ui µ + ui Bν + ui l
a1
2

a2 > >
a2
2
>
>
+ (u>
u l(ui µ + u>
(u0i µ)2 + u>
i l
i Bν) + 2ui µ · ui Bν + 2
i Bν)
a1
a1 i
p
p
a2
|a2 |
M22 + qM22 ν > ν + M22 22 + 2M22 qν > ν + 2M22
(1 + qν > ν)
|a1 |
a1


p
|a2 | 2
M22 1 + qν > ν +
< ∞.
(15)
|a1 |

Thus, using (13) and
p
X
i=1

E(Vi |ν) =

p
X
λ2
i

i=1

n

(1 + δi2 ) = tr(Σ2 )/n + µ>
a,ν Σµa,ν

we get that
√ x̃> Σx̃ − tr(Σ2 )/n − µ>
D
a,ν Σµa,ν
n q
ν −→ N (0, 2)
tr(Σ4 )/n + 2µ0a,ν Σ3 µa,ν

9

(16)

and for a1 x> Σx + 2a2 l> Σx we have
 a22 >
> Σx + 2a l> Σx − a tr(Σ2 )/n + µ> Σµ
a
l Σl
x
2
1
a,ν +
a,ν
√ 1
D
a1
q
n
ν −→ N (0, 2) .

4
3
a21 tr(Σ )/n + 2µ0a,ν Σ µa,ν
Denoting a = (a1 , 2a2 )> and µν = µ + Bν we can rewrite it as
√

"
n a

>

x> Σx
l> Σx

2

!
−a

>

tr(Σ
µ>
ν Σµν + c p
l> Σµν

which implies that the vector

√

)

4

!#
D

ν −→ N

0, a

>

)
3
2c tr(Σ
+ 4µ>
ν Σ µν
p
2l> Σ3 µν

2l> Σ3 µν
l> Σ3 l

! !
a


>
2
tr(Σ )
>
>
>
>
n x Σx − µν Σµν − c p , l Σx − l Σµν
has asymp-

totically multivariate normal distribution conditionally on ν because the vector a is arbitrary.
Taking into account (17), (12) and the fact that ξ, z0 and x are mutually independent we
get


ξ
n
>


√ 
 x Σx
n
>

 l Σx
z0
√
n





1

2
 
)
µ>
Σµν + c tr(Σ
 
ν
p

−
 
l> Σµν
0









D
 ν −→
N



 
  0
0, 
 
  0
0

0

0

0



)
3
2c tr(Σ
+ 4µ>
ν Σ µν
p
2l> Σ3 µν
0

2l> Σ3 µν
l> Σ3 l
0

0
0
1



 .



2
4

The application of the multivariate delta method leads to
√



D
−1
l> Sx − l> Σµν ν −→ N (0, 1)
nσν

(18)

where
tr(Σ2 )
σν = (l Σµν ) + l Σ l + l Σl µν Σµν + c
p
2

>

2

>

3

>



>



The asymptotic distribution does not depend on ν and, thus, it is also the unconditional asymptotic distribution.
Theorem 1 shows that properly normalized bilinear form l> Sx itself can be accurately approximated by a mixture of normal distributions with both mean and variance depending on ν.
Moreover, this central limit theorem delivers the following approximation for the distribution of
l> Sx, namely for large n and p we have
p



2
p−2 σν
−1 >
l Sx|ν ≈ CN p l Σµν ,
,
n

−1 >

(19)

i.e., it has a compound normal distribution with random mean and variance.
The proof of Theorem 1 shows, in particular, that its key point is a stochastic representation
of the product l> Sx which can be presented using a χ2 distributed random variable, a standard
normally distributed random variable, and a random vector which has a location mixture of
normal distributions. Both assumptions (A1) and (A2) guarantee that the asymptotic mean
2 are bounded with probability one and the cop−1 l> Σµν and the asymptotic variance p−2 σν
variance matrix Σ is invertible as the dimension p increases. Note that the case of standard

10

(17)

asymptotics can be easily recovered from our result if we set c → 0.
Although Theorem 1 presents a result similar to the classical central limit theorem, it is not
a proper CLT because of ν which is random and, hence, it provides an asymptotic equivalence
only. In the current settings it is not possible to provide a proper CLT without additional
assumptions imposed on ν. The reason is the finite dimensionality of ν denoted by q which
is fixed and independent of p and n. This assures that the randomness in ν will not vanish
asymptotically. So, in order to justify a classical CLT we need to have an increasing value of
q = q(n) with the following additional assumptions:
(A3) Let

E(ν) = ω and Cov(ν) = Ω. For any vector l it holds that
√ l> ν − l> ω D
q √
−→ z̃0
l> Ωl

(20)

where the distribution of z̃0 may in general depend on l.
(A4) There exists M2 < ∞ such that
|u>
i Bω| ≤ M2 for all i = 1, ..., p uniformly in p.
Moreover, BΩB> has an uniformly bounded spectral norm.
Assumption (A3) assures that the vector of random coefficients ν satisfies itself a sort of concentration property in the sense that any linear combination of its elements asymptotically
concentrates around some random variable. In a special case of z̃0 being standard normally
distributed (A3) will imply a usual CLT for a linear combination of vector ν. As a result,
assumption (A3) is a pretty general and natural condition on ν taking into account that it is
assumed ν ∼ Nq (ω, Ω) (see, e.g., Rao (1965)) in many practical situations. Assumption (A4)
is similar to (A1) and (A2) and it ensures that µ and Bω as well as Σ and BΩB> have the
same behaviour as p → ∞.
Corollary 1 (CLT). Under the assumptions of Theorem 1, assume (A3) and (A4). Let l> Σl =
O(p) and µ> Σµ = O(p). Then it holds that

√ −1  >
D
l Sx̄ − l> Σ(µ + Bω) −→ N (0, 1)
nσ

for p/n = c + o n−1/2 with c ≥ 0 and q/n → γ with γ > 0 as n → ∞ where
σ

2



tr(Σ2 )
>
= (l Σ(Bω + µ)) + l Σ l + l Σl (µ + Bω) Σ(µ + Bω) + c
.
p
>

2

>

3

>

11

Proof. Since


2
(u>
i µa,ν )

=

u>
i µ

+

u>
i Bω

+

u>
i B(ν

− ω) +

a2
u>
i l
a1

2

a2
M22 + M22 + M22 q(ν − ω)> (ν − ω) + M22 22
a



q 1
|a2 |
|a2 |
2M22 1 + 2
+ 2M22 2 +
q(ν − ω)> (ν − ω) < ∞
|a1 |
|a1 |

(A2),(A4)

≤
+

with probability one, we get from the proof of Theorem 1 and assumption (A3) that
√

nσ

−1



√ √>
σν
n l ΣBΩB> Σl
l Sx̄ − l Σ(µ + Bω) =
z0 + √
z̃0 + oP (1) ,
σ
q
σ
>

>



(21)

where z0 ∼ N (0, 1) and z̃0 are independently distributed.
Let λmax (A) denotes the largest eigenvalue of a symmetric positive definite matrix A. First,
we note that for any vector m that satisfies (A2), we get
m> ΣBΩB> Σm ≤ λmax (Σ)2 λmax (UT BΩB> U)

p
X

|m> ui | = O(p),

i=1

where U is the eigenvector matrix of Σ which is a unitary matrix and, consequently, (see, e.g.,
Chapter 8.5 in Lütkepohl (1996))
λmax (UT BΩB> U) = λmax (BΩB> ) = O(1).
As a result, we get l> ΣB(ν − ω) = OP (1) and (µ + Bω)> ΣB(ν − ω) = OP (1).
Furthermore, it holds that
(ν − ω)> B> ΣB(ν − ω) ≤ λmax (Σ)

p
X

p

2
(u>
i B(ν − ω)) = OP (1)

i=1

i=1

≤ OP (1)λmax (BΩB> )
Hence, we have

√

1X >
ui BΩB> ui
q

1
q

p
X

u>
i ui = OP (1).

i=1

l> ΣBΩB> Σl
→0
σ

and, consequently, together with z̃0 = OP (1) and q/n → γ > 0 it ensures that the second
summand in (21) will disappear. At last,
2
σν
a.s.
→1
σ2


for p/n = c + o n−1/2 with c ≥ 0 and q/n → γ with γ > 0 as n → ∞.
Note that the asymptotic regime q/n → γ > 0 ensures that the dimension of the random
coefficients is not increasing much faster than the sample size, i.e., q and n are comparable. It
12

is worth mentioning that due to Corollary 1 the asymptotic distribution of the bilinear form
l> Sx̄ does not depend on the covariance matrix of the random effects Ω. This knowledge may
be further useful in constructing a proper pivotal statistic.

3.2

CLT for the product of inverse sample covariance matrix and sample
mean vector

In this section we consider the distributional properties of the product of the inverse sample
covariance matrix S−1 and the sample mean vector x. Again we prove that proper weighted
bilinear forms involving S−1 and x have asymptotically a normal distribution. This result is
summarized in Theorem 2.
Theorem 2. Assume X ∼ LMN p,n;q (µ, Σ, B; fν ), p < n − 1, with Σ positive definite and let
p/n = c + o(n−1/2 ), c ∈ [0, 1) as n → ∞. Let l be a p-dimensional vector of constants that
satisfies (A2). Then, under (A1) and (A2) it holds that

√ −1 > −1
nσ̃ν l S x −

1 > −1
l Σ µν
1−c



D

−→ N (0, 1)

(22)

where µν = µ + Bν and
2

σ̃ν

=

1
(1 − c)3



>

−1

l Σ

µν

2

>

−1

+l Σ

>

−1

l(1 + µν Σ


µν ) .

Proof. From Theorem 3.4.1 of Gupta and Nagar (2000) and Proposition 2 we get

S−1 ∼ IW p n + p, (n − 1)Σ−1 .
It holds that
l> S−1 x = (n − 1)x> Σ−1 x

l> S−1 x
x> S−1 x
.
x> S−1 x (n − 1)x> Σ−1 x

Since S−1 and x are independently distributed, we get from Theorem 3.2.12 of Muirhead
(1982) that
x> Σ−1 x
ξe = (n − 1) > −1 ∼ χ2n−p
x S x

(23)

and it is independent of x. Moreover, the application of Theorem 3 in Bodnar and Okhrin
(2008) proves that x> S−1 x is independent of l> S−1 x/x> S−1 x for given x. As a result, it is also
independent of x> Σ−1 x · l> S−1 x/x> S−1 x for given x and, consequently,
(n − 1)x> Σ−1 x

l> S−1 x
x> S−1 x
and
.
x> S−1 x
(n − 1)x> Σ−1 x

are independent.

13

From the proof of Theorem 1 of Bodnar and Schmid (2008) we obtain
>

−1

(n − 1)x Σ



> −1
l> S−1 x
> −1
2 x Σ x >
x > −1 x ∼ t n − p + 1; (n − 1)l Σ x; (n − 1)
l Rx l ,
n−p+1
x S x

(24)

where Ra = Σ−1 − Σ−1 aa> Σ−1 /a> Σ−1 a, a ∈ IRp , and the symbol t(k, µ, τ 2 ) denotes a tdistribution with k degrees of freedom, location parameter µ and scale parameter τ 2 .
Combining (23) and (24), we get the stochastic representation of l> S−1 x given by



> Σ−1 x
x
d
· l> Rx l
l> S−1 x = ξe−1 (n − 1) l> Σ−1 x + t0
n−p+1


q
√
t0
−1
> −1
>
>
−1
e
√
l Σ l x Rl x ,
= ξ (n − 1) l Σ x +
n−p+1
s

where (see, Proposition 2)
d

l> Σ−1 x = l> Σ−1 (µ + Bν) +

l> Σ−1 l 0
√
z ,
n 0

e ν, z 0 , and t0 are mutually independent.
with ξe ∼ χ2n−p , z00 ∼ N (0, 1), and t0 ∼ t(n−p+1, 0, 1); ξ,
0
Since Rl ΣRl = Rl , tr(Rl Σ) = p − 1, and Rl ΣΣ−1 l = 0, the application of Corollary 5.1.3a
and Theorem 5.5.1 in Mathai and Provost (1992) leads to
nx> Rl x|ν ∼ χ2p−1 (nδ 2 (ν)) with δ 2 (ν) = (µ + Bν)> Rl (µ + Bν)

(25)

as well as l> Σ−1 x and x> Rl x are independent given ν. Finally, using the stochastic representation of a t-distributed random variable, we get
z000

d

t0 =

q

(26)

ζ
n−p+1

with ζ ∼ χ2n−p+1 , z000 ∼ N (0, 1) and z00 ∼ N (0, 1) being mutually independent. Together with
(25) it yields
> −1

l S

r


z000 √ > −1
p−1
z00 p > −1
−1
> −1
e
x = ξ (n − 1) l Σ (µ + Bν) + √
l Σ l+ √
l Σ l
η
n−p+1
n
n
r


√
p−1
z0
−1
> −1
>
−1
e
= ξ (n − 1) l Σ (µ + Bν) + l Σ l 1 +
η√
,
(27)
n−p+1
n
d

where the last equality in (27) follows from the fact that
z00

+

z000

r

p−1
d
η = z0
n−p+1

>

r
1+

p−1
η
n−p+1

Rl x/(p−1)
with η = nxζ/(n−p+1)
|ν ∼ Fp−1,n−p+1 (nδ 2 (ν)) (non-central F -distribution with p − 1 and
n − p + 1 degrees of freedom and non-centrality parameter nδ 2 (ν)). Moreover, we have that

14

e z0 , and η are mutually independent.
ξe ∼ χ2n−p and z0 ∼ N (0, 1); ξ,
From Lemma 6.4.(b) in Bodnar et al. (2016b) we get


√ 
n 









e
ξ/(n
− p)
1
D
 

2
η
 −  1 + δ (ν)/c  ν −→ N
√
z0 / n
0

 
 
0, 
 
 

2/(1 − c) 0
0
ση2
0
0
0
0

0
0
1
0








for p/n = c + o(n−1/2 ), c ∈ [0, 1) as n → ∞ with
ση2

2
=
c




2
δ 2 (ν)
2
δ 2 (ν)
1+2
+
1+
c
1−c
c

Consequently,


 

e
ξ/(n
− 1)
(1 − c)
√ 
 

n  (p − 1)η/(n − p + 1)  −  (c + δ 2 (ν))/(1 − c)  ν
√
z0 / n
0
 

2(1 − c)
0
0
D
 

2
2
2
−→ N 0, 
0
c ση /(1 − c) 0 
0
0
1
for p/n = c + o(n−1/2 ), c ∈ [0, 1) as n → ∞.
Finally, the application of the delta-method (c.f. DasGupta (2008, Theorem 3.7)) leads to
√


n l> S−1 x −


1 > −1
D
2
l Σ (µ + Bν) ν −→ N (0, σ̃ν
)
1−c

for p/n = c + o(n−1/2 ), c ∈ [0, 1) as n → ∞ with
1
σ̃ν =
(1 − c)3
2

 

2
> −1
> −1
2
2 l Σ (µ + Bν) + l Σ l(1 + δ (ν)) .

Consequently,
√



1 > −1
D
> −1
nσ̃ν l S x −
l Σ (µ + Bν) ν −→ N (0, 1) ,
1−c
−1

where the asymptotic distribution does not depend on ν. Hence, it is also the unconditional
asymptotic distribution.
Again, Theorem 2 shows that the distribution of l> S−1 x can be approximated by a mixture
of normal distributions. Indeed,
p

−1 > −1

l S


x|ν ≈ CN

2
p−1 > −1
p−2 σ̃ν
l Σ µν ,
1−c
n


.

(28)

In the proof of Theorem 2 we can read out that the stochastic representation for the product
15

of the inverse sample covariance matrix and the sample mean vector is presented by using a χ2
distributed random variable, a general skew normally distributed random vector and a standard
t-distributed random variable. This result is itself very useful and allows to generate the values
of l> S−1 x by just simulating three random variables from the standard univariate distributions
and a random vector ν which determines the family of the matrix-variate location mixture of
normal distributions. The assumptions about the boundedness of the quadratic and bilinear
forms involving Σ−1 plays here the same role as in Theorem 1. Note that in this case we need
no assumption either on the Frobenius norm of the covariance matrix or its inverse.
Finally, in Corollary 2 we formulate the CLT for the product of the inverse sample covariance
matrix and the sample mean vector.
Corollary 2 (CLT). Under the assumptions of Theorem 2, assume (A3) and (A4). Let
l> Σ−1 l = O(p) and µ> Σ−1 µ = O(p). Then it holds that
√

nσ̃

−1


l> S−1 x −


1 > −1
D
l Σ (µ + Bω) −→ N (0, 1)
1−c


for p/n = c + o n−1/2 with c ∈ [0, 1) and q/n → γ with γ > 0 as n → ∞ where
σ̃

2

=

1
(1 − c)3



>

−1

l Σ

(µ + Bω)

2

>

−1

+l Σ

>

−1

l(1 + (µ + Bω) Σ


(µ + Bω)) .

Proof. From the proof of Theorem 2 and assumption (A3) it holds that
√

nσ̃

−1


l> S−1 x −

1 > −1
l Σ (µ + Bω)
1−c



√ √ > −1
n l Σ BΩB> Σ−1 l
σ̃ν
1
=
z0 +
z̃0 + oP (1) ,
√
σ̃
1−c q
σ̃

where z0 ∼ N (0, 1) and z̃0 are independently distributed.
Finally, in a similar way like in Corollary 1 we get from (A4) that
(1 − c)−2 γ −1

l> Σ−1 BΩB> Σ−1 l
σ̃ 2
2
σ̃ν
σ̃ 2

→
a.s.

→

0,

(29)

1

(30)


for p/n = c + o n−1/2 with c ∈ [0, 1) and q/n → γ with γ > 0 as n → ∞.

4

Numerical study

In this section we provide a Monte Carlo simulation study to investigate the performance of the
suggested CLTs for the products of the (inverse) sample covariance matrix and the sample mean
vector.
In our simulations we put l = 1p , each element of the vector µ is uniformly distributed
on [−1, 1] while each element of the matrix B is uniformly distributed on [0, 1]. Also, we take
Σ as a diagonal matrix where each diagonal element is uniformly distributed on [0, 1]. It can
be checked that in such a setting the assumptions (A1) and (A2) are satisfied. Indeed, the

16

population covariance matrix satisfies the condition (A1) because the probability of getting
exactly zero eigenvalue equals to zero. On the other hand, the condition (A2) is obviously valid
too because the ith eigenvector of Σ is ui = ei = (0, . . . , 1 , 0, . . . , 0)0 .
ith place

In order to define the distribution for the random vector ν, we consider two special cases.
In the first case we take ν = |ψ|, where ψ ∼ Nq (0, Iq ), i.e., ν has a q-variate truncated normal
distribution. In the second case we put ν ∼ GALq (Iq , 1q , 10), i.e., ν has a q-variate generalized
asymmetric Laplace distribution (c.f., Kozubowski et al. (2013)). Also, we put q = 10.
We compare the results for several values of c ∈ {0.1, 0.5, 0.8, 0.95}. The simulated data
consist of N = 105 independent realizations which are used to fit the corresponding kernel
density estimators with Epanechnikov kernel. The bandwith parameters are determined via
cross-validation for every sample. The asymptotic distributions are simulated using the results
of Theorems 1 and 2. The corresponding algorithm is given next:
a) generate ν = |ψ|, where ψ ∼ Nq (0q , Iq ), or generate ν ∼ GALq (Iq , 1q , 10);
b) generate l> Sx by using the stochastic representation (11) obtained in the proof of Theorem
1, namely
√
ξ >
ξ
l Σ(y + Bν) +
((y + Bν)> Σ(y + Bν)l> Σl − (l> Σ(y + Bν))2 )1/2 z0 ,
l Sx =
n−1
n−1
>

d

where ξ ∼ χ2n−1 , z0 ∼ N (0, 1), y ∼ Np (µ, n1 Σ); ξ, z0 , y, and ν are mutually independent
b’) generate l> S−1 x by using the stochastic representation (27) obtained in the proof of Theorem 2, namely
> −1

l S

r

√
−1
> −1
>
−1
e
x = ξ (n − 1) l Σ (µ + Bν) + l Σ l 1 +
d

p−1
z0
η√
n−p+1
n


,

where ξe ∼ χ2n−p , z0 ∼ N (0, 1), and η ∼ Fp−1,n−p+1 (nδ 2 (ν)) with δ 2 (ν) = (µ+Bν)> Rl (µ+
e z0 and (η, ν) are mutually independent.
Bν), Rl = Σ−1 − Σ−1 ll> Σ−1 /l> Σ−1 l; ξ,
c) compute

√

and

√



−1
l> Sx − l> Σµν
nσν



1 > −1
> −1
nσ̃ν l S x −
l Σ µν
1−c
−1

where
µν
2
σν
2
σ̃ν

= µ + Bν
h
i
2
>
>
2
> 3
= µ>
ν Σµν + c||Σ||F l Σl + (l Σµν ) + l Σ l
 

2
1
> −1
> −1
2
=
2 l Σ µν + l Σ l(1 + δ (ν))
(1 − c)3

−1
−1 > −1 > −1
with δ 2 (ν) = µ>
ν Rl µν , Rl = Σ − Σ ll Σ /l Σ l.

17

d) repeat a)-c) N times.
It is remarkable that for generating l> Sx and l> S−1 x only random variables from the standard distributions are need. Neither the data matrix X nor the sample covariance matrix S are
used.
[ F igures 1 − 8 ]
In Figures 1-4 we present the results of simulations for the asymptotic distribution that
is given in Theorem 1 while the asymptotic distribution as given in Theorem 2 is presented in
Figures 5-8 for different values of c = {0.1, 0.5, 0.8, 0.95}. The suggested asymptotic distributions
are shown as a dashed black line, while the standard normal distribution is a solid black line.
All results demonstrate a good performance of both asymptotic distributions for all considered
values of c. Even in the extreme case c = 0.95 our asymptotic results seem to produce a quite
reasonable approximation. Moreover, we observe a good robustness of our theoretical results
for different distributions of ν. Also, we observe that all asymptotic distributions are slightly
skewed to the right for the finite dimensions. This effect is even more significant in the case
of the generalized asymmetric Laplace distribution. Nevertheless, the skewness disappears with
growing dimension and sample size, i.e., the distribution becomes symmetric one and converges
to its asymptotic counterpart.

5

Summary

In this paper we introduce the family of the matrix-variate location mixture of normal distributions that generalizes a large number of the existing skew normal models. Under the MVLMN we
derive the distributions of the sample mean vector and the sample covariance matrix. Moreover,
we show that they are independently distributed. Furthermore, we derive the CLTs under the
high-dimensional asymptotic regime for the products of the (inverse) sample covariance matrix
and the sample mean vector. In the numerical study, the good finite sample performance of
both asymptotic distributions is documented.

Acknowledgement
The authors are thankful to Professor Niels Richard Hansen, the Associate Editor, and two
anonymous Reviewers for careful reading of the manuscript and for their suggestions which have
improved an earlier version of this paper.

6

Appendix

Proof of Proposition 1.

18

Proof. Straightforward but tedious calculations give
fX (Z)

= C −1

Z
Rq+

∗
∗
fNp,n (µ1>
(Z − Bν ∗ 1>
n )fNq (0,Ω) (ν )dν
n ,Σ⊗In )



Z
(2π)−(np+q)/2
1 ∗> −1 ∗
exp
−
ν
Ω
ν
2
|Ω|1/2 |Σ|n/2 Rq+




1
>
∗ > >
−1
>
∗ >
dν ∗
exp − vec Z − µ1n − Bν 1n (In ⊗ Σ) vec Z − µ1n − Bν 1n
2


1
(2π)−(np+q)/2
> >
−1
>
exp
−
vec(Z
−
µ1
)
(I
⊗
Σ)
vec(Z
−
µ1
)
C −1
n
n
n
2
|Ω|1/2 |Σ|n/2


1
> >
>
exp
vec(Z − µ1>
n ) E DEvec(Z − µ1n )
2


Z
> −1 ∗
i
1h ∗
>
dν ∗
ν − DEvec(Z − µ1>
)
D
ν
−
DEvec(Z
−
µ1
)
exp −
n
n
2
Rq+

= C −1
×
=
×
×



|F|1/2 |D|1/2
>
Φq 0; −DEvec(Z − µ1>
n ), D φpn vec(Z − µ1n ); 0, F
1/2
n/2
|Ω| |Σ|


−1
>
e
= C Φq 0; −DEvec(Z − µ1>
n ), D φpn vec(Z − µ1n ); 0, F

= C −1

> −1
, F = (In ⊗ Σ−1 − E> DE)−1 , and
where D = (nB> Σ−1 B + Ω−1 )−1 , E = 1>
n ⊗B Σ
1/2

1/2

e −1 = C −1 |F| |D| .
C
|Ω|1/2 |Σ|n/2

References
Adcock, C., Eling, M., and Loperfido, N. (2015). Skewed distributions in finance and actuarial
science: a review. The European Journal of Finance, 21:1253–1281.
Amemiya, Y. (1994). On multivariate mixed model analysis. Lecture Notes-Monograph Series,
24:83–95.
Azzalini, A. (2005).
The skew-normal distribution and related multivariate families.
Scandinavian Journal of Statistics, 32:159–188.
Azzalini, A. and Capitanio, A. (1999). Statistical applications of the multivariate skew-normal
distribution. Journal of the Royal Statistical Society: Series B, 61:579–602.
Azzalini, A. and Dalla-Valle, A. (1996). The multivariate skew-normal distribution. Biometrika,
83:715–726.
Bai, J. and Shi, S. (2011). Estimating high dimensional covariance matrices and its applications.
Annals of Economics and Finance, 12:199–215.
Bai, Z. D. and Silverstein, J. W. (2004). CLT for linear spectral statistics of large dimensional
sample covariance matrices. Annals of Probability, 32:553–605.

19

Bartoletti, S. and Loperfido, N. (2010). Modelling air polution data by the skew-normal distribution. Stochastic Environmental Research and Risk Assessment, 24:513–517.
Bodnar, T. and Gupta, A. K. (2011). Estimation of the precision matrix of multivariate elliptically contoured stable distribution. Statistics, 45:131–142.
Bodnar, T., Gupta, A. K., and Parolya, N. (2014a). On the strong convergence of the optimal
linear shrinkage estimator for large dimensional covariance matrix. Journal of Multivariate
Analysis, 132:215–228.
Bodnar, T., Gupta, A. K., and Parolya, N. (2016a). Direct shrinkage estimation of large dimensional precision matrix. Journal of Multivariate Analysis, 146:223–236.
Bodnar, T., Hautsch, N., and Parolya, N. (2016b). Consistent estimation of the high dimensional
efficient frontier. Technical report.
Bodnar, T., Mazur, S., and Okhrin, Y. (2013). On the exact and approximate distributions
of the product of a wishart matrix with a normal vector. Journal of Multivariate Analysis,
125:176–189.
Bodnar, T., Mazur, S., and Okhrin, Y. (2014b). Distribution of the product of singular wishart
matrix and normal vector. Theory of Probability and Mathematical Statistics, 91:1–14.
Bodnar, T., Okhrin, O., and Parolya, N. (2016c). Optimal Shrinkage Estimator for HighDimensional Mean Vector. ArXiv e-prints.
Bodnar, T. and Okhrin, Y. (2008). Properties of the singular, inverse and generalized inverse
partitioned wishart distributions. Journal of Multivariate Analysis, 99:2389–2405.
Bodnar, T. and Okhrin, Y. (2011). On the product of inverse wishart and normal distributions
with applications to discriminant analysis and portfolio theory. Scandinavian Journal of
Statistics, 38:311–331.
Bodnar, T. and Schmid, W. (2008). A test for the weights of the global minimum variance
portfolio in an elliptical model. Metrika, 67(2):127–143.
Bühlmann, P. and van de Geer, S. (2011). Statistics for High-Dimensional Data: Methods,
Theory and Applications. Springer Publishing Company, Incorporated, 1st edition.
Cai, T. T. and Yuan, M. (2012). Adaptive covariance matrix estimation through block thresholding. Ann. Statist., 40:2014–2042.
Cai, T. T. and Zhou, H. H. (2012). Optimal rates of convergence for sparse covariance matrix
estimation. Ann. Statist., 40:2389–2420.
Christiansen, M. and Loperfido, N. (2014). Improved approximation of the sum of random
vectors by the skew normal distribution. Journal of Applied Probability, 51(2):466–482.

20

DasGupta, A. (2008). Asymptotic Theory of Statistics and Probability. Springer Texts in
Statistics. Springer.
De Luca, G. and Loperfido, N. (2015). Modelling multivariate skewness in financial returns: a
SGARCH approach. The European Journal of Finance, 21:1113–1131.
Efron, B. (2006). Minimum volume confidence regions for a multivariate normal mean vector.
Journal of the Royal Statistical Society: Series B, 68:655–670.
Fan, J., Fan, Y., and Lv, J. (2008). High dimensional covariance matrix estimation using a
factor model. Journal of Econometrics, 147:186–197.
Fan, J., Liao, Y., and Mincheva, M. (2013). Large covariance estimation by thresholding principal orthogonal complements. Journal of the Royal Statistical Society: Series B, 75:603–680.
Gupta, A. and Nagar, D. (2000). Matrix Variate Distributions. Chapman & Hall/CRC.
Jorion, P. (1986). Bayes-stein estimation for portfolio analysis.
Quantative Analysis, 21:279–292.

Journal of Financial and

Kotsiuba, I. and Mazur, S. (2015). On the asymptotic and approximate distributions of the
product of an inverse Wishart matrix and a Gaussian random vector. Theory of Probability
and Mathematical Statistics, 93:96–105.
Kozubowski, T. J., Podgórski, K., and Rychlik, I. (2013). Multivariate generalized Laplace
distribution and related random fields. Journal of Multivariate Analysis, 113:59–72.
Liseo, B. and Loperfido, N. (2003). A Bayesian interpretation of the multivariate skew-normal
distribution. Statistics & Probability Letters, 61:395–401.
Liseo, B. and Loperfido, N. (2006). A note on reference priors for the scalar skew-normal
distribution. Journal of Statistical Planning and Inference, 136:373–389.
Loperfido, N. (2010). Canonical transformations of skew-normal variates. TEST, 19:146–165.
Lütkepohl, H. (1996). Handbook of Matrices. New York: John wiley & Sons.
Marčenko, V. A. and Pastur, L. A. (1967). Distribution of eigenvalues for some sets of random
matrices. Sbornik: Mathematics, 1:457–483.
Mathai, A. and Provost, S. B. (1992). Quadratic Forms in Random Variables. Marcel Dekker.
Muirhead, R. J. (1982). Aspects of Multivariate Statistical Theory. Wiley, New York.
Potthoff, R. F. and Roy, S. N. (1964). A generalized multivariate analysis of variance model
useful especially for growth curve problems. Biometrika, 51(3/4):313–326.
Provost, S. and Rudiuk, E. (1996). The exact distribution of indefinite quadratic forms in
noncentral normal vectors. Annals of the Institute of Statistical Mathematics, 48:381–394.

21

Rao, C. R. (1965). The theory of least squares when the parameters are stochastic and its
application to the analysis of growth curves. Biometrika, 52(3/4):447–458.
Silverstein, J. W. (1995). Strong convergence of the empirical distribution of eigenvalues of
large-dimensional random matrices. Journal of Multivariate Analysis, 55:331–339.
Stein, C. (1956). Inadmissibility of the usual estimator of the mean of a multivariate normal distribution. In Neyman, J., editor, Proceedings of the third Berkeley symposium on
mathematical and statistical probability. University of California, Berkley.
Wang, C., Pan, G., Tong, T., and Zhu, L. (2015). Shrinkage estimation of large dimensional
precision matrix using random matrix theory. Statistica Sinica, 25:993–1008.

22

0.5

0.5

0.4
0.3
0.2
0.1
-4

-2

0

2

4

-4

p = 50, n = 500, ν ∼ T N q (0, Iq ).

(b)

0

2

4

p = 100, n = 1000, ν ∼ T N q (0, Iq ).

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

-4

(c)

-2

0.5

0.5

(a)

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

-2

0

2

4

p = 50, n = 500, ν ∼ GALq (1q , Iq , 10).

-4

(d)

-2

0

2

4

p = 100, n = 1000, ν ∼ GALq (1q , Iq , 10).

Figure 1: The kernel density estimator of the asymptotic distribution as given in Theorem 1 for
c = 0.1.

23

0.5

0.5

0.4
0.3
0.2
0.1
-4

-2

0

2

4

-4

p = 250, n = 500, ν ∼ T N q (0, Iq ).

(b)

0

2

4

p = 500, n = 1000, ν ∼ T N q (0, Iq ).

Asymptotic
Standard Normal

0.0

0.0

0.1

0.1

0.2

0.2

0.3

0.3

0.4

0.4

Asymptotic
Standard Normal

-4

(c)

-2

0.5

0.5

(a)

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

-2

0

2

4

−4

p = 250, n = 500, ν ∼ GALq (1q , Iq , 10).

(d)

−2

0

2

4

p = 500, n = 1000, ν ∼ GALq (1q , Iq , 10).

Figure 2: The kernel density estimator of the asymptotic distribution as given in Theorem 1 for
c = 0.5.

24

0.5

0.5

0.3

0.3

0.2

0.2

0.1

0.1

0.0

0.0
-4

-2

0

2

4

−4

p = 400, n = 500, ν ∼ T N q (0, Iq ).

(b)

0

2

4

p = 800, n = 1000, ν ∼ T N q (0, Iq ).

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−4

(c)

−2

0.5

0.5

(a)

Asymptotic
Standard Normal

0.4

0.4

Asymptotic
Standard Normal

−2

0

2

4

−4

p = 400, n = 500, ν ∼ GALq (1q , Iq , 10).

(d)

−2

0

2

4

p = 800, n = 1000, ν ∼ GALq (1q , Iq , 10).

Figure 3: The kernel density estimator of the asymptotic distribution as given in Theorem 1 for
c = 0.8.

25

0.5

0.5

0.4
0.3
0.2
0.1
−4

−2

0

2

4

−4

p = 475, n = 500, ν ∼ T N q (0, Iq ).

(b)

0

2

4

p = 950, n = 1000, ν ∼ T N q (0, Iq ).

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−4

(c)

−2

0.5

0.5

(a)

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−2

0

2

4

−4

p = 475, n = 500, ν ∼ GALq (1q , Iq , 10).

(d)

−2

0

2

4

p = 950, n = 1000, ν ∼ GALq (1q , Iq , 10).

Figure 4: The kernel density estimator of the asymptotic distribution as given in Theorem 1 for
c = 0.95.

26

0.5

0.5

0.4
0.3
0.2
0.1
−4

−2

0

2

4

−4

p = 50, n = 500, ν ∼ T N q (0, Iq ).

(b)

0

2

4

p = 100, n = 1000, ν ∼ T N q (0, Iq ).

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−4

(c)

−2

0.5

0.5

(a)

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−2

0

2

4

p = 50, n = 500, ν ∼ GALq (1q , Iq , 10).

−4

(d)

−2

0

2

4

p = 100, n = 1000, ν ∼ GALq (1q , Iq , 10).

Figure 5: The kernel density estimator of the asymptotic distribution as given in Theorem 2 for
c = 0.1.

27

0.5

0.5

0.4
0.3
0.2
0.1
−4

−2

0

2

4

−4

p = 250, n = 500, ν ∼ T N q (0, Iq ).

(b)

0

2

4

p = 500, n = 1000, ν ∼ T N q (0, Iq ).

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−4

(c)

−2

0.5

0.5

(a)

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−2

0

2

4

−4

p = 250, n = 500, ν ∼ GALq (1q , Iq , 10).

(d)

−2

0

2

4

p = 500, n = 1000, ν ∼ GALq (1q , Iq , 10).

Figure 6: The kernel density estimator of the asymptotic distribution as given in Theorem 2 for
c = 0.5.

28

0.5

0.5

0.4
0.3
0.2
0.1
−4

−2

0

2

4

−4

p = 400, n = 500, ν ∼ T N q (0, Iq ).

(b)

0

2

4

p = 800, n = 1000, ν ∼ T N q (0, Iq ).

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−4

(c)

−2

0.5

0.5

(a)

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−2

0

2

4

−4

p = 400, n = 500, ν ∼ GALq (1q , Iq , 10).

(d)

−2

0

2

4

p = 800, n = 1000, ν ∼ GALq (1q , Iq , 10).

Figure 7: The kernel density estimator of the asymptotic distribution as given in Theorem 2 for
c = 0.8.

29

0.5

0.5

0.4
0.3
0.2
0.1
−4

−2

0

2

4

−4

p = 475, n = 500, ν ∼ T N q (0, Iq ).

(b)

0

2

4

p = 950, n = 1000, ν ∼ T N q (0, Iq ).

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−4

(c)

−2

0.5

0.5

(a)

Asymptotic
Standard Normal

0.0

0.0

0.1

0.2

0.3

0.4

Asymptotic
Standard Normal

−2

0

2

4

−4

p = 475, n = 500, ν ∼ GALq (1q , Iq , 10).

(d)

−2

0

2

4

p = 950, n = 1000, ν ∼ GALq (1q , Iq , 10).

Figure 8: The kernel density estimator of the asymptotic distribution as given in Theorem 2 for
c = 0.95.

30

