data_beans/aux/
ontology.rs

1//! Generic OBO ontology (`cl-basic.obo`, `go-basic.obo`, …) parsed with
2//! `fastobo` and stored as a `petgraph` directed graph (edges point
3//! **child → parent**, tagged with the relation kind).
4//!
5//! This is a thin, reusable wrapper: parse once, then query each term's
6//! ancestors via [`Ontology::ancestors_or_self`] (`is_a` only — the relation a
7//! collapsed-label tree needs) or [`Ontology::ancestors_or_self_with_part_of`]
8//! (`is_a` + `part_of`, the GO "true-path" closure); [`Ontology::edges`] hands
9//! the whole hierarchy to callers that embed or draw it rather than walk it.
10//!
11//! Prefix-agnostic: any term id is kept (`CL:`, `GO:`, …); edges to
12//! unknown/obsolete targets are dropped, as are obsolete terms. The
13//! `{is_inferred="true"}` qualifier blocks fastobo handles natively (a hand
14//! rolled parser missing them was the original prototype's bug).
15
16use anyhow::{Context, Result};
17use fastobo::ast::{EntityFrame, TermClause};
18use petgraph::graph::{DiGraph, NodeIndex};
19use petgraph::visit::EdgeRef;
20use petgraph::Direction::Outgoing;
21use rustc_hash::{FxHashMap, FxHashSet};
22
23/// Ontology relation kind carried on each `child → parent` edge.
24#[derive(Debug, Clone, Copy, PartialEq, Eq)]
25pub enum Rel {
26    /// `is_a` subsumption (subclass).
27    IsA,
28    /// `relationship: part_of` (mereological containment; GO true-path).
29    PartOf,
30}
31
32/// Parsed OBO DAG. Edges are `child → parent` (tagged [`Rel`]), so a node's
33/// ancestors are reachable along `Outgoing` edges and its children along
34/// `Incoming` edges.
35pub struct Ontology {
36    graph: DiGraph<Box<str>, Rel>,
37    idx: FxHashMap<Box<str>, NodeIndex>,
38    names: FxHashMap<Box<str>, Box<str>>,
39    /// `def:` text per term, unescaped, without the xref list.
40    defs: FxHashMap<Box<str>, Box<str>>,
41    /// `namespace:` per term (`biological_process`, …), the header's
42    /// `default-namespace:` where a term names none; obsolete terms too.
43    namespaces: FxHashMap<Box<str>, Box<str>>,
44}
45
46/// One parsed, non-obsolete OBO term gathered in the first pass.
47struct ParsedTerm {
48    id: Box<str>,
49    name: Option<Box<str>>,
50    def: Option<Box<str>>,
51    /// `is_a` parent ids.
52    is_a: Vec<Box<str>>,
53    /// `relationship: part_of` parent ids.
54    part_of: Vec<Box<str>>,
55}
56
57impl Ontology {
58    /// Parse an OBO file into the `is_a` + `part_of` DAG. Obsolete terms and
59    /// edges to unknown/obsolete targets are skipped.
60    pub fn load_obo(path: &str) -> Result<Self> {
61        let doc = fastobo::from_file(path)
62            .with_context(|| format!("failed to parse OBO file: {path}"))?;
63
64        // The header's `default-namespace:`, for terms that name none.
65        let default_namespace: Option<Box<str>> = doc.header().iter().find_map(|c| match c {
66            fastobo::ast::HeaderClause::DefaultNamespace(ns) => Some(ns.to_string().trim().into()),
67            _ => None,
68        });
69        // Pass 1: collect non-obsolete terms; every term's namespace,
70        // obsolete ones' too (an old annotation may still name them).
71        let mut terms: Vec<ParsedTerm> = Vec::new();
72        let mut namespaces: FxHashMap<Box<str>, Box<str>> = FxHashMap::default();
73        for frame in doc.entities() {
74            let EntityFrame::Term(term) = frame else {
75                continue;
76            };
77            // fastobo's `Line<_>` Display can carry a trailing newline / inline
78            // comment, so trim every extracted token down to the bare value.
79            let id: Box<str> = term.id().to_string().trim().into();
80            let mut name: Option<Box<str>> = None;
81            let mut def: Option<Box<str>> = None;
82            let mut namespace: Option<Box<str>> = None;
83            let mut is_a: Vec<Box<str>> = Vec::new();
84            let mut part_of: Vec<Box<str>> = Vec::new();
85            let mut obsolete = false;
86            for line in term.clauses() {
87                match &**line {
88                    TermClause::Name(n) => name = Some(n.to_string().trim().into()),
89                    TermClause::Def(d) => {
90                        let text = d.text().as_str().trim();
91                        if !text.is_empty() {
92                            def = Some(text.into());
93                        }
94                    }
95                    TermClause::Namespace(ns) => {
96                        namespace = Some(ns.to_string().trim().into());
97                    }
98                    TermClause::IsObsolete(b) => obsolete = obsolete || *b,
99                    TermClause::IsA(parent) => is_a.push(parent.to_string().trim().into()),
100                    TermClause::Relationship(rel, target)
101                        if rel.to_string().trim() == "part_of" =>
102                    {
103                        part_of.push(target.to_string().trim().into());
104                    }
105                    _ => {}
106                }
107            }
108            if let Some(ns) = namespace.or_else(|| default_namespace.clone()) {
109                namespaces.insert(id.clone(), ns);
110            }
111            if !obsolete {
112                terms.push(ParsedTerm {
113                    id,
114                    name,
115                    def,
116                    is_a,
117                    part_of,
118                });
119            }
120        }
121
122        // Pass 2: nodes first, then edges (only between known nodes).
123        let mut graph: DiGraph<Box<str>, Rel> = DiGraph::new();
124        let mut idx: FxHashMap<Box<str>, NodeIndex> = FxHashMap::default();
125        let mut names: FxHashMap<Box<str>, Box<str>> = FxHashMap::default();
126        let mut defs: FxHashMap<Box<str>, Box<str>> = FxHashMap::default();
127        for term in &terms {
128            let node = graph.add_node(term.id.clone());
129            idx.insert(term.id.clone(), node);
130            if let Some(n) = &term.name {
131                names.insert(term.id.clone(), n.clone());
132            }
133            if let Some(d) = &term.def {
134                defs.insert(term.id.clone(), d.clone());
135            }
136        }
137        for term in &terms {
138            let child = idx[&term.id];
139            for (parents, rel) in [(&term.is_a, Rel::IsA), (&term.part_of, Rel::PartOf)] {
140                for p in parents {
141                    if let Some(&parent) = idx.get(p) {
142                        graph.add_edge(child, parent, rel);
143                    }
144                }
145            }
146        }
147
148        Ok(Self {
149            graph,
150            idx,
151            names,
152            defs,
153            namespaces,
154        })
155    }
156
157    /// Every (non-obsolete) term id, in no particular order.
158    pub fn ids(&self) -> impl Iterator<Item = &str> + '_ {
159        self.idx.keys().map(|k| &**k)
160    }
161
162    /// The term's `def:` text (`None` if the term is unknown or undefined).
163    #[must_use]
164    pub fn def(&self, id: &str) -> Option<&str> {
165        self.defs.get(id).map(|d| &**d)
166    }
167
168    /// Every `child → parent` edge with its relation kind, in no particular
169    /// order — the hierarchy as a graph, for callers that embed or draw it
170    /// rather than walk it.
171    pub fn edges(&self) -> impl Iterator<Item = (&str, &str, Rel)> + '_ {
172        self.graph.edge_references().map(|e| {
173            (
174                &*self.graph[e.source()],
175                &*self.graph[e.target()],
176                *e.weight(),
177            )
178        })
179    }
180
181    /// Number of (non-obsolete) terms.
182    #[must_use]
183    pub fn len(&self) -> usize {
184        self.graph.node_count()
185    }
186
187    #[must_use]
188    pub fn is_empty(&self) -> bool {
189        self.graph.node_count() == 0
190    }
191
192    #[must_use]
193    pub fn contains(&self, id: &str) -> bool {
194        self.idx.contains_key(id)
195    }
196
197    /// Human-readable term name (`None` if the term is unknown or unnamed).
198    #[must_use]
199    pub fn name(&self, id: &str) -> Option<&str> {
200        self.names.get(id).map(|n| &**n)
201    }
202
203    /// The `namespace:` of `id` (GO's `biological_process`,
204    /// `molecular_function` or `cellular_component`), else the header's
205    /// `default-namespace:`; obsolete terms keep theirs.
206    #[must_use]
207    pub fn namespace(&self, id: &str) -> Option<&str> {
208        self.namespaces.get(id).map(|n| &**n)
209    }
210
211    /// All `is_a` ancestors of `id` plus `id` itself (empty if `id` is unknown).
212    #[must_use]
213    pub fn ancestors_or_self(&self, id: &str) -> FxHashSet<Box<str>> {
214        self.ancestors_impl(id, false)
215    }
216
217    /// All `is_a` + `part_of` ancestors of `id` plus `id` itself — the GO
218    /// "true-path" closure (empty if `id` is unknown).
219    #[must_use]
220    pub fn ancestors_or_self_with_part_of(&self, id: &str) -> FxHashSet<Box<str>> {
221        self.ancestors_impl(id, true)
222    }
223
224    /// Walk ancestors in `NodeIndex` space (Copy — no per-visit string clones),
225    /// following `is_a` edges and, when `with_part_of`, `part_of` edges too;
226    /// then materialize each id exactly once at the boundary.
227    fn ancestors_impl(&self, id: &str, with_part_of: bool) -> FxHashSet<Box<str>> {
228        let Some(&start) = self.idx.get(id) else {
229            return FxHashSet::default();
230        };
231        let mut seen: FxHashSet<NodeIndex> = FxHashSet::default();
232        seen.insert(start);
233        let mut stack = vec![start];
234        while let Some(n) = stack.pop() {
235            for edge in self.graph.edges_directed(n, Outgoing) {
236                let follow = matches!(edge.weight(), Rel::IsA)
237                    || (with_part_of && matches!(edge.weight(), Rel::PartOf));
238                if follow && seen.insert(edge.target()) {
239                    stack.push(edge.target());
240                }
241            }
242        }
243        seen.iter().map(|&n| self.graph[n].clone()).collect()
244    }
245}
246
247#[cfg(test)]
248mod tests {
249    use super::*;
250    use std::io::Write;
251
252    /// Minimal CL-shaped OBO: cell → lymphocyte → T → {CD4, CD8}; B sibling;
253    /// one obsolete term; one `{is_inferred="true"}` qualifier on an is_a line;
254    /// one `part_of` edge (CD8 T part_of an immune-system "compartment") to
255    /// exercise the relation-filtered traversal.
256    fn write_obo() -> tempfile::NamedTempFile {
257        let mut f = tempfile::NamedTempFile::new().unwrap();
258        writeln!(
259            f,
260            "format-version: 1.2\ndefault-namespace: cell_type\n\n\
261             [Term]\nid: CL:0000000\nname: cell\nnamespace: cell\n\n\
262             [Term]\nid: CL:0000542\nname: lymphocyte\nis_a: CL:0000000 ! cell\n\n\
263             [Term]\nid: CL:0000084\nname: T cell\ndef: \"A lymphocyte with a \\\"TCR\\\", made in the thymus.\" [GOC:add]\nis_a: CL:0000542 {{is_inferred=\"true\"}} ! lymphocyte\n\n\
264             [Term]\nid: CL:0000624\nname: CD4 T\nis_a: CL:0000084 ! T cell\n\n\
265             [Term]\nid: CL:0000625\nname: CD8 T\nis_a: CL:0000084 ! T cell\nrelationship: part_of CL:1000000 ! compartment\n\n\
266             [Term]\nid: CL:1000000\nname: immune compartment\n\n\
267             [Term]\nid: CL:0000236\nname: B cell\nis_a: CL:0000542 ! lymphocyte\n\n\
268             [Term]\nid: CL:9999999\nname: dead\nis_obsolete: true\n"
269        )
270        .unwrap();
271        f.flush().unwrap();
272        f
273    }
274
275    #[test]
276    fn parses_and_resolves_ancestry() {
277        let f = write_obo();
278        let onto = Ontology::load_obo(f.path().to_str().unwrap()).unwrap();
279
280        // obsolete term dropped; 7 live terms (incl. the compartment).
281        assert_eq!(onto.len(), 7);
282        assert!(!onto.contains("CL:9999999"));
283        assert_eq!(onto.name("CL:0000084"), Some("T cell"));
284        assert_eq!(onto.namespace("CL:0000000"), Some("cell"));
285        assert_eq!(onto.namespace("CL:0000084"), Some("cell_type"));
286        assert_eq!(onto.namespace("CL:9999999"), Some("cell_type"), "obsolete");
287
288        // The {is_inferred} qualifier must NOT break the edge: CD4/CD8 under T,
289        // T under lymphocyte under cell.
290        let anc = onto.ancestors_or_self("CL:0000624");
291        for a in ["CL:0000624", "CL:0000084", "CL:0000542", "CL:0000000"] {
292            assert!(anc.contains(a), "missing ancestor {a}");
293        }
294        assert!(!anc.contains("CL:0000236"));
295    }
296
297    #[test]
298    fn definitions_are_kept_unescaped_and_the_hierarchy_is_exposed_as_edges() {
299        let f = write_obo();
300        let onto = Ontology::load_obo(f.path().to_str().unwrap()).unwrap();
301        assert_eq!(
302            onto.def("CL:0000084"),
303            Some("A lymphocyte with a \"TCR\", made in the thymus.")
304        );
305        assert_eq!(onto.def("CL:0000000"), None, "no def: line");
306        assert_eq!(onto.def("CL:9999999"), None, "obsolete");
307        let mut edges: Vec<(String, String, Rel)> = onto
308            .edges()
309            .map(|(c, p, r)| (c.to_string(), p.to_string(), r))
310            .collect();
311        edges.sort_by(|a, b| (&a.0, &a.1).cmp(&(&b.0, &b.1)));
312        assert_eq!(edges.len(), 6, "5 is_a + 1 part_of among live terms");
313        assert!(edges.contains(&("CL:0000625".into(), "CL:1000000".into(), Rel::PartOf)));
314        assert!(edges.contains(&("CL:0000084".into(), "CL:0000542".into(), Rel::IsA)));
315        assert!(edges
316            .iter()
317            .all(|(c, p, _)| onto.contains(c) && onto.contains(p)));
318    }
319
320    #[test]
321    fn part_of_only_followed_on_demand() {
322        let f = write_obo();
323        let onto = Ontology::load_obo(f.path().to_str().unwrap()).unwrap();
324
325        // is_a-only walk does NOT cross the part_of edge to the compartment.
326        let isa = onto.ancestors_or_self("CL:0000625");
327        assert!(isa.contains("CL:0000084"), "is_a ancestor missing");
328        assert!(
329            !isa.contains("CL:1000000"),
330            "part_of must not leak into is_a-only walk"
331        );
332
333        // is_a + part_of walk reaches the compartment.
334        let full = onto.ancestors_or_self_with_part_of("CL:0000625");
335        assert!(full.contains("CL:0000084"), "is_a ancestor missing");
336        assert!(full.contains("CL:1000000"), "part_of ancestor missing");
337    }
338}
data_beans/aux/ontology.rs

data_beans/aux/
ontology.rs
