diff --git a/CHANGELOG.md b/CHANGELOG.md index 6afa8011..37e96553 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -252,6 +252,23 @@ All notable changes to PyBNF are documented below. This project adheres to by default. Both surfaces are documented under gradient-based fitting. ### Fixed +- **The PEtab export writes the fit's start point, so a round trip no longer moves the fit + back to a sampled draw (#719).** `nominalValue` is where a PEtab problem states the point + a fit starts from, and since #583 PyBNF reads it: an imported problem's nominal becomes a + `start_point` line. The export direction dropped it twice over. `write_parameter_table` + built each record as four fields and had no `nominalValue` column, so a row that carried + one wrote a file that could not reproduce it; and `_free_parameters_from_conf` built its + free parameters from the `_var` lines alone, never reading the `start_point` lines + beside them, so on the plain "publish my job as PEtab" path there was no value to drop in + the first place. Both are fixed. A PEtab problem now survives import -> export -> import + with its start point intact, and a native conf's `start_point` line exports as that + parameter's `nominalValue`. The column is written only when some parameter declares a + start, and a parameter that declares none writes an empty cell, so a job with no start + point exports exactly the four-column table it did before. A `start_point` outside the + parameter's own box is refused at export, since PEtab has no way to state one and the + re-import rejects exactly that; so is a `start_point` naming a parameter no exported + declaration claims, which was the other way the value could vanish without a word. This closes the gap between the export and the round-trip identity the PEtab + documentation states. - **Tutorial lesson 25 no longer promises rates its one curve cannot settle, and its check no longer passes on the luck of one seed (#703).** The lesson fits a transit-compartment model to a single plasma curve, and its README said island DE "recovers all three" rates. The diff --git a/docs/petab.rst b/docs/petab.rst index 07b561dd..6f116e75 100644 --- a/docs/petab.rst +++ b/docs/petab.rst @@ -101,7 +101,11 @@ following all survive an import and an export: :ref:`start_point ` line, so the imported fit starts from the problem's own published point instead of the box centre; delete the line to start from the centre. A ``nominalValue`` outside the row's own ``lowerBound``/``upperBound`` is a configuration - error rather than a silently relocated start. + error rather than a silently relocated start. The export writes the same fact back: a + ``start_point`` line becomes that parameter's ``nominalValue``, and a parameter with no + declared start writes an empty cell (the column is omitted entirely when the job + declares no start at all). An out-of-box ``start_point``, which PEtab has no way to + state, is refused at export rather than written as a bound. - **Observables and noise** — the ``observables`` table's noise half becomes a per-observable ``(noise model, noise-parameter source)``. Noise may be a fixed value, a data ``_SD`` column, or an estimated parameter, and it can vary by diff --git a/pybnf/petab/export.py b/pybnf/petab/export.py index 4310e89d..d100a9b0 100644 --- a/pybnf/petab/export.py +++ b/pybnf/petab/export.py @@ -73,7 +73,7 @@ from ..objective import _OBJECTIVE_DESUGAR from ..parse import ploop from ..printing import PybnfError -from ..pset import FreeParameter +from ..pset import FreeParameter, OutOfBoundsException from ._bngl import parse_model as parse_bngl_model from ._sbml import parse_model as parse_sbml_model from .conditions import ( @@ -1393,7 +1393,13 @@ def _read_conf_dict(conf_path): def _free_parameters_from_conf(conf): - """Build ``FreeParameter`` objects from the config's ``(keyword, name)`` entries.""" + """Build ``FreeParameter`` objects from the config's ``(keyword, name)`` entries. + + Each one carries the fit's declared start point for that parameter, if it has one, on + ``.value`` -- the source :func:`~pybnf.petab.parameters.petab_parameter_row` writes as + the row's ``nominalValue``, closing the export half of the #583 round trip (#719). + """ + start_points = _start_points_from_conf(conf) free_params = [] for key, value in conf.items(): if not (isinstance(key, tuple) and len(key) == 2 @@ -1414,14 +1420,52 @@ def _free_parameters_from_conf(conf): # is the native ``bounded`` flag, inert for the location families. A one-parameter # unbounded family (exponential/chisquare/rayleigh, #417) carries only p1. p2 = float(value[1]) if len(value) >= 2 else None - free_params.append(FreeParameter(name, keyword, float(value[0]), p2)) + start = start_points.pop(name, None) + try: + free_params.append( + FreeParameter(name, keyword, float(value[0]), p2, value=start)) + except OutOfBoundsException: + # A start point outside the parameter's own box. PEtab cannot express one + # either (a nominalValue outside lowerBound/upperBound is what the importer + # refuses on the way back in), and the bare OutOfBoundsException subclasses + # Exception, so pybnf.main would report it as "an unknown error ... please + # report this bug" on a config the user wrote (#583). + raise PybnfError( + f"start point out of bounds for '{name}'", + f"The config starts '{name}' at {start}, which is outside the box " + f"[{value[0]}, {p2}] its own declaration gives it. A start point is " + f"refused rather than moved.", + "Correct the start point, or widen the parameter's bounds.") if not free_params: raise PybnfError( "No exportable free parameters found in the config (expected one of " f"{sorted(EXPORTABLE_PRIOR_KEYWORDS)}).") + if start_points: + # A start point for a name no exportable free parameter claims. Silently dropping + # it is the failure this whole path exists to remove, and config.py refuses the + # same thing when the job is run. + names = ', '.join(sorted(start_points)) + known = ', '.join(sorted(fp.name for fp in free_params)) + raise PybnfError( + f"start point for unknown parameter(s) {names}", + f"The config declares a start point for {names}, which no exportable free " + f"parameter declaration names. The free parameters being exported are: " + f"{known}.") return free_params +def _start_points_from_conf(conf): + """``{name: theta}`` for every ``start_point = `` line (#583). + + The ``start_point`` half of the two spellings ``Configuration._load_start_point`` + merges. The other, a ``parameter:`` record's ``initial_value:`` field, is not read + here because the exporter does not build free parameters from ``parameter:`` records + at all -- a separate gap (#733), not a second place to populate. + """ + return {key[1]: float(value) for key, value in conf.items() + if isinstance(key, tuple) and len(key) == 2 and key[0] == 'start_point'} + + @dataclass(frozen=True) class _SbmlModelView: """An SBML model presented through the attribute surface the exporter's shared diff --git a/pybnf/petab/parameters.py b/pybnf/petab/parameters.py index 08e5ec3a..de15e264 100644 --- a/pybnf/petab/parameters.py +++ b/pybnf/petab/parameters.py @@ -450,25 +450,39 @@ def _petab_prior_row(fp, parameter_id, dist, stem, is_log, nominal): _PARAMETER_COLUMNS = ['parameterId', 'estimate', 'lowerBound', 'upperBound'] +_NOMINAL_COLUMN = 'nominalValue' _PRIOR_COLUMNS = ['priorDistribution', 'priorParameters'] def write_parameter_table(rows, path): """Write parameter ``rows`` to ``path`` as a PEtab v2 ``parameters.tsv``. - Always writes ``parameterId``/``estimate``/``lowerBound``/``upperBound``. The - prior columns (``priorDistribution``/``priorParameters``) are appended only when - some row carries an explicit prior, so a plain ``uniform_var`` job keeps the - four-column chunk-1 shape (PEtab v2 defaults a prior-less estimated parameter to - uniform-over-bounds). An unbounded location-scale family writes blank bounds; - ``nominalValue`` is optional in PEtab v2 and omitted while unused. + Always writes ``parameterId``/``estimate``/``lowerBound``/``upperBound``. The two + optional chunks are appended only when some row needs them, so a plain + ``uniform_var`` job with no declared start keeps the four-column chunk-1 shape: + + * ``nominalValue`` -- the fit's start point, the reverse of the importer reading it + onto ``FreeParameter.value`` and emitting a ``start_point`` line (#583). It was + omitted while nothing produced it; since #583 it is the point the job starts from, + and dropping it silently moved a re-imported fit back to a sampled draw (#719). A + row without one writes a blank cell (PEtab v2 leaves ``nominalValue`` optional for + an estimated parameter). + * ``priorDistribution``/``priorParameters`` -- an explicit prior (PEtab v2 defaults a + prior-less estimated parameter to uniform-over-bounds). + + An unbounded location-scale family writes blank bounds. """ + has_nominal = any(r.nominal_value is not None for r in rows) has_prior = any(r.prior_distribution is not None for r in rows) - header = _PARAMETER_COLUMNS + (_PRIOR_COLUMNS if has_prior else []) + header = (_PARAMETER_COLUMNS + + ([_NOMINAL_COLUMN] if has_nominal else []) + + (_PRIOR_COLUMNS if has_prior else [])) records = [] for r in rows: rec = [r.parameter_id, 'true' if r.estimate else 'false', num(r.lower_bound), num(r.upper_bound)] + if has_nominal: + rec.append(num(r.nominal_value)) if has_prior: rec += [r.prior_distribution or '', ';'.join(num(p) for p in r.prior_parameters)] diff --git a/tests/test_petab_export.py b/tests/test_petab_export.py index a84624a9..da3ab2fd 100644 --- a/tests/test_petab_export.py +++ b/tests/test_petab_export.py @@ -48,7 +48,10 @@ PetabParameterRow, free_parameter_from_row, petab_parameter_row, + read_parameter_table, + write_parameter_table, ) +from pybnf.printing import PybnfError from pybnf.pset import FreeParameter DEMO_DIR = Path(__file__).resolve().parents[1] / 'examples' / 'demo' @@ -762,6 +765,150 @@ def test_prediction_dependent_noise_exports_verbatim_formula(self, tmp_path): (imp2 / 'imported.conf').read_text() +# --------------------------------------------------------------------------- +# The fit's declared start point survives the export (#719). #583 made a PEtab +# nominalValue the imported fit's start_point; the reverse direction dropped it twice -- +# the conf reader never looked at the start_point lines, and the writer had no +# nominalValue column -- so a round trip silently moved the fit back to a sampled draw. +# --------------------------------------------------------------------------- + +def _demo_job(tmp_path, extra=''): + """A copy of the demo job whose conf has ``extra`` appended (a start_point line).""" + import shutil + job = tmp_path / 'job' + shutil.copytree(DEMO_DIR, job) + conf = job / DEMO_CONF.name + conf.write_text(conf.read_text() + extra) + return conf + + +class TestStartPointRoundTrip: + + # -- the writer, in isolation -------------------------------------------------- + + def test_written_nominal_value_reads_back(self, tmp_path): + # The unit the loss lived in: a row carrying a nominal must survive the TSV. + # write_parameter_table built a four-field record and read_parameter_table found + # no nominalValue column, so the value died between them (#719). + fp = FreeParameter('v1', 'uniform_var', 0.0, 10.0, value=0.5) + row = petab_parameter_row(fp) + path = tmp_path / 'parameters.tsv' + write_parameter_table([row], path) + assert 'nominalValue' in path.read_text().splitlines()[0].split('\t') + back = read_parameter_table(path)[0] + assert back == row # the whole row, not just the nominal + assert free_parameter_from_row(back).value == 0.5 + + def test_no_nominal_keeps_the_four_column_shape(self, tmp_path): + # The column is optional in PEtab v2 and stays absent when nothing declares a + # start, so a job with no start point exports byte-for-byte as it did before. + row = petab_parameter_row(FreeParameter('v1', 'uniform_var', 0.0, 10.0)) + path = tmp_path / 'parameters.tsv' + write_parameter_table([row], path) + assert path.read_text() == 'parameterId\testimate\tlowerBound\tupperBound\nv1\ttrue\t0\t10\n' + + def test_partial_nominals_write_blank_cells(self, tmp_path): + # A start point is partial by design (ADR-0117): declaring one parameter's start + # must not invent a nominal for the others. The column appears once any row needs + # it, and the rows that do not carry one write an empty cell, which reads back None. + rows = [PetabParameterRow('a', True, 0., 10., nominal_value=None), + PetabParameterRow('b', True, 0., 10., nominal_value=4.), + PetabParameterRow('c', True, 0., 10., nominal_value=None)] + path = tmp_path / 'parameters.tsv' + write_parameter_table(rows, path) + assert [r.nominal_value for r in read_parameter_table(path)] == [None, 4.0, None] + + def test_nominal_and_prior_columns_coexist(self, tmp_path): + # Both optional chunks at once: nominalValue precedes the prior pair (PEtab v2's + # own column order), and neither displaces the other's cells. + fp = FreeParameter('k', 'loguniform_var', 1e-3, 1e3, value=1.0) + path = tmp_path / 'parameters.tsv' + write_parameter_table([petab_parameter_row(fp)], path) + header = path.read_text().splitlines()[0].split('\t') + assert header == ['parameterId', 'estimate', 'lowerBound', 'upperBound', + 'nominalValue', 'priorDistribution', 'priorParameters'] + (back,) = read_parameter_table(path) + assert (back.nominal_value, back.prior_distribution) == (1.0, 'log-uniform') + assert back.prior_parameters == pytest.approx((1e-3, 1e3)) + + # -- the conf reader ------------------------------------------------------------ + + def test_conf_start_point_becomes_a_nominal_value(self, tmp_path): + # The plain "publish my job as PEtab" path: a native conf's start_point line is the + # fit's declared start, and it exports as that parameter's nominalValue. The + # parameters it says nothing about keep an empty cell. + export_job(_demo_job(tmp_path, '\nstart_point = v1 0.7\n'), tmp_path / 'out') + nominal = {r['parameterId']: r['nominalValue'] + for r in _tsv_rows(tmp_path / 'out' / 'parameters.tsv')} + assert nominal == {'v1': '0.7', 'v2': '', 'v3': ''} + + def test_conf_without_a_start_point_writes_no_nominal_column(self, tmp_path): + # The unchanged path: the demo conf declares no start, so the exported table keeps + # its four-column shape rather than gaining a column of blanks. + export_job(DEMO_CONF, tmp_path / 'out') + assert 'nominalValue' not in (tmp_path / 'out' / 'parameters.tsv').read_text() + + def test_out_of_box_start_point_is_refused(self, tmp_path): + # PEtab cannot express a nominalValue outside the row's own bounds either -- the + # importer refuses exactly that on the way back in. Refused as a PybnfError, not the + # bare OutOfBoundsException, which pybnf.main reports as "an unknown error" (#583). + with pytest.raises(PybnfError, match='out of bounds'): + export_job(_demo_job(tmp_path, '\nstart_point = v1 42\n'), tmp_path / 'out') + + def test_start_point_for_an_undeclared_parameter_is_refused(self, tmp_path): + # Silently dropping it is the failure this path exists to remove; config.py refuses + # the same line when the job is run. + with pytest.raises(PybnfError, match='unknown parameter'): + export_job(_demo_job(tmp_path, '\nstart_point = nope 1\n'), tmp_path / 'out') + + def test_partial_nominal_problem_is_petab_valid(self, tmp_path): + # The external oracle on the shape the fixtures never produce: an estimated row with + # a blank nominalValue beside rows that carry one. + pytest.importorskip('petab.v2') + export_job(_demo_job(tmp_path, '\nstart_point = v1 0.7\n'), tmp_path / 'out') + assert _petab_validation_errors(tmp_path / 'out' / 'problem.yaml') == [] + + # -- the full round trip, the harm as reported ---------------------------------- + + def test_petab_nominal_survives_import_export_reimport(self, tmp_path): + # fixedsigma_v2 ships nominalValue 0.5/1/3. The first import turns them into + # start_point lines (#583); before #719 the export dropped them and the second + # import produced a conf with none, so the fit began from a sampled draw instead of + # the point the published problem states -- with no warning and no error. + from pybnf.petab.import_ import import_job + imp1, pet2, imp2 = tmp_path / 'imp1', tmp_path / 'pet2', tmp_path / 'imp2' + import_job(FIXEDSIGMA_DIR / 'problem.yaml', imp1) + export_job(imp1 / 'imported.conf', pet2) + import_job(pet2 / 'problem.yaml', imp2) + + nominal = {r['parameterId']: r['nominalValue'] + for r in _tsv_rows(pet2 / 'parameters.tsv')} + assert nominal == {'v1': '0.5', 'v2': '1', 'v3': '3'} + + def starts(conf): + return sorted(line.strip() for line in conf.read_text().splitlines() + if line.startswith('start_point')) + assert starts(imp2 / 'imported.conf') == starts(imp1 / 'imported.conf') + assert starts(imp2 / 'imported.conf') == [ + 'start_point = v1 0.5', 'start_point = v2 1', 'start_point = v3 3'] + + def test_reimported_conf_starts_the_fit_at_the_published_point(self, tmp_path, monkeypatch): + # The user-visible fact behind the tsv cells: the twice-round-tripped conf loads to a + # Configuration whose resolved start point is the published one. Configuration is the + # consumer every start-point optimizer reads (ADR-0117), so this is the end of the chain. + from pybnf import config as config_mod + from pybnf.parse import ploop + from pybnf.petab.import_ import import_job + imp1, pet2, imp2 = tmp_path / 'imp1', tmp_path / 'pet2', tmp_path / 'imp2' + import_job(FIXEDSIGMA_DIR / 'problem.yaml', imp1) + export_job(imp1 / 'imported.conf', pet2) + import_job(pet2 / 'problem.yaml', imp2) + monkeypatch.chdir(imp2) + cfg = config_mod.Configuration( + ploop((imp2 / 'imported.conf').read_text().splitlines(keepends=True))) + assert cfg.start_point == {'v1': 0.5, 'v2': 1.0, 'v3': 3.0} + + # --------------------------------------------------------------------------- # Chunk 5a: the exporter reads the new-era data surface (model: / experiment: / data: / # observable:) as transcription, not the legacy linkage (ADR-0028). The demo twin