Stabilizing Thompson Sampling with Null Hypothesis Bayesian Response-Adaptive Randomization

Abstract

Response-adaptive randomization (RAR) methods can be used to adapt randomization probabilities based on accumulating data, aiming to increase the probability of allocating patients to effective treatments. A popular RAR method is Thompson sampling, which randomizes patients proportionally to the Bayesian posterior probability that each treatment is the most effective. However, its high variability can also increase the risk of assigning patients to inferior treatments and lead to inferential problems such as confidence interval undercoverage. We propose a principled method based on Bayesian hypothesis testing to address these issues: We introduce a null hypothesis postulating equal effectiveness of treatments. Bayesian model averaging then induces shrinkage toward equal randomization probabilities, with the degree of shrinkage controlled by the prior probability of the null hypothesis. Equal randomization and Thompson sampling arise as special cases when the prior probability is set to one or zero, respectively. Simulated and real-world examples illustrate that the method balances highly variable Thompson sampling with static equal randomization. A simulation study demonstrates that the method can mitigate issues with Thompson sampling and has comparable statistical properties to Thompson sampling with common ad hoc modifications such as power transformation and probability capping. We implement the method in the free and open-source R package brar, enabling experimenters to easily perform null hypothesis Bayesian RAR and support more effective randomization of patients.

Figures (7)

• Figure 1: Ternary plot where the color of each point depicts the treatment randomization probability $\pi$ for a particular combination of posterior probabilities of $H_{+}$ (bottom axis), $H_{-}$ (left axis), and $H_{0}$ (right axis). For example, the black asterisk denotes the combination where $\Pr(H_{+} \,\vert\, y) = 0.2$, $\Pr(H_{-} \,\vert\, y) = 0.3$, and $\Pr(H_0 \,\vert\, y) = 0.5$ with corresponding randomization probability $\pi = 45\%$.
• Figure 2: Illustration of spike-and-slab prior for the effect $\theta$. A point prior at 0 is assumed under $H_0$. A normal prior $\theta \sim \mathrm{N}(0, 1)$ with support truncated to the positive or negative side is assumed under $H_{+}$ and $H_{-}$, respectively. These priors are averaged assuming prior hypothesis probabilities $\Pr(H_0) = 0.5$, $\Pr(H_{+}) = 0.25$, and $\Pr(H_{-}) = 0.25$.
• Figure 3: Evolution of Bayesian RAR probabilities for three simulated data sequences (a). In each step, one patient is allocated to the treatment or control group using Thompson sampling assuming a normal prior centered at zero with standard deviation $\tau = 1$. Patient outcomes are then simulated from a normal distribution with a true standard deviation of 1 and assuming a true mean difference $\theta$ between treatment and control group as indicated in the plot panels. Randomization probabilities under the spike-and-slab prior ($\Pr(H_0) > 0$) are computed from the same data sequence to enable comparability. Plot (b) shows the distribution of randomization probabilities based on 1'000 simulations under Bayesian RAR for two different prior standard deviations and $\Pr(H_0)$.
• Figure 4: Illustration of a spike-and-slab prior for a two-dimensional effect $\boldsymbol{\theta} = (\theta_1, \theta_2)^\top$. A point mass prior at $(0, 0)^\top$ is assumed under $H_0$. A normal prior $\theta \sim \mathrm{N}((0, 0)^\top, \boldsymbol{\mathcal{T}})$ with $\boldsymbol{\mathcal{T}}_{ij} = 0.5$ for $i\neq j$ and $\boldsymbol{\mathcal{T}}_{ij} = 1$ for $i = j$, with support truncated to the space of the corresponding hypothesis is assumed under $H_{-}$, $H_{+1}$, and $H_{+2}$. The correlation of 0.5 ensures that all treatments receive equal prior probability $\Pr(H_{-}) = \Pr(H_{+1}) = \Pr(H_{+2}) = \{1 - \Pr(H_0)\}/3$.
• Figure 5: Evolution of Bayesian RAR randomization probabilities (a) and posterior probability of a beneficial ECMO treatment effect (b) for data from the ECMO trial. Plot (c) shows the probability of observing the ECMO allocation sequence as a function of $\Pr(H_0)$.