Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
41,211 characters · 10 sections · 38 citation commands
Portfolio Optimization Using a Consistent Vector-Based MSE Estimation Approach
Decision making with regards to investment in the stock market has become increasingly more complex because of the dynamic nature of the stocks available to investors and the advent of new unconventional and risky options 2019-bodnar-tests. Throughout the years, the portfolio optimization problem has attracted the attention of many signal-processing researchers due to its close relationship to the field. The portfolio optimization problem aims at achieving the maximum possible returns with the least volatility percentage 2020-huni-application. The Economist, Harry Markowitz, introduced the modern portfolio theory, or mean-variance analysis (MVP), in 1952-markowitz-mvp. Other portfolios such as the global MVP and the maximum sharp ratio portfolio (MSRP) have been proposed as improvements of the MVP. Portfolio optimization utilizes the available financial data to reach conclusions regarding the allocation of wealth to each of the available stocks. The most important measurement in portfolio optimization is the data covariance matrix (CM).
CM estimation in the classical signal processing framework relies on asymptotic statistics of a number of observations, $n$, which is assumed to grow largely compared to the population dimension, $p$, i.e., $n/p \rightarrow \infty$ as $n\rightarrow \infty$ 2012-couillet-signal. However, many practical applications, such as finance, bioinformatics and data classification, require an estimate of the CM when the data dimension is large compared to the sample size 2019-ollila-optimal. In such cases, it is well known that the default estimator, i.e., the empirical sample covariance matrix (SCM), is usually ill-conditioned, leading to poor performance.
If the case where $p>n$, the SCM is not invertible; whereas for $p<n$, the SCM is invertible but might be ill-conditioned, which substantially increases estimation error. In other words, for a large $p$, it is not practically guaranteed that the number of observations is sufficient to develop a well-conditioned CM estimator 2004-ledoit-well. Such scenarios have motivated researchers to look into estimation problems in the high-dimensional regime 2012-couillet-signal.
In scenarios with limited data, a regularized SCM (RSCM) estimator of the following general form is widely used 2019-ollila-optimal:
where $\widehat{\hbox{\boldmath$\Sigma$}}$ is the SCM defined in (ref) (further ahead), $\beta$, $\gamma \in \mathbb{R}^{+}$ are the regularization, or shrinkage, parameters. These parameters can be determined based on minimizing the mean-squared error (MSE), which results in oracle shrinkage parameters, $\beta_\text{o}$ and $\gamma_\text{o}$, as follows 2019-ollila-optimal, 2004-ledoit-well:
where $\|.\|_{\text{F}}$ denotes the Forbenius matrix norm. The estimation of $(\beta_\text{o}, \gamma_\text{o})$ based on (ref) depends on the true CM, $\hbox{\boldmath$\Sigma$}$. To circumvent this issue, Ledoit and Wolf 2004-ledoit-well proposed a distribution-free consistent estimator of $(\beta_\text{o}, \gamma_\text{o}) $ in high-dimensional settings. The work in 2019-ollila-optimal assumes that the observations are from unspecified elliptically symmetric distribution. The consistent estimator proposed in 2015-yang-robust uses a hybrid CM estimator based on the Taylor's M-estimator and Ledoit-Wolf shrinkage estimator, which suits a global minimum variance portfolio (GMVP) influenced by outliers. A similar approach based on the M-estimator is proposed in 2020-ollila-shrinking, considering $n>p$ with fully automated selection of the shrinkage parameters. The minimum variance portfolio estimator in 2020-cai-high is based on certain sparsity assumptions imposed on the inverse of the CM. The work presented in 2021-ballal-adaptive proposes a different RSCM estimator by manipulating the expression of the GMVP weights.
In this paper, we propose a single-parameter CM estimator. Instead of minimizing the MSE, as in (ref), we minimize the MSE of the estimation of the sample noise vector. We utilize RMT tools to obtain a consistent estimator of this MSE. The value of the regularization parameter $\gamma$ is selected as the one that minimizes the estimated MSE. By choosing to minimize the MSE of the noise vector's estimation, we consider the inaccuracy of estimating the true mean.
We consider a time series comprising ${\bf y}_1, {\bf y}_2 \cdots, {\bf y}_L$ logarithmic returns of $p$ financial assets over a certain investment period. We assume that the elements of ${\bf y}_t$, ($t=1,2, \cdots, L$) are independent and identically distributed (i.i.d.) and are generated according to the following stochastic model 2012-rubio-performance:
where $\hbox{\boldmath$\mu$}_t \in \mathbb{R}^{n \times 1}$ and $\hbox{\boldmath$\Sigma$}_t \in \mathbb{R}^{p \times p}$ are the mean and the CM of the asset returns over the investment period, and ${\bf x}_t$ is an i.i.d. random noise vector of zero mean and identity CM. For simplicity, we drop the subscript $t$ from $\hbox{\boldmath$\mu$}_t$ and $\hbox{\boldmath$\Sigma$}_t$. For the investment period of interest, we define ${\bf w} \in \mathbb{R}^p$ as the asset holdings vector, also known as the weight vector. The GMVP optimally minimizes the portfolio variance under single-period investment horizon, such that the weight vector is normalized by the outstanding wealth 2012-rubio-performance, i.e.,
where $\mathds{1}_p$ is a column vector of $p$ 1's. The solution of (ref) can be obtained by using the Lagrange-multipliers method, which results in the optimum weights 2015-yang-robust:
The CM in (ref) is unknown and should be estimated. As stated earlier, the SCM estimate does not perform well because it is usually ill-conditioned; hence, we apply the RSCM estimator and (ref) becomes
where $\widehat{\hbox{\boldmath$\Sigma$}}_\text{RSCM}$ is the RSCM which can take the form of (ref), for example. In the following section, we develop a RSCM estimator method and properly set the value of its regularization parameter.
The SCM, $\widehat{\boldsymbol{\Sigma}}$, and the sample mean, $\widehat{\hbox{\boldmath$\mu$}}$, can be estimated from the $n$ past return observations as follows:
We notice that computing $\widehat{\hbox{\boldmath$\Sigma$}}$ using (ref) involves evaluating the sample mean, not the true mean. This can worsen performance, especially for a small number of observations. Subtracting $\widehat{\hbox{\boldmath$\mu$}}$ from both sides of (ref), we obtain
where $\hbox{\boldmath$\delta$} \triangleq \hbox{\boldmath$\mu$} - \widehat{\hbox{\boldmath$\mu$}}$ and $ \hbox{\boldmath$\Delta$} = \hbox{\boldmath$\Sigma$}^{\frac{1}{2}}-\widehat{{\hbox{\boldmath$\Sigma$}}}^ {\frac{1}{2}} $. Eq. (ref) can be viewed as a linear model with bounded uncertainties in both $\hbox{\boldmath$\Sigma$}^\frac{1}{2}$ and $\hbox{\boldmath$\mu$}$ 1998-chandrasekaran-SIAMparameter. We seek an estimate, $\widehat{{\bf x}}_t$ that performs well for any allowed perturbation $(\hbox{\boldmath$\Delta$}, \hbox{\boldmath$\delta$})$ by formulating the following min-max problem 1998-chandrasekaran-SIAMparameter:
A unique solution can exist which takes the form 1998-chandrasekaran-SIAMparameter
$\widehat{{\bf x}}_t $ is a function of $\gamma$, which when properly set leads to the best estimate of ${\bf x}_t$. It is easy to recognize that $(\widehat{\hbox{\boldmath$\Sigma$}} + \gamma \mathbf{I})^{-1}$ can be used as an estimator of the CM inverse, i.e., $\widehat{\hbox{\boldmath$\Sigma$}}_\text{RSCM}^{-1} = \widehat{\hbox{\boldmath$\Sigma$}}_\gamma^{-1}= (\widehat{\hbox{\boldmath$\Sigma$}} + \gamma \mathbf{I})^{-1}$. Such estimator is widely used in the literature, e.g., 2007-guo-regularized,2015-zollanvari-generalized,2017-ballal-bounded,2020-elkhalil-large,1988-carlson-covariance,2003-li-robust,2021-mahadi-low,2021-mahadi-robust,2018-suliman-robust; to name a few. The optimal value of $\gamma$ that estimates $\widehat{\hbox{\boldmath$\Sigma$}}_\gamma^{-1}$ is the one that minimizes the MSE for estimating ${{\bf x}}_t$. That is
We choose the optimal $\gamma_\text{o}$ as follows:
The choice of minimizing the MSE is reasonable because, under certain conditions, the minimization problem in (ref) and the minimum MSE are equivalent 1993-kay-fundamentals, 2008-rubio-generalized. Unlike the other methods, it is remarkable that the uncertainty in estimating the mean is incorporated in (ref). We expect the effect of the uncertainty in the mean estimation to be high when we have a limited number observations. Also, unlike the methods that are based on (ref), when we search for the optimal $\gamma$ that minimizes (ref), we actually estimate the inverse of the CM rather than estimating the CM itself. This is important because we use it in (ref). We obtain the following normalized (by $n$) expression of the MSE (see Appendix (ref)):
We observe that (ref) is expressed in terms of the unknown quantity, $\boldsymbol{\Sigma}$. In this case, using a direct plugin formula, i.e., substituting $\hbox{\boldmath$\Sigma$}$ with $\widehat{\hbox{\boldmath$\Sigma$}}$ results in
However, the estimator in (ref) is an inconsistent estimator in the regime where $n$ and $p$ grow at constant rate 2020-lama-asymptotic. To clarify, Fig. (ref) plots an example of the derived MSE$(\gamma)$(ref) and the plugin estimation method (ref) versus a wide range values of $\gamma$. It is clear that using the plugin strategy does not help obtain the minimum MSE suitably. Instead, as the figure depicts, the plugin estimation method selects an improper $\gamma$ that corresponds to a high MSE.
As an alternative remedy , we seek a consistent estimator of (ref) by leveraging tools from RMT. To this end, we need to first obtain an asymptotic expression of (ref). To do so, the following assumption should hold true.
Assumption 1: As $p$,$n$ $\rightarrow \infty$, $p/n \rightarrow c\in (0,\infty)$.
Assumption 1 leads to the following theorem:
Proof: see Appendix (ref).
Now, we are in a position to reveal the consistent estimator of (ref).
Proof: see Appendix (ref).
Back to Fig. (ref) which compares the derived MSE with the asymptotic formula (ref) and the consistently estimated MSE (ref). It can be seen clearly that the consistent MSE is more suitable to obtain the value of $\gamma$ that minimizes (ref).
A closed form solution for $\gamma$ in (ref) is infeasible, so we rely on using a line search, where we search for $\gamma$ that minimizes (ref) within a predefined range.
As conventionally described in the financial literature, we implement the out-of-sample strategy defined in terms of a rolling window method (see 2015-yang-robust). At a particular day $t$, the training window for CM estimation is formed from the previous $n$ days, i.e., from $t-n$ to $t-1$, to design the portfolio weights, $\widehat{{\bf w}}_{\text{GMVP}}$. The portfolio returns in the following 20 days are computed based on these weights. Next, the window is shifted 20 days forward and the returns for another 20 days are computed. The same procedure is repeated until the end of the data. Finally, the realized risk is computed as the standard deviation of the returns. The following list describes the data from different stock market indices used in our evaluation:
Fig. (ref) shows the annualized realized risk of the aforementioned market indices versus the number of training samples. We compare the proposed vector-based method, VB-MSE, against the elliptical estimators ELL1-RSCM, ELL2-RSCM and ELL3-RSCM 2019-ollila-optimal, 2018-ollila-matlab, the Ledoit-Wolf estimator, LW-RSCM 2004-ledoit-well, 2018-ollila-matlab, the nonlinear estimator Quest 1 2017-ledoit-numerical .
Fig. (ref) (\subref{fig:SP11-12}) plots the result of the S&P 100 index from 2 Jan. 2015 to 30 Dec. 2016. As can be seen from the figure, the performance of the proposed VB-MSE method outperforms all other the methods except at $n= 80$ and $100$, where it is slightly worse than Ell1-RSCM aand Ell3-RSCM. Similarly, VB-MSE has a superior performance in Fig. (ref) (\subref{fig:SP15-16}), which plots the result from 7 Jan. 2014 to 31 Dec. 2015. However, at $n = 20$ and $80$ Ell1-RSCM and Ell3-RSCM perform better. The realized risk for the HSI index is depicted in Fig. (ref) (\subref{fig:HSI14-15}) from 1 Jan. 2016 -- 27 Dec. 2017. The proposed method has a comparable performance to Quest 1, Ell1-RSCM and Ell3-RSCM at $n = 20, 40$ and $60$ but it outperforms all the methods for $ 100 <n\le 340$. The results of the XMI index from 4 Jan. 2016 -- 29 Dec. 2017 and from 10 Jan. 2014 -- 31 Dec. 2015 are shown in Fig. (ref) (\subref{fig:XMI16-17}) and Fig. (ref) (\subref{fig:XMI14-15}), respectively. Overall, in both figures, VB-MSE is the best performing method. Finally, Fig. (ref) (\subref{fig:SP500_15-17}) plots the realized risk of the S&P 500 index from 10 Jan. 2015 -- 31 Dec. 2017. The figure shows clearly that the proposed method outperforms the other methods when $ 200 \le n \le 400$.
From Fig. (ref) (\subref{fig:SP11-12}) -- (\subref{fig:SP500_15-17}), we can conclude that, on average, the proposed VB-MSE method compares favorably to all the benchmark methods tested in this paper. The method is also more consistent over the various datasets.
In this paper, we have proposed a regularized covariance matrix estimator under high-dimensionality settings. The proposed method searches for the optimal regularization parameter based on a consistent estimator of the MSE of the estimated vector. Portfolio optimization results from real financial data show that the proposed method performs reasonably well and outperforms a host of benchmark methods.
For convenience, we write Equation (ref) in matrix form
where ${\bf X} = [{\bf x}_1 {\bf x}_2 \cdots {\bf x}_n]$ with ${\bf x}_i \sim \mathcal{N}(\boldsymbol{0}, {\bf I}_p)$. We need to express the SCM in (ref) in an appropriate matrix form as well, as follows:
where ${\bf B} \in \mathbb{R}^{p \times n}$. It can be immediately recognized from (ref) that ${\bf B}$ is
Also, we can easily verify that
where ${\bf Z} \triangleq \hbox{\boldmath$\Sigma$}^{\frac{1}{2}} {\bf X}$. Finally, we perform the following operations to reach the model of $\widehat{{\hbox{\boldmath$\Sigma$}}}$ at the end:
where ${\bf U}$ and ${\bf T}$ are the matrices of eigenvalue vectors and eigenvalues, respectively, of $( {\bf I}_n - \dfrac{\mathds{1}_n\mathds{1}_n^T}{n})$ obtained using the eigenvalue decomposition. The Gaussian distribution is invariant when multiplying by a unitary matrix; hence, $\widetilde{{\bf X}}$ has the same distribution as ${\bf X}$ 2015-zollanvari-generalized, 2020-elkhalil-large.
The model (ref) is a well-established model is the RMT literature. Based on this model, for $z \in \mathbb{C} - \mathbb{R}^+$ and bounded $\hbox{\boldmath$\Theta$} \in \mathbb{R}^{p \times p}$ the following relations, which will be used throughout the derivations, hold true (under Assumption 1) 2012-rubio-performance:
The MSE in (ref) can be easily obtained from expanding (ref) and computing the resulted terms. The first term results from $\mathbb{E} [ {\bf x}_t{\bf x}_t^T] = {\bf I}_p$. The second term is computed as follows:
where $\tilde{{\bf y}} \triangleq({\bf y} - \widehat{\hbox{\boldmath$\mu$}})= \left( \hbox{\boldmath$\Sigma$}^\frac{1}{2} {\bf x}_t + \hbox{\boldmath$\delta$} \right) $. Using the fact that the expectation and the trace are interchangeable and the cyclic property of traces, we can write
Also, using the eigenvalue decomposition of $\widehat{\hbox{\boldmath$\Sigma$}} $, it is easy to prove that
Hence,
Observing that ${{\bf x}_t}$, $\hbox{\boldmath$\delta$}$ and $\widehat{\hbox{\boldmath$\Sigma$}}$ are independent, and $\hbox{\boldmath$\delta$} \sim \mathcal{N} (\boldsymbol{0}, \frac{\hbox{\boldmath$\Sigma$}}{n})$, we obtain
Finally, the third term is obtained from
$\mathbb{E}[A(\gamma)]$ in (ref) can be directly obtained from (ref) with setting $\hbox{\boldmath$\Theta$} = \hbox{\boldmath$\Sigma$}$ and $z = -\gamma$. The second term, $B(\gamma)$, resulted from adding and subtracting $i\sqrt{\gamma} {\bf I}_p$ with factoring $\widehat{\hbox{\boldmath$\Sigma$}}_\gamma^{-1}$ as follows:
We can further simplify (ref) by noticing that the quantity on the left-hand side is a real quantity, so this implies
Thus, we can express $B(\gamma)$ equivalently as
so we can find $\mathbb{E} [B(\gamma)]$ easily from (ref) with setting $\hbox{\boldmath$\Theta$} = \hbox{\boldmath$\Sigma$}^\frac{1}{2}$, $z = i\sqrt{\gamma}$, and the third term in (ref) resulted.
The MSE($\gamma$) expressed in (ref) converges to a sum of deterministic terms. To find a consistent estimator of (ref), it is sufficient to find a consistent estimator of each of these terms (Theorem 3.2.6 in 1985-takeshi-advanced).
The consistent estimator, $\hat{\delta}_1$, can be derived using (ref) and (ref) in (ref) ( with $\hbox{\boldmath$\Theta$} = {\bf I}_p$),
The consistent estimator, $\hat{\delta}_1$, in (ref) (that satisfies $\hat{\delta} \asymp \delta$) results immediately after rearranging (ref). The derivation of $\hat{\delta}_2$ follows similarly.
To derive the consistent estimator of $\phi \triangleq \dfrac{1}{n}\Bigl[ \hbox{\boldmath$\Sigma$}^2\bigl( \tilde{\delta}_1\hbox{\boldmath$\Sigma$} + \gamma {\bf I} \bigr)^{-2} \Bigr] $ we differentiate $\delta_1$ to obtain
where $\Psi \triangleq \dfrac{1}{n}{\hbox{tr}} [\hbox{\boldmath$\Sigma$} (\tilde{\delta_1}\hbox{\boldmath$\Sigma$} + \gamma{\bf I}_p)^{-2}]$. Hence, we can estimate ${\phi}$ consistently, (i.e., $\hat{\phi} \asymp \phi$) as
where $\widehat{\Psi}$ is the consistent estimator of $\Psi$ can be estimated from (ref) when $\boldsymbol{\Theta} = {\bf I}$, $z=-\gamma$, as follows:
Substituting in ((ref)), we can obtain the second term as