The R Bayesian evidence synthesis Tools (RBesT) have been created to facilitate the use of historical information in clinical trials. Once relevant historical information has been identified, RBesT supports the derivation of informative priors via the Meta-Analytic-Predictive (MAP) approach [1], the evaluation of the trial’s operating characteristics, and the data analysis with the actual trial data. RBesT has been developed for endpoints from a number of well known distributions. Here we consider an example for a normally distributed response variable.

Let’s look at the Crohn’s disease example data [2] (data-set
`crohn`

in RBesT). The primary endpoint is the change from
baseline in Crohn’s Disease Activity Index (CDAI), which is assumed to
be normally distributed. Note that for CDAI, an improved outcome
corresponds to a negative change from baseline.

First from historical studies we get the estimated standard deviation of the response variable as \(\sigma\) = 88, which is used to obtain the standard errors of the effect estimates.

```
<- crohn
dat <- 88
crohn_sigma $y.se <- crohn_sigma/sqrt(dat$n) dat
```

study | n | y | y.se |
---|---|---|---|

Gastr06 | 74 | -51 | 10.23 |

AIMed07 | 166 | -49 | 6.83 |

NEJM07 | 328 | -36 | 4.86 |

Gastr01a | 20 | -47 | 19.68 |

APhTh04 | 25 | -90 | 17.60 |

Gastr01b | 58 | -54 | 11.55 |

The MAP prior can be derived with the function
** gMAP**. The between-trial heterogeneity
parameter \(\tau\) governs how much
information will be shared from the historical trials into the design
and analysis of the future trials. In the normal case with a known
sampling standard deviation \(\sigma\),
the amount of borrowing from historical data depends on the ratio \(\tau/\sigma\). A conservative choice for
the prior on \(\tau\) is a

`HalfNormal(0,`

`)`

`?gMAP`

`set.seed`

```
library(RBesT)
set.seed(689654)
<- gMAP(cbind(y, y.se) ~ 1 | study,
map_mcmc weights=n,data=dat,
family=gaussian,
beta.prior=cbind(0, crohn_sigma),
tau.dist="HalfNormal",tau.prior=cbind(0,crohn_sigma/2))
print(map_mcmc)
```

```
## Generalized Meta Analytic Predictive Prior Analysis
##
## Call: gMAP(formula = cbind(y, y.se) ~ 1 | study, family = gaussian,
## data = dat, weights = n, tau.dist = "HalfNormal", tau.prior = cbind(0,
## crohn_sigma/2), beta.prior = cbind(0, crohn_sigma))
##
## Exchangeability tau strata: 1
## Prediction tau stratum : 1
## Maximal Rhat : 1
## Estimated reference scale : 88
##
## Between-trial heterogeneity of tau prediction stratum
## mean sd 2.5% 50% 97.5%
## 14.50 9.48 1.48 12.60 38.10
##
## MAP Prior MCMC sample
## mean sd 2.5% 50% 97.5%
## -49.9 18.8 -91.6 -48.6 -12.3
```

```
## a graphical representation is also available
<- plot(map_mcmc)
pl
## a number of plots are immediately defined
names(pl)
```

`## [1] "densityThetaStar" "densityThetaStarLink" "forest_model"`

```
## forest plot with model estimates
print(pl$forest_model)
```

Next, the MCMC MAP prior from the previous section is converted to a
parametric representation with the
** automixfit** function. This function fits a
parametric mixture representation using expectation-maximization (EM).
The number of mixture components is chosen automatically using AIC. One
can also specify the number of components for the mixture via

`mixfit`

`automixfit`

```
<- automixfit(map_mcmc)
map print(map)
```

```
## EM for Normal Mixture Model
## Log-Likelihood = -17037.29
##
## Univariate normal mixture
## Reference scale: 88
## Mixture Components:
## comp1 comp2 comp3
## w 0.64809468 0.28699434 0.06491098
## m -47.37713685 -51.26304658 -68.78804534
## s 9.36876807 30.45398192 6.12245218
```

```
## check accuracy of mixture fit
plot(map)$mix
```

The main advantage of using historical information is the possibility
to reduce the number of control patients, as the informative prior is
effectively equivalent to a certain number of control patients. This is
called the effective sample size (ESS) and can be calculated in RBesT
with the ** ess** function.

In the study protocol of this Crohn’s disease data example, the very conservative moment-based ESS of 20 was used to reduce the planned sample size of the control group.

`round(ess(map)) ## default elir method`

`## Using default prior reference scale 88`

`## [1] 39`

`round(ess(map, method="morita"))`

`## Using default prior reference scale 88`

```
## Warning in calc_loc(mix, "mode"): Detected multiple modes.
## The ESS is determined for the largest mode, but ESS concept is ill-defined for multi-modal distributions.
```

`## [1] 78`

`round(ess(map, method="moment"))`

`## Using default prior reference scale 88`

`## [1] 22`

We recommend robustifying [5] the prior with the
** robustify** function, which protects against
type-I error inflation in presence of prior-data conflict. For the
normal case we strongly recommend explicitly choosing the mean of the
robust component. We use \(-50\)
consistent with the mean of the MAP prior. Furthermore, 20% probability
is used for the additional robust (unit-information) mixture component.
The choice of such probability reflects the confidence about the
validitiy of the model assumptions, i.e. the possibility of a
non-exchangable control group to be enrolled per inclusion/exclusion
criteria in the current trial as compared to the historical control
group population. Note that robustification decreases the ESS.

```
## add a 20% non-informative mixture component
<- robustify(map, weight=0.2, mean=-50) map_robust
```

`## Using default prior reference scale 88`

`print(map_robust)`

```
## Univariate normal mixture
## Reference scale: 88
## Mixture Components:
## comp1 comp2 comp3 robust
## w 0.51847574 0.22959547 0.05192878 0.20000000
## m -47.37713685 -51.26304658 -68.78804534 -50.00000000
## s 9.36876807 30.45398192 6.12245218 88.00000000
```

`round(ess(map_robust)) `

`## Using default prior reference scale 88`

`## [1] 28`

Typically, operating characteristics are required to evaluate the
proposed design or compare a few design options. RBesT requires the
input of decision rules via ** decision2S**
function and then calculates the operating characteristics with the

`oc2S`

Consider this 2-arm design of placebo (with an informative prior) against an experimental treatment. The dual-criterion for success is defined as follows: \[ \begin{align*} \textrm{Criterion 1:} & \Pr(\theta_{act} - \theta_{pbo} \lt 0) &> 0.95 \\ \textrm{Criterion 2:} & \Pr(\theta_{act} - \theta_{pbo} \lt -50) &> 0.50. \end{align*} \]

Equivalently, the second criterion requires that the posterior median difference exceeds -50. The dual-criteria account for statistical significance as well as clinical relevance. Note that a negative change from baseline in CDAI corresponds to improvement.

```
## dual decision criteria
## pay attention to "lower.tail" argument and the order of active and pbo
<- decision2S(pc=c(0.95,0.5), qc=c(0,-50), lower.tail=TRUE)
poc print(poc)
```

```
## 2 sample decision function
## Conditions for acceptance:
## P(theta1 - theta2 <= 0) > 0.95
## P(theta1 - theta2 <= -50) > 0.5
## Link: identity
```

For the active group we use the same weakly informative (unit-information) prior as used in the robustification step of the MAP prior. Also, we set up the few design options with different choices of control prior and different sizes of control group.

```
## set up prior for active group
<- mixnorm(c(1,-50,1), sigma=crohn_sigma, param = 'mn')
weak_prior <- 40
n_act <- 20
n_pbo
## four designs
## "b" means a balanced design, 1:1
## "ub" means 40 in active and 20 in placebo
<- oc2S(weak_prior, weak_prior, n_act, n_act, poc,
design_noprior_b sigma1=crohn_sigma, sigma2=crohn_sigma)
<- oc2S(weak_prior, weak_prior, n_act, n_pbo, poc,
design_noprior_ub sigma1=crohn_sigma, sigma2=crohn_sigma)
<- oc2S(weak_prior, map, n_act, n_pbo, poc,
design_nonrob_ub sigma1=crohn_sigma, sigma2=crohn_sigma)
<- oc2S(weak_prior, map_robust, n_act, n_pbo, poc,
design_rob_ub sigma1=crohn_sigma, sigma2=crohn_sigma)
```

The type I can be increased compared to the nominal \(\alpha\) level in case of a conflict between the trial data and the prior. The robustified MAP prior can reduce the type I error inflation in this case to a lower level.

```
# the range for true values
<- seq(-120, -40, by=1)
cfb_truth
<- design_noprior_b(cfb_truth, cfb_truth)
typeI1 <- design_noprior_ub(cfb_truth, cfb_truth)
typeI2 <- design_nonrob_ub(cfb_truth, cfb_truth)
typeI3 <- design_rob_ub(cfb_truth, cfb_truth)
typeI4
<- rbind(data.frame(cfb_truth=cfb_truth, typeI=typeI1,
ocI design="40:40 with non-informative priors"),
data.frame(cfb_truth=cfb_truth, typeI=typeI2,
design="40:20 with non-informative priors"),
data.frame(cfb_truth=cfb_truth, typeI=typeI3,
design="40:20 with non-robust prior for placebo"),
data.frame(cfb_truth=cfb_truth, typeI=typeI4,
design="40:20 with robust prior for placebo")
)qplot(cfb_truth, typeI, data=ocI, colour=design, geom="line", main="Type I Error") +
xlab(expression(paste('True value of change from baseline ', mu[act] == mu[pbo]))) +
ylab('Type I error') +
coord_cartesian(ylim=c(0,0.2)) +
theme(legend.justification=c(1,1),legend.position=c(0.95,0.85))
```

The power shows the gain of using an informative prior for the control arm; i.e. 90% power is reached for smaller \(\delta\) values compared to the design with weakly informative priors for both arms or the balanced design.

```
<- seq(-80,0,by=1)
delta <- summary(map)["mean"]
m <- m + delta # active for 1
cfb_truth1 <- m + 0*delta # pbo for 2
cfb_truth2
<- design_noprior_b(cfb_truth1, cfb_truth2)
power1 <- design_noprior_ub(cfb_truth1, cfb_truth2)
power2 <- design_nonrob_ub(cfb_truth1, cfb_truth2)
power3 <- design_rob_ub(cfb_truth1, cfb_truth2)
power4
<- rbind(data.frame(cfb_truth1=cfb_truth1, cfb_truth2=cfb_truth2,
ocP delta=delta, power=power1,
design="40:40 with non-informative priors"),
data.frame(cfb_truth1=cfb_truth1, cfb_truth2=cfb_truth2,
delta=delta, power=power2,
design="40:20 with non-informative priors"),
data.frame(cfb_truth1=cfb_truth1, cfb_truth2=cfb_truth2,
delta=delta, power=power3,
design="40:20 with non-robust prior for placebo"),
data.frame(cfb_truth1=cfb_truth1, cfb_truth2=cfb_truth2,
delta=delta, power=power4,
design="40:20 with robust prior for placebo")
)
qplot(delta, power, data=ocP, colour=design, geom="line", main="Power") +
xlab('True value of difference (act - pbo)')+ ylab('Power') +
scale_y_continuous(breaks=c(seq(0,1,0.2),0.9)) +
scale_x_continuous(breaks=c(seq(-80,0,20),-70)) +
geom_hline(yintercept=0.9,linetype=2) +
geom_vline(xintercept=-70,linetype=2) +
theme(legend.justification=c(1,1),legend.position=c(0.95,0.85))
```

When the actual trial data are available, the final analysis can be
run with RBesT via the ** postmix** function.
The real data are used for this example, where the trial data led to a
negative conclusion. However, note that

`postmix`

```
## one can either use summary data or individual data. See ?postmix.
<- -29.2
y.act <- 14.0
y.act.se <- 39
n.act
<- -63.1
y.pbo <- 13.9
y.pbo.se <- 20
n.pbo
## first obtain posterior distributions
<- postmix(weak_prior, m=y.act, se=y.act.se)
post_act <- postmix(map_robust, m=y.pbo, se=y.pbo.se)
post_pbo
## then calculate probability for the dual criteria
## and compare to the predefined threshold values
<- pmixdiff(post_act, post_pbo, 0); print(p1) p1
```

`## [1] 0.06169579`

`<- pmixdiff(post_act, post_pbo, -50); print(p2) p2 `

`## [1] 3.623577e-06`

`print(p1>0.95 & p2>0.5)`

`## [1] FALSE`

```
## or we can use the decision function
poc(post_act, post_pbo)
```

`## [1] 0`

[1] Neuenschwander B et. al, *Clin Trials*. 2010;
7(1):5-18

[2] Hueber W. et. al, *Gut*, 2012, 61(12):1693-1700

[3] Kass RE, Wasserman L, *J Amer Statist Assoc*; 1995,
90(431):928-934.

[4] Morita S. et. al, *Biometrics* 2008;64(2):595-602

[5] Schmidli H. et. al, *Biometrics* 2014;70(4):1023-1032

`sessionInfo()`

```
## R version 4.1.0 (2021-05-18)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 20.04.5 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.9.0
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.9.0
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] bayesplot_1.8.1 tidyr_1.1.3 dplyr_1.0.8 ggplot2_3.3.5
## [5] RBesT_1.6-6 knitr_1.33
##
## loaded via a namespace (and not attached):
## [1] Rcpp_1.0.7 mvtnorm_1.1-2 prettyunits_1.1.1
## [4] ps_1.6.0 assertthat_0.2.1 digest_0.6.29
## [7] utf8_1.2.2 V8_3.4.2 plyr_1.8.6
## [10] R6_2.5.1 ggridges_0.5.3 backports_1.2.1
## [13] stats4_4.1.0 evaluate_0.14 highr_0.9
## [16] pillar_1.6.2 rlang_1.0.1 curl_4.3.2
## [19] callr_3.7.0 jquerylib_0.1.4 checkmate_2.0.0
## [22] rmarkdown_2.11 labeling_0.4.2 stringr_1.4.0
## [25] loo_2.4.1 munsell_0.5.0 compiler_4.1.0
## [28] xfun_0.25 rstan_2.21.2 pkgconfig_2.0.3
## [31] pkgbuild_1.2.0 rstantools_2.1.1 htmltools_0.5.2
## [34] tidyselect_1.1.1 tibble_3.1.3 gridExtra_2.3
## [37] codetools_0.2-18 matrixStats_0.60.1 fansi_0.5.0
## [40] crayon_1.4.2 withr_2.4.3 grid_4.1.0
## [43] jsonlite_1.7.2 gtable_0.3.0 lifecycle_1.0.1
## [46] DBI_1.1.2 magrittr_2.0.1 StanHeaders_2.21.0-7
## [49] scales_1.1.1 RcppParallel_5.1.4 cli_3.1.1
## [52] stringi_1.7.3 reshape2_1.4.4 farver_2.1.0
## [55] bslib_0.3.1 ellipsis_0.3.2 generics_0.1.0
## [58] vctrs_0.3.8 Formula_1.2-4 tools_4.1.0
## [61] glue_1.6.1 purrr_0.3.4 processx_3.5.2
## [64] parallel_4.1.0 fastmap_1.1.0 yaml_2.2.1
## [67] inline_0.3.19 colorspace_2.0-2 sass_0.4.0
```