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.
82,559 characters · 24 sections · 36 citation commands
Forecasting Large Realized Covariance Matrices: The Benefits of Factor Models and Shrinkage
\doublespacing
This paper aims to construct models based on economically motivated factor decompositions and shrinkage methods to forecast large-dimensional, and time-varying realized measures of daily covariance matrices of returns on financial assets. Realized measures of a covariance matrix are estimates, based on intraday returns, of the integrated covariance matrix of a multivariate diffusion process. One example of such an estimator used in this paper is the composite realized kernel method recently introduced by lunde2016econometric. Our proposed model is evaluated in terms of its forecasting ability and, more importantly, several performance measures in a conditional mean-variance portfolio allocation problem.
Modeling and forecasting the covariance matrix of financial assets are essential for portfolio allocation and risk management. Moreover, it is an established empirical fact that (conditional) covariance matrices vary considerably over time. A natural way to model such dynamics is to either use multivariate generalizations of the ARCH/GARCH family of models, as proposed by BEEK1995 or EngleDCC2002, or model directly a given realized measure of the covariance matrices by usual multivariate time-series models, as in Bauer2011, Chiriac2011 or golosnoy2012conditional. This is motivated by the close connection between the conditional covariance matrix and the integrated covariance of a multivariate diffusion process.
However, when the number of assets increases, the amount of parameters to be estimated becomes very large. For instance, for a covariance matrix of $N$ assets, there are $N(N+1)/2$ distinct entries to be modeled. If a vector autoregressive specification of order $p$, VAR($p$), is used, then the total number of parameters will be $(N+1)(p+1)/2$. Therefore, the curse of dimensionality precludes the application of the above-referenced methods to moderately large covariance matrices, and most of the previous studies in the literature focused on sets of less than ten assets.
More recently, based on the advances of modern statistical tools to handle high-dimensional models, new alternatives have been proposed in the literature to address a large number of assets. callot2017modeling advocate the use of the least absolute shrinkage and selection operator (LASSO) of tibshirani1996regression and the adaptive LASSO of zou2006adaptive to model the dynamics of the realized covariance of the constituents of the Dow Jones index by a large dimensional VAR model. However, their modeling strategy cannot handle sets of assets much larger than the 30 used in their paper. engle2017large combine nonlinear shrinkage with the DCC-GARCH model and put forward a methodology in which the dynamics of large latent conditional covariance matrices can be modeled. Their approach makes the estimation of multivariate GARCH-type models feasible even in large-dimensional cases. Nevertheless, the authors did not consider realized measures in their models.
In this paper, we propose a method that can be applied to a large set of assets with the difference that we consider realized covariance measures rather than latent ones, as in engle2017large. However, it is important to highlight that our method relies on realized measures of the integrated covariance matrix. Although our methodology can be applied to thousands of assets, illiquid assets may cause severe problems for the integrated covariance matrix estimation. For instance, although the estimator by lunde2016econometric, which is used in this paper, is scalable as it is composed of 2x2 matrix estimates, an illiquid stock will cause a poor estimate of that stock’s covariance with other stocks. Nevertheless, the estimates of covariances between other liquid stocks will not be affected.
We provide empirical evidence that our proposed model drastically improves the forecasting ability and produces portfolios with better performance measures. We apply our model to some of the constituents of the S&P 500 using daily information on the realized covariance matrix. The models based on our framework can improve the forecasts' quality compared to the random walk benchmark by reducing the $\ell_2$-forecast error in 15%. On the other hand, models based on latent volatility are the worst-performing specification regarding forecasting ability. Furthermore, our models deliver a significant reduction in the volatility of minimum variance portfolios, as the annualized standard deviation of returns declines by more than 25% relative to the second-best alternative in the literature and to nearly 40% versus the random walk alternative. We also observe a decline in realized portfolio risk when considering longer investment horizons of up to one month and when we introduce portfolio constraints.
Our model applies a methodology similar to callot2017modeling to the elements of the factor decomposition of realized covariance matrices. Hence, the daily covariance matrix of returns is first written into a covariance matrix for a low-dimensional set of factors plus an idiosyncratic covariance matrix that is (almost) block diagonal. We consider economically motivated factors widely used in the finance literature based on firm-level characteristics such as size, value, and profitability. The dynamics of the variances and covariances of the factors are modeled by a VHAR model estimated with LASSO. To guarantee positive definiteness of the forecasts, we apply the log matrix transformation of chiu1996matrix. The daily factor loadings (“betas”) are computed with high-frequency data. We show that these loadings are time-varying and exhibit a high degree of long-range dependence, and we model their dynamics by a HAR specification. Finally, the dynamics of each block of the idiosyncratic covariance matrix is modeled by a restricted autoregressive model, also estimated by LASSO.
In this paper, we do not address the interesting challenges involved in the construction of realized measures (e.g., how to address microstructure noise), especially in the multivariate case (e.g., asynchronicity in transactions of different assets, which bias covariance estimations toward zero). We focus solely on modeling and forecasting a large realized covariance matrix of returns on hundreds of financial assets estimated elsewhere.
The remainder of the paper is organized as follows. Sections (ref) and (ref) describe the proposed model and forecasting framework. In Section (ref), we present the data and show a descriptive analysis. Forecasting results and portfolio analysis are presented, respectively, in Sections (ref) and (ref). Finally, Section (ref) concludes the paper. Supplementary results are shown in the Appendix.
We base our methodology on the fact that realized covariance matrices are highly persistent over time, which suggests the use of an autoregressive model of a large order $p$, usually larger than 20. Defining $\boldsymbol y_{t}=\textnormal{vech}(\boldsymbol\Sigma_{t})$, where vech is the half-vectorization operator returning a vector with the unique entries of $\boldsymbol\Sigma_{t}$, one possible specification is
where $\boldsymbol\epsilon_{t}$ is zero-mean random noise.
While this model is sensible for a small number of assets, the number of parameters grows quadratically as new assets are added (curse of dimensionality). To circumvent this limitation, callot2017modeling propose using penalized regressions (LASSO) to address the large number of parameters. However, the direct use of penalized regressions is unfeasible for hundreds of assets.
To illustrate, assume that the representation in ((ref)) is used to model the covariance matrix for the constituents of the S&P 500 index\footnote{In our application, we restrict the analysis to stocks that remain in the index for the entire period of our sample. This reduces the number of stocks to 430. Regarding the number of equations and potential predicting variables, our problem is still subject to the curse of dimensionality on the same order of magnitude as the illustration.} with 10 lags, that is, $p=10$. Since each matrix $\boldsymbol\Sigma_{t}$ has $N(N+1)/2$ unique entries, this configuration would result in $125,250$ equations with $10\times 125,250$ variables each, plus constants. In this case, estimation is unfeasible even with LASSO.
To reduce the dimensionality to a manageable one, we propose the use of a factor model as well as economic restrictions based on sector classifications and penalized regressions.
Following the factor model discussed in chamberlain1983arbitrage, the excess return on any asset $i$, denoted $r_{i,t}$, satisfies
where $f_{1,t},\cdots, f_{K,t}$ are the excess returns of $K$ factors, $\beta_{ik,t}$, $k=1,\ldots, K$, are factor loadings for asset $i$, and $\varepsilon_{i,t}$ is the idiosyncratic error term. Note that the factor loadings are allowed to be variable over time. For $N$ assets, the set of equations can be written in matrix form:
where $\boldsymbol B_{t}$ is a $K \times N$ matrix of loadings, $\boldsymbol r_{t}$ is an $N \times 1$ vector of excess returns, and $\boldsymbol\varepsilon_{t}$ is an $N \times 1$ vector of idiosyncratic errors. Throughout, we assume that $\mathbb{E}(\boldsymbol\varepsilon_{t}|\boldsymbol f_{t})=\boldsymbol 0$. The factors used in this work are linear combinations of returns constructed solely with the assets considered here, i.e., long-short stock portfolios where stocks that are part of our sample are sorted on firm characteristics. In matrix form, for all $K$ factors,
where weights are calculated based on accounting and market information, as we describe in Section (ref).
This section describes how we decompose the realized covariance matrix of returns for all assets into a factor covariance matrix and a residual covariance matrix. Let $\boldsymbol\Sigma_{t}$ denote the realized covariance matrix of returns at time $t$, that is, $\boldsymbol \Sigma_{t} =\textnormal{cov}(\boldsymbol r_{t})$. By using equation ((ref)) and the assumption $\mathbb{E}(\boldsymbol \varepsilon_{t}|\boldsymbol f_{t})=\boldsymbol 0$, we have
where $\boldsymbol B_{t}$ is a $K \times N$ matrix of loadings of $N$ assets on $K$ factors, $\boldsymbol\Sigma_{f,t}$ is the $K \times K$ factor covariance matrix and $\boldsymbol\Sigma_{\boldsymbol\varepsilon, t}$ is the $N \times N$ residual covariance matrix, all at time $t$. Since each factor is a linear combination of returns, factor covariance matrices can be obtained by using equation ((ref)) and the known values of $\boldsymbol\Sigma_{t}$, that is,
Factor loadings $\boldsymbol B_{t}$ are calculated using a similar procedure (see Appendix (ref)), and the values of $\boldsymbol\Sigma_{\boldsymbol\varepsilon,t}$ are simply given by the difference $\boldsymbol\Sigma_{t}-\boldsymbol B_{t}'\boldsymbol \Sigma_{f,t}\boldsymbol B_{t}$. It is common to assume that $\boldsymbol\Sigma_{\boldsymbol\varepsilon,t}$ is diagonal, but we will be less restrictive in this work.
With the decomposition achieved in ((ref)), we write the forecasting equation for the complete covariance matrix $\boldsymbol\Sigma_{t}$ by combining separate forecasts of $\boldsymbol B_{t}$, $\boldsymbol\Sigma_{f,t}$, and $\boldsymbol\Sigma_{\boldsymbol\varepsilon, t}$. That is,
Since the number of factors is much smaller than the number of assets, one could propose using an unrestricted VAR model for the factor covariance matrix dynamics, as in equation ((ref)). Despite the reduction in dimensionality achieved by using a factor model, the number of parameters in this configuration is still quite large (note that each equation would have $K(K+1)/2+1$ parameters). To reduce this number to a more manageable one, we follow the heterogeneous autoregressive (HAR) model proposed by corsi2009simple. In this model, the predictors are obtained from the simple average of daily data computed for different horizons (daily, weekly, and monthly). In our case, daily, weekly, and realized covariance matrices of factors are given by, respectively,
where $\boldsymbol y_{f,t}=\textnormal{vech}(\boldsymbol\Sigma_{f,t})$, $\boldsymbol y_{f,t}^{day}=\textnormal{vech}(\boldsymbol\Sigma_{f,t}^{day})$, $\boldsymbol y_{f,t}^{week}=\textnormal{vech}(\boldsymbol\Sigma_{f,t}^{week})$, and $\boldsymbol y_{f,t}^{month}=\textnormal{vech}(\boldsymbol\Sigma_{f,t}^{month})$. Furthermore, the dynamic process for $y_{f,t}$ is defined as
where $\boldsymbol\Phi_{day}$, $\boldsymbol\Phi_{week}$, and $\boldsymbol\Phi_{month}$ are $M \times M$ matrices, where $M=K(K+1)/2$ is the number of unique entries on the factor covariance matrix. $\boldsymbol\omega$ is an $M\times 1$ vector of constants.
Due to the high number of parameters in equation ((ref)), direct estimation with ordinary least squares (OLS) could result in overfitting, harming the precision of model forecasts. LASSO shrinks these estimates by imposing a penalty related to the magnitude of the coefficients, shrinking them towards zero (tibshirani1996regression). This methodology has been shown to almost always provide a higher out-of-sample forecasting precision, while the reduced number of predictors makes interpreting the model easier.
We estimate ((ref)) equation by equation. Consider a sample size of $T$, and let $\boldsymbol{Z}_{t} = \left(1, \boldsymbol{y}_{f,t-1}^{day'},\boldsymbol{y}_{f,t-1}^{week'},\boldsymbol{y}_{f,t-1}^{month'}\right)'\in\mathbb{R}^{3m + 1}$ be the vector of explanatory variables and $\boldsymbol Z=(\boldsymbol Z_{T},\ldots,\boldsymbol Z_1)'$ be the $T \times (3m+1)$ matrix of covariates. Let $\boldsymbol y_{f,i} = (y_{f,T,i},\ldots, y_{f,1,i})'\in\mathbb{R}^T$ be the vector of observations on the $i$th equation of ((ref)), $i=1,\ldots,M$, and $\boldsymbol\epsilon_{f,i}=(\epsilon_{f,T,i},\ldots,\epsilon_{f,1,i})'$ be the corresponding vector of error terms. With $\boldsymbol \gamma_{f,i}=(\omega_{i},\boldsymbol\beta_{i}')'$ being the $3M+1$ vector of true parameters for equation $i$, one can write
where each vector $\boldsymbol\gamma_{f,i}$ is then estimated by minimizing
The penalty parameter $\lambda_{T}$ determines how much penalization is imposed on the size of the coefficients. In our setup, the value of $\lambda_{T}$ is determined by minimizing the Bayesian information criterion (BIC). For equation $i$ and penalty parameter $\lambda$, the BIC is given by
where $\widehat{\boldsymbol\epsilon}_{\lambda,i}$ is the estimated vector of error terms for penalty $\lambda$.
After obtaining $\widehat{\boldsymbol\gamma}_{f,i}$ from equation ((ref)), the one-step-ahead forecast is given by
which can be used to provide $\widehat{\boldsymbol \Sigma}_{f,T+1}$.
We also consider the adaptive version of LASSO (adaLASSO), as in zou2006adaptive. adaLASSO is a two-step procedure that estimates the model via LASSO in the first step and excludes the variables classified as zero from the second step. The second step estimates the parameters with a slightly modified objective function that considers the size of the parameters estimated in the first step. To illustrate, consider the set of indices of the coefficients in the $i$th equation that are different from zero in the first step: $J(\widehat{\boldsymbol\beta}_{i})=\{j \in R^{3M}: \widehat{\beta}_{i,j} \neq 0\}$. adaLASSO estimates the vector of parameters $\boldsymbol\gamma_{f,i}$ by minimizing the following objective function:
As before, the penalty parameter $\lambda_T$ is chosen by minimizing the BIC, and the one-step-ahead forecasts are computed as in equation ((ref)).
Instead of using unconditional loadings on factors, we assume that betas change daily and have long-range dependence. In Appendix (ref), we show the distribution of the estimates of the fractional integration parameter for the beta series. The results show that these series present high persistence, similar to what is observed in realized covariance data. Based on this evidence and to maintain consistency with the rest of our methodology, we forecast each series of betas using HAR models. For each element of $\boldsymbol B_{t}$, we have
where $\beta_{k,i,t}$, the entry $k\times i$ in the matrix $\boldsymbol B_{t}$, is the loading of stock i on factor k at date t. We estimate the coefficients $\Phi_{day}$, $\Phi_{day}$, and $\Phi_{day}$ by OLS.
Since the residual covariance matrix dimension is $N \times N$, the curse of dimensionality remains a concern when forecasting $\boldsymbol\Sigma_{\boldsymbol\varepsilon,t}$. A factor model may not explain all the covariances between assets depending on the number and choice of factors, implying the covariance matrix of the residuals is not diagonal. Instead of imposing the common restriction that $\boldsymbol\Sigma_{\boldsymbol\varepsilon_,t}$ is diagonal (as in fan2008high), we consider a less restrictive assumption that stocks in the same sector may still co-vary even after controlling for the standard factors used in the finance literature.
We assume that $\boldsymbol\Sigma_{\boldsymbol\varepsilon,t}$ is block diagonal, where blocks are defined by industry classification. After controlling for factor exposure, we assume there is no contemporaneous correlation between assets in different sectors. We argue that this assumption is not too strong in Section (ref), where we analyze residual correlations between stocks in different sectors.
We assume that each asset belongs to an industrial sector $s$, where $S$ is the total number of sectors and $S << N$. Without loss of generality, we can order assets so that the residual covariance matrix has sector covariance matrices near its diagonal:
Depending on the selected industry classification, the number of assets in each group $s$ (called $N^{s}$) can still be quite large, as each block has $M^{s}=N^{s}(N^{s}+1)/2 $ unique elements.
Our second assumption is that the dynamics of each block $\boldsymbol\Sigma_{\boldsymbol\varepsilon,t}^{s}$ depend only on the elements of the same block at $t-1$, i.e., only $\boldsymbol\Sigma_{\boldsymbol\varepsilon,t-1}^{s}$ predict values in $\boldsymbol\Sigma_{\boldsymbol\varepsilon,t}^{s}$. To further simplify this model, we also consider an additional restriction that only past variances matter in this prediction. This last assumption relies on the previous evidence in callot2017modeling that LASSO more frequently selects past variances as good predictors of covariance and variance terms than past covariances.
With the notation $\boldsymbol y_{\epsilon,t}^{s}=\textnormal{vech}(\boldsymbol\Sigma_{\boldsymbol\varepsilon,t}^{s})$, we have
where $\boldsymbol\omega_{\boldsymbol\varepsilon}^{s}$ is an $M^{s} \times 1$ vector of intercepts, $\boldsymbol\Phi^{s}$ is an $M^{s} \times N^{s}$ matrix of coefficients, $\boldsymbol\Lambda_{\boldsymbol\varepsilon,t-1}^{s} = \mathsf{diag}(\boldsymbol\Sigma_{\boldsymbol\varepsilon,t-1}^{s})$ is an $N^{s} \times 1$ vector of past variances, and $\boldsymbol u_{\boldsymbol\varepsilon,t}^{s}$ is the vector of errors.
The parameter estimation is done block by block. Each equation in ((ref)) is estimated separately. The procedure is the same as used for the factor covariance matrix model: LASSO/adaLASSO regression equation by equation. We then regroup the one-day ahead forecast for each group, $\widehat{\boldsymbol\Sigma}_{\boldsymbol\varepsilon, T+1}^{s}$, to form the full residual covariance matrix forecast $\widehat{\boldsymbol\Sigma}_{\boldsymbol\varepsilon, T+1}$.
The data consist of daily realized covariance matrices of returns for constituents of the S&P 500 index. These matrices were constructed using 5-minute returns by the composite realized kernel method (discussed in lunde2016econometric and provided by the authors). The full sample comprises all business days between January 2006 and December 2011. We consider companies that remained in the index and had balance sheet data available for the full sample period, resulting in a total of 430 stocks. With these considerations, the dataset consists of 1,495 daily $430 \times 430$ realized covariance matrices of returns.
Our sample spans the 2008 financial crisis and flash crashes from 2010 and 2011. To mitigate the effect of these extreme events, we perform a light cleaning on the realized covariance matrices, as in callot2017modeling. For each day, we verify which (unique) entries are more than four standard deviations (of the series corresponding to that entry) away from their sample average up to then. If more than 25% of the unique entries have extreme values according to this criterion, we flag that specific day. We replace the matrices corresponding to the flagged days by an average of the nearest ten preceding non-flagged matrices.
We construct each factor as the return time series for a long-short stock portfolio derived from individual sorts of underlying stocks on different signals. Since our approach uses the loading matrix $\boldsymbol W$ as an input to calculate factor covariance matrices and factor loadings, we could not use the widely available data on financial factors series of returns (such as on Kenneth French's website). Instead, we construct our versions of the standard factors, ranking our universe of stocks into the different signals to calculate the matrices $\boldsymbol W$ for our sample.
In total, we use seven factors widely used in the finance literature. Besides the market factor, we also consider size (SMB) and value (HML) (fama1993common), gross profitability (novy2013other), investment (lyandres2008new), asset growth (cooper2008asset), and accruals (sloan1996). We report a detailed description of the construction of each factor in Appendix (ref). We use four combinations of factors with one, three, five, and seven factors. They are denoted by 1F (market), 3F (market, size, and value), 5F (3F + gross profitability and investment), and 7F (5F + asset growth and accruals), respectively.
Each stock is classified into one of 10 sectors, following the Standard Industrial Classification (SIC). Table (ref) shows the number of companies from our sample in each sector. Note that some groups are quite large (the group Others, for instance, has more than 100 stocks). This motivates the use of additional restrictions we discussed in Section (ref).
We analyze the residual covariance matrices for the four different combinations of factors. We use the previous notation, $\boldsymbol\Sigma_{\epsilon,t}= \boldsymbol\Sigma_{t}-\boldsymbol B_{t}'\boldsymbol\Sigma_{f,t}\boldsymbol B_{t}$, for each of the 1,495 days in the sample. To analyze whether these matrices are approximately block diagonal, we follow the procedure discussed in ait2016using. First, we transform the covariance matrices into correlation matrices by dividing the entries by the standard deviations. We classify a correlation as significant if it is higher than 0.15 in absolute value in at least 1/3 of the sample. We then plot the significant relations as dots while the remaining points are blank. Sectors in Table (ref) are represented as red squares (in the same order). Figure (ref) shows the results when we apply the criterion to the full covariance matrices $\boldsymbol\Sigma_{t}$ before the factor decomposition. Most correlations are classified as significant.
When we follow this procedure for the residual covariance matrices (Figure (ref)), we obtain plots that are much more sparse than previously. This means that most correlations are not significant according to the criterion. Furthermore, we can see that most of the dots are contained inside the blocks, meaning that the most significant correlations are between stocks in the same sector. Therefore, we interpret the three plots as evidence that we would not lose much information by assuming block diagonality for the residual covariance matrices. Note that the results are robust across different factor configurations.
This section reports the forecasting results obtained using the methodology described in Section (ref). Subsection (ref) shows the forecasting results for $\boldsymbol\Sigma_{f}$, while Subsection (ref) shows the results for the complete covariance matrix $\boldsymbol\Sigma$. Daily forecasts are computed using rolling windows of 1,000 observations in both cases. Since we need 22 days to compute the first monthly regressor, $\boldsymbol\Sigma_{f,t}^{month}$, the out-of-sample forecasts comprise days 1,023 to 1,495 ($T_1$ to $T_2$), totaling 473 daily forecasts.
We evaluate our forecasts by using the $\ell_2$-norm for the vector of errors, that is, $||\widehat{\boldsymbol e}_{T+1}||=||\textnormal{vech}(\widehat{\boldsymbol\Sigma}_{T+1}-\boldsymbol\Sigma_{T+1})||$. We compare different methods by using the average $\ell_2$-forecast error:
In all cases, our reference forecast is a random walk model: $\widehat{\boldsymbol\Sigma}_{f, T+1}=\boldsymbol\Sigma_{f,T}$ and $\widehat{\boldsymbol\Sigma}_{T+1}=\boldsymbol\Sigma_{T}$ for the complete covariance matrix.
Table (ref) shows the forecast results for the factor covariance matrices, following the method described in Section (ref) (we refer to this method as FHAR, from factor HAR, hereafter). However, the results are uniformly stronger when we apply a log-matrix transformation, as proposed by chiu1996matrix. Before estimation, we apply the log-matrix transformation to all data, that is, $\boldsymbol\Omega_{f,t}=\log(\boldsymbol\Sigma_{f,t})$. We then use the FHAR model to compute $\widehat{\boldsymbol\Omega}_{f,t+1}$. Finally, we revert the transformation by applying the exponential-matrix transformation and obtain our forecasts, that is, $\widehat{\boldsymbol\Sigma}_{f,t+1}=\exp(\widehat{\boldsymbol\Omega}_{f,t+1})$. We also report the results without applying the log-matrix transformation.
From Table (ref), we see that using the log-matrix transformation considerably improves the forecasts. After this transformation, the variance and covariance series become smoother than in the original data. This reduces the weights of outliers when fitting the model and improves forecasting precision. Another advantage is that the exponential matrix is positive definite by construction. This consideration is important, as we will need the inverse of the covariance matrix $\boldsymbol\Sigma$ when solving the minimum variance problems in Section (ref). Finally, applying adaLASSO to the best-performing models does not improve the results up to the second decimal place.
Since we estimate all equations daily, we now investigate their evolution over time\footnote{The results for average equation size and change in the number of selected variables obtained with the adaLASSO are similar to the ones obtained via LASSO. These results are readily upon request.}. Figure (ref) shows the average number of variables selected by LASSO over the 473 days for the FHAR Log-matrix models. To illustrate, consider the upper-left panel in this figure. It shows the average model size for the diagonal equations in the FHAR Log-matrix model with three factors. In other words, it displays the mean size for the equations of $\sigma_{market}^2$, $\sigma_{SMB}^2$, and $\sigma_{HML}^2$. Panels on the right show the same results for the covariance equations. From top to bottom, we vary the number of factors\footnote{For the 1F configuration, we forecast only $\sigma_{market}^2$. In this case, the equations estimated by LASSO and adaLASSO have only three predictors, which are, in this case, always selected. For this reason, we do not present the average number of selected variables (constant) or average change in the number of selected variables (zero).}.
In all cases, LASSO reduces the number of selected variables substantially. This reduction is most noticeable for the FHAR Log-matrix model with seven factors, in which the number of selected variables fluctuates around 10 (out of 76 potential predictors) for the variance equations and even less for the covariance equations.
Another feature is that average equation sizes are stable over time. Figure (ref) shows the average change in the number of selected variables from one day to the next. We consider parameters changing from zero to non-zero values or the opposite direction. All results are reported in percentage (relative to the maximum number of variables in each equation). For the three-factor models, variance and covariance equations remain unchanged on most days. Despite the stronger variation in the five- and seven-factor configurations, the percentage changes are smaller than 5% on almost all days. In general, variance equations are slightly more stable. We present results for the residual covariance matrix forecasts in Appendix (ref).
Table (ref) shows the forecast results for the complete covariance matrix (430 stocks). In this case, VHAR denotes the complete methodology described in Section (ref). We compute all models using Log-matrix transformation and estimate these models using LASSO and adaLASSO. Table (ref) also reports alternative models for comparison purposes:
Some alternative models use daily returns data instead of realized covariance matrices as inputs. We compare our results to these models to see whether high-frequency data bring gains in portfolio performance.
Models based on daily returns do not outperform the benchmark, with a particularly bad performance for the EWMA. Despite having performed much better, BEKK-NL, DCC-NL, AFM1-DCC-NL, and IDR-DCC-NL models have almost two times larger errors than the random walk model. Imposing block diagonality on the residual covariance matrices (Blocks 1F, 3F, 5F, and 7F) improves the results, relative to the random walk, by 3%. Our methodology can beat the random walk model by up to 15%. In terms of the average $\ell_2$-forecast error, there is no difference when estimating the models with LASSO or adaLASSO.
Despite the encouraging results obtained when forecasting covariance matrices, one could argue that our measure of forecast precision is too aggregate, especially considering that we have more than 90,000 unique entries being forecasted simultaneously. For instance, these estimates may produce large mistakes in only a few entries but with large implications when applied to practical problems. To evaluate the economic benefits provided by our forecasts, we use them to construct daily investment portfolios\footnote{In this section, we use daily covariance matrix forecasts to build portfolios for the next day. These portfolios are rebalanced daily, using the most recent forecasts. We investigate portfolios formed for longer holding periods in Appendix (ref). Similar results to this section hold even with infrequent rebalancing.}. We construct minimum variance portfolios using the models described in the previous section. We consider the global minimum variance portfolio (without any restrictions on weights), the restricted portfolio (limiting long positions and short-selling), and the long-only portfolios (without any short-selling).
We focus our analysis on variations of the minimum variance problem, as the sole portfolio problem requires only information on the covariance matrix. Other portfolio optimization problems, such as finding tangency portfolios that maximize Sharpe ratios, also require estimates of expected returns for every eligible asset. When we analyze the outcomes, we also focus on risk measures, such as the portfolio standard deviation, as we aim to minimize the portfolio's risk. In finance theory, all efficient portfolios have higher expected returns and Sharpe ratios than the minimum variance portfolio. Minimum variance optimizations that are less precise or biased may, by chance, find portfolios with higher realized returns and higher Sharpe ratios. While we should base our assessment of the models on the risk measures, other measures related to the concentration, trading volume, and even realized return might also provide relevant information.
First, we compute various statistics for indices directly related to our portfolio exercise. Table (ref) shows the standard deviation, kurtosis, skewness, average excess return over the risk-free rate, and Sharpe ratio for the daily series of returns for three S&P 500 indices: minimum volatility, low volatility, and equal weight. We provide these values as reference points to evaluate our results, as these indices are commonly used by practitioners. According to the index provider, the first index “is designed to reflect a managed-volatility equity strategy that seeks to achieve lower total risk, measured by standard deviation, than the S&P 500 while maintaining similar characteristics.” It is based on the Northfield Open Optimizer and is rebalanced every six months. The second index is based on risk measures for each asset individually, focusing on stocks with low volatility.
Consider the problem of an investor at time $t=t_0,\ldots,T-1$ who wishes to construct a minimum variance portfolio to be held in time $t+1$. For this minimization problem, the investor needs to forecast the future covariance matrix, $\hat{\Sigma}_{t+1}$. The optimization problem consists of choosing a vector of weights $\hat{w}_{t+1}$ (dimension $N \times 1$):
We use our VHAR method with log matrix transformation, as it provides the best results in terms of forecasting precision. Estimation is done via LASSO and adaLASSO. We compare our results against the random walk models (RW and Block models), the EWMA, BEKK-NL, DCC-NL, AFM1-DCC-NL, and IDR-DCC-NL models. We evaluate ex-post portfolio performance, using our time $t$ estimated weights, $\widehat{\boldsymbol w}_{t+1}$, with data from $t+1$. In the following, $\widehat{w}_{it}$ is the $i$-th component of $\widehat{\boldsymbol w}_t$. We also define $r_{pt}^*=(r_{pt}-\frac{1}{(T-t_0)}\sum_{t=t_0+1}^{T}r_{pt})$. We show the following statistics:
Table (ref) presents the minimum variance portfolio results. Regarding standard deviations, the VHAR and Block RW models perform better than all others. For these models, the standard deviations go from a high of 8.46 % (one factor and LASSO) to a low of 8.09 % (seven factors and adaLASSO). The lowest value is almost half of the standard deviation for the S&P 500 Minimum Volatility Index in Table (ref). Benchmark models (upper results) also have a smaller standard deviation than the S&P 500 indices (except for the EWMA, which has the worst performance in this sense). Similar results hold for lower partial standard deviations. Among the models that do not use intraday data, the AFM1-DCC-NL model has superior performance.
Some portfolios have extreme short positions. In most cases, the leverage proportion is close to 50%. Average gross leverage is also quite high for many of the alternative models. Since shorting stocks may not be feasible and can be costly, we consider these positions potentially extreme. The values for average maximum and minimum weights are also quite large. Our models appear to generate more balanced portfolios.
Our models are slightly inferior to the IDR-DCC-NL model in terms of average excess and cumulative returns. The VHAR with five factors and LASSO estimation is the one with close performance to the IDR-DCC-NL alternative. Our models have superior performance compared to the RW and Block models and are much better than the EWMA specification. Since the optimization problem focuses on minimizing the portfolio variance, it is not obvious which models should achieve higher returns. For this reason, we do not see the slightly superior performance of the IRD-DCC-NL model as a problem. Moreover, the small difference in realized returns is more than compensated by the substantial reduction in standard deviation (as evidenced by the higher Sharpe ratios for our models).
We now resort to what we consider a more realistic investor problem. In this section, we solve a problem similar to equation ((ref)), except that now we impose two additional restrictions. First, we allow the maximum leverage to be 30% (in some sense, consistent with a 130-30 fund concept in the mutual fund industry). Second, we restrict the maximum weights on individual stocks to be 20% (in absolute value). The problem for an investor at time $t=t_0,\ldots, T-1$ is then given by
In Table (ref), we report the same performance statistics as in the last subsection. Again, our methodology generates lower standard deviations than the competing models. The VHAR with five factors and estimated with LASSO has the smallest standard deviation (12.57). In terms of returns, the Block 5F model is the best-performing alternative. This is the best model also in terms of Sharpe Ratios. The several VHAR specifications come second. One additional interesting result is that the IDR-DCC-NL model now performs better than the DCC-NL and AFM1-DCC-NL models, both in terms of standard deviation as cumulative and average returns.
In this section, we restrict no short-selling. The problem for an investor at time $t=t_0,\ldots, T-1$ is then given by
From Table (ref), we see that the VHAR family outperforms the benchmarks in almost all cases in terms of standard deviation. The VHAR model with seven factors estimated with LASSO is the best-performing specification.
We propose a model to forecast large realized covariance matrices of returns, applying it to some of the constituents of the S&P 500 on a daily basis. To address the curse of dimensionality, we decompose the return covariance matrix using standard firm-level factors (e.g., size, value, and profitability) and use sectoral restrictions in the residual covariance matrix. This restricted model is then estimated using vector heterogeneous autoregressive (VHAR) models with the least absolute shrinkage and selection operator (LASSO). Our methodology improves forecasting precision relative to standard benchmarks and leads to better estimates of minimum variance portfolios.