blrm_exnex {OncoBayes2}R Documentation

Bayesian Logistic Regression Model for N-compounds with EXNEX

Description

Bayesian Logistic Regression Model (BLRM) for N compounds using EXchangability and NonEXchangability (EXNEX) modeling.

Usage

blrm_exnex(
  formula,
  data,
  prior_EX_mu_mean_comp,
  prior_EX_mu_sd_comp,
  prior_EX_tau_mean_comp,
  prior_EX_tau_sd_comp,
  prior_EX_corr_eta_comp,
  prior_EX_mu_mean_inter,
  prior_EX_mu_sd_inter,
  prior_EX_tau_mean_inter,
  prior_EX_tau_sd_inter,
  prior_EX_corr_eta_inter,
  prior_is_EXNEX_inter,
  prior_is_EXNEX_comp,
  prior_EX_prob_comp,
  prior_EX_prob_inter,
  prior_NEX_mu_mean_comp,
  prior_NEX_mu_sd_comp,
  prior_NEX_mu_mean_inter,
  prior_NEX_mu_sd_inter,
  prior_tau_dist,
  iter = getOption("OncoBayes2.MC.iter", 2000),
  warmup = getOption("OncoBayes2.MC.warmup", 1000),
  save_warmup = getOption("OncoBayes2.MC.save_warmup", TRUE),
  thin = getOption("OncoBayes2.MC.thin", 1),
  init = getOption("OncoBayes2.MC.init", 0.5),
  chains = getOption("OncoBayes2.MC.chains", 4),
  cores = getOption("mc.cores", 1L),
  control = getOption("OncoBayes2.MC.control", list()),
  backend = getOption("OncoBayes2.MC.backend", "rstan"),
  prior_PD = FALSE,
  verbose = FALSE
)

## S3 method for class 'blrmfit'
print(x, ..., prob = 0.95, digits = 2)

Arguments

formula

the model formula describing the linear predictors of the model. The lhs of the formula is a two-column matrix which are the number of occured events and the number of times no event occured. The rhs of the formula defines the linear predictors for the marginal models for each drug component, then the interaction model and at last the grouping and optional stratum factors of the models. These elements of the formula are separated by a vertical bar. The marginal models must follow a intercept and slope form while the interaction model must not include an interaction term. See the examples below for an example instantiation.

data

optional data frame containing the variables of the model. If not found in data, the variables are taken from environment(formula).

prior_EX_mu_mean_comp, prior_EX_mu_sd_comp

Mean and sd for the prior on the mean parameters μi=(μαi,μβi)\boldsymbol\mu_i = (\mu_{\alpha i}, \mu_{\beta i}) of each component. Two column matrix (intercept, log-slope) with one row per component.

prior_EX_tau_mean_comp, prior_EX_tau_sd_comp

Prior mean and sd for heterogeniety parameter τsi=(ταsi,τβsi)\boldsymbol\tau_{si} = (\tau_{\alpha s i}, \tau_{\beta s i}) of each stratum. If no differential discounting is required (i.e. if there is only one stratum s=1s = 1), then it is a two-column matrix (intercept, log-slope) with one row per component. Otherwise it is a three-dimensional array whose first dimension indexes the strata, second dimension indexes the components, and third dimension of length two for (intercept, log-slope).

prior_EX_corr_eta_comp

Prior LKJ correlation parameter for each component given as numeric vector. If missing, then a 1 is assumed corresponding to a marginal uniform prior of the correlation.

prior_EX_mu_mean_inter, prior_EX_mu_sd_inter

Prior mean and sd for population mean parameters μηk\mu_{\eta k} of each interaction parameter. Vector of length equal to the number of interactions.

prior_EX_tau_mean_inter, prior_EX_tau_sd_inter

Prior mean and sd for heterogeniety parameter τηsk\tau_{\eta s k} of each stratum. Matrix with one column per interaction and one row per stratum.

prior_EX_corr_eta_inter

Prior LKJ correlation parameter for interaction given as numeric. If missing, then a 1 is assumed corresponding to a marginal uniform prior of the correlations.

prior_is_EXNEX_inter

Defines if non-exchangability is admitted for a given interaction parameter. Logical vector of length equal to the number of interactions. If missing FALSE is assumed for all interactions.

prior_is_EXNEX_comp

Defines if non-exchangability is admitted for a given component. Logical vector of length equal to the number of components. If missing TRUE is assumed for all components.

prior_EX_prob_comp

Prior probability pijp_{ij} for exchangability of each component per group. Matrix with one column per component and one row per group. Values must lie in [0-1] range.

prior_EX_prob_inter

Prior probability pηkjp_{\eta k j} for exchangability of each interaction per group. Matrix with one column per interaction and one row per group. Values must lie in [0-1] range.

prior_NEX_mu_mean_comp, prior_NEX_mu_sd_comp

Prior mean mij\boldsymbol m_{ij} and sd sij=\mboxdiag(Sij)\boldsymbol s_{ij} = \mbox{diag}(\boldsymbol S_{ij}) of each component for non-exchangable case. Two column matrix (intercept, log-slope) with one row per component. If missing set to the same prior as given for the EX part. It is required that the specification be the same across groups j.

prior_NEX_mu_mean_inter, prior_NEX_mu_sd_inter

Prior mean mηkjm_{\eta k j} and sd sηkjs_{\eta k j} for each interaction parameter for non-exchangable case. Vector of length equal to the number of interactions. If missing set to the same prior as given for the EX part.

prior_tau_dist

Defines the distribution used for heterogeniety parameters. Choices are 0=fixed to it's mean, 1=log-normal, 2=truncated normal.

iter

number of iterations (including warmup).

warmup

number of warmup iterations.

save_warmup

save warmup samples (TRUE / FALSE). Only if set to TRUE, then all random variables are saved in the posterior. This substantially increases the storage needs of the posterior.

thin

period of saving samples.

init

positive number to specify uniform range on unconstrained space for random initialization. See stan.

chains

number of Markov chains.

cores

number of cores for parallel sampling of chains.

control

additional sampler parameters for NuTS algorithm.

backend

sets Stan backend to be used. Possible choices are "rstan" (default) or "cmdstanr".

prior_PD

Logical flag (defaults to FALSE) indicating if to sample the prior predictive distribution instead of conditioning on the data.

verbose

Logical flag (defaults to FALSE) controlling if additional output like stan progress is reported.

x

blrmfit object to print

...

not used in this function

prob

central probability mass to report, i.e. the quantiles 0.5-prob/2 and 0.5+prob/2 are displayed. Multiple central widths can be specified.

digits

number of digits to show

Details

blrm_exnex is a flexible function for Bayesian meta-analytic modeling of binomial count data. In particular, it is designed to model counts of the number of observed dose limiting toxicities (DLTs) by dose, for guiding dose-escalation studies in Oncology. To accommodate dose escalation over more than one agent, the dose may consist of combinations of study drugs, with any number of treatment components.

In the simplest case, the aim is to model the probability π\pi that a patient experiences a DLT, by complementing the binomial likelihood with a monotone logistic regression

\mboxlogitπ(d)=logα+βt(d),\mbox{logit}\,\pi(d) = \log\,\alpha + \beta \, t(d),

where β>0\beta > 0. Most typically, dd represents the dose, and t(d)t(d) is an appropriate transformation, such as t(d)=log(d/d)t(d) = \log (d \big / d^*). A joint prior on θ=(logα,logβ)\boldsymbol \theta = (\log\,\alpha, \log\,\beta) completes the model and ensures monotonicity β>0\beta > 0.

Many extensions are possible. The function supports general combination regimens, and also provides framework for Bayesian meta-analysis of dose-toxicity data from multiple historical and concurrent sources.

For an example of a single-agent trial refer to example-single-agent.

Value

The function returns a S3 object of type blrmfit.

Methods (by generic)

Combination of two treatments

For a combination of two treatment components, the basic modeling framework is that the DLT rate π(d1,d2)\pi(d_1,d_2) is comprised of (1) a "no-interaction" baseline model π~(d1,d2)\tilde \pi(d_1,d_2) driven by the single-agent toxicity of each component, and (2) optional interaction terms γ(d1,d2)\gamma(d_1,d_2) representing synergy or antagonism between the drugs. On the log-odds scale,

\mboxlogitπ(d1,d2)=\mboxlogitπ~(d1,d2)+ηγ(d1,d2).\mbox{logit} \,\pi(d_1,d_2) = \mbox{logit} \, \tilde \pi(d_1,d_2) + \eta \, \gamma(d_1,d_2).

The "no interaction" part π~(d1,d2)\tilde \pi(d_1,d_2) represents the probability of a DLT triggered by either treatment component acting independently. That is,

π~(d1,d2)=1(1π1(d1))(1π2(d2)). \tilde \pi(d_1,d_2) = 1- (1 - \pi_1(d_1))(1 - \pi_2(d_2)).

In simple terms, P(no DLT for combination) = P(no DLT for drug 1) * P(no DLT from drug 2). To complete this part, the treatment components can then be modeled with monotone logistic regressions as before.

\mboxlogitπi(di)=logαi+βit(di),\mbox{logit} \, \pi_i(d_i) = \log\, \alpha_i + \beta_i \, t(d_i),

where t(di)t(d_i) is a monotone transformation of the doses, such as t(di)=log(di/di)t(d_i) = \log (d_i \big / d_i^*).

The inclusion of an interaction term γ(d1,d2)\gamma(d_1,d_2) allows DLT rates above or below the "no-interaction" rate. The magnitude of the interaction term may also be made dependent on the doses (or other covariates) through regression. As an example, one could let

γ(d1,d2)=d1d1d1d2.\gamma(d_1, d_2) = \frac{d_1}{d_1^*} \frac{d_1}{d_2^*}.

The specific functional form is specified in the usual notation for a design matrix. The interaction model must respect the constraint that whenever any dose approaches zero, then the interaction term must vanish as well. Therefore, the interaction model must not include an intercept term which would violate this consistency requirement. A dual combination example can be found in example-combo2.

General combinations

The model is extended to general combination treatments consisting of NN components by expressing the probability π\pi on the logit scale as

\mboxlogitπ(d1,,dN)=\mboxlogit(1i=1N(1πi(di)))+k=1Kηkγk(d1,,dN), \mbox{logit} \, \pi(d_1,\ldots,d_N) = \mbox{logit} \Bigl( 1 - \prod_{i = 1}^N ( 1 - \pi_i(d_i) ) \Bigr) + \sum_{k=1}^K \eta_k \, \gamma_k(d_1,\ldots,d_N),

Multiple drug-drug interactions among the NN components are now possible, and are represented through the KK interaction terms γk\gamma_k.

Regression models can be again be specified for each πi\pi_i and γk\gamma_k, such as

\mboxlogitπi(di)=logαi+βit(di) \mbox{logit}\, \pi_i(d_i) = \log\, \alpha_i + \beta_i \, t(d_i)

Interactions for some subset I(k){1,,N}I(k) \subset \{1,\ldots,N \} of the treatment components can be modeled with regression as well, for example on products of doses,

γk(d1,,dN)=iI(k)didi. \gamma_k(d_1,\ldots,d_N) = \prod_{i \in I(k)} \frac{d_i}{d_i^*}.

For example, I(k)={1,2,3}I(k) = \{1,2,3\} results in the three-way interaction term

d1d1d2d2d3d3 \frac{d_1}{d_1^*} \frac{d_2}{d_2^*} \frac{d_3}{d_3^*}

for drugs 1, 2, and 3.

For a triple combination example please refer to example-combo3.

Meta-analytic framework

Information on the toxicity of a drug may be available from multiple studies or sources. Furthermore, one may wish to stratify observations within a single study (for example into groups of patients corresponding to different geographic regions, or multiple dosing dose_info corresponding to different schedules).

blrm_exnex provides tools for robust Bayesian hierarchical modeling to jointly model data from multiple sources. An additional index j=1,,Jj=1, \ldots, J on the parameters and observations denotes the JJ groups. The resulting model allows the DLT rate to differ across the groups. The general NN-component model becomes

\mboxlogitπj(d1,,dN)=\mboxlogit(1i=1N(1πij(di)))+k=1Kηkjγk(d1,,dN), \mbox{logit} \, \pi_j(d_1,\ldots,d_N) = \mbox{logit} \Bigl( 1 - \prod_{i = 1}^N ( 1 - \pi_{ij}(d_i) ) \Bigr) + \sum_{k=1}^K \eta_{kj} \, \gamma_{k}(d_1,\ldots,d_N),

for groups j=1,,Jj = 1,\ldots,J. The component toxicities πij\pi_{ij} and interaction terms γk\gamma_{k} are modelled, as before, through regression. For example, πij\pi_{ij} could be a logistic regression on t(di)=log(di/di)t(d_i) = \log(d_i/d_i^*) with intercept and log-slope θij\boldsymbol \theta_{ij}, and γk\gamma_{k} regressed with coefficient ηkj\eta_{kj} on a product iI(k)(di/di)\prod_{i\in I(k)} (d_i/d_i^*) for some subset I(k)I(k) of components.

Thus, for j=1,,Jj=1,\ldots,J, we now have group-specific parameters θij=(logαij,logβij)\boldsymbol\theta_{ij} = (\log\, \alpha_{ij}, \log\, \beta_{ij}) and νj=(η1j,,ηKj)\boldsymbol\nu_{j} = (\eta_{1j}, \ldots, \eta_{Kj}) for each component i=1,,Ni=1,\ldots,N and interaction k=1,,Kk=1,\ldots,K.

The structure of the prior on (θi1,,θiJ)(\boldsymbol\theta_{i1},\ldots,\boldsymbol\theta_{iJ}) and (ν1,,νJ)(\boldsymbol\nu_{1}, \ldots, \boldsymbol\nu_{J}) determines how much information will be shared across groups jj. Several modeling choices are available in the function.

References

Neuenschwander, B., Roychoudhury, S., & Schmidli, H. (2016). On the use of co-data in clinical trials. Statistics in Biopharmaceutical Research, 8(3), 345-354.

Neuenschwander, B., Wandel, S., Roychoudhury, S., & Bailey, S. (2016). Robust exchangeability designs for early phase clinical trials with multiple strata. Pharmaceutical statistics, 15(2), 123-134.

Neuenschwander, B., Branson, M., & Gsponer, T. (2008). Critical aspects of the Bayesian approach to phase I cancer trials. Statistics in medicine, 27(13), 2420-2439.

Neuenschwander, B., Matano, A., Tang, Z., Roychoudhury, S., Wandel, S. Bailey, Stuart. (2014). A Bayesian Industry Approach to Phase I Combination Trials in Oncology. In Statistical methods in drug combination studies (Vol. 69). CRC Press.

Examples

## Setting up dummy sampling for fast execution of example
## Please use 4 chains and 100x more warmup & iter in practice
.user_mc_options <- options(OncoBayes2.MC.warmup=10, OncoBayes2.MC.iter=20, OncoBayes2.MC.chains=1,
                            OncoBayes2.MC.save_warmup=FALSE)

# fit an example model. See documentation for "combo3" example
example_model("combo3")

# print a summary of the prior
prior_summary(blrmfit, digits = 3)

# print a summary of the posterior (model parameters)
print(blrmfit)

# summary of posterior for DLT rate by dose for observed covariate levels
summ <- summary(blrmfit, interval_prob = c(0, 0.16, 0.33, 1))
print(cbind(hist_combo3, summ))

# summary of posterior for DLT rate by dose for new set of covariate levels
newdata <- expand.grid(
  stratum_id = "BID", group_id = "Combo",
  drug_A = 400, drug_B = 800, drug_C = c(320, 400, 600, 800),
  stringsAsFactors = FALSE
)
summ_pred <- summary(blrmfit, newdata = newdata, interval_prob = c(0, 0.16, 0.33, 1))
print(cbind(newdata, summ_pred))

# update the model after observing additional data
newdata$num_patients <- rep(3, nrow(newdata))
newdata$num_toxicities <- c(0, 1, 2, 2)
library(dplyr)
blrmfit_new <- update(blrmfit,
                      data = rbind(hist_combo3, newdata) %>%
                               arrange(stratum_id, group_id))

# updated posterior summary
summ_upd <- summary(blrmfit_new, newdata = newdata, interval_prob = c(0, 0.16, 0.33, 1))
print(cbind(newdata, summ_upd))
## Recover user set sampling defaults
options(.user_mc_options)


[Package OncoBayes2 version 0.8-9 Index]