Consider allowing JavaScript. Otherwise, you have to be proficient in reading since formulas will not be rendered. Furthermore, the table of contents in the left column for navigation will not be available and codefolding not supported. Sorry for the inconvenience.
Examples in this article were generated with
4.4.1
by the package PowerTOST
.^{1}
More examples are given in the respective vignette.^{2} See also the README on GitHub for an overview and the online manual^{3} for details.
The righthand badges give the respective section’s ‘level’.
Abbreviations are given at the end.
What is Average Bioequivalence with Expanding Limits?
For details about inferential statistics and hypotheses in equivalence see another article.
Definitions:
The concept of Scaled Average Bioequivalence (SABE) for HVDPs is based on the following considerations:
The conventional confidence interval inclusion approach in ABE \[\begin{matrix}\tag{1} \theta_1=1\Delta,\theta_2=\left(1\Delta\right)^{1}\\ H_0:\;\frac{\mu_\text{T}}{\mu_\text{R}}\ni\left\{\theta_1,\theta_2\right\}\;vs\;H_1:\;\theta_1<\frac{\mu_\text{T}}{\mu_\text{R}}<\theta_2\textsf{,}, \end{matrix}\] where \(\small{H_0}\) is the null hypothesis of inequivalence and \(\small{H_1}\) the alternative hypothesis, \(\small{\theta_1}\) and \(\small{\theta_2}\) are the fixed lower and upper limits of the acceptance range, and \(\small{\mu_\text{T}}\) are the geometric least squares means of \(\small{\text{T}}\) and \(\small{\text{R}}\), respectively is in Scaled Average Bioequivalence (SABE)^{6} modified to \[H_0:\;\frac{\mu_\text{T}}{\mu_\text{R}}\Big{/}\sigma_\text{wR}\ni\left\{\theta_{\text{s}_1},\theta_{\text{s}_2}\right\}\;vs\;H_1:\;\theta_{\text{s}_1}<\frac{\mu_\text{T}}{\mu_\text{R}}\Big{/}\sigma_\text{wR}<\theta_{\text{s}_2}\textsf{,}\tag{2}\] where \(\small{\sigma_\text{wR}}\) is the standard deviation of the reference and the scaled limits \(\small{\left\{\theta_{\text{s}_1},\theta_{\text{s}_2}\right\}}\) of the acceptance range depend on conditions given by the agency.
Average Bioequivalence with Expanding Limits (ABEL)^{7} for HVDPs is acceptable in numerous jurisdictions (Fig. 2), whereas directly widening the limits was recommended 2016 – 2022 in the member states of the Gulf Cooperation Council^{8} (Fig. 3). Alas, we are far away from global harmonization. Work on the ICH’s M13C guideline^{9} has not started with July 2024 (expected after finalization of M13B).
Expanding the limits is acceptable for different pharmacokinetic metrics.
PK metric  Jurisdiction 

C_{max}  The EEA,^{10} the WHO,^{11} ^{12} Australia,^{13} Belarus,^{14} the EAC,^{15} the ASEAN states,^{16} the Russian Federation,^{17} the EEU,^{18} Egypt,^{19} New Zealand,^{20} Chile,^{21} Brazil,^{22} members of the GCC,^{23} Mexico.^{24} 
AUC  WHO (if in a 4period full replicate design),^{11} Canada.^{25} 
C_{min}, C_{τ}  The EEA (controlled release products in steady state).^{26} 
_{partial} AUC  The EEA (controlled release products).^{26} 
Based on the switching \(\small{CV_0=30\%}\) we get the switching standard deviation \(\small{s_0=\sqrt{\log_{e}(CV_{0}^{2}+1)}\approx0.2935604\ldots}\), the (rounded) regulatory constant \(\small{k=\frac{\log_{e}1.25}{s_0}\sim0.76}\), and finally the expanded limits \(\small{\left\{\theta_{\text{s}_1},\theta_{\text{s}_2}\right\}=}\) \(\small{100\left(\exp(\mp0.76\cdot s_{\text{wR}})\right)}\). Note that in the guidelines the terminology \(\small{\left\{L,U\right\}}\) is used instead.
In order to apply the methods following conditions have to be fulfilled:
The wording in the guidelines gives the false impression that the methods are straightforward. In fact they are decision schemes (frameworks), which hinge on the estimated standard deviation of the reference treatment \(\small{s_{\text{wR}}}\).
If \(\small{CV_\text{wR}\leq30\%}\) the study has to be assessed for ABE (left branch) or for ABEL (right branch) otherwise.
In the ABELbranch there is an ‘upper cap’ of scaling (\(\small{uc=50\%}\) in all jurisdictions except for Health Canada, where \(\small{uc\approx 57.382\%}\)). Furthermore, the point estimate (\(\small{PE}\)) has to lie within 80.00 – 125.00%.
From March 2016 to August 2022 the GCC recommended a simplified variant with fixed widened limits of 75.00 – 133.33% if \(\small{CV_\text{wR}>30\%}\).^{8} Like in ABEL the \(\small{PE}\) had to lie within 80.00 – 125.00%. As in the FDA’s RSABE there was no ‘upper cap’.
Since the applicability of these approaches hinges on the realized values (\(\small{CV_\text{wR}}\), \(\small{PE}\)) in the particular study – which are random variables and as such unknown beforehand – analytical solutions for power (and hence, the sample size) do not exist. Therefore, extensive simulations of potential combinations have to be employed.
Whereas ABE is bijective (if \(\small{T}\) is equivalent to \(\small{R}\), \(\small{R}\) is also equivalent to \(\small{T}\)), this holds in ABEL only if \(\small{CV_\text{wT}=}\) \(\small{CV_\text{wR}}\). Therefore, switching from \(\small{R}\) to \(\small{T}\) is tolerable if \(\small{CV_\text{wT}<CV_\text{wR}}\) but switching from \(\small{T}\) to \(\small{R}\) can be problematic. Individual Bioequivalence (IBE) would not only compare the averages of products but also their variances and thus, support ‘Switchability’.^{27} However, IBE was never implemented in regulatory practice. For details see this article.
Cave: Under certain conditions the methods may lead to an inflated Type I Error (increased patient’s risk).^{28} ^{29} It is elaborated in another article.
For the evaluation I recommend the package
replicateBE
.^{30}
Interlude 1
Where do these numbers come from?
With the regulatory constant \(\small{k=0.76}\) one gets in the \(\small{\log_{e}}\)scale a linear expansion
from \(\small{CV_\text{wR}=30\%}\) to
\(\small{uc}\) based on \(\small{s_\text{wR}=\sqrt{\log_{e}(CV_\text{wR}^{2})+1}}\)
and \(\small{100\left(\exp(\mp k\cdot
s_\text{wR})\right)}\).
Hence, with \(\small{CV_\text{wR}=30\%}\) one gets
\(\small{\;\;\;\left\{\theta_1,\theta_2\right\}=100\left(\exp(\mp0.76\cdot
0.2935604)\right)=\left\{80.00,\,125.00\right\}}\).
With \(\small{uc=50\%}\) one gets
\(\small{\;\;\;\left\{\theta_{\text{s}_1},\,\theta_{\text{s}_2}\right\}=100\left(\exp(\mp0.76\cdot
0.4723807)\right)=\left\{69.83678,\,143.19102\right\}}\).
For Health Canada with \(\small{uc=57.382\%}\) one gets \(\small{\;\;\;\left\{\theta_{\text{s}_1},\,\theta_{\text{s}_2}\right\}=100\left(\exp(\mp0.76\cdot
0.5335068)\right)=\left\{66.66667,\,150.00000\right\}}\).
As stated above, the regulatory constant \(\small{k}\) is given in all guidelines rounded to two significant digits. Based on the switching \(\small{CV_\text{w0}=30\%}\) one gets \(\small{s_\text{w0}=\sqrt{\log_{e}(CV_\text{w0}^2+1)}\approx 0.2935604}\) and, therefore, the exact value would be \(\small{k=\frac{\log_{e}1.25}{s_\text{w0}}\approx 0.7601283}\). For a rant see Karalis et al.^{31} The wisdom of regulators is unfathomable.
A clinically not relevant \(\small{\Delta}\) of 30% (leading a fixed range of 70.00 – 142.86%) was acceptable in Europe (in exceptional cases even for AUC) if prespecified in the protocol. A replicate design was not required.^{32} ^{33} ^{34}A basic knowledge of R is
required. To run the scripts at least version 1.4.8 (20190829) of
PowerTOST
is required and at least version 1.5.3
(20210118) suggested. Any version of R would
likely do, though the current release of PowerTOST
was only
tested with R 4.3.3 (20240229) and later.
All scripts were run on a Xeon E31245v3 @ 3.40GHz (1/4 cores) 16GB RAM
with R 4.4.1 on Windows 7 build 7601, Service
Pack 1.
Note that in all functions of PowerTOST
the arguments
(say, the assumed T/Rratio theta0
, the assumed coefficient
of variation CV
, etc.) have to be given as ratios and not
in percent.
sampleN.scABEL()
gives balanced sequences
(i.e., an equal number of subjects is allocated to all
sequences). Furthermore, the estimated sample size is the total
number of subjects – and not the number of subjects per sequence like in
some other software packages.
All examples deal with studies where the response variables likely follow a lognormal distribution, i.e., we assume a multiplicative model (ratios instead of differences). We work with \(\small{\log_{e}}\)transformed data in order to allow analysis by the ttest (requiring differences).
It may sound picky but ‘sample size calculation’ (as used in most guidelines and alas, in some publications and textbooks) is sloppy terminology. In order to get prospective power (and hence, a sample size), we need five values:
where
In other words, obtaining a sample size is not an exact calculation like \(\small{2\times2=4}\) but always just an estimation.
“Power Calculation – A guess masquerading as mathematics.
Of note, it is extremely unlikely that all assumptions will be exactly realized in a particular study. Hence, calculating retrospective (a.k.a. estimated, post hoc, a posteriori) power is not only futile but plain nonsense.^{42} For details see another article.
It is a prerequisite that no carryover from one period to the next exists. Only then the comparison of treatments will be unbiased.^{43} Carryover is elaborated in another article.
Except for Health Canada (where a mixedeffects model is required) the recommended evaluation by an ANOVA assumes homoscedasticity (\(\small{s_\text{wT}^2=s_\text{wR}^2}\)), which is – more often than not – wrong.
The sample size cannot be directly estimated, in
ABEL
only
power simulated for an already given sample size.
In contrast to ABE, where power can be exactly calculated and the sample size is determined in an iterative procedure, ABEL is a decision scheme and therefore, simulations must to performed for different combinations of assumed T/Rratios and CVs.
“Power. That which statisticians are always calculating but never have.
Let’s start with PowerTOST
.
library(PowerTOST) # attach it to run the examples
The sample size functions’ defaults are:
Argument  Default  Meaning 

alpha

0.05

Nominal level of the test 
targetpower

0.80

Target (desired) power 
theta0

0.90

Assumed T/Rratio 
theta1

0.80

Lower BElimit in ABE and lower PEconstraint in ABEL 
theta2

1.25

Upper BElimit in ABE and upper PEconstraint in ABEL 
design

"2x3x3"

Treatments × Sequences × Periods 
regulator

"EMA"

Guess… 
nsims

1e05

Number of simulations 
print

TRUE

Output to the console 
details

TRUE

Show regulatory settings and sample size search 
setseed

TRUE

Set a fixed seed (recommended for reproducibility) 
Implemented designs are "2x3x3"
(the three sequence,
three period – partial – replicate design TRRRTRRRT),
"2x2x4"
(two sequence, four period full replicate
designs TRTRRTRT, TRRTRTTR, and TTRRRRTT), and "2x2x3"
(two sequence, three period full replicate designs TRTRTR and
TRRRTT).
For a quick overview of the regulatory limits use the function
scABEL()
– for once in percent according to the guidelines.
Note that Health Canada uses rounding to the first decimal place.
< data.frame(regulator = "EMA", CV = c(30, 50),
df1 L = NA, U = NA,
cap = c("lower", "upper"))
< data.frame(regulator = "HC",
df2 CV = c(30, 57.382),
L = NA, U = NA,
cap = c("lower", "upper"))
for (i in 1:2) {
3:4] < sprintf("%.2f%%", 100*scABEL(df1$CV[i]/100))
df1[i, 3:4] < sprintf("%.1f %%", 100*scABEL(df2$CV[i]/100,
df2[i, regulator = "HC"))
}$CV < sprintf("%.3f%%", df1$CV)
df1$CV < sprintf("%.3f%%", df2$CV)
df2if (packageVersion("PowerTOST") >= "1.5.3") {
< data.frame(regulator = "GCC <2022",
df3 CV = c(30, 50), L = NA, U = NA,
cap = c("lower", "none"))
for (i in 1:2) {
3:4] < sprintf("%.2f%%", 100*scABEL(df3$CV[i]/100,
df3[i, regulator = "GCC"))
}$CV < sprintf("%.3f%%", df3$CV)
df3
}if (packageVersion("PowerTOST") >= "1.53") {
print(df1, row.names = F); print(df2, row.names = F); print(df3, row.names = F)
else {
} print(df1, row.names = F); print(df2, row.names = F)
}
# regulator CV L U cap
# EMA 30.000% 80.00% 125.00% lower
# EMA 50.000% 69.84% 143.19% upper
# regulator CV L U cap
# HC 30.000% 80.0 % 125.0 % lower
# HC 57.382% 66.7 % 150.0 % upper
# regulator CV L U cap
# GCC <2022 30.000% 80.00% 125.00% lower
# GCC <2022 50.000% 75.00% 133.33% none
The sample size functions of PowerTOST
use a
modification of Zhang’s method^{46} for the first guess.
# Note that theta0 = 0.90 is the default
sampleN.scABEL(CV = 0.45, targetpower = 0.80,
design = "2x2x4", details = TRUE)
#
# +++++++++++ scaled (widened) ABEL +++++++++++
# Sample size estimation
# (simulation based on ANOVA evaluation)
# 
# Study design: 2x2x4 (4 period full replicate)
# logtransformed data (multiplicative model)
# 1e+05 studies for each step simulated.
#
# alpha = 0.05, target power = 0.8
# CVw(T) = 0.45; CVw(R) = 0.45
# True ratio = 0.9
# ABE limits / PE constraint = 0.8 ... 1.25
# EMA regulatory settings
#  CVswitch = 0.3
#  cap on scABEL if CVw(R) > 0.5
#  regulatory constant = 0.76
#  pe constraint applied
#
#
# Sample size search
# n power
# 24 0.7539
# 26 0.7846
# 28 0.8112
An alternative to simulating the ‘key statistics’ is by subject data
simulations via the function sampleN.scABEL.sdsims()
.
However, it comes with a price, speed.
# method n power rel.speed
# ‘key statistics’ 28 0.81116 1.0
# subject simulations 28 0.81196 21.5
Throughout the examples I’m referring to studies in a single center – not multiple groups within them or multicenter studies. That’s another pot of tea.
We assume a CV of 0.45, a T/Rratio of 0.90, a target a power of 0.80, and want to perform the study in a two sequence four period full replicate study (TRTRRTRT or TRRTRTTR or TTRRRRTT) for the EMA’s ABEL.
Since theta0 = 0.90
,^{47} targetpower = 0.80
, and
regulator = "EMA"
are defaults of the function, we don’t
have to give them explicitely. As usual in bioequivalence,
alpha = 0.05
is employed (we will assess the study by a
\(\small{100\,(12\,\alpha)=90\%}\)
confidence interval). Hence, you need to specify only the
CV
(assuming \(\small{CV_\text{wT}=CV_\text{wR}}\)) and
design = "2x2x4"
.
To shorten the output, use the argument
details = FALSE
.
sampleN.scABEL(CV = 0.45, design = "2x2x4", details = FALSE)
#
# +++++++++++ scaled (widened) ABEL +++++++++++
# Sample size estimation
# (simulation based on ANOVA evaluation)
# 
# Study design: 2x2x4 (4 period full replicate)
# logtransformed data (multiplicative model)
# 1e+05 studies for each step simulated.
#
# alpha = 0.05, target power = 0.8
# CVw(T) = 0.45; CVw(R) = 0.45
# True ratio = 0.9
# ABE limits / PE constraint = 0.8 ... 1.25
# Regulatory settings: EMA
#
# Sample size
# n power
# 28 0.8112
Sometimes we are not interested in the entire output and want to use
only a part of the results in subsequent calculations. We can suppress
the output by stating the additional argument print = FALSE
and assign the result to a data frame (here x
).
< sampleN.scABEL(CV = 0.45, design = "2x2x4",
x details = FALSE, print = FALSE)
Let’s retrieve the column names of x
:
names(x)
# [1] "Design" "alpha" "CVwT"
# [4] "CVwR" "theta0" "theta1"
# [7] "theta2" "Sample size" "Achieved power"
# [10] "Target power" "nlast"
Now we can access the elements of x
by their names. Note
that double square brackets and double or single quotes
([["…"]]
, [['…']])
have to be used.
"Sample size"]]
x[[# [1] 28
'Achieved power']]
x[[# [1] 0.81116
Although you could access the elements by the number of the
column(s), I don’t recommend that, since in various other functions of
PowerTOST
these numbers are different and hence, difficult
to remember. If you insist in accessing elements by columnnumbers, use
single square brackets […]
.
8:9]
x[# Sample size Achieved power
# 1 28 0.81116
With 28 subjects (14 per sequence) we achieve the power we desire.
What happens if we have one dropout?
power.scABEL(CV = 0.45, design = "2x2x4",
n = x[["Sample size"]]  1)
# Unbalanced design. n(i)=14/13 assumed.
# [1] 0.79848
Below the 0.80 we desire.
Since dropouts are common, it makes sense to include / dose more subjects in order to end up with a number of eligible subjects which is not lower than our initial estimate.
Let us explore that in the next section.
We define two supportive functions:
n
will be rounded up to achieve balance.< function(n, n.seq) {
balance return(as.integer(n.seq * (n %/% n.seq + as.logical(n %% n.seq))))
}
n
and the anticipated droputrate
do.rate
.< function(n, do.rate, n.seq) {
nadj return(as.integer(balance(n / (1  do.rate), n.seq)))
}
In order to come up with a suggestion we have to anticipate a (realistic!) dropout rate. Note that this not the job of the statistician; ask the Principal Investigator.
“It is a capital mistake to theorise before one has data.
The dropoutrate is calculated from the eligible and
dosed subjects
or simply \[\begin{equation}\tag{3}
do.rate=1n_\text{eligible}/n_\text{dosed}
\end{equation}\] Of course, we know it only after the
study was performed.
By substituting \(n_\text{eligible}\) with the estimated sample size \(n\), providing an anticipated dropoutrate and rearrangement to find the adjusted number of dosed subjects \(n_\text{adj}\) we should use \[\begin{equation}\tag{4} n_\text{adj}=\;\upharpoonleft n\,/\,(1do.rate) \end{equation}\] where \(\upharpoonleft\) denotes rounding up to the complete sequence as implemented in the functions above.
An all too common mistake is to increase the estimated sample size \(n\) by the dropoutrate according to \[\begin{equation}\tag{5} n_\text{adj}=\;\upharpoonleft n\times(1+do.rate) \end{equation}\] If you used \(\small{(5)}\) in the past – you are not alone. In a small survey a whopping 29% of respondents reported to use it.^{49} Consider changing your routine.
“There are no routine statistical questions, only questionable statistical routines.
In the following I specified more arguments to make the function more
flexible.
Note that I wrapped the function power.scABEL()
in
suppressMessages()
. Otherwise, the function will throw for
any odd sample size a message telling us that the design is
unbalanced. Well, we know that.
< 0.45 # withinsubject CV
CV < 0.80 # target (desired) power
target < 0.90 # assumed T/Rratio
theta0 < "2x2x4"
design < 0.15 # anticipated dropoutrate 15%
do.rate # might be realively high due
# to the 4 periods
< as.integer(substr(design, 3, 3))
n.seq < scABEL(CV) # expanded limits
lims < sampleN.scABEL(CV = CV, theta0 = theta0,
df targetpower = target,
design = design,
details = FALSE,
print = FALSE)
# calculate the adjusted sample size
< nadj(df[["Sample size"]], do.rate, n.seq)
n.adj # (decreasing) vector of eligible subjects
< n.adj:df[["Sample size"]]
n.elig < paste0("Assumed CV : ",
info
CV,"\nAssumed T/R ratio : ",
theta0,"\nExpanded limits : ",
sprintf("%.4f\u2026%.4f",
1], lims[2]),
lims["\nPE constraints : ",
sprintf("%.4f\u2026%.4f",
0.80, 1.25), # fixed in ABEL
"\nTarget (desired) power : ",
target,"\nAnticipated dropoutrate: ",
do.rate,"\nEstimated sample size : ",
"Sample size"]], " (",
df[["Sample size"]]/n.seq, "/sequence)",
df[["\nAchieved power : ",
signif(df[["Achieved power"]], 4),
"\nAdjusted sample size : ",
" (", n.adj/n.seq, "/sequence)",
n.adj, "\n\n")
# explore the potential outcome for
# an increasing number of dropouts
< signif((n.adj  n.elig) / n.adj, 4)
do.act < data.frame(dosed = n.adj,
df eligible = n.elig,
dropouts = n.adj  n.elig,
do.act = do.act,
power = NA)
for (i in 1:nrow(df)) {
$power[i] < suppressMessages(
dfpower.scABEL(CV = CV,
theta0 = theta0,
design = design,
n = df$eligible[i]))
}cat(info)
print(round(df, 4), row.names = FALSE)
# Assumed CV : 0.45
# Assumed T/R ratio : 0.9
# Expanded limits : 0.7215…1.3859
# PE constraints : 0.8000…1.2500
# Target (desired) power : 0.8
# Anticipated dropoutrate: 0.15
# Estimated sample size : 28 (14/sequence)
# Achieved power : 0.8112
# Adjusted sample size : 34 (17/sequence)
#
# dosed eligible dropouts do.act power
# 34 34 0 0.0000 0.8720
# 34 33 1 0.0294 0.8630
# 34 32 2 0.0588 0.8553
# 34 31 3 0.0882 0.8456
# 34 30 4 0.1176 0.8340
# 34 29 5 0.1471 0.8237
# 34 28 6 0.1765 0.8112
In the worst case (six dropouts) we end up with the originally estimated sample size of 28. Power preserved, mission accomplished. If we have less dropouts, splendid – we gain power.
If we would have adjusted the sample size according to \(\small{(5)}\) we
would have dosed also 34 subjects.
Cave: This might not always be the case… If the anticipated dropout rate
of 15% is realized in the study, we would have also 28 eligible subjects
(power 0.8112). In this example we achieve still more than our target
power but the loss might be relevant in other cases.
As said in the preliminaries, calculating post hoc power is futile.
“There is simple intuition behind results like these: If my car made it to the top of the hill, then it is powerful enough to climb that hill; if it didn’t, then it obviously isn’t powerful enough. Retrospective power is an obvious answer to a rather uninteresting question. A more meaningful question is to ask whether the car is powerful enough to climb a particular hill never climbed before; or whether a different car can climb that new hill. Such questions are prospective, not retrospective.
However, sometimes we are interested in it for planning the next study.
If you give and odd total sample size n
,
power.scABEL()
will try to keep sequences as balanced as
possible and show in a message how that was done.
< 27
n.act signif(power.scABEL(CV = 0.45, n = n.act,
design = "2x2x4"), 6)
# Unbalanced design. n(i)=14/13 assumed.
# [1] 0.79848
Say, our study was more unbalanced. Let us assume that we dosed 34
subjects, the total number of subjects was also 27 but all dropouts
occured in one sequence (unlikely but possible).
Instead of the total sample size n
we can give the number
of subjects of each sequence as a vector (the order is generally^{51} not
relevant, i.e., it does not matter which element refers to
which sequence).
By setting details = TRUE
we can retrieve the components
of the simulations (probability to pass each test).
< "2x2x4"
design < 0.45
CV < 34
n.adj < 27
n.act < n.adj / 2
n.s1 < n.act  n.s1
n.s2 < 0.90
theta0 < suppressMessages(
post.hoc power.scABEL(CV = CV,
n = c(n.s1, n.s2),
theta0 = theta0,
design = design,
details = TRUE))
< power.TOST(CV = CV,
ABE.xact n = c(n.s1, n.s2),
theta0 = theta0,
design = design)
< nchar(as.character(n.adj))
sig.dig < paste0("%", sig.dig, ".0f (%",
fmt ".0f dropouts)")
sig.dig, cat(paste0("Dosed subjects: ", sprintf("%2.0f", n.adj),
"\nEligible : ",
sprintf(fmt, n.act, n.adj  n.act),
"\n Sequence 1 : ",
sprintf(fmt, n.s1, n.adj / 2  n.s1),
"\n Sequence 1 : ",
sprintf(fmt, n.s2, n.adj / 2  n.s2),
"\nPower overall : ",
sprintf("%.5f", post.hoc[1]),
"\n p(ABEL) : ",
sprintf("%.5f", post.hoc[2]),
"\n p(PE) : ",
sprintf("%.5f", post.hoc[3]),
"\n p(ABE) : ",
sprintf("%.5f", post.hoc[4]),
"\n p(ABE) exact: ",
sprintf("%.5f", ABE.xact), "\n"))
# Dosed subjects: 34
# Eligible : 27 ( 7 dropouts)
# Sequence 1 : 17 ( 0 dropouts)
# Sequence 1 : 10 ( 7 dropouts)
# Power overall : 0.77670
# p(ABEL) : 0.77671
# p(PE) : 0.91595
# p(ABE) : 0.37628
# p(ABE) exact: 0.37418
The components of overall power are:
p(ABEL)
is the probability that the confidence interval
is within the expanded / widened limits.p(PE)
is the probability that the point estimate is
within 80.00–125.00%.p(ABE)
is the probability of passing conventional
Average Bioequivalence.power.TOST()
– confirming the simulation’s result.Of course, in a particular study you will provide the numbers in the
n
vector directly.
The CV and the T/Rratio are only assumptions.
Whatever their origin might be (literature, previous studies), they
carry some degree of uncertainty. Hence, believing^{52} that they are the
true ones may be risky.
Some statisticians call that the ‘CarvedinStone’ approach.
Say, we performed a pilot study in 16 subjects and estimated the CV as 0.45.
The \(\small{\alpha}\) confidence interval of the CV is obtained via the \(\small{\chi^2}\)distribution of its error variance \(\small{\sigma^2}\) with \(\small{n2}\) degrees of freedom. \[\eqalign{\tag{6} s_\text{w}^2&=\log_{e}(CV_\text{w}^2+1)\\ L=\frac{(n1)\,s_\text{w}^2}{\chi_{\alpha/2,\,n2}^{2}}&\leq\sigma_\text{w}^2\leq\frac{(n1)\,s_\text{w}^2}{\chi_{1\alpha/2,\,n2}^{2}}=U\\ \left\{lower\;CL,\;upper\;CL\right\}&=\left\{\sqrt{\exp(L)1},\sqrt{\exp(U)1}\right\} }\]
Let’s calculate the 95% confidence interval of the CV to get an idea.
< 16 # pilot study
m < CVCL(CV = 0.45, df = m  2,
ci side = "2sided", alpha = 0.05)
< c(ci, mean(ci))
ci names(ci)[3] < "midpoint"
signif(ci, 4)
# lower CL upper CL midpoint
# 0.3223 0.7629 0.5426
Surprised? Although 0.45 is the best estimate for planning the next study, there is no guarantee that we will get exactly the same outcome. Since the \(\small{\chi^2}\)distribution is skewed to the right, it is more likely that we will face a higher CV than a lower one in the planned study.
If we plan the study based on 0.45, we would opt for 28 subjects like
in the examples before (not adjusted for the dropoutrate).
If the CV will be lower, we loose power (less expansion). But
what if it will be higher? Depends. Since we may expand more, we gain
power. However, if we cross the upper cap of scaling (50% for the EMA),
we will loose power. But how much?
Let’s explore what might happen at the confidence limits of the CV.
< 16
m < CVCL(CV = 0.45, df = m  2,
ci side = "2sided", alpha = 0.05)
< 28
n < data.frame(CV = c(ci[["lower CL"]],
comp 0.45, mean(ci),
"upper CL"]]),
ci[[power = NA)
for (i in 1:nrow(comp)) {
$power[i] < power.scABEL(CV = comp$CV[i],
compdesign = "2x2x4",
n = n)
}1] < signif(comp[, 1], 4)
comp[, 2] < signif(comp[, 2], 6)
comp[, print(comp, row.names = FALSE)
# CV power
# 0.3223 0.73551
# 0.4500 0.81116
# 0.5426 0.79868
# 0.7629 0.60158
Might hurt.
What can we do? The larger the previous study was, the larger the degrees of freedom and hence, the narrower the confidence interval of the CV. In simple terms: The estimate is more certain. On the other hand, it also means that very small pilot studies are practically useless. What happens when we plan the study based on the confidence interval of the CV?
< seq(12, 30, 6)
m < data.frame(n.pilot = m, CV = 0.45,
df l = NA, u = NA,
n.low = NA, n.CV = NA, n.hi = NA)
for (i in 1:nrow(df)) {
3:4] < CVCL(CV = 0.45, df = m[i]  2,
df[i, side = "2sided",
alpha = 0.05)
5] < sampleN.scABEL(CV = df$l[i], design = "2x2x4",
df[i, details = FALSE,
print = FALSE)[["Sample size"]]
6] < sampleN.scABEL(CV = 0.45, design = "2x2x4",
df[i, details = FALSE,
print = FALSE)[["Sample size"]]
7] < sampleN.scABEL(CV = df$u[i], design = "2x2x4",
df[i, details = FALSE,
print = FALSE)[["Sample size"]]
}
3:4] < signif(df[, 3:4], 4)
df[, names(df)[3:4] < c("lower CL", "upper CL")
print(df, row.names = FALSE)
# n.pilot CV lower CL upper CL n.low n.CV n.hi
# 12 0.45 0.3069 0.8744 36 28 56
# 18 0.45 0.3282 0.7300 34 28 42
# 24 0.45 0.3415 0.6685 34 28 38
# 30 0.45 0.3509 0.6334 34 28 34
One leading generic company has an internal rule to perform pilot studies of HVDPs in a four period full replicate design and at least 24 subjects. Makes sense.
Furthermore, we don’t know where the true T/Rratio lies but we can calculate the lower 95% confidence limit of the pilot study’s point estimate to get an idea about a worst case. Say, it was 0.90.
< 16
m < 0.45
CV < 0.90
pe < round(CI.BE(CV = CV, pe = 0.90, n = m,
ci design = "2x2x4"), 4)
if (pe <= 1) {
< ci[["lower"]]
cl else {
} < ci[["upper"]]
cl
}print(cl)
# [1] 0.7515
Explore the impact of a relatively 5% lower CV (less expansion) and a relatively 5% lower T/Rratio on power for the given sample size.
< 28
n < 0.45
CV < 0.90
theta0 < data.frame(CV = c(CV, CV*0.95),
comp1 power = NA)
< data.frame(theta0 = c(theta0, theta0*0.95),
comp2 power = NA)
for (i in 1:2) {
$power[i] < power.scABEL(CV = comp1$CV[i],
comp1theta0 = theta0,
design = "2x2x4",
n = n)
}$power < signif(comp1$power, 5)
comp1for (i in 1:2) {
$power[i] < power.scABEL(CV = CV,
comp2theta0 = comp2$theta0[i],
design = "2x2x4",
n = n)
}$power < signif(comp2$power, 5)
comp2print(comp1, row.names = FALSE)
print(comp2, row.names = FALSE)
# CV power
# 0.4500 0.81116
# 0.4275 0.80095
# theta0 power
# 0.900 0.81116
# 0.855 0.61952
Note the logscale of the xaxis. It demonstrates that power curves are symmetrical around 1 (\(\small{\log_{e}(1)=0}\), where \(\small{\log_{e}(\theta_2)=\left\log_{e}(\theta_1)\right}\)) and we will achieve the same power for \(\small{\theta_0}\) and \(\small{1/\theta_0}\) (e.g., for 0.90 and 1.1111). Contrary to ABE, power is maintained, unless we cross the upper scaling limit, where additionally the PEconstraint becomes increasingly important.
<nitpick>< 0.45
CV < 0.10 # direction unknown
delta < "2x2x4"
design < c(1  delta, 1 / (1 + delta),
theta0s 1 + delta, 1 / (1  delta))
< sampleN.scABEL(CV = CV, theta0 = 1  delta,
n design = design,
details = FALSE,
print = FALSE)[["Sample size"]]
< data.frame(CV = CV, theta0 = theta0s,
comp1 base = c(TRUE, rep(FALSE, 3)),
n = n, power = NA)
for (i in 1:nrow(comp1)) {
$power[i] < power.scABEL(CV = CV,
comp1theta0 = comp1$theta0[i],
design = design, n = n)
}< sampleN.scABEL(CV = CV, theta0 = 1 + delta,
n design = design,
details = FALSE,
print = FALSE)[["Sample size"]]
< data.frame(CV = CV, theta0 = theta0s,
comp2 base = c(FALSE, FALSE, TRUE, FALSE),
n = n, power = NA)
for (i in 1:nrow(comp2)) {
$power[i] < power.scABEL(CV = CV,
comp2theta0 = comp2$theta0[i],
design = design, n = n)
}c(2, 5)] < signif(comp1[, c(2, 5)] , 4)
comp1[, c(2, 5)] < signif(comp2[, c(2, 5)] , 4)
comp2[, print(comp1, row.names = FALSE)
print(comp2, row.names = FALSE)
# CV theta0 base n power
# 0.45 0.9000 TRUE 28 0.8112
# 0.45 0.9091 FALSE 28 0.8388
# 0.45 1.1000 FALSE 28 0.8397
# 0.45 1.1110 FALSE 28 0.8101
# CV theta0 base n power
# 0.45 0.9000 FALSE 26 0.7846
# 0.45 0.9091 FALSE 26 0.8149
# 0.45 1.1000 TRUE 26 0.8140
# 0.45 1.1110 FALSE 26 0.7837
</nitpick>
Essentially this leads to the murky waters of prospective sensitivity
analyses, which is covered in another article.
An appetizer to show the maximum deviations (CV, T/Rratio and
decreased sample size due to dropouts) which give still a minimum
acceptable power of ≥ 0.70:
< 0.45
CV < 0.90
theta0 < 0.80
target < 0.70
minpower < pa.scABE(CV = CV, theta0 = theta0,
pa targetpower = target,
minpower = minpower,
design = "2x2x4")
< 100*(tail(pa$paCV[["CV"]], 1) 
change.CV $plan[["CVwR"]]) /
pa$plan[["CVwR"]]
pa< 100*(head(pa$paGMR$theta0, 1) 
change.theta0 $plan$theta0) /
pa$plan[["theta0"]]
pa< 100*(tail(pa$paN[["N"]], 1) 
change.n $plan[["Sample size"]]) /
pa$plan[["Sample size"]]
pa< data.frame(parameter = c("CV", "theta0", "n"),
comp change = c(change.CV,
change.theta0,
change.n))$change < sprintf("%+.2f%%", comp$change)
compnames(comp)[2] < "rel. change"
print(pa, plotit = FALSE)
print(comp, row.names = FALSE)
# Sample size plan scABE (EMA/ABEL)
# Design alpha CVwT CVwR theta0 theta1 theta2 Sample size
# 2x2x4 0.05 0.45 0.45 0.9 0.8 1.25 28
# Achieved power Target power
# 0.81116 0.8
#
# Power analysis
# CV, theta0 and number of subjects leading to min. acceptable power of ~0.7:
# CV= 0.6629, theta0= 0.8719
# n = 22 (power= 0.7185)
#
# parameter rel. change
# CV +47.32%
# theta0 3.12%
# n 21.43%
Confirms what we have seen above. As expect the method is robust to changes of the CV. The sample size is also not very sensitive; many overrate the impact of dropouts on power.
As we have seen already above, for an ANOVA we have to assume homoscedasticity.
I recommend to perform pilot studies in one of the fully replicated
designs. When you are concerned about dropouts in a four period design
or the bioanalytical method requires large sample volumes, opt for one
the two sequence three period designs (TRTRTR or
TRRRTT).
Contrary to the partial replicate design (TRRRTRRRT) we get estimates
of both \(\small{CV_\text{wT}}\) and
\(\small{CV_\text{wR}}\). Since
pharmaceutical technology improves, it is not uncommon that \(\small{CV_\text{wT}<CV_\text{wR}}\). If
this is the case, we get an incentive in the sample size \(\small{n}\) of the planned pivotal study
(expanding the limits is based on \(\small{CV_\text{wR}}\) but the 90%
CI on the – pooled – \(\small{s_\text{w}^{2}}\)).
\[\eqalign{\tag{7} s_\text{wT}^{2}&=\log_{e}(CV_\text{wT}^{2}+1)\\ s_\text{wR}^{2}&=\log_{e}(CV_\text{wR}^{2}+1)\\ s_\text{w}^{2}&=\left(s_\text{wT}^{2}+s_\text{wR}^{2}\right)/2\\ CV_\text{w}&=\sqrt{\exp(s_\text{w}^{2})1}}\]
Say, we performed two pilot studies. In the partial replicate we estimated the \(\small{CV_\text{wR}}\) with 0.45 (estimation of \(\small{CV_\text{wT}}\) is not possible). In the full replicate we estimated \(\small{CV_\text{wR}}\) with 0.484 and \(\small{CV_\text{wT}}\) with 0.414. Note that the \(\small{CV_\text{w}}\) according to \(\small{(7)}\) is 0.45 as well. How will that impact the sample size \(\small{n}\) of the pivotal 4period full replicate design?
< data.frame(pilot = c("TRRRTRRRT", "TRTRTR"),
comp CVwT = c(0.45, 0.414),
CVwR = c(0.45, 0.484),
CVw = NA,
n = NA, power = NA)
for (i in 1:nrow(comp)) {
4] < signif(
comp[i, mse2CV((CV2mse(comp$CVwT[i]) +
CV2mse(comp$CVwR[i])) / 2), 3)
5:6] < sampleN.scABEL(CV = c(comp$CVwT[i], comp$CVwR[i]),
comp[i, design = "2x2x4", details = FALSE,
print = FALSE)[8:9]
}print(comp, row.names = FALSE)
# pilot CVwT CVwR CVw n power
# TRRRTRRRT 0.450 0.450 0.45 28 0.81116
# TRTRTR 0.414 0.484 0.45 24 0.80193
Since bioanalytics drives study costs to a great extent, we may save ~14%.
Note that when you give CV
as twoelement vector, the
first element has to be \(\small{CV_\text{wT}}\) and the second \(\small{CV_\text{wR}}\).
Although according to the guidelines it is not required to estimate \(\small{CV_\text{wT}}\), its value is ‘nice to know’. Sometimes studies fail only due to the large \(\small{CV_\text{wR}}\) thus inflating the confidence interval. In such a case you have at least ammunation to start an argument.
Even if you plan the pivotal study in a partial replicate design (why on earth?) knowing both \(\small{CV_\text{wT}}\) and \(\small{CV_\text{wR}}\) is useful.
< data.frame(pilot = c("TRRRTRRRT", "TRTRTR"),
comp CVwT = c(0.45, 0.414),
CVwR = c(0.45, 0.484),
CVw = NA,
n = NA, power = NA)
for (i in 1:nrow(comp)) {
4] < signif(
comp[i, mse2CV((CV2mse(comp$CVwT[i]) +
CV2mse(comp$CVwR[i])) / 2), 3)
5:6] < sampleN.scABEL(CV = c(comp$CVwT[i], comp$CVwR[i]),
comp[i, design = "2x3x3", details = FALSE,
print = FALSE)[8:9]
}print(comp, row.names = FALSE)
# pilot CVwT CVwR CVw n power
# TRRRTRRRT 0.450 0.450 0.45 39 0.80588
# TRTRTR 0.414 0.484 0.45 36 0.80973
Again, a smaller sample size is possible and we may still save ~7.7%.
Note that sampleN.scABEL()
is inaccurate for the partial
replicate design if and only if \(\small{CV_\text{wT}>CV_\text{wR}}\).
Let’s reverse the values and compare the results with subject
simulations by the function sampleN.scABEL.sdsims()
.
< data.frame(method = c("key", "subj"), n = NA,
comp power = NA, rel.speed = NA)
< proc.time()[[3]]
st 1, 2:3] < sampleN.scABEL(CV = c(0.484, 0.414),
comp[design = "2x3x3",
details = FALSE,
print = FALSE)[8:9]
< proc.time()[[3]]
et $rel.speed[1] < et  st
comp< proc.time()[[3]]
st 2, 2:3] < sampleN.scABEL.sdsims(CV = c(0.484, 0.414),
comp[design = "2x3x3",
details = FALSE,
print = FALSE)[8:9]
< proc.time()[[3]]
et $rel.speed[2] < et  st
comp$rel.speed < signif(comp$rel.speed /
comp$rel.speed[1], 3)
comp1] < c("\u2018key statistics\u2019",
comp["subject simulations")
print(comp, row.names = FALSE)
# method n power rel.speed
# ‘key statistics’ 45 0.80212 1.0
# subject simulations 48 0.80938 84.6
In such a case, always use the function
sampleN.scABEL.sdsims()
.
For demonstrating bioequivalence the pharmacokinetic metrics C_{max} and AUC_{0–t} are mandatory (in some jurisdictions like the FDA additionally AUC_{0–∞}).
We don’t have to worry about multiplicity issues (inflated Type I Error), since if all tests must pass at level \(\alpha\), we are protected by the intersectionunion principle.^{53} ^{54} We design the study always for the worst case combination, i.e., based on the PK metric requiring the largest sample size. In most jurisdictions wider BE limits are acceptable only for C_{max}. Let’s explore that with different CVs and T/Rratios.
< c("Cmax", "AUC")
metrics < c("ABEL", "ABE")
methods < c(0.45, 0.35)
CV < "2x2x4"
design < c(0.90, 0.925)
theta0 < data.frame(metric = metrics, method = methods,
df CV = CV, theta0 = theta0, n = NA,
power = NA)
1, 5:6] < sampleN.scABEL(CV = CV[1], theta0 = theta0[1],
df[design = design, details = FALSE,
print = FALSE)[8:9]
2, 5:6] < sampleN.TOST(CV = CV[2], theta0 = theta0[2],
df[design = design, print = FALSE)[7:8]
$power < signif(df$power, 5)
df< paste0("Sample size is driven by ",
txt $metric[df$n == max(df$n)], ".\n")
dfprint(df, row.names = FALSE); cat(txt)
# metric method CV theta0 n power
# Cmax ABEL 0.45 0.900 28 0.81116
# AUC ABE 0.35 0.925 36 0.81604
# Sample size is driven by AUC.
Commonly the PK metric evaluated by ABE drives the sample size. Consequently, with 36 subjects the study will be ‘overpowered’ (~0.886) for C_{max}.
Let us assume the same T/Rratios for both metrics. Which are the extreme T/Rratios (largest deviations of T from R) for C_{max} giving still the target power?
< function(x) {
opt power.scABEL(theta0 = x, CV = df$CV[1],
design = design,
n = df$n[2])  target
}< c("Cmax", "AUC")
metrics < c("ABEL", "ABE")
methods < c(0.45, 0.30)
CV < c(0.90, 0.90)
theta0 < "2x2x4" # must be a replicate
design < 0.80
target if (!design %in% c("2x2x4", "2x3x3", "2x2x3"))
stop("design = \"", design, "\" is not supported.")
< data.frame(metric = metrics, method = methods,
df CV = CV, theta0 = theta0, n = NA,
power = NA)
1, 5:6] < sampleN.scABEL(CV = CV[1], theta0 = theta0[1],
df[design = design, details = FALSE,
print = FALSE)[8:9]
2, 5:6] < sampleN.TOST(CV = CV[1], theta0 = theta0[2],
df[design = design, print = FALSE)[7:8]
$power < signif(df$power, 5)
dfif (theta0[1] < 1) {
< uniroot(opt, tol = 1e8,
res interval = c(0.80 + 1e4, theta0[1]))
else {
} < uniroot(opt, tol = 1e8,
res interval = c(theta0[1], 1.25  1e4))
}< unlist(res)
res < c(res[["root"]], 1/res[["root"]])
theta0s < paste0("Target power for ", metrics[1],
txt " and sample size ",
$n[2], "\nachieved for theta0 ",
dfsprintf("%.4f", theta0s[1]), " or ",
sprintf("%.4f", theta0s[2]), ".\n")
print(df, row.names = FALSE); cat(txt)
# metric method CV theta0 n power
# Cmax ABEL 0.45 0.9 28 0.81116
# AUC ABE 0.30 0.9 84 0.80569
# Target power for Cmax and sample size 84
# achieved for theta0 0.8340 or 1.1990.
That means, although we assumed for C_{max} the same T/Rratio as for AUC– with the sample size of 84 required AUC, for C_{max} it can be as low as 0.834 or as high as 1.199, which is a soothing sideeffect.
For HVD(P)s Health Canada accepts ABEL only for AUC, whereas for C_{max} the T/Rratio has to lie within 80.0 – 125.0% (i.e., without assessing its CI).
< c("Cmax", "AUC")
metrics < c("PE", "ABE(L)")
methods < c(0.50, 0.05)
alpha < c(0.65, 0.35)
CV < c(0.90, 0.90)
theta0 < "2x2x4" # must be a replicate
design < 0.90
target if (!design %in% c("2x2x4", "2x3x3", "2x2x3"))
stop("design = \"", design, "\" is not supported.")
< data.frame(design = design, metric = metrics, method = methods,
df CV = CV, theta0 = theta0, n = NA, power = NA)
# alpha = 0.5 leads to a zerowidth confidence
# interval, i.e., only the T/Rratio is relevant
1, 6] < sampleN.TOST(alpha = alpha[1],
df[CV = CV[1],
theta0 = theta0[1],
design = design,
targetpower = target,
print = FALSE)[["Sample size"]]
if (df[1, 6] < 12) df[1, 6] < 12 # acc. to the guidance
1, 7] < power.TOST(alpha = alpha[1],
df[CV = CV[1],
theta0 = theta0[1],
design = design,
n = df[1, 6])
2, 6] < sampleN.scABEL(alpha = alpha[2],
df[CV = CV[2],
theta0 = theta0[2],
design = design,
regulator = "HC",
targetpower = target,
details = FALSE,
print = FALSE)[["Sample size"]]
if (df[2, 6] < 12) df[2, 6] < 12 # acc. to the guidance
2, 7] < power.scABEL(alpha = alpha[2],
df[CV = CV[2],
theta0 = theta0[2],
regulator = "HC",
design = design,
n = df[2, 6])
< which(df$n == max(df$n))
driving < which(df$n != max(df$n))
easy if (driving == 1) {
< power.TOST(alpha = alpha[easy],
easy.pwr CV = CV[easy],
theta0 = theta0[easy],
design = design,
n = df$n[driving])
else {
} < power.scABEL(alpha = alpha[easy],
easy.pwr CV = CV[easy],
theta0 = theta0[easy],
regulator = "HC",
design = design,
n = df$n[driving])
}$power < signif(df$power, 5)
dfnames(df)[5] < "T/R"
< paste0(sprintf("Sample size is driven by %s.\n", df$metric[driving]),
txt sprintf("With a sample size of %.0f power for %s = %.5f.\n",
$n[driving], metrics[easy], easy.pwr))
dfprint(df, row.names = FALSE); cat(txt)
# design metric method CV T/R n power
# 2x2x4 Cmax PE 0.65 0.9 42 0.90058
# 2x2x4 AUC ABE(L) 0.35 0.9 50 0.90897
# Sample size is driven by AUC.
# With a sample size of 50 power for Cmax = 0.91980.
Even for a CV of C_{max} which is substantially higher than the one of AUC, the latter still drives the sample size.
Q: Can we use R in a
regulated environment and is PowerTOST
validated?
A: See this document^{55} about the
acceptability of Base R
and its
Software
Development Life Cycle.^{56}
R
is updated every couple of
months with documented changes^{57} and maintaining a bugtracking system.^{58} I
recommed to use always the latest release.
The authors of PowerTOST
tried to do their best to provide
reliable and valid results. Its ‘NEWS’
documents the development of the package, bugfixes, and introduction of
new methods. Issues are tracked at GitHub (as of
today none is open). So far the package had >118,000 downloads.
Therefore, it is extremely unlikely that bugs were not detected given
its large user base.
Validation of any software (yes, of SAS
as well…)
lies in the hands of the user.^{59} ^{60} ^{61} ^{62}
Execute the script test_ABEL.R
which can be found in the
/tests
subdirectory of the package to reproduce tables
given in the literature.^{63} You will notice some discrepancies: The
authors employed only 10,000 simulations – which is not sufficient for a
stable result (see below). Furthermore, they
reported the minimum sample size which gives at least the target power,
whereas sampleN.scABEL
always rounds up to give balanced
sequences.
Q: Shall we throw away our sample size
tables?
A: Not at all. File them in your archives to collect
dust. Maybe in the future you will be asked by an agency how you arrived
at a sample size. But: Don’t use them any more. What you should not do
(and hopefully haven’t done before): Interpolate. Power and therefore,
the sample size depends in a highly nonlinear fashion on the five
conditions listed above, which makes
interpolation of values given in a table a nontrivial job.
Q: Which of the methods should we use in our
daily practice?
A:
sampleN.scABEL()
/power.scABEL()
for speed
reasons. Only for the partial replicate designs and the –
rare – case of CV_{wT} >
CV_{wR}, use
sampleN.scABEL.sdsims()
/power.scABEL.sdsims()
instead.
Q: I fail to understand your example about
dropouts. We finish the study with 28 eligible subjects as desired. Why
is the dropoutrate ~18% and not the anticipated 15%?
A: That’s due to rounding up the calculated adjusted
sample size (32.94…) to the next even number (34).
If you manage it to dose fractional subjects (I can’t), your dropout
rate would indeed equal the anticipated one: 100 × (1 – 28/32.94…) =
15%. ⬜
Q: Do we have to worry about unbalanced
sequences?
A:
sampleN.scABEL()
/sampleN.scABEL.sdsims()
will
always give the total number of subjects for balanced
sequences.
If you are interested in – irrelevant – post hoc power,^{64} give
the sample size as a vector, i.e.,
power.scABEL(..., n = c(foo, bar, baz)
, where
foo
, bar
, and baz
are the number
of subjects per sequence.
Q: The
default number of simulations in the sample size estimation is
100,000. Why?
A: We found that with this number the simulations are
stable. For the background see another article.
Of course, you can give a larger number in the argument
nsims
. However, you shouldn’t decrease the
number.
Q: How reliable are the results?
A: As stated above an exact
method doesn’t – and will never – exist. We can only compare the empiric
power of the
ABEcomponent to
the exact one obtained by power.TOST()
. For an example see
‘Post hoc Power’ in the section about Dropouts.
Q: Why give links to some references an
HTTPerror 404 ‘Not found’?
A: I checked all URLs in July 2024. Contrary
to us mere mortals who have to maintain a version control of documents,
agencies don’t care. They change the structure of their websites (worst
are the ones of the ANVISA
and the WHO), don’t
establish automatic
redirects, rename or even delete files… Quod
licet Iovi, non licet bovi
If you discover an error, please drop me
a note at [email protected]..
Q: I still have questions. How to proceed?
A: The preferred method is to register at the BEBA Forum and post your question in
the category or (please read the Forum’s Policy
first).
You can contact me at [email protected]. Be
warned – I will charge you for anything beyond most basic
questions.
top of section ↩︎ previous section ↩︎
Abbreviation  Meaning 

(A)BE  (Average) Bioequivalence 
ABEL  Average Bioequivalence with Expanding Limits 
ASEAN  Association of Southeast Asian Nations 
CI  \(\small{(12\,\alpha)}\) Confidence Interval 
CL  Confidence Limit 
\(\small{CV_\text{b}}\)  Betweensubject Coefficient of Variation 
\(\small{CV_\text{w}}\)  (Pooled) withinsubject Coefficient of Variation 
\(\small{CV_\text{wR},CV_\text{wT}}\)  Observed withinsubject Coefficient of Variation of the Reference and Test product 
EAC  East African Community 
EEA  European Economic Area 
EEU  Eurasian Economic Union 
EMA  European Medicines Agency 
EU  European Union 
FDA  (U.S.) Food and Drug Administration 
GCC  Gulf Cooperation Council 
\(\small{H_0}\)  Null hypothesis (inequivalence) 
\(\small{H_1}\)  Alternative hypothesis (equivalence) 
HVD(P)  Highly Variable Drug (Product) 
IBE  Individual Bioequivalence 
ICH  International Council for Harmonisation 
\(\small{k}\)  Regulatory constant (0.76) 
NTID  Narrow Therapeutic Index Drug 
PE  Point Estimate 
\(\small{\text{R}}\)  Reference product 
RSABE  ReferenceScaled Average Bioequivalence (FDA and China’s CDE) 
SABE  Scaled Average Bioequivalence (ABEL and RSABE) 
\(\small{s_\text{bR}^2,s_\text{bT}^2}\)  Observed betweensubject variance of the Reference and Test product 
\(\small{s_\text{wR},s_\text{wT}}\)  Observed withinsubject standard deviation of the Reference and Test product 
\(\small{s_\text{wR}^2,s_\text{wT}^2}\)  Observed withinsubject variance of the Reference and Test product 
\(\small{\text{T}}\)  Test product 
TIE  Type I Error (patient’s risk) 
TIIE  Type II Error (producer’s risk) 
T/R  T/Rratio (observed, PE) or assumed in sample size estimation 
\(\small{uc}\)  Upper cap of expansion in ABEL 
WHO  World Health Organization 
\(\small{\alpha}\)  Nominal level of the test, probability of Type I Error (patient’s risk) 
\(\small{\beta}\)  Probability of Type II Error (producer’s risk) 
\(\small{\mu_\text{T}/\mu_\text{R},\theta_0}\)  True T/Rratio 
\(\small{\pi}\)  (Prospective) power, where \(\small{\pi=1\beta}\) 
\(\small{\sigma_\text{wR}}\)  True withinsubject standard deviation of the Reference product 
\(\small{\theta_1,\theta_2}\)  Fixed lower and upper limits of the acceptance range (generally 80.00 – 125.00%) 
\(\small{\theta_{\text{s}_1},\theta_{\text{s}_2}}\)  Expanded lower and upper limits of the acceptance range \(\small{=100\exp(k\,\cdot s_\text{wR})}\); in the guidelines denoted as \(\small{\left\{L,U\right\}}\) 
Helmut Schütz 2024
R
and PowerTOST
GPL 3.0,
klippy
MIT,
pandoc
GPL 2.0.
1^{st} version March 23, 2021. Rendered July 27, 2024 11:59 CEST
by rmarkdown
via pandoc in 0.95 seconds.
FAX
machines.
Labes D, Schütz H, Lang B. PowerTOST: Power and Sample Size for (Bio)Equivalence Studies. Package version 1.5.6. 20240318. CRAN.↩︎
Schütz H. ReferenceScaled Average Bioequivalence. 20240229. CRAN.↩︎
Labes D, Schütz H, Lang B. Package ‘PowerTOST’. March 18, 2024. CRAN.↩︎
Some gastric resistant formulations of diclofenac are HVDPs, practically all topical formulations are HVDPs, whereas diclofenac itself is not a HVD (\(\small{CV_\text{w}}\) of a solution ~8%).↩︎
An inglorious counterexample: In the approval of dabigatran (the first univalent direct thrombin (IIa) inhibitor) the originator withheld information about severe bleeding events. Although it is highly variable, it is considered an NTID by the FDA and referencescaling not acceptable. The study has to performed in a four period replicate design in order to allow a comparison of \(\small{s_\text{wT}}\) with \(\small{s_\text{wR}}\).↩︎
Tóthfalusi L, Endrényi L, GarcíaArieta A. Evaluation of bioequivalence for highly variable drugs with scaled average bioequivalence. Clin Pharmacokinet. 2009; 48(11): 725–43. doi:10.2165/1131804000000000000000.↩︎
ABEL is one variant of SABE. Another is ReferenceScaled Average Bioequivalence (RSABE), which is recommended by the FDA and China’s CDE. It is covered in another article.↩︎
GCC. Executive Directorate Of Benefits And Risks Evaluation. The GCC Guidelines for Bioequivalence. Version 2.4. 30/3/2016. Internet Archive.↩︎
Final Concept Paper. M13: Bioequivalence for ImmediateRelease Solid Oral Dosage Forms. Geneva. 19 June 2020. Endorsed 10 July 2020. Online.↩︎
EMA, CHMP. Guideline on the Investigation of Bioequivalence. London. 20 January 2010. Online.↩︎
WHO, Essential Medicines and Health Products: Multisource (generic) pharmaceutical products: guidelines on registration requirements to establish interchangeability. WHO Technical Report Series, No. 1003, Annex 6. Geneva. 28 April 2017. Online↩︎
WHO/PQT: medicines. Guidance Document: Application of referencescaled criteria for AUC in bioequivalence studies conducted for submission to PQT/MED. Geneva. 02 July 2021. Online.↩︎
Australian Government, Department of Health and Aged Care, TGA. International scientific guideline: Guideline on the investigation of bioequivalence. Canberra. Effective 16 June 2011. Online.↩︎
Shohin LE, Rozhdestvenkiy DA, Medvedev VYu, Komarow TN, Grebenkin DYu. Russia, Belarus & Kazakhstan. In: Kanfer I, editor. Bioequivalence Requirements in Various Global Jurisdictions. Charm: Springer; 2017. pp. 215–6.↩︎
EAC, Medicines and Food Safety Unit. Compendium of Medicines Evaluation and Registration for Medicine Regulation Harmonization in the East African Community, Part III: EAC Guidelines on Therapeutic Equivalence Requirements. Online.↩︎
ASEAN States Pharmaceutical Product Working Group. ASEAN Guideline for the Conduct of Bioequivalence Studies. Vientiane. March 2015. Online.↩︎
Shohin LE, Rozhdestvenkiy DA, Medvedev VYu, Komarow TN, Grebenkin DYu. Russia, Belarus & Kazakhstan. In: Kanfer I, editor. Bioequivalence Requirements in Various Global Jurisdictions. Charm: Springer; 2017. p. 207.↩︎
Eurasian Economic Commission. Regulations Conducting Bioequivalence Studies within the Framework of the Eurasian Economic Union. 3 November 2016, amended 15 February 2023. Online. [Russian]↩︎
Ministry of Health and Population, The Specialized Scientific Committee for Evaluation of Bioavailability & Bioequivalence Studies. Egyptian Guideline For Conducting Bioequivalence Studies for Marketing Authorization of Generic Products. Cairo. February 2017. Internet Archive.↩︎
New Zealand Medicines and Medical Devices Safety Authority. Guideline on the Regulation of Therapeutic Products in New Zealand. Part 6: Bioequivalence of medicines. Wellington. February 2018. Online.↩︎
Departamento Agencia Nacional de Medicamentos. Instituto de Salud Pública de Chile. Guia para La realización de estudios de biodisponibilidad comparativa en formas farmacéuticas sólidas de administración oral y acción sistémica. Santiago. December 2018. [Spanish]↩︎
ANVISA. Dispõe sobre os critérios para a condução de estudos de biodisponibilidade relativa/bioequivalência (BD/BE) e estudos farmacocinéticos. Brasilia. 17 August 2022. Online. [Portuguese]↩︎
Executive Board of the Health Ministers’ Council for GCC States. The GCC Guidelines for Bioequivalence. Version 3.1. 10 August 2022. Online.↩︎
Secretaría de Salud. Disposiciones para los Estudios de Bioequivalencia para Fármacos Altamente Variables. Ciudad de México. 24 January 2023. Online. [Spanish]↩︎
Health Canada. Guidance Document. Comparative Bioavailability Standards: Formulations Used for Systemic Effects. Ottawa. 08 June 2018. Online.↩︎
EMA, CHMP. Guideline on the pharmacokinetic and clinical evaluation of modified release dosage forms. London. 20 November 2014. Online.↩︎
Midha KK, Rawson MJ, Hubbard JW. Prescribability and switchability of highly variable drugs and drug products. J Contr Rel. 1999; 62(12): 33–40. doi:10.1016/s01683659(99)000504.↩︎
Labes D, Schütz H. Inflation of Type I Error in the Evaluation of Scaled Average Bioequivalence, and a Method for its Control. Pharm Res. 2016: 33(11); 2805–14. doi:10.1007/s1109501620061.↩︎
Schütz H, Labes D, Wolfsegger MJ. Critical Remarks on ReferenceScaled Average Bioequivalence. J Pharm Pharmaceut Sci. 2022; 25: 285–296. doi:10.18433/jpps32892. Open Access.↩︎
Schütz H, Tomashevskiy M, Labes D. replicateBE: Average Bioequivalence with Expanding Limits (ABEL). Package version 1.1.3. 20220502. CRAN.↩︎
Karalis V, Symillides M, Macheras P. On the levelingoff properties of the new bioequivalence limits for highly variable drugs of the EMA guideline. Europ J Pharm Sci. 2011; 44: 497–505. doi:10.1016/j.ejps.2011.09.008.↩︎
Commission of the European Communities, CPMP Working Party on Efficacy of Medicinal Products. Note for Guidance: Investigation of Bioavailability and Bioequivalence. Appendix III: Technical Aspects of Bioequivalence Statistics. Brussels. December 1991.↩︎
Commission of the European Communities, CPMP Working Party on Efficacy of Medicinal Products. Note for Guidance: Investigation of Bioavailability and Bioequivalence. June 1992.↩︎
Blume H, Mutschler E, editors. Bioäquivalenz. Qualitätsbewertung wirkstoffgleicher Fertigarzneimittel. Frankfurt / Main: GoviVerlag; 6. Ergänzungslieferung 1996. [German]↩︎
EMA, CHMP EWP, Therapeutic Subgroup on Pharmacokinetics. Questions & Answers on the Bioavailability and Bioequivalence Guideline. London. 27 July 2006. Online.↩︎
MCC. Registration of Medicines. Biostudies. Pretoria. June 2015. Online.↩︎
League of Arab States, Higher Technical Committee for Arab Pharmaceutical Industry. Harmonised Arab Guideline on Bioequivalence of Generic Pharmaceutical Products. Cairo. March 2014. Online.↩︎
Shohin LE, Rozhdestvenkiy DA, Medvedev VYu, Komarow TN, Grebenkin DYu. Russia, Belarus & Kazakhstan. In: Kanfer I, editor. Bioequivalence Requirements in Various Global Jurisdictions. Charm: Springer; 2017. p. 222.↩︎
Of note, for NTIDs Health Canada’s limits are 90.0 – 112.0% and not 90.00 – 111.11% like in other jurisdictions. Guess the reason.↩︎
That’s contrary to ABE, where \(\small{CV_\text{w}}\) is an assumption as well.↩︎
Senn S. Guernsey McPearson’s Drug Development Dictionary. 21 April 2020. Online.↩︎
Hoenig JM, Heisey DM. The Abuse of Power: The Pervasive Fallacy of Power Calculations for Data Analysis. Am Stat. 2001; 55(1): 19–24. doi:10.1198/000313001300339897. Open Access.↩︎
There is no statistical method to ‘correct’ for unequal carryover. It can only be avoided by design, i.e., a sufficiently long washout between periods. According to the guidelines subjects with predose concentrations > 5% of their C_{max} can by excluded from the comparison if stated in the protocol.↩︎
This is not always the case in ABEL. This issue is elaborated in another article.↩︎
Senn S. Statistical Issues in Drug Development. Chichester: John Wiley; 2^{nd} ed. 2007.↩︎
Zhang P. A Simple Formula for Sample Size Calculation in Equivalence Studies. J Biopharm Stat. 2003; 13(3): 529–38. doi:10.1081/BIP120022772.↩︎
Don’t be tempted to assume a ‘better’ T/Rratio – even if based on a pilot or a previous study. It is a natural property of HVDPs that the T/Rratio varies between studies. Don’t be overly optimistic!↩︎
Doyle AC. The Adventures of Sherlock Holmes. A Scandal in Bohemia. 1892. p. 3.↩︎
Schütz H. Sample Size Estimation in Bioequivalence. Evaluation. 20201023. BEBA Forum. Vienna. 20201023. Online.↩︎
Lenth RV. Two SampleSize Practices that I Don’t Recommend. October 24, 2000. Online.↩︎
The only exception is design = "2x2x3"
(the full replicate with sequences TRTRTR). Then the first element is
for sequence TRT and the second for RTR.↩︎
Quoting my late father: »If you believe, go to church.«↩︎
Berger RL, Hsu JC. Bioequivalence Trials, IntersectionUnion Tests and Equivalence Confidence Sets. Stat Sci. 1996; 11(4): 283–302. JSTOR:2246021.↩︎
Zeng A. The TOST confidence intervals and the coverage probabilities with R simulation. March 14, 2014. Online.↩︎
The R Foundation for Statistical Computing. A Guidance Document for the Use of R in Regulated Clinical Trial Environments. Vienna. October 18, 2021. Online.↩︎
The R Foundation for Statistical Computing. R: Software Development Life Cycle. A Description of R’s Development, Testing, Release and Maintenance Processes. Vienna. October 18, 2021. Online.↩︎
FDA. General Principles of Software Validation; Final Guidance for Industry and FDA Staff. Rockville. January 11, 2002. Download.↩︎
FDA. Statistical Software Clarifying Statement. May 6, 2015. Online.↩︎
WHO. Guidance for organizations performing in vivo bioequivalence studies. Technical Report Series No. 996, Annex 9. Section 4. Geneva. May 2016. Online.↩︎
EMA, GCP IWG. Guideline on computerised systems and electronic data in clinical trials. Amsterdam. 9 March 2023. Online.↩︎
Tóthfalusi L, Endrényi L. Sample Sizes for Designing Bioequivalence Studies for Highly Variable Drugs. J Pharm Pharmacol Sci. 2012; 15(1): 73–84. doi:10.18433/j3z88f. Open Access.↩︎
User ‘BEQool’. ABEL is a framework (decision scheme). BEBA Forum. Vienna. 20240130. Online.↩︎