dmphillippo
commented
1 week ago Hi @rogula, this is currently the intended behaviour.
The relative_effects() function produces (population-average) conditional hazard ratios; when there is a model on the shape parameters with aux_regression you are right that this is not reflected in these conditional hazard ratios, as these are only based on the treatment effects and regression parameters on the rate linear predictor. I don't think that these conditional HRs are particularly interpretable for non-PH models. I also don't think that there is a way to combine these with the model on the shape to get time-varying conditional HRs. If you know of anywhere this has been done please say!
For non-PH models I would suggest using the marignal_effects() function, which will produce time-varying marginal hazard ratios. These are probably what you are after - and they do have a proper interpretation, as ratios of the estimated marginal hazard functions on each treatment.
Hello Dr. Phillippo,
Thanks again for your help with addressing issue #43 (where aux_regression terms were not always included in the predict function). I'd like to report a similar issue with relative_effects, which seems to ignore the ~.trt coefficient in aux_regression and only returns a single set of (non-time varying) HRs for a non-PH survival model.
Similar to before, by adapting your R code (simulated_parametric_survival.R), we can reproduce this error with the following changes:
First, update the simulated data to produce non-proportional hazard (PH) datasets: survdat_AB <- simsurv("weibull", lambdas = weib_sim_pars["scaleAB"], gammas = weib_sim_pars["shapeAB"], x = covdat_AB, betas = betas, maxt = 1, tde = c(trtB = 5)) %>% mutate(status = if_else(runif(nAB) <= cens_rate, 0L, status))
survdat_AC <- simsurv("weibull", lambdas = weib_sim_pars["scaleAC"], gammas = weib_sim_pars["shapeAC"], x = covdat_AC, betas = betas, maxt = 1, tde = c(trtC = 1)) %>% mutate(status = if_else(runif(nAC) <= cens_rate, 0L, status))
Then, when fitting the following model (using the IPD dataset), the model coefficients correctly capture the treatment effect on the auxiliary parameters, but the relative_effects function does not appear to include it and only one set of HRs is returned:
weib_IPD_nph1 <- nma(sim_net_IPD, regression = ~(x1 + x2 + x3)*.trt, likelihood = "weibull", prior_intercept = normal(0, 100), prior_reg = normal(0, 100), prior_trt = normal(0, 100), prior_aux = half_normal(10), aux_regression = ~.trt, QR = TRUE)
tab = relative_effects(weib_IPD_nph1, trt_ref = "A")$summary;