API

Public

Bioequivalence.Bioequivalence — Module

Bioequivalence.jl

This module offers a suite of routines for bioequivalence (BE) analysis.

Bioequivalence.BioequivalenceEndpointOutput — Type

BioequivalenceEndpointOutput

A bioequivalence study output object.

See also: pumas_be.

Fields

endpoint_name the name of the endpoint
data::DataFrame data used for the study
design::NamedTuple number of subjects in each sequence
model statistical models used for the analysis
model_stats statistics for the model
stat_tests p-values
result::DataFrame results for inference in a Dataframe format
result_table::OrderedDict{Symbol,Any} results in an OrderedDict format
input_parameters a copy of the input parameters used

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)

sequence ╲ period │  1   2   3
──────────────────┼───────────
RST               │  9   9   8
RTS               │ 11  11  11
SRT               │ 11  11  11
STR               │ 10  10  10
TRS               │ 11  11  10
TSR               │ 10  10  10

-----------------------------------------------------------------------------------------
Paradigm: 3 formulations                 
Model: Linear model                      
Criteria: Standard ABE                   
Endpoint: Cmax                           
Formulations: Reference(R), Test(S and T)
-----------------------------------------------------------------------------------------
                                             Results(Cmax)    Assessment   Criteria      
-----------------------------------------------------------------------------------------
R             Geometric Marginal Mean        832.4                                       
              Geometric Naive Mean           837.5                                       
-----------------------------------------------------------------------------------------
S             Geometric Marginal Mean        1327                                        
              Geometric Naive Mean           1340                                        
              Geometric Mean S/R Ratio (%)   159.4                                       
              Degrees of Freedom             118                                         
              90% Confidence Interval (%)    [146.1, 173.9]   Fail         CI ⊆ [80, 125]
-----------------------------------------------------------------------------------------
T             Geometric Marginal Mean        1082                                        
              Geometric Naive Mean           1078                                        
              Geometric Mean T/R Ratio (%)   129.9                                       
              Degrees of Freedom             118                                         
              90% Confidence Interval (%)    [119.1, 141.8]   Fail         CI ⊆ [80, 125]
-----------------------------------------------------------------------------------------
Variability   CV (%) | σ̂                     29.69 | 0.2907                              
-----------------------------------------------------------------------------------------

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

endpoint ~ 1 + sequence + formulation + period + sequence & id_within_sequence

Coefficients:
──────────────────────────────────────────────────────────────────────────────────────────────────────────────
                                             Coef.   Std. Error       t  Pr(>|t|)     Lower 95%      Upper 95%
──────────────────────────────────────────────────────────────────────────────────────────────────────────────
(Intercept)                              6.62478      0.0712859   92.93    <1e-99    6.48361       6.76595
sequence: RTS                            0.063206     0.0766071    0.83    0.4110   -0.0884968     0.214909
sequence: SRT                            0.154313     0.0766071    2.01    0.0462    0.00260986    0.306016
sequence: STR                           -0.157545     0.0782601   -2.01    0.0464   -0.312522     -0.00256894
sequence: TRS                            0.139309     0.0773718    1.80    0.0743   -0.0139084     0.292526
sequence: TSR                            0.0515892    0.0782601    0.66    0.5111   -0.103387      0.206566
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.108264     0.0522307    2.07    0.0404    0.00483251    0.211695
period: 3                                0.065041     0.0528643    1.23    0.2210   -0.0396448     0.169727
sequence: RST & id_within_sequence: 2   -0.3262       0.15878     -2.05    0.0421   -0.640626     -0.0117729
sequence: RTS & id_within_sequence: 2    0.331286     0.16         2.07    0.0406    0.0144416     0.64813
sequence: SRT & id_within_sequence: 2   -0.501131     0.16        -3.13    0.0022   -0.817975     -0.184287
sequence: STR & id_within_sequence: 2    0.0716688    0.159198     0.45    0.6534   -0.243587      0.386925
sequence: TRS & id_within_sequence: 2    0.301479     0.160375     1.88    0.0626   -0.0161081     0.619065
sequence: TSR & id_within_sequence: 2    0.102741     0.159198     0.65    0.5199   -0.212515      0.417997
sequence: RST & id_within_sequence: 3   -0.51875      0.15878     -3.27    0.0014   -0.833177     -0.204324
sequence: RTS & id_within_sequence: 3    0.0552092    0.16         0.35    0.7307   -0.261635      0.372053
sequence: SRT & id_within_sequence: 3   -0.507347     0.16        -3.17    0.0019   -0.824191     -0.190502
sequence: STR & id_within_sequence: 3    0.170922     0.159198     1.07    0.2852   -0.144334      0.486178
sequence: TRS & id_within_sequence: 3    0.433913     0.160375     2.71    0.0078    0.116326      0.7515
sequence: TSR & id_within_sequence: 3    0.28292      0.159198     1.78    0.0781   -0.0323359     0.598176
sequence: RST & id_within_sequence: 4   -0.223057     0.15878     -1.40    0.1627   -0.537483      0.0913701
sequence: RTS & id_within_sequence: 4    0.733329     0.16         4.58    <1e-04    0.416485      1.05017
sequence: SRT & id_within_sequence: 4   -0.271254     0.16        -1.70    0.0926   -0.588098      0.0455901
sequence: STR & id_within_sequence: 4   -0.314543     0.159198    -1.98    0.0505   -0.629799      0.00071311
sequence: TRS & id_within_sequence: 4   -0.0758404    0.160375    -0.47    0.6372   -0.393427      0.241746
sequence: TSR & id_within_sequence: 4   -0.696306     0.159198    -4.37    <1e-04   -1.01156      -0.38105
sequence: RST & id_within_sequence: 5    0.0929178    0.15878      0.59    0.5595   -0.221509      0.407344
sequence: RTS & id_within_sequence: 5   -0.396815     0.16        -2.48    0.0145   -0.713659     -0.0799707
sequence: SRT & id_within_sequence: 5    0.590422     0.16         3.69    0.0003    0.273578      0.907266
sequence: STR & id_within_sequence: 5   -0.358415     0.159198    -2.25    0.0262   -0.673671     -0.0431591
sequence: TRS & id_within_sequence: 5    0.275196     0.160375     1.72    0.0888   -0.0423908     0.592783
sequence: TSR & id_within_sequence: 5   -0.31622      0.159198    -1.99    0.0493   -0.631476     -0.000964371
sequence: RST & id_within_sequence: 6   -0.33796      0.15878     -2.13    0.0354   -0.652387     -0.0235336
sequence: RTS & id_within_sequence: 6    0.0150445    0.16         0.09    0.9252   -0.3018        0.331889
sequence: SRT & id_within_sequence: 6   -0.150983     0.16        -0.94    0.3473   -0.467827      0.165862
sequence: STR & id_within_sequence: 6   -0.155991     0.159198    -0.98    0.3292   -0.471246      0.159265
sequence: TRS & id_within_sequence: 6   -0.043903     0.160375    -0.27    0.7848   -0.36149       0.273684
sequence: TSR & id_within_sequence: 6    0.604415     0.159198     3.80    0.0002    0.289159      0.919671
sequence: RST & id_within_sequence: 7   -0.0104784    0.15878     -0.07    0.9475   -0.324905      0.303948
sequence: RTS & id_within_sequence: 7   -0.73212      0.16        -4.58    <1e-04   -1.04896      -0.415276
sequence: SRT & id_within_sequence: 7    0.431845     0.16         2.70    0.0080    0.115001      0.748689
sequence: STR & id_within_sequence: 7   -0.0893864    0.159198    -0.56    0.5755   -0.404642      0.22587
sequence: TRS & id_within_sequence: 7    0.0591428    0.160375     0.37    0.7130   -0.258444      0.376729
sequence: TSR & id_within_sequence: 7   -0.93299      0.159198    -5.86    <1e-07   -1.24825      -0.617734
sequence: RST & id_within_sequence: 8    0.558771     0.15878      3.52    0.0006    0.244345      0.873198
sequence: RTS & id_within_sequence: 8   -0.1594       0.16        -1.00    0.3212   -0.476244      0.157444
sequence: SRT & id_within_sequence: 8   -0.316077     0.16        -1.98    0.0505   -0.632921      0.000767652
sequence: STR & id_within_sequence: 8    0.377635     0.159198     2.37    0.0193    0.0623787     0.692891
sequence: TRS & id_within_sequence: 8   -0.500216     0.160375    -3.12    0.0023   -0.817802     -0.182629
sequence: TSR & id_within_sequence: 8    0.246161     0.159198     1.55    0.1247   -0.0690946     0.561417
sequence: RST & id_within_sequence: 9    0.498742     0.191139     2.61    0.0102    0.120235      0.87725
sequence: RTS & id_within_sequence: 9   -0.0536409    0.16        -0.34    0.7380   -0.370485      0.263203
sequence: SRT & id_within_sequence: 9    0.173038     0.16         1.08    0.2817   -0.143807      0.489882
sequence: STR & id_within_sequence: 9   -0.330719     0.159198    -2.08    0.0399   -0.645975     -0.0154627
sequence: TRS & id_within_sequence: 9    0.409442     0.160375     2.55    0.0120    0.0918552     0.727029
sequence: TSR & id_within_sequence: 9    0.386147     0.159198     2.43    0.0168    0.0708915     0.701403
sequence: RST & id_within_sequence: 10   0.0        NaN          NaN       NaN     NaN           NaN
sequence: RTS & id_within_sequence: 10   0.262141     0.16         1.64    0.1040   -0.0547029     0.578985
sequence: SRT & id_within_sequence: 10  -0.0533425    0.16        -0.33    0.7394   -0.370187      0.263502
sequence: STR & id_within_sequence: 10   0.573641     0.159198     3.60    0.0005    0.258385      0.888897
sequence: TRS & id_within_sequence: 10  -0.343073     0.160375    -2.14    0.0345   -0.66066      -0.0254862
sequence: TSR & id_within_sequence: 10   0.272039     0.159198     1.71    0.0901   -0.0432167     0.587295
sequence: RST & id_within_sequence: 11   0.0        NaN          NaN       NaN     NaN           NaN
sequence: RTS & id_within_sequence: 11  -0.281096     0.16        -1.76    0.0815   -0.59794       0.035748
sequence: SRT & id_within_sequence: 11   0.255426     0.16         1.60    0.1131   -0.0614185     0.57227
sequence: STR & id_within_sequence: 11   0.0        NaN          NaN       NaN     NaN           NaN
sequence: TRS & id_within_sequence: 11  -0.52366      0.19394     -2.70    0.0080   -0.907714     -0.139606
sequence: TSR & id_within_sequence: 11   0.0        NaN          NaN       NaN     NaN           NaN
──────────────────────────────────────────────────────────────────────────────────────────────────────────────

Bioequivalence.EMA_HighlyVariable — Type

EMA_HighlyVariable <: BioequivalenceCriterion

EMA average bioequivalence with expanding limits (ABEL) for highly-variable endpoints (Cmax). When cv_min ≤ CVᵣ ≤ cv_max (defaults 0.30 and 0.50) the acceptance limits expand to exp(±k·σwᵣ) with k = √θ; otherwise standard ABE applies. Point-estimate GMR ∈ [0.80, 1.25] is also required.

Bioequivalence.EMA_NarrowTherapeuticIndex — Type

EMA_NarrowTherapeuticIndex <: BioequivalenceCriterion

EMA narrow-therapeutic-index drug criterion: 90% CI of GMR ⊆ [0.90, 1.11].

Bioequivalence.FDA_HighlyVariable — Type

FDA_HighlyVariable <: BioequivalenceCriterion

FDA reference-scaled average bioequivalence (RSABE) for highly-variable drugs. When CVᵣ ≥ cv_min (default 0.30), Howe's upper bound is required (cb ≤ 0) in addition to the point-estimate GMR ∈ [0.80, 1.25]; otherwise standard ABE.

Bioequivalence.FDA_NarrowTherapeuticIndex — Type

FDA_NarrowTherapeuticIndex <: BioequivalenceCriterion

FDA RSABE for narrow-therapeutic-index drugs. Always applies Howe's upper bound (cb ≤ 0) alongside standard ABE and an upper bound on the within-subject variability ratio σwₜ/σwᵣ.

Bioequivalence.NoAssessmentBioequivalenceCriterion — Type

NoAssessmentBioequivalenceCriterion <: BioequivalenceCriterion

Skip assessment.

Bioequivalence.StandardBioequivalenceCriterion — Type

StandardBioequivalenceCriterion <: BioequivalenceCriterion

Standard average bioequivalence: 90% CI of GMR ⊆ [0.80, 1.25].

Bioequivalence.assess_be — Method

assess_be(
    criterion::Type{<:BioequivalenceCriterion},
    endpoint_output::BioequivalenceEndpointOutput;
    cv_min::Union{Real,Nothing} = nothing,
    cv_max::Union{Real,Nothing} = nothing,
)::BioequivalenceEndpointDecision

Assess a bioequivalence inference output against a regulator criterion and return a BioequivalenceEndpointDecision. assess_be is the assessment stage of the package's lower-level pipeline:

data_pp  = preprocess_be(data; endpoint = :AUC)
output   = run_be(data_pp; ...)
decision = assess_be(EMA_HighlyVariable, output)
set_result_table!(decision)           # optional: populate output.result_table
result_summary_table(output)          # optional: pretty-print the table

pumas_be is the one-call convenience wrapper that does all of the above.

Arguments

criterion: the regulator criterion (a singleton type, not an instance), e.g. EMA_HighlyVariable. See BioequivalenceCriterion subtypes.
endpoint_output: the output of run_be (or pumas_be).
cv_min, cv_max: optional overrides for the switching / cap thresholds on CVᵣ. When nothing, the criterion's regulatory defaults apply (FDA/HV ≈ 0.30, EMA/HV ∈ [0.30, 0.50], others = 0 / Inf). cv_min is consumed only by criteria that switch on CVᵣ (FDA_HighlyVariable, FDA_NarrowTherapeuticIndex, EMA_HighlyVariable); cv_max is consumed only by EMA_HighlyVariable.

Return value

BioequivalenceEndpointDecision with field assessments::BitMatrix — one row per formulation comparison, one column per individual criterion comprising the regulator criterion. The column count depends on the criterion:

StandardBioequivalenceCriterion, EMA_NarrowTherapeuticIndex: 1 column (CI containment).
FDA_HighlyVariable: 2 columns (point estimate, scaled or standard ABE).
EMA_HighlyVariable: 2 columns (point estimate, ABE or ABEL).
FDA_NarrowTherapeuticIndex: 3 columns (standard ABE, Howe RSABE, σ ratio).
NoAssessmentBioequivalenceCriterion: an empty BitMatrix.

A comparison passes the criterion iff all its columns are true.

assess_be does not mutate endpoint_output. To render the assessment into endpoint_output.result_table, call set_result_table!(decision) afterwards.

Bioequivalence.detect_design — Function

detect_design(sequences::AbstractVector) ->
    NamedTuple{(:design, :replicated, :crossover), Tuple{String, ReplicationType, Bool}}
detect_design(
    data::AbstractDataFrame,
    sequence::Union{AbstractString, Symbol} = :sequence
    )::NamedTuple{(:design, :replicated, :supports_rsabe, :supports_abel, :crossover), Tuple{String, ReplicationType, Bool, Bool, Bool}}

Return a named tuple which describes the design. The fields are:

:design - a String describing the design
:num_sequences - an Int of the number of sequences in the design
:num_periods - an Int of the number of periods in the design
:num_formulations - an Int of the number of formulations in the design
:reference_formulation - the Char of the reference formulation
:test_formulations - a Vector of Chars of the test formulations
:replicated - a value of ReplicationType
:supports_rsabe - true if FDA RSABE is supported with this design
:supports_abel - true if EMA ABEL is supported with this design
:crossover - true if a crossover design

Bioequivalence.generate_design — Method

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.geocv2sigma — Method

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.inference_model_type — Method

inference_model_type(be_output::BioequivalenceEndpointOutput)::InferenceModelType

Determines the inference model type used in a bioequivalence study for an endpoint.

Bioequivalence.is_preprocessed — Method

is_preprocessed(data::AbstractDataFrame)

Returns true if data is suspected of having the form returned by preprocess_be.

Bioequivalence.legacy_show — Method

legacy_show(obj::BioequivalenceEndpointOutput)

This can be used to present output of pumas_be in legacy mode, presenting a dataframe obj.result in place of tabular output (out.result_table).

Examples

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

julia> legacy_show(pumas_be(data))
sequence ╲ period │ 1  2  3  4
──────────────────┼───────────
RTTR              │ 8  8  8  8
TRRT              │ 9  9  9  8

Reference scaled using 𝜃 = 0.797

─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
             δ       SE       lnLB     lnUB  trt_margGM  ref_margGM  trt_naiveGM  ref_naiveGM    GMR      LB     UB      CVᵣ     CVₜ  σ_ratio     σ⁺        cb    dof
─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
T - R  0.03568  0.02431  -0.006884  0.07826     7412.14     7152.31      7454.28      7131.48  1.036  0.9931  1.082  0.08022  0.1084     1.35  2.118  0.001813  15.25
─────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────

Bioequivalence.preprocess_be — Method

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

Bioequivalence.pumas_be — Function

pumas_be(
    data::AbstractDataFrame,
    criterion::Type{<:BioequivalenceCriterion} = StandardBioequivalenceCriterion;
    endpoint::Union{AbstractString,Symbol} = "AUC",
    logtransformed::Bool = false,
    reference_scale::Union{Real,Nothing} = nothing,
    cv_min::Union{Real,Nothing} = nothing,
    cv_max::Union{Real,Nothing} = nothing,
    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,
    custom_model::CustomModel = nothing,
    reml::Union{Bool,Nothing} = nothing,
    level_tost::Real = 0.9,
    level_rsabe_stat::Real = 0.95,
    level_var_ratio::Real = 0.9,
    sigdigits::Integer = 4,
    )::BioequivalenceEndpointOutput

This is the main function of the Bioequivalence package for carrying out bioequivalence analysis. It combines the actions of preprocess_be and run_be into one unified action. The return type is a BioequivalenceEndpointOutput object.

Arguments

data: must have id, sequence, period, and an endpoint.
criterion: regulator criterion used to assess bioequivalence (default StandardBioequivalenceCriterion). Use NoAssessmentBioequivalenceCriterion to skip assessment.
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. See notes below.
reference_scale: 𝜃 for reference scale (e.g., FDA/HV (Highly Variable Drug) ≈ 0.797, FDA/NTI (Narrow Therapeutic Index Drug) ≈ 1.11, EMA/HV = 0.5776 = 0.76²). If nothing (default), value is taken from the regulator criterion.
cv_min: switching threshold on CVᵣ below which standard ABE is used instead of reference scaling (FDA/HV = 0.30, EMA/HV = 0.30, others = 0). If nothing (default), the regulator criterion's value is used. Only consumed by criteria that switch on CVᵣ (FDA_HighlyVariable, FDA_NarrowTherapeuticIndex, EMA_HighlyVariable).
cv_max: maximum within-subject CV used for reference scaling (EMA/HV = 0.50, others = Inf). If nothing (default), the regulator criterion's value is used. Only consumed by EMA_HighlyVariable.
custom_model: selects the model fitted to crossover designs. The system's default preference (nothing) is a mixed model (LinearMixedModel, subjects random) for replicated designs and a linear model (LinearModel, subjects fixed) for non-replicated designs. The EMA_HighlyVariable and EMA_NarrowTherapeuticIndex criteria are the exception: they prefer the linear (subjects-fixed) model even on replicated designs. Pass LinearMixedModel to force the mixed model, or LinearModel to force the fixed-effects fit. LinearModel requires homogeneity = true (its residual variance is pooled) and is not supported with the FDA criteria. Passing LinearModel on a non-replicated design is a no-op, since the linear model is already the default there. custom_model also selects the nonparametric model, where the sole accepted value is SignedRankTest: it pools the within-subject differences instead of canceling the period effect between the sequences, and so estimates the same quantity only if that period effect is nil. Passing it for a parallel or paired analysis, or passing any other value on the nonparametric path, raises an ArgumentError.
reml: whether the linear mixed model uses restricted maximum likelihood (REML) or maximum likelihood. Applies only when a LinearMixedModel is fitted; it errors on the fixed-effects (LinearModel), parallel, paired, and nonparametric paths, where REML is meaningless. Defaults to true for the mixed model.
level_tost: applies to the confidence intervals for the GMR.
level_rsabe_stat: applies to the critical boundary for reference-scaled average bioequivalence.
level_var_ratio: applies to the upper bound of the within subject variability ratio (used for NTID with RSABE).
sigdigits: results given with how many significant digits.

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

Description	Formulations	Periods	Sequences	Replicated	Crossover
`R_T`	2	1	2	No	No
`RT`	2	2	1	No	Yes (paired)
`RT_TR`	2	2	2	No	Yes
`RR_RT_TR_TT`	2	2	4	Fully	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	2	2	Fully	No
`RST_RTS_SRT_STR_TRS_TSR`	3	3	6	No	Yes
`ADBC_BACD_CBDA_DCAB`	4	4	4	No	Yes

Note that in most paradigms (all except for nonparametric), the output in includes a Variability section, where CV, CVᵣ, or CVₜ estimates are presented. The meaning varies according to the type of paradigm:

In parallel design, CV is based on a pooled standard deviation from the endpoint, and serves as an estimate for both the within and between subject variability.
In paired design, CV is based on the standard deviation of the differences between formulations. This serves as an estimate for within subject variability.
In crossover design, CV is based on an estimate from the residual of the linear model, and serves as an estimate for the within subject variability.
In crossover design with custom_model = LinearMixedModel, CV is based on an estimate of the residual variability from the mixed model.
In replicated designs, depending on which treatment is replicated, CVᵣ for the reference, and/or CVₜfor the test, estimate the within subject variability for the respective treatment. In these cases, the estimation is via an auxiliary linear model.

Note the the implication homogeneity varies across instances:

In parallel design the nothing default is treated as true which implies an equal variance assumption for the t test. Setting to false spawns an approximate t test that does not assume equality of variance.
In replicated designs (except for criteria = EMA_HighlyVariable), where a mixed model is used, the nothing default is treated as false which implies a mixed model with variance components separate for R and for T. Setting to true implies the same variance component for R and T.
For replicated designs with criteria = EMA_HighlyVariable, the nothing default is instead treated as true (single shared variance component). An explicit homogeneity = false still selects the separate-variance-component model.
In non-replicated crossover designs without a custom_model, a (fixed effects) linear model is used and homogeneity has no effect on the computation.
In non-replicated crossover designs with custom_model = LinearMixedModel, the nothing default is treated as true, and if homogeneity is set to false and error is thrown. This is because some variance parameters are not identifiable in this case.

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)

sequence ╲ period │ 1  2
──────────────────┼─────
RT                │ 9  9
TR                │ 9  9

-----------------------------------------------------------------------------------------
Paradigm: Non replicated crossover design
Model: Linear model
Criteria: Standard ABE
Endpoint: AUC
Formulations: Reference(R), Test(T)
-----------------------------------------------------------------------------------------
                                             Results(AUC)     Assessment   Criteria      
-----------------------------------------------------------------------------------------
R             Geometric Marginal Mean        146.9                                       
              Geometric Naive Mean           146.9                                       
-----------------------------------------------------------------------------------------
T             Geometric Marginal Mean        139.7                                       
              Geometric Naive Mean           139.7                                       
              Geometric Mean T/R Ratio (%)   95.09                                       
              Degrees of Freedom             16                                          
              90% Confidence Interval (%)    [90.76, 99.62]   Pass         CI ⊆ [80, 125]
-----------------------------------------------------------------------------------------
Variability   CV (%) | σ̂                     8.01 | 0.08                                 
-----------------------------------------------------------------------------------------
ANOVA         Formulation (p-value)          0.077                                       
              Sequence (p-value)             2e-5                                        
              Period (p-value)               0.01702                                     
-----------------------------------------------------------------------------------------

julia> output = pumas_be(data, nonparametric = true)

sequence ╲ period │ 1  2
──────────────────┼─────
RT                │ 9  9
TR                │ 9  9

-------------------------------------------------
Paradigm: Nonparametric analysis
Model: Cross-sequence rank sum test (exact)
Endpoint: AUC
Formulations: Reference(R), Test(T)
-------------------------------------------------
                                    Results(AUC) 
-------------------------------------------------
R     Median                        145.3        
-------------------------------------------------
T     Median                        144.4        
-------------------------------------------------
T/R   Median Ratio (%)              94.94        
      90% Confidence Interval (%)   [90.04, 99.6]
-------------------------------------------------

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)

sequence ╲ period │  1   2   3
──────────────────┼───────────
RTT               │ 46  45  43
TRR               │ 47  47  47

-----------------------------------------------------------------------------------------
Paradigm: Replicated crossover (supports ABE and ABEL but not RSABE)
Model: Mixed model with REML (unequal variance)
Criteria: Standard ABE
Endpoint: AUC
Formulations: Reference(R), Test(T)
-----------------------------------------------------------------------------------------
                                             Results(AUC)     Assessment   Criteria      
-----------------------------------------------------------------------------------------
R             Geometric Marginal Mean        102.8                                       
              Geometric Naive Mean           112                                         
-----------------------------------------------------------------------------------------
T             Geometric Marginal Mean        100.1                                       
              Geometric Naive Mean           93.96                                       
              Geometric Mean T/R Ratio (%)   97.38                                       
              Degrees of Freedom             148.7                                       
              90% Confidence Interval (%)    [86.86, 109.2]   Pass         CI ⊆ [80, 125]
-----------------------------------------------------------------------------------------
Variability   CVᵣ (%) | σ̂ᵣ                   42.75 | 0.4097                              
              CVₜ (%) | σ̂ₜ                   69.65 | 0.6289                              
-----------------------------------------------------------------------------------------
ANOVA         Formulation (p-value)          0.7013                                      
              Sequence (p-value)             0.2506                                      
              Period (p-value)               0.01209                                     
-----------------------------------------------------------------------------------------

Bioequivalence.reference_scale_default — Method

reference_scale_default(criterion::Type{<:BioequivalenceCriterion})::Union{Float64,Nothing}

Reference-scaling constant θ baked into the regulator criterion. Returns nothing for criteria that do not use reference scaling (standard ABE, EMA narrow therapeutic index).

Criterion	θ
`FDA_HighlyVariable`	`(ln(1.25) / 0.25)^2 ≈ 0.797` (FDA RSABE)
`FDA_NarrowTherapeuticIndex`	`(ln(1 / 0.9) / 0.1)^2 ≈ 1.110` (FDA NTID)
`EMA_HighlyVariable`	`0.76^2 = 0.5776` (EMA ABEL)
`StandardBioequivalenceCriterion`, `EMA_NarrowTherapeuticIndex`, `NoAssessmentBioequivalenceCriterion`	`nothing`

pumas_be reads this value to drive inference-time σwᵣ-based scaling. Useful for the lower-level workflow `preprocessbe → runbe → assessbewhen the caller wants to mirror the criterion's regulatoryθwithout hard-coding the constant. See alsopumasbe,runbe,assess_be`.

Bioequivalence.reference_scaled_acceptance_bounds — Function

reference_scaled_acceptance_bounds(σwr::Real, θ::Real)
reference_scaled_acceptance_bounds(obj::BioequivalenceEndpointOutput)

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)

sequence ╲ period │ 1  2  3  4
──────────────────┼───────────
RTTR              │ 8  8  8  8
TRRT              │ 9  9  9  8

-----------------------------------------------------------------------------------------
Paradigm: Replicated crossover (supports ABE, ABEL, and RSABE)
Model: Mixed model with REML (unequal variance)
Criteria: Standard ABE
Endpoint: AUC
Formulations: Reference(R), Test(T)
-----------------------------------------------------------------------------------------
                                             Results(AUC)     Assessment   Criteria      
-----------------------------------------------------------------------------------------
R             Geometric Marginal Mean        7152                                        
              Geometric Naive Mean           7131                                        
-----------------------------------------------------------------------------------------
T             Geometric Marginal Mean        7412                                        
              Geometric Naive Mean           7454                                        
              Geometric Mean T/R Ratio (%)   103.6                                       
              Degrees of Freedom             15.25                                       
              90% Confidence Interval (%)    [99.31, 108.2]   Pass         CI ⊆ [80, 125]
-----------------------------------------------------------------------------------------
Variability   CVᵣ (%) | σ̂ᵣ                   8.02 | 0.0801                               
              CVₜ (%) | σ̂ₜ                   10.84 | 0.1081                              
              Variability Ratio (%)          135                                         
-----------------------------------------------------------------------------------------
ANOVA         Formulation (p-value)          0.1624                                      
              Sequence (p-value)             0.3184                                      
              Period (p-value)               0.6651                                      
-----------------------------------------------------------------------------------------

julia> reference_scaled_acceptance_bounds(output)
(0.9310079077618337, 1.0741047327986988)

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

julia> reference_scaled_acceptance_bounds(geocv2sigma(21.42), (log(1 / 0.9) / 0.1)^2)
(0.07362432315197623, 13.582467820258094)

Bioequivalence.result_summary_table — Method

result_summary_table(endpoint_output::BioequivalenceEndpointOutput)

Creates a SummaryTables.Table representing the endpoint output.

Bioequivalence.run_be — Method

run_be(
    data::AbstractDataFrame;
    reference_scale::Union{Real,Nothing} = nothing,
    cv_min::Union{Real,Nothing} = nothing,
    cv_max::Union{Real,Nothing} = nothing,
    wsv_model::Union{Nothing,WSVModel} = nothing,
    nonparametric::Bool = false,
    homogeneity::Union{Bool,Nothing} = nothing,
    level_tost::Real = 0.9,
    level_rsabe_stat::Real = 0.95,
    level_var_ratio::Real = 0.9,
    custom_model::CustomModel = nothing,
    reml::Union{Bool,Nothing} = nothing,
    sigdigits::Integer = 4,
    )::BioequivalenceEndpointOutput

This carries out bioequivalence analysis after it was processed by the preprocess_be function. The return type is a BioequivalenceEndpointOutput object.

Arguments

data: must have the form of data processed by the preprocess_be function. This includes having columns: id, sequence, period, and endpoint.
nonparametric: whether to use a nonparametric analysis (uncommon and usually reserved for Tmax).
homogeneity: whether formulation groups should be modeled with equal variance. See notes below.
reference_scale: 𝜃 for reference scale (FDA/HV (Highly Variable Drug) ≈ 0.797, FDA/NTI (Narrow Therapeutic Index Drug) ≈ 1.11, EMA/HV = 0.5776 = 0.76²). If nothing (default), treated as FDA/HV ≈ 0.797. pumas_be sets this from the regulator criterion.
cv_min: switching threshold on CVᵣ below which standard ABE is used instead of reference scaling. If nothing (default), treated as 0.0 (always scaled). pumas_be sets this from the regulator criterion.
cv_max: maximum within-subject CV for reference scaling. If nothing (default), treated as Inf (no cap). pumas_be sets this from the regulator criterion.
wsv_model: auxiliary linear model used to estimate the within-subject variability σwr (and σwt when applicable) in replicated crossover designs. When nothing (default) the model is auto-selected by within_subject_variability based on the data and homogeneity (the FDA RSABE "dlat" model is used for multi-sequence subsets). Pass additive_sequence_wsv_model to force the EMA SAS scheme. pumas_be sets this from criteria.
custom_model: selects the model fitted to crossover designs. The default preference (nothing) is a mixed model (LinearMixedModel, subjects random) for replicated designs and a linear model (LinearModel, subjects fixed) for non-replicated designs. Unlike pumas_be, run_be is criterion-agnostic and does not apply the EMA exception, so to fit the subjects-fixed (linear) model on a replicated design you must pass LinearModel explicitly, which is what pumas_be does internally for the EMA criteria. LinearModel requires homogeneity = true (its residual variance is pooled). Passing LinearModel on a non-replicated design is a no-op, since the linear model is already the default there. custom_model also selects the nonparametric model, where the sole accepted value is SignedRankTest: it pools the within-subject differences instead of canceling the period effect between the sequences, and so estimates the same quantity only if that period effect is nil. Passing it for a parallel or paired analysis, or passing any other value on the nonparametric path, raises an ArgumentError.
reml: whether the linear mixed model uses restricted maximum likelihood (REML) or maximum likelihood. Applies only when a LinearMixedModel is fitted; it errors on the fixed-effects (LinearModel), parallel, paired, and nonparametric paths, where REML is meaningless. Defaults to true for the mixed model.
level_tost: applies to the confidence intervals for the GMR.
level_rsabe_stat: applies to the critical boundary for reference-scaled average bioequivalence.
level_var_ratio: applies to the upper bound of the within subject variability ratio (used for NTID with RSABE).
sigdigits: results given with how many significant digits.
endpoint_name: the name of the endpoint.

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

Description	Formulations	Periods	Sequences	Replicated	Crossover
`R_T`	2	1	2	No	No
`RT`	2	2	1	No	Yes (paired)
`RT_TR`	2	2	2	No	Yes
`RR_RT_TR_TT`	2	2	4	Fully	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	2	2	Fully	No
`RST_RTS_SRT_STR_TRS_TSR`	3	3	6	No	Yes
`ADBC_BACD_CBDA_DCAB`	4	4	4	No	Yes

In parallel design, CV is based on a pooled standard deviation from the endpoint, and serves as an estimate for both the within and between subject variability.
In paired design, CV is based on the standard deviation of the differences between formulations. This serves as an estimate for within subject variability.
In crossover design, CV is based on an estimate from the residual of the linear model, and serves as an estimate for the within subject variability.
In crossover design with custom_model = LinearMixedModel, CV is based on an estimate of the residual variability from the mixed model.
In replicated designs, depending on which treatment is replicated, CVᵣ for the reference, and/or CVₜfor the test, estimate the within subject variability for the respective treatment. In these cases, the estimation is via an auxiliary linear model.

Note the the implication homogeneity varies across instances:

In parallel design the nothing default is treated as true which implies an equal variance assumption for the t test. Setting to false spawns an approximate t test that does not assume equality of variance.
In replicated designs (except for criteria = EMA_HighlyVariable), where a mixed model is used, the nothing default is treated as false which implies a mixed model with variance components separate for R and for T. Setting to true implies the same variance component for R and T.
For replicated designs with criteria = EMA_HighlyVariable, the nothing default is instead treated as true (single shared variance component). An explicit homogeneity = false still selects the separate-variance-component model.
In non-replicated crossover designs without a custom_model, a (fixed effects) linear model is used and homogeneity has no effect on the computation.
In non-replicated crossover designs with custom_model = LinearMixedModel, the nothing default is treated as true, and if homogeneity is set to false and error is thrown. This is because some variance parameters are not identifiable in this case.

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

julia> output = run_be(pkdata)

sequence ╲ period │ 1  2
──────────────────┼─────
RT                │ 9  9
TR                │ 9  9

───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
                δ        SE        lnLB         lnUB  trt_margGM  ref_margGM  trt_naiveGM  ref_naiveGM       GMR        LB        UB         CV
───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────
T - R  -0.0503868  0.026658  -0.0969286  -0.00384499     139.722     146.943      139.722      146.943  0.950862  0.907621  0.996162  0.0801021
───────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────────

Bioequivalence.set_result_table! — Method

set_result_table!(decision::BioequivalenceEndpointDecision)

Render the assessment carried by decision into the result table of its endpoint output, mutating decision.endpoint_output.result_table in place and returning nothing.

This is the final step of the lower-level workflow preprocess_be → run_be → assess_be → set_result_table!; pumas_be calls it internally. The appropriate per-paradigm renderer is selected from study_paradigm(decision.endpoint_output). After it runs, result_summary_table reflects the assessment.

See also assess_be and result_summary_table.

Bioequivalence.sigma2geocv — Method

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.stat_tests_pvalues — Method

stat_tests_pvalues(data::AbstractDataFrame)

Returns a named tuple of p values with p_formulation, p_period, p_subject_in_group.

Bioequivalence.study_paradigm — Method

study_paradigm(data::AbstractDataFrame, nonparametric::Bool)::EndpointStudyParadigm

Determines the study paradigm of a dataset.

Bioequivalence.study_paradigm — Method

study_paradigm(be_output::BioequivalenceEndpointOutput)::EndpointStudyParadigm

Determines the study paradigm used in a bioequivalence study for an endpoint.

Bioequivalence.within_subject_variability — Method

within_subject_variability(
    data::AbstractDataFrame;
    allownonreplicated::Bool = false,
    homogeneity::Union{Bool,Nothing} = nothing,
    wsv_model::Union{Nothing,WSVModel} = 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. Only observations from subjects with a repeated formulation are used. 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.

The wsv_model argument can override the auxiliary linear model used for the σw estimate. When nothing (default), it is selected automatically based on the data and homogeneity: `formulationpooledwsvmodelwhenhomogeneity = true,periodonlywsvmodelwhen the relevant subset has a single sequence, andperiodbysequencewsvmodelotherwise (the FDA RSABE "dlat" model — see references). Passadditivesequencewsvmodel` to force the EMA SAS scheme.

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.

European Medicines Agency. Questions & Answers: Positions on specific questions addressed to the Pharmacokinetics Working Party. p. 18. https://www.ema.europa.eu/en/documents/scientific-guideline/questions-and-answers-positions-specific-questions-addressed-pharmacokinetics-working-party_en.pdf

U.S. Food and Drug Administration. Draft Guidance on Progesterone. Source of the "dlat" SAS procedure referenced above. https://www.accessdata.fda.gov/drugsatfda_docs/psg/Progesterone_caps_19781_RC02-11.pdf

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.113973       21.0
   2 │ T            missing         missing

Bioequivalence.within_subject_variability_ratio — Function

within_subject_variability_ratio(
    σ::AbstractVector{<:Union{Real, Missing}},
    k::AbstractVector{<:Union{Real, Missing}},
    level_var_ratio::Real = 0.9
    )::Vector{NamedTuple{(:σ_ratio, :σ⁺), ...}}

within_subject_variability_ratio(
    data::AbstractDataFrame,
    level_var_ratio::Real = 0.9
    )::DataFrame

Return the within-subject variability ratio and upper bounds. σ contains the estimated within subject variability. k contains the associated degrees of freedom. level_var_ratio 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::Float64, σ⁺::Float64}}:
 (σ_ratio = 1.349328424089969, σ⁺ = 2.1176303305572595)

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

Private

Bioequivalence.NONPARAMETRIC_DESIGNS — Constant

The designs the rank paradigm accepts, named rather than characterized. They are every two-formulation design the package detects, less Balaam's RR|RT|TR|TT.

Ten of them run a procedure by default: the two parallel routes, where no subject supplies a within-subject difference and the two formulation groups are compared directly; the two single-sequence designs, where formulation is confounded with period and the pooled signed-rank test is all that remains; and the six mirror pairs, where the period effect cancels between the two sequences and Koch's cross-sequence construction applies.

The remaining four hold no usable mirror pair, so their default is refused by _no_mirror_pair_message rather than contaminated. They stay on this list because the pooled analysis is still well defined on them, every subject receiving both formulations, so custom_model = SignedRankTest reaches it. That is the same opt-in the mirror pairs offer, and it costs the same assumption.

Two kinds of design are off the list entirely, and for them custom_model changes nothing:

Balaam's RR|RT|TR|TT. Its RT and TR subjects do form a mirror pair, but its RR and TT subjects supply no difference, so any estimate here rests on part of the study while the medians printed beside it rest on all of it. A block whose two halves describe different subjects is worse than no block.
Anything with more than two formulations. Each pairwise comparison would draw a difference from every sequence, and the block renders a single comparison.

Bioequivalence.BioequivalenceCriterion — Type

BioequivalenceCriterion

An abstract bioequivalence criterion.

Bioequivalence._infer_homogeneity — Method

_infer_homogeneity(m::LinearMixedModel) -> Bool

Heuristic to decide whether the fitted model uses homogeneous residual variance.

In Bioequivalence.jl, heterogeneity is represented by an observation-level random effect grouped by row (see zerocorr(... | row)). If such a term is present we treat the model as heterogeneous.

Bioequivalence._show_html_coeftable — Method

_show_html_coeftable(io::IO, obj::BioequivalenceEndpointOutput)

Write the coefficient table of obj as preformatted HTML.

Used as the fallback when result_table has not been populated. CoefTable provides no text/html method, so its plain-text rendering is escaped and wrapped in a <pre> block.

Bioequivalence.cv_min_default — Method

cv_min_default(criterion::Type{<:BioequivalenceCriterion})::Float64

Criterion's regulatory default cv_min (used when no override is supplied).

Bioequivalence.implied_design — Method

implied_design(data::AbstractDataFrame)::NamedTuple

Sequence and number of subjects per sequence

Bioequivalence.label — Function

label(criterion::Type{<:BioequivalenceCriterion})::String

Short human-readable name shown in result tables (e.g. "FDA RSABE for HV").

Bioequivalence.satterthwaite_confint — Method

satterthwaite_confint(model::LinearMixedModel, δ, SE, level, formulation_str; homogeneity::Bool)

Compute a t-based confidence interval using the Satterthwaite approximation for the denominator degrees of freedom.

homogeneity controls whether Bioequivalence.jl should treat the model as having homogeneous residual variance (true) or heterogeneous residual variance (false). It must be resolved by the caller (the BE inference code decides the default).

Bioequivalence.uses_cv_min — Method

uses_cv_min(criterion::Type{<:BioequivalenceCriterion})::Bool
uses_cv_max(criterion::Type{<:BioequivalenceCriterion})::Bool

Whether the criterion consumes a user-supplied switching/cap threshold on CVᵣ. uses_cv_min reflects assess-time consumption; uses_cv_max reflects either assess-time consumption or the inference-time σwᵣ clamp the criterion relies on. Criteria returning false reject any override other than nothing or the criterion's default (a no-op) — see `validatecvoverrides`.