gwseq-io 0.2.0

Rust library for processing bigWig, bigBed, BAM and HiC files
Documentation
//! Chromosome sizes by genome name.
//!
//! `bundled` (private, so not linked) is a generated table; anything else, and
//! any call with `full`, goes to `api.genome.ucsc.edu` and is memoised for the
//! process lifetime.
//!
//! The bundled table doubles as the place a case-insensitive name is resolved
//! to its canonical spelling, so `get_chr_sizes("SACCER3", true)` fetches
//! `sacCer3`. UCSC assembly names are case-sensitive and several are mixed
//! case, which is why lowercasing the argument alone will not do: that only
//! ever meets a key already lowercase.

use std::collections::HashMap;
use std::sync::OnceLock;

use parking_lot::Mutex;

use crate::error::Result;
use crate::genomic::ChrMap;

mod bundled;

/// Fetched assemblies, kept for the process lifetime.
/// Keyed by `(resolved name, full)` — the two answers differ and both are
/// worth holding.
fn fetch_cache() -> &'static Mutex<HashMap<(String, bool), ChrMap>> {
    static CACHE: OnceLock<Mutex<HashMap<(String, bool), ChrMap>>> = OnceLock::new();
    CACHE.get_or_init(|| Mutex::new(HashMap::new()))
}

/// Chromosome sizes for `genome`.
///
/// `full` includes unplaced contigs and alt loci, and always goes to the
/// network — the bundled tables hold the primary assembly only.
///
/// The result is ordered by chromosome id **as a string**, so `chr10` comes
/// before `chr2`, which is the order callers see in `chr_sizes`.
pub fn get_chr_sizes(genome: &str, full: bool) -> Result<ChrMap> {
    let bundled = bundled::get(genome);
    if !full {
        if let Some((_, sizes)) = bundled {
            return Ok(ChrMap::from_entries(
                sizes.iter().map(|(id, size)| ((*id).to_string(), *size)),
            ));
        }
    }
    // Canonical spelling when we know it, the caller's otherwise.
    let resolved = bundled.map(|(name, _)| name).unwrap_or(genome);
    fetch_cached(resolved, full).map_err(|error| annotate(genome, error))
}

/// The names the crate ships without a network.
pub fn bundled_genomes() -> &'static [&'static str] {
    bundled::GENOME_NAMES
}

/// What a failed fetch is, which depends on why it failed.
///
/// UCSC answers 4xx for an assembly it does not have, and the overwhelmingly
/// common cause of that is a typo — so it is the caller's argument that is
/// wrong, and the message offers the names available with no network at all.
/// Anything else (a refused connection, a 5xx, a timeout) is the network and
/// keeps its own error: telling someone their genome name is wrong when the
/// truth is that their wifi is off sends them to fix the wrong thing.
fn annotate(genome: &str, error: crate::error::Error) -> crate::error::Error {
    let client_error = matches!(
        &error,
        crate::error::Error::Http {
            status: Some(400..=499),
            ..
        }
    );
    if !client_error {
        return error;
    }
    crate::error::Error::invalid(format!(
        "no genome {genome:?} at UCSC ({error}). Bundled, and needing no network: {}",
        bundled::GENOME_NAMES.join(", ")
    ))
}

fn fetch_cached(genome: &str, full: bool) -> Result<ChrMap> {
    let key = (genome.to_string(), full);
    if let Some(hit) = fetch_cache().lock().get(&key) {
        return Ok(hit.clone());
    }
    let fetched = fetch_ucsc(genome, full)?;
    fetch_cache().lock().insert(key, fetched.clone());
    Ok(fetched)
}

/// `https://api.genome.ucsc.edu/list/chromosomes?genome=<genome>`.
///
/// Sorted and filtered exactly as `tools/bundle_genomes.py` sorts and filters
/// the bundled table, so which of the two answered a call shows up only in the
/// latency.
#[cfg(feature = "url")]
fn fetch_ucsc(genome: &str, full: bool) -> Result<ChrMap> {
    use crate::error::Error;

    let url = format!("https://api.genome.ucsc.edu/list/chromosomes?genome={genome}");
    let body = crate::source::http_get_text(&url)?;

    let parsed: serde_json::Value = serde_json::from_str(&body).map_err(|e| Error::Http {
        url: url.clone(),
        status: None,
        message: format!("response was not JSON: {e}"),
    })?;
    let chromosomes = parsed
        .get("chromosomes")
        .and_then(|c| c.as_object())
        .ok_or_else(|| Error::Http {
            url: url.clone(),
            status: None,
            message: format!("no chromosomes for genome {genome:?}"),
        })?;

    let mut entries: Vec<(String, i64)> = chromosomes
        .iter()
        .filter(|(id, _)| full || !id.contains('_'))
        .filter_map(|(id, size)| size.as_i64().map(|size| (id.clone(), size)))
        .collect();
    // By id as a string, which is what puts chr10 before chr2.
    entries.sort_by(|a, b| a.0.cmp(&b.0));

    if entries.is_empty() {
        return Err(Error::Http {
            url,
            status: None,
            message: format!("genome {genome:?} has no chromosomes"),
        });
    }
    Ok(ChrMap::from_entries(entries))
}

#[cfg(not(feature = "url"))]
fn fetch_ucsc(genome: &str, _full: bool) -> Result<ChrMap> {
    Err(crate::error::Error::Unsupported(format!(
        "genome {genome:?} is not bundled and the `url` feature is off"
    )))
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn bundled_genomes_need_no_network() {
        let mm10 = get_chr_sizes("mm10", false).unwrap();
        assert_eq!(mm10.len(), 22);
        assert_eq!(mm10.resolve("chr19").unwrap().size, 61_431_566);
        // Resolution rules apply to the result like any other chromosome map.
        assert_eq!(mm10.resolve("19").unwrap().id, "chr19");
    }

    #[test]
    fn genome_names_are_case_insensitive() {
        let a = get_chr_sizes("sacCer3", false).unwrap();
        let b = get_chr_sizes("SACCER3", false).unwrap();
        assert_eq!(a.names(), b.names());
    }

    #[test]
    fn the_order_is_string_order() {
        let hg38 = get_chr_sizes("hg38", false).unwrap();
        let names = hg38.names();
        assert_eq!(&names[..3], ["chr1", "chr10", "chr11"]);
        // Index follows that order, which is what a reader's chr_sizes shows.
        assert_eq!(hg38.by_index(0).unwrap().id, "chr1");
    }

    /// Which failure an unknown genome is depends on what came back from
    /// UCSC, and the point of the split is that the two are told apart rather
    /// than flattened into one message.
    ///
    /// Against `annotate` rather than against the network. The version that
    /// called `get_chr_sizes` asserted the same thing, but its answer depended
    /// on whether the machine running the suite could reach UCSC — so it made
    /// every `cargo test` issue a request, and offline it waited through three
    /// retries and their backoff to assert the branch it had not meant to take.
    /// A test whose meaning changes with the wifi is two tests wearing one
    /// name.
    #[test]
    fn an_unknown_genome_is_a_bad_argument_and_a_dead_network_is_not() {
        use crate::error::Error;

        // UCSC says 4xx: the name is wrong, and the caller is offered what
        // needs no network at all.
        let refused = Error::Http {
            url: "https://api.genome.ucsc.edu/list/chromosomes?genome=nosuchgenome_zz".into(),
            status: Some(400),
            message: "fetching failed with HTTP 400".into(),
        };
        let message = match annotate("nosuchgenome_zz", refused) {
            Error::InvalidArgument(message) => message,
            other => panic!("a 400 should be the caller's argument: {other:?}"),
        };
        assert!(message.contains("nosuchgenome_zz"), "{message}");
        assert!(message.contains("Bundled"), "{message}");
        assert!(
            message.contains("mm10") && message.contains("sacCer3"),
            "{message}"
        );

        // Anything else is the network, and keeps its own error: telling
        // someone their genome name is wrong when the truth is that their wifi
        // is off sends them to fix the wrong thing.
        for error in [
            Error::Http {
                url: "https://api.genome.ucsc.edu/x".into(),
                status: Some(503),
                message: "fetching failed with HTTP 503".into(),
            },
            Error::Http {
                url: "https://api.genome.ucsc.edu/x".into(),
                status: None,
                message: "fetching failed: dns error".into(),
            },
            Error::io(
                "https://api.genome.ucsc.edu/x",
                std::io::Error::other("down"),
            ),
        ] {
            let kind = format!("{error:?}");
            assert!(
                !matches!(annotate("hg38", error), Error::InvalidArgument(_)),
                "{kind} was flattened into a bad-argument error"
            );
        }
    }

    #[test]
    fn bundled_names_are_listed() {
        assert!(bundled_genomes().contains(&"hg38"));
        assert_eq!(bundled_genomes().len(), 7);
    }
}