From ad0e149b688df9bc86464247d58510051ad925b3 Mon Sep 17 00:00:00 2001 From: sfluegel Date: Wed, 22 Jul 2026 18:00:03 +0200 Subject: [PATCH 1/4] update docs --- README.md | 55 ++++++++++- chebi_utils/extract_properties.py | 154 ++++++++++++++++++++++++++---- chebi_utils/obo_extractor.py | 18 +++- chebi_utils/sample_filters.py | 68 +++++++++---- chebi_utils/sdf_extractor.py | 9 +- 5 files changed, 262 insertions(+), 42 deletions(-) diff --git a/README.md b/README.md index 1b1a532..e79034d 100644 --- a/README.md +++ b/README.md @@ -1,6 +1,12 @@ # python-chebi-utils -Common processing functionality for the ChEBI ontology — download versioned data files, build an ontology graph, extract molecules, assemble labeled datasets, and generate stratified train/val/test splits. +Common processing functionality for the ChEBI ontology — download versioned data files, build an ontology graph, extract molecules, assemble labeled datasets, generate stratified train/validation/test splits, extract first-order-logic molecular properties, and select hierarchy-aware sample subsets. + +> **⚠️ Breaking change in v0.3** +> +> `create_multilabel_splits` now returns the validation split under the key +> `"validation"` instead of `"val"`. Update any code that reads `splits["val"]` +> to use `splits["validation"]`. ## Installation @@ -86,7 +92,7 @@ from chebi_utils import create_multilabel_splits splits = create_multilabel_splits(dataset, train_ratio=0.8, val_ratio=0.1, test_ratio=0.1) train_df = splits["train"] -val_df = splits["val"] +val_df = splits["validation"] # renamed from "val" in v0.3 test_df = splits["test"] ``` @@ -96,6 +102,51 @@ present, `MultilabelStratifiedShuffleSplit` from the `iterative-stratification` package is used; for a single label column, `StratifiedShuffleSplit` from scikit-learn is used. +### Extract molecular properties as first-order-logic facts + +```python +from chebi_utils.extract_properties import mol_to_fol_atoms, get_numerical_facts + +atom_facts, mol_facts = mol_to_fol_atoms(mol, with_rings=True, with_steroids=True) +# atom_facts — dict[str, list] of predicates over atom indices: +# unary (e.g. "c", "charge_p", "has_2_hs", "cip_code_R", "in_ring6", "steroid_3") +# → list[int] of atom indices +# binary (e.g. "has_bond_to", "bSINGLE", "ring6") → list[tuple[int, ...]] +# mol_facts — set[str] of molecule-level predicates that hold for the whole +# molecule (e.g. "net_charge_neutral", "aromatic") + +numerical_facts = get_numerical_facts(mol) +# {"mol_weight": [], "ring_size": [, …]} +``` + +Turns an RDKit `Mol` into a symbolic model suitable for building FOL structures +for reasoning tasks. Facts cover per-atom element, formal charge, hydrogen +counts, and CIP chirality; symmetric bond and bond-stereo relations; ring +membership up to `MAX_RING_SIZE` (8); and steroid-nucleus positions +(`steroid_1` … `steroid_17`) matched against the gonane core via IUPAC steroid +numbering. Ring and steroid extraction can be toggled with `with_rings` and +`with_steroids`. + +### Select hierarchy-aware sample subsets + +```python +from chebi_utils.sample_filters import get_closest_negatives, get_direct_neighbors + +# Nearest negatives: samples that are NOT subclasses of the target but close to +# it in the ontology, expanding outward until min_samples (up to max_samples) is met. +negatives = get_closest_negatives( + samples, graph, target_id="15841", min_samples=25, max_samples=None +) + +# Split samples into positives (descendants of the target) and "direct neighbor" +# negatives (descendants of ALL direct parents of the target, but not the target). +pos_ids, neg_ids = get_direct_neighbors(samples, graph, target_id="15841") +``` + +Useful for constructing balanced positive/negative sets for a given ChEBI class +by leveraging the `is_a` hierarchy. `samples` is a list of ChEBI IDs (as +strings) and `graph` is a graph from `build_chebi_graph`. + ## Running Tests ```bash diff --git a/chebi_utils/extract_properties.py b/chebi_utils/extract_properties.py index 7d2cc84..d2d09c5 100644 --- a/chebi_utils/extract_properties.py +++ b/chebi_utils/extract_properties.py @@ -1,5 +1,8 @@ -# extract basic (and not so basic) properties from molecules. This is used to construct -# FOL structures for reasoning tasks on molecules. +"""Extract basic (and not so basic) properties from molecules. + +These are used to construct first-order-logic (FOL) structures for reasoning +tasks on molecules. +""" import logging @@ -31,12 +34,28 @@ def mol_to_fol_atoms( ) -> tuple[dict[str, list], set[str]]: """Convert an RDKit ``Mol`` into a first-order logic model at the atom level. - Returns ``(atom_extensions, mol_extensions)`` where: - - ``atom_extensions`` is a ``dict[str, list]``: unary predicates map to - ``list[int]`` of atom indices; binary predicates map to - ``list[tuple[int, int]]`` of (left, right) index pairs. - - ``mol_extensions`` is a ``set[str]`` of molecule-level predicate names - that hold for this molecule (e.g. ``net_charge_positive``). + Aggregates atom, bond, ring, and steroid predicates into a single atom-level + extension mapping plus the set of molecule-level predicates. + + Parameters + ---------- + mol : Chem.Mol + The molecule to convert. + with_rings : bool + If ``True`` (default), include ring predicates from :func:`get_rings`. + with_steroids : bool + If ``True`` (default), include steroid-nucleus position predicates from + :func:`get_steroid_positions`. + + Returns + ------- + tuple[dict[str, list], set[str]] + ``(atom_extensions, mol_extensions)`` where: + - ``atom_extensions`` is a ``dict[str, list]``: unary predicates map to + ``list[int]`` of atom indices; binary predicates map to + ``list[tuple[int, int]]`` of (left, right) index pairs. + - ``mol_extensions`` is a ``set[str]`` of molecule-level predicate names + that hold for this molecule (e.g. ``net_charge_positive``). """ atom_extensions: dict[str, list] = {} @@ -55,6 +74,25 @@ def mol_to_fol_atoms( def get_atom_properties(mol: Chem.Mol) -> dict[str, list]: + """Extract per-atom unary predicates. + + For each atom, emits predicates for its element symbol, formal charge + (both a sign predicate ``charge_p``/``charge_n``/``charge0`` and an exact + predicate ``charge{n}``), total hydrogen count (exact ``has_{n}_hs`` and + cumulative ``has_at_least_{n}_hs``), and CIP chirality label + (``cip_code_R``/``cip_code_S``) when assigned. + + Parameters + ---------- + mol : Chem.Mol + The molecule to inspect. + + Returns + ------- + dict[str, list] + Mapping from predicate name to the ``list[int]`` of atom indices for + which the predicate holds. + """ try: Chem.rdCIPLabeler.AssignCIPLabels(mol) except Exception as e: @@ -65,7 +103,6 @@ def get_atom_properties(mol: Chem.Mol) -> dict[str, list]: atom_extensions: dict[str, list] = {} - # For each atom: element symbol, charge, hydrogen counts, chirality for atom in mol.GetAtoms(): atom_idx = atom.GetIdx() atom_symbol = atom.GetSymbol().lower() @@ -93,7 +130,24 @@ def get_atom_properties(mol: Chem.Mol) -> dict[str, list]: def get_bond_properties(mol: Chem.Mol) -> dict[str, list]: - # Bond predicates (symmetric) + """Extract symmetric per-bond binary predicates. + + For each bond, emits a bond-type predicate (e.g. ``bSINGLE``), a generic + ``has_bond_to`` predicate, and a bond-stereo predicate (e.g. + ``bSTEREOE``) when the bond carries stereochemistry. Every predicate is + stored symmetrically, i.e. both ``(left, right)`` and ``(right, left)``. + + Parameters + ---------- + mol : Chem.Mol + The molecule to inspect. + + Returns + ------- + dict[str, list] + Mapping from predicate name to a ``list[tuple[int, int]]`` of ordered + atom-index pairs for which the predicate holds. + """ atom_extensions: dict[str, list] = {} for bond in mol.GetBonds(): left = bond.GetBeginAtomIdx() @@ -110,7 +164,23 @@ def get_bond_properties(mol: Chem.Mol) -> dict[str, list]: def get_molecule_level_properties(mol: Chem.Mol) -> set[str]: - # Molecule-level (global) properties (either true or false for the whole molecule) + """Extract global (molecule-level) predicates. + + Emits exactly one net-charge predicate + (``net_charge_positive``/``net_charge_negative``/``net_charge_neutral``) and + one aromaticity predicate (``aromatic`` if the molecule has at least one + aromatic atom, otherwise ``aliphatic``). + + Parameters + ---------- + mol : Chem.Mol + The molecule to inspect. + + Returns + ------- + set[str] + The set of molecule-level predicate names that hold for this molecule. + """ mol_extensions: set[str] = set() net_charge = Chem.GetFormalCharge(mol) if net_charge > 0: @@ -119,7 +189,6 @@ def get_molecule_level_properties(mol: Chem.Mol) -> set[str]: mol_extensions.add("net_charge_negative") else: mol_extensions.add("net_charge_neutral") - # aliphatic vs aromatic (defined as having at least one aromatic atom) if len(list(mol.GetAromaticAtoms())) > 0: mol_extensions.add("aromatic") else: @@ -128,11 +197,28 @@ def get_molecule_level_properties(mol: Chem.Mol) -> set[str]: def get_rings(mol: Chem.Mol) -> dict[str, list]: - # Rings have two predicates. One for the atom-ring relation and one for the ring itself - # ring{N}(A1, …, AN) – A1…AN form an N-membered ring (all permutations) - # Only for N <= MAX_RING_SIZE. - # in_ring{N}(A) – A belongs to some N-membered ring (N <= MAX_RING_SIZE). - # in_ring(A) – A belongs to some ring of any size. + """Extract ring-membership predicates. + + Emits three kinds of predicate: + + - ``ring{N}(A1, …, AN)`` – ``A1…AN`` form an N-membered ring, stored for all + rotations in both directions. Only for ``N <= MAX_RING_SIZE``. + - ``in_ring{N}(A)`` – ``A`` belongs to some N-membered ring + (``N <= MAX_RING_SIZE``). + - ``in_ring(A)`` – ``A`` belongs to some ring of any size. + + Parameters + ---------- + mol : Chem.Mol + The molecule to inspect. + + Returns + ------- + dict[str, list] + Mapping from predicate name to its extension: ``ring{N}`` maps to a + ``list[tuple[int, ...]]`` of atom-index tuples; ``in_ring`` and + ``in_ring{N}`` map to a sorted ``list[int]`` of atom indices. + """ atom_extensions: dict[str, list] = {} in_ring_atoms: set[int] = set() in_ringN_atoms: dict[int, set[int]] = {} @@ -153,7 +239,23 @@ def get_rings(mol: Chem.Mol) -> dict[str, list]: def get_steroid_positions(mol: Chem.Mol) -> dict[str, list]: - # Steroid nucleus positions (steroid_1 … steroid_17) + """Extract steroid-nucleus position predicates. + + Matches the molecule against the gonane core and, on a match, labels the + ring atoms with their IUPAC steroid position as predicates ``steroid_1`` … + ``steroid_17``. Molecules without a gonane core yield no predicates. + + Parameters + ---------- + mol : Chem.Mol + The molecule to inspect. + + Returns + ------- + dict[str, list] + Mapping from ``steroid_{position}`` to the ``list[int]`` of matched atom + indices. Empty when the molecule has no steroid nucleus. + """ atom_extensions: dict[str, list] = {} steroid_match = mol.GetSubstructMatch(_GONANE_PATTERN, useChirality=False) if steroid_match: @@ -167,8 +269,20 @@ def get_steroid_positions(mol: Chem.Mol) -> dict[str, list]: def get_numerical_facts(mol: Chem.Mol) -> dict[str, list]: """Molecular weight and ring size as numerical values. - Expresses "this molecule has weight ..." and "this molecule has a ring of size ..." - as molecule-integer value relations. + Expresses "this molecule has weight ..." and "this molecule has a ring of + size ..." as molecule-integer value relations. + + Parameters + ---------- + mol : Chem.Mol + The molecule to inspect. + + Returns + ------- + dict[str, list] + ``{"mol_weight": [], "ring_size": [, …]}``. + The ``ring_size`` key is present only when the molecule has at least one + ring. """ atom_extensions: dict[str, list] = {} atom_extensions["mol_weight"] = [round(Descriptors.MolWt(mol))] diff --git a/chebi_utils/obo_extractor.py b/chebi_utils/obo_extractor.py index e893ec5..f804ca8 100644 --- a/chebi_utils/obo_extractor.py +++ b/chebi_utils/obo_extractor.py @@ -141,8 +141,22 @@ def build_chebi_graph(filepath: str | Path, top_class: str | None = "23367") -> def get_hierarchy_subgraph(chebi_graph: nx.DiGraph) -> nx.DiGraph: - """Subgraph of ChEBI including only edges corresponding to hierarchical relations (is_a). - Also removes nodes that are not connected by any is_a edges to other nodes.""" + """Extract the ``is_a`` hierarchy of a ChEBI graph as a subgraph. + + Keeps only edges whose ``relation`` attribute is ``is_a``, and drops nodes + that are not connected to any other node by an ``is_a`` edge. + + Parameters + ---------- + chebi_graph : nx.DiGraph + Full ChEBI ontology graph from :func:`build_chebi_graph`. + + Returns + ------- + nx.DiGraph + Subgraph containing only ``is_a`` edges (child → parent) and the nodes + they connect. + """ return nx.DiGraph( chebi_graph.edge_subgraph( (u, v) for u, v, d in chebi_graph.edges(data=True) if d.get("relation") == "is_a" diff --git a/chebi_utils/sample_filters.py b/chebi_utils/sample_filters.py index d5d6468..1e66853 100644 --- a/chebi_utils/sample_filters.py +++ b/chebi_utils/sample_filters.py @@ -1,4 +1,6 @@ -# functionality for selecting specific sample subsets from the ChEBI dataset +"""Select specific sample subsets from the ChEBI dataset using the hierarchy.""" + +import queue import networkx as nx @@ -8,13 +10,36 @@ def get_closest_negatives( samples: list[str], chebi_graph: nx.DiGraph, target_id: str, min_samples=25, max_samples=None ) -> set[str]: - # from the list of samples, find those that are not subclasses of the target_id, but close - # to it in the hierarchy. - # goal: reach min_samples, but continue collecting samples (until max_samples) if they are - # siblings. + """Find samples close to a target class in the hierarchy but not below it. + + Performs a breadth-first walk outward from ``target_id`` over the undirected + ``is_a`` hierarchy, collecting samples that are subclasses of the visited + neighbours (but not of ``target_id`` itself). The walk aims to reach + ``min_samples``, then keeps collecting further-out (non-sibling) samples only + until ``max_samples`` is reached. + + Parameters + ---------- + samples : list[str] + Candidate ChEBI IDs (as strings) to select from. + chebi_graph : nx.DiGraph + Full ChEBI ontology graph from :func:`build_chebi_graph`. + target_id : str + ChEBI ID whose neighbourhood is searched for negatives. + min_samples : int + Target number of samples to collect before stopping expansion + (default 25). + max_samples : int or None + Hard upper bound on the number of samples collected. ``None`` (default) + means no upper bound. + + Returns + ------- + set[str] + The selected ChEBI IDs. + """ hierarchy_graph = nx.transitive_closure_dag(get_hierarchy_subgraph(chebi_graph)) undirected_graph = get_hierarchy_subgraph(chebi_graph).to_undirected() - import queue q = queue.Queue() q.put(target_id) @@ -48,16 +73,27 @@ def get_direct_neighbors( chebi_graph: nx.DiGraph, target_id: str, ) -> tuple[list[str], list[str]]: - """ - Filter samples and sort into two groups: - positive: sample is a descendant of the target_id - negative: sample is not a descendant of the target_id, but a "direct neighbor" -> a - descendant of all direct parents of the target_id. - - Returns: - pos_ids: list of positive validation molecule IDs - neg_ids: list of negative validation molecule IDs - (empty when target has no siblings) + """Filter samples and sort them into positives and direct-neighbor negatives. + + A sample is *positive* when it is a descendant of ``target_id``, and a + *negative* when it is not a descendant of ``target_id`` but is a "direct + neighbor" — a descendant of all direct parents of ``target_id``. + + Parameters + ---------- + samples : list[str] + Candidate ChEBI IDs (as strings) to filter. + chebi_graph : nx.DiGraph + Full ChEBI ontology graph from :func:`build_chebi_graph`. + target_id : str + ChEBI ID defining the positive class. + + Returns + ------- + tuple[list[str], list[str]] + ``(pos_ids, neg_ids)`` where ``pos_ids`` are positive molecule IDs and + ``neg_ids`` are negative molecule IDs (empty when the target has no + siblings). """ hierarchy_graph = nx.transitive_closure_dag(get_hierarchy_subgraph(chebi_graph)) pos_ids = [str(d) for d in hierarchy_graph.predecessors(target_id) if str(d) in samples] diff --git a/chebi_utils/sdf_extractor.py b/chebi_utils/sdf_extractor.py index 98607e5..22aae48 100644 --- a/chebi_utils/sdf_extractor.py +++ b/chebi_utils/sdf_extractor.py @@ -131,6 +131,13 @@ def extract_molecules(filepath: str | Path) -> pd.DataFrame: DataFrame with one row per molecule. Columns depend on the properties present in the file. Common columns (renamed for convenience): chebi_id, name, inchi, inchikey, smiles, formula, charge, mass, mol. + + Notes + ----- + Records that cannot be parsed into a molecule are dropped. Records that + parse into a valid molecule but contain no atoms are also dropped (e.g. + CHEBI:192499 in v251, cf. + https://github.com/ebi-chebi/ChEBI/issues/4915). """ rows = [] molblocks = [] @@ -160,9 +167,7 @@ def extract_molecules(filepath: str | Path) -> pd.DataFrame: df["mol"] = [_parse_molblock(mb, cid) for mb, cid in zip(molblocks, chebi_ids, strict=False)] df["chebi_id"] = df["chebi_id"].apply(_chebi_id_to_str) - # exclude records without a valid molecule df = df[df["mol"].notna()] - # some molecule records are valid, but have no atoms (e.g. CHEBI:192499 in v251, cf. https://github.com/ebi-chebi/ChEBI/issues/4915) df = df[[mol.GetNumAtoms() > 0 for mol in df["mol"]]] return df From b7d8333fe1a792f5f61b020f5f66d90252242947 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Simon=20Fl=C3=BCgel?= <43573433+sfluegel05@users.noreply.github.com> Date: Wed, 22 Jul 2026 18:10:02 +0200 Subject: [PATCH 2/4] Potential fix for pull request finding Co-authored-by: Copilot Autofix powered by AI <175728472+Copilot@users.noreply.github.com> --- chebi_utils/extract_properties.py | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/chebi_utils/extract_properties.py b/chebi_utils/extract_properties.py index d2d09c5..195b7c4 100644 --- a/chebi_utils/extract_properties.py +++ b/chebi_utils/extract_properties.py @@ -51,9 +51,9 @@ def mol_to_fol_atoms( ------- tuple[dict[str, list], set[str]] ``(atom_extensions, mol_extensions)`` where: - - ``atom_extensions`` is a ``dict[str, list]``: unary predicates map to - ``list[int]`` of atom indices; binary predicates map to - ``list[tuple[int, int]]`` of (left, right) index pairs. + - ``atom_extensions`` maps each predicate name to its extension: unary predicates map + to ``list[int]`` of atom indices; relation predicates map to ``list[tuple[int, ...]]`` + of atom-index tuples (arity depends on the predicate). - ``mol_extensions`` is a ``set[str]`` of molecule-level predicate names that hold for this molecule (e.g. ``net_charge_positive``). """ From ef19276f29a1464ed5b6f2507ac31732578f6425 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Simon=20Fl=C3=BCgel?= <43573433+sfluegel05@users.noreply.github.com> Date: Wed, 22 Jul 2026 18:10:46 +0200 Subject: [PATCH 3/4] Potential fix for pull request finding Co-authored-by: Copilot Autofix powered by AI <175728472+Copilot@users.noreply.github.com> --- chebi_utils/extract_properties.py | 10 +++++----- 1 file changed, 5 insertions(+), 5 deletions(-) diff --git a/chebi_utils/extract_properties.py b/chebi_utils/extract_properties.py index 195b7c4..c96f351 100644 --- a/chebi_utils/extract_properties.py +++ b/chebi_utils/extract_properties.py @@ -76,11 +76,11 @@ def mol_to_fol_atoms( def get_atom_properties(mol: Chem.Mol) -> dict[str, list]: """Extract per-atom unary predicates. - For each atom, emits predicates for its element symbol, formal charge - (both a sign predicate ``charge_p``/``charge_n``/``charge0`` and an exact - predicate ``charge{n}``), total hydrogen count (exact ``has_{n}_hs`` and - cumulative ``has_at_least_{n}_hs``), and CIP chirality label - (``cip_code_R``/``cip_code_S``) when assigned. + For each atom, emits predicates for its element symbol, formal charge (both a sign + predicate ``charge_p``/``charge_n``/``charge0`` and an exact predicate + ``charge{n}`` for positive charges or ``charge_m{n}`` for negative charges), total + hydrogen count (exact ``has_{n}_hs`` and cumulative ``has_at_least_{n}_hs``), and + CIP chirality label (``cip_code_R``/``cip_code_S``) when assigned. Parameters ---------- From 5c2d457861564ce5294eb2642da4ef5865b100db Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Simon=20Fl=C3=BCgel?= <43573433+sfluegel05@users.noreply.github.com> Date: Wed, 22 Jul 2026 18:11:35 +0200 Subject: [PATCH 4/4] Potential fix for pull request finding Co-authored-by: Copilot Autofix powered by AI <175728472+Copilot@users.noreply.github.com> --- chebi_utils/sample_filters.py | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/chebi_utils/sample_filters.py b/chebi_utils/sample_filters.py index 1e66853..7f6a01c 100644 --- a/chebi_utils/sample_filters.py +++ b/chebi_utils/sample_filters.py @@ -14,9 +14,11 @@ def get_closest_negatives( Performs a breadth-first walk outward from ``target_id`` over the undirected ``is_a`` hierarchy, collecting samples that are subclasses of the visited - neighbours (but not of ``target_id`` itself). The walk aims to reach - ``min_samples``, then keeps collecting further-out (non-sibling) samples only - until ``max_samples`` is reached. + neighbours (but not of ``target_id`` itself). The walk expands outward until + it reaches ``min_samples``. If ``min_samples`` is reached while still + processing the immediate neighbours of ``target_id``, that neighbour layer is + completed; otherwise expansion stops as soon as ``min_samples`` is reached. + If ``max_samples`` is set, selection stops once it is reached. Parameters ----------