Negative LCPO

[Felipe Barletta]

Hello INLA group.

I fitted a joint model (Binomial and weibull survival model) with some association structures and some log Conditional predictive ordinates (LCPO) are negative, like table below:

 I expected to find LCPO positive and close values between DIC and WAIC. I have been surprised to find very large differences in several models. 

Someone has idea where is the mistake? For example, below the code to model III is:

m3 <- INLAjoint::joint(
                     formLong = list(Y2 ~ age + weekend + gender + (1  | ID_hosp) #GLMM
                    ),
    formSurv = inla.surv(time = I(los/365.25), event = Y1) ~ gender + age +  weekend, # survival
    dataLong = data,
    dataSurv = data,
    id = "ID_hosp",
    timeVar = "",
    corLong = FALSE, 
    family = "binomial",
    basRisk = "weibullsurv",
    assoc = c("SRE"),
    control = list(int.strategy="eb",dic = TRUE, cpo = TRUE)
)

[Finn Lindgren]

LCPO is log(CPO), so it’s perfectly normal for it to be negative.

Finn

On 13 Dec 2024, at 15:23, Felipe Barletta <felipe.e...@gmail.com> wrote:



Hello INLA group.

I fitted a joint model (Binomial and weibull survival model) with some association structures and some log Conditional predictive ordinates (LCPO) are negative, like table below:

<image.png>

 I expected to find LCPO positive and close values between DIC and WAIC. I have been surprised to find very large differences in several models. 

Someone has idea where is the mistake? For example, below the code to model III is:

m3 <- INLAjoint::joint(
                     formLong = list(Y2 ~ age + weekend + gender + (1  | ID_hosp) #GLMM
                    ),
    formSurv = inla.surv(time = I(los/365.25), event = Y1) ~ gender + age +  weekend, # survival
    dataLong = data,
    dataSurv = data,
    id = "ID_hosp",
    timeVar = "",
    corLong = FALSE, 
    family = "binomial",
    basRisk = "weibullsurv",
    assoc = c("SRE"),
    control = list(int.strategy="eb",dic = TRUE, cpo = TRUE)
)

--
You received this message because you are subscribed to the Google Groups "R-inla discussion group" group.
To unsubscribe from this group and stop receiving emails from it, send an email to r-inla-discussion...@googlegroups.com.
To view this discussion, visit https://groups.google.com/d/msgid/r-inla-discussion-group/CAJz4ZjEvQpNSoiTbZda5yxOqUp0iDkZFFW-PTNK5pe_p-0uv_w%40mail.gmail.com.

[Felipe Barletta]

Is it not -sum(log(CPO))? 

[Finn Lindgren]

As far as I know, no. But one can define scores in both direction.

check by comparing with the calculated individual cpo values…

However, for discrete data, lcpo should either all be negative or all positive, depending on the definition of the aggregated score. For continuous data it’s the log of a density, so can be both positive and negative with either definition…

Finn

On 13 Dec 2024, at 16:03, Felipe Barletta <felipe.e...@gmail.com> wrote:



[Felipe Barletta]

Thanks Finn for your answer. I got now. Actually I calculated -sum(log(CPO)) to each model. 

But I will check like your suggestion. 

Thank you so much! 

[Finn Lindgren]

Then it depends on whether your observation model is discrete or  continuous.

Finn

[Denis Rustand]

Hi Felipe,

To complete Finn's answer and in the specific context of your joint model, in order to account for the time-dependent linear combination of random effects in the survival submodel, we use some fake data that computes this linear combination to share. By default, INLA will include this fake data in the computation of dic, cpo and waic. You want to remove this part by only looking at values that corresponds to your longitudinal and survival submodels. For example for cpo:

CPO_l <- sum(m3$cpo$cpo[which(!is.na(m3$.args$data$Yjoint[[1]]))]) # cpo longitudinal part
CPO_s <- sum(m3$cpo$cpo[which(!is.na(m3$.args$data$Yjoint[[3]]))]) # cpo survival part
CPO_t <- CPO_l + CPO_s # cpo total

When you use "SRE_ind" instead of "SRE", the random effects are just copied and scaled (without the use of fake data) but this cannot be used to share anything time-dependent.

Best,
Denis