# Constrained Markov decision processes for response-adaptive procedures in
clinical trials with binary outcomes

Stef Baas
Stochastic Operations Research Group, University of Twente, 7522 NH Enschede,
The Netherlands
Aleida Braaksma
Stochastic Operations Research Group, University of Twente, 7522 NH Enschede,
The Netherlands Richard J. Boucherie
Stochastic Operations Research Group, University of Twente, 7522 NH Enschede,
The Netherlands
Corresponding author

###### Abstract

A _constrained Markov decision process_ (CMDP) approach is developed for
response-adaptive procedures in clinical trials with binary outcomes. The
resulting CMDP class of Bayesian response-adaptive procedures can be used to
target a certain objective, e.g., patient benefit or power while using
constraints to keep other operating characteristics under control. In the CMDP
approach, the constraints can be formulated under different priors, which can
induce a certain behaviour of the policy under a given statistical hypothesis,
or given that the parameters lie in a specific part of the parameter space. A
solution method is developed to find the optimal policy, as well as a more
efficient method, based on backward recursion, which often yields a near-
optimal solution with an available optimality gap. Three applications are
considered, involving type I error and power constraints, constraints on the
mean squared error, and a constraint on prior robustness. While the CMDP
approach slightly outperforms the _constrained randomized dynamic programming_
(CRDP) procedure known from literature when focussing on type I and II error
and mean squared error, showing the general quality of CRDP, CMDP
significantly outperforms CRDP when the focus is on type I and II error only.

_Keywords_ Bayesian optimisation, Power constraints, Type I error control,
Mean squared error control, Prior misspecification control, adaptive treatment
allocation

##  1 Introduction

The current gold standard for assessing the efficacy of an experimental
treatment is the _randomized controlled trial_ (RCT) (Bhatt, 2010).
Historically, participants enrolled in an RCT are randomized to each
treatment, where the probability of allocation is independent of the history
of outcomes and allocations in the trial (the trial history). This approach
induces a high quality of statistical _operating characteristics_ (OCs) such
as the variance of the treatment effect estimator or power of statistical
tests, making it possible to make statements on how the experimental treatment
will perform for a new patient, not considered in the trial (extrapolation),
which is important in settings with large patient populations (Palmer, 1993).
Furthermore, independence of the trial history and allocations ensures
comparability of treatment groups and robustness to time trends. A downside of
this approach is that a large part of the patients included in the trial might
obtain inferior treatment. For a two-arm trial with equal allocation
probabilities, this could concern half the patients, while for a multi-arm
trial, this percentage can be even higher. The effect of this depends on the
severity and nature of the illness. For instance, the effect is greater when
the illness is very severe. Also, for rare diseases, a large part of the
patients having the disease will be included in the trial, and in this setting
extrapolation of trial results might be considered less important. In such
cases, it could be desirable to make a different trade-off between individual
ethics and collective ethics (Heilig and Weijerb, 2005).

An alternative, aiming to make such a trade-off, is a _response-adaptive_ (RA)
procedure (see, e.g., Rosenberger (1996), Antognini and Giovagnoli (2015), or
Villar and Jacko (2022)). The general idea behind an RA procedure, introduced
in Thompson (1933), is to sequentially base the allocation of the next
arriving participant to a treatment on the current trial history, aiming to
reach a certain objective, e.g., high patient benefit or good statistical OCs.
RA procedures have been proposed where the allocation probability is strictly
between zero and one for all trial histories, or where the allocation is
deterministic for at least one trial history. An RA procedure in the former
class is termed a _response-adaptive randomization_ (RAR) procedure, while an
RA procedure in the latter class is termed a _deterministic response-adaptive_
(DRA) procedure. A clinical trial using an RA procedure will from now on be
termed an RA design, while a design where treatment and outcomes are
independent will be termed a _non-adaptive_ (NA) design. In recent years, RA
procedures have found a growing number of applications in clinical trials
(aBerry and Viele, 2023), where most of these applications involved RAR
procedures. The US FDA has encouraged to consider multiple possibilities for
statistical trial design, including RAR procedures in US FDA (2019), while
this guidance document also highlights the controversies of RA procedures (for
an overview and comment, see Robertson et al. (2023)).

The current paper focuses on optimisation-based _Bayesian RA_ (BRA)
procedures. Optimisation-based BRA procedures focus on the exploration-
exploitation trade-off, which, in a clinical trial, corresponds to the aim of
learning which treatment is best while simultaneously allocating the highest
amount of patients to the best treatment. For a recent literature review
including optimisation-based BRA approaches see Williamson (2020, Chapter 2).
In the literature, BRA procedures can mainly be classified into three
categories:

  * •

Index-based approaches
For an index-based BRA procedure, an index value is independently determined
based on the data collected for each treatment group, and the treatment with
the highest index value has the highest probability of being allocated to the
next patient by the RA procedure. Examples are DRA approaches such as the
Gittins index procedure (Gittins, 2018), Bayes-UCB procedure (Kaufmann et al.,
2012), and RAR procedures such as semi-randomized index-based approaches
(Bather, 1981).

  * •

Thompson sampling (and modifications)
For Thompson sampling (introduced in Thompson (1933)) the treatment with the
highest posterior probability of having the highest expected outcome has the
highest probability of being allocated to the next patient. A modification of
Thompson sampling has been proposed in Thall and Wathen (2007) which induces a
smaller variance in the allocation probabilities.

  * •

Markov decision process (MDP) approaches
MDP approaches, introduced in Bradt et al. (1956), allocate treatment to trial
participants with the aim to maximize the patient benefit up to a fixed
horizon, given a prior distribution on the parameters for the model. The MDP
approach will be the main focus of the current paper.

Out of the three optimisation-based BRA procedures, the MDP approach is the
most computationally intensive, as shown in Villar et al. (2015). As
(modifications of) Thompson sampling can directly be applied in blocked
randomized designs, the second approach is currently the most common approach
when implementing RA procedures in practice (aBerry and Viele, 2023), while
the other two approaches have the potential to show the highest patient
benefit (see, e.g., Villar et al. (2015) or Williamson et al. (2017))

The MDP approach is a natural optimisation method for the setting of a
clinical trial, as every outcome is weighted equally and a finite trial
horizon is taken into account (Hardwick, 1995). A variety of different MDP BRA
procedures have been proposed in the literature, showing the flexibility of
this approach. In .Berry and Eick (1995), an MDP BRA procedure was introduced
which optimises, under a Bayesian model, the expected outcomes of participants
in the trial as well as the expected outcomes of a finite number of patients
which are all allocated to one treatment after the trial is completed. This
procedure hence makes an explicit trade-off between individual and collective
ethics. While RAR procedures can be constructed based on index-based BRA
procedures (Bather, 1981), the randomization component is not taken into
account in the optimisation. Cheng and Berry (2007) introduced an MDP Bayesian
RAR (BRAR) procedure where the randomization was included in the optimisation
in a natural way. Williamson et al. (2021) considered extensions of this
approach, including constraints on the minimum allocations to both treatment
groups, delays in outcomes, and random arrivals. Hardwick and Stout (1991)
introduced a multi-objective MDP BRA procedure, where the two objectives
revolved around patient benefit and best treatment selection. In Merrell et
al. (2022) an MDP BRA procedure was formulated for block-wise allocation,
optimising the block sizes, number of blocks, and treatment allocations for
blocks, where in the objective a trade-off is made between the expected number
of successes and OCs. Yi and Wang (2023) formulated an infinite-horizon MDP
for allocation in a clinical trial with general outcomes and showed that the
allocation ratio converges, ensuring the validity of likelihood-based tests.

The current paper introduces a novel class of constrained MDP BRA procedures,
CMDP procedures for short. The CMDP procedures follow from a modified version
of a constrained Markov decision process (Altman, 1999), where the
modification is that it is possible to use different expectation operators in
each constraint. This modification is needed to optimise the allocation of
treatment while keeping, e.g., type I error and power under control. In MDP
BRA procedures, optimisation is often performed under a vague prior, while the
probability distribution of the outcomes under, for instance, a type I error
constraint can have strong assumptions on the parameters such as equality of
the expectation of the outcomes, hence different prior distributions are
needed to formulate Bayesian type I and power constraints. A CMDP procedure
can be made to satisfy requirements on OCs for the trial using constraints
while optimising patient benefit. In comparison to previous methods from the
literature, the CMDP approach can be considered a more natural way to directly
impose a desired behaviour for an RA procedure while ensuring that, given the
imposed constraints, the obtained policy is optimal, which makes it possible
to reach a higher patient benefit. An additional advantage over penalized
methods is that adding constraints will have an effect that is known
beforehand, namely it will shrink the feasible region for the policies,
whereas the effect of adding a penalty term to the objective is less clear
beforehand.

Throughout the current paper we make use of an efficient implementation of
backward and forward recursion following Jacko (2019), where we make use of a
conservation law for the states, use a storage mapping function to store
values efficiently, and overwrite elements of the value function not used
further in the algorithm. This makes it possible to compute the MDP BRA
procedures considered in this paper in a relatively short amount of time.
Furthermore, this makes it possible to compute the values of clinical trial
OCs directly instead of approximating them by simulation, avoiding Monte Carlo
error.

The current paper is structured as follows. Section 2 introduces the model of
a binary two-arm clinical trial where outcomes are collected using a response-
adaptive procedure, as well as relevant operating characteristics. Section 3
introduces the class of CMDP procedures, as well as an algorithmic method for
determining CMDP procedures. In section 4 CMDP procedures are constructed that
optimise patient benefit under a restriction on the power and type I error, a
restriction on the mean squared error, and also on a restriction on robustness
to prior misspecification. Section 5 concludes the paper and gives directions
for future research.

##  2 Model and operating characteristics

###  2.1 Two-arm Response adaptive design with binary outcomes

We consider a trial, in which there are two treatments (arms) with unknown
outcome distributions, the control and developmental treatment. Trial
participants are sequentially allocated to a treatment using a response-
adaptive (RA) procedure which, given the current trial history, determines the
probability that the next participant obtains a given treatment. After
allocation, the (binary) outcome of the participant, sampled from a Bernoulli
distribution, becomes available and is added to the trial history before
allocating the next participant.

We now make the above formal. Let
$\text{\boldmath$\mathbf{\theta}$}=(\theta_{\text{C}},\theta_{\text{D}})\in[0,1]^{2}$
be a tuple of (unknown) success probabilities, where C denotes the control
treatment and D denotes the developmental treatment. In the following, the
same convention (i.e., first C then D) will be used to construct tuples from
variables for the control and developmental treatment. For a fixed trial size
$n\in\mathbb{N}$, let
$\bm{Y}=(Y_{a,t})_{a\in\\{\text{C},\text{D}\\},t\in\\{1,\dots,n\\}}$ be a
sequence of independent Bernoulli random variables, where
$\mathbb{P}_{\text{\boldmath$\mathbf{\theta}$}}(Y_{a,t}=1)=\theta_{a}$ for
$a\in\\{\text{C},\text{D}\\}$. The random variable $Y_{a,t}$ denotes a
potential outcome for trial participant $t$ under treatment $a$. Let
$\mathcal{H}=\bigcup_{t=0}^{\infty}\mathcal{H}_{t}$ where
$\mathcal{H}_{0}=\\{()\\}$ only contains the empty tuple and
$\mathcal{H}_{t}=\\{(a_{1},y_{1},a_{2},y_{2},\dots,a_{t},y_{t}):y_{w}\in\\{0,1\\},\,a_{w}\in\\{\text{C},\text{D}\\},\;\forall
w\leq t\\}$ is the set of possible trial histories of outcomes and actions up
to participant $t$. An RA procedure $\pi:\mathcal{H}\mapsto[0,1]$ maps a trial
history to the probability that the subsequent participant is allocated to the
control treatment. Trial participants are allocated sequentially, after which
the outcome for that trial participant is observed before allocating the next
participant, i.e., letting $\bm{H}_{0}=()$ we recursively define the realised
trial history as
$\bm{H}_{t}=(A_{1},Y_{A_{1},1},A_{2},Y_{A_{2},2},\dots,A_{t},Y_{A_{t},t})$,
where each $A_{t}$ is an independently drawn Bernoulli random variable with
$\mathbb{P}(A_{t}=1)=\pi(\bm{H}_{t-1})$ for $t=1,\dots,n$.

Let $S_{a,t}$ be the recorded number of successes, and $N_{a,t}$ denote the
number of allocations for arm $a$ up to time $t$, i.e.,

$S_{a,t}=\sum_{u=1}^{t}Y_{a,u}\mathbb{I}(A_{u}=a),\indent
N_{a,t}=\sum_{u=1}^{t}\mathbb{I}(A_{u}=a),\indent\forall
a\in\\{\text{C},\text{D}\\},\;t\in\\{1,\dots,n\\},~{}$

where $\mathbb{I}$ denotes the indicator function. Let
$\bm{X}_{t}=(\bm{S}_{t},\bm{N}_{t})$ be a state variable containing the
successes and allocations for each arm up to time $t$. Letting
$\bm{X}_{0}=((0,0),(0,0))$, the state space for
$\bm{X}=(\bm{X}_{t})_{t=0}^{n}$ is $\mathcal{X}=\cup_{t=0}^{n}\mathcal{X}_{t}$
where for all $t$:

$\mathcal{X}_{t}=\\{((x_{11},x_{12}),(x_{21},x_{22})):x_{ij}\in\\{0,\dots,t\\},\,x_{1j}\leq
x_{2j},\;x_{21}+x_{22}=t,\;\forall i,j\in\\{1,2\\}\\}.$

Let $\partial\bm{s}_{\text{C}}=((1,0),(1,0))$ and
$\partial\bm{f}_{\text{C}}=((0,0),(1,0))$ be the change in $\bm{X}_{t}$ after
a success and failure for the control arm, and let
$\partial\bm{s}_{\text{D}}$, $\partial\bm{f}_{D}$ be defined similarly.
Letting
$q(\bm{x}_{t+1},\bm{x}_{t},a)=\mathbb{P}_{\text{\boldmath$\mathbf{\theta}$}}(\bm{X}_{t+1}=\bm{x}_{t+1}\mid\bm{X}_{t}=\bm{x}_{t},A_{t}=a)$,
we have for  all $t\in\\{0,\dots,n\\}$,
$\bm{x}_{t}\in\mathcal{X}_{t},\bm{x}_{t+1}\in\mathcal{X}_{t+1}$,
$a\in\\{\text{C},\text{D}\\}$, and
$\text{\boldmath$\mathbf{\theta}$}\in[0,1]^{2}$

$q(\bm{x}_{t+1},\bm{x}_{t},a)=\begin{cases}\theta_{a},\indent&\text{
if~{}$\bm{x}_{t+1}=\bm{x}_{t}+\partial\bm{s}_{a}$,}\\\
(1-\theta_{a}),\indent&\text{
if~{}$\bm{x}_{t+1}=\bm{x}_{t}+\partial\bm{f}_{a}$,}\end{cases}$ (1)

where the addition of tuples is understood to be element-wise. If the RA
procedure $\pi$ is a function of $\bm{X}_{t}$, i.e, can be written as a
function $\pi:\mathcal{X}\mapsto[0,1]$, the process $\bm{X}=(\bm{X}_{t})_{t}$
is a Markov process with transition structure

$\mathbb{P}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}(\bm{X}_{t+1}=\bm{x}_{t+1}\mid\bm{X}_{t}=\bm{x}_{t})=\sum_{a\in\\{\text{C},\text{D}\\}}q(\bm{x}_{t+1},\bm{x}_{t},a)\pi(\bm{x}_{t})^{\iota_{a}}(1-\pi(\bm{x}_{t}))^{1-\iota_{a}}$

for all $\bm{x}_{t}\in\mathcal{X}_{t},\bm{x}_{t+1}\in\mathcal{X}_{t+1}$, where
$\mathbb{P}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}$ denotes the probability
measure on states induced by the RA procedure $\pi$ and (1) and
$\iota_{a}=\mathbb{I}(a=\text{C})$ for all $a\in\\{\text{C},\text{D}\\}$. In
this case, $\pi$ is a Markov RA procedure (Yi, 2013).

###  2.2 Definition of operating characteristics

In this section, we define important OCs that can be calculated based on the
model introduced at the start of Section 2.

First, at decision epoch $t=n$, denoted the trial horizon, we test for a
treatment effect $\theta_{\text{D}}-\theta_{\text{C}}$, i.e., we test

$H_{0}:\theta_{\text{D}}-\theta_{\text{C}}=0\indent\text{v.s.}\indent
H_{1}:\theta_{\text{D}}-\theta_{\text{C}}\neq 0.$

The test of choice is often Fisher’s exact test which, for a significance
level $\alpha\in(0,1)$, rejects when $\mathrm{T}(\bm{X}_{n})\leq\alpha$,
where, letting ${s}_{a}(\bm{x}_{t})$ and ${n}_{a}(\bm{x}_{t})$ be the number
of successes and allocations encoded in $\bm{x}_{t}$ and
$s(\bm{x})=s_{\text{C}}(\bm{x})+s_{\text{D}}(\bm{x})$ for all
$\bm{x}\in\mathcal{X}$,

$\indent\mathrm{T}(\bm{X}_{n})=\sum_{\begin{subarray}{c}\bm{x}^{\prime}_{n}\in\mathcal{X}_{n}\\\
P(\bm{x}^{\prime}_{n})\leq
P(\bm{x}_{n})\end{subarray}}P(\bm{x}^{\prime}_{n}),\indent\text{ and }\indent
P(\bm{x}_{n})=\frac{\binom{n_{\text{C}}(\bm{x}_{n})}{s_{\text{C}}(\bm{x}_{n})}\binom{n_{\text{D}}(\bm{x}_{n})}{s_{\text{D}}(\bm{x}_{n})}}{\binom{n}{s(\bm{x}_{n})}}\indent\forall\bm{X}_{n}\in\mathcal{X}_{n}.$

This test is exact, i.e., the type I error is bounded by (and as close as
possible up to a discreteness error to) the significance level $\alpha$ when
$\pi(\bm{x})=1/2$ for all $\bm{x}\in\mathcal{X}$ (Agresti, 1992).

According to the above testing situation, we can define four OCs, which depend
on the parameters $\mathbf{\theta}$, where we let
$\mathbb{E}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}$ denote the expectation
w.r.t. $\mathbb{P}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}$.

  * •

Patient Benefit:
The patient benefit is calculated as

$\mathbb{E}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}[n_{\text{C}}(\bm{X}_{n})/n]\cdot\mathbb{I}(\theta_{\text{C}}>\theta_{\text{D}})+\mathbb{E}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}[1-n_{\text{C}}(\bm{X}_{n})/n]\cdot\mathbb{I}(\theta_{\text{D}}>\theta_{\text{C}})+1/2\cdot\mathbb{I}(\theta_{\text{C}}=\theta_{\text{D}}).$

This OC represents the patient benefit in the trial given the parameters,
which we want to be high.

  * •

Rejection rate (RR):
This OC equals
$\mathbb{P}_{\text{\boldmath$\mathbf{\theta}$}}^{\pi}(\mathrm{T}(\bm{X}_{n})\leq\alpha)$,
i.e., the probability of rejecting $H_{0}$ at the end of the trial. If
$\theta_{\text{C}}=\theta_{\text{D}}$, this probability equals the type I
error and we want the rejection rate to be less than $\alpha$. If
$\theta_{\text{C}}\neq\theta_{\text{D}}$, this OC equals the power and we want
the rejection rate to be high.

  * •

Bias:
This OC equals
$\mathbb{E}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}[\hat{\theta}_{\text{D}}(\bm{X}_{n})-\hat{\theta}_{\text{C}}(\bm{X}_{n})]-\theta_{\text{D}}-\theta_{\text{C}},$
where for all $\bm{x}_{n}\in\mathcal{X}_{n}$, $a\in\\{\text{C},\text{D}\\}$

$\hat{\theta}_{a}(\bm{x}_{n})=\begin{cases}s_{a}(\bm{x}_{n})/n_{a}(\bm{x}_{n}),\indent&\text{if~{}$\min_{a}n_{a}(\bm{x}_{n})>0$,}\\\
(s_{a}(\bm{x}_{n})+1)/(n_{a}(\bm{x}_{n})+2),&\text{else.}\end{cases}$

The above adjustment to the maximum likelihood estimator
$s_{a}(\bm{x}_{n})/n_{a}(\bm{x}_{n})$ is made in order to provide an estimate
when either of the treatment groups contains zero observations. The OC equals
the expected error in the estimated treatment effect at the end of the trial,
and we want the bias to be as close to zero as possible.

  * •

Mean squared error (MSE):
This OC equals
$\mathbb{E}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}[(\hat{\theta}_{\text{D}}(\bm{X}_{n})-\hat{\theta}_{\text{C}}(\bm{X}_{n})-(\theta_{\text{D}}-\theta_{\text{C}}))^{2}]$
and corresponds to the quality of the treatment effect estimate
$\hat{\theta}_{\text{D}}(\bm{X}_{n})-\hat{\theta}_{\text{C}}(\bm{X}_{n})$, we
want this OC to be low.

##  3 Constrained Markov decision processes

In this section, we introduce the class of CMDP procedures. The class of CMDP
procedures is based on constrained optimisation. Trial participants are
sequentially allocated treatment, where the probability of allocating a
treatment is determined based on a trade-off between exploration and
exploitation, while the resulting distribution over states should also be such
that certain constraints are satisfied. It is allowed to formulate the
constraint under a different transition structure. This makes it possible to
define, e.g., type I error or power constraints, where optimisation can be
performed under a model that does not assume any prior information on
$\mathbf{\theta}$, while the kernel under the null hypothesis would restrict
$\theta_{\text{C}}=\theta_{\text{D}}$, and the kernel under the alternative
hypothesis could, e.g., assume $|\theta_{\text{D}}-\theta_{\text{C}}|\geq\xi$
for some $\xi\in(0,1).$

###  3.1 Formulation of the optimisation problem

We determine a CMDP procedure according to a Bayesian optimisation problem,
where the expected number of successes is maximized under a prior predictive
distribution, given a set of policy constraints. We assume an independent
$\text{Beta}(\tilde{s}_{a,0},\tilde{f}_{a,0})$ prior $\Pi$ for each arm
$a\in\\{\text{C},\text{D}\\}$, where $\text{Beta}(\cdot,\cdot)$ denotes the
Beta distribution, and $\tilde{s}_{a,0},\tilde{f}_{a,0}\in(0,\infty)$ are a
prior number of successes and failures. The
$\text{Beta}(\tilde{s}_{a,0},\tilde{f}_{a,0})$ prior is a conjugate prior for
the Bernoulli distribution, which brings the advantage that the posterior
distribution for the success probability is known in closed form. As in Cheng
and Berry (2007), for a CMDP procedure, the actions $\delta_{t}\in[1-p,p]$
correspond to the probability of allocating the next participant to the
control arm, where $p\in[1/2,1]$ is the degree of randomization. Let
$\tilde{s}_{a}(\bm{x}_{t})=s_{a}(\bm{x}_{t})+\tilde{s}_{a,0}$ and
$\tilde{n}_{a}(\bm{x}_{t})=n_{a}(\bm{x}_{t})+\tilde{s}_{a,0}+\tilde{f}_{a,0}$
for each arm $a\in\\{\text{C},\text{D}\\}$. The transition dynamics (1) can be
made independent of $\mathbf{\theta}$ by taking the expectation of the
transition probabilities with respect to the prior, in which case, following
(Williamson et al., 2017, Section 2.3) and  letting
$q(\bm{x}_{t+1},\bm{x}_{t},\delta_{t})=\mathbb{P}(\bm{X}_{t+1}=\bm{x}_{t+1}\mid\bm{X}_{t}=\bm{x}_{t},\delta_{t})$:

$q(\bm{x}_{t+1},\bm{x}_{t},\delta_{t})=\begin{cases}\delta_{t}\cdot\tilde{s}_{\text{C}}(\bm{x}_{t})/\tilde{n}_{\text{C}}(\bm{x}_{t}),\indent&\text{
if~{}$\bm{x}_{t+1}={\bm{x}}_{t}+\partial\bm{s}_{\text{C}}$,}\\\
\delta_{t}\cdot(1-\tilde{s}_{\text{C}}(\bm{x}_{t})/\tilde{n}_{\text{C}}(\bm{x}_{t})),\indent&\text{
if~{}${\bm{x}}_{t+1}={\bm{x}}_{t}+\partial\bm{f}_{\text{C}}$,}\\\
(1-\delta_{t})\cdot\tilde{s}_{\text{D}}(\bm{x}_{t})/\tilde{n}_{\text{D}}(\bm{x}_{t}),\indent&\text{
if~{}${\bm{x}}_{t+1}={\bm{x}}_{t}+\partial\bm{s}_{\text{D}},$}\\\
(1-\delta_{t})\cdot(1-\tilde{s}_{\text{D}}(\bm{x}_{t})/\tilde{n}_{\text{D}}(\bm{x}_{t})),\indent&\text{
if~{}${\bm{x}}_{t+1}={\bm{x}}_{t}+\partial\bm{f}_{\text{D}}.$}\\\ \end{cases}$
(2)

For all $\bm{x}_{n}\in\mathcal{X}_{n}$, $t<n$,
$\delta_{t},\delta_{n}\in[1-p,p]$, and $\bm{x}_{t}\in\mathcal{X}_{t}$ let

$r({\bm{x}}_{t},\delta_{t})=\delta_{t}\tilde{s}_{\text{C}}(\bm{x}_{t})/\tilde{n}_{\text{C}}(\bm{x}_{t})+(1-\delta_{t})\tilde{s}_{\text{D}}(\bm{x}_{t})/\tilde{n}_{\text{D}}(\bm{x}_{t}),\indent
r({\bm{x}}_{n},\delta_{n})=0,$

be the posterior mean rewards after choosing allocation probability
$\delta_{t}$ for the control treatment in state ${\bm{x}}_{t}$, for $t\leq n$.
A CMDP procedure now maximizes the total expected reward, i.e., the expected
sum of $r(\bm{X}_{t},\delta_{t})$ over decision epochs.

Constraints are enforced ensuring good OCs for the resulting procedure. The
constraints are defined in terms of total expected reward under an alternative
prior for $\mathbf{\theta}$, i.e., defining priors $\Pi_{c}$ for
$\mathbf{\theta}$ for $c$ in a finite countable index set $\mathcal{C}$,
transition probabilities
$\mathbb{P}_{c}(\bm{X}_{t+1}=\bm{x}_{t+1}\mid\bm{X}_{t}=\bm{x}_{t},\,\delta_{t})$
can be defined according to (1) by integrating out $\mathbf{\theta}$ w.r.t.
$\Pi_{c}$. Defining rewards $r_{c}(\bm{x}_{t},\delta_{t})$, and values
$V_{c}\in\mathbb{R}$ for all
$\bm{x}_{t+1},\bm{x}_{t}\in\mathcal{X},\;\delta_{t}\in[1-p,p]$ for
$c\in\mathcal{C}\subseteq\mathbb{N}_{0}$, a CMDP procedure is defined as the
solution of

$\displaystyle\underset{\pi}{\max}\,\mathbb{E}^{\pi}\left[\sum_{t=0}^{n}r({\bm{X}}_{t},\delta_{t})\right]$
(3a) s.t.
$\displaystyle\mathbb{E}^{\pi}_{c}\left[\sum_{t=0}^{n}r_{c}({\bm{X}}_{t},\delta_{t})\right]\leq
V_{c}\indent\forall c\in\mathcal{C}.$ (3b)

Problem (3) will be referred to as the CMDP problem. As the Constraints (3b)
are specified under expectation operators different from the expectation
operator in the Objective (3a), Problem (3) is a generalization of a finite
horizon constrained Markov decision problem as defined in Altman (1999). We
now give three examples of CMDP procedures:

.

######  Example 1 (Control of power and type I error).

A downside of RA procedures is that reliable statistical inference is
complicated due to an imbalance and variance over treatment group sizes, as
well as dependence in the outcomes induced by the RA procedure, which can lead
to a reduction of power and/or type I error inflation in comparison to a NA
design (Robertson et al., 2023). To test a certain hypothesis on the
difference $\theta_{\text{D}}-\theta_{\text{C}}$ after collecting data using
the CMDP procedure, the constraints (3b) can be specified such that the
Bayesian expected power and type I error remain under control. This idea has
similarities with Bayesian sample size determination, where Bayesian expected
power and type I error are used to determine the sample size of a clinical
trial based on historical data (Brown et al., 1987).

We propose to control the type I error by defining a prior $\Pi_{0}$ with
support
$\Theta_{0}=\\{\text{\boldmath$\mathbf{\theta}$}\in[0,1]^{2}:\theta_{\text{C}}=\theta_{\text{D}}\\}$,
agreeing with the null hypothesis $H_{0}$. In this evaluation, we let
$\Pi_{0}$ correspond to a $\text{Beta}(\tilde{s}_{0},\tilde{f}_{0})$ prior on
$\theta_{\text{C}}$ with the additional restriction
$\theta_{\text{C}}=\theta_{\text{D}},$ where
$\tilde{s}_{0},\tilde{f}_{0}\in(0,\infty)$. As in (2), we integrate out the
parameter $\mathbf{\theta}$ in (1) but now using the prior $\Pi_{0}$. Letting
$\tilde{s}(\bm{x}_{t})=s(\bm{x}_{t})+\tilde{s}_{0}$,
$n(\bm{x}_{t})=n_{\text{C}}(\bm{x}_{t})+n_{\text{D}}(\bm{x}_{t})$,
$\tilde{n}(\bm{x}_{t})=n(\bm{x}_{t})+\tilde{s}_{0}+\tilde{f}_{0}$, and
$q_{0}(\bm{x}_{t+1},\bm{x}_{t},\delta_{t})=\mathbb{P}_{0}(\bm{X}_{t+1}=\bm{x}_{t+1}\mid\bm{X}_{t}=\bm{x}_{t},\delta_{t})$
we obtain for all
${\bm{x}}_{t}\in{\mathcal{X}_{t}},{\bm{x}}_{t+1}\in{\mathcal{X}_{t+1}}$

$q_{0}({\bm{x}}_{t+1}{\bm{x}}_{t},\delta_{t})=\begin{cases}\delta_{t}\cdot\tilde{s}(\bm{x}_{t})/\tilde{n}(\bm{x}_{t}),\indent&\text{
if~{}$\bm{x}_{t+1}=\bm{x}_{t}+\partial\bm{s}_{\text{C}}$,}\\\
\delta_{t}\cdot(1-\tilde{s}(\bm{x}_{t})/\tilde{n}(\bm{x}_{t})),\indent&\text{
if~{}$\bm{x}_{t+1}=\bm{x}_{t}+\partial\bm{f}_{\text{C}}$,}\\\
(1-\delta_{t})\cdot\tilde{s}(\bm{x}_{t})/\tilde{n}(\bm{x}_{t}),\indent&\text{
if~{}$\bm{x}_{t+1}=\bm{x}_{t}+\partial\bm{s}_{\text{D}}$,}\\\
(1-\delta_{t})\cdot(1-\tilde{s}(\bm{x}_{t})/\tilde{n}(\bm{x}_{t})),\indent&\text{
if~{}$\bm{x}_{t+1}=\bm{x}_{t}+\partial\bm{f}_{\text{D}}.$}\\\ \end{cases}$

We enforce the type I error constraint

$\mathbb{P}^{\pi}_{0}(\mathrm{T}(\bm{X}_{n})\leq\alpha)=\int_{[0,1]^{2}}\mathbb{P}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}(\mathrm{T}(\bm{X}_{n})\leq\alpha)\Pi_{0}(d\text{\boldmath$\mathbf{\theta}$})\leq\alpha^{*}.$
(4)

We see from (4) that the constraint enforces that the weighted average of the
frequentist type I error
$\mathbb{P}^{\pi}_{\text{\boldmath$\mathbf{\theta}$}}(\mathrm{T}(\bm{X}_{n})\leq\alpha)$
is bounded by $\alpha^{*}$, hence to ensure
$\mathbb{P}^{\pi}_{0}(\mathrm{T}(\bm{X}_{n})\leq\alpha\mid\text{\boldmath$\mathbf{\theta}$})\leq\alpha$
for all $\mathbf{\theta}$, it is necessary that $\alpha^{*}\leq\alpha$ above.

Similarly, we can define a prior $\Pi_{1}$ with support in $[0,1]^{2}$, and,
given a required maximum type II error $\beta\in(0,1)$, the power constraint
becomes

$\mathbb{P}^{\pi}_{1}(\mathrm{T}(\bm{X}_{n})\leq\alpha)\geq 1-\beta.$

The resulting optimisation problem, the solution of which is denoted as CMDP-T
procedure, can be written as

$\displaystyle\underset{\pi}{\max}\,\mathbb{E}^{\pi}\left[\sum_{t=0}^{n}r(\bm{X}_{t},\delta_{t})\right]$
(5a) s.t.
$\displaystyle\mathbb{P}^{\pi}_{0}(\mathrm{T}(\bm{X}_{n})\leq\alpha)\leq\alpha^{*},$
(5b) $\displaystyle\mathbb{P}_{1}^{\pi}(\mathrm{T}(\bm{X}_{n})\leq\alpha)\geq
1-\beta.$ (5c)

Note that (5) corresponds to (3) with $r_{c}(\bm{x}_{t},\delta_{t})=0$ for
$\bm{x}_{t}\in\mathcal{X}_{t}$ and $t<n$, $V_{0}=\alpha^{*}$,
$V_{1}=-(1-\beta)$ and
$r_{c}(\bm{x}_{n},\delta_{n})=(-1)^{\mathbb{I}(c=1)}\mathbb{I}(\mathrm{T}(\bm{X}_{n})\leq\alpha)$
for $\bm{x}_{n}\in\mathcal{X}_{n}$.

.

######  Example 2 (Control of estimation error).

Next to difficulties in testing, using an RA procedure can also lead to large
errors when estimating treatment effects (Bowden and Trippa, 2017). To improve
the estimate of the treatment effect $\theta_{\text{D}}-\theta_{\text{C}}$
after collecting data using the CMDP procedure, we specify the constraints
(3b) such that the MSE remains under control.

We optimise the successes incurred under the
$\text{Beta}(\tilde{s}_{a,0},\tilde{f}_{a,0})$ prior introduced in Section 2,
under the constraint that the posterior MSE is small. In order to control the
shape of the MSE curve over $[0,1]^{2}$, we discretize $[0,1]^{2}$ into a
collection of disjoint two-dimensional intervals

$\mathcal{S}=\\{[\theta_{\text{C}}^{\ell},\theta_{\text{C}}^{u})\times[\theta_{\text{D}}^{\ell},\theta_{\text{D}}^{u}):0\leq\theta_{a}^{\ell}\leq\theta_{a}^{u}\leq
1\;\forall a\in\\{\text{C},\text{D}\\}\\}$

such that $\cup_{\sigma\in\mathcal{S}}\sigma=[0,1)^{2}$ and define
$\Pi_{\sigma}$ to be a $\text{Beta}(\tilde{s}^{E}_{a,0},\tilde{f}^{E}_{a,0})$
prior truncated to $\sigma$ for all $\sigma\in\mathcal{S}$,
$a\in\\{\text{C},\text{D}\\}.$ By the law of total expectation, letting
$\Pi_{\sigma}(\text{\boldmath$\mathbf{\theta}$}\mid\bm{X}_{n})$ denote the
posterior distribution under prior $\Pi_{\sigma}$, the posterior MSE for a
policy $\pi$ over $\sigma$ in $\mathcal{S}$ can be expressed as

$\mathbb{E}_{\sigma}^{\pi}[\textstyle\int_{\sigma}(\hat{\theta}_{\text{D}}(\bm{X}_{n})-\hat{\theta}_{\text{C}}(\bm{X}_{n})-(\theta_{\text{D}}-\theta_{\text{C}}))^{2}d\Pi_{\sigma}(\text{\boldmath$\mathbf{\theta}$}\mid\bm{X}_{n})].$
(6)

The expression in (6) is the expectation of a function of the final state
$\bm{X}_{n}$, which does not depend on the policy as all allocations have been
realised, hence a bound on this quantity can be written as (3b). To see
whether policies with high patient benefit, type I error control and high
power, as well as low MSE can be found, we furthermore add the power and type
I error constraints as in the CMDP-T formulation (5). The resulting
optimisation problem, the solution of which is denoted as CMDP-E RA procedure,
can be written as

$\displaystyle\underset{\pi}{\max}\,\mathbb{E}^{\pi}\left[\sum_{t=0}^{n}r(\bm{X}_{t},\delta_{t})\right]$
(7a) s.t.
$\displaystyle\mathbb{P}^{\pi}_{0}(\mathrm{T}(\bm{X}_{n})\leq\alpha)\leq\alpha^{*},$
(7b) $\displaystyle\mathbb{P}_{1}^{\pi}(\mathrm{T}(\bm{X}_{n})\leq\alpha)\geq
1-\beta,$ (7c)
$\displaystyle\mathbb{E}_{\sigma}^{\pi}[\textstyle\int_{\sigma}(\hat{\theta}_{\text{D}}(\bm{X}_{n})-\hat{\theta}_{\text{C}}(\bm{X}_{n})-(\theta_{\text{D}}-\theta_{\text{C}}))^{2}d\Pi_{\sigma}(\text{\boldmath$\mathbf{\theta}$}\mid\bm{X}_{n})]\leq
V_{\sigma},\indent\forall\sigma\in\mathcal{S}.$ (7d)

The priors under the constraints (7d) are designed in such a way that the
constraints reflect the average behaviour of the policy on a specific part of
the parameter space $[0,1]^{2}$.

.

######  Example 3 (Robustness to prior misspecification).

In clinical trials, historical data and expert opinion are often available on
one or both treatments. In Bayesian RA procedures, historical data can be used
to construct a prior distribution on the success probabilities, which may be
leveraged to increase patient benefit. However, if the actual treatment
effects have low probability mass under the prior, i.e., when there is prior
misspecification, the number of participants allocated to the optimal arm can
be lower than under a less informative prior. The constraints (3b) can be
specified such that the RA procedure is robust against prior misspecification,
by controlling the number of treatment failures $1-Y_{A_{t},t}$ under a less
informative prior distribution.

We assume that the incorporation of historical data and expert opinion leads
to an independent $\text{Beta}(\tilde{s}_{a,0},\tilde{f}_{a,0})$ prior on the
success probability for each arm $a$ in the trial, corresponding to the prior
$\Pi$ used for optimisation in (3). Let $\Pi_{\text{LI}}$ be a prior on the
success probabilities that is less informative, e.g., a uniform prior, let

$r_{\text{LI}}(\bm{x}_{t},\delta_{t})=-\delta_{t}\mathbb{E}_{\text{LI}}[\theta_{\text{C}}\mid\bm{X}_{t}=\bm{x}_{t}]+(1-\delta_{t})\mathbb{E}_{\text{LI}}[\theta_{\text{D}}\mid\bm{X}_{t}=\bm{x}_{t}],\;\;\;\;\forall
t<n,$

and $r_{\text{LI}}(\bm{x}_{n},\delta_{n})=0$ for all $t\in\\{0,\dots,n\\}$,
$\bm{x}_{t}\in\mathcal{X}_{t}$, $\bm{x}_{n}\in\mathcal{X}_{n}$,
$\delta_{t},\delta_{n}\in[1-p,p]$. We optimise the successes incurred under
$\Pi$ under the constraint that the total expected number of successes under
$\Pi_{\text{LI}}$ is within a percentage of the maximum under that prior,
i.e., letting
$v_{\text{LI}}=\max_{\pi}\mathbb{E}_{\text{LI}}^{\pi}\left[\sum_{t=0}^{n}r_{\text{LI}}(\bm{X}_{t},\delta_{t})\right]$
and $\xi\in(0,1)$, the resulting optimisation problem, the solution of which
is denoted as CMDP-R RA procedure, can be written as

$\displaystyle\underset{\pi}{\max}\,\mathbb{E}^{\pi}\left[\sum_{t=0}^{n}r(\bm{X}_{t},\delta_{t})\right]$
(8a) s.t.
$\displaystyle-\mathbb{E}_{\text{LI}}^{\pi}\left[\sum_{t=0}^{n}r_{\text{LI}}(\bm{X}_{t},\delta_{t})\right]\leq-\xi\cdot
v_{\text{LI}},$ (8b)

which is of the form (3) with $V_{\text{LI}}=-\xi\cdot v_{\text{LI}}$. Note
that the unconstrained version of (8) would correspond to the RA procedure in
Cheng and Berry (2007) when not considering participants outside the trial.

###  3.2 Properties of the constrained Markov decision process problem

This section provides theoretical properties of the CMDP problem (3) and
introduces a computational approach to obtain an optimal policy $\pi^{*}$ for
(3) which involves solving a linear program (LP). Section 3.3 provides a
computationally more efficient procedure using backward recursion, which finds
a feasible, but possibly suboptimal solution for (3), for which the (relative)
optimality gap of the approximation can be calculated.

We first show that problem (3) can be rewritten to a standard finite-horizon
constrained Markov decision process. First, we need the following lemma.

.

######  Lemma 4.

For all $c\in\mathcal{C}$, $t\leq n$, and states
$\bm{x}_{t}\in\mathcal{X}_{t}$

$\mathbb{P}^{\pi}(\bm{X}_{t}=\bm{x}_{t})=g_{t}^{\pi}(\bm{x}_{t})q(\bm{x}_{t}),\;\;\;\;\mathbb{P}_{c}^{\pi}(\bm{X}_{t}=\bm{x}_{t})=g_{t}^{\pi}(\bm{x}_{t})q_{c}(\bm{x}_{t}),$
(9)

where

$\displaystyle q(\bm{x}_{t})$
$\displaystyle=\int_{[0,1]^{2}}\prod_{a\in\\{\text{C},\text{D}\\}}\theta_{a}^{s_{a}(\bm{x}_{t})}(1-\theta_{a})^{n_{a}(\bm{x}_{t})-s_{a}(\bm{x}_{t})}\Pi(d\text{\boldmath$\mathbf{\theta}$}),$
(10) $\displaystyle q_{c}(\bm{x}_{t})$
$\displaystyle=\int_{[0,1]^{2}}\prod_{a\in\\{\text{C},\text{D}\\}}\theta_{a}^{s_{a}(\bm{x}_{t})}(1-\theta_{a})^{n_{a}(\bm{x}_{t})-s_{a}(\bm{x}_{t})}\Pi_{c}(d\text{\boldmath$\mathbf{\theta}$}),$

and where $g_{t}^{\pi}$ is defined recursively by

$\displaystyle g_{0}^{\pi}(\bm{x}_{0})$ $\displaystyle=1,$
$\displaystyle\indent g_{t}^{\pi}(\bm{x}_{t})$
$\displaystyle=\sum_{\begin{subarray}{c}a\in\\{\text{C},\text{D}\\}\\\
\partial\bm{x}_{a}\in\\{\partial\bm{s}_{a},\partial{\bm{f}}_{a}\\}\end{subarray}}g_{t-1}^{\pi}(\bm{x}_{t}-\partial\bm{x}_{a})\pi(\bm{x}_{t}-\partial\bm{x}_{a})^{\mathbb{I}(a=\text{C})}(1-\pi(\bm{x}_{t}-\partial\bm{x}_{a}))^{\mathbb{I}(a=\text{D})},$

for all $\bm{x}_{t}\in\mathcal{X}_{t}$, $t\in\mathbb{N}$.

###### Proof.

From (Yi, 2013, Equation (1)), we have

$\mathbb{P}_{\text{\boldmath$\mathbf{\theta}$}}^{\pi}(\bm{X}_{t}=\bm{x}_{t})=g_{t}^{\pi}(\bm{x}_{t})\prod_{a\in\\{\text{C},\text{D}\\}}\theta_{a}^{s_{a,t}}(1-\theta_{a})^{n_{a,t}-s_{a,t}}\indent\forall
t.$ (11)

The statement of the lemma follows as

$\displaystyle\mathbb{P}^{\pi}(\bm{X}_{t}=\bm{X}_{t})=\int_{[0,1]^{2}}\mathbb{P}_{\text{\boldmath$\mathbf{\theta}$}}^{\pi}(\bm{X}_{t}=\bm{x}_{t})\Pi(d\text{\boldmath$\mathbf{\theta}$}),$
$\displaystyle\mathbb{P}_{c}^{\pi}(\bm{X}_{t}=\bm{X}_{t})=\int_{[0,1]^{2}}\mathbb{P}_{\text{\boldmath$\mathbf{\theta}$}}^{\pi}(\bm{X}_{t}=\bm{x}_{t})\Pi_{c}(d\text{\boldmath$\mathbf{\theta}$}).$

∎

In Theorem 5 below, Lemma 4 is used to rewrite the constraints (3b) under the
same expectation operator as the Objective (3a) using a change of measure.

.

######  Theorem 5.

If $q({\bm{x}})=0$ implies $q_{c}(\bm{x})=0$ for all $c$ and
$\bm{x}\in\mathcal{X}$, problem (3) can be rewritten as

$\displaystyle\underset{\pi}{\max}\,\mathbb{E}^{\pi}\left[\sum_{t=0}^{n}r({\bm{X}}_{t},\delta_{t})\right]$
(12a) s.t.
$\displaystyle\mathbb{E}^{\pi}\left[\sum_{t=0}^{n}\tilde{r}_{c}({\bm{X}}_{t},\delta_{t})\right]\leq
V_{c}\indent\forall c\in\mathcal{C},$ (12b)

where

$\tilde{r}_{c}({\bm{x}}_{t},\delta_{t})=\begin{cases}r_{c}(\bm{x}_{t},\delta_{t})q_{c}(\bm{x}_{t})/q(\bm{x}_{t}),\indent&\text{if~{}$q(\bm{x}_{t})>0$},\\\
0,&\text{else.}\end{cases}$

###### Proof.

The result follows as, by Lemma 4,

$\displaystyle\mathbb{E}^{\pi}_{c}\left[\sum_{t=0}^{n}r_{c}(\bm{X}_{t},\delta_{t})\right]=\sum_{t=0}^{n}\sum_{{\bm{x}}_{t}\in\mathcal{X}_{t}}\mathbb{P}^{\pi}_{c}({\bm{X}}_{t}={\bm{x}}_{t})r_{c}({\bm{x}}_{t},\pi(\bm{x}_{t}))$
$\displaystyle\stackrel{{\scriptstyle\eqref{expression_likelihood}}}{{=}}\sum_{t=0}^{n}\sum_{\begin{subarray}{c}{\bm{x}}_{t}\in\mathcal{X}_{t}\\\
q({\bm{x}}_{t})>0\end{subarray}}g_{t}^{\pi}({\bm{x}}_{t})q_{c}({\bm{x}}_{t})r_{c}({\bm{x}}_{t},\pi({\bm{x}}_{t}))=\sum_{t=0}^{n}\sum_{\begin{subarray}{c}{\bm{x}}_{t}\in\mathcal{X}_{t}\\\
q({\bm{x}}_{t})>0\end{subarray}}g_{t}^{\pi}({\bm{x}}_{t})q({\bm{x}}_{t})\tilde{r}_{c}(\bm{x}_{t},\pi(\bm{x}_{t}))$
$\displaystyle\stackrel{{\scriptstyle\eqref{expression_likelihood}}}{{=}}\sum_{t=0}^{n}\sum_{\begin{subarray}{c}{\bm{x}}_{t}\in{\mathcal{X}_{t}}\end{subarray}}\mathbb{P}^{\pi}({\bm{X}}_{t}={\bm{x}}_{t})\tilde{r}_{c}(\bm{x}_{t},\pi(\bm{x}_{t}))=\mathbb{E}^{\pi}\left[\sum_{t=0}^{n}\tilde{r}_{c}(\bm{X}_{t},\delta_{t})\right].$

∎

The next result follows from Altman (1999, Theorem 2.1 and Theorem 3.8), by
reformulating the finite horizon constrained MDP (3) to a discounted infinite
horizon constrained MDP.

.

######  Theorem 6.

If $q({\bm{x}})=0$ implies $q_{c}(\bm{x})=0$, then for any feasible policy
which depends on the full history of states and actions, there is a feasible
randomized Markov policy $\pi^{\prime}:{\mathcal{X}}\mapsto[1-p,p]$ inducing
at least the same total expected reward in (3), i.e., under Definition 2.2 in
Altman (1999), the class of randomized Markov policies dominates the class of
history-dependent policies for the CMDP. Furthermore, there exists an optimal
Markov policy such that the number of random actions given the current state
is at most $|\mathcal{C}|$.

.

Observe that randomization, a property that is usually desired in clinical
trials, arises naturally for CMDP procedures according to Theorem 6, i.e.,
without restricting the action space. Note that this is not the case for the
MDP BRA procedures introduced in Cheng and Berry (2007) or Williamson et al.
(2017), where setting $p=1.00$ leads to a DRA (Cheng and Berry, 2007, Theorem
1), while randomized optimal actions are enforced when setting $p<1$.
Furthermore, for the model introduced in Section 2, the condition that
$q({\bm{x}})=0$ implies $q_{c}(\bm{x})=0$ is only violated when $\Pi$
concentrates on $\\{0,1\\}^{2}$ as, following (10), $q(\bm{x})$ is only zero
if either $\theta_{\text{C}}^{s_{\text{C}}(\bm{x})}$,
$(1-\theta_{\text{C}})^{n_{\text{C}}(\bm{x})-s_{\text{C}}(\bm{x})},~{}$
$\theta_{\text{D}}^{s_{\text{D}}(\bm{x})},$ or
$(1-\theta_{\text{D}})^{n_{\text{D}}(\bm{x})-s_{\text{D}}(\bm{x})}$ is zero
with probability one under the prior $\Pi$. When defining a CMDP procedure,
such measures will usually not be considered for $\Pi$.

We now turn to solving (12), and hence (3), by reformulating (12) as an LP.
Let $\mathcal{X}_{<n}=\cup_{t=0}^{n-1}\mathcal{X}_{t}$,
$d_{<n}=|\mathcal{X}_{<n}|$, $d_{n}=|\mathcal{X}_{n}|$, $d_{\leq
n}=|\mathcal{X}|$, and $d=d_{n}+2d_{<n}$ be the number of state-action pairs
in the CMDP. Let

$\bm{A}=\begin{bmatrix}\bm{I}_{d_{<n}}&\bm{I}_{d_{<n}}&\mathcal{O}_{d_{<n},d_{n}}\\\
\mathcal{O}_{d_{n},d_{<n}}&\mathcal{O}_{d_{n},d_{<n}}&\bm{I}_{d_{n}}\end{bmatrix}-\begin{bmatrix}\bm{P}_{C}&\bm{P}_{D}&\mathcal{O}_{d_{\leq
n},d_{n}}\end{bmatrix}\in\mathbb{R}^{d_{\leq n}\times d},$

where $\bm{I}_{k},\mathcal{O}_{k,\ell}$ are the identity and zero matrix in
$\mathbb{R}^{k\times k},\mathbb{R}^{k\times\ell}$ for all
$k,\ell\in\mathbb{N}$, and $\bm{P}_{C},\bm{P}_{D}\in\mathbb{R}^{d_{\leq
n}\times d_{<n}}$ are matrices such that

$\displaystyle[\bm{P}_{C}]_{i(\bm{x}),i(\bm{x}^{\prime})}$
$\displaystyle=\mathbb{P}(\bm{X}_{t+1}=\bm{x}\mid\bm{X}_{t}=\bm{x}^{\prime},\,\delta_{t}=p),\indent$
$\displaystyle\forall\bm{x}\in\mathcal{X},\bm{x}^{\prime}\in\mathcal{X}_{<n},$
$\displaystyle[\bm{P}_{D}]_{i(\bm{x}),i(\bm{x}^{\prime})}$
$\displaystyle=\mathbb{P}(\bm{X}_{t+1}=\bm{x}\mid\bm{X}_{t}=\bm{x}^{\prime},\,\delta_{t}=1-p),\indent$
$\displaystyle\forall\bm{x}\in\mathcal{X},\bm{x}^{\prime}\in\mathcal{X}_{<n},$

where $i:\mathcal{X}\mapsto\\{1,\dots,d_{\leq n}\\}$ is a storage mapping
function (Jacko, 2019) for the states such that
$i(\mathcal{X}_{<n})=\\{1,\dots,d_{<n}\\}$. Let $\bm{r}\in\mathbb{R}^{d}$ be a
vector defined as

$r_{j}=\sum_{\bm{x}\in\mathcal{X}}r(\bm{x},p)\mathbb{I}(j=i(\bm{x}))+r(\bm{x},1-p)\mathbb{I}(j=i(\bm{x})+d_{<n}).$

Let $\tilde{\bm{r}}_{c}$ be a similarly defined vector for all
$c\in\mathcal{C}$, and $\bm{b}$ be a vector denoting the initial distribution,
defined such that $b_{j}=\mathbb{I}(j=i(\bm{x}_{0}))$ for all $j.$ The next
theorem states that the optimal policy for (12), hence (3), can be found by
solving a linear program.

.

######  Theorem 7.

If (3) is feasible and $q({\bm{x}})=0$ implies $q_{c}(\bm{x})=0$, an optimal
policy $\pi^{*}$ for (3) is given by

$\pi^{*}(\bm{x})=\begin{cases}\frac{p\mu_{i(\bm{x})}+(1-p)\mu_{i(\bm{x})+d_{<n}}}{(\mu_{i(\bm{x})}+\mu_{i(\bm{x})+d_{<n}})},\indent&\text{if~{}$\mu_{i(\bm{x})}+\mu_{i(\bm{x})+d_{<n}}>0$,}\\\
1/2,&\text{else,}\end{cases}$

where the vector $\mathbf{\mu}$ is the solution to

$\displaystyle\max_{\text{\boldmath$\mathbf{\mu}$}\in\mathbb{R}^{d}}\text{\boldmath$\mathbf{\mu}$}^{\top}\bm{r}$
(13a) s.t. $\displaystyle\bm{A}\text{\boldmath$\mathbf{\mu}$}=\bm{b},$ (13b)
$\displaystyle\text{\boldmath$\mathbf{\mu}$}^{\top}\tilde{\bm{r}}_{c}\leq
V_{c}\indent\forall c\in\mathcal{C},$ (13c)
$\displaystyle\text{\boldmath$\mathbf{\mu}$}\geq\bm{0}.$ (13d)

###### Proof.

The result can be shown by writing the dual linear program formulation of the
Markov decision process (12) without constraints (Puterman, 2014), consisting
of (13a), (13b) and (13d), where $\delta_{t}\in[1-p,p]$ can be viewed as a
Randomization over actions $a_{t}\in\\{1-p,p\\}$. The statement of the theorem
follows by adding the constraints in (12b) as (13c), which follow as
$\mathbf{\mu}$ is the probability vector for the state-action pairs in the
CMDP, a fact that is also used in determining the optimal policy $\pi^{*}$. ∎

###  3.3 Proposed solution method

Solving (13) can be computationally heavy for large values of $d$. The next
theorem provides a way to find a solution to (3) using backward recursion,
which is computationally more tractable. Conditions are given under which the
obtained solution is optimal.

.

######  Theorem 8.

Let $V$ be the value of (3). If (3) is feasible and $q({\bm{x}})=0$ implies
$q_{c}(\bm{x})=0$, we have

$V=\min_{\text{\boldmath$\mathbf{\lambda}$}\in\mathbb{R}^{C}_{+}}L(\text{\boldmath$\mathbf{\lambda}$}),\;\;L(\text{\boldmath$\mathbf{\lambda}$})=\max_{\pi}\mathbb{E}^{\pi}\left[\sum_{t=0}^{n}{r}(\bm{X}_{t},\delta_{t})+\sum_{c}\lambda_{c}\left(V_{c}/n-\tilde{r}_{c}(\bm{X}_{t},\delta_{t})\right)\right].$
(14)

If (14) yields an optimal point
$(\hat{\text{\boldmath$\mathbf{\lambda}$}},\hat{\pi})$ such that

$\sum_{c}\hat{\lambda}_{c}\mathbb{E}^{\hat{\pi}}\left[\sum_{t=0}^{n}\left(V_{c}/n-\tilde{r}_{c}(\bm{X}_{t},\delta_{t})\right)\right]=0,$
(15)

i.e., $(\hat{\text{\boldmath$\mathbf{\lambda}$}},\hat{\pi})$ satisfies the
Karush-Kuhn Tucker (KKT) conditions (Kuhn and Tucker, 2014), then $\hat{\pi}$
is an optimal solution to (3).

###### Proof.

Using Lagrange multipliers, letting $\tilde{\bm{R}},\;\tilde{\bm{V}}$ be the
horizontal concatenation of the vectors $\tilde{\bm{r}}_{c}$ and constants
$V_{c}$, we reformulate (13) to

$\displaystyle\max_{\text{\boldmath$\mathbf{\mu}$}\in\mathbb{R}^{d}}\min_{\text{\boldmath$\mathbf{\lambda}$}\geq\bm{0}}\,\text{\boldmath$\mathbf{\mu}$}^{\top}\bm{r}+\text{\boldmath$\mathbf{\lambda}$}^{\top}(\tilde{\bm{V}}-\tilde{\bm{R}}\text{\boldmath$\mathbf{\mu}$})$
(16a) s.t. $\displaystyle\bm{A}\text{\boldmath$\mathbf{\mu}$}=\bm{b},$ (16b)
$\displaystyle\text{\boldmath$\mathbf{\mu}$}\geq\bm{0}.$ (16c)

By (Altman, 1999, Lemma 9.2), we have that (16) is equivalent to

$\displaystyle\min_{\text{\boldmath$\mathbf{\lambda}$}\geq\bm{0}}\max_{\text{\boldmath$\mathbf{\mu}$}\in\mathbb{R}^{d}}\,\text{\boldmath$\mathbf{\mu}$}^{\top}\bm{r}+\text{\boldmath$\mathbf{\lambda}$}^{\top}(\tilde{\bm{V}}-\tilde{\bm{R}}\text{\boldmath$\mathbf{\mu}$})$
(17a) s.t. $\displaystyle\bm{A}\text{\boldmath$\mathbf{\mu}$}=\bm{b},$ (17b)
$\displaystyle\text{\boldmath$\mathbf{\mu}$}\geq\bm{0}.$ (17c)

Choosing $\mathbf{\lambda}$ fixed in (17) results in the LP that returns
$L(\text{\boldmath$\mathbf{\lambda}$})$, from which (14) follows. The
conditions for an optimal solution follow from Kuhn and Tucker (2014). ∎

The function $L$ is a convex function in $\mathbf{\lambda}$, as it is the
maximum over a set of affine functions in $\mathbf{\lambda}$ (one for each
$\pi$). Backward recursion (Puterman, 2014) can be performed to find the value
$L(\text{\boldmath$\mathbf{\lambda}$})$ and a maximizer
$\pi^{*}_{\text{\boldmath$\mathbf{\lambda}$}}$ of (14) (letting
$\pi^{*}_{\text{\boldmath$\mathbf{\lambda}$}}(\bm{x}_{t})=1/2$ in case of
ties). Hence, without loss of generality, the set of randomized Markov
policies optimized over in (14) can be restricted to the set $\mathcal{P}$ of
policies $\pi:\mathcal{X}\mapsto\\{1-p,1/2,p\\}$ when determining $V$. The
optimal policy $\pi^{*}_{\text{\boldmath$\mathbf{\lambda}$}}$ in this set can
be used to determine a subgradient $\nabla
L(\text{\boldmath$\mathbf{\lambda}$})~{}$ of $L$ such that $\nabla
L(\text{\boldmath$\mathbf{\lambda}$})_{c}=V_{c}-\sum_{t=0}^{n}\tilde{r}_{c}(\bm{X}_{t},\pi_{\text{\boldmath$\mathbf{\lambda}$}}^{*}(\bm{X}_{t}))$
for all $c\in\mathcal{C}$. We propose Algorithm 1 to find $V$ up to a given
numerical precision, which corresponds to a cutting plane method (Kelley,
1960). In case of no feasible solution, the algorithm returns $f^{*}=-\infty.$
Otherwise, as $L$ is the maximum over a finite set of affine functions in
$\mathbf{\lambda}$, one for each $\pi\in\mathcal{P}$, the algorithm finds the
minimizer after a finite number of iterations. However, due to the (possibly)
large amount of affine functions included in the objective, this could take a
long time, while a small error $\epsilon_{\text{tol}}$ may usually be
sufficient. If the policy found under backward recursion is infeasible, the
elements of $\mathbf{\lambda}$ corresponding to the infeasible constraints are
multiplied by a factor $(1+\varphi)$ until the policy becomes feasible.
Letting $\hat{\pi}$ be the resulting policy, the relative optimality gap can
be found as
$(V-\mathbb{E}^{\hat{\pi}}[\sum_{t=0}^{n}r(\bm{X}_{t},\delta_{t})])/V$. If
(15) is satisfied, it can be verified that $\hat{\pi}$ is the unique optimiser
of (3) by verifying that $\hat{\pi}(\bm{x})\in\\{1-p,p\\}$ for all
$\bm{x}\in\mathcal{X}$.

Algorithm 1 Cutting plane algorithm for calculating V

1:Inputs:

2: $L,\nabla L,|\mathcal{C}|,\epsilon_{\text{tol}}$;

3:Initialize:

4: Set $\bm{A}=\bm{0}^{\top}_{|\mathcal{C}|+1}$, $b=0$,
$\epsilon=\infty,\;f^{*}=\infty$;

5:while $\epsilon>\epsilon_{\text{tol}}$ and $f^{*}\neq-\infty$ do

6: Solve  $\displaystyle\min_{\bm{x}\in\mathbb{R}^{|\mathcal{C}|+1}}x_{1}$
s.t. $\displaystyle\bm{A}\bm{x}\leq\bm{b}$ $\displaystyle\bm{x}\geq\bm{0};$

7: Set
$\underaccent{\bar}{f}=x_{1},\;\text{\boldmath$\mathbf{\lambda}$}^{*}=[x_{2},\dots,x_{|\mathcal{C}|+1}],\;f^{*}=\min(f^{*},L(\text{\boldmath$\mathbf{\lambda}$}^{*})),\;\epsilon=f^{*}-\underaccent{\bar}{f},\newline
\;\bm{g}~{}=~{}\nabla L(\text{\boldmath$\mathbf{\lambda}$}^{*});$

8: Redefine $\bm{A}=[\bm{A};[-1,\bm{g}^{\top}]]$,
$\bm{b}=[\bm{b};-L(\text{\boldmath$\mathbf{\lambda}$}^{*})+\bm{g}^{\top}\text{\boldmath$\mathbf{\lambda}$}^{*}]$
($[\bm{v};\bm{w}]$ denotes concatenation);

9: if $L(\text{\boldmath$\mathbf{\lambda}$}^{*})\leq 0$ then

10: $f^{*}=-\infty$, $\text{\boldmath$\mathbf{\lambda}$}^{*}=\emptyset$;

11: end if

12:end while

13:Outputs: $f^{*},\text{\boldmath$\mathbf{\lambda}$}^{*}$;

##  4 Applications

In this section, we evaluate the performance of CMDP procedures in three
applications. First, we evaluate the performance of CMDP-T from Example 1.
Second, we evaluate the performance of CMDP-E from Example 2. Third, we
evaluate the performance of CMDP-R from Example 3. All CMDP procedures are
evaluated according to the OCs introduced in Section 2.2, which are calculated
using forward recursion, instead of simulation, to determine the distribution
of $\bm{X}_{n}$.

We compare the performance of the CMDP procedures with three other Markov RA
procedures known from literature (Williamson et al., 2017), namely:

  * •

Equal randomization (ER):
$\pi(\bm{x})=1/2$ for all $\bm{x}\in\mathcal{X}$.

  * •

Dynamic Programming (DP):
This is the RA procedure $\pi$ found from the unconstrained version of (3)
when taking $p=1.00$.

  * •

Constrained Randomized Dynamic Programming (CRDP):
This is the RA procedure proposed in (Williamson et al., 2017) that follows
from the unconstrained version of (3) when taking $p=0.9$ and adding a penalty
equal to $-n$ to the objective whenever $\min_{a}n_{a}(\bm{x}_{n})<0.15n$.

The code for calculation of the CMDP policies and the OCs uses an efficient
implementation of backward and forward recursion based on Jacko (2019), i.e.,
uses the conservation law for the states, a storage mapping function, and
overwrites elements of the value function in backward recursion and
probability distribution in forward recursion that are not used further in the
algorithm. All experiments are performed in Julia 1.9.0 on a laptop with
Intel® CoreTM i7-9750H CPU with 2.60GHz clock speed and 16 GB RAM. The
absolute tolerance $\epsilon_{\text{tol}}$ for Algorithm 1 is set to $10^{-9}$
for each evaluation. In order to obtain a feasible policy in each considered
scenario, we set $\varphi=0.01$. For the evaluation, we let
$\theta_{\text{D}}\in\\{0.00,0.01\dots,1.00\\}$ while fixing
$\theta_{\text{C}}=0.5.$ In the numerical evaluation, we set the significance
level $\alpha~{}$ to $0.1$, in agreement with the comparison in Williamson et
al. (2017). For all CMDP procedures, we set $p=0.95$, while in the appendices
results are also shown for $p=1.00$. From now on, RA procedures based on a
CMDP problem and solved using Algorithm 1 will also be denoted CMDP RA
procedures.

###  4.1 Application 1: Control of power and type I error

First, we evaluate the CMDP-T procedure (Example 1) for $n=75$. We set the
prior $\Pi_{1}$ for the power constraint equal to $\Pi$. After evaluating
several choices based on the OCs, we set $\alpha^{*}=0.05$, $\beta=0.4$. The
prior parameters $\tilde{s}_{a,0},\tilde{f}_{a,0}$ for $\Pi_{0}$ are set to 1,
as are the parameters for the prior $\Pi$, hence both $\Pi_{0}$ and $\Pi$ are
uniform priors on their respective supports. Computation of the CMDP-T policy
took 7.44 seconds and resulted in a relative optimality gap of $1.34\cdot
10^{-4}$ and total expected reward 46.9 under the prior $\Pi$ in comparison to
$0.165$ seconds for CRDP with total expected reward $45.3$. Figure 1 shows the
patient benefit, RR, bias, and MSE for the policies ER, DP, CRDP and CMPD-T
for different values of $\theta_{\text{D}}$ and $\theta_{\text{C}}=0.5$. The
same results for policies ER, DP, and CRDP can also be found in Williamson et
al. (2017), where the sign of the bias is reversed because we chose to vary
the first success rate used in the treatment effect instead of the second one.
Figure 1 shows that the CMDP-T procedure performs very well in terms of
patient benefit, having performance in-between CRDP and DP for all evaluated
parameter values. Furthermore, the CMDP-T procedure has a similar performance
in terms of power, while having a worse performance in terms of bias and MSE
when compared to CRDP. It is hence seen for this scenario that a trade-off was
made between patient benefit and OCs bias and MSE for CMDP-T, keeping the
rejection rate roughly the same. The bias and MSE of CMDP-T are however not as
extreme as for DP, and closer to the bias and MSE of CRDP.

Second, we evaluate the results for CMDP-T for $n=200$. After evaluating
several choices based on the OCs, we set $\alpha^{*}=0.07$, and $\beta=0.23$.
Again, the prior parameters $\tilde{s}_{a,0},\tilde{f}_{a,0}$ are set to 1, as
are the parameters for $\Pi$, hence both priors correspond to uniform
distributions on their respective supports. Computation of the CMDP-T policy
took 296.3 seconds, and resulted in a relative optimality gap of $5.09\cdot
10^{-5}$ and total expected reward 128 under the prior $\Pi$ in comparison to
$1.30$ seconds for CRDP with total expected reward $122$. Figure 2 shows the
results. These results were not shown in Williamson et al. (2017), where the
maximum trial size considered was $n=100$ participants as the focus was on
rare disease trials with a small amount of participants. The figure shows that
in terms of patient benefit, the CRDP procedure ends up around a value of
$0.85$ for $\theta_{\text{D}}$ close to zero and one, which is (roughly) the
maximal value this RA procedure can attain, due to the penalty incurred when
either $N_{\text{C},n}/n$ or $N_{\text{D},n}/n$ are lower than $0.15$. The DP
and CMDP-T procedures end up at higher values, around $0.95$ (maximum for
CMDP-T) and $1.00$ (maximum for DP) for $\theta_{\text{D}}$ close to zero and
one. The power plots show that the power quickly grows to 1.00 in
$|\theta_{\text{D}}-\theta_{\text{C}}|$ for policies ER, CRDP and CMDP-T where
a power of $80\%$ is roughly attained when
$|\theta_{\text{D}}-\theta_{\text{C}}|\geq 0.25$, while for the DP policy, the
power is slightly higher than in Figure 1, with highly irregular behaviour
around high values of $\theta_{\text{D}}.$

In conclusion, for $n=200$ the CMDP-T procedure outperforms CRDP in terms of
patient benefit and power, while CMDP-T significantly outperforms ER in terms
of patient benefit and has slightly lower power. The bias and MSE for DP have
not changed much in comparison with $n=75,$ which is possibly due to the
algorithm allocating all trial participants to one treatment after a certain
time point. The bias and MSE for ER, CRDP and CMDP-T have decreased
significantly in comparison with Figure 1, where the same ordering as for
$n=75$ is seen in the OCs. Figure 2 shows again that a trade-off was made
between patient benefit and OCs bias and MSE for the CMDP-T procedure.

Figure 1: Patient benefit (subfigure a), rejection rate (subfigure b), bias
(subfigure c), and mean squared error (subfigure d) vs. $\theta_{\text{D}}$
for $\theta_{\text{C}}=0.5$, $n=75$, $\alpha^{*}=0.05$, $\beta=0.4$, and RA
procedures ER (solid), DP (dashed), CRDP (dotted) and CMDP-T (dash-dotted)
Figure 2: Patient benefit (subfigure a), rejection rate (subfigure b), bias
(subfigure c), and mean squared error (subfigure d) vs. $\theta_{\text{D}}$
for $\theta_{\text{C}}=0.5$, $n=200$, $\alpha^{*}=0.07$, $\beta=0.23$, and RA
procedures ER (solid), DP (dashed), CRDP (dotted) and CMDP-T (dash-dotted)

Figures 7 and 8 in Appendix A show the OCs according to the comparison
described above when $p=1.00$ for CMDP-T instead of $0.95$, i.e., when a DRA
CMDP-T procedure is used. Figures 7 and 8 show that for the DRA CMDP-T
procedure, the patient benefit increases and the type I error and power remain
under control, while the bias and MSE are worse than for CMDP-T with $p=0.95$
and take on values similar to those for DP when
$\theta_{\text{D}}\gg\theta_{\text{C}}.$ Hence, the randomization incorporated
in the CMDP-T procedure with $p=0.95$ mitigates a large part of the bias and
MSE. Figures 9 \- 12 in Appendix A show the OCs according to the comparison
above when $p\in\\{0.95,1.00\\}$, $n\in\\{75,200\\}$, for several values of
$\theta_{\text{D}}$, and $\theta_{\text{C}}\in\\{0.25,0.75\\}$. The behaviour
of the power remains similar when varying $p$ and $\theta_{\text{C}}$. When
$\theta_{\text{C}}=0.25$, a uniform outperformance of CRDP by CMDP-T is no
longer seen in terms of power as the rejection rate for CRDP is higher for
$\theta_{\text{D}}\approx 0.25$. Note, however, that the rejection rate for
CRDP is only significantly higher than that for CMDP-T for parameter values
close to the null, where the rejection rate is low anyway.

###  4.2 Application 2: Control of estimation error

Two CMDP-E policies, calculated using (7) in Example 2 and denoted CMDP-E1 and
CMDP-E2, are evaluated. For policy CMDP-E1, the values of
$\alpha^{*},\beta,V_{\sigma}$ are chosen such that the MSE is similar to that
of ER, while there is a gain in patient benefit. For CMDP-E1 we set
$\mathcal{S}=\\{[0,1]^{2}\\}$, i.e., $\mathcal{S}$ only contains the unit
square, and we calculate the value of (6) realised under the ER RA procedure
and set $V_{\sigma}$ to $\xi_{1}$ times this value for $\xi_{1}\in[1,\infty)$,
we set $\alpha^{*}=\beta=1$, i.e., the type I and II error constraints are
automatically satisfied, hence we focus solely on MSE for the policy CMDP-E1.
For CMDP-E2 we take
$[\theta_{a}^{\ell},\theta_{a}^{u})\in\\{[0,0.25),[0.25,0.5),[0.5,0.75),[0.75,0.9),[0.9,1.0)\\}$
for all $a\in\\{\text{C},\text{D}\\}$ in order to construct $\mathcal{S}$,
dividing the unit square in blocks with a surface $0.25^{2}$, as well as some
additional blocks for $\theta_{\text{C}}\geq 0.9$ or $\theta_{\text{D}}\geq
0.9$ where the largest imbalances in treatment group sizes occur (for RA
procedures inducing high patient benefit). We set $V_{\sigma}$ equal to
$\xi_{2}$ times the value of (6) realised under CRDP for
$\xi_{2}\in[1,\infty).$

First, we evaluate the results for CMDP-E for $n=75$. After evaluating several
choices based on the OCs, we set $\xi_{1}=1.05$. The prior parameters
$\tilde{s}_{a,0},\tilde{f}_{a,0}$ for $\Pi_{0}$ are set to 1, as are the
parameters for the prior $\Pi$. Computation of the CMDP-E1 policy took 28.2
seconds, and resulted in a relative optimality gap of $8.18\cdot 10^{-7}$ and
total expected reward 41.3 under the prior $\Pi$, while computation of the
CMDP-E2 policy, where we set $\xi_{2}=1.00,$ $\alpha^{*}=0.05,$ $\beta=0.4$,
took 708 seconds, and resulted in a relative optimality gap of $4.59\cdot
10^{-4}$ and total expected reward 45.4 under the prior $\Pi$ in comparison to
$0.165$ seconds for CRDP with total expected reward $45.3$. Figure 3 shows the
results. Figure 3 shows that the CMDP-E1 procedure induces higher patient
benefit than the ER RA procedure, while having similar power and MSE. The
CMDP-E1 procedure induces higher patient benefit when $\theta_{\text{D}}\leq
0.5$ in comparison to when $\theta_{\text{D}}\geq 0.5.$ This behaviour could
be explained by the fact that allocating more patients to the arm with the
highest variance (in the spirit of Neyman allocation (Rosenberger and Hu,
2004)), which is always the control arm in this evaluation, maximizes power
and hence induces good OCs. Figure 3 shows that the CMDP-E2 procedure induces
similar patient benefit, power and MSE to CRDP, while showing slightly higher
bias. This result is curious, as the two policies are found using two
completely different procedures, but could be explained by the fact that the
parameters $\ell$ and $p$ for CRDP are tuned specifically to balance patient
benefit, power, bias, and MSE, while the constraints for CMDP-E2 are based on
the attained MSE for CRDP.

Second, we evaluate the results for $n=200$. For CMDP-E1, we set
$\xi_{1}=1.1$. The prior parameters $\tilde{s}_{a,0},\tilde{f}_{a,0}$ for
$\Pi_{0}$ are set to 1, as are the parameters for the prior $\Pi$. Computation
of the CMDP-E1 policy took 137 seconds, and resulted in a relative optimality
gap of $3.98\cdot 10^{-7}$ and total expected reward 113 under the prior
$\Pi$, while computation of the CMDP-E2 policy, where we set
$\xi_{2}=1.05,\;$$\alpha^{*}=0.07$, $\beta=0.753$, took 190$\cdot 10^{2}$
seconds, and resulted in a relative optimality gap of $2.75\cdot 10^{-4}$ and
total expected reward 123 under the prior $\Pi$ in comparison to $1.30$
seconds for CRDP with total expected reward $122$. Figure 4 shows the results.
Figure 4 shows that the CMDP-E1 procedure again induces higher patient benefit
than the ER RA procedure while having similar power and MSE, while inducing
higher patient benefit when $\theta_{\text{D}}\leq 0.5$ in comparison to when
$\theta_{\text{D}}\geq 0.5.$ Figure 4 shows that the CMDP-E2 procedure induces
similar power and MSE to CRDP, while showing slightly higher bias and slightly
higher patient benefit.

In conclusion, for $n=200$ the CMDP-E2 procedure outperforms CRDP in terms of
patient benefit, while showing similar MSE and power and slightly higher bias,
which can possibly be mitigated by bias reduction techniques (Bowden and
Trippa, 2017). Such procedures can however also lead to higher variance of the
treatment effect estimator. The outperformance of CRDP by the CMDP-E1
procedure in terms of patient benefit is not as high as in the previous
application, indicating that CRDP is a well-performing policy for a vast range
of parameter values and OCs, while a CMDP procedure such as CMDP-T may
outperform CRDP for specific settings/OCs.

Figures 13 and 14 in Appendix A show the OCs according to the comparison
described above when $p=1.00$ instead of $0.95$ for CMDP-E1 and CMDP-E2, i.e.,
when a DRA version of CMDP-E1 and CMDP-E2 is used. Figures 13 and 14 show very
similar measures to the CMDP-E1 and CMDP-E2 RA procedures for $p=0.95$, in
contrast to the comparison with $p=0.95$ vs. $p=1.00$ for CMDP-T, where large
differences in bias and MSE were seen. The RAR procedures CMDP-E1 and CMDP-E2
($p=0.95$) are preferred over the DRA procedures CMDP-E1 and CMDP-E2
($p=1.00$) as the patient benefit is similar, while using an RAR policy brings
advantages, e.g., in terms of selection bias mitigation. Figures 15 \- 18 in
Appendix A show the OCs according to the comparison above when
$p\in\\{0.95,1.00\\}$, $n\in\\{75,200\\}$, and
$\theta_{\text{C}}\in\\{0.25,0.75\\}$. Figures 15 \- 18 show that the
difference in the power and patient benefit between CRDP and CMDP-E2 remains
small when varying $p$ and $\theta_{\text{C}}$, while larger differences are
seen between the bias and MSE of CRDP and CMDP-E2 for
$\theta_{C}\in\\{0.25,0.75\\}$ and $n=200$.

Figure 3: Patient benefit (subfigure a), rejection rate (subfigure b), bias
(subfigure c), and mean squared error (subfigure d) vs. $\theta_{\text{D}}$
for $\theta_{\text{C}}=0.5$, $n=75$, and RA procedures ER (solid), CMDP-E2
(dashed), CRDP (dotted) and CMDP-E1 (dash-dotted) Figure 4: Patient benefit
(subfigure a), rejection rate (subfigure b), bias (subfigure c), and mean
squared error (subfigure d) vs. $\theta_{\text{D}}$ for
$\theta_{\text{C}}=0.5$, $n=200$, and RA procedures ER (solid), CMDP-E2
(dashed), CRDP (dotted) and CMDP-E1 (dash-dotted)

###  4.3 Application 3: Robustness to prior misspecification

First, we evaluate CMDP-R (Example 3) for $n=200$ and
$\tilde{s}_{\text{C},0}=3$, $\tilde{f}_{\text{C},0}=7$,
$\tilde{s}_{\text{D},0}=6$, $\tilde{f}_{\text{D},0}=4$, hence under the
informative prior the control and development treatment have a prior
probability of success of 0.3 and 0.6, with a prior sample size of 10. The
prior $\Pi_{\text{LI}}$ corresponds to an independent uniform prior for both
success probabilities. We do not consider CRDP in this section as this RA
procedure was not constructed with sensitivity to prior misspecification in
mind, hence, we will also not present results for $n=75$ as in Williamson et
al. (2017). In fact, to the best of our knowledge, we are the first to
introduce a BRA procedure focusing on robustness against prior
misspecification.

To investigate different settings of robustness to the prior for the CMDP-R,
we compare policies found under (8) for $\xi=0.00,\;0.990,\;0.999$ and $1.00$
according to the OCs, which took 4.84, 123, 1260, and 3.34 seconds to compute,
with optimality gaps 0.00, $4.87\cdot 10^{-3}$, $4.10\cdot 10^{-3}$, and 0.00,
and optimal values 118, 117, 117, and 116. Figure 5 shows the patient benefit,
RR, bias, and MSE for the four different CMDP-R policies, where we now set
$\theta_{\text{C}}=0.3$ and again vary $\theta_{\text{D}}$ to investigate
prior robustness. Figure 5 shows that the risk under prior misspecification
decreases with $\xi$, as the patient benefit increases for
$\theta_{\text{D}}<0.3$ in $\xi.$ The ordering of patient benefit flips after
$\theta_{\text{D}}$ passes the value $0.3$, i.e., depending on whether the
prior rightly specifies $\theta_{\text{D}}\geq\theta_{\text{C}}$, and the
patient benefit curve becomes symmetric when $\xi=1.0.$ As a higher amount of
allocations to the optimal treatment means a larger imbalance in treatment
group sizes, the reverse ordering in quality, flipped after
$\theta_{\text{D}}$ passes $0.3$, is seen for the power and MSE, where this
behaviour is seen to occur more drastically for the power, where the ordering
flips exactly at $\theta_{\text{D}}=0.3.$ Comparing the bias curves, the
curves become more symmetrical around $0.3$ when $\xi$ increases. This is
because the robustness constraint enforces that the policy depends less on the
informative prior when $\xi$ increases, hence the observations tend to guide
the policy more when $\xi~{}$ is close to $1,$ inducing a change in sign
around $\theta_{\text{D}}=0.3$ and a symmetric behaviour for the policy in
both arms, leading to a more symmetric looking bias.

Second, we evaluate the results for $n=200$ and $\tilde{s}_{\text{C},0}=30$,
$\tilde{f}_{\text{C},0}=70$, $\tilde{s}_{\text{D},0}=60$,
$\tilde{f}_{\text{D},0}=40$, hence under the informative prior the control and
developmental treatment have a prior probability of success of 0.3 and 0.6,
with a prior sample size of 100. The less informative prior was again chosen
as a uniform prior for both success probabilities. We compare policies found
for (8) for $\xi=0.00,\;0.900,\;0.990$ and $1.00$, which took 3.00, 165, 171,
and 2.95 seconds to compute, with optimality gaps 0.00, $5.60\cdot 10^{-3}$,
$3.09\cdot 10^{-3}$, and 0.00, and optimal values 117, 117, 117, and 115.
Figure 6 shows the patient benefit, RR, bias, and MSE for the four different
CMDP-R policies, where we again set $\theta_{\text{C}}=0.3$ to investigate
prior robustness. Figure 6 again shows that the risk under prior
misspecification decreases with $\xi$, where the ordering of patient benefit
flips after $\theta_{\text{D}}$ passes the value $0.3$ and patient benefit
curve becomes symmetric when $\xi=1.0$. The risk of prior misspecification is
higher in this case, as the patient benefit stays around $5\%$ even for
$\theta_{\text{D}}=0.2,$ for $\xi=0.0$. The flip in the ordering for
$\theta_{\text{D}}<0.3$ vs. $\theta_{\text{D}}\geq 0.3$ is seen for the curves
corresponding to $\xi>0$, while the curve for $\xi=0.0$ shows a different
behaviour. Due to the large amount of prior certainty, the policy exploits
knowledge on $\theta_{\text{D}}>\theta_{\text{C}}$ for a large part of the
parameter space, effectively acting as a fixed randomization policy with a
probability of $95\%~{}$ of allocating to the control arm. This is why there
is low bias and roughly constant MSE for $\theta_{\text{D}}>0.3$. As the
treatment groups are extremely unbalanced, the power is low and the MSE is
high for such values of $\theta_{\text{D}}.$ As exploration is seen to kick in
for CMDP-R with $\xi=0.0$ when $\theta_{\text{D}}\ll 0.3$, the treatment
groups become more balanced when $\theta_{\text{D}}$ is close to zero,
resulting in an increase in power.

Note that, rewriting (14), we have

$L(\lambda)=\max_{\pi}\mathbb{E}^{\pi}\left[\sum_{t=0}^{n}{r}(\bm{X}_{t},\delta_{t})\right]+\lambda\mathbb{E}^{\pi}_{\text{LI}}\left[\sum_{t=0}^{n}{r}_{\text{LI}}(\bm{X}_{t},\delta_{t})\right]-\lambda
v_{\text{LI}}\xi.~{}$

The two expectations above can be rewritten to the expected sum of successes
under a mixture prior of $\Pi$ and $\Pi_{\text{LI}}$, so an alternative method
would be to formulate a mixture prior and optimise the sum of rewards under
this mixture prior. The difference with the CMDP-R approach is that (8) gives
a robustness guarantee of the resulting policy, something that might be of
value when designing a clinical trial. Furthermore, one might easily increase
the number of constraints in (8), while adding more priors to the mixture
would quickly make the procedure intractable as all weights need to be
elicited beforehand.

Figure 5: Patient benefit (subfigure a), rejection rate (subfigure b), bias
(subfigure c), and mean squared error (subfigure d) vs. $\theta_{\text{D}}$
for $\theta_{\text{C}}=0.3$, $n=200$, $\tilde{s}_{1,0}=3$,
$\tilde{f}_{1,0}=7$, $\tilde{s}_{2,0}=6$, $\tilde{f}_{2,0}=4$ and the CMDP-R
procedure (8) with $\xi=0.00,0.990,0.999,1.00$ denoted by the solid, dashed,
dotted and dash-dotted lines respectively.  Figure 6: Patient benefit
(subfigure a), rejection rate (subfigure b), bias (subfigure c), and mean
squared error (subfigure d) vs. $\theta_{\text{D}}$ for
$\theta_{\text{C}}=0.3$, $n=200$, $\tilde{s}_{1,0}=30$, $\tilde{f}_{1,0}=70$,
$\tilde{s}_{2,0}=60$, $\tilde{f}_{2,0}=40$ and the CMDP-R procedure (8) with
$\xi=0.00,0.990,0.999,1.00$ denoted by the solid, dashed, dotted and dash-
dotted lines respectively.

##  5 Discussion

The current paper introduced the class of constrained Markov decision process
response-adaptive procedures. The constraints for the Markov decision process
facilitate the construction of response-adaptive procedures with good OCs
while prioritizing one operating characteristic in the trial (e.g., patient
benefit) by putting it in the objective.

Three applications were presented in the current paper. In the first
application, type I error and power constraints were formulated, and the
resulting randomized response-adaptive procedure, CMDP-T, was shown to yield
comparable power to the constrained randomized dynamic programming RA
procedure introduced in Williamson et al. (2017) while inducing higher patient
benefit. Furthermore, it was shown that the deterministic version of CMDP-T
also has the desired power properties, while inducing higher patient benefit
than the randomized CMDP-T procedure and constrained randomized dynamic
programming. A deterministic procedure might however have certain non-
desirable properties, e.g., allocation biases, but it is nonetheless curious
that a deterministic policy with good power properties could be constructed. A
topic of future research for CMDP-T could be to evaluate the performance of
CMDP-T for other tests such as tests for superiority or non-inferiority.

In the second application, we constructed two constrained Markov decision
process procedures by formulating additional constraints on the MSE. The
resulting procedures, CMDP-E1 and CMDP-E2, showed low MSE in addition to a
similar power and type I error behaviour as seen for CMDP-T. Here, the priors
under the constraints were designed in such a way that the constraint
reflected the average behaviour of the policy on a specific part of the
parameter space. The CMDP-E1 procedure showed higher patient benefit than the
non-adaptive equal allocation procedure while having similar MSE and power.
The CMDP-E2 procedure showed OCs similar to constrained randomized dynamic
programming, where only for a trial horizon of 200 participants, slight
outperformance in patient benefit by the constrained approach was seen,
indicating the good performance and general applicability of the constrained
randomized dynamic programming procedure.

In the third application, a constrained Markov decision process procedure,
CMDP-R, was developed for prior robustness, where the objective is to maximize
patient benefit under an informative prior, with the restriction that the
patient benefit should also be sufficiently high under an uninformative prior
on the success probabilities. It was shown that the constraint parameter $\xi$
could be tuned such that the desired robustness to misspecification of prior
information could be attained. To the best of our knowledge, ours is the first
approach to incorporate robustness to prior misspecification directly into the
design of the response-adaptive procedure.

Other applications might also be of interest to explore. For instance, in line
with recent literature on higher moments and tail bounds of the regret
distribution for multi-armed bandits (Fan and Glynn, 2023), it could be
possible to define constraints on moments and tails of the number of
suboptimal decisions made by the Markov decision process procedure, ensuring
that the performance of the response-adaptive procedure in terms of patient
benefit is consistent over different sampled data sets. The constraints could
also be used to constrain the allocations to both treatment groups to a
minimum amount with high probability, as is the purpose of the penalty term in
constrained randomized dynamic programming, and differences in OCs between
such a constrained approach and constrained randomized dynamic programming
could be investigated.

In application 2 (where the focus was on the mean squared error), it was seen
that in order to obtain a good performance in all OCs, the parameters in the
constraints had to be chosen in a specific manner. This is mainly because the
constraints are Bayesian, hence they only guarantee a certain behaviour when
averaged over a prior. In order to obtain good frequentist OCs in a more
straightforward or automatic manner, it might be of interest to explore
approaches more in line with robust optimisation, or to add a second layer of
optimisation, guaranteeing good OCs for all possible parameter configurations.

In this paper, we proposed an algorithm based on Lagrange multipliers and
backward recursion to calculate CMDP RA procedures. In our applications, we
found that this algorithm yields a policy with a small relative optimality gap
in all cases. However, there is no guarantee that this will always be the
case. Namely, in settings where the amount of constraints is larger in
comparison to the number of state-action pairs, this method might result in a
highly suboptimal policy. Future research could focus on alternative solution
methods that approximate the optimal solution to the CMPD problem in different
ways, e.g., by using the linear programming formulation more directly.

For the sake of simplicity and tractability, the current paper considers a
trial with two arms and binary outcomes. It would be interesting to extend the
approach to more general settings, such as multiple trial stages, outcomes,
and arms. As the computational effort for calculating constrained Markov
decision process policies would quickly increase in such settings, it would be
interesting to investigate whether approaches such as approximate dynamic
programming or reinforcement learning could be (successfully) applied. It
would furthermore be interesting to consider whether the constrained Markov
decision process framework and results such as Lemma 4 could be extended to
different outcome types. One could think of a generalization of the model in
Section 2 along the lines in Yi and Wang (2023), where a clinical trial with
general outcome types was described using a Markov decision process.
Furthermore, the assumption underlying the applicability of response-adaptive
procedures is that the inter-arrival times of trial participants are longer
than the follow-up times. If this is not the case, the constrained Markov
decision process needs to be able to deal with delayed responses. In
Williamson et al. (2021), the robustness of constrained randomized dynamic
programming with respect to delayed responses was demonstrated. It would be
interesting to investigate whether the same result can be shown for the
constrained Markov decision process procedure. Lastly, as blocked RA
procedures are often preferred in clinical trials (Merrell et al., 2022), it
would furthermore be interesting to incorporate blocked allocation in the
constrained Markov decision approach.

The current paper provides theoretical results and a method for finding an
optimal policy for a class of finite-horizon constrained Markov decision
processes that is more general than usually considered in the literature, as a
different expectation operator, corresponding to a different prior belief, can
be used in each constraint. Such optimisation problems arise naturally in
settings where multiple hypotheses are considered for the treatment effect,
e.g., when considering type I error control and power constraints. It is
possible to also apply Lemma 4 to multi-objective Markov decision processes,
where each objective is formulated under a different prior. In principle,
Markov decision process response-adaptive procedures follow from a partially
observable Markov decision framework, and hence it would be interesting to
investigate whether there are other settings involving partial observability
where the solution approach can be applied.

## References

  * Bhatt (2010) Arun Bhatt.  Evolution of clinical research: a history before and beyond James Lind.  _Perspectives in clinical research_ , 1(1):6–10, 2010.
  * Palmer (1993) C. R. Palmer.  Ethics and statistical methodology in clinical trials.  _Journal of Medical Ethics_ , 19(4):219–222, 1993\.  doi:10.1136/jme.19.4.219.
  * Heilig and Weijerb (2005) Charles M. Heilig and Charles Weijerb.  A critical history of individual and collective ethics in the lineage of Lellouch and Schwartz.  _Clinical Trials_ , 2(3):244–253, 2005.  doi:10.1191/1740774505cn084oa.
  * Rosenberger (1996) William F. Rosenberger.  New directions in adaptive designs.  _Statistical Science_ , 11(2):137 – 149, 1996\.  doi:10.1214/ss/1038425657.
  * Antognini and Giovagnoli (2015) Alessandro Baldi Antognini and Alessandra Giovagnoli.  _Adaptive designs for sequential treatment allocation_.  CRC Press, Boca Raton, first edition, 2015.
  * Villar and Jacko (2022) Sofia S. Villar and Peter Jacko.  _Bandit Procedures for Designing Patient-Centric Clinical Trials_ , pages 365–389.  Springer International Publishing, New York, Cham, 2022.
  * Thompson (1933) William R. Thompson.  On the likelihood that one unknown probability exceeds another in view of the evidence of two samples.  _Biometrika_ , 25(3–4):285–294, 1933.  doi:10.1093/biomet/25.3-4.285.
  * aBerry and Viele (2023) Scott M. aBerry and Kert Viele.  Comment: Response Adaptive Randomization in Practice.  _Statistical Science_ , 38(2):229 – 232, 2023\.  doi:10.1214/23-STS865F.
  * US FDA (2019) US FDA.  Adaptive design clinical trials for drugs and biologics. Draft guidance for industry., 2019.  URL https://www.fda.gov/media/78495/download.  accessed 9-9-2023.
  * Robertson et al. (2023) David S. Robertson, Kim May Lee, Boryana C. López-Kolkovska, and Sofía S. Villar.  Response-Adaptive Randomization in Clinical Trials: From Myths to Practical Considerations.  _Statistical Science_ , 38(2):185 – 208, 2023\.  doi:10.1214/22-STS865.
  * Williamson (2020) S. Faye Williamson.  _Bayesian bandit models for the design of clinical trials_.  2020\.  Doctoral thesis, University of Lancaster.
  * Gittins (2018) J. C. Gittins.  Bandit Processes and Dynamic Allocation Indices.  _Journal of the Royal Statistical Society: Series B (Methodological)_ , 41(2):148–164, 2018.  doi:10.1111/j.2517-6161.1979.tb01068.x.
  * Kaufmann et al. (2012) Emilie Kaufmann, Olivier Cappe, and Aurelien Garivier.  On Bayesian Upper Confidence Bounds for Bandit Problems.  In Neil D. Lawrence and Mark Girolami, editors, _Proceedings of the Fifteenth International Conference on Artificial Intelligence and Statistics_ , volume 22 of _Proceedings of Machine Learning Research_ , pages 592–600, 2012.
  * Bather (1981) J. A. Bather.  Randomized Allocation of Treatments in Sequential Experiments.  _Journal of the Royal Statistical Society: Series B (Methodological)_ , 43(3):265–283, 1981.  doi:10.1111/j.2517-6161.1981.tb01672.x.
  * Thall and Wathen (2007) Peter F. Thall and J. Kyle Wathen.  Practical Bayesian adaptive randomisation in clinical trials.  _European Journal of Cancer_ , 43(5):859–866, 2007.  doi:10.1016/j.ejca.2007.01.006.
  * Bradt et al. (1956) R. N. Bradt, S. M. Johnson, and S. Karlin.  On Sequential Designs for Maximizing the Sum of $n$ Observations.  _The Annals of Mathematical Statistics_ , 27(4):1060–1074, 1956.  doi:10.1214/aoms/1177728073.
  * Villar et al. (2015) Sofía S. Villar, Jack Bowden, and James Wason.  Multi-armed Bandit Models for the Optimal Design of Clinical Trials: Benefits and Challenges.  _Statistical Science_ , 30(2):199–215, 2015.  doi:10.1214/14-STS504.
  * Williamson et al. (2017) S. Faye Williamson, Peter Jacko, Sofía S. Villar, and Thomas Jaki.  A Bayesian adaptive design for clinical trials in rare diseases.  _Computational Statistics & Data Analysis_, 113:136–153, 2017.  doi:10.1016/j.csda.2016.09.006.
  * Hardwick (1995) Janis P. Hardwick.  A modified bandit as an approach to ethical allocation in clinical trials.  _Lecture Notes-Monograph Series_ , 25:65–87, 1995.
  * .Berry and Eick (1995) Donald A. .Berry and Stephen G. Eick.  Adaptive assignment versus balanced randomization in clinical trials: A decision analysis.  _Statistics in Medicine_ , 14(3):231–246, 1995\.  doi:10.1002/sim.4780140302.
  * Cheng and Berry (2007) Yi Cheng and Donald A. Berry.  Optimal adaptive randomized designs for clinical trials.  _Biometrika_ , 94(3):673–689, 2007.  doi:10.1093/biomet/asm049.
  * Williamson et al. (2021) S. Faye Williamson, Peter Jacko, and Thomas Jaki.  Generalisations of a Bayesian decision-theoretic randomisation procedure and the impact of delayed responses.  _Computational Statistics & Data Analysis_, 174:1–26, 2021.  doi:10.1016/j.csda.2021.107407.
  * Hardwick and Stout (1991) Janis P. Hardwick and Quentin F. Stout.  Bandit strategies for ethical sequential allocation.  _Computing Science and Statistics_ , 23:421–424, 1991.
  * Merrell et al. (2022) David Merrell, Thevaa Chandereng, and Yeonhee Park.  A Markov decision process for response-adaptive randomization in clinical trials.  _Computational Statistics & Data Analysis_, 178:1–11, 2022.  doi:10.1016/j.csda.2022.107599.
  * Yi and Wang (2023) Yanqing Yi and Xikui Wang.  A Markov decision process for response adaptive designs.  _Econometrics and Statistics_ , 25:125–133, 2023.  doi:10.1016/j.ecosta.2021.10.015.
  * Altman (1999) Eitan Altman.  _Constrained Markov decision processes_.  Routledge, New York, first edition, 1999.
  * Jacko (2019) P. Jacko.  BinaryBandit: An Efficient Julia Package for Optimization and Evaluation of the Finite-Horizon Bandit Problem with Binary Responses.  Management Science Working Paper 2019:4, 2019.  URL https://eprints.lancs.ac.uk/id/eprint/136340/1/Jacko2019_binarybandit_wp.pdf.  Accessed 13-7-2023.
  * Yi (2013) Yanqing Yi.  Exact statistical power for response adaptive designs.  _Computational Statistics & Data Analysis_, 58:201–209, 2013.  doi:10.1016/j.csda.2012.09.003.
  * Agresti (1992) Alan Agresti.  A Survey of Exact Inference for Contingency Tables.  _Statistical Science_ , 7(1):131 – 153, 1992\.  doi:10.1214/ss/1177011454.
  * Brown et al. (1987) B. W. Brown, J. Herson, E. N. Atkinson, and M. E. Rozell.  Projection from previous studies: A Bayesian and frequentist compromise.  _Controlled Clinical Trials_ , 8(1):29–44, 1987\.  doi:10.1016/0197-2456(87)90023-7.
  * Bowden and Trippa (2017) Jack Bowden and Lorenzo Trippa.  Unbiased estimation for response adaptive clinical trials.  _Statistical Methods in Medical Research_ , 26(5):2376–2388, 2017.  doi:10.1177/0962280215597716.
  * Puterman (2014) Martin L. Puterman.  _Markov decision processes: discrete stochastic dynamic programming_.  John Wiley & Sons, Hoboken, 2014.
  * Kuhn and Tucker (2014) Harold W. Kuhn and Albert W. Tucker.  Nonlinear programming.  In _Traces and emergence of nonlinear programming_ , pages 247–258. Springer Basel, Basel, first edition, 2014.
  * Kelley (1960) J. E. Kelley.  The cutting-plane method for solving convex programs.  _Journal of the Society for Industrial and Applied Mathematics_ , 8(4):703–712, 1960.  doi:10.1137/010805.
  * Rosenberger and Hu (2004) William F. Rosenberger and Feifang Hu.  Maximizing power and minimizing treatment failures in clinical trials.  _Clinical Trials_ , 1(2):141–147, 2004.  doi:10.1191/1740774504cn016oa.
  * Fan and Glynn (2023) Lin Fan and Peter W. Glynn.  The fragility of optimized bandit algorithms.  https://arxiv.org/abs/2109.13595, 2023.  accessed 29-9-2023.

##  Appendix A Additional plots

###  A.1 Application 1: Control of power and type I error

Figure 7: Patient benefit (subfigure a), rejection rate (subfigure b), bias
(subfigure c), and mean squared error (subfigure d) vs. $\theta_{\text{D}}$
for $\theta_{\text{C}}=0.5$, $n=75$, $\alpha^{*}=0.05$, $\beta=0.4$, and RA
procedures ER (solid), DP (dashed), CRDP (dotted) and CMDP-T with $p=1.00$
(dash-dotted)  Figure 8: Patient benefit (subfigure a), rejection rate
(subfigure b), bias (subfigure c), and mean squared error (subfigure d) vs.
$\theta_{\text{D}}$ for $\theta_{\text{C}}=0.5$, $n=200$, $\alpha^{*}=0.07$,
$\beta=0.23$, and RA procedures ER (solid), DP (dashed), CRDP (dotted) and
CMDP-T with $p=1.00$ (dash-dotted)  Figure 9: Patient benefit (subfigure a),
rejection rate (subfigure b), bias (subfigure c), and mean squared error
(subfigure d) vs. $\theta_{\text{D}}$ for $\theta_{\text{C}}=0.25$, $n=75$,
$\alpha^{*}=0.05$, $\beta=0.4$, and RA procedures ER (solid), DP (dashed,
long), CRDP (dotted), CMDP-T with $p=0.95$ (dash-dotted), CMDP-T with $p=1.00$
(dashed, short) Figure 10: Patient benefit (subfigure a), rejection rate
(subfigure b), bias (subfigure c), and mean squared error (subfigure d) vs.
$\theta_{\text{D}}$ for $\theta_{\text{C}}=0.25$, $n=200$, $\alpha^{*}=0.07$,
$\beta=0.23$, and RA procedures ER (solid), DP (dashed, long), CRDP (dotted),
CMDP-T with $p=0.95$ (dash-dotted), CMDP-T with $p=1.00$ (dashed, short)
Figure 11: Patient benefit (subfigure a), rejection rate (subfigure b), bias
(subfigure c), and mean squared error (subfigure d) vs. $\theta_{\text{D}}$
for $\theta_{\text{C}}=0.75$, $n=75$, $\alpha^{*}=0.05$, $\beta=0.4$, and RA
procedures ER (solid), DP (dashed, long), CRDP (dotted), CMDP-T with $p=0.95$
(dash-dotted), CMDP-T with $p=1.00$ (dashed, short) Figure 12: Patient benefit
(subfigure a), rejection rate (subfigure b), bias (subfigure c), and mean
squared error (subfigure d) vs. $\theta_{\text{D}}$ for
$\theta_{\text{C}}=0.75$, $n=200$, $\alpha^{*}=0.07$, $\beta=0.23$, and RA
procedures ER (solid), DP (dashed, long), CRDP (dotted), CMDP-T with $p=0.95$
(dash-dotted), CMDP-T with $p=1.00$ (dashed, short)

###  A.2 Application 2: Control of estimation error

Figure 13: Patient benefit (subfigure a), rejection rate (subfigure b), bias
(subfigure c), and mean squared error (subfigure d) vs. $\theta_{\text{D}}$
for $\theta_{\text{C}}=0.5$, $n=75$, and RA procedures ER (solid), CMDP-E2
with $p=1.00$ (dashed), CRDP (dotted) and CMDP-E1 with $p=1.00$ (dash-dotted)
Figure 14: Patient benefit (subfigure a), rejection rate (subfigure b), bias
(subfigure c), and mean squared error (subfigure d) vs. $\theta_{\text{D}}$
for $\theta_{\text{C}}=0.5$, $n=200$, and RA procedures ER (solid), CMDP-E2
with $p=1.00$ (dashed), CRDP (dotted) and CMDP-E1 with $p=1.00$ (dash-dotted)
Figure 15: Patient benefit (subfigure a), rejection rate (subfigure b), bias
(subfigure c), and mean squared error (subfigure d) vs. $\theta_{\text{D}}$
for $\theta_{\text{C}}=0.25$, $n=75$, and RA procedures ER (solid), CMDP-E2
with $p=0.95$ (dashed, moderate), CMDP-E2 with $p=1.00$ (dashed, short), CRDP
(dotted), CMDP-E1 with $p=0.95$ (dash dot), CMDP-E2 with $p=1.00$ (dashed,
long) Figure 16: Patient benefit (subfigure a), rejection rate (subfigure b),
bias (subfigure c), and mean squared error (subfigure d) vs.
$\theta_{\text{D}}$ for $\theta_{\text{C}}=0.25$, $n=200$, and RA procedures
ER (solid), CMDP-E2 with $p=0.95$ (dashed, moderate), CMDP-E2 with $p=1.00$
(dashed, short), CRDP (dotted), CMDP-E1 with $p=0.95$ (dash dot), CMDP-E2 with
$p=1.00$ (dashed, long) Figure 17: Patient benefit (subfigure a), rejection
rate (subfigure b), bias (subfigure c), and mean squared error (subfigure d)
vs. $\theta_{\text{D}}$ for $\theta_{\text{C}}=0.75$, $n=75$, and RA
procedures ER (solid), CMDP-E2 with $p=0.95$ (dashed, moderate), CMDP-E2 with
$p=1.00$ (dashed, short), CRDP (dotted), CMDP-E1 with $p=0.95$ (dash dot),
CMDP-E2 with $p=1.00$ (dashed, long) Figure 18: Patient benefit (subfigure a),
rejection rate (subfigure b), bias (subfigure c), and mean squared error
(subfigure d) vs. $\theta_{\text{D}}$ for $\theta_{\text{C}}=0.75$, $n=200$,
and RA procedures ER (solid), CMDP-E2 with $p=0.95$ (dashed, moderate),
CMDP-E2 with $p=1.00$ (dashed, short), CRDP (dotted), CMDP-E1 with $p=0.95$
(dash dot), CMDP-E2 with $p=1.00$ (dashed, long)