Tutorials › Population Pharmacokinetics › Testing Covariate Effects on Between-Subject Variability

Covariate Selection Strategies and Their Statistics · Tutorial 298 of 599

Testing Covariate Effects on Between-Subject Variability

Compare covariate-dependent variance models with likelihood-based tests, interpret variance effects, and distinguish regular tests from boundary cases.

Advanced 12 min read

What You'll Learn

  • Specify a covariate effect on the variance of a population PK random effect.
  • Interpret a coefficient in a log-variance model as a variance and standard-deviation ratio.
  • Form a likelihood-ratio test for a nested covariate variance model.
  • Determine when the likelihood-ratio statistic has a chi-square or boundary-mixture reference distribution.
  • Recognize weak identifiability from sparse information, flat likelihoods, and unstable variance estimates.
  • Compare likelihood-ratio, Wald, and bootstrap approaches without confusing their assumptions.

Why Test a Covariate Effect on Variability?

A covariate can help explain differences in a typical pharmacokinetic parameter, but it can also predict how much subjects vary around that typical value. For example, a population may have more between-subject variability in clearance among subjects with a particular characteristic, even if their covariate-adjusted typical clearance is unchanged. Testing these are distinct questions: a covariate effect on the typical value concerns the location of the random-effects distribution, whereas an effect on between-subject variability concerns its dispersion.

Here the target is a covariate effect on the variance of a random effect. The central statistical issue is not only how to compare the fitted models, but also whether the null hypothesis lies in the interior of the parameter space. That distinction determines whether the usual chi-square likelihood-ratio reference distribution is justified.

A Covariate Model for Random-Effect Variance

Let \(\eta_i\) be a subject-specific random effect for a PK parameter in subject \(i\), and let \(z_i\) be a measured covariate. For example, a clearance model might include \(\eta_i\) on the log scale. A convenient variance model is

$$ \eta_i \mid z_i \sim N\left(0,\omega^2\exp(\gamma z_i)\right). $$

Here \(\omega^2\) is the variance when \(z_i=0\), and \(\gamma\) is the covariate coefficient for log variance. Centering \(z_i\) at a meaningful reference value makes \(\omega^2\) the variance at that value. The model assumes a zero conditional mean for \(\eta_i\), conditional normality, positive variance, and the specified exponential relationship between the covariate and variance. It does not by itself assert that the covariate changes the typical PK parameter.

Definition: In the log-variance model \(\omega_i^2=\omega^2\exp(\gamma z_i)\), a one-unit increase in \(z_i\) multiplies the random-effect variance by \(\exp(\gamma)\). It multiplies the random-effect standard deviation by \(\exp(\gamma/2)\). The null hypothesis of no variance covariate effect is \(H_0:\gamma=0\); the two-sided alternative is \(H_a:\gamma\ne 0\).

The exponential form guarantees a positive variance for every finite \(z_i\). A categorical covariate can instead define separate variances by group, or be represented through indicator variables in a log-variance model. In either case, the hypothesis must match the actual parameterization: a joint test of several coefficients is not a one-degree-of-freedom test.

Worked Example: Interpreting a Log-Variance Coefficient

Suppose a hypothetical fitted model for a clearance random effect uses a centered covariate \(z_i\), measured in 10-unit increments, and estimates \(\hat{\gamma}=-0.25\). The reference variance is estimated as \(\hat{\omega}^2=0.16\) on the log-clearance scale. Interpret the coefficient and compare subjects one covariate unit apart.

The fitted variance ratio for a one-unit increase is

$$ \exp(\hat{\gamma})=\exp(-0.25)=0.7788. $$

Thus the estimated variance is about 0.779 times as large per 10-unit covariate increase, a reduction of \(1-0.7788=0.2212\), or 22.1%. The estimated standard-deviation ratio is

$$ \exp(\hat{\gamma}/2)=\exp(-0.125)=0.8825. $$

At the reference covariate value, the estimated standard deviation is \(\sqrt{0.16}=0.400\). One covariate unit higher, the fitted variance is \(0.16(0.7788)=0.1246\), and the standard deviation is \(\sqrt{0.1246}=0.3530\). These are estimates of variability in the random effect, not a 22.1% decrease in typical clearance. They describe a conditional model relationship and do not establish a causal effect.

Likelihood-Ratio Comparison

To test a covariate variance effect, fit a null model with \(\gamma=0\) and an alternative model with \(\gamma\) freely estimated. Keep the structural PK model, residual-error model, other covariate effects, and other random-effects assumptions the same. The null model must be a special case of the alternative, and both fits must use the same observations and likelihood definition.

Definition: The likelihood-ratio statistic is \(D=2(\ell_A-\ell_0)\), where \(\ell_A\) and \(\ell_0\) are the maximized log likelihoods under the alternative and null models. Larger values of \(D\) indicate greater improvement in fit under the alternative.

For a regular one-parameter test, the null model must be correctly specified, the parameter \(\gamma\) must be identifiable under \(H_0\), the true parameter must be an interior point of its allowed parameter space, and standard likelihood regularity conditions must hold. Subject-level observations are ordinarily assumed independent conditional on modeled covariates and population parameters, with the random-effects distribution and observation model correctly specified. The likelihood must account for the random effects, rather than treating estimated subject-specific effects as error-free observations.

Under these conditions, the likelihood-ratio statistic has an approximate \(\chi^2_1\) distribution under \(H_0\), where the subscript is the difference in the number of free parameters. The basis is the local quadratic approximation to the log likelihood: near the null estimate, twice the log-likelihood improvement behaves like the squared standardized estimate of the one added parameter. That squared standard-normal quantity has a \(\chi^2_1\) distribution. For \(q\) jointly tested, identifiable parameters under regular conditions, the corresponding reference distribution is \(\chi^2_q\).

Worked Example: Testing a Continuous Variance Covariate

In a hypothetical population PK analysis, the null model has maximized log likelihood \(\ell_0=-812.6\), and the alternative model, which adds one log-variance coefficient, has \(\ell_A=-810.5\). Assume both models are fitted by maximum likelihood to the same data, the alternative is nested, and the estimated baseline variance is positive and not near its boundary. Test the covariate effect at \(\alpha=0.05\).

The likelihood-ratio statistic is

$$ D=2\left[(-810.5)-(-812.6)\right] =2(2.1)=4.20. $$

Because the null value \(\gamma=0\) is interior and the baseline variance is positive, the regular one-parameter approximation applies. The p-value is

$$ P(\chi^2_1\ge 4.20)\approx 0.0404. $$

Since \(0.0404<0.05\), reject \(H_0\) at the 5% significance level. Under the stated model and regularity assumptions, the data provide evidence that the covariate predicts between-subject variance. This result alone does not show that the effect is clinically important, that it improves prediction in new subjects, or that the covariate causes variability to change.

When the Null Is on a Boundary

A variance is constrained to be nonnegative. This matters when the null hypothesis sets a variance component itself to zero. Under that null, zero is on the boundary of the parameter space, and the usual regular chi-square result need not apply. A familiar special case is testing whether one nonnegative variance component is zero, with the other nuisance parameters identifiable and regular. Under those conditions, the asymptotic null distribution is a 50:50 mixture of a point mass at zero and \(\chi^2_1\).

The mixture can be understood through the local likelihood approximation. The unconstrained standardized estimate would be approximately normal. Constraining the variance component to be nonnegative projects a negative estimate onto zero, giving a likelihood-ratio statistic of zero for that half of the limiting normal outcomes. For the other half, the squared positive part has the \(\chi^2_1\) tail. Thus, for observed \(D>0\), the upper-tail p-value is one-half the ordinary \(\chi^2_1\) upper-tail probability.

Worked Example: A Single Variance Component at Zero

Suppose a hypothetical analysis compares a model with no between-subject random effect for one parameter against a model with one nonnegative variance component. The likelihood-ratio statistic is \(D=3.00\). Assume the single-component boundary result applies: nuisance parameters are identifiable under the null, and no additional parameters become unidentified when the variance is zero.

The appropriate p-value is the mixture tail:

$$ P(D\ge 3.00\mid H_0) =\tfrac{1}{2}P(\chi^2_1\ge 3.00) =\tfrac{1}{2}(0.0832646) =0.0416323. $$

The p-value is approximately 0.0416. At \(\alpha=0.05\), this is evidence against the null model with zero variance, under the stated assumptions. It is not a valid general correction for every test involving a variance model.

In particular, testing \(H_0:\gamma=0\) in \(\omega_i^2=\omega^2\exp(\gamma z_i)\) is ordinarily a regular test if \(\omega^2>0\): \(\gamma\) is identifiable, and zero is an interior value for \(\gamma\). Do not use the 50:50 mixture merely because the model contains a variance. If instead the null sets \(\omega^2=0\), the covariate coefficient \(\gamma\) may cease to be identifiable, since all covariate values imply zero variance. The simple mixture result then may not apply. More complicated boundary constraints or unidentified parameters require a problem-specific null distribution; a parametric bootstrap under the null is one practical way to calibrate the statistic when its asymptotic reference law is uncertain.

Identifiability, Diagnostics, and Alternatives

A variance covariate effect is identifiable only when the data contain information to distinguish its pattern from other sources of variation. There must be meaningful covariate variation across subjects. The individual concentration-time data, together with the PK model, must also contain enough information to estimate the relevant random-effect distribution. If subjects have very sparse observations, differences in estimated variability can be difficult to separate from residual error or structural-model inadequacy. A large nominal sample size does not guarantee that a variance coefficient is estimable if the observations are weakly informative about individual deviations.

Inspect whether the alternative fit converged reliably, whether the estimated variance is near zero, and whether the log likelihood changes materially across reasonable starting values. A nearly singular information matrix, a very large standard error for \(\hat{\gamma}\), a flat profile likelihood, or strong dependence between variance parameters signals weak identification. Compare the fitted variance pattern with the range and distribution of the covariate; a fitted relationship driven by a few extreme covariate values should not be extrapolated casually. Check that the structural and residual-error assumptions are plausible, since misspecification can be absorbed by random-effect variance.

The likelihood-ratio test is generally preferable to a Wald test when the likelihood is sufficiently well behaved but the coefficient estimate is skewed or the standard error is unreliable. A Wald test divides \(\hat{\gamma}\) by its estimated standard error and uses a normal approximation; it can be especially misleading near boundaries or with a flat likelihood. A score test evaluates the likelihood slope at the null and can avoid fitting the full alternative, but still requires an appropriate reference distribution and reliable information under the null. A parametric bootstrap can estimate the null distribution by repeatedly simulating and refitting both models, but it depends on a credible null model and successful fitting. AIC and BIC compare penalized fit for model selection; they do not supply the same null-hypothesis p-value as a likelihood-ratio test. When many covariates are screened, the testing threshold and multiplicity strategy also matter.

Common pitfall: “This is a variance model, so the likelihood-ratio test must use a 50:50 mixture.” That reasoning is incomplete. The mixture applies to particular tests where a nonnegative variance component is zero under the null and the required boundary conditions hold. A test of a log-variance covariate coefficient with a positive baseline variance is usually regular and uses \(\chi^2_1\), not the mixture.

The main limitation is that evidence of a variance relationship depends on the chosen random-effects distribution, covariate function, structural model, and residual-error model. The likelihood-ratio result is an asymptotic approximation, not a guarantee of reliable calibration in small or weakly informative samples. Report the variance parameterization, the tested hypotheses, the fitting method, the reference distribution, and the estimated variance ratio with its units or covariate scale. A small p-value should be accompanied by an assessment of magnitude and stability, not treated as proof that the covariate must be retained.

Check Your Understanding

Answer each question, then compare with the solutions below.

  1. In \(\omega_i^2=\omega^2\exp(\gamma z_i)\), what is the variance ratio for a one-unit increase in \(z_i\) when \(\gamma=0.30\)? What is the standard-deviation ratio?
  2. A test adds two identifiable variance-covariate coefficients, with a regular null and a nested alternative. What is the usual asymptotic reference distribution?
  3. Why is a test of \(\gamma=0\) usually regular when the baseline variance is positive?
  4. For a single variance component tested against zero, with the stated 50:50 mixture assumptions and \(D=3.00\), what is the p-value to four decimal places?
  5. Name two signs that the covariate variance effect may be weakly identified.

Solutions. (1) The variance ratio is \(\exp(0.30)=1.350\); the standard-deviation ratio is \(\exp(0.15)=1.162\). (2) Use \(\chi^2_2\), because two parameters are added and the test is regular. (3) With positive baseline variance, \(\gamma=0\) is an interior value and remains identifiable: it specifies equal variance across covariate values, rather than eliminating the random effect. (4) The p-value is \(\tfrac12 P(\chi^2_1\ge3.00)=0.0416323\), or 0.0416. (5) Examples include a flat profile likelihood, a very large coefficient standard error, unstable estimates across starting values, or a nearly singular information matrix.

Clinical Biostats