use anyhow::{bail, Context, Result};
use clap::Parser;
use std::collections::HashMap;
use std::path::PathBuf;
use std::time::Instant;
use uno_core::im::ImModel;
use uno_core::{
match_features, match_features_topk, score, DigestConfig, FeatureMatrix, MatchConfig,
PeptideDb, ProteinConfig, RtModel, ScoreConfig,
};
mod im_cli;
mod koth_ingest;
mod quant_cli;
mod quant_export;
mod rt_cli;
#[cfg(test)]
mod rt_command_tests {
use super::*;
#[test]
fn rt_commands_do_not_require_search_inputs() {
for argv in [
vec!["uno", "rt", "inspect", "custom.urtc"],
vec![
"uno",
"rt",
"train",
"--data",
"labels.tsv",
"--out",
"custom.urtc",
"--units",
"minutes",
"--method",
"C18",
],
vec![
"uno",
"rt",
"evaluate",
"--model",
"custom.urtc",
"--data",
"heldout.tsv",
"--units",
"minutes",
],
] {
assert!(Cli::try_parse_from(argv).is_ok());
}
assert!(Cli::try_parse_from(["uno", "--features", "run.tsv"]).is_err());
assert!(Cli::try_parse_from(["uno", "--fasta", "db.fasta"]).is_err());
}
#[test]
fn im_defaults_to_bundled_with_explicit_disable_and_override() {
let base = ["uno", "--features", "run.tsv", "--fasta", "db.fasta"];
let cli = Cli::try_parse_from(base).unwrap();
assert!(load_models(&cli).unwrap().im.unwrap().cnn_info().is_some());
let cli = Cli::try_parse_from(base.into_iter().chain(["--no-im"])).unwrap();
assert!(load_models(&cli).unwrap().im.is_none());
assert!(Cli::try_parse_from(base.into_iter().chain([
"--no-im",
"--im-model",
"custom.uimc"
]))
.is_err());
let cli =
Cli::try_parse_from(base.into_iter().chain(["--im-model", "missing.uimc"])).unwrap();
assert!(load_models(&cli).is_err());
}
#[test]
fn long_and_batch_inputs_are_explicit_and_exclusive() {
for args in [
vec!["uno", "--matrix-long", "bundle", "--fasta", "db"],
vec![
"uno",
"--matrix-long",
"bundle",
"--fasta",
"db",
"--matrix-batch-out",
"out",
],
] {
assert!(Cli::try_parse_from(args).is_ok());
}
assert!(Cli::try_parse_from([
"uno",
"--matrix-long",
"bundle",
"--features",
"f",
"--fasta",
"db"
])
.is_err());
assert!(Cli::try_parse_from([
"uno",
"--matrix-long",
"bundle",
"--matrix-batch-out",
"out",
"--out-proteins",
"p",
"--fasta",
"db"
])
.is_err());
}
#[test]
fn joint_mode_requires_long_input_and_cannot_be_a_batch_search() {
assert!(Cli::try_parse_from([
"uno",
"--matrix-long",
"b",
"--matrix-joint",
"--fasta",
"db"
])
.is_ok());
assert!(
Cli::try_parse_from(["uno", "--matrix", "b", "--matrix-joint", "--fasta", "db"])
.is_err()
);
assert!(Cli::try_parse_from([
"uno",
"--matrix-long",
"b",
"--matrix-joint",
"--matrix-batch-out",
"o",
"--fasta",
"db"
])
.is_err());
}
#[test]
fn rt_sources_are_mutually_exclusive_and_default_needs_no_flag() {
assert!(
Cli::try_parse_from(["uno", "--features", "run.tsv", "--fasta", "db.fasta"]).is_ok()
);
assert!(Cli::try_parse_from([
"uno",
"--features",
"run.tsv",
"--fasta",
"db.fasta",
"--rt-model",
"a",
"--rt-library",
"b"
])
.is_err());
assert!(Cli::try_parse_from([
"uno",
"rt",
"train",
"--data",
"d",
"--out",
"o",
"--units",
"irt",
"--method",
"m",
"--init",
"a",
"--from-scratch"
])
.is_err());
}
}
#[derive(Parser, Debug)]
#[command(
name = "uno",
version,
args_conflicts_with_subcommands = true,
subcommand_negates_reqs = true
)]
struct Cli {
#[command(subcommand)]
command: Option<rt_cli::Command>,
#[arg(long, required_unless_present_any = ["raw", "matrix", "matrix_long"], conflicts_with_all = ["raw", "matrix", "matrix_long"])]
features: Option<PathBuf>,
#[arg(long)]
raw: Option<PathBuf>,
#[arg(long, conflicts_with_all = ["features", "raw"])]
matrix: Option<PathBuf>,
#[arg(long, conflicts_with_all = ["features", "raw", "matrix"])]
matrix_long: Option<PathBuf>,
#[arg(long, requires="matrix_long", conflicts_with_all=["features","raw","matrix","matrix_batch_out","matrix_no_feature_quality","matrix_legacy_values","matrix_prefer_isotope_support","dump_candidate_peptides"])]
matrix_joint: bool,
#[arg(long, requires="matrix_long", conflicts_with_all=["out_proteins","out_peptides","out_lfq","out_abundance","dump_candidate_peptides","matrix_no_feature_quality","matrix_legacy_values"])]
matrix_batch_out: Option<PathBuf>,
#[arg(long, requires = "matrix_long")]
matrix_no_feature_quality: bool,
#[arg(long, requires = "matrix_long")]
matrix_legacy_values: bool,
#[arg(long, requires="matrix_long", conflicts_with_all=["matrix_no_feature_quality","matrix_batch_out"])]
matrix_prefer_isotope_support: bool,
#[arg(long)]
matrix_min_runs: Option<u32>,
#[arg(long, required = true)]
fasta: Option<PathBuf>,
#[arg(long, conflicts_with = "rt_library")]
rt_model: Option<PathBuf>,
#[arg(long, conflicts_with = "no_im")]
im_model: Option<PathBuf>,
#[arg(long)]
no_im: bool,
#[arg(long)]
config: Option<PathBuf>,
#[arg(long)]
rt_library: Option<PathBuf>,
#[arg(long)]
dump_candidate_peptides: Option<PathBuf>,
#[arg(long)]
out_proteins: Option<PathBuf>,
#[arg(long)]
out_peptides: Option<PathBuf>,
#[arg(long, conflicts_with_all = ["matrix", "matrix_long", "matrix_batch_out"])]
out_quant_evidence: Option<PathBuf>,
#[arg(long, conflicts_with_all = ["matrix", "matrix_long", "matrix_batch_out"])]
out_peptide_scores: Option<PathBuf>,
#[arg(long)]
out_lfq: Option<PathBuf>,
#[arg(long)]
out_abundance: Option<PathBuf>,
#[arg(long)]
profile_reassign: bool,
#[arg(long)]
seed: Option<u64>,
}
#[derive(Debug, Default, serde::Deserialize)]
#[serde(default)]
struct Config {
digest: DigestConfig,
matching: MatchConfig,
scoring: ScoreConfig,
protein: ProteinConfig,
koth: Option<koth_ms::config::KothConfig>,
}
impl Config {
fn apply_cli_seed(&mut self, seed: Option<u64>) {
if let Some(seed) = seed {
self.digest.decoy_seed = seed;
self.scoring.seed = seed;
}
}
}
fn main() -> Result<()> {
env_logger::Builder::from_env(env_logger::Env::default().default_filter_or("info")).init();
let cli = Cli::parse();
if let Some(command) = cli.command {
return rt_cli::run(command);
}
let mut config: Config = match &cli.config {
Some(p) => {
let s =
std::fs::read_to_string(p).with_context(|| format!("reading {}", p.display()))?;
toml::from_str(&s).with_context(|| "parsing config")?
}
None => Config::default(),
};
config.apply_cli_seed(cli.seed);
if let Some(path) = &cli.out_peptide_scores {
config.scoring.diagnostic_peptides_path = Some(path.to_string_lossy().into_owned());
}
if let Some(out) = &cli.matrix_batch_out {
return search_batch(&cli, &config, out);
}
let matrix_min_runs = cli
.matrix_min_runs
.unwrap_or(if cli.matrix_joint { 1 } else { 2 });
let t = Instant::now();
let mut matrix = if let Some(path) = &cli.matrix_long {
let bundle = uno_core::long_input::Bundle::open(path)?;
if cli.matrix_joint {
bundle.load_joint(matrix_min_runs)?
} else {
bundle.load_consensus_with_representative(
matrix_min_runs,
!cli.matrix_no_feature_quality,
cli.matrix_legacy_values,
cli.matrix_prefer_isotope_support,
)?
}
} else {
match (&cli.features, &cli.raw, &cli.matrix) {
(Some(features), None, None) => FeatureMatrix::load(features)?,
(None, Some(raw), None) => {
let koth_cfg = match &config.koth {
Some(c) => c.clone(),
None => koth_ingest::default_config_for(raw)?,
};
if koth_ingest::is_bruker(raw) {
log::info!("raw input is a Bruker .d: using the bundled timsTOF koth profile");
if config.matching.min_isotopes == 0
&& config.matching.min_combined_score == 0.0
{
config.matching.min_isotopes = 4;
config.matching.min_combined_score = 0.85;
log::info!(
"Bruker defaults: matching.min_isotopes=4, matching.min_combined_score=0.85 \
(timsTOF mass accuracy is ~3 ppm, so the search must be restricted to \
well-characterised features; override in [matching])"
);
}
}
koth_ingest::features_from_raw(raw, &koth_cfg)
.with_context(|| format!("finding features in {}", raw.display()))?
}
(None, None, Some(mp)) => {
let m = FeatureMatrix::load_matrix(mp, matrix_min_runs)
.with_context(|| format!("reading aligned matrix {}", mp.display()))?;
log::info!(
"aligned matrix: {} runs, {} consensus features (min_runs={})",
m.n_samples(),
m.len(),
matrix_min_runs
);
m
}
_ => bail!("provide exactly one of --features, --raw, or --matrix"),
}
};
uno_core::joint::restrict_observations(&mut matrix, &config.matching);
log::info!(
"loaded {} features from input ({:.1}s)",
matrix.len(),
t.elapsed().as_secs_f64()
);
let t = Instant::now();
let db = PeptideDb::build(
cli.fasta
.as_ref()
.context("--fasta is required for search")?,
&config.digest,
)?;
log::info!(
"built peptide db: {} targets + {} decoys ({:.1}s)",
db.n_targets,
db.n_decoys,
t.elapsed().as_secs_f64()
);
let models = load_models(&cli)?;
search_loaded(&cli, &config, &matrix, &db, &models, None)
}
struct Models {
rt: Option<RtModel>,
im: Option<ImModel>,
library: Option<HashMap<String, f64>>,
}
fn load_models(cli: &Cli) -> Result<Models> {
if cli.dump_candidate_peptides.is_some() {
return Ok(Models {
rt: None,
im: None,
library: None,
});
}
let rt_library: Option<HashMap<String, f64>> = match &cli.rt_library {
Some(p) => {
let s =
std::fs::read_to_string(p).with_context(|| format!("reading {}", p.display()))?;
let mut map = HashMap::new();
for line in s.lines() {
let mut it = line.split('\t');
if let (Some(seq), Some(rt)) = (it.next(), it.next()) {
if let Ok(rt) = rt.parse::<f64>() {
map.insert(seq.to_string(), rt);
}
}
}
log::info!(
"loaded RT library with {} peptides from {}",
map.len(),
p.display()
);
Some(map)
}
None => None,
};
let rt_model = match &cli.rt_model {
Some(p) => {
let m = RtModel::load(p)?;
log::info!("loaded RT model from {}", p.display());
Some(m)
}
None if rt_library.is_none() => {
log::info!("using bundled RT CNN");
Some(RtModel::bundled()?)
}
None => None,
};
let im_model = match &cli.im_model {
Some(p) => {
let m = ImModel::load(p)?;
log::info!("loaded IM model from {}", p.display());
Some(m)
}
None if cli.no_im => None,
None => {
log::info!("using bundled IM CNN (soft scoring after RT shortlisting)");
Some(ImModel::bundled()?)
}
};
Ok(Models {
rt: rt_model,
im: im_model,
library: rt_library,
})
}
fn search_loaded(
cli: &Cli,
config: &Config,
matrix: &FeatureMatrix,
db: &PeptideDb,
models: &Models,
dir: Option<&std::path::Path>,
) -> Result<()> {
let out_proteins = dir
.map(|p| p.join("proteins.tsv"))
.or_else(|| cli.out_proteins.clone());
let out_peptides = dir
.map(|p| p.join("peptides.tsv"))
.or_else(|| cli.out_peptides.clone());
let out_lfq = dir
.map(|p| p.join("lfq.tsv"))
.or_else(|| cli.out_lfq.clone());
let out_abundance = dir
.map(|p| p.join("abundance.tsv"))
.or_else(|| cli.out_abundance.clone());
if let Some(path) = &cli.dump_candidate_peptides {
let t = Instant::now();
let pfms = match_features(matrix, db, &config.matching);
log::info!(
"generated {} PFMs at {} ppm, offsets {:?} ({:.1}s)",
pfms.len(),
config.matching.ppm_tolerance,
config.matching.isotope_offsets,
t.elapsed().as_secs_f64()
);
use std::io::Write;
let mut seen = std::collections::HashSet::new();
let mut w = std::io::BufWriter::new(std::fs::File::create(path)?);
for p in &pfms {
let seq = db.peptides[p.peptide_idx as usize].sequence.as_ref();
if seen.insert(seq.to_string()) {
writeln!(w, "{}", seq)?;
}
}
log::info!(
"dumped {} unique candidate peptides to {}",
seen.len(),
path.display()
);
return Ok(());
}
let t = Instant::now();
let pfms = match_features_topk(matrix, db, &config.matching, config.scoring.top_k);
let n_decoy = pfms.iter().filter(|p| p.decoy).count();
log::info!(
"generated {} reduced PFMs ({} target, {} decoy) at {} ppm, offsets {:?}, top-{} ({:.1}s)",
pfms.len(),
pfms.len() - n_decoy,
n_decoy,
config.matching.ppm_tolerance,
config.matching.isotope_offsets,
config.scoring.top_k,
t.elapsed().as_secs_f64()
);
let t = Instant::now();
let result = score(
matrix,
db,
&pfms,
models.rt.as_ref(),
models.im.as_ref(),
models.library.as_ref(),
&config.scoring,
);
log::info!("scored in {:.1}s", t.elapsed().as_secs_f64());
let proteins = result.proteins;
let n_prot_1pct = proteins
.iter()
.filter(|p| !p.decoy && p.q_value <= 0.01)
.count();
let n_prot_decoy_1pct = proteins
.iter()
.filter(|p| p.decoy && p.q_value <= 0.01)
.count();
log::info!(
"proteins: {} target at 1% FDR ({} decoy proteins below 1%)",
n_prot_1pct,
n_prot_decoy_1pct
);
if let Some(path) = &out_proteins {
use std::io::Write;
let confident_all: std::collections::HashSet<String> = proteins
.iter()
.filter(|p| !p.decoy)
.map(|p| p.accession.clone())
.collect();
let (groups, _) = uno_core::group::build_groups(db, &result.peptides, &confident_all);
let members = uno_core::group::members_by_leading(&groups);
let mut w = std::io::BufWriter::new(std::fs::File::create(path)?);
writeln!(
w,
"accession\tdecoy\tn_group_members\tgroup_members\tobserved_peptides\ttheoretical_peptides\tcoverage\tbinomial_score\tmean_rt_diff\tq_value"
)?;
for p in &proteins {
let mem = members.get(&(p.decoy, p.accession.clone()));
let n_mem = mem.map(|m| m.len()).unwrap_or(1);
let mem_str = mem
.map(|m| m.join(";"))
.unwrap_or_else(|| p.accession.clone());
writeln!(
w,
"{}\t{}\t{}\t{}\t{}\t{}\t{:.4}\t{:.3}\t{:.3}\t{:.5}",
p.accession,
p.decoy as u8,
n_mem,
mem_str,
p.observed_peptides,
p.theoretical_peptides,
p.coverage,
p.binomial_score,
p.mean_rt_diff,
p.q_value
)?;
}
log::info!("wrote {} proteins to {}", proteins.len(), path.display());
}
if let Some(path) = &out_peptides {
use std::io::Write;
let mut w = std::io::BufWriter::new(std::fs::File::create(path)?);
write!(w, "peptide\tdecoy\tproteins\tcharge\tmass\trt")?;
for rn in &matrix.run_names {
write!(w, "\t{rn}")?;
}
writeln!(w)?;
for hit in &result.peptides {
let pep = &db.peptides[hit.peptide_idx as usize];
let f = &matrix.features[hit.feature_idx as usize];
let prots = pep
.proteins
.iter()
.map(|p| p.as_ref())
.collect::<Vec<_>>()
.join(";");
write!(
w,
"{}\t{}\t{}\t{}\t{:.5}\t{:.4}",
pep.sequence.as_ref(),
hit.decoy as u8,
prots,
f.charge,
f.mass,
f.rt
)?;
for v in &f.intensities {
write!(w, "\t{:.5e}", v)?;
}
writeln!(w)?;
}
log::info!(
"wrote {} RT-gated peptide-intensity rows to {}",
result.peptides.len(),
path.display()
);
}
if let Some(path) = &cli.out_quant_evidence {
quant_export::write(path, matrix, db, &result.peptides)?;
}
if out_lfq.is_some() || out_abundance.is_some() {
use std::io::Write;
let prot_score: std::collections::HashMap<(bool, &str), f64> = proteins
.iter()
.map(|p| ((p.decoy, p.accession.as_str()), p.binomial_score))
.collect();
let mut quant_peptides =
uno_core::quant::resolve_feature_competition(&result.peptides, db, &prot_score);
log::info!(
"feature competition: {} -> {} peptide-features after 1-feature-1-protein claim",
result.peptides.len(),
quant_peptides.len()
);
if cli.profile_reassign {
let pfms = match_features(matrix, db, &config.matching);
let candidates: Vec<(u32, u32, bool)> = pfms
.iter()
.filter(|p| p.isotope_offset == 0)
.filter(|p| {
db.peptides[p.peptide_idx as usize]
.proteins
.iter()
.any(|pr| prot_score.contains_key(&(p.decoy, pr.as_ref())))
})
.map(|p| (p.feature_idx, p.peptide_idx, p.decoy))
.collect();
let before = quant_peptides.len();
quant_peptides = uno_core::quant::profile_reassign(
matrix,
db,
&quant_peptides,
&candidates,
0.5,
0.25,
);
log::info!(
"profile-reassign: {} candidates, {} -> {} peptide-features",
candidates.len(),
before,
quant_peptides.len()
);
}
if let Some(path) = &out_lfq {
let mut w = std::io::BufWriter::new(std::fs::File::create(path)?);
write!(w, "peptide\tcharge\tdecoy\tproteins\tmass\trt")?;
for rn in &matrix.run_names {
write!(w, "\t{rn}")?;
}
writeln!(w)?;
for h in &quant_peptides {
let pep = &db.peptides[h.peptide_idx as usize];
let f = &matrix.features[h.feature_idx as usize];
let prots = pep
.proteins
.iter()
.map(|p| p.as_ref())
.collect::<Vec<_>>()
.join(";");
write!(
w,
"{}\t{}\t{}\t{}\t{:.5}\t{:.4}",
pep.sequence.as_ref(),
f.charge,
h.decoy as u8,
prots,
f.mass,
f.rt
)?;
for v in &f.intensities {
write!(w, "\t{:.5e}", v)?;
}
writeln!(w)?;
}
log::info!(
"wrote {} precursor LFQ rows to {}",
quant_peptides.len(),
path.display()
);
}
if let Some(path) = &out_abundance {
let id_confident: std::collections::HashSet<String> = proteins
.iter()
.filter(|p| !p.decoy && p.q_value <= 0.01)
.map(|p| p.accession.clone())
.collect();
let (groups, pep_group) =
uno_core::group::build_groups(db, &quant_peptides, &id_confident);
let members = uno_core::group::members_by_leading(&groups);
let n_runs = matrix.run_names.len();
let ab = uno_core::quant::group_abundance(matrix, &quant_peptides, &pep_group, n_runs);
let mut w = std::io::BufWriter::new(std::fs::File::create(path)?);
write!(w, "accession\tn_group_members\tgroup_members\tn_peptides")?;
for rn in &matrix.run_names {
write!(w, "\t{rn}")?;
}
writeln!(w)?;
#[allow(clippy::type_complexity)]
let mut rows: Vec<(&String, &(Vec<f64>, std::collections::HashSet<u32>))> =
ab.iter().collect();
rows.sort_by(|a, b| a.0.cmp(b.0));
for (leading, (ints, peps)) in rows {
let mem = members.get(&(false, leading.to_string()));
let n_mem = mem.map(|m| m.len()).unwrap_or(1);
let mem_str = mem.map(|m| m.join(";")).unwrap_or_else(|| leading.clone());
write!(w, "{}\t{}\t{}\t{}", leading, n_mem, mem_str, peps.len())?;
for v in ints {
write!(w, "\t{:.5e}", v)?;
}
writeln!(w)?;
}
log::info!("wrote {} protein groups to {}", ab.len(), path.display());
}
}
Ok(())
}
#[cfg(test)]
mod seed_config_tests {
use super::*;
#[test]
fn cli_seed_keeps_legacy_joint_behavior() {
let cli = Cli::try_parse_from([
"uno",
"--features",
"run.tsv",
"--fasta",
"db.fasta",
"--seed",
"9",
])
.unwrap();
let mut cfg: Config =
toml::from_str("[digest]\ndecoy_seed = 2\n[scoring]\nseed = 3\n").unwrap();
cfg.apply_cli_seed(cli.seed);
assert_eq!(cfg.digest.decoy_seed, 9);
assert_eq!(cfg.scoring.seed, 9);
assert_eq!(cfg.scoring.effective_rt_seed(), 9);
assert_eq!(cfg.scoring.effective_cv_seed(), 9);
assert_eq!(cfg.scoring.effective_classifier_seed(), 9);
}
#[test]
fn explicit_stage_seeds_survive_cli_override() {
let mut cfg: Config = toml::from_str(
"[digest]\ndecoy_seed = 2\n[scoring]\nseed = 3\nrt_seed = 0\ncv_seed = 5\nclassifier_seed = 6\n",
).unwrap();
cfg.apply_cli_seed(Some(9));
assert_eq!(cfg.digest.decoy_seed, 9);
assert_eq!(cfg.scoring.seed, 9);
assert_eq!(cfg.scoring.effective_rt_seed(), 0);
assert_eq!(cfg.scoring.effective_cv_seed(), 5);
assert_eq!(cfg.scoring.effective_classifier_seed(), 6);
}
#[test]
fn absent_cli_seed_preserves_independent_config_seeds() {
let mut cfg: Config =
toml::from_str("[digest]\ndecoy_seed = 2\n[scoring]\nseed = 3\ncv_seed = 5\n").unwrap();
cfg.apply_cli_seed(None);
assert_eq!(cfg.digest.decoy_seed, 2);
assert_eq!(cfg.scoring.effective_rt_seed(), 3);
assert_eq!(cfg.scoring.effective_cv_seed(), 5);
assert_eq!(cfg.scoring.effective_classifier_seed(), 3);
}
}
fn search_batch(cli: &Cli, config: &Config, out: &std::path::Path) -> Result<()> {
use std::io::Write;
let bundle = uno_core::long_input::Bundle::open(cli.matrix_long.as_ref().unwrap())?;
std::fs::create_dir_all(out)?;
anyhow::ensure!(
!out.join("batch_manifest.json").exists(),
"batch output already complete; select a new output directory"
);
let db = PeptideDb::build(cli.fasta.as_ref().unwrap(), &config.digest)?;
let models = load_models(cli)?;
let mut runs = Vec::new();
let mut union: std::collections::BTreeMap<String, Vec<usize>> =
std::collections::BTreeMap::new();
uno_core::rt_cache::with_memory_cache(|| -> Result<()> {
for r in &bundle.runs {
let started = Instant::now();
let matrix = bundle.load_run(r.run_id)?;
let dir = out.join(format!("run{:03}", r.run_id));
std::fs::create_dir_all(&dir)?;
search_loaded(cli, config, &matrix, &db, &models, Some(&dir))?;
let mut reader = csv::ReaderBuilder::new()
.delimiter(b'\t')
.from_path(dir.join("proteins.tsv"))?;
let header = reader.headers()?.clone();
let col = |name| header.iter().position(|h| h == name).unwrap();
for row in reader.records() {
let row = row?;
if &row[col("decoy")] == "0" && row[col("q_value")].parse::<f64>()? <= 0.01 {
for acc in row[col("group_members")].split(';') {
let ids = union.entry(acc.to_owned()).or_default();
if !ids.contains(&r.run_id) {
ids.push(r.run_id);
}
}
}
}
runs.push(serde_json::json!({"run_id":r.run_id,"name":r.name,"seconds":started.elapsed().as_secs_f64(),"n_features":matrix.len(),"directory":dir}));
}
Ok(())
})?;
let mut writer =
std::io::BufWriter::new(std::fs::File::create(out.join("accession_union.tsv"))?);
writeln!(writer, "accession\trun_ids")?;
for (acc, ids) in &union {
writeln!(
writer,
"{}\t{}",
acc,
ids.iter()
.map(|i| i.to_string())
.collect::<Vec<_>>()
.join(";")
)?;
}
writer.flush()?;
std::fs::write(
out.join("batch_manifest.json"),
serde_json::to_vec_pretty(
&serde_json::json!({"runs":runs,"union_accessions":union.len(),"interpretation":"Union of per-run calls at q <= 0.01; no calibrated cohort q-value is claimed."}),
)?,
)?;
Ok(())
}