Sequential design with early stopping (restricted action set) - risk based
Warning: package 'cmdstanr' was built under R version 4.4.3source("./R/util.R")
source("./R/data.R")
log_info("Called simulation-results 7 notebook")
toks <- unlist(tstrsplit(getwd(), "/"))
if(toks[length(toks)] == "roadmap-sim"){
prefix_cfg <- "./etc/sim07/"
prefix_stan <- "./stan"
prefix_fig <- "./fig"
} else {
prefix_cfg <- "../etc/sim07/"
prefix_stan <- "../stan"
prefix_fig <- "../fig"
}# Each input file corresponds to the results from a single simulation
# scenario/configuration.
# Load all the files into a single list.
# files of interest
sim_lab <- "sim07-05"
flist <- list.files(paste0("data/", sim_lab), pattern = "sim07")
toks <- list()
l <- list()
i <- 1
for(i in 1:length(flist)){
l[[i]] <- qs::qread(file.path(paste0("data/", sim_lab), flist[i]))
toks[[i]] <- unlist(tstrsplit(flist[i], "[-.]"))
}Data is generated based on the empirical distributions obtained mostly from the PIANO study.
We simulate silo membership from a multinomial distribution with probabilities 0.3, 0.5 and 0.2 for early, late and chronic.
Conditional on silo membership, we generate domain level allocations and entry indicators that are independent of the design but parameterised to silo specific probabilities. In contrast to the previous simulations, we have simplified the model by ignoring joint level heterogeneity. Adjustment for joint will be reintroduced in the actual analysis model.
For the early silo, we assume:
Similalry, for the late acute silo, we assume:
For the chronic silo, we assume:
After the design independent allocations are made, we finalise the participant level domain entry and allocation based on the known design dependencies.
As the trial progresses, decisions may be made which would lead to some allocations being shut off. For example, if revision is deemed superiori to DAIR, then subsequent allocations within the late acute silo would direct all participants to revision with a split between one-stage and two-stage continuing to be based on the clinical preferences.
The true log-odds of response by subject is calculated as the sum of their parameters indexed by the levels generated above. Treatment success is simulated as a bernoulli random variable with probability equal to the inverse logit transform of the log-odds from the linear predictor. To speed up the model, we aggregate number of successes and number of trials by covariate group which gives the analogous binomial random variable representation.
We adopt a revised model where we allow for silo-specific effects in the surgical domain.
If, for whatever reason, the early or chronic silo show variation in the surgical domain effects, adjustment for silo is inadequate to account for this, and we end up with a polluted version of the surgical domain parameters. The remainder of the model is analogous to the previous versions. The model is used to compute unit level risk (probability) and assesses decisions based on risk difference for the intervention comparisons by domain.
For the simulations, we have a single, multivariable logistic regression model with a linear predictor that incorporates all domains and is specified as follows:
\[ \begin{aligned} Y &\sim \text{Binomial}(n, p) \\ \text{logit}(p) &= \begin{cases} \mu + \beta_{s} + \beta_{p} + \beta_{d1[k_{d1}, s]} + \beta_{d4[k_{d4}]} & \text{(dair)} \\ \mu + \beta_{s} + \beta_{d1[k_{d1}, s]} + \beta_{d2[k_{d2}]} + \beta_{d4[k_{d4}]} & \text{(one-stage)} \\ \mu + \beta_{s} + \beta_{p} + \beta_{d1[k_{d1}, s]} + \beta_{d3[k_{d3}]} + \beta_{d4[k_{d4}]} & \text{(two-stage)} \end{cases} \end{aligned} \]
for each distinct covariate grouping where:
The way that terms enter the model is somewhat convoluted and understanding the dependency implications and consequently interpretations is fairly challenging. For the actual trial analysis, the model will be revised to an equivalent Bernoulli likelihood on unit level data but extended to account for joint, site, randomisation period (if required) and prognostic covariates. Bar the surgical domain, for which ‘by silo’ deviations are implicit in the existing parameterisation, no further interactions are included. Given that heterogeneity by site of infection (joint) is of interest (primarily in the surgical domain) the analysis plan will include specification to account for this.
At each interim, we assess the posterior and if a decision threshold is met, we act. For example, if a superiority decision is reached in one of the domains for which this decision type is relevant, then we consider that domain dealt with and all subsequent participants are assigned to receive the superior intervention. We can (and presently do) continue to update the posterior inference for the comparison that has stopped in subsequent interim analyses until we get to the point where all questions have been answered in all domains, at which point the trial will stop.
Superiority and non-inferiority are applicable to some domains and not others, however, we define reference and threshold values for all domains, just in case. Decisions are made with respect to the average within any `random-effect’ terms that might exist within the model.
For the superiority decision, a reference value of 0 was used and the probability thresholds were:
There was no particular reason for the choice of these thresholds bar the fact that they lead to preferred operating characteristics and give nominal control of the type-i assertion probabilities.
For the futility decision (in relation to superiority) a reference value of 0.05 was used and the probability thresholds were:
The above means, for example, that if the probability that the risk difference is greater than 0.05 is less than in the surgical domain comparison, then we say the superiority goal is futile.
For the ni decision, a reference value of -0.05 was used and the probability thresholds were:
The above means, for example, that if the probability that the risk difference is greater than -0.05 is greater than 0.945, in the antibiotic duration domain, then we will say the intervention is non-inferior.
The futility decision (in relation to non-inferiority) has a reference value of 0 and the probability thresholds were:
This means, for example, that if the probability that the risk difference is greater than 0 is less than 0.1, in the antibiotic duration domain, then we say the non-inferiority goal is futile.
The priors were as follows:
The prior predictive distribution is consistent with a mean probability across of response over the covariate groupings that covers the entire support on the probability scale, is centred around 0.6 and has first and third quarters at 0.4 to 0.8.
For this set of simulations, the number of simulated trials per scenario was 10000, the simulation label is sim07-05.
# Cumulative probability of decisions:
# Traverse the list of simulation results and for each one summarise the
# cumulative probability of each decision type.
i <- 8
d_tbl_1 <- data.table()
# For each scenario that was simulated
for(i in 1:length(l)){
# extract the decision matrix - sim, analysis, quantity, domain level decision
d_dec_1 <- copy(l[[i]]$d_decision)
# config for scenario
l_cfg <- copy(l[[i]]$cfg)
# number of enrolments at each interim (interim sample size sequence)
d_N <- data.table(ia = seq_along(l_cfg$N_pt), N = cumsum(l_cfg$N_pt))
# long version of decisions
d_dec_2 <- melt(d_dec_1,
id.vars = c("sim", "ia", "quant"),
variable.name = "domain")
# Should be right, but just in case...
if(any(is.na(d_dec_2$value))){
message("Some of the decision values are NA in index ", i, " file ", flist[i])
d_dec_2[is.na(value), value := FALSE]
}
d_dec_2[, domain := as.numeric(gsub("d", "", domain))]
# Domains 1, 3 and 4 will stop for superiority or futility for superiority.
# Domaain 2 will stop for NI or futility for NI.
# No other stopping rules apply and so we only evaluate the operating
# characteristics on these, i.e. we do not care about the results for the
# cumualative probability of ni for domain 1, 3 and 4 because we would never
# stop for this.
d_dec_2 <- rbind(
d_dec_2[domain %in% c(1, 3, 4) & quant %in% c("sup", "fut_sup")],
d_dec_2[domain %in% c(2) & quant %in% c("ni", "fut_ni")]
)
# compute the cumulative instances of a decision being made by sim, each
# decision type and by parameter
d_dec_2[, value := as.logical(cumsum(value)>0), keyby = .(sim, quant, domain)]
d_dec_2 <- merge(d_dec_2, d_N, by = "ia")
# cumulative proportion for which each decision quantity has been met by
# analysis and domain
d_dec_cumprob <- d_dec_2[, .(pr_val = mean(value)), keyby = .(ia, N, quant, domain)]
d_tbl_1 <- rbind(
d_tbl_1,
cbind(scenario = i, desc = l_cfg$desc, d_dec_cumprob)
)
}# Expected sample size
# Similar to above but focus on expected sample size
# Traverse the list of simulation results and for each one summarise the
# sample sizes at which stopping for a domain occurs for any reason.
# All we are trying to get to is the expected sample size by domain and
# are not really interested in what decision was made. The cumulative prob
# of each decision type is computed previously.
i <- 2
d_tbl_2 <- data.table()
for(i in 1:length(l)){
# extract the decision matrix - sim, analysis, quantity, domain level decision
d_dec_1 <- copy(l[[i]]$d_decision)
# config for scenario
l_cfg <- copy(l[[i]]$cfg)
# interim looks
d_N <- data.table(ia = seq_along(l_cfg$N_pt), N = cumsum(l_cfg$N_pt))
# long version of decision
d_dec_2 <- melt(d_dec_1,
id.vars = c("sim", "ia", "quant"),
variable.name = "domain")
# Should be right, but just in case...
if(any(is.na(d_dec_2$value))){
message("Some of the decision values are NA in index ", i, " file ", flist[i])
d_dec_2[is.na(value), value := FALSE]
}
d_dec_2[, domain := as.numeric(gsub("d", "", domain))]
d_dec_2 <- merge(d_dec_2, d_N, by = "ia")
# First instance of any decision rule being hit by sim and domain.
# Domains 1, 3 and 4 will stop for superiority or futility for superiority.
# Domaain 2 will stop for NI or futility for NI.
# Sometimes a decision rule will not be hit for a domain and it will continue
# to the max sample size. We will deal with that in a minute.
d_dec_stop <- rbind(
d_dec_2[
domain %in% c(1, 3, 4) & value == T & quant %in% c("sup", "fut_sup"),
.SD[1], keyby = .(sim, domain)],
d_dec_2[
domain %in% c(2) & value == T & quant %in% c("ni", "fut_ni"),
.SD[1], keyby = .(sim, domain)]
)
setnames(d_dec_stop, "N", "N_stopped")
setnames(d_dec_stop, "value", "stopped_early")
# Add in any rows for which no early stopping happened
d_dec_stop <- merge(d_dec_stop[, .SD, .SDcols = !c("ia")],
# all combinations of sim and domain
# which with leave non-stoppers with NA
unique(d_dec_2[, .(sim, domain)]),
by = c("sim", "domain"), all.y = T)
# If domain or trial not stopped then record as having run to the
# maximum sample size with no decision made.
d_dec_stop[is.na(N_stopped), N_stopped := max(d_N$N)]
d_dec_stop[is.na(stopped_early), stopped_early := F]
# So now we know where every domain was stopped and the reason for stopping
# within each sim. Great.
d_dec_stop[is.na(quant), quant := "na"]
d_tbl_2 <- rbind(
d_tbl_2,
cbind(
scenario = i,
desc = l_cfg$desc,
d_dec_stop
)
)
}# Distribution of posterior means for parameters of interest.
# Some simulated trials will have stopped prior to the maximum sample size and
# these will have NA for their posterior means. If you were to summarise these
# posterior means, they would be conditional on the trial having 'survived'
# until the relevant interim. This means that you have missing data at later
# interims, which creates a selection bias in that your selection of sims at any
# given interim are not a random sample, but rather a sample conditioned on the
# stopping rules.
# If you do not account for this in some way then a summary can be either
# optimistic or pessimistic depending on how the stopping rules interact
# with the data. Here we try to account for this missingness by imputing the
# missing posterior means with locf within each simulation.
# Note that this is really only a partial 'fix' to get a sense of whether
# our estimates is representative of the parameter values we used to simulate
# the data.
i <- 1
d_tbl_3 <- data.table()
for(i in 1:length(l)){
# config for scenario
l_cfg <- copy(l[[i]]$cfg)
# params
d_pars <- copy(l[[i]]$d_post_smry_1)
d_pars <- d_pars[par %in% c("lor", "rd")]
# interim looks
d_N <- data.table(ia = seq_along(l_cfg$N_pt), N = cumsum(l_cfg$N_pt))
d_pars <- dcast(d_pars, sim + ia + domain ~ par, value.var = c("mu", "se"))
# locf
d_pars[, `:=`(mu_lor = nafill(mu_lor, type = "locf"),
mu_rd = nafill(mu_rd, type = "locf"),
se_lor = nafill(se_lor, type = "locf"),
se_rd = nafill(se_rd, type = "locf")
),
keyby = .(sim, domain)]
#
d_pars <- merge(d_pars, d_N, by = "ia")
d_tbl_3 <- rbind(
d_tbl_3,
cbind(
scenario = i, desc = l_cfg$desc,
d_pars[, .(ia, sim, domain, N, mu_lor, mu_rd, se_lor, se_rd)]
)
)
}#
i <- 1
d_N_by_arm <- data.table()
for(i in 1:length(l)){
# extract the decision matrix - sim, analysis, quantity, domain level decision
d_dec_1 <- copy(l[[i]]$d_decision)
# long version of decisions
d_dec_2 <- melt(
d_dec_1, id.vars = c("sim", "ia", "quant"), variable.name = "domain")
# Should be right, but just in case...
if(any(is.na(d_dec_2$value))){
message("Some of the decision values are NA in index ", i, " file ", flist[i])
d_dec_2[is.na(value), value := FALSE]
}
d_dec_2[, domain := as.numeric(gsub("d", "", domain))]
# Domains 1, 3 and 4 will stop for superiority or futility for superiority.
# Domaain 2 will stop for NI or futility for NI.
# No other stopping rules apply and so we only evaluate the operating
# characteristics on these, i.e. we do not care about the results for the
# cumualative probability of ni for domain 1, 3 and 4 because we would never
# stop for this.
d_dec_2 <- rbind(
d_dec_2[domain %in% c(1, 3, 4) & quant %in% c("sup", "fut_sup")],
d_dec_2[domain %in% c(2) & quant %in% c("ni", "fut_ni")]
)
# config for scenario
l_cfg <- copy(l[[i]]$cfg)
# interim looks
d_enrolment <- data.table(ia = seq_along(l_cfg$N_pt), N_enrol = cumsum(l_cfg$N_pt))
# observed trial data sets
d_tru <- copy(l[[i]]$d_all)
d_tru[, `:=`(
d1 = as.integer(d1),
d2 = as.integer(d2),
d3 = as.integer(d3),
d4 = as.integer(d4))]
# late acute, surgical arms
# *** silo isn't included in the keyby because we want to average across
# the entire trial population, not silo specific sample sizes.
d_1_tmp <- d_tru[
s == 2, .(N = sum(N), domain = 1), keyby = .(sim, ia, d1)]
d_1_tmp[, N := cumsum(N), keyby = .(sim, d1)]
setnames(d_1_tmp, old = "d1", "arm")
d_1_tmp <- merge(
d_1_tmp,
CJ(sim = 1:max(d_tru$sim), ia = 1:5, arm = 1:3, domain = 1),
by = c("sim", "ia", "arm", "domain"),
all.y = T)
d_1_tmp[, N := nafill(N, type = "locf"), keyby = .(sim, arm)]
setkey(d_1_tmp, sim, ia, arm)
# one-stage revision, ab duration randomised arms
d_2_tmp <- d_tru[
d1 == 2 & d2 %in% 2:3,
.(N = sum(N), domain = 2), keyby = .(sim, ia, d2)]
d_2_tmp[, N := cumsum(N), keyby = .(sim, d2)]
setnames(d_2_tmp, old = "d2", "arm")
d_2_tmp <- merge(
d_2_tmp,
CJ(sim = 1:max(d_tru$sim), ia = 1:5, arm = 2:3, domain = 2),
by = c("sim", "ia", "arm", "domain"),
all.y = T)
d_2_tmp[ia == 1 & is.na(N), N := 0]
d_2_tmp[, N := nafill(N, type = "locf"), keyby = .(sim, arm)]
setkey(d_2_tmp, sim, ia, arm)
# two-stage revision, extended prophy randomised arms
d_3_tmp <- d_tru[
d1 == 3 & d3 %in% 2:3,
.(N = sum(N), domain = 3), keyby = .(sim, ia, d3)]
d_3_tmp[, N := cumsum(N), keyby = .(sim, d3)]
setnames(d_3_tmp, old = "d3", "arm")
d_3_tmp <- merge(
d_3_tmp,
CJ(sim = 1:max(d_tru$sim), ia = 1:5, arm = 2:3, domain = 3),
by = c("sim", "ia", "arm", "domain"),
all.y = T)
d_3_tmp[ia == 1 & is.na(N), N := 0]
d_3_tmp[, N := nafill(N, type = "locf"), keyby = .(sim, arm)]
setkey(d_3_tmp, sim, ia, arm)
# ab choice randomised arms
d_4_tmp <- d_tru[
d4 %in% 2:3, .(N = sum(N), domain = 4), keyby = .(sim, ia, d4)]
d_4_tmp[, N := cumsum(N), keyby = .(sim, d4)]
setnames(d_4_tmp, old = "d4", "arm")
d_4_tmp <- merge(
d_4_tmp,
CJ(sim = 1:max(d_tru$sim), ia = 1:5, arm = 2:3, domain = 4),
by = c("sim", "ia", "arm", "domain"),
all.y = T)
d_4_tmp[ia == 1 & is.na(N), N := 0]
d_4_tmp[, N := nafill(N, type = "locf"), keyby = .(sim, arm)]
setkey(d_4_tmp, sim, ia, arm)
d_arm_tmp <- rbind(d_1_tmp, d_2_tmp, d_3_tmp, d_4_tmp)
d_arm_tmp[, `:=`(
scenario = i, desc = l_cfg$desc
)]
d_arm_tmp <- merge(
d_arm_tmp,
d_enrolment,
by = "ia")
d_N_by_arm <- rbind(d_N_by_arm, d_arm_tmp)
}i <- 1
d_tbl_4 <- data.table()
for(i in 1:length(l)){
# extract the decision matrix - sim, analysis, quantity, domain level decision
d_dec_1 <- copy(l[[i]]$d_decision)
# config for scenario
l_cfg <- copy(l[[i]]$cfg)
# interim looks
d_enrolment <- data.table(ia = seq_along(l_cfg$N_pt), N_enrol = cumsum(l_cfg$N_pt))
# observed data
d_all <- copy(l[[i]]$d_all)
d_all <- merge(d_all, d_enrolment , by = "ia")
# long version of decision
d_dec_2 <- melt(
d_dec_1, id.vars = c("sim", "ia", "quant"), variable.name = "domain")
# Should be right, but just in case...
if(any(is.na(d_dec_2$value))){
message("Some of the decision values are NA in index ", i, " file ", flist[i])
d_dec_2[is.na(value), value := FALSE]
}
d_dec_2[, domain := as.numeric(gsub("d", "", domain))]
d_dec_2 <- merge(d_dec_2, d_N, by = "ia")
# First instance of any decision rule being hit by sim and domain.
# Domains 1, 3 and 4 will stop for superiority or futility for superiority.
# Domaain 2 will stop for NI or futility for NI.
# Sometimes a decision rule will not be hit for a domain and it will continue
# to the max sample size. We will deal with that in a minute.
d_dec_stop <- rbind(
d_dec_2[
domain %in% c(1, 3, 4) & value == T & quant %in% c("sup", "fut_sup"),
.SD[1], keyby = .(sim, domain)],
d_dec_2[
domain %in% c(2) & value == T & quant %in% c("ni", "fut_ni"),
.SD[1], keyby = .(sim, domain)]
)
setnames(d_dec_stop, "N", "N_stopped")
setnames(d_dec_stop, "value", "stopped_early")
# Add in any rows for which no early stopping happened
d_dec_stop <- merge(
d_dec_stop[, .SD, .SDcols = !c("ia")],
# all combinations of sim and domain
unique(d_dec_2[, .(sim, domain)]),
by = c("sim", "domain"), all.y = T)
# If domain or trial not stopped then record as having run to the
# maximum sample size with no decision made.
d_dec_stop[is.na(N_stopped), N_stopped := max(d_N$N)]
d_dec_stop[is.na(stopped_early), stopped_early := F]
# So now we know where every domain was stopped and the reason for stopping
# within each sim. Great.
d_dec_stop[is.na(quant), quant := "null"]
d_domain_1 <- d_all[
s == 2 & d1 %in% 1:3, .(sim, ia, s, pref, arm = d1, y, N, N_enrol)]
d_domain_1 <- merge(d_domain_1, d_dec_stop[domain == 1], by = c("sim"))
d_domain_1 <- d_domain_1[N_enrol <= N_stopped, ]
d_domain_1 <- d_domain_1[, .(domain = 1, y = sum(y), N = sum(N)), keyby = .(arm)]
d_domain_1[, p_hat := y / N]
d_domain_2 <- d_all[
d1 == 2 & d2 %in% 2:3, .(sim, ia, s, pref, arm = d2, y, N, N_enrol)]
d_domain_2 <- merge(d_domain_2, d_dec_stop[domain == 2], by = c("sim"))
d_domain_2 <- d_domain_2[N_enrol <= N_stopped, ]
d_domain_2 <- d_domain_2[, .(domain = 2, y = sum(y), N = sum(N)), keyby = .(arm)]
d_domain_2[, p_hat := y / N]
d_domain_3 <- d_all[
d1 == 3 & d3 %in% 2:3, .(sim, ia, s, pref, arm = d3, y, N, N_enrol)]
d_domain_3 <- merge(d_domain_3, d_dec_stop[domain == 3], by = c("sim"))
d_domain_3 <- d_domain_3[N_enrol <= N_stopped, ]
d_domain_3 <- d_domain_3[, .(domain = 3, y = sum(y), N = sum(N)), keyby = .(arm)]
d_domain_3[, p_hat := y / N]
d_domain_4 <- d_all[
d4 %in% 2:3, .(sim, ia, s, pref, arm = d4, y, N, N_enrol)]
d_domain_4 <- merge(d_domain_4, d_dec_stop[domain == 4], by = c("sim"))
d_domain_4 <- d_domain_4[N_enrol <= N_stopped, ]
d_domain_4 <- d_domain_4[, .(domain = 4, y = sum(y), N = sum(N)), keyby = .(arm)]
d_domain_4[, p_hat := y / N]
d_domain <- rbind(
d_domain_1, d_domain_2, d_domain_3, d_domain_4
)
d_tbl_4 <- rbind(
d_tbl_4,
cbind(scenario = i, desc = l_cfg$desc, d_domain)
)
}Table 1 shows the cumulative probability of a superiority decision across each of the scenarios simulated (the same information is shown in Figure 1). Operating characteristics are shown only for the relevant domains and the futility of a superiority decision is included in parentheses.
d_tbl_1_cur <- d_tbl_1[quant %in% c("sup", "fut_sup") & domain %in% c(1, 3, 4), .SD]
setorderv(d_tbl_1_cur, cols = c("scenario", "domain", "ia", "quant"),
order = c(1, 1, 1, -1))
d_tbl_1_cur <- dcast(d_tbl_1_cur, scenario + desc + domain ~ quant + N, value.var = "pr_val")
d_tbl_1_cur <- d_tbl_1_cur[, .SD, .SDcols = !c("scenario")]
d_tbl_1_cur[, domain := factor(domain,
levels = c(1, 3, 4),
labels = c("Surgical", "AB Ext-proph", "AB Choice"))]
g_tbl <- d_tbl_1_cur |>
gt(groupname_col = "desc") |>
gt::text_transform(
locations = cells_row_groups(),
fn = function(x) {
lapply(x, function(x) {
gt::md(paste0("*", x, "*"))
})
}
) |>
cols_align(
columns = 1,
align = "left"
) |>
cols_align(
columns = 2:ncol(d_tbl_1_cur),
align = "center"
) |>
cols_merge(
columns = c("sup_500", "fut_sup_500"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("sup_1000", "fut_sup_1000"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("sup_1500", "fut_sup_1500"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("sup_2000", "fut_sup_2000"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("sup_2500", "fut_sup_2500"
),
pattern = "<<{1}>><< ({2})>>"
) |>
tab_spanner(
label = md("Cumulative probability of **superiority** (futility) decision"),
columns = 3:ncol(d_tbl_1_cur)
) |>
cols_label(
domain = "Domain",
sup_500 = html("500"),
sup_1000 = html("1000"),
sup_1500 = html("1500"),
sup_2000 = html("2000"),
sup_2500 = html("2500")
) |>
tab_options(
table.font.size = "70%"
) |>
fmt_number(decimals = 3, drop_trailing_zeros = TRUE)
g_tbl| Domain |
Cumulative probability of superiority (futility) decision
|
||||
|---|---|---|---|---|---|
| 500 | 1000 | 1500 | 2000 | 2500 | |
| RD = 0 in all domains +silo specific d1 | |||||
| Surgical | 0.021 (0.538) | 0.037 (0.709) | 0.051 (0.79) | 0.06 (0.839) | 0.066 (0.869) |
| AB Ext-proph | 0.04 (0.543) | 0.062 (0.685) | 0.076 (0.755) | 0.085 (0.797) | 0.092 (0.824) |
| AB Choice | 0.01 (0.653) | 0.016 (0.83) | 0.021 (0.905) | 0.024 (0.94) | 0.027 (0.96) |
| RD = 0.12: surgical revision (one and two-stage) +silo specific d1 | |||||
| Surgical | 0.318 (0.068) | 0.557 (0.091) | 0.703 (0.099) | 0.795 (0.104) | 0.842 (0.105) |
| AB Ext-proph | 0.042 (0.545) | 0.067 (0.727) | 0.085 (0.818) | 0.095 (0.858) | 0.1 (0.881) |
| AB Choice | 0.009 (0.661) | 0.014 (0.837) | 0.02 (0.911) | 0.024 (0.946) | 0.025 (0.963) |
| RD = 0.12: surgical revision (one-stage only) +silo specific d1 | |||||
| Surgical | 0.043 (0.393) | 0.085 (0.539) | 0.117 (0.618) | 0.144 (0.671) | 0.167 (0.707) |
| AB Ext-proph | 0.037 (0.542) | 0.059 (0.69) | 0.074 (0.764) | 0.086 (0.809) | 0.092 (0.837) |
| AB Choice | 0.011 (0.647) | 0.019 (0.825) | 0.024 (0.901) | 0.028 (0.94) | 0.03 (0.957) |
| RD = 0.12: surgical revision (two-stage only) +silo specific d1 | |||||
| Surgical | 0.201 (0.132) | 0.379 (0.177) | 0.509 (0.202) | 0.602 (0.215) | 0.664 (0.223) |
| AB Ext-proph | 0.038 (0.545) | 0.063 (0.722) | 0.079 (0.808) | 0.087 (0.854) | 0.093 (0.88) |
| AB Choice | 0.01 (0.654) | 0.018 (0.832) | 0.022 (0.904) | 0.026 (0.941) | 0.028 (0.959) |
| RD = 0.12: abx duration 6wk effect +silo specific d1 | |||||
| Surgical | 0.021 (0.518) | 0.036 (0.691) | 0.049 (0.777) | 0.06 (0.827) | 0.066 (0.861) |
| AB Ext-proph | 0.038 (0.55) | 0.061 (0.692) | 0.073 (0.764) | 0.082 (0.807) | 0.087 (0.835) |
| AB Choice | 0.008 (0.654) | 0.016 (0.829) | 0.021 (0.906) | 0.025 (0.941) | 0.027 (0.958) |
| RD = 0.12: ext-proph 12wk effect +silo specific d1 | |||||
| Surgical | 0.028 (0.486) | 0.047 (0.66) | 0.061 (0.748) | 0.073 (0.801) | 0.08 (0.837) |
| AB Ext-proph | 0.385 (0.088) | 0.607 (0.106) | 0.72 (0.113) | 0.785 (0.118) | 0.822 (0.119) |
| AB Choice | 0.009 (0.664) | 0.017 (0.831) | 0.022 (0.905) | 0.026 (0.942) | 0.029 (0.958) |
| RD = 0.12: abx choice rif effect +silo specific d1 | |||||
| Surgical | 0.019 (0.538) | 0.035 (0.71) | 0.045 (0.795) | 0.055 (0.842) | 0.062 (0.871) |
| AB Ext-proph | 0.042 (0.545) | 0.062 (0.688) | 0.075 (0.758) | 0.083 (0.805) | 0.09 (0.831) |
| AB Choice | 0.458 (0.037) | 0.802 (0.04) | 0.925 (0.041) | 0.954 (0.042) | 0.96 (0.042) |
| RD = 0.12: all domains +silo specific d1 | |||||
| Surgical | 0.38 (0.046) | 0.618 (0.063) | 0.753 (0.07) | 0.828 (0.074) | 0.872 (0.077) |
| AB Ext-proph | 0.434 (0.074) | 0.74 (0.089) | 0.861 (0.093) | 0.896 (0.094) | 0.907 (0.094) |
| AB Choice | 0.473 (0.035) | 0.815 (0.04) | 0.924 (0.041) | 0.953 (0.041) | 0.96 (0.041) |
| RD = -0.05: abx duration 6wk effect +silo specific d1 | |||||
| Surgical | 0.019 (0.533) | 0.036 (0.695) | 0.048 (0.778) | 0.059 (0.826) | 0.066 (0.859) |
| AB Ext-proph | 0.035 (0.547) | 0.057 (0.684) | 0.069 (0.755) | 0.078 (0.801) | 0.085 (0.832) |
| AB Choice | 0.011 (0.651) | 0.016 (0.835) | 0.019 (0.91) | 0.023 (0.944) | 0.025 (0.961) |
| RD = 0.12 surgical, RD = 0.08 abx choice +silo specific d1 | |||||
| Surgical | 0.322 (0.07) | 0.566 (0.091) | 0.709 (0.101) | 0.792 (0.107) | 0.842 (0.11) |
| AB Ext-proph | 0.038 (0.56) | 0.062 (0.743) | 0.077 (0.832) | 0.086 (0.873) | 0.091 (0.895) |
| AB Choice | 0.482 (0.032) | 0.828 (0.035) | 0.936 (0.035) | 0.963 (0.036) | 0.968 (0.036) |
d_fig <- d_tbl_1[quant %in% c("sup", "fut_sup") & domain %in% c(1, 3, 4), .SD]
setorderv(d_fig, cols = c("scenario", "domain", "ia", "quant"),
order = c(1, 1, 1, -1))
d_fig[, domain := factor(domain,
levels = c(1, 3, 4),
labels = c("Surgical", "AB Ext-proph", "AB Choice"))]
d_fig[, quant := factor(quant,
levels = c("sup", "fut_sup"),
labels = c("Superiority", "Futility"))]
d_fig[, desc := factor(desc,
levels = unique(d_fig$desc),
labels = unique(d_fig$desc))]
ggplot(d_fig, aes(x = N, y = pr_val, group = quant, col = quant)) +
geom_line(lwd = 0.25) +
scale_y_continuous("") +
scale_color_discrete("") +
# facet_grid(desc ~ domain) +
# facet_grid2(desc ~ domain, render_empty = FALSE)
ggh4x::facet_grid2(desc ~ domain ,
labeller = labeller(desc = label_wrap_gen(35)),
scales = "free",
axes = "y",
independent = "y") +
# theme_minimal() +
theme(
legend.position = "bottom",
strip.text.y.right = element_text(angle = 0,
hjust = 0,
vjust = 0.2,
size = 4),
axis.ticks = element_blank(),
strip.text.x.top = element_text(size = 4),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "grey",
linewidth = 0.1,
linetype = 1),
axis.text.y = element_text(size = 4),
axis.text.x = element_text(size = 4),
axis.title.x = element_text(size = 5),
legend.text = element_text(size = 4)
)
Table 2 shows the cumulative probability of a non-inferiority decision with futility shown in parentheses (the same information is shown in Figure 2). The results are only shown for the domains for which non-inferiority is evaluated.
The “Null effect in all domains” is actually a bit of a misnomer as the true null with regards to the NI decision would be at the NI margin, whereas the scenario refers to the setting where all effects are set to zero. Thus an inflation over the usual type-i assertion probability is to be expected.
d_tbl_1_cur <- d_tbl_1[quant %in% c("ni", "fut_ni") & domain %in% c(2), .SD]
setorderv(d_tbl_1_cur, cols = c("scenario", "domain", "ia", "quant"),
order = c(1, 1, 1, -1))
d_tbl_1_cur <- dcast(d_tbl_1_cur, scenario + desc + domain ~ quant + N, value.var = "pr_val")
d_tbl_1_cur <- d_tbl_1_cur[, .SD, .SDcols = !c("scenario")]
d_tbl_1_cur[, domain := factor(domain,
levels = c(2),
labels = c("AB Duration"))]
g_tbl <- d_tbl_1_cur |>
gt(groupname_col = "desc") |>
gt::text_transform(
locations = cells_row_groups(),
fn = function(x) {
lapply(x, function(x) {
gt::md(paste0("*", x, "*"))
})
}
) |>
cols_align(
columns = 1,
align = "left"
) |>
cols_align(
columns = 2:ncol(d_tbl_1_cur),
align = "center"
) |>
cols_merge(
columns = c("ni_500", "fut_ni_500"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("ni_1000", "fut_ni_1000"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("ni_1500", "fut_ni_1500"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("ni_2000", "fut_ni_2000"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("ni_2500", "fut_ni_2500"
),
pattern = "<<{1}>><< ({2})>>"
) |>
tab_spanner(
label = md("Cumulative probability of **NI** (futility) decision"),
columns = 3:ncol(d_tbl_1_cur)
) |>
cols_label(
domain = "Domain",
ni_500 = html("500"),
ni_1000 = html("1000"),
ni_1500 = html("1500"),
ni_2000 = html("2000"),
ni_2500 = html("2500")
) |>
tab_options(
table.font.size = "70%"
) |>
fmt_number(decimals = 3, drop_trailing_zeros = TRUE)
g_tbl| Domain |
Cumulative probability of NI (futility) decision
|
||||
|---|---|---|---|---|---|
| 500 | 1000 | 1500 | 2000 | 2500 | |
| RD = 0 in all domains +silo specific d1 | |||||
| AB Duration | 0.112 (0.088) | 0.181 (0.131) | 0.231 (0.164) | 0.275 (0.185) | 0.309 (0.203) |
| RD = 0.12: surgical revision (one and two-stage) +silo specific d1 | |||||
| AB Duration | 0.115 (0.088) | 0.203 (0.145) | 0.274 (0.183) | 0.336 (0.212) | 0.387 (0.231) |
| RD = 0.12: surgical revision (one-stage only) +silo specific d1 | |||||
| AB Duration | 0.11 (0.092) | 0.186 (0.14) | 0.243 (0.173) | 0.284 (0.194) | 0.32 (0.212) |
| RD = 0.12: surgical revision (two-stage only) +silo specific d1 | |||||
| AB Duration | 0.104 (0.093) | 0.185 (0.141) | 0.247 (0.177) | 0.305 (0.203) | 0.353 (0.222) |
| RD = 0.12: abx duration 6wk effect +silo specific d1 | |||||
| AB Duration | 0.38 (0.012) | 0.595 (0.016) | 0.717 (0.017) | 0.796 (0.018) | 0.853 (0.018) |
| RD = 0.12: ext-proph 12wk effect +silo specific d1 | |||||
| AB Duration | 0.113 (0.087) | 0.193 (0.132) | 0.241 (0.162) | 0.282 (0.184) | 0.315 (0.202) |
| RD = 0.12: abx choice rif effect +silo specific d1 | |||||
| AB Duration | 0.107 (0.09) | 0.176 (0.138) | 0.227 (0.17) | 0.268 (0.193) | 0.305 (0.21) |
| RD = 0.12: all domains +silo specific d1 | |||||
| AB Duration | 0.417 (0.007) | 0.728 (0.009) | 0.878 (0.009) | 0.946 (0.01) | 0.97 (0.01) |
| RD = -0.05: abx duration 6wk effect +silo specific d1 | |||||
| AB Duration | 0.058 (0.166) | 0.09 (0.26) | 0.109 (0.324) | 0.124 (0.374) | 0.135 (0.414) |
| RD = 0.12 surgical, RD = 0.08 abx choice +silo specific d1 | |||||
| AB Duration | 0.114 (0.091) | 0.206 (0.142) | 0.281 (0.179) | 0.344 (0.204) | 0.4 (0.223) |
d_fig <- d_tbl_1[quant %in% c("ni", "fut_ni") & domain %in% c(2), .SD]
setorderv(d_fig, cols = c("scenario", "domain", "ia", "quant"),
order = c(1, 1, 1, -1))
d_fig[, domain := factor(domain,
levels = c(2),
labels = c("AB Duration"))]
d_fig[, quant := factor(quant,
levels = c("ni", "fut_ni"),
labels = c("NI", "Futility for NI"))]
d_fig[, desc := factor(desc,
levels = unique(d_fig$desc),
labels = unique(d_fig$desc))]
ggplot(d_fig, aes(x = N, y = pr_val, group = quant, col = quant)) +
geom_line(lwd = 0.25) +
scale_y_continuous("", breaks = seq(0, 1, by = 0.2)) +
scale_color_discrete("") +
# facet_grid(desc ~ domain, space = "free") +
ggh4x::facet_grid2(desc ~ domain ,
labeller = labeller(desc = label_wrap_gen(35)),
scales = "free",
axes = "y",
independent = "y") +
ggh4x::force_panelsizes(cols = unit(c(4), "cm"), TRUE) +
# facet_manual(
# . ~ desc, design = matrix(1:17, ncol = 1),
# widths = unit(3, "cm")
# ) +
# theme_minimal() +
theme(
legend.position = "bottom",
strip.text.y.right = element_text(angle = 0,
hjust = 0,
vjust = 0.2,
size = 4),
axis.ticks = element_blank(),
strip.text.x.top = element_text(size = 4),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "grey",
linewidth = 0.1,
linetype = 1),
axis.text.y = element_text(size = 4),
axis.text.x = element_text(size = 4),
axis.title.x = element_text(size = 5),
legend.text = element_text(size = 4)
)
Table 3 shows the mean number of enrolments into the platform until a stopping rule is hit for each of the domains by scenario and stopping rule. In parentheses are the percentage of sims for each decision type and those that were not stopped early (last row) run to the maximum sample size (2500).
Figure 3 shows the number of participants entering into each of the randomised comparisons by domain and scenario. The expected cumulative sample sizes (Figure 3) are calculated by extracting the number of participants entering into each randomised comparison (i.e. restricted to the relevant strata) and taking the cumulative sum of these over time. For example, the domain 1 (surgical) expected values are based on the participants in the late acute silo that receive randomised surgical intervention. Similarly, the domain 2 (antibiotic duration) expected values are based on the participants across all silos that received one-stage revision and then receive one of the randomised treatment allocations (i.e. are not in the non-randomised treatment group for any reason).
If a decision was made in a domain, then subsequent enrolments would be assigned to the remaining arms and so the allocation would be expected to deviate away from 1:1. For example, if a superiority decision was made for domain 1 (and the trial was still ongoing) then all subsequent participants are assigned to the superior intervention. Finally, if a decision was made for all research questions, the trial stops early. To avoid any bias in the calculation, we use LOCF to propagate the expectations forward through to the maximum number of enrolments.
d_tbl_2_cur <- d_tbl_2[
, .(N_mu = mean(N_stopped),
pct_sims = sprintf("%.0f%%", 100*.N/l[[1]]$cfg$nsim)),
keyby = .(scenario, desc, domain, quant)]
d_tbl_2_cur <- dcast(d_tbl_2_cur, scenario + desc + quant ~ domain, value.var = list("N_mu", "pct_sims"))
d_tbl_2_cur[, quant := factor(
quant,
levels = c("sup", "fut_sup", "ni", "fut_ni", "na"),
labels = c("superiority", "futility (sup)", "ni", "futility (ni)", "-"))]
d_tbl_2_cur <- d_tbl_2_cur[order(scenario, quant)]
d_tbl_2_cur <- d_tbl_2_cur[, .SD, .SDcols = !c("scenario")]
setcolorder(d_tbl_2_cur, c("desc", "quant",
"N_mu_1", "pct_sims_1",
"N_mu_2", "pct_sims_2",
"N_mu_3", "pct_sims_3",
"N_mu_4", "pct_sims_4"
))
g_tbl <- d_tbl_2_cur |>
gt(groupname_col = "desc") |>
gt::text_transform(
locations = cells_row_groups(),
fn = function(x) {
lapply(x, function(x) {
gt::md(paste0("*", x, "*"))
})
}
) |>
cols_align(
columns = 1:2,
align = "left"
) |>
cols_align(
columns = 3:ncol(d_tbl_2_cur),
align = "center"
) |>
cols_merge(
columns = c("N_mu_1", "pct_sims_1"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("N_mu_2", "pct_sims_2"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("N_mu_3", "pct_sims_3"
),
pattern = "<<{1}>><< ({2})>>"
) |>
cols_merge(
columns = c("N_mu_4", "pct_sims_4"
),
pattern = "<<{1}>><< ({2})>>"
) |>
tab_spanner(
label = html("Expected number of total enrolments to hit stopping rule by domain"),
columns = 3:ncol(d_tbl_2_cur)
) |>
cols_label(
N_mu_1 = "Surgical",
N_mu_2 = "AB Duration",
N_mu_3 = "AB Ext-proph",
N_mu_4 = "AB choice"
) |>
tab_options(
table.font.size = "70%"
) |>
fmt_number(decimals = 0, drop_trailing_zeros = TRUE)
g_tbl| quant |
Expected number of total enrolments to hit stopping rule by domain
|
|||
|---|---|---|---|---|
| Surgical | AB Duration | AB Ext-proph | AB choice | |
| RD = 0 in all domains +silo specific d1 | ||||
| superiority | 1,220 (7%) | 1,061 (9%) | 1,137 (3%) | |
| futility (sup) | 839 (86%) | 811 (82%) | 766 (96%) | |
| ni | 1,208 (31%) | |||
| futility (ni) | 1,101 (20%) | |||
| - | 2,500 (7%) | 2,500 (49%) | 2,500 (9%) | 2,500 (1%) |
| RD = 0.12: surgical revision (one and two-stage) +silo specific d1 | ||||
| superiority | 1,085 (84%) | 1,035 (10%) | 1,127 (2%) | |
| futility (sup) | 775 (10%) | 823 (88%) | 757 (96%) | |
| ni | 1,300 (39%) | |||
| futility (ni) | 1,135 (23%) | |||
| - | 2,500 (6%) | 2,500 (38%) | 2,500 (3%) | 2,500 (1%) |
| RD = 0.12: surgical revision (one-stage only) +silo specific d1 | ||||
| superiority | 1,315 (16%) | 1,092 (9%) | 1,141 (3%) | |
| futility (sup) | 924 (70%) | 821 (83%) | 767 (95%) | |
| ni | 1,212 (32%) | |||
| futility (ni) | 1,090 (21%) | |||
| - | 2,500 (13%) | 2,500 (47%) | 2,500 (8%) | 2,500 (2%) |
| RD = 0.12: surgical revision (two-stage only) +silo specific d1 | ||||
| superiority | 1,221 (66%) | 1,055 (9%) | 1,100 (3%) | |
| futility (sup) | 866 (22%) | 832 (88%) | 761 (96%) | |
| ni | 1,308 (35%) | |||
| futility (ni) | 1,117 (22%) | |||
| - | 2,500 (12%) | 2,500 (42%) | 2,500 (3%) | 2,500 (2%) |
| RD = 0.12: abx duration 6wk effect +silo specific d1 | ||||
| superiority | 1,245 (7%) | 1,047 (9%) | 1,203 (3%) | |
| futility (sup) | 865 (86%) | 813 (83%) | 761 (96%) | |
| ni | 1,041 (85%) | |||
| futility (ni) | 771 (2%) | |||
| - | 2,500 (7%) | 2,500 (13%) | 2,500 (8%) | 2,500 (2%) |
| RD = 0.12: ext-proph 12wk effect +silo specific d1 | ||||
| superiority | 1,177 (8%) | 976 (82%) | 1,200 (3%) | |
| futility (sup) | 884 (83%) | 718 (12%) | 754 (96%) | |
| ni | 1,186 (32%) | |||
| futility (ni) | 1,099 (20%) | |||
| - | 2,500 (9%) | 2,500 (48%) | 2,500 (6%) | 2,500 (1%) |
| RD = 0.12: abx choice rif effect +silo specific d1 | ||||
| superiority | 1,244 (6%) | 1,032 (9%) | 862 (96%) | |
| futility (sup) | 840 (87%) | 817 (83%) | 571 (4%) | |
| ni | 1,225 (31%) | |||
| futility (ni) | 1,094 (21%) | |||
| - | 2,500 (7%) | 2,500 (48%) | 2,500 (8%) | 2,500 (0%) |
| RD = 0.12: all domains +silo specific d1 | ||||
| superiority | 1,019 (87%) | 881 (90%) | 849 (96%) | |
| futility (sup) | 850 (8%) | 640 (9%) | 603 (4%) | |
| ni | 970 (97%) | |||
| futility (ni) | 670 (1%) | |||
| - | 2,500 (5%) | 2,500 (2%) | 2,500 (0%) | 2,500 (0%) |
| RD = -0.05: abx duration 6wk effect +silo specific d1 | ||||
| superiority | 1,254 (6%) | 1,091 (8%) | 1,136 (2%) | |
| futility (sup) | 848 (86%) | 824 (83%) | 762 (96%) | |
| ni | 1,094 (14%) | |||
| futility (ni) | 1,143 (41%) | |||
| - | 2,500 (8%) | 2,500 (45%) | 2,500 (8%) | 2,500 (1%) |
| RD = 0.12 surgical, RD = 0.08 abx choice +silo specific d1 | ||||
| superiority | 1,074 (83%) | 1,032 (9%) | 838 (96%) | |
| futility (sup) | 815 (11%) | 812 (89%) | 571 (4%) | |
| ni | 1,318 (40%) | |||
| futility (ni) | 1,120 (22%) | |||
| - | 2,500 (6%) | 2,500 (38%) | 2,500 (2%) | 2,500 (0%) |
d_fig_1 <- d_N_by_arm[, .(mu_N = mean(N)), keyby = .(desc, domain, arm, N_enrol)]
d_fig_1[, desc := factor(
desc,
levels = unique(d_N_by_arm$desc))]
d_fig_1[
, domain := factor(
domain,
levels = 1:4,
labels = c("Surgical", "AB Duration", "AB Ext-proph", "AB Choice"))]
# additional arm to show the total N for revision(1) and revision(2) in the
# surgical domain
d_fig_2 <- d_N_by_arm[
domain == 1 & arm %in% 2:3,
.(N = sum(N), arm = 4),
keyby = .(ia, sim, domain, scenario, desc, N_enrol)]
d_fig_2 <- d_fig_2[, .(mu_N = mean(N)), keyby = .(desc, domain, arm, N_enrol)]
d_fig_2[, desc := factor(
desc,
levels = unique(d_N_by_arm$desc))]
d_fig_2[
, domain := factor(
domain,
levels = 1:4,
labels = c("Surgical", "AB Duration", "AB Ext-proph", "AB Choice"))]
d_fig_1[, arm := factor(arm)]
ggplot(d_fig_1, aes(x = N_enrol, y = mu_N, col = arm)) +
geom_line(lwd = 0.2) +
geom_point(size = 0.4) +
geom_line(
data = d_fig_2,
aes(x = N_enrol, y = mu_N), lwd = 0.2, lty = 2, inherit.aes = F) +
scale_x_continuous("") +
scale_y_continuous("Expected number of participants") +
scale_color_discrete("Treatment arm: ") +
ggh4x::facet_grid2(
desc ~ domain,
labeller = labeller(desc = label_wrap_gen(35)),
scales = "free_y", independent = "y") +
# theme_minimal() +
theme(
legend.position = "bottom",
strip.text.y.right = element_text(angle = 0,
hjust = 0,
vjust = 0.2,
size = 4),
axis.ticks = element_blank(),
strip.text.x.top = element_text(size = 4),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "grey",
linewidth = 0.1,
linetype = 1),
axis.title.y=element_text(size = 5),
axis.text.y = element_text(size = 4),
axis.text.x = element_text(size = 4),
axis.title.x = element_text(size = 5),
legend.text = element_text(size = 4)
)
Figure 4 shows the median value and 95% quantiles from the posterior mean odds ratios obtained from the simulations by each domain and scenario. The estimates are unconditional, in that they ignore the stopping rule and propagate estimates forward with LOCF so that we are working from a random sample (albeit with static imputation) rather than a dependent sample of simulations.
Figure 5 shows the median value and 95% quantiles for the posterior mean risk differences obtained from the simulations for each domain and scenario. The transformation to a risk difference suggests that the odds ratios amount to effects in the order of 10-20% on the risk scale, dependent on the silo, domain etc.
The AB duration domain has high variance in the distribution of posterior means that we anticipate to observe. That is, when there is no true effect, we might still see posterior means as large as \(\pm 0.2\) on the absolute risk scale.
# ggplot2::theme_update(text = element_text(size = 8))
# ggplot2::theme_update(legend.position = "bottom")
# # ggplot2::theme_update(legend.title = element_blank())
# ggplot2::theme_update(axis.text.x = element_text(size = 8))
# ggplot2::theme_update(axis.text.y = element_text(size = 8))
d_fig <- d_tbl_3[,
.(or = median(exp(mu_lor)),
q_025 = quantile(exp(mu_lor), prob = 0.025),
q_975 = quantile(exp(mu_lor), prob = 0.975)),
keyby = .(scenario, desc, ia, domain, N)]
# setorderv(d_fig, cols = "scenario", order = -1L)
d_fig[, desc := factor(desc, levels = unique(d_fig$desc))]
d_fig[, domain := factor(domain,
levels = 1:4,
labels = c("Surgical", "AB Duration", "AB Ext-proph", "AB Choice"))]
# d_tbl_3[scenario == 8, range(lor)]
ggplot(data = d_fig,
aes(x = N, y = or)) +
geom_line(lwd = 0.25) +
geom_errorbar(aes(ymin = q_025, ymax = q_975), lwd = 0.25, width = 50) +
scale_x_continuous("") +
scale_y_continuous("OR") +
ggh4x::facet_grid2(desc ~ domain ,
labeller = labeller(desc = label_wrap_gen(35)),
scales = "free",
axes = "y",
independent = "y") +
# theme_minimal() +
theme(text = element_text(size = 6),
strip.text.y.right = element_text(angle = 0,
hjust = 0,
vjust = 0.2,
size = 4),
strip.text.x = element_text(angle = 0, size = 5),
axis.ticks = element_blank(),
axis.text.x = element_text(angle = 45, vjust = 1, hjust=1, size = 5),
axis.text.y = element_text(size = 5),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "grey",
linewidth = 0.1,
linetype = 1))
d_fig <- d_tbl_3[,
.(rd = mean(mu_rd),
q_025 = quantile(mu_rd, prob = 0.025),
q_975 = quantile(mu_rd, prob = 0.975)),
keyby = .(scenario, desc, ia, domain, N)]
# setorderv(d_fig, cols = "scenario", order = -1L)
d_fig[, desc := factor(desc, levels = unique(d_fig$desc))]
d_fig[, domain := factor(
domain,
levels = 1:4,
labels = c("Surgical", "AB Duration", "AB Ext-proph", "AB Choice"))]
# d_tbl_3[scenario == 8, range(lor)]
ggplot(data = d_fig,
aes(x = N, y = rd)) +
geom_line(lwd = 0.25) +
geom_errorbar(aes(ymin = q_025, ymax = q_975), lwd = 0.25, width = 50) +
scale_x_continuous("") +
scale_y_continuous("Risk difference") +
ggh4x::facet_grid2(desc ~ domain ,
labeller = labeller(desc = label_wrap_gen(35)),
scales = "free",
axes = "y",
independent = "y") +
# theme_minimal() +
theme(text = element_text(size = 6),
strip.text.y.right = element_text(angle = 0,
hjust = 0,
vjust = 0.2,
size = 4),
strip.text.x = element_text(angle = 0, size = 5),
axis.ticks = element_blank(),
axis.text.x = element_text(angle = 45, vjust = 1, hjust=1, size = 5),
axis.text.y = element_text(size = 5),
panel.grid.minor = element_blank(),
panel.grid.major = element_line(color = "grey",
linewidth = 0.1,
linetype = 1))
d_tbl <- d_tbl_3[,
.(rd = mean(mu_rd),
q_025 = quantile(mu_rd, prob = 0.025),
q_975 = quantile(mu_rd, prob = 0.975)),
keyby = .(scenario, desc, ia, domain, N)]
# setorderv(d_fig, cols = "scenario", order = -1L)
d_tbl[, desc := factor(desc, levels = unique(d_fig$desc))]
d_tbl[, domain := factor(
domain,
levels = 1:4,
labels = c("Surgical", "AB Duration", "AB Ext-proph", "AB Choice"))]
d_tbl <- dcast(d_tbl, desc + domain ~ N, value.var = list("rd", "q_025", "q_975"))
ci_names <- function(x = 500){
paste0(c("rd_","q_025_", "q_975_"), x)
}
setcolorder(
d_tbl,
c("desc", "domain",
ci_names(500), ci_names(1000), ci_names(1500), ci_names(2000), ci_names(2500)))
g_tbl <- d_tbl |>
gt(groupname_col = "desc") |>
gt::text_transform(
locations = cells_row_groups(),
fn = function(x) {
lapply(x, function(x) {
gt::md(paste0("*", x, "*"))
})
}
) |>
cols_align(
columns = 1,
align = "left"
) |>
cols_align(
columns = 2:ncol(d_tbl),
align = "center"
) |>
cols_align(
columns = 3:ncol(d_tbl_2_cur),
align = "right"
) |>
cols_merge(
columns = c("rd_500", "q_025_500", "q_975_500"),
pattern = "<<{1}>><< ({2}, {3})>>"
) |>
cols_merge(
columns = c("rd_1000", "q_025_1000", "q_975_1000"),
pattern = "<<{1}>><< ({2}, {3})>>"
) |>
cols_merge(
columns = c("rd_1500", "q_025_1500", "q_975_1500"),
pattern = "<<{1}>><< ({2}, {3})>>"
) |>
cols_merge(
columns = c("rd_2000", "q_025_2000", "q_975_2000"),
pattern = "<<{1}>><< ({2}, {3})>>"
) |>
cols_merge(
columns = c("rd_2500", "q_025_2500", "q_975_2500"),
pattern = "<<{1}>><< ({2}, {3})>>"
) |>
tab_spanner(
label = html("Risk difference (expectation of posterior means and 95% interval"),
columns = 2:ncol(d_tbl)
) |>
cols_label(
rd_500 = "500",
rd_1000 = "1000",
rd_1500 = "1500",
rd_2000 = "2000",
rd_2500 = "2500"
) |>
tab_options(
table.font.size = "70%"
) |>
fmt_number(decimals = 2, drop_trailing_zeros = TRUE)
g_tbl
Risk difference (expectation of posterior means and 95% interval
|
|||||
|---|---|---|---|---|---|
| domain | 500 | 1000 | 1500 | 2000 | 2500 |
| RD = 0 in all domains +silo specific d1 | |||||
| Surgical | 0 (−0.14, 0.14) | 0 (−0.12, 0.11) | −0.01 (−0.12, 0.1) | −0.01 (−0.11, 0.1) | −0.01 (−0.11, 0.1) |
| AB Duration | 0 (−0.22, 0.22) | 0 (−0.21, 0.2) | 0 (−0.2, 0.2) | 0 (−0.2, 0.2) | 0 (−0.2, 0.2) |
| AB Ext-proph | 0 (−0.16, 0.16) | 0 (−0.13, 0.13) | 0 (−0.13, 0.13) | 0 (−0.13, 0.12) | −0.01 (−0.13, 0.12) |
| AB Choice | 0 (−0.11, 0.11) | 0 (−0.09, 0.08) | 0 (−0.08, 0.07) | 0 (−0.08, 0.06) | 0 (−0.08, 0.06) |
| RD = 0.12: surgical revision (one and two-stage) +silo specific d1 | |||||
| Surgical | 0.11 (−0.03, 0.24) | 0.11 (−0.01, 0.23) | 0.12 (0, 0.22) | 0.12 (0, 0.22) | 0.12 (0, 0.22) |
| AB Duration | 0 (−0.21, 0.21) | 0 (−0.19, 0.19) | 0 (−0.18, 0.19) | 0 (−0.18, 0.18) | 0 (−0.18, 0.18) |
| AB Ext-proph | 0 (−0.15, 0.15) | 0 (−0.13, 0.12) | 0 (−0.12, 0.11) | −0.01 (−0.12, 0.11) | −0.01 (−0.12, 0.11) |
| AB Choice | 0 (−0.1, 0.1) | 0 (−0.08, 0.08) | 0 (−0.08, 0.07) | 0 (−0.08, 0.06) | 0 (−0.08, 0.06) |
| RD = 0.12: surgical revision (one-stage only) +silo specific d1 | |||||
| Surgical | 0.02 (−0.12, 0.16) | 0.02 (−0.1, 0.14) | 0.02 (−0.09, 0.13) | 0.02 (−0.09, 0.13) | 0.02 (−0.09, 0.13) |
| AB Duration | 0 (−0.22, 0.22) | 0 (−0.2, 0.2) | 0 (−0.19, 0.19) | 0 (−0.19, 0.19) | 0 (−0.19, 0.19) |
| AB Ext-proph | 0 (−0.16, 0.15) | 0 (−0.13, 0.13) | 0 (−0.13, 0.12) | −0.01 (−0.13, 0.12) | −0.01 (−0.12, 0.11) |
| AB Choice | 0 (−0.11, 0.11) | 0 (−0.08, 0.08) | 0 (−0.08, 0.07) | 0 (−0.08, 0.06) | 0 (−0.08, 0.06) |
| RD = 0.12: surgical revision (two-stage only) +silo specific d1 | |||||
| Surgical | 0.08 (−0.06, 0.22) | 0.09 (−0.03, 0.2) | 0.09 (−0.03, 0.2) | 0.09 (−0.02, 0.2) | 0.09 (−0.02, 0.2) |
| AB Duration | 0 (−0.22, 0.22) | 0 (−0.2, 0.2) | 0 (−0.19, 0.19) | 0 (−0.19, 0.19) | 0 (−0.19, 0.19) |
| AB Ext-proph | 0 (−0.15, 0.15) | 0 (−0.13, 0.12) | 0 (−0.12, 0.11) | −0.01 (−0.12, 0.1) | −0.01 (−0.11, 0.1) |
| AB Choice | 0 (−0.11, 0.11) | 0 (−0.09, 0.08) | 0 (−0.08, 0.07) | 0 (−0.08, 0.06) | 0 (−0.08, 0.06) |
| RD = 0.12: abx duration 6wk effect +silo specific d1 | |||||
| Surgical | 0 (−0.14, 0.14) | 0 (−0.13, 0.12) | −0.01 (−0.12, 0.11) | −0.01 (−0.12, 0.11) | −0.01 (−0.12, 0.11) |
| AB Duration | 0.1 (−0.12, 0.31) | 0.11 (−0.08, 0.3) | 0.12 (−0.06, 0.3) | 0.12 (−0.04, 0.3) | 0.13 (−0.03, 0.3) |
| AB Ext-proph | 0 (−0.15, 0.16) | 0 (−0.14, 0.13) | 0 (−0.13, 0.13) | −0.01 (−0.13, 0.13) | −0.01 (−0.13, 0.13) |
| AB Choice | 0 (−0.11, 0.11) | 0 (−0.09, 0.08) | 0 (−0.09, 0.07) | 0 (−0.08, 0.06) | 0 (−0.08, 0.06) |
| RD = 0.12: ext-proph 12wk effect +silo specific d1 | |||||
| Surgical | 0.01 (−0.13, 0.15) | 0 (−0.11, 0.12) | 0 (−0.11, 0.11) | 0 (−0.11, 0.11) | 0 (−0.11, 0.1) |
| AB Duration | 0 (−0.22, 0.22) | 0 (−0.2, 0.2) | 0 (−0.19, 0.2) | 0 (−0.19, 0.19) | 0 (−0.19, 0.19) |
| AB Ext-proph | 0.11 (−0.04, 0.26) | 0.12 (−0.02, 0.24) | 0.12 (−0.01, 0.24) | 0.12 (−0.01, 0.24) | 0.12 (−0.01, 0.24) |
| AB Choice | 0 (−0.11, 0.11) | 0 (−0.09, 0.08) | 0 (−0.08, 0.07) | 0 (−0.08, 0.06) | 0 (−0.08, 0.06) |
| RD = 0.12: abx choice rif effect +silo specific d1 | |||||
| Surgical | 0 (−0.14, 0.13) | 0 (−0.12, 0.11) | −0.01 (−0.12, 0.1) | −0.01 (−0.11, 0.09) | −0.01 (−0.11, 0.09) |
| AB Duration | 0 (−0.22, 0.22) | 0 (−0.2, 0.2) | 0 (−0.19, 0.19) | 0 (−0.19, 0.19) | 0 (−0.19, 0.19) |
| AB Ext-proph | 0 (−0.16, 0.15) | 0 (−0.13, 0.13) | 0 (−0.13, 0.12) | 0 (−0.12, 0.12) | −0.01 (−0.12, 0.12) |
| AB Choice | 0.12 (0.01, 0.22) | 0.12 (0.04, 0.2) | 0.12 (0.05, 0.2) | 0.12 (0.05, 0.2) | 0.12 (0.06, 0.2) |
| RD = 0.12: all domains +silo specific d1 | |||||
| Surgical | 0.12 (−0.01, 0.24) | 0.12 (0, 0.23) | 0.12 (0.01, 0.22) | 0.12 (0.01, 0.22) | 0.12 (0.01, 0.22) |
| AB Duration | 0.1 (−0.09, 0.29) | 0.11 (−0.04, 0.28) | 0.12 (−0.01, 0.28) | 0.12 (0, 0.28) | 0.13 (0.01, 0.28) |
| AB Ext-proph | 0.11 (−0.02, 0.25) | 0.12 (0, 0.23) | 0.12 (0.01, 0.23) | 0.12 (0.01, 0.23) | 0.12 (0.01, 0.23) |
| AB Choice | 0.12 (0.02, 0.22) | 0.12 (0.04, 0.2) | 0.12 (0.05, 0.19) | 0.12 (0.05, 0.19) | 0.12 (0.05, 0.19) |
| RD = -0.05: abx duration 6wk effect +silo specific d1 | |||||
| Surgical | 0 (−0.14, 0.14) | 0 (−0.12, 0.11) | −0.01 (−0.12, 0.1) | −0.01 (−0.11, 0.1) | −0.01 (−0.11, 0.1) |
| AB Duration | −0.04 (−0.26, 0.18) | −0.05 (−0.25, 0.16) | −0.05 (−0.24, 0.16) | −0.05 (−0.24, 0.16) | −0.05 (−0.24, 0.16) |
| AB Ext-proph | 0 (−0.15, 0.15) | 0 (−0.13, 0.13) | 0 (−0.13, 0.12) | −0.01 (−0.13, 0.12) | −0.01 (−0.13, 0.12) |
| AB Choice | 0 (−0.11, 0.11) | 0 (−0.09, 0.08) | 0 (−0.08, 0.07) | 0 (−0.08, 0.06) | 0 (−0.08, 0.06) |
| RD = 0.12 surgical, RD = 0.08 abx choice +silo specific d1 | |||||
| Surgical | 0.11 (−0.03, 0.24) | 0.11 (0, 0.22) | 0.11 (0, 0.22) | 0.11 (0, 0.22) | 0.12 (0.01, 0.21) |
| AB Duration | 0 (−0.2, 0.2) | 0 (−0.18, 0.18) | 0 (−0.17, 0.17) | 0 (−0.17, 0.17) | 0 (−0.17, 0.17) |
| AB Ext-proph | 0 (−0.15, 0.14) | 0 (−0.12, 0.11) | −0.01 (−0.11, 0.1) | −0.01 (−0.11, 0.1) | −0.01 (−0.11, 0.1) |
| AB Choice | 0.12 (0.02, 0.22) | 0.12 (0.04, 0.2) | 0.12 (0.05, 0.2) | 0.12 (0.06, 0.2) | 0.12 (0.06, 0.2) |
Table 5 shows the empirical proportion having treatment success within the randomised comparisons subsets averaged over the simulations. For example, the values for the surgical domain shows the empirical proportion having treatment success within the late acute silo by treatment arm (dair, rev(1), rev(2)).
d_tbl_4_cur <- dcast(d_tbl_4, scenario + desc + domain ~ arm, value.var = "p_hat")
d_tbl_4_cur[, domain := factor(
domain,
levels = c(1, 2, 3, 4),
labels = c("Surgical", "AB Duration", "AB Ext-proph", "AB Choice"))]
d_tbl_4_cur <- d_tbl_4_cur[, .(desc, domain, `1`, `2`, `3`)]
g_tbl <- d_tbl_4_cur |>
gt(groupname_col = "desc") |>
gt::text_transform(
locations = cells_row_groups(),
fn = function(x) {
lapply(x, function(x) {
gt::md(paste0("*", x, "*"))
})
}
) |>
cols_align(
columns = 1,
align = "left"
) |>
cols_align(
columns = 2:ncol(d_tbl_4_cur),
align = "center"
) |>
cols_label(
domain = "Domain"
) |>
tab_spanner(
label = html("Empirical risk by domain and treatment arm"),
columns = 2:ncol(d_tbl_4_cur)
) |>
tab_options(
table.font.size = "70%"
) |>
fmt_number(decimals = 2, drop_trailing_zeros = TRUE) |>
sub_missing(
columns = everything(),
rows = everything(),
missing_text = "-"
)
g_tbl
Empirical risk by domain and treatment arm
|
|||
|---|---|---|---|
| Domain | 1 | 2 | 3 |
| RD = 0 in all domains +silo specific d1 | |||
| Surgical | 0.59 | 0.62 | 0.57 |
| AB Duration | - | 0.65 | 0.65 |
| AB Ext-proph | - | 0.6 | 0.6 |
| AB Choice | - | 0.62 | 0.62 |
| RD = 0.12: surgical revision (one and two-stage) +silo specific d1 | |||
| Surgical | 0.56 | 0.71 | 0.67 |
| AB Duration | - | 0.7 | 0.7 |
| AB Ext-proph | - | 0.65 | 0.65 |
| AB Choice | - | 0.64 | 0.64 |
| RD = 0.12: surgical revision (one-stage only) +silo specific d1 | |||
| Surgical | 0.56 | 0.68 | 0.55 |
| AB Duration | - | 0.68 | 0.68 |
| AB Ext-proph | - | 0.59 | 0.59 |
| AB Choice | - | 0.61 | 0.61 |
| RD = 0.12: surgical revision (two-stage only) +silo specific d1 | |||
| Surgical | 0.56 | 0.6 | 0.68 |
| AB Duration | - | 0.63 | 0.63 |
| AB Ext-proph | - | 0.66 | 0.66 |
| AB Choice | - | 0.63 | 0.63 |
| RD = 0.12: abx duration 6wk effect +silo specific d1 | |||
| Surgical | 0.59 | 0.66 | 0.57 |
| AB Duration | - | 0.63 | 0.75 |
| AB Ext-proph | - | 0.6 | 0.6 |
| AB Choice | - | 0.62 | 0.62 |
| RD = 0.12: ext-proph 12wk effect +silo specific d1 | |||
| Surgical | 0.58 | 0.62 | 0.62 |
| AB Duration | - | 0.65 | 0.65 |
| AB Ext-proph | - | 0.58 | 0.7 |
| AB Choice | - | 0.63 | 0.63 |
| RD = 0.12: abx choice rif effect +silo specific d1 | |||
| Surgical | 0.62 | 0.65 | 0.61 |
| AB Duration | - | 0.69 | 0.69 |
| AB Ext-proph | - | 0.63 | 0.63 |
| AB Choice | - | 0.6 | 0.72 |
| RD = 0.12: all domains +silo specific d1 | |||
| Surgical | 0.6 | 0.79 | 0.76 |
| AB Duration | - | 0.71 | 0.83 |
| AB Ext-proph | - | 0.66 | 0.78 |
| AB Choice | - | 0.64 | 0.76 |
| RD = -0.05: abx duration 6wk effect +silo specific d1 | |||
| Surgical | 0.59 | 0.62 | 0.57 |
| AB Duration | - | 0.67 | 0.62 |
| AB Ext-proph | - | 0.6 | 0.6 |
| AB Choice | - | 0.62 | 0.62 |
| RD = 0.12 surgical, RD = 0.08 abx choice +silo specific d1 | |||
| Surgical | 0.6 | 0.75 | 0.71 |
| AB Duration | - | 0.74 | 0.74 |
| AB Ext-proph | - | 0.69 | 0.69 |
| AB Choice | - | 0.62 | 0.74 |