fit_PLAT {BayesMoFo} | R Documentation |
A function to fit the stochastic mortality model "PLAT".
Description
Carry out Bayesian estimation of the stochastic mortality model referred to as the "PLAT" model as proposed by Plat (2009).
Usage
fit_PLAT(
death,
expo,
share_alpha = FALSE,
share_gamma = FALSE,
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 |
expo |
expo data that has been formatted through the function |
share_alpha |
a logical value indicating if |
share_gamma |
a logical value indicating if the cohort parameter |
n_iter |
number of iterations to run. Default is |
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 |
forecast |
a logical value indicating if forecast is to be performed (default is |
h |
a numeric value giving the number of years to forecast. Default is |
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="log"
, 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}+k^1_{t,p} + k^2_{t,p}(\bar{x}-x) + k^3_{t,p}(\bar{x}-x)^+ +\gamma_{c,p} ,
where c=t-x
is the cohort index, \bar{x}
is the mean of x
, (\bar{x}-x)^+ = \text{max}(\bar{x}-x,0)
,
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}+k^1_{t,p} + k^2_{t,p}(\bar{x}-x) + k^3_{t,p}(\bar{x}-x)^+ +\gamma_{c,p} ,
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}+k^1_{t,p} + k^2_{t,p}(\bar{x}-x) + k^3_{t,p}(\bar{x}-x)^+ +\gamma_{c,p} ,
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_{t} k^1_{t,p} = \sum_{t} k^2_{t,p} = \sum_{t} k^3_{t,p} = 0, \sum_{c} \gamma_{c,p} = \sum_{c} c\gamma_{c,p} = \sum_{c} c^2\gamma_{c,p} = 0 \text{\ \ \ for each stratum }p.
If share_alpha=TRUE
, then the additive age-specific parameter is the same across all strata p
, i. e.
a_{x}+k^1_{t,p} + k^2_{t,p}(\bar{x}-x) + k^3_{t,p}(\bar{x}-x)^+ +\gamma_{c,p} .
If share_gamma=TRUE
, then the cohort parameter is the same across all strata p
, i. e.
a_{x,p}+k^1_{t,p} + k^2_{t,p}(\bar{x}-x) + k^3_{t,p}(\bar{x}-x)^+ +\gamma_{c} .
If forecast=TRUE
, then the following time series models are fitted on k^1_{t,p}
as follows:
k^1_{t,p} = \eta^1_1+\eta^1_2 t +\rho (k^1_{t-1,p}-(\eta^1_1+\eta^1_2 (t-1))) + \epsilon^1_{t,p} \text{ for }p=1,\ldots,P \text{ and } t=1,\ldots,T,
where \epsilon^1_{t,p}\sim N(0,(\sigma^1)_k^2)
for t=1,\ldots,T
, while \eta^1_1,\eta^1_2,\rho^1,(\sigma^1)_k^2
are additional parameters to be estimated.
Similarly for k^2_{t,p}
and k^3_{t,p}
.
Also,
\gamma_{c,p} = \gamma_{c-1,p}+ \rho_\gamma (\gamma_{c-1,p}-\gamma_{c-2,p}) + \epsilon^\gamma_{c,p} \text{ for }p=1,\ldots,P \text{ and } c=3,\ldots,C,
\gamma_2=\gamma_1+\frac{1}{\sqrt{1-\rho_\gamma^2}}\epsilon^\gamma_{2,p},
\gamma_1=100\epsilon^\gamma_{1,p}
where \epsilon^\gamma_{c,p}\sim N(0,\sigma_\gamma^2)
for c=1,\ldots,C
, while \rho_\gamma,\sigma_\gamma^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
Plat R. (2009). On Stochastic Mortality Modeling. Insurance: Mathematics and Economics, 45(3), 393–404. doi:10.1016/j.insmatheco.2009.08.006
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. doi: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. doi: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_PLAT_result<-fit_PLAT(death=death,expo=expo,n_iter=50,n.adapt=50)
head(fit_PLAT_result)
#if sharing the cohorts only (poisson family)
fit_PLAT_result2<-fit_PLAT(death=death,expo=expo,n_iter=1000,family="poisson",share_gamma=TRUE)
head(fit_PLAT_result2)