Wald and Score Tests in Nonlinear Mixed-Effects Models

Three Routes to Testing a Covariate Effect

A covariate effect in a nonlinear mixed-effects model may be assessed in three standard ways: a Wald test, a score test, or a likelihood-ratio test. Each asks whether the data are compatible with a specified null effect, but they use different parts of the likelihood. A Wald test uses the estimated effect and its uncertainty from the fitted alternative model. A score test evaluates the likelihood's slope at the null model. A likelihood-ratio test compares the maximized likelihoods of the null and alternative models.

These distinctions matter in population pharmacokinetics because model parameters are estimated jointly. An estimated covariate effect can be correlated with typical PK parameters, random-effects parameters, or residual-error parameters. The tests can therefore give different results in finite samples, even when their large-sample reference distributions agree.

Definition: Let \(\psi\) denote the covariate-effect parameter or vector of covariate-effect parameters of interest, and let \(\lambda\) collect the remaining model parameters. The null hypothesis specifies \(q\) restrictions, commonly \(H_0:\psi=0\), and the alternative allows those parameters to vary. A test statistic is compared with an asymptotic chi-square distribution with \(q\) degrees of freedom when the relevant regularity conditions hold.

For example, if clearance for individual \(i\) is modeled as \(CL_i=CL_{\mathrm{typ}}\exp(\beta z_i+\eta_{CL,i})\), then testing whether \(z_i\) is associated with clearance means testing \(H_0:\beta=0\) against \(H_a:\beta\ne 0\). Here \(z_i\) is the covariate, \(\beta\) is its coefficient, and \(\eta_{CL,i}\) is an individual random effect. Under this parameterization, \(\exp(\beta)\) is the multiplicative change in typical clearance for a one-unit increase in \(z_i\), with other model terms held fixed.

Likelihood and Regularity Conditions

Write \(\ell(\psi,\lambda)\) for the log likelihood. In a nonlinear mixed-effects model, the population likelihood is ordinarily a marginal likelihood: it accounts for within-individual observations by integrating over the random effects. In schematic form, for independent individuals with data \(y_i\) and covariates \(x_i\),

\ell(\theta)=\sum_{i=1}^{n}\log\left\{\int f(y_i\mid b_i,x_i;\theta)\,g(b_i;\theta)\,db_i\right\}, \qquad \theta=(\psi,\lambda).

The tests below use this likelihood or a stated approximation to it. Their usual chi-square calibration is an asymptotic result, not a guarantee for every fitted model. It requires, in essence, independent individuals; a correctly specified and identifiable model; a sufficiently smooth likelihood; a true parameter in the interior of the allowed parameter space; and nonsingular information for the parameters being tested after accounting for nuisance parameters. The effective sample size for this approximation is generally the number of independent individuals, not the number of concentration records.

A covariate coefficient equal to zero is often an interior value, but a variance equal to zero is a boundary value. Testing a variance component at zero therefore does not generally satisfy the usual conditions for a standard chi-square test. Similarly, near-nonidentifiability, very sparse covariate groups, or a poorly determined likelihood can make the standard approximations unreliable.

Wald Test: Estimate Divided by Its Uncertainty

The Wald test asks how far the estimated covariate effect lies from its null value, measured in estimated standard errors. For a scalar parameter with null value \(\psi_0\), let \(\hat{\psi}\) be its estimate in the alternative model and let \(SE(\hat{\psi})\) be its standard error from the estimated covariance matrix. The statistic is:

W=\frac{(\hat{\psi}-\psi_0)^2}{SE(\hat{\psi})^2}.

For a \(q\)-dimensional effect, let \(\widehat{\mathrm{Var}}(\hat{\psi})\) be the estimated covariance matrix. Then:

W=(\hat{\psi}-\psi_0)^{T} \left\{\widehat{\mathrm{Var}}(\hat{\psi})\right\}^{-1} (\hat{\psi}-\psi_0).

Under \(H_0\) and the regularity conditions, \(W\) is approximately distributed as \(\chi^2_q\). The covariance matrix matters: a joint Wald test accounts for correlations among the estimated effects rather than testing each coefficient separately.

Worked Example: Wald Test for a Clearance Covariate

Suppose a hypothetical alternative model estimates \(\hat{\beta}=0.180\) for a covariate measured in units of 10 kg, with \(SE(\hat{\beta})=0.0750\). Test \(H_0:\beta=0\) at \(\alpha=0.05\), using the one-degree-of-freedom chi-square approximation.

The Wald statistic is:

W=\frac{(0.180-0)^2}{0.0750^2} =\frac{0.0324}{0.005625} =5.76.

For \(df=1\), the approximate p-value is \(0.0164\). Since \(0.0164<0.05\), the Wald test rejects \(H_0\) at the 5% significance level. A corresponding 95% Wald interval is \(0.180\pm1.96(0.0750)=[0.0330,0.327]\). On the clearance-ratio scale, exponentiating the endpoints gives approximately \([1.034,1.387]\). Under the fitted model, the estimated multiplicative change in typical clearance is \(\exp(0.180)=1.197\), or about 19.7% per 10 kg increase. This is an estimated association conditional on the model, not by itself evidence that the covariate causes the change.

Wald inference is convenient once the alternative model has been fitted and a reliable covariance estimate is available. Its weakness is that it relies on a local, approximately quadratic representation around the alternative estimate. If the likelihood is asymmetric, the estimate is near a parameter boundary, or the standard error is unstable, that local approximation may be poor.

Score Test: Slope at the Null Model

The score test asks whether the likelihood has a sufficiently strong slope in the direction of the alternative when evaluated at the null fit. It requires estimating the null model and evaluating derivatives with respect to the tested effect there; it does not require maximizing the likelihood over the alternative effect. Nuisance parameters must still be estimated under the null.

Let \(U_\psi\) be the score vector, the derivative of \(\ell\) with respect to \(\psi\), evaluated at the null fit. Let \(I\) be the information matrix, partitioned into tested and nuisance parameters. The nuisance-adjusted, or efficient, score and information are:

U_{\psi\cdot\lambda} =U_\psi-I_{\psi\lambda}I_{\lambda\lambda}^{-1}U_\lambda, \qquad I_{\psi\cdot\lambda} =I_{\psi\psi}-I_{\psi\lambda}I_{\lambda\lambda}^{-1}I_{\lambda\psi}.

At the constrained maximum-likelihood estimate under \(H_0\), the nuisance score \(U_\lambda\) is zero, so the efficient score equals \(U_\psi\). The adjustment in the information is still important: it accounts for uncertainty in nuisance parameters that are correlated with the effect being tested. The score statistic is:

S=U_{\psi\cdot\lambda}^{T} I_{\psi\cdot\lambda}^{-1} U_{\psi\cdot\lambda}.

The information may be computed from the expected information or an appropriate observed-information version, provided the chosen method is used consistently and is supported by the fitting approach. Under the null and regularity conditions, \(S\) is approximately \(\chi^2_q\).

Worked Example: Score Test Without an Alternative Fit

In a hypothetical analysis, the null model has been fitted with the covariate effect fixed at zero. At that constrained fit, suppose the nuisance-adjusted score for the scalar effect is \(U_{\psi\cdot\lambda}=3.20\), and its efficient information is \(I_{\psi\cdot\lambda}=2.50\). These are illustrative derivative and information values in the model's parameterization.

The score statistic is:

S=\frac{(3.20)^2}{2.50} =\frac{10.24}{2.50} =4.096.

With \(df=1\), the approximate p-value is \(0.0430\). At \(\alpha=0.05\), this score test rejects the null. The result uses the likelihood's slope and curvature at the null fit; it does not establish that a fully fitted alternative model would produce exactly the same estimate or p-value.

Likelihood-Ratio Test: Compare the Two Fits

The likelihood-ratio test (LRT) compares the maximum log likelihood under the alternative with that under the null. The null model must be nested within the alternative, meaning it can be obtained by imposing the null restrictions. Let \(\hat{\ell}_{A}\) and \(\hat{\ell}_{0}\) be the maximized log likelihoods under the alternative and null, respectively. The statistic is:

LR=2(\hat{\ell}_{A}-\hat{\ell}_{0}).

When the alternative contains \(q\) additional free parameters and the regularity conditions hold, \(LR\) is approximately \(\chi^2_q\) under \(H_0\). Both models must be fitted using comparable likelihood definitions and data. The test uses the improvement in fit after allowing the effect to vary, rather than relying only on the alternative estimate's reported standard error.

Worked Example: Likelihood-Ratio Test for the Same Effect

Suppose the hypothetical null and alternative fits for the clearance covariate give \(\hat{\ell}_{0}=-812.40\) and \(\hat{\ell}_{A}=-809.85\). The alternative adds one covariate coefficient, so \(q=1\).

The statistic is:

$$ LR=2[-809.85-(-812.40)] =2(2.55) =5.10. $$

For \(df=1\), the approximate p-value is \(0.0239\). Since \(0.0239<0.05\), the LRT rejects \(H_0\) at the 5% level. This illustrative result and the Wald result above both reject, although their p-values differ. The difference is plausible because the Wald test uses a local approximation around the alternative estimate, whereas the LRT uses the fitted likelihood improvement between models.

Why the Three Tests Can Agree Asymptotically

Their large-sample relationship can be seen by considering the profile log likelihood \(\ell_p(\psi)\), which maximizes over nuisance parameters for each fixed \(\psi\). Near the null value \(\psi_0\), a quadratic approximation is:

\ell_p(\psi_0+\delta)-\ell_p(\psi_0) \approx U_p\delta-\frac{1}{2}I_p\delta^2,

where \(U_p\) is the profile score at the null and \(I_p\) is the profile information. Maximizing the right-hand side gives \(\hat{\delta}\approx U_p/I_p\). Substituting this value into twice the log-likelihood difference gives \(U_p^2/I_p\), the scalar score statistic. The Wald statistic using this local approximation is also \(\hat{\delta}^2 I_p\), which equals \(U_p^2/I_p\). Thus all three statistics have the same leading quadratic form near the null. For multiple parameters, the corresponding result uses vectors and information matrices. Their agreement is asymptotic; finite-sample equality is not expected.

Worked Example: Joint Wald Test for a Two-Contrast Covariate

Suppose a hypothetical three-level categorical covariate is represented by two coefficients for effects on log clearance. In the alternative fit, let \(\hat{\psi}=(0.30,-0.10)^T\), with estimated covariance matrix:

\widehat{\mathrm{Var}}(\hat{\psi})= \begin{pmatrix} 0.0400 & 0.0100\\ 0.0100 & 0.0225 \end{pmatrix}.

The determinant is \(0.0400(0.0225)-0.0100^2=0.000800\), so the inverse covariance matrix is:

\left\{\widehat{\mathrm{Var}}(\hat{\psi})\right\}^{-1} = \begin{pmatrix} 28.125 & -12.500\\ -12.500 & 50.000 \end{pmatrix}.

For the joint null \(H_0:\psi=(0,0)^T\), the Wald statistic is:

W=(0.30,-0.10) \begin{pmatrix} 28.125 & -12.500\\ -12.500 & 50.000 \end{pmatrix} \begin{pmatrix} 0.30\\ -0.10 \end{pmatrix} =3.781.

The reference distribution has \(df=2\), giving an approximate p-value of \(0.151\). The joint test does not reject at \(\alpha=0.05\). Testing the two coefficients separately would not answer the same joint question; the covariance between their estimates contributes to the result.

Choosing and Diagnosing a Test

Test	Fit required	Main information used	Practical concern
Wald	Alternative model	Estimate and covariance matrix	Can be sensitive to poor standard errors or a non-quadratic likelihood
Score	Null model, plus derivatives for the tested effect	Score and nuisance-adjusted information at the null	Reliable derivatives and information may be difficult to obtain
Likelihood-ratio	Both nested models	Difference in maximized log likelihoods	Requires comparable, successfully optimized fits

The score test can be useful when fitting many alternative models is costly, because it needs only the null fit and derivatives in the tested direction. It does not eliminate the need for a sound null fit or correct nuisance adjustment. The Wald test is often readily available from a fitted model, but a reported covariance matrix can be unreliable when parameters are highly correlated or weakly identified. The LRT directly evaluates the improvement in likelihood, but requires both fits and valid optimization of each.

Several diagnostics should prompt caution rather than automatic reliance on the nominal p-value:

Large changes across tests. Material disagreement can indicate that the likelihood is not well approximated by a quadratic near the null, or that an estimate or covariance is unstable.
Extreme covariance or standard errors. Very large standard errors, high parameter correlations, or a nearly singular covariance matrix can undermine Wald calculations and point to weak identification.
Unstable derivative calculations. Score tests depend on derivatives and information at the null; check that the software can calculate these quantities reliably for the fitted model.
Boundary or sparse-data situations. A boundary parameter, very small covariate subgroup, or nearly empty category can invalidate standard chi-square calibration or make it a poor approximation.
Failed or incomparable fits. An LRT is not meaningful if the fits use different records or likelihood definitions, or if either model has not been adequately optimized.

When these concerns are substantial, report the limitation and investigate the fit, parameterization, and data support before treating a nominal threshold as decisive. The tests quantify evidence against a specified model restriction; they do not establish causality, clinical importance, or predictive value on their own.

Common Mistake

A common mistake is to treat a Wald p-value, a score p-value, and an LRT p-value as interchangeable because all three are compared with a chi-square distribution. The correction is that they are different statistics, based on different parts of the likelihood, and their chi-square calibration is approximate under stated conditions. In particular, a Wald test is not simply an LRT with a different name, and a score test does not require fitting the alternative model. Report which test was used, the restrictions tested, the degrees of freedom, and any important fitting or regularity concerns.

Key takeaway: Wald tests use the fitted effect and its uncertainty; score tests use the null-fit slope and information; likelihood-ratio tests compare the maximized fits. Their large-sample similarity does not remove the need to check parameter uncertainty, model fit, and the assumptions behind chi-square inference.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions.

A scalar covariate estimate is \(0.12\), with standard error \(0.06\). What is its Wald statistic for testing a zero effect, and what is its approximate one-degree-of-freedom p-value?
What model must be fitted for a score test, and why does the test still need information about nuisance parameters?
Two nested models have maximized log likelihoods \(-500.0\) and \(-497.5\), with the latter model adding two parameters. Calculate the LRT statistic and state its reference degrees of freedom.
Why can the Wald, score, and likelihood-ratio tests produce different p-values in a finite sample?

Solutions. (1) \(W=(0.12/0.06)^2=4.00\), which has an approximate \(\chi^2_1\) p-value of \(0.0455\). (2) Fit the null model, then evaluate the score and information for the tested effect there. Nuisance parameters may be correlated with that effect, so the information must adjust for their uncertainty. (3) \(LR=2[-497.5-(-500.0)]=5.00\), with \(df=2\). (4) They use different local or global likelihood information and rely on approximations that need not be equally accurate in finite samples; likelihood asymmetry, uncertain covariance estimates, or weak identification can increase their differences.