stan_jm joint model with non-linear association between biomarker value and hazard

[JonNDCN]

Summary:

Can non-linear functions of the biomarker value / slope be specified? For example if a biomarker was associated with an increased hazard at high or low values

Reproducible Steps:

If logbili increased the hazard of death if it was very low or very high, how might one alter the code below to allow for this?

mod1 <- stan_jm(formulaLong = logBili ~ sex + trt + year + (year | id), 
                dataLong = pbcLong,
                formulaEvent = survival::Surv(futimeYears, death) ~ sex + trt, 
                dataEvent = pbcSurv,
                time_var = "year",
                chains = 1, refresh = 2000, seed = 12345)

RStanARM Version:

2.21.1

R Version:

4.0.2

Many thanks

[sambrilleman]

Hi @JonNDCN, sorry for the delay in getting back to you.

So there is no explicit association structure for non-linear effects, but you can potentially do something very crude using a quadratic effect of the biomarker. There isn't a "quadratic" association structure, but there is one that is designed to make interactions between different biomarkers in a multivariate joint model. But even with a single biomarker model you can use that association structure to "trick" the model into interacting the biomarker with itself, thereby forming a quadratic effect. This would look something like:

mod1 <- stan_jm(formulaLong = logBili ~ sex + trt + year + (year | id), 
                dataLong = pbcLong,
                formulaEvent = survival::Surv(futimeYears, death) ~ sex + trt, 
                dataEvent = pbcSurv,
                assoc = c("etavalue", "etavalue_etavalue(1)"),
                time_var = "year",
                chains = 1, refresh = 2000, seed = 12345)

There are now two terms forming the association structure: i) the usual etavalue (current value of the biomarker) and ii) etavalue interacted with itself (i.e. a quadratic effect of etavalue).

Note that you would probably want to centre the biomarker data before fitting that type of model I guess.

There isn't really any other non-linear association structures I can think of. You can however interact etavalue with covariates. So if you have some external covariate data that might help identify those "high" or "low" biomarker individuals then you could have etavalue only have an effect (or have a different effect) in certain individuals by interacting it with dummy covariates or something -- see the etavalue_data() association structure in the docs.

Hope that helps a bit? Even if not exactly what you were hoping for (e.g. no cubic spline effects of the biomarker or anything like that unfortunately).

[JonNDCN]

I wouldn’t consider that a delay!

Thank you for your very helpful answer.

Do you imagine splines or fractional polynomial functions of the biomarker may be supported in the future?

Thanks again for your help.

[sambrilleman]

Hi @JonNDCN.

Sadly, I doubt that an addition like that will happen.

I don't have much time to work on stan_jm anymore aside from bug fixes etc, no real time to add larger new features. It would be great if someone from the Stan community was able to contribute an addition like this, but to be honest, whilst the codebase for the new stan_surv function is simpler and cleaner and therefore I foresee the possibility of users contributing features, the stan_jm codebase is much more complicated, so it is harder for someone to get across the details and contribute an addition like this. Although, never say never! 😆

Splines or fractional polynomials would need to be implemented in the actual .stan file, and calculated using a custom function definition in the Stan code (not in R) since the basis terms would need to be calculated at each iteration of the sampling (based on the estimated current value of the biomarker). So it wouldn't be as easy as just calling the splines2 package or something in R, I don't think. It's probably not too crazy, but would be a little bit of work.