High-dimensional portfolios (R) — corpcor · quadprog¶
Covariance shrinkage and sparse asset selection on 48 stocks¶
The R cross-check for the high-dimensional portfolio notebook. corpcor::cov.shrink implements analytic (Ledoit-Wolf-type, Schäfer-Strimmer) covariance shrinkage; quadprog::solve.QP solves the long-only minimum-variance quadratic program, whose no-short constraint acts as implicit shrinkage and selects a sparse set of assets (Jagannathan-Ma, 2003). Same 48 stocks, weekly, as the Risk & Asset Allocation section (stocks_weekly.csv). The task is the minimum-variance portfolio $\min_w w'\Sigma w$ s.t. $\mathbf 1'w=1$; the point is that with 48 assets and a one-year window, the sample $\Sigma$ is near-singular and the naive optimiser explodes.
options(repr.plot.width=13, repr.plot.height=4.6)
.libPaths(c("C:/Users/user/R/win-library/4.6", .libPaths()))
suppressMessages({library(corpcor); library(quadprog)})
d<-read.csv("stocks_weekly.csv", row.names=1); R<-as.matrix(d)/100; N<-ncol(R); tickers<-colnames(R)
# Low tol so that a badly-conditioned (but still invertible) sample S is not rejected -- that ill
# conditioning IS the phenomenon under study. It is only legitimate while T > N: a covariance from T
# observations has rank <= T-1, so at T <= 48 the matrix is exactly singular and the minimum-variance
# portfolio does not exist. Every window used below therefore exceeds N.
minvar<-function(S){ w<-solve(S, rep(1,N), tol=1e-30); w/sum(w) }
ann_vol<-function(r) sd(r)*sqrt(52)*100
ann_sh<-function(r) mean(r)/sd(r)*sqrt(52)
cat(sprintf("%d weeks x %d stocks; covariance has N(N+1)/2 = %d params vs a 52-week window.\n", nrow(R), N, N*(N+1)/2))
312 weeks x 48 stocks; covariance has N(N+1)/2 = 1176 params vs a 52-week window.
1. The blow-up — sample Markowitz in high dimensions¶
The sample minimum-variance portfolio on a single one-year window takes enormous offsetting long/short positions — double-digit gross leverage — that are pure estimation noise. The right panel traces out-of-sample volatility as the estimation window shrinks toward $N=48$: the sample portfolio's risk explodes, while cov.shrink stays low. The high-dimensional curse the factor notebook lacked.
w_s<-minvar(cov(R[1:52,]))
options(repr.plot.width=13.5, repr.plot.height=4.6); par(mfrow=c(1,2), mar=c(4,4,3,1))
o<-order(w_s); barplot(w_s[o], col=ifelse(w_s[o]<0,"#c53030","#2b6cb0"), main=sprintf("Sample min-var weights, 1-yr (gross lev %.0fx)", sum(abs(w_s))), xlab="stock (sorted)", ylab="weight")
roll_vol<-function(win, shrink){ r<-numeric(0)
for(k in win:(nrow(R)-1)){ X<-R[(k-win+1):k,]; S<-if(shrink) as.matrix(cov.shrink(X,verbose=FALSE)) else cov(X); r<-c(r, minvar(S)%*%R[k+1,]) }
ann_vol(r) }
wins<-seq(N+2,200,20); vs<-sapply(wins,roll_vol,shrink=FALSE); vl<-sapply(wins,roll_vol,shrink=TRUE)
plot(wins, vs, type="b", pch=19, col="#c53030", lwd=2, ylim=range(c(vs,vl)), xlab="estimation window (weeks)", ylab="OOS annualized vol (%)", main="The N/T curse: sample risk explodes")
lines(wins, vl, type="b", pch=19, col="#2b6cb0", lwd=2); abline(v=N, lty=3)
legend("topright", c("sample covariance","cov.shrink (Ledoit-Wolf)"), col=c("#c53030","#2b6cb0"), lwd=2, pch=19, bty="n")
par(mfrow=c(1,1))
cat(sprintf("At a 1-year window the sample min-var portfolio carries %.0fx gross leverage; its OOS vol climbs from %.0f%%\n", sum(abs(w_s)), vs[length(vs)]))
cat(sprintf("at the longest window to %.0f%% as T approaches N=48, while shrinkage holds near %.0f%% throughout. The mechanism\n", vs[1], mean(vl)))
cat(sprintf("is the conditioning: kappa(Sigma) is %.0e at a 150-week window against %.0e at 52 weeks. Windows shorter than\n", kappa(cov(R[1:150,]),exact=TRUE), kappa(cov(R[1:52,]),exact=TRUE)))
cat("N are not plotted at all -- there the sample covariance has no inverse, so the portfolio is undefined rather\n")
cat("than merely noisy, and any number drawn there would be an artefact of floating point.\n")
At a 1-year window the sample min-var portfolio carries 14x gross leverage; its OOS vol climbs from 13%
at the longest window to 61% as T approaches N=48, while shrinkage holds near 14% throughout. The mechanism
is the conditioning: kappa(Sigma) is 9e+02 at a 150-week window against 5e+04 at 52 weeks. Windows shorter than
N are not plotted at all -- there the sample covariance has no inverse, so the portfolio is undefined rather
than merely noisy, and any number drawn there would be an artefact of floating point.
2. The horse race¶
A rolling backtest (one-year window, weekly) of five portfolios: sample min-variance, Ledoit-Wolf (cov.shrink), ridge (L2) covariance, long-only sparse (L1) via quadprog, and equal-weight 1/N. Metrics a risk desk cares about: out-of-sample volatility (the objective), Sharpe, gross leverage, turnover — plus how many stocks the sparse portfolio holds. (The Sharpe ratio is mean excess return per unit of volatility, $\text{SR}=\bar r/\sigma_r$, annualized by $\sqrt{52}$ for weekly data. It is deliberately *not the objective here — these portfolios minimise variance and never estimate expected returns at all — so volatility is the column to judge them on.)*
longonly<-function(S){ D<-2*(S+1e-8*diag(N)); A<-cbind(rep(1,N),diag(N)); b<-c(1,rep(0,N))
solve.QP(D, rep(0,N), A, b, meq=1)$solution }
ridge<-function(S) minvar(S + 0.1*sum(diag(S))/N*diag(N))
win<-52; nm<-c("sample","Ledoit-Wolf","ridge (L2)","long-only (sparse L1)","equal 1/N")
ret<-setNames(lapply(nm,function(x)numeric(0)),nm); lev<-setNames(lapply(nm,function(x)numeric(0)),nm)
W<-setNames(lapply(nm,function(x)list()),nm); nz<-numeric(0)
for(k in win:(nrow(R)-1)){ X<-R[(k-win+1):k,]; S<-cov(X); nxt<-R[k+1,]
ws<-list(minvar(S), minvar(as.matrix(cov.shrink(X,verbose=FALSE))), ridge(S), longonly(S), rep(1/N,N))
for(j in seq_along(nm)){ w<-ws[[j]]; ret[[j]]<-c(ret[[j]], w%*%nxt); lev[[j]]<-c(lev[[j]], sum(abs(w))); W[[j]][[length(W[[j]])+1]]<-w }
nz<-c(nz, sum(ws[[4]]>1e-4)) }
turn<-function(wl){ M<-do.call(rbind,wl); mean(rowSums(abs(diff(M)))) }
tab<-data.frame(OOS_vol=sapply(ret,ann_vol), OOS_Sharpe=sapply(ret,ann_sh), gross_lev=sapply(lev,mean), turnover=sapply(W,turn))
print(round(tab,2))
cat(sprintf("\nlong-only sparse portfolio holds on average %.0f of %d stocks.\n", mean(nz), N))
cat(sprintf("cov.shrink and the ridge covariance land on the same volatility to two decimals (%.4f vs %.4f) without being\n",
ann_vol(ret[["Ledoit-Wolf"]]), ann_vol(ret[["ridge (L2)"]])))
cat(sprintf("the same portfolio -- their weekly returns correlate %.3f but differ week to week. Both are linear shrinkage\n",
cor(ret[["Ledoit-Wolf"]], ret[["ridge (L2)"]])))
cat("toward a diagonal target, so landing together is unsurprising; sklearn's Ledoit-Wolf, with a different target,\n")
cat("reaches 15.5% in the Python notebook. The choice of target matters more than the choice of package.\n")
gap<-ann_sh(ret[["long-only (sparse L1)"]])-ann_sh(ret[["equal 1/N"]])
set.seed(0); T<-length(ret[["equal 1/N"]])
bs<-replicate(2000,{i<-sample(T,T,replace=TRUE); ann_sh(ret[["long-only (sparse L1)"]][i])-ann_sh(ret[["equal 1/N"]][i])})
cat(sprintf("\nSparse minus 1/N on Sharpe: %+.3f, 95%% bootstrap CI [%+.3f, %+.3f] -- straddles zero, so on risk-adjusted\n",
gap, quantile(bs,0.025), quantile(bs,0.975)))
cat(sprintf("return the two are indistinguishable. On volatility, the objective, the gap is %+.2f pp and that one is real.\n",
ann_vol(ret[["long-only (sparse L1)"]])-ann_vol(ret[["equal 1/N"]])))
OOS_vol OOS_Sharpe gross_lev turnover sample 46.00 -0.07 15.07 7.86 Ledoit-Wolf 16.26 0.32 2.35 0.40 ridge (L2) 16.26 0.20 2.64 0.48 long-only (sparse L1) 15.01 0.71 1.00 0.20 equal 1/N 20.49 0.69 1.00 0.00
long-only sparse portfolio holds on average 12 of 48 stocks.
cov.shrink and the ridge covariance land on the same volatility to two decimals (16.2554 vs 16.2616) without being
the same portfolio -- their weekly returns correlate 0.979 but differ week to week. Both are linear shrinkage
toward a diagonal target, so landing together is unsurprising; sklearn's Ledoit-Wolf, with a different target,
reaches 15.5% in the Python notebook. The choice of target matters more than the choice of package.
Sparse minus 1/N on Sharpe: +0.020, 95% bootstrap CI [-0.598, +0.602] -- straddles zero, so on risk-adjusted
return the two are indistinguishable. On volatility, the objective, the gap is -5.48 pp and that one is real.
3. Reading the race, and the sparse portfolio¶
The sample portfolio has by far the highest OOS volatility and leverage — the optimiser maximised estimation error. Shrinkage and constraints collapse both to sane levels; the long-only sparse portfolio reaches the lowest volatility holding only a dozen stocks (asset selection) at the lowest turnover. Naive 1/N remains a strong Sharpe benchmark (DeMiguel et al.) — and the bootstrap in the previous cell shows that is not a figure of speech: its Sharpe is statistically indistinguishable from the optimised portfolios'. Where those portfolios do win is on volatility, the quantity they actually minimise.
options(repr.plot.width=15, repr.plot.height=4.4); par(mfrow=c(1,3), mar=c(7,4,3,1)); cols<-c("#c53030","#2b6cb0","#2f855a","#6b46c1","#a0aec0")
barplot(tab$OOS_vol, names.arg=nm, las=2, col=cols, ylab="OOS annualized vol (%)", main="OOS volatility (objective)", cex.names=0.7)
barplot(tab$gross_lev, names.arg=nm, las=2, col=cols, ylab="avg gross leverage", main="Leverage", cex.names=0.7)
barplot(tab$OOS_Sharpe, names.arg=nm, las=2, col=cols, ylab="OOS Sharpe", main="Risk-adjusted return", cex.names=0.7)
par(mfrow=c(1,1))
# sparse holdings latest year
wl<-longonly(cov(R[(nrow(R)-51):nrow(R),])); held<-which(wl>1e-4); ho<-held[order(-wl[held])]
options(repr.plot.width=9, repr.plot.height=4.2); par(mar=c(6,4,3,1))
barplot(wl[ho], names.arg=tickers[ho], las=2, col="#6b46c1", ylab="weight", main=sprintf("Sparse long-only portfolio (latest year): %d/%d stocks held", length(held), N), cex.names=0.8)
cat(sprintf("Shrinkage/constraints cut OOS vol and leverage ~an order of magnitude vs sample; the sparse portfolio holds\n"))
cat(sprintf("~%d of %d names -- asset selection by L1, the portfolio counterpart of factor selection.\n", length(held), N))
Shrinkage/constraints cut OOS vol and leverage ~an order of magnitude vs sample; the sparse portfolio holds
~15 of 48 names -- asset selection by L1, the portfolio counterpart of factor selection.
Summary¶
corpcor and quadprog reproduce the Python result: give a minimum-variance optimiser more assets than data and it maximises estimation error — the sample portfolio ran to double-digit leverage and multiples of the sensible out-of-sample volatility. Ledoit-Wolf shrinkage and the long-only sparse (L1) portfolio both tame it to ~sane volatility with 1-2× leverage, the sparse one selecting ~12 of 48 stocks (Jagannathan-Ma), while naive 1/N stays a strong Sharpe benchmark (DeMiguel et al.) — statistically indistinguishable here, on a paired bootstrap, from every optimised portfolio. The optimisers' real and verifiable win is on volatility, which is what they minimise. Shrinkage's payoff scales with dimensionality, and portfolio construction is where it is indispensable — the frequentist twin of Shrinkage Estimation of Mean and Covariance, Bayesian Estimation & Estimation Risk and The Black–Litterman Model, which treat the same $\Sigma^{-1}$ instability with a prior, and via long-only/L1 asset selection, of the Variable Selection arc. Package mirror of hdim_portfolios_python.ipynb.