Girosi King 2007 - Lee-Carter Model Estimation and Forecasting

Summary

Lee-Carter models log mortality $m_{at}={\alpha }_a+{\beta }_a{\gamma }_t+{ϵ}_{at}$. With identification constraints ${\sum }_t{\gamma }_t=0$ and ${\sum }_a{\beta }_a^2=1$ (Girosi-King convention; Lee-Carter use ${\sum }_a{\beta }_a=1$), ${\alpha }_a$ is the time average and $(\beta ,\gamma )$ come from the first singular vectors of the centred matrix. The index ${\gamma }_t$ is forecast as a random walk with drift, so forecasts are a straight line through the first and last $\hat{\gamma }$.

Overview

For one country, $m_{at}$ is log mortality in age group $a=\mathrm{1..}A$ and year $t=\mathrm{1..}T$. The model is a rank-one product, a special case of log-linear contingency-table models with independence of age and time effects.

Main Content

Lee-Carter model (Eqs. 1-3)

$$m_{at}={\alpha }_a+{\beta }_a{\gamma }_t+{ϵ}_{at}.$$

Non-unique: invariant to $\beta \to c\beta ,\gamma \to \gamma /c$ and ${\alpha }_a\to {\alpha }_a-{\beta }_{a}c,{\gamma }_t\to {\gamma }_t+c$. Constraints ${\sum }_t{\gamma }_t=0$ and ${\sum }_a{\beta }_a^2=1$ imply ${\alpha }_a={\bar{m}}_a$ (time average). Centred model ${\tilde{m}}_{at}=m_{at}-{\bar{m}}_a\sim N({\beta }_a{\gamma }_t,{\sigma }^2)$ uses $A+T$ parameters for $A\times T$ cells. Interpretation: ${\beta }_a$ is the sensitivity of age group $a$ to the index; independence assumption = no age$\times$time interaction (the age shape $\beta$ is constant over time).

Estimation by SVD (Sec. 1.2)

Maximum likelihood under Gaussian errors is solved by the SVD $\tilde{m}=BL{U}^{\prime }$: $\hat{\beta }$ is the first left singular vector (unit length) and ${\hat{\gamma }}_t={\hat{\beta }}^{\prime }{\tilde{m}}_t$. Equivalently $\hat{\beta }$ is the leading eigenvector of $C=\tilde{m}{\tilde{m}}^{\prime }$. Likelihood has many equivalent maxima, so SVD avoids optimiser trouble. Lee and Carter also propose a second-stage re-estimation of ${\gamma }_t$ (see Girosi King 2007 - Second-Stage Reestimation and Wilmoth Estimators).

Forecast model for the index (Eqs. 4-9)

${\hat{\gamma }}_t={\hat{\gamma }}_{t-1}+\theta +{\xi }_t$, ${\xi }_t\sim N(0,{\sigma }_{rw}^2)$. MLE of the drift is $\hat{\theta }=({\hat{\gamma }}_T-{\hat{\gamma }}_1)/(T-1)$ (depends only on first and last values), $\text{Var}[\hat{\theta }]={\sigma }_{rw}^2/(T-1)$. Forecast: ${\hat{\gamma }}_{T+\Delta t}={\hat{\gamma }}_T+\Delta t\hat{\theta }+{\sum }_{l=1}^{\Delta t}{\xi }_{T+l-1}$, conditional SE grows as $\sqrt{\Delta t}$. Point forecast of log mortality: ${\mu }_{T+\Delta t}=\bar{m}+\hat{\beta }[{\hat{\gamma }}_T+\Delta t\hat{\theta }]$, linear in $\Delta t$ (Eq. 8). Example: forecast for 2030 from data 1950-2000: ${\hat{\mu }}_{2030}=\bar{m}+\hat{\beta }[{\hat{\gamma }}_{2000}+30({\hat{\gamma }}_{2000}-{\hat{\gamma }}_{1950})/50]$.

Other ARIMA models for ${\gamma }_t$ are possible, but random walk with drift “has been used almost exclusively”. Lee and Miller (2001) note the evidence is against constant ${\beta }_a$ and recommend starting from 1950 (Tuljapurkar, Li and Boe 2000); Booth, Maindonald and Smith (2002) show this fails in some countries; Carter and Prskawetz (2000) extend to structural shifts.

Examples

Rank-1 fit on synthetic data

Setup: $A=17$ age groups, $T=50$ years, truth $\beta \propto |\sin (a/4)|+0.3$, $\gamma$ falling linearly from 4 to -4 plus noise. Solution: centre by age, SVD, take the first singular pair. Result: the first component explains 97% of centred variance; forecast mean profile $\bar{m}+\hat{\beta }[{\hat{\gamma }}_T+\Delta t\hat{\theta }]$.

import numpy as np
def lee_carter(m):                       # m: A x T log-mortality
    mbar = m.mean(1); mt = m - mbar[:, None]
    U, S, Vt = np.linalg.svd(mt, full_matrices=False)
    beta = U[:, 0] * np.sign(U[:, 0].sum()); g = beta @ mt
    theta = (g[-1] - g[0]) / (m.shape[1] - 1)
    return mbar, beta, g, theta, S[0]**2 / (S**2).sum()
def forecast(mbar, beta, g, theta, h): return mbar + beta * (g[-1] + h * theta)

PyMC-style Bayesian version (untested). E111 replaces the Gaussian-on-log-rate with counts: $d_{at}\sim \text{NegBin}(E_{at}{e}^{{\alpha }_a+{\beta }_a{\gamma }_t},\varphi )$, ${\gamma }_t$ a random walk with drift, zero-sum identification on $\gamma$ and unit-norm or simplex constraint on $\beta$ (the E111 text notes an unidentified level shown failing). This follows Wilmoth’s Poisson/MLE suggestion (Girosi King 2007 - Second-Stage Reestimation and Wilmoth Estimators).

Connections

Girosi King 2007 - Lee-Carter as a Random Walk with Drift, Factor Analysis and PPCA, Concept Note - Abridged Life Tables Years of Life Lost and the Counterfactual Path.

See Also

Girosi King 2007 - Understanding the Lee-Carter Method - Overview