rsomics_limma_array_weights/
lib.rs

1//! limma arrayWeights: per-array quality weights estimated by REML.
2//!
3//! Clean-room reimplementation of limma's `arrayWeights`. The method is Ritchie et al. (2006),
4//! BMC Bioinformatics 7:261, doi:10.1186/1471-2105-7-261, with the REML scoring adapted from
5//! Smyth (2002), J. Comput. Graph. Statist. 11:836. No limma (GPL) source was consulted; the
6//! implementation follows the published method and is validated black-box against the binary.
7
8mod matrix;
9mod reml;
10
11use std::io::Write;
12use std::path::Path;
13
14use rsomics_common::{Result, RsomicsError};
15
16pub use matrix::{read_design, read_expr};
17
18pub struct Options<'a> {
19    pub expr: &'a Path,
20    pub design: &'a Path,
21    pub prior_n: f64,
22    pub maxiter: usize,
23    pub tol: f64,
24}
25
26pub struct Results {
27    pub samples: Vec<String>,
28    pub weights: Vec<f64>,
29    pub n_genes: usize,
30    pub iters: usize,
31    pub converged: bool,
32}
33
34pub fn run(opts: &Options) -> Result<Results> {
35    let expr = read_expr(opts.expr)?;
36    let design = read_design(opts.design)?;
37    if design.x.len() != expr.samples.len() {
38        return Err(RsomicsError::InvalidInput(format!(
39            "design has {} rows, expression has {} samples",
40            design.x.len(),
41            expr.samples.len()
42        )));
43    }
44
45    let fit = reml::array_weights(&expr.y, &design.x, opts.prior_n, opts.maxiter, opts.tol)?;
46
47    Ok(Results {
48        samples: expr.samples,
49        weights: fit.weights,
50        n_genes: expr.y.len(),
51        iters: fit.iters,
52        converged: fit.converged,
53    })
54}
55
56pub fn write_weights(res: &Results, out: &mut dyn Write) -> Result<()> {
57    let mut buf = ryu::Buffer::new();
58    writeln!(out, "sample\tweight").map_err(RsomicsError::Io)?;
59    for (s, w) in res.samples.iter().zip(&res.weights) {
60        writeln!(out, "{s}\t{}", buf.format(*w)).map_err(RsomicsError::Io)?;
61    }
62    Ok(())
63}
rsomics_limma_array_weights/lib.rs

rsomics_limma_array_weights/
lib.rs
