Why Assertions Matter in an Executable QSP Model
An executable QSP model can produce plausible-looking trajectories while using inconsistent units, entering an impossible state, or evaluating a mechanism outside its intended domain. Runtime assertions make selected requirements explicit and check them during execution. They do not prove that a model is biologically correct, but they can reveal violations that would otherwise be difficult to notice.
An assertion is a condition that the implementation expects to hold at a specified point in a run. For example, a concentration must have units of amount per volume, a drug amount may be required to remain nonnegative, and a bioavailability parameter must lie between zero and one. The condition should be stated in the model’s terms, checked at a meaningful point, and linked to a response when it fails. A bare check such as “value is not negative” is insufficient if the variable, units, tolerance, and failure handling are unclear.
Dimensional Assertions
A dimensional assertion checks that quantities combined in an equation have compatible physical dimensions. In a pharmacometric model, units often include amount, volume, time, mass, and concentration. The equation must be coherent both in its definition and in its implementation: matching units in a written model cannot protect against a variable accidentally stored in a different unit in code.
Consider a central-compartment amount \(A\) measured in mg, volume \(V\) in L, concentration \(C=A/V\) in mg/L, clearance \(CL\) in L/h, and input rate \(R\) in mg/h. The amount balance is dimensionally consistent because \(CL\,C\) has units of amount per time.
The derivative \(dA/dt\) must have units mg/h, as do both terms on the right-hand side. A useful implementation can attach units to variables or store unit metadata and check equations during model setup or execution. If the software uses plain numeric values, explicit conversion and validation at input boundaries are especially important. A numeric value of 3 might represent 3 L/h or 3 mL/h; the number alone does not determine its meaning.
Worked Example: Checking a Central-Compartment Rate
Suppose a hypothetical model has \(A=10\) mg, \(V=5\) L, \(CL=3\) L/h, and \(R=8\) mg/h. Check the concentration and amount rate, then consider what a dimensional assertion should catch if clearance is mistakenly supplied in mg/h.
First calculate concentration:
The elimination rate is \(CL\,C=(3\ \mathrm{L/h})(2\ \mathrm{mg/L})=6\) mg/h. Thus:
The units agree, and the calculated rate is positive. If the supplied clearance were labeled mg/h instead of L/h, multiplying it by concentration would give units mg2/(L·h), not mg/h. A unit-aware assertion should reject that multiplication or the parameter at input. Merely checking that the computed rate is a finite number would not detect the dimensional defect.
Dimensional checks also help with parameter transformations and derived quantities. For example, an exponent must be dimensionless: if \(k\) is a first-order rate constant in h\(^{-1}\), then \(kt\) is dimensionless when time \(t\) is in hours. A common failure is to pass time in minutes while retaining a rate expressed per hour. The resulting value may be numerically smooth and finite, yet represent a different exposure duration. Assertions should therefore check units and conversion conventions where quantities enter the model, not only whether later outputs appear reasonable.
Positivity and Boundary Conditions
Many QSP states represent amounts, concentrations, cell counts, or fractions and have restricted domains. For a state \(x\) required to be nonnegative, the admissible set is \(x\geq 0\). For multiple states, the admissible region may be the nonnegative orthant, where every component is nonnegative. Whether that region is preserved depends on the model’s differential equations and inputs, not on the desired interpretation alone.
For a continuous-time model \(dx/dt=f(x,t)\), a standard sufficient condition for the nonnegative orthant to be forward invariant is that the vector field points into or along the region at every boundary face. Specifically, for each component \(i\), require \(f_i(x,t)\geq 0\) whenever \(x_i=0\) and all other components are nonnegative. This condition is paired with regularity assumptions, such as local Lipschitz continuity in the state, that give a unique local solution. Inputs and parameters must also satisfy the assumptions used to establish the boundary condition. For time-varying systems, the boundary condition must hold throughout the time interval of interest.
The boundary test is not a claim that the derivative must be strictly positive at zero. It may equal zero, and a zero derivative at one instant does not by itself settle whether a trajectory can later leave the region. The invariance conclusion follows from the boundary condition together with the regularity assumptions, not from an argument that any crossing would require a strictly negative derivative at the crossing.
A short comparison argument explains the result. Write \(x^+\) for the vector formed by replacing negative components of \(x\) with zero, and define \(x_i^-=\max(0,-x_i)\). Let \(V(x)=\frac12\sum_i(x_i^-)^2\), which is zero exactly when all components are nonnegative. On a bounded part of a trajectory, local Lipschitz continuity supplies a constant \(L\) such that \(|f_i(x,t)-f_i(x^+,t)|\leq L\|x-x^+\|\). If \(x_i<0\), then \(x_i^+=0\), so the boundary condition gives \(f_i(x^+,t)\geq0\), and consequently \(f_i(x,t)\geq-L\|x^-\|\). At differentiable times, only negative components contribute to the derivative of \(V\), yielding
where \(m\) is the number of state components. Gronwall’s inequality implies \(V(t)\leq V(0)\exp(2L\sqrt{m}\,t)\). If the initial state is nonnegative, then \(V(0)=0\), so \(V(t)=0\) throughout the interval. This argument does not require a strictly positive derivative at the boundary. It does require the stated regularity and boundary conditions; it is not a justification for ignoring negative values produced by a numerical algorithm.
Worked Example: Testing a Receptor-Fraction Boundary
Consider a hypothetical receptor occupancy \(P\), constrained to \([0,1]\), with
Assume \(C\geq0\), \(k_{\mathrm{on}}\geq0\), and \(k_{\mathrm{off}}\geq0\). At the lower boundary \(P=0\), the rate is \(k_{\mathrm{on}}C\geq0\), so the vector field does not point below zero. At the upper boundary \(P=1\), the rate is \(-k_{\mathrm{off}}\leq0\), so it does not point above one. The same inward-boundary logic applies to the interval \([0,1]\), with the sign reversed at its upper face.
For a numerical boundary check, take \(C=5\) nM, \(k_{\mathrm{on}}=0.1\) nM\(^{-1}\)h\(^{-1}\), and \(k_{\mathrm{off}}=0.2\) h\(^{-1}\). At \(P=0\), \(dP/dt=(0.1)(5)=0.5\) h\(^{-1}\). At \(P=1\), \(dP/dt=-0.2\) h\(^{-1}\). Both signs are consistent with the allowed interval. If a run uses negative concentration or negative rate constants, this check no longer establishes invariance: the input assumptions themselves have been violated and should be asserted separately.
Parameter Ranges and Boundary Assertions
Parameter assertions encode declared domains, not universal truths about every possible model. For example, a fraction parameter \(F\) may be required to satisfy \(0\leq F\leq1\), a volume \(V\) to satisfy \(V>0\), and a first-order rate constant \(k\) to satisfy \(k\geq0\). Some models use signed parameters or transformed parameterizations, so the valid range must match the model definition and the intended use. A check should be applied after any transformation if it is the transformed quantity that enters the mechanism.
Boundary assertions test behavior where a mechanism reaches a limiting state or input: zero amount, full occupancy, zero concentration, or a parameter at an allowed endpoint. They should evaluate the model’s predicted rate or update at that boundary and compare its direction with the specified domain. Such tests complement conservation and limiting-case checks: they focus on whether the executable model respects its admissible region at its edges.
A runtime assertion should state how much numerical tolerance is allowed. For a state expected to be nonnegative, a small negative value can arise from floating-point arithmetic or a solver step, while a materially negative value may indicate a defect. The tolerance should be defined in the state’s units and scaled appropriately to its magnitude. Do not silently replace every negative value with zero: clipping can hide an incorrect equation or solver behavior, change mass balance, and make an invalid run appear acceptable. Record the failed check and investigate before deciding whether a documented correction is justified.
Worked Example: Separating a Model Boundary Failure From a Numerical Overshoot
Consider a hypothetical amount \(A\geq0\) governed by \(dA/dt=-kA\), with \(k=0.8\) h\(^{-1}\), no input, and initial amount \(A(0)=1\) mg. The continuous model has \(dA/dt=0\) at \(A=0\), so its boundary condition is satisfied. Its exact solution at \(t=2\) h is:
Now apply one explicit Euler step of length \(h=2\) h from \(A=1\) mg. The update is \(A_{\mathrm{new}}=A+h(-kA)=1+(2)(-0.8)(1)=-0.6\) mg. A nonnegativity assertion fails even though the continuous model’s boundary condition is valid. This result points first to a numerical-step or implementation issue, not evidence that the continuous equations permit negative amounts. The assertion has done its job by exposing an inadmissible computed state; the next step is to investigate the calculation and its assumptions, not to silently clip the result.
Interpreting Failures and Choosing Checks
An assertion failure is a diagnostic signal, not a root-cause label. The investigation should distinguish at least four possibilities: a coding error (such as a sign or index mistake), an input or parameter outside its declared domain, a numerical artifact, or a model assumption that is not appropriate for the situation. A unit mismatch often points to a coding or data-interface problem; negative states may arise from the model equations, a numerical update, or an invalid input. The assertion itself usually cannot decide among these explanations.
A practical set of assertions checks model inputs and parameters before integration, state domains and finite values during integration, and boundary rates or updates at specified edges. Include units and variable names in failure messages, and preserve enough information to reproduce the failing condition. A check that fails should normally stop or clearly mark the run as invalid rather than allow downstream predictions to be interpreted as if they passed.
Assertions have limits. Checking a few boundary points does not prove that every possible trajectory is valid, and passing dimensional checks does not establish biological plausibility. A boundary condition may hold only under assumptions about inputs, rates, or state dependencies. A violated assumption should be documented and assessed for the intended use, rather than relabeled automatically as a software defect. Likewise, a numerical tolerance that is too loose can conceal a material failure, while one that is too strict can report harmless rounding as an error.
Compared with assertions, a test of a known analytical solution checks numerical outputs against a reference trajectory, while a conservation check tests a specified balance across states. Boundary assertions answer a narrower question: does the implementation respect declared domains and directional behavior at their edges? These checks are complementary. None replaces careful model specification or the interpretation of whether its assumptions are suitable.
Common Pitfalls
- Inferring invariance from a single derivative value. A zero derivative at the boundary is not, by itself, a complete argument. Use the inward-pointing boundary condition with the required regularity assumptions and an invariance argument.
- Treating every negative state as proof of a wrong mechanism. A numerical update can leave the admissible region even when the continuous equations satisfy the boundary condition. Check both the equations and the numerical implementation.
- Checking only that values are finite. Finite values can still have incompatible units or lie outside the permitted parameter range.
- Clipping failed states without recording the event. Silent clipping can mask a defect and alter model balances. Report and investigate failures before applying any explicit, justified correction.
- Assuming the assertion identifies the cause. A failed condition says what requirement was violated, not whether the implementation, input, numerical method, or modeling assumptions are responsible.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions.
- A state amount is stored in mg, and its volume is stored in L. What are the units of the concentration computed as amount divided by volume?
- For \(dX/dt=u-kX\), where \(X\geq0\), \(u\) is a nonnegative input rate, and \(k\geq0\), evaluate the boundary condition at \(X=0\).
- A continuous model satisfies its nonnegative boundary condition, but a computed state is \(-0.02\) mg after one update. Does that alone show the model equations are wrong? Name two issues to investigate.
- A fraction parameter \(F\) is specified to lie in \([0,1]\), but the supplied value is \(1.04\). What does the assertion establish, and what does it not establish?
Solutions. (1) Concentration has units mg/L. (2) At \(X=0\), \(dX/dt=u\geq0\), so the vector field does not point below the boundary under the stated input assumption. (3) No. Investigate the numerical update and implementation, and check whether inputs or other assumptions were valid. (4) The declared parameter-range condition is violated. The assertion does not determine whether the cause is a coding or data-interface error, a deliberately out-of-range scenario, or an unsuitable model-domain assumption.