Modelling USA Age-Cohort Mortality: A Comparison of Multi-Factor Affine Mortality Models ^†

Abstract

Affine mortality models are well suited for theoretical and practical application in pricing and risk management of mortality risk. They produce consistent, closed-form stochastic survival curves allowing for the efficient valuation of mortality-linked claims. We model USA age-cohort mortality data using five multi-factor affine mortality models. We focus on three-factor models and compare four Gaussian models along with a model based on the Cox–Ingersoll–Ross (CIR) process, allowing for Gamma-distributed mortality rates. We compare and assess the Gaussian Arbitrage-Free Nelson–Siegel (AFNS) mortality model, which incorporates level, slope and curvature factors, and the canonical Gaussian factor model, both with and without correlations in the factor dynamics. We show that for USA mortality data, the probability of negative mortality rates in the Gaussian models is small. Models are estimated using discrete time versions of the models with age-cohort data capturing variability in cohort mortality curves. Poisson variation in mortality data is included in the model estimation using the Kalman filter through the measurement equation. We consider models incorporating factor dependence to capture the effects of age-dependence in the mortality curves. The analysis demonstrates that the Gaussian independent-factor AFNS model performs well compared to the other affine models in explaining and forecasting USA age-cohort mortality data.

1. Introduction

Longevity risk is the now well-recognised risk that the overall survival probability of a reference population is higher than expected (Cairns et al. 2006a). Improvements in mortality experienced in recent decades show significant volatility, as evidenced by the impact of COVID-19. Life insurance companies and pension funds, as the holders of substantial longevity risk, are significantly impacted by mortality trends and uncertainty (Blake et al. 2014). The Continuous Mortality Investigation (2018) highlighted the potential impact of the underestimation of longevity risk, estimating that an extra $450 billion would be required in pension payments per year (The Joint Forum 2013). Mortality uncertainty requires increase capital requirements for insurers (Barrieu et al. 2012). The quantification of longevity risk is fundamental to the successful operation of life insurance companies and pension funds.

Insurers and pension funds can manage the risk of mortality uncertainty by transferring longevity risk to counterparties using longevity swaps and to capital markets using securitization (The Joint Forum 2013). Capital markets are expected to be increasingly important in longevity risk management, arising from the potential for a more efficient and effective approach (Blake et al. 2018). Longevity-linked securities, including longevity bonds (Blake and Burrows 2001), longevity swaps (Dowd et al. 2006), and q-forwards (Coughlan et al. 2007), have been proposed. Risk quantification and fair pricing are challenging aspects of these developments. The Life and Longevity Markets Association (2010) acknowledges that mortality improvements and uncertainty are important inputs in pricing longevity risk, but modelling and forecasting improvement trends and uncertainty continue to attract significant research interest.

Stochastic mortality models capture the stochastic evolution of mortality rates at an aggregate population level. Models include the single-factor Lee–Carter model (Lee and Carter 1992) and numerous extensions and improvements. Cairns et al. (2009) provide a comprehensive summary and comparison of these extensions. The Cairns–Blake–Dowd model (Cairns et al. 2006b), a two-factor model, and the age-period-cohort model (Renshaw and Haberman 2006) are popular models. Many of these models do not have closed-form solutions for survival curves, requiring simulation to compute future expected survival probabilities for applications.

Continuous-time affine mortality models, based on mathematical finance models for interest rate and credit risk, were introduced in Milevsky and Promislow (2001), Dahl (2004), and Cairns et al. (2006a), amongst others. Continuous-time mortality models apply diffusion processes to the dynamics of mortality intensities. They are designed to be incorporated into the modelling of longevity risk using consistent model frameworks for mortality and financial risks with applications in the valuation of longevity-linked securities (Jevtic et al. 2013). A general framework for continuous time mortality models for multiple populations is provided in Jevtić and Regis (2019) with applications in UK and Dutch mortality data.

We confine our comparison of mortality models to the continuous-time arbitrage-free modelling framework. The Lee–Carter model is a popular, simple, one-factor model (which is affine in the log mortality rates). A comparison of the Lee–Carter model with three continuous time models, which include two single factor affine mortality models, is found in Novokreshchenova (2016). The single-factor affine mortality models in that study have lower mean absolute prediction error for UK and Australian mortality data compared with the Lee–Carter model. We consider three factor affine mortality models that provide better fit and prediction compared to single-factor models. SriDaran et al. (2022) consider extensions to the Lee–Carter mortality model in the form of generalized age-period-cohort mortality models using regularization to select factors from a large number of factors. The empirical results show that a combination of level, slope, and curvature factors can explain the dynamics of many of the countries in the Human Mortality Database. Although these models differ from ours, the results motivate our use of affine mortality models and the inclusion of an affine mortality model with level, slope, and curvature factors.

Affine mortality models apply concepts underlying Affine Term Structure Models (ATSMs) for interest rate modelling, as in Duffie and Kan (1996) and Dai and Singleton (2000), to mortality rates. Affine mortality models are similar to interest-rate models, provide an integrated pricing framework (Barrieu et al. 2012), and allow the derivation of closed-form solutions for survival probabilities (Dahl 2004). Affine processes provide flexibility and analytical tractability. We develop models with consistency between the dynamics and the functional form of the survival curve. We also impose consistency in the cross-sectional survival probabilities through the arbitrage-free assumption applied to these probabilities. Models that satisfy a consistency requirement (Björk and Christensen 1999) have stable parameters and ensure consistency between the dynamics of mortality rates and the functional form for the survival curve. Projected survival curves are consistent with the dynamics of mortality rates as discussed in De Rossi (2004) and demonstrated empirically in Blackburn and Sherris (2013). Survival curves are exponential affine functions of factors, with factor loadings determining how the risk factors impact differing ages through time.

Pricing of longevity-linked cash flows requires risk-adjusted survival probabilities for a cohort in a reference population (Xu et al. 2020a). Affine mortality models have been predominantly considered for capturing the mortality dynamics of a single cohort, as in Dahl and Møller (2006), Biffis (2005) and Luciano et al. (2008). Affine models have been calibrated to age-period mortality data for a reference population, as in Schrager (2006). Gaussian affine mortality models fit age-period mortality data well, although less so at older ages Blackburn and Sherris (2013). Blackburn and Sherris (2013) shows how three-factor affine mortality models in canonical form capture older-age mortality variation better than two-factor models. Gaussian models do not capture mortality heterogeneity at older ages (Pitacco 2016). Alai et al. (2019) show that the Gamma distribution fits mortality intensities well, which is consistent with mortality heterogeneity, suggesting non-Gaussian models may improve model fit at older ages. Jevtić and Regis (2021) develop square-root latent factor affine mortality models and show how these models provide a good fit to UK mortality data. Cohort effects have been observed in age-period data for many countries, as discussed, for example, in Willets (2004), Cairns et al. (2009) and Gallop (2008).

This motivates our modelling of age-cohort mortality data. We focus on a single-cohort mortality curve. Extensions to multiple-cohort affine age-cohort mortality models are found in Jevtic et al. (2013), Chang and Sherris (2018), Jevtić and Regis (2019), Xu et al. (2020a), and Jevtić and Regis (2021). We model the older ages of a single cohort using age-cohort mortality data. We develop mortality models in a risk-adjusted framework and, through a change of measure, mortality dynamics that can be calibrated to real-world or historical data. We estimate and compare five continuous-time affine cohort mortality models using age-cohort mortality data from the USA for males from ages 50 to 100 for cohorts with complete data born from 1883 to 1915. Blackburn and Sherris (2013) shows how three-factor models perform well in explaining mortality variations at the ages we consider, so we focus on three-factor affine mortality models.

The five models we consider and compare include a canonical affine mortality model; an Arbitrage-Free Nelson–Siegel (AFNS) mortality model (Christensen et al. 2011) with identifiable factors of level, slope, and curvature of the mortality curve; and a mortality model based on the Cox–Ingersoll–Ross (CIR) model (Cox et al. 1985; Jevtić and Regis 2021), allowing for a gamma distribution for mortality rates. We also investigate the impact of incorporating factor dependence to capture correlations for the Gaussian mortality models, giving another two models. We capture cohort effects directly using age-cohort data to calibrate and assess the model survival curve fit and forecasting performance. The continuous dynamics of mortality rates are discretized for estimation and implementation. Maximum likelihood with a Kalman filter are used to estimate model parameters. We estimate parameters for the factor dynamics in the real world or historical measure and the cross sectional survival curve parameters in the Q, or pricing measure, in the Kalman filter measurement equation.

Our contributions to mortality modelling are to provide empirical support for Gaussian affine mortality models based on USA data, to provide a detailed comparison of a number of multi-factor affine mortality models for the first time, and to investigate how modelling age-cohort data with affine mortality models can capture mortality dynamics. We show that, empirically, the independent-factor AFNS mortality model performs well. It better captures the variation in cohort mortality rates in USA data and produces a better fit at older ages than the independent-factor canonical Blackburn–Sherris model. Incorporating dependence in the factors for the Blackburn–Sherris mortality model improves in-sample model fit and out-of-sample forecasting performance. We also show that, for the 1916 birth cohort, the independent-factor AFNS mortality model has better predictive performance compared to the other models. Negative mortality rates, a potential limitation of Gaussian mortally models, have very low empirical probabilities for the AFNS mortality models. The CIR mortality model has the best in-sample model fit but performs poorly in predictive performance for the 1916 birth cohort.

This paper is structured as follows. Section 2 summarizes the framework for affine mortality models and specifies the structure of the continuous-time cohort mortality models. Section 3 describes the US mortality data we use for calibration. Section 4 outlines the estimation methodology using the Kalman filter and provides an analysis of the estimation results and model comparison. In Section 5, we estimate the out-of-sample expected survival probabilities for the latest cohort with full mortality data and assess the out-of-sample forecasting ability of the different models. In Section 6, we discuss the implications of the results for mortality modelling as well as further research. Section 7 concludes the paper with a summary and major findings.

2. Affine Mortality Models

Our model is initially developed in an arbitrage-free valuation framework. This is then linked to the historical dynamics through a change in measure. We first outline the continuous-time model framework of affine mortality models, then introduce the AFNS and the CIR mortality models. We discuss the incorporation of factor dependence in the affine models. We initially derive the mortality models in a financial modelling setting with a risk-neutral pricing measure. This allows the derivation of an expression for the risk-neutral survival probabilities that is consistent with the dynamics assumed for the latent factors driving mortality rates. We then give the link between the risk-neutral dynamics and the real-world dynamics, which allows us to calibrate the models to historical mortality data while preserving cross sectional consistency in survival curves. The arbitrage-free model assumption ensures consistency between the dynamics of the mortality rates and the functional form for the cross-sectional survivor curve as well as ensuring consistency in the arbitrage-free cross-sectional survival probabilities

2.1. Mortality Model Framework

The models are formally defined based on a filtered probability space $\left(\Omega ,\mathcal{F},\mathbb{F},P\right)$, where $\Omega$ is the set of possible states of nature and $\mathbb{F}={\left\{{\mathcal{F}}_t\right\}}_{0\le t\le T}$, where ${\mathcal{F}}_t={\mathcal{H}}_{t}V{\mathcal{M}}_t$ is the combined filtration for both the term structure of interest rates and mortality, assumed to satisfy the conditions of right continuity; ${\mathcal{H}}_t$ is the filtration generated by the term structure of interest rates up to time t; and ${\mathcal{M}}_t$ is the filtration containing all the information generated by the evolution of the mortality rates and survival curves for mortality up to time t.

With an incomplete market for longevity risk, there exists no unique risk-neutral measure Q for pricing mortality linked cash flows (Xu et al. 2020a). In pricing longevity risk and longevity-linked financial products, the risk-neutral measure Q is defined in terms of a zero-coupon longevity bond (Cairns et al. 2006a; Bauer et al. 2008; Blackburn and Sherris 2013). The real-world measure P reflects the best estimate of mortality, which is estimated from historical mortality data (Bauer et al. 2008).

We use $\tilde{S}\left(t,T,x\right)$ to denote the survival probability of an individual aged x at time t surviving to time T ($T\ge t$). The price of a zero-coupon bond paying $1 at time T is denoted as $P\left(t,T\right)$. The price of a longevity bond at time t that pays the amount $\tilde{S}\left(t,T,x\right)$ at time T is denoted by $\tilde{P}\left(t,T,x\right)$. Under the measure Q, equivalent to the real-world measure P for all t, T, and x, for an arbitrage-free financial market, we have

$$P\left(t,T\right)=E^Q\left[exp\left(-{\int }_t^T{r}_{u}du\right)|{\mathcal{H}}_t\right],$$

(1)

$$\tilde{P}\left(t,T,x\right)=E^Q\left[exp\left(-{\int }_t^T{r}_{u}du\right)\tilde{S}\left(t,T,x\right)|{\mathcal{F}}_t\right].$$

(2)

Assuming independence between interest rates and mortality, the price of a zero-coupon longevity bond can be written as:

$$\tilde{P}\left(t,T,x\right)=E^Q\left[exp\left(-{\int }_t^T{r}_{u}du\right)|{\mathcal{H}}_t\right]E^Q\left[\tilde{S}\left(t,T,x\right)|{\mathcal{M}}_t\right]=P\left(t,T\right)S\left(t,T,x\right),$$

(3)

where $S\left(t,T,x\right)$ is the risk-neutral survival probability used for pricing survival contingent claims.

The instantaneous mortality intensity for any age x at time t is defined as an affine function of latent factors, with

$${\mu }_x\left(t\right)={\rho }_1^{\prime }{X}_t,$$

(4)

where ${\rho }_1\in {\mathbb{R}}^n$, and $X_t\in {\mathbb{R}}^n$ is a vector of n latent factors that are assumed to drive the mortality intensity.

We consider a single cohort so that the subscript x indicates the age of the cohort at time t and individuals age as time increases in the cohort. We will later fix the initial age of a cohort at 50 for cohorts with birth years ranging from 1883 to 1915, and the age-cohort mortality rates for each cohort will be from ages 50 to 100. The historical age-cohort mortality rates are single-cohort curves differing by the birth year of the cohort.

The dynamics of the latent factors $X_t$ are given by the following system of stochastic differential equations (SDEs) under the risk-neutral measure Q (Duffie and Kan 1996; Christensen et al. 2011):

$$d{X}_t=K^Q\left[{\theta }^Q-X_t\right]dt+\mathsf{\Sigma }D\left(X_t,t\right)d{W}_t^Q,$$

(5)

where $K^Q\in {\mathbb{R}}^{n\times n}$ is the mean reversion matrix, ${\theta }^Q\in {\mathbb{R}}^n$ is the long-term mean, $\mathsf{\Sigma }\in {\mathbb{R}}^{n\times n}$ is the volatility matrix, $W_t^Q\in {\mathbb{R}}^n$ is a standard Brownian motion, and $D\left(X_t,t\right)$ is a diagonal matrix with the ith diagonal entry as $\sqrt{{\alpha }^i\left(t\right)+{\beta }_1^i\left(t\right)x_t^1+\dots +{\beta }_n^i\left(t\right)x_t^n}$. $\alpha$ and $\beta$ are bounded continuous functions. Instantaneous changes in the latent factors are modelled as stationary processes.

Under these dynamics, the risk-neutral survival probabilities for a specified age x at time t to time T are given by (Blackburn and Sherris 2013):

$$S\left(t,T,x\right)=E^Q\left[exp\left(-{\int }_x^{x+T-t}\mu \left(s\right)ds\right)\right]=exp\left(B{\left(t,T\right)}^{\prime }{X}_t+A\left(t,T\right)\right),$$

(6)

where $B\left(t,T\right)$ and $A\left(t,T\right)$ are the solutions to the following system of ordinary differential equations (ODEs):

$$\frac{dB\left(t,T\right)}{dt}={\rho }_1+{\left(K^Q\right)}^{\prime }B\left(t,T\right),$$

(7)

$$\frac{dA\left(t,T\right)}{dt}=-B{\left(t,T\right)}^{\prime }{K}^Q{\theta }^Q-\frac{1}{2}\sum _{j=1}^3{\left({\mathsf{\Sigma }}^{\prime }B\left(t,T\right)B{\left(t,T\right)}^{\prime }\mathsf{\Sigma }\right)}_{j,j},$$

(8)

with boundary conditions $B\left(T,T\right)=A\left(T,T\right)=0$.

The average force of mortality over the duration $\left(T-t\right)$ for age x is affine in the latent factors and is defined as (Blackburn and Sherris 2013; Xu et al. 2020a):

$$\overline{\mu }\left(t,T,x\right)=-\frac{1}{T-t}log\left[S\left(t,T,x\right)\right]=-\frac{B{\left(t,T\right)}^{\prime }}{T-t}{X}_t-\frac{A\left(t,T\right)}{T-t}.$$

(9)

2.2. Multi-Factor Affine Cohort Mortality Models

We present the dynamics of the latent factors in three-factor affine models. For the different models, we specify the dynamics for Equations (4) and (5) and give the solutions to the ODEs for Equations (7) and (8).

2.2.1. Independent Factor Models with Gaussian Processes

We first consider three-factor affine mortality models with factors following Gaussian processes. We consider the three-factor independent model in Blackburn and Sherris (2013) (hereafter the Blackburn–Sherris model). We also consider an affine mortality model based on the AFNS interest-rate term structure model with the survival curve driven by factors for level ($L_t$), slope ($S_t$), and curvature ($C_t$).

Table 1. Affine Mortality Models—Independent Factor Model Specifications.

Model	Factors ${\mathit{X}}_{\mathit{t}}$	${\mathit{\rho }}_1$	${\mathit{K}}^{\mathit{Q}}$	${\mathit{K}}^{\mathit{P}}$	$\mathsf{\Sigma }$
Blackburn–Sherris Model	$\left(\begin{array}{c}{X}_t^1\\ X_t^2\\ X_t^3\end{array}\right)$	$\left(\begin{array}{c}1\\ 1\\ 1\end{array}\right)$	$\left(\begin{array}{ccc}{\delta }_1& 0& 0\\ 0& {\delta }_2& 0\\ 0& 0& {\delta }_3\end{array}\right)$	$\left(\begin{array}{ccc}{k}_1^P& 0& 0\\ 0& k_2^P& 0\\ 0& 0& k_3^P\end{array}\right)$	$\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ 0& {\sigma }_{22}& 0\\ 0& 0& {\sigma }_{33}\end{array}\right)$
AFNS Model	$\left(\begin{array}{c}{L}_t\\ S_t\\ C_t\end{array}\right)$	$\left(\begin{array}{c}1\\ 1\\ 0\end{array}\right)$	$\left(\begin{array}{ccc}0& 0& 0\\ 0& \delta & -\delta \\ 0& 0& \delta \end{array}\right)$	$\left(\begin{array}{ccc}{k}_1^P& 0& 0\\ 0& k_2^P& 0\\ 0& 0& k_3^P\end{array}\right)$

summarises the assumptions of these models. We show both the Q and P dynamics and will cover the P dynamics in Section 2.3. The independent factor models have an identity matrix for $D\left(X_t,t\right)$ in Equation (5). The long-term mean ${\theta }^Q$ in the models is assumed to be a vector of zeros for both models, as explained in Blackburn and Sherris (2013).

The model dynamics are given as follows:

• The independent Blackburn–Sherris model has instantaneous mortality rate given by

  $$\mu \left(t\right)=X_t^1+X_t^2+X_t^3,$$

  (10)

  with ${\rho }_1={\left(1,1,1\right)}^T$ and $X_t=\left(X_t^1,X_t^2,X_t^3\right)$ in Equation (4).
  The dynamics of the state variables $X_t$ have the following form under the risk-neutral measure Q

  $$\begin{array}{c}\hfill \left(\begin{array}{c}d{X}_t^1\\ d{X}_t^2\\ d{X}_t^3\end{array}\right)=-\left(\begin{array}{ccc}{\delta }_{11}& 0& 0\\ 0& {\delta }_{22}& 0\\ 0& 0& {\delta }_{33}\end{array}\right)\left(\begin{array}{c}{X}_t^1\\ X_t^2\\ X_t^3\end{array}\right)dt+\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ 0& {\sigma }_{22}& 0\\ 0& 0& {\sigma }_{33}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,Q}\\ d{W}_t^{2,Q}\\ d{W}_t^{3,Q}\end{array}\right).\end{array}$$

  (11)
• The independent AFNS mortality model has an instantaneous mortality rate given by

  $$\mu \left(t\right)=L_t+S_t,$$

  (12)

  with ${\rho }_1={\left(1,1,0\right)}^T$ and $X_t=\left(L_t,S_t,C_t\right)$ in Equation (4).
  The dynamics of the factors under the Q-measure are given by:

  $$\begin{array}{c}\hfill \left(\begin{array}{c}d{L}_t\\ d{S}_t\\ d{C}_t\end{array}\right)=-\left(\begin{array}{ccc}0& 0& 0\\ 0& \delta & -\delta \\ 0& 0& \delta \end{array}\right)\left(\begin{array}{c}{L}_t\\ S_t\\ C_t\end{array}\right)dt+\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ 0& {\sigma }_{22}& 0\\ 0& 0& {\sigma }_{33}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,Q}\\ d{W}_t^{2,Q}\\ d{W}_t^{3,Q}\end{array}\right).\end{array}$$

  (13)

The solutions for the survival curve require $B\left(t,T\right)$, the factor loadings, and $A\left(t,T\right)$ from Equations (7) and (8), which can be explicitly solved. The results are as follows:

• The independent Blackburn–Sherris model (Blackburn and Sherris 2013)

  $$B^j\left(t,T\right)=-\frac{1-e^{-{\delta }_{jj}\left(T-t\right)}}{{\delta }_{jj}},j=1,2,3,$$

  (14)

  $$A\left(t,T\right)=\frac{1}{2}\sum _{j=1}^3\frac{{\sigma }_{jj}^2}{{\delta }_{jj}^3}\left[\frac{1}{2}\left(1-e^{-2{\delta }_{jj}\left(T-t\right)}\right)-2\left(1-e^{-{\delta }_{jj}\left(T-t\right)}\right)+{\delta }_{jj}\left(T-t\right)\right].$$

  (15)
• The independent AFNS model (Christensen et al. 2011)

  $$\begin{array}{cc}\hfill B^1\left(t,T\right)& =-\left(T-t\right),B^2\left(t,T\right)=-\frac{1-e^{-\delta \left(T-t\right)}}{\delta },\hfill \\ \hfill B^3\left(t,T\right)& =\left(T-t\right)e^{-\delta \left(T-t\right)}-\frac{1-e^{-\delta \left(T-t\right)}}{\delta },\hfill \end{array}$$

  (16)

  $$\begin{array}{cc}\hfill \frac{A\left(t,T\right)}{T-t}& ={\sigma }_{11}^2\frac{\left(T-t\right)}{6}+{\sigma }_{22}^2\left[\frac{1}{2{\delta }^2}-\frac{1}{{\delta }^3}\frac{1-e^{-\delta \left(T-t\right)}}{T-t}+\frac{1}{4{\delta }^3}\frac{1-e^{-2\delta \left(T-t\right)}}{T-t}\right]+\hfill \\ & {\sigma }_{33}^2\left[\frac{1}{2{\delta }^2}+\frac{1}{{\delta }^2}{e}^{-\delta \left(T-t\right)}-\frac{1}{4\delta }\left(T-t\right)e^{-2\delta \left(T-t\right)}-\frac{3}{4{\delta }^2}{e}^{-2\delta \left(T-t\right)}\right.\hfill \\ & \left.-\frac{2}{{\delta }^3}\frac{1-e^{-\delta \left(T-t\right)}}{T-t}+\frac{5}{8{\delta }^3}\frac{1-e^{-2\delta \left(T-t\right)}}{T-t}\right].\hfill \end{array}$$

  (17)

The independent Blackburn–Sherris model factor loadings all have the same functional form but differ in the values of the fitted ${\delta }_{jj}$. This parameter results in different factor impacts across ages for the cohort mortality curve. The independent AFNS model’s factor loadings have direct interpretation. $B^1\left(t,T\right)$, the level factor, is constant, so it has the same impact on all ages in the cohort mortality curve. $B^2\left(t,T\right)$, the slope factor, is increasing, so it impacts mortality rates at older ages more than at younger ages, resulting in changes in slope. $B^3\left(t,T\right)$, the curvature factor, is decreasing and, along with the other factors, produces curvature in the dynamics of the survival curve through time.

The selection of ${\rho }_1$ and the structure of the mean reversion matrix $K^Q$ in Equation (13) ensures that the factor loadings $-\frac{B\left(t,T\right)}{T-t}$ (Equation (9)) of the AFNS model have a consistent Nelson–Siegel structure through time with the latent factors interpreted as level, slope, and curvature factors driving changes in the mortality curves (Diebold and Li 2006; Diebold and Rudebusch 2013). Since the AFNS model is based on the ATSM in Duffie and Kan (1996), this model maintains an arbitrage-free affine structure when applied to available longevity linked market prices, making it suitable for financial and pricing applications (Christensen et al. 2011).

Björk and Christensen (1999) argues that the Nelson–Siegel model does not satisfy the consistency requirement (proposed by Björk and Christensen 1999). However, Diebold and Rudebusch (2013) explain that not meeting the consistency requirement is met by the AFNS model through the yield-adjustment term, which is $-\frac{A\left(t,T\right)}{T-t}$.

2.2.2. Dependent Factor Models with Gaussian Processes

Dependent factor models can improve how the models capture the correlation in mortality rate dynamics between different ages. We set out the factor dynamics with solutions to the dependent factor models for our Gaussian models. The Blackburn–Sherris model and the AFNS model assume that the volatility matrix $\mathsf{\Sigma }$ is lower-triangular, allowing correlated shocks in the models. Correlation can also be incorporated through $K^Q$; however, for the AFNS model, the structure of $K^Q$ has to be the same as for the independent factor model to preserve the Nelson–Siegel structure for the model factors.

The risk-neutral dynamics of the dependent factor models are

• The dependent Blackburn–Sherris model

  $$\begin{array}{c}\hfill \left(\begin{array}{c}d{X}_t^1\\ d{X}_t^2\\ d{X}_t^3\end{array}\right)=-\left(\begin{array}{ccc}{\delta }_{11}& 0& 0\\ {\delta }_{21}& {\delta }_{22}& 0\\ {\delta }_{31}& {\delta }_{32}& {\delta }_{33}\end{array}\right)\left(\begin{array}{c}{X}_t^1\\ X_t^2\\ X_t^3\end{array}\right)dt+\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ {\sigma }_{21}& {\sigma }_{22}& 0\\ {\sigma }_{31}& {\sigma }_{32}& {\sigma }_{33}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,Q}\\ d{W}_t^{2,Q}\\ d{W}_t^{3,Q}\end{array}\right).\end{array}$$

  (18)
• The dependent AFNS model

  $$\begin{array}{c}\hfill \left(\begin{array}{c}d{L}_t\\ d{S}_t\\ d{C}_t\end{array}\right)=-\left(\begin{array}{ccc}0& 0& 0\\ 0& \delta & -\delta \\ 0& 0& \delta \end{array}\right)\left(\begin{array}{c}{L}_t\\ S_t\\ C_t\end{array}\right)dt+\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ {\sigma }_{21}& {\sigma }_{22}& 0\\ {\sigma }_{31}& {\sigma }_{32}& {\sigma }_{33}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,Q}\\ d{W}_t^{2,Q}\\ d{W}_t^{3,Q}\end{array}\right).\end{array}$$

  (19)

The factor loadings $B\left(t,T\right)$ and $A\left(t,T\right)$ for the dependent Blackburn–Sherris model are provided in Appendix A. For the dependent AFNS model, the factor loadings $B\left(t,T\right)$ are the same as for the independent factor model and an explicit expression for $A\left(t,T\right)$ is given in Christensen et al. (2011).

2.2.3. The Cox–Ingersoll–Ross Mortality Model

To date, the Gaussian models allow mortality rates to become negative, although, as we will show later, this effect is empirically small for the models we consider. We also consider a multi-factor affine mortality model with each factor following a square-root process, based on the Cox–Ingersoll–Ross model (CIR) (Cox et al. 1985) frequently used as affine term structure models for interest rates and credit risk to avoid negative mortality rates.

The Cox–Ingersoll–Ross model (CIR) model can also be seen as a model that can capture the effect of mortality heterogeneity. For the CIR mortality model, mortality rates follow a non-central Chi-square distribution and are asymptotically Gamma distributed (Cox et al. 1985).

Following Chen and Scott (2003) and Geyer and Pichler (1999) for interest rates, we define the instantaneous mortality intensity to be affine with:

$${\mu }_x^i\left(t\right)={\rho }_1^{\prime }{X}_t=X_t^1+X_t^2+X_t^3,$$

(20)

where $X_t=\left(X_t^1,X_t^2,X_t^3\right)$ are the state variables that are driving the mortality intensity and ${\rho }_1$ is assumed to be ${\left(1,1,1\right)}^T$.

The factor dynamics driving the mortality survival curve are then given by the following system of SDEs under the risk-neutral measure Q:

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill \left(\begin{array}{c}d{X}_t^1\\ d{X}_t^2\\ d{X}_t^3\end{array}\right)& =-\left(\begin{array}{ccc}{\delta }_{11}& 0& 0\\ 0& {\delta }_{22}& 0\\ 0& 0& {\delta }_{33}\end{array}\right)\left[\left(\begin{array}{c}{\theta }_1^Q\\ {\theta }_2^Q\\ {\theta }_3^Q\end{array}\right)-\left(\begin{array}{c}{X}_t^1\\ X_t^2\\ X_t^3\end{array}\right)\right]dt\hfill \\ & +\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ 0& {\sigma }_{22}& 0\\ 0& 0& {\sigma }_{33}\end{array}\right)\left(\begin{array}{ccc}\sqrt{X_t^1}& 0& 0\\ 0& \sqrt{X_t^2}& 0\\ 0& 0& \sqrt{X_t^3}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,Q}\\ d{W}_t^{2,Q}\\ d{W}_t^{3,Q}\end{array}\right).\hfill \end{array}\end{array}$$

(21)

Each factor follows a single-factor CIR square-root process. The matrix $D\left(X_t,t\right)$ in Equation (5) is defined as a diagonal matrix with the j-th element on the diagonal as $\sqrt{X_t^j}$ ($j=1,2,3$).

The explicit expressions for $B\left(t,T\right)$ and $A\left(t,T\right)$ are

$$B^j\left(t,T\right)=-\frac{2\left(e^{{\gamma }_j\left(T-t\right)}-1\right)}{\left({\delta }_{jj}+{\gamma }_j\right)\left(e^{{\gamma }_j\left(T-t\right)}-1\right)+2{\gamma }_j},j=1,2,3,$$

(22)

$$A\left(t,T\right)=\sum _{j=1}^3\frac{2{\delta }_{jj}{\theta }_j^Q}{{\sigma }_{jj}^2}ln\left[\frac{2{\gamma }_{j}exp\left(\frac{\left({\delta }_{jj}+{\gamma }_j\right)\left(T-t\right)}{2}\right)}{\left({\delta }_{jj}+{\gamma }_j\right)\left(e^{{\gamma }_j\left(T-t\right)}-1\right)+2{\gamma }_j}\right],$$

(23)

with ${\gamma }_j=\sqrt{{\delta }_{jj}^2+2{\sigma }_{jj}^2},j=1,2,3$ (Duan and Simonato 1999; Chen and Scott 2003; Geyer and Pichler 1999).

2.3. Real-World Dynamics and Change of Measure

We develop the affine mortality models under the risk-neutral measure Q. Since we fit the models using historical data, we change this measure to the real-world measure P. From Girsanov’s theorem, the relationship between the dynamics under the P and Q measure is given by:

$$d{W}_t^Q=d{W}_t^P+{\Lambda }_{t}dt,$$

(24)

where ${\Lambda }_t$ is the risk premium vector.

To specify the structure of the risk premium of longevity risk, we adopt the essentially affine model proposed by Duffee (2002). The essentially affine model removes the strong link between the factor loadings $-\frac{B\left(t,T\right)}{T-t}$ and the drift term under the real-world measure (Blackburn and Sherris 2013) and preserves the affine dynamics under the P-measure (Christensen et al. 2011).

The form of the risk premium is (Duffee 2002):

$${\Lambda }_t=\left\{\begin{array}{cc}{\lambda }^0+{\lambda }^1{X}_t,\hfill & \mathrm{for}\mathrm{models}\mathrm{with}\mathrm{Gaussian}\mathrm{processes};\hfill \\ D\left(X_t,t\right){\lambda }^0,\hfill & \mathrm{for}\mathrm{the}\mathrm{CIR}\mathrm{model}.\hfill \end{array}\right.$$

(25)

where ${\Lambda }_t\in {\mathbb{R}}^{n\times 1}$, ${\lambda }^0\in {\mathbb{R}}^{n\times 1}$ and ${\lambda }^1\in {\mathbb{R}}^{n\times n}$.

With these assumptions, the SDEs for factors under the measure P can be written as:

$$d{X}_t=\left\{\begin{array}{cc}{K}^P\left[{\theta }^P-X_t\right]dt+\mathsf{\Sigma }d{W}_t^P,\hfill & \mathrm{for}\mathrm{models}\mathrm{with}\mathrm{Gaussian}\mathrm{processes};\hfill \\ K^P\left[{\theta }^P-X_t\right]dt+\mathsf{\Sigma }D\left(X_t,t\right)d{W}_t^P,\hfill & \mathrm{for}\mathrm{the}\mathrm{CIR}\mathrm{model}.\hfill \end{array}\right.$$

(26)

The forms of $K^P$ and ${\theta }^P$ are derived in Appendix B. In the essentially affine model, we are free to choose the mean reversion matrix $K^P$ and the mean vector ${\theta }^P$.

Under the real-world measure P, the dynamics of the factors in each model that we will estimate from historical mortality data are:

• The independent Blackburn–Sherris model (Blackburn and Sherris 2013)

  $$\begin{array}{c}\hfill \left(\begin{array}{c}d{X}_t^1\\ d{X}_t^2\\ d{X}_t^3\end{array}\right)=-\left(\begin{array}{ccc}{k}_{11}^P& 0& 0\\ 0& k_{22}^P& 0\\ 0& 0& k_{33}^P\end{array}\right)\left(\begin{array}{c}{X}_t^1\\ X_t^2\\ X_t^3\end{array}\right)dt+\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ 0& {\sigma }_{22}& 0\\ 0& 0& {\sigma }_{33}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,P}\\ d{W}_t^{2,P}\\ d{W}_t^{3,P}\end{array}\right).\end{array}$$

  (27)
• The independent AFNS model (Christensen et al. 2011)

  $$\begin{array}{c}\hfill \left(\begin{array}{c}d{L}_t\\ d{S}_t\\ d{C}_t\end{array}\right)=-\left(\begin{array}{ccc}{k}_{11}^P& 0& 0\\ 0& k_{22}^P& 0\\ 0& 0& k_{33}^P\end{array}\right)\left(\begin{array}{c}{L}_t\\ S_t\\ C_t\end{array}\right)dt+\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ 0& {\sigma }_{22}& 0\\ 0& 0& {\sigma }_{33}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,P}\\ d{W}_t^{2,P}\\ d{W}_t^{3,P}\end{array}\right).\end{array}$$

  (28)
• The dependent Blackburn–Sherris model

  $$\begin{array}{c}\hfill \left(\begin{array}{c}d{X}_t^1\\ d{X}_t^2\\ d{X}_t^3\end{array}\right)=-\left(\begin{array}{ccc}{k}_{11}^P& 0& 0\\ 0& k_{22}^P& 0\\ 0& 0& k_{33}^P\end{array}\right)\left(\begin{array}{c}{X}_t^1\\ X_t^2\\ X_t^3\end{array}\right)dt+\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ {\sigma }_{21}& {\sigma }_{22}& 0\\ {\sigma }_{31}& {\sigma }_{32}& {\sigma }_{33}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,P}\\ d{W}_t^{2,P}\\ d{W}_t^{3,P}\end{array}\right).\end{array}$$

  (29)
• The dependent AFNS model

  $$\begin{array}{c}\hfill \left(\begin{array}{c}d{L}_t\\ d{S}_t\\ d{C}_t\end{array}\right)=-\left(\begin{array}{ccc}{k}_{11}^P& 0& 0\\ 0& k_{22}^P& 0\\ 0& 0& k_{33}^P\end{array}\right)\left(\begin{array}{c}{L}_t\\ S_t\\ C_t\end{array}\right)dt+\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ {\sigma }_{21}& {\sigma }_{22}& 0\\ {\sigma }_{31}& {\sigma }_{32}& {\sigma }_{33}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,P}\\ d{W}_t^{2,P}\\ d{W}_t^{3,P}\end{array}\right).\end{array}$$

  (30)
• The CIR model

  $$\begin{array}{c}\hfill \begin{array}{cc}\hfill \left(\begin{array}{c}d{X}_t^1\\ d{X}_t^2\\ d{X}_t^3\end{array}\right)& =-\left(\begin{array}{ccc}{k}_{11}^P& 0& 0\\ 0& k_{22}^P& 0\\ 0& 0& k_{33}^P\end{array}\right)\left[\left(\begin{array}{c}{\theta }_1^P\\ {\theta }_2^P\\ {\theta }_3^P\end{array}\right)-\left(\begin{array}{c}{X}_t^1\\ X_t^2\\ X_t^3\end{array}\right)\right]dt\hfill \\ & +\left(\begin{array}{ccc}{\sigma }_{11}& 0& 0\\ 0& {\sigma }_{22}& 0\\ 0& 0& {\sigma }_{33}\end{array}\right)\left(\begin{array}{ccc}\sqrt{X_t^1}& 0& 0\\ 0& \sqrt{X_t^2}& 0\\ 0& 0& \sqrt{X_t^3}\end{array}\right)\left(\begin{array}{c}d{W}_t^{1,P}\\ d{W}_t^{2,P}\\ d{W}_t^{3,P}\end{array}\right).\hfill \end{array}\end{array}$$

  (31)

3. Mortality Data

We use USA age-cohort mortality data, as opposed to age-period data, from the Human Mortality Database (2018) (HMD) to calibrate and compare the mortality models. Our focus is on actuarial and financial applications, which generally require age-cohort mortality models rather than age-period models (Blackburn and Sherris 2013; Xu et al. 2020a; Chang and Sherris 2018). Age-period models are often modified to include a cohort factor to capture cohort effects. Our approach is to model age-cohort data directly.

We use mortality data for males from ages 50 to 100 for the cohorts born from 1883 to 1915, since these are the complete cohort data we have available. We are interested in the older ages for post-retirement applications and require full mortality-rate data for all ages for each complete cohort. Cohort death rates are derived from the age-period mortality rates by using the diagonals.

Historical survival probabilities, $S^i\left(x;t,T\right)$, are derived from the data along with historical average forces of mortality ${\overline{\mu }}^i\left(x;t,T\right)$ over the period $\tau =T-t$ for each cohort i aged x at time t, using:

$$S^i\left(x;t,T\right)=\prod _{s=1}^{T-t}\left[1-q^i\left(x+s-1,t+s-1\right)\right],$$

(32)

$${\overline{\mu }}^i\left(x;t,T\right)=-\frac{1}{T-t}log\left[S^i\left(x;t,T\right)\right],$$

(33)

where $q^i\left(x,t\right)$ is the one-year death probability for an individual aged x at time t in cohort i.

Figure 1 gives the average force of mortality for cohorts born between 1883 and 1915, aged 50 to 100. Mortality improvement across cohorts is seen from the downward trend of the average force of mortality at each age. The rate of mortality improvement differs by age. The average force of mortality of each cohort grows exponentially in each cohort.

Figure 2 gives a principal component analysis (PCA) for the change in mortality intensity for all cohorts as they age. The first three principal components are able to explain approximately 90% of the total variance. This is consistent with the results in Blackburn and Sherris (2013) and supports our choice of three-factor affine mortality models.

4. Model Assessment and Comparison

We use the Kalman filter to estimate the model parameters for all the models. We then compare the fitted models using a number of model selection criteria. When using the Kalman filter, we estimate the parameters of the latent factor dynamics in the real-world or historical measure and estimate the Q measure parameters in the measurement equation that links the cross-sectional survival probabilities to the latent factors.

4.1. Parameter Estimation

We follow Christensen et al. (2011) and Blackburn and Sherris (2013) and use the Kalman filter with maximum-likelihood estimation to estimate the parameters in the affine mortality models. Our models capture the volatility of the underlying mortality rates but observed deaths are used to estimate historical mortality rates from the data. These historical mortality rates also include Poisson variation based on the number of individuals in each age. Our model needs to include this. We do so by including an exponentially increasing Poisson variation term in the measurement equation of the Kalman filter (Xu et al. 2020b).

The estimation process is as follows:

• Represent the affine mortality models in the state space form, which consists of two components, the measurement equation and the state transition Equation (Xu et al. 2020a; Shumway and Stoffer 2017).
  The measurement equation describes the affine relationship between the average force of mortality and the state variables (Xu et al. 2020a; Durbin and Koopman 2012). Based on Blackburn and Sherris (2013) and Xu et al. (2020a), the measurement equation in terms of the average forces of mortality is

  $$\overline{\mu }\left(t,T\right)=-\frac{B{\left(t,T\right)}^{\prime }}{T-t}{X}_t-\frac{A\left(t,T\right)}{T-t}+{\epsilon }_t,{\epsilon }_t\sim N\left(0,H\right),$$

  (34)

  where the measurement error ${\epsilon }_t$ is independently and identically distributed noise with the covariance matrix of the measurement error, H, being diagonal. We use the Q or pricing measure parameters in the measurement equation to ensure cross-sectional consistency in the survival probabilities.
  To capture the increasing nature of the Poisson variation, the parametric form assumed for the diagonal of the covariance matrix H is

  $$H\left(t,T\right)=\frac{1}{T-t}\sum _{i=1}^{T-t}\left[r_c+r_1{e}^{r_{2}i}\right],$$

  (35)

  where the values of $r_c$, $r_1$ and $r_2$ are estimated from the data.
  The state transition equation represents the unobserved dynamics of the state variables (Xu et al. 2020a; Durbin and Koopman 2012) and is given by:

  $$X_t=exp\left(-K^P\right)X_{t-1}+{\eta }_t,{\eta }_t\sim N\left(0,R\right),$$

  (36)

  where ${\eta }_t$ is the transition error vector with diagonal matrix R, the covariance matrix of the transition error. We use the real-world or historical parameters for the latent factor dynamics in the state transition equation to capture the real-world dynamics.
  The matrix R has the following structure:

  $$R={\int }_{t-1}^t{e}^{-K^P\left(t-s\right)}\mathsf{\Sigma }{\mathsf{\Sigma }}^{\prime }{e}^{-{\left(K^P\right)}^{\prime }\left(t-s\right)}ds.$$

  (37)
• Use the Kalman filter to evaluate the likelihood function of affine mortality models and to extract the values of the state variables. The information available at time t is denoted by $Y_t=\left(y_1,\dots ,y_t\right)$, and the model parameters are given by $\psi$.
  In the forecasting step, we use the state update $X_{t-1}$ and its mean square error ${\mathsf{\Sigma }}_{t-1}$ obtained at $t-1$ to obtain

  $$X_{t|t-1}=E\left[X_t|Y_{t-1}\right]=\Phi \left(\psi \right)X_{t-1},$$

  (38)

  $${\mathsf{\Sigma }}_{t|t-1}=\Phi \left(\psi \right){\mathsf{\Sigma }}_{t-1}\Phi {\left(\psi \right)}^{\prime }+R\left(\psi \right),$$

  (39)

  where $\Phi =exp\left(-K^P\right)$.
  In the update step, the information at time t, $Y_t$, is used to update the forecasts $X_{t|t-1}$, and we obtain:

  $$X_t=E\left[X_t|Y_t\right]=X_{t|t-1}+{\mathsf{\Sigma }}_{t|t-1}B{\left(\psi \right)}^{\prime }{F}_t^{-1}{\nu }_t,$$

  (40)

  $${\mathsf{\Sigma }}_t={\mathsf{\Sigma }}_{t|t-1}-{\mathsf{\Sigma }}_{t|t-1}B{\left(\psi \right)}^{\prime }{F}_t^{-1}B\left(\psi \right){\mathsf{\Sigma }}_{t|t-1},$$

  (41)

  where

  $${\nu }_t=y_t-E\left[y_t|Y_{t-1}\right]=y_t-A\left(\psi \right)-B\left(\psi \right)X_{t|t-1},$$

  (42)

  $$F_t=cov\left({\nu }_t\right)=B\left(\psi \right){\mathsf{\Sigma }}_{t|t-1}B{\left(\psi \right)}^{\prime }+H\left(\psi \right).$$

  (43)
• Evaluate the following log-likelihood function with the values obtained in the previous step:

  $$logL\left(y_1,\dots ,y_t;\psi \right)=\sum _{t=1}^T\left(-\frac{N}{2}log\left(2\pi \right)-\frac{1}{2}log\left|F_t\right|-\frac{1}{2}{\nu }_t^{\prime }{F}_t{\nu }_t\right),$$

  (44)

  where N is the number of observed average forces of mortality.
  The log-likelihood function is maximized with respect to $\psi$ to obtain the optimal parameter set. For the CIR mortality model we use quasi-maximum likelihood estimation. Jevtić and Regis (2021) estimate similar models using the Kalman filter and quasi-maximum likelihood estimation. The estimation uses the Gaussian Kalman filter with the moments from the CIR model in the likelihood for parameter estimation.

4.2. Model Parameter Estimation Results

Parameter estimates for each model, along with the standard errors are summarized in

Table 2. Estimated Parameters.

	The Blackburn–Sherris Model		The AFNS Model		The CIR Model
	Independent Factor	Dependent Factor	Independent Factor	Dependent Factor
${\delta }_{11}$ (AFNS: $\delta$)	−0.01106 (0.00123)	−0.20183 ($2.944\times {10}^{-5}$)	−0.08348 ($1.580\times {10}^{-4}$)	−0.04725 ($8.096\times {10}^{-5}$)	−0.09652 ($2.513\times {10}^{-4}$)
${\delta }_{21}$	-	0.56206 ($4.091\times {10}^{-5}$)	-	-	-
${\delta }_{22}$	0.07484 (0.00432)	−0.07092 ($1.766\times {10}^{-5}$)	-	-	0.12627 ($2.412\times {10}^{-3}$)
${\delta }_{31}$	-	0.24075 ($1.555\times {10}^{-5}$)	-	-	-
${\delta }_{32}$	-	0.80809 ($4.102\times {10}^{-5}$)	-	-	-
${\delta }_{33}$	−0.06883 ($2.452\times {10}^{-4}$)	0.77825 ($1.461\times {10}^{-5}$)	-	-	−0.11153 ($3.060\times {10}^{-4}$)
$k_{11}^P$	0.38753	−0.04248	0.18793	0.01810	0.00077
$k_{22}^P$	0.13910	0.01869	0.01361	0.02002	0.59402
$k_{33}^P$	0.00718	0.01827	0.02701	0.04972	0.06842
${\sigma }_{11}$	0.00782	$7.557\times {10}^{-11}$	$9.593\times {10}^{-4}$	0.00400	0.00265 ($6.894\times {10}^{-5}$)
${\sigma }_{21}$	-	0.01110	-	−0.00387	-
${\sigma }_{22}$	0.00125	$3.370\times {10}^{-11}$	$1.120\times {10}^{-4}$	0.00091	0.02848 ($1.250\times {10}^{-3}$)
${\sigma }_{31}$	-	−0.01190	-	−0.00183	-
${\sigma }_{32}$	-	0.00047	-	0.00123	-
${\sigma }_{33}$	$5.409\times {10}^{-4}$	0.00029	$3.549\times {10}^{-5}$	0.00023	0.01360 ($9.494\times {10}^{-5}$)
$r_1$	$1.071\times {10}^{-11}$	$4.337\times {10}^{-8}$	$1.422\times {10}^{-10}$	$6.272\times {10}^{-8}$	$5.498\times {10}^{-10}$
$r_2$	0.37797	0.11375	0.17784	0.10742	$6.646\times {10}^{-7}$
$r_c$	$4.360\times {10}^{-8}$	$5.705\times {10}^{-8}$	$4.963\times {10}^{-7}$	$4.636\times {10}^{-13}$	$3.410\times {10}^{-7}$
The CIR Model
${\theta }_1^Q$	0.00080	${\theta }_2^Q$	0.01010	${\theta }_3^Q$	0.00137
${\theta }_1^P$	0.00697	${\theta }_2^P$	0.00415	${\theta }_3^P$	0.00356

. Risk neutral model parameter values are reported along with the historical or real world parameter values. Transitions are based on the historical parameters to reflect the dynamics driving the underlying factors and the measurement equation is based on measurement equation to ensure the arbitrage-free structure of the survival curves.

The $\delta$ parameters in $K^Q$ determine the impact of each of the factors and the significance of the factor loadings for mortality rates at different ages. We see that the diagonal components are largely negative. The ${\delta }^{\prime }$s in the dependent-factor AFNS model are smaller in absolute value than for the independent-factor AFNS model. For the two AFNS models, ${\delta }^{\prime }$s are both negative, impacting the sensitivity for the slope factor at older ages.

The mean reversion $K^P$ parameters, giving the speed of reversion to long-term means, vary across the models. There are different rates of mean reversion in the real world measure. There are negative correlations between factors in the two dependent-factor models, which only impacts the adjustment term $A\left(t,T\right)$.

All the long-term mean $\theta$ parameters in the CIR model are positive. This ensures the factors are positive and positive mortality rates in the model. The second factor, $X^2$, has the largest mean reversion speed, $k_{22}^P$, and largest volatility, ${\sigma }_{22}$. This factor has more impact on the short term (Geyer and Pichler 1999). The mean reversion rate $k_{11}^P$ and the volatility ${\sigma }_{11}$ of $X^1$ are lowest, compared with the other two factors, so the first factor has less impact on mortality dynamics and is less volatile.

As seen in the structure of the measurement error matrix H in Equation (35), the measurement errors are age-dependent and exponentially increasing with age. For the $r_2$, the scalar in the exponential function in matrix H, the independent-factor Blackburn–Sherris model has the largest $r_2$, and the largest estimated measurement error volatility. For the CIR model, the values of all the parameters in matrix H are negligible. These smaller measurement errors reflects in a better in-sample model fit.

4.3. Assessing Model Goodness-of-Fit

Table 3. Comparison of Affine Mortality Models.

	The Blackburn–Sherris Model		The AFNS Model		The CIR Model
	Independent Factor	Dependent Factor	Independent Factor	Dependent Factor
Log Likelihood	9896.419	9938.696	9665.801	9887.878	10,045.70
RMSE	0.00250	$7.601\times {10}^{-4}$	$6.856\times {10}^{-4}$	$9.160\times {10}^{-4}$	$5.227\times {10}^{-4}$
No. of Parameters	12	18	10	13	18
AIC	−19,570.837	−19,643.392	−19,113.602	−19,551.757	−19,857.40
BIC	−18,968.292	−19,008.277	−18,521.914	−18,943.783	−19,222.29
Probability of Negative Mortality	0.02700	$1.011\times {10}^{-32}$	$1.722\times {10}^{-31}$	$4.34\times {10}^{-14}$	-

shows, for each model, the Root Mean Square Error (RMSE), the Akaike information criterion (AIC), and the Bayesian information criterion (BIC). Since the models with Gaussian processes allow for negative mortality, we show the probabilities of negative mortality for these models.

The CIR model has the highest log-likelihood and the smallest RMSE. Although it has more parameters, the AIC and the BIC of the CIR model indicate this is a better model. As noted earlier, the CIR model precludes the probability of negative mortality.

We note that the Gaussian models perform well, particularly the dependent-factor models. The dependent factor Blackburn–Sherris model and the AFNS models all have low probabilities of negative mortality rates. The dependent-factor Blackburn–Sherris model has the largest log-likelihood and better AIC and BIC than the other Gaussian models.

4.4. Factors and Factor Loadings

It is interesting to consider the factors and the factor loadings for the cohort survival curves for the independent AFNS mortality model, where the factors for level, slope, and curvature have a direct interpretation, and the CIR mortality model which, based on the criteria used, is the best-performing model.

Figure 3 and Figure 4 show the fitted values of the factors and factor loadings of the independent AFNS model. For the factor loading $B_1$, the impact of the level factor L is constant across all ages. We see an increase in the factor level for all ages for the cohorts born around 1900 onwards. This increase in the level factor is offset by the reduction in the slope factor. The interaction is impacted by the factor loadings for these factors and the resulting changes in mortality curves for cohorts born later.

The slope factor loading $B_2$ increases exponentially, so it impacts older ages more than younger ages. Mortality rates at older ages are more sensitive to the slope factor S. For cohorts born after 1900, corresponding to the rise in the level factor L, there is a decline in the slope factor S. These factors interact to fit the changes in the historical age-cohort survival curves.

The factor loading $B_3$ is negative and decreasing across all ages. As a result, the convexity of the survival curve at older ages decreases faster than at younger ages. For cohorts born after 1900, the decline in C results in mortality rate curves that are less convex across age. As a result, mortality improvement at older ages is larger than for the younger ages in the age-cohort curves.

The adjustment term A in the survival curve is negative and decreasing.

The estimated latent factors and factor loadings for the CIR mortality model are shown in Figure 5 and Figure 6. As expected, the factors and factor loadings differ from those of the independent AFNS mortality model, reflecting the different model assumptions.

The first factor $X^1$ is relatively stable for cohorts born before 1900, then increases, followed by a moderate decline for cohorts born around 1910 and after. The factor loading, $B_1$, is positive and increases with age. The interaction between the factor dynamics and the factor loading for the first factor $X^1$ results in older ages in the age-cohort mortality curve being affected more by the first factor than for younger ages. Since $X^1$ is the largest factor throughout the time period and the factor loadings are largest for this factor, it is the dominant factor impacting changes in mortality in the CIR model.

The second factor $X^2$ has a downward trend over time. The factor loading $B_2$ is decreasing with age and smaller than $B_1$. As a result, younger ages have mortality improvement from the second factor, but the size of improvement from this factor is smaller than for the first factor. As a result, it has a much smaller overall impact compared to $X^1$.

The third factor $X^3$ is relatively constant for cohorts born before 1900, decreasing afterwards, with a reduction in the rate of decrease for cohorts born after around 1910. The factor loading $B_3$ has a convex shape, is positive, and impacts older ages more than younger ages. This results in curvature changes in the age-cohort survival curve.

The adjustment term A, which produces the consistency in the age-cohort survival curves, is negative and, following a small increase for younger ages, is then decreasing.

The combination of the factor dynamics for all the factor dynamics along with the factor loadings determine the shape and dynamics of the age-cohort survival curves. These dynamics are estimated from the age-cohort mortality data and allow an analysis of the underlying trends and changes in the age-cohort survival curve over time. The more complex dynamics than single factor models explain the variation at older ages better.

4.5. Residual Analysis

The residuals or the affine mortality models are shown in Figure 7. The residuals shown are the differences between the average force of mortality from the historical mortality data and those determined from the fitted mortality models. Plots are on the same scale on the z-axis, except for the independent Blackburn–Sherris model, which has large residuals at older ages, reflecting a poorer fit at these ages.

Excluding the independent Blackburn–Sherris model, all of the other affine age-cohort mortality models show similar residuals. With three factors for level, slope, and curvature, we see that the independent AFNS model (Figure 7c) fits well. The AFNS model reduces the magnitude of residuals at the older ages compared to the independent Blackburn–Sherris model without adding additional parameters. This shows how, by using factors for level, slope, and curvature factors, the mortality model can capture variation in mortality curves, especially at older ages.

We see how factor dependence in the Blackburn–Sherris model produces a better fitting model, reducing the size of residuals, and accounting better for mortality variation at older ages. This is seen by comparing the dependent Blackburn–Sherris model (Figure 7b) with the independent model (Figure 7a). Although including dependence in the factors of the Blackburn–Sherris model improves the residuals of the model fit, we do not see this for the independent AFNS model in Figure 7c. The residuals in the dependent AFNS model (Figure 7d) are larger, particularly at older ages. This highlights how the independent AFNS model captures mortality variability better than the dependent AFNS model, in contrast to the dependent Blackburn–Sherris model.

The CIR model has lower residuals at older ages and residuals similar to those of the independent AFNS model in Figure 7c. The CIR model has slightly smaller residuals at older ages and ages younger than 60, compared with the other models.

We see a hump shape running diagonally across the cohorts in all of the residual plots. Diagonal effects in an age-cohort model correspond to period effects. The residuals highlight a period mortality effect that impacts all of the cohorts around the year 1970. This corresponds to when period mortality improvement trends in age-period models showed a change to a higher level of improvement. We also see that for later cohorts at older ages, the residuals are higher, which would be indicative of a slowing of mortality-improvement rates in recent years.

4.6. In-Sample Analysis

We now consider an in-sample model performance analysis. We compare estimated cohort survival probabilities from the fitted mortality models with the cohort survival probabilities from the historical data. Figure 8 summarizes the in-sample model fit using the Mean Absolute Percentage Error (MAPE) for each age across all cohorts. The scale of the percentage error is different above and below age 85 because of the significant differences in the size of the errors.

Figure 8a shows the MAPE for the affine mortality models with Gaussian processes. Below age 85, all models have similar performance, and the differences between the percentage errors of the different mortality models are relatively small. What is more interesting is the fit above age 85. The independent Blackburn–Sherris model has significantly larger percentage errors. The other affine age-cohort Gaussian mortality models, the dependent Blackburn–Sherris model and the independent AFNS model, have similar and improved model fit at these older ages.

We compare the MAPE of these better performing Gaussian mortality models with the CIR age-cohort mortality model in Figure 8b. We see that the CIR mortality model and the independent AFNS mortality model are similar. The latter has only slightly larger percentage errors at most ages. The dependent Blackburn–Sherris model is similar to the CIR model below age 75, but the percentage errors increase after age 75 with values as high as 10%.

Based on MAPE, we see that both the independent AFNS mortality model and the CIR mortality model have similar and satisfactory performance for the historical age-cohort data for complete cohorts in the USA historical mortality.

5. Forecasts of Survival Probabilities

We compare the predictive performance of the affine mortality models using an out-of-sample forecast with the fitted parameter values estimated from the cohorts born 1883 to 1915 for ages 50 to 100 used to forecast the survival curve of the cohort born in 1916. For this cohort we have full historical mortality data.

Following Christensen et al. (2011), who use optimal forecasts to predict yields to maturity, we use optimal forecasts, also referred to as best-estimate forecasts, to project average forces of mortality and survival probabilities.

At time t, the average force of mortality over $\tau =\left(T+1\right)-\left(t+1\right)$ periods at time $t+1$ for cohort i, ${\overline{\mu }}^i\left(t+1,T+1\right)$, is

$${\overline{\mu }}^i\left(t+1,T+1\right)=-\frac{B{\left(t,T\right)}^{\prime }}{T-t}E\left[X_{t+1}|X_t\right]-\frac{A\left(t,T\right)}{T-t},$$

(45)

where $B\left(t,T\right)$ and $A\left(t,T\right)$ only depend on $\tau =T-t$.

The forecasts of survival probabilities are then

$$S\left(t+1,T+1\right)=exp\left(B{\left(t,T\right)}^{\prime }E\left[X_{t+1}|X_t\right]+A\left(t,T\right)\right).$$

(46)

Since the factor dynamics under measure P in the independent Blackburn–Sherris model and the three-factor independent AFNS model are the same, the conditional expectations of state variables for these two models are as follows:

$$E\left[X_{t+1}^1|X_t^1\right]=e^{-k_{11}^P}{X}_t^1,E\left[X_{t+1}^2|X_t^2\right]=e^{-k_{22}^P}{X}_t^2,E\left[X_{t+1}^3|X_t^3\right]=e^{-k_{33}^P}{X}_t^3.$$

(47)

For the independent AFNS model, the conditional mean has the same structure but with $X_t=\left(L_t,S_t,C_t\right)$.

The SDEs describing the P-dynamics of the dependent Blackburn–Sherris model and the dependent AFNS model are the same as for the independent factor model in Equation (47).

The conditional mean of the CIR model is given in Geyer and Pichler (1999), so for the mortality model:

$$\begin{array}{c}E\left[X_{t+1}^1|X_t^1\right]=e^{-k_{11}^P}{X}_t^1+{\theta }_1^P\left(1-e^{-k_{11}^P}\right),E\left[X_{t+1}^2|X_t^2\right]=e^{-k_{22}^P}{X}_t^2+{\theta }_2^P\left(1-e^{-k_{22}^P}\right),\\ E\left[X_{t+1}^3|X_t^3\right]=e^{-k_{33}^P}{X}_t^3+{\theta }_3^P\left(1-e^{-k_{33}^P}\right).\end{array}$$

(48)

The RMSE under each mortality model for projecting the 1916 cohort survival curve are shown in

Table 4. RMSE by Comparing the Actual and Best-Estimate Survival Probabilities of the 1916 Cohort.

	The Blackburn–Sherris Model		The AFNS Model		The CIR Model
	Independent	Dependent	Independent	Dependent
RMSE	0.03197	0.00726	0.00668	0.00754	0.01835

. We see that the independent AFNS mortality model performs best on these criteria. The dependent AFNS mortality model and the dependent Blackburn–Sherris models perform similarly. In contrast, the independent Blackburn–Sherris mortality model shows the poorest performance. Although the CIR mortality model has reasonable RMSE, it is outperformed by the AFNS mortality models.

To illustrate these differences, Figure 9 shows the survival probabilities for the different mortality models using the best estimate forecasts compared to the actual survival probabilities from the historical mortality data. The models produce reasonable survival curve fits, except the independent Blackburn–Sherris mortality model and the CIR mortality model, which both underestimate the survival rates of the 1916 cohort.

Figure 10 shows the absolute percentage errors across ages for all the mortality models. This confirms the better forecasting performance of the independent AFNS mortality model. The dependent Blackburn–Sherris mortality models performs better than the independent Blackburn–Sherris mortality model, showing the benefit in forecasting of including correlations between the factors in this model. The level, slope, and curvature structure of the factors in the AFNS mortality model reduce the need for including dependence in the factors as compared to the Blackburn–Sherris mortality model.

6. Discussion

We have shown how affine mortality models can be applied to age-cohort mortality data using USA mortality data. We show that three-factor versions of these models fit the data well; the Gaussian models produce low probabilities of negative mortality rates; and the AFNS model, with level, slope, and curvature as factors, provides reliable forecasts of full age-cohort survival curves. We model age-cohort data since age-cohort survival curves are required for practical actuarial applications. Most other research uses age-period mortality, often with a cohort effect adjustment. Our analysis has used USA mortality data since the population is large and reflects many developments in mortality improvement expected to be found in other developed countries.

We develop our models in an arbitrage-free pricing framework assuming independence between interest rates and mortality. Correlation between interest rates and mortality rates can be incorporated into the models by an appropriate change in measure. Calibrating the models and estimating the correlation between interest rates and mortality rates would then require mortality-linked security prices that capture both interest rate and mortality risk.

Our results provide interesting insights into age-cohort mortality curve dynamics. Although we have focused on USA data, estimating and comparing age-cohort affine mortality models for other countries will provide a deeper understanding of the models and their ability to capture differing mortality dynamics.

There are a number of directions for which this research provides a foundation. We use the Kalman filter with maximum likelihood to estimate parameters. An area of research that would improve the model estimation is the development of efficient numerical estimation processes for the Kalman filter along with developing code that can be used by other researchers to estimate and implement the models.

We use only historical data for complete cohorts in our age-cohort mortality curve modelling. Incorporating incomplete cohorts into the estimation will allow the use of more recent age-cohort data in the model estimation and forecasting. Age-period data use the latest calendar year mortality data, whereas age-cohort data for later calendar years are incomplete.

We do not focus on parameter stability, and although this has been considered in (Blackburn and Sherris 2013), assessing this for several other countries would provide additional support for the empirical performance of the multi-factor affine models.

Although we include dependence between the factors by incorporating correlation, we do not include age effects in the factors. In the affine framework, it is possible to capture age dependence in the factors to improve the model fit at older ages with potentially fewer factors in the model. It is also possible to capture age correlations and cohort correlations more effectively using age-dependent factors.

An important issue to also consider is the extension of the models to incorporate jump events such as the COVID-19 pandemic. These events not only have an impact on all ages to a greater or lesser extent in a particular number of periods, but can also have longer lasting impacts, as is occurring with long COVID-19. The affine models we consider do not account for larger jumps that impact several periods across all ages, as would be the case for a pandemic. The models can incorporate jumps such as COVID-19 into the affine mortality framework.

These extensions are all topics that current CEPAR actuarial research is investigating.

7. Conclusions

We have applied several continuous-time affine mortality models to age-cohort survival curve data using USA historical age-cohort data to fit and assess the mortality models. We provide a comprehensive analysis and comparison of the models. We outline, compare, and assess several independent-factor and dependent-factor affine mortality models with Gaussian processes, including the Blackburn–Sherris mortality model (Blackburn and Sherris 2013; Christensen et al. 2011), as well as an affine mortality model with square-root processes (the CIR mortality model). The CIR mortality model precludes negative mortality rates that can occur in the Gaussian models. The CIR latent factors and the mortality intensity have non-central Chi-square distributions, which can better reflect mortality heterogeneity. We also assess the performance of an AFNS mortality model with interpretable latent stochastic factors for level, slope, and curvature of the survival curve.

We incorporate dependence in the Blackburn–Sherris mortality model by including correlation between the factors and show how this improves in-sample model fit and out-of-sample forecasting performance. The CIR mortality model shows the best in-sample model fit across a range of criteria, including model residuals. The superior in-sample performance of the CIR mortality model is likely to reflect the more realistic assumption of Gamma-distributed mortality rates compared to the Gaussian models.

We find that the independent-factor AFNS mortality model performs well. It better captures the variation in age-cohort mortality rates in USA data and produces a better fit at older ages than the independent-factor Blackburn–Sherris model. Negative mortality rates have very low probability in the AFNS mortality models.

We project a complete cohort based on the fitted models to assess forecasting performance of the models. We forecast the 1916 cohort in the USA data, for which we have a complete cohort of mortality data. The independent-factor AFNS mortality model shows better predictive performance compared to the other models.

Based on our detailed assessment and comparison of these affine age-cohort mortality models, we show that, for USA age-cohort mortality data, the independent AFNS model provides satisfactory model fit and satisfactory predictive performance. The model is parsimonious relative to other models and can be readily estimated using the Kalman filter, allowing for Poisson mortality variation in the measurement equation. The model has intuitive factor interpretations in terms of level, slope, and curvature for the dynamics of the mortality survival curve and is well suited for financial and insurance applications, including pricing and longevity risk management.