Joint model for Bidirectional survival data using FastJM
jmcsB.RdThe function fits joint model for survival data with two events. It utilizes the FastJM package for obtaining the model parameter estimates.
Arguments
- dtlong
longitudinal data
- dtsurv
survival data with two event status along with event time
- longm
longitudinal model e.g. list(serBilir~drug * year,serBilir ~ drug * year)
- survm
survival model e.g. list(Surv(years,status2)~drug,Surv(time_2,status_2)~drug+age)
- rd
random effect component e.g. list(~year|id,~year|id)
- id
ID variable
- timeVar
time variable
- BIGdata
logical argument TRUE or FALSE
- samplesize
samplesize for bigdata
References
Li, Shanpeng, et al. "Efficient Algorithms and Implementation of a Semiparametric Joint Model for Longitudinal and Competing Risk Data: With Applications to Massive Biobank Data." Computational and Mathematical Methods in Medicine 2022 (2022).
Bhattacharjee, A., Rajbongshi, B. K., & Vishwakarma, G. K. (2024). jmBIG: enhancing dynamic risk prediction and personalized medicine through joint modeling of longitudinal and survival data in big routinely collected data. BMC Medical Research Methodology, 24(1), 172.
Examples
library(FastJM)
#> Warning: package 'FastJM' was built under R version 4.3.3
#> Loading required package: statmod
#> Loading required package: MASS
#>
#> Attaching package: 'MASS'
#> The following object is masked from 'package:JMbayes2':
#>
#> area
#> The following object is masked from 'package:GLMMadaptive':
#>
#> negative.binomial
#>
#> Attaching package: 'FastJM'
#> The following objects are masked from 'package:JMbayes2':
#>
#> fixef, ranef
#> The following objects are masked from 'package:GLMMadaptive':
#>
#> fixef, ranef
#> The following objects are masked from 'package:nlme':
#>
#> fixef, ranef
library(JMbayes2)
st_pbcid<-function(){
new_pbcid<-pbc2.id
new_pbcid$time_2<-rexp(n=nrow(pbc2.id),1/10)
cen_time<-runif(nrow(pbc2.id),min(new_pbcid$time_2),max(new_pbcid$time_2))
status_2<-ifelse(new_pbcid$time_2<cen_time,1,0)
new_pbcid$status_2<-status_2
new_pbcid$time_2<-ifelse(new_pbcid$time_2<cen_time,new_pbcid$time_2,cen_time)
new_pbcid$time_2<-ifelse(new_pbcid$time_2<new_pbcid$years,new_pbcid$years,new_pbcid$time_2)
new_pbcid
}
new_pbc2id<-st_pbcid()
pbc2$status_2<-rep(new_pbc2id$status_2,times=data.frame(table(pbc2$id))$Freq)
pbc2$time_2<-rep(new_pbc2id$time_2,times=data.frame(table(pbc2$id))$Freq)
pbc2_new<-pbc2[pbc2$id%in%c(1:50),]
new_pbc2id<-new_pbc2id[new_pbc2id$id%in%c(1:50),]
model_jmcs<-jmcsB(dtlong=pbc2_new,dtsurv=new_pbc2id,
longm=list(serBilir~drug*year,
serBilir~drug*year),
survm=list(Surv(years,status2)~drug,
Surv(time_2,status_2)~drug+age),
rd=list(~1|id,~1|id),
id='id',timeVar='year')
model_jmcs
#>
#> Joint model using FastJM
#> =========================================================
#> Summary of first joint model with Event 1:
#> -------------------------------------------
#> Longitudinal process:
#> Estimate SE Zvalue Pvalue
#> (Intercept) 4.747 1.268 3.744 0.000
#> drugD-penicil -0.303 2.308 -0.131 0.896
#> year 0.354 0.037 9.458 0.000
#> drugD-penicil:year 0.132 0.049 2.697 0.007
#> sigma^2 8.948 0.282 31.753 0.000
#>
#> Survival process for event1:
#> Estimate SE Zvalue Pvalue
#> drugD-penicil_1 0.172 0.608 0.283 0.777
#>
#> Association parameters :
#> Estimate SE Zvalue Pvalue
#> 1 0.25 0.064 3.91 0
#>
#> Variance Covariance matrix of Random effects:
#> [,1]
#> [1,] 34.56277
#>
#> Summary of second joint model with Event 2:
#> --------------------------------------------
#> Longitudinal process:
#> Estimate SE Zvalue Pvalue
#> (Intercept) 3.799 1.203 3.158 0.002
#> drugD-penicil 0.737 2.240 0.329 0.742
#> year 0.323 0.038 8.564 0.000
#> drugD-penicil:year 0.150 0.051 2.938 0.003
#> sigma^2 8.998 0.282 31.929 0.000
#>
#> Survival process for event2:
#> Estimate SE Zvalue Pvalue
#> drugD-penicil_1 1.094 0.555 1.970 0.777
#> age_1 0.018 0.024 0.757 0.777
#>
#> Association parameters :
#> Estimate SE Zvalue Pvalue
#> 1 0.195 0.054 3.052 0.002
#>
#> Variance Covariance matrix of Random effects:
#> [,1]
#> [1,] 32.898