Bioequivalence Analysis (BE)

Bioequivalence.jl is a package for performing bioequivalence analysis.

The full API is available in the next section and provides the signatures and examples for using all the available functionality.

Quickstart

In order to use Bioequivalence.jl, load the package through:

using Bioequivalence

A bioequivalence study is an instance of the type BioequivalenceStudy and can be constructed through the pumas_be function. Alternatively, one can use preprocess_be and run_be to first prepare the data and run the analysis in separate steps.

What is bioequivalence?

Bioequivalence is a concept in pharmacokinetics that captures the idea that various pharmaceutical products administrated in a similar manner (e.g., same molar dose of the same active ingredient, route of administration) can be expected to have, for all intents and purposes, the same effect on individuals in a defined population.

Clinical studies collect data which can be analyzed such as through noncompartmental analysis to obtain insightful descriptives about the contentration curve also known as pharmacokinetic endpoints. These endpoints relate to the rate (e.g., maximum concentration, time of peak concentration) and extent of absorption (e.g., area under the curve). Bioequivalence relies on the study design and pharmacokinetic endpoints from clinical trials or simulation models to make a determination about the expected effects of formulations.

Three major types of bioequivalence studies are regularly used:

  1. Average (ABE): are the mean values of the distributions of the pharmacokinetic endpoints for the reference and the test formulations similar enough? The concept is the most popular with a rise in adoption in the early 1990's by the United States and the European Union. It is considered to be the easiest criterion for a new formulation to achieve bioequivalence. It is required by most regulatories agencies for the product to be approved under a bioequivalence process.
  2. Population (PBE): are the distributions of the pharmacokinetic endpoints for the reference and the test formulations similar enough? In this case, it is not longer comparing just the expected value of the distributions but the full distribution. PBE is especially important for determining prescribability or the decision to assign a patient one of the formulations as part of a treatment for the first time.
  3. Individual (IBE): are the distributions of the pharmacokinetic endpoints for the reference and the test formulations similar enough across a large proportion of the intended population? IBE is particularly relevant for switchability or the decision to substitute an ongoing regimen (change formulation) without detrimental effects to the patient.

PBE and IBE can be assessed through three different methods:

  1. Constant scaling: the regulatory agency provides a value to be used in determining PBE or IBE.
  2. Reference scaling: the estimated total variance of the reference formulation is used for determining PBE or IBE.
  3. Mixed scaling: use reference scaling when the estimated total variance of the reference formulation is greater than that of the test formulation and the constant scaling otherwise.
Note

One argument for using the mixed scaling is that if the estimate of the total variance of the test formulation is greater than that of the reference, then bioequivalence would be very conservative.

d55dd31 (Typos fixes on BE docs)

Info

The reference scaling method is mostly used when working with highly variable drugs (HVD), those with intrasubject variability > 30%, and narrow therapeutic index drugs (NTI), those drugs where small differences in dose or blood concentration may lead to serious therapeutic failures and/or adverse drug reactions that are life-threatening or result in persistent or significant disability or incapacity.

Designs

There are three major categories of bioequivalence study designs.

Nonparametric for endpoints such as time of maximum concentration which typically do not have (nor can easily transformed) a normal-like distribution.

Parallel designs which are typically used when crossover designs are not feasible.

Crossover (replicated and nonreplicated) designs which are the most commonly used by the industry.

Designs are fully characterized by:

  • Subjects: participants in the study (each is assigned to a sequence)
  • Formulations: the different formulations being compared (i.e., reference and additional test formulations)
  • Periods: each dosing period at which each subject is administrated a formulation based on the sequence it has been assigned to
  • Sequences: a dosing regimen which establishes what formulation is given at each period
Warning

The periods should be spaced enough such that there are no carryover effects from dosings in the previous periods.

Nonparametric analysis (i.e., x formulations, y sequences, z periods)

With a nonparametric design, the bioequivalence analysis is to perform a a Wilcoxon signed rank test of the null hypothesis that the distribution of the reference formulation and the distribution of an alternative formulation have the same median.

When there are no tied ranks and ≤ 50 samples, or tied ranks and ≤ 15 samples, it will perform an exact signed rank test. Otherwise an approximate signed rank test is carried out.

Since the study design can be inferred from the data argument (i.e., based on sequences, formulations, and periods), the inferred study design approach will be automatically selected. One can manually overwrite the method for the nonparametric option by selecting nonparametric = true in pumas_be.

Parallel design (i.e., x formulations, y periods, z sequences, e.g. R|S|T)

Perform a Welch's t-test (i.e., unequal variance two-sample t-test) of the null hypothesis that the distribution of the reference formulation and the distribution of an alternative formulation have equal means. The value of degrees of freedom of the test is determined via the Welch-Satterthwaite formula:

\[ ν_{χ'} ≈ \frac{\left(\sum_{i=1}^n k_i s_i^2\right)^2}{\sum_{i=1}^n \frac{(k_i s_i^2)^2}{ν_i}}\]

Crossover designs

Crossover designs are divided into two categories:

  1. Nonreplicated
  2. Replicated

Replicated crossover designs are those with subjects receiving the same formulation multiple times. A key feature of replicated designs is that it allows to estimate within-subject variances per formulation which is key for highly variable drugs and narrow therapeutic index drugs.

Analysis of crossover designs

Nonreplicated crossover designs are analyzed through a linear model

\[\log\left(endpoint\right) = β₀ + β₁ formulation + β₂ period + β₃ id + ε\]

where $β_j, j ∈ \{1, 2, 3\}$, are vectors for features where formulation and id use indicators and period uses contrast coding.

Replicated crossover designs are by default fitted with a linear mixed model that allows for separate variance components for each formulation (homogeneity = false):

log(endpoint) ~ formulation + sequence + period + (formulation + 0 | id) + zerocorr(formulation + 0 | row)

which is preferred by the U.S. Food and Drug Administration (FDA). The European Medicines Agency's (EMA) labels this model "Method C". EMA's preferred model for replicated crossover design (labeled "Method B") assumes equal variances for all formulations (homogeneity = true), i.e.

log(endpoint) ~ formulation + sequence + period + (1 | id)
Tip

The default estimator uses restricted maximum likelihood (REML) as the objective. One can request the maximum likelihood objective function through passing reml = false argument to pumas_be/run_be.

Supported designs

DescriptionTreatmentsPeriodsSequencesReplicatedCrossover
R|T212NoNo
RT|TR222NoYes
RR|RT|TR|TT224NoYes
RTR|TRT232FullyYes
RTR|TRR232PartiallyYes
RTT|TRR232FullyYes
RRT|RTR|TRR233PartiallyYes
RTRT|TRTR242FullyYes
RRTT|TTRR242FullyYes
RTTR|TRRT242FullyYes
RRTT|RTTR|TRRT|TTRR244FullyYes
RTRT|RTTR|TRRT|TRTR244FullyYes
RR|TT2>12FullyNo
RST|RTS|SRT|STR|TRS|TSR336NoYes
ADBC|BACD|CBDA|DCAB444NoYes
Note

In the United States, there is a minimum requirement of at least 12 evaluable subjects for any bioequivalence study.

Validation

Each design has been tested using various sources including:

  • Chow, Shein-Chung, and Jen-pei Liu. 2009. Design and Analysis of Bioavailability and Bioequivalence Studies. 3rd ed. Chapman & Hall/CRC Biostatistics Series 27. Boca Raton: CRC Press. DOI: 10.1201/9781420011678.
  • Fuglsang, Anders, Helmut Schütz, and Detlew Labes. 2015. "Reference Datasets for Bioequivalence Trials in a Two-Group Parallel Design." The AAPS Journal 17 (2): 400–404. DOI: 10.1208/s12248-014-9704-6.
  • Patterson, Scott D, and Byron Jones. 2017. Bioequivalence and Statistics in Clinical Pharmacology. 2nd ed. Chapman & Hall/CRC Biostatistics Series. DOI: 10.1201/9781315374161.
  • Schütz, Helmut, Detlew Labes, and Anders Fuglsang. 2014. "Reference Datasets for 2-Treatment, 2-Sequence, 2-Period Bioequivalence Studies." The AAPS Journal 16 (6): 1292–97. DOI: 10.1208/s12248-014-9661-0.
  • Schütz H, Tomashevskiy M, Labes D, Shitova A, González-de la Parra M, Fuglsang A. 2020. Reference Datasets for Studies in a Replicate Design Intended for Average Bioequivalence with Expanding Limits. AAPS J. 22(2): Article 44. DOI: 10.1208/s12248-020-0427-6.

API

Public

Bioequivalence.EMA_CMAXConstant
EMA_CMAX = [GeometricMeanRatioBounds(), AverageBioequivalenceWithExpandingLimits(0.76, CV_min = 0.3, CV_max = 0.5)]

Criteria for the European Medicines Agency average bioequivalence for Cmax endpoints.

Bioequivalence.EMA_NarrowTherapeuticIndexConstant
EMA_NarrowTherapeuticIndex = [AverageBioequivalenceWithExpandingLimits(0, lower_bound = 0.9, upper_bound = 1.11)]

Criteria for the European Medicines Agency narrow therapeutic index drugs.

Bioequivalence.FDAConstant
FDA = [GeometricMeanRatioBounds(), AverageBioequivalenceWithExpandingLimits((log(1.25) / 0.25)^2)]

Criteria for the U.S. Food and Drug Administration (non-narrow therapeutic index drugs):

If the within subject variability of the reference formulation is less than 30%, the unscaled confidence interval should pass the standard critertia: LB ≥ 80% ^ UB ≤ 125%. Otherwise, the 95% upper confidence bound for (T̄ₜ - Ȳᵣ)² - 𝜃 * σwᵣ² ≤ 0 (numbers should be kept to a minimum of four significant figures for comparison) with 𝜃 ≈ 0.7967 and 80% ≤ GMR ≤ 125%.

Bioequivalence.FDA_NarrowTherapeuticIndexConstant
FDA_NarrowTherapeuticIndex = [
    AverageBioequivalenceWithExpandingLimits(0),
    AverageBioequivalenceWithExpandingLimits((log(1 / 0.9) / 0.1)^2, CV_min = 0, CV_max = 0.2142),
    WithinSubjectVariabilityRatio()
    ]

Criteria for the U.S. Food and Drug Administration Narrow Therapeutic Index drug:

a) The 95% upper confidence bound for (T̄ₜ - Ȳᵣ)² - 𝜃 * σwᵣ² ≤ 0 (numbers should be kept to a minimum of four significant figures for comparison). 𝜃 ≈ 1.11.

b) Regular unscaled bioequivalence limits of [80.00%, 125.00%] should be passed.

c) The proposed requirement for the upper limit of the 90% equal-tails confidence interval for σwₜ / σwᵣ ≤ 2.5.

Bioequivalence.AverageBioequivalenceWithExpandingLimitsType
AverageBioequivalenceWithExpandingLimits(
    θ::Real = -1;
    CV_min::Real = 0.3,
    CV_max::Real = Inf,
    lower_bound::Real = 0.8,
    upper_bound::Real = 1.25
    )

Criteria for the average bioequivalence acceptable bounds. If θ < 0, use the θ from the BioequivalenceStudy if available. If θ > 0, overwrite BioequivalenceStudy parameters in assessment. If θ == 0, use unscaled confidence interval with specified interval.

See also assess_be

Bioequivalence.BioequivalenceStudyType
BioequivalenceStudy

Return a bioequivalence study.

See also: pumas_be.

Fields

  • data::DataFrame data used for the study
  • data_stats::NamedTuple
    • total::Int refers to the number of observations the data passed to the function had.
    • used_for_analysis::Int refers to the number of observations used for fitting the model (e.g., drop missing values)
    • treatment::DataFrame gives a DataFrame with the summary statistics of the statistical model's response by treatment
    • sequence::DataFrame gives a DataFrame with the summary statistics of the statistical model's response by sequence
    • period::DataFrame gives a DataFrame with the summary statistics of the statistical model's response by period
  • design::NamedTuple number of subjects in each sequence
  • model statistical models used for the analysis
  • model_stats statistics for the model
  • result::DataFrame results for inference

Examples

julia> data = dataset(joinpath("bioequivalence", "RST_RTS_SRT_STR_TRS_TSR", "PJ2017_4_5"))
186×5 DataFrame
 Row │ id     sequence  period  AUC      Cmax
     │ Int64  String3   Int64   Int64?   Int64?
─────┼───────────────────────────────────────────
   1 │     1  SRT            1     7260     1633
   2 │     1  SRT            2     6463     1366
   3 │     1  SRT            3     8759     2141
   4 │     2  RTS            1     3457      776
   5 │     2  RTS            2     6556     2387
   6 │     2  RTS            3     4081     1355
   7 │     4  TSR            1     4006     1326
   8 │     4  TSR            2     4879     1028
   9 │     4  TSR            3     3817     1052
  10 │     5  STR            1     4250      945
  11 │     5  STR            2     3487     1041
  ⋮  │   ⋮       ⋮        ⋮        ⋮        ⋮
 177 │    61  RTS            3     3779     1144
 178 │    62  SRT            1     5787     1461
 179 │    62  SRT            2     7069     1995
 180 │    62  SRT            3     6530     1236
 181 │    63  TRS            1     2204      495
 182 │    63  TRS            2     2927      770
 183 │    63  TRS            3  missing  missing
 184 │    67  RST            1     4045     1025
 185 │    67  RST            2     7865     2668
 186 │    67  RST            3  missing  missing
                                 165 rows omitted

julia> output = pumas_be(data, endpoint = :Cmax)
Design: RST|RTS|SRT|STR|TRS|TSR

Sequences: RST|RTS|SRT|STR|TRS|TSR (6)
Periods: 1:3 (3)
Subjects per Sequence: (RST = 9, RTS = 11, SRT = 11, STR = 10, TRS = 11, TSR = 10)

Average Bioequivalence
───────────────────────────────────────────────────────────────────────────────────
              δ         SE      lnLB      lnUB      GMR       LB       UB        CV
───────────────────────────────────────────────────────────────────────────────────
S - R  0.466225  0.0525563  0.379094  0.553357  1.59397  1.46096  1.73908  0.296903
T - R  0.261828  0.0525354  0.174731  0.348925  1.2993   1.19093  1.41754  0.296903
───────────────────────────────────────────────────────────────────────────────────

julia> output.data_stats.treatment
3×10 DataFrame
 Row │ formulation  exp_mean  mean     std       min      q25      median   q75      max      n
     │ Cat…         Float64   Float64  Float64   Float64  Float64  Float64  Float64  Float64  Int64
─────┼──────────────────────────────────────────────────────────────────────────────────────────────
   1 │ R             837.478  6.7304   0.466938  5.89715  6.3257   6.66568  7.09589  7.71334     62
   2 │ S            1339.96   7.20039  0.419893  6.09131  6.92952  7.20117  7.5251   8.02027     61
   3 │ T            1078.08   6.98294  0.473224  5.75574  6.62539  7.00851  7.29641  7.90286     61

julia> output.data_stats.sequence
6×10 DataFrame
 Row │ sequence  exp_mean  mean     std       min      q25      median   q75      max      n
     │ Cat…      Float64   Float64  Float64   Float64  Float64  Float64  Float64  Float64  Int64
─────┼───────────────────────────────────────────────────────────────────────────────────────────
   1 │ RST        997.282  6.90503  0.493722  5.90263  6.6385   6.92755  7.17642  7.88908     26
   2 │ RTS       1084.03   6.98844  0.499403  6.09131  6.4677   7.05618  7.32449  7.77779     33
   3 │ SRT       1187.43   7.07954  0.482725  5.89715  6.63068  7.11964  7.45124  7.90286     33
   4 │ STR        869.299  6.76769  0.404857  6.04501  6.4758   6.8663   6.95607  7.71913     30
   5 │ TRS       1180.5    7.07369  0.466733  6.20456  6.71254  7.11698  7.39368  7.96797     32
   6 │ TSR       1071.51   6.97682  0.556906  5.75574  6.69448  7.02452  7.28049  8.02027     30

julia> output.data_stats.period
3×10 DataFrame
 Row │ period  exp_mean  mean     std       min      q25      median   q75      max      n
     │ Cat…    Float64   Float64  Float64   Float64  Float64  Float64  Float64  Float64  Int64
─────┼─────────────────────────────────────────────────────────────────────────────────────────
   1 │ 1        1009.89  6.9176   0.498246  5.75574  6.46653  6.97018  7.28186  7.72356     62
   2 │ 2        1108.56  7.01082  0.460876  5.89715  6.63035  7.06641  7.31235  8.02027     62
   3 │ 3        1076.82  6.98176  0.516518  6.03787  6.63167  6.99805  7.30986  7.96797     60

julia> output.model
StatsModels.TableRegressionModel{GLM.LinearModel{GLM.LmResp{Vector{Float64}}, GLM.DensePredChol{Float64, LinearAlgebra.CholeskyPivoted{Float64, Matrix{Float64}, Vector{Int64}}}}, Matrix{Float64}}

endpoint ~ 1 + formulation + period + id

Coefficients:
─────────────────────────────────────────────────────────────────────────────────
                      Coef.  Std. Error      t  Pr(>|t|)    Lower 95%   Upper 95%
─────────────────────────────────────────────────────────────────────────────────
(Intercept)      7.18626      0.170519   42.14    <1e-72   6.84859      7.52394
formulation: S   0.466225     0.0525563   8.87    <1e-14   0.362149     0.570301
formulation: T   0.261828     0.0525354   4.98    <1e-05   0.157794     0.365862
period: 2        0.0504953    0.0302861   1.67    0.0981  -0.00947942   0.11047
period: 3        0.00727282   0.0306504   0.24    0.8128  -0.0534233    0.067969
id: 2           -0.214448     0.237319   -0.90    0.3680  -0.684403     0.255508
id: 4           -0.401034     0.237319   -1.69    0.0937  -0.87099      0.0689214
id: 5           -0.606075     0.237319   -2.55    0.0119  -1.07603     -0.136119
id: 6           -0.356888     0.237319   -1.50    0.1353  -0.826843     0.113068
id: 7           -0.349386     0.237319   -1.47    0.1436  -0.819342     0.12057
id: 8           -0.850534     0.237319   -3.58    0.0005  -1.32049     -0.380578
id: 9           -0.237702     0.237319   -1.00    0.3186  -0.707658     0.232254
id: 10          -0.589593     0.237319   -2.48    0.0144  -1.05955     -0.119637
id: 11          -0.109224     0.237319   -0.46    0.6462  -0.57918      0.360732
id: 12          -0.0629286    0.237319   -0.27    0.7913  -0.532884     0.407027
id: 13          -0.829916     0.237319   -3.50    0.0007  -1.29987     -0.35996
id: 14          -0.169207     0.237319   -0.71    0.4773  -0.639163     0.300749
id: 15          -0.85675      0.237319   -3.61    0.0005  -1.32671     -0.386794
id: 16          -0.385301     0.237319   -1.62    0.1071  -0.855257     0.084655
id: 17           0.069506     0.237319    0.29    0.7701  -0.40045      0.539462
id: 18          -0.490339     0.237319   -2.07    0.0410  -0.960295    -0.0203832
id: 19          -0.440248     0.237319   -1.86    0.0661  -0.910203     0.0297083
id: 20           0.292819     0.237319    1.23    0.2197  -0.177137     0.762775
id: 21          -1.02247      0.237319   -4.31    <1e-04  -1.49242     -0.552511
id: 22          -0.975804     0.237319   -4.11    <1e-04  -1.44576     -0.505848
id: 23          -1.14843      0.237319   -4.84    <1e-05  -1.61839     -0.678477
id: 24          -0.620657     0.237319   -2.62    0.0101  -1.09061     -0.150702
id: 25           0.241019     0.237319    1.02    0.3119  -0.228937     0.710974
id: 26          -0.768347     0.237319   -3.24    0.0016  -1.2383      -0.298391
id: 27          -0.837325     0.237319   -3.53    0.0006  -1.30728     -0.367369
id: 28          -0.726773     0.237319   -3.06    0.0027  -1.19673     -0.256817
id: 29          -0.0892112    0.237319   -0.38    0.7077  -0.559167     0.380745
id: 30          -1.01968      0.237319   -4.30    <1e-04  -1.48963     -0.549721
id: 31          -0.425466     0.237319   -1.79    0.0756  -0.895421     0.0444903
id: 32          -0.40831      0.237319   -1.72    0.0880  -0.878266     0.0616457
id: 33          -0.410798     0.237319   -1.73    0.0861  -0.880754     0.0591575
id: 34          -0.817252     0.237319   -3.44    0.0008  -1.28721     -0.347296
id: 35          -0.500386     0.237319   -2.11    0.0371  -0.970342    -0.0304301
id: 36           0.152288     0.237319    0.64    0.5223  -0.317668     0.622244
id: 37          -0.750648     0.237319   -3.16    0.0020  -1.2206      -0.280692
id: 39          -1.38512      0.237319   -5.84    <1e-07  -1.85507     -0.915161
id: 40          -1.17263      0.237319   -4.94    <1e-05  -1.64259     -0.702674
id: 41           0.0824414    0.237319    0.35    0.7289  -0.387514     0.552397
id: 42          -0.305264     0.237319   -1.29    0.2009  -0.77522      0.164691
id: 43          -0.59991      0.237319   -2.53    0.0128  -1.06987     -0.129954
id: 44          -0.841676     0.237319   -3.55    0.0006  -1.31163     -0.371721
id: 45          -0.66548      0.237319   -2.80    0.0059  -1.13544     -0.195524
id: 46          -0.864623     0.237319   -3.64    0.0004  -1.33458     -0.394667
id: 47          -0.283627     0.237319   -1.20    0.2344  -0.753583     0.186329
id: 48          -0.205965     0.237319   -0.87    0.3872  -0.675921     0.26399
id: 49           0.0450347    0.237319    0.19    0.8498  -0.424921     0.514991
id: 50          -0.514194     0.237319   -2.17    0.0323  -0.98415     -0.0442386
id: 51          -0.176366     0.237319   -0.74    0.4589  -0.646322     0.29359
id: 52          -0.99198      0.237319   -4.18    <1e-04  -1.46194     -0.522024
id: 53          -0.494151     0.237319   -2.08    0.0395  -0.964107    -0.0241951
id: 54          -0.0659794    0.237319   -0.28    0.7815  -0.535935     0.403976
id: 55          -0.0876202    0.237319   -0.37    0.7126  -0.557576     0.382336
id: 56          -0.180088     0.237319   -0.76    0.4495  -0.650043     0.289868
id: 57          -0.402746     0.237319   -1.70    0.0923  -0.872702     0.06721
id: 58           0.0550554    0.237319    0.23    0.8169  -0.4149       0.525011
id: 59          -0.178369     0.237319   -0.75    0.4538  -0.648325     0.291587
id: 60          -0.70748      0.237319   -2.98    0.0035  -1.17744     -0.237524
id: 61          -0.721606     0.237319   -3.04    0.0029  -1.19156     -0.25165
id: 62          -0.0939777    0.237319   -0.40    0.6928  -0.563934     0.375978
id: 63          -0.888067     0.266187   -3.34    0.0011  -1.41519     -0.360944
id: 67          -0.00497377   0.266231   -0.02    0.9851  -0.532183     0.522235
─────────────────────────────────────────────────────────────────────────────────
Bioequivalence.GeometricMeanRatioBoundsType
GeometricMeanRatioBounds(lower_bound::Real = 0.8, upper_bound::Real = 1.25)

Criteria for the point estimate. Checks that the geometric mean ratio (GMR) lies within the accepted range.

Bioequivalence.ReferenceScaledAverageBioequivalanceType
ReferenceScaledAverageBioequivalance(
    data::AbstractDataFrame,
    𝜃::Real,
    σwᵣ::Real = -1;
    k::Union{Real, Missing} = -1,
    level::Real = 0.9,
    level_y::Real = 0.95,
    userepeatedobsonly::Bool = true
    ) -> ReferenceScaledAverageBioequivalance

Reference-scaled average bioequivalence.

Arguments

  • data::AbstractDataFrame: dataset adhering to the schema from preprocess_be.
  • 𝜃::Real: parameter for reference-scaling.
  • σwᵣ::Real: estimate for the within-subject variability. If negative, the function will compute it.
  • k::Real: degrees of freedom of the within-subject variability. If negative, the function will compute it.
  • level::Real: confidence level for the confidence interval.
  • level_y::Real: confidence level for the upper bound of (Ȳₜ - Ȳᵣ)² - 𝜃 * σwᵣ².
  • userepeatedobsonly: controls the behavior of within_subject_variability if used to estimate σwᵣ and k.

See also preprocess_be and within_subject_variability

Reference

Food and Drug Administration (2021). Bioequivalence Studies With Pharmacokinetic Endpoints for Drugs Submitted Under an ANDA Guidance for Industry. https://www.fda.gov/media/87219/download.

Examples

julia> data = dataset(joinpath("bioequivalence", "RTTR_TRRT", "SLTGSF2020_DS16"));

julia> pkdata = preprocess_be(data, endpoint = :PK);

julia> rsabe = ReferenceScaledAverageBioequivalance(pkdata, 0.76)
Critical boundary: -0.0416
Regulatory parameter: 0.76
95.0% upper confidence bound with 36.0 degrees of freedom
Bioequivalence.WithinSubjectVariabilityRatioType
WithinSubjectVariabilityRatio(upper_bound::Real = 2.5)

Requirement for the upper limit of the 90% equal-tails confidence interval for within-subject variability test/reference ratio is less than or equal to the upper limit.

Bioequivalence.assess_beMethod
assess_be(criteria::AbstractVector{<:BioequivalenceCriterion}, result::BioequivalenceStudy)

Return the average bioequivalence conclusions based on the specified criteria.

Examples

julia> data = dataset(joinpath("bioequivalence", "RTTR_TRRT", "SLTGSF2020_DS16"))
152×4 DataFrame
 Row │ id     sequence  period  PK
     │ Int64  String7   Int64   Float64
─────┼──────────────────────────────────
   1 │     1  RTTR           1   0.2813
   2 │     1  RTTR           2   0.2947
   3 │     1  RTTR           3   0.5471
   4 │     1  RTTR           4   0.651
   5 │     2  RTTR           1   0.2024
   6 │     2  RTTR           2   0.1782
   7 │     2  RTTR           3   0.2076
   8 │     2  RTTR           4   0.3604
   9 │     3  TRRT           1   0.4332
  10 │     3  TRRT           2   0.3131
  11 │     3  TRRT           3   0.2132
  ⋮  │   ⋮       ⋮        ⋮        ⋮
 143 │    38  RTTR           3   2.6522
 144 │    38  RTTR           4   1.2808
 145 │    39  RTTR           1   0.3847
 146 │    39  RTTR           2   0.2991
 147 │    39  RTTR           3   0.3353
 148 │    39  RTTR           4   0.3483
 149 │    40  RTTR           1   1.0289
 150 │    40  RTTR           2   0.753
 151 │    40  RTTR           3   0.7894
 152 │    40  RTTR           4   0.6129
                        131 rows omitted

julia> output = pumas_be(data, endpoint = :PK)
Design: RTTR|TRRT

Sequences: RTTR|TRRT (2)
Periods: 1:4 (4)
Subjects per Sequence: (RTTR = 20, TRRT = 18)
Reference scaled using 𝜃 = 0.797

Average Bioequivalence
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
             δ       SE     lnLB     lnUB     GMR      LB      UB     CVᵣ     CVₜ  σ_ratio     σ⁺        cb    dof
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
T - R  -0.2378  0.07738  -0.3665  -0.1091  0.7883  0.6931  0.8966  0.4972  0.5141    1.031  1.361  -0.04805  86.56
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────

julia> assess_be(StandardBioequivalenceCriterion, output)
Formulation: T - R -- Fail
  LB ≥ 0.8% ^ UB ≤ 1.25% -- Fail

julia> assess_be(FDA, output)
Formulation: T - R -- Fail
  0.8 ≤ GMR ≤ 1.25 -- Fail
  Reference Scaled Average Bioequivalance w/ 𝜃: 0.796689, CV_min: 0.3, CV_max: Inf -- Pass

julia> assess_be(FDA_NarrowTherapeuticIndex, output)
Formulation: T - R -- Fail
  LB ≥ 0.8% ^ UB ≤ 1.25% -- Fail
  Reference Scaled Average Bioequivalance w/ 𝜃: 1.110084, CV_min: 0.0, CV_max: 0.2142 -- Fail
  Upper bound of the within-subject variability ratio ≤ 2.5 -- Pass
Bioequivalence.detect_designFunction
detect_design(sequences::AbstractVector) ->
    NamedTuple{(:design, :replicated, :crossover), Tuple{String, Replicated, Bool}}
detect_design(
    data::AbstractDataFrame,
    sequence::Union{AbstractString, Symbol} = :sequence
    ) -> NamedTuple{(:design, :replicated, :crossover), Tuple{String, Replicated, Bool}}

Return the design class (Parallel, Nonreplicatedcrossover, Replicatedcrossover) and design.

Bioequivalence.generate_designMethod
generate_design(
    sequences::AbstractVector{<:AbstractString},
    amt::Union{Number,AbstractVector{<:Number}},
    subjects_per_sequence::Union{<:Integer,AbstractVector{<:Integer}},
    ) -> DataFrame

Returns a DataFrame with id, sequence, period, amt, evid, cmt, and time. It can be used to quickly set up data for Pumas, NCA, and Bioequivalence. In order to add covariates, use innerjoin to join the result of this function with another DataFrame with covariates.

Examples

julia> skeleton = generate_design(["RT", "TR"], [0, 50], 10)
40×8 DataFrame
 Row │ id     sequence  period  formulation  amt    time   evid   cmt
     │ Int64  Cat…      Int64   Char         Int64  Int64  Int64  Int64
─────┼──────────────────────────────────────────────────────────────────
   1 │     1  RT             1  R                0      0      4      1
   2 │     1  RT             2  T               50      0      4      1
   3 │     2  RT             1  R                0      0      4      1
   4 │     2  RT             2  T               50      0      4      1
   5 │     3  RT             1  R                0      0      4      1
   6 │     3  RT             2  T               50      0      4      1
   7 │     4  RT             1  R                0      0      4      1
   8 │     4  RT             2  T               50      0      4      1
   9 │     5  RT             1  R                0      0      4      1
  10 │     5  RT             2  T               50      0      4      1
  11 │     6  RT             1  R                0      0      4      1
  ⋮  │   ⋮       ⋮        ⋮          ⋮         ⋮      ⋮      ⋮      ⋮
  31 │    16  TR             1  T               50      0      4      1
  32 │    16  TR             2  R                0      0      4      1
  33 │    17  TR             1  T               50      0      4      1
  34 │    17  TR             2  R                0      0      4      1
  35 │    18  TR             1  T               50      0      4      1
  36 │    18  TR             2  R                0      0      4      1
  37 │    19  TR             1  T               50      0      4      1
  38 │    19  TR             2  R                0      0      4      1
  39 │    20  TR             1  T               50      0      4      1
  40 │    20  TR             2  R                0      0      4      1
                                                         19 rows omitted

julia> skeleton = generate_design(["RTRT", "TRTR"], [50, 75], [12, 10])
88×8 DataFrame
 Row │ id     sequence  period  formulation  amt    time   evid   cmt
     │ Int64  Cat…      Int64   Char         Int64  Int64  Int64  Int64
─────┼──────────────────────────────────────────────────────────────────
   1 │     1  RTRT           1  R               50      0      4      1
   2 │     1  RTRT           2  T               75      0      4      1
   3 │     1  RTRT           3  R               50      0      4      1
   4 │     1  RTRT           4  T               75      0      4      1
   5 │     2  RTRT           1  R               50      0      4      1
   6 │     2  RTRT           2  T               75      0      4      1
   7 │     2  RTRT           3  R               50      0      4      1
   8 │     2  RTRT           4  T               75      0      4      1
   9 │     3  RTRT           1  R               50      0      4      1
  10 │     3  RTRT           2  T               75      0      4      1
  11 │     3  RTRT           3  R               50      0      4      1
  ⋮  │   ⋮       ⋮        ⋮          ⋮         ⋮      ⋮      ⋮      ⋮
  79 │    20  TRTR           3  T               75      0      4      1
  80 │    20  TRTR           4  R               50      0      4      1
  81 │    21  TRTR           1  T               75      0      4      1
  82 │    21  TRTR           2  R               50      0      4      1
  83 │    21  TRTR           3  T               75      0      4      1
  84 │    21  TRTR           4  R               50      0      4      1
  85 │    22  TRTR           1  T               75      0      4      1
  86 │    22  TRTR           2  R               50      0      4      1
  87 │    22  TRTR           3  T               75      0      4      1
  88 │    22  TRTR           4  R               50      0      4      1
                                                         67 rows omitted
Bioequivalence.geocv2sigmaMethod
geocv2sigma(CV::Union{Real,Missing}) = √(log(1 + CV^2))

Transform the the coefficient of variation (CV) to the σ parameter of a log-normal distribution.

Examples

julia> geocv2sigma(0.30)
0.29356037920852396
Bioequivalence.preprocess_beMethod
preprocess_be(
    data::AbstractDataFrame,
    id::Union{AbstractString, Symbol} = :id,
    sequence::Union{AbstractString, Symbol} = :sequence,
    period::Union{AbstractString, Symbol} = :period,
    endpoint::Union{AbstractString, Symbol} = :AUC,
    logtransformed::Bool = false
    ) -> DataFrame

Return the standardized dataset with id, sequence, period, formulation, and endpoint. The standardized dataset:

  • selects only relevant variables
  • applies the natural log transformation to the endpoint if in natural scale
  • renames variables to their canonical name for the functions
  • drops missing observations
  • computes the formulation based on the sequence and period variables
  • converts the variables into the appropiate format for analysis (e.g., factors)

The sequence/formulation take values RT, RST, or ABCD based on alphabetical order and number of formulations.

Examples

julia> data = dataset(joinpath("bioequivalence", "RT_TR", "SLF2014_1"))
36×4 DataFrame
 Row │ id     sequence  period  AUC
     │ Int64  String3   Int64   Float64
─────┼──────────────────────────────────
   1 │     1  RT             1   181.09
   2 │     1  RT             2   210.14
   3 │     2  RT             1   114.48
   4 │     2  RT             2    98.72
   5 │     3  TR             1   225.95
   6 │     3  TR             2   241.09
   7 │     4  RT             1   176.91
   8 │     4  RT             2   186.65
   9 │     5  TR             1   147.01
  10 │     5  TR             2   139.56
  11 │     6  TR             1    97.53
  ⋮  │   ⋮       ⋮        ⋮        ⋮
  27 │    14  TR             1   179.96
  28 │    14  TR             2   181.09
  29 │    15  TR             1   173.86
  30 │    15  TR             2   206.66
  31 │    16  RT             1   144.0
  32 │    16  RT             2   143.25
  33 │    17  RT             1   185.1
  34 │    17  RT             2   192.22
  35 │    18  TR             1   117.99
  36 │    18  TR             2   125.5
                         15 rows omitted

julia> preprocess_be(data)
36×5 DataFrame
 Row │ id    sequence  period  formulation  endpoint
     │ Cat…  Cat…      Cat…    Cat…         Float64
─────┼───────────────────────────────────────────────
   1 │ 1     RT        1       R             5.19899
   2 │ 1     RT        2       T             5.34777
   3 │ 2     RT        1       R             4.7404
   4 │ 2     RT        2       T             4.59229
   5 │ 3     TR        1       T             5.42031
   6 │ 3     TR        2       R             5.48517
   7 │ 4     RT        1       R             5.17564
   8 │ 4     RT        2       T             5.22924
   9 │ 5     TR        1       T             4.9905
  10 │ 5     TR        2       R             4.93849
  11 │ 6     TR        1       T             4.58016
  ⋮  │  ⋮       ⋮        ⋮          ⋮          ⋮
  27 │ 14    TR        1       T             5.19273
  28 │ 14    TR        2       R             5.19899
  29 │ 15    TR        1       T             5.15825
  30 │ 15    TR        2       R             5.33107
  31 │ 16    RT        1       R             4.96981
  32 │ 16    RT        2       T             4.96459
  33 │ 17    RT        1       R             5.2209
  34 │ 17    RT        2       T             5.25864
  35 │ 18    TR        1       T             4.7706
  36 │ 18    TR        2       R             4.83231
                                      15 rows omitted
Bioequivalence.pumas_beMethod
pumas_be(
    data::AbstractDataFrame;
    endpoint::Union{AbstractString,Symbol} = "AUC",
    logtransformed::Bool = false,
    reference_scale::Real = (log(1.25) / 0.25)^2,
    cv_max::Real = Inf,
    id::Union{AbstractString,Symbol} = "id",
    sequence::Union{AbstractString,Symbol} = "sequence",
    period::Union{AbstractString,Symbol} = "period",
    nonparametric::Bool = occursin(r"(?i)tmax", string(endpoint)),
    homogeneity::Union{Bool, Nothing} = nothing,
    userepeatedobsonly::Bool = true,
    reml::Bool = true,
    level::Real = 0.9,
    alpha::Real = 0.1,
    level_y::Real = 0.95,
    sigdigits::Integer = 4
    ) -> BioequivalenceStudy

BioequivalenceStudy constructor. See also: BioequivalenceStudy, preprocess_be, and run_be.

Arguments

  • data: must have id, sequence, period, and an endpoint.
  • endpoint: which variable is the endpoint?
  • logtransformed: has the endpoint been log transformed?
  • id: which variable is the subject identifier?
  • sequence: which variable is the sequence?
  • period: which variable is the period?
  • nonparametric: whether to use a nonparametric (default if endpoint includes tmax ignoring case) or parametric model.
  • homogeneity: whether formulation groups should be modeled with equal variance
  • userepeatedobsonly: whether estimating the within subject variability should only repeated observations
  • reference_scale: 𝜃 for reference scale (e.g., FDA ≈ 0.797, FDA/NTI ≈ 1.11, EMA = 0.76)
  • cv_max: maximum within subject variability for reference scaling (FDA = Inf, FDA/NTI = 0.2142, EMA = 0.5)
  • reml: whether the linear mixed model should use restricted maximum likelihood or maximum likelihood
  • level: applies to the confidence intervals for the GMR
  • alpha: applies to the upper bound of the within subject variability ratio
  • level_y: applies to the critical boundary for reference-scaled average bioequivalence
  • sigdigits: results given with how many significant digits.

Current designs include: nonparametric, parallel, and various crossover designs | Description | Treatments | Periods | Sequences | Replicated | Crossover | |––––––––––––-|––––––|––––-|–––––-|––––––|–––––-| | R|T | 2 | 1 | 2 | No | No | | RT|TR | 2 | 2 | 2 | No | Yes | | RR|RT|TR|TT | 2 | 2 | 4 | No | Yes | | RTR|TRT | 2 | 3 | 2 | Fully | Yes | | RTR|TRR | 2 | 3 | 2 | Partially | Yes | | RTT|TRR | 2 | 3 | 2 | Fully | Yes | | RRT|RTR|TRR | 2 | 3 | 3 | Partially | Yes | | RTRT|TRTR | 2 | 4 | 2 | Fully | Yes | | RRTT|TTRR | 2 | 4 | 2 | Fully | Yes | | RTTR|TRRT | 2 | 4 | 2 | Fully | Yes | | RRTT|RTTR|TRRT|TTRR | 2 | 4 | 4 | Fully | Yes | | RTRT|RTTR|TRRT|TRTR | 2 | 4 | 4 | Fully | Yes | | RR|TT | 2 | >1 | 2 | Fully | No | | RST|RTS|SRT|STR|TRS|TSR | 3 | 3 | 6 | No | Yes | | ADBC|BACD|CBDA|DCAB | 4 | 4 | 4 | No | Yes |

Examples

julia> data = dataset(joinpath("bioequivalence", "RT_TR", "SLF2014_1"))
36×4 DataFrame
 Row │ id     sequence  period  AUC
     │ Int64  String3   Int64   Float64
─────┼──────────────────────────────────
   1 │     1  RT             1   181.09
   2 │     1  RT             2   210.14
   3 │     2  RT             1   114.48
   4 │     2  RT             2    98.72
   5 │     3  TR             1   225.95
   6 │     3  TR             2   241.09
   7 │     4  RT             1   176.91
   8 │     4  RT             2   186.65
   9 │     5  TR             1   147.01
  10 │     5  TR             2   139.56
  11 │     6  TR             1    97.53
  ⋮  │   ⋮       ⋮        ⋮        ⋮
  27 │    14  TR             1   179.96
  28 │    14  TR             2   181.09
  29 │    15  TR             1   173.86
  30 │    15  TR             2   206.66
  31 │    16  RT             1   144.0
  32 │    16  RT             2   143.25
  33 │    17  RT             1   185.1
  34 │    17  RT             2   192.22
  35 │    18  TR             1   117.99
  36 │    18  TR             2   125.5
                         15 rows omitted
julia> output = pumas_be(data)
Design: RT|TR

Sequences: RT|TR (2)
Periods: 1:2 (2)
Subjects per Sequence: (RT = 9, TR = 9)

Average Bioequivalence
─────────────────────────────────────────────────────────────────────────────────────────────
                δ        SE        lnLB         lnUB       GMR        LB        UB         CV
─────────────────────────────────────────────────────────────────────────────────────────────
T - R  -0.0503868  0.026658  -0.0969286  -0.00384499  0.950862  0.907621  0.996162  0.0801021
─────────────────────────────────────────────────────────────────────────────────────────────
julia> output = pumas_be(data, nonparametric = true)
Design: RT|TR

Sequences: RT|TR (2)
Periods: 1:2 (2)
Subjects per Sequence: (RT = 9, TR = 9)

Average Bioequivalence
────────────────────────────────────────
          lnLB       lnUB     LB      UB
────────────────────────────────────────
T - R  -0.1064  -0.001216  0.899  0.9988
────────────────────────────────────────

julia> data = dataset(joinpath("bioequivalence", "RTT_TRR", "PJ2017_4_1"))
285×5 DataFrame
 Row │ id     sequence  period  AUC       Cmax
     │ Int64  String3   Int64   Float64?  Float64?
─────┼─────────────────────────────────────────────
   1 │   101  TRR            1     12.26     0.511
   2 │   101  TRR            2     16.19     0.688
   3 │   101  TRR            3     11.34     0.533
   4 │   102  TRR            1    397.98    13.27
   5 │   102  TRR            2    267.63     7.933
   6 │   102  TRR            3    487.55    12.952
   7 │   103  TRR            1    243.81    16.771
   8 │   103  TRR            2    141.7      6.926
   9 │   103  TRR            3    198.44     9.257
  10 │   109  TRR            1    182.52     8.816
  11 │   109  TRR            2    112.34     4.921
  ⋮  │   ⋮       ⋮        ⋮        ⋮         ⋮
 276 │   186  RTT            3     87.63     4.87
 277 │   190  RTT            1     82.78     3.88
 278 │   190  RTT            2    164.56     7.37
 279 │   190  RTT            3    213.98     7.01
 280 │   191  RTT            1     98.86     4.59
 281 │   191  RTT            2     99.02     2.96
 282 │   191  RTT            3     75.48     2.38
 283 │   194  RTT            1     21.29     1.51
 284 │   194  RTT            2     46.3      2.74
 285 │   194  RTT            3     15.41     1.41
                                   264 rows omitted
julia> output = pumas_be(data)
Design: RTT|TRR

Sequences: RTT|TRR (2)
Periods: 1:3 (3)
Subjects per Sequence: (RTT = 46, TRR = 48)
Reference scaled using 𝜃 = 0.797

Average Bioequivalence
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
              δ       SE     lnLB     lnUB     GMR      LB     UB     CVᵣ     CVₜ  σ_ratio     σ⁺        cb    dof
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
T - R  -0.02654  0.06906  -0.1409  0.08777  0.9738  0.8686  1.092  0.4275  0.6965    1.536  1.983  -0.09135  148.7
──────────────────────────────────────────────────────────────────────────────────────────────────────────────────
Bioequivalence.reference_scaled_acceptance_boundsFunction
reference_scaled_acceptance_bounds(σwr::Real, θ::Real)
reference_scaled_acceptance_bounds(obj::BioequivalenceStudy)

Return the lower and upper bounds based on the reference within subject variability and regulatory parameter.

Reference: European Medicines Agency (2010). "GUIDELINE ON THE INVESTIGATION OF BIOEQUIVALENCE": 4.1.10 Highly variable drugs or drug products. Doc. Ref.: CPMP/EWP/QWP/1401/98 Rev. 1/ Corr **

Examples

julia> data = dataset(joinpath("bioequivalence", "RTTR_TRRT", "PJ2017_4_3"));

julia> output = pumas_be(data)
Design: RTTR|TRRT

Sequences: RTTR|TRRT (2)
Periods: 1:4 (4)
Subjects per Sequence: (RTTR = 8, TRRT = 9)
Reference scaled using 𝜃 = 0.797

Average Bioequivalence
───────────────────────────────────────────────────────────────────────────────────────────────────────────────────
             δ       SE       lnLB     lnUB    GMR      LB     UB      CVᵣ     CVₜ  σ_ratio     σ⁺        cb    dof
───────────────────────────────────────────────────────────────────────────────────────────────────────────────────
T - R  0.03568  0.02431  -0.006884  0.07826  1.036  0.9931  1.082  0.08022  0.1084     1.35  2.118  0.001829  15.25
───────────────────────────────────────────────────────────────────────────────────────────────────────────────────

julia> reference_scaled_acceptance_bounds(output)
(0.9381852023891446, 1.0658876280007832)

julia> reference_scaled_acceptance_bounds(geocv2sigma(30), 0.76)
(0.13774539477808087, 7.259770837428591)

julia> reference_scaled_acceptance_bounds(geocv2sigma(21.42), (log(1 / 0.9) / 0.1)^2)
(0.06401589407631314, 15.62112057371101)
Bioequivalence.run_beMethod
run_be(
    data::AbstractDataFrame;
    reference_scale::Real = (log(1.25) / 0.25)^2,
    cv_max::Real = Inf,
    nonparametric::Bool = false,
    homogeneity::Union{Bool, Nothing} = nothing,
    userepeatedobsonly::Bool = true,
    level::Real = 0.9,
    alpha::Real = 0.1,
    level_y::Real = 0.95,
    reml::Bool = true,
    sigdigits::Integer = 4
    ) -> BioequivalenceStudy

BioequivalenceStudy constructor.

See also: BioequivalenceStudy.

Arguments

  • data: must have id, sequence, period, and an endpoint.
  • nonparametric: whether to use a nonparametric analysis (uncommon and usually reserved for Tmax)
  • homogeneity: whether formulation groups should be modeled with equal variance
  • userepeatedobsonly: whether estimating the within subject variability should only repeated observations
  • reference_scale: 𝜃 for reference scale (FDA ≈ 0.797, FDA/NTI ≈ 1.11, EMA = 0.76)
  • cv_max: maximum within subject variability for reference scaling (FDA = Inf, FDA/NTI = 0.2142, EMA = 0.5)
  • reml: whether the linear mixed model should use restricted maximum likelihood or maximum likelihood.
  • level: applies to the confidence intervals for the GMR.
  • alpha: applies to the upper bound of the within subject variability ratio
  • level_y: applies to the critical boundary for reference-scaled average bioequivalence
  • sigdigits: results given with how many significant digits.

Current designs include: nonparametric, parallel, and various crossover designs

DescriptionTreatmentsPeriodsSequencesReplicatedCrossover
RT212No
RTTR222No
RRRTTRTT22
RTRTRT232Fully
RTRTRR232Partially
RTTTRR232Fully
RRTRTRTRR233
RTRTTRTR242Fully
RRTTTTRR242Fully
RTTRTRRT242Fully
RRTTRTTRTRRTTTRR24
RTRTRTTRTRRTTRTR24
RRTT2>12Fully
RSTRTSSRTSTRTRSTSR
ADBCBACDCBDADCAB44

Examples

julia> data = dataset(joinpath("bioequivalence", "RT_TR", "SLF2014_1"))
36×4 DataFrame
 Row │ id     sequence  period  AUC
     │ Int64  String3   Int64   Float64
─────┼──────────────────────────────────
   1 │     1  RT             1   181.09
   2 │     1  RT             2   210.14
   3 │     2  RT             1   114.48
   4 │     2  RT             2    98.72
   5 │     3  TR             1   225.95
   6 │     3  TR             2   241.09
   7 │     4  RT             1   176.91
   8 │     4  RT             2   186.65
   9 │     5  TR             1   147.01
  10 │     5  TR             2   139.56
  11 │     6  TR             1    97.53
  ⋮  │   ⋮       ⋮        ⋮        ⋮
  27 │    14  TR             1   179.96
  28 │    14  TR             2   181.09
  29 │    15  TR             1   173.86
  30 │    15  TR             2   206.66
  31 │    16  RT             1   144.0
  32 │    16  RT             2   143.25
  33 │    17  RT             1   185.1
  34 │    17  RT             2   192.22
  35 │    18  TR             1   117.99
  36 │    18  TR             2   125.5
                         15 rows omitted

julia> pkdata = preprocess_be(data)
36×5 DataFrame
 Row │ id    sequence  period  formulation  endpoint
     │ Cat…  Cat…      Cat…    Cat…         Float64
─────┼───────────────────────────────────────────────
   1 │ 1     RT        1       R             5.19899
   2 │ 1     RT        2       T             5.34777
   3 │ 2     RT        1       R             4.7404
   4 │ 2     RT        2       T             4.59229
   5 │ 3     TR        1       T             5.42031
   6 │ 3     TR        2       R             5.48517
   7 │ 4     RT        1       R             5.17564
   8 │ 4     RT        2       T             5.22924
   9 │ 5     TR        1       T             4.9905
  10 │ 5     TR        2       R             4.93849
  11 │ 6     TR        1       T             4.58016
  ⋮  │  ⋮       ⋮        ⋮          ⋮          ⋮
  27 │ 14    TR        1       T             5.19273
  28 │ 14    TR        2       R             5.19899
  29 │ 15    TR        1       T             5.15825
  30 │ 15    TR        2       R             5.33107
  31 │ 16    RT        1       R             4.96981
  32 │ 16    RT        2       T             4.96459
  33 │ 17    RT        1       R             5.2209
  34 │ 17    RT        2       T             5.25864
  35 │ 18    TR        1       T             4.7706
  36 │ 18    TR        2       R             4.83231
                                      15 rows omitted

julia> output = run_be(pkdata)
Design: RT|TR

Sequences: RT|TR (2)
Periods: 1:2 (2)
Subjects per Sequence: (RT = 9, TR = 9)

Average Bioequivalence
─────────────────────────────────────────────────────────────────────────────────────────────
                δ        SE        lnLB         lnUB       GMR        LB        UB         CV
─────────────────────────────────────────────────────────────────────────────────────────────
T - R  -0.0503868  0.026658  -0.0969286  -0.00384499  0.950862  0.907621  0.996162  0.0801021
─────────────────────────────────────────────────────────────────────────────────────────────
Bioequivalence.sigma2geocvMethod
sigma2geocv(σw::Union{Real,Missing}) = √(exp(σw^2) - 1)

Transform the σ parameter of a log-normal distribution to the coefficient of variation (CV).

Examples

julia> sigma2geocv(0.294)
0.3004689459216001
Bioequivalence.update_rsabe_theta!Method
update_rsabe_theta!(obj::ReferenceScaledAverageBioequivalance, 𝜃::Real) -> ReferenceScaledAverageBioequivalance

Modifies in-place the 𝜃 of the ReferenceScaledAverageBioequivalance and updates the associated critical boundary.

julia> data = dataset(joinpath("bioequivalence", "RTTR_TRRT", "SLTGSF2020_DS16"));

julia> pkdata = preprocess_be(data, endpoint = :PK);

julia> rsabe = ReferenceScaledAverageBioequivalance(pkdata, 0.76)
Critical boundary: -0.0416
Regulatory parameter: 0.76
95.0% upper confidence bound with 36.0 degrees of freedom

julia> update_rsabe_theta!(rsabe, 1.11)
Critical boundary: -0.1019
Regulatory parameter: 1.11
95.0% upper confidence bound with 36.0 degrees of freedom
Bioequivalence.within_subject_variabilityMethod
within_subject_variability(
    data::AbstractDataFrame;
    userepeatedobsonly::Bool = true
    allownonreplicated::Bool = false,
    homogeneity::Union{Bool,Nothing} = nothing,
)

Return the within-subject variability estimate and degrees of freedom. It uses a linear model of the log-transformed PK response with fixed effects for subject ID and period. If userepeatedobsonly, data used is a sample with only observations from repeated subject/treatment. If allownonreplicated, the function returns a NamedTuple of missings instead of throwing an error when the passed formulation is not replicated in any sequences. This option is only relevant for partially replicated designs.

The homogeneity argument specifies if equal variances are assumed between formulations. When set to true, multiple formulations are expected in data. When set to false, only a single formulation should be available in data. In this case, design must be replicated. The default value of homogeneity is false for replicated designs and true otherwise.

References:

Schütz H, Tomashevskiy M, Labes D, Shitova A, González-de la Parra M, Fuglsang A. 2020. Reference Datasets for Studies in a Replicate Design Intended for Average Bioequivalence with Expanding Limits. AAPS J. 22(2): Article 44. DOI: 10.1208/s12248-020-0427-6.

Examples

julia> data = dataset(joinpath("bioequivalence", "RTTR_TRRT", "PJ2017_4_3"));

julia> pkdata = preprocess_be(data);

julia> combine(groupby(pkdata, :formulation), within_subject_variability)
2×3 DataFrame
 Row │ formulation  σw         k
     │ Cat…         Float64    Float64
─────┼─────────────────────────────────
   1 │ R            0.0800952     15.0
   2 │ T            0.108075      14.0

julia> data = dataset(joinpath("bioequivalence", "RRT_RTR_TRR", "SLTGSF2020_DS02"));

julia> pkdata = preprocess_be(data, endpoint = :PK);

julia> combine(groupby(pkdata, :formulation), t -> within_subject_variability(t; allownonreplicated=true))
2×3 DataFrame
 Row │ formulation  σw              k
     │ Cat…         Float64?        Float64?
─────┼────────────────────────────────────────
   1 │ R                  0.111361       22.0
   2 │ T            missing         missing
Bioequivalence.within_subject_variability_ratioFunction
within_subject_variability_ratio(
    σ::AbstractVector{<:Union{Real, Missing}},
    k::AbstractVector{<:Union{Real, Missing}},
    level::Real = 0.95
    ) -> Vector{NamedTuple{(:σ_ratio, :σ⁺), ...}}

within_subject_variability_ratio(
    data::AbstractDataFrame,
    level::Real = 0.95
    ) -> DataFrame

Return the within-subject variability ratio and upper bounds. σ contains the estimated within subject variability. k contains the associated degrees of freedom. level determines the confidence level for the upper bound.

Examples

julia> data = dataset(joinpath("bioequivalence", "RTTR_TRRT", "PJ2017_4_3"));

julia> pkdata = preprocess_be(data);

julia> wsv_estimates = combine(groupby(pkdata, :formulation), within_subject_variability)
2×3 DataFrame
 Row │ formulation  σw         k
     │ Cat…         Float64    Float64
─────┼─────────────────────────────────
   1 │ R            0.0800952     15.0
   2 │ T            0.108075      14.0

julia> within_subject_variability_ratio(wsv_estimates[!, :σw], wsv_estimates[!, :k])
1-element Vector{NamedTuple{(:σ_ratio, :σ⁺), Tuple{Float64, Float64}}}:
 (σ_ratio = 1.3493284240899683, σ⁺ = 2.117630330557258)

julia> within_subject_variability_ratio(wsv_estimates)
2×5 DataFrame
 Row │ formulation  σw         k        σ_ratio        σ⁺
     │ Cat…         Float64    Float64  Float64?       Float64?
─────┼───────────────────────────────────────────────────────────────
   1 │ R            0.0800952     15.0  missing        missing
   2 │ T            0.108075      14.0        1.34933        2.11763