Why Forward Inclusion Uses Nested Model Comparisons
A covariate may explain systematic differences in a pharmacokinetic or pharmacodynamic parameter across individuals. Forward inclusion asks a specific question at each step: does adding this candidate effect improve the fit enough to justify the additional parameter, given the effects already in the model? The likelihood-ratio test (LRT) answers that question by comparing two fitted, nested models.
For subject \(i\) at observation time \(t_{ij}\), a typical NLME model has response \(Y_{ij}=f(t_{ij},\theta_i)+\epsilon_{ij}\), where \(f\) is the structural model, \(\theta_i\) contains individual parameters, and \(\epsilon_{ij}\) follows a specified residual-error model. Individual parameters may depend on population parameters, covariates, and random effects. For example, a positive clearance parameter may be written
Here \(x_i\) is a measured, possibly centered covariate, \(\beta\) is its fixed-effect coefficient, and \(\eta_{CL,i}\) is an individual random effect. The reduced model sets \(\beta=0\); the expanded model estimates \(\beta\). This parameterization makes the covariate multiplier positive. The exact covariate model can differ, but for an LRT the reduced model must be obtained from the expanded model by fixing the added effect or effects to their null values.
Assumptions and the Likelihood-Ratio Statistic
Let \(\ell_0\) be the maximized log-likelihood for the reduced model and \(\ell_1\) the maximized log-likelihood for the expanded model. The usual LRT relies on the following conditions: both models use the same response observations and likelihood definition; the models are nested; the expanded model is successfully optimized; and the null value of the added coefficient is an interior point of its parameter space. The likelihood and random-effects structure must also be appropriate enough for the model comparison to be meaningful. Under the null and standard large-sample regularity conditions, the statistic is approximately chi-square distributed.
where \(q\) is the number of freely estimated coefficients added in the expanded model. The p-value is the upper-tail probability \(P(\chi^2_q\ge D)\). Because the expanded model contains the reduced model, its maximized log-likelihood should be at least as large, so \(D\) should be nonnegative, apart from small numerical discrepancies. The approximation is not an exact finite-sample result.
Many NLME programs report an objective function equal to \(-2\ell\), possibly plus a constant shared by the compared models. If the reported values are \(OFV_0\) and \(OFV_1\), calculate \(D=OFV_0-OFV_1\): a lower objective function for the expanded model gives a positive statistic. Do not subtract the values in the opposite order.
Worked Example: Testing a Weight Effect on Clearance
Hypothetical data. Suppose 120 subjects have usable concentration observations. The current model has no weight effect on clearance, and the expanded model adds one coefficient using \(x_i=(\text{weight}_i-70\text{ kg})/10\text{ kg}\). Both models use the same observations, structural model, random-effects structure, and residual-error model, and both optimizations converge. Their maximized log-likelihoods are \(\ell_0=-412.80\) and \(\ell_1=-409.60\).
The expanded model adds one coefficient, so \(q=1\). The statistic is
For one degree of freedom, \(P(\chi^2_1\ge6.40)\approx0.0114\). At a prespecified \(\alpha=0.05\), this comparison provides evidence against \(H_0:\beta=0\) in favor of the model with a weight effect. It does not establish that weight causes a clearance difference.
Suppose the fitted coefficient is \(\hat\beta=0.18\), and the expanded model estimates \(CL_{\mathrm{pop}}=5.05\) L/h at the reference weight of 70 kg. A 10 kg higher weight corresponds to \(x=1\), so the modeled typical clearance ratio is \(\exp(0.18)\approx1.197\), or about 19.7% higher than at 70 kg, holding the other model features fixed. The reduced model estimated typical clearance at 5.20 L/h; that estimate and other parameters can shift when the covariate is included. The LRT tests the added coefficient, not whether every changed parameter is individually significant.
Conducting a Forward Search
Forward inclusion begins with a defined base model, often containing the structural, random-effects, and residual-error components considered necessary before covariate selection. At each iteration, fit a separate expanded model for each remaining candidate, adding that candidate to the current model. Compare each expanded model with the same current model, using the same data and likelihood specification. Among candidates meeting the prespecified inclusion rule, add the strongest eligible candidate, then repeat the comparisons from this updated model. If no candidate meets the rule, stop.
Record its maximized log-likelihood or objective function, parameter estimates, and convergence diagnostics.
For each remaining candidate, fit a model that differs from the current model only by the specified added effect.
Calculate the statistic with the correct sign, use the appropriate added-parameter count, and apply the prespecified inclusion rule.
Add an eligible effect, refit the resulting model, record changes in estimates and diagnostics, and test remaining candidates conditional on the updated model.
For transparency, keep a record of the current model at each step, the candidate tested, both likelihood values, the statistic, the reference degrees of freedom, the p-value, the decision, and relevant parameter estimates. In particular, report the estimate of each accepted covariate coefficient and the population parameter it modifies. Centering and scaling should be stated, because they determine the reference value and units of a coefficient.
Worked Example: A Candidate’s Evidence Changes After an Addition
Hypothetical data. A base model is fitted to the same observations for 150 subjects. Its log-likelihood is \(-500.0\). Two one-coefficient candidates, age and sex, are tested separately against this base model. Adding age gives \(-494.9\), while adding sex gives \(-497.4\).
For age, \(D=2\{-494.9-(-500.0)\}=10.2\), with \(q=1\) and \(p\approx0.0014\). For sex, \(D=2\{-497.4-(-500.0)\}=5.2\), also with \(q=1\) and \(p\approx0.0226\). With an inclusion rule of \(p<0.05\), both are eligible at this step; age has the larger LRT statistic and is selected first.
Suppose age is centered at 50 years and the accepted age model estimates \(CL_{\mathrm{pop}}=4.62\) L/h at age 50, with \(\hat\beta_{\mathrm{age}}=-0.12\) per 10 years. The model’s typical clearance multiplier for a 10-year increase is \(\exp(-0.12)\approx0.887\), about 11.3% lower. These are model-based associations, not causal effects.
Now test sex conditional on age. Suppose the age model has log-likelihood \(-494.9\), and the age-plus-sex model has \(-493.7\). Then \(D=2\{-493.7-(-494.9)\}=2.4\), \(q=1\), and \(p\approx0.121\). Sex no longer meets the same inclusion rule. Its initial comparison asked whether sex improved the base model; the second asks whether it adds evidence beyond age. These are different questions. Suppose adding sex changes the age coefficient from \(-0.12\) to \(-0.11\) and estimates a sex coefficient of \(0.08\). Record those changes, but do not include sex under the stated rule.
Justification and Practical Checks
The statistic follows from comparing how well the two models can maximize the same likelihood. The likelihood ratio is \(L_1/L_0\), where \(L_1\) and \(L_0\) are the maximized likelihoods. Nesting means the expanded model can reproduce the reduced model by setting its added coefficient to the null value; therefore its maximum cannot be smaller. Taking twice the logarithm gives the displayed statistic. Under \(H_0\), regularity conditions yield the approximate chi-square reference distribution, with one degree of freedom for each independent coefficient added.
This justification also shows why the comparison must be controlled. If the models use different subjects, different response data, or different likelihood definitions, their likelihoods do not provide the intended nested comparison. A change to the random-effects or residual-error structure at the same time as adding a covariate also prevents attributing the likelihood improvement specifically to that covariate.
Before interpreting a result, check that the expanded-model objective function is no larger than the reduced-model value, the optimizer reports successful convergence, and parameter estimates are finite and plausible. Review gradients, boundary estimates, standard errors, and sensitivity to starting values where available. An apparent likelihood improvement from a failed or local optimization is not reliable evidence. If the added parameter is on a boundary under the null, such as a variance component constrained to be nonnegative, the standard chi-square calibration may fail; that is not the ordinary fixed covariate-coefficient test described here.
Forward-search p-values are conditional on the model reached at that step. Repeatedly testing candidates and selecting among them makes the sequence data-adaptive, so the nominal p-values do not by themselves provide an overall error rate for the entire search. The selected model should be interpreted with that limitation in mind and supported by pharmacological plausibility, covariate support, and model diagnostics.
| Approach | What it compares | Practical distinction |
|---|---|---|
| Forward LRT | Nested current and expanded models | Tests the incremental likelihood improvement from an addition, conditional on included effects. |
| Wald test | An estimate relative to its estimated uncertainty within a fitted model | Can be convenient, but relies on a local standard-error approximation and is not the same refitted nested-model comparison. |
| AIC comparison | Fit balanced against a parameter-count penalty | Ranks candidate models by a criterion; it is not an LRT p-value or a test at a specified significance level. |
These approaches answer related but different questions. The LRT is useful when the scientific question is whether a specified nested addition improves fit under a chosen inclusion rule. It does not guarantee that the selected model is uniquely correct, that estimates are unbiased after selection, or that a statistically supported covariate is clinically important.
Common Mistakes and Limitations
- Using the wrong subtraction order. A negative statistic usually signals that likelihoods or objective functions were subtracted in the wrong direction, or that optimization did not reach the expected nested-model maximum. Check the fits rather than taking an absolute value.
- Comparing models fitted to different records. Missingness can silently change the observations. Refit the reduced model on exactly the records used by the expanded model before making the comparison.
- Changing several model features at once. If covariates, random effects, and residual error all change together, the LRT cannot isolate the contribution of the candidate covariate.
- Treating a nonsignificant result as proof of no effect. It means the comparison did not meet the chosen rule; limited data or weak covariate support may also explain the result.
- Reading the final p-value as if the candidate had been specified in advance. Forward selection uses the observed data to choose additions, so nominal sequential p-values do not account for the full search process.
Check Your Understanding
Answer each question, then check the worked solutions below.
- A reduced model has log-likelihood \(-220.0\), and a nested one-coefficient model has log-likelihood \(-217.5\). Calculate the LRT statistic. At the 0.05 level, is there evidence against the null coefficient value?
- Why is it invalid to compare a reduced model fitted to 100 subjects with an expanded model fitted to 94 subjects using the ordinary LRT described here?
- A candidate has \(p=0.03\) when added to the base model but \(p=0.14\) when added after another covariate. Explain the difference in the questions being tested.
- What should be checked if the expanded model reports a worse objective function than its nested reduced model?
Solutions. (1) \(D=2\{-217.5-(-220.0)\}=5.0\). With \(q=1\), \(p\approx0.0253\), so there is evidence against the null at \(\alpha=0.05\). (2) The likelihoods are based on different response records, so they do not constitute the required same-data nested comparison. Refit both models to the same subjects. (3) The first test asks whether the candidate improves the base model; the second asks whether it improves a model already containing the other covariate. Conditional evidence can differ. (4) Check the subtraction convention, ensure the models are nested and use the same data and likelihood, and inspect convergence and optimization. Do not use the absolute difference as a substitute for diagnosing the problem.