This function fits joint latent class models for multivariate longitudinal markers and competing causes of event. It is a multivariate extension of the Jointlcmm function. It defines each longitudinal dimension as a latent process (mp in mpjlcmm is for multivariate processes), possibly measured by sereval continuous markers (Gaussian or curvilinear). For the moment, theses processes are assumed independent given the latent classes. The (optional) survival part handles competing risks, right censoring and left truncation. The specification of the function is similar to other estimating functions of the package.
Usage
mpjlcmm(
longitudinal,
subject,
classmb,
ng,
survival,
hazard = "Weibull",
hazardtype = "Specific",
hazardnodes = NULL,
TimeDepVar = NULL,
data,
B,
convB = 1e-04,
convL = 1e-04,
convG = 1e-04,
maxiter = 100,
nsim = 100,
prior,
logscale = FALSE,
subset = NULL,
na.action = 1,
posfix = NULL,
partialH = FALSE,
verbose = FALSE,
nproc = 1,
clustertype = NULL
)Arguments
- longitudinal
list of longitudinal models of type hlme, lcmm or multlcmm. Each model defines the structure of one latent process.
- subject
name of the covariate representing the grouping structure (called subject identifier)
- classmb
optional one-sided formula describing the covariates in the class-membership multinomial logistic model
- ng
number of latent classes considered
- survival
two-sided formula object specifying the survival part of the model
- hazard
optional family of hazard function assumed for the survival model (Weibull, piecewise or splines)
- hazardtype
optional indicator for the type of baseline risk function (Specific, PH or Common)
- hazardnodes
optional vector containing interior nodes if
splinesorpiecewiseis specified for the baseline hazard function inhazard- TimeDepVar
optional vector specifying the name of the time-dependent covariate in the survival model
- data
data frame containing all the variables used in the model
- B
optional specification for the initial values of the parameters. Three options are allowed: (1) a vector of initial values is entered (the order in which the parameters are included is detailed in
detailssection). (2) nothing is specified. Initial values are extracted from the models specified inlongitudinal, and default initial values are chosen for the survival part (3) when ng>1, a mpjlcmm object is entered. It should correspond to the exact same structure of model but with ng=1. The program will automatically generate initial values from this model. Note that due to possible local maxima, theBvector should be specified and several different starting points should be tried. This can be done automatically using gridsearch function.- convB
optional threshold for the convergence criterion based on the parameter stability
- convL
optional threshold for the convergence criterion based on the log-likelihood stability
- convG
optional threshold for the convergence criterion based on the derivatives
- maxiter
optional maximum number of iterations for the Marquardt iterative algorithm
- nsim
optional number of points for the predicted survival curves and predicted baseline risk curves
- prior
optional name of a covariate containing a prior information about the latent class membership
- logscale
optional boolean indicating whether an exponential (logscale=TRUE) or a square (logscale=FALSE -by default) transformation is used to ensure positivity of parameters in the baseline risk functions
- subset
a specification of the rows to be used: defaults to all rows. This can be any valid indexing vector for the rows of data or if that is not supplied, a data frame made up of the variable used in formula.
- na.action
Integer indicating how NAs are managed. The default is 1 for 'na.omit'. The alternative is 2 for 'na.fail'. Other options such as 'na.pass' or 'na.exclude' are not implemented in the current version.
- posfix
Optional vector specifying the indices in vector B of the parameters that should not be estimated. Default to NULL, all parameters are estimated
- partialH
optional logical for Piecewise and Splines baseline risk functions and Splines link functions only. Indicates whether the parameters of the baseline risk or link functions can be dropped from the Hessian matrix to define convergence criteria (can solve non convergence due to estimates at the boundary of the parameter space - usually 0).
- verbose
logical indicating whether information about computation should be reported. Default to FALSE.
- nproc
the number cores for parallel computation. Default to 1 (sequential mode).
- clustertype
optional character indicating the type of cluster for parallel computation.
Examples
if (FALSE) { # \dontrun{
paquid$age65 <- (paquid$age-65)/10
##############################################################################
### EXAMPLE 1 : ###
###two outcomes measuring the same latent process along with dementia onset###
##############################################################################
## multlcmm model for MMSE and BVRT for 1 class
mult1 <- multlcmm(MMSE+BVRT~age65+I(age65^2)+CEP+male,random=~age65+I(age65^2),
subject="ID",link=c("5-quant-splines","4-quant-splines"),data=paquid)
summary(mult1)
## joint model for 1 class
m1S <- mpjlcmm(longitudinal=list(mult1),subject="ID",ng=1,data=paquid,
survival=Surv(age_init,agedem,dem)~1)
summary(m1S)
##### joint model for 2 classes #####
## specify longitudinal model for 2 classes, without estimation
mult2 <- multlcmm(MMSE+BVRT~age65+I(age65^2)+CEP+male,random=~age65+I(age65^2),
subject="ID",link=c("5-quant-splines","4-quant-splines"),ng=2,
mixture=~age65+I(age65^2),data=paquid,B=random(mult1),maxiter=0)
## estimation of the associated joint model
m2S <- mpjlcmm(longitudinal=list(mult2),subject="ID",ng=2,data=paquid,
survival=Surv(age_init,agedem,dem)~1)
## estimation by a grid search with 50 replicates, initial values
## randomly generated from m1S
m2S_b <- gridsearch(mpjlcmm(longitudinal=list(mult2),subject="ID",ng=2,
data=paquid,survival=Surv(age_init,agedem,dem)~1), minit=m1S, rep=50, maxiter=30)
##### joint model for 3 classes #####
mult3 <- multlcmm(MMSE+BVRT~age65+I(age65^2)+CEP+male,random=~age65+I(age65^2),
subject="ID",link=c("5-quant-splines","4-quant-splines"),ng=3,
mixture=~age65+I(age65^2),data=paquid,B=random(mult1),maxiter=0)
m3S <- mpjlcmm(longitudinal=list(mult3),subject="ID",ng=3,data=paquid,
survival=Surv(age_init,agedem,dem)~1)
m3S_b <- gridsearch(mpjlcmm(longitudinal=list(mult3),subject="ID",ng=3,
data=paquid,survival=Surv(age_init,agedem,dem)~1), minit=m1S, rep=50, maxiter=30)
##### summary of the models #####
summarytable(m1S,m2S,m2S_b,m3S,m3S_b)
##### post-fit #####
## update longitudinal models :
mod2 <- update(m2S)
mult2_post <- mod2[[1]]
## -> use the available functions for multlcmm on the mult2_post object
## fit of the longitudinal trajectories
par(mfrow=c(2,2))
plot(mult2_post,"fit","age65",marg=TRUE,shades=TRUE,outcome=1)
plot(mult2_post,"fit","age65",marg=TRUE,shades=TRUE,outcome=2)
plot(mult2_post,"fit","age65",marg=FALSE,shades=TRUE,outcome=1)
plot(mult2_post,"fit","age65",marg=FALSE,shades=TRUE,outcome=2)
## predicted trajectories
dpred <- data.frame(age65=seq(0,3,0.1),male=0,CEP=0)
predL <- predictL(mult2_post,newdata=dpred,var.time="age65",confint=TRUE)
plot(predL,shades=TRUE) # in the latent process scale
predY <- predictY(mult2_post,newdata=dpred,var.time="age65",draws=TRUE)
plot(predY,shades=TRUE,ylim=c(0,30),main="MMSE") #in the 0-30 scale for MMSE
plot(predY,shades=TRUE,ylim=c(0,15),outcome=2,main="BVRT") #in 0-15 for BVRT
## baseline hazard and survival curves :
plot(m2S,"hazard")
plot(m2S,"survival")
## posteriori probabilities and classification :
postprob(m2S)
####################################################################################
### EXAMPLE 2 : ###
### two latent processes measured each by one outcome along with dementia onset ###
####################################################################################
## define the two longitudinal models
mMMSE1 <- lcmm(MMSE~age65+I(age65^2)+CEP,random=~age65+I(age65^2),subject="ID",
link="5-quant-splines",data=paquid)
mCESD1 <- lcmm(CESD~age65+I(age65^2)+male,random=~age65+I(age65^2),subject="ID",
link="5-quant-splines",data=paquid)
## joint estimation
mm1S <- mpjlcmm(longitudinal=list(mMMSE1,mCESD1),subject="ID",ng=1,data=paquid,
survival=Surv(age_init,agedem,dem)~CEP+male)
## with 2 latent classes
mMMSE2 <- lcmm(MMSE~age65+I(age65^2)+CEP,random=~age65+I(age65^2),subject="ID",
link="5-quant-splines",data=paquid,ng=2,mixture=~age65+I(age65^2),
B=random(mMMSE1),maxiter=0)
mCESD2 <- lcmm(CESD~age65+I(age65^2)+male,random=~age65+I(age65^2),subject="ID",
link="5-quant-splines",data=paquid,ng=2,mixture=~age65+I(age65^2),
B=random(mCESD1),maxiter=0)
mm2S <- mpjlcmm(longitudinal=list(mMMSE2,mCESD2),subject="ID",ng=2,data=paquid,
survival=Surv(age_init,agedem,dem)~CEP+mixture(male),classmb=~CEP+male)
mm2S_b <- gridsearch(mpjlcmm(longitudinal=list(mMMSE2,mCESD2),subject="ID",ng=2,
data=paquid,survival=Surv(age_init,agedem,dem)~CEP+mixture(male),
classmb=~CEP+male),minit=mm1S,rep=50,maxiter=50)
summarytable(mm1S,mm2S,mm2S_b)
mod1_biv <- update(mm1S)
mod2_biv <- update(mm2S)
## -> use post-fit functions as in exemple 1
} # }