Why the Estimation Method Matters
Reference scaling uses the variability of repeated reference-product measurements within the same subject. This quantity is not the total variability among all observations: it excludes between-subject differences and must be separated from treatment, period, and other fixed effects. In a replicated crossover model, restricted maximum likelihood (REML) estimates the variance components after accounting for those fixed effects.
The distinction matters because the estimated within-subject reference variance can affect an EMA expanded limit or the FDA linearized reference-scaled criterion. A variance estimate that is too large or too small can therefore change a scaling calculation. REML is a model-based estimation method, not a correction that guarantees a suitable model or a stable estimate.
Model and Restricted Likelihood
A general linear mixed model for log-scale responses is \(y=X\beta+Zu+\varepsilon\). Here \(y\) is the vector of observed responses, \(X\beta\) represents fixed effects such as sequence, period, and formulation, \(u\) represents random effects, and \(\varepsilon\) is residual variation. The covariance matrix of \(y\) is \(V(\theta)\), determined by variance components \(\theta\). The model must specify which covariance component represents the within-subject reference variation. In an appropriate replicate design, repeated reference responses within a subject identify that component.
The model assumes that the random effects and errors have the stated mean-zero distributions and covariance structure, commonly multivariate normal distributions, and that different subjects are independent. It also assumes the fixed-effects design is correctly specified and the variance components are meaningful for the design. These assumptions are specific to the fitted model; “REML” alone does not establish that they hold.
Ordinary maximum likelihood estimates fixed effects and variance components using the likelihood of all observations. REML instead constructs a likelihood from linear combinations of the data that eliminate the fixed effects. For a full-rank fixed-effects design matrix \(X\), its log likelihood, up to constants, is
where
REML chooses the allowable variance-component values that maximize \(\ell_R(\theta)\), subject to the components being nonnegative. In general, numerical optimization is required; the estimate of \(\sigma_{WR}^2\) is one component of the fitted covariance model, not necessarily a simple variance of all reference observations.
Why Fixed-Effect Adjustment Is Built In
A useful derivation comes from a balanced, simplified replicate setting. Suppose each of \(n\) subjects has two reference responses, and fixed period effects have been adjusted so that the replicate contrast for subject \(i\) is \(D_i\). The subject-specific intercept cancels in the difference. If the two reference residuals are independent with variance \(\sigma_{WR}^2\), then
where \(\delta\) is any remaining common mean contrast and \(e_{i1},e_{i2}\) are independent with variance \(\sigma_{WR}^2\). Estimating the unknown mean contrast removes one degree of freedom. The restricted likelihood for these contrasts is proportional to the likelihood of their deviations from \(\bar D\). Maximizing it with respect to \(\sigma_{WR}^2\) gives
The denominator contains two factors: \(n-1\) because the contrast mean is estimated, and 2 because each difference has variance \(2\sigma_{WR}^2\). This result is exact for the simplified assumptions just stated. A fuller crossover model may include additional fixed effects, random effects, unequal replication, or missing observations; its REML estimate then comes from the full model likelihood and need not equal this shortcut.
Worked Example: Estimating Variance From Adjusted Replicate Contrasts
Suppose a hypothetical balanced study has six subjects with two reference administrations each. After the model's period adjustment, the six replicate contrasts, in log units, are \(-0.20,-0.10,-0.05,0.05,0.10,0.20\). Assume independent subjects, equal within-subject reference variance, and the normal contrast model above. The contrasts sum to zero, so \(\bar D=0\).
The restricted-likelihood estimate is
The estimated within-subject reference standard deviation is \(\sqrt{0.01050}=0.1025\) log units. The variance, \(0.01050\) log-units squared, is the quantity used in variance-based criteria; the standard deviation is often needed for limits expressed using a standard deviation. This estimate is not the variance of the six raw reference responses, because that would include subject-to-subject differences and potentially period effects.
From the Estimate to a Scaling Calculation
The REML estimate is an input to the applicable regulatory calculation, not itself a bioequivalence decision. For example, EMA expanded limits use the within-subject reference standard deviation when it exceeds the specified scaling threshold, with a cap on the resulting limits. In the following hypothetical calculation, assume the relevant endpoint and analysis meet the applicable requirements, and use the estimated variance from the adjusted contrasts below.
Worked Example: Checking an EMA Expansion Threshold
In a second hypothetical balanced study, six adjusted replicate-reference contrasts are \(-0.60,-0.35,-0.10,0.10,0.35,0.60\) log units. Their mean is zero. Assume the same contrast model and independent, normally distributed contrasts.
Their sum of squared deviations is \(2(0.60^2+0.35^2+0.10^2)=0.9850\). Thus,
Because \(0.3138\) exceeds the EMA scaling threshold \(0.294\), the expansion calculation is triggered. Before the cap, the upper limit is \(\exp(0.760\times0.294)=1.250\) after using the capped standard-deviation value \(0.294\); the corresponding lower limit is its reciprocal, approximately \(0.800\). The expanded limits therefore reach the cap, \(0.800\) to \(1.250\). This is the limit calculation only; the study's treatment estimate and its confidence interval must still satisfy the applicable bioequivalence requirements.
For the FDA linearized criterion, \(\widehat{\sigma}_{WR}^2\) enters a quantity involving both the squared treatment-reference difference and the scaled reference variance. For illustration, let the hypothetical estimated log-scale treatment-reference difference be \(d=0.120\), and let \(\theta=\{\ln(1.25)/0.25\}^2\), approximately \(0.7967\). Using the variance estimate \(0.09850\) from the preceding example gives the plug-in quantity
Worked Example: Calculating the FDA Criterion Component
The point estimate of the criterion component is
The negative value shows that this plug-in component is below zero. It is not by itself the FDA decision: the applicable procedure evaluates an upper confidence bound for the criterion, and the point-estimate constraint and other requirements also matter. Reporting only this negative value as proof of bioequivalence would omit the uncertainty calculation.
Interpreting the Fit and Choosing an Estimator
In a full mixed model, REML commonly finds variance components by numerical optimization. The estimate is conditional on the chosen fixed-effects design and covariance structure. Check that the software has converged, that the reported component is specifically the within-subject reference variance, and that the fitted model corresponds to the replicate design and intended analysis. A zero estimate can occur when the likelihood is maximized at the nonnegative boundary; it means the data and model provide no positive estimate beyond that boundary, not that reference responses are known to be identical.
The simplified contrast calculation offers a useful check in a balanced setting, but it is not a universal replacement for fitting the specified model. Warning signs include a component reported for the wrong grouping factor, convergence warnings, a boundary estimate that conflicts with the observed replicate differences, or strong discrepancies between model-based fitted variation and the replicate contrasts. Such signs call for review of the model, data coding, and assumptions; they do not identify a remedy on their own.
| Approach | What it does | Key distinction |
|---|---|---|
| REML | Estimates variance components from likelihoods of contrasts that remove fixed effects | Generally accounts for fixed-effect estimation when estimating variance |
| Maximum likelihood | Maximizes the likelihood of the full data for fixed effects and variance components | Can give more downward-biased variance estimates in finite samples because it does not make the same fixed-effect adjustment |
| Method of moments or balanced contrast formula | Equates observed contrast variation to its model expectation | Can be transparent and useful in simple balanced designs, but may not extend directly to complex or unbalanced data |
REML is not automatically superior in every respect. Its estimate can still be imprecise in a small sample, sensitive to model specification, or at the boundary. Nor does REML remove the need for an appropriate endpoint, design, or regulatory analysis. Its main role here is to estimate a variance component while accounting for the fixed effects in the specified model.
Common Mistake
Wrong reasoning: “The sample variance of all reference observations is the within-subject reference variance, so I can insert it directly into a scaled criterion.” This includes differences between subjects and may also include period effects; it does not isolate within-subject reference variation.
Correction: Use a model or valid replicate contrasts that separate the within-subject reference component from fixed effects and other sources of variation. Then use the variance or standard deviation required by the particular regulatory formula, and complete that formula's full decision procedure.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions.
- Why does the balanced replicate-contrast formula divide the sum of squared deviations by \(2(n-1)\), rather than by \(n\)?
- A REML fit reports \(\widehat{\sigma}_{WR}^2=0.0400\). What is the corresponding within-subject reference standard deviation, in log units?
- Why is the variance of all raw reference responses generally not an estimate of \(\sigma_{WR}^2\)?
- A plug-in FDA criterion component is negative. Does that alone establish bioequivalence? Explain.
- What does a REML estimate of zero imply, and what does it not imply?
Solutions. (1) The replicate difference has variance \(2\sigma_{WR}^2\), giving the factor 2, and estimating the common contrast mean uses one degree of freedom, giving \(n-1\). (2) \(\widehat{\sigma}_{WR}=\sqrt{0.0400}=0.200\) log units. (3) Raw reference responses can include between-subject variation and fixed period effects, whereas \(\sigma_{WR}^2\) is the within-subject reference component. (4) No. The FDA procedure requires the specified upper confidence bound and satisfaction of the other applicable constraints; the point component alone is insufficient. (5) The fitted likelihood is maximized at the nonnegative boundary for that component. It does not establish that true reference variability is exactly zero.