simulate.fitStMoMo: Simulate future sample paths from a Stochastic Mortality Model

Description

Simulate future sample paths from a Stochastic Mortality Model. The period indexes $\kappa_t^{(i)}, i = 1,..N,$ are modelled using ether a Multivariate Random Walk with Drift (MRWD) or $N$ independent ARIMA$(p, d, q)$ models. The cohort index $\gamma_{t-x}$ is modelled using an ARIMA$(p, d, q)$. By default an ARIMA$(1, 1, 0)$ with a constant is used.

Usage

# S3 method for fitStMoMo
simulate(object, nsim = 1000, seed = NULL, h = 50,
  oxt = NULL, gc.order = c(1, 1, 0), gc.include.constant = TRUE,
  jumpchoice = c("fit", "actual"), kt.method = c("mrwd", "iarima"),
  kt.order = NULL, kt.include.constant = TRUE, kt.lookback = NULL,
  gc.lookback = NULL, ...)

Arguments

object

an object of class "fitStMoMo" with the fitted parameters of a stochastic mortality model.

nsim

number of sample paths to simulate.

seed

either NULL or an integer that will be used in a call to set.seed before simulating the time series. The default, NULL will not change the random generator state.

number of years ahead to forecast.

oxt

optional matrix/vector or scalar of known offset to be added in the simulations. This can be used to specify any a priori known component to be added to the simulated predictor.

gc.order

a specification of the ARIMA model for the cohort effect: the three components $(p, d, q)$ are the AR order, the degree of differencing, and the MA. The default is an ARIMA$(1, 1, 0)$.

gc.include.constant

a logical value indicating if the ARIMA model should include a constant value. The default is TRUE.

jumpchoice

option to select the jump-off rates, i.e. the rates from the final year of observation, to use in projections of mortality rates. "fit"(default) uses the fitted rates and "actual" uses the actual rates from the final year.

kt.method

optional forecasting method for the period index. The alternatives are "mrwd"(default) and "iarima". See details.

kt.order

an optional matrix with one row per period index specifying the ARIMA models: for the ith row (ith period index) the three components $(p, d, q)$ are the AR order, the degree of differencing, and the MA order. If absent the arima models are fitted using auto.arima. This argument is only used when kt.method is "iarima".

kt.include.constant

an optional vector of logical values indicating if the ARIMA model for the ith period index should include a constant value. The default is TRUE. This argument is only used when kt.method is "iarima".

kt.lookback

optional argument to specify the look-back window to use in the estimation of the time series model for the period indexes. By default all the estimated values are used. If kt.lookback is provided then the last kt.lookback years of $\kappa_t^{(i)}, i = 1,..N,$ are used.

gc.lookback

optional argument to specify the look-back window to use in the estimation of the ARIMA model for the cohort effect. By default all the estimated values are used in estimating the ARIMA model. If gc.lookback is provided then the last gc.lookback years of $\gamma_{t-x}$ are used.

...

other arguments.

Value

A list of class "simStMoMo" with components:

rates

a three dimensional array with the future simulated rates.

ages

vector of ages corresponding to the first dimension of rates.

years

vector of years for which a simulations has been produced. This corresponds to the second dimension of rates.

kt.s

information on the simulated paths of the period indexes of the model. This is a list with the model fitted to $\kappa_t$; the simulated paths (sim); and the years for which simulations were produced. If the mortality model does not have any age-period terms (i.e. $N=0$) this is set to NULL.

gc.s

information on the simulated paths of the cohort index of the model. This is a list with the model fitted to $\gamma_c$; the simulated paths (sim); and the cohorts for which simulations were produced. If the mortality model does not have a cohort effect this is set to NULL.

oxt.s

a three dimensional array with the offset used in the simulations.

fitted

a three dimensional array with the in-sample rates of the model for the years for which the mortality model was fitted.

jumpchoice

Jump-off method used in the simulation.

kt.method

method used in the modelling of the period index.

model

the model fit from which the simulations were produced.

Details

If kt.method is "mrwd", fitting and simulation of the time series model for the period indexes is done with a Multivariate Random Walk with Drift using the function mrwd.

If kt.method is "iarima", fitting and simulation of the time series model for the period indexes is done with $N$ independent arima models using the function iarima. See this latter function for details on input arguments kt.order and kt.include.constant.

Fitting and simulation of the ARIMA model for the cohort index is done with function Arima from package forecast. See the latter function for further details on input arguments gc.order and gc.include.constant.

Note that in some cases simulations of the cohort effects may be needed for a horizon longer than h. This is the case when in the fitted model the most recent cohorts have been zero weighted. The simulated cohorts can be seen in gc.s$cohorts.

Examples

Run this code

# NOT RUN {
#Lee-Carter
LCfit <- fit(lc(), data = EWMaleData, ages.fit = 55:89)
LCsim.mrwd <- simulate(LCfit, nsim = 100)
LCsim.iarima <- simulate(LCfit, nsim = 100, kt.method = "iarima", 
                         kt.order = c(1, 1, 2))

par(mfrow=c(2, 2))
plot(LCfit$years, LCfit$kt[1, ], xlim = range(LCfit$years, LCsim.mrwd$kt.s$years),
     ylim = range(LCfit$kt, LCsim.mrwd$kt.s$sim), type = "l", 
     xlab = "year", ylab = "kt", 
     main = "Lee-Carter: Simulated paths of the period index kt (mrwd)")
matlines(LCsim.mrwd$kt.s$years, LCsim.mrwd$kt.s$sim[1, , ], type = "l", lty = 1)

plot(LCfit$years, (LCfit$Dxt / LCfit$Ext)["65", ], 
     xlim = range(LCfit$years, LCsim.mrwd$years),
     ylim = range((LCfit$Dxt / LCfit$Ext)["65", ], LCsim.mrwd$rates["65", , ]), 
     type = "l", xlab = "year", ylab = "rate", 
     main = "Lee-Carter: Simulated mortality rates at age 65")
matlines(LCsim.mrwd$years, LCsim.mrwd$rates["65", , ], type = "l", lty = 1)

plot(LCfit$years, LCfit$kt[1, ], xlim = range(LCfit$years, LCsim.iarima$kt.s$years),
     ylim = range(LCfit$kt, LCsim.iarima$kt.s$sim), type = "l", 
     xlab = "year", ylab = "kt", 
     main = "Lee-Carter: Simulated paths of the period index kt (ARIMA(1, 1, 2))")
matlines(LCsim.iarima$kt.s$years, LCsim.iarima$kt.s$sim[1, , ], type = "l", lty = 1)

plot(LCfit$years, (LCfit$Dxt / LCfit$Ext)["65", ], 
     xlim = range(LCfit$years, LCsim.iarima$years),
     ylim = range((LCfit$Dxt / LCfit$Ext)["65", ], LCsim.iarima$rates["65", , ]), 
     type = "l", xlab = "year", ylab = "rate", 
     main = "Lee-Carter: Simulated mortality rates at age 65 (ARIMA(1, 1, 2))")
matlines(LCsim.iarima$years, LCsim.iarima$rates["65", , ], type = "l", lty = 1)

#APC
par(mfrow=c(1, 3))
wxt <- genWeightMat(55:89,  EWMaleData$years, clip = 3)
APCfit <- fit(apc(), data = EWMaleData, ages.fit = 55:89, wxt = wxt)
APCsim <- simulate(APCfit, nsim = 100, gc.order = c(1, 1, 0))

plot(APCfit$years, APCfit$kt[1, ], 
     xlim = range(APCfit$years, APCsim$kt.s$years),
     ylim = range(APCfit$kt, APCsim$kt.s$sim), type = "l",
     xlab = "year", ylab = "kt",
     main = "APC: Simulated paths of the period index kt")
matlines(APCsim$kt.s$years, APCsim$kt.s$sim[1, , ], type = "l", lty = 1)

plot(APCfit$cohorts, APCfit$gc, 
     xlim = range(APCfit$cohorts, APCsim$gc.s$cohorts),
     ylim = range(APCfit$gc, APCsim$gc.s$sim, na.rm = TRUE), type = "l",
     xlab = "year", ylab = "kt", 
     main = "APC: Simulated paths of the cohort index (ARIMA(1,1,0))")
matlines(APCsim$gc.s$cohorts, APCsim$gc.s$sim, type = "l", lty = 1)

plot(APCfit$years, (APCfit$Dxt / APCfit$Ext)["65", ], 
     xlim = range(APCfit$years, APCsim$years),
     ylim = range((APCfit$Dxt/APCfit$Ext)["65", ], APCsim$rates["65", , ]), 
     type = "l", xlab = "year", ylab = "rate", 
     main = "APC: Simulated of mortality rates at age 65")
matlines(APCsim$years, APCsim$rates["65", , ], type = "l", lty = 1)

#Compare LC and APC
library(fanplot)
par(mfrow=c(1, 1))
plot(LCfit$years, (LCfit$Dxt / LCfit$Ext)["65", ], 
     xlim = range(LCfit$years, LCsim.mrwd$years),
     ylim = range((LCfit$Dxt / LCfit$Ext)["65", ], LCsim.mrwd$rates["65", , ], 
     APCsim$rates["65", , ]), type = "l", xlab = "year", ylab = "rate", 
     main = "Fan chart of mortality rates at age 65 (LC vs. APC)")
fan(t(LCsim.mrwd$rates["65", , ]), start = LCsim.mrwd$years[1], 
    probs = c(2.5, 10, 25, 50, 75, 90, 97.5), n.fan = 4,
    fan.col = colorRampPalette(c(rgb(1, 0, 0), rgb(1, 1, 1))), ln = NULL)
fan(t(APCsim$rates["65", 1:(length(APCsim$years) - 3), ]), 
    start = APCsim$years[1], probs = c(2.5, 10, 25, 50, 75, 90, 97.5), 
    n.fan = 4, fan.col = colorRampPalette(c(rgb(0, 0, 1), rgb(1, 1, 1))), 
    ln = NULL) 
# }