Identifiability Of Random Effects On PD Parameters

Why a Small Variance Estimate Can Be Misleading

A population PD model may estimate a random-effect variance close to zero even when subjects truly differ. The estimate might indicate genuinely little between-subject variability, but it can also arise because the observed concentrations do not make individual PD parameters distinguishable, responses barely move across the observed range, or residual noise obscures the differences. These explanations have different implications: collecting more informative observations may resolve weak identification, whereas a well-supported near-zero estimate may justify a simpler model.

The key question is not only whether an estimation method returns a variance estimate. It is whether the design and model can recover a known variance when that variance is present. This tutorial develops that question using sensitivity and simulation-based parameter recovery.

Definition: A random-effect variance is structurally identifiable if distinct values of that variance imply distinct probability distributions for the observations under the model. It is practically identifiable if the available design and amount of data provide enough information to estimate it with useful precision. A variance may be structurally identifiable but practically weakly identified.

Where Information About a PD Variance Comes From

Let \(Y_{ij}\) be the response for subject \(i\) at observation \(j\), and let \(C_{ij}\) be the corresponding drug concentration. A basic population model is

Y_{ij}=f(C_{ij},\theta_i)+\epsilon_{ij}, \qquad \theta_i=\theta_{\mathrm{pop}}\exp(\eta_i), \qquad \eta_i\sim N(0,\omega^2).

Here \(\theta_{\mathrm{pop}}\) is the typical PD parameter, \(\eta_i\) is the subject-specific random effect, and \(\omega^2\) is its variance on the log scale. The exponential parameterization ensures a positive individual parameter. For this model, assume random effects are independent across subjects, residual errors have mean zero and variance \(\sigma^2\), and residual errors are independent across observations conditional on the individual parameters. The concentration values are treated as known for the illustration; uncertainty in estimated concentrations or PK parameters can further reduce information about PD variability.

A local approximation shows why the concentration range matters. Write \(g_{ij}\) for the derivative of the predicted response with respect to the subject’s random effect, evaluated at \(\eta_i=0\). For small random effects, a first-order Taylor expansion gives

f(C_{ij},\theta_{\mathrm{pop}}\exp(\eta_i)) \approx f(C_{ij},\theta_{\mathrm{pop}})+g_{ij}\eta_i.

Because the same \(\eta_i\) affects all of subject \(i\)’s responses, this approximation implies an additional response variance of about \(g_{ij}^2\omega^2\), and a covariance between two responses from that subject of about \(g_{ij}g_{ik}\omega^2\). If the response is insensitive to the parameter at the observed concentrations, \(g_{ij}\) is small: even a nonzero \(\omega^2\) then creates little observable variation. Repeated measurements can help distinguish a shared subject effect from independent residual error, but only if those measurements are informative about the parameter.

Local sensitivity result: When the first-order approximation is adequate, the contribution of a random-effect variance to response variation scales with the squared response sensitivity. A variance is difficult to estimate when the responses are insensitive to the parameter, when residual error is large relative to the between-subject signal, or when the data do not distinguish the random effect from other model features.

Concentration Range, Response Range, and Parameter Confounding

Consider an Emax response model with a baseline \(E_0\), maximum drug effect \(E_{\max}\), and concentration for half-maximum effect \(EC_{50}\):

f(C)=E_0+\frac{E_{\max}C}{EC_{50}+C}.

For a random effect on \(E_{\max}\), using \(\log(E_{\max,i})=\log(E_{\max})+\eta_i\), the derivative with respect to \(\eta_i\) at the typical value is \(E_{\max}C/(EC_{50}+C)\). If all observed concentrations are far below \(EC_{50}\), this derivative is small, so variation in \(E_{\max}\) produces only small response differences. For a random effect on \(EC_{50}\), the corresponding derivative is \(-E_{\max}EC_{50}C/(EC_{50}+C)^2\). Its magnitude is greatest around the transition region, near \(C=EC_{50}\); observations only at very low or very high concentrations give relatively little information about individual \(EC_{50}\) values.

Limited ranges can also make random effects on different parameters difficult to separate. At concentrations far below \(EC_{50}\), the drug effect is approximately \(E_{\max}C/EC_{50}\). The response then primarily reflects the ratio \(E_{\max}/EC_{50}\), not each parameter separately. Variation in either parameter can produce similar changes in the observed responses. A narrow response range, a measurement ceiling or floor, or high residual error can further hide those changes.

Worked Example: Checking Sensitivity at Low Concentration

Suppose the typical Emax model has \(E_{\max}=8\) response units and \(EC_{50}=10\) concentration units. All observations are near \(C=0.2\). Calculate the predicted drug effect and the local sensitivities to log-\(E_{\max}\) and log-\(EC_{50}\).

The predicted drug effect is

\frac{8(0.2)}{10+0.2}=\frac{1.6}{10.2}=0.1569\text{ response units}.

The sensitivity to log-\(E_{\max}\) is the predicted drug effect, \(0.1569\) response units per log unit. The magnitude of the sensitivity to log-\(EC_{50}\) is

\frac{8(10)(0.2)}{(10+0.2)^2} =\frac{16}{104.04}=0.1538\text{ response units per log unit}.

Both sensitivities are small compared with, for example, a residual standard deviation of \(0.5\) response units. They are also nearly equal in magnitude, with opposite signs: at these low concentrations, increasing \(E_{\max}\) and decreasing \(EC_{50}\) have similar effects on the response. The design therefore offers weak information about either random-effect variance and limited separation of the two effects. This calculation is a local diagnostic, not a proof that the fitted model is unidentifiable.

Use Parameter-Recovery Simulations

A parameter-recovery simulation asks whether the planned design and estimation procedure can recover known parameter values. Specify the model and design, simulate many datasets using chosen values of \(\omega^2\) and the other parameters, fit each dataset using the intended estimation procedure, and summarize the estimates. A useful assessment includes average or median bias, the frequency of estimates at the boundary, interval coverage, interval width, and correlations between estimates of the variance and other parameters. The simulation should use realistic concentration schedules, sample sizes, residual error, missingness, and the same fitting settings intended for the analysis.

Recovery under one nonzero variance is not enough. Include a near-zero or zero variance scenario as well as one or more scientifically meaningful nonzero scenarios. If the method repeatedly estimates a substantial known variance as nearly zero, a low estimate in the actual analysis cannot reliably distinguish negligible variability from weak identification. If meaningful variances are recovered accurately and precisely under the same design, a near-zero estimate is more informative, though it still does not prove that the true variance is exactly zero.

Worked Example: Comparing Sparse and Informative Designs

Consider a hypothetical study of 80 subjects. The true log-scale random-effect variance on \(E_{\max}\) is \(\omega^2=0.09\), so the log-scale standard deviation is \(\sqrt{0.09}=0.30\). Residual errors have standard deviation \(0.5\) response units. Compare a sparse design with two observations per subject, at \(C=0.1\) and \(0.3\), with a design having six observations per subject at \(C=0.5, 2, 5, 10, 20,\) and \(30\). Use \(E_{\max}=8\) and \(EC_{50}=5\).

In the sparse design, at \(C=0.3\), the typical drug effect is

\frac{8(0.3)}{5+0.3}=\frac{2.4}{5.3}=0.4528\text{ response units}.

The local sensitivity to log-\(E_{\max}\) is also \(0.4528\). Multiplying by the true log-scale standard deviation gives an approximate between-subject response standard deviation of \(0.4528(0.30)=0.1358\) response units at that concentration. This is small relative to the residual standard deviation of \(0.5\), and both concentrations are below \(EC_{50}\).

Now suppose 200 datasets are simulated and fitted for each design. An illustrative simulation might find that, in the sparse design, the median estimated variance is \(0.018\), corresponding to a log-scale standard deviation of \(\sqrt{0.018}=0.1342\), and 61% of estimates are at the zero boundary. In the six-concentration design, the median estimated variance might be \(0.082\), corresponding to a standard deviation of \(\sqrt{0.082}=0.2864\), with 8% of estimates at the boundary. These hypothetical results show that the sparse design tends to under-recover the known variance, while the wider concentration range gives substantially better recovery. The simulation outputs are illustrative, not empirical findings about a real study.

Distinguishing Weak Identification from Small Variability

The fitted estimate alone cannot make this distinction. Instead, compare how the analysis behaves when the true variance is set to zero and when it is set to meaningful nonzero values. Examine the distribution of estimates and intervals under each scenario, including how often the estimator returns zero. Then locate the actual estimate within those simulated distributions. This is a design-specific calibration: it asks whether the observed result would be unsurprising even if a meaningful variance were present.

Worked Example: Interpreting an Estimate Near Zero

Suppose an analysis estimates \(\omega^2=0.004\), or a log-scale standard deviation of \(\sqrt{0.004}=0.0632\). In 200 hypothetical simulations under the actual sparse design, suppose 126 fits produce an estimate at or below \(0.004\) when the true variance is \(0.09\). Under a zero-variance scenario, suppose 12 of 200 fits produce an estimate at or below \(0.004\).

Under the nonzero scenario, the simulated frequency is \(126/200=0.63\), or 63%. Thus an estimate this small is common even when the true variance is \(0.09\); it provides little evidence for distinguishing that variance from zero under this design. The frequency under the zero scenario is \(12/200=0.06\), or 6%. These frequencies are simulation summaries conditional on the specified scenarios, not probabilities that either scenario is true.

If the same \(0.004\) estimate occurred under a richer design, and simulations showed that meaningful variances were almost always recovered while estimates this small were uncommon under those nonzero scenarios, the result would be more supportive of small variability. It would still be evidence about the specified model and design, not proof that the population variance is exactly zero.

Diagnostics, Alternatives, and Limitations

Warning signs of weak practical identification include frequent boundary estimates across simulated datasets, broad or asymmetric uncertainty intervals, strong dependence of \(\omega^2\) estimates on starting values, and pronounced trade-offs with other parameter estimates. Also check whether the observed concentrations cover the relevant PD transition and whether responses span a useful range rather than clustering near baseline or a ceiling. A successful optimizer termination does not establish that the variance is well identified.

Simulation-based recovery is especially useful because it can reproduce the full design and fitting procedure, including boundary behavior. It is computationally demanding and its conclusions depend on the assumed model, parameter scenarios, and simulation settings. It cannot reveal an omitted biological mechanism simply by recovering parameters from data generated under the same assumed model.

Profile likelihood offers a complementary view: a flat or highly asymmetric likelihood over plausible variance values signals weak information. A standard error alone can be misleading for a variance near the zero boundary, where symmetric, approximately normal uncertainty is often a poor description. Model comparison or a simpler model with the variance fixed at zero can be useful, but a zero estimate should not be treated as proof that the random effect is absent. These approaches answer related but distinct questions; none substitutes for checking whether the design can recover scientifically meaningful variability.

Do not equate a boundary estimate with no variability. Boundary estimates can be common when the data are weakly informative.
Do not simulate only the fitted estimate. Recovery at zero and at plausible nonzero values is needed for comparison.
Do not confuse response noise with individual variation. Repeated responses can help, but only when the PD parameter affects predictions at the sampled concentrations.
Do not overgeneralize a recovery result. It applies to the simulated design, model, and parameter scenarios.

Key takeaway: A small estimated PD random-effect variance is interpretable only in light of the information in the observed concentration and response ranges. Parameter-recovery simulations show whether the design and fitting procedure could have detected meaningful variability if it were present.

Check Your Understanding

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

What is the difference between structural and practical identifiability of a random-effect variance?
For an Emax model, why can observations only at concentrations far below \(EC_{50}\) make random effects on \(E_{\max}\) and \(EC_{50}\) difficult to distinguish?
A recovery simulation uses a true log-scale variance of \(0.16\). What is the corresponding log-scale standard deviation?
In 100 simulations with a meaningful nonzero variance, 70 estimates are at the boundary. What does this say about interpreting a boundary estimate in the actual analysis?

Solutions.

Structural identifiability means distinct variance values imply distinct observation distributions under the model. Practical identifiability means the available design and data estimate that variance with useful precision. A structurally identifiable variance can still be practically weakly identified.
At low concentrations, the drug effect is approximately \(E_{\max}C/EC_{50}\). Increasing \(E_{\max}\) or decreasing \(EC_{50}\) can therefore produce similar response changes, and the effect may be small relative to residual error.
The standard deviation is \(\sqrt{0.16}=0.40\) on the log scale.
The boundary estimate is not, by itself, evidence that variability is absent. Since \(70/100=70\%\) of fits hit the boundary despite meaningful true variability, the design and fitting procedure often fail to recover that variability; a small estimate is therefore ambiguous.