Why Review Mechanistic Model Code?
A mechanistic model is expressed both as a scientific specification and as an implementation that a computer can execute. Review asks whether those two representations agree. A model can be well documented and still have a coding error: a state may be updated using the wrong parameter, a concentration may be substituted where an amount is required, or an equation may be implemented with a missing factor. Such errors can produce plausible-looking outputs, so review should examine the model’s meaning term by term rather than relying on whether a simulation runs.
This tutorial focuses on peer review of model implementation: equation-to-code correspondence, state and parameter mapping, units, and readability. It is not a substitute for evaluating whether the scientific model is appropriate for a particular use. Earlier material on validation and qualification addresses evidence about predictive performance; here the question is narrower: does the code implement the specified model as intended?
Set the Review Scope Before Reading Expressions
A reviewer needs a reference against which to compare the code. That reference may include the model equations, state and parameter definitions, units, reaction or process descriptions, and documented assumptions. Names and units should agree with the shared model documentation; if the reference is incomplete or inconsistent, the appropriate finding may be a request for clarification rather than a claim that the code is wrong.
For a continuous-time model with states \(x_1,\ldots,x_m\) and parameters \(\theta_1,\ldots,\theta_q\), the specification commonly gives a rate equation for each state:
The implementation must assign each declared state to the corresponding computed quantity and implement each \(f_i\) with the intended variables, parameters, signs, and units. A reviewer should also establish the model boundary: which equations and implementation files are in scope, what behavior is specified, and which questions are outside this review. This prevents a code comment from silently becoming a scientific redesign.
Locate the documented equations, state definitions, parameter meanings, units, and assumptions. Record ambiguities before judging code against them.
For each equation, identify which implementation variable represents each mathematical state and which parameter each code symbol represents. Check that mappings are explicit and consistent.
Check the implemented operation, sign, factors, and dimensional units against the specified rate. Work a simple value through the expression if it helps expose a discrepancy.
State the location or expression, the expected specification, the observed difference, its possible effect, and an actionable request. Separate confirmed defects from questions requiring author input.
Trace Every Equation into the Implementation
A useful review proceeds one equation at a time. For each term, ask what process it represents, which state or parameter it uses, whether its sign matches the direction of change, and whether any multiplier or denominator is present in both representations. Then check that the full expression contributes to the correct state derivative. Matching variable names alone is not enough: a symbol can be named correctly but used in the wrong place.
Dimensional analysis supplies a compact justification for many checks. If a state \(C\) is a concentration measured in nM and time is measured in hours, then \(dC/dt\) has units nM/h. Every term added to that derivative must have units nM/h. Addition of quantities with incompatible dimensions is not meaningful, even if the programming language accepts the expression.
This check is necessary but not sufficient. A term can have the right units and still use the wrong biological state or parameter. Conversely, a codebase may use unitless internal values; in that case the reviewer still needs the declared conversion and mapping that explain how those values represent the specified quantities.
Worked Example: Detecting a Missing State Factor
Suppose the specification defines a hypothetical concentration state \(C\), in nM, with synthesis rate \(k_{\mathrm{syn}}=2.0\) nM/h and first-order degradation constant \(k_{\mathrm{deg}}=0.25\) h\(^{-1}\):
The implementation computes the derivative as “synthesis rate minus degradation constant,” without multiplying the latter by \(C\). At \(C=8.0\) nM, the specified derivative is
The implemented expression would instead subtract \(0.25\) h\(^{-1}\) directly from \(2.0\) nM/h. Those terms have different units, so the expression is dimensionally inconsistent; it also omits the state dependence required by the specification. A focused review comment would identify the derivative expression, quote the specified degradation term \(k_{\mathrm{deg}}C\), and request that the implementation include the state factor or explain a documented alternative. The discrepancy is established by comparison with the stated equation, not by whether a simulation appears reasonable.
Check State and Parameter Mapping
A mapping links the mathematical model’s symbols to the variables and parameters used by the implementation. Reviewers should check more than spelling. Confirm that each code variable represents the documented quantity, compartment, and form: for example, free concentration is not total concentration, and an amount is not a concentration. Check that parameters with similar names are not interchanged and that units and any conversions are consistent with the documented definitions.
For a reaction, a reliable check follows the species affected by the reaction and verifies that each state equation receives the appropriate contribution. For reversible binding of free ligand \(L\) and free receptor \(R\) to form complex \(LR\), mass action gives the forward rate \(k_{\mathrm{on}}LR\) and reverse rate \(k_{\mathrm{off}}LR\). The forward rate consumes one \(L\) and one \(R\), and forms one \(LR\); the reverse rate has the opposite effects. This sign and stoichiometry pattern is a useful reference when inspecting the three state equations.
Worked Example: Finding a Free-versus-Total Receptor Mismatch
Consider hypothetical states \(L\), free ligand; \(R\), free receptor; and \(LR\), bound complex, all in nM. The specified complex formation rate is \(k_{\mathrm{on}}LR\), where \(k_{\mathrm{on}}=0.10\) nM\(^{-1}\)h\(^{-1}\). Suppose the implementation uses \(R_{\mathrm{tot}}=R+LR\) in place of free \(R\). At \(L=2.0\) nM, \(R=3.0\) nM, and \(LR=1.0\) nM, the specified formation rate is
The implementation instead uses \(R_{\mathrm{tot}}=3.0+1.0=4.0\) nM and computes
The rate units are valid in both expressions, so a units-only review would miss the mapping error. The specification calls for free receptor, whereas the code uses total receptor; in this hypothetical state, the formation rate is 0.20 nM/h higher than specified. A useful comment names the expected free-receptor state, identifies the total-receptor expression, and asks the author to correct the mapping or document and justify a different model definition.
Review Units and Readability Together
Units are part of the model definition, not decoration. For every state and parameter, check the stated unit and time basis, then verify that each equation is dimensionally coherent. Pay particular attention to conversions between amount and concentration, between time units, and between quantities in different compartments. A conversion can be correct in isolation but wrong at the point where it is applied if the code’s state mapping is unclear.
Readability affects whether another scientist can perform these checks. Names should distinguish quantities with different meanings, and compound expressions should make the model structure visible. Avoid inferring biological meaning from a convenient code name when the documentation does not support it. If an expression is difficult to follow, a review can request a clearer name, a short explanation, or a decomposition into named intermediate quantities. These changes support inspection; they do not by themselves establish scientific correctness.
Worked Example: Reviewing a Time-Unit Conversion
Suppose a documented decay constant is \(k=0.020\) h\(^{-1}\), but the implementation advances a model whose time variable is in minutes and directly multiplies the state by \(k\) per minute. For first-order decay over one minute, the consistent rate constant is
For a hypothetical concentration of 10 nM, the initial decay rate should therefore be \(0.0003333\ \mathrm{min}^{-1}\times10\ \mathrm{nM}=0.003333\) nM/min in magnitude. Using \(0.020\) as though it were per minute gives \(0.200\) nM/min, 60 times the intended rate. A review comment should state the documented time unit, the implementation time unit, and the conversion required. It should also ask that the mapping be made explicit so the conversion is not applied twice elsewhere.
Writing Findings and Recognizing Review Limits
A strong review comment is specific and verifiable. It avoids vague judgments such as “this seems wrong” and does not prescribe a scientific change that the specification does not support. A useful structure is: location (where the expression appears), evidence (what the specification says and what the code does), impact (why the difference matters), and request (what correction or clarification is needed). When the evidence is incomplete, phrase the finding as a question and identify the missing definition.
| Review question | Evidence to inspect | Possible finding |
|---|---|---|
| Does each equation match its implementation? | Specified terms, signs, factors, and the computed state derivative | A term is missing, has the wrong sign, or contributes to the wrong state |
| Are states and parameters mapped correctly? | Documented definitions, code variables, compartments, and parameter units | Free and total quantities are confused, or a parameter is used for a different quantity |
| Are units coherent? | Units of every state, parameter, time basis, and conversion | Terms added together have incompatible units or a conversion is absent |
| Can another reviewer follow the implementation? | Names, structure, comments, and consistency with model documentation | Meaning is ambiguous and needs a clearer name or explicit explanation |
A code review is not a guarantee that the model is scientifically valid, that all possible defects have been found, or that every output is fit for a decision. It is also distinct from checking numerical solver choices or event handling, and from automated checks of behavior across code changes. Those questions require their own review or testing activities. For this tutorial’s scope, a finding should stay tied to implementation fidelity and interpretability.
A common mistake is to reason, “The code runs and produces plausible curves, so the implementation must match the model.” The conclusion does not follow: multiple coding errors can produce smooth, plausible output, and a numerically executable expression can still violate the specified units or state mapping. The correction is to trace each specified equation and mapping directly into the implementation, documenting the evidence for each finding.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions below.
- A state \(A\) is an amount in \(\mu\)mol and volume \(V\) is in L. What units does \(A/V\) have, and what mapping question should a reviewer ask if the equation specifies concentration?
- A specified rate is \(kC\), with \(k\) in h\(^{-1}\) and \(C\) in nM. What are the rate units, and why is subtracting \(k\) alone from a term in nM/h inconsistent?
- A code expression uses total receptor where the specification names free receptor. If both quantities have units of nM, can dimensional analysis alone detect the discrepancy?
- What elements should a reviewer include in a useful comment about a missing factor in a state derivative?
Solutions. (1) \(A/V\) has units \(\mu\)mol/L. The reviewer should check whether the implementation uses this concentration or the amount \(A\), and whether any volume conversion is specified. (2) \(kC\) has units \((\mathrm{h}^{-1})(\mathrm{nM})=\mathrm{nM/h}\). The constant \(k\) alone has units h\(^{-1}\), so it cannot be subtracted from a concentration rate. (3) No. Both expressions can be dimensionally coherent; the reviewer must check the state definitions and mapping. (4) Identify the code location, the specified equation and missing factor, the likely consequence or violated units, and a concrete correction or request for clarification.