View source: R/MortalityLaw_main.R
MortalityLaw | R Documentation |
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 more than 27 parametric models have been
implemented; check the 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 the
availableLF
function.
MortalityLaw(x, Dx = NULL, Ex = NULL, mx = NULL, qx = NULL,
law = NULL,
opt.method = "LF2",
parS = NULL,
fit.this.x = x,
custom.law = NULL,
show = FALSE, ...)
x |
Vector of ages at the beginning of the age interval. |
Dx |
Object containing death counts. An element of the |
Ex |
Exposure in the period. |
mx |
Life table 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 used. e.g.
|
opt.method |
How would you like to find the parameters? Specify the
function to be optimize. Available options: the Poisson likelihood function
|
parS |
Starting parameters used in the optimization process (optional). |
fit.this.x |
Select the ages to be considered in model fitting.
By default |
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: |
... |
Arguments to be passed to or from other methods. |
Depending on the complexity of the model, one of following optimization strategies is employed:
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
PORT routines: provides unconstrained optimization and optimization
subject to box constraints for complicated functions. For details check
nlminb
Levenberg-Marquardt algorithm: damped least-squares method.
For details check nls.lm
The output is of the "MortalityLaw"
class with the components:
input |
List with arguments provided in input. Saved for convenience. |
info |
Brief 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 |
Resultant optimization object useful for checking the convergence etc. |
Marius D. Pascariu
availableLaws
availableLF
LifeTable
ReadHMD
# 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')
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')
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)
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 and build a life table.
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)
LifeTable(x = x, qx = fitted(M4))
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.