Five risk based allocation methods, backtested in Python

Two academic papers sit underneath this. The first is Marcos Lopez de Prado's Hierarchical Risk Parity (2016), which I covered in an earlier video using a Monte Carlo experiment rather than market history. The second is Thierry Raffinot's Hierarchical Clustering based Asset Allocation (2018), which adds a family of tree based methods that lean on correlation instead of variance. This time everything runs on real data, one long sample of eleven asset class ETFs, so the question becomes simple. Do the clean theoretical claims survive contact with a single messy history.

These are risk based, not return seeking

One framing point before the code. None of the five methods forecasts returns. There is no expected return vector and no efficient frontier point being chased. Even the Minimum-Variance portfolio here is the pure risk minimiser, the far left tip of the frontier, which needs no return estimate at all. So the fair yardstick is risk adjusted, Sharpe and drawdown, not raw return. Equal-Weight is included only as the benchmark that is famously hard to beat.

The two interesting rows are Naive Risk Parity and HCAA. Naive Risk Parity uses variance but is blind to correlation, so it happily piles into whichever bucket looks calmest on its own. HCAA is the mirror image. It reads the correlation structure to build the tree, then ignores variance entirely and splits capital in half at every fork. HRP sits between them, using the tree from correlation and then tilting each split by inverse variance.

The data

Eleven ETFs spanning global equity, a full ladder of US Treasuries, international sovereigns, emerging markets and gold. This is Raffinot's universe, chosen for replication rather than because it is the allocation I would run for a client. In practice the universe is the single most consequential choice you make, and it is driven by the client's investment policy statement, home bias included.

import numpy as np, pandas as pd, matplotlib.pyplot as plt
import scipy.cluster.hierarchy as sch
from scipy.spatial.distance import squareform
from scipy.optimize import minimize

NAMES = {"SPY":"US Large Cap","IWM":"US Small Cap","EZU":"Eurozone Equity",
         "EWU":"UK Equity","SHY":"US 1_3y Tsy","IEI":"US 3_7y Tsy","IEF":"US 7_10y Tsy",
         "TLT":"US 20y+ Tsy","BWX":"Intl Sovereign","EEM":"EM Equity","GLD":"Gold"}
TICKERS = list(NAMES)

px   = pd.read_csv('prices.csv', index_col=0, parse_dates=True)[TICKERS].rename(columns=NAMES)
rets = px.pct_change().dropna()
print(f'{px.shape[1]} assets, {px.shape[0]} trading days: '
      f'{px.index[0].date()} to {px.index[-1].date()}')

The five allocators

This is the heart of the script. Every allocator takes a covariance or correlation matrix estimated from a rolling window and returns a weight vector. The HRP block is a Python 3 port of de Prado's appendix, cluster the assets by a correlation distance, quasi diagonalise, then split risk by inverse variance down the tree. HCAA reuses the tree idea but splits capital in half at each node. Minimum-Variance is a plain long only quadratic programme, and Naive Risk Parity is just the inverse variance portfolio, correlation ignored.

# ---------- HRP: de Prado 2016 appendix ----------
def getIVP(cov):
    ivp = 1./np.diag(cov.values if hasattr(cov,'values') else cov)
    return ivp/ivp.sum()

def getClusterVar(cov, items):
    c = cov.loc[items, items]; w = getIVP(c).reshape(-1,1)
    return float((w.T @ c.values @ w)[0,0])

def getQuasiDiag(link):
    link = link.astype(int); s = pd.Series([link[-1,0], link[-1,1]]); n = link[-1,3]
    while s.max() >= n:
        s.index = range(0, s.shape[0]*2, 2)
        df0 = s[s >= n]; i = df0.index; j = df0.values - n
        s[i] = link[j,0]; df0 = pd.Series(link[j,1], index=i+1)
        s = pd.concat([s, df0]).sort_index(); s.index = range(s.shape[0])
    return s.tolist()

def getRecBipart(cov, sortIx):
    w = pd.Series(1.0, index=sortIx); clusters = [sortIx]
    while clusters:
        clusters = [c[j:k] for c in clusters for j,k in ((0,len(c)//2),(len(c)//2,len(c))) if len(c) > 1]
        for i in range(0, len(clusters), 2):
            c0, c1 = clusters[i], clusters[i+1]
            v0, v1 = getClusterVar(cov,c0), getClusterVar(cov,c1)
            a = 1 - v0/(v0+v1)
            w[c0] *= a; w[c1] *= 1-a
    return w

def correlDist(corr): return ((1-corr)/2.)**0.5

def getHRP(cov, corr):
    dist = correlDist(corr)
    link = sch.linkage(squareform(dist.values, checks=False), 'single')
    sortIx = [corr.index[i] for i in getQuasiDiag(link)]
    return getRecBipart(cov, sortIx).reindex(cov.index)

# ---------- HCAA: Raffinot 2018, equal split down the dendrogram ----------
def getHCAA(corr, method='average'):
    dist = np.sqrt(2*(1-corr))
    link = sch.linkage(squareform(dist.values, checks=False), method=method)
    tree = sch.to_tree(link); w = pd.Series(0.0, index=corr.index)
    def rec(node, wt):
        if node.is_leaf(): w.iloc[node.id] = wt
        else: rec(node.left, wt/2); rec(node.right, wt/2)
    rec(tree, 1.0); return w

# ---------- Minimum-Variance: Markowitz risk minimiser, long only ----------
def getMinVar(cov):
    n = len(cov); S = cov.values*252; x0 = np.ones(n)/n
    r = minimize(lambda w: w @ S @ w, x0, method='SLSQP', bounds=[(0,1)]*n,
                 constraints=[{'type':'eq','fun':lambda w: w.sum()-1}])
    return pd.Series(r.x, index=cov.index)

def getEW(cols): return pd.Series(1/len(cols), index=cols)

# ---------- Naive Risk Parity = inverse variance, correlation ignored ----------
def getNaiveRP(cov): return pd.Series(getIVP(cov), index=cov.index)

The backtest engine

One walk forward loop drives every method. At each rebalance date it re estimates the covariance from the trailing window, asks the allocator for weights, charges a transaction cost proportional to turnover, and holds until the next rebalance. Nothing from the future ever touches a weight. That is the whole point, the weights on any given day are built only from data available on that day.

LOOKBACK = 756   # trading days: 252 = 1y, 756 = 3y, 1260 = 5y
COST_BPS = 10    # one way cost per unit of turnover

def backtest(rets, alloc_fn, lookback=LOOKBACK, rebal='ME'):
    dates = rets.resample(rebal).last().index
    w_prev = pd.Series(0.0, index=rets.columns); port = []; turn = {}; sspw = []
    cur_w = None
    for t in rets.index:
        if cur_w is not None:
            port.append((t, float((cur_w * rets.loc[t]).sum())))
        if t in dates:
            win = rets.loc[:t].tail(lookback)
            if len(win) >= lookback:
                cur_w = alloc_fn(win).reindex(rets.columns).fillna(0)
                turn[t] = float((cur_w - w_prev).abs().sum())
                sspw.append(float((cur_w**2).sum())); w_prev = cur_w
    return pd.Series(dict(port)).sort_index(), pd.Series(turn), float(np.mean(sspw))

allocs = {
  'Equal-Weight': lambda w: getEW(w.columns),
  'Naive RP':     lambda w: getNaiveRP(w.cov()),
  'Min-Variance': lambda w: getMinVar(w.cov()),
  'HRP':          lambda w: getHRP(w.cov(), w.corr()),
  'HCAA':         lambda w: getHCAA(w.corr()),
}

Two summary numbers travel with every run. SSPW is the sum of squared weights, a Herfindahl style concentration score. Its reciprocal, the effective N, reads as the number of genuinely independent positions the portfolio holds. A portfolio of eleven names with an effective N near one is really a single bet wearing a diversified costume.

Results on historical volatility

Estimation window of three years, monthly rebalance, net of ten basis points of cost. The methods split cleanly into two risk regimes.

Equal-Weight and HCAA live in the high risk regime. Both hold a lot of equity, so both earn roughly six percent a year with nine to ten percent volatility and drawdowns past twenty percent. HCAA gets there for a specific reason, it never looks at variance, so it does not shrink the volatile equity sleeves. It is well spread across names, yet concentrated in risk.

The other three cluster in the low risk regime, one and a half to one and seven tenths percent return on under two percent volatility, because they load heavily on the calm short Treasuries. Inside that group the ranking is the real story:

• Minimum-Variance wins on raw Sharpe at 1.01, but its adjusted Sharpe falls to 0.79 once you penalise skew and fat tails, and its annual turnover is by far the highest of the field.
• HRP posts the best adjusted Sharpe, 0.86, and the shallowest drawdown, 7.4 percent, at a fraction of Minimum-Variance's trading. That is exactly de Prado's claim, comparable risk with far less turnover.
• Naive Risk Parity looks strong on Sharpe, 0.85, but its concentration gives it away. SSPW of 0.71 is an effective N near 1.4, so it is essentially one bond bucket. HRP reaches similar risk while staying genuinely spread.

Turnover is where HRP earns its keep. Its allocation barely moves from one three year snapshot to the next, so you spend little on rebalancing. HCAA is the opposite. Because it keys off the correlation structure, and correlations reorder violently in stress, it churns far more, more than half the book a year here.

The biases, and why the knobs exist

A single backtest on a single history is the easiest number in finance to fool yourself with. Three biases do most of the damage, and the design of this script is a direct response to each.

Overfitting, the selection bias

Notice that the winner changes if you move the estimation window from three years to five, or switch the rebalance from monthly to annual. If you quietly try every combination and report only the flattering one, you have selected a result, not discovered one. That is why the notebook exposes the window and the rebalance frequency as sliders rather than hiding a single lucky setting. The instability is the finding. De Prado sidestepped it in the first video by using Monte Carlo, thousands of synthetic histories, precisely so no one path could be cherry picked.

Look ahead, the past is not the future

Every weight above is built from trailing variance, an implicit bet that the last three years describe the next month. The honest partial fix is to stop trusting history and forecast instead, which is the GARCH section below.

Survivorship

The cleaner move is to include losers, not only the assets that happened to do well. There is a subtle catch, though. The moment you fix a universe of eleven ETFs, you have already excluded most of the investable world, so some survivorship creeps in through the back door of universe selection. Honesty here is about disclosing the choice, not pretending it is neutral.

GARCH to forecast volatility

Trailing variance treats last year and last week as equally informative and forgets that shocks cluster. A GARCH(1,1) model fixes both. It carries a persistent variance state and lets yesterday's surprise push tomorrow's forecast.

Read it left to right. ω is a long run variance floor, α weights yesterday's squared return shock, and β says how much of yesterday's variance persists into today. The (1,1) means one lag of each, the shock and the variance. That order is the workhorse for a reason. Hansen and Lunde (2005) compared hundreds of specifications and found nothing reliably beats GARCH(1,1) out of sample, so it is the honest default rather than a tuned choice.

GARCH gives you volatilities, but the allocators need a covariance matrix. The standard move is to separate the two, forecast each asset's volatility with its own GARCH, keep the sample correlation, and rebuild the covariance as

This matters most for HCAA, which consumes correlation directly, so the split between volatility and correlation is not optional bookkeeping, it is the interface the method needs.

from arch import arch_model
_gcache = {}

def garch_cov(win, corr):
    key = (win.index[-1], len(win))
    if key in _gcache: return _gcache[key]
    vols = {}
    for a in win.columns:
        r = win[a].dropna().values * 100.0        # scale up for numerical stability
        try:
            res = arch_model(r, mean='Zero', vol='Garch', p=1, q=1,
                             dist='normal', rescale=False).fit(disp='off')
            fc = res.forecast(horizon=1, reindex=False)
            vols[a] = float(np.sqrt(fc.variance.values[-1,0])) / 100.0
        except Exception:
            vols[a] = float(win[a].std())          # fall back to sample vol
    d = pd.Series(vols).reindex(win.columns); D = np.diag(d.values)
    Sig = pd.DataFrame(D @ corr.values @ D, index=win.columns, columns=win.columns)
    _gcache[key] = Sig; return Sig

g_allocs = {
  'Equal-Weight': lambda w: getEW(w.columns),
  'Naive RP':     lambda w: getNaiveRP(garch_cov(w, w.corr())),
  'Min-Variance': lambda w: getMinVar(garch_cov(w, w.corr())),
  'HRP':          lambda w: getHRP(garch_cov(w, w.corr()), w.corr()),
  'HCAA':         lambda w: getHCAA(w.corr()),
}

Swapping historical variance for a one step GARCH forecast pulls Minimum-Variance and HRP almost on top of each other, 0.83 against 0.82 on Sharpe, and trims a little off the low risk trio overall. It does not crown a new champion. That is the correct lesson. GARCH is not a magic return engine, it is a way to stop pretending the trailing window is the future. The argument for using it is the bias it attacks, not a higher number on this particular sample.

What actually held up

Three things survive the whole exercise. First, there is no single winner, the best method flips with the window and the rebalance frequency, and that instability is the honest headline rather than an inconvenience to hide. Second, HRP does deliver on its central promise, similar risk to Minimum-Variance with a fraction of the turnover and the best tail adjusted Sharpe, which is a real edge once trading costs are on the table. Third, the methods that skip a piece of the risk picture pay for it, Naive Risk Parity concentrates because it ignores correlation, and HCAA rides equity risk and churns because it ignores variance.

If you play with the sliders you will find, as I did, that some of the paper's claims hold on real data and some do not. That is the value of running it yourself rather than trusting one clean chart. The first video reproduces the de Prado result with Monte Carlo, and the papers are linked below if you want the full derivations.

Part 1, HRP with Monte Carlo: the reproduction in Python
More paper walkthroughs and notes: the papers page

Related

• MVO Portfolio Optimization, mean variance from scratch, the starting point for all of these.
• Hierarchical Risk Parity, clustering the correlation matrix before allocating.