Tutorials › Pharmacometric Workflows and Reproducible Modeling › Reviewing Residual Error Implementation

Quality Control and Independent Review · Tutorial 333 of 399

Reviewing Residual Error Implementation

Learn to translate residual-error equations into model code and verify that their implied means, variances, and units match the stated specification.

Intermediate 12 min read

What You'll Learn

  • Distinguish additive, proportional, and combined residual-error models by how their variance changes with the prediction.
  • Check whether a residual draw is centered at zero and whether its scale has the correct units.
  • Verify the absolute-value scale required when a proportional-error standard deviation is written explicitly.
  • Derive the combined variance and identify the independence assumption behind it.
  • Compare component-wise simulation code with an equivalent direct variance implementation.
  • Use code-level checks to locate common residual-error translation mistakes.

Why Residual-Error Code Needs Its Own Review

A model can have the intended structural equations and random-effects implementation yet still produce the wrong observation distribution if residual error is translated incorrectly. A misplaced scale, a nonzero-centered draw, or an unintended dependence between error components can change the conditional mean or variance without changing the model’s apparent structure. Reviewing residual-error implementation therefore means checking the equation, the code, and the distribution implied by that code against one another.

The notation here uses \(F\) for the model prediction conditional on the relevant model parameters and subject effects, and \(Y\) for the observation. The residual is \(Y-F\). The review concerns the observation model conditional on \(F\); it does not assess whether the structural model predicts well.

Definition: A residual-error model specifies the conditional distribution of \(Y\) given the model prediction \(F\). Its conditional mean is \(E(Y\mid F)\), and its conditional variance is \(\operatorname{Var}(Y\mid F)\). An implementation is consistent with a stated residual-error model when its code implies the stated conditional distribution, including its mean, variance, units, and required dependence assumptions.

Three Standard Forms and Their Variance Behavior

In an additive model, error has a constant scale in the units of the observation. Let \(\sigma_a\geq 0\) denote its standard deviation. The standard form is:

$$ Y=F+\epsilon_a,\qquad \epsilon_a\mid F\sim N(0,\sigma_a^2). $$

It follows that \(E(Y\mid F)=F\) and \(\operatorname{Var}(Y\mid F)=\sigma_a^2\). Thus the conditional standard deviation does not change with \(F\). If code draws \(Z\sim N(0,1)\), an equivalent implementation is \(Y=F+\sigma_a Z\). A code expression that adds a centered normal draw with standard deviation \(\sigma_a\) has the intended form; adding a draw with variance \(\sigma_a\), rather than \(\sigma_a^2\), does not.

In a proportional model, the residual scale changes in proportion to the magnitude of the prediction. Write the proportional parameter as \(\sigma_p\geq 0\), dimensionless when \(Y\) and \(F\) have the same units:

$$ Y=F(1+\epsilon_p),\qquad \epsilon_p\sim N(0,\sigma_p^2). $$

Equivalently, \(Y=F+\epsilon\), where \(\epsilon\mid F\sim N(0,F^2\sigma_p^2)\). Therefore:

$$ E(Y\mid F)=F,\qquad \operatorname{Var}(Y\mid F)=F^2\sigma_p^2,\qquad SD(Y\mid F)=|F|\sigma_p. $$

The absolute value matters. When predictions are nonnegative, as is usual for concentration predictions, the standard deviation can be written \(F\sigma_p\). For a prediction that might be negative, a standard deviation written explicitly must instead be \(|F|\sigma_p\). The equation \(Y=F+F\sigma_p Z\), with \(Z\sim N(0,1)\), is still valid for signed \(F\): multiplying a centered normal draw by a negative coefficient reverses its sign but leaves its variance \(F^2\sigma_p^2\). It is not valid to pass a negative \(F\sigma_p\) as a normal standard deviation.

A combined additive and proportional model includes both scales. In the usual form, the two error components are conditionally independent, centered normal variables:

$$ Y=F+\epsilon_a+\epsilon_p,\qquad \epsilon_a\mid F\sim N(0,\sigma_a^2),\qquad \epsilon_p\mid F\sim N(0,F^2\sigma_p^2). $$

Under those assumptions, the conditional mean is \(F\), and the conditional variance is the sum of the component variances:

$$ \operatorname{Var}(Y\mid F)=\sigma_a^2+F^2\sigma_p^2,\qquad SD(Y\mid F)=\sqrt{\sigma_a^2+F^2\sigma_p^2}. $$

The additive part dominates near \(F=0\), while the proportional part increasingly contributes as \(|F|\) grows. The combined model is not generally implemented by adding the two standard deviations: independent variances, not standard deviations, add.

A Short Justification for the Combined Variance

For any two random variables, the variance of their sum equals the sum of their variances plus twice their covariance. Conditional on \(F\), the combined residual is \(\epsilon_a+\epsilon_p\). If the two components are independent, their covariance is zero. Consequently:

$$ \begin{aligned} \operatorname{Var}(Y\mid F) &=\operatorname{Var}(\epsilon_a+\epsilon_p\mid F)\\ &=\operatorname{Var}(\epsilon_a\mid F)+\operatorname{Var}(\epsilon_p\mid F) +2\operatorname{Cov}(\epsilon_a,\epsilon_p\mid F)\\ &=\sigma_a^2+F^2\sigma_p^2. \end{aligned} $$

Independence is sufficient for the covariance to be zero; zero covariance alone is enough for this variance identity, but does not in general make the joint distribution independent. If the components are correlated, a covariance term must be included. If they are jointly normal, the sum is still normal, with that covariance included in its variance.

Worked Example: Checking an Additive Implementation

Suppose a hypothetical observation model states that concentrations follow an additive residual error with \(\sigma_a=2.4\) mg/L. At a record where \(F=18\) mg/L, code generates \(Z\sim N(0,1)\) and sets \(Y=F+2.4Z\). Check the implied mean, standard deviation, and variance.

The added error is \(\epsilon_a=2.4Z\). Since \(E(Z)=0\) and \(\operatorname{Var}(Z)=1\), its mean is \(2.4(0)=0\) mg/L and its variance is \((2.4)^2(1)=5.76\) (mg/L)\(^2\). Thus:

$$ E(Y\mid F)=18+0=18\ \text{mg/L},\qquad SD(Y\mid F)=2.4\ \text{mg/L}. $$

The result matches the stated additive model: the prediction is the conditional mean, and the error standard deviation is constant in concentration units. This conclusion assumes the code’s \(Z\) is a standard normal draw for the observation and is not given a nonzero mean or a scale of its own.

Worked Example: Checking a Proportional Scale, Including Its Sign

Consider a hypothetical proportional model with \(\sigma_p=0.15\). For a nonnegative prediction \(F=80\) mg/L, code implements \(Y=F+F\sigma_p Z\), where \(Z\sim N(0,1)\). Find the conditional standard deviation and variance, then determine the correct scale if a prediction were \(F=-5\) mg/L.

At \(F=80\) mg/L, the multiplier on \(Z\) is \(80(0.15)=12\) mg/L. The error is centered at zero, and:

$$ SD(Y\mid F)=|80|(0.15)=12\ \text{mg/L},\qquad \operatorname{Var}(Y\mid F)=80^2(0.15)^2=144\ (\text{mg/L})^2. $$

For \(F=-5\) mg/L, the proportional residual standard deviation is \(|-5|(0.15)=0.75\) mg/L, and its variance is \((-5)^2(0.15)^2=0.5625\) (mg/L)\(^2\). The expression \(F\sigma_p=-0.75\) mg/L is a signed coefficient, not a valid standard deviation. The code form \(F\sigma_p Z\) remains a valid centered residual draw, because its variance is \((F\sigma_p)^2\). If the code instead requests a normal draw with a standard deviation argument, it must use \(|F|\sigma_p\).

Worked Example: Verifying a Combined Model

A hypothetical concentration model specifies independent additive and proportional components with \(\sigma_a=2.0\) mg/L and \(\sigma_p=0.10\). At \(F=30\) mg/L, code generates independent \(Z_a,Z_p\sim N(0,1)\) and sets \(Y=F+2.0Z_a+0.10FZ_p\). Compute the implied conditional variance and standard deviation.

The additive component has variance \(2.0^2=4.00\) (mg/L)\(^2\). The proportional component has standard deviation \(|30|(0.10)=3.0\) mg/L and variance \(30^2(0.10)^2=9.00\) (mg/L)\(^2\). Because the two standard normal draws are independent, the components are independent and their covariance is zero. Therefore:

$$ \operatorname{Var}(Y\mid F)=4.00+9.00=13.00\ (\text{mg/L})^2,\qquad SD(Y\mid F)=\sqrt{13.00}=3.606\ \text{mg/L}. $$

A code expression using one standard normal draw, \(Y=F+\sqrt{13.00}Z\), would give the same conditional distribution at this fixed \(F\). It would not, however, represent two separately generated component draws. Adding the component standard deviations to obtain \(2.0+3.0=5.0\) mg/L would be wrong under the stated independence assumption.

A Practical Code-Review Sequence

A reviewer should translate the code into a conditional distribution rather than relying on the name assigned to an error parameter. The following sequence makes the comparison explicit.

1
Write the stated equation.
Record the prediction term, each residual component, the distribution of each component, and any stated dependence assumption.
2
Identify the actual random draw.
Determine whether the code draws a standard normal variate, a normal variate with an explicit standard deviation, or a variance-scaled quantity. Check its centering and scale.
3
Derive the implied moments.
Calculate \(E(Y\mid F)\) and \(\operatorname{Var}(Y\mid F)\) from the code. For multiple components, include their covariance unless independence or another zero-covariance condition is justified.
4
Check units and parameter constraints.
An additive standard deviation has the units of \(Y\); a proportional standard deviation parameter is dimensionless when \(Y\) and \(F\) share units. Scales and standard deviations must be nonnegative.
5
Compare across prediction values.
Evaluate the implied variance at \(F=0\), at a typical positive prediction, and—if the model can produce them—at negative predictions. Confirm that the behavior matches the stated error form.

The principal alternatives are not interchangeable. An additive model has constant variance; a proportional model has variance increasing as \(F^2\); a combined model has both a nonzero additive variance at \(F=0\) and a prediction-dependent contribution. A combined model can also be coded directly as \(Y=F+\sqrt{\sigma_a^2+F^2\sigma_p^2}Z\), with \(Z\sim N(0,1)\). At a fixed \(F\), this has the same normal conditional distribution as the independent-component form. The two implementations need not produce identical simulated values when separate components are generated, but their conditional distributions agree under the stated assumptions.

Diagnostics, Pitfalls, and Limits

Code-level warning signs include a residual draw added without a zero-centering check, a scale multiplied by the prediction twice, a proportional parameter treated as an observation-scale quantity, or an explicit normal standard deviation that can become negative. Another common mismatch is using \(\sigma_a+\sigma_pF\) as a combined standard deviation while the specification calls for independent additive and proportional components. That expression does not yield the stated variance sum and can even become negative.

A particularly important mistake is to write “the proportional standard deviation is \(F\sigma_p\)” without qualifying the sign. That is only a nonnegative standard deviation when \(F\geq 0\). The general conditional standard deviation is \(|F|\sigma_p\). Alternatively, the signed coefficient \(F\sigma_p\) may multiply a centered standard normal draw; the resulting variance is still \(F^2\sigma_p^2\).

The variance calculation alone does not verify every feature of an observation model. For example, an incorrectly centered error can preserve its variance while shifting \(E(Y\mid F)\); reusing a random draw can induce dependence across records; and correlated components change the combined variance. A normal error model also permits simulated values outside physically plausible ranges. Whether that behavior is acceptable depends on the model specification; changing or truncating the distribution is a modeling decision, not a code-review correction to make silently.

A useful review record states the written equation, the corresponding code expression, the derived conditional mean and variance, and any assumptions needed for equivalence. This makes it possible to distinguish a true implementation discrepancy from two different code expressions that imply the same conditional distribution.

Key takeaway: For every residual-error implementation, derive \(E(Y\mid F)\) and \(\operatorname{Var}(Y\mid F)\) from the code. Check centering, units, nonnegative scales, the absolute-value interpretation of proportional standard deviation, and the independence or covariance assumptions behind combined error.

Check Your Understanding

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

  1. A model uses \(Y=F+3Z\), where \(Z\sim N(0,1)\). What are the conditional mean and standard deviation, and what units must the coefficient 3 have?
  2. For \(F=-4\) and \(\sigma_p=0.25\), what are the proportional residual standard deviation and variance? Is \(F\sigma_p\) a valid standard deviation?
  3. A combined model has \(\sigma_a=1.5\) mg/L and \(\sigma_p=0.20\), with independent components. At \(F=10\) mg/L, calculate its conditional variance and standard deviation.
  4. Why can two independent residual components not generally be combined by adding their standard deviations?
  5. What additional term is required in the combined variance if the additive and proportional components have nonzero conditional covariance?

Solutions. (1) Since \(E(Z)=0\) and \(\operatorname{Var}(Z)=1\), \(E(Y\mid F)=F\) and \(SD(Y\mid F)=3\); the coefficient has the units of \(Y\). (2) The standard deviation is \(|-4|(0.25)=1\), and the variance is \((-4)^2(0.25)^2=1\) in squared observation units. \(F\sigma_p=-1\) is not a valid standard deviation, though it may multiply a standard normal draw. (3) The variance is \(1.5^2+10^2(0.20)^2=2.25+4.00=6.25\) (mg/L)\(^2\), so the standard deviation is \(2.5\) mg/L. (4) Independent variances add, so the combined standard deviation is the square root of the sum of squared component standard deviations, not their sum. (5) Add \(2\operatorname{Cov}(\epsilon_a,\epsilon_p\mid F)\) to \(\sigma_a^2+F^2\sigma_p^2\).

Clinical Biostats