Jointlcmm error in chol.default(varcov)

[wy2020hku]

Dear Cécile and Viviane,

I rencently met an error using the jointlcmm function. My data have 100,000 individuals with a median of 9 annual measurements. My code is:

mod <- gridsearch(rep=8, maxiter=200, minit=initmod, cl=16, Jointlcmm(
  fixed = bmi ~ ns(time, knots = c(4,8)),
  mixture = ~ ns(time, knots = c(4,8)),
  random = ~ time,
  subject = "patient_pssn",
  ng = 2,
  data = data,
  survival = Surv(first_cvd_fu_period2, first_cvd_after) ~ hba1c_bl + dm_duration + 
    age +
    smoking_bl +
    systolic_bp_bl +
    diastolic_bp_bl +
    ldl_c_bl +
    bmi_bl +
    egfr_bl +
    charlson_index + ht_b4 +  
insulin_b4 +
    dm_drug_b4 +
    ht_drug_b4 +
    lipid_drug_b4,
  logscale = FALSE,
  B = initmod))

And the error is:

Error in chol.default(varcov) : 
  the leading minor of order 1 is not positive definite
Calls: gridsearch -> do.call -> Jointlcmm -> chol -> chol.default
Execution halted

I got this error in one subset of my dataset but didn't in other subsets, however, I have tried to modify this subset (changing the cut off of different subset) but I always get this error in some of the datasets. Such problem didn't appear in my other models with different dependent variable in the fixed model but all other parameters similar. May I know what this error means in jointlcmm model or if you have any experience to deal with it?

Moreover, May I ask more about a previous problem #175 : why can we consider the model with a stabilization of score test results as satisfying the conditional independence? Although the decline of score test statistics become slower or unobservable, it may still be possible to reach insignificant bound when we have many classes (like 200 or 1000 classes). I noticed that this point was also mentioned in the vignettes for JointLcmm and sincerely hope to know the reason. Thank you in advance for your help!

Best, Yuan

[VivianePhilipps]

Hi,

this means that the variance-covariance matrix of the random effects provided in initmod is not positive definite. Does the one class model converged correctly with meaningful estimates?

Viviane

[wy2020hku]

Hi Viviane,

The one class model details are:

Joint latent class model for quantitative outcome and competing risks fitted by maximum likelihood method

Jointlcmm(fixed = egfr ~ ns(time, knots = c(4, 8)), random = ~time, subject = "patient_pssn", ng = 1, survival = Surv(first_cvd_fu_period2, first_cvd_after) ~ hba1c_bl + dm_duration + age + sex + smoking_bl + systolic_bp_bl + diastolic_bp_bl + ldl_c_bl + bmi_bl + egfr_bl + charlson_index + ht_b4 + insulin_b4 + dm_drug_b4 + ht_drug_b4 + lipid_drug_b4, data = data, maxiter = 5000, logscale = FALSE, nproc = 16)

Statistical Model: Dataset: data Number of subjects: 30319 Number of observations: 306917 Number of latent classes: 1 Number of parameters: 26
Event 1: Number of events: 5702 Weibull baseline risk function

Iteration process: Convergence criteria satisfied Number of iterations: 279 Convergence criteria: parameters= 1.6e-06 : likelihood= 4.4e-05 : second derivatives= 1.1e-05

Goodness-of-fit statistics: maximum log-likelihood: -1113492.85
AIC: 2227037.69
BIC: 2227254

Maximum Likelihood Estimates:

Parameters in the proportional hazard model:

                          coef      Se     Wald p-value

event1 +/-sqrt(Weibull1) 0.02285 0.00295 7.747 0.00000 event1 +/-sqrt(Weibull2) 1.07529 0.00691 155.597 0.00000 hba1c_bl 0.06440 0.00974 6.613 0.00000 dm_duration 0.01847 0.00223 8.294 0.00000 age 0.03081 0.00355 8.689 0.00000 sex 0.43776 0.02924 14.969 0.00000 smoking_bl 0.33204 0.04463 7.441 0.00000 systolic_bp_bl 0.00610 0.00105 5.797 0.00000 diastolic_bp_bl -0.00585 0.00181 -3.233 0.00123 ldl_c_bl 0.02821 0.01319 2.139 0.03242 bmi_bl 0.03051 0.00347 8.785 0.00000 egfr_bl -0.00415 0.00168 -2.473 0.01342 charlson_index 0.08905 0.01814 4.910 0.00000 ht_b4 0.18721 0.04974 3.763 0.00017 insulin_b4 0.26637 0.05998 4.441 0.00001 dm_drug_b4 0.04311 0.03404 1.266 0.20537 ht_drug_b4 0.20893 0.05243 3.985 0.00007 lipid_drug_b4 0.14016 0.03254 4.308 0.00002

Fixed effects in the longitudinal model:

           coef      Se     Wald p-value

intercept 77.93431 0.07298 1067.900 0.00000 ns(...)1 -9.68677 0.08114 -119.378 0.00000 ns(...)2 -13.58961 0.15100 -89.997 0.00000 ns(...)3 -13.65894 0.10670 -128.015 0.00000

Variance-covariance matrix of the random-effects: intercept time intercept 131.47566
time 0.75764 2.39834

                         coef      Se

Residual standard error 6.36045 0.00908

I didn't observe any problem in this stage, can you find anything special that might be related to the error in the following models?

Best, Yuan

[VivianePhilipps]

Maybe it is related to the parameters' variance. Could you please show : cbind(estimates(initmod), sqrt(diag(vcov(initmod))))

[wy2020hku]

Thanks Viviane, the results are:

> cbind(estimates(initmod), sqrt(diag(vcov(initmod))))
                                    [,1]        [,2]
event1 +/-sqrt(Weibull1)     0.022861326 0.003025871
event1 +/-sqrt(Weibull2)     1.075279709 0.006909724
hba1c_bl                     0.064404761 0.009324639
dm_duration                  0.018464348 0.002440469
age                          0.030821261 0.003452958
sex                          0.437680037 0.029820493
smoking_bl                   0.332044985 0.044580043
systolic_bp_bl               0.006100327 0.001043560
diastolic_bp_bl             -0.005838289 0.001885677
ldl_c_bl                     0.028023407 0.015737764
bmi_bl                       0.030507497 0.003463963
egfr_bl                     -0.004163676 0.001553711
charlson_index               0.088951048 0.015773130
ht_b4                        0.186924632 0.051585194
insulin_b4                   0.266390575 0.058951345
dm_drug_b4                   0.042543993 0.061290931
ht_drug_b4                   0.209112871 0.053736952
lipid_drug_b4                0.140004042 0.034089621
intercept                   77.934308247 0.072978982
ns(time, knots = c(4, 8))1  -9.686673650 0.081185856
ns(time, knots = c(4, 8))2 -13.589387049 0.151068233
ns(time, knots = c(4, 8))3 -13.658781053 0.106745470
cholesky 1                  11.466305367 0.052364446
cholesky 2                   0.066045491 0.010631928
cholesky 3                   1.547251910 0.008201156
stderr                       6.360448470 0.009076073

Best, Yuan