Age-period-cohort model

Description

Fits and forecasts mortality rates using age-period-cohort model.

Usage

APCS(
  x,
  M,
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

The age-period-cohort (APC) model is specified as

ln(m_{x,t}) = \alpha_x + \kappa_t + \gamma_{t-x} + \epsilon_{x,t}.

The model is estimated by Newton updating scheme and is forecasted by ARIMA applied to \kappa_t and \gamma_c. Constraints include sum of \kappa_t is zero and and sum of \gamma_c is zero. It can be applied to whole age range.

Value

An object of class APCS with associated S3 methods coef, forecast, plot, residuals, and simulate (nsim for setting number of simulations; seed for initialising random number generator).

References

Bray, I. (2002). Application of Markov chain Monte Carlo methods to projecting cancer incidence and mortality. Journal of the Royal Statistical Society Series C, 51(2), 151-164.

Examples

x <- 60:89
a <- c(-4.8499,-4.7676,-4.6719,-4.5722,-4.4847,-4.3841,-4.2813,-4.1863,-4.0861,-3.9962,
-3.8885,-3.7896,-3.6853,-3.5737,-3.4728,-3.3718,-3.2586,-3.1474,-3.0371,-2.9206,
-2.7998,-2.6845,-2.5653,-2.4581,-2.3367,-2.2159,-2.1017,-1.9941,-1.8821, -1.7697)
b <- rep(1/30,30)
k <- c(12.11,10.69,11.18,9.64,9.35,8.21,6.89,5.74,4.56,3.6,
3.27,2.04,1.11,-0.44,-1.05,-1.03,-1.84,-2.9,-4.03,-4.12,
-5.18,-5.64,-6,-6.51,-6.91,-6.9,-8.32,-8.53,-9.69,-9.31)
set.seed(123)
M <- exp(outer(k,b)+matrix(a,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.035))
fit <- APCS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

Common age effect model

Description

Fits and forecasts mortality rates of two populations using common age effect model.

Usage

CAES(
  x,
  M1,
  M2,
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

The common age effect (CAE) model is specified as

ln(m_{x,t,i}) = \alpha_{x,i} + \beta_x \kappa_{t,i} + \epsilon_{x,t,i}.

The model is estimated by Newton updating scheme and is forecasted by ARIMA applied to \kappa_{t,i}. Constraints include sum of \beta_x is one and sum of \kappa_{t,i} is zero. It can be applied to whole age range.

Value

An object of class CAES with associated S3 methods coef, forecast, plot (which = 1 gives parameter estimates; which = 2 gives residuals and forecasts), and residuals.

References

Kleinow, T. (2015). A common age effect model for the mortality of multiple populations. Insurance: Mathematics and Economics, 63(C), 147-152.

Examples

x <- 60:89
a1 <- c(-5.18,-5.12,-4.98,-4.92,-4.82,-4.73,-4.66,-4.53,-4.45,-4.35,
-4.26,-4.17,-4.05,-3.95,-3.84,-3.73,-3.65,-3.52,-3.40,-3.29,
-3.14,-3.02,-2.88,-2.76,-2.64,-2.49,-2.37,-2.25,-2.12,-2.00)
a2 <- c(-4.78,-4.68,-4.57,-4.49,-4.39,-4.29,-4.19,-4.10,-4.00,-3.89,
-3.80,-3.69,-3.60,-3.49,-3.39,-3.29,-3.17,-3.07,-2.96,-2.85,
-2.71,-2.62,-2.49,-2.37,-2.26,-2.14,-2.04,-1.91,-1.82,-1.72)
b <- c(0.0381,0.0340,0.0420,0.0389,0.0423,0.0414,0.0406,0.0393,0.0415,0.0400,
0.0411,0.0362,0.0387,0.0381,0.0384,0.0385,0.0356,0.0314,0.0317,0.0337,
0.0316,0.0298,0.0284,0.0270,0.0248,0.0262,0.0205,0.0215,0.0142,0.0145)
k1 <- c(8.68,8.34,7.99,6.87,8.18,5.73,4.83,5.20,2.74,3.22,
2.99,1.59,1.67,-0.65,-0.39,-1.07,-0.95,-2.78,-3.46,-2.45,
-4.12,-4.66,-4.98,-4.58,-6.30,-4.39,-5.56,-6.52,-8.26,-6.92)
k2 <- c(11.81,11.01,10.59,10.40,9.75,8.15,6.07,6.45,4.60,4.57,
4.15,1.49,1.77,-1.08,-1.44,-0.96,-1.66,-2.25,-4.67,-4.62,
-4.38,-6.37,-6.27,-6.91,-8.22,-7.35,-8.39,-7.87,-9.72,-8.65)
set.seed(123)
M1 <- exp(outer(k1,b)+matrix(a1,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.07))
M2 <- exp(outer(k2,b)+matrix(a2,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.07))
fit <- CAES(x=x,M1=M1,M2=M2,curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

CBD with cohort model

Description

Fits and forecasts mortality rates using CBD with cohort model.

Usage

CBDCS(
  x,
  M,
  curve = c("gompertz", "makeham", "perks", "weibull", "beard", "martinelle", "thatcher",
    "gompertz2", "makeham2", "perks2", "weibull2", "beard2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

The CBD with cohort (M6) model is specified as

ln(m_{x,t}) = \kappa_{1,t} + \kappa_{2,t} (x-\bar{x}) + \gamma_{t-x} + \epsilon_{x,t}.

The model is estimated by Newton updating scheme and is forecasted by ARIMA applied to \kappa_{1,t}, \kappa_{2,t}, and \gamma_c. Constraints include sum of \gamma_c is zero and sum of c\gamma_c is zero. It is designed for ages 50-90.

Value

An object of class CBDCS with associated S3 methods coef, forecast, plot, residuals, and simulate (nsim for setting number of simulations; seed for initialising random number generator).

References

Cairns, A.J.G., Blake, D., Dowd, K., Coughlan, G.D., Epstein, D., Ong, A., and Balevich, I. (2009). A quantitative comparison of stochastic mortality models using data from England and Wales and the United States. North American Actuarial Journal, 13(1), 1-35.

Examples

x <- 60:89
k1 <- -2.97-0.0245*(0:29)
k2 <- 0.101+0.000345*(0:29)
set.seed(123)
M <- exp(matrix(k1,nrow=30,ncol=30,byrow=FALSE)+outer(k2,(x-mean(x)))+rnorm(900,0,0.035))
fit <- CBDCS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

CBD with curvature and cohort model

Description

Fits and forecasts mortality rates using CBD with curvature and cohort model.

Usage

CBDQCS(
  x,
  M,
  curve = c("gompertz", "makeham", "perks", "weibull", "beard", "martinelle", "thatcher",
    "gompertz2", "makeham2", "perks2", "weibull2", "beard2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

The CBD with curvature and cohort (M7) model is specified as

ln(m_{x,t}) = \kappa_{1,t} + \kappa_{2,t} (x-\bar{x}) + \kappa_{3,t} ((x-\bar{x})^2-\sigma^2) + \gamma_{t-x} + \epsilon_{x,t}.

The model is estimated by Newton updating scheme and is forecasted by ARIMA applied to \kappa_{1,t}, \kappa_{2,t}, \kappa_{3,t}, and \gamma_c. Constraints include sum of \gamma_c is zero, sum of c\gamma_c is zero, and sum of c^{2}\gamma_c is zero. It is designed for ages 50-90.

Value

An object of class CBDQCS with associated S3 methods coef, forecast, plot, residuals, and simulate (nsim for setting number of simulations; seed for initialising random number generator).

References

Cairns, A.J.G., Blake, D., Dowd, K., Coughlan, G.D., Epstein, D., Ong, A., and Balevich, I. (2009). A quantitative comparison of stochastic mortality models using data from England and Wales and the United States. North American Actuarial Journal, 13(1), 1-35.

Examples

x <- 60:89
k1 <- -2.97-0.0245*(0:29)
k2 <- 0.101+0.000345*(0:29)
set.seed(123)
M <- exp(matrix(k1,nrow=30,ncol=30,byrow=FALSE)+outer(k2,(x-mean(x)))+rnorm(900,0,0.035))
fit <- CBDQCS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

CBD model

Description

Fits and forecasts mortality rates using CBD model.

Usage

CBDS(
  x,
  M,
  curve = c("gompertz", "makeham", "perks", "weibull", "beard", "martinelle", "thatcher",
    "gompertz2", "makeham2", "perks2", "weibull2", "beard2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

The CBD (M5) model is specified as

ln(m_{x,t}) = \kappa_{1,t} + \kappa_{2,t} (x-\bar{x}) + \epsilon_{x,t}.

The model is estimated by regression and is forecasted by ARIMA applied to \kappa_{1,t} and \kappa_{2,t}. It is designed for ages 50-90.

Value

An object of class CBDS with associated S3 methods coef, forecast, plot, residuals, and simulate (nsim for setting number of simulations; seed for initialising random number generator).

References

Cairns, A.J.G., Blake, D., and Dowd, K. (2006). A two-factor model for stochastic mortality with parameter uncertainty: Theory and calibration. Journal of Risk and Insurance, 73(4), 687-718.

Examples

x <- 60:89
k1 <- -2.97-0.0245*(0:29)
k2 <- 0.101+0.000345*(0:29)
set.seed(123)
M <- exp(matrix(k1,nrow=30,ncol=30,byrow=FALSE)+outer(k2,(x-mean(x)))+rnorm(900,0,0.035))
fit <- CBDS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

Mortality curve calculation

Description

Calculates mortality rates given mortality curve parameter values.

Usage

CC(
  x,
  par,
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher")
)

Arguments

Details

See MC() for more details of different mortality curves.

Value

Calculated mortality rates based on selected mortality curve and parameter values.

Examples

CC(x=60:89,par=c(0.0000082,0.10771),curve="gompertz")

Common factor model with global age pattern

Description

Fits and forecasts mortality rates of two populations using common factor model with global age pattern.

Usage

CFM2S(
  x,
  M1,
  M2,
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

The common factor model with global age pattern is specified as

ln(m_{x,t,i}) = \alpha_{x,i} + B_x K_t + \beta_x \kappa_{t,i} + \epsilon_{x,t,i}.

The model is estimated by Newton updating scheme and is forecasted by ARIMA applied to K_t and \kappa_{t,i}. Constraints include sum of B_x is one, sum of K_t is zero, sum of \beta_x is one, and sum of \kappa_{t,i} is zero. It can be applied to whole age range.

Value

An object of class CFM2S with associated S3 methods coef, forecast, plot (which = 1 gives parameter estimates; which = 2 gives residuals and forecasts), and residuals.

References

Li, J., Wang, M., Liu, J., and Leung, J.W.Y. (2026). Financial valuation of retirement village via stochastic modelling of disability prevalence rates. ASTIN Bulletin, 56(2), 447-473.

Examples

x <- 60:89
a1 <- c(-5.18,-5.12,-4.98,-4.92,-4.82,-4.73,-4.66,-4.53,-4.45,-4.35,
-4.26,-4.17,-4.05,-3.95,-3.84,-3.73,-3.65,-3.52,-3.40,-3.29,
-3.14,-3.02,-2.88,-2.76,-2.64,-2.49,-2.37,-2.25,-2.12,-2.00)
a2 <- c(-4.78,-4.68,-4.57,-4.49,-4.39,-4.29,-4.19,-4.10,-4.00,-3.89,
-3.80,-3.69,-3.60,-3.49,-3.39,-3.29,-3.17,-3.07,-2.96,-2.85,
-2.71,-2.62,-2.49,-2.37,-2.26,-2.14,-2.04,-1.91,-1.82,-1.72)
B <- c(0.0381,0.0340,0.0420,0.0389,0.0423,0.0414,0.0406,0.0393,0.0415,0.0400,
0.0411,0.0362,0.0387,0.0381,0.0384,0.0385,0.0356,0.0314,0.0317,0.0337,
0.0316,0.0298,0.0284,0.0270,0.0248,0.0262,0.0205,0.0215,0.0142,0.0145)
K <- c(9.66,9.89,10.66,9.83,9.52,7.39,7.64,6.36,2.32,4.18,
2.91,-0.61,0.28,-0.38,-1.79,-3.34,-1.74,-3.50,-4.28,-4.77,
-4.98,-7.13,-5.09,-6.41,-5.56,-5.65,-6.12,-5.64,-7.35,-6.28)
b <- c(-0.00010,0.01195,0.03030,0.02170,0.03690,0.02365,0.02280,0.03850,0.05845,0.04415,
0.04185,0.05175,0.03670,0.04195,0.04090,0.02775,0.04990,0.02865,0.03935,0.03820,
0.04000,0.02790,0.03705,0.03370,0.02940,0.02850,0.03400,0.02310,0.02675,0.03430)
k1 <- c(-1.24,-1.38,-3.48,-2.51,-1.32,-1.90,-3.42,-0.94,0.24,-0.48,
-0.26,2.70,1.39,-0.46,1.74,2.53,0.90,1.43,0.76,2.48,
0.74,2.32,0.42,1.69,-0.64,1.30,0.19,-0.69,-1.11,-1.01)
k2 <- c(2.35,0.62,-0.38,0.12,0.00,0.80,-1.39,0.38,2.47,0.40,
0.76,3.06,1.42,-0.73,0.79,1.94,0.12,0.60,-0.43,0.29,
0.17,0.98,-1.01,-0.13,-2.46,-1.24,-1.65,-2.48,-2.32,-3.06)
set.seed(123)
M1 <- exp(outer(k1,b)+outer(K,B)+matrix(a1,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.07))
M2 <- exp(outer(k2,b)+outer(K,B)+matrix(a2,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.07))
fit <- CFM2S(x=x,M1=M1,M2=M2,curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

Common factor model

Description

Fits and forecasts mortality rates of two populations using common factor model.

Usage

CFMS(
  x,
  M1,
  M2,
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

The common factor model (CFM) is specified as

ln(m_{x,t,i}) = \alpha_{x,i} + B_x K_t + \beta_{x,i} \kappa_{t,i} + \epsilon_{x,t,i}.

The model is estimated by Newton updating scheme and is forecasted by ARIMA applied to K_t and \kappa_{t,i}. Constraints include sum of B_x is one, sum of K_t is zero, sum of \beta_{x,i} is one, and sum of \kappa_{t,i} is zero. It can be applied to whole age range.

Value

An object of class CFMS with associated S3 methods coef, forecast, plot (which = 1 gives parameter estimates; which = 2 gives residuals and forecasts), and residuals.

References

Li, N. and Lee, R. (2005). Coherent mortality forecasts for a group of populations: An extension of the Lee-Carter method. Demography, 42(3), 575-594.

Examples

x <- 60:89
a1 <- c(-5.18,-5.12,-4.98,-4.92,-4.82,-4.73,-4.66,-4.53,-4.45,-4.35,
-4.26,-4.17,-4.05,-3.95,-3.84,-3.73,-3.65,-3.52,-3.40,-3.29,
-3.14,-3.02,-2.88,-2.76,-2.64,-2.49,-2.37,-2.25,-2.12,-2.00)
a2 <- c(-4.78,-4.68,-4.57,-4.49,-4.39,-4.29,-4.19,-4.10,-4.00,-3.89,
-3.80,-3.69,-3.60,-3.49,-3.39,-3.29,-3.17,-3.07,-2.96,-2.85,
-2.71,-2.62,-2.49,-2.37,-2.26,-2.14,-2.04,-1.91,-1.82,-1.72)
B <- c(0.0381,0.0340,0.0420,0.0389,0.0423,0.0414,0.0406,0.0393,0.0415,0.0400,
0.0411,0.0362,0.0387,0.0381,0.0384,0.0385,0.0356,0.0314,0.0317,0.0337,
0.0316,0.0298,0.0284,0.0270,0.0248,0.0262,0.0205,0.0215,0.0142,0.0145)
K <- c(9.66,9.89,10.66,9.83,9.52,7.39,7.64,6.36,2.32,4.18,
2.91,-0.61,0.28,-0.38,-1.79,-3.34,-1.74,-3.50,-4.28,-4.77,
-4.98,-7.13,-5.09,-6.41,-5.56,-5.65,-6.12,-5.64,-7.35,-6.28)
b1 <- c(0.0012,-0.0033,0.0523,0.0161,0.0529,0.0220,0.0312,0.0437,0.0709,0.0444,
0.0398,0.0361,0.0403,0.0396,0.0506,0.0315,0.0428,0.0261,0.0384,0.0388,
0.0300,0.0269,0.0275,0.0256,0.0239,0.0421,0.0314,0.0284,0.0174,0.0314)
k1 <- c(-1.24,-1.38,-3.48,-2.51,-1.32,-1.90,-3.42,-0.94,0.24,-0.48,
-0.26,2.70,1.39,-0.46,1.74,2.53,0.90,1.43,0.76,2.48,
0.74,2.32,0.42,1.69,-0.64,1.30,0.19,-0.69,-1.11,-1.01)
b2 <- c(-0.0014,0.0272,0.0083,0.0273,0.0209,0.0253,0.0144,0.0333,0.0460,0.0439,
0.0439,0.0674,0.0331,0.0443,0.0312,0.0240,0.0570,0.0312,0.0403,0.0376,
0.0500,0.0289,0.0466,0.0418,0.0349,0.0149,0.0366,0.0178,0.0361,0.0372)
k2 <- c(2.35,0.62,-0.38,0.12,0.00,0.80,-1.39,0.38,2.47,0.40,
0.76,3.06,1.42,-0.73,0.79,1.94,0.12,0.60,-0.43,0.29,
0.17,0.98,-1.01,-0.13,-2.46,-1.24,-1.65,-2.48,-2.32,-3.06)
set.seed(123)
M1 <- exp(outer(k1,b1)+outer(K,B)+matrix(a1,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.07))
M2 <- exp(outer(k2,b2)+outer(K,B)+matrix(a2,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.07))
fit <- CFMS(x=x,M1=M1,M2=M2,curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

Mortality ensemble interval

Description

Generates ensemble interval forecast of mortality rates.

Usage

ENI(..., width = 0.95, method = 1, wm, nsim = 10, seed = 123)

Arguments

Details

Ensemble interval forecast is constructed by combining interval forecasts from individual stochastic mortality models using different methods including simple averaging, weighted averaging, envelope, and interior trimming. See LCS(), RHS(), APCS(), CBDS(), CBDCS(), CBDQCS(), and STARS() for more details of different stochastic mortality models.

Value

An object of class ENI with associated S3 methods forecast and plot.

References

Li, J., Wang, M., Liu, J., and Tickle, L. (2025). Ensemble interval forecasts of mortality. Scandinavian Actuarial Journal, 2025(6), 598-616.

Examples

x <- 60:69
a <- c(-4.8499,-4.7676,-4.6719,-4.5722,-4.4847,-4.3841,-4.2813,-4.1863,-4.0861,-3.9962)
b <- c(0.0801,0.0909,0.0948,0.0951,0.0965,0.1014,0.1042,0.1141,0.1110,0.1118)
k <- c(12.11,10.69,11.18,9.64,9.35,8.21,6.89,5.74,4.56,3.60,
3.27,2.04,1.11,-0.44,-1.05,-1.03,-1.84,-2.90,-4.03,-4.12,
-5.18,-5.64,-6.00,-6.51,-6.91,-6.90,-8.32,-8.53,-9.69,-9.31)
set.seed(123)
M <- exp(outer(k,b)+matrix(a,nrow=30,ncol=10,byrow=TRUE)+rnorm(300,0,0.035))
fit1 <- LCS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
fit2 <- RHS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
fit3 <- APCS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
fit4 <- CBDS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
fit <- ENI(fit1,fit2,fit3,fit4)
forecast::forecast(fit)
plot(fit)

Mortality ensemble

Description

Fits and forecasts mortality rates using mortality ensemble.

Usage

ENS(
  x,
  M,
  wm = rep(1/7, 7),
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

Ensemble forecast is obtained as a weighted average of forecasts from individual stochastic mortality models. See LCS(), RHS(), APCS(), CBDS(), CBDCS(), CBDQCS(), and STARS() for more details of different stochastic mortality models.

Value

An object of class ENS with associated S3 methods forecast and plot.

References

Li, J. (2023). A model stacking approach for forecasting mortality. North American Actuarial Journal, 27(3), 530-545.

Examples

x <- 60:69
a <- c(-4.8499,-4.7676,-4.6719,-4.5722,-4.4847,-4.3841,-4.2813,-4.1863,-4.0861,-3.9962)
b <- c(0.0801,0.0909,0.0948,0.0951,0.0965,0.1014,0.1042,0.1141,0.1110,0.1118)
k <- c(12.11,10.69,11.18,9.64,9.35,8.21,6.89,5.74,4.56,3.60,
3.27,2.04,1.11,-0.44,-1.05,-1.03,-1.84,-2.90,-4.03,-4.12,
-5.18,-5.64,-6.00,-6.51,-6.91,-6.90,-8.32,-8.53,-9.69,-9.31)
set.seed(123)
M <- exp(outer(k,b)+matrix(a,nrow=30,ncol=10,byrow=TRUE)+rnorm(300,0,0.035))
fit <- ENS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
forecast::forecast(fit)
plot(fit)

Mortality model fitting

Description

Fits and forecasts mortality rates using different stochastic mortality models with different estimation methods.

Usage

FCS(
  x,
  M,
  model = c("LC", "RH", "APC", "CBD", "CBDC", "CBDQC", "STAR"),
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

See LCS(), RHS(), APCS(), CBDS(), CBDCS(), CBDQCS(), and STARS() for more details of different stochastic mortality models.

Value

An object of class based on selected stochastic mortality model with associated S3 methods coef, forecast, plot, and residuals.

Examples

x <- 60:89
a <- c(-4.8499,-4.7676,-4.6719,-4.5722,-4.4847,-4.3841,-4.2813,-4.1863,-4.0861,-3.9962,
-3.8885,-3.7896,-3.6853,-3.5737,-3.4728,-3.3718,-3.2586,-3.1474,-3.0371,-2.9206,
-2.7998,-2.6845,-2.5653,-2.4581,-2.3367,-2.2159,-2.1017,-1.9941,-1.8821, -1.7697)
b <- c(0.0283,0.0321,0.0335,0.0336,0.0341,0.0358,0.0368,0.0403,0.0392,0.0395,
0.0396,0.0399,0.0397,0.0386,0.039,0.0375,0.0367,0.0368,0.035,0.0354,
0.0336,0.0323,0.0313,0.0295,0.0282,0.0265,0.024,0.0226,0.0219,0.0183)
k <- c(12.11,10.69,11.18,9.64,9.35,8.21,6.89,5.74,4.56,3.6,
3.27,2.04,1.11,-0.44,-1.05,-1.03,-1.84,-2.9,-4.03,-4.12,
-5.18,-5.64,-6,-6.51,-6.91,-6.9,-8.32,-8.53,-9.69,-9.31)
set.seed(123)
M <- exp(outer(k,b)+matrix(a,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.035))
fit <- FCS(x=x,M=M,model="LC",curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

Lee-Carter model

Description

Fits and forecasts mortality rates using Lee-Carter model.

Usage

LCS(
  x,
  M,
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

The Lee-Carter (LC) model is specified as

ln(m_{x,t}) = \alpha_x + \beta_x \kappa_t + \epsilon_{x,t}.

The model is estimated by singular value decomposition and is forecasted by ARIMA applied to \kappa_t. Constraints include sum of \beta_x is one and sum of \kappa_t is zero. It can be applied to whole age range.

Value

An object of class LCS with associated S3 methods coef, forecast, plot, residuals, and simulate (nsim for setting number of simulations; seed for initialising random number generator).

References

Lee, R.D. and Carter, L.R. (1992). Modeling and forecasting U.S. mortality. Journal of the American Statistical Association, 87(419), 659-671.

Examples

x <- 60:89
a <- c(-4.8499,-4.7676,-4.6719,-4.5722,-4.4847,-4.3841,-4.2813,-4.1863,-4.0861,-3.9962,
-3.8885,-3.7896,-3.6853,-3.5737,-3.4728,-3.3718,-3.2586,-3.1474,-3.0371,-2.9206,
-2.7998,-2.6845,-2.5653,-2.4581,-2.3367,-2.2159,-2.1017,-1.9941,-1.8821, -1.7697)
b <- c(0.0283,0.0321,0.0335,0.0336,0.0341,0.0358,0.0368,0.0403,0.0392,0.0395,
0.0396,0.0399,0.0397,0.0386,0.039,0.0375,0.0367,0.0368,0.035,0.0354,
0.0336,0.0323,0.0313,0.0295,0.0282,0.0265,0.024,0.0226,0.0219,0.0183)
k <- c(12.11,10.69,11.18,9.64,9.35,8.21,6.89,5.74,4.56,3.6,
3.27,2.04,1.11,-0.44,-1.05,-1.03,-1.84,-2.9,-4.03,-4.12,
-5.18,-5.64,-6,-6.51,-6.91,-6.9,-8.32,-8.53,-9.69,-9.31)
set.seed(123)
M <- exp(outer(k,b)+matrix(a,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.035))
fit <- LCS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

Mortality curve fitting

Description

Fits parametric mortality curves to observed mortality rates using multiple optimisation strategies.

Usage

MC(
  x,
  m,
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  w = rep(1, length(x))
)

Arguments

Details

User can choose one of the following 14 mortality curves (with their suitable age ranges in brackets):

• Gompertz (1825) (ages 30-90) - m_x = B e^{Cx}

• Makeham (1860) (ages 20-90) - m_x = A + B e^{Cx}

• Oppermann (1870) (ages 0-20) - m_x = \frac{A}{\sqrt{x+1}} + B + C\sqrt{x+1}

• Thiele (1871) (ages 0-90) - m_x = A_1 e^{-B_1 x} + A_2 e^{-0.5 B_2 (x-C)^2} + A_3 e^{B_3 x}

• Wittstein & Bumsted (1883) (ages 0-90) - m_x = \frac{A^{-(Bx)^N}}{B} + A^{-(M-x)^N}

• Perks (1932) (ages 20-100+) - m_x = \frac{A+BC^x}{1+DC^x}

• Weibull (1939) (ages 50-100+) - m_x = Bx^C

• Van der Maen (1943) (ages 80-100+) - m_x = A + Bx + Cx^2 + \frac{I}{N-x}

• Beard (1971) (ages 50-100+) - m_x = \frac{A e^{Cx}}{1+B e^{Cx}}

• Heligman & Pollard (1980) (ages 0-100+) - m_x = A^{(x+B)^C} + D e^{-E(ln(x)-ln(F))^2} + \frac{GH^x}{1+GH^x}

• Rogers & Planck (1983) (ages 0-100+) - m_x = A_0 + A_1 e^{-Ax} + A_2 e^{-B(x-U)-e^{-C(x-U)}} + A_3 e^{Dx}

• Siler (1983) (ages 0-90) - m_x = A_1 e^{-B_1 x} + A_2 + A_3 e^{B_3 x}

• Martinelle (1987) (ages 20-100+) - m_x = \frac{A+B e^{Cx}}{1+D e^{Cx}} + E e^{Cx}

• Thatcher (1999) (ages 20-100+) - m_x = A + \frac{B e^{Cx}}{1+D e^{Cx}}

The mortality curves are fitted by employing multiple optimisation strategies simultaneously (including PORT routines, Nelder-Mead method, and Levenberg-Marquardt algorithm) and selecting one with smallest weighted least squares on ln(m_x). Constrained parameterisations offer traditional interpretability while unconstrained parameterisations provide increased flexibility.

Value

An object of class based on selected mortality curve with associated S3 methods coef, fitted, predict, plot, deviance, and residuals.

References

Gompertz, B. (1825). On the nature of the function expressive of the law of human mortality, and on a new mode of determining the value of life contingencies. Philosophical Transactions of the Royal Society of London, 115(1825), 513-583.

Makeham, W.M. (1860). On the law of mortality and the construction of annuity tables. Journal of the Institute of Actuaries, 8(6), 301-310.

Oppermann, L.H.F. (1870). On the graduation of life tables, with special application to the rate of mortality in infancy and childhood. The Institute of Actuaries.

Thiele, T.N. (1871). On a mathematical formula to express the rate of mortality throughout the whole of life, tested by a series of observations made use of by the Danish Life Insurance Company of 1871. Journal of the Institute of Actuaries and Assurance Magazine, 16(5), 313-329.

Wittstein, T. and Bumsted, D.A. (1883). The mathematical law of mortality. Journal of the Institute of Actuaries and Assurance Magazine, 24(3), 153-173.

Perks, W. (1932). On some experiments in the graduation of mortality statistics. Journal of the Institute of Actuaries, 63(1), 12-57.

Weibull, W. (1951). A statistical distribution function of wide applicability. Journal of Applied Mechanics, 18(3), 293-297.

Beard, R.E. (1971). Some aspects of theories of mortality, cause of death analysis, forecasting and stochastic processes. Biological Aspects of Demography, 57-68.

Heligman, L. and Pollard, J.H. (1980). The age pattern of mortality. Journal of the Institute of Actuaries, 107(1), 49-80.

Rogers, A. and Planck, F. (1983). Model: A general program for estimating parametrized model schedules of fertility, mortality, migration, and marital and labor force status transitions. IIASA Working Paper WP-83-102.

Siler, W. (1983). Parameters of mortality in human populations with widely varying life spans. Statistics in Medicine, 2(3), 373-380.

Martinelle, S. (1987). A generalized Perks formula for old-age mortality. Statistiska Centralbyran.

Thatcher, A.R. (1999). The long-term pattern of adult mortality and the highest attained age. Journal of the Royal Statistical Society Series A, 162(1), 5-43.

Tabeau, E. (2001). A review of demographic forecasting models for mortality. Forecasting Mortality in Developed Countries. European Studies of Population, Vol 9.

Examples

x <- 60:89
set.seed(123); m <- 0.0000082*exp(0.10771*c(60:89)+rnorm(30,0,0.1))
fit <- MC(x=x,m=m,curve="gompertz")
coef(fit)
fitted(fit)
predict(fit,62.5)
plot(fit)
deviance(fit)
residuals(fit)

Renshaw-Haberman model

Description

Fits and forecasts mortality rates using Renshaw-Haberman model.

Usage

RHS(
  x,
  M,
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  h = 10,
  jumpoff = 1
)

Arguments

Details

The Renshaw-Haberman (RH) model is specified as

ln(m_{x,t}) = \alpha_x + \beta_x \kappa_t + \gamma_{t-x} + \epsilon_{x,t}.

The model is estimated by Newton updating scheme and is forecasted by ARIMA applied to \kappa_t and \gamma_c. Constraints include sum of \beta_x is one, sum of \kappa_t is zero, and sum of \gamma_c is zero. It can be applied to whole age range.

Value

An object of class RHS with associated S3 methods coef, forecast, plot, residuals, and simulate (nsim for setting number of simulations; seed for initialising random number generator).

References

Haberman, S. and Renshaw, A. (2011). A comparative study of parametric mortality projection models. Insurance: Mathematics and Economics, 48(1), 35-55.

Examples

x <- 60:89
a <- c(-4.8499,-4.7676,-4.6719,-4.5722,-4.4847,-4.3841,-4.2813,-4.1863,-4.0861,-3.9962,
-3.8885,-3.7896,-3.6853,-3.5737,-3.4728,-3.3718,-3.2586,-3.1474,-3.0371,-2.9206,
-2.7998,-2.6845,-2.5653,-2.4581,-2.3367,-2.2159,-2.1017,-1.9941,-1.8821, -1.7697)
b <- c(0.0283,0.0321,0.0335,0.0336,0.0341,0.0358,0.0368,0.0403,0.0392,0.0395,
0.0396,0.0399,0.0397,0.0386,0.039,0.0375,0.0367,0.0368,0.035,0.0354,
0.0336,0.0323,0.0313,0.0295,0.0282,0.0265,0.024,0.0226,0.0219,0.0183)
k <- c(12.11,10.69,11.18,9.64,9.35,8.21,6.89,5.74,4.56,3.6,
3.27,2.04,1.11,-0.44,-1.05,-1.03,-1.84,-2.9,-4.03,-4.12,
-5.18,-5.64,-6,-6.51,-6.91,-6.9,-8.32,-8.53,-9.69,-9.31)
set.seed(123)
M <- exp(outer(k,b)+matrix(a,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.035))
fit <- RHS(x=x,M=M,curve="makeham",h=30,jumpoff=2)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

Spatial-temporal autoregressive model

Description

Fits and forecasts mortality rates using spatial-temporal autoregressive model.

Usage

STARS(
  x,
  M,
  curve = c("gompertz", "makeham", "oppermann", "thiele", "wittsteinbumsted", "perks",
    "weibull", "vandermaen", "beard", "heligmanpollard", "rogersplanck", "siler",
    "martinelle", "thatcher", "gompertz2", "makeham2", "oppermann2", "thiele2",
    "wittsteinbumsted2", "perks2", "weibull2", "vandermaen2", "beard2",
    "heligmanpollard2", "rogersplanck2", "siler2", "martinelle2", "thatcher2"),
  h = 10
)

Arguments

Details

The spatial-temporal autoregressive (STAR) model is specified as

ln(M_t) = \mu + R ln(M_{t-1}) + E_t,

where \mu contains intercept parameters and R is banded matrix with non-zero entries in main diagonal and two subdiagonals below that.

The model is estimated by constrained regression and is forecasted via its autoregressive nature. Constraints in R include: all non-zero entries are between zero and one and sum of each row is one. It can be applied to whole age range.

Value

An object of class STARS with associated S3 methods coef, forecast, plot, residuals, and simulate (nsim for setting number of simulations; seed for initialising random number generator).

References

Li, H. and Lu, Y. (2017). Coherent forecasting of mortality rates: A sparse vector-autoregression approach. ASTIN Bulletin, 47(2), 563-600.

Examples

x <- 60:89
a <- c(-4.8499,-4.7676,-4.6719,-4.5722,-4.4847,-4.3841,-4.2813,-4.1863,-4.0861,-3.9962,
-3.8885,-3.7896,-3.6853,-3.5737,-3.4728,-3.3718,-3.2586,-3.1474,-3.0371,-2.9206,
-2.7998,-2.6845,-2.5653,-2.4581,-2.3367,-2.2159,-2.1017,-1.9941,-1.8821, -1.7697)
b <- c(0.0283,0.0321,0.0335,0.0336,0.0341,0.0358,0.0368,0.0403,0.0392,0.0395,
0.0396,0.0399,0.0397,0.0386,0.039,0.0375,0.0367,0.0368,0.035,0.0354,
0.0336,0.0323,0.0313,0.0295,0.0282,0.0265,0.024,0.0226,0.0219,0.0183)
k <- c(12.11,10.69,11.18,9.64,9.35,8.21,6.89,5.74,4.56,3.6,
3.27,2.04,1.11,-0.44,-1.05,-1.03,-1.84,-2.9,-4.03,-4.12,
-5.18,-5.64,-6,-6.51,-6.91,-6.9,-8.32,-8.53,-9.69,-9.31)
set.seed(123)
M <- exp(outer(k,b)+matrix(a,nrow=30,ncol=30,byrow=TRUE)+rnorm(900,0,0.035))
fit <- STARS(x=x,M=M,curve="makeham",h=30)
coef(fit)
forecast::forecast(fit)
plot(fit)
residuals(fit)

Common S3 methods for mortality curves

Description

Common S3 methods for objects returned by MC():

coef(object) gives estimated parameter values.

fitted(object) gives fitted values.

predict(object,x) gives estimated mortality rate(s) for given age(s) x.

plot(object) plots data and fitted mortality curve.

deviance(object) gives weighted sum of squared residuals on ln(m_x).

residuals(object) gives residuals on ln(m_x).

Common S3 methods for stochastic mortality models

Description

Common S3 methods for objects returned by LCS(), RHS(), APCS(), CBDS(), CBDCS(), CBDQCS(), STARS(), ENS(), ENI(), CFMS(), CFM2S(), and CAES():

coef(object) gives estimated parameter values.

forecast::forecast(object) gives mortality forecasts smoothed by selected mortality curve (or ensemble interval forecasts for ENI()).

plot(object) plots estimated parameter values, standardised residuals heatmap, data, and mortality forecasts smoothed by selected mortality curve (or ensemble interval forecasts for ENI()).

residuals(object) gives standardised residuals on ln(m_{x,t}).

simulate(object,nsim,seed) gives simulated future mortality rates (unsmoothed).