EconBase
← Back to paper

Robustifying Markowitz

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.

66,616 characters · 15 sections · 82 citation commands

Rendered from LaTeX for readability, not typeset faithfully. Citation keys are highlighted; maths is left as source; figures, tables and equation environments are summarised rather than reproduced; unrecognised commands are greyed out so nothing is silently dropped. Email addresses are removed.

Robustifying Markowitz

\sloppy

abstractMarkowitz mean-variance portfolios with sample mean and covariance as input parameters feature numerous issues in practice. They perform poorly out of sample due to estimation error, they experience extreme weights together with high sensitivity to change in input parameters. The heavy-tail characteristics of financial time series are in fact the cause for these erratic fluctuations of weights that consequently create substantial transaction costs. In robustifying the weights we present a toolbox for stabilizing costs and weights for global minimum Markowitz portfolios. Utilizing a projected gradient descent (PGD) technique, we avoid the estimation and inversion of the covariance operator as a whole and concentrate on robust estimation of the gradient descent increment. Using modern tools of robust statistics we construct a computationally efficient estimator with almost Gaussian properties based on median-of-means uniformly over weights. This robustified Markowitz approach is confirmed by empirical studies on equity markets. We demonstrate that robustified portfolios reach the lowest turnover compared to shrinkage-based and constrained portfolios while preserving or slightly improving out-of-sample performance.

Introduction

The cornerstone mean-variance portfolio theory proposed by Markowitz1952 plays a significant role in research and practice. Efficient mean-variance portfolios (MV) experience a number of attractive properties and have a simple and straightforward analytical solution with only two input parameters: the expected mean and covariance matrix of asset returns. Mean-variance analysis is naturally connected to the Capital Asset Pricing Model (CAPM), a standard tool in asset pricing.

Despite its simplicity and theoretical appeal, implementation of mean-variance portfolios is often impractical. The traditional approach to use the sample moments as input parameters leads to extreme negative and positive weights, and extensive literature documents poor out-of-sample performance of such plug-in approach, see frost1986empirical, Frost29, Best1991, Chopra1993, broadie1993computing, litterman2004modern, merton1980estimating. {The problem might be seen as an inverse problem, and it} has high sensitivity to even small perturbations of the input estimates: the mean and {the} covariance matrix. It is possibly surprising that the MV portfolios are more sensitive to changes in the mean estimate, jagannathan2003risk spell this out explicitly by writing that the error of mean estimation is so large that nothing is lost when one ignores the mean at all, and michaud1989markowitz describes the influence of the mean error as “error-maximization”. {Following the majority of research on this topic, } we focus here on global minimum variance portfolios (GMV), which only depend on the covariance.

However, even with the mean left out of the equation, traditional policies suffer from extreme instability, which means that the portfolio weights fluctuate significantly over time. Drastic changes in the portfolio composition lead to increasing management and transaction costs and consequently reducing the popularity of MV policies among investors. In order to improve upon the stability of portfolio weights, one has to resort to alternative, robust estimation techniques. A robust estimator is one that performs well even when the observations do not follow the standard (normality) assumptions, have heavy tails, or are even subject to contamination. Although in case of normal distributions, sample moment estimators are asymptotically optimal MLEs, they are not necessarily the best choice when the data deviates from normality huber2004robust. This is of particular importance in financial applications, where it is well known that the data is not only non-Gaussian, but also exhibits heavy tails.

To tackle the problem of heavy tails, demiguel2009portfolio construct a portfolio optimization procedure based on M- and S-estimation technique and analyze the stability of the estimator analytically; they also demonstrate empirically that their approach reduces portfolio turnover, whereas it slightly improves the out-of-sample performance. fan2019robust {construct} an elementwise covariance estimator through an M-estimation procedure with Huber loss, providing statistical high-probability guarantees. Robust portfolio optimization problem has gained significant attention, xidonas2020robust categorize 148 {studies} conducted during the last 25 years and focused on this topic.

Failures of MV portfolios become even more pronounced with a growing investment universe, especially for cases when a sample size is less than the number of assets. Evidence was investigated by kan2007optimal, Bai2009, Karoui2012, and Chen2016. To overcome this curse of dimensionality, structured covariance matrix estimators are proposed for asset return data. fan2008high considered estimators based on factor models with observable factors. stock2002forecasting, bai2012statistical, fan2013large studied covariance matrix estimators based on latent factor models. ledoit2003improved, ledoit2004well, ledoit2004honey proposed linear and ledoit2017nonlinear non-linear shrinkage of sample eigenvalues. These estimators are commonly based on the sample covariance matrix, and sub-Gaussian tail assumptions are required to guarantee sharp non-asymptotic bounds (see Section 4.7 in vershynin2018high).

The goal of our robustifying Markowitz approach is stabilizing portfolio weights by solving two problems at the same time: How to optimize the GMV portfolio when the dimension is possibly higher than the sample size and the distribution of the returns has heavy tails? Moreover, even if the returns are not heavy-tailed, how can one avoid the usual Gaussian assumption in the theoretical analysis?

Our theoretical and algorithmic contributions dwell on some recent breakthroughs in statistical literature with regard to robust estimation. lugosi2019sub constructed a multivariate mean estimator based on the idea of median-of-means that dates back to nemirovsky1983problem. The remarkable property of their estimator is that it pertains favorable properties of the Gaussian sample mean: it allows deviation bounds with high probability without much loss in the accuracy. Their only condition is that the second moment of each component of the random vector is bounded, which is the minimal possible condition to have a square-root convergence, even on average. However, their original estimator was not computationally tractable and in the past years the problem has attracted a lot of attention. hopkins2018sub first proposed an estimator with polynomial computational (efficient) complexity, and subsequent research led to nearly linear-time algorithms depersin2019robust, hopkins2020robust, thus making practical applications possible. As for the covariance estimation, mendelson2020robust proposed an abstract algorithm that achieves performance of Gaussian sample covariance estimator under four bounded moments assumption. So far, it remains an open question whether such performance can be achieved with a polynomial algorithm, with some conjecturing that the answer is no, cherapanamjeri2020algorithms.

We bypass these algorithmic problems appearing in robust estimation of the covariance matrix. In fact, our approach does not require estimating the covariance operator directly. It is based on a simple iterative gradient descent, that requires estimating only the action of the covariance operator on a current approximation at each step.

Our contribution to robustifying Markowitz is threefold:

itemize• Based on the algorithm from hopkins2020robust, we introduce a robust and computationally tractable algorithm that achieves nearly Gaussian performance under only four moments assumption on the distribution of the return vector. This means, in particular, that the estimator works with almost any distribution with four bounded moments as good as it works with Gaussian data. • We provide theoretical guarantees for our method in two cases. In the first case, we assume that the covariance matrix is well-conditioned, which means that the objective of our optimization problem enjoys the strong convexity property, and the convergence guarantees are provided even for the mean-variance objective. However, in that case we require that the dimension of the investment universe ($N$) is much smaller than the size of the sample ($T$). Moreover, the assumption that the covariance is well-conditioned is impractical due to the presence of strong factors in financial panel data, and we provide this result merely out of theoretical curiosity. In the general case where we have no control over the small eigenvalues of the covariance matrix, we only consider the GMV objective. However, we can take advantage of possibly small effective rank of the covariance matrix, which allows the dimension $N$ to be a lot larger than the sample size $T$. • In our empirical study, we compare behavior of the proposed portfolio estimator to the traditional portfolio benchmarks on equity data for two cases: when the size of the sample $T$ is comparable with the dimension $N$, and when $T<N$. For the first case, we consider the S&P100 data, and for the second case, we take the constituents of the Russell3000 index. In both cases, we conduct study for daily data over the course of 12 months and 24 months, which corresponds to $T = 252$ and $T = 500$. We demonstrate that our approach enjoys more stable weights than the traditional portfolios, while preserving (or slightly improving) their out-of-sample performance.

Let us also recall a well known hypothesis of Green1992 that extreme portfolio weights appear not entirely due to high estimation errors, but rather due to the population optimal portfolios themselves having extreme weights and being poorly diversified. Specifically, they show that asset returns generated by a model with a single dominant factor result in excessive short and long positions. This {leads} to {the study of} restricted portfolio policies. In a seminal work, jagannathan2003risk consider portfolios with non-negative constraints. Despite considering a lesser class of portfolios, they demonstrate a better out-of-sample performance. Furthermore, fan2012vast introduce Gross Exposure Constraints, which work similarly but allow negative allocation weights. Contrary to these ideas, we demonstrate that applying our robust procedure leads to desirable properties of weights without any constraints enforced a priori, which contradicts the original hypothesis of Green and Hollifield.

Notation

Throughout the paper we write $ a \lesssim b $ and $ b \gtrsim a $ if there is a constant $C$ such that $ a \leq C b $. If we have both $ a \lesssim b $ and $ b \lesssim a $, we write $ a\sim b $.

For a vector $ x \in \mathbb{R}^{d} $, we denote by $ \| x \| = \sqrt{x_1^2 + \dots + x_d^2} $ its Euclidean norm. If $ A $ is a matrix, we denote $ \| A \| = \sup_{u, v \in \mathbb{S}^{d - 1}} u^{\top} A v $ its spectral norm, where $\mathbb{S}^{d - 1}$ is a sphere in $\mathbb{R}^{d}$. If $ A \in \mathbb{R}^{d \times d} $ is symmetric, we write $ \lambda_1(A) \geq \lambda_2(A) \geq \dots \geq \lambda_{d}(A) $ to denote {its} eigenvalues in descending order. We also denote $ \lambda_{\max}(A) = \lambda_{1}(A) $ and $ \lambda_{\min}(A) = \lambda_{d}(A) $ --- its largest and smallest eigenvalues, respectively. In particular, we have that $ \| A \| = \max\{ | \lambda_{\max}(A)|, |\lambda_{\min}(A)|\} $. We say that a symmetric matrix $ A $ is positive semi-definite (PSD) if $ v^{\top} A v \geq 0 $ for all $v$. We also write $ A \preceq B $ and $ B \succeq A $ if $ B - A $ is a PSD {matrix}. Furthermore, $I$ denotes the identity matrix whose dimensions are clear from the context. We denote $ \mathbf{1} = (1, \dots, 1)^{\top} $ of dimension $N$, so that $ w^{\top} \mathbf{1} = \sum_{i = 1}^{N} w_i $.

For a PSD matrix $A \in \mathbb{R}^{d \times d}$ we {denote its effective rank by \[ \mathbf{r}(A) = \operatorname{Tr}(A) / \lambda_{\max}(A) = \sum_{j = 1}^{d} \lambda_{j}(A) / \lambda_{\max}(A). \] The effective rank is clearly always smaller than both the dimension $d$ and the matrix rank of $A$. This quantity plays an important role in covariance estimation problems. In particular, it was shown by koltchinskii2017concentration that in the Gaussian case, the performance of the sample covariance matrix is governed by the effective rank of the covariance matrix and is not sensitive to a potentially larger dimension of the ambient space.}

Mean-Variance and Global Minimum Variance portfolios

Suppose we have an opportunity to invest into $N$ assets and $ r_{1}, \dots, r_{N} $ denote their log-returns. Let $ X = (r_1, \dots, r_{N})^{\top} $ be the multivariate return vector with mean $ \mu $ and covariance $\Sigma$. Then a portfolio with allocation weights $w = (w_1, \dots, w_N)^{\top}$ has returns with expectation $ \mu^{\top} w $ and variance $ w^{\top} \Sigma w $.

One of the fundamental portfolio policies, the mean-variance portfolio (MV), is based on maximizing the utility

equation*[equation* omitted — 160 chars of source]

which takes as input the mean $ \mu $ and the covariance operator $\Sigma$. Moreover, $ \gamma $ is a fixed risk aversion parameter provided by the investor. The quadratic term in the above expression represents the variance of the portfolio return $ \operatorname{Var}(w^{\top} X) $, and the linear term is its mean ${\mathsf{E}} (w^{\top} X)$.

{Some researchers} often discard the dependence on mean and concentrate on an alternative portfolio policy that minimizes the risk measure

equation[equation omitted — 134 chars of source]

which corresponds to finding a global minimum variance portfolio (GMV). The quantity $ w^{\top} \Sigma w = \operatorname{Var}(w^{\top} X) $ is often regarded as risk of a portfolio allocation $w$ in the financial literature.

Suppose we have an i.i.d. sample $ X_1, \dots, X_T $ that comes from a distribution with mean $\mu$ and covariance $\Sigma$. Our goal is to construct portfolio allocation weights $\hat{w}$ that are as close to optimum as possible. Below we provide theoretical high-probability guarantees in terms of the gap between the estimator and the population optimal solution, that is, \[ R(\hat{w}; \Sigma) - \min_{w^{\top}\mathbf{1}} R({w}; \Sigma) \] in the GMV case, or \[ \max_{w^{\top}\mathbf{1}} M_{\gamma}({w}; \mu, \Sigma) - M_{\gamma}(\hat{w}; \mu, \Sigma) \] for the mean-variance portfolio. Notice that both of these entities are non-negative.

We analyze two different situations, focusing on high-dimensional non-asymptotic bounds. Firstly, we consider the hypothetical situation where the covariance matrix $\Sigma$ is well-conditioned, that is when the condition number $\kappa \overset{\mathrm{def}}{=} \lambda_{\max}(\Sigma) / \lambda_{\min}(\Sigma)$ is constant. Such situation is unlikely in practice, and we present the following result mainly out of theoretical curiosity. {

theorem*[A simplified statement] Suppose that $ \Sigma $ is well-conditioned, i.e., its condition number is bounded by an absolute constant. Fix $\delta \in (0, 1)$. There is a computationally efficient estimator $\hat{w}_{\delta}$ that satisfies, with probability at least $1 - \delta$, \[ \max_{w^{\top}\mathbf{1} = 1} M_{\gamma}(w; \mu, \Sigma) - M_{\gamma}(\hat{w}_{\delta}; \mu, \Sigma) \lesssim \|\Sigma\|\frac{N \log N + \log(1/\delta)}{T} \,, \] even when the distribution is non-Gaussian and has heavy tails.

} The above result demonstrates that the MV portfolio can be robustly estimated as long as the ratio $ N \log N / T$ remains small. It is known that for convergence to the optimal risk one has to have that $ N / T = o(1) $. For instance, karoui2013realized considers a “large $N$, large $T$” situation where $N/T$ converges to some constant $ \gamma \in (0, 1)$, and shows that there is a constant gap between the realized risks of the empirical and the optimal solutions. More recently, bartl2021monte studied a similar portfolio optimization setup in a well-conditioned case. Although their algorithm is robust with respect to heavy-tailed data and achieves similar rates of convergence, their estimator is not computationally feasible.

Secondly, we consider the case where $\Sigma$ is allowed to be ill-conditioned. {In this case, the condition number $\kappa$ can be large and the effective rank $\mathbf{r}(\Sigma) = \operatorname{Tr}(\Sigma)/\|\Sigma\|$ can be much smaller than $N$.} This case corresponds to a less regular optimization problem and we provide slower convergence guarantees with respect to the number of observations. {

theorem*[A simplified statement] There is a computationally efficient estimator $\hat{w}_{\delta}$, such that, with probability at least $1 - \delta$, \[ R(\hat{w}_{\delta}; \Sigma) - \min_{w^{\top}\mathbf{1} = 1} R(w; \Sigma) \lesssim \|\Sigma\|\sqrt{\frac{\mathbf{r}(\Sigma) \log \mathbf{r}(\Sigma) + \log(1/\delta)}{T}} \, , \] even when the distribution is non-Gaussian and has heavy tails.

} This result suggests that the GMV portfolio converges to optimum as long as $ \mathbf{r}(\Sigma) $ is much smaller than $ T $, which is a rather adequate assumption. For example, for the S&P100 dataset, we evaluate that $ \mathbf{r}(\Sigma) \approx 3 $ and for the Russell3000 constituents, $ \mathbf{r}(\Sigma) \approx 7 $. Moreover, implying that effective rank is also a measure of effective dimensionality of the investment universe, this measure also reflects the conclusions of widely used three and five Fama-French factor models, who proved that a cross-sectional variation in average stock returns can be explained with only three fama1992cross extended later to five observable factors, fama2015five. Notice also that the condition number $\kappa = \lambda_{\max}(\Sigma) / \lambda_{\min}(\Sigma)$ is bounded from below by $ N / \mathbf{r}(\Sigma) $, hence the covariance matrix is indeed ill-conditioned in these two applications.

Recent advances in robust statistics

The covariance matrix and the mean are not known in practice and must be estimated based on the observed log-returns.

In an idealized situation where $ X_i \sim \mathcal{N}(\mu, \Sigma) $ are Gaussian, we have that the standard empirical mean estimator $ \hat{\mu} = T^{-1} \sum_{i = 1}^{T} X_i $ provides one with optimal high probability deviation bounds. In particular, for {any} $ \delta \in (0, 1) $, we have that, with probability at least $ 1- \delta$,

equation[equation omitted — 167 chars of source]

{See Example 5.7 in boucheron2013concentration for this derivation. The sharp deviation term $\sqrt{\frac{2\| \Sigma \| \log(1/\delta)}{T}}$ is very specific to the Gaussian assumption and could not be expected for less regular distributions. In particular, } here the dependence on the confidence level is logarithmic and additive, in the sense that the bound separates into the strong term scaled with $\sqrt{\operatorname{Tr}(\Sigma)}$ and corresponding to the error on average, and the weak term that is scaled with $\sqrt{\|\Sigma\|}$. The weak term can potentially be a lot smaller than the strong one, even for very small values of $\delta$.

Similarly, koltchinskii2017concentration proved that in the case of i.i.d. zero mean Gaussian observations, the sample covariance $ \hat{\Sigma} = T^{-1} \sum_{i = 1}^{T} X_iX_i^{\top}$ satisfies the following deviation bound. For any $\delta \in (0, 1)$, with probability at least $1-\delta$,

equation[equation omitted — 144 chars of source]

whenever $T \gtrsim \mathbf{r}(\Sigma) + \log(1/\delta)$. For the version of this inequality with explicit constants we refer to zhivotovskiy2021dimension.

Having these performance bounds in mind, one is interested if the same bounds can be achieved under milder assumptions in a computationally efficient manner. lugosi2019sub developed an estimator that matches the bound (ref) {up to multiplicative constant factors} under the only assumption that the covariance matrix exists (i.e., {the} two moments assumption). Loosely speaking, they propose to control the deviations of the median-of-means of the projections $ X_i^{\top} v $ uniformly in all directions $ v \in \mathbb{R}^N $. Based on this bound, {they came up with an estimator, which is however not practical}. Further developments were made in lugosi2021robust; see also lugosi2019mean for a thorough review on this topic.

hopkins2018sub first proposed an estimator that can be computed in polynomial time, and the time complexity was subsequently improved to nearly linear cherapanamjeri2019fast, depersin2019robust, hopkins2020robust. An alternative method called spectral sample reweighing was developed in the context of robust estimation with outliers. Given data points $ \{ x_i\}_{i = 1, \dots, k} $ the goal is to reweigh the points $x_i$ with some weights $ u_{i} \in [0, 1]$ and find a center $ \nu \in \mathbb{R}^{N} $ such that the largest eigenvalue of the weighted covariance $ \sum_{i} u_i (x_i - \nu)(x_i - \nu)^{\top} $ is small. \citet*{hopkins2020robust} develop an algorithm that does this in nearly linear time; see also diakonikolas2017being, zhu2020robust. More importantly for us, \citet*{hopkins2020robust} establish a direct connection between the sample reweighing and the method developed in lugosi2019sub, which makes this approach applicable in the heavy-tailed {setup} as well. We discuss this connection in greater detail in Section (ref).

{The problem of robust covariance estimation is more challenging. mendelson2020robust construct an abstract estimator that matches the bound (ref) up to some logarithmic factors and a more recent paper abdalla2022covariance removes the remaining logarithmic factor in the bound. Unfortunately, it is still not known whether such performance can be achieved with a computationally efficient algorithm. Existing computationally efficient implementations are showing suboptimal statistical guarantees ke2019user, ostrovskii2019affine, cherapanamjeri2020algorithms and sometimes require additional assumptions such as the so-called SoS-hypercontractivity that are hard to verify in the non-Gaussian situation. Moreover, cherapanamjeri2020algorithms conjecture that as long the median-of-means approach is used, it is algorithmically hard to robustly estimate the sample covariance matrix in the presence of heavy-tailed data.}

After this short excursion to some recent results in robust estimation, let us now come back to our portfolio optimization problem. Since we cannot get our hands on a robust covariance estimator, we take another route by observing that both MV and GMV are convex optimization problems.

Gradient descent for portfolio optimization

Our goal is to avoid the estimation of the whole covariance matrix, but rather resort to the estimation of the action of this operator $ \Sigma w $ on some limited set of vectors $w$. We will use a procedure based on projected gradient descent (PGD), which is a standard convex optimization method. For instance, if we want to minimize the GMV objective with known $ \Sigma$, the following sequence of approximations converges to an optimal solution (which is not necessarily unique): we start with arbitrary initial vector $w_0$ and then take the update steps,

equation[equation omitted — 133 chars of source]

where $\Pi_{1}$ is the orthogonal projector onto the restricted (convex) set $ \{w: \; w^{\top}\mathbf{1} = 1\} $, which can be explicitly defined by the mapping,

equation[equation omitted — 127 chars of source]

It is straightforward to see that $ R(\,\cdot\,; \Sigma)$ is convex (since the covariance operator is positive semi-definite), and $\| \Sigma\|$-smooth {(we say that a continuously differentiable function $f$ is $\beta$-smooth if the gradient $\nabla f$ is $\beta$-Lipschitz, that is, $\|\nabla f(x) - \nabla f(y)\| \le \beta \|x - y\|$)}. By Theorem {3.7} from bubeck2014convex the sequence (ref) converges to a minimum at a rate $1/s$ as long as $ \eta \leq 1/ \| \Sigma \|$. Moreover, in the case where $ \Sigma$ is non-degenerate, the objective becomes strongly convex, and the sequence converges at a faster exponential rate under the same requirement on the step size, see Theorem 3.10 in bubeck2014convex.

The case of MV portfolio is similar, only this time we need to maximize a concave function instead of minimizing a convex one. If we replace $\nabla_{w} R(w_{s-1}; \Sigma)$ with $ -\nabla_{w} M_{\gamma}(w_{s-1}; \Sigma)$ in (ref), then by the same reasons, the sequence converges to the maximum of $M_{\gamma}$ as long as $ \eta \leq (\gamma \| \Sigma \|)^{-1}$, with exponential rate when $\Sigma$ is non-degenerate.

The PGD iterations require computation of the following gradients, \[ \nabla_{w}R(w; \Sigma) = \Sigma w \qquad \text{or} \qquad \nabla_{w}M_{\gamma}(w; \Sigma) = \mu - \gamma \Sigma w. \] where the mean $\mu$ and covariance $\Sigma$ are typically replaced with their empirical counterparts that are calculated using given historical observations $X_1, \dots, X_T$. As discussed in the previous section, there is a practical robust mean estimator in hopkins2020robust with all desired properties. Since such {an} estimator is not available for the covariance operator, we instead produce an estimator $\hat{a}_{\delta}(w)$ for the PGD increment that estimates $ \Sigma w$ for each $w$ separately, and plug it into the update steps of PGD.

To see how it can be done, suppose for a moment that the expectation vanishes. Then, we can represent this product as a mean of a random vector as follows, \[ \Sigma w = {\mathsf{E} XX^{\top} w\, }. \] We therefore can apply the robust mean algorithm to the vectors {$X_iX_i^{\top}w$} and obtain a robust estimator of $\Sigma w$. However, we need to take additional care to ensure that the estimator is an appropriate approximation uniformly in all directions. For this, we slightly adjust the procedure in the spirit of mendelson2020robust. In the latter work, the only assumption used is the equivalence of the fourth and the second moments in all directions {sometimes called the bounded kurtosis assumption.}

assumption[Bounded kurtosis] The return vectors $ X_1, \dots, X_T $ are i.i.d. observations of a random vector $X $, that has mean $\mu$, covariance $\Sigma$, and satisfies for all $ u \in \mathbb{R}^{N} $, \begin{equation} \mathsf{E}^{1/4} |u^{\top}(X - \mu)|^{4} \leq K \mathsf{E}^{1/2} |u^{\top}(X - \mu)|^{2}, \end{equation} where $ K \geq 1$ is some fixed constant. For the rest of the paper, we ignore the explicit dependence on $K$ in our bounds treating $K$ as an absolute constant.

All our results will be stated under this assumption. Importantly, no other assumptions are made on the covariance matrix $\Sigma$ throughout the paper. This makes our setup quite different from traditional factor models where one assumes some specific structure of $\Sigma$. The assumption is clearly satisfied for Gaussian random vectors with $K$ being an absolute constant. However, we even cover some heavy-tailed distributions. {For some heavy-tailed examples where Assumption (ref) is satisfied we refer to mendelson2020robust}. In the context of covariance estimation, a robust estimator is one that has Gaussian deviation bounds (i.e., as in (ref)) but only requires the underlying distribution to follow the bounded kurtosis assumption. The next step is to provide a robust estimator of $ \Sigma w $ that works simultaneously in all directions. Recall that $\mathbf{r}(\Sigma)$ denotes the effective rank $ \mathbf{r}(\Sigma) = \operatorname{Tr}(\Sigma) / \| \Sigma\| \, . $

propositionSuppose that Assumption (ref) holds. Fix $\delta \in (0, 1)$ There is a computationally efficient estimator $ \hat{a}_{\delta}(w) $ that depends on direction $w$ and $T$ i.i.d. observations and such that, with probability at least $1-\delta$, \[ \| \hat{a}_{\delta}(w) - \Sigma w\| \lesssim \| \Sigma \| \sqrt{\frac{\mathbf{r}(\Sigma) \log \mathbf{r}(\Sigma) + \log(1/\delta)}{T}} \| w\|, \] uniformly for all vectors $w$.

We postpone the proof and detailed description of the estimation algorithm until Section (ref). {We note that its computational complexity is $\widetilde{O}(TN + T^{1/2} N \log(1/\delta)^{2})$, where the notation $ \widetilde{O} $ suppresses some multiplicative logarithmic factors.}

remarkWe expect that the proposed estimator is robust to finite number of outliers as well, as typically happens with estimators based on lugosi2019mean. To avoid technically involved proofs, we do not check this rigorously, instead we demonstrate robustness to outliers in the simulation study in Section (ref).
remarkFor technical reasons, the estimator $ \hat{a}_{\delta}(w) $ depends on a {norm-truncation} parameter $ D $ that needs to be of order $ \left(\frac{T \mathbf{r}(\Sigma)}{\log \mathbf{r}(\Sigma)}\right)^{1/4} $, which is unknown in general. {It appears that since $D$ increases with $T$, in many natural situations this truncation parameter is of a much larger order than {$\max\limits_i\|X_i\|$} and can be mostly ignored in practice. For more details see Section (ref).}
remarkWe point out that when one has access to some covariance estimator $\hat{\Sigma}$, one can simply take a family of estimators $ \hat{a}(w) = \hat{\Sigma} w $. For instance, in the Gaussian case, taking the standard empirical covariance estimator would yield thanks to (ref), \[ \| \hat{a}(w) - \Sigma w\| \lesssim \| \Sigma \| \sqrt{\frac{\mathbf{r}(\Sigma) + \log(1/\delta)}{T}} \| w\| \ , \] with probability at least $1 - \delta$, uniformly for all $w$. The estimator of Proposition (ref) achieves the same rate of convergence under minimal distributional assumptions.

Now we can plug in the estimator of $\Sigma w$ (appearing in Proposition (ref)) into the update rule (ref). To be precise, in the case of GMV optimization, our updates look as follows, \[ w_{s} = \Pi_{1} \left[ w_{s-1} - \eta \hat{a}_{\delta}(w_{s-1}) \right], \qquad s = 1, 2, \dots \] Naturally, the error may accumulate with each update, and we need to carefully analyze how the resulting solution differs from the optimum, to which the sequence (ref) converges. We analyze this update rule in two separate cases.

First, we consider the case of a well-conditioned matrix $ \Sigma $, meaning that the ratio of its maximal and minimal eigenvalues is bounded by a constant. This means that the problem of maximizing the MV utility is a strongly-convex optimization problem, so that the gradient descent sequence enjoys exponential convergence rate and, as we show below, the error of estimation does not accumulate. However, in that situation the effective rank $\mathbf{r}(\Sigma) = \operatorname{Tr}(\Sigma) / \| \Sigma \|$ is of order $N$, so the convergence only works in the case where $ N /T $ is small. Moreover, in typical applications, the covariance matrix is ill-conditioned, which is one of the reasons the MV portfolio performs so poorly. For instance, this can be checked through evaluation of the effective rank: for the S&P100 dataset we estimate $ \mathbf{r}(\Sigma) \approx 3 $ and for the Russell3000 set we estimate that $ \mathbf{r}(\Sigma) \approx 7 $, in both cases much smaller than the dimension $N$. This brings us to the second part of our GD analysis, where we only consider the case of GMV optimization with ill-conditioned covariance matrix that has small effective rank. This scenario corresponds to non-strongly convex optimization and has weaker convergence rate. However, it enjoys dimension-free bounds, meaning that the convergence is guaranteed as long as the number of observations is much larger than $\mathbf{r}(\Sigma)$, regardless of how high the total number of assets is. We also point out that in this case, one has to stop after an appropriate number of steps to avoid overfitting.

Well-conditioned case

For maximizing the MV utility $ M_{\gamma}(w; \mu, \Sigma) $, we consider the following updates,

equation[equation omitted — 158 chars of source]

where $ \hat{\mu} $ is some estimator of mean $ \mu $, and $ \hat{a}(w) $ is some family of estimators for the action of covariance operator $ \Sigma w $. We first show a deterministic result that controls the convergence through the errors of estimators $ \hat{\mu} $ and $ \hat{a}(w) $.

lemmaDenote, $ {w}^* = \arg\max_{w^{\top}\mathbf{1} = 1} M_{\gamma}(w; \mu, {\Sigma})$. Suppose that we have {an} access to an estimator $\hat{\mu}$ satisfying $$ \| \hat{\mu} - \mu\| \leq \Delta_{\mu}, $$ and {an} access to a family of estimators $\hat{a}(w)$ satisfying uniformly for all $ w\in \mathbb{R}^{N} $, $$ \| \hat{a}(w) - \Sigma w \| \leq \Delta_{\Sigma} \| w \| \, . $$ Let $ \lambda_{\max}$, $\lambda_{\min}$ denote the maximal and minimal eigenvalues of $\Sigma$, respectively. Assume that $ \eta \leq 1/(\gamma \lambda_{\max}) $ and $ \Delta_{\Sigma} < \lambda_{\min} $. Then, the sequence (ref) satisfies \[ \| w_{s} - w^{*} \| < \left(1 - {\gamma\eta(\lambda_{\min} -\Delta_{\Sigma})}\right)^{s} \| w_{0} - w^{*} \| + \frac{ {\gamma^{-1}} \Delta_{\mu} + \Delta_{\Sigma}\| w^{*}\|}{\lambda_{\min} - \Delta_\Sigma} \, . \]

We now apply this lemma to the case where we use $ \hat{\mu}_{\delta} $ from hopkins2020robust and $ \hat{a}_{\delta}(w) $ from Proposition (ref). {Theorem D.3 from hopkins2020robust shows that their estimator satisfies the sub-Gaussian bound (ref). That is,} with probability at least $1-\delta$,

equation[equation omitted — 147 chars of source]

Furthermore, by Proposition (ref), with probability at least $1-\delta$, simultaneously for all $w \in \mathbb{R}^N$,

equation[equation omitted — 182 chars of source]

{Substituting the two error terms (ref), (ref) into Lemma (ref), we arrive at the following result. We postpone the derivations to Section (ref).}

theoremSuppose, we are given independent $ X_{1}, \dots, X_T $ that have mean $\mu$ and covariance $\Sigma$, and the distribution satisfies the bounded kurtosis assumption (ref). Let $ \kappa = \lambda_{\max}(\Sigma) / \lambda_{\min}(\Sigma) $ {denote} the condition number. There is an absolute constant $ C>0$, such that the following holds. {If for $\delta \in (0, 1)$, the following holds} \begin{equation} T \geq C \kappa^{2} \left(\mathbf{r}(\Sigma)\log\mathbf{r}(\Sigma) + \log(1/\delta)\right), \end{equation} {then} there is an estimator $ \hat{w}_{\delta} $ such that, with probability at least $1-\delta$, \[ \max_{w^{\top}\mathbf{1} = 1} M_{\gamma}(w; \Sigma, \mu) - M_{\gamma}(\hat{w}_{\delta}; \Sigma, \mu) \lesssim \kappa^{2} \left({\gamma^{-1} + \gamma} \| \Sigma\| \|w^{*}\|^{2} \right) \frac{\mathbf{r}(\Sigma)\log\mathbf{r}(\Sigma) + \log(1/\delta)}{T} \, . \]

The above result has a number of favorable properties:

itemize• The estimator only requires $\mathcal{O}(\log T)$ gradient descent updates. In addition, the amount of steps only has to be sufficiently large, i.e., there is no danger of overfitting by running the gradient descent for too long; • the bound scales with $1/T $ when all the other parameters are fixed. In the optimization {and statistics} literature, this is regarded as a fast rate convergence. This rate is typical for strongly convex stochastic optimization problems; • The value {$\| w \|^{2}$} is often considered as a diversification measure of an allocation strategy, see strongin2000beating. For instance, for the EW portfolio {this} value is $1/N$. One may expect it to be very small for the optimal portfolio.

However, the dependence on the condition number of the covariance matrix outweighs some of these useful properties. It is straightforward to verify that $ \kappa \mathbf{r}(\Sigma) \geq N $. Hence, the above result only works in the setting, where the number of observations $T$ is greater than the dimension $N$. The remaining term $ \kappa $ may further worsen the bound, so our result is rather limited to well conditioned covariance matrices. Unfortunately, it is rarely the case in practice: for our two datasets we estimate that $ \mathbf{r}(\Sigma) \approx 3 $ for S&P100 and $ \mathbf{r}(\Sigma) \approx 7 $ for Russell3000. The naive lower bound $\kappa \geq N/\mathbf{r}(\Sigma)$ yields that $ \kappa \geq 27 $ for S&P100 and $ \kappa > 250 $ for Russell3000. Therefore, our result does not contradict a commonly accepted evidence that MV portfolios perform poorly even when $T$ is moderately larger than $N$ ao2019approaching.

Ill-conditioned case

We now consider the case where we have no control over the condition number of $\Sigma$ and it can even be degenerate. We will state the bound in the regime where only the effective rank $ \mathbf{r}(\Sigma) $ has to be small, and no requirements on the total dimension $N$ are needed. In this section, we only consider the GMV portfolio.

For minimizing the GMV risk $ R(w; \Sigma) $, we consider the following updates,

equation*[equation* omitted — 113 chars of source]

where $ \hat{a}(w) $ is some family of estimators for the action of covariance operator $ \Sigma w $.

Similarly to the previous section, we first show a deterministic result that controls the convergence through the error of this estimator.

lemmaDenote, $ w^{*} = \arg\min_{w^{\top}\mathbf{1} = 1} R(w; \Sigma)$. Suppose that we have {an} access to a family of estimators $\hat{a}(w)$ satisfying uniformly for all $ w\in \mathbb{R}^{d} $, $$ \| \hat{a}(w) - \Sigma w \|_{2} \leq \Delta_{\Sigma} \| w \|_2 \, . $$ Assume that $ \eta \leq 1/\lambda_{\max} $ and let the number of {weights updates $S$ satisfy $ S \eta \Delta_{\Sigma} \leq 1$}. Then, \[ R(w_{S}; \Sigma) - R(w^{*}; \Sigma) \lesssim \max\{ \| w_{0} - w^{*}\|, \| w^{*}\|\}^{2} \left( \frac{1}{\eta S} + \eta \Delta_{\Sigma}^{2} S \right) . \] {In particular, }for the optimal choice $ S \sim 1 / (\eta \Delta_{\Sigma}) $, we have \[ R(w_{s}; \Sigma) - R(w^{*}; \Sigma) \lesssim \max\{ \| w_{0} - w^{*}\|, \| w^{*}\|\}^{2} \Delta_{\Sigma} . \] When the solution to $\min_{w^{\top}\mathbf{1} = 1} R(w; \Sigma)$ is not unique, it is sufficient to pick any $ w^{*} \in \mathrm{Argmin}_{w^{\top} 1 = 1} R(w; \Sigma) $ and the result still holds.

Once again, we plug our estimator $ \hat{a}_{\delta}(w) $ into the update rule. In addition, we take the initial approximation to be {an} EW portfolio. Namely, our sequence {is} as follows

equation[equation omitted — 166 chars of source]
corollarySuppose that Assumption (ref) holds. Take $ \eta = 1/\lambda_{\max} $ and $ S \sim T (\mathbf{r}(\Sigma) \log \mathbf{r}(\Sigma) + \log(1/\delta))$, and set $ \hat{w}_{\delta} = w_{S} $. Then, with probability at least $1- \delta$, \[ R(\hat{w}_{\delta}; \Sigma) - R(w^{*}; \Sigma) \lesssim \| \Sigma \| \| w^{*}\| ^{2} \sqrt{\frac{\mathbf{r}(\Sigma) \log \mathbf{r}(\Sigma) + \log(1/\delta)}{T}} \, . \]
proofSimply substitute the bound (ref) into Lemma (ref). We also notice that $ \| w^{*} - w_0\|^{2} \leq 2\| w^{*}\|^{2} + 2\|w_0\|^{2} $, and that $ \| w^{*}\|^{2} \geq (\mathbf{1}^{\top} w^{*})^{2} / N = \| w_{0}\|^{2} $.
remarkWe remark that the scaling value $ \| \Sigma \| \| w^{*}\|^{2} $ is only an upper bound on the optimal risk $ R(w^{*}; \Sigma) = \frac{1}{2} (w^{*})^{\top} \Sigma w^{*}$ and we cannot guarantee a ratio-type bound of the form $ R(\hat{w}_{\delta};\Sigma) = (1 + o(1))R(w^{*}; \Sigma) $. However, this is not uncommon. For instance, fan2012vast {show} that a portfolio with GEC {(Gross Exposure Constraints)} constraints $ \sum_{i} |w_i| \leq C $ satisfies, \[ R(\hat{w}; \Sigma) - R(w^{*}; \Sigma) \leq (1 + C)^{2} \max_{ij} | \hat{\Sigma}_{ij} - \Sigma_{ij}|, \] where one typically has a bound $\max_{ij} | \hat{\Sigma}_{ij} - \Sigma_{ij}| \lesssim \sqrt{(\log N) / T}$. Our bound {is more beneficial when the optimal portfolio} is well-diversified (i.e., $\| w^{*}\|^{2} \sim 1/N$), even though we do not impose any restrictions on the selected portfolio.

Evaluation of empirical results

To test the performance of our approach, we apply it to two data sets of stocks. The first data set consists of 81 constituents of S&P100 index (as on January 1, 2021) and covers time span from January 2, 2000 to December 31, 2020 summing up to 5282 daily log-returns. These 81 stocks have a continuous return time-series over the period of our study. The second data set consists of 600 random constituents of Russell3000 index as on January 1, 2021, period of analysis is limited by 11 years: from January 2, 2010 to December 31, 2020. The length of analyzed time series is 2768 observations.

For the portfolio construction, we employ a rolling-window approach with monthly rebalancing. Specifically, we choose an estimation window of length $T$ days starting on date $T + 1$, for each rebalancing period $l$ (${l = 1, \ldots, L}$, with $L$ the number of rebalancing periods) we use the data in the previous $T$ days to estimate the parameters required to implement a particular strategy. The input parameters are estimated using daily returns of the most recent 12 and 24 months, corresponding roughly to 252 (500) daily returns of past data (with the length of estimation windows $T=252$ and $T = 500$). Results for the window length 24 months are discussed in Section (ref) and summary tables for the length window 12 months are provided in Section (ref). Thus, the out-of-sample period for the S&P100 data set starts on January 2, 2002 with 4781 observations, which corresponds to the number of rebalancing periods $L = 216$, and for the Russell3000 data set --- on January 3, 2012 with 2265 out-of-sample observations corresponding to $L = 108$. The source for both data sets is Thompson Reuters.

Benchmark portfolios

Here we present the empirical results for GMV portfolio and evaluate its relative performance. The allocation rules included into the empirical study with corresponding reference and abbreviations are listed in Table (ref).

\paragraph{Equally weighted (EW).} demiguel2009portfolio argue that a naive allocation strategy with weights $w_i = 1/N$ is hard to outperform in practice. It is often used as a benchmark for comparative analysis.

\paragraph{Sample-based Global minimum portfolio (GMV).} This is the most straightforward way to GMV optimization. The sample covariance matrix $\hat{\Sigma}$ is plugged into the objective in (ref). We should note that this strategy is only included for S&P100 data set, since for the Russell3000 we have $N > T$, and the sample covariance matrix is not invertible.

\paragraph{Global minimum portfolio with short-sale constraint (GMV_long).} This portfolio is a sample-based GMV portfolio with only long positions allowed. This means that GMV objective corresponding to the empirical covariance is optimized subject to the constraints $ w_j \geq 0 $.

\paragraph{Global minimum portfolio with linear shrinkage estimator (GMV_lin).} ledoit2004well propose an asymptotically optimal convex linear combination of the sample covariance matrix $\hat{\Sigma}$ with the identity matrix. Optimality is meant with respect to a quadratic loss function, asymptotically, as the number of observations and the number of assets go to infinity together. ledoit2004well use as a covariance matrix estimator as a convex linear combination of the sample covariance matrix and the identity matrix (shrinkage target) as follows:

equation*[equation* omitted — 75 chars of source]

where $\rho$ is the shrinkage intensity parameter and $\hat{\Sigma}$ is the sample covariance matrix. \href{https://www.econ.uzh.ch/en/people/faculty/wolf/publications.html}{Their R code} is used in this horse race exercise. \paragraph{Global minimum portfolio with non-linear shrinkage estimator (GMV_nlin).}

ledoit2017nonlinear use the spectral decomposition for the empirical covariance

equation*[equation* omitted — 82 chars of source]

where $\widehat{D} \overset{\mathrm{def}}{=} \operatorname{diag}\left(\widehat{d} \left(\lambda_{ 1}\right), \ldots, \widehat{d}\left(\lambda_{N}\right)\right)$, where $\lambda_{1}, \dots, \lambda_{N}$ are the sample eigenvalues, and $\hat{d}$ is some nonlinear cutoff threshold based on value of $N/T$ and the magnitude of the eigenvalues $\lambda_j$.

table[table omitted — 812 chars of source]

Performance measures

We report the following five out-of-sample performance measures for each benchmark portfolio rule.

itemize• Turnover (TO) The main practical objective of the introduced methodology is stabilizing of portfolio weights, aiming at reduction of transaction costs. To assess the impact of potential trading costs associated with portfolio rebalancing, we calculate two measures for turnover. First, following demiguel2009optimal and demiguel2009portfolio, we present Turnover, which is defined as an average sum of the absolute value of the rebalancing trades across the $N$ assets of the investment universe and over the $L$ rebalancing months (ref). \begin{equation} TO=L^{-1} \sum_{l=1}^{L} \sum_{j=1}^{N}\left|\hat{w}_{ j, l+1}-\hat{w}_{j, l+}\right|\, . \end{equation} where $\hat{w}_{j,l}$ and $\hat{w}_{j,l+1}$ are the weights assigned to the asset $j$ for rebalancing periods $l$ and $l + 1$ and $\hat{w}_{j,l+}$ denotes its weight just before rebalancing at ${l + 1}$. Thus, one accounts for the price change over the period, as one needs to execute trades to rebalance the portfolio towards the ${w_l}$ target. High turnover will imply significant transaction costs; consequently, the lower TO of a strategy, the less its performance would be harmed by non-zero transaction costs. Whereas, for market makers, the reduced turnover is not a crucial parameter in a decision-making process, for buy-side investors, we believe, that a lower turnover could be an important advantage. • Target Turnover (TTO) Further, following petukhina2021investing we also calculate a target turnover, which is constructed as follows: \begin{equation*} {TTO} = L^{-1}\sum_{l=1}^{L}{\sum_{j=1}^N}\vert \hat{w}_{j,l+1}-\hat{w}_{j,l}\vert \, . \end{equation*} In contrast to equation (ref) this definition of turnover implies by construction a value of zero for the EW portfolio. We provide this measure to focus on modifications of the target portfolio weights due to active management decisions and cleaned from the influence of assets’ price dynamics. We calculate the following four performance measures for both, gross and net, returns series. Gross returns are raw returns that do not take transaction costs into account and net returns with subtracted fees. When the portfolio is rebalanced at time $l+1$, it gives rise to a trade in each asset of magnitude $\left|\hat{\mathrm{w}}_{k, j, l+1}-\hat{\mathrm{w}}_{k, j, l^{+}}\right|$. Thus, to calculate the net return we reduce the return by the cost of such a trade over all assets, given by $$tc = c\sum_{j=1}^N\left|\hat{\mathrm{w}}_{k, j, l+1}-\hat{\mathrm{w}}_{k, j, l^{+}}\right|,$$ where $c$ is the proportional transaction cost. In our empirical study we set a proportional transactions cost of 50 basis points per transaction. The same level was accepted in e.g., demiguel2009optimal, fays2011. • Cumulative wealth (CW) CW generated by each benchmark strategy with initial investment $W_0 = 1 $ is computed as follows for gross returns series: \begin{equation*} W_{l+1} = W_{l} + {\hat{w}_{l}^{\top}X_{l+1}}. \end{equation*} CW for net returns series: \begin{equation*} W_{l+1} = W_{l} + {\hat{w}_{l}^{\top}X_{l+1}} - tc. \end{equation*} • Standard Deviation (SD) We compute $SD$ of out-of-sample daily returns series. • Sharpe ratio (SR) \begin{equation*} SR = \frac{{AV}}{{SD}}, \end{equation*} where $AV$ and $SD$ are average out-of-sample returns and their standard deviations for each strategy are calculated based on daily out-of-sample returns series. We also test, whether GMV_{robust} delivers a better risk measured by SD, and risk-adjusted performance, measured by SR, than other portfolios at a level that is statistically significant. A two-sided $p$-values for the null hypothesis of equal SD and SR are obtained by HAC method described in ledoit2011robust for SD and in ledoit2008robust for SR. • \emph{Calmar Ratio (CR)} We also compare \emph{(CR)}, metric for risk-adjusted return, which is widely used by practitioners due to its asymmetry property: focus on the maximum drawdown not on volatility like in $SR$. \begin{equation*} CR = \frac{{252\times AV}}{{D_{P}}}, \end{equation*} where $AV$ is the average out-of-sample return multiplied by 252 to annulize and $D_P$ is the maximum drawdown of the portfolio returns.

Since the focus of this research is the reduction of portfolio weights' fluctuations, following ledoit2017nonlinear we also compute the following five characteristics of weights' vectors $\hat w_{t}$ averaged through number of rebalancing periods. Thus, we calculate minimum weight (min) for every benchmark strategy as follows:

equation*[equation* omitted — 104 chars of source]

We similarly compute maximum weight (max), standard deviation (sd), and range of weights (max-min).

In addition, we provide MAD from EW portfolio (mad-ew), which is defined as:

equation*[equation* omitted — 150 chars of source]

Empirical study

Discussion of S&P100 data set results

First, we discuss portfolio weight stability, since it is the main focus of the research. Figure (ref) demonstrates the dynamics of weights for S&P100. It can be observed that weights of plug-in GMV portfolio are characterized by a lot of extremes in comparison with all other policies. The least dispersed weights are observed for, introduced in this paper, GMV_{robust} approach. This visual result is confirmed by descriptive statistics of portfolio weights reported in the Table (ref). It can be found that the average range of weights for GMV_{robust} 0.07 is the lowest one, what is almost four times less than the range of GMV_lin, GMV_nlin and GMV_long and five times less than plug-in GMV policy; $mad-ew$ also is the lowest for GMV_{robust} strategy, pointing out the more balanced distribution of weights. Table (ref) reports these results, which can be summarized as follows.

The two main characteristics of interest for this research would be Turnover and Target Turnover. As expected, the best performing policy in this dimension is GMV_long with imposed non-negative constraints; it requires on average almost 13% (TO) of trading volume to rebalance the portfolio. GMV_{long} is followed by EW with 20% and GMV_{robust} with 25%. The highest turnover is reached by GMV with almost 71% of portfolio value to rebalance the portfolio to target weights $\hat{w_t}$. The ranking of strategies via $TTO$ stays unchanged. The sample conducted t-test, confirms significant difference of GMV_{robust} $TO$ and $TTO$ in comparison with other benchmarks.

Naturally, for investors, Cumulative Wealth ($CW$) is of high interest as a measure of performance for the period considered. The gross $CW$ for portfolio benchmarks does not vary considerably, ranging from 256% for GMV_{long} portfolio to 273% for EW. However, the inclusion of transaction fees changes the results: the best performing EW portfolio gains 249%, followed by GMV_{robust} with CW 242%. The worst performing portfolio is GMV with 188% of CW. The evolvement of gross CW of all benchmark strategies for the considered period is plotted in Figure (ref).

Shrinkage estimators exhibit the lowest risk, measured by SD and tests conducted point out the significance of this difference with GMV_robust. In terms of risk-adjusted performance, the winning strategies for Gross returns are shrinkage portfolios and GMV with SR 5.71% and 5.65% respectively. However, for net returns they lose half of their performance, ending up with 2.3%, whereas SR of our approach stays more stable changing from 4.26% to 3.2%. In addition, conducted tests of difference significance support the hypothesis of equal SRs, although it can be noticed that the difference is becoming more prominent for net returns.

Finally, GMV_{robust} dominates other benchmarks in terms of Calmar ratio. The drawdowns often determine whether a buy-side investor can keep an investment or will have to unwind and thus miss subsequent recoveries, what makes it an important risk measure for buy-side investors and managers.

Discussion of Russell3000 data set results

Outcomes of weights stability analysis are consistent with ones described for S&P100 data set. GMV_{robust} weights are characterized by harmonized weights without extreme short or long values. It is visible in the Figure (ref) and in the Table (ref): GMV_{robust} $MAD-EW$ and $max-min$ are the lowest in comparison with benchmark portfolios (excluding EW).

In terms of accumulated wealth, GMV_{robust} for large portfolios performs very close to shrinkage estimators, Table (ref) summarizes investment performance characteristics. Thus, GMV_{robust} gains in the end of the period 210 % of initial value while GMV_lin and GMV_nlin 209% and 201.9%. But considering non-zero trading fees would change the rank drastically, making the best performing benchmark Index itself, followed by GMV_{robust}. The findings for other risk measures stay in line with findings for S&P100:

itemize• the less risky portfolios with the lowest SD are shrinkage portfolios, and the difference is significant in relation to GMV_{robust}. • Shrinkage estimators demonstrate the highest SR for gross returns, but for the net returns GMV_{robust} overperforms, although difference tests do not prove the significance of the result. • GMV_{robust} dominates other strategies in terms of CR for both net and gross returns.

Thus, according to outcomes of empirical experiments we can claim that GMV_{robust} portfolio policy achieves its goal and substantially reduces fluctuations of weights, leading to the lowest level of accumulated transaction costs. The risk adjusted performance is equal or slightly lower than shrinkage benchmarks and higher than constrained rules. Our approach is of specific interest for invetors and managers who focus on the worst drawdowm. This conclusion stays robust for small and large portfolios.

table[table omitted — 1,039 chars of source]
table[table omitted — 487 chars of source]
table[table omitted — 1,086 chars of source]
figure[figure omitted — 924 chars of source]
table[table omitted — 503 chars of source]

Conclusion and discussion

“Robustifying Markowitz” has seen many attempts that are mostly based on robustifying the original simple inversion formula for exact determination of optimal GMV weights. In bypassing this “error maximizing” technique, we have presented a tool fixing the portfolio weights in a low cost re-balancing ballpark. Using modern results from robust statistics, we have constructed an algorithm that provides a computationally effective estimator for GMV Markowitz portfolios. We have shown that it suffices to utilize a PGD procedure to optimize the portfolio weights without estimating the covariance operator itself. The focus on just the PGD updates significantly distinguishes our approach from the previous techniques. We have successfully derived almost Gaussian properties of this estimator in nice ($N/T$ small) and not-so-nice ($N/T$ big) condition cases.

The weights developed with the robustified approach are less sensitive to deviations of the asset-return distribution from normality than those of the traditional minimum-variance policy. Empirical studies confirm that the proposed policies are indeed more stable and cost-reducing. The stability of the proposed portfolios makes them a feasible alternative to traditional portfolios.

The proposed toolbox improves the stability properties of weights, leading to better investment characteristics of allocation policies. Although inferior at the risk level measured by SD, our algorithm reaches a superior risk-adjusted performance (Sharpe and Calmar ratios) for net returns due to a substantial reduction of trading volume measured by turnover. Finally, these performance results are confirmed across small and large portfolios. Even for dimensions of portfolio size larger than the length of the estimation window (e.g., the Russell3000 data), the above claim pertains.

Acknowledgments

Wolfgang H\"ardle and Alla Petukhina gratefully acknowledge the financial support of the European Union's Horizon 2020 research and innovation program “FIN-TECH: A Financial supervision and Technology compliance training programme" under the grant agreement No 825215 (Topic: ICT-35-2018, Type of action: CSA), the European Cooperation in Science & Technology COST Action grant CA19130 - Fintech and Artificial Intelligence in Finance - Towards a transparent financial industry, the Deutsche Forschungsgemeinschaft's IRTG 1792 grant, Wolfgang H\"ardle - the Yushan Scholar Program of Taiwan, the Czech Science Foundation's grant no. 19-28231X / CAS: XDA 23020303. Nikita Zhivotovskiy was funded in part by ETH Foundations of Data Science (ETH-FDS). We wish to thank three anonymous referees for their thoughtful comments and efforts toward improving our manuscript. We also gratefully acknowledge the comments and discussions from Valerio Poti, Natalie Packham, and Jörg Osterrieder, as well as participants of COST FinAI Annual Meeting in Bucharest, Research seminars "Mathematics" at HTW Berlin and the Department of Economics and Management, University of Pavia. We appreciate the editorial assistance of Elie Tamer and Terry Liu.