R/MC.R

In demofit: Parametric Mortality Curve Fitting and Mortality Forecasting Tools

#' Mortality curve fitting
#'
#' Fits parametric mortality curves to observed mortality rates using multiple optimisation strategies.
#'
#' @param x vector of ages.
#' @param m vector of mortality rates.
#' @param curve name of mortality curve (including gompertz, makeham, oppermann, thiele, wittsteinbumsted, perks, weibull, vandermaen, beard, heligmanpollard, rogersplanck, siler, martinelle, thatcher, gompertz2, makeham2, oppermann2, thiele2, wittsteinbumsted2, perks2, weibull2, vandermaen2, beard2, heligmanpollard2, rogersplanck2, siler2, martinelle2, thatcher2, where first 14 curves' parameters are unconstrained and last 14 curves' parameters are generally restricted to be positive).
#' @param w vector of weights (default = 1).
#'
#' @details
#' User can choose one of the following 14 mortality curves (with their suitable age ranges in brackets):
#'
#' \itemize{
#' \item Gompertz (1825) (ages 30-90) - \eqn{m_x = B e^{Cx}}
#' \item Makeham (1860) (ages 20-90) - \eqn{m_x = A + B e^{Cx}}
#' \item Oppermann (1870) (ages 0-20) - \eqn{m_x = \frac{A}{\sqrt{x+1}} + B + C\sqrt{x+1}}
#' \item Thiele (1871) (ages 0-90) - \eqn{m_x = A_1 e^{-B_1 x} + A_2 e^{-0.5 B_2 (x-C)^2} + A_3 e^{B_3 x}}
#' \item Wittstein & Bumsted (1883) (ages 0-90) - \eqn{m_x = \frac{A^{-(Bx)^N}}{B} + A^{-(M-x)^N}}
#' \item Perks (1932) (ages 20-100+) - \eqn{m_x = \frac{A+BC^x}{1+DC^x}}
#' \item Weibull (1939) (ages 50-100+) - \eqn{m_x = Bx^C}
#' \item Van der Maen (1943) (ages 80-100+) - \eqn{m_x = A + Bx + Cx^2 + \frac{I}{N-x}}
#' \item Beard (1971) (ages 50-100+) - \eqn{m_x = \frac{A e^{Cx}}{1+B e^{Cx}}}
#' \item Heligman & Pollard (1980) (ages 0-100+) - \eqn{m_x = A^{(x+B)^C} + D e^{-E(ln(x)-ln(F))^2} + \frac{GH^x}{1+GH^x}}
#' \item Rogers & Planck (1983) (ages 0-100+) - \eqn{m_x = A_0 + A_1 e^{-Ax} + A_2 e^{-B(x-U)-e^{-C(x-U)}} + A_3 e^{Dx}}
#' \item Siler (1983) (ages 0-90) - \eqn{m_x = A_1 e^{-B_1 x} + A_2 + A_3 e^{B_3 x}}
#' \item Martinelle (1987) (ages 20-100+) - \eqn{m_x = \frac{A+B e^{Cx}}{1+D e^{Cx}} + E e^{Cx}}
#' \item Thatcher (1999) (ages 20-100+) - \eqn{m_x = A + \frac{B e^{Cx}}{1+D e^{Cx}}}
#' }
#'
#' The mortality curves are fitted by employing multiple optimisation strategies simultaneously (including PORT routines, Nelder-Mead method, and Levenberg-Marquardt algorithm) and selecting one with smallest weighted least squares on \eqn{ln(m_x)}. Constrained parameterisations offer traditional interpretability while unconstrained parameterisations provide increased flexibility.
#'
#' @return 
#' An object of class based on selected mortality curve with associated S3 methods coef, fitted, predict, plot, deviance, and residuals.
#'
#' @references
#' Gompertz, B. (1825). On the nature of the function expressive of the law of human mortality, and on a new mode of determining the value of life contingencies. Philosophical Transactions of the Royal Society of London, 115(1825), 513-583.
#'
#' Makeham, W.M. (1860). On the law of mortality and the construction of annuity tables. Journal of the Institute of Actuaries, 8(6), 301-310.
#'
#' Oppermann, L.H.F. (1870). On the graduation of life tables, with special application to the rate of mortality in infancy and childhood. The Institute of Actuaries.
#'
#' Thiele, T.N. (1871). On a mathematical formula to express the rate of mortality throughout the whole of life, tested by a series of observations made use of by the Danish Life Insurance Company of 1871. Journal of the Institute of Actuaries and Assurance Magazine, 16(5), 313-329.
#'
#' Wittstein, T. and Bumsted, D.A. (1883). The mathematical law of mortality. Journal of the Institute of Actuaries and Assurance Magazine, 24(3), 153-173.
#'
#' Perks, W. (1932). On some experiments in the graduation of mortality statistics. Journal of the Institute of Actuaries, 63(1), 12-57.
#'
#' Weibull, W. (1951). A statistical distribution function of wide applicability. Journal of Applied Mechanics, 18(3), 293-297.
#'
#' Beard, R.E. (1971). Some aspects of theories of mortality, cause of death analysis, forecasting and stochastic processes. Biological Aspects of Demography, 57-68.
#'
#' Heligman, L. and Pollard, J.H. (1980). The age pattern of mortality. Journal of the Institute of Actuaries, 107(1), 49-80.
#'
#' Rogers, A. and Planck, F. (1983). Model: A general program for estimating parametrized model schedules of fertility, mortality, migration, and marital and labor force status transitions. IIASA Working Paper WP-83-102.
#'
#' Siler, W. (1983). Parameters of mortality in human populations with widely varying life spans. Statistics in Medicine, 2(3), 373-380.
#'
#' Martinelle, S. (1987). A generalized Perks formula for old-age mortality. Statistiska Centralbyran.
#'
#' Thatcher, A.R. (1999). The long-term pattern of adult mortality and the highest attained age. Journal of the Royal Statistical Society Series A, 162(1), 5-43.
#'
#' Tabeau, E. (2001). A review of demographic forecasting models for mortality. Forecasting Mortality in Developed Countries. European Studies of Population, Vol 9.
#'
#' @examples
#' x <- 60:89
#' set.seed(123); m <- 0.0000082*exp(0.10771*c(60:89)+rnorm(30,0,0.1))
#' fit <- MC(x=x,m=m,curve="gompertz")
#' coef(fit)
#' fitted(fit)
#' predict(fit,62.5)
#' plot(fit)
#' deviance(fit)
#' residuals(fit)
#'
#' @export
MC <- function(x,m,curve=c("gompertz","makeham","oppermann","thiele","wittsteinbumsted","perks","weibull","vandermaen","beard","heligmanpollard","rogersplanck","siler","martinelle","thatcher","gompertz2","makeham2","oppermann2","thiele2","wittsteinbumsted2","perks2","weibull2","vandermaen2","beard2","heligmanpollard2","rogersplanck2","siler2","martinelle2","thatcher2"),w=rep(1,length(x))) {
curve <- tryCatch(match.arg(curve),error = function(e) { stop("invalid curve choice") })
tryCatch({
fit <- switch(curve,
gompertz = fit1(x,m,w),
makeham = fit2(x,m,w),
oppermann = fit3(x,m,w),
thiele = fit4(x,m,w),
wittsteinbumsted = fit5(x,m,w),
perks = fit6(x,m,w),
weibull = fit7(x,m,w),
vandermaen = fit8(x,m,w),
beard = fit9(x,m,w),
heligmanpollard = fit10(x,m,w),
rogersplanck = fit11(x,m,w),
siler = fit12(x,m,w),
martinelle = fit13(x,m,w),
thatcher = fit14(x,m,w),
gompertz2 = cfit1(x,m,w),
makeham2 = cfit2(x,m,w),
oppermann2 = cfit3(x,m,w),
thiele2 = cfit4(x,m,w),
wittsteinbumsted2 = cfit5(x,m,w),
perks2 = cfit6(x,m,w),
weibull2 = cfit7(x,m,w),
vandermaen2 = cfit8(x,m,w),
beard2 = cfit9(x,m,w),
heligmanpollard2 = cfit10(x,m,w),
rogersplanck2 = cfit11(x,m,w),
siler2 = cfit12(x,m,w),
martinelle2 = cfit13(x,m,w),
thatcher2 = cfit14(x,m,w)
)
invisible(fit)
}, error = function(e) { stop(paste0("curve fitting is unsuccessful - please make sure the data and age range are suitable for the curve\n",e$message),call.=FALSE) })
}