Computing Epsilon Shrinkage from Conditional Residuals

Why Measure Epsilon Shrinkage?

An individual conditional prediction uses a subject’s data to estimate that subject’s random effects. As a result, the fitted prediction can follow the observed data more closely than a population prediction does. The remaining conditional residuals may have less spread than the residual-error model would predict. Epsilon shrinkage summarizes this reduction in spread on a standardized residual scale.

The quantity is useful as a descriptive warning: conditional residuals with little spread may provide limited visual information about residual model adequacy. It does not, by itself, establish that the residual model is wrong or that any particular subject’s estimates are unreliable. The introductory material on conditional distributions behind empirical Bayes estimates distinguishes individual conditional predictions from population predictions; here the focus is specifically on measuring the spread left in residuals after conditioning on individual estimates.

Definition: For observation $j$ from subject $i$, let $y_{ij}$ be the observation, $f_{ij}(\hat\eta_i)$ its individual conditional prediction, and $V_{ij}(\hat\eta_i)$ the modeled residual variance at that prediction. Define the individual conditional weighted residual (IWRES) as $$ z_{ij}=\frac{y_{ij}-f_{ij}(\hat\eta_i)}{\sqrt{V_{ij}(\hat\eta_i)}}. $$ For $N$ observations, let $s_z$ be the sample standard deviation of all $N$ IWRES values, centered on their sample mean and using denominator $N-1$. The standard SD-based epsilon-shrinkage estimate is $$ \widehat{ES}=1-s_z, \qquad \widehat{ES}_{\%}=100(1-s_z)\%. $$

This definition uses a standard-deviation scale: a spread of 1 corresponds to zero estimated epsilon shrinkage. Positive values indicate a spread below 1; zero indicates a spread of 1; and a negative value indicates a spread above 1. Do not replace the formula with $1-s_z^2$: that is a different, variance-scale summary. State which convention is used when reporting the result.

Model the Residual Error Before Standardizing

The numerator of an IWRES is the difference between an observation and its individual conditional prediction. Its denominator must be the residual standard deviation implied by the specified error model at that prediction. Since this denominator can vary across observations, raw residuals generally cannot be compared directly or used in place of IWRES.

Proportional error: On the observation scale, write $y_{ij}=f_{ij}(\eta_i)(1+\epsilon_{ij})$, where $\epsilon_{ij}$ has mean 0 and standard deviation $\sigma_p$. Conditional on the individual prediction $f_{ij}$, the residual standard deviation is $\sigma_p|f_{ij}|$. Thus, $$ z_{ij}=\frac{y_{ij}-f_{ij}(\hat\eta_i)}{\sigma_p|f_{ij}(\hat\eta_i)|}. $$ Combined additive and proportional error: If the residual is the sum of independent mean-zero additive and proportional errors, with standard deviations $\sigma_a$ and $\sigma_p|f_{ij}|$, respectively, then $$ V_{ij}=\sigma_a^2+\sigma_p^2f_{ij}(\hat\eta_i)^2, \qquad z_{ij}=\frac{y_{ij}-f_{ij}(\hat\eta_i)}{\sqrt{\sigma_a^2+\sigma_p^2f_{ij}(\hat\eta_i)^2}}. $$

These expressions assume the residual error has mean zero and the stated variance, the additive and proportional components in the combined model are independent, and the variance parameters are expressed in units consistent with the observations and predictions. The conditional prediction and variance must use the same individual estimates and model. If the residual distribution is not adequately represented by the assumed model, the IWRES spread may reflect model misspecification as well as conditional shrinkage.

Why Conditional Fitting Can Reduce Residual Spread

With correctly specified residual variance and known model parameters, a residual divided by its standard deviation has unit variance under the model. Conditioning on individual estimates changes the situation: those estimates are informed by the same observations whose residuals are then calculated. Some of the observed variation can be absorbed into the individual fit, leaving less variation in the residuals.

A simple linear least-squares calculation illustrates the mechanism. Suppose a vector of $N$ observations follows a linear model with design matrix $X$, $p$ fitted coefficients, and independent errors with variance $\sigma^2$. The fitted residual vector is $e=(I-H)y=(I-H)\epsilon$, where $H=X(X^\mathsf{T}X)^{-1}X^\mathsf{T}$ is the projection matrix. Since $I-H$ has rank $N-p$, $$ E(e^\mathsf{T}e)=\sigma^2\operatorname{tr}(I-H)=\sigma^2(N-p). $$

Consequently, residual sum of squares per observation is expected to be $\sigma^2(1-p/N)$, below $\sigma^2$ when fitted coefficients absorb variation. This is an illustration, not an exact finite-sample formula for nonlinear mixed-effects models: individual conditional estimation, nonlinear predictions, estimated variance parameters, and dependence among observations all affect the residual spread. It explains why a spread below 1 can arise even when conditional fitting is behaving as expected.

Worked Example: Epsilon Shrinkage Under Proportional Error

Consider four hypothetical observations with individual predictions of 10, 20, 15, and 25 concentration units. Suppose the proportional residual standard deviation is $\sigma_p=0.10$, and the individual conditional residuals, in the same units, are 0.5, -1.0, 0.75, and -2.5. Calculate epsilon shrinkage.

For proportional error, the residual standard deviations are $0.10|f_{ij}|$, giving 1, 2, 1.5, and 2.5 concentration units. Dividing each conditional residual by its corresponding standard deviation gives IWRES values $0.5,-0.5,0.5,-1.0$. Their mean is $$ \bar z=\frac{0.5-0.5+0.5-1.0}{4}=-0.125. $$

The sum of squared deviations from the mean is $1.6875$. Therefore, $$ s_z=\sqrt{\frac{1.6875}{4-1}}=\sqrt{0.5625}=0.75, \qquad \widehat{ES}=1-0.75=0.25. $$ The estimated epsilon shrinkage is 25%. This means the pooled IWRES sample standard deviation is 0.75 rather than the reference value 1; it does not mean that 25% of observations have been shrunk or that each subject has the same degree of shrinkage.

Worked Example: Epsilon Shrinkage Under Combined Error

For four hypothetical observations, suppose the individual predictions are 10, 20, 30, and 40 concentration units, with $\sigma_a=2$ concentration units and $\sigma_p=0.10$. Let the conditional residuals be $\sqrt{5}$, $-\sqrt{8}$, $0.5\sqrt{13}$, and $-0.5\sqrt{20}$ concentration units. Calculate epsilon shrinkage using the combined error model.

The conditional residual variances are $2^2+(0.10f)^2$. At predictions 10, 20, 30, and 40, these variances are 5, 8, 13, and 20 concentration-units squared. The corresponding residual standard deviations are $\sqrt{5},\sqrt{8},\sqrt{13},\sqrt{20}$. Dividing each residual by its matching standard deviation gives the IWRES values $1,-1,0.5,-0.5$.

Their mean is 0, and their squared deviations sum to $1+1+0.25+0.25=2.5$. Thus, $$ s_z=\sqrt{\frac{2.5}{4-1}}=\sqrt{0.8333}=0.9129, \qquad \widehat{ES}=1-0.9129=0.0871. $$ The estimated epsilon shrinkage is 8.71%. Using only the proportional component in the denominator would give different IWRES values and a different summary; the additive component matters, particularly when predictions are small.

Interpreting the Summary and Checking Its Calculation

Epsilon shrinkage is a pooled descriptive statistic, so the set of residuals included in the calculation must be stated. In the definition above, all observations are pooled, each contributes one IWRES, and the sample standard deviation is centered on the pooled mean. A large number of observations from one subject can therefore contribute more to the pooled result than a small number from another subject. The statistic does not give a subject-by-subject measure.

Worked Example: A Negative Estimate Is Not a Calculation Error

Suppose four hypothetical standardized conditional residuals are $-1.5,-0.5,0.5,1.5$. Their mean is 0 and their squared deviations sum to $2.25+0.25+0.25+2.25=5$. Therefore, $$ s_z=\sqrt{\frac{5}{4-1}}=1.291, \qquad \widehat{ES}=1-1.291=-0.291. $$ The estimated epsilon shrinkage is $-29.1\%$. The residual spread is above the unit reference, so this estimate indicates no shrinkage by this measure; it flags excess spread instead. It should not be reset to zero. Possible explanations include residual-model misspecification, influential observations, or ordinary sampling variability.

Before interpreting a value, check the calculation and the residual behavior:

Verify the scale. Confirm that each residual was divided by the model-implied residual standard deviation at that observation. Under proportional error, the denominator changes with the prediction; under combined error, include both variance components.
Inspect the mean as well as the spread. Epsilon shrinkage uses the standard deviation around the sample mean. A nonzero mean can indicate bias or a systematic prediction error that the shrinkage value alone does not describe.
Look for patterns across prediction, time, and subject. A pooled spread can hide residual trends, changing variability, or a few influential observations. Plotting and stratifying residuals can expose problems that one number cannot.
Check what the software calls the residual. IWRES use individual conditional predictions. Residuals standardized using population predictions, often given a different label such as CWRES, are not interchangeable with IWRES for this calculation. Software conventions can differ, so verify the definition rather than relying only on a label.
Consider the number and distribution of observations. A small or unbalanced dataset can make a pooled SD unstable or dominated by a subset of subjects. The estimate is a descriptive sample quantity, not a precise correction for the residual model.

Alternatives, Limitations, and a Common Pitfall

The main practical alternative is to examine other residual definitions. Population-based weighted residuals can be useful for assessing population-level predictions, while individual conditional residuals describe what remains after subject-specific fitting. Raw residuals retain the observation’s units but are not on a comparable scale when residual standard deviations vary. These quantities answer different diagnostic questions; substituting one for another changes the interpretation of the spread.

Epsilon shrinkage compresses an entire pooled residual distribution to one number. It does not reveal why the spread differs from 1, distinguish variance misspecification from conditional fitting, or establish that residuals are independent or normally distributed. It is also not a test with a universal cutoff: interpretation depends on the model, data design, residual definition, and amount of information per subject. Examine residual plots and model assumptions alongside the summary.

Common mistake: A reader sees a residual SD of 0.75 and concludes that the residual error variance should be multiplied by 0.75. This confuses a standardized residual spread with an adjustment to the model’s residual variance. The correct calculation is $1-0.75=0.25$, or 25% SD-based epsilon shrinkage; it is a descriptive measure, not an automatic instruction to alter the fitted error model.

Key takeaway: Compute epsilon shrinkage from the pooled sample SD of individual conditional residuals standardized by the correct conditional residual SD. Report the residual definition and scaling, retain negative estimates, and treat the result as a diagnostic summary rather than a stand-alone verdict on model adequacy.

Check Your Understanding

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

Under proportional error with $\sigma_p=0.08$, a prediction is 25 concentration units and its conditional residual is 1 concentration unit. What is its IWRES?
For combined error, explain why the residual standard deviation is not generally $\sigma_a+\sigma_p|f|$.
A pooled set of IWRES has sample SD 0.60. Calculate SD-based epsilon shrinkage and state what the result describes.
Why should a negative epsilon-shrinkage estimate not be changed to zero?
What distinction matters when choosing between individual conditional residuals and residuals based on population predictions?

Solutions. (1) The proportional residual SD is $0.08(25)=2$ concentration units, so the IWRES is $1/2=0.5$. (2) Independent variances add: the variance is $\sigma_a^2+\sigma_p^2f^2$, so the SD is its square root, not the sum of the component SDs. (3) $1-0.60=0.40$, or 40%; the pooled standardized conditional residuals have SD 0.60. (4) A negative estimate records a spread above 1 and may signal excess variability or sampling variation; clipping it hides that information. (5) Individual conditional residuals use subject-specific estimates, whereas population-based residuals use population predictions; they represent different conditioning and are not interchangeable.