Creating a life table from StMoMo Forecasts

[ibo2]

Dear all,

We have a problem, with which we hope one of you are able to help out with. We want to create a life expectancy table using either the 'lifecontingencies' or 'life.expectanies' function in the demography package. However, the life.expectancies package takes only demogdata objects, and we have StMoMo objects from our forecast.

If anyone has had this problem, or has an idea on how to do solve it, then contact us.

[pmilloss]

Hi,

check the code below.

Pietro

`library("StMoMo") library(demography)

LCfit <- fit(object = lc(), data = EWMaleData, ages.fit = 55:89, years.fit = 1961 : 2011)

LCfor <- forecast(LCfit, h = 50)

forecast.rates <- LCfor$rates

forecast.as.demogdata <- demogdata(data = forecast.rates, pop = matrix(10 ^ 5, nrow = nrow(forecast.rates), ncol = ncol(forecast.rates)), ages = 55 : 89, years = 2012 : 2061, type = "mortality", label = "UK", name = "male")

life.expectancy(forecast.as.demogdata) `

[spedygiorgio]

@ibo2, the issue is handled in the vignette on mortality modeling (currently in development version) that is going to be submitted to cran in these days

[spedygiorgio]

@ibo2, the new version of lifecontingencies package has been released on cran

[andersmaxDK]

@pmilloss, do you have a solution for this problem when using the bootstrap-command in your package?

I have bootstrapped Danish data for 1947 to 2016, and I now have 5.000 elements of nBootParameters. The rates are not given directly from nBootparameters, so I am trying to create a loop that calculate the mortality rates for each age for every year of my data, and then I will repeat this process 5.000 times.

Do you have an easier solution than this?

[pmilloss]

@andersmaxDK ,

there is no simple way to extract coefficients or rates from the boostrap.fitStMoMo, definitely something we may want to add in the package in one of the next releases.

For the moment you can use a calculation like the one below (which works for the LC model). In your case, it depends on the model you used, but once you have the parameters you can just reconstruct the rates.

Regards,

Pietro

library("StMoMo")

LCfit <- fit(lc(), data = EWMaleData)

number of bootsrap samples

nBoot <- 5

LCResBoot <- bootstrap(LCfit, nBoot = nBoot, type = "residual")

just for the first bootstrap sample

ax1 <- LCResBoot$bootParameters[[1]]$ax bx1 <- as.numeric(LCResBoot$bootParameters[[1]]$bx) kt1 <- as.numeric(LCResBoot$bootParameters[[1]]$kt) mxt1 <- exp(ax1 + outer(bx1, kt1))

for all bootstrap samples

axBoot <- sapply(LCResBoot$bootParameters, function(x) x$ax) bxBoot <- sapply(LCResBoot$bootParameters, function(x) x$bx) ktBoot <- sapply(LCResBoot$bootParameters, function(x) x$kt)

mxtBoot <- array(data = NA, dim = c(nrow(axBoot), nrow(ktBoot), nBoot))

for (i in 1:nBoot) { mxtBoot[, , i] <- exp(axBoot[, i] + outer(bxBoot[, i], ktBoot[, i])) }

age 65

matplot(mxtBoot[66, , ], type = "l", col = "black")