Mapping Simulated States to Model Observations

From Simulated States to Measured Quantities

A mechanistic model evolves biological states over time, but experiments usually report quantities defined by measurement procedures. A model may track several molecular pools separately while an assay reports their total, or track drug amount while a laboratory reports concentration. The observation function connects these two descriptions: it takes the simulated state at a specified time and returns the quantity the model predicts would be measured.

This mapping is conceptually separate from the ODE system. The ODE determines how the state changes; the observation function interprets the state. Keeping those roles separate makes it easier to preserve the model’s biological balances and to change the measurement definition without rewriting the dynamics.

Definition: Let $x(t)$ be the vector of model states produced by an ODE system, and let $\theta$ denote its parameters. An observation function for measurement $k$ is a specified mapping $h_k$ from the relevant state and parameter values to a predicted observable: $$ \hat{y}_k(t)=h_k(x(t),\theta). $$ The symbol $\hat{y}_k(t)$ denotes the model-predicted observation, not an additional ODE state. The mapping may also use known quantities such as a sampling volume or a measurement-specific covariate when required.

For example, an ODE might track free and bound receptor as two states, while an assay reports total receptor. The observation function adds the two states at the time of measurement. It does not replace them with a total state or add the total to the ODE vector. If the model has $m$ states, its state vector remains $m$-dimensional unless the biological model itself calls for another dynamic state.

Define the Measurement Before Coding the Map

A useful observation function begins with a precise description of what the assay reports. Identify which states contribute, whether the output is a sum, ratio, or unit conversion, and what units the result should have. If only some model states represent the assayed material, include only those states. A mathematically convenient combination is not automatically the quantity measured in the experiment.

For a sum of compatible pools $x_1,\ldots,x_r$, the predicted total is

\hat{y}_{\mathrm{total}}(t)=\sum_{j=1}^{r}x_j(t).

This sum has the same units as its component states. Adding quantities with incompatible units, such as an amount and a concentration, is not meaningful unless one is first converted to the units of the other.

For a ratio, both numerator and denominator must correspond to the intended measurement, and the denominator must be nonzero. For instance, the fraction of a receptor pool that is bound may be calculated as bound receptor divided by total receptor. A ratio of quantities with the same units is dimensionless; that does not make it a percentage until it is multiplied by 100.

Key distinction: A state is a quantity whose time evolution is specified by the ODE system. An observation is a quantity calculated from one or more states at a time of interest. Use an observation function for a derived measurement; add a new state only when the model represents a quantity with its own dynamic evolution.

Worked Example: Summing Receptor Pools

Suppose a hypothetical model tracks free receptor $R$ and ligand-bound receptor $RL$, both in nM. At a sampling time, the simulated values are $R=42\ \mathrm{nM}$ and $RL=18\ \mathrm{nM}$. An assay reports total receptor across these two pools.

The observation function is $$ \hat{y}_{\mathrm{total}}=R+RL. $$

Substituting the simulated values gives $$ \hat{y}_{\mathrm{total}}=42\ \mathrm{nM}+18\ \mathrm{nM}=60\ \mathrm{nM}. $$

Thus, the model predicts a total receptor concentration of $60\ \mathrm{nM}$ at that time. The sum is appropriate because both states represent receptor and have the same concentration units. The ODE still tracks $R$ and $RL$ separately; the calculated total is not written back into either state.

Ratios and Fractions

Ratios often express occupancy, the fraction of a total pool in a particular form, or the relative abundance of two species. The numerator and denominator should be stated explicitly, rather than inferred from a variable name. For receptor occupancy, if $RL$ is bound receptor and $R$ is free receptor, then the total receptor pool is $R+RL$, and the bound fraction is

\hat{f}_{\mathrm{bound}}(t)=\frac{RL(t)}{R(t)+RL(t)}.

This function is defined when $R(t)+RL(t)\ne 0$. If both pools are zero, the bound fraction is undefined, not automatically zero. If the total is positive and the pools are nonnegative, the result lies between zero and one. Those conditions follow from the meaning of the pools; they should not be silently assumed if the model permits states outside that range.

Worked Example: Calculating Receptor Occupancy

Use the hypothetical receptor values $R=42\ \mathrm{nM}$ and $RL=18\ \mathrm{nM}$. Calculate the bound fraction and express it as a percentage.

The total receptor concentration is $42+18=60\ \mathrm{nM}$, so $$ \hat{f}_{\mathrm{bound}}=\frac{18\ \mathrm{nM}}{60\ \mathrm{nM}}=0.300. $$

The units cancel because numerator and denominator are both concentrations. As a percentage, the predicted bound fraction is $0.300\times100=30.0\%$. This interpretation requires a positive total receptor concentration; here the denominator is $60\ \mathrm{nM}$, so the ratio is defined.

A ratio can also be biologically meaningful when it is not constrained to the interval from zero to one. For example, a concentration ratio between two distinct species may exceed one. The range depends on the definition of the numerator, denominator, and underlying states; do not impose fraction-specific bounds on every ratio.

Concentrations and Sampling Times

Many models track an amount because amount is convenient for mass balance, while an assay measures concentration. For an amount $A(t)$ in a volume $V(t)$, the concentration observation is

\hat{C}(t)=\frac{A(t)}{V(t)}.

The units follow from the division. For example, amount in nmol divided by volume in L gives nmol/L. If a volume changes over time, use the volume corresponding to the same measurement time. If the reported units differ, include the conversion factor explicitly and verify that the final units match the assay output.

An ODE solver may return states at internal integration times that do not coincide with experimental sampling times. The observation function should be evaluated using the model state at each requested sampling time. A solver that provides continuous output can evaluate its numerical solution there; alternatively, a numerical method may interpolate between saved states. The observation is then calculated from the evaluated state, not from whichever saved time happens to be nearest.

Worked Example: Sampling a Concentration from an Amount

A hypothetical model tracks drug amount $A$ in a central compartment. At a requested sampling time of 3 h, the solver’s continuous numerical solution gives $A(3\ \mathrm{h})=72\ \mathrm{nmol}$. The compartment volume at that time is $V(3\ \mathrm{h})=3.0\ \mathrm{L}$. Calculate the predicted concentration.

The observation function is $\hat{C}(t)=A(t)/V(t)$, so $$ \hat{C}(3\ \mathrm{h})= \frac{72\ \mathrm{nmol}}{3.0\ \mathrm{L}} =24\ \mathrm{nmol/L}. $$

The model therefore predicts $24\ \mathrm{nmol/L}$ at 3 h. The time identifies where the state and volume are evaluated; it is not a factor in the unit conversion. The calculation uses the numerical state at the requested time and leaves the amount state unchanged.

Why the Mapping Does Not Change the ODE

The separation can be stated directly. Suppose the state satisfies

\frac{dx}{dt}=f(t,x,\theta), \qquad x(t_0)=x_0.

Once this system has been solved, the observation at time $t$ is obtained by applying $h$ to the resulting state. Calculating $h(x(t),\theta)$ does not alter the right-hand side $f$, the initial condition $x_0$, or the solution $x(t)$. This is why an algebraic total can be reported without replacing the component states that participate in the ODE.

For a differentiable observation function with no explicit time dependence, the chain rule gives a useful justification for this distinction:

\frac{d}{dt}h(x(t),\theta) = \nabla_x h(x(t),\theta)^{\mathsf T}f(t,x(t),\theta).

This derivative describes how the derived observation changes because the underlying states change. It does not make $h$ an independent state: its time course is determined by $x(t)$. If a measured quantity has its own delayed or lagged response, then a static observation map may not represent that biology; a dynamic state or other model structure may be needed.

A practical implementation should therefore read the required state elements, calculate the observation, and return the result separately. It should not overwrite entries in the state vector, append algebraic outputs as if they were dynamically integrated, or alter the state while generating predictions at multiple sampling times.

Checks, Alternatives, and Limitations

Before using a mapping, check that its ingredients match the measurement definition. Confirm the included states, their units, and the intended sampling time. For ratios, verify the denominator is nonzero where the observation is requested. For amount-to-concentration conversions, verify the volume basis and unit conversion. As an additional check, evaluate the map at a simple state with an interpretable result, such as equal pools for a fraction or a known amount and volume.

There are two common alternatives to an algebraic observation function, but each represents a different modeling choice. One can choose state variables that directly correspond to measured quantities, provided those states still support the required mechanistic dynamics and balances. Or one can introduce a dynamic measurement process, for example when an assay signal has its own lag. Neither choice is interchangeable with a static map: a directly measured state changes the state representation, and a dynamic measurement process adds modeled evolution. For an instantaneous sum, ratio, or concentration conversion, an observation function is usually the simplest representation.

A static mapping also has limits. It does not by itself represent assay noise, bias, detection limits, sample processing, or a time lag between the biological state and the recorded value. Those features require additional specifications in the broader model or analysis. The observation function should be no more elaborate than the measurement process requires, but it should not omit a transformation that defines what the assay reports.

Common pitfall: If an assay reports total receptor, it may seem convenient to replace the free-receptor state with the total and then add the bound state to it. That counts the bound pool twice. Define the reported total as the sum of the distinct contributing states, $R+RL$, and retain both states unchanged in the ODE system.

Other warning signs include outputs with unexpected units, ratios computed with the wrong total, and predictions evaluated at a nearby saved time rather than the requested sample time. These are reasons to revisit the observation definition and implementation. A correctly specified map is transparent: each output can be traced to its contributing states, operation, units, and time of evaluation.

Check Your Understanding

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

A model tracks free and bound ligand in $\mathrm{nM}$, with values $7.5\ \mathrm{nM}$ and $2.5\ \mathrm{nM}$. What observation function gives total ligand, and what is its predicted value?
A hypothetical occupancy measurement is $RL/(R+RL)$, with $R=24\ \mathrm{nM}$ and $RL=6\ \mathrm{nM}$. Calculate the fraction and percentage, and state whether the ratio is defined.
A model predicts an amount of $45\ \mathrm{mg}$ in a volume of $1.5\ \mathrm{L}$ at the sample time. What concentration does the observation function return, in $\mathrm{mg/L}$?
Why should an algebraic total generally not be written back into the state vector after it is calculated?

Solutions. (1) Use $\hat{y}_{\mathrm{total}}=x_{\mathrm{free}}+x_{\mathrm{bound}}$. The prediction is $7.5+2.5=10.0\ \mathrm{nM}$. (2) The total is $24+6=30\ \mathrm{nM}$, so the fraction is $6/30=0.200$, or $20.0\%$. It is defined because the denominator is positive. (3) The predicted concentration is $45\ \mathrm{mg}/1.5\ \mathrm{L}=30\ \mathrm{mg/L}$. (4) The total is derived from existing states and has no independent ODE evolution in this model. Writing it back can overwrite a component or otherwise change the state that determines the dynamics; calculate and return it separately instead.