Learn R Programming

MortalityLaws (version 1.4.0)

MortalityLaw: Fit Mortality Laws

Description

Fit parametric mortality models given a set of input data which can be represented by death counts and mid-interval population estimates (Dx, Ex) or age-specific death rates (mx) or death probabilities (qx). Using the argument law one can specify the model to be fitted. So far 27 parametric model have been implemented; check availableLaws function to learn about the available options. The models can be fitted under the maximum likelihood methodology or by selecting a loss function to be optimised. See the implemented loss function by running availableLF function.

Usage

MortalityLaw(x, Dx = NULL, Ex = NULL, mx = NULL, qx = NULL, 
                law = NULL, 
                opt.method = "poissonL", 
                parS = NULL, 
                fit.this.x = x,
                scale.x = FALSE,
                custom.law = NULL, 
                show = TRUE)

Arguments

x

Vector of ages at the beginning of the age interval.

Dx

Object containing death counts. An element of the Dx object, represents the number of deaths during the year to persons aged x to x+n.

Ex

Exposure in the period. Ex can be approximated by the mid-year population aged x to x+n.

mx

Death rate in age interval [x, x+n).

qx

Probability of dying in age interval [x, x+n).

law

The name of the mortality law/model to be fitted. e.g. gompertz, makeham, ... To investigate all the possible options, see availableLaws function.

opt.method

How would you like to find the parameters? Specify the function to be optimize. Available options: the Poisson likelihood function poissonL; the Binomial likelihood function -binomialL; and other 6 loss functions. For more details, check availableLF function.

parS

Starting parameters used in optimization process (optional).

fit.this.x

Select the ages to be considered in model fitting. By default fit.this.x = x. One may want exclude from the fitting procedure say the advance ages were the data is sparse.

scale.x

Logical. Scale down "x" vector so that is begins with a small value. This is useful in order to obtain meaningful estimates and sometimes a better fit. Default: FALSE. Method: new.x = x - min(x) + 1.

custom.law

Allows you to fit a model that is not defined in the package. Accepts as input a function.

show

Choose whether to display a progress bar during the fitting process. Logical. Default: TRUE.

Value

The output is of "MortalityLaw" class with the components:

input

List with arguments provided in input. Saved for convenience

info

Short information about the model

coefficients

Estimated coefficients

fitted.values

Fitted values of the selected model

residuals

Deviance residuals

goodness.of.fit

List containing goodness of fit measures like AIC, BIC and log-Likelihood

opt.diagnosis

Resulted optimization object useful for checking the convergence etc.

stats

List containing statistical measures like: parameter correlation, standard errors, degrees of freedom, deviance, gradient matrix, QR decomposition, covariance matrix etc.

Details

Depending on the complexity of the model, one of following optimization strategies are employed:

  1. Nelder-Mead method: approximates a local optimum of a problem with n variables when the objective function varies smoothly and is unimodal. For details see optim

  2. PORT routines: provides unconstrained optimization and optimization subject to box constraints for complicated functions. For details check nlminb

  3. Levenberg-Marquardt algorithm: damped least-squares method. For details check nls.lm

Examples

Run this code
# NOT RUN {
# Example 1: ---
# Fit Makeham Model for Year of 1950.

x  <- 45:75
Dx <- ahmd$Dx[paste(x), "1950"]
Ex <- ahmd$Ex[paste(x), "1950"]

M1 <- MortalityLaw(x   = x, 
                   Dx  = Dx, 
                   Ex  = Ex, 
                   law = 'makeham',
                   scale.x = TRUE)
M1
ls(M1)
coef(M1)
summary(M1)
fitted(M1)
predict(M1, x = 45:95)
plot(M1)


# Example 2: ---
# We can fit the same model using a different data format 
# and a different optimization method.
x  <- 45:75
mx <- ahmd$mx[paste(x), ]
M2 <- MortalityLaw(x   = x, 
                   mx  = mx, 
                   law = 'makeham', 
                   opt.method = 'LF1', 
                   scale.x = TRUE)
M2
fitted(M2)
predict(M2, x = 55:90)

# Example 3: ---
# Now let's fit a mortality law that is not defined 
# in the package, say a reparameterized Gompertz in 
# terms of modal age at death
# hx = b*exp(b*(x-m)) (here b and m are the parameters to be estimated)

# A function with 'x' and 'par' as input has to be defined, which returns at least
# an object called 'hx' (hazard rate).
my_gompertz <- function(x, par = c(b = 0.13, M = 45)){
  hx  <- with(as.list(par), b*exp(b*(x - M)) )
  return(as.list(environment()))
}

M3 <- MortalityLaw(x  = x,
                   Dx = Dx, 
                   Ex = Ex, 
                   custom.law = my_gompertz,
                   scale.x = TRUE) 
summary(M3)
# predict M3 for different ages
predict(M3, x = 85:130)


# Example 4: ---
# Fit Heligman-Pollard model for a single 
# year in the dataset between age 0 and 100.

x  <- 0:100
mx <- ahmd$mx[paste(x), "1950"] # select data
M4 <- MortalityLaw(x   = x, 
                   mx  = mx, 
                   law = 'HP', 
                   opt.method = 'LF2')
M4
plot(M4)
# }

Run the code above in your browser using DataLab