Tutorials › Stochastic Differential Equations in Pharmacokinetics › A Model-Choice Workshop with Pharmacokinetic Scenarios

Choosing Between ODEs, SDEs and Mixed-Effects Models · Tutorial 3 of 100

A Model-Choice Workshop with Pharmacokinetic Scenarios

Learn to match a pharmacokinetic model to where variability arises: in fixed individual parameters, across patients, or within a patient over time.

Beginner 12 min read

What You'll Learn

  • Distinguish a fixed individual parameter from a population distribution of parameters.
  • Identify when a stochastic state represents evolving within-patient dynamics.
  • Classify PK scenarios and justify the simplest suitable model.
  • Calculate population concentration summaries under a simple mixed-effects model.
  • Interpret a multiplicative SDE example and its limitations.
  • Recognize when the available sampling design cannot distinguish competing explanations.

Start With the Source of Variation

A pharmacokinetic model should represent the source of variation relevant to the question. If one patient's kinetics are fixed and the aim is to describe their concentration trajectory, a deterministic ordinary differential equation (ODE) may be sufficient. If patients differ in fixed kinetic parameters, a mixed-effects model can represent that population heterogeneity. If a patient's latent PK state itself changes stochastically over time, a stochastic differential equation (SDE) may be appropriate.

These choices are about the model for the underlying state, not simply whether observed concentrations scatter around a curve. Sparse measurements can make distinct sources of variation difficult to distinguish, as discussed in the earlier material on sparse sampling and variability sources. A model choice therefore needs both a scientific justification and a check against what the data can identify.

Definition: A deterministic ODE model specifies a state trajectory from its initial condition and parameters; conditional on those inputs, the trajectory is fixed. A mixed-effects model allows parameters to vary across individuals, commonly by assigning them a population distribution, while each individual's parameter is fixed over the modeled course. An SDE model includes a stochastic term in the state dynamics, so the state can evolve differently across realizations even with the same initial condition and fixed parameters.

For example, a one-compartment elimination model can be written \(dC_i(t)/dt=-k_iC_i(t)\), where \(C_i(t)\) is the concentration for individual \(i\) and \(k_i\) is the elimination-rate constant. With a fixed \(k_i\) and initial concentration \(C_i(0)=C_{i0}\), its solution is \(C_i(t)=C_{i0}e^{-k_it}\). The model class depends on how \(k_i\) and the evolving state are understood: a single fixed value gives an ODE; a population distribution of fixed \(k_i\)'s gives a mixed-effects model; a stochastic state evolution calls for an SDE.

A Decision Rule With Explicit Assumptions

Use the following questions in order. They do not replace scientific judgment; they make the judgment explicit.

1
Is the target one individual's expected or reference trajectory?
If a fixed-parameter trajectory answers the question, begin with an ODE. Specify the initial condition, dosing input, and parameters.
2
Does the question concern differences among individuals?
If individuals have different but course-stable kinetic parameters drawn from a population distribution, use a mixed-effects model. The distribution describes the population; a realized parameter describes one individual's kinetics.
3
Is there a scientific reason to model state fluctuations within an individual?
If the latent concentration or another PK state is hypothesized to evolve stochastically, an SDE can represent that evolving state. A random-looking set of observations alone is not enough to establish this explanation.

In an illustrative SDE, \(dC(t)=-kC(t)\,dt+\sigma C(t)\,dW_t\), \(W_t\) denotes standard Brownian motion, and \(\sigma\) controls the magnitude of the multiplicative stochastic term. The equation represents randomness in the state dynamics. It is different in purpose from assigning different fixed \(k_i\)'s to different patients. Neither model, by itself, specifies how assay measurement error should be modeled.

Result: Suppose that, conditional on fixed parameters \(\theta_i\), an individual's state follows a deterministic ODE with a unique solution on the time interval of interest. Then at each fixed time \(t\), the conditional variance of that state is zero: \(\operatorname{Var}(C_i(t)\mid\theta_i,C_i(0))=0\).

Proof. Under the stated uniqueness assumption, fixing \(\theta_i\) and \(C_i(0)\) fixes one solution value \(c(t;\theta_i,C_i(0))\). Thus the conditional random variable \(C_i(t)\) is equal to that constant with probability one. By the definition of variance,

$$ \operatorname{Var}(C_i(t)\mid\theta_i,C_i(0)) =E\!\left[(C_i(t)-E[C_i(t)\mid\theta_i,C_i(0)])^2\mid\theta_i,C_i(0)\right] =0. $$

This result does not say that a population of patients has zero variance. If their parameters differ, their state values can differ across patients even though each conditional trajectory is deterministic. \(\square\)

Worked Scenarios

Worked Example: One Patient With Fixed Elimination

A patient receives an IV dose that establishes \(C(0)=10\) mg/L. For this illustration, assume the elimination rate is fixed at \(k=0.15\) per hour and there are no time-varying inputs after the dose. Choose a model and calculate the predicted concentration at 4 and 8 hours.

The question concerns one patient's trajectory, and the assumptions specify a fixed rate with no within-patient stochastic state term. The one-compartment ODE is \(dC/dt=-0.15C\), with \(C(0)=10\). Its solution is \(C(t)=10e^{-0.15t}\).

$$ C(4)=10e^{-0.60}=5.488116\ldots\ \text{mg/L}, \qquad C(8)=10e^{-1.20}=3.011942\ldots\ \text{mg/L}. $$

The ODE is appropriate for the stated target: a reference trajectory conditional on the specified rate. It does not claim that all patients share this rate, or that real concentrations are measured without error. If the clinical question were instead about variation among patients, the fixed-\(k\) assumption would need reconsideration.

Worked Example: Two Patient Types in a Population

Suppose patients have the same \(C(0)=10\) mg/L but differ in a fixed elimination rate. In a deliberately simple population model, \(K_i=0.10\) per hour with probability \(1/2\), and \(K_i=0.20\) per hour with probability \(1/2\). Calculate the population mean and variance of concentration at 5 hours, assuming each patient follows \(C_i(t)=10e^{-K_it}\).

The rate is random across patients but fixed for a given patient. This is population heterogeneity, so a mixed-effects formulation is suitable. At 5 hours the two possible concentrations are

$$ c_1=10e^{-0.10(5)}=10e^{-0.5}=6.065306597\ldots, \qquad c_2=10e^{-0.20(5)}=10e^{-1}=3.678794412\ldots\ \text{mg/L}. $$

Each value occurs with probability \(1/2\). Therefore,

$$ E[C_i(5)] =\frac{c_1+c_2}{2} =\frac{6.065306597\ldots+3.678794412\ldots}{2} =4.872050504\ldots\ \text{mg/L}. $$

For two equally likely values, the variance is the average squared distance from their mean. Equivalently, it is one quarter of the squared difference:

$$ \operatorname{Var}(C_i(5)) =\frac{(c_1-4.872050504\ldots)^2+(c_2-4.872050504\ldots)^2}{2} =\left(\frac{c_1-c_2}{2}\right)^2 =(1.193256093\ldots)^2 =1.423860\ldots\ (\text{mg/L})^2. $$

The variance here is across patients. Conditional on either rate, each patient's modeled concentration at 5 hours is fixed. Thus the population spread does not require stochastic within-patient dynamics.

Worked Example: A Stochastically Evolving Concentration State

Consider the model \(dC(t)=-0.20C(t)\,dt+0.10C(t)\,dW_t\), with \(C(0)=10\) mg/L. The question is specifically whether a latent concentration state may fluctuate within a patient, not merely whether patients have different fixed rates. Explain the model choice and calculate the mean and variance at 2 hours under this model.

The diffusion term makes this an SDE. For this geometric SDE, the solution is \(C(t)=C(0)\exp[(-k-\sigma^2/2)t+\sigma W_t]\). Since \(W_t\sim N(0,t)\), \(C(t)\) is lognormal. Its mean and variance are

$$ E[C(t)]=C(0)e^{-kt}, \qquad \operatorname{Var}(C(t))=C(0)^2e^{-2kt}\left(e^{\sigma^2t}-1\right). $$

At \(t=2\), \(k=0.20\), and \(\sigma=0.10\), the mean is

$$ E[C(2)]=10e^{-0.4}=6.703200\ldots\ \text{mg/L}. $$

The variance is

$$ \operatorname{Var}(C(2)) =100e^{-0.8}(e^{0.02}-1) =100(0.449328\ldots)(0.0202013\ldots) =0.9077\ldots\ (\text{mg/L})^2. $$

With \(\sigma=0\), this model reduces to the deterministic ODE and the variance is zero. With \(\sigma=0.10\), it assigns nonzero conditional variation to the state at a fixed time. The SDE is justified only if such within-patient state evolution is scientifically relevant and supported; it should not be selected merely because observations do not fall exactly on an ODE curve.

How the Model Choices Differ

A compact comparison helps prevent a change of model class from being used to answer the wrong question.

ModelWhat varies?Typical targetKey limitation
ODENothing stochastic in the modeled state, conditional on inputs and parametersA reference trajectory for specified kineticsDoes not represent population parameter variation or stochastic state evolution on its own
Mixed-effects modelParameters across individuals; each individual's realized parameter is fixed during the modeled coursePopulation distribution and individual differencesDoes not by itself represent within-individual stochastic state evolution
SDEThe state evolves stochastically over timeWithin-individual random dynamics in a latent stateNeeds a mechanistic rationale and enough information to distinguish it from alternatives
Result: If \(\theta_i\) is a random individual parameter and \(C_i(t)\) has finite second moment, then $$ \operatorname{Var}(C_i(t)) =E[\operatorname{Var}(C_i(t)\mid\theta_i)] +\operatorname{Var}(E[C_i(t)\mid\theta_i]). $$ For a deterministic ODE conditional on \(\theta_i\), the first term is zero, so the population variance is the variance of the individual-specific deterministic predictions.

Derivation. Let \(m(\theta_i)=E[C_i(t)\mid\theta_i]\). Add and subtract this conditional mean inside the squared deviation from the population mean:

$$ C_i(t)-E[C_i(t)] =(C_i(t)-m(\theta_i))+(m(\theta_i)-E[C_i(t)]). $$

Squaring and taking expectations gives two squared terms and a cross term. The cross term has expectation zero, because

$$ E[(C_i(t)-m(\theta_i))(m(\theta_i)-E[C_i(t)])] =E[(m(\theta_i)-E[C_i(t)])E[C_i(t)-m(\theta_i)\mid\theta_i]] =0. $$

The first squared term has expectation \(E[\operatorname{Var}(C_i(t)\mid\theta_i)]\); the second is \(\operatorname{Var}(m(\theta_i))\). This proves the stated decomposition. If the conditional state is deterministic, its conditional variance vanishes by the previous result. \(\square\)

Diagnostics, Pitfalls, and Limits

A useful diagnostic is to ask whether the model's source of variation matches the scientific story. For a proposed ODE, inspect whether persistent patient-to-patient differences remain after accounting for covariates. For a mixed-effects model, check whether its individual parameters are being interpreted as stable over the modeled period. For an SDE, ask whether the stochastic term is intended to represent changing latent dynamics rather than fixed differences among patients.

Sparse sampling can limit these checks: several models may produce similar predictions at the few observed times. A model that fits those observations is not automatically evidence for its assumed source of variation. Consider the sampling schedule, the clinical mechanism, repeated measurements, and whether parameters can be estimated with useful precision. These are diagnostics for model adequacy, not automatic tests that uniquely identify a model class.

  • Pitfall: treating all scatter as evidence for an SDE. The wrong reasoning is that concentrations vary, so the state must fluctuate randomly. The correction is to ask whether the variation is within-patient state evolution, between-patient parameter heterogeneity, or another feature of the observation process.
  • Pitfall: calling a random effect a stochastic trajectory. A rate drawn once for each patient can create different trajectories across patients without ongoing stochastic changes within any one trajectory.
  • Pitfall: choosing the most complex model by default. An SDE adds assumptions and parameters. Prefer an ODE or mixed-effects model if it addresses the question and is adequate; add stochastic state dynamics only when they are scientifically meaningful.

These models are not mutually exclusive in every application: population parameters and within-patient stochastic dynamics can both matter. This workshop focuses on identifying the primary feature motivating the choice. The model should also state separately how measured concentrations relate to the latent state; the details of that distinction are treated in the later material on process noise and observation error.

Check Your Understanding

For each question, choose or interpret a model, then compare your reasoning with the brief solutions.

  1. Every patient follows a one-compartment curve with their own elimination rate, and each rate stays constant during follow-up. Which model class represents this population structure?
  2. A single patient's latent concentration is hypothesized to receive continuous stochastic perturbations during follow-up. Which model class is a natural candidate, and why is an ODE alone insufficient for that hypothesis?
  3. A deterministic ODE gives a smooth curve, while observed concentrations scatter around it. Does this fact alone establish that the latent state follows an SDE?
  4. In a mixed-effects model, what does the population distribution of a rate parameter describe, and what does one patient's realized rate describe?

Solutions. (1) A mixed-effects model: rates vary across patients but are fixed within a patient's modeled course. (2) An SDE, because the hypothesis concerns stochastic state evolution; a deterministic ODE has a fixed trajectory conditional on its inputs and parameters. (3) No. The scatter alone does not identify its source; assess the scientific mechanism and competing explanations. (4) The distribution describes between-patient heterogeneity; a patient's realized rate describes that individual's kinetics and is fixed in the stated model.

Clinical Biostats