fit_M1U: A function to fit the stochastic mortality model M1U.

View source: R/fit_M1U.R

fit_M1UR Documentation

A function to fit the stochastic mortality model M1U.

Description

Carry out Bayesian estimation of the stochastic mortality model M1U.

Usage

fit_M1U(
  death,
  expo,
  n_iter = 10000,
  family = "nb",
  n.chain = 1,
  thin = 1,
  n.adapt = 1000,
  forecast = FALSE,
  h = 5,
  quiet = FALSE
)

Arguments

death

death data that has been formatted through the function preparedata_fn.

expo

expo data that has been formatted through the function preparedata_fn.

n_iter

number of iterations to run. Default is n_iter=10000.

family

a string of characters that defines the family function associated with the mortality model. "poisson" would assume that deaths follow a Poisson distribution and use a log link; "binomial" would assume that deaths follow a Binomial distribution and a logit link; "nb" (default) would assume that deaths follow a Negative-Binomial distribution and a log link.

n.chain

number of parallel chains for the model.

thin

thinning interval for monitoring purposes.

n.adapt

the number of iterations for adaptation. See ?rjags::adapt for details.

forecast

a logical value indicating if forecast is to be performed (default is FALSE). See below for details.

h

a numeric value giving the number of years to forecast. Default is h=5.

quiet

if TRUE then messages generated during compilation will be suppressed, as well as the progress bar during adaptation.

Details

The model can be described mathematically as follows: If family="poisson", then

d_{x,t,p} \sim \text{Poisson}(E^c_{x,t,p} m_{x,t,p}) ,

\log(m_{x,t,p})=a_{x,p}+b_xk_t ,

where d_{x,t,p} represents the number of deaths at age x in year t of stratum p, while E^c_{x,t,p} and m_{x,t,p} represents respectively the corresponding central exposed to risk and central mortality rate at age x in year t of stratum p. Similarly, if family="nb", then a negative binomial distribution is fitted, i.e.

d_{x,t,p} \sim \text{Negative-Binomial}(\phi,\frac{\phi}{\phi+E^c_{x,t,p} m_{x,t,p}}) ,

\log(m_{x,t,p})=a_{x,p}+b_xk_t ,

where \phi is the overdispersion parameter. See Wong et al. (2018). But if family="binomial", then

d_{x,t,p} \sim \text{Binomial}(E^0_{x,t,p} , q_{x,t,p}) ,

\text{logit}(q_{x,t,p})=a_{x,p}+b_xk_t ,

where q_{x,t,p} represents the initial mortality rate at age x in year t of stratum p, while E^0_{x,t,p}\approx E^c_{x,t,p}+\frac{1}{2}d_{x,t,p} is the corresponding initial exposed to risk. Constraints used are:

\sum_x b_x = 1, \sum_t k_t = 0 .

If forecast=TRUE, then the following time series models are fitted on k_{t,p} as follows:

k_{t,p} = \eta_1+\eta_2 t +\rho (k_{t-1,p}-(\eta_1+\eta_2 (t-1))) + \epsilon_{t,p} \text{ for }p=1,\ldots,P \text{ and } t=1,\ldots,T,

where \epsilon_{t,p}\sim N(0,\sigma_k^2) for t=1,\ldots,T, while \eta_1,\eta_2,\rho,\sigma_k^2 are additional parameters to be estimated. Note that the forecasting models are inspired by Wong et al. (2023).

Value

A list with components:

post_sample

An mcmc.list object containing the posterior samples generated.

param

A vector of character strings describing the names of model parameters.

death

The death data that was used.

expo

The expo data that was used.

family

The family function used.

forecast

A logical value indicating if forecast has been performed.

h

The forecast horizon used.

References

Jackie S. T. Wong, Jonathan J. Forster, and Peter W. F. Smith. (2018). Bayesian mortality forecasting with overdispersion, Insurance: Mathematics and Economics, Volume 2018, Issue 3, 206-221. \Sexpr[results=rd]{tools:::Rd_expr_doi("https://doi.org/10.1016/j.insmatheco.2017.09.023")}

Jackie S. T. Wong, Jonathan J. Forster, and Peter W. F. Smith. (2023). Bayesian model comparison for mortality forecasting, Journal of the Royal Statistical Society Series C: Applied Statistics, Volume 72, Issue 3, 566–586. \Sexpr[results=rd]{tools:::Rd_expr_doi("https://doi.org/10.1093/jrsssc/qlad021")}

Examples

#load and prepare mortality data
data("dxt_array_product");data("Ext_array_product")
death<-preparedata_fn(dxt_array_product,strat_name = c("ACI","DB","SCI"),ages=35:65)
expo<-preparedata_fn(Ext_array_product,strat_name = c("ACI","DB","SCI"),ages=35:65)

#fit the model (negative-binomial family)
#NOTE: This is a toy example, please run it longer in practice.
fit_M1U_result<-fit_M1U(death=death,expo=expo,n_iter=50,n.adapt=50)
head(fit_M1U_result)

BayesMoFo documentation built on Aug. 11, 2025, 1:07 a.m.