Core

Data model

Core data types: ligands, atom mappings, transformations, and networks.

Every type here is a frozen dataclass that validates its own invariants on construction. That choice is deliberate and load-bearing: the soft-core repair, the scorers, and the planner all assume a mapping is well formed, and the cheapest place to guarantee that is at the boundary. An AtomMapping that exists is a valid one.

The central type is AtomMapping, which promotes the {"sc1", "sc2", "cc1", "cc2"} dictionary contract used by amberstudio’s BuildEdges into a real type with enforced invariants. AtomMapping.from_contract() and AtomMapping.to_contract() round-trip that dictionary exactly, which is the seam that lets BuildEdges(mapping_method=...) call this package through a small shim.

rbfenetmap.core.models.EDGE_SEPARATOR = '~'

Separator between the two endpoints of an edge in file names and CLI arguments. Ligand names are forbidden from containing it so parse_edge_key is unambiguous.

class rbfenetmap.core.models.AtomMapping(cc1, cc2, sc1, sc2, n_atoms_1, n_atoms_2, method='unknown')

Bases: object

A common-core / soft-core partition of two molecules.

cc1[i] and cc2[i] are corresponding atoms – the pairing is positional, not by sorted order of cc2. sc1 and sc2 are the unmapped (soft-core) atoms on each side.

Parameters:

• cc1 (tuple[int, ...]) – The common core on each side, paired positionally. cc1 is sorted ascending.

• cc2 (tuple[int, ...]) – The common core on each side, paired positionally. cc1 is sorted ascending.

• sc1 (tuple[int, ...]) – The soft-core atoms on each side, sorted ascending.

• sc2 (tuple[int, ...]) – The soft-core atoms on each side, sorted ascending.

• n_atoms_1 (int) – Atom counts of the two molecules, so the partition can be checked for completeness without holding a reference to the molecules.

• n_atoms_2 (int) – Atom counts of the two molecules, so the partition can be checked for completeness without holding a reference to the molecules.

• method (str) – Name of the plugin that produced the mapping.

Raises:

ValueError – If any invariant below is violated.

Notes

The enforced invariants are:

• no duplicate indices within any of the four tuples;

• every index lies in range(n_atoms_k);

• sc_k and cc_k are disjoint and together cover range(n_atoms_k) (every atom is either transformed or held in common – there is no third state);

• len(cc1) == len(cc2);

• cc2 contains no duplicates, i.e. the correspondence is injective.

The last two are worth stating explicitly because together they are exactly Amber’s linear-scaling constraint, len(TI1) - len(SC1) == len(TI2) - len(SC2). Since len(TI_k) - len(SC_k) is just len(cc_k), an AtomMapping that exists already satisfies it. The Amber exporter re-checks anyway, but only hand-authored maps can ever trip it.

classmethod from_contract(contract, *, n_atoms_1, n_atoms_2, method='unknown')

Build from the {"sc1", "sc2", "cc1", "cc2"} dictionary contract.

This is the amberstudio BuildEdges interchange format. cc1/cc2 are taken as already paired positionally and are not re-sorted independently – doing so would silently scramble the correspondence.

Parameters:

• contract (Mapping[str, Sequence[int]]) – Must contain the keys sc1, sc2, cc1, cc2.

• n_atoms_1 (int) – Atom counts of the two molecules.

• n_atoms_2 (int) – Atom counts of the two molecules.

• method (str, optional) – Name recorded as the producing method.

Return type:

AtomMapping

classmethod from_core_pairs(core, *, n_atoms_1, n_atoms_2, method='unknown')

Build from a {idx1: idx2} correspondence, inferring the soft-core.

Every atom not appearing in core becomes soft-core on its side. This is the usual entry point for a mapper, which naturally produces a correspondence rather than a four-way partition.

Parameters:

• core (Mapping[int, int] or Sequence[tuple[int, int]]) – The common-core correspondence from molecule 1 to molecule 2.

• n_atoms_1 (int) – Atom counts of the two molecules.

• n_atoms_2 (int) – Atom counts of the two molecules.

• method (str, optional) – Name recorded as the producing method.

Return type:

AtomMapping

to_contract()

Return the {"sc1", "sc2", "cc1", "cc2"} dictionary contract.

Byte-for-byte what amberstudio’s cartograph_mapping_method returns, so a shim on the amberstudio side can hand this straight to BuildEdges. That shim also absorbs BuildEdges’ two unused parmed.Structure positionals – which is why ParmEd is not a dependency of this package.

Return type:

dict[str, tuple[int, …]]

property forward: dict[int, int]

The common-core correspondence as {idx1: idx2}.

property reverse: dict[int, int]

The common-core correspondence as {idx2: idx1}.

property n_common_core: int

Number of mapped atom pairs.

property n_softcore_1: int

Soft-core atom count on side 1, hydrogens included.

property n_softcore_2: int

Soft-core atom count on side 2, hydrogens included.

swapped()

Return the mapping with the two sides exchanged.

Used when a transformation is reoriented. The positional pairing is preserved by re-sorting on the new side 1.

Return type:

AtomMapping

class rbfenetmap.core.models.EdgeKind(*values)

Bases: str, Enum

Which alchemical experiment an edge stands for.

RBFE is the relative calculation this package was built to plan: one common core held fixed while two soft-core regions are grown and shrunk. CBFE is a counterpoised binding free energy – two absolute calculations run simultaneously in opposite directions, one ligand decoupling as the other couples.

A CBFE edge needs no correspondence between the two ligands, because neither molecule is being morphed into the other. That is what makes it useful to a network planner: it is available between any two ligands, including the pairs an MCS search cannot relate at all, so it can join subnetworks that RBFE alone leaves disconnected.

A str enum for the same reason as RejectionReason – it serializes to JSON as itself and compares equal to a plain string.

class rbfenetmap.core.models.EdgeScore(total=inf, feasible=False, descriptors=<factory>, contributions=<factory>, rejections=(), scorer='unknown')

Bases: object

The cost of a candidate transformation. Lower is better.

Parameters:

• total (float) – The scalar cost, math.inf when the edge is infeasible.

• feasible (bool) – Whether the edge may be selected at all.

• descriptors (Mapping[str, float]) – Raw descriptor values, as produced by rbfenetmap.core.descriptors.compute_descriptors().

• contributions (Mapping[str, float]) – Weighted terms, summing to total for feasible edges.

• rejections (tuple[RejectionReason, ...]) – Why the edge is infeasible; empty when it is not.

• scorer (str) – Name of the scoring plugin.

Notes

Feasibility and cost are kept strictly separate. Rejection is structural and originates only in the mapper, the repair, or validation – never from a weighted sum crossing a threshold. A merely bad edge has a large finite total and stays in the candidate pool where the planner can still use it if the alternative is a disconnected network. This is why infeasible candidates are retained on Network.candidates rather than dropped: they are the audit trail that explains a disconnection.

classmethod rejected(*reasons, scorer='unknown', **kwargs)

Build an infeasible score carrying reasons.

Parameters:

• reasons (RejectionReason)

• scorer (str)

• kwargs (Any)

Return type:

EdgeScore

class rbfenetmap.core.models.IntermediateRecord(source, target, generator, accepted=False, names=(), rejection=None, trace=())

Bases: object

One attempt to invent a ligand bridging a pair, successful or not.

Parameters:

• source (str) – The gap the generator was asked to bridge.

• target (str) – The gap the generator was asked to bridge.

• generator (str) – Registered name of the generator plugin.

• accepted (bool, optional) – Whether any molecule survived posing and was added to the ligand set.

• names (tuple[str, ...], optional) – Names of the ligands the attempt contributed.

• rejection (str, optional) – Why nothing was contributed. A plain string, not a RejectionReason: that enum is the vocabulary of edge feasibility, and overloading it would make core_geometry_mismatch mean two different things depending on where it was read from.

• trace (tuple[str, ...], optional) – Human-readable log of what the generator and the poser did.

Notes

Retained for the same reason rejected candidates are retained on Network. Without a record per gap attempted, a network where generation ran and found nothing is indistinguishable from one where generation was never enabled – and those two call for opposite responses from the user.

class rbfenetmap.core.models.Ligand(name, mol, charge, source=None, metadata=<factory>, provenance=None)

Bases: object

A network vertex: one molecule with a single 3D conformer.

Parameters:

• name (str) – The vertex identifier, also the file-name token used by exporters. Restricted to [A-Za-z0-9_.+-]+ so it is filesystem-safe and cannot contain the edge separator ~.

• mol (rdkit.Chem.Mol) – The molecule. Must carry explicit hydrogens and exactly one 3D conformer.

• charge (int) – Net formal charge, cached at construction so scorers never need to import RDKit. Use from_mol() to compute it.

• source (pathlib.Path, optional) – Where the molecule was read from, for diagnostics.

• metadata (Mapping[str, Any], optional) – Free-form annotations carried through to exported networks.

• provenance (LigandProvenance, optional) – Set only on a ligand this package invented; None for every molecule read from an input file. Appended last and defaulted so every existing positional construction – including the one in rbfenetmap.io.networkio – keeps working untouched.

Raises:

ValueError – If the name is malformed, the molecule is empty, it does not have exactly one 3D conformer, or any atom carries implicit hydrogens.

Notes

The implicit-hydrogen check is not pedantry. Mappers address atoms positionally by index, and every invariant in AtomMapping is stated over range(mol.GetNumAtoms()). An implicit hydrogen is an atom that participates in the chemistry but has no index, so a mapping over a molecule with implicit Hs is silently incomplete: the soft-core region it describes omits atoms that a downstream engine will nonetheless have to transform. Requiring Chem.AddHs up front makes that impossible rather than merely unlikely.

A synthetic ligand is this class with provenance set, not a subclass. The network loader constructs Ligand directly, so a subclass would silently downgrade to a plain ligand on every round-trip and every isinstance check downstream would become round-trip-fragile. A field survives that, and survives the dataclasses.replace() calls in rbfenetmap.core.align as well.

classmethod from_mol(mol, name, *, source=None, **metadata)

Build a Ligand, computing the net charge from mol.

Parameters:

• mol (rdkit.Chem.Mol) – Molecule with explicit hydrogens and one 3D conformer.

• name (str) – Vertex identifier.

• source (pathlib.Path, optional) – Origin of the molecule.

• **metadata – Stored on the ligand and carried into exports.

Return type:

Ligand

classmethod synthesized(mol, name, provenance, *, source=None, metadata=None)

Build a ligand this package invented, recording how.

Parameters:

• mol (rdkit.Chem.Mol) – Molecule with explicit hydrogens and the single posed 3D conformer.

• name (str) – Vertex identifier, normally from intermediate_name().

• provenance (LigandProvenance) – Where the molecule came from.

• source (pathlib.Path, optional) – Origin, if the molecule was also written somewhere.

• metadata (Mapping[str, Any], optional) – Free-form annotations.

Returns:

With synthetic true.

Return type:

Ligand

Notes

A separate constructor rather than a provenance= keyword on from_mol(), because from_mol collects **metadata: a user annotating a ligand with from_mol(mol, name, provenance="ChEMBL") would silently be setting this field instead of storing their note. Splitting the two constructors makes the collision impossible rather than merely documented.

property synthetic: bool

Whether this ligand was invented by the planner rather than supplied.

property n_atoms: int

Total atom count, including hydrogens.

property heavy_indices: tuple[int, ...]

Indices of the non-hydrogen atoms.

property n_heavy: int

Number of non-hydrogen atoms.

property atom_names: tuple[str, ...]

Per-atom names, used by the Amber exporter to build masks.

Prefers an existing name from the Tripos mol2 property or PDB residue info; falls back to element + 1-based index (C1, N2, …), which is unique by construction.

class rbfenetmap.core.models.LigandProvenance(kind, generator, parents, pose_method, pose_rmsd, detail=<factory>)

Bases: object

Where a ligand came from, when it was not read from a file.

Parameters:

• kind (str) – What sort of construction produced the ligand. "intermediate" is the only value this package writes today; the field is a string rather than an enum so a third-party generator can record its own kind without patching this module.

• generator (str) – Registered name of the plugin that proposed the molecule.

• parents (tuple[str, ...]) – Names of the real ligands it was derived from, sorted. These are the endpoints of the gap the intermediate was invented to bridge.

• pose_method (str) – How the conformer was produced – "parent_atom_map" when the generator handed over a complete correspondence, "mcs_fallback" when one had to be recovered. The weaker method is named rather than hidden precisely so it is visible in a report.

• pose_rmsd (float) – In-place RMSD, in angstroms, of the posed atoms against the parent coordinates they were taken from. Measured with the same core_rmsd() the feasibility gate uses, so it is directly comparable to softcore.core_rmsd_threshold.

• detail (Mapping[str, Any], optional) – Free-form annotations from the generator and the poser.

Raises:

ValueError – If kind or generator is empty, parents is empty, or pose_rmsd is negative.

Notes

A bare synthetic: bool cannot answer the first question anyone asks of an invented molecule – from what, by what, and how good is the pose. Recording the answer at construction is also the only place it is knowable: by the time the network reaches an exporter the generator is long gone.

class rbfenetmap.core.models.Network(ligands, edges=(), candidates=(), planner='unknown', options=None, unmet_constraints=(), intermediates=())

Bases: object

A planned perturbation network.

Parameters:

• ligands (Mapping[str, Ligand]) – The vertices, in insertion order.

• edges (tuple[Transformation, ...]) – The selected transformations.

• candidates (tuple[Transformation, ...]) – Every transformation that was scored, feasible or not. Retained as an audit trail: when the planner reports a disconnection, this is what explains it.

• planner (str) – Name of the planning plugin.

• options (NetworkOptions, optional) – The options the network was planned under.

• unmet_constraints (tuple[str, ...]) – Best-effort constraints that could not be satisfied (for example a requested edges_per_ligand the candidate pool could not support). Hard conflicts raise NetworkPlanError instead of landing here.

• intermediates (tuple[IntermediateRecord, ...]) – One record per gap intermediate generation was attempted on, whether or not it produced anything. Appended last and defaulted, so every existing construction site is unaffected and a network planned with generation off carries an empty tuple that serializes to nothing at all.

to_networkx()

Return the selected edges as an undirected networkx.Graph.

Nodes carry synthetic; each edge carries transformation, weight (the score total), and kind (the EdgeKind value as a plain string) attributes. kind is duplicated out of the transformation so consumers that only style or filter edges – the SVG renderer, the GraphML exporter – never have to reach back through the object.

synthetic is a plain bool rather than the whole LigandProvenance for the same reason and one more: GraphML types every attribute, so a nested mapping would not survive the trip at all, and a graph exported without it would show an invented vertex as an ordinary ligand – which is precisely the mistake this feature can make that costs somebody a simulation.

Return type:

nx.Graph

validate(*, require_connected=True)

Check the network is structurally sound.

Parameters:

require_connected (bool, optional) – Whether to require that the selected edges span every ligand.

Raises:

ValueError – If an endpoint is unknown, a self-loop is present, an unordered pair appears twice, or (when required) the network is disconnected.

Return type:

None

property rejected: tuple[Transformation, ...]

Candidates that were found infeasible.

property synthetic_ligands: tuple[Ligand, ...]

Vertices this package invented rather than read from an input file.

Derived from the ligands themselves rather than tracked separately: a second registry of which names are synthetic is a second source of truth that can disagree with the first, and the disagreement would surface as a wrong export.

property rbfe_edges: tuple[Transformation, ...]

Selected edges that are relative transformations.

property cbfe_edges: tuple[Transformation, ...]

Selected edges that are counterpoised (paired absolute) calculations.

The two properties partition edges; they are separate because the two kinds are set up and run by different machinery downstream, so almost every consumer wants one or the other rather than the mixed list.

class rbfenetmap.core.models.RejectionReason(*values)

Bases: str, Enum

Why a candidate transformation is infeasible.

A str enum so the values serialize to JSON as themselves and compare equal to plain strings, which keeps the exported network readable without a custom decoder.

class rbfenetmap.core.models.SoftcoreRepair(applied=False, n_fragments_before=(0, 0), n_fragments_after=(0, 0), demoted_1=(), demoted_2=(), iterations=0, rejection=None, trace=())

Bases: object

Outcome of the soft-core connectivity repair for one transformation.

Parameters:

• applied (bool) – Whether any atom was demoted from the common core.

• n_fragments_before (tuple[int, int]) – Soft-core connected-component counts per side, before and after repair. After a successful repair both entries are 0 or 1 – the constraint is at most one region, and an empty soft-core (zero regions) is perfectly legal.

• n_fragments_after (tuple[int, int]) – Soft-core connected-component counts per side, before and after repair. After a successful repair both entries are 0 or 1 – the constraint is at most one region, and an empty soft-core (zero regions) is perfectly legal.

• demoted_1 (tuple[int, ...]) – Atoms moved from the common core into the soft-core, per side.

• demoted_2 (tuple[int, ...]) – Atoms moved from the common core into the soft-core, per side.

• iterations (int) – Repair loop iterations consumed.

• rejection (RejectionReason, optional) – Set when the repair concluded the edge is infeasible.

• trace (tuple[str, ...]) – Human-readable log of each repair step, surfaced by rbfenet inspect. This is the only window into why a given edge grew the soft-core it did, so it is retained even on success.

property succeeded: bool

True when the repair produced at most one soft-core region per side.

property n_demoted: int

Total atoms demoted across both sides.

class rbfenetmap.core.models.Transformation(source, target, mapping, repair=<factory>, score=<factory>, kind=EdgeKind.RBFE)

Bases: object

A candidate or selected network edge: an alchemical transformation.

Parameters:

• source (str) – Ligand names. Must differ.

• target (str) – Ligand names. Must differ.

• mapping (AtomMapping) – The common-core / soft-core partition, after repair.

• repair (SoftcoreRepair) – What the repair did to get there.

• score (EdgeScore) – The cost and feasibility verdict.

• kind (EdgeKind) – Which alchemical experiment this edge stands for. Defaults to EdgeKind.RBFE, so every existing construction site keeps its meaning.

property key: str

The directed "source~target" key.

property unordered_key: tuple[str, str]

The endpoints as a sorted pair.

Selection is undirected – the free energy of a transformation is antisymmetric, so a -> b and b -> a are the same experiment. Direction only matters when writing files and assigning timask/scmask. Keying the candidate pool by this rather than by key is what keeps the two from being double-counted.

property feasible: bool

Whether this edge may be selected.

reversed()

Return the transformation with its direction flipped.

The mapping, fragment counts, and demoted-atom lists are all swapped along with the endpoints, so the result stays internally consistent rather than describing the old direction under new labels.

The repair trace is free-form text written during the repair and cannot be rewritten, so a note is prepended recording that its “side 1” and “side 2” refer to the original orientation. Without it a reader comparing the trace against the edge’s reported soft-core sizes sees them transposed and reasonably concludes one of the two is wrong.

Return type:

Transformation

rbfenetmap.core.models.edge_key(source, target)

Return the canonical "source~target" key for a directed edge.

Parameters:

• source (str)

• target (str)

Return type:

str

rbfenetmap.core.models.orient_edge(edge, ligands, direction)

Return edge oriented according to direction.

Parameters:

• edge (Transformation) – The selected edge, in whatever orientation selection left it.

• ligands (Mapping[str, Ligand]) – Must contain both endpoints; consulted only by "heavier_second".

• direction (str) – A EdgeDirection value.

Return type:

Transformation

Notes

Orientation is applied once, after selection, because selection itself is undirected: the free energy of a transformation is antisymmetric, so a -> b and b -> a are the same experiment. Direction only starts to matter when files are written and timask/scmask are assigned.

"fewer_softcore_first" starts from the side that has less to grow, so the transformation builds outward into the larger ligand. For a CBFE edge every atom is soft-core, so that rule degenerates to “smaller ligand first”. That is still the convention one wants – the source is the molecule being decoupled from the site – but it is arrived at by a different route than the rationale above describes, which is worth knowing before touching this.

This lives in the data model rather than in the planner because it is a property of an edge, not of a selection strategy: post-planning surgery orients the edges it adds by exactly the same rule, and a second copy of it would be free to drift.

rbfenetmap.core.models.parse_edge_key(key)

Split a "source~target" key back into its endpoints.

Parameters:

key (str) – An edge key such as "lig_a~lig_b".

Returns:

The (source, target) ligand names.

Return type:

tuple[str, str]

Raises:

ValueError – If key does not contain exactly one separator, or either side is empty. Ligand names cannot contain ~ (see Ligand), so a key with two separators is user error rather than an ambiguous name.

Options

User-tunable options for mapping, repair, scoring, and network selection.

All frozen dataclasses. Conflicting knobs are rejected here, at construction, rather than deep inside the planner: a user who asks for a connected 12-ligand network with 8 edges should be told immediately, not after the mapping stage has burned several minutes.

rbfenetmap.core.options.COMPAT_LEVELS: tuple[Literal['v0.4'], ...] = ('v0.4',)

Released behaviours a run can be pinned to. Versioned rather than a single legacy flag: “legacy” stops meaning anything the moment there are two of them, and the whole point of the mechanism is to still be unambiguous several releases from now.

rbfenetmap.core.options.DESIGN_CRITERIA: tuple[Literal['none', 'a_optimal', 'd_optimal'], ...] = ('none', 'a_optimal', 'd_optimal')

Statistical design criteria, plus the "none" that means “select on cost alone”. "none" is a member rather than an Optional because every other selection knob in this class spells “off” as a value, and one knob that spells it as None would be a second convention for the same idea.

rbfenetmap.core.options.CONSISTENCY_SCOPES: tuple[Literal['pairwise', 'component', 'graph'], ...] = ('pairwise', 'component', 'graph')

How widely one core – and therefore one soft-core – is required per ligand, ordered from loosest to strictest. A scope’s position is what callers test against, rather than a chain of equality checks that would drift as scopes are added.

class rbfenetmap.core.options.AlignmentOptions(method='mcs', reference=None, min_mcs_atoms=3, max_matches=200)

Bases: object

Controls the optional pre-alignment of a ligand set into a common frame.

Parameters:

• method ({"mcs", "o3a"}) – "mcs" (default) fits each ligand onto an already-aligned neighbour through their maximum common substructure, which gives an auditable set of atoms and a residual RMSD that means something. "o3a" uses RDKit’s Open3DAlign, which needs no shared substructure and is the fallback for a set too diverse for an MCS to bite on.

• reference (str, optional) – Name of the ligand whose frame everything else is brought into. None picks the ligand with the most heavy atoms, ties broken by name. In a congeneric series the largest ligand usually contains the shared scaffold, so its substructure overlap with every partner is large, and the rule costs no MCS searches to evaluate.

• min_mcs_atoms (int) – Refuse to fit on fewer corresponding atoms than this. A ligand that cannot clear the bar is left in its own frame and reported, rather than moved on the strength of an overlap too small to determine where it should go.

• max_matches (int) – Cap on the substructure embeddings enumerated while resolving a symmetric overlap.

Raises:

ValueError – If the method is unknown, min_mcs_atoms is below three, or max_matches is not positive.

Notes

There is deliberately no "none" method. Alignment is either requested or not requested; a do-nothing member would be a second way to express “off” that every caller downstream would then have to test for.

max_matches defaults well below MappingOptions’ 1000 because the jobs are not comparable. Mapping is choosing the common core an alchemical transformation will actually run, once per candidate edge; this is choosing a frame, once per ligand, and the answer is a rigid motion that a few hundred embeddings pin down as well as a thousand would.

class rbfenetmap.core.options.CorePruningPolicy(demote_element_mismatch=False, demote_degree_mismatch=False, demote_formal_charge_mismatch=True, demote_aromaticity_mismatch=False, demote_ring_membership_mismatch=False, demote_light_element_swap=True)

Bases: object

Which mapped atom pairs to demote before the connectivity repair runs.

Generalizes BuildEdges._classify_softcore_method0/1/2 from a three-way string into independent flags. The named presets reproduce the original three methods.

Parameters:

• demote_element_mismatch (bool) – Demote pairs whose atomic numbers differ (the MCSS-E2 behaviour).

• demote_degree_mismatch (bool) – Demote pairs whose heavy-atom connectivity differs (the MCSS-E behaviour).

• demote_formal_charge_mismatch (bool) – Demote pairs whose formal charges differ. On by default: a charge that changes across the core/soft-core boundary is a common source of unphysical setups.

• demote_aromaticity_mismatch (bool) – Demote pairs disagreeing on aromaticity or ring membership. Off by default; geometry mappers already filter ring/non-ring pairs upstream.

• demote_ring_membership_mismatch (bool) – Demote pairs disagreeing on aromaticity or ring membership. Off by default; geometry mappers already filter ring/non-ring pairs upstream.

• demote_light_element_swap (bool) – Demote across a hydrogen/heavy-atom pairing, taking the attached branch with it.

classmethod preset(name)

Return a named preset: "mcss", "mcss-e", or "mcss-e2".

Parameters:

name (str)

Return type:

CorePruningPolicy

class rbfenetmap.core.options.MappingOptions(timeout=60, ring_matches_ring_only=True, complete_rings_only=True, match_valences=False, match_chiral_tag=False, max_matches=1000, match_selection='fewest_fragments', distance_threshold=2.0, core_pruning=<factory>)

Bases: object

Controls how a mapper proposes an atom correspondence.

Parameters:

• timeout (int) –

  Seconds allowed for a single MCS search.

  This is the memory knob as much as the time knob. FindMCS allocates monotonically while it searches and frees nothing until it returns, at roughly 40 MB per second on drug-sized ligands, so peak usage is about 40 MB/s * timeout * jobs. The default of 60 with jobs=8 is therefore some 20 GB of search structures before a single candidate is retained. Raise it knowingly.

• ring_matches_ring_only (bool) – RDKit FindMCS settings, mirroring the values BuildEdges._find_mcs uses.

• complete_rings_only (bool) – RDKit FindMCS settings, mirroring the values BuildEdges._find_mcs uses.

• match_valences (bool) – Further FindMCS settings.

• match_chiral_tag (bool) – Further FindMCS settings.

• max_matches (int) – Cap on substructure embeddings enumerated when resolving a symmetric core.

• match_selection ({"fewest_fragments", "best_rmsd", "first"}) – How to choose among those embeddings. See the note below.

• distance_threshold (float) – Geometric cutoff, in angstroms, for the geometry-based mappers.

• core_pruning (CorePruningPolicy) – Pre-repair demotions applied to the raw correspondence.

Notes

match_selection defaults to "fewest_fragments" rather than to first-match for a specific reason. BuildEdges._find_mcs calls the singular GetSubstructMatch on each molecule independently and zips the two results together. For any symmetric substructure – a para-substituted ring being the everyday case – the two matches can correspond to different orientations, and the zip then pairs atoms that sit on opposite sides of the ring. The mapping is topologically valid, so nothing complains until the geometry check much later. Enumerating embeddings and picking one by an explicit criterion removes the coin flip.

class rbfenetmap.core.options.NetworkOptions(pair_strategy='all_unordered_pairs', hub=None, explicit_pairs=(), n_edges=None, edges_per_ligand=2, min_cycle_coverage=1.0, forced_edges=(), banned_edges=(), require_connected=True, edge_direction='fewer_softcore_first', prefilter='none', prefilter_k=8, prefilter_min_tanimoto=0.4, selection_objective='uniform_redundancy', cycle_coverage_mode='node', max_cycle_size=None, max_diameter=None, n_redundancy=2, hub_selection='most_partners', pair_evaluation='eager', adaptive_initial_neighbors=3, adaptive_batch_size=32, show_progress=False, jobs=1, consistency='pairwise', cbfe_mode='off', cbfe_base_cost=8.0, cbfe_atom_weight=0.05, cluster_by='none', cluster_bridges=2, design='none', design_candidate_factor=3.0, design_refine=False, design_total_ns=None, design_lambda_min=12, design_lambda_max=24, softcore=<factory>, intermediates=<factory>, compat=None)

Bases: object

Controls candidate generation and final edge selection.

Parameters:

• pair_strategy (PairStrategy) – How candidate pairs are enumerated before scoring.

• hub (str, optional) – Ligand to place at the centre of a star network, or to bias the MST toward.

• explicit_pairs (tuple[str, ...]) – "a~b" specifications used by the explicit strategy.

• n_edges (int, optional) – Cap on the total number of selected edges.

• edges_per_ligand (int) – Target minimum degree for every ligand. Best-effort.

• min_cycle_coverage (float) – Target fraction of ligands lying on at least one cycle. Best-effort. Cycles are what make a network’s free energies checkable against themselves, so this is the knob that buys statistical confidence rather than raw coverage.

• forced_edges (tuple[str, ...]) – "a~b" specifications, normalized to unordered pairs.

• banned_edges (tuple[str, ...]) – "a~b" specifications, normalized to unordered pairs.

• require_connected (bool) – Whether the selected network must span every ligand.

• edge_direction (EdgeDirection) – How each selected edge is oriented once selection is done.

• prefilter ({"none", "fingerprint"}) – Optional similarity prefilter applied before mapping.

• prefilter_k (int) – Neighbours retained per ligand by the prefilter.

• prefilter_min_tanimoto (float) – Similarity floor for the prefilter.

• selection_objective ({"uniform_redundancy", "connectivity_then_cycles"}) – Whether redundancy first tries to raise degree targets uniformly, or instead focuses on putting as many ligands as possible on at least one cycle after the spanning network has been built.

• cycle_coverage_mode ({"node", "edge"}) – What min_cycle_coverage is a fraction of. "node" (the default, and LOMAP’s rule) measures the ligands that lie on at least one cycle. "edge" measures the selected edges that lie on one, which is FEP+’s stated invariant and is exactly 2-edge-connectivity: at coverage 1.0 the network has no bridges at all. The edge form is strictly the harder target – every bridge has covered endpoints as soon as something else puts them on a cycle – so it is opt-in rather than a correction to the node form.

• max_cycle_size (int, optional) – Maximum cycle length allowed when adding redundancy edges to improve cycle coverage. None permits any cycle size.

• max_diameter (int, optional) –

  Target upper bound on the network’s diameter – the longest shortest path between any two ligands, counted in edges. Statistical error accumulates along a path, so LOMAP caps it at 6 and FEP+ below 5. None (the default) imposes no bound.

  Best-effort, like the other redundancy targets: selection here is additive, so the bound is approached by buying shortcut edges rather than, as LOMAP does, by refusing to remove one. A pool with no shortcut left to sell warns and records the shortfall instead of raising.

• n_redundancy (int) – Number of spanning trees the redundant-mst planner overlays. Ignored by every other planner. Konnektor defaults to 2; the paper that introduced the topology uses 3.

• hub_selection ({"most_partners", "min_total_cost"}) – How the star planner picks a hub when none is named. "most_partners" (the default) ranks by feasible partner count and only breaks ties on cost, so it never compares cost across ligands of differing connectivity. "min_total_cost" ranks by summed cost to every other ligand, charging an unreachable partner the worst cost in the pool; that is LOMAP’s pick_lead and HiMap’s ref_lig_gen, which sum a similarity matrix in which an unrelatable pair scores zero. OpenEye’s own documentation calls hub choice the dominant factor in a star map’s performance, which is why it is a knob rather than a constant.

• pair_evaluation ({"eager", "adaptive"}) – Whether to map every candidate before planning, or evaluate fingerprint-ranked batches until the requested network targets are met.

• adaptive_initial_neighbors (int) – Fingerprint-nearest neighbours evaluated per ligand in the first adaptive batch.

• adaptive_batch_size (int) – Maximum number of additional pairs evaluated in each adaptive expansion.

• show_progress (bool) – Write pair-evaluation progress to stderr. Disabled by default for library use; the CLI enables it automatically on interactive terminals.

• jobs (int) – Worker processes used for mapping and scoring.

• consistency ({"pairwise", "component", "graph"}) –

  How widely a ligand is required to hold one common core – and therefore, since the core and the soft-core are a strict partition of the ligand’s atoms, one soft-core. That corollary is the point of the knob: the Amber scmask is the soft-core, so a ligand with one core across its edges has one scmask across them too.

  A ladder from loosest to strictest:

  • "pairwise" (default) – no requirement. Each edge is mapped independently and holds the largest core its own pair supports, which is the cheapest transformation for that pair. A ligand on three edges holds three cores, and three soft-cores.

  • "component" – one core per ligand within each connected component of the RBFE-only selected subgraph. Costs nothing to configure and is the scope that works on a set whose scaffolds do not all map to each other, because such a set fragments into exactly those components anyway.

  • "graph" – one core per ligand across all of its selected RBFE edges. The strongest form, and the one that fails first: the shared core is an intersection over the whole network, so one chemically distant ligand shrinks everyone’s.

  In every case the surviving core is the intersection of the pairwise ones, the rest is demoted, and the soft-core repair re-runs on what remains, iterated to a fixed point.

  Two limits apply to all of them. CBFE edges are exempt – a counterpoised edge has no common core by construction, so reading one as “this ligand’s core is empty here” would erase the core of every ligand a bridge touches. And the pass never re-selects; it can leave a selected edge infeasible, which is raised rather than absorbed. See rbfenetmap.core.consistency.

• cbfe_mode ({"off", "bridge", "cycles", "all"}) –

  How freely the planner may spend counterpoised (CBFE) edges. A CBFE edge needs no atom mapping, so it is available between any two ligands – including the pairs an MCS search cannot relate – at the price of two absolute calculations.

  • "off" – never. Every edge is RBFE.

  • "bridge" – only to join subnetworks the feasible RBFE pool leaves disconnected. This is the mode that turns a hard connectivity failure into a planned network.

  • "cycles" – everything "bridge" does, and additionally to put ligands on a cycle when no RBFE candidate can.

  • "all" – the whole network is CBFE. Mapping is skipped entirely.

  The modes form a strict ladder, so raising the setting only ever adds possibilities.

  Eligibility is a gate applied before cost competition, and this is the point most easily misread. cbfe_base_cost decides which CBFE edge is chosen among the ones the mode makes eligible, and orders RBFE against CBFE inside cycle closure. It never lets a CBFE edge outbid a feasible RBFE edge inside an already connected component: under "bridge" a CBFE edge that does not join two components is not in the pool at all, at any price.

• cbfe_base_cost (float) – Fixed cost of a CBFE edge, on the same scale as the scorer’s edge totals. The default sits at the linear scorer’s charge-change ceiling – a CBFE edge costs about what the most expensive thing that can happen to a still-feasible RBFE edge costs – so CBFE never wins on price alone, only on availability.

• cbfe_atom_weight (float) – Added to cbfe_base_cost for each heavy atom summed over both ligands. A counterpoised calculation decouples both molecules in full, so its expense scales with how much there is to decouple.

• cluster_by ({"none", "charge", "scaffold", "fingerprint"}) –

  Partition the ligands and plan each cluster as its own subnetwork, joined to the others by a few deliberately chosen edges. "none" (default) is the unpartitioned behaviour.

  This is an edge-budget knob, not a feasibility one. The precision floor of an RBFE network goes as n ln n, and that is superlinear, so sum_i n_i ln n_i is strictly smaller than n ln n for any real partition: planning five balanced clusters of twenty to the floor costs roughly 190 edges where one hundred ligands cost 460. Nothing here changes which edges are feasible – cross-cluster mappings are as available as they ever were, they are simply not worth buying in quantity.

  • "charge" – net formal charge classes. Exact, thresholdless, and it isolates the transformation the scorer already penalises hardest.

  • "scaffold" – the Bemis-Murcko framework, which is what “series” usually means.

  • "fingerprint" – average-linkage hierarchical clustering on Tanimoto distance, for a set with neither a clean charge split nor a shared framework.

• cluster_bridges (int) –

  Edges spent joining each pair of clusters that gets joined, when cluster_by is set. Ignored otherwise.

  The default of 2 is deliberate and is the reason this is not simply 1. Two edges between the same two clusters put the crossing itself on a cycle, since the paths inside each cluster close the loop. Cross-cluster edges are the least similar and therefore the least trustworthy edges in the network, so applying the every-edge-in-a-cycle invariant precisely there buys more per edge than anywhere else. Setting 1 gives the minimal spanning join and leaves each crossing unchecked.

• design ({"none", "a_optimal", "d_optimal"}) –

  Statistical design criterion the optimal planner minimises. "none" (default) selects on cost alone, which is what every release up to v0.4 did.

  The criterion is evaluated on the network’s Fisher information matrix, which for a set of relative measurements is the weighted graph Laplacian – see rbfenetmap.core.design. a_optimal minimises the summed variance of the estimates; d_optimal minimises the volume of their joint confidence ellipsoid, which because the Laplacian’s pseudo-determinant counts spanning trees produces a markedly more cyclic network at the same edge count. Prefer d_optimal when a cycle-closure correction will be applied downstream, a_optimal otherwise.

  This is an objective, so it is meaningful only to a planner that optimises it. Naming it alongside a planner that does not is refused rather than ignored; see check_design_support().

• design_candidate_factor (float) – The design’s candidate pool is capped at design_candidate_factor * n_ligands edges – the M = 3m of Xu’s Appendix-H heuristic. Raising it widens the search at quadratic cost in criterion evaluations; the published 1.10x bound is measured at the default of 3.0.

• design_refine (bool) – Run a Fedorov exchange pass after the heuristic, swapping edges one at a time while the criterion improves. Off by default: the heuristic is already within a published 1.10x of the optimum and the refinement costs far more criterion evaluations.

• design_total_ns (float, optional) – Total simulation budget, in nanoseconds, to distribute A-optimally across the selected edges. None (default) emits no allocation at all. Set it and the Amber exporter writes a per-edge lambda-window and nanosecond budget into each .runconfig. Static, computed once from the predicted variances – the iterative refit against measured variances needs a round trip through the MD engine and is not part of this.

• design_lambda_min (int) – Bounds on the per-edge lambda-window count the allocation is mapped onto. The defaults, 12 and 24, bracket what an Amber RBFE edge normally runs at.

• design_lambda_max (int) – Bounds on the per-edge lambda-window count the allocation is mapped onto. The defaults, 12 and 24, bracket what an Amber RBFE edge normally runs at.

• softcore (SoftcorePolicy) – Feasibility policy handed to the repair.

• intermediates (IntermediateOptions) –

  Whether the pipeline may invent ligands to bridge pairs no mapping can relate, and how many. Off by default.

  This is the only knob in the class that changes the vertex set rather than the edge set, which is why it sits between max_softcore_atoms and cbfe_mode in the precedence table: it widens the pool the planner is handed, and it does so before the planner runs, so a gap an intermediate closed is simply not a gap by the time CBFE eligibility is evaluated. That ordering is the whole of “stay relative, then fall back to counterpoised” – there is no precedence flag behind it.

• compat (str, optional) – The released behaviour this run was pinned to, or None. Set by preset(); recorded so a planned network states which behaviour produced it. Purely a label – it changes nothing on its own, because preset() has already written the values it stands for.

Raises:

ValueError – If an edge appears in both the forced and banned sets, if a knob is out of range, or if compat names an unknown level.

classmethod preset(level, **overrides)

Return the options a released version of the package planned with.

Parameters:

• level (str) – A member of COMPAT_LEVELS.

• **overrides – Applied on top of the pinned values. Intended for the settings that describe this run rather than this behaviour – the ligand-specific intent (hub, forced_edges, banned_edges, explicit_pairs) and the operational knobs (jobs, show_progress). Overriding an algorithmic knob is permitted here and rejected at the CLI, where the user’s intent is unambiguous enough to call it a contradiction.

Returns:

With compat set to level.

Return type:

NetworkOptions

Raises:

ValueError – If level is unknown.

Notes

Every value below is written out literally, and that is the entire point. Building this from the dataclass defaults would be shorter and would defeat the mechanism: the moment a later release moves a default, the preset would move with it and silently stop reproducing the version it names. These numbers are a record of what v0.4.0 did, not a view onto what the current code does, so they must be edited only to fix a transcription error – never to track a new default.

The pinned surface is the algorithmic one. Ligand-specific intent is not pinned: banning an edge or naming a hub is a statement about one ligand set, not about a version’s behaviour, so those stay available alongside a compat level.

property forced_pairs: frozenset[tuple[str, str]]

Forced edges as unordered endpoint pairs.

property banned_pairs: frozenset[tuple[str, str]]

Banned edges as unordered endpoint pairs.

property cbfe_bridges_components: bool

Whether CBFE edges may join otherwise-disconnected subnetworks.

property cbfe_closes_cycles: bool

Whether CBFE edges may be spent putting ligands on a cycle.

property generates_intermediates: bool

Whether the pipeline may invent ligands for this run.

intermediate_headroom(n_ligands)

Return how many ligands may be invented before n_edges runs out.

Parameters:

n_ligands (int) – Size of the real ligand set.

Returns:

None when n_edges is unset, so nothing constrains generation here. Otherwise n_edges - (n_ligands - 1), possibly zero or negative.

Return type:

int or None

Notes

The budget is spent, never inflated. Every invented ligand is another vertex, so a spanning network over the augmented set needs one more edge than it did before. Quietly raising n_edges to pay for a molecule the user never asked for would be exactly the silent over-spend check_edge_budget() refuses to make in the other direction. When the headroom runs out, generation stops and says so on unmet_constraints; it does not raise, because unlike a spanning tree that cannot fit, an intermediate that cannot fit leaves a perfectly valid network.

check_edge_budget(n_ligands)

Verify n_edges can support a spanning network over n_ligands.

Raises:

ValueError – If n_edges is below n_ligands - 1 while connectivity is required.

Parameters:

n_ligands (int)

Return type:

None

Notes

This is the single most likely knob conflict, and it is a hard error rather than a silent override in either direction. Trimming the spanning tree to honour n_edges would produce a disconnected network the user explicitly forbade; quietly raising n_edges would ignore a budget the user explicitly set. Only the user can say which they meant.

With cbfe_mode bridging enabled this becomes the only way a spanning network can fail. A CBFE edge exists between every pair, so the candidate pool can no longer be too sparse to connect the ligands – only an edge budget below n_ligands - 1, or a ban on every bridging pair, can prevent it.

class rbfenetmap.core.options.SoftcorePolicy(ring_policy='ring_system', max_softcore_atoms=12, max_softcore_fraction=0.6, min_core_atoms=4, min_mcs_fraction=0.35, core_rmsd_threshold=2.0, charge_change_policy='penalize', max_iterations=None, core_pruning=<factory>)

Bases: object

Controls the soft-core connectivity repair and the feasibility budget.

Parameters:

• ring_policy ({"ring_system", "none"}) – "ring_system" (default) never leaves a ring half soft-core: touching any ring atom absorbs the whole ring, and fused systems cascade. "none" permits half-broken rings, which is what a deliberate ring-opening study needs.

• max_softcore_atoms (int) – Reject when either side’s heavy soft-core exceeds this. The single most effective knob for controlling how aggressive the repair is allowed to be.

• max_softcore_fraction (float) – Reject when either side’s heavy soft-core exceeds this fraction of that molecule. Catches the small-ligand case that an absolute count misses.

• min_core_atoms (int) – Reject when fewer than this many heavy atoms remain in the common core.

• min_mcs_fraction (float) – Reject before repair when the core covers less than this fraction of the smaller molecule. Cheap scaffold-hop filter that avoids wasting Steiner work.

• core_rmsd_threshold (float) – Reject when the mapped core’s in-place RMSD exceeds this, in angstroms.

• charge_change_policy ({"allow", "penalize", "reject"}) – How to treat a net formal charge change across the edge.

• max_iterations (int, optional) – Repair loop bound. None derives it from the molecule sizes, which is already a proven upper bound – this is only a backstop.

• core_pruning (CorePruningPolicy)

rbfenetmap.core.options.normalize_edge_specs(specs)

Turn "a~b" edge specifications into a set of unordered endpoint pairs.

Selection is undirected – the free energy of a transformation is antisymmetric, so a -> b and b -> a name the same experiment. Normalizing forced and banned edges to unordered pairs here means a user who writes --banned-edge b~a gets the ban they intended rather than one that silently misses.

Parameters:

specs (tuple[str, ...] | list[str] | None)

Return type:

frozenset[tuple[str, str]]

Soft-core repair

Soft-core connectivity repair.

This module enforces the constraints the whole package is organised around: a transformation has at most one connected soft-core region per side, and each region attaches to the common core through exactly one bond. A mapper is free to return a correspondence whose unmapped atoms fall into several disconnected pieces – that is the normal outcome for, say, a benzene to para-xylene transformation, where two hydrogens on opposite sides of the ring both disappear. Such a partition cannot be run as a single alchemical transformation, so it must either be repaired or rejected.

The repair works by demoting common-core atoms into the soft-core until the pieces join up. Choosing which atoms to demote is a node-weighted Steiner tree problem: the soft-core fragments are the terminals, the bond graph is the network, and the cost of recruiting an atom is how much soft-core that recruitment ultimately drags in.

A second kind of demotion has nothing to do with fragmentation. The mapper may keep an atom in the common core that is reachable from the rest of the core only through the soft-core – a terminal methyl retained while the carbon joining it goes soft. The soft-core is then a bridge, core -- soft-core -- core, and attaches by two bonds however connected it is. Such a stranded fragment is absorbed into the soft-core, which restores the single-attachment invariant at the cost of a few atoms. It is not a closure rule: it is checked and applied in the same loop as Steiner recruitment, because bridging two regions can strand a fragment that was not stranded before.

Three closure rules run to a fixpoint after every recruitment:

(a) whole-ring – a ring is never left half soft-core. Touching one ring atom absorbs the ring. Fused systems then cascade on their own, because absorbing one ring pulls in the atoms it shares with its neighbours, which makes those neighbours intersected-but-not-contained on the next sweep.

(b) hydrogen-follows-parent – a hydrogen joins the soft-core when its heavy parent does. Deliberately one-way: a soft-core hydrogen whose parent is common core stays put as a one-atom region. That asymmetry is not an oversight. R-H -> R-CH3 is the single most common transformation in the field, and its soft-core on the R-H side is exactly one hydrogen attached to a core carbon. A two-way rule would demote that carbon, then its ring, and destroy the edge.

(c) mapped-partner – demoting an atom demotes whatever it is mapped to. This is the only rule that couples the two molecules, and it is why the repair is genuinely joint: fixing a fragmentation on side 1 can create a new one on side 2, which the loop must then fix in turn.

class rbfenetmap.core.softcore.RepairContext(graph_1, graph_2, rings_1, rings_2, hydrogen_parent_1, hydrogen_parent_2, forward, reverse, heavy_1, heavy_2, policy, n_atoms_1=0, n_atoms_2=0, _cost_cache=<factory>)

Bases: object

Everything the repair needs about one candidate pair, precomputed.

Built once per edge by build(). Holding the bond graphs, ring lists, hydrogen parentage, and demotion costs here keeps them out of the repair loop, which would otherwise recompute them on every iteration.

Parameters:

• graph_1 (Graph)

• graph_2 (Graph)

• rings_1 (tuple[frozenset[int], ...])

• rings_2 (tuple[frozenset[int], ...])

• hydrogen_parent_1 (dict[int, int])

• hydrogen_parent_2 (dict[int, int])

• forward (dict[int, int])

• reverse (dict[int, int])

• heavy_1 (frozenset[int])

• heavy_2 (frozenset[int])

• policy (SoftcorePolicy)

• n_atoms_1 (int)

• n_atoms_2 (int)

• _cost_cache (dict[tuple[int, int], float])

classmethod build(source, target, mapping, policy)

Precompute the graphs, rings, hydrogen parentage, and correspondence.

Raises:

RepairError – If the mapping’s atom counts disagree with the ligands. This is a programming error rather than a chemistry one – a mapping built for a different pair.

Parameters:

• source (Ligand)

• target (Ligand)

• mapping (AtomMapping)

• policy (SoftcorePolicy)

Return type:

RepairContext

side(side)

Return (graph, rings, hydrogen_parents, heavy_indices) for side.

Parameters:

side (int)

Return type:

tuple[Graph, tuple[frozenset[int], …], dict[int, int], frozenset[int]]

n_heavy(side)

Heavy-atom count for side.

Parameters:

side (int)

Return type:

int

demote_cost(atom, side)

Cost of demoting atom, measured as the closure it triggers.

Parameters:

• atom (int) – The common-core atom under consideration.

• side (int) – 1 or 2.

Returns:

Total atoms across both molecules that demoting this one atom ultimately pulls into the soft-core.

Return type:

float

Notes

Using the closure size as the cost is what makes the Steiner search behave chemically without any hand-tuned table of per-element weights. A hydrogen costs about 2 (itself and its partner). A peripheral heavy atom costs a little more. An aromatic carbon costs its entire fused ring system, plus every attached hydrogen, plus all of their partners on the other side – so the solver routes around rings whenever an acyclic path exists, and only pays for a ring when there is no alternative.

The cost is measured in isolation, from an empty soft-core, so it slightly overestimates once part of the closure is already soft-core. It is a search heuristic, not an accounting of the final result, and computing it once keeps the repair loop cheap.

cost_fn(side)

Return a one-argument cost callable bound to side.

Parameters:

side (int)

Return type:

Callable[[int], float]

rbfenetmap.core.softcore.detect_fragments(graph, softcore)

Return the connected components of the soft-core, largest first.

An empty soft-core yields an empty list. Zero regions is legal – the constraint is at most one region, and a transformation that only reorders a common core has none.

Parameters:

• graph (Graph)

• softcore (Iterable[int])

Return type:

list[set[int]]

rbfenetmap.core.softcore.joint_closure(softcore_1, softcore_2, context)

Apply the three closure rules to a fixpoint over both sides.

Parameters:

• softcore_1 (set[int]) – Current soft-core atom sets. Not modified in place.

• softcore_2 (set[int]) – Current soft-core atom sets. Not modified in place.

• context (RepairContext) – Precomputed graphs, rings, hydrogen parentage, and correspondence.

Returns:

The closed soft-core sets.

Return type:

tuple[set[int], set[int]]

Notes

Terminates because every rule only ever adds atoms, and the atom sets are finite.

rbfenetmap.core.softcore.precheck_mapping(source, target, mapping, policy)

Cheap rejections applied before the repair runs.

Parameters:

• source (Ligand) – The two ligands.

• target (Ligand) – The two ligands.

• mapping (AtomMapping) – The mapper’s raw output.

• policy (SoftcorePolicy) – Feasibility thresholds.

Return type:

RejectionReason or None

Notes

Ordering matters for cost, not just for message quality. A scaffold hop whose common core covers a tenth of either molecule will certainly fail the soft-core budget, but only after the Steiner solver has done real work on a large fragmented soft-core. Catching it on the MCS fraction first skips that entirely.

rbfenetmap.core.softcore.repair_softcore_connectivity(source, target, mapping, policy=None)

Repair the mapping so each side has at most one connected soft-core region.

Parameters:

• source (Ligand) – The two ligands.

• target (Ligand) – The two ligands.

• mapping (AtomMapping) – The mapper’s output, whose soft-core may be fragmented.

• policy (SoftcorePolicy, optional) – Feasibility thresholds and the ring policy. Defaults are used if omitted.

Returns:

The repaired mapping and a record of what was done. On rejection the original mapping is returned unchanged alongside a repair carrying the RejectionReason: an edge that will not be used should not be silently mutated, and the caller may still want to show the user the mapping that failed.

Return type:

tuple[AtomMapping, SoftcoreRepair]

Raises:

rbfenetmap.core.exceptions.RepairError – Only for malformed input – a mapping that does not describe these molecules, or a molecule whose bond graph is disconnected. An edge that genuinely cannot be repaired is not an error; it comes back as a rejection.

Notes

The loop terminates in at most n_atoms_1 + n_atoms_2 iterations. Both soft-core sets grow monotonically within finite atom sets, and any iteration that does not return adds at least one atom: a bridge joining two or more fragments must contain an atom that is not already soft-core, and an iteration entered only because something was stranded absorbs a non-empty component.

rbfenetmap.core.softcore.softcore_attachment_edges(graph, softcore)

Return bonds crossing from the soft-core to the common core.

Each tuple is oriented (softcore_atom, common_core_atom). Counting edges, rather than distinct common-core atoms, expresses the alchemical topology rule directly: one soft-core region must be a singly attached substituent, not a bridge or ring path.

Parameters:

• graph (Graph)

• softcore (Iterable[int])

Return type:

list[tuple[int, int]]

Pre-repair demotion of mapped atom pairs that should not be common core.

An MCS is a topological answer to a topological question, and it will happily map a carbon onto a nitrogen or a CH2 onto a CH3 because the graphs match. Whether such a pair may be held fixed through an alchemical transformation is a different question, and the answer depends on what the user is willing to accept.

This module generalizes BuildEdges._classify_softcore_method0/1/2 from a three-way method string into independent flags on CorePruningPolicy, with presets reproducing the original three. It runs before rbfenetmap.core.softcore, so anything demoted here becomes part of the fragmentation problem the repair then has to solve.

rbfenetmap.core.coreprune.choose_softcore_branch(graph, core, center)

Return the atoms of the least-conserved branches hanging off center.

Partitions the graph into the branches reachable from center through acyclic bonds, counts how many common-core atoms each branch holds, and returns the atoms of every branch below the maximum.

Parameters:

• graph (networkx.Graph) – Bond graph of the molecule.

• core (set[int]) – Atom indices currently in the common core.

• center (int) – The atom whose substituents are being classified.

Returns:

Atoms to demote to the soft-core.

Return type:

set[int]

Notes

This is a corrected port of BuildEdges._choose_sc_atoms. The original partitions across every bond partner via _atoms_beyond_bond, which blocks only the central atom. For a ring atom that traversal wraps around the ring and returns, so each “branch” contains nearly the whole molecule and the branches overlap almost completely. Demoting all-but-the-largest then demotes almost everything.

Restricting the partition to acyclic bonds keeps the branches genuinely disjoint. A ring atom simply yields fewer branches – possibly none, in which case nothing is demoted here and the whole-ring closure rule in rbfenetmap.core.softcore deals with the ring.

rbfenetmap.core.coreprune.prune_core(source, target, mapping, policy=None)

Demote mapped pairs that the policy says cannot be held in common.

Parameters:

• source (Ligand) – The two ligands.

• target (Ligand) – The two ligands.

• mapping (AtomMapping) – The mapper’s raw correspondence.

• policy (CorePruningPolicy, optional) – Which demotion rules to apply. Defaults are used if omitted.

Returns:

A mapping with the offending pairs moved into the soft-core. Returned unchanged when nothing is demoted, so the common case costs almost nothing.

Return type:

AtomMapping

Notes

A light-element swap (a hydrogen mapped onto a heavy atom) is handled specially: demoting just the pair would leave the heavy atom’s substituents attached to a soft-core atom while remaining common core themselves. The branch-selection helper takes the less-conserved substituents along, which is what _classify_softcore_method1/2 intends.

Molecular graphs

Molecular graph utilities, including the node-weighted Steiner tree solver.

Everything here operates on networkx.Graph objects whose nodes are atom indices. Only mol_to_graph(), ring_systems(), and hydrogen_parents() touch RDKit; the rest – crucially node_weighted_steiner(), which is the heart of the soft-core repair – is pure graph theory and is unit-testable against hand-built graphs with no chemistry involved.

This generalizes cartograph._connected_subsets and the ParmEd-based traversal helpers in BuildEdges (_atoms_beyond_bond, _partition_across_atom) onto a single graph representation.

rbfenetmap.core.molgraph.component_beyond_bond(graph, keep, start)

Return every node reachable from start without passing through keep.

The graph-native equivalent of BuildEdges._atoms_beyond_bond.

Note that for a ring bond this does not partition the molecule: the traversal wraps around the ring and comes back, so the result contains almost everything. Callers that want genuine branches must use acyclic_branches() instead – see the note there.

Parameters:

• graph (Graph)

• keep (int)

• start (int)

Return type:

set[int]

rbfenetmap.core.molgraph.connected_components_of(graph, nodes)

Return the connected components of the subgraph induced by nodes.

Sorted largest first, then by smallest member, so the result is deterministic.

This is the rdkit-free generalization of cartograph._connected_subsets.

Parameters:

• graph (Graph)

• nodes (Iterable[int])

Return type:

list[set[int]]

rbfenetmap.core.molgraph.hydrogen_parents(mol)

Map each terminal hydrogen index to its heavy-atom neighbour.

Bridging hydrogens (degree > 1) are excluded: they have no single parent, and the hydrogen-follows-parent rule is not well defined for them.

Parameters:

mol (Chem.Mol)

Return type:

dict[int, int]

rbfenetmap.core.molgraph.mol_to_graph(mol)

Return the bond graph of mol.

Nodes are atom indices carrying element (atomic number), is_ring, and degree attributes; edges are bonds carrying in_ring.

Parameters:

mol (rdkit.Chem.Mol)

Return type:

networkx.Graph

rbfenetmap.core.molgraph.node_weighted_steiner(graph, terminals, cost)

Find a cheap set of nodes connecting every terminal fragment.

Parameters:

• graph (networkx.Graph) – The bond graph.

• terminals (Sequence[set[int]]) – Disjoint node sets to be joined. Fewer than two means there is nothing to do.

• cost (Callable[[int], float]) – Cost of recruiting a node. Terminal nodes are never charged for.

Returns:

The nodes to recruit, and whether the result came from the approximate solver. The flag is propagated into the repair trace: an approximate bridge is still valid, but it is not guaranteed reproducible across networkx versions, and a user comparing two runs deserves to know which is which.

Return type:

tuple[set[int], bool]

Raises:

ValueError – If the terminals cannot be connected at all – i.e. the molecule itself is disconnected, which no valid ligand should be.

Notes

Two terminals is a shortest-path problem and is solved exactly by Dijkstra on the node-split digraph.

Three or more terminals is NP-hard, and is solved by iterative cheapest merge: repeatedly find the cheapest node-weighted path joining any two of the current components, recruit its interior nodes, and merge. This is the classic greedy Steiner heuristic, costing O(k^2) Dijkstra runs for k fragments – trivial at molecular scale, where k is rarely above five.

An earlier version enumerated candidate subsets exhaustively for small instances, claiming exactness. That was a mistake: “small” was bounded at 25 candidate nodes, and C(25, 12) is 5.2 million subsets, so real ligands (20-32 candidates, 3-4 fragments) hung rather than solving. Exponential search is not viable here even at molecular size.

Ties are broken by sorted node order everywhere, so repeated runs on the same input give the same answer regardless of dictionary or set iteration order.

rbfenetmap.core.molgraph.ring_systems(mol)

Return the SSSR rings of mol as frozensets of atom indices.

Individual rings, not merged fused systems. That is deliberate: the whole-ring closure rule in rbfenetmap.core.softcore iterates to a fixpoint, so absorbing one ring of a fused system pulls in the shared atoms, which makes the next ring intersected-but-not-contained, which absorbs it in turn. Fused systems therefore cascade automatically and need no special case – matching the conservative behaviour of cartograph._filter_fused_rings without duplicating its logic.

Parameters:

mol (Chem.Mol)

Return type:

tuple[frozenset[int], …]

rbfenetmap.core.molgraph.acyclic_branches(graph, center)

Partition the graph into the branches hanging off center by acyclic bonds only.

Parameters:

• graph (networkx.Graph) – The bond graph.

• center (int) – The atom whose substituents are being separated.

Returns:

{neighbour: nodes_in_that_branch}, one entry per neighbour reached through a bond that is not in a ring. Branches are guaranteed disjoint.

Return type:

dict[int, set[int]]

Notes

This deliberately differs from BuildEdges._partition_across_atom, which traverses across every bond partner while blocking only the central atom. For a ring atom that traversal wraps around the ring, so each “branch” contains nearly the whole molecule and the branches overlap almost completely. The downstream heuristic then demotes all-but-the-largest branch, which is to say almost the entire molecule.

Restricting the partition to acyclic bonds keeps the branches genuinely disjoint. A ring atom simply has fewer branches (possibly none), and the whole-ring closure rule in rbfenetmap.core.softcore handles the ring itself.

rbfenetmap.core.molgraph.stranded_components(graph, removed)

Return the nodes that removed cuts off from the graph’s main body.

Deleting removed may split the remainder into several pieces. The largest is taken to be the main body and every other piece is “stranded” – reachable from the rest of the graph only by passing through removed.

Determinism comes from connected_components_of(), which orders by size and then by smallest member, so an exact tie in size resolves the same way on every run and in every networkx version. That matters: the caller uses this to decide which atoms move into a soft-core, and a report that reshuffled between runs would not be diffable.

Parameters:

• graph (networkx.Graph)

• removed (Iterable[int]) – Nodes to delete before looking for components.

Returns:

The union of every component but the largest. Empty when removed leaves the remainder connected, which is the ordinary case.

Return type:

set[int]

Rigid-body superposition and core-geometry metrics.

NumPy only – no RDKit, no ParmEd. Callers pass coordinate arrays extracted from whatever molecule representation they hold.

The core RMSD computed here is a genuine quality signal for a mapping, not just a diagnostic. Two molecules can share a large maximum common substructure whose atoms sit in completely different places once the ligands are posed in the binding site; such a mapping is topologically defensible and physically useless. A high core RMSD is how that shows up.

rbfenetmap.core.kabsch.apply_transform(coords, rotation, translation)

Return coords rotated and translated by a rigid_transform() result.

Parameters:

• coords (numpy.ndarray) – (n, 3) coordinates. Need not be the points the transform was fitted on.

• rotation (numpy.ndarray) – A (3, 3) rotation matrix.

• translation (numpy.ndarray) – A length-3 translation vector.

Returns:

The transformed (n, 3) coordinates.

Return type:

numpy.ndarray

rbfenetmap.core.kabsch.core_rmsd(mobile, reference, *, superpose_first=False)

RMSD between corresponding points.

Parameters:

• mobile (numpy.ndarray) – (n, 3) coordinate arrays of corresponding atoms.

• reference (numpy.ndarray) – (n, 3) coordinate arrays of corresponding atoms.

• superpose_first (bool, optional) – Whether to rigidly superpose before measuring. Default False.

Returns:

The RMSD, or 0.0 for an empty core.

Return type:

float

Notes

The default is not to superpose. Ligands are normally supplied already posed in a common binding-site frame, and in that frame the in-place deviation of the mapped core is the physically meaningful quantity: it says whether the mapping pairs atoms that actually occupy the same region of the pocket. Superposing first would discard precisely that information and reward a mapping that is self-consistent but misplaced. Pass superpose_first=True only when the inputs are not co-posed.

Superposing here measures a single pair in isolation and throws the transform away. To bring a whole set of ligands into a common frame before planning – the case where the inputs were prepared separately, for instance converted to mol2 from independent Amber topologies – use rbfenetmap.core.align instead.

rbfenetmap.core.kabsch.kabsch_rotation(mobile, reference)

Return the rotation matrix best superposing mobile onto reference.

Both arrays must already be centred on their own centroids.

Parameters:

• mobile (numpy.ndarray) – (n, 3) centred coordinate arrays.

• reference (numpy.ndarray) – (n, 3) centred coordinate arrays.

Returns:

A (3, 3) proper rotation matrix (determinant +1).

Return type:

numpy.ndarray

Notes

The determinant correction is what keeps this a rotation rather than a rotoinversion. Without it, a near-planar or otherwise degenerate core can superpose onto its own mirror image, giving a flatteringly low RMSD for a mapping that is in fact chirally wrong.

rbfenetmap.core.kabsch.pair_distances(mobile, reference)

Per-pair distances between corresponding points.

Returns:

A length-n array. Useful for finding the single worst-placed mapped pair, which is often more diagnostic than the aggregate RMSD.

Return type:

numpy.ndarray

Parameters:

• mobile (ndarray)

• reference (ndarray)

rbfenetmap.core.kabsch.rigid_transform(mobile, reference)

Return the (rotation, translation) best superposing mobile onto reference.

Parameters:

• mobile (numpy.ndarray) – (n, 3) coordinate arrays of corresponding points.

• reference (numpy.ndarray) – (n, 3) coordinate arrays of corresponding points.

Returns:

• rotation (numpy.ndarray) – A (3, 3) proper rotation matrix.

• translation (numpy.ndarray) – A length-3 vector, such that mobile @ rotation.T + translation is the superposed result.

Raises:

ValueError – If the two arrays do not have the same shape.

Return type:

tuple[ndarray, ndarray]

Notes

This is the piece superpose() cannot give you, and the reason it exists. The fit is computed over a subset of corresponding atoms – typically a common core – but the transform it yields is a property of the whole rigid body, so it can be applied to every atom of the mobile molecule, including the ones that took no part in the fit. Superposing coordinate arrays in isolation re-centres only the points it was handed, which is correct for measuring an RMSD and useless for moving a molecule.

rbfenetmap.core.kabsch.superpose(mobile, reference)

Return mobile rigidly superposed onto reference.

Parameters:

• mobile (numpy.ndarray) – (n, 3) coordinate arrays of corresponding points.

• reference (numpy.ndarray) – (n, 3) coordinate arrays of corresponding points.

Returns:

The transformed mobile coordinates.

Return type:

numpy.ndarray

See also

rigid_transform

Returns the transform itself, for applying to further points.

Bring a ligand set into a common frame.

The rest of this package assumes ligands arrive co-posed, and measures the mapped core’s RMSD in place so that a mapping pairing atoms in different parts of the pocket is caught rather than flattered. That assumption breaks for structures prepared separately – ligands set up individually for ABFE runs, then written to mol2 from their own Amber topologies, each sitting wherever its own simulation box put it. The conformations are real bound poses; only the frames disagree. Left alone, every candidate edge between them is rejected for geometry, which is the correct answer to the wrong question.

This module answers the right one. It fits each ligand onto an already-aligned neighbour and moves it there rigidly, so the in-place core RMSD downstream measures conformational difference rather than an accident of where each box was centred.

What it cannot do is worth stating as plainly. Rigid alignment recovers a common frame, never a common conformation. Each independently relaxed structure keeps its own ring puckers, exocyclic torsions, and bond-length noise, so a residual core RMSD survives alignment and should. The per-ligand records returned here are what let a caller tell the two apart.

class rbfenetmap.core.align.AlignmentResult(ligands, reference, records)

Bases: object

The aligned ligands, plus an account of how each one got there.

Parameters:

• ligands (tuple[Ligand, ...]) – In the order they were supplied, not the order they were aligned in.

• reference (str) – The ligand whose frame the set now shares.

• records (tuple[LigandAlignment, ...]) – One per ligand, in the same order as ligands.

property failures: tuple[LigandAlignment, ...]

Records for ligands left in their own frame.

property median_rmsd: float

Median post-fit RMSD over the ligands that were actually moved.

Returns:

0.0 when nothing was moved. This is the number to quote when advising on core_rmsd_threshold: it is what alignment achieved, not what it hoped for.

Return type:

float

class rbfenetmap.core.align.LigandAlignment(name, reference, method, n_fit_atoms, rmsd, ok, note='')

Bases: object

What alignment did to one ligand, and how well it worked.

Parameters:

• name (str) – The ligand this record describes.

• reference (str, optional) – The already-aligned ligand it was fitted onto. None for the root, which defines the frame and is therefore never moved.

• method ({"mcs", "o3a", "reference", "none"}) – How the transform was obtained. "reference" marks the root; "none" marks a ligand left in its own frame because no usable fit was found.

• n_fit_atoms (int) – How many corresponding atoms the transform was fitted on. Small values deserve suspicion even when the RMSD looks good.

• rmsd (float) – RMSD over those atoms after the fit, in angstroms.

• ok (bool) – Whether the ligand was moved into the common frame.

• note (str) – Why not, when ok is false, or a warning about the fit when it is true.

as_metadata()

Return a JSON-safe record for metadata.

Returns:

Plain types only, so it survives the network JSON round trip unchanged.

Return type:

dict

rbfenetmap.core.align.align_ligands(ligands, *, options=None, mapping_options=None)

Bring a ligand set into a common frame.

Parameters:

• ligands (Sequence[Ligand] or Mapping[str, Ligand]) – The set to align. Order is preserved in the result.

• options (AlignmentOptions, optional) – Method, reference, and the fit thresholds. Defaults to AlignmentOptions.

• mapping_options (MappingOptions, optional) – Supplies the FindMCS settings, so that alignment maximises the same substructure the mapper will later work from.

Return type:

AlignmentResult

Raises:

ValueError – If the set is empty, or options.reference names a ligand that is not in it.

Notes

A ligand that cannot be fitted is left in its own frame, recorded with ok=False, and logged – not raised over. This package’s established line is that an infeasible edge is data rather than an error, and an unalignable ligand is the same kind of fact: its edges will be rejected for geometry, which is now an honest answer rather than a mystery, and the rest of the set still gets a usable network.

rbfenetmap.core.align.choose_reference(ligands, requested=None)

Pick the ligand whose frame the set should adopt.

Parameters:

• ligands (Mapping[str, Ligand])

• requested (str, optional) – An explicit choice. Validated against ligands.

Return type:

str

Raises:

ValueError – If ligands is empty, or requested names a ligand that was not loaded.

Notes

The automatic rule is the ligand with the most heavy atoms, ties broken by name. In a congeneric series the largest member usually contains the shared scaffold, so its overlap with every partner is large. The name in the sort key is not decoration: it is what makes the choice reproducible when the same set is supplied in a different order.

Pipeline

The pipeline: map, repair, score, bridge, plan.

build_network() is the package’s main entry point. Everything the CLI does, and everything an embedding program needs, goes through here.

The stage that most shapes the result is the second one. A mapper is allowed to return a fragmented soft-core; the repair either fixes it or rejects the edge. Rejection is a normal outcome recorded on the candidate, never an exception – one impossible pair among several hundred must not abort a run.

The fourth stage, augment_with_intermediates(), is the newest and the only one that changes the vertex set. It is off by default and runs between scoring and planning, where it can see which pairs the first three stages could not relate and hand exactly those to a generator. Why there rather than per-pair, inside build_candidate(), or as a second pass over an augmented ligand set is argued in that function’s docstring; the short version is that whether an intermediate is worth making is a question about the network, and only a stage that runs once over the settled pool can answer it.

class rbfenetmap.core.pipeline.AugmentationResult(ligands, candidates, records=(), unmet_constraints=())

Bases: object

What intermediate generation added to the pipeline’s inputs.

Parameters:

• ligands (Mapping[str, Ligand]) – The real ligands followed by every invented one, in acceptance order. The same mapping object as the input when generation was off, so the default path pays nothing for the stage existing.

• candidates (tuple[Transformation, ...]) – The original pool plus the sub-edges of accepted proposals – infeasible ones included, for the same reason the pipeline keeps every other rejection.

• records (tuple[IntermediateRecord, ...], optional) – One per gap attempted, in the order attempted.

• unmet_constraints (tuple[str, ...], optional) – Best-effort budgets generation could not satisfy, phrased for the user and merged into the planned network’s own list.

Notes

A result object rather than a mutated network because generation runs before the planner: there is no network yet to mutate, and building one only to replan over it would discard the rejections that justified the intermediates in the first place.

property synthetic_names: tuple[str, ...]

Names of the invented vertices, in acceptance order.

rbfenetmap.core.pipeline.augment_with_intermediates(ligands, candidates, generator, mapper, scorer, mapping_options, network_options)

Invent ligands for the gaps no mapping could cross, and score the new sub-edges.

The fifth stage, between scoring and planning: map -> repair -> score -> bridge -> plan.

Parameters:

• ligands (Mapping[str, Ligand]) – The real ligands.

• candidates (Sequence[Transformation]) – The settled candidate pool, feasible and not. Its rejections are what define the gaps.

• generator (AbstractIntermediateGenerator) – Proposes molecules. It never decides feasibility and never prices anything.

• mapper (AbstractMapper or AbstractScorer) – The user’s own, unchanged. Every proposed sub-edge is put through build_candidate() with them and with the run’s SoftcorePolicy.

• scorer (AbstractMapper or AbstractScorer) – The user’s own, unchanged. Every proposed sub-edge is put through build_candidate() with them and with the run’s SoftcorePolicy.

• mapping_options (MappingOptions)

• network_options (NetworkOptions) – Supplies intermediates, the bans, the forced pairs, and the n_edges headroom.

Returns:

With the input mapping and pool returned unchanged when intermediates.mode == "off".

Return type:

AugmentationResult

Notes

The generator proposes; the existing feasibility machinery is the sole judge. Nothing here fabricates a Transformation the way make_cbfe_transformation() legitimately does, and the difference is not stylistic: a counterpoised edge has no geometry to check, while an intermediate edge is nothing but geometry. A badly posed molecule has to come back as an ordinary core_geometry_mismatch, which is exactly what routing through build_candidate() makes happen.

A proposal is accepted or dropped whole. If the surviving feasible sub-edges do not connect the two ends of the gap, the molecules go too – there is no such thing as a partially useful intermediate, and an orphan synthetic vertex would be a ligand nobody can compute a free energy for.

The edge budget is spent, not inflated. See intermediate_headroom(). Running out is recorded on unmet_constraints rather than raised, because it leaves a perfectly valid network.

This stage is deliberately serial. Posing consumes a fixed seed per molecule and naming is content-addressed, but the order in which gaps consume the shared budget is not commutative – so the pool is evaluated under jobs while the augmentation over it is not, and the output is identical at any jobs.

rbfenetmap.core.pipeline.build_candidate(source, target, mapper, scorer, mapping_options, network_options)

Map, repair, and score a single candidate pair.

Parameters:

• source (Ligand)

• target (Ligand)

• mapper (AbstractMapper)

• scorer (AbstractScorer)

• mapping_options (MappingOptions)

• network_options (NetworkOptions) – Contributes the SoftcorePolicy.

Returns:

Always a transformation, never an exception. A pair that cannot be mapped or repaired comes back marked infeasible with the reason attached, so it stays visible in the audit trail and can explain a later disconnection.

Return type:

Transformation

rbfenetmap.core.pipeline.build_network(ligands, *, mapper='mcss-e2', scorer='linear', planner='mst', mapping_options=None, network_options=None, progress_callback=None)

Plan a perturbation network over ligands.

Parameters:

• ligands (Sequence[Ligand] or Mapping[str, Ligand]) – The vertices. Names must be unique.

• mapper (AbstractMapper or AbstractScorer or AbstractNetworkPlanner or str) – Plugin instances, or names to look up in the built-in registries. Under network_options.cbfe_mode == "all" the mapper and scorer are never used – a counterpoised edge has no core to map and a closed-form cost – and a mapper name is not even resolved, so an unavailable optional mapper is not an error there.

• scorer (AbstractMapper or AbstractScorer or AbstractNetworkPlanner or str) – Plugin instances, or names to look up in the built-in registries. Under network_options.cbfe_mode == "all" the mapper and scorer are never used – a counterpoised edge has no core to map and a closed-form cost – and a mapper name is not even resolved, so an unavailable optional mapper is not an error there.

• planner (AbstractMapper or AbstractScorer or AbstractNetworkPlanner or str) – Plugin instances, or names to look up in the built-in registries. Under network_options.cbfe_mode == "all" the mapper and scorer are never used – a counterpoised edge has no core to map and a closed-form cost – and a mapper name is not even resolved, so an unavailable optional mapper is not an error there.

• mapping_options (MappingOptions, optional)

• network_options (NetworkOptions, optional)

• progress_callback (callable, optional) –

  Called with the number of candidate pairs just finished, repeatedly, as the mapping stage proceeds. Divide by the pair count to get a fraction.

  For a caller with a progress bar of its own, and the reason this exists rather than show_progress: that knob writes to stderr, which suits a terminal and suits nothing else. Mapping is where the time goes – quadratic in the ligand count, and over a thousand pairs by fifty ligands – so an embedding program that cannot report on it has nothing to show for minutes at a stretch.

  Two limits worth knowing. Under pair_evaluation="adaptive" the loop stops as soon as the targets are met, so the increments sum to at most the pair count rather than exactly it. And the sub-edges of a generated intermediate are not counted: they are a small, bounded stage that runs after this one.

Returns:

With edges selected and candidates holding everything scored.

Return type:

Network

Raises:

• ValueError – If fewer than two ligands are supplied, or two share a name.

• rbfenetmap.core.exceptions.NetworkPlanError – If the user’s constraints cannot be satisfied.

Examples

>>> network = build_network(ligands, mapper="cartograph")
>>> len(network.edges)
11

rbfenetmap.core.pipeline.evaluate_pairs(ligands, pairs, mapper, scorer, mapping_options, network_options, *, progress_callback=None)

Map, repair, and score every pair.

Parallelised over network_options.jobs threads, which keeps the immutable ligand and scorer mappings shared by reference; Python process pools cannot serialize them.

Threads win only on the native part of the work – FindMCS and the substructure search. Core selection, pruning, soft-core repair, and descriptors are pure Python and hold the GIL, so scaling is sublinear and flattens well below the core count.

Parameters:

• ligands (Mapping[str, Ligand])

• pairs (Sequence[tuple[str, str]])

• mapper (AbstractMapper)

• scorer (AbstractScorer)

• mapping_options (MappingOptions)

• network_options (NetworkOptions)

• progress_callback (Callable[[int], None] | None)

Return type:

list[Transformation]

rbfenetmap.core.pipeline.evaluate_pairs_adaptively(ligands, pairs, mapper, scorer, planner, mapping_options, network_options)

Evaluate pairs adaptively and plan over what was evaluated.

Parameters:

• ligands (Mapping[str, Ligand])

• pairs (Sequence[tuple[str, str]]) – Every pair the loop is allowed to reach for, in any order.

• mapper (AbstractMapper)

• scorer (AbstractScorer)

• planner (AbstractNetworkPlanner) – Drives the loop as well as producing the result: the probe plans are what tell the loop whether any target is still unmet.

• mapping_options (MappingOptions)

• network_options (NetworkOptions)

Return type:

Network

Notes

A thin wrapper over _adaptive_candidate_pool(), which does the work. It stays public and keeps returning a Network because that is its released signature; handing back a candidate list instead would be a real API break for anything that calls it directly.

Intermediate generation is deliberately not run here. It belongs after the loop settles and before planning, which is build_network()’s job – and a caller who reaches for this function directly is asking for adaptive evaluation of the pairs it was given, not for new vertices it never mentioned.

The final plan runs outside the loop’s warning suppression, so any genuinely unmet best-effort target is visible exactly once. For a required but impossible connection it raises with diagnostics, after every component-bridging possibility was attempted.

rbfenetmap.core.pipeline.feasible_graph(names, candidates)

Build the undirected graph of candidates that passed feasibility checks.

Parameters:

• names (Sequence[str]) – Every ligand, so an isolated one is a node of its own rather than absent.

• candidates (Sequence[Transformation]) – Scored candidates, feasible or not.

Return type:

networkx.Graph

Notes

Public because two stages now need exactly this graph and they must agree on it: the adaptive loop decides which pairs still cross a component boundary, and augment_with_intermediates() decides which gaps are worth offering to a generator. Two private copies that drifted apart would show up as a generator being offered a gap that no longer exists.

Candidate pair generation and prefiltering.

Decides which transformations are worth mapping and scoring at all. For anything past a couple of dozen ligands the all-pairs set is dominated by pairs no one would consider, and mapping is the expensive stage, so a cheap similarity prefilter pays for itself many times over.

The prefilter carries one obligation, discharged by reconnect_pairs(): it must not be allowed to disconnect the candidate pool behind the user’s back.

rbfenetmap.core.pairs.expand_pairs(names, strategy='all_unordered_pairs', *, hub=None, explicit=())

Enumerate candidate pairs under strategy.

A port of BuildEdges._expand_edges.

Parameters:

• names (Sequence[str]) – Ligand names, in input order.

• strategy (PairStrategy) – "all_unordered_pairs", "all_pairs", "star", "linear", or "explicit".

• hub (str, optional) – Required by "star".

• explicit (Sequence[str]) – "a~b" specifications, required by "explicit".

Returns:

Ordered pairs, deduplicated.

Return type:

list[tuple[str, str]]

Raises:

ValueError – For an unknown strategy, a missing or unknown hub, an unknown ligand in explicit, or a strategy that yields no pairs at all.

rbfenetmap.core.pairs.fingerprint_pair_similarities(ligands, pairs)

Return Morgan/Tanimoto similarity for each requested pair.

This is deliberately mapping-free and therefore cheap enough to rank an all-pairs pool before any MCS searches are launched.

Parameters:

• ligands (Mapping[str, Ligand])

• pairs (Sequence[tuple[str, str]])

Return type:

dict[tuple[str, str], float]

rbfenetmap.core.pairs.fingerprint_prefilter(ligands, pairs, *, top_k=8, min_similarity=0.4)

Keep only the pairs most likely to yield a usable transformation.

For each ligand, retains its top_k most similar partners, plus every pair above min_similarity. Similarity is Morgan/Tanimoto.

Parameters:

• ligands (Mapping[str, Ligand])

• pairs (Sequence[tuple[str, str]]) – Candidate pairs to filter.

• top_k (int, optional) – Neighbours retained per ligand.

• min_similarity (float, optional) – Tanimoto floor for unconditional retention.

Returns:

The surviving pairs, in the input order.

Return type:

list[tuple[str, str]]

Notes

Callers must follow this with reconnect_pairs(). Retaining each ligand’s nearest neighbours says nothing about whether the resulting graph is connected: a series containing two distinct chemical families will happily split into two components, each internally well connected, and the planner would then report a disconnection the user never asked for.

rbfenetmap.core.pairs.generate_candidate_pairs(ligands, options)

Produce the pairs to map and score, applying strategy, prefilter, and forcing.

Returns:

The candidate pairs and any pairs restored by the reconnection pass.

Return type:

tuple[list[tuple[str, str]], list[tuple[str, str]]]

Parameters:

• ligands (Mapping[str, Ligand])

• options (NetworkOptions)

rbfenetmap.core.pairs.reconnect_pairs(names, pairs, all_pairs, ligands=None)

Add pairs back until the candidate graph spans every ligand.

Parameters:

• names (Sequence[str]) – Every ligand name.

• pairs (Sequence[tuple[str, str]]) – The prefiltered pairs.

• all_pairs (Sequence[tuple[str, str]]) – The unfiltered pairs, from which bridges may be restored.

• ligands (Mapping[str, Ligand], optional) – Used to rank restoration candidates by similarity. Without it, restoration is arbitrary but deterministic.

Returns:

The reconnected pair list, and the pairs that had to be restored – reported so the user can see the prefilter was overridden rather than silently corrected.

Return type:

tuple[list[tuple[str, str]], list[tuple[str, str]]]

Notes

Mandatory after fingerprint_prefilter(). Prefiltering is an optimisation, and an optimisation that changes the answer – here, by making a connected network impossible – is a bug. Restoring the best available bridge keeps the prefilter honest: it may reorder the work, but it cannot remove an outcome.

Edge descriptors: the single place where scoring inputs are computed.

Every number a scorer sees originates here. That centralisation is deliberate. Scorers receive a plain Mapping[str, float] and never import RDKit, which means re-scoring a network under new weights requires no remapping, a scorer can be tested against hand-written dictionaries, and a third-party scorer cannot quietly grow a dependency on how the mapping was produced.

Descriptors are raw, unnormalised, and unweighted. Turning them into a cost – deciding that eight soft-core atoms is “one unit of bad” – is the scorer’s job, not this module’s.

rbfenetmap.core.descriptors.DESCRIPTOR_NAMES = ('n_core_heavy', 'n_softcore_heavy_1', 'n_softcore_heavy_2', 'n_softcore_max_heavy', 'softcore_asymmetry', 'n_heavy_1', 'n_heavy_2', 'heavy_atom_delta', 'charge_delta', 'n_rings_1', 'n_rings_2', 'ring_delta', 'n_ring_atoms_in_softcore', 'mcs_fraction', 'core_rmsd', 'core_max_pair_distance', 'rotatable_delta', 'logp_delta', 'n_fragments_before_1', 'n_fragments_before_2', 'n_demoted_atoms')

Every key compute_descriptors() produces. Scorers use this to validate weights.

rbfenetmap.core.descriptors.compute_descriptors(source, target, mapping, repair=None)

Compute every scoring descriptor for one candidate transformation.

Parameters:

• source (Ligand) – The two ligands.

• target (Ligand) – The two ligands.

• mapping (AtomMapping) – The mapping after repair, so soft-core sizes reflect what will actually run.

• repair (SoftcoreRepair, optional) – The repair record, contributing the fragmentation and demotion counts.

Returns:

Keyed by DESCRIPTOR_NAMES. Every value is a float, including the counts, so scorers never have to think about integer division.

Return type:

dict[str, float]

Notes

Soft-core sizes are counted in heavy atoms. Hydrogens follow their parent heavy atom into the soft-core automatically, so including them would mostly measure how hydrogenated a substituent is rather than how large the perturbation is – a -CH3 to -CF3 change would look like a shrinking soft-core.

core_rmsd is measured in place, without superposition. Ligands are normally supplied already posed in a common binding-site frame, and in that frame the deviation of mapped core atoms says whether the mapping pairs atoms that occupy the same part of the pocket. Superposing first would hide exactly that.

The maximum-common-substructure search, in one place.

FindMCS has a dozen switches and the answer changes with every one of them. Two callers configuring it separately will drift, and the drift is invisible: an aligner that maximises one substructure while the geometry gate measures a different one produces a run whose alignment report looks healthy next to edges rejected for core_geometry_mismatch, with nothing on screen to reconcile the two.

So the settings live here and both MCSSMapper and rbfenetmap.core.align call in. This module holds no policy of its own – every switch comes from the caller’s MappingOptions.

rbfenetmap.core.mcs.mcs_embeddings(mol_1, mol_2, pattern, options)

Return every embedding of pattern in each molecule.

Parameters:

• mol_1 (rdkit.Chem.Mol)

• mol_2 (rdkit.Chem.Mol)

• pattern (rdkit.Chem.Mol) – A query molecule, normally from mcs_query().

• options (MappingOptions) – Supplies max_matches.

Returns:

The embeddings in mol_1 and in mol_2. Either can be empty.

Return type:

tuple[tuple[tuple[int, …], …], tuple[tuple[int, …], …]]

Notes

uniquify=False is what makes the symmetry visible. With it set, RDKit collapses embeddings related by an automorphism of the query and returns one representative – so a para-substituted ring offers a single embedding and the flip that pairs atoms across the ring from one another can never be examined, let alone rejected.

rbfenetmap.core.mcs.mcs_query(mol_1, mol_2, options, *, match_elements=False)

Return the MCS of two molecules as a query molecule.

Parameters:

• mol_1 (rdkit.Chem.Mol)

• mol_2 (rdkit.Chem.Mol)

• options (MappingOptions) – Supplies timeout, the four FindMCS comparison flags, and the ring settings.

• match_elements (bool, optional) – Require paired atoms to be the same element. Default False, which is what a mapper wants. See the note below before switching it on or off.

Returns:

None when the molecules share no substructure, or when the SMARTS RDKit produces cannot be parsed back into a query. Both are ordinary outcomes for a caller deciding what to do about a pair, so neither raises here.

Return type:

rdkit.Chem.Mol or None

Notes

bondCompare=CompareAny is not laxness, and is not negotiable per caller. Ligands routinely arrive from force-field topologies whose bond orders are approximate – an Amber mol2 recording a carbonyl as C-O single is the everyday case, and rbfenetmap.io.loaders._prepare() deliberately preserves it rather than “correcting” it. A strict bond compare would silently shrink the common substructure on exactly those inputs.

The atom comparison is a different matter, and is the one setting callers legitimately disagree about. A mapper can afford CompareAny because everything downstream of it is a safety net: prune_core() demotes element mismatches out of the core, and the geometry gate catches whatever survives. A caller with no such net – alignment, where the correspondence is the answer and nothing revisits it – cannot. Left permissive there, FindMCS will pair a methoxy oxygen with a methyl carbon and an amide nitrogen with a hydrogen, which superposes eleven atoms to a convincing fraction of an angstrom while putting the scaffold several angstroms wrong. It is precisely the failure a good-looking RMSD hides.

Counterpoised binding free energy (CBFE) edges.

A CBFE edge is two absolute calculations run simultaneously in opposite directions: one ligand decoupling from the site as the other couples into it. It yields the same relative quantity an RBFE edge does, by a different route, and the difference that matters to a network planner is that neither molecule is morphed into the other. There is no common core to find, no soft-core region to repair, and therefore no way for a CBFE edge to be infeasible. It exists between every pair of ligands.

That makes it the natural repair for the failure mode this package hits on real ligand sets: an MCS-based candidate pool that comes back in several disconnected pieces, with no RBFE edge able to cross between them. CBFE edges can cross, and the cost model here is what keeps them from being used anywhere else.

Cost

cbfe_base_cost + cbfe_atom_weight * (n_heavy_1 + n_heavy_2), deliberately on the same scale as the scorer’s edge totals so the two kinds can be compared by the planner. The default base, 8.0, is the linear scorer’s charge-change ceiling: a CBFE edge is priced at roughly the most expensive thing that can happen to a still-feasible RBFE edge. Realistic totals land in ~[9, 13] against ~0.3 for a good RBFE edge and ~5-6 for a bad one, so CBFE never wins on price – only on being available where nothing else is. The per-atom term exists because a counterpoised calculation decouples both molecules in full, so its expense really does scale with how much there is to decouple.

Bridge selection

Choosing which CBFE edges join the subnetworks is a separate question from cost, and is answered by select_cbfe_bridges(): a maximum-merit spanning forest over the component quotient graph, where merit combines pairwise similarity with how well connected each endpoint is inside its own subnetwork. See bridge_rank_key() for why both terms are there.

The forest sweep itself is select_bridges(), which takes the partition as an argument instead of deriving it. Connected components are only one interesting partition of a ligand set – a chemical clustering is another, and joining those groups is the same problem with the same ranking. Keeping the partition a parameter is what lets clustered planning in rbfenetmap.core.clustering reuse this machinery rather than grow a second copy of it that would drift.

rbfenetmap.core.cbfe.BRIDGE_CENTRALITY_WEIGHT = 0.5

How much endpoint centrality is worth relative to similarity when ranking bridges. Both inputs are on [0, 1], so this reads as “being maximally well connected inside both subnetworks is worth as much as a 0.5 jump in Tanimoto similarity”. Kept a module constant rather than a user option: it is a tie-break heuristic, and the option surface is already three knobs wide. It is the obvious thing to promote if it ever needs tuning.

rbfenetmap.core.cbfe.CBFE_METHOD = 'cbfe'

Recorded as method on a CBFE edge. It names the reason the mapping is empty rather than leaving it as "unknown", which would be indistinguishable from a mapper that failed.

rbfenetmap.core.cbfe.CBFE_SCORER = 'cbfe'

Recorded as scorer. No scoring plugin is involved: the cost is closed-form, so attributing it to the configured scorer would be a lie that shows up in the exported JSON.

rbfenetmap.core.cbfe.bridge_rank_key(pair, *, similarity, centrality, cost)

Rank one candidate bridge. Lower sorts better.

Merit is a sum of two terms both on [0, 1]:

merit = similarity + BRIDGE_CENTRALITY_WEIGHT * mean(centrality of the endpoints)

Similarity is there because a bridge between chemically close ligands is the one most likely to give a trustworthy number even though it is being run as two absolute calculations. Centrality is there because a bridge landing on a hub propagates through the subnetwork – it participates in cycles and shares its endpoints with many RBFE edges – whereas one landing on a leaf leaves a dangling path that nothing checks.

A sum rather than a product so the trade-off stays legible: the weight states exactly what one term is worth in units of the other, which a product cannot do.

Parameters:

• pair (tuple[str, str]) – Sorted endpoint pair.

• similarity (Mapping[tuple[str, str], float]) – Fingerprint similarity per pair.

• centrality (Mapping[str, float]) – Per-ligand centrality from component_centrality().

• cost (Mapping[tuple[str, str], float]) – CBFE cost per pair, used only to break merit ties toward the cheaper edge.

Returns:

Ends in the pair itself, so the ordering is total and the selection is reproducible regardless of dictionary iteration order.

Return type:

tuple[float, float, tuple[str, str]]

rbfenetmap.core.cbfe.build_cbfe_pool(ligands, options, *, exclude=())

Return every eligible CBFE pair mapped to its cost.

Costs, not Transformation objects. The pool is quadratic in the ligand count and only a handful of its entries are ever selected, so materializing the rest would build tens of thousands of objects – each carrying full soft-core index tuples over every atom – to throw them away. The planner builds the transformation for a pair once it has decided to keep it.

Parameters:

• ligands (Mapping[str, Ligand])

• options (NetworkOptions) – banned_edges is honoured here; a ban applies to CBFE exactly as it does to RBFE, since it expresses “do not run this pair” rather than “do not map it”.

• exclude (Iterable[tuple[str, str]], optional) – Unordered pairs to omit – in practice the pairs that already have a feasible RBFE candidate. Filtering them here, once, is what keeps a pair from being offered as both kinds and tripping validate()’s duplicate check much later.

Returns:

Keyed by sorted endpoint pair.

Return type:

dict[tuple[str, str], float]

Notes

The pool spans all unordered pairs and deliberately ignores both pair_strategy and the fingerprint prefilter. Those exist to hold down the number of MCS searches, and a CBFE edge runs none. Applying the prefilter here would be worse than merely unnecessary: it drops dissimilar pairs, and dissimilar pairs are precisely the ones that need a bridge, because their similar neighbours already had RBFE edges available.

rbfenetmap.core.cbfe.cbfe_cost(source, target, options)

Return the cost of the CBFE edge between two ligands.

Parameters:

• source (Ligand)

• target (Ligand)

• options (NetworkOptions) – Supplies cbfe_base_cost and cbfe_atom_weight.

Returns:

Always finite – a CBFE edge cannot be infeasible.

Return type:

float

rbfenetmap.core.cbfe.component_centrality(graph, components)

Score each ligand by how well connected it is within its own subnetwork.

Degree divided by the largest degree in the same component, so the result is on [0, 1] and comparable across components of very different sizes.

Max-normalising rather than using networkx.degree_centrality() (degree over size - 1) matters for the singleton case. A ligand nothing could be mapped to forms a component of one with degree 0, and degree_centrality would score it 0.0 – the worst possible entry point. But a singleton has exactly one way into the network, and ranking its only option last is backwards. Here it scores 1.0.

Parameters:

• graph (networkx.Graph) – The feasible RBFE graph, before any bridges are added.

• components (Sequence[set[str]]) – Its connected components.

Return type:

dict[str, float]

rbfenetmap.core.cbfe.make_cbfe_transformation(source, target, options)

Build the CBFE edge between two ligands.

The mapping is the empty-core partition: nothing in common, every atom on both sides soft-core. That is not a placeholder standing in for a mapping nobody computed – it is a literal description of the experiment, in which both molecules are decoupled in full. It also satisfies every AtomMapping invariant, so the edge is a first-class transformation rather than a special case downstream.

Parameters:

• source (Ligand)

• target (Ligand)

• options (NetworkOptions)

Returns:

Always feasible, with CBFE as its kind.

Return type:

Transformation

rbfenetmap.core.cbfe.select_bridges(partition, ligands, pool, *, graph=None, n_per_pair=1)

Choose the edges that join the groups of partition into one network.

A maximum-merit spanning forest over the quotient graph of the partition, by the same union-find sweep the planner’s Kruskal pass uses. With g groups exactly g - 1 group pairs must be joined, and picking which g - 1 is precisely a spanning selection – ranking each group pair separately would produce C(g, 2) winners and still leave the same problem to solve.

The partition is a parameter rather than something derived here, and that is what makes the function reusable. select_cbfe_bridges() passes the connected components of the feasible RBFE graph, which is the “these pieces cannot reach each other” case. Clustered planning passes the clusterer’s partition, which is the “these pieces can reach each other but the budget is better spent inside them” case. Both want the same thing – the most trustworthy few edges crossing a boundary – and neither wants to reimplement the ranking.

Parameters:

• partition (Mapping[str, int]) – Group index per ligand name. Pairs whose endpoints are missing from it are ignored.

• ligands (Mapping[str, Ligand])

• pool (Mapping[tuple[str, str], float]) – Eligible pairs and their costs. CBFE costs from build_cbfe_pool() for the connectivity case, RBFE edge costs for the clustered case.

• graph (networkx.Graph, optional) – The graph the endpoints live in, used only for component_centrality(). Not mutated. Without it, centrality contributes nothing and bridges rank on similarity and cost alone – which is the right degradation, since “how well connected is this ligand” has no answer without a graph to ask it of.

• n_per_pair (int, optional) – Edges to take across each joined group pair. 1 gives the minimal spanning selection. 2 puts the crossing itself on a cycle: two edges between the same two groups, plus the paths inside each group, form a loop through both of them – which applies the every-edge-in-a-cycle invariant precisely to the edges that most need checking.

Returns:

Sorted endpoint pairs, in selection order. Shorter than n_per_pair * (g - 1) whenever the pool cannot supply that many crossings – every remaining pair having been banned, or simply absent. The caller reports that; this function does not raise, because the planner has a much better diagnostic to hand than anything available here.

Return type:

list[tuple[str, str]]

Notes

Centrality is computed once, on the graph as given, and not refreshed as bridges are added. Refreshing would make the outcome depend on selection order and would also contradict the intent: the question is how well connected a ligand is inside the group it came from, not how connected it became by being chosen.

rbfenetmap.core.cbfe.select_cbfe_bridges(graph, ligands, pool)

Choose the CBFE edges that join graph into one component.

A thin caller of select_bridges() over the partition induced by graph’s own connected components: with c components, exactly c - 1 bridges. One crossing per joined pair, because a second CBFE edge between the same two components would buy a cycle at the price of two more absolute calculations – and the mode that pays for counterpoised cycle coverage is cbfe_mode="cycles", which the planner applies afterwards on the merits of each ligand rather than blindly per component pair.

Parameters:

• graph (networkx.Graph) – The feasible RBFE graph. Not mutated.

• ligands (Mapping[str, Ligand])

• pool (Mapping[tuple[str, str], float]) – Eligible CBFE pairs and their costs, from build_cbfe_pool().

Returns:

Sorted endpoint pairs, in selection order.

Return type:

list[tuple[str, str]]

Diagnostics and cost

Network-level metrics over an already-planned network.

Everything here is read-only and pure: a Network goes in, numbers come out, and nothing is ever selected, rejected, or re-planned. That is what makes it safe to call from a report renderer and from the CLI without either of them being able to change what the other sees.

The distinction from rbfenet inspect is the unit of analysis. inspect answers questions about one edge – its mapping, its soft-core, why it was rejected. diagnose answers questions about the network: how long the longest comparison path is, how much of it survives a failed edge, whether the edge budget is anywhere near the statistical floor. A reviewer’s second question about a planned network, after “why isn’t X connected to Y”, is always one of these, and until now the HTML report was a picture with no numbers beside it.

The metric set follows Konnektor’s, with two deliberate departures.

Every seed is mandatory. failure_robustness() is the only Monte-Carlo function in the package, and tests/test_softcore.py already asserts that planning is deterministic. A defaulted seed is a defaulted seed until someone forgets it, so it is a required argument: an unseeded run is not possible to write by accident.

The n ln n edge-budget floor is advice, not a warning. Pitman et al., JCIM 2023, 63, 1776-1793 derive k_min ~ n ln n edges, below which precision degrades worse as n grows – at 40 ligands that is 148 edges where the default edges_per_ligand=2 buys about 40. Routing that through warnings.warn() would fire it on essentially every run this package has ever planned, and a warning that always fires is a warning nobody reads. It belongs in a report the user asked for, which is edge_budget_advice().

class rbfenetmap.core.diagnostics.DegreeSummary(degrees, minimum, maximum, mean)

Bases: object

Per-ligand edge counts and their extremes.

Parameters:

• degrees (dict[str, int]) – Ligand name to number of selected edges touching it.

• minimum (int) – The extremes. 0 for a network with no ligands.

• maximum (int) – The extremes. 0 for a network with no ligands.

• mean (float) – Mean degree, which is 2 * n_edges / n_ligands.

property isolated: tuple[str, ...]

Ligands no selected edge touches, sorted.

The single most actionable line in a diagnostic report: a ligand at degree zero has no measured free energy at all, whatever the rest of the network looks like.

class rbfenetmap.core.diagnostics.EdgeBudgetAdvice(n_ligands, n_edges, recommended, shortfall)

Bases: object

How the planned edge count compares with the published precision floor.

Parameters:

• n_ligands (int) – What was planned.

• n_edges (int) – What was planned.

• recommended (int) – ceil(n * ln n), the floor Pitman 2023 derives.

• shortfall (int) – recommended - n_edges, clamped at zero.

Notes

Advisory, and deliberately not a warning – see this module’s docstring. The floor is also a floor for precision, not for correctness: a network below it is a perfectly valid network whose free energies simply carry more statistical uncertainty than a denser one over the same ligands would, and the gap widens as the series grows.

property message: str

One line stating the comparison, suitable for a report or the plan summary.

class rbfenetmap.core.diagnostics.FailureRobustness(connected_fraction, mean_ligands_retained, failure_rate, n_repeats, seed)

Bases: object

What survives when edges fail, measured by Monte-Carlo removal.

Parameters:

• connected_fraction (float) – Fraction of trials in which the surviving network still spans every ligand.

• mean_ligands_retained (float) – Mean size of the largest surviving connected component, in ligands. The natural companion to the fraction above: a network that stays connected 40% of the time but keeps 95% of its ligands the rest of the time is in a very different position from one that shatters.

• failure_rate (float) – The per-edge failure probability the trials used.

• n_repeats (int) – How many trials were run.

• seed (int) – The seed they were run with, carried so a reported figure can be reproduced without going back to the command line that produced it.

rbfenetmap.core.diagnostics.count_cycles(network, max_length=4)

Count the simple cycles of at most max_length ligands.

Parameters:

• network (Network)

• max_length (int, optional) – Longest cycle counted. The default of 4 follows cinnabar’s convention.

Return type:

int

Raises:

ValueError – If max_length is below three, which no cycle can be.

Notes

The bound is not a performance nicety, it is what makes the function terminate in useful time: the simple cycles of a dense graph are exponential in the node count, and a redundant network over fifty ligands has enough of them to hang a report. Short cycles are also the ones that matter – a cycle-closure residual over twenty edges localises nothing.

rbfenetmap.core.diagnostics.degree_summary(network)

Summarize how many edges each ligand carries.

Parameters:

network (Network)

Return type:

DegreeSummary

rbfenetmap.core.diagnostics.diameter(network)

Longest shortest path between two ligands, in edges.

Parameters:

network (Network)

Returns:

None when the network is disconnected, because the diameter is then infinite rather than large and reporting a number for it would be a lie. 0 for a single ligand.

Return type:

int or None

Notes

Computed with usebounds=True: the bounded form (the FastLomap optimisation, arXiv:2304.04713) replaces the all-pairs sweep with a handful of BFS runs, which is what keeps this affordable in a report over a hundred ligands.

rbfenetmap.core.diagnostics.edge_budget_advice(n_ligands, n_edges)

Compare an edge count against the n ln n precision floor.

Parameters:

• n_ligands (int)

• n_edges (int)

Return type:

EdgeBudgetAdvice

Notes

Takes two integers rather than a Network on purpose: the most useful moment to ask this is before planning, when the only thing that exists is a ligand count and a budget.

rbfenetmap.core.diagnostics.failure_robustness(network, *, failure_rate=0.05, n_repeats=100, seed)

Estimate how much of network survives independent edge failures.

An alchemical edge fails for reasons a planner cannot see – a sampling problem, a crashed run, a pose that turns out wrong. This asks what the network looks like afterwards: remove each edge independently with probability failure_rate, and see whether the rest still hangs together.

Parameters:

• network (Network)

• failure_rate (float, optional) – Independent per-edge failure probability, in [0, 1].

• n_repeats (int, optional) – Number of Monte-Carlo trials.

• seed (int) – Required, not optional. This is the only stochastic function in the package, and everything around it asserts determinism. A default here would be a default right up until someone left it off, and a diagnostic number that changes between two runs of the same command is worse than no number.

Return type:

FailureRobustness

Raises:

ValueError – If failure_rate is outside [0, 1] or n_repeats is not positive.

Notes

Edge failures are treated as independent, which is optimistic: in practice the ligand that breaks one edge tends to break its neighbours too, so read the result as an upper bound on robustness rather than an estimate of it.

rbfenetmap.core.diagnostics.network_cost(network)

Total cost of the selected edges, on the scorer’s scale.

Parameters:

network (Network)

Return type:

float

rbfenetmap.core.diagnostics.network_efficiency(network)

Mean cost per selected edge, on the scorer’s scale.

Parameters:

network (Network)

Returns:

0.0 for a network with no edges, rather than a division error: an empty network is a legitimate thing to hand a report renderer.

Return type:

float

Notes

Useful only between networks over the same ligands. Comparing it across ligand sets compares two arbitrary difficulty scales, and comparing it against network_cost() compares a mean against a sum – a denser network is expected to have the higher total and the lower mean at the same time.

rbfenetmap.core.diagnostics.summarize(network, *, seed=0, failure_rate=0.05, n_repeats=100, max_cycle_length=4)

Run every diagnostic over network and return the results together.

Parameters:

• network (Network)

• seed (int, optional) – Passed to failure_robustness(). Defaulted here and nowhere else: this is the presentation layer, where every call is one of many and the caller is asking for a report rather than for a number, so a stable default is what keeps two runs of rbfenet diagnose on the same file agreeing with each other.

• failure_rate (float, int, optional) – Passed to failure_robustness().

• n_repeats (float, int, optional) – Passed to failure_robustness().

• max_cycle_length (int, optional) – Passed to count_cycles().

Returns:

Keys n_ligands, n_edges, n_rbfe, n_cbfe, cost, efficiency, n_cycles, max_cycle_length, degrees (a DegreeSummary), diameter, robustness (a FailureRobustness), and budget (an EdgeBudgetAdvice).

Return type:

dict[str, Any]

rbfenetmap.core.diagnostics.summarize_json(network, **kwargs)

Return summarize() flattened into JSON-ready primitives.

Parameters:

• network (Network)

• **kwargs – Passed through to summarize().

Returns:

The same keys summarize() produces, except that degrees, robustness and budget – which are dataclasses there – arrive as the nested degree, robustness and edge_budget objects, and the cost summary from rbfenetmap.core.cost.network_cost_summary() is merged in as gpu_hours and price.

Return type:

dict

Notes

This is the shape rbfenet diagnose --format json emits, and it is a function rather than a literal inside that command because it is now emitted from two places. A second consumer building its own dictionary would be a second answer to “how did this network do”, differing from the first the moment either gained a field – and these numbers are meant to be compared across runs, which only works while every run reports them the same way.

Turning a planner cost into machine time and money, for reporting only.

The scorer’s edge total is a difficulty number on an arbitrary scale. It orders edges, which is all selection needs, and it is meaningless to anyone deciding whether a network fits in the cluster allocation they have this month. This module supplies the other translation: how many GPU-hours the planned network will take, and what that costs.

Nothing here feeds selection. cbfe_base_cost is untouched, no planner reads a CostModel, and switching --cost-units cannot move a single edge. That is a deliberate boundary rather than an omission: a wall-clock price is a constant per edge kind and would make CBFE edges uniformly unaffordable, which is precisely the confusion the “eligibility is a gate, not a price” rule in Network selection exists to prevent. Feeding cost into selection is a later phase’s job, where variance-weighted edge costs give it a principled basis.

The defaults come from measurements rather than estimates: Tsai et al., JCIM 2026, 66, 1626-1636 (10.1021/acs.jcim.5c02204) Table 1 reports 3.97 GPU-hours for an RBFE edge at 12 lambda windows and 12.81 for a counterpoised one at 25, a ratio of 3.2. The price per GPU-hour follows Pitman et al., JCIM 2023, 63, 1776-1793, at $0.40.

rbfenetmap.core.cost.COST_UNITS: tuple[Literal['score', 'gpu_hours'], ...] = ('score', 'gpu_hours')

The units a cost report can be expressed in. "score" is the scorer’s own arbitrary difficulty scale, which orders edges; "gpu_hours" is machine time, which orders budgets. Neither is a conversion of the other – they answer different questions.

class rbfenetmap.core.cost.CostModel(rbfe_gpu_hours=3.97, cbfe_multiplier=3.2, price_per_gpu_hour=0.4)

Bases: object

Per-edge machine cost, in GPU-hours and in currency.

Parameters:

• rbfe_gpu_hours (float) – GPU-hours for one relative edge. The default, 3.97, is Tsai 2026’s measurement at 12 lambda windows.

• cbfe_multiplier (float) –

  What a counterpoised edge costs relative to a relative one. The default, 3.2, is the same paper’s 12.81 / 3.97 at 25 lambda windows against 12.

  Expressed as a multiplier rather than as a second absolute figure on purpose: a user who changes their lambda schedule or their hardware moves both numbers together, and the ratio is the part that survives. Overriding rbfe_gpu_hours alone then keeps the CBFE figure honest.

• price_per_gpu_hour (float) – Currency per GPU-hour. The default, 0.40, is the rate Pitman 2023 costs its networks at.

Raises:

ValueError – If any field is negative, or cbfe_multiplier is below one. A counterpoised edge is two absolute calculations; it cannot be cheaper than the relative edge it replaces, and a multiplier below one would silently invert every report.

property cbfe_gpu_hours: float

GPU-hours for one counterpoised edge.

edge_gpu_hours(edge)

GPU-hours for one edge, chosen by its EdgeKind.

Parameters:

edge (Transformation)

Return type:

float

network_gpu_hours(network)

GPU-hours for every selected edge of network.

Parameters:

network (Network)

Return type:

float

network_price(network)

Currency cost of running network, at price_per_gpu_hour.

Parameters:

network (Network)

Return type:

float

rbfenetmap.core.cost.network_cost_summary(network, *, model=None)

Summarize what network costs, in both the scorer’s units and machine time.

Parameters:

• network (Network)

• model (CostModel, optional) – Defaults to the published figures.

Returns:

score (summed edge totals), gpu_hours, and price. All three are reported together so a caller choosing between --cost-units never has to run this twice, and so a reader comparing two networks sees the difficulty figure and the wall-clock figure move independently – which they do, because difficulty varies per edge and machine time does not.

Return type:

dict[str, float]

Partitioning a ligand set into clusters, so the network can be planned per cluster.

Why a partition changes the edge budget

Pitman et al. (JCIM 2023, 63, 1776-1793) put the precision floor of an RBFE network at k_min ~ n ln n edges: below it, precision degrades worse as the set grows. That floor is superlinear, and superlinear costs are exactly the ones a partition beats. For clusters of sizes n_1 ... n_d summing to n:

sum_i n_i ln n_i  <  n ln n

with equality only for a single cluster. Planning each cluster to the floor and joining the clusters with a handful of bridges therefore buys the same per-cluster precision for n ln(n/d)-ish edges: 100 ligands in five balanced clusters need roughly 190 edges rather than 460. Even a badly imbalanced split saves 30-50%, because the term that dominates is the largest cluster and it is still smaller than the whole set.

What a cluster is not

Not a feasibility statement. Nothing here consults the soft-core budget, a mapping, or a rejection: clustering is a selection-level objective, in the sense the planning notes insist on, and it shapes which of the feasible edges are worth spending on rather than which edges exist. Two ligands in different clusters may well have a perfectly good RBFE mapping between them; the point is that paying for many such edges buys less than paying for edges inside a cluster, because the within-cluster edges are the ones a cycle can check.

The three clusterers

cluster_by_charge() is exact and free – net formal charge is a property of the molecule, not of a similarity threshold, and a charge-changing edge is the one alchemical transformation the package already penalises hardest. cluster_by_scaffold() groups on the Bemis-Murcko framework, which is what a medicinal chemist means by “series”. cluster_by_fingerprint() is the general fallback for a set with neither a clean charge split nor a shared framework.

Every clusterer returns {ligand name: cluster index}, and the indices are canonicalised by _label_groups() so that the same ligand set always yields the same numbering regardless of dictionary order or of how the underlying grouping key sorted.

rbfenetmap.core.clustering.CLUSTER_METHODS: tuple[str, ...] = ('none', 'charge', 'scaffold', 'fingerprint')

The clustering methods assign_clusters() understands, "none" included so the option surface has a single vocabulary. "none" is a real member rather than a sentinel the caller tests for separately, which keeps the planner’s dispatch honest.

rbfenetmap.core.clustering.DEFAULT_FINGERPRINT_CUTOFF = 0.6

Average-linkage cut, as a Tanimoto distance, used when the caller names neither a cluster count nor a cutoff. 0.6 is one minus 0.4, and 0.4 is already this package’s stated notion of “similar enough to be worth mapping” – it is the default prefilter_min_tanimoto. Reusing that number rather than inventing a second one means a user who has tuned the prefilter for their chemistry has tuned this in the same direction, and there is only one similarity threshold in the package to reason about.

rbfenetmap.core.clustering.assign_clusters(ligands, method, **kwargs)

Dispatch to a clusterer by name.

Parameters:

• ligands (Mapping[str, Ligand])

• method ({"none", "charge", "scaffold", "fingerprint"}) – "none" puts every ligand in cluster 0, which is the partition that makes every downstream clustering step a no-op. Returning that rather than raising means a caller can hand the option through unconditionally.

• **kwargs – Forwarded to the chosen clusterer. Only "fingerprint" accepts any.

Returns:

Cluster index per ligand name.

Return type:

dict[str, int]

Raises:

ValueError – If method is unknown, or a keyword is passed to a clusterer that takes none.

rbfenetmap.core.clustering.cluster_by_charge(ligands)

Cluster on net formal charge.

The one clusterer with no threshold in it. Charge is a property of the molecule rather than of a similarity measure, and a net charge change is the transformation this package already treats as the most expensive thing that can happen to a still-feasible edge – so grouping by charge concentrates the edge budget on the edges whose free energies are most trustworthy, and forces the charge-crossing edges to be few and deliberately chosen.

Parameters:

ligands (Mapping[str, Ligand])

Returns:

Cluster index per ligand name. The clusters are exactly the charge classes.

Return type:

dict[str, int]

rbfenetmap.core.clustering.cluster_by_fingerprint(ligands, *, n_clusters=None, cutoff=None)

Cluster by average-linkage hierarchical clustering on Tanimoto distance.

Uses scipy.cluster.hierarchy.linkage() on the condensed 1 - Tanimoto matrix, then fcluster() to cut it. scipy rather than sklearn: scipy is already a core dependency of this package and sklearn is not, and average linkage on Tanimoto distance is the standard chemoinformatics choice – the density methods the neighbouring tools reach for (HDBSCAN in Konnektor, DBSCAN in HiMap) would buy a noise label this package has no use for, since every ligand must land in some cluster to be planned at all.

Average linkage specifically, rather than single or complete: single linkage chains a congeneric series into one cluster through a string of near-duplicates, and complete linkage refuses to admit a ligand that is far from any member, which splits a real series on its most-substituted compound.

Parameters:

• ligands (Mapping[str, Ligand])

• n_clusters (int, optional) – Cut the dendrogram to exactly this many clusters. Takes precedence over cutoff.

• cutoff (float, optional) – Cut at this Tanimoto distance; ligands merged below it share a cluster. Defaults to DEFAULT_FINGERPRINT_CUTOFF when n_clusters is also unset.

Returns:

Cluster index per ligand name.

Return type:

dict[str, int]

Raises:

ValueError – If n_clusters is not positive, or cutoff is outside [0, 1].

rbfenetmap.core.clustering.cluster_by_scaffold(ligands)

Cluster on the Bemis-Murcko scaffold.

RDKit’s MurckoScaffold.GetScaffoldForMol strips every side chain, leaving the ring systems and the linkers between them – which is close to what a medicinal chemist means by “series”. Ligands sharing a framework are the ones an MCS search relates cleanly, so a scaffold cluster is usually also a well-connected RBFE subnetwork, and the edges the partition sacrifices are the scaffold hops that were the least trustworthy anyway.

An acyclic ligand has an empty scaffold. Those are grouped together under that empty scaffold rather than each becoming a singleton: “has no ring system” is a genuine shared property here, and scattering them into singletons would produce a partition with as many bridges as ligands.

Parameters:

ligands (Mapping[str, Ligand])

Returns:

Cluster index per ligand name.

Return type:

dict[str, int]

rbfenetmap.core.clustering.cluster_edge_budget(partition)

Report the precision floor of partition against the unclustered floor.

Both numbers are the Pitman n ln n floor: one evaluated over the whole set, one summed over the clusters. The ratio is the saving clustering buys, and reporting it rather than asserting it is deliberate – a partition into one cluster, or into n singletons, saves nothing, and a user who has picked a clusterer that does that should be able to see it.

Parameters:

partition (Mapping[str, int]) – Cluster index per ligand name.

Returns:

n_ligands, n_clusters, sizes, clustered_floor (sum_i n_i ln n_i), unclustered_floor (n ln n), and saving, the fraction of the unclustered floor the partition avoids. saving is 0.0 for a set too small to have a floor at all.

Return type:

dict[str, Any]

rbfenetmap.core.clustering.cluster_sizes(partition)

Return the cluster sizes of partition, largest first.

Parameters:

partition (Mapping[str, int])

Return type:

list[int]

Optimal design

Optimal experimental design over a perturbation network.

The one fact this whole module rests on: the Fisher information matrix of a network of relative free energy measurements is the weighted graph Laplacian.

An RBFE edge measures \(\Delta G_j - \Delta G_i\) with variance \(\sigma_{ij}^2\). Stacking those observations and forming \(X^T \Sigma^{-1} X\) gives

\[F_{ii} = \sum_{k \ne i} \sigma_{ik}^{-2}, \qquad F_{ij} = -\sigma_{ij}^{-2}\]

which is exactly the Laplacian of the graph with edge weights \(\sigma_{ij}^{-2}\). DiffNet, HiMap, Yang’s MLE and cinnabar’s network analysis are all this same object, so one implementation serves selection, allocation, and analysis rather than three.

The covariance of the estimated free energies is \(C = F^{-1}\), and the two classical criteria are

• A-optimal – minimise \(\operatorname{tr} C\), the total variance of the estimates;

• D-optimal – minimise \(\ln \det C\), the volume of the joint confidence ellipsoid.

Singularity, and why it is not a problem

\(F\) is singular for an RBFE-only network: no relative measurement can pin the absolute offset, so the all-ones vector is always in the null space. NetBFE regularises by restraining the mean,

\[F^*(\omega) = F + \omega m^{-2} \mathbb{1}\mathbb{1}^T ,\]

and takes \(\omega \to \infty\) through the bordered system. Because \(F \mathbb{1} = 0\), the two terms commute and

\[\left(F + \tfrac{\omega}{m} P\right)^{-1} = F^{+} + \tfrac{m}{\omega} P \;\xrightarrow[\omega \to \infty]{}\; F^{+},\]

where \(P = \mathbb{1}\mathbb{1}^T / m\) projects onto the null space. So the limit is just the Moore-Penrose pseudo-inverse, and – the part that makes the problem well posed – the optimal design does not depend on \(\omega\) at all. The regulariser only fixes the unidentifiable offset; it never trades against the criterion.

D-optimality is spanning trees

The pseudo-determinant of a graph Laplacian is \(m\) times the number of spanning trees (Kirchhoff’s matrix-tree theorem, weighted). Minimising \(\ln \det C\) therefore maximises the weighted spanning-tree count, which is why D-optimal designs come out markedly more cyclic than A-optimal ones at the same edge count. That is Pitman’s argument for preferring D-optimality whenever a cycle-closure correction will be applied downstream.

What this buys, and what it does not

Precision, and only precision. Over five TYK2 iterations NetBFE’s \(\operatorname{tr} C\) fell monotonically from 1.08 to 0.78 while the RMSE against experiment rose from 0.84 to 0.91. An optimal design makes the numbers it produces more reproducible; it says nothing about whether the force field they come from is right.

rbfenetmap.core.design.DESIGN_CRITERIA: tuple[str, ...] = ('a_optimal', 'd_optimal')

The criteria criterion_value() understands, lowest-is-best in both cases.

rbfenetmap.core.design.a_optimal_criterion(fisher)

Return \(\operatorname{tr} C\), or inf for a disconnected network.

inf rather than a large number: a ligand no edge reaches genuinely has unbounded variance, and any finite stand-in would let a disconnected design win a comparison against a connected one that happened to be poor.

Parameters:

fisher (ndarray)

Return type:

float

rbfenetmap.core.design.a_optimal_gradient(nodes, edges, fisher)

Return \(-\partial \operatorname{tr} C / \partial w_e\) for every edge.

Returns:

\(u_e^T C^2 u_e\), one per edge, where \(u_e = e_i - e_j\).

Return type:

numpy.ndarray

Parameters:

• nodes (Sequence[str])

• edges (Sequence[tuple[str, str]])

• fisher (ndarray)

Notes

From \(\mathrm{d}C = -C \,\mathrm{d}F\, C\) and \(\partial F/\partial w_e = u_e u_e^T\), so \(\partial \operatorname{tr} C/\partial w_e = -\operatorname{tr}(C u_e u_e^T C) = -u_e^T C^2 u_e\). Positive by construction, so more weight on any edge always lowers the total variance – what the design chooses is where the next unit of weight helps most.

rbfenetmap.core.design.allocate_effort(nodes, edges, sigmas, *, total=1.0, iterations=200, tolerance=1e-09)

Distribute a fixed sampling budget over the edges, A-optimally.

Parameters:

• nodes (Sequence[str])

• edges (Sequence[tuple[str, str]])

• sigmas (Sequence[float]) – Predicted standard deviation of each edge at unit effort.

• total (float, optional) – Budget to divide, in whatever unit the caller wants back (nanoseconds, GPU-hours, or 1.0 for fractions).

• iterations (int, optional) – Cap on the multiplicative iteration.

• tolerance (float, optional) – Stop once the largest relative change in an allocation falls below this.

Returns:

Effort per edge, summing to total. Keys are the edges as given.

Return type:

dict[tuple[str, str], float]

Raises:

ValueError – If total is not positive, or the network is disconnected – an unreachable ligand’s variance is unbounded, so no finite budget makes the criterion finite and there is nothing to optimise.

Notes

The model is the usual one: variance falls as \(1/t\), so an edge given effort \(t_e\) contributes Fisher weight \(w_e = t_e / v_e\) with \(v_e = \sigma_e^2\) its unit-effort variance. Minimising \(\operatorname{tr} C\) over \(\sum t_e = T\) is convex, and this uses Titterington’s multiplicative algorithm – at each step scale every allocation by the square root of its directional derivative \(g_e = u_e^T C^2 u_e / v_e\) (see a_optimal_gradient()) and renormalise. The exponent of one half is the standard choice for A-optimality, and it is what makes the update scale-invariant: \(g_e\) is homogeneous of degree \(-2\) in the effort, so \(t_e \sqrt{g_e}\) is of degree zero and the iteration cannot drift with the budget.

The stationary point satisfies \(g_e = \operatorname{tr}(C) / T\) for every edge, which is worth stating in words: at the optimum every edge returns the same variance reduction per nanosecond. Any other allocation has an edge worth moving time to.

Published payoff is roughly a twofold variance reduction at equal total cost. This is the static first pass – it predicts variances from the descriptors rather than measuring them, so it cannot refit against what the simulations actually produced.

rbfenetmap.core.design.criterion_value(fisher, criterion)

Evaluate criterion on fisher. Lower is better for both.

Raises:

ValueError – If criterion is not in DESIGN_CRITERIA.

Parameters:

• fisher (ndarray)

• criterion (str)

Return type:

float

rbfenetmap.core.design.covariance(fisher)

Return \(C = F^{+}\), the mean-restrained covariance of the estimates.

Notes

The pseudo-inverse is the \(\omega \to \infty\) limit of the bordered system, as the module docstring derives – there is no approximation here and no regularisation parameter to choose.

Parameters:

fisher (ndarray)

Return type:

ndarray

rbfenetmap.core.design.d_optimal_criterion(fisher)

Return \(\ln \det C\) on the identifiable subspace, or inf if disconnected.

Computed as \(-\sum \ln \lambda_i\) over the non-null eigenvalues – the log of the pseudo-determinant, negated. Summing logs rather than taking a determinant is not a micro-optimisation: for a 40-ligand network the product of the eigenvalues underflows long before the log of it would.

Parameters:

fisher (ndarray)

Return type:

float

rbfenetmap.core.design.effective_resistances(nodes, edges, fisher)

Return each edge’s effective resistance under the current design.

Parameters:

• nodes (Sequence[str])

• edges (Sequence[tuple[str, str]])

• fisher (numpy.ndarray) – The Fisher matrix the resistances are measured against.

Returns:

R_e = C_ii + C_jj - 2 C_ij, one per edge.

Return type:

numpy.ndarray

Notes

This is the D-optimal gradient, up to sign: \(\partial \ln\det C / \partial w_e = -R_e\). It also has the direct reading that gives it its name – an edge with a large effective resistance sits where the network is weakest, carrying information no parallel path supplies – and it sums to \(n - 1\) over any design, which is Foster’s theorem and a cheap check that the matrix was built right.

Not the A-optimal gradient, which is the same quadratic form taken against \(C^2\) rather than \(C\); see a_optimal_gradient(). Confusing the two is easy and self-consistent enough to survive a smoke test, because both are positive and both are largest on the same sort of edge – but only one of them is homogeneous of the right degree, so an allocation built on the wrong one drifts instead of converging.

rbfenetmap.core.design.fisher_information(nodes, edges, sigmas)

Return the Fisher information matrix of a network, i.e. its weighted Laplacian.

Parameters:

• nodes (Sequence[str]) – Ligand names, in the order the matrix rows and columns take.

• edges (Sequence[tuple[str, str]]) – Endpoint pairs. Repeated pairs accumulate, which is correct: two independent measurements of the same transformation add their information.

• sigmas (Sequence[float]) – Predicted standard deviation per edge, in the same order. Must be positive.

Returns:

Symmetric (len(nodes), len(nodes)) array.

Return type:

numpy.ndarray

Raises:

ValueError – If the lengths disagree, an endpoint is not in nodes, or a sigma is non-positive. A zero sigma is a claim of infinite information from one edge, which would make every criterion below read as zero regardless of the rest of the network – silently the best possible design.

rbfenetmap.core.design.summarize(nodes, edges, sigmas)

Return both criteria for a design, for reporting.

Returns:

{"a_optimal": tr(C), "d_optimal": ln det(C)}.

Return type:

Mapping[str, float]

Parameters:

• nodes (Sequence[str])

• edges (Sequence[tuple[str, str]])

• sigmas (Sequence[float])

Graph-wide core consistency: one core per ligand, not one per edge.

The default, consistency="pairwise", maps each edge independently. A ligand sitting on three edges therefore holds three different common cores, one per partner, and nothing requires them to agree. That is the right default – each edge gets the largest core its own pair supports, which is the cheapest transformation for that pair – but it means the series has no single shared scaffold. Whether an atom is “in the core” is a question that can only be answered per edge.

consistency="graph" answers it per ligand. Each ligand keeps the intersection of the cores it holds across all of its selected RBFE edges; everything else is demoted to soft-core and the repair is re-run on what remains. The result is a genuine common core for the whole (connected) network rather than a merely pairwise-compatible one, which is what makes a group of ligands share a scaffold in the sense a per-cluster Amber setup wants.

Why it iterates

Intersecting is not a single pass. Demoting an atom on one side drops its partner on the other, which shrinks that ligand’s core, which changes its intersection; and the soft-core repair may demote further atoms still to keep the soft-core in one connected piece. Cores only ever shrink, so the iteration is monotone in a finite set and terminates – it is run to a fixed point rather than applied once.

What it does not do

It does not re-select. Shrinking a core makes an edge dearer, and in principle a different network would be optimal under the reduced cores; recomputing selection here would mean re-mapping the whole candidate pool under a constraint that depends on which edges were selected, which is circular. This is a post-selection refinement of the mappings of the edges that were chosen, and the costs it recomputes are reported honestly rather than fed back into selection.

CBFE edges are ignored. A counterpoised edge has no common core by construction, so reading one as “this ligand’s core is empty here” would erase the core of every ligand a bridge touches – an artefact of the bridge, not a statement about the scaffold.

rbfenetmap.core.consistency.apply_graph_consistency(network, *, scorer='linear', policy=None, scope='graph')

Rewrite network’s selected edges onto one core per ligand.

Parameters:

• network (Network) – A planned network. Its candidate pool, planner, and edge selection are unchanged; only the mappings, repairs, and costs of the selected RBFE edges are rewritten.

• scorer (AbstractScorer or str, optional) – Used to re-cost the reduced edges. A string is looked up in the scorer registry. Pass the same scorer the network was planned with – costs computed by two different scorers are not comparable, and the returned network holds a mixture of neither.

• policy (SoftcorePolicy, optional) – Feasibility policy for the re-run repair. Defaults to the network’s own.

• scope (str, optional) – How widely one core is required per ligand: "graph" (default) over all of its selected RBFE edges, or "component" only within its connected component of the RBFE-only selected subgraph. See consistency_groups().

Returns:

With graph-consistent mappings and recomputed costs.

Return type:

Network

Raises:

rbfenetmap.core.exceptions.NetworkPlanError – If any selected edge becomes infeasible under the reduced core.

Notes

The failure mode is a hard error rather than a per-edge rejection, and that is the one design decision here worth arguing about. Elsewhere in the package an infeasible edge is recorded and kept, because it is a candidate nobody has to run. These are selected edges: a network handed back containing an edge marked infeasible is a network that cannot be run, and quietly reverting the offending edges to their pairwise cores would hand back something that is not graph-consistent while claiming to be – the exact failure --consistency graph was reported for in the first place.

A raise here also carries real information: it means these ligands do not share a core large enough to run on, which is a fact about the series, and the message names the edges and reasons so the user can loosen a threshold, drop a ligand, or plan the subsets separately.

rbfenetmap.core.consistency.consistency_groups(network, scope)

Partition the ligands into the sets that must each share one core.

Parameters:

• network (Network)

• scope (str) – A member of CONSISTENCY_SCOPES other than "pairwise", which asks for no consistency at all and never reaches here.

Returns:

Ligand name to group index, or None for "graph" – one group covering everything, which is the same thing said without building a dictionary the callers would then have to check every lookup against.

Return type:

Mapping[str, int] or None

Notes

"component" groups by the connected components of the RBFE-only selected subgraph, and each of those three words is load-bearing. Selected, because the candidate pool holds edges the planner rejected. RBFE-only, because a CBFE bridge joins components without relating any atoms, so counting it would merge two groups that share no scaffold and hand the intersection a pair of ligands with nothing in common. And components rather than clusters, because that partition is already implied by the pool: a set whose scaffolds cannot be mapped to each other is several components.

rbfenetmap.core.consistency.graph_consistent_cores(network, *, policy=None, scope='graph')

Return the atoms each ligand keeps in the graph-consistent core.

Parameters:

• network (Network) – A planned network. Only its selected RBFE edges are read.

• policy (SoftcorePolicy, optional) – Used for the repair run between intersection passes. Defaults to the network’s own policy, or to library defaults if the network carries no options.

• scope (str, optional) – "graph" (default) or "component". See consistency_groups().

Returns:

Atom indices, per ligand name. A ligand with no selected RBFE edge is absent: it is under no constraint, because consistency is a statement about atoms shared across edges and it has none.

Return type:

dict[str, frozenset[int]]

Warning

This is not a feasibility statement. Unlike apply_graph_consistency(), it never raises, and the cores it reports can be ones no runnable network could use: when the repair rejects an edge it returns that edge’s mapping unchanged, so a core that survived only because its repair failed is reported here exactly like one that survived on merit. Read it as “the atoms these ligands have in common”, not as “the core your edges will run with”. Call apply_graph_consistency() for the latter – it is the one that checks.

Notes

Exposed separately from apply_graph_consistency() because the surviving core is the answer to “do these ligands share a scaffold at all, and how big is it?”, which is worth asking without rewriting a network to find out.

rbfenetmap.core.consistency.maybe_apply_graph_consistency(network, options, *, scorer='linear')

Apply apply_graph_consistency() when options asks for it.

A single gate, called on every path out of the pipeline, so that a consistency scope cannot be honoured on one route and silently dropped on another – which is the shape of the bug the option had before it did anything at all.

"pairwise" is the only scope that does nothing, and it is tested for by name rather than by position in the ladder: a scope added later should have to state that it is a no-op, not inherit it from being listed first.

Parameters:

• network (Network)

• options (NetworkOptions)

• scorer (AbstractScorer | str)

Return type:

Network

Surgery and diagnostics

Post-planning network surgery: add a ligand, drop an edge, join two networks.

Nobody plans once. A campaign gains compounds in batches, loses edges when a run fails to converge, and grows by joining a new series onto one that is already running. Re-planning from scratch each time is the wrong answer to all three: it discards the mappings that were already computed, and – worse – it silently reshuffles edges that are already set up, queued, or finished, so the network you get back is not the network you were running.

Everything here therefore edits. Network is frozen, so each function returns a new one, leaving the input untouched. The edges that were already there keep their identity, their mappings, and their costs; only the requested change and its consequences are new.

Invariants

Each function validates its result before returning it, so a surgery that would produce an inconsistent network fails rather than handing one back. Beyond that:

• Nothing here re-scores an existing edge. Costs are comparable across a surgery only because the untouched edges are literally the same objects.

• Connectivity is protected by default. delete_edge() refuses an edge whose removal would split the network, and names the two sides. Deleting one anyway is available, spelled out, and recorded on unmet_constraints.

• A CBFE edge is never spent to satisfy a degree target, exactly as in the planner. It is used only where a relative edge cannot reach at all – appending a ligand nothing maps to, or bridging two components – and only when cbfe_mode allows it.

cyclize_around_component deserves a note: Konnektor declares it and raises NotImplementedError. It is implemented here because a ligand that lies on no cycle has a free energy nothing checks, and after a deletion or an append that is exactly the ligand you have.

rbfenetmap.core.surgery.append_ligand(network, ligand, *, n_edges=2, mapper='mcss-e2', scorer='linear', mapping_options=None)

Add ligand to network, connecting it with its n_edges cheapest edges.

Parameters:

• network (Network) – The network to extend. Not modified.

• ligand (Ligand) – The new vertex. Its name must not already be in the network.

• n_edges (int, optional) – How many edges to attach it by. Two is the default because one leaves the new ligand hanging off a bridge, where its free energy is checked by nothing; the second edge is what puts it on a cycle.

• mapper (AbstractMapper or AbstractScorer or str, optional) – Used to evaluate the new ligand against every existing one. Pass the same ones the network was planned with – costs from two different scorers are not comparable.

• scorer (AbstractMapper or AbstractScorer or str, optional) – Used to evaluate the new ligand against every existing one. Pass the same ones the network was planned with – costs from two different scorers are not comparable.

• mapping_options (MappingOptions, optional)

Returns:

With the new ligand, its new edges, and the new candidates appended to the audit trail. Existing edges are untouched.

Return type:

Network

Raises:

• ValueError – If the ligand’s name is already present.

• rbfenetmap.core.exceptions.NetworkPlanError – If nothing in the network maps to the new ligand and cbfe_mode does not permit a counterpoised edge. The message lists why each pair was rejected.

Notes

Only the new ligand’s pairs are mapped: appending to an n-ligand network costs n mappings rather than the n(n+1)/2 a re-plan would.

A shortfall – fewer feasible partners than n_edges – is best-effort and lands on unmet_constraints, matching how the planner treats edges_per_ligand. Having no feasible partner is not a shortfall but a failure, because the result would be a disconnected network the caller did not ask for.

rbfenetmap.core.surgery.concatenate_networks(a, b, *, n_bridges=2, mapper='mcss-e2', scorer='linear', mapping_options=None)

Join two disjoint networks with n_bridges new edges.

Parameters:

• a (Network) – Neither is modified. Their ligand sets must be disjoint.

• b (Network) – Neither is modified. Their ligand sets must be disjoint.

• n_bridges (int, optional) – How many edges to build across the join. Two by default, for the same reason append_ligand() attaches two: a single bridge is checked by nothing, whereas two put the join itself on a cycle and make the relative offset between the two halves verifiable.

• mapper (AbstractMapper or AbstractScorer or str, optional)

• scorer (AbstractMapper or AbstractScorer or str, optional)

• mapping_options (MappingOptions, optional)

Returns:

The two ligand sets, both edge sets, and the new bridges.

Return type:

Network

Raises:

• ValueError – If the two share a ligand – merge_networks() handles that case.

• rbfenetmap.core.exceptions.NetworkPlanError – If no cross pair is feasible and cbfe_mode does not permit a counterpoised bridge.

Notes

Every cross pair is evaluated, which is len(a) * len(b) mappings. That is the honest cost of finding the best join rather than a plausible one, and it is still far below re-planning the union.

Bridges after the first are chosen to land on ligands the earlier bridges did not already use. Two bridges sharing an endpoint make that one ligand a single point of failure for the whole join, which is most of what the second bridge was bought to avoid.

rbfenetmap.core.surgery.cyclize_around_component(network, component=None, *, max_cycle_size=None, mapper=None, scorer='linear', mapping_options=None)

Add edges until every ligand in component lies on a cycle.

Parameters:

• network (Network) – Not modified.

• component (Iterable[str], optional) – The ligands to put on cycles. None means every ligand in the network. A connected component’s name is the usual argument – after a concatenate_networks() or a deletion, it is the newly attached or newly exposed part that has ligands hanging off bridges.

• max_cycle_size (int, optional) – Ignore closures longer than this. None takes the network’s own max_cycle_size.

• mapper (AbstractMapper or str, optional) – When given, pairs inside component that were never evaluated are mapped now. None (the default) restricts the search to the candidates the network already carries, which costs nothing and is usually enough: the pool from the original plan holds far more feasible pairs than the plan selected.

• scorer (AbstractScorer or str, optional)

• mapping_options (MappingOptions, optional)

Returns:

With the added edges. A ligand that could not be put on a cycle is recorded on unmet_constraints rather than being an error, matching how the planner reports a cycle-coverage shortfall.

Return type:

Network

Notes

A ligand on no cycle has a free energy nothing checks: every path to it runs through a bridge, so an error on that bridge moves the ligand’s number and shows up nowhere. That is why this exists as its own operation rather than as a re-plan – after an append or a deletion, one or two ligands are in exactly that state and the rest of the network is fine.

Konnektor declares the same operation and raises NotImplementedError for it.

Candidates are ranked by how many new ligands they put on a cycle, then by kind, then by cycle length and cost – the planner’s ranking, including its preference for a relative edge over a counterpoised one that would buy the same coverage.

rbfenetmap.core.surgery.delete_edge(network, pair, *, must_stay_connected=True)

Remove one edge from network.

Parameters:

• network (Network) – Not modified.

• pair (str or tuple[str, str]) – The edge, as "a~b" or as a pair of names. Direction is irrelevant – selection is undirected, so "b~a" names the same edge.

• must_stay_connected (bool, optional) – Refuse the deletion if it would split the network. On by default.

Returns:

Without that edge. The remaining edges keep their identity and costs.

Return type:

Network

Raises:

ValueError – If the edge is not in the network, or it is a bridge and must_stay_connected is set. The bridge message names the two groups the edge is the only link between, because “that would disconnect the network” alone does not say what to add instead.

Notes

This is the failure path of a campaign: an edge whose λ windows will not converge, or whose setup turns out to be wrong. Deleting it does not re-plan – see rbfenetmap.core.replanning.replan_after_diagnostics() for the version that refills the gap from the candidate pool.

rbfenetmap.core.surgery.merge_networks(a, b)

Merge two networks that share at least one ligand.

Parameters:

• a (Network) – Neither is modified. They must have at least one ligand name in common, and any shared name must denote the same molecule in the same atom order.

• b (Network) – Neither is modified. They must have at least one ligand name in common, and any shared name must denote the same molecule in the same atom order.

Returns:

The union of both ligand sets and both edge sets. A pair selected by both networks keeps the cheaper of the two edges. a’s options are carried through.

Return type:

Network

Raises:

ValueError – If the two share no ligand – which is concatenate_networks()’ job, and the message says so – or if a shared name denotes different molecules.

Notes

Sharing a ligand is what makes the result comparable rather than merely combined: free energies from the two networks are on the same scale only through a path that joins them, and a shared vertex is that path. Two networks with several shared ligands also gain cycles through them for free, which is why no bridging is done here.

The result is not re-planned, so it may exceed edges_per_ligand around the shared ligands. That is deliberate: the extra edges already exist and dropping them would discard work.

Diagnostics-driven replanning: prune the edges the analysis distrusts, refill the gaps.

An RBFE campaign is a loop – plan, run, analyse, replan – and this module is the return leg. It takes a per-edge diagnostic from the analysis stage, drops the edges that diagnostic condemns, and hands the pruned pool back to the planner so the network is rebuilt around what is left.

The Lagrange Multiplier Index

FE-ToolKit’s edgembar fits the whole network at once under the constraint that every cycle closes. Each edge’s Lagrange multiplier measures how hard that constraint had to pull on it: a large Lagrange Multiplier Index (LMI) means the edge disagrees with the consensus its cycles impose, which is the signature of a poorly converged or badly set up transformation. The CBFE paper does exactly this by hand on BACE1 and BRD4 – inspect the worst edges, drop them, re-run.

What LMI pruning is and is not worth

Pruning high-LMI edges substantially reduces cycle-closure error and leaves MUE and RMSE against experiment essentially unchanged. That is not a disappointing result, it is the correct interpretation of the quantity: hysteresis is a sampling diagnostic, not an accuracy predictor. A network can close every cycle perfectly and still sit a kcal/mol off the experimental values, because a systematic error in the force field or the protonation state moves every edge in a cycle the same way and cancels exactly where hysteresis would have shown it.

So use this to find edges that are internally inconsistent and worth re-running or replacing. Do not use it, and do not report it, as a route to better agreement with experiment.

Ingesting the diagnostic

load_edge_lmi() reads a small, explicit format of this package’s own – a mapping of "source~target" to a number, as a dict or a JSON file. It does not read edgembar’s on-disk output. Writing a parser against a format that could not be verified against a real file would be guesswork dressed up as an integration; extracting the multipliers from an edgembar analysis and writing this JSON is a short script on the user’s side today, and a first-class reader is a follow-up that needs a real file to develop against.

rbfenetmap.core.replanning.cycle_closure_errors(network, values)

Sum a per-edge quantity around each independent cycle of the network.

Parameters:

• network (Network)

• values (Mapping) – Per-edge quantities keyed directionally, as "source~target" or as a (source, target) tuple: the value is the quantity for that direction, typically the computed ΔΔG. The reverse direction is filled in as its negative, so only one orientation of each edge need be supplied.

Returns:

Keyed by the cycle’s ligand names in traversal order; the value is the signed sum around it, which should be zero and is not.

Return type:

dict[tuple[str, …], float]

Notes

Uses a cycle basis, not every cycle in the graph: the sums around a basis determine the sums around all the rest, and the number of cycles in a dense network is exponential. A cycle with an edge missing from values is dropped rather than summed over what is present, since a partial sum around a loop is not a closure error.

This is here to make the claim in the module docstring checkable on your own data – prune, replan, re-run, and watch these shrink while the errors against experiment do not.

rbfenetmap.core.replanning.load_edge_lmi(source)

Read per-edge Lagrange Multiplier Indices into a keyed mapping.

Parameters:

source (Mapping or str or pathlib.Path) –

Either a mapping already in memory, or the path to a JSON file. Accepted shapes:

• {"lig_a~lig_b": 0.42, ...} – the usual one;

• {("lig_a", "lig_b"): 0.42, ...} – in memory only, since JSON has no tuple keys;

• {"edges": {"lig_a~lig_b": 0.42, ...}} – a wrapper, so a file may carry other analysis output alongside.

Returns:

Keyed by unordered endpoint pair, because an LMI is a property of the transformation and the transformation is undirected: an analysis that reports b~a describes the edge the planner selected as a~b.

Return type:

dict[tuple[str, str], float]

Raises:

ValueError – If the document is not one of the shapes above, a key is not an edge, a value is not a number, or the same unordered pair appears twice with different values – which means the analysis and the network disagree about what the edges are, and picking one silently is the failure this package refuses everywhere else.

Notes

This is deliberately a small format of this package’s own, not edgembar’s. See the module docstring: a parser written against a file format nobody could check would be a guess presented as an integration.

rbfenetmap.core.replanning.lmi_threshold(values, *, quantile=0.9)

Return the quantile cut point of values, linearly interpolated.

Parameters:

• values (Sequence[float])

• quantile (float, optional) – In [0, 1]. 0.9 cuts at the worst tenth.

Returns:

The cut. Edges strictly above it are the ones pruned, so a network whose LMIs are all equal loses none of them however low the quantile.

Return type:

float

Raises:

ValueError – If values is empty or quantile lies outside [0, 1].

Notes

Interpolated rather than nearest-rank, and on a small network the difference decides whether anything is pruned at all. With five edges, nearest-rank at 0.9 lands exactly on the largest observed value, and “strictly above the largest value” is nothing – so the worst edge in the network would survive a request to prune the worst tenth of it. Interpolation puts the cut strictly between the two largest values instead, where the answer does not depend on a tie at the boundary.

rbfenetmap.core.replanning.replan_after_diagnostics(network, lmi, *, threshold=None, quantile=0.9, max_pruned=None, require_complete=True, keep_existing=True, planner='mst')

Prune the high-LMI edges and re-plan the network without them.

Parameters:

• network (Network) – A network whose candidate pool is intact – rbfenet plan writes it, and it is what makes this cheap. Not modified.

• lmi (Mapping[tuple[str, str], float])

• threshold (float | None) – Passed to select_high_lmi_edges().

• quantile (float) – Passed to select_high_lmi_edges().

• max_pruned (int | None) – Passed to select_high_lmi_edges().

• require_complete (bool) – Passed to select_high_lmi_edges().

• keep_existing (bool, optional) – Hold the surviving edges of the current network in place, so the replan changes only the gaps. On by default, and it is the difference between a replan you can act on and one you cannot: those edges are set up, queued, or already finished, and a selection pass free to reshuffle them hands back a network that is not the one being run. Turn it off for a clean re-selection over the pruned pool – the right choice before anything has been submitted.

• planner (str, optional) – Planner plugin used for the re-plan. The default re-runs the one this package plans with.

Returns:

The replanned network and the pairs that were pruned. The pruned pairs are also on the returned network’s options.banned_edges.

Return type:

tuple[Network, tuple[tuple[str, str], …]]

Raises:

rbfenetmap.core.exceptions.NetworkPlanError – If the pool cannot support a network once the pruned pairs are banned. The planner’s usual diagnostics apply: it names the components and the rejected candidates that would have bridged them.

Notes

Pruning is expressed as a ban, and the re-plan is the ordinary planner. That is the whole design. Deleting the edges and patching the holes by hand would need a second, parallel selection strategy that would drift from the real one; banning them and re-running means the replanned network satisfies exactly the same guarantees the first one did – spanning, degree targets, cycle coverage, the CBFE eligibility ladder – with a smaller pool. Nothing new is mapped: the replacements come from the candidates the original run already scored and did not select.

A pruned edge is banned rather than merely dropped because a re-plan over an unmodified pool would simply select it again – it was, after all, the cheapest edge there.

keep_existing is expressed the same way, as forced edges. Both halves of the request therefore travel to the planner as ordinary constraints, and the planner resolves them with the machinery it already has – including refusing outright if a surviving edge and a pruned one cannot both be honoured.

rbfenetmap.core.replanning.select_high_lmi_edges(network, lmi, *, threshold=None, quantile=0.9, max_pruned=None, require_complete=True)

Return the selected edges whose LMI exceeds the cut, worst first.

Parameters:

• network (Network) – Only its selected edges are considered; a candidate that was never run has no diagnostic to read.

• lmi (Mapping[tuple[str, str], float]) – From load_edge_lmi().

• threshold (float, optional) – Absolute cut. Edges with an LMI strictly greater than this are selected. When omitted the cut is taken from quantile.

• quantile (float, optional) – Used only when threshold is None.

• max_pruned (int, optional) – Keep at most this many, taking the worst. A guard for the case where the whole network scores badly: pruning half of it is a statement that the run failed, not a repair, and it should be a deliberate act rather than a quantile’s side effect.

• require_complete (bool, optional) – Raise if any selected edge has no LMI value. On by default: treating a missing value as zero would silently exempt exactly the edges an analysis failed to produce a number for, which are not the edges one wants to trust by default.

Returns:

Unordered pairs, ordered by descending LMI.

Return type:

tuple[tuple[str, str], …]

Raises:

ValueError – If require_complete is set and an edge is missing a value, or if no selected edge has one at all.

Notes

Forced edges are never returned. A user who pinned an edge has asserted it must be in the network, and a diagnostic does not override that – it is reported as skipped instead, since “your forced edge is the worst edge in the network” is worth reading.

Intermediate ligands

What an intermediate generator proposes, and how a proposal becomes a ligand.

A generator invents a molecule that sits between two ligands the mapper cannot relate cheaply, turning one hard edge into two easier ones. This module holds the vocabulary it speaks in – ProposedMolecule, ProposedLink, IntermediateProposal – plus the naming rule and the one function that turns a proposal into a real Ligand.

These types live here rather than in rbfenetmap.core.models on purpose: they carry RDKit molecules and an atom-level correspondence, and models is deliberately free of chemistry so the data model stays importable and testable without one.

Three constraints the types enforce, and why

A proposed molecule carries no conformer. Posing is centralised in rbfenetmap.core.posing, for the same reason rbfenetmap.core.descriptors centralises scoring inputs: a generator that poses its own output makes the quality of every intermediate depend on which generator produced it, and makes the pose unauditable. A generator that knows where its atoms belong says so through a complete ProposedMolecule.parent_atom_map, which is strictly more useful than coordinates – the poser can act on a correspondence, and can tell you which parent every atom came from afterwards.

A link’s hint is advisory and can never become an total. A generator may know that one of its proposals is more promising than another, and that is worth recording; it is not a cost. Letting a hint reach a score would put a second scoring system in the package on a different scale, which is the dual of the rule that a scorer must not invent rejections.

A proposal’s rejection is a plain string. RejectionReason is the vocabulary of edge feasibility. Reusing it here would make core_geometry_mismatch mean two different things depending on which object it was read from.

Naming

intermediate_name() is content-addressed – int_{a}_{b}_{sha1(smiles)[:6]} with the parents sorted – rather than counter-based. The same intermediate is often reachable from the same gap by two routes, and a hash makes that a natural dedupe instead of a pair of near-identical ligands; it also keeps a run reproducible under jobs > 1, where a counter’s value depends on which worker finished first.

rbfenetmap.core.intermediates.INTERMEDIATE_KIND = 'intermediate'

The kind this package writes.

rbfenetmap.core.intermediates.INTERMEDIATE_NAME_PREFIX = 'int_'

Prefix every generated ligand name carries. Matches Ligand’s name pattern and contains no EDGE_SEPARATOR, so a generated name survives edge_key() and parse_edge_key() unchanged.

rbfenetmap.core.intermediates.INTERMEDIATE_MODES: tuple[str, ...] = ('off', 'bridge', 'gaps')

How freely the pipeline may invent ligands.

There is deliberately no "all". For NetworkOptions.cbfe_mode that member means “every edge is counterpoised”, which is a coherent request; “every pair gets an intermediate” is not one, because most pairs already have a perfectly good direct edge and inventing a molecule for them adds two calculations to avoid nothing.

class rbfenetmap.core.intermediates.IntermediateOptions(mode='off', generator='pairmap', max_intermediates=None, max_gaps=None, max_molecules=4, seed=61453, max_pose_attempts=10, pose_rmsd_factor=0.5, min_link_score=0.2, max_dist=3, max_cycle=4, max_subgraph_dist=4, beta=0.1)

Bases: object

Whether to invent ligands, how many, and how hard to try posing them.

Parameters:

• mode ({"off", "bridge", "gaps"}, optional) –

  Which gaps are offered to the generator.

  • "off" (default) – none. No generator is even constructed, so a run that does not ask for intermediates never imports one.

  • "bridge" – only pairs whose endpoints fall in different components of the feasible candidate graph. This is the mode that turns a hard connectivity failure into a planned network, and it is the analogue of cbfe_mode="bridge".

  • "gaps" – everything "bridge" does, and additionally infeasible pairs inside a component. Those are already reachable by some path, so an intermediate there buys accuracy rather than connectivity.

  Forced pairs with no feasible mapping are offered under both non-off modes: the user demanded that comparison, and an intermediate is the only way to keep it relative.

• generator (str, optional) – Registered name of the generator plugin to construct. Resolved lazily, and only when mode is not "off".

• max_intermediates (int, optional) – Cap on how many ligands one run may invent in total. None means only the per-gap cap and the edge budget constrain it.

• max_gaps (int, optional) – Cap on how many gaps are offered to the generator, taken in decreasing fingerprint similarity. None offers every gap. This is the knob that keeps generation off the tail of an O(n^2) rejection list, where the pairs are least similar and least likely to be bridgeable anyway.

• max_molecules (int, optional) – Cap on how many molecules a generator may propose for one gap. A generator that enumerates every single-substituent swap on a heavily decorated pair can produce dozens; each one costs an embedding and a minimisation, so the cap is the knob that keeps generation from dominating a run.

• seed (int, optional) – Base RDKit random seed for posing, so a run is reproducible.

• max_pose_attempts (int, optional) – Embedding attempts spent on one molecule before giving up on it.

• pose_rmsd_factor (float, optional) – Fraction of SoftcorePolicy.core_rmsd_threshold an accepted pose must stay under. See POSE_RMSD_FACTOR.

• min_link_score (float, optional) – Lowest link score a generator may consider worth proposing, on a (0, 1] similarity scale where 1 is “no atoms change at all”.

• max_dist (int, optional) – Longest source-to-target path, in links, the generator may propose. At least 2: a one-link path is the direct transformation that was already rejected.

• max_cycle (int, optional) – Largest cycle the generator may build to give a proposed link a second, independent route. Cycles are what turn a chain of intermediates into a network with a closure error to check.

• max_subgraph_dist (int, optional) – How far from either parent, in links, a molecule may sit and still be considered for the subnetwork. Bounds the search, and must be at least max_dist.

• beta (float, optional) – Decay rate of the exponential link score, in inverse heavy atoms. The published default of 0.1 is the same constant LOMAP’s similarity uses, and it is what makes min_link_score 0.2 mean “at most about sixteen heavy atoms change”.

Raises:

ValueError – If mode is unknown, or any budget is out of range.

Notes

Reachable as intermediates, nested the way SoftcorePolicy is. It was deliberately left off NetworkOptions while nothing read it, because a field on the serialized options block that no stage consumes is a knob that lies. The pipeline now consumes every one of these, so the nesting is what makes a planned network state the settings that invented its vertices.

The last five fields keep the names and the published defaults of the PairMap constants (MIN_SCORE, MAX_DIST, MAX_CYCLE, MAX_SUBGRAPH_DIST, beta) from Furui et al., J. Chem. Inf. Model. 2025, 65, 705-721 (doi:10.1021/acs.jcim.4c01634), so a reader can grep them against the paper. They live on the shared options object rather than on the generator because a plugin’s parameters have to survive serialization to make an invented ligand reproducible, and describe_parameters() is a report, not a record. A generator that does not search a subnetwork simply ignores them.

property enabled: bool

Whether generation runs at all.

property bridges_components: bool

Whether cross-component gaps are offered to the generator.

property fills_internal_gaps: bool

Whether infeasible pairs inside one component are offered too.

class rbfenetmap.core.intermediates.IntermediateProposal(source, target, generator, molecules=(), links=(), rejection=None, trace=())

Bases: object

A generator’s complete answer for one gap.

Parameters:

• source (str) – The gap the generator was asked about.

• target (str) – The gap the generator was asked about.

• generator (str) – Registered plugin name.

• molecules (tuple[ProposedMolecule, ...], optional) – What it suggests inserting. Empty is a legitimate answer.

• links (tuple[ProposedLink, ...], optional) – Sub-edges it expects to become feasible.

• rejection (str, optional) – Why it proposed nothing. A plain str, never a RejectionReason.

• trace (tuple[str, ...], optional) – Human-readable log of what it tried.

Raises:

ValueError – If the two endpoints are the same ligand.

property proposed: bool

Whether the generator suggested anything at all.

class rbfenetmap.core.intermediates.ProposedLink(source, target, hint=None, detail=<factory>)

Bases: object

A sub-edge the generator expects its molecule to make possible.

Parameters:

• source (str) – Endpoints. One is normally the proposed molecule, the other a parent.

• target (str) – Endpoints. One is normally the proposed molecule, the other a parent.

• hint (float, optional) – How promising the generator believes the link to be, lower being better.

• detail (Mapping[str, Any], optional) – Free-form annotations.

Notes

The hint can never become an total. Nothing in this package reads it as a cost, and nothing should: a generator’s opinion about its own output is on a scale only that generator knows, so promoting it to a score would mean two incomparable numbers competing inside the planner’s edge ordering. The link exists to say which edges are worth evaluating; the scorer says what they cost.

class rbfenetmap.core.intermediates.ProposedMolecule(mol, parents, parent_atom_map=<factory>, hint=None, detail=<factory>)

Bases: object

One molecule a generator suggests inserting into a gap.

Parameters:

• mol (rdkit.Chem.Mol) – The molecule, without a conformer. Any conformers present are stripped at construction rather than rejected, so a generator that happened to build its molecule from a posed parent does not have to remember to clear them – but the coordinates are discarded either way, because posing is centralised.

• parents (tuple[str, ...]) – Names of the ligands it was derived from, sorted at construction.

• parent_atom_map (Mapping[str, Mapping[int, int]]) – {parent name: {proposed atom index: parent atom index}}. May be empty, or cover only some parents; the poser falls back to an MCS search for whatever is missing and records that it had to.

• hint (float, optional) – The generator’s own ordering preference among its proposals, lower being more promising. Advisory only – see the module docstring.

• detail (Mapping[str, Any], optional) – Free-form annotations carried into the ligand’s provenance.

Raises:

ValueError – If parents is empty, or parent_atom_map names a parent that is not in parents.

rbfenetmap.core.intermediates.describe_intermediate_attempts(records)

Summarise what generation was asked to do and what came of it.

Parameters:

records (Sequence[IntermediateRecord]) – One per gap attempted, in the order they were attempted.

Returns:

A short paragraph naming each gap offered to the generator and why it was refused, or the empty string when records is empty.

Return type:

str

Notes

The paragraph exists to discharge the same obligation the CBFE branch of _describe_disconnection() discharges: “disconnected” alone tells a user nothing they can act on, and a user who switched generation on and still got a disconnection needs to know whether their gaps were never offered, offered and declined, or bridged by molecules that failed the geometry gate. Those three call for entirely different responses – raise max_gaps, change generator, loosen core_rmsd_threshold – and only this record distinguishes them.

rbfenetmap.core.intermediates.intermediate_name(parents, mol)

Return the content-addressed name for an invented molecule.

Parameters:

• parents (Sequence[str]) – Names of the ligands it bridges. Sorted here, so the caller need not.

• mol (rdkit.Chem.Mol) – The molecule. Only its canonical SMILES is used, so a name is stable across a re-embedding and across atom reordering.

Returns:

int_{a}_{b}_{sha1(canonical_smiles)[:6]}.

Return type:

str

Raises:

ValueError – If parents is empty, or the resulting name would not be a legal ligand name.

Notes

Content-addressed rather than counter-based, and the parents are part of the address rather than only the structure. The hash makes the same molecule proposed twice for the same gap collapse to one ligand – which is the common case, since a gap is usually reachable from either end. The parent tokens stay in because two gaps that happen to want the same molecule want it for different reasons, and a name that hid that would leave two different provenances competing for one vertex.

Hydrogens are suppressed before canonicalisation. A generator that hands over a molecule with explicit hydrogens and one that leaves them implicit are proposing the same thing, and they must not get different names for it.

rbfenetmap.core.intermediates.reserve_intermediate_names(names, *, enabled)

Refuse user ligand names that would collide with generated ones.

Parameters:

• names (Iterable[str]) – The user’s ligand names.

• enabled (bool) – Whether intermediate generation is switched on for this run.

Raises:

ValueError – If enabled and any name starts with INTERMEDIATE_NAME_PREFIX.

Return type:

None

Notes

The prefix is reserved only when the feature is on, which is the whole point of taking enabled rather than reading it from a global. A user with a ligand honestly called int_3 is running a plain network plan that cannot possibly generate a conflicting name, and refusing it would be a compatibility break bought for nothing.

rbfenetmap.core.intermediates.synthesize_ligand(proposed, parents, *, generator, softcore=None, options=None, mapping_options=None)

Pose proposed against its parents and wrap the result in a ligand.

Parameters:

• proposed (ProposedMolecule) – The conformer-free molecule a generator suggested.

• parents (Mapping[str, Ligand]) – The real ligands, keyed by name. Must contain every name in ProposedMolecule.parents.

• generator (str) – Registered name of the generator, recorded in the provenance.

• softcore (SoftcorePolicy, optional) – Supplies core_rmsd_threshold. Defaults are used when omitted.

• options (IntermediateOptions, optional) – Supplies the seed, the attempt budget, and the RMSD factor.

• mapping_options (MappingOptions, optional) – Only consulted by the poser’s MCS fallback.

Returns:

• ligand (Ligand or None) – None when the pose was rejected; the PoseResult says why.

• result (PoseResult) – Always returned, successful or not, because the trace is what a user needs when an intermediate does not appear.

Raises:

KeyError – If parents does not contain a name the proposal claims. That is a caller bug rather than a chemistry outcome, so it raises where posing failures do not.

Return type:

tuple[Ligand | None, PoseResult]

Notes

The parents are handed to the poser in sorted-name order, which is also the order intermediate_name() uses. That matters because donors are consumed in priority order: an unordered iteration would let two runs of the same input take scaffold coordinates from different parents and produce two slightly different poses under the same name.

Posing an invented molecule into its parents’ binding-site frame.

Everything else in this package assumes its ligands arrive already posed in a shared frame – the in-place core_rmsd() gate depends on it. A molecule the planner invents has no such pose, and has no crystal structure to inherit one from, so it has to be built. This module is the only place that happens.

Neither of the two obvious tools is right

AllChem.ConstrainedEmbed pins a scaffold and re-embeds everything else from scratch. That throws away the one thing an intermediate uniquely has: every heavy atom of a hybrid molecule has a specific corresponding atom in a parent whose position in the pocket is already known. Re-embedding from the scaffold outwards lets ETKDG place an ortho substituent where the parent had it meta – a valid conformer of the right molecule in the wrong place, which is exactly the failure the geometry gate is least able to explain to a user.

rbfenetmap.core.align is the other near-miss. It recovers a common frame for molecules that already have good conformers. An intermediate has none until this module makes one, so there is nothing for it to align.

The algorithm

1. Add explicit hydrogens before embedding. Ligand forbids implicit hydrogens, and AddHs(addCoords=True) after the fact places them by rule rather than by geometry – fine for a picture, not for a starting structure. RDKit appends the new atoms, so donor maps built on the heavy molecule stay valid.

2. Take the correspondence from the generator. It built the molecule by editing a parent and therefore knows which atom came from where. Re-deriving it with GetSubstructMatch re-introduces precisely the symmetry coin-flip that match_selection="fewest_fragments" exists to remove. Only when a generator declines to supply one do we fall back to an MCS search, and the result is tagged "mcs_fallback" so the weaker provenance is visible in a report.

3. Seed a coordMap from the donors’ heavy coordinates and embed with useRandomCoords. Heavy only, because a coordMap naming every atom leaves the distance-bounds smoothing nothing to solve and RDKit aborts on the degenerate system; the mapped hydrogens are restrained in the next step instead.

4. Restore exactness with a force field. coordMap applies distance-bounds constraints, not exact placements; an embedding that satisfies every bound can still sit an angstrom from where it was asked to. Fixed extra points plus distance constraints, then a minimisation, is what actually pins the mapped atoms.

5. Rigid Kabsch fit onto the donor coordinates over the mapped set. Minimisation moves everything a little, and the fit removes the net drift without distorting the relaxed geometry the minimiser just produced.

6. Gate on the in-place core_rmsd – the same measurement, with the same superpose_first=False, that the pipeline will judge the resulting edges by.

Why the gate is not the safety net

You do not need to trust this module. The in-place core RMSD check on the A~M and M~B sub-edges is already a complete test of whether M is posed in the parents’ frame: a badly posed intermediate comes back as an ordinary CORE_GEOMETRY_MISMATCH and the proposal is dropped for failing to close the gap it was invented for. The job here is to make that check pass often, not to make it unnecessary.

Failures are data

Every way this can fail returns a PoseResult carrying a PoseRejection, never an exception. One molecule that will not embed among hundreds must not stop a run, for the same reason one impossible pair does not stop the mapping stage.

rbfenetmap.core.posing.POSE_RMSD_FACTOR = 0.5

Fraction of SoftcorePolicy.core_rmsd_threshold an accepted pose must stay under.

Named rather than inlined because it encodes a judgement, not a tolerance: an intermediate posed at exactly the threshold would be accepted here and then rejected by the very next stage, having spent an embedding and a minimisation to learn it. Half the budget leaves room for the mapping the pipeline will find to differ slightly from the one the generator handed over.

class rbfenetmap.core.posing.PoseDonor(name, mol, atom_map=None)

Bases: object

One parent, and which of its atoms lend their coordinates.

Parameters:

• name (str) – The parent ligand’s name, recorded in the trace and the provenance.

• mol (rdkit.Chem.Mol) – The parent molecule, already posed in the binding-site frame.

• atom_map (Mapping[int, int], optional) – {intermediate atom index: parent atom index}. None asks pose_intermediate() to recover a correspondence by MCS instead, which is strictly weaker and is reported as such.

Notes

Donors are consumed in order and the first to claim an intermediate atom keeps it. An intermediate is a hybrid, so its two parents will generally both map some of the same scaffold atoms; taking the earlier donor’s coordinates rather than averaging is deliberate, since the average of two poses is a pose neither parent has.

class rbfenetmap.core.posing.PoseRejection(*values)

Bases: str, Enum

Why a molecule could not be posed in its parents’ frame.

A str enum so a value drops straight into rejection, which is a plain string by design – RejectionReason is the vocabulary of edge feasibility and is deliberately not reused here.

INVALID_MOLECULE

The proposed molecule does not sanitize. A generator bug, recorded rather than raised so the rest of the run continues.

CHARGE_MISMATCH

The intermediate’s net formal charge differs from a parent’s. Refused outright: an intermediate exists to split one hard edge into two easier ones, and a charge change would instead split it into two harder ones.

STEREO_UNDEFINED

An unassigned stereocentre or double-bond configuration. Refused because nobody downstream can parameterise it – the ambiguity would be resolved arbitrarily by whichever tool touched the molecule first.

NO_DONOR_ATOMS

No correspondence to any parent could be established, so there is nothing to pose against.

EMBED_FAILED

EmbedMolecule could not satisfy the distance bounds within the attempt budget.

FORCEFIELD_FAILED

Neither MMFF nor UFF could be set up for the molecule.

POSE_RMSD_EXCEEDED

A conformer was produced, but it sits too far from the donor coordinates to be worth handing to the feasibility stage.

class rbfenetmap.core.posing.PoseResult(mol=None, rmsd=inf, method='parent_atom_map', rejection=None, attempts=0, trace=(), detail=<factory>)

Bases: object

The outcome of one posing attempt.

Parameters:

• mol (rdkit.Chem.Mol, optional) – The posed molecule, with explicit hydrogens and exactly one 3D conformer, or None when rejection is set.

• rmsd (float) – In-place RMSD of the mapped atoms against their donor coordinates. Reported even on a PoseRejection.POSE_RMSD_EXCEEDED rejection, because how badly the pose missed is the diagnostic.

• method (str) – "parent_atom_map" or "mcs_fallback".

• rejection (str, optional) – A PoseRejection value, or None on success.

• attempts (int) – Embedding attempts spent.

• trace (tuple[str, ...]) – Human-readable log, carried into the IntermediateRecord.

• detail (Mapping[str, Any]) – Structured extras, notably n_mapped and the donors consulted.

property posed: bool

Whether a usable conformer was produced.

rbfenetmap.core.posing.pose_intermediate(mol, donors, *, core_rmsd_threshold=2.0, rmsd_factor=0.5, seed=61453, max_attempts=10, mapping_options=None)

Give mol a conformer sitting in the frame its donors occupy.

Parameters:

• mol (rdkit.Chem.Mol) – The proposed molecule. Carries no conformer; one is added here.

• donors (Sequence[PoseDonor]) – The posed parents whose coordinates the pose is built from, in priority order.

• core_rmsd_threshold (float, optional) – The pipeline’s SoftcorePolicy.core_rmsd_threshold. The accept/reject bar is this times rmsd_factor.

• rmsd_factor (float, optional) – Fraction of the threshold an accepted pose must stay under. Default POSE_RMSD_FACTOR.

• seed (int, optional) – Base RDKit random seed. Attempt k uses seed + k, so retries explore rather than repeat.

• max_attempts (int, optional) – Embedding attempts before giving up.

• mapping_options (MappingOptions, optional) – Only consulted for the MCS fallback. Defaults are used when omitted.

Returns:

With PoseResult.mol set on success, or a PoseRejection value in PoseResult.rejection on any failure. Nothing here raises.

Return type:

PoseResult

Examples

>>> result = pose_intermediate(mol, [PoseDonor("lig_a", parent.mol, atom_map)])
>>> result.posed, round(result.rmsd, 3)
(True, 0.041)

Plugin registry

Lazy plugin registry.

Adapted from pharmaforge.core.pluginregistry, keeping that package’s conventions so the two are recognisably the same mechanism.

The design point worth preserving: a PluginSpec describes a plugin without importing it. Registration is pure metadata, and the implementation module is imported only when PluginRegistry.create() is called. That is what lets rbfenet plugins list every backend – including ones whose dependencies are absent – without importing RDKit, kartograf, or anything else, and what lets the whole test suite run with no optional dependency installed.

class rbfenetmap.core.pluginregistry.PluginRegistry

Bases: object

Store plugin metadata and construct plugins on demand.

The registry tracks registered plugins, which is not the same as active ones. Registration says a plugin exists; activation says the current configuration intends to use it. Keeping the two apart lets a caller enumerate everything installed while still restricting a given run to a chosen subset.

Initialize an empty registry.

register(spec)

Register spec.

Raises:

PluginError – If a plugin with the same (kind, name) is already registered. Silently overwriting would make plugin behaviour depend on import order.

Parameters:

spec (PluginSpec)

Return type:

None

get_spec(name, kind)

Return the spec for (kind, name).

Raises:

PluginError – If no such plugin is registered. The message lists the registered names for that kind, since the usual cause is a typo.

Parameters:

• name (str)

• kind (str)

Return type:

PluginSpec

activate(name, kind)

Mark (kind, name) active, registering nothing new.

Parameters:

• name (str)

• kind (str)

Return type:

None

deactivate(name, kind)

Mark (kind, name) inactive.

Parameters:

• name (str)

• kind (str)

Return type:

None

is_active(name, kind)

Whether (kind, name) is active.

Parameters:

• name (str)

• kind (str)

Return type:

bool

create(name, kind, **kwargs)

Import and instantiate a plugin.

This is the only method that imports anything.

Raises:

PluginError – If the plugin is unknown, its module or class cannot be imported, or its declared requirements are missing. A missing requirement is reported before the import is attempted, so the user sees "needs kartograf" rather than a raw ModuleNotFoundError.

Parameters:

• name (str)

• kind (str)

• kwargs (Any)

Return type:

Any

list_plugins(kind=None, *, active_only=False)

Return registered specs, optionally filtered by kind or activity.

Parameters:

• kind (str | None)

• active_only (bool)

Return type:

tuple[PluginSpec, …]

class rbfenetmap.core.pluginregistry.PluginSpec(name, kind, target, description='', requires=())

Bases: object

Describe a plugin without importing its implementation.

Parameters:

• name (str) – Unique name within a given plugin kind.

• kind (str) – The plugin category: "mapper", "scorer", "planner", or "exporter".

• target (str) – Import target in "package.module:ClassName" form.

• description (str, optional) – Short human-readable summary.

• requires (tuple[str, ...], optional) – Top-level modules the plugin’s backend needs. Probed with importlib.util.find_spec(), so availability can be reported without importing anything.

property missing_requirements: tuple[str, ...]

Required modules that cannot be located in this environment.

property available: bool

Whether every required module can be located.

Contracts

The mapper contract: propose an atom correspondence between two ligands.

class rbfenetmap.core.meta.mappers.AbstractMapper

Bases: ABC

Produce a common-core / soft-core partition for one candidate pair.

A mapper is responsible only for the correspondence. It does not need to produce a connected soft-core region: repairing fragmentation is the job of rbfenetmap.core.softcore.repair_softcore_connectivity(), which runs afterwards on every mapper’s output. Trying to enforce connectivity inside a mapper duplicates that logic and makes mappers harder to write and compare.

name

The registered plugin name, used in diagnostics and recorded on the mapping.

Type:

str

abstractmethod map_pair(source, target, options)

Return the atom correspondence between source and target.

Parameters:

• source (Ligand) – The two ligands, each with explicit hydrogens and one 3D conformer.

• target (Ligand) – The two ligands, each with explicit hydrogens and one 3D conformer.

• options (MappingOptions) – Search settings and the pre-repair core-pruning policy.

Returns:

A validated mapping. Construction enforces the contract, so an implementation that builds one via from_core_pairs() cannot return something malformed.

Return type:

AtomMapping

Raises:

rbfenetmap.core.exceptions.MappingError – If no correspondence can be produced at all. The caller converts this into a MAPPER_FAILED rejection rather than letting it abort the whole run – one impossible pair among hundreds should not stop the planning.

supports_pair(source, target)

Whether this mapper can handle the pair at all.

Cheap pre-check, called before map_pair(). The default accepts everything.

Parameters:

• source (Ligand)

• target (Ligand)

Return type:

bool

The scorer contract: reduce edge descriptors to a scalar cost.

class rbfenetmap.core.meta.scorers.AbstractScorer

Bases: ABC

Turn precomputed edge descriptors into a cost. Lower is better.

Notes

A scorer receives a plain Mapping[str, float] and nothing else – no molecules, no mapping object, no RDKit. Descriptors are computed once, centrally, by rbfenetmap.core.descriptors.compute_descriptors().

That narrow interface buys three things. Re-scoring a network under different weights costs nothing, because no mapping has to be recomputed. A scorer can be tested against hand-written dictionaries, with no chemistry in the test at all. And a third-party scorer cannot accidentally reach past its inputs and reintroduce a dependency on how the mapping was produced.

A scorer must not invent rejections. Feasibility is decided upstream by the mapper and the repair; rejections is passed in so the scorer can propagate it into the returned EdgeScore, not so it can add to it. A scorer that wants to express “this edge is terrible” returns a large finite cost – which leaves the planner free to use it anyway if the alternative is a disconnected network.

abstractmethod score_edge(descriptors, *, rejections)

Return the cost of an edge described by descriptors.

Parameters:

• descriptors (Mapping[str, float]) – Raw descriptor values from rbfenetmap.core.descriptors.compute_descriptors().

• rejections (Sequence[RejectionReason]) – Structural rejections already determined upstream. Non-empty means the implementation must return EdgeScore.rejected(*rejections, scorer=self.name).

Return type:

EdgeScore

describe_weights()

Return the scorer’s tunable weights, for display by rbfenet score.

Return type:

Mapping[str, float]

The planner contract: select the final edge set from scored candidates.

class rbfenetmap.core.meta.planners.AbstractNetworkPlanner

Bases: ABC

Choose which candidate transformations make up the network.

Notes

A planner selects; it does not judge feasibility. Candidates arrive already scored, with infeasible ones marked. An implementation must filter on feasible and must still place every candidate – feasible or not – on candidates.

Retaining the infeasible ones is what makes a disconnected result explicable. When the planner has to report that two groups of ligands cannot be joined, the rejected candidates that span the gap, and their reasons, are the actionable part of the message; without them the user gets “disconnected” and no idea what to loosen.

supports_cbfe: ClassVar[bool] = False

Whether this planner knows how to place counterpoised (CBFE) edges. Both of the modes that need planner cooperation – bridge and cycles – are expressed as decisions about where an edge goes, which only a planner that reasons about components and cycles can make. all needs nothing from the planner, because the pipeline hands it a pool that is already entirely CBFE.

supports_design: ClassVar[bool] = False

Whether this planner optimises a statistical design criterion. design names an objective, not a filter, so there is nothing a planner can do with it halfway: a planner that does not optimise the criterion would ignore the flag entirely, and a knob that is silently ignored is worse than one that is absent.

check_design_support(options)

Raise if options names a design criterion this planner cannot optimise.

Parameters:

options (NetworkOptions)

Raises:

rbfenetmap.core.exceptions.NetworkPlanError –

Return type:

None

Notes

The counterpart to check_cbfe_support(), and refused for the same reason. --design a_optimal under the mst planner would produce a perfectly ordinary minimum-spanning-tree network, with nothing anywhere to connect the result to the flag the user set – the --consistency graph failure mode this package already has one instance of and does not want a second.

check_cbfe_support(options)

Raise if options asks for CBFE placement this planner cannot do.

Parameters:

options (NetworkOptions)

Raises:

rbfenetmap.core.exceptions.NetworkPlanError –

Return type:

None

Notes

Called rather than silently ignored. A user who passes --cbfe bridge and a planner that cannot honour it would otherwise get a disconnected network, or a connectivity error, with nothing to connect either outcome to the flag they set – the same failure mode the scorers refuse for an unknown weight name.

abstractmethod plan(ligands, candidates, options)

Select edges and return the planned network.

Parameters:

• ligands (Mapping[str, Ligand]) – Every vertex, including any the planner ends up unable to connect.

• candidates (Sequence[Transformation]) – Scored candidates, feasible and infeasible.

• options (NetworkOptions) – The user’s selection knobs.

Return type:

Network

Raises:

rbfenetmap.core.exceptions.NetworkPlanError – If the constraints are unsatisfiable: a forced edge that is infeasible, an n_edges too small to span, or a disconnected candidate pool while require_connected is set. Constraints that are merely tight – an edges_per_ligand the pool cannot support – are recorded on unmet_constraints instead.

The exporter contract: serialize a planned network for a downstream program.

This is the seam the package hangs its “hooks to other programs” on. An exporter is how a network reaches Amber, a workflow engine, a viewer, or anything else, without any of those systems’ concerns leaking back into the core.

class rbfenetmap.core.meta.exporters.AbstractExporter

Bases: ABC

Write a planned network out in some downstream format.

name

The registered plugin name.

Type:

str

default_suffix

Extension used when destination names a file with none.

Type:

str

abstractmethod export(network, destination, **options)

Write network to destination.

Parameters:

• network (Network) – The planned network.

• destination (pathlib.Path) – A file or a directory, depending on the exporter. Exporters that emit one file per edge take a directory.

• **options – Exporter-specific settings, passed through from --exporter-opt.

Returns:

Every path written, so a caller can report or clean up.

Return type:

tuple[pathlib.Path, …]

Raises:

rbfenetmap.core.exceptions.ExporterError – If the network cannot be represented in the target format.

validate(network)

Check network can be exported, without writing anything.

Format-specific constraints that the core does not enforce belong here – the Amber exporter’s atom-name uniqueness requirement being the motivating case.

Called early by rbfenet plan --validate-exporter so a constraint that would only surface at export time is caught before the expensive mapping work runs, rather than after.

Raises:

rbfenetmap.core.exceptions.ExporterError – If the network violates a format constraint.

Parameters:

network (Network)

Return type:

None

The intermediate-generator contract: invent a molecule to bridge a gap.

class rbfenetmap.core.meta.intermediates.AbstractIntermediateGenerator

Bases: ABC

Propose molecules that split one hard transformation into two easier ones.

The fifth plugin kind, and the only one whose output changes the ligand set rather than the network over it. That is why the contract is narrow: a generator proposes, and nothing else. It does not pose its molecules – rbfenetmap.core.posing does, once, for everyone – it does not decide whether the resulting edges are feasible, and it does not price them. Each of those already has an owner, and a generator that took any of them over would make an intermediate’s quality depend on which generator happened to invent it.

name

The registered plugin name, recorded on every ligand the generator’s proposals become.

Type:

str

Notes

A generator that cannot help with a pair returns an IntermediateProposal with no molecules and a rejection string. It does not raise: an intermediate is an optimisation, and one gap that cannot be bridged is an ordinary outcome rather than an impossible request. This is the same rule that keeps a rejected edge out of rbfenetmap.core.exceptions.

abstractmethod propose(source, target, options, mapping_options)

Suggest molecules bridging the gap between source and target.

Parameters:

• source (Ligand) – The two real ligands, each with explicit hydrogens and one 3D conformer.

• target (Ligand) – The two real ligands, each with explicit hydrogens and one 3D conformer.

• options (IntermediateOptions) – How many molecules may be proposed, and the posing budget that will be spent on them.

• mapping_options (MappingOptions) – The same settings the mappers run under, so a generator that needs an MCS finds the one the pipeline would have found.

Returns:

Possibly empty, with rejection set to say why. Every proposed molecule must carry no conformer; a ProposedMolecule strips any it is given, so this is enforced rather than merely asked for.

Return type:

IntermediateProposal

supports_pair(source, target)

Whether this generator can attempt the pair at all.

Cheap pre-check, called before propose(). The default accepts everything.

Parameters:

• source (Ligand)

• target (Ligand)

Return type:

bool

describe_parameters()

Return the generator’s own settings, for the run record.

Returns:

JSON-friendly values. The default is empty.

Return type:

Mapping[str, Any]

Notes

Generators are the plugin kind most likely to carry knobs of their own – how far to search, which transformations to consider – and those knobs change what molecules a run invents. Reporting them alongside the network is what makes an invented ligand reproducible by someone who was not there when it was invented.

Validation and errors

Validation of the common-core / soft-core mapping contract.

A parmed-free port of BuildEdges._validate_mapping_result. Kept in its own module rather than inlined into __post_init__() so it can be called directly on a mapper’s raw output during debugging, and so the error messages live in one place where they can be kept specific.

Every message names the offending indices. A mapping is rejected at construction, so a vague message here becomes a vague failure hundreds of lines from the actual mistake.

rbfenetmap.core.validate.validate_mapping(mapping)

Enforce every invariant of the mapping contract.

Parameters:

mapping (AtomMapping) – The mapping to check.

Raises:

ValueError – With a message naming the specific offending atom indices.

Return type:

None

Notes

The partition check (sc_k and cc_k disjoint and jointly covering every atom) is the one that catches the largest class of real mapper bugs. An atom that is in neither set has no defined behaviour under the transformation: it is neither held fixed nor alchemically transformed. Mappers that forget hydrogens, or that build the soft-core from a stale copy of the core, fail exactly here.

Exception hierarchy for rbfenetmap.

Every error the package raises deliberately derives from RBFENetworkMapError, so a caller embedding this in a larger workflow can catch one type and know the failure came from network planning rather than from RDKit, NumPy, or its own code.

Note the distinction this hierarchy encodes, which matters throughout the package: a rejected edge is not an error. An edge whose soft-core cannot be repaired within the configured budget is a normal, expected outcome recorded as a RejectionReason on the transformation’s score. These exceptions are for situations where the caller asked for something impossible or inconsistent – a forced edge that cannot exist, a mapping that violates its own invariants, a plugin that is not installed.

exception rbfenetmap.core.exceptions.RBFENetworkMapError

Bases: Exception

Base class for every error raised by rbfenetmap.

exception rbfenetmap.core.exceptions.MappingError

Bases: RBFENetworkMapError

An atom mapping could not be produced, or violates the mapping contract.

exception rbfenetmap.core.exceptions.RepairError

Bases: RBFENetworkMapError

The soft-core repair could not run.

Raised for malformed input to the repair (for example a mapping whose indices do not match the molecules). A repair that runs correctly and concludes the edge is infeasible returns a RejectionReason instead – that is an answer, not an error.

exception rbfenetmap.core.exceptions.NetworkPlanError(message, rejected=())

Bases: RBFENetworkMapError

The requested network cannot be planned.

Raised when user constraints are unsatisfiable rather than merely tight: a forced edge that is infeasible, n_edges too small to span the ligands, or a candidate pool that is disconnected while require_connected is set.

Parameters:

• message (str) – The actionable description, which is what a bare str() of the error yields.

• rejected (Sequence, optional) – The infeasible candidates the planner had to work around, when it has them.

Return type:

None

Notes

The payload exists because the pattern of rejections often diagnoses the run better than any single one of them does. A caller that wants to say “these all failed the same way, and here is what that usually means” should not have to re-derive the pool or, worse, parse it back out of the message it just formatted.

exception rbfenetmap.core.exceptions.ExporterError

Bases: RBFENetworkMapError

A network could not be serialized for a downstream program.

exception rbfenetmap.core.exceptions.PluginError

Bases: RBFENetworkMapError

A plugin is unknown, duplicated, or its backend is not importable.