Joint Modeling of Longitudinal and Survival Data

Shared parameter models for the joint modeling of longitudinal and time-to-event data.

JM: Joint Models for Longitudinal and Survival Data using Maximum Likelihood

Travis-CI Build Status CRAN status Download counter Research software impact


This repository contains the source files for the R package JM. This package fits joint models for longitudinal and time-to-event data using maximum likelihood. These models are applicable in mainly two settings. First, when focus is on the survival outcome and we wish to account for the effect of an endogenous (aka internal) time-dependent covariates measured with error. Second, when focus is on the longitudinal outcome and we wish to correct for nonrandom dropout.

The basic joint-model-fitting function of the package is jointModel(). This accepts as main arguments a linear mixed model fitted by function lme() from the nlme package and a Cox model fitted using function coxph() from the survival package.

Basic Features

  • It can fit joint models for a single continuous longitudinal outcome and a time-to-event outcome.

  • For the survival outcome a relative risk models is assumed. The method argument of jointModel() can be used to define the type of baseline hazard function. Options are a B-spline approximation, a piecewise-constant function, the Weibull hazard and a completely unspecified function (i.e., a discrete function with point masses at the unique event times).

  • The user has now the option to define custom transformation functions for the terms of the longitudinal submodel that enter into the linear predictor of the survival submodel (arguments derivForm, parameterization). For example, the current value of the longitudinal outcomes, the velocity of the longitudinal outcome (slope), the area under the longitudinal profile. From the aforementioned options, in each model up to two terms can be included. In addition, using argument InterFact interactions terms can be considered.

Dynamic predictions

  • Function survfitJM() computes dynamic survival probabilities.

  • Function predict() computes dynamic predictions for the longitudinal outcome.

  • Function aucJM() calculates time-dependent AUCs for joint models, and function rocJM() calculates the corresponding time-dependent sensitivities and specifies.

  • Function prederrJM() calculates prediction errors for joint models.


Changes in version: JM_1.4-7

  • Changed the definition of JM:::rmvt() to work under the shifted method (noted by Graeme Hickey).
  • Resolved issue with non-ordered id numbers in the longitudinal submodel.


Changes in version: JM_1.4-6

  • Updated aucJM() method.


Changes in version: JM_1.4-5

  • Fixed a bug in fixef.jointModel().


Changes in version: JM_1.4-4

  • Fixed a bug in


Changes in version: JM_1.4-3

  • Small updates.


Changes in version: JM_1.4-2

  • Small updates.


Changes in version: JM_1.3-0

  • Several minor improvements.


Changes in version: JM_1.2-0

  • the new generic function aucJM() calculates time-dependent AUCs for joint models.

  • an updated version of function dynCJM() calculates a dynamic discrimination index (weighted average of time-dependent AUCs) for joint models.

  • the new generic function prederrJM() calculates prediction errors for joint models.

  • survfitJM() is now a generic function with a method for 'jointModel' objects.

  • new versions of functions ins() and ibs() with updated '' argument, and makepredictcall() methods.


Changes in version: JM_1.1-0

  • a small bug has been corrected in the plot() method for 'jointModel' objects, when method = "piecewise-PH-aGH" or method "piecewise-PH-GH" was used.


Changes in version: JM_1.0-1

  • use of globalVariables() in source code.

  • a small bug has been corrected in the plot() method for 'jointModel' objects, when a random intercepts linear mixed model was used


Changes in version: JM_1.0-0

  • This is the version of the package related to the book: Rizopoulos, D. (2012). Joint Models for Longitudinal and Time-to-Event Data: with Applications in R. Boca Raton: Chapman & Hall/CRC.

  • functions dns(), dbs(), ins() and ibs() calculate numerically derivative and integrals for functions ns() and bs(), respectively.

  • a coef() method has been added for objects of class 'summary.jointModel'.


Changes in version: JM_0.9-2

  • prediction.jointModel() can also compute now prediction intervals.

  • rocJM() has a new logical argument 'directionSmaller' denoting whether smaller values for the longitudinal outcome are associated with higher risk for an event.

  • fitted.jointModel() has the new option 'Slope' for the 'type' argument that returns the fitted values corresponding to the slope term when parameterization %in% c('slope', 'both') in jointModel().


Changes in version: JM_0.9-1

  • minor bug fixes.


Changes in version: JM_0.9-0

  • jointModel() can now also handle exogenous time-dependent covariates when method = "spline-PH-aGH".

  • jointModel() can now also handle competing risks settings when method = "spline-PH-aGH".

  • new function crLong() expands a data frame in the long format in the competing risks setting.

  • predict() method now calculates marginal and subject-specific predictions for the longitudinal outcome.


Changes in version: JM_0.8-4

  • method ranef() has now the extra argument 'type' the specifies whether to compute the mean (default) or the mode of the posterior distribution of the random effects.

  • the anova() method now also produces marginal Wald tests when a single joint model is provided.

  • plot.survfitJM() produces a more informative plot when argument 'include.y' is set to TRUE.

  • a bug has been corrected in residuals.jointModel() that it did not work when 'MI = TRUE', and 'parameterization = "slope"' in jointModel().


Changes in version: JM_0.8-3

  • the default method is now the Weibull model under the relative risk parameterization using the pseudo-adaptive Gauss-Hermite rule.

  • the plot() method has a new logical argument called 'return', which if set to TRUE the values use to create the plot are returned.

  • a typo in the code creating the scaling for the pseudo-adaptive Gauss-Hermite points has been corrected. This was primarily affecting the standard errors in the longitudinal submodel. The point estimates may also slightly change in some datasets.


Changes in version: JM_0.8-2

  • a bug was corrected in the internal function ModelMats().


Changes in version: JM_0.8-1

  • the new function xtable.jointModel() in conjunction with the xtable package can be used to produce a LaTeX table with the results of joint modeling analysis.

  • the new function simulateJM() and the simulate() method for objects of class 'jointModel' can be used to simulate data from a joint model.


Changes in version: JM_0.8-0

  • the new argument 'interFact' added in jointModel() allows the specification of interaction terms between the longitudinal outcome and baseline covariates.

  • for all joint models fitted in JM there is now the option to use a pseudo adaptive Gauss-Hermite rule. This is much faster than the default option and produces results of equal or better quality.

  • a predict() method has been added. Currently this only calculates fitted average longitudinal evolutions based on the information provided in the 'newdata' argument.

  • a new algorithm for calculating the starting values has been implemented. In most of the cases these will be closer to the MLEs than in the previous version.

  • some small changes have been made in the default Gauss-Hermite quadrature rule. This will result in minor changes in parameter estimates, standard errors and log-likelihood value compared to the previous version.

  • a bug has been corrected in the code used to specify the design matrix for the random effects in the longitudinal outcome, that did not allow this matrix not to be a subset of the design matrix of the fixed effects.


Changes in version: JM_0.7-0

  • the new function rocJM() has been added that calculates time-dependent ROC curves and the corresponding AUCs for joint models.

  • methods "weibull-AFT-GH", "weibull-PH-GH", "piecewise-PH-GH", and "spline-PH-GH" support now the true slope parameterization. This is invoked be specifying the 'parameterization' and 'derivForm' arguments accordingly.


Changes in version: JM_0.6-2

  • a small bug was corrected in summary.jointModel().


Changes in version: JM_0.6-1

  • jointModel() has now the extra argument 'scaleWB' that allows to fix the scale parameter for the Weibull baseline hazard to a specific value.


Changes in version: JM_0.6-0

  • method = "spline-PH-GH" allows now to include stratification factors for which different spline coefficients are estimated. By default the knots positions are the same across strata -- this can be changed by either directly specifying the knots or by setting the control argument 'equal.strata.knots' to FALSE.

  • the new function wald.strata() can be used to test for equality of the spline coefficients among strata.

  • a confint() method has been introduced for 'jointModel' objects.

  • jointModel() has now the extra argument 'lag' that allows for lagged effects in the time-dependent covariate represented by the linear mixed model.

  • a bug was corrected in joint models with piecewise constant baseline risk function. In particular, the 'xi' parameters were reported as double their actual value.


Changes in version: JM_0.5-0

  • function dynC() has been added that calculates a dynamic concordance index for joint models.

  • method = "ch-GH" has been replaced by method = "spline-PH-GH" that fits a relative risk model with a B-spline-approximated baseline risk function.

  • method = "ph-GH" that fits a relative risk with an unspecified baseline risk function has been renamed to method = "Cox-PH-GH".


Changes in version: JM_0.4-0

  • function survfitJM() has been added that calculates predictions of subject-specific probabilities of survival given a history of longitudinal responses.

  • the multiple-imputation residuals now work also for joint models with piecewise constant baseline risk functions.

  • faster optimization algorithms have implemented for 'method = "weibull-PH-GH"' and 'method = "piecewise-PH-GH".


Changes in version: JM_0.3-0

  • the Weibull model is now available under both the relative risk and accelerated failure time parameterizations.

  • a number of enhancements have been implemented in the functions that compute the MI-based residuals.

  • new more robust algorithms have been written for the numerical approximation of integrals; this will lead to some discrepancies in the results, especially in the survival part, compared to the previous versions of the package.


Changes in version: JM_0.2-1

  • changes in e-mail addresses in .Rd files.


Changes in version: JM_0.2-0

  • the jointModel method for the residuals generic has further options: (i) MI residuals for fixed and random visit times for the longitudinal process, and (ii) martingale, Cox-Snell, and AFT residuals for the survival process.

  • Function weibull.frailty() is introduced (along with supporting methods) for fitting multivariate survival data using the Weibull model with Gamma multiplicative frailties under maximum likelihood.

  • several typos have been corrected in .Rd files.


Changes in version: JM_0.1-1

  • corrected some typos in .Rd files.

Reference manual

It appears you don't have a PDF plugin for this browser. You can click here to download the reference manual.


1.4-8 by Dimitris Rizopoulos, 4 years ago

Browse source code at

Authors: Dimitris Rizopoulos <[email protected]>

Documentation:   PDF Manual  

Task views: Survival Analysis

GPL (>= 2) license

Depends on MASS, nlme, splines, survival

Enhances xtable

Imported by winRatioAnalysis.

Depended on by SEchart.

Suggested by insight, joineRML.

See at CRAN