1use std::sync::Arc;
2
3use log::{debug, info, warn};
4
5use crate::convert::try_open_or_convert;
6use crate::sparse_io::COLUMN_SEP;
7use crate::sparse_io_vector::{data_name_of, SparseIoVec};
8use legume_numeric::matrix::common_io::{self, read_lines};
9use rustc_hash::FxHashSet;
10
11use crate::aux::feature_names::FeatureNameKind;
12use crate::sparse_io_vector::{ColumnAlignment, RowAlignment};
13
14const MULTIMODAL_HINT_DISJOINTNESS_FRACTION: f64 = 0.5;
21
22#[derive(Default)]
24pub struct ReadSharedRowsArgs {
25 pub data_files: Vec<Box<str>>,
26 pub batch_files: Option<Vec<Box<str>>>,
27 pub preload: bool,
28 pub feature_kind: Option<FeatureNameKind>,
33 pub row_alignment: RowAlignment,
41 pub column_alignment: ColumnAlignment,
50 pub qc: Option<crate::qc_lib::QcConfig>,
57 pub keep_empty_barcodes: bool,
68 pub qc_exempt_files: Option<Vec<bool>>,
77 pub qc_block_size: Option<usize>,
79 pub qc_report_out: Option<Box<str>>,
81 pub per_file_feature_suffix: Option<Vec<Box<str>>>,
88 pub per_file_barcode_suffix: Option<Vec<Option<Box<str>>>>,
98}
99
100pub struct SparseDataWithBatch {
102 pub data: SparseIoVec,
103 pub batch: Vec<Box<str>>,
104 pub output_keep_idx: Option<Vec<usize>>,
108}
109
110pub fn read_data_on_shared_rows(args: ReadSharedRowsArgs) -> anyhow::Result<SparseDataWithBatch> {
117 type OpenedBackend = Box<dyn crate::sparse_io::SparseIo<IndexIter = Vec<usize>>>;
120 let mut opened: Vec<(Box<str>, OpenedBackend)> = Vec::with_capacity(args.data_files.len());
121 for data_file in args.data_files.iter() {
122 info!("Importing data file: {}", data_file);
123 let mut data = try_open_or_convert(data_file)?;
124 if args.preload || crate::sparse_io::auto_preload_enabled() {
125 data.preload_columns()?;
126 }
127 opened.push((data_file.clone(), data));
128 }
129
130 let mut data_vec = SparseIoVec::new()
137 .with_row_alignment(args.row_alignment)
138 .expect("with_row_alignment on empty SparseIoVec")
139 .with_column_alignment(args.column_alignment)
140 .expect("with_column_alignment on empty SparseIoVec");
141
142 if let Some(suffix) = args.per_file_feature_suffix.clone() {
146 anyhow::ensure!(
147 suffix.len() == args.data_files.len(),
148 "per_file_feature_suffix has {} entries but {} data files were given",
149 suffix.len(),
150 args.data_files.len(),
151 );
152 data_vec = data_vec
153 .with_per_backend_row_suffix(suffix)
154 .expect("with_per_backend_row_suffix on empty SparseIoVec");
155 }
156
157 if let Some(sfx) = args.per_file_barcode_suffix.as_ref() {
160 anyhow::ensure!(
161 sfx.len() == args.data_files.len(),
162 "per_file_barcode_suffix has {} entries but {} data files were given",
163 sfx.len(),
164 args.data_files.len(),
165 );
166 }
167
168 use crate::aux::feature_names::FeatureNameKind;
169
170 let needs_names = args.feature_kind.is_none()
173 || args
174 .feature_kind
175 .as_ref()
176 .is_some_and(|k| k.needs_global_pass());
177 let mut file_ends: Vec<usize> = Vec::with_capacity(opened.len());
182 let all_names: Option<Vec<Box<str>>> = if needs_names {
183 let mut acc: Vec<Box<str>> = Vec::new();
184 for (_, d) in opened.iter() {
185 acc.extend(d.row_names()?);
186 file_ends.push(acc.len());
187 }
188 Some(acc)
189 } else {
190 None
191 };
192
193 let kind_was_auto = args.feature_kind.is_none();
194 let resolved_kind: FeatureNameKind = match args.feature_kind.clone() {
195 Some(k) => k,
196 None => {
197 let names = all_names.as_ref().expect("peeked when auto");
204 let per_file: Vec<FeatureNameKind> = std::iter::once(0)
208 .chain(file_ends.iter().copied())
209 .zip(file_ends.iter().copied())
210 .map(|(start, end)| FeatureNameKind::auto_detect(&names[start..end]))
211 .collect();
212 let k = FeatureNameKind::reconcile(&per_file);
213 debug!(
214 "Row alignment: auto-detected feature name kind → {:?} (per file: {:?}; {} rows)",
215 k,
216 per_file,
217 names.len()
218 );
219 k
220 }
221 };
222
223 match &resolved_kind {
228 FeatureNameKind::Mixed => {
229 let names = all_names.as_ref().expect("peeked for Mixed").clone();
230 debug!(
231 "Row alignment: building MIXED-kind canonical map over {} names \
232 across {} file(s)",
233 names.len(),
234 opened.len()
235 );
236 let canon = crate::aux::feature_names::build_mixed_kind_canonicalizer(&names);
237 data_vec = data_vec
238 .with_row_canonicalizer(move |name| canon(name))
239 .expect("with_row_canonicalizer on empty SparseIoVec");
240 }
241 FeatureNameKind::Locus {
242 merge_overlapping: true,
243 } => {
244 let names = all_names
245 .as_ref()
246 .expect("peeked for Locus merge_overlapping")
247 .clone();
248 debug!(
249 "Row alignment: building locus-overlap canonical map over {} names \
250 across {} file(s)",
251 names.len(),
252 opened.len()
253 );
254 let canon = crate::aux::feature_names::build_locus_overlap_canonicalizer(&names);
255 data_vec = data_vec
256 .with_row_canonicalizer(move |name| canon(name))
257 .expect("with_row_canonicalizer on empty SparseIoVec");
258 }
259 kind => {
260 if let Some(canon) = kind.clone().into_canonicalizer() {
261 debug!(
262 "Row alignment: applying {:?} canonicalizer across {} file(s)",
263 kind,
264 opened.len()
265 );
266 data_vec = data_vec
269 .with_row_canonicalizer(move |name| canon(name))
270 .expect("with_row_canonicalizer on empty SparseIoVec");
271 }
272 }
273 }
274 info!(
275 "Row alignment: {:?} · {:?} canon{} · {} file(s)",
276 args.row_alignment,
277 resolved_kind,
278 if kind_was_auto { " (auto)" } else { "" },
279 opened.len(),
280 );
281 let n_files = opened.len();
282 for (file_idx, (data_file, data)) in opened.into_iter().enumerate() {
283 let data_name = match args.column_alignment {
285 ColumnAlignment::Disjoint => data_name_of(&data_file, n_files)?,
286 ColumnAlignment::Union => None,
287 };
288 let barcode_suffix: Option<&str> = args
290 .per_file_barcode_suffix
291 .as_ref()
292 .and_then(|v| v[file_idx].as_deref());
293 data_vec.push_with_barcode_suffix(Arc::from(data), data_name, barcode_suffix)?;
294 }
295
296 let intersection_size = data_vec.num_rows();
300 for j in 0..data_vec.len() {
301 let backend_rows = data_vec[j].num_rows().unwrap_or(0);
302 if backend_rows != intersection_size {
303 info!(
304 "Backend {} has {} rows; using {} shared rows for fitting",
305 j, backend_rows, intersection_size
306 );
307 }
308 }
309
310 if args.column_alignment == ColumnAlignment::Disjoint && data_vec.len() >= 2 {
316 maybe_warn_multimodal_pattern(&data_vec);
317 }
318
319 let n_cells = data_vec.num_columns();
321 let mut batch_membership: Vec<Box<str>> = match args.column_alignment {
322 ColumnAlignment::Disjoint => {
323 resolve_batch_disjoint(&args.data_files, &data_vec, args.batch_files.as_deref())?
324 }
325 ColumnAlignment::Union => {
326 resolve_batch_union(&data_vec, args.batch_files.as_deref(), n_cells)?
327 }
328 };
329
330 if batch_membership.len() != data_vec.num_columns() {
331 return Err(anyhow::anyhow!(
332 "# batch membership {} != # of columns {}",
333 batch_membership.len(),
334 data_vec.num_columns()
335 ));
336 }
337
338 if !args.keep_empty_barcodes {
343 if let Some(flags) = args.qc_exempt_files.as_ref() {
344 anyhow::ensure!(
345 flags.len() == args.data_files.len(),
346 "qc_exempt_files has {} entries for {} data files",
347 flags.len(),
348 args.data_files.len(),
349 );
350 }
351 if let Some(keep) = empty_barcode_keep(&data_vec, args.qc_exempt_files.as_deref()) {
352 data_vec.mask_columns(&keep)?;
353 batch_membership = crate::qc_lib::filter_by_keep(&batch_membership, &keep);
354 }
355 }
356
357 let output_keep_idx = if let Some(cfg) = args.qc.as_ref() {
361 let mut exempt: Option<Vec<bool>> = None;
365 if let Some(flags) = args.qc_exempt_files.as_ref() {
366 anyhow::ensure!(
367 flags.len() == args.data_files.len(),
368 "qc_exempt_files has {} entries for {} data files",
369 flags.len(),
370 args.data_files.len(),
371 );
372 if flags.iter().any(|&f| f) {
373 exempt = Some(
374 (0..data_vec.num_columns())
375 .map(|c| data_vec.column_source(c).is_some_and(|b| flags[b]))
376 .collect(),
377 );
378 }
379 }
380 let report = crate::qc_lib::compute_qc_exempting(
381 &data_vec,
382 cfg,
383 args.qc_block_size,
384 exempt.as_deref(),
385 )?;
386 if let Some(path) = args.qc_report_out.as_deref() {
387 crate::qc_lib::write_qc_report(path, &data_vec.column_names()?, &report)?;
388 }
389 let n_near_empty = report.near_empty.iter().filter(|&&e| e).count();
390 info!(
391 "QC: dropped {}/{} cells from training, {} near-empty masked at output, {}/{} features dropped",
392 report.n_cells_dropped,
393 report.train_keep.len(),
394 n_near_empty,
395 report.n_features_dropped,
396 report.feature_keep.len(),
397 );
398 if report.n_features_dropped > 0 {
400 data_vec.mask_rows(&report.feature_keep)?;
401 }
402 let keep_idx = report.output_keep_idx();
405 if report.n_cells_dropped > 0 {
406 data_vec.mask_columns(&report.train_keep)?;
407 batch_membership = crate::qc_lib::filter_by_keep(&batch_membership, &report.train_keep);
408 }
409 Some(keep_idx)
410 } else {
411 None
412 };
413
414 Ok(SparseDataWithBatch {
415 data: data_vec,
416 batch: batch_membership,
417 output_keep_idx,
418 })
419}
420
421fn empty_barcode_keep(data_vec: &SparseIoVec, exempt: Option<&[bool]>) -> Option<Vec<bool>> {
428 let mut missing_indptr: Vec<usize> = Vec::new();
429 let cutoffs: Vec<Option<u64>> = (0..data_vec.len())
430 .map(|b| {
431 if exempt.is_some_and(|e| e[b]) {
432 return None;
433 }
434 let backend = &data_vec[b];
435 let ncol = backend.num_columns().unwrap_or(0);
436 let nnz: Option<Vec<f32>> = (0..ncol)
437 .map(|c| backend.column_nnz(c).map(|x| x as f32))
438 .collect();
439 let Some(nnz) = nnz else {
440 missing_indptr.push(b);
441 return None;
442 };
443 crate::qc::suggest_nnz_cutoff(&nnz).map(|c| c as u64)
444 })
445 .collect();
446 if !missing_indptr.is_empty() {
447 warn!(
448 "Empty-barcode gate: file index(es) {} have no resident column indptr; \
449 skipping cell call for those backends",
450 missing_indptr
451 .iter()
452 .map(ToString::to_string)
453 .collect::<Vec<_>>()
454 .join(", "),
455 );
456 }
457 if cutoffs.iter().all(Option::is_none) {
458 return None;
459 }
460
461 let keep: Vec<bool> = (0..data_vec.num_columns())
462 .map(|c| {
463 data_vec.column_locations(c).iter().any(|loc| {
464 let b = loc.backend as usize;
465 cutoffs[b].is_none_or(|cut| {
466 data_vec[b]
467 .column_nnz(loc.local_col as usize)
468 .is_none_or(|x| x >= cut)
469 })
470 })
471 })
472 .collect();
473 let n_drop = keep.iter().filter(|&&k| !k).count();
474 info!(
475 "Empty-barcode gate: {} / {} columns called empty (per-file nnz cutoffs: {})",
476 n_drop,
477 keep.len(),
478 cutoffs
479 .iter()
480 .map(|c| c.map_or("none".to_string(), |x| x.to_string()))
481 .collect::<Vec<_>>()
482 .join(", "),
483 );
484 (n_drop > 0 && n_drop < keep.len()).then_some(keep)
486}
487
488fn maybe_warn_multimodal_pattern(data_vec: &SparseIoVec) {
492 let n_backends = data_vec.len();
493 if n_backends < 2 {
494 return;
495 }
496 let intersection = data_vec.num_rows_in_at_least(n_backends);
497 let min_backend_rows = (0..n_backends)
498 .map(|j| data_vec[j].num_rows().unwrap_or(0))
499 .min()
500 .unwrap_or(0);
501 if min_backend_rows == 0 {
502 return;
503 }
504 let disjointness = 1.0_f64 - (intersection as f64) / (min_backend_rows as f64);
505 if disjointness < MULTIMODAL_HINT_DISJOINTNESS_FRACTION {
506 return;
507 }
508
509 let mut shared: Option<FxHashSet<Box<str>>> = None;
513 for j in 0..n_backends {
514 let names = match data_vec[j].column_names() {
515 Ok(n) => n,
516 Err(_) => return, };
518 let set: FxHashSet<Box<str>> = names.into_iter().collect();
519 match shared.as_mut() {
520 None => shared = Some(set),
521 Some(prev) => {
522 prev.retain(|k| set.contains(k));
523 if prev.is_empty() {
524 return;
525 }
526 }
527 }
528 }
529 let shared_count = shared.map(|s| s.len()).unwrap_or(0);
530 if shared_count == 0 {
531 return;
532 }
533
534 warn!(
535 "Inputs look multi-modal-shaped (feature-axis disjointness {:.0}% across {} \
536 backends) and {} barcode(s) overlap across files. To glue cells across \
537 modalities, pass `--multiome` (or the equivalent ColumnAlignment::Union). \
538 Continuing with default Disjoint stacking — cells with shared barcodes \
539 will be treated as distinct.",
540 disjointness * 100.0,
541 n_backends,
542 shared_count
543 );
544}
545
546fn resolve_batch_disjoint(
548 data_files: &[Box<str>],
549 data_vec: &SparseIoVec,
550 batch_files: Option<&[Box<str>]>,
551) -> anyhow::Result<Vec<Box<str>>> {
552 let mut batch_membership: Vec<Box<str>> = Vec::with_capacity(data_vec.num_columns());
553
554 if let Some(batch_files) = batch_files {
555 if batch_files.len() != data_files.len() {
556 return Err(anyhow::anyhow!("# batch files != # of data files"));
557 }
558 let column_counts = data_vec.num_columns_by_data()?;
561 for (file_idx, batch_file) in batch_files.iter().enumerate() {
562 info!("Reading batch file: {}", batch_file);
563 let labels = read_lines(batch_file)?;
564 if labels.len() != column_counts[file_idx] {
565 return Err(anyhow::anyhow!(
566 "batch file {} has {} lines but data file {} has {} columns",
567 batch_file,
568 labels.len(),
569 data_files[file_idx],
570 column_counts[file_idx]
571 ));
572 }
573 batch_membership.extend(labels);
574 }
575 } else {
576 let column_counts = data_vec.num_columns_by_data()?;
577 let column_names = data_vec.column_names()?;
578 let mut col_start = 0usize;
579
580 for (file_idx, &ncols) in column_counts.iter().enumerate() {
581 let data_file = data_files[file_idx].clone();
582 let (_dir, file_base, _ext) = common_io::dir_base_ext(&data_file)?;
583 let col_end = col_start + ncols;
584 let file_columns = &column_names[col_start..col_end];
585
586 let appended_suffix = data_name_of(&data_file, data_files.len())?
587 .map(|name| format!("{COLUMN_SEP}{name}").into_boxed_str());
588 let (tags, used_embedded) = infer_batch_from_columns(
589 file_columns,
590 file_base.as_ref(),
591 appended_suffix.as_deref(),
592 );
593 if used_embedded {
594 info!(
595 "File {}: using embedded batch from column names (file '{}')",
596 file_idx, file_base
597 );
598 } else {
599 info!(
600 "File {}: using file name '{}' as batch",
601 file_idx, file_base
602 );
603 }
604 batch_membership.extend(tags);
605 col_start = col_end;
606 }
607 }
608 Ok(batch_membership)
609}
610
611fn resolve_batch_union(
625 data_vec: &SparseIoVec,
626 batch_files: Option<&[Box<str>]>,
627 n_cells: usize,
628) -> anyhow::Result<Vec<Box<str>>> {
629 if let Some(batch_files) = batch_files {
630 if batch_files.len() != 1 {
631 return Err(anyhow::anyhow!(
632 "Under ColumnAlignment::Union, --batch-files must have exactly one \
633 file listing one label per unified cell (got {} files for {} \
634 unified cells). A cell shared across modalities cannot carry two \
635 batch labels.",
636 batch_files.len(),
637 n_cells
638 ));
639 }
640 info!("Reading unified batch file: {}", batch_files[0]);
641 let labels: Vec<Box<str>> = read_lines(&batch_files[0])?;
642 if labels.len() != n_cells {
643 return Err(anyhow::anyhow!(
644 "Unified batch file {} has {} lines but data has {} unified cells",
645 batch_files[0],
646 labels.len(),
647 n_cells
648 ));
649 }
650 return Ok(labels);
651 }
652
653 let unified_names = data_vec.column_names()?;
661 let (per_cell_tags, used_embedded) = infer_batch_from_columns(&unified_names, "all", None);
664 if used_embedded {
665 info!(
666 "Union mode: per-cell @batch tag taken from unified barcodes ({} cells)",
667 n_cells
668 );
669 } else {
670 info!(
671 "No --batch-files and no embedded @batch tags — falling back to single \
672 batch 'all' (Union mode: per-file batch fallback is ambiguous)."
673 );
674 }
675 Ok(per_cell_tags)
676}
677
678fn infer_batch_from_columns(
693 file_columns: &[Box<str>],
694 file_base: &str,
695 appended_suffix: Option<&str>,
696) -> (Vec<Box<str>>, bool) {
697 fn raw_of<'a>(name: &'a str, suffix: Option<&str>) -> &'a str {
698 match suffix {
699 Some(sfx) => name.strip_suffix(sfx).unwrap_or(name),
700 None => name,
701 }
702 }
703
704 let mut used_embedded = false;
705 let tags = file_columns
706 .iter()
707 .map(|col_name| {
708 let raw = raw_of(col_name.as_ref(), appended_suffix);
709 match raw.rsplit_once('@') {
710 Some((_, embedded)) => {
711 used_embedded = true;
712 embedded.into()
713 }
714 None => file_base.into(),
715 }
716 })
717 .collect();
718 (tags, used_embedded)
719}
720
721#[cfg(test)]
722#[path = "data_loading_tests.rs"]
723mod data_loading_tests;
724
725#[cfg(test)]
726mod tests {
727 use super::*;
728
729 fn cols(v: &[&str]) -> Vec<Box<str>> {
730 v.iter()
731 .map(|s| (*s).to_string().into_boxed_str())
732 .collect()
733 }
734
735 #[test]
736 fn embedded_donor_survives_push_suffix() {
737 let names = cols(&[
741 "ACGT-1@donorA@mix",
742 "ACGT-2@donorB@mix",
743 "ACGT-3@donorA@mix",
744 "ACGT-4@donorB@mix",
745 ]);
746 let (tags, used_embedded) = infer_batch_from_columns(&names, "mix", Some("@mix"));
747 assert!(used_embedded);
748 assert_eq!(
749 tags.iter().map(|b| b.as_ref()).collect::<Vec<_>>(),
750 vec!["donorA", "donorB", "donorA", "donorB"]
751 );
752 }
753
754 #[test]
755 fn no_embedded_batch_falls_back_to_file_base() {
756 let names = cols(&["AAAA@s1", "CCCC@s1"]);
758 let (tags, used_embedded) = infer_batch_from_columns(&names, "s1", Some("@s1"));
759 assert!(!used_embedded);
760 assert_eq!(
761 tags.iter().map(|b| b.as_ref()).collect::<Vec<_>>(),
762 vec!["s1", "s1"]
763 );
764 }
765
766 #[test]
767 fn single_file_embedded_batch() {
768 let names = cols(&["ACGT-1@donorA", "ACGT-2@donorB"]);
770 let (tags, used_embedded) = infer_batch_from_columns(&names, "only", None);
771 assert!(used_embedded);
772 assert_eq!(
773 tags.iter().map(|b| b.as_ref()).collect::<Vec<_>>(),
774 vec!["donorA", "donorB"]
775 );
776 }
777
778 #[test]
779 fn single_file_no_embedded_batch() {
780 let names = cols(&["AAAA", "CCCC"]);
781 let (tags, used_embedded) = infer_batch_from_columns(&names, "only", None);
782 assert!(!used_embedded);
783 assert_eq!(
784 tags.iter().map(|b| b.as_ref()).collect::<Vec<_>>(),
785 vec!["only", "only"]
786 );
787 }
788
789 #[test]
790 fn each_name_is_read_for_its_own_tag() {
791 for (names, want) in [
794 (
795 cols(&["AAAA@donorA", "CCCC", "GGGG@donorB"]),
796 ["donorA", "all", "donorB"],
797 ),
798 (
799 cols(&["CCCC", "AAAA@donorA", "GGGG@donorB"]),
800 ["all", "donorA", "donorB"],
801 ),
802 ] {
803 let (tags, used_embedded) = infer_batch_from_columns(&names, "all", None);
804 assert!(used_embedded);
805 assert_eq!(tags.iter().map(|b| b.as_ref()).collect::<Vec<_>>(), want);
806 }
807 }
808
809 #[test]
810 fn empty_file_columns() {
811 let names: Vec<Box<str>> = vec![];
812 let (tags, used_embedded) = infer_batch_from_columns(&names, "x", Some("@x"));
813 assert!(!used_embedded);
814 assert!(tags.is_empty());
815 }
816}