1 Survival Analysis with Frailty model.

Why Frailty model: When Cox fails (extension of the Cox proportional Hazard model) and for recurrent events . * Model variation : shared, nested ,joint, additve

1.1 Parametric Frailty Model : Kidney dataset

##   id time status age sex disease frail
## 1  1    8      1  28   1   Other   2.3
## 2  1   16      1  28   1   Other   2.3
## 3  2   23      1  48   2      GN   1.9
## 4  2   13      0  48   2      GN   1.9
## 5  3   22      1  32   1   Other   1.2
## 6  3   28      1  32   1   Other   1.2
## 
## Frailty distribution: gamma 
## Baseline hazard distribution: Exponential 
## Loglikelihood: -333.248 
## 
##        ESTIMATE SE    p-val    
## theta   0.301   0.156          
## lambda  0.025   0.014          
## sex    -1.485   0.396 <.001 ***
## age     0.005   0.011 0.657    
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Kendall's Tau: 0.131
##   low    up 
## 0.104 0.492
## Gamma frailty model with Exponential baseline
##  id frailty
##  1  1.325  
##  2  1.207  
##  3  1.109  
##  4  0.632  
##  5  1.199  
##  6  1.065  
##  7  1.386  
##  8  0.723  
##  9  1.002  
##  10 0.604  
##  11 0.929  
##  12 0.999  
##  13 1.258  
##  14 0.684  
##  15 0.634  
##  16 1.072  
##  17 0.874  
##  18 0.86   
##  19 0.713  
##  20 1.031  
##  21 0.205  
##  22 0.704  
##  23 1.353  
##  24 1.103  
##  25 1.077  
##  26 0.759  
##  27 1.053  
##  28 1.403  
##  29 1.266  
##  30 1.188  
##  31 1.344  
##  32 1.176  
##  33 1.134  
##  34 0.919  
##  35 1.288  
##  36 0.871  
##  37 1.097  
##  38 0.756
## 
## Frailty distribution: inverse Gaussian 
## Baseline hazard distribution: Exponential 
## Loglikelihood: -333.85 
## 
##        ESTIMATE SE    p-val    
## theta   0.375   0.259          
## lambda  0.022   0.013          
## sex    -1.310   0.370 <.001 ***
## age     0.004   0.011 0.686    
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Kendall's Tau: 0.125

1.2 Semiparametric Frailty : CGD

##   tstart tstop status id    sex   treat
## 1      0   219      1  1 female  rIFN-g
## 2    219   373      1  1 female  rIFN-g
## 3    373   414      0  1 female  rIFN-g
## 4      0     8      1  2   male placebo
## 5      8    26      1  2   male placebo
## 6     26   152      1  2   male placebo
## Call: 
## emfrail(formula = Surv(tstart, tstop, status) ~ sex + treat + 
##     cluster(id), data = cgd)
## 
## Regression coefficients:
##               coef exp(coef) se(coef) adj. se      z    p
## sexfemale   -0.227     0.797    0.396   0.396 -0.575 0.57
## treatrIFN-g -1.052     0.349    0.310   0.310 -3.389 0.00
## Estimated distribution: gamma / left truncation: FALSE 
## 
## Fit summary:
## Commenges-Andersen test for heterogeneity: p-val  0.00172 
## no-frailty Log-likelihood: -331.997 
## Log-likelihood: -326.619 
## LRT: 1/2 * pchisq(10.8), p-val 0.00052
## 
## Frailty summary:
##                    estimate lower 95% upper 95%
## Var[Z]                0.821     0.231     1.854
## Kendall's tau         0.291     0.104     0.481
## Median concordance    0.289     0.101     0.491
## E[logZ]              -0.464    -1.164    -0.120
## Var[logZ]             1.241     0.260     4.341
## theta                 1.218     0.539     4.326
## Confidence intervals based on the likelihood function

## Warning: Removed 2 row(s) containing missing values (geom_path).

1.3 Frailty

1.3.1 Shared

## 
## Be patient. The program is computing ... 
## The program took 0.46 seconds
## Call:
## frailtyPenal(formula = Surv(time, event) ~ cluster(id) + as.factor(dukes) + 
##     as.factor(charlson) + sex + chemo, data = readmission, cross.validation = TRUE, 
##     n.knots = 10, kappa = 100)
## 
## 
##   Shared Gamma Frailty model parameter estimates  
##   using a Penalized Likelihood on the hazard function 
## 
##                   coef exp(coef) SE coef (H) SE coef (HIH)        z          p
## dukesC        0.294307  1.342196    0.160245      0.160245  1.83660 6.6269e-02
## dukesD        1.050124  2.858007    0.193975      0.193975  5.41370 6.1734e-08
## charlson1-2   0.450914  1.569746    0.257823      0.257823  1.74892 8.0304e-02
## charlson3     0.407673  1.503316    0.136585      0.136585  2.98475 2.8381e-03
## sexFemale    -0.536102  0.585024    0.138560      0.138560 -3.86910 1.0924e-04
## chemoTreated -0.209464  0.811019    0.142754      0.142754 -1.46731 1.4229e-01
## 
##            chisq df global p
## dukes    30.0693  2 2.95e-07
## charlson 10.9292  2 4.23e-03
## 
##     Frailty parameter, Theta: 0.660088 (SE (H): 0.141368 ) p = 1.5113e-06 
##  
##       penalized marginal log-likelihood = -3243.13
##       Convergence criteria: 
##       parameters = 4.87e-07 likelihood = 3.94e-06 gradient = 2.17e-11 
## 
##       LCV = the approximate likelihood cross-validation criterion
##             in the semi parametrical case     = 3.78877 
## 
##       n= 861
##       n events= 458  n groups= 403
##       number of iterations:  14 
## 
##       Exact number of knots used:  10 
##       Best smoothing parameter estimated by
##       an approximated Cross validation:  2622783, DoF:  11.00
## 
## Be patient. The program is computing ... 
## The program took 0.03 seconds
## Call:
## frailtyPenal(formula = Surv(time, status) ~ cluster(id) + sex + 
##     age, data = kidney, n.knots = 10, kappa = 1000)
## 
## 
##   Shared Gamma Frailty model parameter estimates  
##   using a Penalized Likelihood on the hazard function 
## 
##            coef exp(coef) SE coef (H) SE coef (HIH)         z          p
## sex -1.64968062  0.192111   0.4750662     0.4750662 -3.472528 0.00051558
## age  0.00597124  1.005989   0.0120986     0.0120986  0.493548 0.62163000
## 
##     Frailty parameter, Theta: 0.481138 (SE (H): 0.239159 ) p = 0.022121 
##  
##       penalized marginal log-likelihood = -323.01
##       Convergence criteria: 
##       parameters = 1.34e-06 likelihood = 1.04e-07 gradient = 1.03e-12 
## 
##       LCV = the approximate likelihood cross-validation criterion
##             in the semi parametrical case     = 4.44717 
## 
##       n= 76
##       n events= 58  n groups= 38
##       number of iterations:  23 
## 
##       Exact number of knots used:  10 
##       Value of the smoothing parameter:  1000, DoF:  12.00

1.3.2 Nested Model

## 
## Be patient. The program is computing ... 
## The program took 1.26 seconds
## Call:
## frailtyPenal(formula = Surv(t1, t2, event) ~ cluster(group) + 
##     subcluster(subgroup) + cov1 + cov2, data = dataNested, n.knots = 8, 
##     kappa = 10000)
## 
##       left truncated structure used
## 
##  Nested Frailty model parameter estimates using  
##  using a Penalized Likelihood on the hazard functions 
## 
##      coef                 exp(coef)           SE coef (H)        
## cov1 "-0.558268748017781" "0.572198827190525" "0.13322900468541" 
## cov2 "1.27396901855314"   "3.57501373663203"  "0.148808655939855"
##      SE coef (HIH)       z                   p           
## cov1 "0.133146544498935" "-4.19029436822713" "2.7859e-05"
## cov2 "0.14816220138121"  "8.56112173385969"  "< 1e-16"   
## 
##   Frailty parameters: 
##    alpha  (group effect):  0.373203  (SE(H): 0.107021 ) p = 0.00024405 
##    eta (subgroup effect):  0.109654  (SE(H): 0.0726317 ) p = 0.065557 
##  
##     penalized marginal log-likelihood = -1020.88
##     Convergence criteria: 
##     parameters = 0.000443 likelihood = 0.000782 gradient = 3.9e-07 
## 
##     LCV = the approximate likelihood cross-validation criterion
##           in the semi parametrical case     = 2.57851 
## 
##     n= 400
##     n events= 287  n groups= 20
##     number of iterations:  15 
##    Number of nodes for the Gauss-Laguerre quadrature:  20 
##     Exact number of knots used:  8 
##     Value of the smoothing parameter:  10000, DoF:  6.75

1.3.3 JOint Frailty model

dependence between recurrent events and a time to event data

## 
## Be patient. The program is computing ... 
## The program took 10.43 seconds
## Call:
## frailtyPenal(formula = Surv(time, event) ~ cluster(id) + dukes + 
##     sex + chemo + terminal(death), formula.terminalEvent = ~dukes + 
##     sex + chemo, data = readmission, jointGeneral = TRUE, n.knots = 8, 
##     kappa = c(9.55e+09, 1.41e+12))
## 
## 
##   General Joint gamma frailty model for recurrent and a terminal event processes 
##   using a Penalized Likelihood on the hazard function 
## 
## Recurrences:
## ------------- 
##              coef                 exp(coef)           SE coef (H)        
## dukesC       "0.412429890112696"  "1.51048363875758"  "0.163798227716201"
## dukesD       "1.47220451638146"   "4.35883367735597"  "0.194265279041877"
## sexFemale    "-0.52675550058018"  "0.590517799570723" "0.14177346252013" 
## chemoTreated "-0.154786774986889" "0.856597806090863" "0.145815012919338"
##              SE coef (HIH)       z                   p           
## dukesC       "0.163609079006325" "2.5179142403621"   "0.011805"  
## dukesD       "0.194112486424715" "7.57832034444047"  "3.4972e-14"
## sexFemale    "0.141742540582018" "-3.71547320081422" "0.00020282"
## chemoTreated "0.145658372142173" "-1.06152838372352" "0.28845"   
## 
##         chisq df global p
## dukes 63.7708  2 1.42e-14
## 
## Terminal event:
## ---------------- 
##              coef                 exp(coef)           SE coef (H)        
## dukesC       "1.59253789308488"   "4.91620992447146"  "0.344961906921893"
## dukesD       "3.99593644924482"   "54.3767378438515"  "0.360716188206594"
## sexFemale    "-0.291660551949043" "0.747022068125473" "0.234927978483215"
## chemoTreated "0.993711829977274"  "2.70124243958545"  "0.245168397230133"
##              SE coef (HIH)       z                  p           
## dukesC       "0.344728065565823" "4.6165616003667"  "3.9015e-06"
## dukesD       "0.360417942715141" "11.0777851948142" "< 1e-16"   
## sexFemale    "0.234904185839602" "-1.241489216534"  "0.21443"   
## chemoTreated "0.245014404567108" "4.05318075740612" "5.0526e-05"
## 
##       chisq              df  global p 
## dukes "140.668443624939" "2" "< 1e-16"
## 
##  Frailty parameters: 
##    theta (variance of u, association between recurrences and terminal event): 0.764919 (SE (H): 0.134567 ) p = 6.5683e-09 
##    eta (variance of v, intra-subject correlation): 0.00736244 (SE (H): 0.00051788 ) p = < 1e-16 
##  
## Penalized marginal log-likelihood = -3861.49
##    Convergence criteria: 
##    parameters = 0.000284 likelihood = 6.77e-05 gradient = 3.43e-08 
## 
## Likelihood Cross-Validation (LCV) criterion in the semi parametrical case:
##    approximate LCV = 4.51958 
## 
##    n observations= 861  n subjects= 403
##    n recurrent events= 458
##    n terminal events= 109
##    n censored events= 403
##    number of iterations:  10 
##    Number of nodes for the Gauss-Laguerre quadrature:  20 
## 
##    Exact number of knots used:  8 
##    Value of the smoothing parameters:  9.55e+09 1.41e+12

1.3.4 Additive frailty model

heterogeneity of the baseline risk heterogeneity of treatment effects across trials

## 
## Be patient. The program is computing ... 
## The program took 4.34 seconds
## Call:
## additivePenal(formula = Surv(t1, t2, event) ~ cluster(group) + 
##     var1 + var2 + slope(var1), data = dataAdditive, correlation = TRUE, 
##     cross.validation = TRUE, n.knots = 10, kappa = 1)
## 
##       left truncated structure used
## 
##   Additive gaussian frailty model parameter estimates  
##   using a Penalized Likelihood on the hazard function 
## 
##            coef exp(coef) SE coef (H) SE coef (HIH)         z       p
## var1 -0.0270881  0.973276   0.1063351     0.1056913 -0.254743 0.79892
## var2  0.0312746  1.031769   0.0880684     0.0879292  0.355118 0.72250
## 
##     Covariance (between the two frailty terms,  
##         the intercept and the slope): -0.0459305 (SE: 0.0502047 ) 
## 
##     Corresponding correlation between the two frailty terms : -0.93785 
##     Variance for random intercept: 0.105869 (SE (H): 0.0763206 ) 
##     Variance for random slope: 0.022655 (SE (H): 0.0463304 ) 
##  
##     penalized marginal log-likelihood = -2413.85
##     Convergence criteria: 
##     parameters = 1.84e-12 likelihood = 1.25e-09 gradient = 1.63e-09 
## 
##     LCV = the approximate likelihood cross-validation criterion
##     in the semi parametrical case     = 3.03398 
## 
##     n= 800
##     n events= 540  n groups= 8
##     number of iterations:  11 
## 
##     Exact number of knots used:  10 
##     Smoothing parameter estimated by Cross validation:  591492, DoF:  8.19