# Conditional Power and Friends: The Why and How of (Un)planned, Unblinded
Sample Size Recalculations in Confirmatory Trials

Kevin Kunzmann
MRC Biostatistics Unit
University of Cambridge
Michael J. Grayling
Population Health Sciences Institute
Newcastle University
Kim May Lee
Pragmatic Clinical Trials Unit
Queen Mary University of London
David S. Robertson
MRC Biostatisics Unit
University of Cambridge
Kaspar Rufibach
Methods, Collaboration, and Outreach Group (MCO)
Department of Biostatistics
F. Hoffmann-La Roche, Basel
James M. S. Wason
Population Health Sciences Institute
Newcastle University
and
MRC Biostatistics Unit
University of Cambridge  MRC Biostatistics Unit, University of Cambridge,
Cambridge Institute of Public Health, Forvie Site, Robinson Way, Cambridge
Biomedical Campus, Cambridge CB2 0SR, United Kingdom, kevin.kunzmann@mrc-
bsu.cam.ac.uk

###### Abstract

Adapting the final sample size of a trial to the evidence accruing during the
trial is a natural way to address planning uncertainty. Designs with adaptive
sample size need to account for their optional stopping to guarantee strict
type-I error-rate control. A variety of different methods to maintain type-I
error-rate control after unplanned changes of the initial sample size have
been proposed in the literature. This makes interim analyses for the purpose
of sample size recalculation feasible in a regulatory context. Since the
sample size is usually determined via an argument based on the power of the
trial, an interim analysis raises the question of how the final sample size
should be determined conditional on the accrued information. Conditional power
is a concept often put forward in this context. Since it depends on the
unknown effect size, we take a strict estimation perspective and compare
assumed conditional power, observed conditional power, and predictive power
with respect to their properties as estimators of the unknown conditional
power. We then demonstrate that pre-planning an interim analysis using
methodology for unplanned interim analyses is ineffective and naturally leads
to the concept of optimal two-stage designs. We conclude that unplanned design
adaptations should only be conducted as reaction to trial-external new
evidence, operational needs to violate the originally chosen design, or post
hoc changes in the objective criterion. Finally, we show that commonly
discussed sample size recalculation rules can lead to paradoxical outcomes and
propose two alternative ways of reacting to newly emerging trial-external
evidence.

Keywords: Adaptive Design - Conditional Power - Interim Analysis - Optimal
Design - Predictive Power - Sample Size Recalculation

##  1 Introduction

The planning phase of confirmatory clinical trials is typically characterized
by substantial uncertainty about the magnitude of the parameters underlying
the hypothesis of interest. Often the alternative hypothesis is the
superiority of a new intervention arm over a treatment-as-usual control arm
and is formalized via the average treatment effect, e.g., the difference of
means between the two randomly allocated treatment arms. Clinical data
collection is expensive and time consuming leading to a strong economic
incentive to reach the study goals with as little data as possible. The
conventional statistical criteria to determine the sample size of a trial are
a one-sided type-I error rate of 2.5% and a power of 80% or 90%. Since power
is a function of the unknown effect size the initial design must be specified
under substantial uncertainty about the magnitude of the effect which the
trial hopes to detect eventually. Mainly, two ways of addressing this
challenge have been put forward in the literature.

Firstly, the so-called ‘hybrid’ approach to sample size derivation takes a
Bayesian view on determining the initial sample size of a clinical trial ? .
The Bayesian component requires the specification of an informative prior on
the parameters of interest. This then allows reasoning about the ‘expected
power’ of a trial as a function of its sample size and, consequently,
determination of the sample size such that the expected power exceeds a target
threshold. This concept is a straight-forward extension of the usual practice
of computing the power under a fixed point alternative, which can easily be
recovered as a special case when considering a point prior. The advantage lies
in the fact that the a priori information about the magnitude of the effect
size is faithfully reflected in the sample size derivation. Since the actual
analysis is still frequentist, type-I error rate control is not compromised.

Secondly, authors have proposed to apply the concept of adaptive design
changes to recalculate the initial sample size of a design based on data
observed within the trial itself ? . The rationale behind this approach is to
use the accruing evidence about the unknown parameters driving the sample size
derivation to ‘correct’ the sample size mid-trial. During the interim
analysis, the accrued data can either be unblinded or not. Methods for the
latter ‘blinded interim analysis’ are particularly useful if relevant nuisance
parameters such as the variance are also unknown ? ? . We focus on the
unblinded case which allows for a more precise interim assessment of the
effect size than methods that retain the blind. The maximal type-I error rate
constraint is then usually protected by applying the conditional error
principle ? ? or an equivalent formulation via $\alpha$-spending or
combination functions ? . Often the sample size of the current trial is
adjusted such that the conditional power given the data observed up to the
interim analysis again exceeds the initial threshold for unconditional power ?
. Any recalculation based on conditional power arguments must address the
problem that conditional power, just as unconditional power, depends on the
unknown underlying parameters. It must thus be estimated to inform a sample
size recalculation. Yet, precise estimation of the relevant parameters before
the conclusion of a trial is hard since only a fraction of the final sample
size is available. Three approaches to estimating conditional power have been
discussed in the literature.

Firstly, ‘conditional power’ is the probability to reject the null hypothesis
given the interim data as a function of the unknown parameters. In a slight
abuse of terminology, the same term is often used to refer to the assumed
conditional power that is obtained by plugging in the point alternative used
for the initial sample size derivation. Evidently, this assumed conditional
makes no use of the accrued trial data since the effect size is kept fixed.
Secondly, authors have proposed to use the ‘observed conditional power’
instead, which replaces the unknown parameters with their maximum-likelihood
estimates. The latter ‘plug-in’ approach is often criticized for ‘… [using]
the interim estimate of the effect […] twice …’ Bauer et al., p. 330 ? and ? .
Thirdly, conditional power can be evaluated as Bayesian expected power
conditional on the observed interim data, i.e., by averaging conditional power
as function of the unknown parameters with respect to the posterior density
after conducting the interim analysis. Within the hybrid Bayesian framework
this is usually referred to as ‘predictive power’ ? ? .

The purpose of this paper is to argue that pre-planned sample size adaptations
using methodology intended for unplanned interim analyses to react to within-
trial interim data are inefficient and unnecessary if the original design was
planned optimally. We then discuss cases where an unplanned design adaptation
might still be warranted and propose two criteria that allow such an
adaptation that is consistent with the original design. We first review
assumed conditional-, observed conditional-, and predictive power from an
estimation perspective. For sake of simplicity, we use a simple single-arm
one-stage design to illustrate core characteristics. We then discuss the risks
of constructing an ‘adaptive’ two-stage design by naïvely applying a
conditional-power-based sample size recalculation rule and the conditional
error principle. The drawbacks of this naïve approach directly lead to the
concept of optimal-two stage designs and we build on ideas presented in ? ? ,
and ? to derive optimal designs for the situation studied in this manuscript.
Without loss of generality, we focus on minimal expected sample size as
optimality criterion. Optimal two-stage designs, by definition, cannot be made
more efficient by sample size recalculation. The optimality of such designs
does, however, depend on the initial trial-external evidence that feeds into
the planning assumptions. This trial-external evidence might change during an
ongoing study and may thus mandate a sample size recalculation. In Section 5,
we discuss when an unplanned sample size recalculation is reasonable, how such
a recalculation interacts with the concept of optimal two-stage designs, and
the practical implementation.

In the following, we assume that the interest lies in testing a new treatment
for efficacy. For the sake of simplicity, we consider the case of a single-arm
trial. All considerations can easily be extended to the more practically
relevant two-arm case. We further assume that the individually observed
outcomes of the study participants $X_{i},i=1,\ldots,n$ are $iid$ and that
their distribution has finite first moment $\theta$ and unit variance. Again,
all considerations can be extended to the case of generic known variances and,
at least approximately, to the case of unknown variance. A suitable test
statistic for the null hypothesis of interest ${\mathcal{H}_{0}:\theta\leq 0}$
is

$\displaystyle Z_{n}=\frac{1}{\sqrt{n}}\sum_{i=1}^{n}X_{i}\ .$ (1)

Invoking the central limit theorem,
$Z_{n}\stackrel{{\scriptstyle\cdot}}{{\sim}}\mathcal{N}(\sqrt{n}\,\theta,1)$
and, on the boundary of the null hypothesis,
$Z_{n}\stackrel{{\scriptstyle\cdot}}{{\sim}}\mathcal{N}(0,1)$.

Further assuming a maximal permissible type-I error rate of $\alpha=0.025$
(which we use for the remainder of the paper), the critical value for a
single-stage fixed-sample-size design is the $1-\alpha$ quantile of the
standard normal distribution, i.e., approximately $c=1.96$. The test then
rejects $\mathcal{H}_{0}$ if and only if $Z_{n}>c$ after the outcomes of $n$
subjects have been observed. The required size of the trial $n$, given a
maximal allowable type-I error rate, is usually determined by some form of
restriction on the (minimal) statistical power of the test. Approaches to
defining such a power constraint under different assumptions were reviewed in
? .

##  2 Monitoring power

No matter which rationale was used to derive the overall sample size $n$,
after observing $m$, $0<m<n$, outcomes, an independent data monitoring
committee might be interested to learn about the prospects of eventually
rejecting the null hypothesis. Both the rejection probability conditional on
the presence of an effect (conditional power) or the joint probability of
rejection and presence of an effect (conditional probability of success) could
be of interest. Here, we focus on the former due to the central role of power
(and conditional power) in the planning of clinical trials. This effectively
constitutes an estimation problem where the estimand is the _conditional_
probability of rejecting the null hypothesis given $Z_{m}=z_{m}$ and
$\theta>0$. Here the conditioning on a positive effect is crucial to justify
the notion of ‘power’ in contrast to the unconditional ‘probability of
success’ ? .

Since $Z_{n}$ and $Z_{m}$ are (asymptotically) jointly normal, the conditional
distribution of $Z_{n}$ given $Z_{m}$ and $\theta$ is again normal and given
by

$\displaystyle\mathcal{L}_{\theta}\big{[}\,Z_{n}\mid
Z_{m}=z_{m}\,\big{]}\stackrel{{\scriptstyle\cdot}}{{=}}\mathcal{N}\Big{(}\sqrt{n}\,\theta+\sqrt{\tau}\big{(}z_{m}-\sqrt{m}\,\theta\big{)},\
1-\tau\Big{)}$ (2)

where $\tau:=m/n$ is the ‘information fraction’ at the interim analysis.

The probability of rejecting the null hypothesis at the end of the trial as a
function of $\theta$ is referred to as the ‘conditional power’ in the
literature. It is defined as

$\displaystyle\operatorname{CP}(z_{m},c,\theta):$
$\displaystyle=\operatorname{Pr}_{\theta}[\,Z_{n}>c\,|\,Z_{m}=z_{m}\,]$ (3)
$\displaystyle=1-\Phi\left(\frac{c-\sqrt{n}\,\theta-\sqrt{\tau}\,z_{m}+\sqrt{\tau}\sqrt{m}\,\theta}{\sqrt{1-\tau}}\right)$
(4)
$\displaystyle=1-\Phi\left(\frac{c/\sqrt{m}+\sqrt{\tau}\,\big{(}\theta-\widehat{\theta}_{m}\big{)}-\theta/\sqrt{\tau}}{\sqrt{(n-m)/(n\,m)}}\right)\
.$ (5)

Ignoring the implicit dependency on the time point of the interim analysis,
$m$, conditional power is thus a function of the observed effect
$\widehat{\theta}_{m}:=z_{m}/\sqrt{m}$ (or equivalently the observed test
statistic $z_{m}$), the critical value for the test decision $c$, and the
standardized effect size $\theta$. Since $\theta$ is unknown, so is
$\operatorname{CP}(z_{m},c,\theta)$ and it thus cannot be evaluated directly
upon observing $Z_{m}=z_{m}$. Yet, as a function of the unknown quantity
$\theta$, $\operatorname{CP}(z_{m},c,\theta)$ can be estimated from observed
data. The estimation perspective on ‘evaluating’ conditional power is less
commonly taken in the literature but provides a consistent framework to
compare the characteristics of different methods ? .

Firstly, conditional power can be estimated based on a fixed point alternative
$\theta_{1}$. In a slight abuse of terminology, this quantity is often also
referred to as ‘conditional power’. To clearly distinguish it from
$\operatorname{CP}(z_{m},c,\theta)$ we denote it ‘assumed conditional power’

$\displaystyle\operatorname{ACP}(z_{m},c):=\operatorname{CP}\big{(}\,z_{m},c,\theta_{1}\,\big{)}\
.$ (6)

Secondly, the observed effect $\widehat{\theta}_{m}$ can be used as a plug-in
estimator for the unknown effect size ? . This quantity is sometimes referred
to as ‘observed conditional power’ in the literature and is defined as

$\displaystyle\operatorname{OCP}(z_{m},c):=\operatorname{CP}\big{(}\,z_{m},c,\widehat{\theta}_{m}\,\big{)}\
.$ (7)

Thirdly, a Bayesian approach can be taken if one is willing to quantify the
uncertainty about the unknown parameter $\theta$ by modeling it as the
realization of a random variable $\Theta\sim\varphi(\cdot)$ where
$\varphi(\theta)$ is the prior probability density function evaluated at the
parameter value $\theta$. Our definition of this so-called ‘predictive power’

$\displaystyle\operatorname{PP}_{\varphi}(z_{m},c):$
$\displaystyle=\operatorname{Pr}_{\varphi}[\,Z_{n}>c\,|\,\Theta>0,Z_{m}=z_{m}\,]$
(8)
$\displaystyle=\int_{0}^{\infty}\operatorname{CP}\big{(}z_{m},c,\theta\big{)}\,\varphi\big{(}\theta\,|\,Z_{m}=z_{m},\Theta\geq
0\big{)}\operatorname{d}\theta$ (9)

differs slightly from the one proposed by ? in that we condition on a positive
effect size $\Theta>0$. This is more consistent with the notion of frequentist
power as discussed in ? , although the difference is only of practical
relevance when a substantial fraction of the a priori probability mass is
concentrated on the null hypothesis.

We now compare $\operatorname{ACP},\operatorname{OCP}$, and
$\operatorname{PP}$ by means of a concrete example. To this end we assume that
the available prior information can be summarized in a truncated normal prior
with density

$\displaystyle\varphi(\theta):=\boldsymbol{1}_{[-0.5,1]}(\theta)\,\frac{\displaystyle\phi\left(\frac{\theta-0.4}{0.2}\right)}{\displaystyle\Phi\left(\frac{1-0.4}{0.2}\right)-\Phi\left(\frac{-0.5-0.4}{0.2}\right)}\
,$ (10)

(see Figure 1)

Figure 1: Assumed prior density function. The gray area indicates the null
hypothesis, $\mathcal{H}_{0}$, of no effect.

and that all positive effect sizes are clinically relevant. A truncated normal
prior is convenient since it allows the analytic computation of the posterior
distribution under a normal likelihood. It is also the maximum entropy
distribution on a compact interval given mean and standard deviation and thus
‘least informative’ given a lower and upper boundaries on plausible effects as
well as the location and vagueness (variance) of the prior. Following ? the
required sample size $n$ is then determined by requiring a minimal expected
power of $1-\beta=80\%$, where expected power is defined as

$\displaystyle\operatorname{EP}_{\varphi}(n,c):=\operatorname{Pr}\big{[}\,Z_{n}>c\,|\,\Theta\geq
0\,\big{]}=\int\operatorname{Pr}_{\theta}\big{[}\,Z_{n}>c\,\big{]}\,\varphi(\theta\,|\,\Theta\geq
0)\,\operatorname{d}\theta\ .$ (11)

Under these assumptions, the required sample size is $n=79$. No point
alternative is necessary to derive this sample size. We thus use the prior
mean rounded to the first decimal point to evaluate $\operatorname{ACP}$ and
define $\theta_{1}:=0.4$.

The three estimators are depicted in Figure 2(a) as functions of the observed
interim outcome, where the interim time-point $m=26$.

(a) Estimates as function of the interim data.

(b) Histograms of the sampling distributions.

(c) Bias, mean absolute error (MAE) and mean squared error (MSE).

Figure 2: Properties of $\operatorname{ACP},\operatorname{OCP}$, and
$\operatorname{PP}$ as estimators of the unknown conditional power at $m=26$
an overall sample size $n=79$.

The most sensitive measure is $\operatorname{OCP}$ while $\operatorname{ACP}$
is the least sensitive (see also spread of the sampling distributions in
Figure 2(b)). This is due to $\operatorname{ACP}$ assuming a fixed effect size
and thus only being affected by direct changes in $z_{m}$. Both
$\operatorname{OCP}$ and $\operatorname{PP}$, however, are also indirectly
affected by updating the belief about the effect size with the interim
results. The plug-in estimate $\operatorname{OCP}$ assumes that the observed
effect is the true effect and does not quantify the uncertainty around this
value in any way. $\operatorname{PP}$, on the other hand, invokes Bayes’
theorem and then integrates $\operatorname{CP}(z_{m},c,\theta)$ with respect
to the obtained posterior distribution. The degree of adaptation of
$\operatorname{PP}$ thus depends on the vagueness of the chosen prior. In
fact, as the prior approaches a point mass at $\theta_{1}$,
$\operatorname{PP}$ converges to $\operatorname{ACP}$. $\operatorname{ACP}$
can thus be seen as a mere special case of $\operatorname{PP}$ with a point
prior on $\theta_{1}>0$.

A fundamental difference between $\operatorname{ACP}$ and $\operatorname{PP}$
on the one hand and $\operatorname{OCP}$ on the other hand is that
$\operatorname{OCP}$ is the only estimator that does not condition on
$\Theta$. In contrast, $\operatorname{ACP}$ implicitly conditions on
$\Theta=\theta_{1}>0$ and $\operatorname{PP}$ on ${\Theta>0}$. The practical
implication of this can be seen by looking at the sampling distributions of
the three estimators under different assumptions about the value of $\theta$
(see Figure 2(b)). Since $\operatorname{OCP}$ does not condition on
$\Theta>0$, the sampling distribution is much more left-skewed for small
effects than for the other two estimators. Also, the negligence of the
sampling variation of $\widehat{\theta}_{m}$, which $\operatorname{OCP}$ plugs
in the expression for conditional power, leads to a larger variance for
intermediate values of $\theta$ and the characteristic U-shape ? . This high
variance of $\operatorname{OCP}$ directly translates to a relatively large
mean absolute error (MAE) and mean squared error (MSE) when estimating
conditional power, see Figure 2(c). $\operatorname{PP}$ is the posterior
expectation of conditional power with respect to the chosen prior and thus
minimizes the quadratic Bayes risk, i.e. the MSE. This leads to relative good
precision for parameter values with high a priori likelihood (around
$\theta=0.4$). If one wanted to minimize the MAE directly, the posterior
median would be optimal. The principle, however, remains the same and the
posterior mean is more consistent with the derivation of the initial sample
size using expected power. $\operatorname{ACP}$ is clearly the best in terms
of bias and MAE/MSE for values close to or above $\theta_{1}$, but its
performance as estimator of the conditional power quickly deteriorates for
small effect sizes.

$\operatorname{ACP}$ is merely an extreme special case of $\operatorname{PP}$.
$\operatorname{OCP}$ is both hard to justify theoretically and has inferior
precision for parameter values with a high a priori likelihood.
$\operatorname{PP}$ is thus the natural choice for monitoring power during the
course of a trial. Depending on the prior chose, it’s properties can be either
more similar to $\operatorname{ACP}$ (low prior variance) or
$\operatorname{OCP}$ (high prior variance).

##  3 Naïve unplanned sample size adaptations

Monitoring the predictive power of an ongoing trial naturally raises the
question of whether this information can be leveraged to improve the operating
characteristics of a trial. Typically, sample size recalculation is considered
in the context of group-sequential designs at the pre-specified interim
analysis ? ? ? . From a statistical perspective, there is no reason why a
sample size recalculation should not be performed for a study that initially
was planned using a single-stage design. For the sake of simplicity the
following considerations thus assume a simple single-stage design with fixed
$n$ and $c$ as starting point and are restricted to a single interim analysis.
All results can be generalized to the case of more complex starting designs
and multiple interim analyses.

A method for sample size recalculation that is often discussed in the
literature is to derive a new sample size $n^{\prime}$ and a new critical
value $c^{\prime}$ such that some estimate of conditional power exceeds a
threshold $1-\beta_{cond.}$. Different choices for the estimator of
conditional power and the choice of $\beta_{cond.}$ have been discussed in the
literature ? ? ? . Strict overall type-I error rate control can be maintained
by invoking the conditional error principle ? ? ? , i.e., by limiting the
maximal type-I error rate of the new design to the maximal conditional type-I
error rate of the original design given the observed $z_{m}$. Often, trial
protocols leave the exact choice of $1-\beta_{cond.}$ open since it can be
chosen ad hoc without compromising strict type-I-error rate control. The basic
concept of readjusting the sample size to achieve the desired conditional
power is, however, a common approach, see for instance, ? and ? .

In accordance with the discussion in the previous section and to be consistent
with the derivation of the initial sample size of $n=79$ in the example
considered earlier, we use predictive power as an estimator of conditional
power. Furthermore, we set $\beta_{cond.}=\beta=0.2$. For given $Z_{m}=z_{m}$,
the recalculation rule then corresponds to solving the optimization problem

$\displaystyle\operatornamewithlimits{argmin:\ }_{n^{\prime},c^{\prime}}$
$\displaystyle n^{\prime}$ (12) $\displaystyle\operatorname{subject\ to:\ }$
$\displaystyle\operatorname{CP}(m,n^{\prime},z_{m},c^{\prime},0)$
$\displaystyle\ \leq\ \operatorname{CP}_{n}(m,n,z_{m},c,0)$ (13)
$\displaystyle\operatorname{PP}_{\varphi}(m,n^{\prime},z_{m},c^{\prime})$
$\displaystyle\ \geq\ 1-\beta_{cond.}$ (14) $\displaystyle n^{\prime}$
$\displaystyle\ \geq\ n_{min}>m$ (15) $\displaystyle n^{\prime}$
$\displaystyle\ \leq\ n_{max}\ .$ (16)

Here, we indicate the dependency on the final sample size and the interim time
point explicitly by redefining

$\displaystyle\operatorname{CP}(m,n,z_{m},c,\theta):$
$\displaystyle=\operatorname{Pr}_{\theta}[\,Z_{n}>c\,|\,Z_{m}=z_{m}\,]$ (17)
$\displaystyle\operatorname{PP}_{\varphi}(m,n,z_{m},c):$
$\displaystyle=\operatorname{Pr}_{\varphi}[\,Z_{n}>c\,|\,\Theta>0,Z_{m}=z_{m}\,]\
.$ (18)

Constraint (13) implements the conditional error principle by limiting the
conditional (type-I) error rate under the new design to the conditional
(type-I) error rate under the original design. The minimal sample size
constraint (15) allows $n_{min}>m$ for cases where a minimal sample size upon
rejection of the null hypothesis is deemed necessary. This could e.g. be of
interest if a minimal precision of the final maximum likelihood estimator is
sought. Note that the trial cannot stop at $n^{\prime}=m$ without immediately
accepting the null hypothesis since the conditional error under the new design
in the case of early acceptance would be $1$ but is always smaller than $1$
for the original single-stage design. Imposing a maximal sample size of
$n_{max}$ via constraint (16) is a practical necessity since the recalculated
sample size could otherwise tend to infinity as $z_{m}$ approaches negative
infinity. In cases where this constraint prevents a solution (low observed
effect), the trial is usually stopped early for futility declaring that the
null hypothesis cannot be rejected. Alternatively, one could continue with the
maximal sample size accepting the fact that the conditional power constraint
is not met. This would, however, raise the question as to why a hard
constraint on conditional power was imposed in the first place.

If the recalculation rule defined in (12)–(16) is made mandatory in the study
protocol, it yields two functions $n(z_{m})$ and $c(z_{m})$ as the point-wise
solution of the problem which jointly define an ‘adaptive design’ for any
$z_{m}$. Here, ‘mandatory’ means that it is decided a priori to always (for
all $z_{m}$) recalculate the final sample size and critical value in the above
specified way after a fixed number of outcomes $m$ has been observed. The term
‘adaptive design’ is slightly misleading in this context. The design itself is
pre-specified (and thus not adapted or changed) but only $n$ and $c$ are
adaptive since they vary as functions of $z_{m}$. We would like to emphasize
that such a detailed specification of the sample size recalculation rule is
not common in practice since most protocols are focused on guaranteeing
control of the maximal type-I-error rate control and maintaining maximal
flexibility during the interim analysis. Still, if a sample-size reassessment
was prescribed in the study protocol, there must have been some rationale as
to the goal of that reassessment and we merely investigate the consequences
for a design’s unconditional properties when applying such a rationale in a
mandatory way. The described approach should thus be seen as a hypothetical
example highlighting the ineffectiveness that can arise from a naïve
implementation of a recalculation based on a fixed threshold for conditional
power. The corresponding sample size and critical value functions for $m=26$
are depicted in Figure 3.

Figure 3: Original single-stage design; naïvely adapted design with sample
size $n^{\prime}(z_{m})$ and critical value $c^{\prime}(z_{m})$ functions
defined by the conditional optimization problem (12)–(16) and
${n=79},{m=26},{n_{max}=160},{n_{min}=30}$ and $1-\beta_{cond.}=0.8$; the
optimal design is the solution of (19)–(21).

The mandatory application of the adaptation rule results in a two-stage
design. The final sample size is thus a random variable $n(Z_{m})$. This
implies that the objective criterion of ‘minimal sample size’ is no longer
meaningful. Instead, one could use a weighted sum,
$\operatorname{E}[\,n(Z_{m})\,]+\eta\sqrt{\operatorname{V}[\,n(Z_{m})\,]}$, of
the expected value and standard deviation of the sample size as the objective
criterion. The rationale for adding a penalty depending on the standard
deviation of the sample size is that it might incur additional costs due to
more complicated logistics (e.g., on demand production of additional drug
doses). For the single-stage design this proposal always reduces to the fixed
sample size $n$ since the standard deviation is 0. Whether or not the two-
stage design is considered better then depends on the choice of $\eta$, i.e,
the relative weight of the the standard deviation in the objective. In the
particular situation considered here, the expected sample size is $49.5$ and
its standard deviation is $29.7$. Since the original $n$ was $79$, the two-
stage design would be considered ‘better’ than the single-stage design for
$\eta>0.99$.

This comparison, however, ignores the fact that the original unconditional
constraint of and expected power of more than $80\%$ is no longer fulfilled by
the design with mandatory sample size adaptation. Both expected power
($72.1\%$ _versus_ $80.0\%$) and maximal type-I error rate ($2.1\%$ _versus_
$2.5\%$) are lower than for the original design. The difference in operating
characteristics is mainly due to the fact that the new design implicitly
introduces a binding early-stopping boundary for futility when the
recalculation problem does not yield a $n^{\prime}\leq n_{max}$. This effect
is less pronounced when the same recalculation rule is applied to an initial
two-stage group-sequential design which already defines binding early stopping
boundaries for futility and efficacy.

Similar to suggestions in the literature on binding futility stopping, one
could modify the initial nominal $\alpha$ and the conditional power threshold
$\beta_{cond.}$ until the recalculated design again fulfills the required
operating characteristics ? . This approach to define a fully pre-specified
two-stage design is, however, unnecessarily complicated in that it uses
methods like the conditional error principle or a combination function, which
are originally intended for _unplanned_ design adaptations.

##  4 Optimal pre-planned sample size adaptations

Evidently, the naïve application of the predictive power adaptation rule
outlined earlier is needlessly conservative. Fixing the operating
characteristics by tuning the nominal $\alpha$ and $\beta_{cond.}$ does not
address the deeper issue that the sample size and critical value functions are
still derived based on a merely heuristic recalculation criterion and that
they still implicitly depend on an original design that is never actually
realized via the conditional error constraint (13).

Instead, the interim time-point $m$, the sample size function $n(z_{m})$, the
critical value function $c(z_{m})$, and the early futility- and efficacy
boundaries ($n(z_{m})=m$ and $c(z_{m})=\pm\infty$) can be optimized directly.
For the sample size function $n(z_{m})$ ? discussed this approach by fitting a
polynomial of degree 4 but they did not optimize any of the other parameters.
In a more general form, ? approached the problem using variational methods. An
R implementation that uses cubic splines for both the $n$ and $c$ functions is
available via the package adoptr ? .

In the following we restrict the considerations to the case where only
expected sample size is of interest and the variation is ignored, i.e.
$\eta=0$. The corresponding optimal two-stage design is the solution of

$\displaystyle\operatornamewithlimits{argmin:\ }_{m,n(\cdot),c(\cdot)}$
$\displaystyle\operatorname{E}_{\varphi}[\,n(Z_{m})\,]$ (19)
$\displaystyle\operatorname{subject\ to:\ }$
$\displaystyle\operatorname{Pr}_{0}[\,Z_{n(Z_{m})}>c(Z_{m})\,]$
$\displaystyle\leq\alpha$ (20)
$\displaystyle\operatorname{EP}_{\varphi}\big{(}n(\cdot),c(\cdot)\big{)}$
$\displaystyle\geq 1-\beta$ (21)

and can be derived numerically using the R package adoptr ? . Here,
$\operatorname{E}_{\varphi}[\,n(Z_{m})\,]$ is the expected sample size under
the prior $\varphi$ and expected power needs to be redefined to account for
the fact that $n$ and $c$ are now functions of the interim results

$\displaystyle\operatorname{EP}_{\varphi}\big{(}n(\cdot),c(\cdot)\big{)}:=\operatorname{Pr}_{\varphi}\big{[}\,Z_{n(Z_{m})}>c(Z_{m})\,|\,\Theta\geq
0\,\big{]}\ .$ (22)

Figure 3 shows the sample size-, critical value-, and the predictive power
function in the example situation discussed earlier together with the single-
stage design and the naïve adaptation based on predictive power discussed in
Section 3. The optimal two-stage design complies with both the predictive
power and the maximal type-I error rate constraints and is thus comparable
with the original single-stage design in terms of sample size. The expected
sample size is much lower ($56.4$ _versus_ $79$) at the cost of a non-zero
standard deviation of the expected sample size ($28.5$ _versus_ $0$). Since it
is entirely pre-specified, the optimal two-stage design does not need to fall
back to the conditional error principle to control the maximal type-I error
rate. Instead, it achieves the desired design characteristics in an optimal
way since they are incorporated as constraints to problem (19)–(21). The
optimal design fully exhausts the allowable maximal type-I error rate and
complies with the expected power constraint. This comes at the cost of a
predictive power that can drop as low as 40% close to the futility boundary as
shown in the right panel of Figure 3. The sample size function of the optimal
design is also characteristically different from the naïve design. The latter
is convex on the continuation region whereas the optimal shape has a mode
close to the early-futility boundary. A recalculation based on exceeding a
fixed lower threshold for some estimator of conditional power always results
in a convex sample size function and can thus never be fully optimal
irrespective of how the nominal values of $\alpha$ and $\beta_{cond.}$ are
chosen.

The fact that the optimal two-stage design exhibits a monotonously increasing
predictive power rather than a constant predictive power of, say, 80%
indicates that any recalculation rule keeping predictive power close to a
fixed target value must be inefficient in terms of minimizing expected sample
size. The mandatory application of the heuristic recalculation rule introduced
in Section 3 would change the optimal design for almost every value of $z_{m}$
and thus increase expected sample size of the resulting design. It may thus
appear as if any ‘recalculation’ of an optimal design’s sample size using
methods for unplanned interim analyses was not only unnecessary but even
counterproductive. This is, however, only the case if the initial design was
optimal and the planning assumptions still hold.

##  5 Consistent unplanned recalculation

The preceding section made the case that the direct optimization of all design
parameters is superior to the mandatory application of methods for unplanned
design adaptations during the planning phase of a trial. By definition, no
optimal design ever needs to be changed as a reaction to trial internal data
during the course of the study. However, there are still situations in which
an unplanned sample size recalculation is warranted. The optimality of a
design typically depends on the chosen prior density $\varphi$ and the
objective function. At any point in time, the prior density encodes all trial-
external evidence about the effect size. This might change over the course of
a trial. For instance, another study might publish new results that trigger a
reassessment of the considerations leading to the choice of $\varphi$.

In the following, we discuss two approaches of reacting to such trial-external
events in a consistent way if the original design minimizes expected sample
size. Here, ‘consistent’ means that the original optimal design is invariant
under the recalculation rule unless the planning assumptions (prior) are
changed or the time-point of the interim analysis is different from the
optimal one. The considerations are restricted to cases where the result of an
adaptation is a final sample size and no further interim analyses are planned
for the remainder of the trial. This means that the unplanned adaptation
replaces the pre-planned interim analysis although it might occur at a
different point in time. All methods can be generalized to an unplanned
interim analysis in the first stage with a subsequent second interim analysis.
This might be necessary when the unplanned interim analysis has to be
conducted shortly after the start of a trial.

To obtain a consistent sample size recalculation rule, the conditional error
principle can be applied not only to the maximal type-I error rate but also to
the (average) type-II error rate. Let $m^{\prime}$ be the sample size of the
unplanned interim analysis and $\psi$ be the probability density function of
the revised prior. We propose to recalculate the new sample size, $n^{\prime}$
and the new critical value $c^{\prime}$ as the solution of the point-wise
optimization problem

$\displaystyle\operatornamewithlimits{argmin:\ }_{n^{\prime},c^{\prime}}$
$\displaystyle n^{\prime}$ (23) $\displaystyle\operatorname{subject\ to:\ }$
$\displaystyle\operatorname{CP}(m^{\prime},n^{\prime},z_{m^{\prime}},c^{\prime},0)$
$\displaystyle\leq\operatorname{Pr}_{0}[\,Z_{n(Z_{m})}>c(Z_{m})\,|\,Z_{m^{\prime}}=z_{m^{\prime}}\,]$
(24)
$\displaystyle\operatorname{PP}_{\psi}(m^{\prime},n^{\prime},z_{m^{\prime}},c^{\prime})$
$\displaystyle\geq$
$\displaystyle\operatorname{Pr}_{\varphi}[\,Z_{n(Z_{m})}>c(Z_{m})\,|\,\Theta>0,Z_{m^{\prime}}=z_{m^{\prime}}\,]\
.$ (25)

Here, constraint (25) means that the modified design under the new prior has
at least as much predictive power as the original design under the original
prior. Constraint (24) implements the conditional error principle and thus
ensures strict type-I error rate control of the procedure. The crucial
difference to the procedure discussed in Section 5 is the fact that the lower
boundary on predictive power varies with $z_{m}$.

In cases where the interim analysis is conducted at the original $m$, the
problem simplifies to

$\displaystyle\operatornamewithlimits{argmin:\ }_{n^{\prime},c^{\prime}}$
$\displaystyle n^{\prime}$ (26) $\displaystyle\operatorname{subject\ to:\ }$
$\displaystyle\operatorname{CP}(m,n^{\prime},z_{m},c^{\prime},0)$
$\displaystyle\leq\operatorname{CP}(m,n(z_{m}),z_{m},c(z_{m}),0)$ (27)
$\displaystyle\operatorname{PP}_{\psi}(m,n^{\prime},z_{m},c^{\prime})$
$\displaystyle\geq\operatorname{PP}_{\varphi}(m,,n(z_{m}),z_{m},c(z_{m}))\ .$
(28)

Since $\operatorname{CP}$ is monotone in $n^{\prime}$, the solution is already
uniquely defined by the constraints, if it exists. It follows directly that,
for $m^{\prime}=m,\psi=\varphi$, the original design ${n^{\prime}=n(z_{m})}$,
${c^{\prime}=c(z_{m})}$ fulfills both constraints with equality. The original
optimal design is thus indeed invariant under the proposed recalculation rule
if neither the planning prior not the interim time-point are changed.

Problem (24)–(25) is in no way the only conceivable means of deriving a
consistently recalculated sample size. For the sake of simplicity, we consider
the case of $m^{\prime}=m$ (i.e. problem (26)–(28)). Depending on how much the
prior changes, the required sample size adjustment to meet the predictive
power constraint can be quite drastic. In practice, it might thus be
preferable to relax the predictive power constraint in a principled way. This
can be done via the Lagrangian mechanism. Let

$\displaystyle f(n^{\prime},c^{\prime}):$ $\displaystyle=n^{\prime}$ (29)
$\displaystyle g(n^{\prime},c^{\prime}):$
$\displaystyle=\operatorname{CP}(m,n^{\prime},z_{m},c^{\prime},0)-\operatorname{CP}(m,n(z_{m}),z_{m},c(z_{m}),0)$
(30) $\displaystyle h_{\psi}(n^{\prime},c^{\prime}):$
$\displaystyle=\operatorname{PP}_{\varphi}(m,,n(z_{m}),z_{m},c(z_{m}))-\operatorname{PP}_{\psi}(m,n^{\prime},z_{m},c^{\prime})$
(31)

As discussed earlier, the original designs $n$ and $c$ trivially are the
unique solution of

$\displaystyle\operatornamewithlimits{argmin:\ }_{n^{\prime},c^{\prime}}$
$\displaystyle f(n^{\prime},c^{\prime})$ (32)
$\displaystyle\operatorname{subject\ to:\ }$ $\displaystyle
g(n^{\prime},c^{\prime})$ $\displaystyle\leq 0$ (33) $\displaystyle
h_{\varphi}(n^{\prime},c^{\prime})$ $\displaystyle\leq 0\ .$ (34)

This problem is equivalent to solving the unconstrained problem

$\displaystyle\operatornamewithlimits{argmin:\
}_{n^{\prime},c^{\prime}}f(n^{\prime},c^{\prime})-\lambda_{g}g(n^{\prime},c^{\prime})-\lambda_{h}h_{\varphi}(n^{\prime},c^{\prime})$
(35)

for suitable choice of the Lagrange multipliers $\lambda_{g}$ and
$\lambda_{h}$. Since the original $n,c$ are optimal, $\lambda_{g}$ and
$\lambda_{h}$ can be obtained directly as solution of

$\displaystyle 0$ $\displaystyle=\nabla f(n,c)-\lambda_{g}\nabla
g(n,c)-\lambda_{h}\nabla h_{\varphi}(n,c)$ (36)
$\displaystyle\Leftrightarrow(1,0)^{\top}$ $\displaystyle=\big{(}\nabla
g(n,c),\nabla h_{\varphi}(n,c)\big{)}\cdot(\lambda_{g},\lambda_{h})^{\top}$
(37) $\displaystyle\Leftrightarrow(\lambda_{g},\lambda_{h})^{\top}$
$\displaystyle=\big{(}\nabla g(n,c),\nabla
h_{\varphi}(n,c)\big{)}^{-1}\cdot(1,0)^{\top}\ .$ (38)

Here, $\nabla$ denotes the gradient as column vector. The Lagrange multiplier
corresponding to the predictive power constraint, $\lambda_{h}$, encodes the
trade-off between predictive power and sample size of the original, optimal
design at $z_{m}$.

A recalculated sample size under a different prior $\psi$ can then be obtained
by keeping $\lambda_{h}$ fixed and minimizing

$\displaystyle\operatornamewithlimits{argmin:\ }_{n^{\prime},c^{\prime}}$
$\displaystyle n^{\prime}$
$\displaystyle-\lambda_{h}\,h_{\psi}(n^{\prime},c^{\prime})$ (39)
$\displaystyle\operatorname{subject\ to:\ }$ $\displaystyle
g(n^{\prime},c^{\prime})$ $\displaystyle\leq 0\ .$ (40)

This effectively relaxes the predictive power constraint by allowing the same
trade-off between sample size and predictive power that was implied by the
original optimal design. The conditional type-I error rate constraint is not
relaxed to comply with the conditional error principle. This approach is
closely related to the formulation of the optimal-design problem given in ?
and the utility perspective on sample size derivation presented in ? . In
contrast to the situation analyzed by ? , however, $\lambda_{h}$ does depend
on $z_{m}$. This is due to the fact that expected sample size is formed with
respect to the unconditional prior while expected power uses the prior
conditional on $\Theta>0$.

We illustrate both approaches by example of the optimal two-stage design that
minimizes the objective $\operatorname{E}_{\varphi}[\,n(Z_{m})\,]$ subject to
a maximal type-I error rate of $0.025$ and a minimal expected power of $0.8$,
see Figure 3. This design has a pre-planned interim analysis after $m=35$
individuals. Now assume that the prior $\varphi$ needs to be revised due to
new trial-external information at the pre-planned interim analysis. We
consider members of the class of shifted truncated Normal distributions

$\displaystyle\psi_{\mu}(\theta):=\frac{\displaystyle\boldsymbol{1}_{[-0.5,1]}(\theta)\,\phi\left(\frac{\theta-\mu}{0.2}\right)}{\displaystyle\Phi\left(\frac{1-\mu}{0.2}\right)-\Phi\left(\frac{-0.5-\mu}{0.2}\right)}$
(41)

as revised priors, i.e. the original prior $\varphi=\psi_{0.4}$. Figure 4(a)
shows the behavior of the recalculated sample size and critical value in this
scenario for observed effects of $z_{m}/\sqrt{m}=0.25,0.3,0.35$ and varying
mean of the revised prior if the originally planned interim time-point is
maintained.

(a) Interim analysis at the originally planned interim sample size of $m=35$
for varying observed effect.

(b) Earlier-than-planned interim analysis at an observed effect of $0.3$ and
varying interim time-point.

(c) Interim analysis at the originally planned interim sample size of $m=35$
for varying observed effect and the $\lambda$ approach.

Figure 4: Comparison of consistent unplanned recalculation schemes.

Clearly, the original design is indeed invariant if the prior remains
unchanged. The sample size is increased if the new prior is more conservative
than $\varphi$ and vice versa.

Of course, it could also be necessary to conduct the interim analysis earlier
or later than originally anticipated, at $m^{\prime}\neq m$. In this case, the
design is no longer invariant since the earlier or later time point already
implies a deviation from the original optimal design even if the prior remains
the same. Figure 4(b) shows the resulting adaptations for earlier-than-planned
time points of the interim analysis while the observed effect is kept fixed at
$0.3$. The formula for the conditional error in this situation is given in the
Appendix.

Figure 4(c) shows the recalculated sample size for $m^{\prime}=m=35$ and the
$\lambda$-approach. The original design is again (approximately) invariant
under recalculation. Compared to the fixed conditional type-II error approach,
predictive power under the $\lambda$-approach is more flexible and might drop
below the predictive power under the original design, see Figure 4(c), right
panel. This flexibility in terms of predictive power leads to a much smaller
deviation from the original design’s sample size. The $\lambda$-approach
transitions smoothly from an increase of the original sample size due to a
slightly more conservative prior to a decrease under very conservative priors.
This behavior leads to less extreme adjustments and is thus more realistic in
practice. After all, a slight loss in (predictive) power as compared to the
original design might be tolerable if the sample size can be kept at a
moderate level in exchange. Of course, it is still possible to stop the trial
early for futility at any point without compromising type one error rate
control if the predictive power is considered too low to warrant a
continuation.

##  6 Discussion

Monitoring the power of the remainder of an ongoing trial constitutes an
estimation problem since the conditional power is a function of the unknown
treatment effect. An often overlooked aspect of this estimation problem is
that it is conditional on a positive treatment effect since the notion of
power is conditional on there being an non-null effect ? . In many practical
situations where the prior mass is concentrated on non-null effects the effect
of conditioning is negligible but its importance increases with the vagueness
of a the prior. This lack of conditioning on a positive treatment effect is
one of the reasons observed conditional power performs rather poorly as an
estimator of the unknown conditional power. Assumed conditional power for
$\theta_{1}>0$ does reflect this conditionality. It is, however, merely a
special case of predictive power with a point prior on the effect size. In
situations where the effect size is fairy well known, a point prior might
constitute an acceptable, simpler approximation. In any other case, the
Bayesian predictive power allows a priori information to be incorporated in a
more fine-grained way and guarantees optimal mean-squared-error performance.
If an unconditional measure for the probability to reject the null hypothesis
is sought, a conditional version of the probability of success (or assurance)
can be derived in an analogue way ? .

Great care should be taken when predictive power or another estimator of
conditional power is used to modify an ongoing trial’s sample size. Altering a
simple starting design with a seemingly intuitive recalculation rule can
actually lead to overall _less effective_ trials. Techniques originally
intended for _unplanned_ adaptations of trials should not be used for a _pre-
planned_ sample size recalculation. Instead, an optimal two-stage design
should directly be derived for the objective criterion of interest. By
definition, no trial-internal event can then justify a sample size
recalculation. Only trial external events, such as a change in the objective
criterion, the emergence of new trial-external evidence, or unforeseen
deviations from the planned interim analysis time-point may require a
reassessment of the sample size of an optimal design. A generic recalculation
rule might lead to the paradoxical situation that even an optimal starting
design is always modified during the interim analysis \- even if the planning
assumptions remain unaltered throughout the trial. This is clearly ineffective
and, consequently, an optimal starting design should be invariant under a
sensible adaptation rule if the planning assumptions remain unchanged. For
instance, this minimal consistency property is not fulfilled when
recalculating a design’s sample size based on a fixed threshold for its
minimal conditional power since designs optimizing expected sample size tend
to have variable conditional power.

We propose two consistent approaches to adjusting an ongoing optimal design to
newly emerging trial-external data. The first method applies the conditional
error principle to the (average) type-II error rate in a similar way to its
use in controlling the maximal type-I error rate. The method is easy to
implement and consistent in the above defined way. We would like to stress
that this is by no means the only way of conducting a consistent sample size
recalculation. The alternative $\lambda$-approach allows a trade-off between
predictive power and sample size leading to smaller sample size deviations
from the original design. This comes at the risk of potentially loosing
predictive power. Both methods differ by the quantity that is held fixed
during the recalculation. For the first method, the predictive power under the
old design and the old prior is maintained for the new design under the new
prior whereas the $\lambda$-approach only keeps the trade-off between
predictive power and sample size fixed but allows the actual predictive power
to deviate form the original design.

These considerations show how absolutely crucial the initial planning stage of
a trial is. Adaptive methods should not be taken as an excuse to start with a
sub-optimal design and rely on a later sample size recalculation. Ideally, all
uncertainty is quantified during the planning phase to the best possible
extend and integrated in the initial design via a planning prior. This does
not mean that a complex design is always the best choice. If the variability
of the final sample size and the operational burden of conducting interim
analyses is penalized strong enough in the objective criterion, a simple one-
stage design might very well perform better than more complex alternatives.
Also, optimal-two stage designs require optimizing over function-spaces to
find the optimal sample size and critical value functions. Obtaining a stable
solution might thus be hard in practice and group-sequential designs are a
viable approximation since optimizing the stage-wise sample sizes and stopping
boundaries only requires optimizing over a small set of real parameters. It is
well-known, that optimal group-sequential designs approximate the performance
of optimal two-stage designs with variable sample sizes sufficiently well ? .
Still, the same principle considerations apply: An optimal (group-sequential)
design only needs to be revised if either the objective function or the
underlying planning assumptions, i.e., the prior on the effect size, changes.

Adaptations of initially sub-optimal designs are another potential application
of the new recalculation scheme. Clearly, the single-stage design with $n=79$
is not minimizing _expected_ sample size. To switch to a near-optimal adaptive
design mid-trial, one could derive the optimal two-stage design given the
current best prior post-hoc and then use the $\lambda$-approach to recalculate
the initial design while controlling the maximal conditional type-I error rate
in the usual way. The resulting design will not be optimal from an
unconditional perspective since the conditional error function of the initial,
sub-optimal design has to be respected but the power-sample size trade-off
during the interim analysis will be similar to the optimal design due to the
choice of $\lambda_{h}$.

A limitation of the proposed recalculation methods is the fact that they rely
on the conditional error principle for strict type-I-error rate control. The
conditional error principle can be difficult to extend to cases with nuisance
parameters ? . In practice, a simple plug-in approach similar to the approach
in blinded sample size reassessment might be viable but have not yet been
investigated more thoroughly. Alternatively, any other form of pre-specifying
a combination function for the stage-wise p values might be used to control
the type-I-error rate. This combination function should then also be optimized
over during the planning stage to avoid inconsistencies. In particular, the
optimal combination function of a two stage-design minimizing expected sample
size can be approximated with an inverse normal combination test ? .

### Data Availibility Statement

Data sharing is not applicable to this article as no datasets were generated
and the research is entirely theoretical. The source code is available at
github.com and zenodo.org.

## Bibliography

  * [1] Spiegelhalter DJ, Freedman LS, Parmar MK. Bayesian approaches to randomized trials. Journal of the Royal Statistical Society: Series A (Statistics in Society) 1994; 157(3): 357–387.
  * [2] Bauer P, Bretz F, Dragalin V, König F, Wassmer G. Twenty-five years of confirmatory adaptive designs: Opportunities and pitfalls. Statistics in Medicine 2016; 35(3): 325–47.
  * [3] Birkett MA, Day SJ. Internal pilot studies for estimating sample size. Statistics in Medicine 1994; 13(23-24): 2455–63.
  * [4] Müller HH, Schäfer H. A general statistical principle for changing a design any time during the course of a trial. Statistics in Medicine 2004; 23(16): 2497–508.
  * [5] Brannath W, Gutjahr G, Bauer P. Probabilistic foundation of confirmatory adaptive designs. Journal of the American Statistical Association 2012; 107(498): 824–32.
  * [6] Proschan MA, Hunsberger SA. Designed extension of studies based on conditional power. Biometrics 1995; 51(4): 1315–24.
  * [7] Bauer P, Koenig F. The reassessment of trial perspectives from interim data–A critical view. Statistics in Medicine 2006; 25(1): 23–36.
  * [8] Brannath W, Bauer P. Optimal conditional error functions for the control of conditional power. Biometrics 2004; 60(3): 715–23.
  * [9] Pilz M, Kunzmann K, Herrmann C, Rauch G, Kieser M. A variational approach to optimal two-stage designs. Statistics in Medicine 2019; 38(21): 4159–71.
  * [10] Kunzmann K, Pilz M, Herrmann C, Rauch G, Kieser M. The adoptr package: Adaptive optimal designs for clinical trials in R. Journal of Statistical Software 2020; accepted.
  * [11] Kunzmann K, Grayling MJ, Lee KM, Robertson DS, Rufibach K, Wason J. A review of Bayesian perspectives on sample size derivation for confirmatory trials. arXiv preprint arXiv:2006.15715 2020\.
  * [12] Müller HH, Schäfer H. Adaptive group sequential designs for clinical trials: Combining the advantages of adaptive and of classical group sequential approaches. Biometrics 2001; 57(3): 886–91.
  * [13] Bhatt DL, Stone GW, Mahaffey KW, et al. Effect of platelet inhibition with cangrelor during PCI on ischemic events. New England Journal of Medicine 2013; 368(14): 1303–13.
  * [14] Mehta C, Gao P, Bhatt DL, Harrington RA, Skerjanec S, Ware JH. Optimizing trial design: Sequential, adaptive, and enrichment strategies. Circulation 2009; 119(4): 597–605.
  * [15] Jennison C, Turnbull BW. Adaptive sample size modification in clinical trials: Start small then ask for more?. Statistics in Medicine 2015; 34(29): 3793–810.
  * [16] Wassmer G, Brannath W. Group sequential and confirmatory adaptive designs in clinical trials.  Springer .  2016\.
  * [17] Gutjahr G, Brannath W, Bauer P. An approach to the conditional error rate principle with nuisance parameters. Biometrics 2011; 67(3): 1039–46.

## Appendix

The conditional error for an earlier-than-planned interim analysis
($m^{\prime}<m$) is

$\displaystyle\operatorname{Pr}_{\theta}\big{[}\,Z_{n(Z_{m})}>c(Z_{m})\,|\,Z_{m^{\prime}}=z_{m^{\prime}}\,\big{]}$
(42) $\displaystyle=$
$\displaystyle\int\operatorname{CP}_{\theta}\big{(}m,n(z_{m}),z_{m},c(z_{m}),\theta\big{)}\
p_{\theta}(z_{m}\,|\,z_{m^{\prime}})\operatorname{d}z_{m}$ (43)

where $p_{\theta}(z_{m}\,|\,z_{m^{\prime}})$ denotes the probability density
function of the outcome at $m$ given the partial observed outcome up to
$m^{\prime}$ (see equation (2)).


# The hadronic nucleus-nucleus cross section and the nucleon size

Govert Nijs  Center for Theoretical Physics, Massachusetts Institute of
Technology, Cambridge, MA 02139, USA  Wilke van der Schee  Theoretical Physics
Department, CERN, CH-1211 Genève 23, Switzerland

###### Abstract

Even though the total hadronic nucleus-nucleus cross section is among the most
fundamental observables, it has only recently been measured precisely for
lead-lead collisions at the LHC. This measurement implies the nucleon width
should be below 0.7 fm, which is in contradiction with all known state-of-the-
art Bayesian estimates. We study the implications of the smaller nucleon width
on quark-gluon plasma properties such as the bulk viscosity. The smaller
nucleon width dramatically improves the description of several triple-
differential observables.

††preprint: CERN-TH-2022-103//MIT-CTP/5445

Introduction - The understanding of the creation quark-gluon plasma (QGP) as
created at colliders such as the Large Hadron Collider (LHC) in Geneva
requires the understanding of several stages of the collision of heavy ions
Heinz and Snellings (2013); Busza  _et al._ (2018). The first stage is far-
from-equilibrium and involves an initial condition together with an
understanding how this evolves towards a hydrodynamic QGP Berges  _et al._
(2021). Secondly, there is a hydrodynamic stage, which involves understanding
the temperature-dependent transport coefficients such as the shear viscosity
Romatschke and Romatschke (2019). Lastly, the QGP undergoes particlization
into a gas of interacting hadrons which can then be detected experimentally.

Relatively little is known about the initial stage and many analyses use a
phenomenological parameterization called the TRENTo model Moreland  _et al._
(2015). Here a heavy ion is composed of a superposition of nucleons with
Gaussian distributions of energy of width $w$, while fixing the interaction
rate to fit the hadronic nucleon-nucleon cross section $\sigma_{NN}$
determined from $pp$ collisions. Hence, interactions of nucleons with smaller
$w$ depend strongly on the impact parameter, while large $w$ nucleons do not
always interact even if they collide at zero impact parameter. Interestingly,
all recent state-of-the-art global analyses of a wide variety of experimental
data have preferred a large nucleon width in fm of $0.98\pm 0.18$ Moreland
_et al._ (2020), $0.96\pm 0.05$ Bernhard  _et al._ (2019), $0.94\pm 0.18$ Nijs
_et al._ (2021a), $1.05\pm 0.13$ Everett  _et al._ (2021) or $0.82\pm 0.23$
Nijs and van der Schee (2021). Since $w$ equals the Gaussian width this
implies that the resulting energy profile is then much larger than the proton
charge radius of $0.841\,$fm Pohl  _et al._ (2010), or the radius implied by
$\sigma_{NN}=68\,$mb in the black disk approximation ($r=0.74\,$fm).

In this Letter we show that the recent ALICE measurement of the PbPb total
hadronic cross section $\sigma_{AA}{}$ of $7.67\pm 0.24\,$b at $\sqrt{s_{\rm
NN}}=5.02\,$TeV ALICE (2022) implies a nucleon width smaller than
approximately $0.7\,$fm, which is smaller than the width from all quoted
Bayesian estimates. This measurement hence raises two important questions.
First, why did the Bayesian probability estimates not result in the correct
nucleon width? Second, what are the implications of this smaller nucleon
width? Part of the answer to the first question must be an inaccurate estimate
of the systematic uncertainty covariance matrix. Here we note that the
covariance matrix does not only include the systematic and statistical
experimental and theoretical uncertainties, but it is essential to also
include correlations or anti-correlations between observables. For the second
question, the smaller nucleon width implies a larger bulk viscosity. Finally,
we will show the improved analysis implies a better description of
statistically difficult triple-differential observables.

The initial condition and the cross section -  All computations of
$\sigma_{AA}$ start with the Monte Carlo Glauber model Alver  _et al._ (2008);
Loizides  _et al._ (2015, 2018); d’Enterria and Loizides (2021). In this model
nucleons are placed randomly according to a Woods-Saxon distribution,
$\rho(r)=\rho_{0}/(1+e^{(r-R)/a})$, specified by the half-density radius $R$
and the diffusivity $a$ with potentially the requirement of a minimum nucleon-
nucleon distance $d_{\rm min}$. In the current implementation the parameters
$R$ and $a$ follow the point density distributions that are determined
separately for protons and neutrons using electron-ion scattering experiments
Klos  _et al._ (2007); Tarbert  _et al._ (2014); Loizides  _et al._ (2018).
$\sigma_{AA}$ is then determined by the condition that at least a single
nucleon-nucleon interaction occurs.

In the black disk approximation nucleons interact if their distance $d$
satisfies $d<\sqrt{\sigma_{NN}/\pi}$, with $\sigma_{NN}$ the nucleon-nucleon
cross section as determined from $pp$ collisions at the same collision energy.
In Pythia 8 a normalized overlap function is specified as a function of the
impact parameter $b$ as $T_{pp}(b)\propto\exp(-(b/w)^{m})$ with $w$ the
nucleon width and $m=1.85$ for the Monash tune Skands  _et al._ (2014);
d’Enterria and Loizides (2021). In this work we will use the TRENTo model
Moreland  _et al._ (2015), which uses a Gaussian overlap function ($m=2$) and
the nucleon width $w$ as a parameter. Nucleons then interact with probability
$P(b)=1-\exp(-\sigma_{gg}\int d^{2}x_{T}T_{A}(x_{T})T_{B}(x_{T}))$, whereby
$\sigma_{gg}$ is determined by $\sigma_{NN}{}$. In this work we use nucleons
that are composed of smaller constituents (as in Moreland  _et al._ (2020));
we verified that this does not affect $\sigma_{AA}$.

It is no surprise that the (traditionally used) black-disk approximation
produces the smallest $\sigma_{AA}$. Indeed, all models by construction have
an equal nucleon-nucleon cross section, but wider overlap functions allow a
nucleon-nucleon interaction to occur more easily for a nucleus-nucleus
collision at relatively large impact parameter.

Figure 1: We show the dependence of the PbPb cross section on the nucleon-
nucleon cross section $\sigma_{\text{NN}}$, the nucleon width $w$, the
centrality normalization $\text{cent}_{\text{norm}}$ and the minimal nucleon-
nucleon distance $d_{\text{min}}$. The parameters other than the one varied
are kept fixed at the value indicated by the ALICE data point ALICE (2022).

In Fig. 1 we see that $\sigma_{AA}$ can increase up to 23% for $w$ as large as
$1.2\,$fm. Perhaps surprisingly, the dependence on $\sigma_{NN}$ is fairly
mild, and in this work we will keep the measured value of $61.2$ and
$67.6\,$mb for 2.76 and 5.02 TeV collisions respectively Loizides  _et al._
(2018). We also observe an almost negligible dependence on $d_{\rm min}$. The
cross section depends linearly on the centrality normalization Nijs and van
der Schee (2021), e.g. on which events to count as a collision (as in the
models above), or experimentally on how many collisions are recorded. Both
theoretically and experimentally this contains an uncertainty, which we
include as a separate parameter in Fig. 1 (bottom-left). Motivated by ALICE
(2022) we give this parameter a prior probability distribution for our
Bayesian analysis of a Gaussian with unit width. We note that this parameter
propagates into a significant uncertainty mainly for peripheral spectra and
multiplicities, but also for more central elliptic flow coefficients.

Crucially, $\sigma_{AA}$ only depends sensitively on $w$, $\sigma_{NN}$ and
the centrality normalization, whereby the latter two are well constrained
experimentally. A measurement of $\sigma_{AA}$ hence provides robust
constraints on $w$ with only weak theoretical modelling uncertainties. The
recently measured value of $7.67\pm 0.26\,$b implies $w\approx\text{0.4 --
0.5}\,$fm (see Fig. 1), which as noted in the introduction is in direct
contradiction with all state-of-the-art global analyses of heavy ion
collisions so far.

Figure 2: We show posterior distributions for the nucleon width $w$, initial
state fluctuation parameter $\sigma_{\text{fluct}}$, free streaming time
$\tau_{\text{fs}}$ and shear relaxation time $\tau_{\pi}$. We show the results
from four different fits, where we turn weighting on and off, and do or do not
fit to the cross section $\sigma_{\text{AA}}$.

Implications for QGP properties -  Fig. 2 shows the posterior parameter
analysis for the width $w$, the TRENTo fluctuation parameter $\sigma_{\rm
fluct}$, the free streaming time $\tau_{\text{fs}}$ and the second order hydro
relaxation time $\tau_{\pi}$ when including (green dotted) or excluding
(purple dot-dashed) $\sigma_{AA}$. This is the result of a global analysis of
653 experimental data points Aamodt  _et al._ (2011); Adam  _et al._ (2016a,
b); Acharya  _et al._ (2019, 2020); Abelev  _et al._ (2013, 2014a); Acharya
_et al._ (2018a); Adam  _et al._ (2016c); Acharya  _et al._ (2018b); ALICE
(2022); Khachatryan  _et al._ (2016) (see also the Supplemental Material) with
a 21-dimensional parameter space within the publicly available _Trajectum_ 1.3
framework tra  similar to Nijs  _et al._ (2021a, b); Nijs and van der Schee
(2021), but full details will be presented elsewhere Nijs and van der Schee .
The analysis with $\sigma_{AA}{}$ also includes the proton-lead cross section
collisions at $\sqrt{s_{\rm NN}}=5.02\,$TeV of $\sigma_{p\text{Pb}}=2.06\pm
0.08\,$b as measured earlier by CMS Khachatryan  _et al._ (2016). It is
striking that adding two experimental observations ($\sigma_{AA}$ for PbPb
ALICE (2022) and $p$Pb at 5.02 TeV Abelev  _et al._ (2014b); Khachatryan  _et
al._ (2016)) has such a strong effect on these parameters.

Nevertheless, we see that the nucleon width is still not quite compatible with
the value required from the $\sigma_{AA}$ measurement (the green dotted $w$ is
peaked around $0.7\,$fm). This is consistent with the previous Bayesian
analyses that strongly prefer larger widths and hence a Bayesian updated
estimate now reduces the width, but not all the way to be compatible. Given
the theoretical robustness of the $\sigma_{AA}$ comparison this remaining
discrepancy has to be taken seriously. A crude way would be to fix
$w=0.45\,$fm by hand or to use a prior distribution close to this value. A
more quantitative approach is to give higher weights to observables we have a
higher trust in on physical grounds, which is what we will attempt in this
work.

Figure 3: We show the temperature dependent specific shear and bulk viscosities $\eta/s$ and $\zeta/s$ using the posterior distributions from four different global analyses. This includes (blue) or excludes (red) $\sigma_{AA}$, and includes weighting (top) or does not use weighting (bottom). Unweighted viscosities are more precise (especially for $\zeta/s$), but do not necessarily give the correct nucleon width when considering $\sigma_{AA}{}$. |  $\sigma_{\text{PbPb}}$ [b] |  $\sigma_{p\text{Pb}}$ [b]
---|---|---
$\sigma_{\text{AA}}$ & weights | $8.02\pm 0.19$ | $2.20\pm 0.06$
weights | $8.95\pm 0.36$ | $2.48\pm 0.10$
$\sigma_{\text{AA}}$ | $8.19\pm 0.19$ | $2.25\pm 0.06$
neither | $8.83\pm 0.29$ | $2.45\pm 0.09$
ALICE/CMS | $7.67\pm 0.24$ | $2.06\pm 0.08$
Table 1: Posterior values for the PbPb and $p$Pb cross sections for the four
different fits, compared to the ALICE ALICE (2022) (PbPb) and CMS Khachatryan
_et al._ (2016) ($p$Pb) values. Theoretical emulation uncertainty is
negligible and the uncertainty comes almost entirely from the posterior
uncertainty on the nucleon width (see Fig. 2). Figure 4: We show the
correlation between elliptic flow $v_{2}\\{2\\}^{2}$ (left),
$v_{3}\\{2\\}^{2}$ (middle) and $v_{4}\\{2\\}^{2}$ (right) with mean
transverse momentum $\langle p_{T}\rangle$ for weighted fit including (solid)
and without (dotted) the cross section $\sigma_{AA}{}$ and also the unweighted
fit including $\sigma_{AA}{}$ (dashed) for PbPb collisions. We show
statistical uncertainties in gray and systematic uncertainties from the
posterior as a band (first case only). We also include ATLAS data with
systematic (boxes) and statistical uncertainties ATLAS (2022).

Alternatively, we could say that the inconsistency of the nucleon width before
and after $\sigma_{AA}{}$ implies that we underestimated the systematic
uncertainties. By far the best way would be to update those uncertainties,
importantly including all (anti-)correlations between all observables. In the
current analysis data points from nearby bins are treated as correlated (see
also Bernhard (2018)), although even here there is a certain level of
arbitrariness. Moreover, there are theoretical uncertainties from e.g.
particlization that are even in principle hard to quantify. Inaccurate
estimates of correlations can cause the affected class of observables to be
given too much or too little relative importance in the Bayesian analysis. The
weighting is a way to modify the relative importance of observable classes and
hence correct for this effect. This interpretation of extra uncertainty is
more applicable to the theoretical uncertainty (where model uncertainties are
hard to estimate), but we stress that even experimentally it is often not
clear how all systematic uncertainties are correlated among observables.

As mentioned we highly trust $\sigma_{AA}$ on grounds that it is relatively
model independent theoretically and in fact only strongly depends on the
nucleon width. Weaker but similar arguments can be made for integrated
unidentified observables, such as integrated multiplicities, mean transverse
momentum and integrated ansiotropic flow coefficients $v_{n}\\{k\\}$. Particle
identified observables are theoretically more difficult to model and here we
assign a weight $1/2$. We define the weight $\omega$ to mean that we multiply
the difference in an observable between theory and experiment by $\omega$.
Note that this preserves the correlation matrix. Also transverse momentum
($p_{T}$) differential observables are more model-dependent, especially at
larger $p_{T}$, and similar arguments can be made for observables in
peripheral centrality classes. We hence chose to weight $p_{T}$-differential
observables by an extra factor 1/2 and as well as an extra factor of
$(2.5-p_{T}[\text{GeV}])/1.5$ if $p_{T}>1\,$GeV. We also multiply the weight
for any observable by $(50-c[\%])/50$ if the centrality class $c$ is beyond
50%.

The posteriors including the weights are shown in Fig. 2 as blue solid (with
$\sigma_{AA}$) and red dashed (without $\sigma_{AA}$) curves. Indeed by using
lower weights for more model-dependent observables, we see the nucleon width
is in agreement with the estimate from $\sigma_{AA}$. Also clear is that as
expected the lower weights decrease the precision of the posterior parameter
determination, which is especially clear for $\tau_{\pi}$ and
$\tau_{\text{fs}}$. As shown in Table 1 both including $\sigma_{AA}{}$ and
including weights improves the theoretical postdiction of $\sigma_{AA}{}$ so
that the agreement is within 1.1 (weighted) or 1.7 (unweighted) standard
deviations from the ALICE result for PbPb and 1.5 and 1.9 standard deviations
from the CMS result for $p$Pb.

In Fig. 3 we show the temperature dependent specific shear and bulk
viscosities, $\eta/s$ and $\zeta/s$ including (blue) and without $\sigma_{AA}$
(red). Given the smaller nucleon width we expect a larger bulk and shear
viscosity to reduce the average radial and elliptic flow that is induced by
the larger radial gradient. This is indeed what is observed.

The nucleon width and observables -  After showing the $\sigma_{AA}$ updated
QGP properties an important question remains. Why do analyses without
$\sigma_{AA}$ favor a large nucleon width, and is there an inconsistency given
that the nucleon width is small? To answer this question Table 2 shows per
class of observables the average discrepancy with the experimental result.
Naturally, including the $\sigma_{AA}{}$ dramatically improves the agreement
with the $\sigma_{AA}{}$ postdiction. What is however perhaps surprising is
that the fit of the other observable classes only worsens mildly, on average
worsening from 0.83 to 0.89 standard deviations for the weighted case.
Virtually all observables get slightly worse, with the notable exception of
the mean $p_{T}$ fluctuations. Naturally, in the weighted case the observables
with a lower weight have a worse agreement, which was precisely the motivation
for introducing the weights. We note that these deviations do not trivially
translate into a $\chi^{2}$ value, since many observables are highly
correlated (see also Bernhard (2018)). A more complete overview of the match
of all observables is presented in the Supplemental Material, where indeed by
eye it is difficult to see the difference between the fit with and without
$\sigma_{AA}{}$.

|  $\sigma_{\text{AA}}$ & $\omega$ | $\omega$ | $\sigma_{\text{AA}}$ | neither | $\bar{\omega}$
---|---|---|---|---|---
$dN_{\text{ch}}/d\eta$ | 0.55 | 0.60 | 1.23 | 1.22 | 1.00
$dN_{\pi^{\pm},k^{\pm},p^{\pm}}/dy$ | 0.76 | 0.70 | 0.60 | 0.57 | 0.48
$dE_{T}/d\eta$ | 1.59 | 1.51 | 0.82 | 0.77 | 0.48
$\langle p_{T}\rangle_{\text{ch},\pi^{\pm},K^{\pm},p^{\pm}}$ | 0.66 | 0.60 | 0.88 | 0.72 | 0.46
$\delta p_{T}/\langle p_{T}\rangle$ | 0.56 | 0.62 | 0.51 | 0.58 | 0.49
$v_{n}\\{k\\}$ | 0.58 | 0.51 | 0.54 | 0.49 | 1.00
$d^{2}N_{\pi^{\pm}}/dy\,dp_{T}$ | 1.19 | 1.07 | 0.86 | 0.92 | 0.20
$d^{2}N_{K^{\pm}}/dy\,dp_{T}$ | 1.41 | 1.27 | 0.79 | 0.73 | 0.20
$d^{2}N_{p^{\pm}}/dy\,dp_{T}$ | 1.35 | 1.21 | 0.73 | 0.67 | 0.25
$v_{2}^{\pi^{\pm}}(p_{T})$ | 0.81 | 0.74 | 0.46 | 0.44 | 0.19
$v_{2}^{K^{\pm}}(p_{T})$ | 0.92 | 0.89 | 0.55 | 0.55 | 0.19
$v_{2}^{p^{\pm}}(p_{T})$ | 0.49 | 0.47 | 0.34 | 0.35 | 0.25
$v_{3}^{\pi^{\pm}}(p_{T})$ | 0.65 | 0.57 | 0.69 | 0.57 | 0.24
average | 0.89 | 0.83 | 0.69 | 0.66 |
$\sigma_{AA}$ | 1.13 | 3.80 | 1.53 | 3.40 | 1.00
Table 2: Average number of standard deviations from experimental data for
different classes of observables for the four fits presented in Fig. 2 and the
average weight $\bar{\omega}$ per observable class when used in the weighted
analysis. Uncertainties include experimental uncertainty and theoretical
uncertainty from the emulation (the latter is dominant for the $v_{n}$
classes). Including the cross section $\sigma_{AA}{}$ in the fit strongly
improves the agreement with $\sigma_{AA}{}$ but leads to only a mild worsening
for the other observables.

An excellent test of our model is provided by the triple-differential
observable $\rho(v_{n}\\{2\\}^{2},\langle p_{T}\rangle)$, which measures the
correlation between anisotropic flow and the mean $p_{T}$ Bozek (2016); Bozek
and Mehrabpour (2020); Giacalone  _et al._ (2021). This observable is
statistically expensive to compute (we simulate 625k hydro events for 20
parameter settings from the posterior) and hence cannot be included in the
Bayesian fit. It moreover sensitively depends on the precise experimental
procedure, including cuts on pseudorapidity, cuts on transverse momentum and
the method to select centrality bins ATLAS (2022). Nevertheless, the
observable is conjectured to be sensitive only to the hydrodynamic initial
conditions, and in particular the nucleon width Giacalone  _et al._ (2022).

Fig. 4 presents $\rho(v_{n}\\{2\\}^{2},\langle p_{T}\rangle)$ as compared with
ATLAS data ATLAS (2022). Due to the expensive nature of this analysis we only
include systematic uncertainty from the posterior for the weighted case
including $\sigma_{AA}{}$ and show Maximum a Posteriori (MAP) results for the
unweighted with $\sigma_{AA}{}$ and weighed without $\sigma_{AA}{}$ cases.
Clearly including $\sigma_{AA}{}$ dramatically improves the description of
this observable, which can at least in part be attributed to the smaller
nucleon width Giacalone  _et al._ (2022). A comparison with ALICE data Acharya
_et al._ (2021) is included in the Supplemental Material. Curiously the
systematic uncertainty for XeXe collisions is significantly larger than for
PbPb collisions (see also Supplemental Material).

Discussion -  A question deserving further study is in what way our weighting
procedure realistically captures the theoretical and experimental uncertainty.
Indeed, systematic off-sets on average observables are within one standard
deviation from the experimental results and a naive $\chi^{2}$ would say that
the uncertainties are accurate (see also Table 2 and the Supplemental
Material). This, however, ignores the fact that the 653 data points are highly
correlated, which is difficult to fully take into account. Also, while
observables are on average one standard deviation away from the experimental
results, the deviations are not Gaussian. Instead most observables are well
within one standard deviation, while a small number deviates significantly.
This is a further indication that a naive $\chi^{2}$ should not be trusted. We
hence believe the weighting leads to more physically realistic results, in
particular giving a physically realistic $\sigma_{AA}$, nucleon width and bulk
viscosity.

On a superficial level the results in this Letter show that the newly measured
$\sigma_{AA}$ improves estimates on the nucleon width, and subsequently the
transport coefficients. We wish to caution here, however, that this fact
should make us rethink the growing popularity of Bayesian analyses. Indeed,
prior to this Letter all such analyses ruled out nucleon widths smaller than
about 0.8 fm, which in light of this new analysis was not warranted. The
question we should ask is if the data really convincingly implied such large
nucleon widths? Small changes in many data points (often as small as 0.1
standard deviation) quickly add up to large Bayes factors that can
artificially constrain posterior distributions.

In principle, given accurate estimates of the correlated uncertainties of the
data and theoretical model the Bayes posterior is accurate. Note, however,
that even underestimating uncertainties by 10% is enough to cause an otherwise
good fit to now be off by about 0.1 standard deviation, which, if this occurs
for many data points, will add up, and the Bayesian analysis will try to
compensate for this elsewhere. In a sense this means that the robustness of
the Bayesian posterior depends sensitively on the accuracy of the uncertainty.
We stress that the hardest part is to accurately estimate the theoretical
model uncertainty. Currently this only includes emulation uncertainty, which
can be sizeable but is not a physical uncertainty. Instead, many modelling
choices such as our viscous particlization scheme (see Bernhard (2018)) or the
particulars of the hadronic afterburner have uncertainties that are not
included in the analysis. In fact, the full modeling uncertainty is difficult
to quantify, though we note that the versatility of our 21-dimensional model
attempts to include a wide scope of theoretical uncertainties.

It is perhaps curious that an early Bayesian study found a smaller width of
$0.48\pm 0.1$ Bernhard  _et al._ (2016). As noted in Bernhard (2018), however,
this is mostly due to having a simplistic initial stage where no radial flow
is created during the first $0.4\,$fm$/c$ and due to ignoring bulk viscous
corrections during particlization. This earlier value is hence mostly a
feature of the model being simplified, albeit in hindsight it got a reasonable
nucleon width estimate (but for the wrong reason).

Another relevant comment not included in this study can be that the nucleon
width could depend on their position within the nucleus. It is not so
unreasonable that nucleons in the skin of the nucleus could be smaller than
nucleons in the center. The size of nucleons in the center would not
significantly affect the nucleus-nucleus cross section.

A smaller nucleon width is consistent with other analyses, in particular
including the IP-Glasma model Schenke  _et al._ (2012); Mäntysaari and Schenke
(2016a, b). Indeed, the IP-Glasma model is based on HERA data, which extracts
a transverse two-gluon radius from $J/\psi$ scattering of $0.50\pm 0.03\,$fm
Caldwell and Kowalski (2010). Ref. Giacalone  _et al._ (2022) argues that the
characteristic nucleon length scale should be of order $0.5\,$fm or smaller,
although this can refer to the nucleon width or the constituent width. Since
$\sigma_{AA}{}$ is uniquely sensitive to the nucleon width the combination of
these measurements leads to a unique characterization of the nucleon spatial
profile.

Acknowledgements -  We thank Somadutta Bhatta, Andrea Dainese, Giuliano
Giacalone, Constantin Loizides, Aleksas Mazeliauskas, Rosi Reed, Mike Sas, Urs
Wiedemann and You Zhou for interesting discussions.

## References

  * Heinz and Snellings (2013) U. Heinz and R. Snellings, Ann. Rev. Nucl. Part. Sci. 63, 123 (2013), arXiv:1301.2826 [nucl-th] .
  * Busza  _et al._ (2018) W. Busza, K. Rajagopal, and W. van der Schee, Ann. Rev. Nucl. Part. Sci. 68, 339 (2018), arXiv:1802.04801 [hep-ph] .
  * Berges  _et al._ (2021) J. Berges, M. P. Heller, A. Mazeliauskas, and R. Venugopalan, Rev. Mod. Phys. 93, 035003 (2021), arXiv:2005.12299 [hep-th] .
  * Romatschke and Romatschke (2019) P. Romatschke and U. Romatschke,  _Relativistic Fluid Dynamics In and Out of Equilibrium_, Cambridge Monographs on Mathematical Physics (Cambridge University Press, 2019) arXiv:1712.05815 [nucl-th] .
  * Moreland  _et al._ (2015) J. S. Moreland, J. E. Bernhard, and S. A. Bass, Phys. Rev. C 92, 011901 (2015), arXiv:1412.4708 [nucl-th] .
  * Moreland  _et al._ (2020) J. S. Moreland, J. E. Bernhard, and S. A. Bass, Phys. Rev. C 101, 024911 (2020), arXiv:1808.02106 [nucl-th] .
  * Bernhard  _et al._ (2019) J. E. Bernhard, J. S. Moreland, and S. A. Bass, Nature Phys. 15, 1113 (2019).
  * Nijs  _et al._ (2021a) G. Nijs, W. van der Schee, U. Gürsoy, and R. Snellings, Phys. Rev. Lett. 126, 202301 (2021a), arXiv:2010.15130 [nucl-th] .
  * Everett  _et al._ (2021) D. Everett _et al._ (JETSCAPE), Phys. Rev. C 103, 054904 (2021), arXiv:2011.01430 [hep-ph] .
  * Nijs and van der Schee (2021) G. Nijs and W. van der Schee, (2021), arXiv:2110.13153 [nucl-th] .
  * Pohl  _et al._ (2010) R. Pohl _et al._ , Nature 466, 213 (2010).
  * ALICE (2022) ALICE, (2022), arXiv:2204.10148 [nucl-ex] .
  * Alver  _et al._ (2008) B. Alver, M. Baker, C. Loizides, and P. Steinberg, (2008), arXiv:0805.4411 [nucl-ex] .
  * Loizides  _et al._ (2015) C. Loizides, J. Nagle, and P. Steinberg, SoftwareX 1-2, 13 (2015), arXiv:1408.2549 [nucl-ex] .
  * Loizides  _et al._ (2018) C. Loizides, J. Kamin, and D. d’Enterria, Phys. Rev. C 97, 054910 (2018), [Erratum: Phys.Rev.C 99, 019901 (2019)], arXiv:1710.07098 [nucl-ex] .
  * d’Enterria and Loizides (2021) D. d’Enterria and C. Loizides, Ann. Rev. Nucl. Part. Sci. 71, 315 (2021), arXiv:2011.14909 [hep-ph] .
  * Klos  _et al._ (2007) B. Klos _et al._ , Phys. Rev. C 76, 014311 (2007), arXiv:nucl-ex/0702016 .
  * Tarbert  _et al._ (2014) C. M. Tarbert _et al._ , Phys. Rev. Lett. 112, 242502 (2014), arXiv:1311.0168 [nucl-ex] .
  * Skands  _et al._ (2014) P. Skands, S. Carrazza, and J. Rojo, Eur. Phys. J. C 74, 3024 (2014), arXiv:1404.5630 [hep-ph] .
  * Aamodt  _et al._ (2011) K. Aamodt _et al._ (ALICE), Phys. Rev. Lett. 106, 032301 (2011), arXiv:1012.1657 [nucl-ex] .
  * Adam  _et al._ (2016a) J. Adam _et al._ (ALICE), Phys. Rev. Lett. 116, 222302 (2016a), arXiv:1512.06104 [nucl-ex] .
  * Adam  _et al._ (2016b) J. Adam _et al._ (ALICE), Phys. Rev. C 94, 034903 (2016b), arXiv:1603.04775 [nucl-ex] .
  * Acharya  _et al._ (2019) S. Acharya _et al._ (ALICE), Phys. Lett. B 793, 420 (2019), arXiv:1805.05212 [nucl-ex] .
  * Acharya  _et al._ (2020) S. Acharya _et al._ (ALICE), Phys. Rev. C 101, 044907 (2020), arXiv:1910.07678 [nucl-ex] .
  * Abelev  _et al._ (2013) B. Abelev _et al._ (ALICE), Phys. Rev. C 88, 044910 (2013), arXiv:1303.0737 [hep-ex] .
  * Abelev  _et al._ (2014a) B. B. Abelev _et al._ (ALICE), Eur. Phys. J. C 74, 3077 (2014a), arXiv:1407.5530 [nucl-ex] .
  * Acharya  _et al._ (2018a) S. Acharya _et al._ (ALICE), JHEP 07, 103 (2018a), arXiv:1804.02944 [nucl-ex] .
  * Adam  _et al._ (2016c) J. Adam _et al._ (ALICE), JHEP 09, 164 (2016c), arXiv:1606.06057 [nucl-ex] .
  * Acharya  _et al._ (2018b) S. Acharya _et al._ (ALICE), JHEP 09, 006 (2018b), arXiv:1805.04390 [nucl-ex] .
  * Khachatryan  _et al._ (2016) V. Khachatryan _et al._ (CMS), Phys. Lett. B 759, 641 (2016), arXiv:1509.03893 [hep-ex] .
  * (31) https://sites.google.com/view/govertnijs/trajectum.
  * Nijs  _et al._ (2021b) G. Nijs, W. van der Schee, U. Gürsoy, and R. Snellings, Phys. Rev. C 103, 054909 (2021b), arXiv:2010.15134 [nucl-th] .
  * (33) G. Nijs and W. van der Schee, arXiv:220x.xxxxx [nucl-th] .
  * Abelev  _et al._ (2014b) B. B. Abelev _et al._ (ALICE), JINST 9, P11003 (2014b), arXiv:1405.1849 [nucl-ex] .
  * ATLAS (2022) ATLAS, (2022), arXiv:2205.00039 [nucl-ex] .
  * Bernhard (2018) J. E. Bernhard,  _Bayesian parameter estimation for relativistic heavy-ion collisions_ , Ph.D. thesis, Duke U. (2018), arXiv:1804.06469 [nucl-th] .
  * Bozek (2016) P. Bozek, Phys. Rev. C 93, 044908 (2016), arXiv:1601.04513 [nucl-th] .
  * Bozek and Mehrabpour (2020) P. Bozek and H. Mehrabpour, Phys. Rev. C 101, 064902 (2020), arXiv:2002.08832 [nucl-th] .
  * Giacalone  _et al._ (2021) G. Giacalone, F. G. Gardim, J. Noronha-Hostler, and J.-Y. Ollitrault, Phys. Rev. C 103, 024909 (2021), arXiv:2004.01765 [nucl-th] .
  * Giacalone  _et al._ (2022) G. Giacalone, B. Schenke, and C. Shen, Phys. Rev. Lett. 128, 042301 (2022), arXiv:2111.02908 [nucl-th] .
  * Acharya  _et al._ (2021) S. Acharya _et al._ (ALICE), (2021), arXiv:2111.06106 [nucl-ex] .
  * Bernhard  _et al._ (2016) J. E. Bernhard, J. S. Moreland, S. A. Bass, J. Liu, and U. Heinz, Phys. Rev. C 94, 024907 (2016), arXiv:1605.03954 [nucl-th] .
  * Schenke  _et al._ (2012) B. Schenke, P. Tribedy, and R. Venugopalan, Phys. Rev. C 86, 034908 (2012), arXiv:1206.6805 [hep-ph] .
  * Mäntysaari and Schenke (2016a) H. Mäntysaari and B. Schenke, Phys. Rev. Lett. 117, 052301 (2016a), arXiv:1603.04349 [hep-ph] .
  * Mäntysaari and Schenke (2016b) H. Mäntysaari and B. Schenke, Phys. Rev. D 94, 034042 (2016b), arXiv:1607.01711 [hep-ph] .
  * Caldwell and Kowalski (2010) A. Caldwell and H. Kowalski, Phys. Rev. C 81, 025203 (2010).
  * Aaboud  _et al._ (2020) M. Aaboud _et al._ (ATLAS), JHEP 01, 051 (2020), arXiv:1904.04808 [nucl-ex] .

Figure 5: Similar to Fig. 4 we show the comparison with available ALICE data
for PbPb (top) and XeXe (bottom) Acharya  _et al._ (2021). It is important
that ALICE uses different $p_{T}$ cuts and selects centrality bins based on
forward multiplicity instead of forward transverse energy. After taking this
into account the weighted and unweighted analyses including $\sigma_{AA}{}$
show satisfactory agreement. Figure 6: Similar to Fig. 4, we show the
comparison to ATLAS data for XeXe ATLAS (2022). We note that the
pseudorapidity cuts used in the theory prediction are different from those
used in ATLAS (2022) (see also main text). It is interesting that the
systematic uncertainty is significantly larger than for the PbPb case. This is
not due to statistical uncertainty (the gray statistical uncertainties are
much smaller), but presumably due to the deformed nature of the Xe nucleus.

##  I Supplemental Material

The $\rho(v_{n}\\{2\\}^{2},\langle p_{T}\rangle)$ observables depend
sensitively on details such as kinematic cuts on the particles used to compute
$\rho$ itself, and in addition to this also on which particles are used for
the determination of centrality (offline tracks), and whether these offline
tracks are weighted by their transverse energy $E_{T}$ or not. To illustrate
this, in Fig. 5 (top) we show the results from our analysis for PbPb compared
to ALICE data Acharya  _et al._ (2021), which uses different kinematic cuts
compared to ATLAS ATLAS (2022). ALICE also uses the forward V0 detectors which
measure particle multiplicity, whereas ATLAS uses the forward FCal that
measures the total transverse energy $E_{T}$. One can see that the results
from both of our fits which include the cross section are in good agreement
with the ALICE data, while they are in slight tension with the ATLAS data,
especially in the semi-central region (10–30%). Here we note that the ALICE
and ATLAS data points are largely consistent despite their different
experimental cuts (both give $\rho(v_{n}\\{2\\}^{2},\langle p_{T}\rangle)$ in
the range 0.19–0.23 in this centrality region), whereas for _Trajectum_ we
find a significantly larger $\rho(v_{n}\\{2\\}^{2},\langle p_{T}\rangle)$ for
the ATLAS experimental cuts.

Indeed, in Fig. 4 we find for $\rho(v_{2}\\{2\\}^{2},\langle p_{T}\rangle)$ a
value in the 10–20% centrality class of 0.263, whereas for the ALICE
comparison we compute only 0.197. In _Trajectum_ this difference can be
explained by the lower $p_{T}$ cut (lowering $\rho_{2}$ by 14%) and by the
$N_{\rm ch}$ versus $E_{T}$ centrality selection (lowering $\rho_{2}$ by
another 21%). We note that ATLAS also compared an $N_{\rm ch}$ versus $E_{T}$
comparison within its own dataset (Fig. 5 in ATLAS (2022)), but here it is
essential that $N_{\rm trk}$ is at mid-rapidity and furthermore provides an
inferior centrality selector than selecting on $E_{T}$ Aaboud  _et al._
(2020). This complicates the comparison, as within Trajectum we find that a
mid-rapidity centrality selection increases $\rho_{2}$ by 13%, which competes
with the decrease due to $N_{\rm ch}$ versus $E_{T}$ selection such that in
the end $\rho_{2}$ is about 7% lower for the mid-rapidity centrality selector
(somewhat in tension with Fig. 5 of ATLAS (2022), but see also Aaboud  _et
al._ (2020)). Finally, on the other hand, the ATLAS measurement depends
sensitively on the $\eta$ cut (Fig. 4 in ATLAS (2022)), by approximately the
tension we find in Fig. 4, whereas in Trajectum we do not find any $\eta$-cut
dependence in $\rho_{2}$ (we however present the calculation in Fig. 4 for
$|\eta|<1$, since this is much cheaper to compute when simulating all 20
parameter settings).

In Fig. 5 (bottom) and 6, we also show $\rho(v_{n}\\{2\\}^{2},\langle
p_{T}\rangle)$ for XeXe collisions at $\sqrt{s_{\rm NN}}=5.44\,$TeV. The
normalization of the profile at $5.44\,$TeV is determined by assuming a power-
law extrapoliation from the fits at 2.76 and $5.02\,$TeV. We compare results
with respectively ALICE Acharya  _et al._ (2021) and ATLAS data ATLAS (2022).
The comparison shows good agreement, even though in Trajectum for numerical
convenience we use a smaller rapidity cut than ATLAS uses. We verified,
however, for a single parameter setting that our boost invariant model does
not depend on the rapidity cut. The ATLAS results depend significantly on the
rapidity cut and since our MCMC fit is done using mid-rapidity experimental
data from ALICE it is preferable to compare with mid-rapidity ATLAS data. For
XeXe collisions these are however unfortunately not available and hence the
comparison as shown. In the XeXe simulations we use a deformed Woods-Saxon
distribution, where the radius $R$ is multiplied by
$1+\beta_{2}Y_{2}^{0}(\theta)+\beta_{3}Y_{3}^{0}(\theta)+\beta_{4}Y_{4}^{0}(\theta)$,
with $Y_{n}^{0}$ spherical harmonics, $\beta_{2}=0.162$, $\beta_{3}=0$ and
$\beta_{4}=-0.003$. Comparing to the PbPb results, the systematic errors are
considerably larger for XeXe, indicating that XeXe is more sensitive to the
remaining uncertainty in the posterior parameters.

In Fig. 7 and 8 we show the experimental data points used together with the
posterior distributions from the weighted fit including the new
$\sigma_{AA}{}$ (odd rows) and without the $\sigma_{AA}{}$ mesaurement (even
rows). In order of appearance (though 2.76 TeV always appears before 5.02 TeV
collisions) we include charged particle multiplicity $dN_{\text{ch}}/d\eta$ at
2.76 Aamodt  _et al._ (2011) and 5.02 TeV Adam  _et al._ (2016a) and
transverse energy $dE_{T}/d\eta$ at 2.76 TeV Adam  _et al._ (2016b). We also
include identified yields $dN/dy$ and mean $p_{T}$ for pions, kaons and
protons at 2.76 TeV Abelev  _et al._ (2013) and 5.02 TeV Acharya  _et al._
(2019, 2020). We then include identified transverse momentum spectra for
pions, kaons and protons at 2.76 TeV Abelev  _et al._ (2013) and 5.02 TeV
Acharya  _et al._ (2020) in seven coarse grained $p_{T}$-bins with bin
boundaries at $(0.25,0.5,0.75,1.0,1.4,1.8,2.2,3.0)\,\text{GeV}$. Next we show
$p_{T}$ fluctuations, which are only available at 2.76 TeV Abelev  _et al._
(2014a). Finally we show $p_{T}$ integrated charged particle anisotropic flow
at 2.76 and 5.02 TeV Acharya  _et al._ (2018a), as well as the
$p_{T}$-differential identified anisotropic flow coefficient
$v_{2}\\{2\\}(p_{T})$ for pions, kaons and protons, and $v_{3}\\{2\\}(p_{T})$
for pions at 2.76 TeV Adam  _et al._ (2016c) and 5.02 TeV Acharya  _et al._
(2018b).

Figure 7: We show all our 653 data points including a band that is taken from
the posterior distribution of our parameters. Odd rows show the results
including weights $\sigma_{AA}{}$ (hence having small nucleon width), whereas
even rows show the results without $\sigma_{AA}{}$ (also weighted). Often it
is hard to notice the difference, but we stress that the odd rows fit
$\sigma_{AA}{}$ much better at the cost of a slightly worse fit for the rest
of the data points as presented here (continues on next page).  Figure 8:
Continued from previous page.