Worked Population PD Model-Building Analysis

From a Concentration-Response Model to a Population Model

A population pharmacodynamic (PD) analysis must describe both the average concentration-response relationship and the variation among individuals. A curve that fits the typical subject may still misrepresent a subgroup, while a covariate effect that improves fit may be implausible or poorly supported by the data. A useful model-building analysis therefore makes its assumptions explicit, estimates between-subject variability, evaluates prespecified covariates, and checks predictions against observed responses at the population level.

The worked analysis below uses entirely hypothetical data. Forty-eight subjects, with ages 40, 50, or 60 years, each have continuous responses measured at concentrations $C=0,\ 0.5,\ 2,$ and $8$ mg/L. The concentration is treated as known. The analysis considers age as a covariate on maximum drug effect and sex as a covariate on baseline response. These covariates are prespecified on biological grounds; model comparisons do not replace that rationale.

Model and assumptions: For subject $i$ at observation $j$, let $Y_{ij}$ be the continuous response at concentration $C_{ij}$. The candidate model is $$ Y_{ij}=E_0+E_{\max}\exp\left\{\beta_{\text{age}}\frac{\text{age}_i-50}{10}+\eta_i\right\} \frac{C_{ij}}{EC_{50}+C_{ij}}+\epsilon_{ij}, $$ where $E_0$ is baseline response, $E_{\max}$ is the median-subject maximum effect at age 50, $EC_{50}$ is the concentration producing half of that subject's maximum effect above baseline, and $\beta_{\text{age}}$ is the age effect per decade on the log scale of $E_{\max}$. The random effect satisfies $\eta_i\sim N(0,\omega^2)$; observation errors satisfy $\epsilon_{ij}\sim N(0,\sigma^2)$, independently across observations and independently of $\eta_i$. The model assumes a common $EC_{50}$, no age effect on $E_0$, and no concentration error or unmodeled time dependence. The concentration-response function is the Emax model.

The lognormal random effect ensures each subject's $E_{\max}$ is positive. Because $E(\eta_i)=0$, $E_{\max}$ is the median individual maximum effect at the reference age, not the arithmetic mean across subjects. That distinction matters whenever the analysis compares model predictions with arithmetic subgroup means.

What Is Estimated and How?

The fixed effects are $E_0$, $E_{\max}$, $EC_{50}$, and any included covariate coefficients. The variance components are $\omega^2$, the between-subject variance on the log-$E_{\max}$ scale, and $\sigma^2$, the within-subject observation-error variance. For a given subject, the repeated observations share $\eta_i$, so they are not independent after integrating over that subject's random effect.

In a likelihood-based fit, the contribution from subject $i$ is obtained by integrating the conditional response density over the random-effect distribution:

L_i=\int\left\{\prod_j f(Y_{ij}\mid\eta_i)\right\} \phi(\eta_i;0,\omega^2)\,d\eta_i, \qquad L=\prod_i L_i.

Here $f$ is the normal density implied by the observation model and $\phi$ is the normal density of $\eta_i$. The integral accounts for uncertainty about each subject's unobserved random effect. It is not equivalent to fitting a separate ordinary regression to each subject's data and treating the resulting effects as known.

Worked Example: Interpreting the Random-Effects Estimate

Suppose the fitted age model estimates $\omega^2=0.0900$. The standard deviation of $\eta_i$ is $\omega=\sqrt{0.0900}=0.300$ on the log scale. The corresponding coefficient of variation in individual $E_{\max}$ values is $$ \sqrt{\exp(\omega^2)-1} =\sqrt{\exp(0.0900)-1} =\sqrt{1.0942-1} =0.3069. $$

Thus, under the lognormal random-effects assumption, the between-subject coefficient of variation in $E_{\max}$, conditional on age, is about 30.7%. This describes variation in the individual maximum-effect parameter, not a 30.7% standard deviation in observed responses: response variability also depends on concentration and includes the separate observation-error component.

Evaluating Prespecified Covariates

The analysis starts with a base model that includes the random effect on $E_{\max}$, but sets $\beta_{\text{age}}=0$ and includes no sex effect. The age model adds the prespecified age coefficient. To isolate the sex question, a separate sex model adds a sex coefficient on $E_0$ to the base model. Each comparison is between nested models fitted to the same observations with the same likelihood and random-effects structure, apart from the covariate being tested.

For a nested comparison, the likelihood-ratio statistic is $2(\ell_1-\ell_0)$, where $\ell_1$ and $\ell_0$ are the maximized log likelihoods of the larger and smaller models. Under the null hypothesis that the added coefficient is zero, and under the usual regularity conditions, the statistic is compared with a chi-square distribution with degrees of freedom equal to the number of added parameters. Here each comparison adds one coefficient, so $df=1$. This approximation can be unreliable for variance parameters on a boundary; the tests here concern fixed covariate coefficients, not whether $\omega^2=0$.

Worked Example: Testing the Prespecified Age Effect

Suppose the hypothetical base-model estimates are $E_0=80.1$ response units, $E_{\max}=39.6$ response units, $EC_{50}=2.10$ mg/L, and $\omega^2=0.100$. The age model estimates $E_0=80.0$ response units, $E_{\max}=40.0$ response units, $EC_{50}=2.00$ mg/L, $\omega^2=0.0900$, and $\beta_{\text{age}}=-0.120$ per decade. Its observation-error standard deviation is $\sigma=7.80$ response units. Suppose the fitted log likelihoods give a likelihood-ratio statistic of 15.20 for adding age.

The hypotheses are $H_0:\beta_{\text{age}}=0$ and $H_a:\beta_{\text{age}}\ne0$. With one added parameter, $$ P(\chi^2_1\ge 15.20)=9.67034\times10^{-5}. $$ At $\alpha=0.05$, the data provide evidence against $H_0$ under the model and test assumptions. The estimated coefficient also has a direct interpretation: a 10-year increase in age multiplies the median $E_{\max}$ by $$ \exp(-0.120)=0.8869. $$ That is an estimated 11.3% lower median $E_{\max}$ per decade in this dataset, conditional on the model. It is not evidence that changing an individual's age causes a change in drug response.

The statistical result alone does not settle whether the model should be retained. The direction and size of the estimate should be considered alongside biological plausibility, data support across ages, and subgroup prediction checks. Those checks assess whether the covariate improves predictions where the data are informative, rather than merely improving a single fit statistic.

Worked Example: Comparing the Prespecified Sex Effect

Suppose adding sex on $E_0$ to the base model estimates a difference of $+1.10$ response units for females relative to males, with a likelihood-ratio statistic of 0.400. The comparison tests $H_0:\beta_{\text{sex}}=0$ against a two-sided alternative and adds one parameter: $$ P(\chi^2_1\ge0.400)=0.5271. $$

At $\alpha=0.05$, this comparison does not provide sufficient evidence to distinguish the sex-effect model from the base model. It does not show that the true sex effect is exactly zero. The estimate is small relative to the residual response scale, but uncertainty and the study's ability to detect effects still matter. In this hypothetical analysis, the sex term is not retained in the final parameterization because its estimated contribution is limited and the likelihood comparison does not support its addition. The decision should be reported with that qualification, not as proof of equivalence.

Typical Predictions and Population Means

A prediction at $\eta_i=0$ is a typical-subject prediction, because zero is the median of the normal random-effect distribution and corresponds to the median individual $E_{\max}$. A population arithmetic mean instead averages over the random effects. For this model, the identity $E[\exp(\eta_i)]=\exp(\omega^2/2)$ gives the marginal mean response at age $a$ and concentration $C$:

E[Y\mid a,C]=E_0+E_{\max} \exp\left\{\beta_{\text{age}}\frac{a-50}{10}+\frac{\omega^2}{2}\right\} \frac{C}{EC_{50}+C}.

This expression assumes the stated lognormal random-effects distribution and the model's mean-zero observation error. It also shows why the typical-subject curve is not generally the population mean curve: for positive $\omega^2$, the factor $\exp(\omega^2/2)$ exceeds one. This derivation uses the normal moment-generating function: if $\eta_i\sim N(0,\omega^2)$, then $E[\exp(\eta_i)]=\exp(0+\omega^2/2)$.

Worked Example: Checking Age-Subgroup Predictions

Use the age-model estimates $E_0=80.0$ response units, $E_{\max}=40.0$ response units, $EC_{50}=2.00$ mg/L, $\beta_{\text{age}}=-0.120$ per decade, and $\omega^2=0.0900$. At $C=2.00$ mg/L, the concentration fraction is $2/(2+2)=0.500$. For age 40, the population mean prediction is $$ 80.0+40.0\exp\{(-0.120)(-1)+0.0900/2\}(0.500) =80.0+20.0\exp(0.165) =103.6\ \text{response units}. $$ For age 60, it is $$ 80.0+40.0\exp\{(-0.120)(1)+0.0900/2\}(0.500) =80.0+20.0\exp(-0.0750) =98.55\ \text{response units}. $$

The corresponding typical-subject predictions, obtained by setting $\eta_i=0$, are $80.0+20.0\exp(0.120)=102.6$ response units at age 40 and $80.0+20.0\exp(-0.120)=97.74$ response units at age 60. These are different targets. Comparing the typical predictions directly with arithmetic subgroup means would fail to account for the estimated between-subject variability.

For a prediction check, suppose the hypothetical arithmetic subgroup means at $C=2.00$ mg/L are 104.0 response units for the 16 subjects aged 40 and 97.9 response units for the 16 subjects aged 60. The population predictions differ from these observed means by $103.6-104.0=-0.400$ and $98.55-97.9=0.650$ response units, respectively. Both subgroup means are close to their corresponding model predictions in this example. Such agreement is reassuring but does not establish that the model is correct or that it will predict well in other populations.

Putting the Analysis Together

The full hypothetical dataset also includes subjects aged 50, for whom the population mean prediction at $C=2.00$ mg/L is $$ 80.0+20.0\exp(0.0450)=100.9\ \text{response units}. $$ At $C=0$, the model predicts $E_0=80.0$ response units for every age because age was placed on $E_{\max}$, not baseline. At $C=8.00$ mg/L, the fraction is $8/(2+8)=0.800$, so the age-40 population mean is $80.0+32.0\exp(0.165)=117.7$ response units. A model-building check should compare predictions and observations across the concentration range and within covariate subgroups, not only at a single concentration. The earlier material on population and individual prediction checks and covariate-stratified visual predictive checks provides the relevant prediction-checking frameworks.

A concise report of this hypothetical result could state: “A one-compartment-style Emax concentration-response model with a lognormal subject-specific random effect on $E_{\max}$ was fitted to continuous responses. The estimated typical $E_{\max}$ at age 50 was 40.0 response units, $EC_{50}$ was 2.00 mg/L, and $E_0$ was 80.0 response units. Between-subject variance on log $E_{\max}$ was 0.0900, corresponding to a 30.7% coefficient of variation. Adding age on log $E_{\max}$ improved fit (likelihood-ratio statistic 15.20, $df=1$, $p=9.67\times10^{-5}$); the estimated age coefficient was -0.120 per decade. Adding sex on baseline was not supported in the separate comparison (likelihood-ratio statistic 0.400, $df=1$, $p=0.527$). At 2.00 mg/L, marginal predictions were 103.6 response units at age 40 and 98.55 response units at age 60, close to the hypothetical subgroup means.” This reports the reference value, units, variability scale, and prediction type so that the parameterization is interpretable.

Diagnostics, Alternatives, and Limitations

A covariate model can fail even when its likelihood improves. Check that the concentration-response shape is compatible with the data, that the estimated parameters are stable, and that subgroup predictions do not show systematic discrepancies across the observed concentration range. Assess whether each age or sex subgroup has enough observations and adequate concentration coverage. A coefficient estimated mainly by extrapolation, or one that changes sharply when a few subjects are removed, should be treated cautiously. Simulation checks can help reveal instability, but they cannot compensate for weak data support or incorrect structural assumptions.

The analysis also depends on the chosen random-effects structure. A random effect on $E_{\max}$ permits individual differences in maximal effect, but does not represent individual variation in $EC_{50}$ or baseline. The common-$EC_{50}$, Gaussian-error, and age-effect assumptions may be inadequate for other data. Alternative approaches include evaluating a different structural concentration-response model, placing a random effect on a different parameter when supported, or using a different observation distribution when responses require it. These are not interchangeable fixes: each changes the scientific assumptions and should be supported by the design and data.

A common mistake is to set $\eta_i=0$, calculate a typical-subject prediction, and compare it with the arithmetic mean observed in a subgroup. The two quantities differ whenever the random effect is nonzero and acts nonlinearly. For the lognormal $E_{\max}$ model, the correction is to integrate over $\eta_i$, which multiplies the effect term by $\exp(\omega^2/2)$. Another mistake is to interpret a nonsignificant covariate test as proof of no effect. The test instead indicates that the data and model did not provide sufficient evidence for an improvement under that comparison; the estimate, uncertainty, design, and scientific context remain relevant.

Key takeaway: A population PD model report must identify which predictions are typical-subject predictions and which are marginal population means. Estimate random-effect variability, test prespecified covariates in a defined nested comparison, and check predictions on the same scale as the observations being summarized.

Check Your Understanding

Answer each question using the model and hypothetical estimates above. Brief solutions follow the questions.

What does $E_{\max}=40.0$ response units represent when the random effect is lognormal?
At concentration 2.00 mg/L, why is the population mean prediction at age 40 not the same as the prediction with $\eta_i=0$?
Interpret $\exp(-0.120)$ in the age model.
What conclusion is supported by the sex-effect likelihood-ratio result $0.400$ with $df=1$?
Name one limitation of the age model that should be checked before extrapolating predictions.

Solutions. (1) It is the median individual maximum effect at age 50, conditional on the model, not the population arithmetic mean of $E_{\max}$. (2) The population mean integrates over $\eta_i$ and includes the factor $\exp(\omega^2/2)=\exp(0.0450)$; the $\eta_i=0$ prediction does not. (3) Each 10-year increase in age multiplies median $E_{\max}$ by 0.8869, an estimated 11.3% decrease per decade. (4) The comparison does not provide sufficient evidence at $\alpha=0.05$ to add the sex effect; it does not prove the effect is zero. (5) For example, the model assumes a common $EC_{50}$ across ages and does not include age-related change in baseline; extrapolation beyond the observed age range is especially uncertain.