Fix: InChI-first shattered 461 metal compounds; formula column dropped every R group - #296
Merged
samseaver merged 6 commits intoSep 23, 2026
Conversation
Review of ModelSEED#295 found compounds losing their R groups. Confirmed, and the cause was this bundle's formula column, not Marvin and not the structures. WHAT WAS WRONG. The first cut wrote Marvin's getFormula() into the formula column and deferred Print_Structure_Formula_Charge.py as a follow-up step. Marvin omits wildcard atoms from a formula; this repository renders them as R, in one line of that script: formula = re.sub(r'\*', 'R', formula) So SMILE rows whose formula contains R fell from 23.4's 8,712 to 8. Stearoyl-ACPs went from C32H60N3O9PR2S to C32H60N3O9PS, and Update_Compound_Structures_Formulas_Charge.py carried that into 6,052 compound records -- cpd00049 "carboxylic acid" became CHO2, a generic compound quietly ceasing to be generic. THE STRUCTURES WERE NEVER AFFECTED: 8,727 SMILE structures carry a `*` in both bundles, identically. Only the column was wrong. That is why it survived the whole validation suite -- coverage, compound sets, InChI agreement, InChIKey consistency and stereochemical fidelity all inspect structures, and the structures were always right. A column-only defect was invisible to every check. A formula-column comparison against 23.4 is now part of the report. THE FIX. formula and charge come from Print_Structure_Formula_Charge's own parse_structure -- imported, not reimplemented -- computed per row from that row's structure string. Deferring it was the error; the cascade consumes this column directly, so the bundle has to be correct as written. Two further cases had to be reconciled to make that hold everywhere: parse_structure's substitution only fires when RDKit produced the formula, since RDKit renders a dummy atom as `*`. OpenBabel and Marvin both omit them, so structures like ISOCITHASE-P's `*OP(=O)(=O)=O` -- an invalid valence only OpenBabel will read -- still lost their R. The convention is now enforced as an invariant: a SMILE structure carrying `*` gets an R in its formula. This bundle violates it 0 times; 23.4 violated it 23 times. Counting wildcards to enforce that has its own trap: a molblock R atom, which is how RDKit writes every dummy atom and therefore how every structure arrives through the import bridge, reads back from Marvin as symbol "R#", not "R". A wildcard set without R# counts zero on a molecule that plainly has them. RESULT. R-bearing SMILE rows 8,735 against 23.4's 8,712; per source ChEBI 2,071=2,071, KEGG 964=964, Rhea 30=30, MetaCyc 5,670 vs 5,647. Formula-column agreement with 23.4 rises from 68,466 rows to 76,479. Coverage, structures, InChI, InChIKey and stereochemistry are unchanged -- 141,715 rows, still exactly 23.4's count, 0 non-standard InChI, 0 InChIKey disagreements. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Self-review of the numbers in the previous commit, not a data change --
the bundles are untouched and were already correct.
Counting R groups with `"R" in formula` also matches the elements Ru, Rb, Rh,
Re and Rn. That inflated every figure by 8 and, worse, made the regression look
milder than it was: the eight rows that appeared to survive with an R group
were ruthenium and rubidium compounds, not R groups at all.
SMILE rows with an R group was reported actually
this bundle 8,735 8,727
marvin 23.4 8,712 8,704
as merged in ModelSEED#295 8 0
So ModelSEED#295 did not leave 8 R groups standing; it left none. Per source the
corrected counts are ChEBI 2,071 = 2,071, KEGG 957 = 957, Rhea 30 = 30, and
MetaCyc 5,669 against 5,646 -- the +23 delta, and every other conclusion, is
unchanged.
Re-verified against the committed bundles while correcting this:
coverage 141,715 rows and 53,127 compounds, both exactly 23.4's
invariant 0 violations in either direction; 23.4 has 23
vs merged ModelSEED#295 9,185 formulas and 84 charges changed, 0 structures,
0 rows added or removed
columns one generated_on, tool_version and ph across all four files
reproducibility re-running a source rewrites it byte for byte
Both the report and the script docstring now state the rule (`R(?![a-z])`) next
to the numbers, so the next person to check does not repeat it.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
freiburgermsu
force-pushed
the
marvin-26.1-protonation-r-group-fix
branch
from
September 22, 2026 06:35
c129c91 to
ee55974
Compare
The manifest is what a maintainer reads before touching one of these bundles, and it described the tautomer gap and the RDKit import bridge but said nothing about the column that actually broke. Each 26.1 protonation entry now records that formula and charge come from Print_Structure_Formula_Charge.parse_structure, that substituting Marvin's getFormula() is what dropped the R group from 8,704 rows, and the invariant the file satisfies: a SMILE structure carrying `*` carries an R in its formula, 8,727 times, with 0 violations. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
Review of ModelSEED#296 flagged 74 formula rows regressing on metals. Chasing it found the formula was a symptom. The structures were broken. source Cl[Sn](c1ccccc1)(c1ccccc1)c1ccccc1 intact 23.4 Cl[Sn](c1ccccc1)(c1ccccc1)c1ccccc1 correct ModelSEED#295 [Cl-].[SnH3+].[c]1ccccc1.[c]1ccccc1.[c]1ccccc1 INCHI DISCONNECTS METAL-LIGAND BONDS BY DESIGN. Triphenyltin chloride's InChI is `3C6H5.ClH.Sn` -- five components before any tool touches it. INPUT SELECTION took InChI-first unconditionally, inherited from Run_Marvin_pKas.py where it is correct, so Marvin was handed an already-shattered molecule and protonated the pieces. 487 compounds have an InChI more fragmented than their SMILES -- cobalamins, Ni/Fe/Mg porphyrins, molybdenum cofactors, organotins, iron-sulfur clusters. 461 shipped fragmented. THE RULE is now InChI-first UNLESS the InChI is the more fragmented of the two. Comparing fragment counts rather than screening for metals keeps it general. 493 compounds take the SMILES route on that test. Two subtleties, both of which cost a regeneration to find: The InChI is DEMOTED, not reordered. Sharing one ladder lets a writable-but-wrong rung beat a correct one: on KEGG C18384 the SMILES yields the right dative-bonded magnesium propionate, which Marvin's SMILES writer refuses, while the disconnected InChI yields a writable three-fragment answer and wins. The InChI ROW is written from the InChI-derived molecule. Each column carries what its own representation can express, which is what 23.4 did: CPD-18407 ships a connected 8-iron cluster in SMILE and disconnected `InChI=1S/C.8Fe.6HS.3S/...` in InChI. This is also load bearing for stability. Asking Marvin to write an InChI for a CONNECTED metal cluster aborts the JVM with `free(): double free detected in tcache 2` from InChINativeGenerateInChICall -- a native fault no handler can catch, and one the previous regeneration died on. RESULT, measured against 23.4 rather than against the previous run: fragmented vs 23.4 461 -> 16 -> 8 R rows 8,727 vs 8,704, 0 invariant violations rows 141,715, exactly 23.4's ChEBI/KEGG/Rhea 0 compounds more fragmented than 23.4 The remaining 8 are not input selection. Marvin 26.1 breaks metal coordination bonds during protonation itself, from either representation: ferrocene splits into iron and two cyclopentadienyls, C12862 sheds both ammines. Confirmed by protonating the connected SMILES directly. Recorded, not worked around. WHY THIS SHIPPED. ModelSEED#296 claimed "0 structures changed" -- true only against ModelSEED#295, and both were wrong. Comparing a fix to a broken baseline instead of to 23.4 is the check that would have caught it, and is now the first table in the report. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
freiburgermsu
changed the title
Correction to the previous commit. It routed the InChI row of the 493
disconnected-InChI compounds to the source string unprotonated, on the premise
that 23.4 left them alone. That premise was checked on a handful of
iron-sulfur clusters and generalised without measuring it.
artifact compounds with a 23.4 InChI row 493
23.4 InChI == source (protonation a no-op) 234
23.4 InChI protonated, differs from source 259
23.4 protonated this form UNIFORMLY. The no-op was an outcome for 47% of them,
not a rule, and passing the source through silently un-protonates the other 259
-- chlorophylls and cobalamins lose a /p-2 layer they should carry.
So the InChI row is built from the InChI-derived molecule and protonated, which
is 23.4's method. The SMILE row still comes from the connected SMILES, which is
what fixed the shattering, and Marvin is still never asked to write an InChI
for a connected metal cluster, which is what aborted the JVM.
The cost is now stated rather than optimised away: 26.1 reads the detached
`4Fe.4S` as four FREE sulfide ions and protonates them to H2S at pH 7, so ChEBI
33722 ships Fe4S4 as H8Fe4S4 and 136511 ships MnO2 as H4MnO2. Both bundles are
internally consistent -- each formula matches its own structure -- and the
difference is which protonation state the engine assigns a ligand InChI
detached from its metal.
Measured against 23.4, unchanged by this correction:
fragmented vs 23.4 8 (461 in ModelSEED#295)
R rows 8,727 vs 8,704, 0 invariant violations
rows 141,715, exactly 23.4's
ChEBI/KEGG/Rhea 0 compounds more fragmented than 23.4
Against ModelSEED#295, 9,491 formula rows move: 8,388 restore agreement with 23.4, 1,028
agree with neither, 75 break it (54 metal-bearing). Those 75 are the same class
as the 20.3% net-charge delta -- 26.1 protonating what 23.4 did not -- and are
recorded in the report rather than suppressed. The pass-through scored 60
instead of 75 only by diverging from 23.4's method.
Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
.gitignore covers .DS_Store but not the AppleDouble ._* sidecar, so a `git add Biochemistry/Structures` picks it up. Not part of this change. Co-Authored-By: Claude Opus 5 (1M context) <noreply@anthropic.com>
freiburgermsu
added a commit
to freiburgermsu/ModelSEEDDatabase
that referenced
this pull request
Sep 25, 2026
…n-invariant tree The paper now describes the database as #7 (1ad34ad, on top of the merged ModelSEED#296) leaves it, so the numbers this branch introduced were re-measured there: M05 gold-graded transport that translocates a proton 2.1% -> 2.7% M13 reactions satisfying mass and charge balance 58% -> 69% balanced and carrying an annotated role 14,041 -> 15,128 of which outside every current template 8,089 -> 8,826 The balance share moves because #7 repairs the reaction refresh: Rebuild_Stoichiometry.py had never refreshed the charge embedded in each reagent, and Rebalance_Reactions.py rejected its own `save` flag, so the proton-adjustment step had not been persisting. Between dev (41b20c2) and 1ad34ad, 5,428 live reactions go from imbalanced to balanced: 3,130 by a proton or water adjustment the fixed step now writes, 969 through a participant whose record the re-pick changed, and 1,329 whose flag had been stale against unchanged records and stoichiometry; 25 go the other way, every one through a changed participant. That is the 58% -> 69%. The template union (8,597), the ATP-coupled share of gold transport (96%), dGPredictor's commit rates on and off transport (1.4% / 17.3%) and the role-annotation share (41%) did not move. Measured with Papers/NAR_Update_2026/analysis/review_transport_and_llm.py --live and the reaction records of that tree. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
freiburgermsu
added a commit
to freiburgermsu/ModelSEEDDatabase
that referenced
this pull request
Sep 25, 2026
…nt tree The paper now describes the database as #7 (1ad34ad, stacked on ModelSEED#296) leaves it. Its reaction refresh retires 19 more reactions, three of them among those no thermodynamic source directs, so the coverage figures move by the same small amount in M06 and S04: reactions no predictor directs 23,729 -> 23,726 of which the ensemble directs 22,902 -> 22,899 (96.5%, unchanged) eQuilibrator/ensemble co-commitments 8,085 -> 8,084 (94.7% agree, unchanged) reverse calls 5.1% -> 5.2% of the ensemble's calls Re-verified on that tree and unchanged: accuracy against measurement at tau (165 of 215, kappa 0.50; eQuilibrator 54 of 55, kappa 0.92), the sensitivity range from a bare sign test (68.1%, kappa 0.33) to RT ln 1000 (65.2%, kappa 0.38) with the maximum at tau (76.7%, kappa 0.50) and eQuilibrator never below 98.8%, the 85.3% forward share, and the one-sided errors (119 of 120 measured forward, 46 of 95 measured reverse). Measured with Papers/NAR_Update_2026/analysis/review_transport_and_llm.py --live. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
freiburgermsu
added a commit
to freiburgermsu/ModelSEEDDatabase
that referenced
this pull request
Sep 25, 2026
…the figures there The paper now describes the database as #7 (1ad34ad, stacked on ModelSEED#296) leaves it. That tree's reaction refresh retires 19 more reactions (48,403 -> 48,384 live) and rebalances the reactions against the re-picked compound records, and seven more compounds carry a structure (36,900 -> 36,907). Every count this branch touched was re-measured there, and the three figures were regenerated from the same tree with the analysis scripts as they stand on this branch: M09 MetaCyc 26,143/18,896 -> 26,140/18,893; KEGG 10,843/5,294 -> 10,840/5,291 structures 36,900 -> 36,907; reactions added since 2020 12,261 -> 12,259 M11 48,403 -> 48,384 reactions; Marvin 79.9% -> 79.8% of reactions; no ionizable site 1,238 -> 1,237; dGPredictor median 17.56 -> 17.57 M12 directed 24,674 -> 24,658; undirected 23,729 -> 23,726; dGPredictor-only 3,890 -> 3,877; graded 27,355 -> 27,338 (2,485 / 15,513 / 9,340); cascade 17,404 / 3,877 / 3,377 M13 atom mappings 26,256 -> 26,242 (19,809 clean, 6,433 salvaged) S03 six rows of the grade table and its total; unpaired 4,690 -> 4,686 Fig 1 caption 12,261 -> 12,259, twice The anchor set (365 distinct measured reactions), the growth rates (34% / 34%), the Rhea intake (8,409), the central-metabolism count (41) and the LLM accuracy at tau (165 of 215; 54 of 55) did not move. The response plan keeps its counts as measured at 41b20c2 and now says so at the top. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
freiburgermsu
added a commit
to freiburgermsu/ModelSEEDDatabase
that referenced
this pull request
Sep 25, 2026
Rebuilt on the protonation-invariant tree (#7, 1ad34ad, stacked on ModelSEED#296). The sixty reaction shards this writes are the same sixty that #7's reaction refresh rewrites, so the branch as first pushed (a4039a2 and 5da5c1e on dev at 41b20c2) conflicted with #7 in every one of them. The script is the one reviewed there, unchanged except for two counts in its docstring (below); its output was regenerated against #7's records. This commit replaces those two. WHAT THIS DOES. Not one of the 12,261 reactions added since the 2020 release (rxn48576-rxn60859) carried a pathway annotation: the alias file stops at rxn48568 and the `pathways` field with it, because the step that produced them, Scripts/Archived_Perl_Scripts/Compile_External_Pathways.pl, was never ported to Python. Build_Reaction_Pathways.py joins ModelSEED -> MetaCyc/KEGG alias -> the committed pathway tables, walks the parent closure so class-level annotation is present, and writes both the alias file and the `pathways` field, for reactions that carry no annotation only. On #7's tree: 2,835 reactions gain an annotation (2,285 post-2020, 550 earlier), the alias file gains 17,311 rows (121,444 -> 138,755), 24,569 of 56,012 records are annotated and 31,586 remain, most of them the Rhea intake. Multi-source records list MetaCyc before KEGG as every existing one does; 81 newly annotated records carry both. On dev the counts were 2,820 / 2,283 / 17,109 / 24,554 / 31,601 / 69; the difference tracks the 89 MetaCyc and KEGG reaction-alias rows that #7's alias refresh adds. WHY ADDITIVE. The committed pathway tables are a smaller snapshot than the distributions that produced the original file. The docstring had said a full rebuild reproduces 40,279 of the 121,444 rows and drops rxn00001's own annotation; the script's --check reports 104,203 reproduced and 17,241 dropped, on dev and on this tree alike, with 601 reactions losing every row and rxn00001 keeping all six of its own. The docstring now says that. The conclusion is unchanged: a regeneration would silently lose annotation that cannot be rebuilt from a clean checkout, so existing rows are never touched and --check keeps the remaining count visible. EC numbers are still not addressed: there is no committed EC table to join against, and Rhea's assignments need Data/rhea.rdf, which is not in the repository. VERIFIED on this tree: all 121,444 pre-existing alias rows preserved; exactly 2,835 records changed and `pathways` is the only field that differs on any of them; no existing value overwritten; every TSV keeps its line count, and the JSON shards differ only where a `null` became that record's list, so no shard is reformatted. Co-Authored-By: Claude Fable 5.1 <noreply@anthropic.com>
This was referenced Sep 25, 2026
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Follow-up to #295, which is already merged. Review of it surfaced two independent defects, both of which shipped. This PR fixes both.
1. InChI-first shattered 461 metal compounds
Review flagged 74 formula rows regressing on metals. The formula was a symptom; the structures were broken.
InChI disconnects metal–ligand bonds by design. Triphenyltin chloride's InChI is
3C6H5.ClH.Sn— five components before any tool touches it.INPUT SELECTIONtook InChI-first unconditionally (inherited fromRun_Marvin_pKas.py, where it is correct), so Marvin was handed an already-shattered molecule and protonated the pieces.487 compounds have an InChI more fragmented than their SMILES — cobalamins, Ni/Fe/Mg porphyrins, molybdenum cofactors, organotins, iron–sulfur clusters. 461 shipped fragmented.
The rule is now InChI-first unless the InChI is the more fragmented of the two. Comparing fragment counts rather than screening for metals keeps it general; 493 compounds take the SMILES route on that test.
Two subtleties, each of which cost a regeneration to find:
C18384the SMILES yields the right dative-bonded magnesium propionate, which Marvin's SMILES writer refuses, while the disconnected InChI yields a writable three-fragment answer and wins.CPD-18407ships a connected 8-iron cluster in SMILE and disconnectedInChI=1S/C.8Fe.6HS.3S/...in InChI. This is also load-bearing for stability: asking Marvin to write an InChI for a connected metal cluster aborts the JVM withfree(): double free detected in tcache 2fromInChINativeGenerateInChICall, a native fault no handler can catch, and one a previous regeneration died on.The remaining 8 are not input selection: Marvin 26.1 breaks metal coordination bonds during protonation itself, from either representation (ferrocene splits into iron plus two cyclopentadienyls;
C12862sheds both ammines). Confirmed by protonating the connected SMILES directly. Recorded in the report, not worked around.2. The formula column was dropping every R group
#295 wrote Marvin's
getFormula()into the formula column and deferredPrint_Structure_Formula_Charge.pyto a follow-up step. Marvin omits wildcard atoms; this repo renders them as R, in one line of that script:SMILE rows with an R group fell from 23.4's 8,704 to zero, and
Update_Compound_Structures_Formulas_Charge.pycarried that into 6,052 compound records —cpd00049"carboxylic acid" becameCHO2.formulaandchargenow come fromparse_structure— this repo's own function, imported rather than reimplemented, per row. Two further cases: the*→R substitution only fires on the RDKit path (OpenBabel omits dummy atoms), so the convention is enforced as an invariant — a SMILE structure carrying*gets an R in its formula; and a molblockRatom reads back from Marvin as symbolR#, so a wildcard set withoutR#silently counts zero.Result, measured against 23.4
ChEBI, KEGG and Rhea have zero compounds more fragmented than 23.4; all 8 residuals are in MetaCyc.
Why both of these shipped
#296 originally claimed "0 structures changed" — true only against #295, and both were wrong. Comparing a fix to a broken baseline rather than to 23.4 is the check that would have caught the metal defect, and the R-group defect survived because every check in the suite inspected structures, which on that axis were correct. Both checks — fragmentation vs 23.4, and the formula column vs 23.4 — are now in the report.
Generated by Claude Code