Component-dose network meta-analysis with dose-dependent interactions

This vignette shows the basic workflow for cdtmbnma. The example is the sacubitril and valsartan sitting systolic blood-pressure dose plane. Each row is one arm. d_sac and d_val are component dose columns; zero means absent.

library(cdtmbnma)
sv <- sacval_example()
head(sv)
#>   study                    arm d_sac d_val   n      y    sd       se
#> 1  Izzo           Izzo_placebo     0     0  57  -7.00 22.40 2.966952
#> 2  Izzo      Izzo_valsartan320     0   320 143 -16.10 22.20 1.856457
#> 3  Izzo  Izzo_freeSac50_Val320    50   320 132 -19.35 14.62 1.272508
#> 4  Izzo Izzo_freeSac100_Val320   100   320 141 -21.26 14.52 1.222805
#> 5  Izzo Izzo_freeSac200_Val320   200   320 144 -23.64 14.81 1.234167
#> 6  Izzo Izzo_freeSac400_Val320   400   320 143 -20.95 14.73 1.231784
#>                               sd_source
#> 1            imputed (per-study median)
#> 2                    reported/NMA-table
#> 3 derived from LS-mean SE (sensitivity)
#> 4 derived from LS-mean SE (sensitivity)
#> 5 derived from LS-mean SE (sensitivity)
#> 6 derived from LS-mean SE (sensitivity)

Building the single-timepoint design

cdt_data() turns the arm-level table into the model design. Studies with an all-zero arm use that arm as the within-study reference. Studies without an all-zero arm use their lowest-total-dose arm as the baseline and issue a warning.

d <- cdt_data(
  sv,
  study = "study",
  components = c("d_sac", "d_val"),
  outcome = "continuous",
  y = "y",
  se = "se",
  dstar = c(d_sac = 200, d_val = 320)
)
#> Warning: Study 'Cheung' has no all-zero arm; using its lowest-total-dose arm as
#> baseline.
#> Warning: Study 'Huo' has no all-zero arm; using its lowest-total-dose arm as
#> baseline.
#> Warning: Study 'Rakugi' has no all-zero arm; using its lowest-total-dose arm as
#> baseline.
#> Warning: Study 'Ruilope' has no all-zero arm; using its lowest-total-dose arm
#> as baseline.
d
#> <cdt_data>
#>   arms:       18 across 5 studies
#>   components: d_sac, d_val
#>   outcome:    continuous
#>   interaction pairs: d_sac x d_val

Fitting the bilinear interaction model

cdt_fit() compiles the Stan model on first use and samples. The bilinear interaction parameter is interpreted at the dstar corner: here, sacubitril 200 mg and valsartan 320 mg.

fit <- cdt_fit(
  d,
  interaction = "bilinear",
  chains = 4,
  iter_warmup = 1000,
  iter_sampling = 1000,
  seed = 1
)
#> Running MCMC with 4 parallel chains...
#> Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
#> Chain 1 Exception: normal_lpdf: Location parameter[2] is -inf, but must be finite! (in '/tmp/RtmpKcYyi0/model-23f73c17463052.stan', line 97, column 4 to column 23)
#> Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
#> Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
#> Chain 1
#> Chain 1 finished in 0.5 seconds.
#> Chain 2 finished in 0.6 seconds.
#> Chain 3 finished in 0.5 seconds.
#> Chain 4 finished in 0.5 seconds.
#> 
#> All 4 chains finished successfully.
#> Mean chain execution time: 0.5 seconds.
#> Total execution time: 0.7 seconds.
summary(fit)
#>                     variable       mean         sd       2.5%      97.5%
#> 1                emax: d_sac -10.958955  4.4289041 -20.981215  -3.790150
#> 2                emax: d_val -16.797890  4.5391811 -25.698727  -7.735944
#> 3 interaction: d_sac x d_val   1.682457  2.0916998  -2.358932   5.732784
#> 4                      omega   2.685903  0.7441266   1.455292   4.363215
#> 5                ED50: d_sac  72.783504 56.4619799  10.421806 211.643737
#> 6                ED50: d_val 134.180238 76.7936804  31.804358 322.117577
#>        rhat ess_bulk
#> 1 1.0013255 2301.701
#> 2 1.0014045 1881.681
#> 3 1.0013801 2279.402
#> 4 0.9998781 1710.026
#> 5 1.0007031 3259.655
#> 6 1.0016402 2740.434

The marginal component curves vary one component at a time while holding the other at zero.

plot(fit)
#> Warning: Dropping 'draws_df' class as required metadata was removed.
#> Warning: Dropping 'draws_df' class as required metadata was removed.

Prediction at an unobserved dose pair

predict() returns the relative effect against the all-zero reference at any component-dose combination, including dose pairs not directly observed in a trial.

predict(fit, newdata = data.frame(d_sac = 100, d_val = 160))
#> Warning: Dropping 'draws_df' class as required metadata was removed.
#> Warning: Dropping 'draws_df' class as required metadata was removed.
#> Warning: Dropping 'draws_df' class as required metadata was removed.
#>   d_sac d_val      mean       sd      2.5%     97.5%
#> 1   100   160 -15.40289 2.383413 -19.99665 -10.81706

Additive benchmark

The additive model is fitted by setting interaction = "none".

fit_add <- cdt_fit(d, interaction = "none")

If the loo package is installed, the pointwise log likelihood stored in Stan’s generated quantities can be used for information-criterion comparisons.

loo::loo_compare(
  loo::loo(fit$fit$draws("log_lik")),
  loo::loo(fit_add$fit$draws("log_lik"))
)

Stage-one longitudinal model

The package also contains a narrower two-component longitudinal interface. It is useful when repeated arm means over time are available and an exponential approach-to-asymptote is reasonable. The example below builds a small design but does not fit it during vignette building.

long_dat <- data.frame(
  study = rep("trial1", 6),
  arm = rep(c("placebo", "A", "AB"), each = 2),
  week = rep(c(4, 8), 3),
  dose_A = rep(c(0, 10, 10), each = 2),
  dose_B = rep(c(0, 0, 20), each = 2),
  y = c(0, 0, -1, -2, -2, -3),
  se = rep(1, 6)
)

time_design <- cdt_time_data(
  long_dat,
  study = "study",
  arm = "arm",
  time = "week",
  components = c("dose_A", "dose_B"),
  y = "y",
  se = "se",
  dstar = c(dose_A = 10, dose_B = 20)
)
time_design
#> <cdt_time_data>
#>   observations: 6
#>   arms:         3 across 1 studies
#>   components:   dose_A, dose_B
#>   model:        two-component exponential time-course with bilinear interaction
fit_time <- cdt_time_fit(time_design)
predict(fit_time, data.frame(dose_A = 10, dose_B = 20))