- Disclaimer
- Introduction
- Methods
- Tested designs
- Data structure
- Notes on the methods
- Applicability, outlook
- Cross-validation
- Contributors
- Session Information

**Program offered for Use without any Guarantees and Absolutely No Warranty. No Liability is accepted for any Loss and Risk to Public Health Resulting from Use of this R-Code.**

The library provides data sets (internal `.rda`

and in CSV-format in `/extdata/`

) given by Schütz *et al*.^{1} which support users in a black-box performance qualification (PQ) of their software installations.

The methods given in the EMA’s Q&A document for reference-scaling according to the EMA’s Guideline on the Investigation of Bioequivalence are implemented. Potential influence of outliers on the variability of the reference can be assessed by box plots of studentized and standardized residuals as suggested at a joint EGA/EMA workshop.

`method.A()`

A linear model of log-transformed PK responses and effects

*sequence*, *subject(sequence)*, *period*, *treatment*

where all effects are fixed (*i.e.*, ANOVA). Estimated via function `lm()`

of library `stats`

.

`method.B()`

A linear model of log-transformed PK responses and effects

*sequence*, *subject(sequence)*, *period*, *treatment*

where *subject(sequence)* is a random effect and all others are fixed.

Two options are provided

- Estimated via function
`lmer()`

of library`lmerTest`

. Uses Satterthwaite’s degrees of freedom.

`method.B(..., option = 1)`

- Estimated via function
`lme()`

of library`nlme`

. Uses degrees of freedom equivalent to SAS’`DDFM=CONTAIN`

and Phoenix/WinNonlin’s`DF Residual`

. Implicitly preferred according to the Q&A document.

`method.B(..., option = 2)`

, the default.

Details about the reference data sets:

Both the test and the reference treatments are administered at least once.

`TRTR | RTRT`

`TRRT | RTTR`

`TTRR | RRTT`

`TRTR | RTRT | TRRT | RTTR`

^{3}

`TRRT | RTTR | TTRR | RRTT`

^{4}

The test treatment is administered at least once to ½ of the subjects and the reference treatment at least once to the respective other ½ of the subjects.

`TRT | RTR`

`TRR | RTT`

The test and reference treatments are administered once to ½ of the subjects (for the estimation of the CI). *I.e.*, the first group of subjects follow a conventional 2×2×2 trial. In the second group the test and reference treatments are administered at least once to ¼ of the subjects, repectively (for the estimation of *CV _{wT}* and

`TR | RT | TT | RR`

^{5}

The test treatment is administered once and the reference treatment at least once.

`TRR | RTR | RRT`

`TRR | RTR`

^{6}

Columns must have the headers `subject`

, `period`

, `sequence`

, `treatment`

, `PK`

, and/or `logPK`

.^{7} Any order of columns is acceptable. Uppercase and mixed case headers will be internally converted to lowercase headers.

variable | format |
---|---|

`subject` |
Integer numbers or any combination of alphanumerics (A-Z, a-z, -, _, #, 0-9) |

`period` |
Integer numbers |

`sequence` |
Numbers or any other literal sequences not listed in the tested designs are not accepted (e.g., `ABAB` ). |

`treatment` |
The Test treatment must be coded `T` and the Reference `R` (both uppercase). |

Relevant data are used for the estimation of *CV _{wR}* (and

If a subject drops out from the study in a higher period, data of repeated administrations will still be used for the estimation of *CV _{w}*, although data of the other treatment might be missing. Examples for the estimation of

`.`

):`RTRT`

→ `RTR.`

`TRRT`

→ `TRR.`

`RRTT`

→ `RRT.`

or `RR..`

`RRT`

→ `RR..`

If a subject drops out from the study in a higher period, data with at least one administration of the Test and Reference will be used in the assessment of BE. Examples (missings denoted by `.`

):

`TRTR`

→ `TRT.`

or → `TR..`

`RTRT`

→ `RTR.`

or → `RT..`

`TRRT`

→ `TRR.`

or → `TR..`

`RTTR`

→ `RTT.`

or → `RT..`

`TTRR`

→ `TTR.`

`TRT`

→ `TR.`

`RTR`

→ `RT.`

`TRR`

→ `TR.`

`RTR`

→ `RT.`

`RTT`

→ `RT.`

Study started with 16 subjects. In sequence `RTRT`

one dropout in the 2^{nd} period and one in the 4^{th} period. In sequence `TRTR`

one dropout in the 3^{rd} period and one in the 4^{th}.

`1 RTR. 5 RTRT 9 TRTR 13 RTRT`

`2 RTRT 6 TR.. 10 TRTR 14 TRT.`

`3 RTRT 7 RTRT 11 RTRT 15 TRTR`

`4 TRTR 8 R... 12 TRTR 16 TRTR`

We obtain these data sets:

purpose | included subjects | excluded |
---|---|---|

Estimation of CV_{wR} |
13 who received 2 treatments `R` |
6, 8, 14 |

Estimation of CV_{wT} |
12 who received 2 treatments `T` |
1, 6, 8 |

Assessment of BE | 15 who received ≥1 treatment `T` and `R` |
8 |

The EMA proposed a linear model of log-transformed PK responses of the reference treatment

*sequence*, *subject(sequence)*, *period*

where all effects are fixed. Estimated via function `lm()`

of library `stats`

:

For informational purposes in full replicate designs (required by the WHO for reference-scaling of AUC) the same model is run with `data = data[data$treatment = "T", ]`

.

Special conditions for the sample size in three period full replicate designs:

The question raised asks if it is possible to use a design where subjects are randomised to receive treatments in the order of TRT or RTR.

[…] if a 3-period replicate design, where treatments are given in the order TRT or RTR, is to be used to justify widening of a confidence interval for C_{max}then it is considered that at least 12 patients would need to provide data from the RTR arm. This implies a study with at least 24 patients in total would be required if equal number of subjects are allocated to the 2 treatment sequences.— Q&A document (June 2015 and later revisions)

If less than twelve subjects remain in sequence `RTR`

of a `TRT | RTR`

design (and in analogy in sequence `TRR`

of a `TRR | RTT`

design), the user is notified about the uncertain estimate of *CV _{wR}*. However, in a sufficiently powered study such a case is extremely unlikely.

The EMA’s models assumes equal [*sic*] intra-subject variances of Test and Reference (like in 2×2×2 trials) – even if proven false in one of the full replicate designs (were *both* *CV _{wT}* and

The nested structure *subject(sequence)* of the methods leads to an over-specificed model.^{8} The simple model

*sequence*, *subject*, *period*, *treatment*

gives identical estimates of the residual variance and the treatment effect and hence, its confidence interval.

The same holds true for the EMA’s model to estimate *CV _{wR}*. The simple model

gives an identical estimate of the residual variance.

Reference-scaling is acceptable for C_{max} (immediate release products^{9}) and C_{max}, C_{max,ss}, C_{τ,ss}, _{partial}AUC (modified release products^{10}) if high variability of the reference product (*CV _{wR}* >30%) was demonstrated.

The intention to widen the limits has to be stated in the protocol and – contrary to the FDA’s RSABE – a clinical justification provided.

Those HVDP for which a wider difference in C_{max}is considered clinically irrelevant based on a sound clinical justification can be assessed with a widened acceptance range. The request for widened interval must be prospectively specified in the protocol.— BE Guideline (January 2010)

The limits can be expanded based on *CV _{wR}*.

*CV*≤ 30%_{wR}

Lower cap,*i.e.*, no scaling, conventional limits:

\(\left[ {L,U} \right] = {80.00 - 125.00\%}\)- 30% <
*CV*≤ 50%_{wR}

Expanded limits based on \(s_{wR}\):

\(\left[ {L,U} \right] = 100{e^{ \mp 0.760 \cdot {s_{wR}}}}\) *CV*> 50%_{wR}

Upper cap,*i.e.*, applying \(s^*_{wR}=\sqrt{log(0.50^2+1)}\) in the expansion formula or

\(\left[ {L,U} \right] = {69.84 - 143.19\%}\).

In reference-scaling a so-called *mixed* (a.k.a. *aggregate*) criterion is applied. In order to pass BE,

- the 90% confidence interval has to lie entirely within the acceptance range \(\left[ {L,U} \right]\)
*and* - the point estimate has to lie within 80.00 – 125.00%.

In order to avoid discontinuities due to double rounding, expanded limits are calculated in full numeric precision and only the confidence interval is rounded according to the GL.

```
# Calculate limits with library PowerTOST
CV <- c(29, 30, 40, 50, 51)
df <- data.frame(CV = CV, L = NA, U = NA, cap = "",
stringsAsFactors = FALSE)
for (i in seq_along(CV)) {
df[i, 2:3] <- sprintf("%.8f", PowerTOST:::scABEL(CV[i]/100)*100)
}
df$cap[df$CV <= 30] <- "lower"
df$cap[df$CV >= 50] <- "upper"
names(df)[1:3] <- c("CV(%)", "L(%)", "U(%)")
print(df, row.names = FALSE)
#> CV(%) L(%) U(%) cap
#> 29 80.00000000 125.00000000 lower
#> 30 80.00000000 125.00000000 lower
#> 40 74.61770240 134.01645559
#> 50 69.83678198 143.19101936 upper
#> 51 69.83678198 143.19101936 upper
```

The SAS code provided by the EMA in the Q&A document does not specify how the degrees of freedom should be calculated in ‘Method B’. Hence, the default in `PROC MIXED`

, namely `DDFM=CONTAIN`

is applied, *i.e.*, `method.B(..., option = 2`

). For incomplete data (missing periods) Satterthwaite’s approximation of the degrees of freedom, *i.e.*, `method.B(..., option = 1)`

might be a better choice – if stated as such in the SAP.

The EMA seemingly prefers ‘Method A’:

A simple linear mixed model, which assumes identical within-subject variability (Method B), may be acceptable as long as results obtained with the two methods do not lead to different regulatory decisions. However, in borderline cases […] additional analysis using Method A might be required.— Q&A document (January 2011 and later revisions)

The half-width of the confidence interval in log-scale allows a comparison of methods (B *v.s.* A) where a higher value *might* point towards a more conservative decision.^{11} In the provided example data sets – with one exception – the conclusion of BE (based on the mixed criterion) agrees between ‘Method A’ and ‘Method B’.

However, for the highly incomplete data set 14 ‘Method B’ was *liberal* (passing by ANOVA but failing by the mixed effects model):

```
# Compare Method B acc. to the GL with Method A for all reference data sets.
ds <- substr(grep("rds", unname(unlist(data(package = "replicateBE"))),
value = TRUE), start = 1, stop = 5)
for (i in seq_along(ds)) {
A <- method.A(print=FALSE, details = TRUE, data = eval(parse(text = ds[i])))$BE
B <- method.B(print=FALSE, details = TRUE, data = eval(parse(text = ds[i])))$BE
r <- paste0("A ", A, ", B ", B, " \u2013 ")
cat(paste0(ds[i], ":"), r)
if (A == B) {
cat("Methods agree.\n")
} else {
if (A == "fail" & B == "pass") {
cat("Method A is conservative.\n")
} else {
cat("Method B is conservative.\n")
}
}
}
#> rds01: A pass, B pass – Methods agree.
#> rds02: A pass, B pass – Methods agree.
#> rds03: A pass, B pass – Methods agree.
#> rds04: A fail, B fail – Methods agree.
#> rds05: A pass, B pass – Methods agree.
#> rds06: A pass, B pass – Methods agree.
#> rds07: A pass, B pass – Methods agree.
#> rds08: A pass, B pass – Methods agree.
#> rds09: A pass, B pass – Methods agree.
#> rds10: A pass, B pass – Methods agree.
#> rds11: A pass, B pass – Methods agree.
#> rds12: A fail, B fail – Methods agree.
#> rds13: A fail, B fail – Methods agree.
#> rds14: A pass, B fail – Method B is conservative.
#> rds15: A fail, B fail – Methods agree.
#> rds16: A fail, B fail – Methods agree.
#> rds17: A fail, B fail – Methods agree.
#> rds18: A fail, B fail – Methods agree.
#> rds19: A fail, B fail – Methods agree.
#> rds20: A fail, B fail – Methods agree.
#> rds21: A fail, B fail – Methods agree.
#> rds22: A pass, B pass – Methods agree.
#> rds23: A pass, B pass – Methods agree.
#> rds24: A pass, B pass – Methods agree.
#> rds25: A pass, B pass – Methods agree.
#> rds26: A fail, B fail – Methods agree.
#> rds27: A pass, B pass – Methods agree.
#> rds28: A pass, B pass – Methods agree.
```

Exploring data set 14:

```
A <- method.A(print = FALSE, details = TRUE, data = rds14)
B2 <- method.B(print = FALSE, details = TRUE, data = rds14, option = 2)
B1 <- method.B(print = FALSE, details = TRUE, data = rds14, option = 1)
# Rounding of CI according to the GL
A[15:19] <- round(A[15:19], 2) # all effects fixed
B1[15:19] <- round(B1[15:19], 2) # Satterthwaite's df
B2[15:19] <- round(B2[15:19], 2) # df acc. to Q&A
df <- rbind(A[c(2, 15:23)], B1[c(2, 15:23)], B2[c(2, 15:23)])
names(df)[c(1, 10)] <- c("Meth.", "hw")
df$hw <- signif(df$hw, 5)
print(df[order(df$BE, df$hw, decreasing = c(FALSE, TRUE)), ],
row.names = FALSE)
#> Meth. EL.lo(%) EL.hi(%) CI.lo(%) CI.hi(%) PE(%) CI GMR BE hw
#> B-1 69.84 143.19 69.21 121.27 91.62 fail pass fail 0.28043
#> B-2 69.84 143.19 69.21 121.28 91.62 fail pass fail 0.28046
#> A 69.84 143.19 69.99 123.17 92.85 pass pass pass 0.28261
```

Both variants of ‘Method B’ are more conservative than ‘Method A’. Before rounding the confidence interval, `option = 2`

with 192 degrees of freedom would be more conservative (lower CL 69.21029) than `option = 1`

with 197.44 degrees of freedom (lower CL 69.21286). However, given the incompleteness of this data set (four missings in period 2, twelve in period 3, and 19 in period 4), Satterthwaite’s degrees of freedom probably are the better choice.

It is an open issue how outliers should be handled.

The applicant should justify that the calculated intra-subject variability is a reliable estimate and that it is not the result of outliers.

Box plots were ‘suggested’ by the author as a joke [*sic*] at the EGA/EMA workshop, being aware of their nonparametric nature and the EMA’s reluctance towards robust methods. Unfortunately, this *joke* was included in the Q&A document.

[…] a study could be acceptable if the bioequivalence requirements are met both including the outlier subject (using the scaled average bioequivalence approach and the within-subject CV with this subject) and after exclusion of the outlier (using the within-subject CV without this subject).

An outlier test is not an expectation of the medicines agencies but outliers could be shown by a box plot. This would allow the medicines agencies to compare the data between them.— Q&A

With the additional argument `ola = TRUE`

in `method.A()`

and `method.B()`

an outlier analysis is performed, where the default `fence = 2`

^{12} is used.

Results slightly differ depending on software’s algorithms to calculate the median and quartiles. Example with the methods implemented in R:

```
### Compare different types with some random data
x <- rnorm(48)
p <- c(25, 50, 75)/100
q <- matrix(data = "", nrow = 9, ncol = 4,
dimnames=list(paste("type", 1:9),
c("1st quart.", "median", "3rd quart.",
"software / default")))
for (i in 1:9) {
q[i, 1:3] <- sprintf("%.5f", quantile(x, prob = p, type = i))
}
q[c(2, 4:8), 4] <- c("SAS, STaTa", "SciPy", "",
"Phoenix, Minitab, SPSS",
"R, S, MATLAB, Excel", "Maple")
print(as.data.frame(q))
#> 1st quart. median 3rd quart. software / default
#> type 1 -0.67025 -0.19369 0.46120
#> type 2 -0.65502 -0.17713 0.47472 SAS, STaTa
#> type 3 -0.67025 -0.19369 0.46120
#> type 4 -0.67025 -0.19369 0.46120 SciPy
#> type 5 -0.65502 -0.17713 0.47472
#> type 6 -0.66263 -0.17713 0.48148 Phoenix, Minitab, SPSS
#> type 7 -0.64740 -0.17713 0.46796 R, S, MATLAB, Excel
#> type 8 -0.65755 -0.17713 0.47698 Maple
#> type 9 -0.65692 -0.17713 0.47641
```

- Box plots of studentized
^{13}and standarized^{14}model residuals are constructed.^{15} - Potential outliers are flagged based on the argument
`fence`

provided by the user. - With the additional argument
`verbose = TRUE`

detailed information is shown in the console.

Example for the reference data set 01:

```
Outlier analysis
(externally) studentized residuals
Limits (2×IQR whiskers): -1.717435, 1.877877
Outliers:
subject sequence stud.res
45 RTRT -6.656940
52 RTRT 3.453122
standarized (internally studentized) residuals
Limits (2×IQR whiskers): -1.69433, 1.845333
Outliers:
subject sequence stand.res
45 RTRT -5.246293
52 RTRT 3.214663
```

If based on studentized residuals outliers are detected, additionally to the scaled limits based on the complete reference data, tighter limits are calculated based on *CV _{wR}* after exclusion of outliers and BE assessed with the new limits. Standardized residuals are shown for informational purposes only and are not used for exclusion of outliers.

Output for the reference data set 01 (re-ordered for clarity):

```
CVwR : 46.96% (reference-scaling applicable)
swR : 0.44645
Expanded limits : 71.23% ... 140.40% [100exp(±0.760·swR)]
Assessment based on original CVwR 46.96%
────────────────────────────────────────
Confidence interval: 107.11% ... 124.89% pass
Point estimate : 115.66% pass
Mixed (CI & PE) : pass
╟────────┼─────────────────────┼───────■────────◊─────────■───────────────╢
Outlier fence : 2×IQR of studentized residuals.
Recalculation due to presence of 2 outliers (subj. 45|52)
─────────────────────────────────────────────────────────
CVwR (outl. excl.) : 32.16% (reference-scaling applicable)
swR (recalculated) : 0.31374
Expanded limits : 78.79% ... 126.93% [100exp(±0.760·swR)]
Assessment based on recalculated CVwR 32.16%
────────────────────────────────────────────
Confidence interval: pass
Point estimate : pass
Mixed (CI & PE) : pass
╟┼─────────────────────┼───────■────────◊─────────■─╢
```

Note that the PE and its CI are not affected since the entire data are used and therefore, these values not reported in the second analysis (only the conclusion of the assessement).

The line plot is given for informational puposes (its resolution is only ~0.5%). The filled squares ` ■ `

are the lower and upper 90% confidence limits, the rhombus ` ◊ `

the point estimate, the vertical lines ` │ `

at 100% and the PE restriction (80.00 – 125.00%), and the double vertical lines ` ║ `

the expanded limits. The PE and CI take presedence over other symbols. In this case the upper limit of the PE restriction is not visible.

Since both analyses arrive at the same conclusion, the study should be acceptable according to the Q&A document.

The EMA’s approach of reference-scaling for highly variable drugs / drug products is currently recommended in other jurisdictions as well (*e.g.*, the WHO; ASEAN States, Australia, Brazil, Egypt, the Russian Federation, the Eurasian Economic Union, the East African Community, New Zealand).

The WHO accepts reference-scaling for AUC (four period full replicate studies are mandatory in order to assess the variability associated with each product). It is not evident how this assessment should be done. In Population Bioequivalence (PBE) and Individual Bioequivalence (IBE) the *s _{wT}*/

Results of reference data sets agree with ones obtained in SAS (9.3, 9.4) and Phoenix/WinNonlin (6.4 – 8.1).

- Helmut Schütz (Author)
- Michael Tomashevskiy (Contributor)
- Detlew Labes (Contributor)

Inspect this information for reproducibility. Of particular importance are the versions of R and the packages used to create this workflow. It is considered good practice to record this information with every analysis.

```
options(width = 80)
devtools::session_info()
#> - Session info ---------------------------------------------------------------
#> setting value
#> version R version 3.6.0 (2019-04-26)
#> os Windows 7 x64 SP 1
#> system x86_64, mingw32
#> ui RTerm
#> language EN
#> collate C
#> ctype German_Germany.1252
#> tz Europe/Vienna
#> date 2019-06-14
#>
#> - Packages -------------------------------------------------------------------
#> package * version date lib source
#> assertthat 0.2.1 2019-03-21 [2] CRAN (R 3.6.0)
#> backports 1.1.4 2019-04-10 [2] CRAN (R 3.6.0)
#> boot 1.3-22 2019-04-02 [2] CRAN (R 3.6.0)
#> callr 3.2.0 2019-03-15 [2] CRAN (R 3.6.0)
#> cellranger 1.1.0 2016-07-27 [2] CRAN (R 3.6.0)
#> cli 1.1.0 2019-03-19 [2] CRAN (R 3.6.0)
#> colorspace 1.4-1 2019-03-18 [2] CRAN (R 3.6.0)
#> crayon 1.3.4 2017-09-16 [2] CRAN (R 3.6.0)
#> cubature 2.0.3 2018-12-18 [2] CRAN (R 3.6.0)
#> desc 1.2.0 2018-05-01 [2] CRAN (R 3.6.0)
#> devtools 2.0.2 2019-04-08 [2] CRAN (R 3.6.0)
#> digest 0.6.19 2019-05-20 [2] CRAN (R 3.6.0)
#> dplyr 0.8.1 2019-05-14 [2] CRAN (R 3.6.0)
#> evaluate 0.14 2019-05-28 [2] CRAN (R 3.6.0)
#> fs 1.3.1 2019-05-06 [2] CRAN (R 3.6.0)
#> ggplot2 3.1.1 2019-04-07 [2] CRAN (R 3.6.0)
#> glue 1.3.1 2019-03-12 [2] CRAN (R 3.6.0)
#> gtable 0.3.0 2019-03-25 [2] CRAN (R 3.6.0)
#> htmltools 0.3.6 2017-04-28 [2] CRAN (R 3.6.0)
#> knitr 1.23 2019-05-18 [2] CRAN (R 3.6.0)
#> lattice 0.20-38 2018-11-04 [2] CRAN (R 3.6.0)
#> lazyeval 0.2.2 2019-03-15 [2] CRAN (R 3.6.0)
#> lme4 1.1-21 2019-03-05 [2] CRAN (R 3.6.0)
#> lmerTest 3.1-0 2019-02-11 [2] CRAN (R 3.6.0)
#> magrittr 1.5 2014-11-22 [2] CRAN (R 3.6.0)
#> MASS 7.3-51.4 2019-03-31 [2] CRAN (R 3.6.0)
#> Matrix 1.2-17 2019-03-22 [2] CRAN (R 3.6.0)
#> memoise 1.1.0 2017-04-21 [2] CRAN (R 3.6.0)
#> minqa 1.2.4 2014-10-09 [2] CRAN (R 3.6.0)
#> munsell 0.5.0 2018-06-12 [2] CRAN (R 3.6.0)
#> mvtnorm 1.0-10 2019-03-05 [2] CRAN (R 3.6.0)
#> nlme 3.1-140 2019-05-12 [2] CRAN (R 3.6.0)
#> nloptr 1.2.1 2018-10-03 [2] CRAN (R 3.6.0)
#> numDeriv 2016.8-1.1 2019-06-06 [2] CRAN (R 3.6.0)
#> pillar 1.4.1 2019-05-28 [2] CRAN (R 3.6.0)
#> pkgbuild 1.0.3 2019-03-20 [2] CRAN (R 3.6.0)
#> pkgconfig 2.0.2 2018-08-16 [2] CRAN (R 3.6.0)
#> pkgload 1.0.2 2018-10-29 [2] CRAN (R 3.6.0)
#> plyr 1.8.4 2016-06-08 [2] CRAN (R 3.6.0)
#> PowerTOST 1.4.7.9000 2019-05-14 [2] local
#> prettyunits 1.0.2 2015-07-13 [2] CRAN (R 3.6.0)
#> processx 3.3.1 2019-05-08 [2] CRAN (R 3.6.0)
#> ps 1.3.0 2018-12-21 [2] CRAN (R 3.6.0)
#> purrr 0.3.2 2019-03-15 [2] CRAN (R 3.6.0)
#> R6 2.4.0 2019-02-14 [2] CRAN (R 3.6.0)
#> Rcpp 1.0.1 2019-03-17 [2] CRAN (R 3.6.0)
#> readxl 1.3.1 2019-03-13 [2] CRAN (R 3.6.0)
#> remotes 2.0.4 2019-04-10 [2] CRAN (R 3.6.0)
#> replicateBE * 1.0.8 2019-06-14 [1] local
#> rlang 0.3.4 2019-04-07 [2] CRAN (R 3.6.0)
#> rmarkdown 1.13 2019-05-22 [2] CRAN (R 3.6.0)
#> rprojroot 1.3-2 2018-01-03 [2] CRAN (R 3.6.0)
#> scales 1.0.0 2018-08-09 [2] CRAN (R 3.6.0)
#> sessioninfo 1.1.1 2018-11-05 [2] CRAN (R 3.6.0)
#> stringi 1.4.3 2019-03-12 [2] CRAN (R 3.6.0)
#> stringr 1.4.0 2019-02-10 [2] CRAN (R 3.6.0)
#> TeachingDemos 2.10 2016-02-12 [2] CRAN (R 3.6.0)
#> testthat 2.1.1 2019-04-23 [2] CRAN (R 3.6.0)
#> tibble 2.1.3 2019-06-06 [2] CRAN (R 3.6.0)
#> tidyselect 0.2.5 2018-10-11 [2] CRAN (R 3.6.0)
#> usethis 1.5.0 2019-04-07 [2] CRAN (R 3.6.0)
#> withr 2.1.2 2018-03-15 [2] CRAN (R 3.6.0)
#> xfun 0.7 2019-05-14 [2] CRAN (R 3.6.0)
#> yaml 2.2.0 2018-07-25 [2] CRAN (R 3.6.0)
#>
#> [1] C:/Users/HS/AppData/Local/Temp/Rtmp4qQvuK/Rinst1f983ebf72c0
#> [2] D:/Program Files/R/R-3.6.0/library
```

Schütz H, Tomashevskiy M, Labes D, Shitova A, González-de la Parra M, Fuglsang A.

*Reference Datasets for Studies in a Replicate Design intended for Average Bioequivalence with Expanding Limits.*Manuscript in preparation 2019.↩Gulf Cooperation Council: Bahrain, Kuwait, Oman, Qatar, Saudi Arabia, and the United Arab Emirates.↩

Confounded effects (design

*not recommended*).↩Confounded effects (design

*not recommended*).↩Balaam’s design (

*not recommended*due to its poor power characteristics).↩Extra-reference design; biased in the presence of period effects (design

*not recommended*).↩Napierian logarithm (base

*e*). The decadic logarithm (base 10) is not supported).↩Contradics the law of simplicity. Such a nesting is superfluous since in BE trials subjects are uniquely coded. If, say, subject

`1`

is allocated to sequence`TRTR`

there is not yet ‘another’ subject`1`

allocated to sequence`RTRT`

. This explains the many lines in SAS`PROC GML`

given with`.`

and in Phoenix/WinNonlin as`not estimable`

.↩Guideline on the pharmacokinetic and clinical evaluation of modified release dosage forms↩

Of course, only if the point estimates are identical.↩

The fences are given by the lowest datum still within

*m*×IQR of the lower quartile, and the highest datum still within*m*×IQR of the upper quartile, where IQR is the interquartile range (difference between the 3^{rd}and 1^{st}quartiles). Data outside fences are considered outliers. Decreasing the multiplier*m*to*e.g.*, 1.5 might result in many outliers whereas increasing the multiplier in only a few.

Different methods exist to calculate quartiles (nine are available in R, where the default is`type = 7`

). R’s default is used by S, MATLAB, and Excel. Phoenix/WinNonlin, Minitab, and SPSS use`type = 6`

, ScyPy uses`type = 4`

, whereas the default in SAS and Stata is`type = 2`

(though others are available as well).↩*Externally*studentized: \(\widehat{\sigma}_{(i)}^2={1 \over n-m-1}\sum_{\begin{smallmatrix}j = 1\\j \ne i\end{smallmatrix}}^n \widehat{\varepsilon\,}_j^{\,2}\)↩*Internally*studentized: \(\widehat{\sigma}^2={1 \over n-m}\sum_{j=1}^n \widehat{\varepsilon\,}_j^{\,2}\)↩Both are available in SAS and R, whereas only the latter in

*e.g.*, Phoenix/WinNonlin. In general the former are slightly more restrictive. Which one will be used has to be stated in the SAP.↩