Optimal shrinkage-based portfolio selection in high dimensions

Taras Bodnara , Yarema Okhrinb and Nestor Parolyac,∗

arXiv:1611.01958v3 [q-fin.ST] 9 May 2017

c

1

a

Department of Mathematics, Stockholm University, Roslagsvägen 101, SE-10691 Stockholm, Sweden

b

Department of Statistics, University of Augsburg, Universitätsstr. 16, D-86159 Augsburg, Germany

Institute of Statistics, Leibniz University Hannover, Königsworther Platz 1, D-30167 Hannover, Germany

Abstract
In this paper we estimate the mean-variance (MV) portfolio in the high-dimensional case using
the recent results from the theory of random matrices. We construct a linear shrinkage estimator which is distribution-free and is optimal in the sense of maximizing with probability 1 the
asymptotic out-of-sample expected utility, i.e., mean-variance objective function. Its asymptotic
properties are investigated when the number of assets p together with the sample size n tend to infinity such that p/n → c ∈ (0, +∞). The results are obtained under weak assumptions imposed on
the distribution of the asset returns, namely the existence of the fourth moments is only required.
Thereafter we perform numerical and empirical studies where the small- and large-sample behavior
of the derived estimator is investigated. The suggested estimator shows significant improvements
over the naive diversification and it is robust to the deviations from normality.

JEL Classification: G11, C13, C14, C58, C65
Keywords: expected utility portfolio, large-dimensional asymptotics, covariance matrix estimation,
random matrix theory.

1

Introduction

In the seminal paper of Markowitz (1952) the author suggests to determine the optimal composition of
a portfolio of financial assets by minimizing the portfolio variance assuming that the expected portfolio
return attains some prespecified fixed value. By varying this value we obtain the whole efficient frontier
in the mean-standard deviation space. Despite of its simplicity, this approach justifies the advantages
of diversification and is a standard technique and benchmark in asset management. Equivalently (see,
Tobin (1958), Bodnar et al. (2013)) we can obtain the same portfolios by maximizing the expected
quadratic utility (EU) with the optimization problem given by
w 0 µn −

γ 0
w Σn w → max subject to w0 1p = 1 ,
2

(1.1)

where w = (ω1 , . . . , ωp )0 is the vector of portfolio weights, 1p is the p-dimensional vector of ones,
µn and Σn are the p-dimensional mean vector and the p × p covariance matrix of asset returns,
respectively. The quantity γ > 0 determines the investor’s behavior towards risk. It must be noted
1∗

Corresponding author. E-mail address: parolya@statistik.uni-hannover.de

1

that the maximization of the mean-variance objective function (1.1) is equivalent to the maximization
of the exponential utility (CARA) function under the assumption of normality of the asset returns.
In this case γ equals the investor’s absolute risk aversion coefficient (see, e.g., Pratt (1964)).
The solution of the optimization problem (1.1) is well known and it is given by
wEU = wGM V + γ −1 Qn µn ,

(1.2)

where
Qn =

Σn−1

−

0 −1
Σ−1
n 1p 1p Σn

(1.3)

10p Σ−1
n 1p

and
wGM V =

Σ−1
n 1p
0
1p Σ−1
n 1p

(1.4)

is the vector of the weights of the global minimum variance portfolio.
In practice, however, the above mentioned approach of constructing an optimal portfolio frequently
shows poor out-of-sample performance in terms of various performance measures. Even naive portfolio
strategies, e.g., equally weighted portfolio (see, DeMiguel et al. (2009)), often outperform the meanvariance strategy. One of the reasons is the estimation risk. The unknown parameters µn and Σn
have to be estimated using historical data on asset returns. This results in the ”plug-in” estimator of
the EU portfolio (1.2) which is a traditional and simple way to evaluate the portfolio in practice. This
estimator is constructed by replacing in (1.2) the mean vector µn and the covariance matrix Σn with
their sample counterparts. Okhrin and Schmid (2006) derive the expected return and the variance of
the sample portfolio weights under the assumption that the asset returns follow a multivariate normal
distribution, while Bodnar and Schmid (2011) obtain the exact finite-sample distribution. Recently,
Bodnar et al. (2016) extended these results to the case n < p.
The estimation of the parameters has a negativ impact on the performance of the asset allocation
strategy. This was noted in a series of papers with Merton (1980), Best and Grauer (1991), Chopra and
Ziemba (1993) among others. Several approaches have arisen to reduce the consequences of estimation
risk. One strand of research opted for the Bayesian framework and using appropriate priors take the
estimation risk into account already while building the portfolio. The second strand relied on the
shrinkage techniques and is related to the method exploited in this paper. A straightforward way
to improve the properties of the estimators for µn and Σn is to use the shrinkage approach (see,
Jorion (1986), Ledoit and Wolf (2004)). Alternatively, one may apply the shrinkage estimation to the
portfolio weights directly. Golosnoy and Okhrin (2007) consider the multivariate shrinkage estimator
by shrinking the portfolios with and without the riskless asset to an arbitrary static portfolio. A similar
technique is used by Frahm and Memmel (2010), who constructed a feasible shrinkage estimator for
the GMV portfolio which dominates the traditional one. At last, Bodnar et al. (2017) suggest a
shrinkage estimator for the GMV portfolio which is feasible even for the singular sample covariance
matrix.
An important issue nowadays is, however, the asset allocation for large portfolios. The sample
estimators work well only in the case when the number of assets p is fixed and substantially smaller
than the sample size n. This case is known as the standard asymptotics in statistics (see, Le Cam and
2

Lo Yang (2000)). Under this asymptotics the traditional sample estimator is a consistent estimator
for the EU portfolio. But what happens when the dimension p and the sample size n are comparable
of size, say p = 900 and n = 1000? Technically, here we are in the situation when both the number
of assets p and the sample size n tend to infinity. In the case when p/n tends to some concentration
ratio c > 0 this asymptotics is known as high-dimensional asymptotics or “Kolmogorov” asymptotics
(see, e.g., Bai and Silverstein (2010)). If c is close to one the sample covariance matrix tends to
be close to a singular one and when c > 1 it becomes singular. The sample estimator of the EU
portfolio behaves badly in this case both from the theoretical and practical points of view (see, e.g.,
El Karoui (2010)). The sample covariance matrix is very unstable and tends to under- or overestimate
the true parameters for c smaller but close to 1 (see, Bai and Shi (2011)). For c > 1 the inverse
sample covariance does not exist and the portfolio cannot be constructed in the traditional way.
Taking the above mentioned information into account the aim of the paper is to construct a feasible
and simple shrinkage estimator of the EU portfolio which is optimal in an asymptotic sense and is
additionally distribution-free. The estimator is developed using the fast growing branch of probability
theory, namely random matrix theory. The main result of this theory is proved by Marčenko and
Pastur (1967) and further extended under very general conditions by Silverstein (1995). Now it is
called Marc̆enko-Pastur equation. Its importance arises in many areas of science because it shows how
the true covariance matrix and its sample estimator are connected asymptotically. Knowing this we
can build suitable estimators for high-dimensional quantities which depend on Σn . In our case this
refers to the shrinkage intensities. Note however, that the optimal shrinkage intensity depends again
on the unknown characteristics of the asset returns. To overcome this problem we derive consistent
estimators for specific functions (quadratic and bilinear forms) of the inverse sample covariance matrix
and succeed to provide consistent estimators for the optimal shrinkage intensities.
The rest of paper is organized as follows. In the next section, we construct a shrinkage estimator
for the optimal portfolio weights obtained by shrinking the mean-variance weights to an arbitrary
target portfolio. The oracle shrinkage intensity and the corresponding feasible bona-fide estimators
for c < 1 and c > 1 are established as well. The derived results are evaluated in Section 3 in extensive
simulation and empirical studies. All proofs are moved to in the Appendix.

2

Optimal shrinkage estimator of mean-variance portfolio

Let Yn = (y1 , y2 , ..., yn ) be the p × n data matrix which consist of n vectors of the returns on p ≡ p(n)
assets. Let E(yi ) = µn and Cov(yi ) = Σn for i ∈ 1, ..., n. We assume that p/n → c ∈ (0, +∞) as
n → ∞. This type of limiting behavior is known as ”the large dimensional asymptotics” or ”the
Kolmogorov asymptotics”. In this case the traditional sample estimators perform poor or even very
poor and tend to over/underestimate the unknown parameters of the asset returns, e.g., the mean
vector and the covariance matrix.
Throughout the paper it is assumed that there exists a p × n random matrix Xn which consists
of independent and identically distributed (i.i.d.) real random variables with zero mean and unit
variance such that
1

Yn = µn 10n + Σn2 Xn .

(2.1)

It must be noted that the observation matrix Yn has the dependent rows but independent columns.
3

Broadly speaking, this means that we allow arbitrary cross-sectional correlations of the asset returns
but assume their independence over time. Although this assumption looks quite restrictive for financial
applications, there exist stronger results from random matrix theory which show that the model can be
extended to (weakly) depending variables by demanding more complicated conditions on the elements
of Yn (see, Bai and Zhou (2008)) or by controlling the number of dependent entries as dimension
increases (see, Hui and Pan (2010), Friesen et al. (2013), Wei et al. (2016)). Nevertheless, this will
only make the proofs more technical, but leave the results unchanged. For that reason we assume
independent asset returns over time only to simplify the proofs of the main theorems.
The three assumptions which are used throughout the paper are
(A1) The covariance matrix of the asset returns Σn is a nonrandom p-dimensional positive definite
matrix.
(A2) The elements of the matrix Xn have uniformly bounded 4 + ε moments for some ε > 0.
−1
0
(A3) There exist Ml , Mu ∈ (0, +∞) such that Ml ≤ 10p Σ−1
n 1p , µn Σn µn ≤ Mu .

All of these regularity assumptions are general enough to fit many real world situations. The
assumption (A1) together with (2.1) are usual for financial and statistical problems and impose no
strong restrictions. The assumption (A2) is a technical one and can be relaxed for practical purposes
(see, section with simulations). The assumption (A3) requires that the quantities which are used in
the calculations are finite. This assumption is quite general and imposes no additional constrains
neither on the mean vector µn , like its Euclidean norm is bounded, nor on the covariance matrix Σn ,
like its eigenvalues lie in the compact interval. The last point allows us to assume a factor model for
the data matrix Yn which implies that the largest eigenvalue of Σn is of order p (c.f. Fan et al. (2008),
Fan et al. (2012), Fan et al. (2013)). Finally, assumption (A3) ensures that 10p Σ−1
n µn is finite as well
which follows directly from the Cauchy-Schwarz inequality. This implies that the variance and the
expected return of the global minimum variance portfolio are uniformly bounded in p.
The sample covariance matrix is given by
Sn =

1
1
1
1 1
1
Yn (In − 1n 10n )Yn0 = Σn2 Xn (In − 1n 10n )X0n Σn2 ,
n
n
n
n

(2.2)

where the symbol In stands for the n-dimensional identity matrix. The sample mean vector becomes
ȳn =

2.1

1
1
1
Yn 1n = µn + Σn2 x̄n with x̄n = Xn 1n .
n
n

(2.3)

Oracle estimator. Case c < 1

In this section we consider the optimal shrinkage estimator for the EU portfolio weights presented in
the introduction by optimizing the shrinkage parameter α and fixing some target portfolio b.
The resulting estimator for c < 1 is given by
ŵGSE = αn ŵS + (1 − αn )b with b0 1p = 1 and b0 Σn b ≤ ∞ ,

(2.4)

where the vector ŵS is the sample estimator of the EU portfolio given in (1.2), namely
ŵS =

S−1
n 1p
+ γ −1 Q̂n ȳn
10p S−1
1
n
p

(2.5)

4

with
Q̂n = S−1
n −

0 −1
S−1
n 1p 1p Sn

10p S−1
n 1p

The target portfolio b ∈

.

(2.6)

Rp is a given nonrandom (and further random independent of Yn) vector

with b0 1p = 1 and uniformly bounded variance. No assumption is imposed on the shrinkage intensity
αn which is the object of our interest.
The aim is now to find the optimal shrinkage intensity for a given nonrandom target portfolio b
which maximizes the mean-variance objective function (1.1). The maximization problem is given by
0
U = ŵGSE
(αn )µn −

γ 0
ŵ
(αn )Σn ŵGSE (αn ) −→ max with respect to αn .
2 GSE

(2.7)

The expected utility function (2.7) can be rewritten as
U = αn ŵS0 µn + (1 − αn )b0 µn −
with respect to αn .


γ 2 0
α ŵ Σn ŵS + 2αn (1 − αn )b0 Σn ŵS + (1 − αn )2 b0 Σn b → max
2 n S
(2.8)

Next, taking the derivative of U with respect to αn and setting it equal to zero we get
∂U
= (ŵS − b)0 µn − γ (αn ŵS0 Σn ŵS + (1 − 2αn )b0 Σn ŵS − (1 − αn )b0 Σn b) .
∂αn

From the last equation it is easy to find the optimal shrinkage intensity αn∗ given by
αn∗ = γ −1

(ŵS − b)0 (µn − γΣn b)
.
(ŵS − b)0 Σn (ŵS − b)

(2.9)

To ensure that αn∗ is the unique maximizer of (2.7) the second derivative of U must be negative, which
is always fulfilled. Indeed, it follows from the positive definitiveness of the matrix Σn , namely
∂2U
= −γ(ŵS − b)0 Σn (ŵS − b) < 0 .
∂αn2

(2.10)

In the next theorem we show the asymptotic properties of the optimal shrinkage intensity αn∗ under
large dimensional asymptotics.
Theorem 2.1. Assume (A1)-(A3). Then it holds that
a.s.

|αn∗ − α∗ | −→ 0 for

p
→ c ∈ (0, 1) as n → ∞
n

with
(RGM V
α∗ = γ −1

γ −1
+ γ(Vb − VGM V ) +
s
1−c


, (2.11)

s
c
γ −1
−2
+ 1−c (Rb − RGM V ) + γ
+
+
V
b
(1 − c)3 (1 − c)3


− Rb ) 1 +


1
VGM V − 2 VGM V
1−c

1
1−c



10 Σ−1
1
n µn
and VGM V = 0 −1 are the expected return and the variance of the true
−1
0
1 Σn 1
1 Σn 1
global minimum variance portfolio, while s = µ0n Qn µn is the slope parameter of the efficient frontier.
where RGM V =

The quantities Rb = b0 µn and Vb = b0 Σn b denote the expected return and the variance of the target
5

portfolio b.
Next we assess the performance of the classical estimator of the portfolio weights ŵS and the
optimal shrinkage weights ŵGSE . As a measure of performance we consider the relative increase in
the utility of the portfolio return compared to the portfolio based on true parameters of asset returns.
The results are summarized in the following corollary.
Corollary 2.1. (a) Let UEU and US be the expected quadratic utilities for the true EU portfolio and its
traditional estimator. Then under the assumptions of Theorem 2.1, the relative loss of the traditional
estimator of the EU portfolio is given by
UEU −
RS =
UEU
p
n

for

γ
US a.s. 2

−→



1
1−c



− 1 · VGM V + γ −1 − 21 −
RGM V −

γ −1
2

1
(1−c)

+

1
2(1−c)3



·s+

γ −1
2

·

c
(1−c)

· s − γ2 VGM V

(2.12)

→ c ∈ (0, 1) as n → ∞.

(b) Let UGSE be the expected quadratic utility for optimal shrinkage estimator of the EU portfolio.
Under the assumptions of Theorem 2.1, the relative loss of the optimal shrinkage estimator is given
by
RGSE =

2.2

UEU − UGSE a.s.
−→ (α∗ )2 RS + (1 − α∗ )2 Rb
UEU

for

p
→ c ∈ (0, 1) as n → ∞.
n

(2.13)

Oracle estimator. Case c > 1.

Here, similarly as in Bodnar et al. (2017), we will use the generalized inverse of the sample covariance
matrix Sn . Particularly, we use the following generalized inverse of the sample covariance matrix Sn
S∗n = Σ−1/2
n



1
Xn X0n − x̄n x̄0n
n

+

Σn−1/2 ,

(2.14)

where 0 +0 denotes the Moore-Penrose inverse. It can be shown that S∗n is a generalized inverse of Sn
satisfying S∗n Sn S∗n = S∗n and Sn S∗n Sn = Sn . However, S∗n is not exactly equal to the Moore-Penrose
inverse because it does not satisfy the conditions (S∗n Sn )0 = S∗n Sn and (Sn S∗n )0 = Sn S∗n . In case c < 1
the generalized inverse S∗n coincides with the usual inverse S−1
n . Moreover, if Σn is a multiple of
∗
identity matrix then S∗n is equal to the Moore-Penrose inverse S+
n . In this section, Sn is used only

to determine an oracle estimator for the weights of the MV portolio. The bona fide estimator is
constructed in the next section.
Thus, the oracle estimator for c > 1 is given by
∗
ŵGSE
= αn+ ŵS ∗ + (1 − αn+ )b with b0 1p = 1 and b0 Σn b < ∞ ,

(2.15)

where the vector ŵS ∗ is the sample estimator of the EU portfolio given in (1.2), namely
ŵS ∗ =

S∗n 1p
+ γ −1 Q̂∗n ȳn
10p S∗n 1p

(2.16)

with
Q̂∗n

=

S∗n

S∗n 1p 10p S∗n
− 0 ∗
.
1p Sn 1p

(2.17)
6

Again, the shrinkage intensity αn+ is the object of our interest. In order to save place we skip the
optimization procedure for αn+ as it is only slightly different from the case c < 1. Thus, the optimal
shrinkage intensity αn+ in case c > 1 is given by
αn+ = γ −1

(ŵS ∗ − b)0 (µn − γΣn b)
.
(ŵS ∗ − b)0 Σn (ŵS ∗ − b)

(2.18)

In the next theorem we find the asymptotic equivalent quantity for αn+ in the case p/n → c ∈
(1, +∞) as n → ∞.
Theorem 2.2. Assume (A1)-(A3). Then it holds that
a.s.

αn+ − α+ −→ 0 for

p
→ c ∈ (1, +∞) as n → ∞
n

with
(RGM V
α+ = γ −1


− Rb ) 1 +

1
c(c − 1)


c2
VGM V − 2 VGM V +
(c − 1)


+ γ(Vb − VGM V ) +

γ −1
c(c−1) (Rb


− RGM V ) +

γ −1
s
c(c − 1)

γ −2
(s + c) + Vb
(c − 1)3

,

(2.19)

1
10 Σ−1
n µn
and VGM V = 0 −1 are the expected return and the variance of the true
−1
0
1 Σn 1
1 Σn 1
global minimum variance portfolio, while s = µ0n Σn µn is the slope parameter of the efficient frontier.
where RGM V =

The quantities Rb = b0 µn and Vb = b0 Σn b denote the expected return and the variance of the target
portfolio b.
Similarly as for the case c < 1 we provide here the expression for the relative losses.
Corollary 2.2. (a) Let UEU and US be the expected quadratic utilities for the true EU portfolio and its
traditional estimator. Then under the assumptions of Theorem 2.2, the relative loss of the traditional
estimator of the EU portfolio is given by
UEU −
RS =
UEU
p
n

for

γ
US a.s. 2

−→



c2
c−1



− 1 · VGM V + γ −1 − 12 −
RGM V −

γ −1
2

1
c(c−1)

+

1
2(c−1)3



· s − γ2 VGM V

·s+

γ −1
2

·

c
(1−c)3

(2.20)

→ c ∈ (1, +∞) as n → ∞.

(b) Let UGSE be the expected quadratic utility for the optimal shrinkage estimator of the EU portfolio.
Under the assumptions of Theorem 2.2, the relative loss of the optimal shrinkage estimator is given
by
RGSE =

2.3

UEU − UGSE a.s.
−→ (α+ )2 RS + (1 − α+ )2 Rbn
UEU

for

p
→ c ∈ (1, +∞) as n → ∞. (2.21)
n

Estimation of unknown parameters. Bona fide estimator

The limiting shrinkage intensity α∗ is not feasible in practice, since the quantities it depends on, i.e.
RGM V , VGM V , s, Rb and Vb , are unknown. In this subsection we derive the corresponding consistent
estimators. The results are summarized in two propositions dealing with the cases c ∈ (0, 1) and
c ∈ (1, ∞) respectively. The statements follow directly from the proofs of Theorems 2.1 and 2.2.
7

Proposition 2.1. The consistent estimators of RGM V , VGM V , s, Rb and Vb under large dimensional
asymptotics p/n → c < 1 as n → ∞ are given by
a.s.

R̂c = R̂GM V −→ RGM V
1
a.s.
V̂c =
V̂GM V −→ VGM V
1 − p/n

(2.22)
(2.23)

a.s.

ŝc = (1 − p/n)ŝ − p/n −→ s

(2.24)

a.s.

0

R̂b = b ȳn −→ Rb

(2.25)

a.s.

0

V̂b = b Sn b −→ Vb ,

(2.26)

where R̂GM V , V̂GM V and ŝ are traditional plug-in estimators.
Using Proposition 2.1 we can immediately construct a bona-fide estimator for expected utility
portfolio weights in case c < 1. It holds that
ŵBF GSE = α
b

∗

S−1
n 1p
+ γ −1 Q̂n ȳn
−1
0
1p Sn 1p

!
+ (1 − α
b∗ )b

(2.27)

with

1
γ −1
(R̂c − R̂b ) 1 +
+ γ(V̂b − V̂c ) +
ŝc
1 − p/n
1 − p/n
∗
−1


α
b =γ


1
p/n
ŝc
γ −1
V̂c − 2 V̂c + 1−p/n
+
+ V̂b
(R̂b − R̂c ) + γ −2
1 − p/n
(1 − p/n)3 (1 − p/n)3


(2.28)

where R̂c , V̂c , ŝc , R̂b and V̂b are given above in (2.1). The expression (2.27) is the optimal shrinkage
estimator for a given target portfolio b in sense that the shrinkage intensity α
b∗ tends almost surely to

its optimal value α∗ for p/n → c ∈ (0, 1) as n → ∞.

The situation is more complex in case c > 1. Here we can present only oracle estimators for the
unknown quantities RGM V , VGM V and s.
Proposition 2.2. The consistent oracle estimators of RGM V , VGM V , s under large dimensional
asymptotics p/n → c > 1 as n → ∞ are given by
a.s.

R̂co = R̂GM V −→ RGM V
1
a.s.
V̂GM V −→ VGM V
V̂co =
p/n(p/n − 1)
a.s.

ŝoc = p/n[(p/n − 1)ŝ − 1] −→ s ,
where R̂GM V , V̂GM V and ŝ are the traditional plug-in estimators.
Here, the quantities from Proposition 2.2 are not the bona fide estimators, since the matrix S∗n
depends on the unknown quantities. Thus, we propose a reasonable approximation using the applica2
tion of the Moore-Penrose inverse S+
n . It is easy to verify that in case of Σn = σ Ip for any σ > 0 the

considered approximation becomes the exact one.
Taking into account the above discussion and the result of Theorem 2.2, the bona fide estimators

8

of the quantities RGM V , VGM V and x in case c > 1 is approximated by
R̂c+ ≈
V̂c+ ≈

ȳn0 S+
n 1p
for c ∈ (1, +∞)
+
0
1p Sn 1p
1
1
for c ∈ (1, +∞)
p/n(p/n − 1) 10p S+
n 1p

ŝ+
≈ p/n[(p/n − 1)ȳn0 Q+
c
n ȳn − 1] for

(2.29)
(2.30)

c ∈ (1, +∞) .

(2.31)

The application of (2.29) leads to the bona fide optimal shrinkage estimator of the GMV portfolio in
case c > 1 expressed as
+
ŵBF
b+
GSE = α

S+
n 1p
+ γ −1 Q̂+
n ȳn
0
1p S+
n 1p

!
+ (1 − α
b+ )bn ,

(2.32)

with
γ −1
ŝ+
c
p/n(p/n
−
1)
α
b+ = γ −1
,



(p/n)2 +
γ −2
γ −1
+
+
+
ŝc + p/n + V̂b
V̂c − 2 V̂c + p/n(p/n−1) (R̂b − R̂c ) +
p/n − 1
(p/n − 1)3

(R̂c+ − R̂b ) 1 +

1
p/n(p/n − 1)



+ γ(V̂b − V̂c+ ) +

+
where R̂b and V̂b are given in (2.25) and (2.26), respectively; Q+
n = Sn −

Moore-Penrose pseudo-inverse of the sample covariance matrix Sn .

3

(2.33)

0 +
S+
n 11 Sn
and S+
n is the
10 S+
n1

Simulation and empirical studies

In this section we illustrate the performance and the advantages of the derived results using simulated
and real data. Particularly we address the estimation precision of the shrinkage coefficient and compare
the traditional estimator with the asymptotic intensity and its consistent estimator.

3.1

Simulation study

For simulation purposes we select the structure of the spectrum of the covariance matrix and of the
mean vector to make it consistent with the characteristics of the empirical data. Particularly, for
each dimension p we select the expected returns equally spread on the interval -0.3 to 0.3, capturing
a typical spectrum of daily returns measured in percent. The covariance matrix has a strong impact
on the properties of the shrinkage intensity and for this reason we consider several structures of its
spectra. Replicating the properties of empirical data we generate covariance matrices with eigenvalues
satisfying the equation λi = 0.1eδc·(i−1)/p for i = 1, ..., p. Thus the smallest eigenvalue is 0.1 and by
selecting appropriate values for c we control the largest eigenvalue and thus the condition index of
the covariance matrix. Large condition indices imply ill-conditioned covariance matrices, with the
eigenvalues very sensitive to changes of the elements. We choose δ to attain the condition indices of
150, 1000 and 8000. The target portfolio weights are set equal to the weights of the equally weighted
portfolio, i.e. bi = 1/p for i = 1, ..., p.
First we assess the general behavior of the oracle and bona-fide shrinkage intensities as functions
of c. The oracle shrinkage intensities are computed using expressions in (2.11) and (2.19) for the cases
c < 1 and c > 1, respectively. The parameters these expressions depend on are computed using the
9

true mean vector and the true covariance matrix. For the bona-fide shrinkage intensities we estimate
these parameters consistently and thus use expressions in (2.28) and (2.33). The results are illustrated
for different condition indices in Figure 1. We observe that in all cases the shrinkage intensity falls to
zero if c → 1− and increases with c if c > 1. Thus if c is small the shrinkage estimator puts higher
weight on the classical estimator of the portfolio weights, due to lower estimation risk. If c tends to 1
the system becomes unstable because of eigenvalues which are close to zero. In this case the portfolio
weights collapse to the target portfolio weights. With c further increasing the shrinkage intensity
increases too, implying that the pseudo-inverse covariance matrix can be evaluated in a proper way.
The fraction of the EU portfolio increases with c in this case. Furthermore, the discrepancy between
the oracle and bona-fide estimators naturally increases with higher condition indices. However, the
shrinkage intensities are smaller for lower condition indices. This can be explained by the specific
impact of the condition indices on the portfolio characteristics VGM V and s. Simulations showed
that VGM V increases in the condition index and s decreases. Thus their complex interaction in the
expression for the limit of αn+ results in the decreasing behavior of the shrinkage intensity as a function
of the condition index.
In a similar fashion we analyze the relative losses of portfolios based on the traditional, the oracle
and the bona-fide estimators. As a benchmark, we take the equally weighted portfolio which is also
the target portfolio of the shrinkage estimator. The relative losses as functions of c for fixed p = 150
are plotted in Figure 2. For c < 1 the losses of the traditional estimator show explosive behavior and
are comparable to the shrinkage-based estimators only for very small values of c. Thus the traditional
estimator is reliable only if the sample size is at least three times larger than the dimension. The
performance of the two shrinkage-based estimators is similar and stable over the whole range of c
excluding a small neighborhood of one. It is important that the shrinkage-based portfolio clearly
dominates the equally weighted benchmark. However, the losses are not monotonous and attain the
maximum for c between 0.5 and 0.8. For c between 0.9 and 1.2 the shrinkage-based estimators lead to
discontinuous and unstable losses, disadvising the investors from using samples sizes comparable with
the dimension. For c > 1 the situation is similar with an extremely small difference between oracleand bona-fide-based losses also in the limit. The traditional estimators attains the minimal loss for c
between two and three and increase steadily for larger intensities.
The behavior of losses as functions of the dimension p is illustrated in Figure 3. The fraction c
is set to 0.2, 0.5, 0.8 and 2 from left to right, while the condition index equals 150, 1000 and 8000
from top to bottom. From financial perspective it is important to note that the traditional estimator
outperforms the equally weighted portfolio only for small values of c (in the particular setup for
c < 0.3), thus when the classical estimators are stable and robust. For larger c the losses of the
traditional estimator increase dramatically and are comparable with the remaining estimators only
for c > 1. As before the oracle and the bona-fide estimators show similar results and clearly beat the
benchmarks. Furthermore, the performance is stable for a wide range of dimensions. Note that all
losses are decreasing functions of the dimension. Several exceptions can be found for the shrinkagebased strategies, but these can be due to poor validity of the asymptotic results for relatively small
dimensions.

10

3.2

Empirical study

The data sets used in this study covers daily data on 445 S&P500 constituents available for the
whole period from 01.01.2004 till 10.10.2014. The investor allocates his wealth to the constituents
for prespecified values of p and c with daily reallocation. For simplicity we neglect the transaction
costs in the below discussion. For illustration purposes we consider the first p = 50 and 100 assets
in alphabetic ordering and c = 0.2, 0.5, 0.8 and 2. The risk aversion coefficient equals 10. In
contrary to the simulation study the true parameters are not known for real data. Thus the oracle
shrinkage coefficient is computed using the plug-in and not the consistent estimators of the portfolio
characteristics.
The time series of estimated shrinkage coefficients (oracle and bona fide) is given in Figure 4. We
observe that for small values of c and thus a low estimation risk the shrinkage intensities are close
to one. The behavior is very stable, but mimics the periods of high and low volatility of financial
markets, particularly with spikes in 2009 and 2011. Thus high volatility on financial markets causes
higher shrinkage coefficients and thus larger fraction of the mean-variance portfolio. This can be
justified by stronger effects of diversification during turmoil periods. With larger c the certainty in
the classical portfolio diminishes. This results in lower and more volatile the shrinkage intensities.
The investors prefers equally weighted portfolio due large estimation risk. The bona fide estimator is
smaller than the oracle counterpart for moderate values of c, but becomes larger if c tends to one. For
c > 2 the shrinkage coefficients are unstable, with a better behavior for larger portfolio consisting of
100 assets.
To illustrate the economic advantages of the derived theoretical results we consider the certainty
equivalent, the Sharpe ratio, the Value-at-risk and the expected shortfall (both at 0.01 and 0.05 levels)
as performance measures of real-data trading strategies. To robustify the procedure we randomly draw
1000 portfolios of size p = 50, 100, 300 from the population of size 455. For each of the portfolios we
implement dynamic trading strategies for the last 200 days. We compute the optimal portfolio weights
using the most recent n observations obtained for c = 0.2, 0.5, 0.8 and 2. The realized returns for each
portfolio are used to compute the above mentioned performance measures. The averages over all 1000
portfolios of given size are summarized in Tables ??-2. For small c values, e.g. c = 0.2 the traditional
estimator performs well and dominates all the alternatives. This coincides with the conclusions in
the simulation study. For c = 0.5 the trading strategies based on the oracle and bona-fide estimators
become superior. For c = 0.8 the evidence is mixed if we take the benchmarks into account. For some
parameter constellations and measures the equally weighted portfolio is dominating, while for other
the shrinkage-based strategies stay superior. For c = 2 the equally weighted portfolio dominate in
the majority of the cases. The traditional estimator obviously fails and provides frequently useless
results. Additionally to the evidence from the simulation study in the previous section, we conclude
that the shrinkage-based strategies are economically advantageous too if c < 1 and not close to
one. This conclusion is robust with respect to the portfolio composition and the chosen performance
measure. For large c the shrinkage-based portfolios still improve the traditional estimator, but cannot
outperform the benchmark.

11

4

Summary

In this paper we consider the portfolio selection in high-dimensional framework. Particularly, we
assume that the number of assets p and the sample size n tend to infinity, but their ratio p/n tends
to constant c. Note that the c maybe larger than one, implying that we have more assets than
observations. Because of the large estimation risk we suggest a shrinkage-based estimator of the
portfolio weights, which shrinks the mean-variance portfolio to the equally weighted portfolio. For the
established shrinkage intensity we derive the limiting value. It depends on c and on the characteristics
of the efficient frontier. Unfortunately the result is only an oracle value and is not feasible in practice,
since it depends on unknown quantities. Thus we suggest a bona-fide estimator which overcomes this
problem. From the technical point of view we rely on the theory of random matrices and work with
the asymptotic behavior of linear and quadratic forms in mean vector and (pseudo)-inverse covariance
matrix. In extensive simulation and empirical studies we evaluate the performance of established
results with artificial and real data. If the sample size is smaller than or comparable to the dimension,
then the performance is poor and some simple alternatives, e.g. equally weighted portfolio, should be
favored.

5

Appendix A: Proofs

Here the proofs of the theorems are given. Recall that the sample mean vector and the sample
covariance matrix are given by
ȳn =

1
1
1
Yn 1n = µn + Σn2 x̄n with x̄n = Xn 1n
n
n

(5.1)

Sn =

1
1
1
1
1
1
Yn (I − 110 )Yn0 = Σn2 Vn Σn2 with Vn = Xn (I − 110 )X0n ,
n
n
n
n

(5.2)

and

respectively. Later on, we also make use of Ṽn defined by
Ṽn =

1
Xn X0n
n

(5.3)

and the formula for the 1-rank update of usual inverse given by (c.f., Horn and Johnsohn (1985))
Vn−1 = (Ṽn − x̄n x̄0n )−1 = Ṽn−1 +

Ṽn−1 x̄n x̄0n Ṽn−1
1 − x̄0n Ṽn−1 x̄n

(5.4)

as well as the formula for the 1-rank update of Moore-Penrose inverse (see,Meyer (1973)) expressed as
Vn+ =



Ṽ0 − x̄n x̄0n

= Ṽn+ −

+

Ṽn+ x̄n x̄0n (Ṽn+ )2 + (Ṽn+ )2 x̄n x̄0n (Ṽn+ )
x̄0n (Ṽn+ )3 x̄n
+
Ṽn+ x̄n x̄0n Ṽn+ .
x̄0n (Ṽn+ )2 x̄n
(x̄0n (Ṽn+ )2 x̄n )2

(5.5)

First, we present an important lemma which is a special case of Theorem 1 in Rubio and Mestre
(2011).

12

Lemma 5.1. Assume (A2). Let a nonrandom p × p-dimensional matrix Θp and a nonrandom n × ndimensional matrix Θn possess a uniformly bounded trace norms (sum of singular values). Then it
holds that


a.s.
tr Θp (Ṽn − zIp )−1 − m(z)tr (Θp ) −→ 0

a.s.
tr Θn (1/nX0n Xn − zIn )−1 − m(z)tr (Θn ) −→ 0

(5.6)
(5.7)

for p/n −→ c ∈ (0, +∞) as n → ∞, where
m(z) = (x(z) − z)−1 and m(z) = −

1−c
+ cm(z)
z

(5.8)

with
x(z) =


p
1
1 − c + z + (1 − c + z)2 − 4z .
2

(5.9)

Proof of Lemma 5.1: The application of Theorem 1 in Rubio and Mestre (2011) leads to (5.6) where
x(z) is a unique solution in

C+ of the following equation

1 − x(z)
c
=
.
x(z)
x(z) − z

(5.10)

The two solutions of (5.10) are given by
x1,2 (z) =


p
1
1 − c + z ± (1 − c + z)2 − 4z .
2

(5.11)

In order to decide which of two solutions is feasible, we note that x1,2 (z) is the Stieltjes transform
√
with a positive imaginary part. Thus, without loss of generality, we can take z = 1 + c + i2 c and get
 

n
√ 
√ o √ 
√ 
√
√
1
Im{x1,2 (z)} = Im
2 + i2 c ± i2 2c
= Im 1 + i c(1 ± 2) = c 1 ± 2 , (5.12)
2
which is positive only if the sign ” + ” is chosen. Hence, the solution is given by
x(z) =


p
1
1 − c + z + (1 − c + z)2 − 4z .
2

(5.13)

The second assertion of the lemma follows directly from Bai and Silverstein (2010).

Second, we will need the following technical lemmas.
Lemma 5.2. Assume (A2). Let θ and ξ be universal nonrandom vectors with bounded Euclidean

13

norms. Then it holds that
ξ 0 Ṽn−1 θ − (1 − c)−1 ξ 0 θ

a.s.

−→ 0 ,

(5.14)

a.s.

x̄0n Ṽn−1 x̄n −→ c ,
x̄0n Ṽn−1 θ
ξ 0 Ṽn−2 θ − (1 − c)−3 ξ 0 θ

(5.15)

a.s.

−→ 0 ,

(5.16)

a.s.

−→ 0 ,
a.s.

x̄0n Ṽn−2 x̄n −→

(5.17)

c
,
(1 − c)

(5.18)

a.s.

x̄0n Ṽn−2 θ −→ 0

(5.19)

for p/n −→ c ∈ (0, 1) as n → ∞ .
Proof of Lemma 5.2: Since the trace norm of θξ 0 is uniformly bounded, i.e.
||θξ 0 ||tr ≤

√

θ0 θ

p
ξ 0 ξ < ∞,

we get from Lemma 5.1 that
a.s.

|tr((Ṽn − zIp )−1 θξ 0 ) − m(z)tr(θξ 0 )| −→ 0 for p/n → c < 1 as n → ∞
Furthermore, the application of m(z) → (1 − c)−1 as z → 0 leads to
a.s.

|ξ 0 Ṽn−1 θ − (1 − c)−1 ξ 0 θ| −→ 0 for p/n → c < 1 as n → ∞ ,
which proves (5.14).
For deriving (5.15) we consider
"

#
−1

0
1
1
1
1
1
n
n
√ Xn
x̄0n Ṽn−1 x̄n = lim tr √ X0n
Xn X0n − zIp
n
n
n
n
z→0+
"


−1 
#
1n 10n
1 0
1n 10n
= lim tr
+ ztr
X Xn − zIn
,
n
n n
n
z→0+
where the last equality follows from the Woodbury formula (e.g., Horn and Johnsohn (1985)). The
application of Lemma 5.1 leads to

tr

1n 10n
n

"


+ ztr

1 0
X Xn − zIn
n n

−1 

1n 10n
n

#

a.s.

−→ [1 + (c − 1) + czm(z)] tr



1n 10n
n



for p/n −→ c < 1 as n → ∞ where m(z) is given by (5.8). Setting z → 0+ and taking into account
1
we get
lim m(z) =
+
1−c
z→0
a.s.

x̄0n Ṽn−1 x̄n −→ 1 + c − 1 = c for

p
−→ c ∈ (0, 1) as n → ∞ .
n

The result (5.16) was derived in Pan (2014) (see, p. 673 of this reference).
Next, we prove (5.17). It holds that


−1
∂
0 −2
0
ξ Ṽn θ =
tr Ṽn − zIp
θξ
∂z

=
z=0

14

∂
ζn (z)
∂z

z=0

where ζn (z) = tr




−1
Ṽn − zI
θξ 0 . From Lemma 5.1 ζn (z) tends a.s. to m(z)ξ 0 θ as n → ∞.

Furthermore,
∂
m(z)
∂z

=
z=0

∂
1
∂z x(z) − z

=−
z=0

x0 (z) − 1
(x(z) − z)2

=−

1
2



1− √

z=0

1+c−z
(1−c+z)2 −4z

(x(z) −


−1
=

z)2

1
. (5.20)
(1 − c)3

z=0

Consequently,
a.s.

−3 0 −1
|ξ 0 S−2
n θ − (1 − c) ξ Σn θ| −→ 0 for p/n → c < 1 as n → ∞ .

Let ηn (z) =

x̄0n (Ṽn

x̄0n Ṽn−2 x̄n =

−

zI)−1 x̄

∂
ηn (z)
∂z


n

and Θn =


1n 10n
. Then
n

,
z=0

where
#
−1
1
1
0
√ Xn Θn
Xn Xn − zIp
n
n

 a.s.
= tr(Θn ) + ztr (1/nX0n Xn − zIn )−1 Θn −→ 1 + zm(z) = c + czm(z)
"

1
ηn (z) = tr √ X0n
n

for



p
→ c ∈ (0, 1) as n → ∞. Hence,
n
a.s.

x̄0n Ṽn−2 x̄n −→ cm(0) + cz

∂
m(z)
∂z

=
z=0

c
1−c

p
→ c ∈ (0, 1) as n → ∞.
n
Finally, we get

for

x̄0n Ṽn−2 θ
for

 
−1 
∂
0
tr x̄n Ṽn − zIp
=
θ
∂z

a.s.

−→ 0
z=0

p
→ c ∈ (0, 1) as n → ∞.
n

Lemma 5.3. Assume (A2). Let θ and ξ be universal nonrandom vectors with bounded Euclidean
norms. Then it holds that
ξ 0 Vn−1 θ − (1 − c)−1 ξ 0 θ

a.s.

−→ 0 ,

(5.21)

c
a.s.
x̄0n Vn−1 x̄n −→
,
1−c
a.s.
x̄0n Vn−1 θ −→ 0 ,
ξ 0 Vn−2 θ − (1 − c)−3 ξ 0 θ
x̄0n Vn−2 x̄n

(5.22)
(5.23)

a.s.

−→ 0 ,
a.s.

−→

(5.24)

c
,
(1 − c)3

(5.25)

a.s.

x̄0n Vn−2 θ −→ 0

(5.26)
15

for p/n −→ c ∈ (0, 1) as n → ∞.
Proof of Lemma 5.3: From (5.4) we obtain
ξ 0 Vn−1 θ = ξ 0 Ṽn−1 θ +

ξ 0 Ṽn−1 x̄n x̄0n Ṽn−1 θ a.s.
−→ (1 − c)−1 ξ 0 θ
1 − x̄0n Ṽn−1 x̄n

for p/n −→ c ∈ (0, 1) as n → ∞ following (5.14)-(5.16). Similarly, we get (5.22) and (5.23).
In case of (5.23), we get
ξ 0 Vn−2 θ = ξ 0 Ṽn−2 θ +
+ x̄0n Ṽn−2 x̄n

ξ 0 Ṽn−2 x̄n x̄0n Ṽn−1 θ ξ 0 Ṽn−1 x̄n x̄0n Ṽn−2 θ
+
1 − x̄0n Ṽn−1 x̄n
1 − x̄0n Ṽn−1 x̄n

ξ 0 Ṽn−1 x̄n x̄0n Ṽn−1 θ a.s.
−→ (1 − c)−1 ξ 0 θ
(1 − x̄0n Ṽn−1 x̄n )2

for p/n −→ c ∈ (0, 1) as n → ∞. Similarly,
x̄0n Vn−2 x̄n = x̄0n Ṽn−2 x̄n +
+ x̄0n Ṽn−2 x̄n

x̄0n Ṽn−2 x̄n x̄0n Ṽn−1 x̄n x̄0n Ṽn−1 x̄n x̄0n Ṽn−2 x̄n
+
1 − x̄0n Ṽn−1 x̄n
1 − x̄0n Ṽn−1 x̄n

x̄0n Ṽn−1 x̄n x̄0n Ṽn−1 x̄n
x̄0n Ṽn−2 x̄n
c
a.s.
=
−→
−1
−1
0
2
0
2
(1
−
c)3
(1 − x̄n Ṽn x̄n )
(1 − x̄n Ṽn x̄n )

and
x̄0n Vn−2 θ = x̄0n Ṽn−2 θ +
+ x̄0n Ṽn−2 x̄n

x̄0n Ṽn−2 x̄n x̄0n Ṽn−1 θ x̄0n Ṽn−1 x̄n x̄0n Ṽn−2 θ
+
1 − x̄0n Ṽn−1 x̄n
1 − x̄0n Ṽn−1 x̄n

x̄0n Ṽn−1 x̄n x̄0n Ṽn−1 θ a.s.
−→ 0
(1 − x̄0n Ṽn−1 x̄n )2

for p/n −→ c ∈ (0, 1) as n → ∞.

Lemma 5.4. Assume (A2). Let θ and ξ be universal nonrandom vectors with bounded Euclidean
norms and let Pn = Vn−1 −

−1
−1
Vn
ηη 0 Vn
−1
η 0 Vn
η

where η is a universal nonrandom vectors with bounded Eu-

clidean norm. Then it holds that


ξ 0 ηη 0 θ
a.s.
0
−1
0
ξ Pn θ −→ (1 − c)
ξθ−
,
η0η
c
a.s.
x̄0n Pn x̄n −→
,
1−c
a.s.
x̄0n Pn θ −→ 0 ,


ξ 0 ηη 0 θ
a.s.
0 2
−3
0
ξ Pn θ −→ (1 − c)
ξθ−
,
η0η
c
a.s.
x̄0n P2n x̄n −→
,
(1 − c)3

(5.27)
(5.28)
(5.29)
(5.30)
(5.31)

a.s.

x̄0n P2n θ −→ 0

(5.32)

for p/n −→ c ∈ (0, 1) as n → ∞.
Proof of Lemma 5.4: It holds that
0

ξ Pn θ = ξ

0

Vn−1 θ



ξ 0 ηη 0 θ
ξ 0 Vn−1 ηη 0 Vn−1 θ a.s.
−1
0
−
−→ (1 − c)
ξθ−
η0η
η 0 Vn−1 η
16

for p/n −→ c ∈ (0, 1) as n → ∞ following (5.21). Similarly, we get
x̄0n Pn x̄n = x̄0n Vn−1 x̄n −

x̄0n Vn−1 ηη 0 Vn−1 x̄n a.s. c
−→
1−c
η 0 Vn−1 η

and
x̄0n Pn θ = x̄0n Vn−1 θ −

x̄0n Vn−1 ηη 0 Vn−1 θ a.s.
−→ 0
η 0 Vn−1 η

for p/n −→ c ∈ (0, 1) as n → ∞.
The rest of the proof follows from the equality
P2n = Vn−2 −

−1
0 −1
Vn−2 ηη 0 Vn−1 Vn−1 ηη 0 Vn−2
0 −2 Vn ηη Vn
−
+
η
V
η
n
η 0 Vn−1 η
η 0 Vn−1 η
(η 0 Vn−1 η)2

and Lemma 5.3.

Proof of Theorem 2.1: The optimal shrinkage intensity can be rewritten in the following way
∗
αn

=

=

γ −1

ŵS0 (µn − γΣn b) − b0 (µn − γΣn b)
ŵS0 Σn ŵS − 2b0 Σn ŵS + b0 Σn b

γ −1

10 S−1
n (µn − γΣn b)
+ γ −1 ȳn0 Q̂n (µn − γΣn b) − b0 (µn − γΣn b)
10 S−1
n 1
−1
10 S−1
ȳ0 Q̂ Σ S−1 1
b0 Σn S−1
n Σn Sn 1
n 1
−1 n n n n
−2 ȳ0 Q̂ Σ Q̂ ȳ − 2
+
2γ
+
γ
− 2γ −1 b0 Σn Q̂n ȳn + b0 Σn b
n n n n
n
−1
−1
(10 Sn 1)2
10 Sn 1
10 S−1
n 1
1/2

1/2

(5.33)

1/2

and ȳn = µn + Σn x̄n and Sn = Σn Vn Σn .
−1/2

From Assumption (A3), we get that the following vectors Σn
and

1/2
Σn b

−1/2

1, Σn

−1/2

(µn − γΣn b), Σn

possess bounded Euclidean norms. As a result, the application of Lemma 5.2 leads to
a.s.

0 −1/2 −1 −1/2
10 S−1
Vn Σn 1 −→ (1 − c)−1 10 Σ−1
n 1 = 1 Σn
n 1,
a.s.

0 −1/2 −1 −1/2
10 S−1
Vn Σn (µn − γΣn b) −→ (1 − c)−1 10 Σ−1
n (µn − γΣn b) = 1 Σn
n (µn − γΣn b),
a.s.

0
−1/2 −1 −1/2
−1
b0 Σn S−1
Vn Σn 1 −→ (1 − c)−1 b0 Σn Σ−1
n 1 = b Σn Σn
n 1 = (1 − c) ,
a.s.

−1
0 −1/2 −2 −1/2
10 S−1
Vn Σn 1 −→ (1 − c)−3 10 Σ−1
n Σn Sn 1 = 1 Σn
n 1

for p/n −→ c ∈ (0, 1) as n → ∞.
Finally, from Lemma 5.3 and 5.4 as well as by using the equalities
Q̂n = Σ−1/2
Pn Σ−1/2
and Pn Vn−1 = Vn−2 −
n
n

Vn−1 ηη 0 Vn−2
η 0 Vn−1 η

17

µn ,

−1/2

with η = Σn

1 we obtain

ȳn0 Q̂n (µn − γΣn b) = µ0n Σ−1/2
Pn Σ−1/2
(µn − γΣn b) + x̄0n Pn Σ−1/2
(µn − γΣn b)
n
n
n
a.s.

−→ (1 − c)−1 µ0n Qn (µn − γΣn b),
0
−1/2
ȳn0 Q̂n Σn S−1
Pn Vn−1 Σ−1/2
1 + x̄0n Pn Vn−1 Σ−1/2
1
n 1 = µn Σn
n
n
a.s.

−→ (1 − c)−3 µ0n Qn 1 = 0,
ȳn0 Q̂n Σn Q̂n ȳn = x̄0n P2n x̄n + 2µ0n Σn−1/2 P2n x̄n + µ0n Σn−1/2 P2n Σn−1/2 µn
c
a.s.
+ (1 − c)−3 µ0n Qn µn ,
−→
(1 − c)3
0 1/2
−1/2
b0 Σn Q̂n ȳn = b0 Σ1/2
µn
n Pn x̄n + b Σn Pn Σn


−1
0
1 Σn µn
a.s.
0
−1
b µn − 0 −1
−→ (1 − c)
1 Σn 1

for p/n −→ c ∈ (0, 1) as n → ∞.
a.s.
Substituting the above results into the expression of the shrinkage intensity, we get αn∗ −→ α∗ ,
where

α∗ = γ −1

10 Σ−1
γ −1 0
n (µn − γΣn b)
+
µ Qn (µn − γΣn b) − b0 (µn − γΣn b)
−1
0
1−c n
1 Σn 1




10 Σ−1
1
c
1
1
1
n µn
γ −1
0µ −
0 Q µ +
−2
−
2
b
µ
+
γ
−
2
+ b0 Σn b
n n
n
1−c
1 − c 10 Σ−1
(1 − c)3 n
(1 − c)3
10 Σ−1
10 Σ−1
n 1
n 1
n 1
(5.34)

10 Σ−1
1
n µn
,
R
=
, s = µ0n Qn µn , Rb = b0 µn and making some
GM
V
−1
0
0
1 Σn 1
1 Σ−1
n 1
technical manipulations we get the statement of Theorem 2.1.
Denoting now VGM V =

For the proof of Theorem 2.2 we need several results about the properties of Moore-Penrose inverse
which are summarized in the following three lemmas.
Lemma 5.5. Assume (A2). Let θ and ξ be universal nonrandom vectors with bounded Euclidean
norms. Then it holds that
a.s.

ξ 0 Ṽn+ θ −→ c−1 (c − 1)−1 ξ 0 θ ,
ξ

0

(Ṽn+ )2 θ

x̄0 Ṽn+ x̄n
x̄0 (Ṽn+ )2 x̄n
x̄0 (Ṽn+ )3 x̄n
x̄0 (Ṽn+ )4 x̄n

a.s.

(5.35)

−3 0

−→ (c − 1)

ξ θ,

(5.36)

=

1,
1
a.s.
−→
,
c−1
c
a.s.
,
−→
(c − 1)3
c(c + 1)
a.s.
−→
,
(c − 1)5

(5.37)
(5.38)
(5.39)
(5.40)

a.s.

x̄0 Ṽn+ θ −→ 0 ,

(5.41)

a.s.

x̄0 (Ṽn+ )2 θ −→ 0 ,

(5.42)

a.s.

x̄0 (Ṽn+ )3 θ −→ 0

(5.43)

18

for p/n −→ c ∈ (1, +∞) as n → ∞ .
Proof of Lemma 5.5: It holds that


+

Ṽ =

1
Xn X0n
n

+

1
= √ Xn
n



1 0
X Xn
n n

−2

1
√ X0n
n

and, similarly,
1
(Ṽ ) = √ Xn
n
+ i



1 0
X Xn
n n

−(i+1)

1
√ X0n for i = 2, 3, 4.
n

Let Θ = θξ 0 . It holds that
"

−2
1 0
1
0 +
ξ Ṽn θ = tr √ Xn
X Xn
n n
n
"
−3

1
1 0
0
+ 2
X Xn
ξ (Ṽn ) θ = tr √ Xn
n n
n

#
#
"

−1
1 0
1 0
∂
1 0
1
√ Xn Θ =
√ Xn Θ
,
tr √ Xn
X Xn − zIn
∂z
n n
n
n
n
z=0
"
#
#
−1

1 0
1
1 0
1 0
1 ∂2
√ Xn Θ =
√ Xn Θ
tr √ Xn
X Xn − zIn
2 ∂z 2
n n
n
n
n

z=0

The application of Woodbury formula (matrix inversion lemma, see, e.g., Horn and Johnsohn
(1985)),
1
√ Xn
n



1 0
X Xn − zIn
n n

−1

1
√ X0n = Ip + z
n



1
Xn X0n − zIp
n

−1
(5.44)

leads to
ξ

0

Ṽn+ θ

=

ξ 0 (Ṽn+ )2 θ =

−1 #
1
0
Xn Xn − zIp
Θ
,
n
z=0
"
−1 #
2
1 ∂
1
ztr
Xn X0n − zIp
Θ
2
2 ∂z
n
∂
ztr
∂z

"

.

z=0

From the proof of Lemma 5.2 we know that the matrix Θ possesses the bounded trace norm. Then
the application of Lemma 5.1 leads to
a.s.

ξ 0 Ṽn+ θ −→
a.s.

ξ 0 (Ṽn+ )2 θ −→

z
∂
∂z x(z) − z

ξ0 θ ,
z=0

1 ∂2
z
2
2 ∂z x(z) − z

ξ0 θ
z=0

for p/n → c > 1 as n → ∞, where x(z) is given in (5.9).
Let us make the following notations
θ(z) =

z
x(z) − z

and φ(z) =

x(z) − zx0 (z)
.
z2

Then the first and the second derivatives of θ(z) are given by
θ0 (z) = θ2 (z)φ(z)

00

and θ (z) = 2θ(z)θ0 (z)φ(z) + θ2 (z)φ0 (z) .

19

(5.45)

.

Using L’Hopital’s rule, we get
θ(0) = lim θ(z) = lim
z→0+

z→0+

φ(0) = lim φ(z) = lim
z→0+

z→0+

z
1
c−1
1

=−
= lim
= 
, (5.46)
0
+
1
1+c
x(z) − z z→0 (x (z) − 1)
c
1−
−1
2
|1 − c|

c
x(z) − zx0 (z)
1
1
−2c
00
=
= − lim x (z) = − lim
,
2
3/2
2
+
+
z
2 z→0
2 z→0 ((1 − c + z) − 4z)
(c − 1)3
(5.47)

and
00

2(x(z) − zx0 (z)) + z 2 x (z)
lim φ (z) = − lim
z2
z→0+
z→0+
00
2φ(z) + x (z)
000
= − lim
= − lim (2φ0 (z) + x (z)) ,
z
z→0+
z→0+
0

(5.48)

which implies
φ0 (0) = lim φ0 (z) = −
z→0+

1
6c(z − c − 1)
2c(c + 1)
1
000
=
lim x (z) = − lim
.
3 z→0+
3 z→0+ ((1 − c + z)2 − 4z)5/2
(c − 1)5

Combining (5.45), (5.46), (5.47), and (5.49), we get
0

0

θ (0) = lim θ (z) = θ2 (0)φ(0) =
z→0+

1
c(c − 1)

and
00

00

θ (0) = lim θ (z) = 2θ3 (0)φ2 (0) + θ2 (0)φ0 (0) =
z→0+

2
.
(c − 1)3

Hence,
a.s.

ξ 0 Ṽn+ θ −→
a.s.

ξ 0 (Ṽn+ )2 θ −→

1
ξ 0 Σ−1
n θ for p/n → c > 1 as n → ∞,
c(c − 1)
1
ξ 0 Σ−1
n θ for p/n → c > 1 as n → ∞ .
(c − 1)3

Taking into account that
x̄0n Ṽn+ x̄n =

1 0 0
1
1 X Xn (X0n Xn )−2 X0n Xn 1n = 10n 1n = 1 .
n n n
n

we get (5.37). Similarly, using
x̄0n (Ṽn+ )i x̄n = 1/n10n (1/nX0n Xn )−(i−1) 1n for i = 2, 3, 4
we get
1/n10n (1/nX0n Xn )−1 1n =
1/n10n (1/nX0n Xn )−2 1n =
1/n10n (1/nX0n Xn )−2 1n =

a.s.

lim tr[(1/nX0n Xn − zI)−1 Θ] −→ m(0) ,

z→0+

∂
a.s.
tr[(1/nX0n Xn − zI)−1 Θ] −→ m0 (0) ,
∂z
1
∂2
a.s. 1
lim
tr[(1/nX0n Xn − zI)−1 Θ] −→ m00 (0)
2
+
2 z→0 ∂z
2
lim

z→0+

20

(5.49)

for p/n → c > 1 as n → ∞, where Θn = 1/n1n 10n .
The application of the equality (c.f., Bai and Silverstein (2010))

m(z) = − z −

c
1 + m(z)

−1
.

leads to2
1
m(0) ≡ lim m(z) = (1 + m(0)) ,
+
c
z→0
which is linear in m(0) and results
m(0) =

1
.
c−1

(5.50)

For the next one we investigate the first derivative of m(z) with respect to z, namely


c
m (z) = z −
1 + m(z)
0

−2 
1+

cm0 (z)
(1 + m(z))2


.

Now for z → 0+ we obtain
m0 (0) =



c
1 + m(0)

−2 
1+

cm0 (0)
(1 + m(0))2




=

1 + m(0)
c

2
+

m0 (0)
,
c

which is again a linear equation in m0 (0). Thus, using (5.50) we have


1+


(1 − c−1 )m0 (0) = 

1 2
c−1

c

As a result, we get
c
.
(c − 1)3

m0 (0) =

Finally, we calculate

00
m (z) = −2 z −

−3 
2
c
cm0 (z)
1+
1 + m(z)
(1 + m(z))2

−2  00

cm (z)((1 + m(z))2 ) − 2c(m0 (z))2 (1 + m(z))
c
+ z−
1 + m(z)
(1 + m(z))4

which leads to
(m0 (0))2
m (0) = 2
+ (m(0))2
m(0)
00



cm00 (0)
2c(m0 (0))2
−
(1 + m(0))2 (1 + m(0))3



2
Note that m(z) is bounded as z → 0+ because for any contour C = {x + iy|y ∈ [0, y0 ]} the Stieltjes transform m(z)
is analytic on C and, thus, there exists δ > 0 s.t. sup |m(k) (z)| < δ k+1 k! (see, e.g., Bai and Silverstein (2004), p. 585).
z∈C

21

and, consequently,
m00 (0) = 2

c(c + 1)
.
(c − 1)5

As a result, we obtain (5.38)-(5.40).
For (5.41) we consider
"
x̄0 Ṽn+ θ = tr Ip + z



1
Xn X0n − zIp
n

−1

#


−1
θx̄0 = x̄0 θ + zx̄0 Ṽn − zIp
θ.


−1
Because of (5.16), it holds that x̄0 Ṽn − zIp
θ is uniformly bounded as z −→ 0. Moreover,
a.s.

x̄0 θ −→ 0 as p −→ ∞ following Kolmogorov’s strong law of large numbers (c.f., Sen and Singer (1993,
a.s.

Theorem 2.3.10), since θ has a bounded Euclidean norm. Hence, x̄0 Ṽn+ θ −→ 0 for p/n −→ c ∈
(1, +∞) as n → ∞.
Finally, in the case of (5.42) and (5.43), we get
#
"
−1
∂
1
a.s.
x̄0 (Ṽn+ )2 θ =
θx̄0
−→ 0 ,
ztr
Xn X0n − zIp
∂z
n
z=0
"
#
−1
2
1
1
∂
a.s.
x̄0 (Ṽn+ )3 θ =
−→ 0.
ztr
Xn X0n − zIp
θx̄0
2
2 ∂z
n
z=0

for p/n −→ c ∈ (1, +∞) as n → ∞.

Lemma 5.6. Assume (A2). Let θ and ξ be universal nonrandom vectors with bounded Euclidean
norms. Then it holds that
a.s.

ξ 0 Vn+ θ −→ c−1 (c − 1)−1 ξ 0 θ ,
1
a.s.
x̄0n Vn+ x̄n −→
,
c−1
a.s.
x̄0n Vn+ θ −→ 0 ,

(5.51)
(5.52)
(5.53)

a.s.

ξ 0 (Vn+ )2 θ −→ (c − 1)−3 ξ 0 θ ,
c
a.s.
,
x̄0n (Vn+ )2 x̄n −→
(c − 1)3

(5.54)
(5.55)

a.s.

x̄0n (Vn+ )2 θ −→ 0

(5.56)

for p/n −→ c ∈ (1, +∞) as n → ∞.
Proof of Lemma 5.6: From (5.5) we get
ξ 0 Vn+ θ = ξ 0 Ṽn+ θ −
+

ξ 0 Ṽn+ x̄n x̄0n (Ṽn+ )2 θ + ξ 0 (Ṽn+ )2 x̄n x̄0n Ṽn+ θ
x̄0n (Ṽn+ )2 x̄n

x̄0n (Ṽn+ )3 x̄n 0 +
a.s.
ξ Ṽn x̄n x̄0n Ṽn+ θ −→ c−1 (c − 1)−1 ξ 0 θ
+ 2
0
2
(x̄n (Ṽn ) x̄n )

22

for p/n −→ c ∈ (1, +∞) as n → ∞ following (5.35)-(5.37). Similarly, we get
x̄0n Vn+ x̄n = x̄0n Ṽn+ x̄n −
+

x̄0n Ṽn+ x̄n x̄0n (Ṽn+ )2 x̄n + x̄0n (Ṽn+ )2 x̄n x̄0n Ṽn+ x̄n
x̄0n (Ṽn+ )2 x̄n

x̄0n (Ṽn+ )3 x̄n 0 +
1
a.s.
x̄n Ṽn x̄n x̄0n Ṽn+ x̄n −→
+ 2
0
2
c−1
(x̄n (Ṽn ) x̄n )

and
x̄0n Vn+ θ = x̄0n Ṽn+ θ −
+

x̄0n Ṽn+ x̄n x̄0n (Ṽn+ )2 θ + x̄0n (Ṽn+ )2 x̄n x̄0n Ṽn+ θ
x̄0n (Ṽn+ )2 x̄n

x̄0n (Ṽn+ )3 x̄n 0 +
a.s.
x̄n Ṽn x̄n x̄0n Ṽn+ θ −→ 0
(x̄0n (Ṽn+ )2 x̄n )2

for p/n −→ c ∈ (1, +∞) as n → ∞.
Now, we consider the equality
h
i2
(Ṽn − x̄n x̄0n )+

2
x̄0n (Ṽn+ )3 x̄n
Ṽn+ x̄n x̄0n (Ṽn+ )2 + (Ṽn+ )2 x̄n x̄0n (Ṽn+ )
+
0
+
+
Ṽ
x̄
x̄
Ṽ
n
n
n
n
x̄0n (Ṽn+ )2 x̄n
(x̄0n (Ṽn+ )2 x̄n )2
2
 0
(Ṽn+ )2 x̄n x̄0n Ṽn+ + Ṽn+ x̄n x̄0n (Ṽn+ )2
x̄n (Ṽn+ )3 x̄n
+
0
+
Ṽ
x̄
x̄
Ṽ
−
(Ṽn+ )2 +
n n n n
(x̄0n (Ṽn+ )2 x̄n )2
x̄0n (Ṽn+ )2 x̄n


=
=

Ṽn+ −

−

2(Ṽn+ )2 x̄n x̄0n (Ṽn+ )2 + Ṽn+ x̄n x̄0n (Ṽn+ )3 + (Ṽn+ )3 x̄n x̄0n Ṽn+
x̄0n (Ṽn+ )2 x̄n

+

x̄0n (Ṽn+ )3 x̄n
(Ṽn+ x̄n x̄0n (Ṽn+ )2 + (Ṽn+ )2 x̄n x̄0n Ṽn+ )
(x̄0n (Ṽn+ )2 x̄n )2

=

(Ṽn+ )2 +

+
−

(x̄0 (Ṽ+ )3 x̄n )2 +
x̄0n (Ṽn+ )4 x̄n
Ṽn+ x̄n x̄0n Ṽn+ − n n+
Ṽn x̄n x̄0n Ṽn+
+ 2
0
2
(x̄n (Ṽn ) x̄n )
(x̄0n (Ṽn )2 x̄n )3
i
x̄0n (Ṽn+ )3 x̄n h + 2
0
+
+
0
+ 2
(
Ṽ
)
x̄
x̄
Ṽ
+
Ṽ
x̄
x̄
(
Ṽ
)
n
n
n
n
n
n
n
n
(x̄0n (Ṽn+ )2 x̄n )2
(Ṽn+ )2 x̄n x̄0n (Ṽn+ )2 + Ṽn+ x̄n x̄0n (Ṽn+ )3 + (Ṽn+ )3 x̄n x̄0n Ṽn+
x̄0n (Ṽn+ )2 x̄n

Hence,
x̄0n (Ṽn+ )4 x̄n 0 +
(x̄0n (Ṽn+ )3 x̄n )2 0 +
0
+
ξ
Ṽ
x̄
x̄
Ṽ
θ
−
ξ Ṽn x̄n x̄0n Ṽn+ θ
n
n
n n
(x̄0n (Ṽn+ )2 x̄n )2
(x̄0n (Ṽn+ )2 x̄n )3
i
x̄0n (Ṽn+ )3 x̄n h 0 + 2
0
+
0 +
0
+ 2
ξ
(
Ṽ
)
x̄
x̄
Ṽ
θ
+
ξ
Ṽ
x̄
x̄
(
Ṽ
)
θ
n n n
n
n n n
n
(x̄0n (Ṽn+ )2 x̄n )2

ξ 0 (Vn+ )2 θ = ξ 0 (Ṽn+ )2 θ +
+
−

ξ 0 (Ṽn+ )2 x̄n x̄0n (Ṽn+ )2 θ + ξ 0 Ṽn+ x̄n x̄0n (Ṽn+ )3 θ + ξ 0 (Ṽn+ )3 x̄n x̄0n Ṽn+ θ a.s.
−→ (c − 1)−3 ξ 0 θ ,
+ 2
0
x̄n (Ṽn ) x̄n

x̄0n (Vn+ )2 x̄n = x̄0n (Ṽn+ )2 x̄n +

x̄0n (Ṽn+ )4 x̄n
(x̄0n (Ṽn+ )3 x̄n )2
x̄0n (Ṽn+ )3 x̄n
−
+
2
(x̄0n (Ṽn+ )2 x̄n )2 (x̄0n (Ṽn+ )2 x̄n )3
x̄0n (Ṽn+ )2 x̄n

− x̄0n (Ṽn+ )2 x̄n − 2
=

x̄0n (Ṽn+ )3 x̄n
x̄0n (Ṽn+ )2 x̄n

x̄0n (Ṽn+ )4 x̄n
(x̄0n (Ṽn+ )3 x̄n )2 a.s.
c
−
−→
,
+ 2
+ 2
0
2
0
3
(c − 1)3
(x̄n (Ṽn ) x̄n )
(x̄n (Ṽn ) x̄n )

23

and
x̄0n (Ṽn+ )4 x̄n 0 +
(x̄0n (Ṽn+ )3 x̄n )2 0 +
x̄
Ṽ
θ
−
x̄n Ṽn θ
n
n
(x̄0n (Ṽn+ )2 x̄n )2
(x̄0n (Ṽn+ )2 x̄n )3
i
x̄0n (Ṽn+ )3 x̄n h 0
0
+
0
+ 2
+ 2
x̄
Ṽ
θ
+
x̄
(
Ṽ
)
θ
x̄
x̄
(
Ṽ
)
n
n n
n
n
n
n
(x̄0n (Ṽn+ )2 x̄n )2

x̄0n (Vn+ )2 θ = x̄0n (Ṽn+ )2 θ +
+
−

x̄0n (Ṽn+ )2 x̄n x̄0n (Ṽn+ )2 θ + x̄0n (Ṽn+ )3 θ + x̄0n (Ṽn+ )3 x̄n x̄0n Ṽn+ θ a.s.
−→ 0
x̄0n (Ṽn+ )2 x̄n

for p/n −→ c ∈ (1, +∞) as n → ∞.

Lemma 5.7. Assume (A2). Let θ and ξ be universal nonrandom vectors with bounded Euclidean
+

0

+

Vn ηη Vn
+
norms and let P+
n = Vn − η 0 V+ η where η is a universal nonrandom vectors with bounded Euclidean
n

norm. Then it holds that
ξ

0

P+
nθ

a.s.

−1

−→ c

−1

(c − 1)



ξ 0 ηη 0 θ
0
,
ξθ−
η0η

1
,
c−1
a.s.
x̄0n P+
n θ −→ 0 ,
a.s.

x̄0n P+
n x̄n −→

(5.58)
(5.59)



ξ 0 ηη 0 θ
a.s.
0
+ 2
−3
0
ξ (Pn ) θ −→ (c − 1)
ξθ−
,
η0η
c
a.s.
2
x̄0n (P+
,
n ) x̄n −→
(c − 1)3
a.s.

2
x̄0n (P+
n ) θ −→ 0

Proof of Lemma 5.7: It holds that
ξ

P+
nθ

= ξ

0

Vn+ θ



ξ 0 Vn+ ηη 0 Vn+ θ a.s. −1
ξ 0 ηη 0 θ
−1
0
−
−→ c (c − 1)
ξθ−
η0η
η 0 Vn+ η

for p/n −→ c ∈ (1, +∞) as n → ∞ following (5.35). Similarly, we get
0
+
x̄0n P+
n x̄n = x̄n Vn x̄n −

x̄0n Vn+ ηη 0 Vn+ x̄n a.s. 1
−→
c−1
η 0 Vn+ η

and
0
+
x̄0n P+
n θ = x̄n Vn θ −

x̄0n Vn+ ηη 0 Vn+ θ a.s.
−→ 0
η 0 Vn+ η

for p/n −→ c ∈ (1, +∞) as n → ∞.
The rest of the proof follows from the equality
2
+ 2
(P+
n ) = (Vn ) −

(5.60)
(5.61)
(5.62)

for p/n −→ c ∈ (1, +∞) as n → ∞.

0

(5.57)

+
0 +
(Vn+ )2 ηη 0 Vn+ Vn+ ηη 0 (Vn+ )2
0
+ 2 Vn ηη Vn
−
+
η
(V
)
η
n
η 0 Vn+ η
η 0 Vn+ η
(η 0 Vn+ η)2

and Lemma 5.6.

24

Proof of Theorem 2.2: In case of c > 1, the optimal shrinkage intensity is given by
+
αn

=

=

γ −1

ŵS0 ∗ (µn − γΣn b) − b0 (µn − γΣn b)
ŵS0 ∗ Σn ŵS ∗ − 2b0 Σn ŵS ∗ + b0 Σn b

γ −1

10 S∗n (µn − γΣn b)
+ γ −1 ȳn0 Q̂∗n (µn − γΣn b) − b0 (µn − γΣn b)
10 S∗n 1
,
10 S∗n Σn S∗n 1
ȳn0 Q̂∗n Σn S∗n 1
b0 Σn S∗n 1
−1
−2
−1
0
∗
0
0
∗
∗
+ 2γ
+ γ ȳn Q̂n Σn Q̂n ȳn − 2 0 ∗
− 2γ b Σn Q̂n ȳn + b Σn b
(10 S∗n 1)2
10 S∗n 1
1 Sn 1
−1/2

1/2

where ȳn = µn + Σn x̄n and S∗n = Σn

−1/2

Vn+ Σn

.
−1/2

From Assumption (A3), we get that the following vectors Σn
and

1/2
Σn b

−1/2

1, Σn

−1/2

(µn − γΣn b), Σn

µn ,

possess bounded Euclidean norms. As a result, the application of Lemma 5.2 leads to
a.s.

10 S∗n 1 = 10 Σ−1/2
Vn+ Σ−1/2
1 −→ c−1 (c − 1)−1 10 Σ−1
n
n
n 1,
a.s.

10 S∗n (µn − γΣn b) = 10 Σ−1/2
Vn+ Σ−1/2
(µn − γΣn b) −→ c−1 (c − 1)−1 10 Σ−1
n
n
n (µn − γΣn b),
a.s.

−1
−1
b0 Σn S∗n 1 = b0 Σn Σ−1/2
Vn+ Σ−1/2
1 −→ c−1 (c − 1)−1 b0 Σn Σ−1
n
n
n 1 = c (c − 1) ,
a.s.

10 S∗n Σn S∗n 1 = 10 Σ−1/2
(Vn+ )2 Σ−1/2
1 −→ (c − 1)−3 10 Σ−1
n
n
n 1
for p/n −→ c ∈ (1, +∞) as n → ∞.
Finally, from Lemma 5.3 and 5.4 as well as by using the equalities
−1/2
+
+ 2
Q̂∗n = Σ−1/2
P+
and P+
n
n Σn
n Vn = (Vn ) −

−1/2

with η = Σn

Vn+ ηη 0 (Vn+ )2
η 0 Vn+ η

1 we obtain

−1/2
P+
(µn − γΣn b) + x̄0n Pn Σ−1/2
(µn − γΣn b)
ȳn0 Q̂∗n (µn − γΣn b) = µ0n Σ−1/2
n
n Σn
n
a.s.

−→ c−1 (c − 1)−1 µ0n Qn (µn − γΣn b),
+ −1/2
+ −1/2
ȳn0 Q̂∗n Σn S∗n 1 = µ0n Σ−1/2
P+
1 + x̄0n P+
1
n
n Vn Σn
n Vn Σn
a.s.

−→ (c − 1)−3 µ0n Qn 1 = 0,
2
0
−1/2
2
0
−1/2
2 −1/2
ȳn0 Q̂∗n Σn Q̂∗n ȳn = x̄0n (P+
(P+
(P+
µn
n ) x̄n + 2µn Σn
n ) x̄n + µn Σn
n ) Σn
c
a.s.
+ (c − 1)−3 µ0n Qn µn ,
−→
(c − 1)3
+
0 1/2 + −1/2
b0 Σn Q̂∗n ȳn = b0 Σ1/2
µn
n Pn x̄n + b Σn Pn Σn


−1
0
1 Σn µn
a.s. −1
−1
0
−→ c (c − 1)
b µn − 0 −1
1 Σn 1

for p/n −→ c ∈ (1, +∞) as n → ∞.
a.s.
Hence, αn+ −→ α+ , where
α+ = γ −1

10 Σ−1
γ −1
n (µn − γΣn b)
+
µ0 Qn (µn − γΣn b) − b0 (µn − γΣn b)
−1
0
c(c − 1) n
1 Σn 1


10 Σ−1
c2
1
γ −2
1
n µn
γ −1
0 Q µ + c) − 2
0µ −
+
(µ
−
2
b
+ b 0 Σn b
n n
n
c(c−1)
(c − 1) 10 Σ−1
(c − 1)3 n
10 Σ−1
10 Σ−1
n 1
n 1
n 1

This completes the proof of Theorem 2.2.

References
Bai, J. and Shi, S. (2011). Estimating high dimensional covariance matrices and its applications.
Annals of Economics and Finance, 12:199–215.
25

Bai, Z. and Zhou, W. (2008). Large sample covariance matrices without independence structures in
columns. Statistica Sinica, 18(2):425–442.
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.
Bai, Z. D. and Silverstein, J. W. (2010). Spectral Analysis of Large Dimensional Random Matrices.
Springer, New York.
Best, M. J. and Grauer, R. R. (1991). On the sensitivity of mean-variance-efficient portfolios to changes
in asset means: some analytical and computational results. The Review of Financial Studies, 4:315–
342.
Bodnar, T., Mazur, S., and Podgórski, K. (2016). Singular inverse Wishart distribution with application to portfolio theory. Journal of Multivariate Analysis, 143:314–326.
Bodnar, T., Parolya, N., and Schmid, W. (2013). On the equivalence of quadratic optimization
problems commonly used in portfolio theory. European Journal of Operational Research, 229:637–
644.
Bodnar, T., Parolya, N., and Schmid, W. (2017). Estimation of the global minimum variance portfolio
in high dimensions. European Journal of Operational Research, under revision.
Bodnar, T. and Schmid, W. (2011). On the exact distribution of the estimated expected utility
portfolio weights: Theory and applications. Statistics & Risk Modeling, 28:319–342.
Chopra, V. K. and Ziemba, W. T. (1993). The effect of errors in means, variances, and covariances
on optimal portfolio choice. The Review of Financial Studies, 19:6–11.
DeMiguel, V., L., G., and U., R. (2009). Optimal versus naive diversification: How inefficient is the
1/N portfolio strategy? The Review of Financial Studies, 22:1915–1953.
El Karoui, N. (2010). High-dimensionality effects in the Markowitz problem and other quadratic
programs with linear constraints: Risk underestimation. Annals of Statistics, 38:3487–3566.
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.
Fan, J., Zhang, J., and Yu, K. (2012). Vast portfolio selection with gross-exposure constraints. Journal
of American Statistical Association, 107:592–606.
Frahm, G. and Memmel, C. (2010). Dominating estimators for minimum-variance portfolios. Journal
of Econometrics, 159:289–302.
Friesen, O., Löwe, M., and Stolz, M. (2013). Gaussian fluctuations for sample covariance matrices
with dependent data. Journal of Multivariate Analysis, 114:270–287.
Golosnoy, V. and Okhrin, Y. (2007). Multivariate shrinkage for optimal portfolio weights. European
Journal of Finance, 13:441–458.
Horn, R. A. and Johnsohn, C. R. (1985). Matrix Analysis. Cambridge University Press, Cambridge.
Hui, J. and Pan, G. (2010). Limiting spectral distribution for large sample covariance matrices with
m-dependent elements. Communications in Statistics - Theory and Methods, 39(6):935–941.

26

Jorion, P. (1986).

Bayes-Stein estimation for portfolio analysis.

The Journal of Financial and

Quantitative Analysis, 21:279–292.
Le Cam, L. and Lo Yang, G. (2000). Asymptotics in statistics: some basic concepts. Springer, New
York.
Ledoit, O. and Wolf, M. (2004). A well-conditioned estimator for large-dimensional covariance matrices. Journal of Multivariate Analysis, 88:365–411.
Marčenko, V. A. and Pastur, L. A. (1967). Distribution of eigenvalues for some sets of random
matrices. Sbornik: Mathematics, 1:457–483.
Markowitz, H. (1952). Portfolio selection. The Journal of Finance, 7:77–91.
Merton, C. R. (1980). On estimating the expected return on the market: An exploratory investigation.
Journal of Financial Economics, 8:323–361.
Meyer, C. R. (1973).

Generalized inversion of modified matrices.

SIAM Journal of Applied

Mathematics, 24:315–323.
Okhrin, Y. and Schmid, W. (2006).

Distributional properties of portfolio weights.

Journal of

Econometrics, 134:235–256.
Pan, G. M. (2014). Comparison between two types of large sample covariance matrices. Ann. Inst. H.
Poincaré Probab. Statist., 50:655–677.
Pratt, J. W. (1964). Risk aversion in the small and in the large. Econometrica, 32:122–136.
Rubio, F. and Mestre, X. (2011). Spectral convergence for a general class of random matrices. Statistics
& Probability Letters, 81:592–602.
Silverstein, J. W. (1995). Strong convergence of the empirical distribution of eigenvalues of largedimensional random matrices. Journal of Multivariate Analysis, 55:331–339.
Tobin, J. (1958). Liquidity preference as behavior towards risk. Review of Economic Studies, 25:65–86.
Wei, M., Yang, G., and Yang, L. (2016). The limiting spectral distribution for large sample covariance
matrices with unbounded m-dependent entries. Communications in Statistics - Theory and Methods,
45(22):6651–6662.

6

Appendix B

Tables and Figures

27

1.0
0.8
0.6
0.4

α*
0.0

0.2

oracle, CN=150
bona fide, CN=150
oracle, CN=1000
bona fide, CN=1000
oracle, CN=8000
bona fide, CN=8000
0.0

0.5

1.0

1.5

2.0

c

Figure 1: The asymptotic optimal shrinkage intensity as a function of c (left) and of the sample size
n (right).

28

3.5

3.5

2.0

2.5

3.0

traditional
equally weighted
bona fide
oracle shrinkage

0.0

0.5

1.0

1.5

relative loss

2.0
1.5
0.0

0.5

1.0

relative loss

2.5

3.0

traditional
equally weighted
bona fide
oracle shrinkage

0.0

0.2

0.4

0.6

0.8

1

2

3

5

6

7

3.5
2.0

2.5

3.0

traditional
equally weighted
bona fide
oracle shrinkage

0.0

0.5

1.0

1.5

relative loss

2.0
1.5
0.0

0.5

1.0

relative loss

2.5

3.0

traditional
equally weighted
bona fide
oracle shrinkage

0.0

0.2

0.4

0.6

0.8

1

2

3

4

5

6

7

5

c

3

4

traditional
equally weighted
bona fide
oracle shrinkage

1
0

0

1

2

3

relative loss

4

traditional
equally weighted
bona fide
oracle shrinkage

2

5

c

relative loss

4
c

3.5

c

0.0

0.2

0.4

0.6

0.8

1

c

2

3

4

5

6

7

c

Figure 2: The relative losses for the portfolios based on optimal shrinkage estimator, the traditional
estimator, the bona-fide estimator and the equally weighted portfolio as a function of c for different
values of c, fixed p = 150 and different values of the condition index (150, 1000, 8000 from top to
bottom).

29

4

3

2

1

0

4

relative loss

relative loss

3

2

1

0

5

4

3

150

150

100

100

150

100

200

200

200

250

p

250

p

250

300

300

300

400

400

350

400

traditional
equally weighted
bona fide
oracle shrinkage

350

traditional
equally weighted
bona fide
oracle shrinkage

350

traditional
equally weighted
bona fide
oracle shrinkage

150

150

100

100

150

100

200

200

200

250

p

250

p

250

300

300

300

400

400

350

400

traditional
equally weighted
bona fide
oracle shrinkage

350

traditional
equally weighted
bona fide
oracle shrinkage

350

traditional
equally weighted
bona fide
oracle shrinkage

150

150

100

100

150

100

200

200

200

250

p

250

p

250

300

300

300

400

400

350

400

traditional
equally weighted
bona fide
oracle shrinkage

350

traditional
equally weighted
bona fide
oracle shrinkage

350

traditional
equally weighted
bona fide
oracle shrinkage

p

p

150

150

100

100

150

100

200

200

200

p

250

p

250

p

250

300

300

300

400

400

350

400

traditional
equally weighted
bona fide
oracle shrinkage

350

traditional
equally weighted
bona fide
oracle shrinkage

350

traditional
equally weighted
bona fide
oracle shrinkage

Figure 3: The relative losses for the portfolios based in the optimal shrinkage estimator, the traditional estimator, the bona-fide estimator and the
equally weighted portfolio as a function of the dimension p, for different values of the condition index (150, 1000, 8000 from top to bottom) and different
values of c (0.2, 0.5, 0.8, 2 from left to right).

p

1

0

2

relative loss
2
0

relative loss

6
5
4
3
2
1
0
6
5
4
3
2
1
0
10
8
6

relative loss
2

4

relative loss
relative loss

3.0
2.5
2.0
1.5
1.0
0.5
0.0
3.0
2.5
2.0
1.5
1.0
0.5
0.0
10
8
6

relative loss

0

4

relative loss
relative loss

3.0
2.5
2.0
1.5
1.0
0.5
0.0
3.0
2.5
2.0
1.5
1.0
0.5
0.0
10
8
6
4

relative loss
relative loss

2
0

30

1.0

1.0

0.8
2007

2008

2009

2010

2011

2012

2013

2014

2004

2005

2006

2007

2008

2009

2010

2011

2012

2013

shrinkage

0.6

0.8

oracle shrinkage
bona fide

0.0

0.2

0.4

0.6
0.4
0.0

0.2

2007

2008

2009

2010

2011

2012

2013

2014

2004

1.0

1.0

2006

2007

2008

2009

2010

2011

2012

2013

2014

shrinkage

0.0

0.2

0.4
0.0

0.2

0.6

0.8

oracle shrinkage
bona fide

0.4

0.6

0.8

oracle shrinkage
bona fide

2006

2006

2007

2008

2009

2010

2011

2012

2013

2014

2004

1.0

1.0

2004

2007

2008

2009

2010

2011

2012

2013

2014

0.6

0.8

oracle shrinkage
bona fide

0.2
0.0

0.0

0.2

0.4

shrinkage

0.6

0.8

oracle shrinkage
bona fide

2006

0.4

shrinkage

0.8

oracle shrinkage
bona fide

2004

shrinkage

0.6

shrinkage

0.2
0.0
2006

1.0

1.0

2004 2005

shrinkage

oracle shrinkage
bona fide

0.4

0.6
0.4
0.0

0.2

shrinkage

0.8

oracle shrinkage
bona fide

2004

2006

2007

2008

2009

2010

2011

2012

2013

2014

2004

2006

2007

2008

2009

2010

2011

2012

2013

2014

Figure 4: The oracle and the bona fide shrinkage coefficients for the first 50 assets (left) and 100 assets
(right) and c = 0.2, 0.5, 0.8 and 2 (from top to bottom).

31

p = 50 and c = 0.2
SR
average
median
average median
-1.77657
-1.72203
0.04145 0.04230
-1.66834
-1.63833
0.03892 0.03985
-1.64857
-1.62718
0.03766 0.03823
-2.45857
-2.45519
0.01821 0.01823
p = 50 and c = 0.5
CE
SR
average
median
average median
-3.19607
-2.94847
0.04238 0.04301
-2.22557
-2.14978
0.03586 0.03566
-2.07126
-2.02931
0.03314 0.03279
-2.47466
-2.45833
0.01786 0.01769
p = 50 and c = 0.8
CE
SR
average
median
average median
-19.03932 -17.92892 0.02148 0.02255
-2.39935
-2.39654
0.02056 0.02099
-2.62398
-2.61298
0.02007 0.01994
-2.46909
-2.45949
0.01789 0.01741
p = 50 and c = 2
CE
SR
average
median
average median
-24.65195 -24.01919 -0.0137 -0.01339
-9.45345
-9.26547 -0.01346 -0.01487
-12.25876 -12.05119 -0.01428 -0.01449
-2.46299
-2.45453
0.01766 0.01833
CE

trad
oracle
bona fide
equal

trad
oracle
bona fide
equal

trad
oracle
bona fide
equal

trad
oracle
bona fide
equal

VaR
α = 0.05 α = 0.01
-0.95454 -1.44646
-0.95813 -1.45896
-0.96505 -1.47793
-1.24743 -2.21880

ES
α = 0.05 α = 0.01
-1.29504 -1.75337
-1.30192 -1.74567
-1.31451 -1.75670
-1.83318 -2.44254

VaR
α = 0.05 α = 0.01
-1.21000 -1.87739
-1.09712 -1.70294
-1.09301 -1.76240
-1.25117 -2.23102

ES
α = 0.05 α = 0.01
-1.67650 -2.33106
-1.50856 -2.04258
-1.53212 -2.07562
-1.83897 -2.44971

VaR
α = 0.05 α = 0.01
-2.88101 -4.89415
-1.2319 -2.07394
-1.26838 -2.05094
-1.25167 -2.22567

ES
α = 0.05 α = 0.01
-4.32690 -6.59821
-1.77196 -2.37602
-1.78578 -2.41087
-1.83918 -2.45105

VaR
α = 0.05 α = 0.01
-3.52458 -6.01948
-2.23232 -3.63413
-2.51985 -4.17245
-1.25125 -2.22469

ES
α = 0.05 α = 0.01
-5.34400 -8.22013
-3.24817 -4.83833
-3.72291 -5.58275
-1.8364 -2.45083

Table 1: Certainty equivalent, Sharpe ratio, Value-at-Risk and expected shortfall averaged over 1000
randomly chosen portfolios of size 50 and different values of c. The risk aversion is set to 10. The
trading period is fixed at 200 days from 01.01.2014 to 10.10.2014.

32

CE
trad
oracle
bona fide
equal

average
-1.50263
-1.43215
-1.41780
-2.39851

median
-1.42892
-1.38207
-1.37610
-2.39802

trad
oracle
bona fide
equal

average
-2.50439
-1.92386
-1.82298
-2.39526

trad
oracle
bona fide
equal

CE
average
median
-13.89420 -12.93480
-2.30811
-2.30721
-2.47852
-2.46092
-2.39781
-2.39383

trad
oracle
bona fide
equal

CE
average
median
-63.10889 -61.07481
-26.12505 -25.61899
-30.11004 -29.40366
-2.39220
-2.38997

CE
median
-2.46847
-1.90458
-1.80993
-2.39368

p = 100 and c = 0.2
SR
average median
0.07636 0.07835
0.07162 0.07262
0.06988 0.07078
0.01826 0.01799
p = 100 and c = 0.5
SR
average median
0.04978 0.04968
0.04506 0.04521
0.04266 0.04253
0.01839 0.01842
p = 100 and c = 0.8
SR
average median
0.02525 0.02473
0.02189 0.02106
0.02245 0.02133
0.01916 0.01888
p = 100 and c = 2
SR
average median
-0.00440 -0.00529
-0.00769 -0.00701
-0.00782 -0.00740
0.01807 0.01815

VaR
α = 0.05 α = 0.01
-0.83652 -1.23564
-0.84756 -1.26211
-0.85319 -1.27461
-1.23748 -2.22931

α = 0.05
-1.12160
-1.14206
-1.15202
-1.82561

ES
α = 0.01
-1.51789
-1.53457
-1.54414
-2.42345

VaR
α = 0.05 α = 0.01
-1.12221 -1.71939
-1.03271 -1.58796
-1.02794 -1.62191
-1.24098 -2.22600

α = 0.05
-1.53858
-1.41473
-1.42840
-1.82481

ES
α = 0.01
-2.13992
-1.92871
-1.93419
-2.42005

VaR
α = 0.05 α = 0.01
-2.50742 -3.98555
-1.21885 -2.10329
-1.23968 -2.05522
-1.23787 -2.2234

α = 0.05
-3.55934
-1.76802
-1.76667
-1.82392

ES
α = 0.01
-5.13121
-2.33794
-2.34685
-2.42085

VaR
α = 0.05 α = 0.01
-5.35257 -9.20368
-3.56410 -5.88731
-3.82082 -6.28483
-1.23619 -2.22278

α = 0.05
-8.03425
-5.19438
-5.53534
-1.82301

ES
α = 0.01
-12.07135
-7.61019
-8.06647
-2.41900

Table 2: Certainty equivalent, Sharpe ratio, Value-at-Risk and expected shortfall averaged over 1000
randomly chosen portfolios of size 100 and different values of c. The risk aversion is set to 10. The
trading period is fixed at 200 days from 01.01.2014 to 10.10.2014.

33

CE
trad
oracle
bona fide
equal

average
-1.46983
-1.45115
-1.44714
-2.35315

median
-1.46470
-1.44582
-1.44190
-2.35454

trad
oracle
bona fide
equal

average
-1.78646
-1.50546
-1.45753
-2.35470

trad
oracle
bona fide
equal

average
-9.05407
-2.20594
-2.43972
-2.35623

trad
oracle
bona fide
equal

CE
average
median
-86.06159 -85.4061
-38.24696 -38.00040
-43.66032 -43.36126
-2.35435
-2.35629

CE
median
-1.76039
-1.49837
-1.45435
-2.35703
CE
median
-8.92031
-2.20383
-2.43447
-2.35311

p = 300 and c = 0.2
SR
average median
0.09563 0.09559
0.09383 0.09378
0.09339 0.09339
0.01825 0.01824
p = 300 and c = 0.5
SR
average median
0.07088 0.07094
0.06469 0.06509
0.06268 0.06317
0.01825 0.01830
p = 300 and c = 0.8
SR
average median
0.01743 0.01827
0.01976 0.02006
0.01989 0.02001
0.01830 0.01828
p = 300 and c = 2
SR
average median
-0.03174 -0.03060
-0.02245 -0.02107
-0.02666 -0.02532
0.01830 0.01820

VaR
α = 0.05 α = 0.01
-0.77747 -1.32552
-0.77448 -1.31440
-0.77396 -1.31201
-1.23726 -2.24940

α = 0.05
-1.19412
-1.18989
-1.18914
-1.81895

ES
α = 0.01
-2.00830
-2.01176
-2.01262
-2.40111

VaR
α = 0.05 α = 0.01
-0.89048 -1.29435
-0.85834 -1.24359
-0.85873 -1.25050
-1.23801 -2.24995

α = 0.05
-1.16944
-1.12347
-1.12740
-1.81966

ES
α = 0.01
-1.52942
-1.45169
-1.45647
-2.40180

VaR
α = 0.05 α = 0.01
-2.15178 -3.32989
-1.25841 -2.03914
-1.29504 -1.98745
-1.23702 -2.25272

α = 0.05
-2.97405
-1.76796
-1.76910
-1.81999

ES
α = 0.01
-4.14909
-2.28501
-2.30087
-2.40182

VaR
α = 0.05 α = 0.01
-6.81331 -10.16486
-4.54427 -6.70359
-4.86611 -7.18398
-1.23762 -2.25067

α = 0.05
-9.16014
-6.05765
-6.49152
-1.81970

ES
α = 0.01
-12.42224
-8.18999
-8.76704
-2.40095

Table 3: Certainty equivalent, Sharpe ratio, Value-at-Risk and expected shortfall averaged over 1000
randomly chosen portfolios of size 300 and different values of c. The risk aversion is set to 10. The
trading period is fixed at 200 days from 01.01.2014 to 10.10.2014.

34

