Case study: principal components of the yield curve

Five yields, one covariance matrix

The data are daily changes in the yields of government bonds maturing in 1, 2, 5, 10 and 30 years. They cover 1,000 trading days, about four years.

Changes are in basis points (bp), hundredths of a percentage point (case study 1.6.4). A yield going from 4.00% to 4.05% has risen 5 bp.

The numbers are invented for illustration, generated with the seed np.random.default_rng(36). Each day is built from three hidden moves plus noise:

• a level move of typical size (standard deviation) 5 bp, shifting every yield the same way: the 2-year a little more, the 30-year less;
• a slope move of typical size 2 bp, tilting the curve: the short end one way, the long end the other;
• a curvature move of typical size 1 bp, bending the middle against the two ends;
• and independent noise of 0.5 bp at each maturity.

The analysis below does not use these recipes. It recovers the three shapes from the data alone, which is what a desk must do with real yields. Case study 3.8.4 shows that the singular value decomposition finds the same components without forming the covariance matrix.

Put the changes in a $1000 \times 5$ matrix: row $t$ holds day $t$'s five changes. On the first day every yield fell, by between 4.53 and 6.41 bp. On the second, the 1-year fell 3.21 bp while the 30-year rose 1.25 bp.

Their sample covariance matrix $S$ is np.cov, as in case study 3.2.5. Its entries are in bp² per day:

The daily volatilities, the square roots of the diagonal, are 5.30, 5.30, 5.03, 4.60 and 4.49 bp. The correlations of neighbouring maturities are high, 0.968 between the 1-year and the 2-year, and fall with distance: 0.67 between the 1-year and the 30-year. The curve moves mostly together, but not rigidly.

Principal components

$S$ is symmetric, so lesson 3.6.3 applies: eigh gives five real eigenvalues and five orthonormal eigenvectors. Each eigenvector $\mathbf{q}_k$ is a shape of curve move, with one number, its loading, for each maturity. These shapes are the principal components of the yield changes.

On any day, the size of the move along $\mathbf{q}_k$ is the coordinate $z_k = \mathbf{q}_k \cdot \Delta\mathbf{y}$, where $\Delta\mathbf{y}$ is that day's vector of changes (lesson 3.6.3). It is called the day's score on component $k$.

The score is a portfolio of yield changes with weights $\mathbf{q}_k$, so its variance is $\mathbf{q}_k^\top S\,\mathbf{q}_k$ (case study 3.2.5). Since $S\mathbf{q}_k = \lambda_k\mathbf{q}_k$ and $\mathbf{q}_k \cdot \mathbf{q}_k = 1$, that is exactly $\lambda_k$.

The total is $\operatorname{tr} S = 122.824$ bp², the sum of the five daily variances.

Reading the loadings

An eigenvector's sign is arbitrary (lesson 3.6.1), so fix one before reading it. Here each is flipped to make its 30-year loading positive.

• Level. Every loading is positive: the whole curve moves up or down together. The 2-year moves most and the 30-year least, so a level move is not quite parallel.
• Slope. Negative at the short end, positive at the long end, crossing zero near 5 years. A positive score steepens the curve; a negative one flattens it.
• Curvature. Positive at the ends and negative in the middle: the 5-year moves against the 1-year and the 30-year.

How many components?

One component explains 90.08% of the variance, two explain 98.13% and three 99.58%. After the third, the eigenvalues level off: 0.272 and 0.238 bp², close to $0.5^2 = 0.25$, the variance of the noise put in at each maturity. The last two components are that noise.

So three numbers a day describe the five yield changes. A common rule keeps components until the running total passes 99%, or until the eigenvalues stop falling and flatten out; here both say three.

Each day in three numbers

The scores $z_k$ are the coordinates of each day's change along the components. They behave as the theory says: their sample variances are the eigenvalues, and they are uncorrelated, because the eigenvectors are perpendicular (lesson 3.6.3).

np.cov subtracts each column's mean before it multiplies, as case study 3.1.4 centred return vectors. So the scores are computed from the centred changes, and rebuilding a day adds the mean back: $\Delta\mathbf{y} \approx \bar{\mathbf{y}} + z_1\mathbf{q}_1 + z_2\mathbf{q}_2 + z_3\mathbf{q}_3$, with $\bar{\mathbf{y}}$ the mean change.

Take the day with the largest move, day 808 counting from 0:

Its scores are 34.87 on level and 6.03 on slope: a large sell-off that also steepened the curve. The level score is 3.32 times the level's daily volatility of 10.52 bp. Level alone misses by up to 4.01 bp, at the 30-year; three components come within 0.37 bp at every maturity.

Over all 1,000 days, the share of the squared changes left out is 9.92% after one component, 1.87% after two and 0.42% after three. Each is one minus the running total of shares.

A position's exposure to each component

The 10-year position of the hook has a DV01 of $10,000: it loses about that much when the 10-year yield rises 1 bp (case study 2.2.4). Write its DV01 at each maturity as a vector, $\mathbf{d} = (0, 0, 0, 10000, 0)$, in dollars per bp. A day's change in value is then about $-\mathbf{d} \cdot \Delta\mathbf{y}$.

Split the day's change along the components, $\Delta\mathbf{y} - \bar{\mathbf{y}} = z_1\mathbf{q}_1 + \cdots + z_5\mathbf{q}_5$ (lesson 3.6.3). Then $\mathbf{d} \cdot (\Delta\mathbf{y} - \bar{\mathbf{y}})$ splits too, with one term $z_k\,(\mathbf{d} \cdot \mathbf{q}_k)$ per component.

For the 10-year position, each $e_k$ is 10,000 times the 10-year loading: 4212.6 on level, 3554.6 on slope and $-2171.9$ on curvature. Its daily volatility, $\sqrt{\mathbf{d}^\top S\mathbf{d}}$, is $45,983. Level carries 92.86% of the variance, slope 5.91% and curvature 0.40%.

Hedging the level move

Sell 2-year bonds with a DV01 of $h$ dollars: the hedged position is $\mathbf{d} = (0, -h, 0, 10000, 0)$.

A DV01-matched hedge sells $h = 10{,}000$, so that a parallel shift of 1 bp at every maturity costs nothing. But the level component moves the 2-year more than the 10-year, by 0.4937 against 0.4213. So this hedge sells too much: its level exposure is $-724.2$.

A level-neutral hedge sets $e_1 = 10000 \times 0.4213 - h \times 0.4937 = 0$. So $h = 10000 \times \frac{0.4213}{0.4937}$: 8533.1 with the unrounded loadings, whose ratio is 0.8533.

What is left after the level hedge is almost all slope, 91.29% of the variance. The 2-year and the 10-year load on slope with opposite signs, so long one and short the other is a bet on the slope.

Case study 3.5.5's minimum-variance hedge, $h = 10000 \times \frac{\operatorname{Cov}(\Delta y_{10}, \Delta y_2)}{\operatorname{Var}(\Delta y_2)} = 7788.8$, leaves a little less, $20,196 a day, by trading some level exposure for less slope exposure.

To remove slope too, add a second hedge: one instrument for each component, found by solving a $2 \times 2$ system (lesson 3.3.1). With the 2-year and the 30-year, the system says to sell $3,329.1 and $6,933.9 of DV01. The daily volatility falls to $9,519, and curvature is now the largest part, 61.96% of what is left.

On the desk

import numpy as np

SHAPES = np.array([[1.0, 1.05, 1.0, 0.9, 0.8], [-1.0, -0.6, 0.0, 0.5, 1.0], [0.6, -0.2, -0.8, -0.3, 0.7]])
_rng = np.random.default_rng(36)
_moves = _rng.standard_normal((1000, 3)) * np.array([5.0, 2.0, 1.0])
CHANGES = _moves @ SHAPES + 0.5 * _rng.standard_normal((1000, 5))   # the case study's data, bp per day
TEN_YEAR = np.array([0, 0, 0, 10_000.0, 0])                         # $10,000 per bp at 10 years


def _raises(fn, what):
    try:
        fn()
    except ValueError:
        return
    raise AssertionError(f"expected ValueError for {what}")


def test_variances():
    """The case study's data give variances 110.638, 9.889 and 1.787 bp², largest first"""
    variances, _, _ = curve_pca(CHANGES)
    expect_close(variances, [110.638, 9.889, 1.787], what="variances", rel=2e-4)


def test_all_five_descending():
    """With k = 5, all five variances come back, largest first"""
    variances, loadings, _ = curve_pca(CHANGES, k=5)
    expect_equal(len(variances), 5, what="number of variances")
    expect_equal(bool(np.all(np.diff(variances) <= 0)), True, what="variances in decreasing order")
    expect_equal(loadings.shape, (5, 5), what="shape of loadings")


def test_shares_of_the_total():
    """Shares are of the total variance, 122.824, not of the three kept: 90.08%, 8.05%, 1.46%"""
    _, _, shares = curve_pca(CHANGES)
    expect_close(shares, [0.90078, 0.08051, 0.01455], what="shares", rel=1e-3)


def test_level_loading():
    """The first loading is the level: (0.4693, 0.4937, 0.4705, 0.4213, 0.3705)"""
    _, loadings, _ = curve_pca(CHANGES)
    expect_close(loadings[:, 0], [0.4693, 0.4937, 0.4705, 0.4213, 0.3705], what="level loading", abs_tol=1e-4)


def test_signs():
    """Each loading has a positive 30-year entry, so the slope runs from −0.5731 to 0.6644"""
    _, loadings, _ = curve_pca(CHANGES)
    expect_equal(bool(np.all(loadings[-1] > 0)), True, what="30-year entries all positive")
    expect_close(loadings[[0, -1], 1], [-0.5731, 0.6644], what="slope loading at 1 and 30 years", abs_tol=1e-4)


def test_orthonormal():
    """The loadings are orthonormal: their columns give LᵀL = I"""
    _, loadings, _ = curve_pca(CHANGES)
    expect_close(loadings.T @ loadings, np.eye(3), what="LᵀL", abs_tol=1e-10)


def test_centred():
    """Adding 3 bp to every change leaves the components unchanged: the data are centred first"""
    variances, loadings, _ = curve_pca(CHANGES + 3.0)
    expect_close(variances, [110.638, 9.889, 1.787], what="variances of the shifted data", rel=2e-4)
    expect_close(loadings[:, 0], [0.4693, 0.4937, 0.4705, 0.4213, 0.3705], what="level loading of the shifted data", abs_tol=1e-4)


def test_small_example():
    """Changes (2, 2), (−2, −2), (1, −1), (−1, 1) give variances 16/3 and 4/3, shares 0.8 and 0.2"""
    X = np.array([[2, 2], [-2, -2], [1, -1], [-1, 1]])
    variances, loadings, shares = curve_pca(X, k=2)
    expect_close(variances, [16 / 3, 4 / 3], what="variances")
    expect_close(shares, [0.8, 0.2], what="shares")
    expect_close(loadings, np.array([[1, -1], [1, 1]]) / np.sqrt(2), what="loadings (last entry positive)")


def test_bad_inputs():
    """A missing value, k = 0, k = 6, a single day or a 1-D array: ValueError"""
    gap = CHANGES.copy()
    gap[10, 2] = np.nan
    _raises(lambda: curve_pca(gap), "a NaN in the data")
    _raises(lambda: curve_pca(CHANGES, k=0), "k = 0")
    _raises(lambda: curve_pca(CHANGES, k=6), "k = 6 with 5 maturities")
    _raises(lambda: curve_pca(CHANGES[:1]), "a single day")
    _raises(lambda: curve_pca(CHANGES[:, 0]), "a 1-D array")


def test_level_hedge():
    """Hedging the 10-year's level exposure with the 2-year sells $8,533.1 of DV01"""
    _, loadings, _ = curve_pca(CHANGES)
    h = neutral_hedge(TEN_YEAR, [1], loadings[:, :1])
    expect_close(h, [8533.1], what="2-year DV01 to sell", rel=1e-5)


def test_level_and_slope_hedge():
    """Hedging level and slope with the 2-year and 30-year sells $3,329.1 and $6,933.9 of DV01"""
    _, loadings, _ = curve_pca(CHANGES)
    h = neutral_hedge(TEN_YEAR, [1, 4], loadings[:, :2])
    expect_close(h, [3329.1, 6933.9], what="2-year and 30-year DV01 to sell", rel=1e-4)
    hedged = TEN_YEAR.copy()
    hedged[[1, 4]] -= h
    expect_close(hedged @ loadings[:, :2], [0.0, 0.0], what="exposures after the hedge", abs_tol=1e-6)


def test_hedge_errors():
    """Two hedges for one component, or the same maturity twice: ValueError"""
    _, loadings, _ = curve_pca(CHANGES)
    _raises(lambda: neutral_hedge(TEN_YEAR, [1, 4], loadings[:, :1]), "two hedges for one component")
    _raises(lambda: neutral_hedge(TEN_YEAR, [1, 1], loadings[:, :2]), "the 2-year used twice")
    _raises(lambda: neutral_hedge(TEN_YEAR[:4], [1], loadings[:, :1]), "a dv01 with 4 entries for 5 maturities")