ciri_toolkit/
circ_catalog.rs

1//! External circRNA catalog parsing for multi-sample second-pass runs.
2//!
3//! The catalog is intentionally narrower than CIRI `.out`: `ciri-merge` writes
4//! BED6 sites that define cohort-level BSJ candidate boundaries, while each
5//! sample still rescans its own BAM/SAM to decide whether those candidates have
6//! read support. Keeping this as a typed module avoids leaking BED coordinate
7//! conventions into Scan2's Java-compatible candidate payload.
8
9use anyhow::{bail, Context, Result};
10use std::collections::HashSet;
11use std::fs::File;
12use std::io::{BufRead, BufReader};
13
14/// One external BSJ candidate from a BED6 circRNA catalog.
15#[derive(Debug, Clone, Eq, PartialEq, Hash)]
16pub struct CircCatalogRecord {
17    /// Reference sequence name.
18    pub chrom: String,
19    /// CIRI-style 1-based BSJ start coordinate.
20    pub start: i32,
21    /// CIRI-style 1-based BSJ end coordinate.
22    pub end: i32,
23    /// Strand label carried by BED6 (`+`, `-`, or `NA` when unknown).
24    pub strand: String,
25}
26
27/// Reads a BED6 circRNA catalog produced by `ciri-merge`.
28///
29/// BED uses 0-based half-open coordinates, while Scan2 candidate records use the
30/// same 1-based coordinates as CIRI `.out`. The conversion is performed exactly
31/// once here so downstream code can work in CIRI coordinates.
32pub fn read_circ_bed6(path: &str) -> Result<Vec<CircCatalogRecord>> {
33    let file = File::open(path).with_context(|| format!("open circ catalog {path}"))?;
34    let reader = BufReader::with_capacity(1024 * 1024, file);
35    let mut records = Vec::new();
36    let mut seen = HashSet::new();
37    for (line_idx, line_result) in reader.lines().enumerate() {
38        let line = line_result?;
39        let trimmed = line.trim();
40        if trimmed.is_empty()
41            || trimmed.starts_with('#')
42            || trimmed.starts_with("track")
43            || trimmed.starts_with("browser")
44        {
45            continue;
46        }
47        let fields: Vec<&str> = trimmed.split_whitespace().collect();
48        if fields.len() < 6 {
49            bail!(
50                "{}:{} has {} columns, expected BED6",
51                path,
52                line_idx + 1,
53                fields.len()
54            );
55        }
56        let chrom_start = fields[1]
57            .parse::<i64>()
58            .with_context(|| format!("parse chromStart at {}:{}", path, line_idx + 1))?;
59        let chrom_end = fields[2]
60            .parse::<i64>()
61            .with_context(|| format!("parse chromEnd at {}:{}", path, line_idx + 1))?;
62        if chrom_start < 0 || chrom_end <= chrom_start {
63            bail!(
64                "{}:{} has invalid BED interval {}-{}",
65                path,
66                line_idx + 1,
67                chrom_start,
68                chrom_end
69            );
70        }
71        let start = chrom_start + 1;
72        if start > i32::MAX as i64 || chrom_end > i32::MAX as i64 {
73            bail!(
74                "{}:{} exceeds supported i32 coordinate range",
75                path,
76                line_idx + 1
77            );
78        }
79        let record = CircCatalogRecord {
80            chrom: fields[0].to_string(),
81            start: start as i32,
82            end: chrom_end as i32,
83            strand: normalize_bed_strand(fields[5]).to_string(),
84        };
85        if seen.insert((
86            record.chrom.clone(),
87            record.start,
88            record.end,
89            record.strand.clone(),
90        )) {
91            records.push(record);
92        }
93    }
94    Ok(records)
95}
96
97/// Normalizes BED strand labels without inventing biological direction.
98fn normalize_bed_strand(raw: &str) -> &str {
99    match raw {
100        "+" | "-" => raw,
101        _ => "NA",
102    }
103}
104
105#[cfg(test)]
106mod tests {
107    use super::*;
108    use std::io::Write;
109    use tempfile::NamedTempFile;
110
111    #[test]
112    fn bed6_catalog_converts_to_ciri_coordinates() {
113        let mut file = NamedTempFile::new().unwrap();
114        writeln!(file, "chr1\t9\t20\tchr1:10|20:+\t7\t+").unwrap();
115        writeln!(file, "chr1\t9\t20\tduplicate\t7\t+").unwrap();
116        writeln!(file, "chr2\t29\t40\tunknown\t3\t.").unwrap();
117
118        let records = read_circ_bed6(file.path().to_str().unwrap()).unwrap();
119
120        assert_eq!(records.len(), 2);
121        assert_eq!(
122            records[0],
123            CircCatalogRecord {
124                chrom: "chr1".to_string(),
125                start: 10,
126                end: 20,
127                strand: "+".to_string(),
128            }
129        );
130        assert_eq!(records[1].strand, "NA");
131    }
132}
ciri_toolkit/circ_catalog.rs

ciri_toolkit/
circ_catalog.rs
