Evaluate BNGL parameter expressions instead of dropping them

[wshlavacek]

A BNGL parameters block may give a parameter a value that is an expression over other parameters, such as kon koff/(KdNAV). That is ordinary BNGL rather than an edge case. Across 303 models drawn from the BioNetGen model collections, 1934 of 9323 parameter declarations are expression-valued.

Right now get_parameter_value raises NotImplementedError for those, and get_free_parameter_ids_with_values skips them without saying anything, so a PEtab problem built from a BNGL model quietly loses parameters. To reproduce on main:

from petab.v1.models.bngl_model import BnglModel, parse_bngl

text = """
begin parameters
  NA    6.022e23
  V     1e-12
  Kd    5.0
  koff  0.1
  kon   koff/(Kd*NA*V)
end parameters
"""
m = BnglModel(parse_bngl(text), model_id="demo")
print(list(m.get_parameter_ids()))
print(m.get_free_parameter_ids_with_values())

kon is listed as a parameter but is missing from the values, and nothing is said about it.

Resolving these needs no BNG2.pl and no reaction network, because a parameters block is only arithmetic over other parameters. This adds a parser and evaluator for that arithmetic and points BnglModel at it. Expressions are tokenized and parsed rather than passed to eval.

Getting the arithmetic right

BNGL arithmetic is not Python arithmetic, and the differences give a wrong number rather than an error. I checked each rule against BNG2.pl 2.9.3 by putting the expression in a parameters block and running writeNET({evaluate_expressions=>1}), which is the only export path that prints numbers instead of copying the source text back out. Three cases are easy to get wrong:

-2^2       is 4, because unary minus binds tighter than the power operator
2^3^2      is 64, because the power operator groups from the left
rint(2.5)  is 3, because rint is floor(x + 0.5)

The function table and the operator precedence follow Perl2/Expression.pm in the BioNetGen source. So _pi and _e take no arguments and are written _pi(), floor and ceil are rejected because they are commented out there as unsupported, and the comparison operators, the logical operators and if are supported because real models use them. A bare log is rejected, because BNGL's natural logarithm is ln and quietly accepting log would turn a typo into a plausible wrong number. TFUN is left out on purpose, because it reads a data file while a simulation runs and so is not a constant.

Partial resolution

get_free_parameter_ids_with_values now resolves what it can and warns about the rest by name. Failing the whole block on one bad definition would lose more than the current bug does. A real model can carry a construct we do not evaluate, and taking its other sixteen parameters down with it is a worse outcome than one silent drop.

Tests

The table of expressions and their expected values is checked from both sides. test_bng_verified_table pins the evaluator against it and needs no BioNetGen, so it runs in continuous integration. test_table_still_matches_bng2pl rebuilds the same values from a real BNG2.pl in a single run where one is available, and skips otherwise, so the table cannot go stale unnoticed if BioNetGen changes. I checked that the second test can actually fail by putting a wrong value in the table and watching it report the mismatch.

test_expression_valued_parameter_is_not_evaluated pinned the old behaviour, so I updated and renamed it.

The test suite passes apart from test_combine_archive and test_data_overview, which also fail on main here because of optional packages I do not have installed.

Credit

The evaluator started as a contribution from @arpitjain099 in lanl/PyBNF#673, which is where the parser design and many of the test cases come from. The BNG2.pl corrections and the partial resolution came out of review there.

[wshlavacek]

The Windows 3.14 job failed while downloading the test corpus rather than in
a test. It could not fetch Chattaraj_2021.bngl from the content delivery
network, reporting that the connection was forcibly closed by the remote
host. The other five jobs were cancelled by fail-fast rather than failing on
their own, so nothing here got as far as running a test.

This looks unrelated to the change. The same step failed on main on 26 August
at commit 847b16f, on a different file, and that same commit passed on the
days either side.

Locally the corpus fetches all 21 files with checksums verified, and
test_bngl_corpus.py passes.

Could someone re-run the job when you have a moment? I do not have the rights
to do it myself.

One thing worth looking at separately: the fetch step has no retry, so a
single transient network error fails the whole build.

[wshlavacek]

@dweindl this is ready for review whenever you have time. I do not have the
rights to assign a reviewer directly.

It will also need a CI re-run, for the corpus download reason in the comment
above.

[codecov-commenter]

Codecov Report

❌ Patch coverage is 88.74459% with 26 lines in your changes missing coverage. Please review.
✅ Project coverage is 76.48%. Comparing base (847b16f) to head (8f5b11a).
⚠️ Report is 6 commits behind head on main.

Files with missing lines	Patch %	Lines
petab/v1/models/bngl_model.py	88.74%	20 Missing and 6 partials ⚠️

Additional details and impacted files

@@            Coverage Diff             @@
##             main     #517      +/-   ##
==========================================
+ Coverage   76.01%   76.48%   +0.47%     
==========================================
  Files          67       67              
  Lines        7524     7744     +220     
  Branches     1341     1383      +42     
==========================================
+ Hits         5719     5923     +204     
- Misses       1302     1314      +12     
- Partials      503      507       +4     

☔ View full report in Codecov by Harness.
📢 Have feedback on the report? Share it here.

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[dweindl]

Thanks for the contribution. I already triggered a CI re-run. I should be able to review within a week.

[arpitjain099]

@dweindl sounds good.

[dweindl]

Thanks again for your work on PEtab+BNGL. I currently have two concerns regarding these changes related to maintainability and semantics:

• Maintainability / duplication:

  Together with Add support for BNGL models (cont.) #508 (Add support for BNGL models #501), this is adding a substantial amount of generic, PEtab-independent BNGL parsing/evaluation code to this repo. Would you be open to spinning this out into a small, separate Python library instead? Not a BNG2.pl replacement, just the parsing/evaluation functionality (block reading, expression arithmetic (ideally BNGL->sympy conversion), etc.) as a standalone, reusable package. libpetab-python would then depend on it rather than growing its own BNGL implementation.

  I think that would be much more maintainable long-term than having this logic duplicated between here and PyBNF (and potentially other consumers), and it would let the BNGL-parsing code be used and get tested independent of PEtab.

• Semantics, regarding the changes to BnglModel specifically:

  I think we have a problem there if get_parameter_value/get_free_parameter_ids_with_values only evaluates kon = koff/(Kd*NA*V) using the BNGL file's own defaults for koff, Kd, etc. In an PEtab problem, koff could be present in the PEtab parameter table (fixed override or estimated) and its effective value can differ from the file default. Then, kon needs to be recomputed from that value, not the file's. Model has no way to ask "evaluate this against these parameter-table values," only "give me the model's own value," so this can only be correct as long as nothing kon depends on is ever overridden or estimated, which isn't knowable at this point in the code.

  SBML has the same gap but fails safe: SbmlModel.get_free_parameter_ids_with_values only resolves an InitialAssignment when it's self-contained; anything referencing another parameter is left out entirely rather than evaluated against file defaults. I'd rather BnglModel do the same for now: eagerly evaluate only expressions with no parameter references, and leave anything referencing another parameter out of get_free_parameter_ids_with_values (and get_parameter_value).

  The proper fix is out of scope of this PR as it affects all Model classes: Model needs to expose the raw expression for a derived parameter, so it can be evaluated after the parameter table has been applied (see "Pre-initialization" at https://petab.readthedocs.io/en/latest/v2/documentation_data_format.html#v2-initialization-semantics). Tracked as part of Improve model abstractions #518.

[wshlavacek]

Thanks for the review. You were right about the semantics, and the change is in.

A parameter now gets a value only when its right hand side refers to no other
parameter, which covers a plain number and arithmetic over numbers and built-in
functions such as 2*_pi(). A parameter defined in terms of another parameter is
left out of get_free_parameter_ids_with_values, and get_parameter_value raises
for it.

I checked your reasoning against BNG2.pl, and
BioNetGen turns out to be a worse case than SBML rather than an exception to it.
Take kon = koff/Kd with koff 0.1 and Kd 5. If a consumer applies the value this
branch used to report, 0.02, and then applies a fitted koff of 10, BNG2.pl
writes "kon 0.02 # Constant". The expression is gone and the number is wrong by
a factor of a hundred, with nothing said. Leave kon alone and the same model
writes "kon 2 # ConstantExpression". In SBML the initial assignment would have
won and protected the user. In BioNetGen it does not.

Measured against released 0.9.0, which is what the version range resolves to
today, the change takes nothing away. Released 0.9.0 raises NotImplementedError
for every parameter whose right hand side is not a bare number, including
arithmetic over constants, and skips it in the free value list. Over 1047
BioNetGen models the branch reports 326 values that 0.9.0 does not, loses none,
and changes none. What the conservative rule gives up is relative only to the
first commit on this branch, which never shipped.

Three details worth knowing:

1. get_parameter_value now raises ValueError rather than NotImplementedError.
   create_parameter_df catches only ValueError, and NotImplementedError is not a
   subclass of it, so on 0.9.0 the error escapes that helper instead of leaving the
   nominal value empty. The change fixes that.

2. A parameter that refers to a name the parameters block does not declare is
   still reported in a warning, because BNG2.pl rejects such a model too. That
   keeps a broken model separate from a parameter that is only waiting for the
   parameter table, and waiting is silent as it is for SBML. The warning names at
   most five parameters and then counts the rest, because without a cap it ran to
   2865 characters on a real model that leaves a placeholder on 37 of its
   parameters.

3. In parameter_mapping, an empty value in the condition table for a derived
   parameter now reaches the "Not sure how to handle NaN in condition table" error,
   where before it was filled in from the model file. That is the trade you asked
   for.

I also corrected something I had written in the commit message. SbmlModel
matches this for get_free_parameter_ids_with_values, but its get_parameter_value
still returns the value written in the file, so BnglModel is now stricter than
SBML rather than the same.

Something else came up while measuring, which I think wants its own issue rather
than a change here. A derived parameter marked estimate=1 in the parameter table
never reaches the simulator, because _apply_parameter_table skips a row that is
not already in the mapping and the mapping is seeded from
get_free_parameter_ids_with_values. On a parameter table generated for one real
model that was 30 rows out of 46, all varied by the optimizer with no effect on
the objective, with no warning and no lint error. SBML does the same thing on main today for a parameter set by an initial
assignment, and released 0.9.0 dropped strictly more BioNetGen parameters than
this branch does. It looks like further evidence for issue #518.

A question about scope. evaluate_bngl_parameters and
evaluate_bngl_parameters_partial resolve a whole parameters block in dependency
order, and BnglModel no longer calls either of them. They are what a consumer
needs once the parameter table has been applied, which is the case in #518. If
they stay, they should go into all, because neither they nor the two
exception classes are exported today, and petab.v2.models.bngl_model re-exports
with a star import, so none of them can be reached from v2 at all. If you would
rather keep the surface small, I can cut the branch back to
evaluate_bngl_expression and bring the rest back with #518. Tell me which you
prefer.

On the duplication, I agree that one home for this code is ideal. I propose
that PEtab-dev owns the repository and the name on PyPI, and
we seed it and look after the BioNetGen side. The contents would be the reader
for the model file and the evaluator for parameter expressions, together with
the table of expressions and the values BNG2.pl computes for them, which is what
pins the arithmetic down. It needs only the standard library and it never calls
BNG2.pl. libpetab-python and PyBNF would both depend on it. On sympy, a
conversion to sympy is what #518 needs, and I would put it behind an optional
install so that reading a model file does not require sympy for everyone.

I would rather not hold this pull request while the package is sorted out. I can
open an issue for the move and do the extraction as separate work.