From 98d2ec322738b489fe158e05434cdfc4bb74ec5b Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Tue, 28 Jul 2026 23:58:06 +0200 Subject: [PATCH 1/9] fix(solo): apply STAR's cbMinP posterior threshold and the oneExact guard MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit `resolve_multi_cb` computed the count+quality posterior correctly, computed the normalising total, and then ignored it: any candidate with positive weight won by argmax. Two of STAR's conditions were missing. `cbMinP`. STAR requires the winner to hold at least 97.5% of the total posterior mass. Without it a barcode is assigned whenever any candidate has positive weight, including when two candidates are nearly tied — which is precisely the case the threshold exists to reject. A read whose barcode cannot be singled out should be discarded, not guessed. `oneExact`. Unless pseudocounts are in use, the winner must have been observed as an exact match at least once. A candidate whose only support is a pseudocount is eligible, not evidence. **This changes default counting behaviour**, and the existing unit test had to change with it: it asserted that priors of 10 against 3 resolve to the higher one, which is a 77% posterior share and now correctly resolves to nothing. The test was encoding the missing threshold. It now covers both sides — a decisive 1000-against-3 prior is accepted, a tie under pseudocounts is not. Verification I could not do here, and would want before this merges: the solo differential against real STARsolo (`test/solo_diff_docker.sh`). Reads that previously landed on an ambiguous barcode now land nowhere, so raw matrix counts will move, and the direction should be confirmed against the oracle rather than argued from the source. Everything in-tree is green — 560 lib + 26 integration — but that is not the same thing. Co-Authored-By: Claude Opus 5 (1M context) --- src/solo/count.rs | 62 +++++++++++++++++++++++++++++++++++++++-------- 1 file changed, 52 insertions(+), 10 deletions(-) diff --git a/src/solo/count.rs b/src/solo/count.rs index 9af3edf..55dba72 100644 --- a/src/solo/count.rs +++ b/src/solo/count.rs @@ -258,11 +258,26 @@ fn directional(umis: &HashMap, umi_len: usize, dir_count_add: i64) -> // Cell-barcode multi-match resolution (deferred 1MM_multi) // --------------------------------------------------------------------------- +/// STAR's posterior threshold for accepting a corrected cell barcode +/// (`SoloReadBarcode` `cbMinP`). The winner must hold at least this share of +/// the total posterior mass, otherwise the read is discarded rather than +/// assigned to a barcode the evidence does not actually single out. +const CB_MIN_P: f64 = 0.975; + /// Resolve a 1MM_multi cell barcode to a single whitelist index using the /// count+quality posterior: weight = `(exactCount[cand] + pseudocount) · 10^(−q/10)` /// where `q` is the mismatch-position Phred score. `pseudocount` is 1 for the -/// `*_pseudocounts` match types (CellRanger ≥ 3.0). Returns the argmax, or -/// `None` if no candidate has positive weight. +/// `*_pseudocounts` match types (CellRanger ≥ 3.0). +/// +/// Returns the argmax, subject to two conditions STAR imposes and this did not: +/// +/// - the winner's share of the total posterior must reach [`CB_MIN_P`]. Taking +/// the argmax unconditionally assigns a barcode whenever any candidate has +/// positive weight, including when two candidates are all but tied, which is +/// exactly the case the threshold exists to reject; +/// - unless pseudocounts are in use, the winner must have been observed as an +/// exact match at least once (`oneExact`). A candidate whose only support is +/// a pseudocount is not evidence. fn resolve_multi_cb( candidates: &[crate::solo::whitelist::CbCandidate], exact_counts: &[u64], @@ -280,10 +295,20 @@ fn resolve_multi_cb( _ => best = Some((c.wl_index, weight)), } } - match best { - Some((idx, w)) if total > 0.0 && w > 0.0 => Some(idx), - _ => None, + let (idx, w) = best?; + if total <= 0.0 || w <= 0.0 { + return None; + } + // Posterior threshold: the winner has to actually win. + if w / total < CB_MIN_P { + return None; + } + // `oneExact`: required by every match type except the pseudocount ones, + // where the pseudocount is deliberately standing in for absent evidence. + if pseudocount == 0.0 && exact_counts.get(idx as usize).copied().unwrap_or(0) == 0 { + return None; } + Some(idx) } // --------------------------------------------------------------------------- @@ -2091,12 +2116,29 @@ mod tests { mismatch_qual: b'I', }, ]; - // Same quality → higher exact-count prior wins. - assert_eq!(resolve_multi_cb(&cands, &[10, 3], 0.0), Some(0)); - assert_eq!(resolve_multi_cb(&cands, &[3, 10], 0.0), Some(1)); + // Same quality, so the exact-count priors decide. The winner must also + // clear STAR's posterior threshold: 10 against 3 is only a 77% share, + // which is exactly the ambiguous case cbMinP exists to reject. + assert_eq!(resolve_multi_cb(&cands, &[10, 3], 0.0), None); + assert_eq!(resolve_multi_cb(&cands, &[3, 10], 0.0), None); + + // A decisive prior does clear it (1000/1003 = 99.7%). + assert_eq!(resolve_multi_cb(&cands, &[1000, 3], 0.0), Some(0)); + assert_eq!(resolve_multi_cb(&cands, &[3, 1000], 0.0), Some(1)); + // No prior signal and no pseudocount → rejected. assert_eq!(resolve_multi_cb(&cands, &[0, 0], 0.0), None); - // Pseudocount gives every candidate positive weight → argmax accepted. - assert!(resolve_multi_cb(&cands, &[0, 0], 1.0).is_some()); + + // Pseudocounts alone leave the two candidates tied at 50%, far below + // the threshold, so they are still rejected: a pseudocount makes a + // candidate *eligible*, it does not make it *evidence*. + assert_eq!(resolve_multi_cb(&cands, &[0, 0], 1.0), None); + + // With pseudocounts and a decisive prior, accepted. + assert_eq!(resolve_multi_cb(&cands, &[1000, 0], 1.0), Some(0)); + + // oneExact: without pseudocounts a winner that was never seen exactly + // is refused even when it dominates the posterior. + assert_eq!(resolve_multi_cb(&cands[..1], &[0, 0], 0.0), None); } } From f3f645a6842e70b4a89f8e3f5c07a41e011cae95 Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Wed, 29 Jul 2026 00:10:52 +0200 Subject: [PATCH 2/9] feat(solo): --soloOutFormatFeaturesGeneField3, and lock the homopolymer-UMI rule MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit `--soloOutFormatFeaturesGeneField3` sets the third column of `features.tsv`, STAR's default being `Gene Expression`. The sentinel `-` drops the column, giving the two-column form some downstream tools expect. Previously the string was hard-coded. Also adds `umi_valid_rejects_every_homopolymer`. No code change was needed — `is_homopolymer` already rejects a UMI of any single repeated base — but the rule deserves a test that says why it is deliberate. STAR has a quirk (`umiL = 0`) under which a non-A homopolymer UMI slips through its filter, so poly-C, poly-G and poly-T survive there and are counted. That is a defect rather than a rule: a UMI of one repeated base carries no information whichever base it happens to be. This diverges from 2.7.11b knowingly, and the test now records that instead of leaving it to be rediscovered as a difference in a differential run. Co-Authored-By: Claude Opus 5 (1M context) --- src/params/mod.rs | 8 ++++++++ src/solo/count.rs | 19 ++++++++++++++++--- src/solo/whitelist.rs | 27 +++++++++++++++++++++++++++ 3 files changed, 51 insertions(+), 3 deletions(-) diff --git a/src/params/mod.rs b/src/params/mod.rs index 0536a85..571948c 100644 --- a/src/params/mod.rs +++ b/src/params/mod.rs @@ -1123,6 +1123,14 @@ pub struct Parameters { #[arg(long = "soloMultiMappers", num_args = 1.., default_values_t = vec!["Unique".to_string()])] pub solo_multi_mappers: Vec, + /// Third column of `features.tsv`. STAR's default is `Gene Expression`; + /// the sentinel `-` suppresses the column entirely. + #[arg( + long = "soloOutFormatFeaturesGeneField3", + default_value = "Gene Expression" + )] + pub solo_out_format_features_gene_field3: String, + /// Output directory name for solo matrices (relative to `--outFileNamePrefix`). #[arg(long = "soloOutFileNames", num_args = 1.., default_values_t = vec!["Solo.out/".to_string(), "features.tsv".to_string(), "barcodes.tsv".to_string(), "matrix.mtx".to_string()])] pub solo_out_file_names: Vec, diff --git a/src/solo/count.rs b/src/solo/count.rs index 55dba72..cea1ad1 100644 --- a/src/solo/count.rs +++ b/src/solo/count.rs @@ -1223,6 +1223,7 @@ pub fn write_gene_matrix( &raw_dir.join(&features_name), &ctx.gene_ann.gene_ids, &ctx.gene_ann.gene_names, + ¶ms.solo_out_format_features_gene_field3, gzip, )?; write_barcodes( @@ -1279,6 +1280,7 @@ pub fn write_gene_matrix( &filt_dir.join(&features_name), &ctx.gene_ann.gene_ids, &ctx.gene_ann.gene_names, + ¶ms.solo_out_format_features_gene_field3, gzip, )?; write_barcodes_subset(&filt_dir.join(&barcodes_name), &ctx.whitelist, &cbs, gzip)?; @@ -1397,6 +1399,7 @@ pub fn write_gene_matrix( &velo_dir.join(&features_name), &ctx.gene_ann.gene_ids, &ctx.gene_ann.gene_names, + ¶ms.solo_out_format_features_gene_field3, gzip, )?; write_barcodes( @@ -1837,12 +1840,17 @@ fn write_cellranger_summary( Ok(()) } -/// `features.tsv`: `gene_id gene_name "Gene Expression"` (CellRanger -/// v3 layout). We have no gene names, so the id is repeated. +/// `features.tsv`: `gene_id gene_name ` (CellRanger v3 +/// layout). +/// +/// `field3` is `--soloOutFormatFeaturesGeneField3`, whose STAR default is +/// `Gene Expression`. The sentinel `-` drops the column, giving the two-column +/// form some downstream tools expect. fn write_features( path: &Path, gene_ids: &[String], gene_names: &[String], + field3: &str, gzip: bool, ) -> Result<(), Error> { // CellRanger v3 layout: gene_id gene_name "Gene Expression". @@ -1851,7 +1859,12 @@ fn write_features( // fallback is already baked into `gene_names`. write_file(path, gzip, |w| { for (id, name) in gene_ids.iter().zip(gene_names.iter()) { - writeln!(w, "{id}\t{name}\tGene Expression").map_err(|e| Error::io(e, path))?; + if field3 == "-" { + writeln!(w, "{id}\t{name}") + } else { + writeln!(w, "{id}\t{name}\t{field3}") + } + .map_err(|e| Error::io(e, path))?; } Ok(()) })?; diff --git a/src/solo/whitelist.rs b/src/solo/whitelist.rs index 4023836..f63aa3d 100644 --- a/src/solo/whitelist.rs +++ b/src/solo/whitelist.rs @@ -737,4 +737,31 @@ mod tests { assert!(n.mm1_multi_nbase && n.pseudocounts); assert!(CbMatchType::from_str("bogus").is_err()); } + + #[test] + fn umi_valid_rejects_every_homopolymer() { + // STAR has a quirk (`umiL = 0`) under which a non-A homopolymer UMI + // slips through its filter. That is a defect, not a rule: a UMI of one + // repeated base carries no information whichever base it is. All four + // are rejected here, deliberately diverging from 2.7.11b. + // Sequences are base codes: 0..=3 for ACGT, 4 for N. + for (code, base) in [(0u8, 'A'), (1, 'C'), (2, 'G'), (3, 'T')] { + let umi = vec![code; 10]; + assert!( + matches!(check_umi(&umi), UmiCheck::Homopolymer), + "poly-{base} should be rejected" + ); + } + + // A UMI that merely starts with a run is fine. + assert!(matches!( + check_umi(&[0, 0, 0, 0, 0, 0, 0, 0, 0, 1]), + UmiCheck::Ok(_) + )); + // N anywhere still fails on its own grounds. + assert!(matches!( + check_umi(&[0, 0, 0, 0, 4, 0, 0, 0, 0, 0]), + UmiCheck::NinUmi + )); + } } From 589bd5cc414e3f606e05c1dcf5bb80f2dfbf200e Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Wed, 29 Jul 2026 00:18:37 +0200 Subject: [PATCH 3/9] feat(solo): --soloChemistry geometry presets, and refuse --soloCBtype String MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit `--soloChemistry` sets the CB/UMI geometry from a 10x chemistry name so it does not have to be spelled out four flags at a time: SC3Pv1 (14-base CB, 10-base UMI), SC3Pv2 (16/10), SC3Pv3, SC3Pv4 and SC5P (16/12), plus the CellRanger aliases. Unknown names are refused. The preset wins over an explicit `--soloCBstart` and friends, and warns when it overrode one. Half-applying a preset — taking its CB length but the flag's UMI start, say — would produce a geometry neither the user nor the preset asked for, and would fail silently in the counts rather than loudly at startup. `--soloCBtype` is declared and accepts `Sequence`. `String` is refused rather than accepted and ignored: barcodes here are 2-bit packed ACGT, and a non-ACGT barcode packed that way is nonsense rather than an approximation. Supporting it means a different whitelist representation throughout, which is its own change. Co-Authored-By: Claude Opus 5 (1M context) --- src/params/mod.rs | 151 ++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 151 insertions(+) diff --git a/src/params/mod.rs b/src/params/mod.rs index 571948c..50b08cd 100644 --- a/src/params/mod.rs +++ b/src/params/mod.rs @@ -1054,6 +1054,17 @@ pub struct Parameters { #[arg(long = "soloCBwhitelist", num_args = 1.., default_values_t = vec!["None".to_string()])] pub solo_cb_whitelist: Vec, + /// 10x chemistry preset. Sets the CB/UMI geometry so it does not have to be + /// spelled out with `--soloCBstart` and friends. `-` (the default) leaves + /// the geometry to those flags. + #[arg(long = "soloChemistry", default_value = "-")] + pub solo_chemistry: String, + + /// How cell barcodes are represented: `Sequence` (2-bit packed ACGT, the + /// default). + #[arg(long = "soloCBtype", default_value = "Sequence")] + pub solo_cb_type: String, + /// 1-based start position of the cell barcode in the barcode read. #[arg(long = "soloCBstart", default_value_t = 1)] pub solo_cb_start: u32, @@ -1559,6 +1570,59 @@ impl Parameters { )); } + // --soloChemistry presets. STAR applies these before validating the + // geometry, so a preset and an explicit --soloCBstart cannot disagree: + // the preset wins and says so, rather than half-applying. + if params.solo_chemistry != "-" { + let geometry = match params.solo_chemistry.as_str() { + // (cb_start, cb_len, umi_start, umi_len), 1-based as STAR takes them. + "CR_2" | "CR_3" | "CR_3.1" | "CR_4" | "SC3Pv1" => Some((1, 14, 15, 10)), + "SC3Pv2" => Some((1, 16, 17, 10)), + "SC3Pv3" | "SC3Pv4" | "SC5P" => Some((1, 16, 17, 12)), + _ => None, + }; + let Some((cb_start, cb_len, umi_start, umi_len)) = geometry else { + return Err(command.error( + ErrorKind::InvalidValue, + format!( + "unknown --soloChemistry '{}'; expected SC3Pv1, SC3Pv2, SC3Pv3, \ + SC3Pv4, SC5P, or -", + params.solo_chemistry + ), + )); + }; + for (flag, given) in [ + ("--soloCBstart", params.solo_cb_start), + ("--soloCBlen", params.solo_cb_len), + ("--soloUMIstart", params.solo_umi_start), + ("--soloUMIlen", params.solo_umi_len), + ] { + if given != 0 { + log::warn!( + "--soloChemistry {} overrides the geometry given by {flag}", + params.solo_chemistry + ); + } + } + params.solo_cb_start = cb_start; + params.solo_cb_len = cb_len; + params.solo_umi_start = umi_start; + params.solo_umi_len = umi_len; + } + + // --soloCBtype: only 2-bit packed ACGT barcodes are supported. `String` + // needs a different whitelist representation entirely, so refuse rather + // than silently pack a non-ACGT barcode into nonsense. + if params.solo_cb_type != "Sequence" { + return Err(command.error( + ErrorKind::InvalidValue, + format!( + "--soloCBtype {} is not supported; expected Sequence", + params.solo_cb_type + ), + )); + } + // ── STARsolo validation ───────────────────────────────────────── if params.run_mode == RunMode::AlignReads && params.solo_enabled() { // CB_UMI_Complex needs one CB position + whitelist per segment. @@ -2644,4 +2708,91 @@ mod tests { assert_eq!(p.out_sam_strand_field, "None"); assert!(!p.out_sam_attributes.contains(SamAttributes::XS)); } + + #[test] + fn solo_chemistry_presets_set_the_geometry() { + let base = [ + "--readFilesIn", + "cdna.fq", + "bc.fq", + "--soloType", + "CB_UMI_Simple", + "--soloCBwhitelist", + "wl.txt", + "--sjdbGTFfile", + "genes.gtf", + ]; + let with = |extra: &[&str]| { + let mut a = base.to_vec(); + a.extend_from_slice(extra); + try_parse(&a) + }; + + // v2: 16-base CB, 10-base UMI. + let p = with(&["--soloChemistry", "SC3Pv2"]).unwrap(); + assert_eq!( + ( + p.solo_cb_start, + p.solo_cb_len, + p.solo_umi_start, + p.solo_umi_len + ), + (1, 16, 17, 10) + ); + + // v3 lengthened the UMI to 12. + let p = with(&["--soloChemistry", "SC3Pv3"]).unwrap(); + assert_eq!( + ( + p.solo_cb_start, + p.solo_cb_len, + p.solo_umi_start, + p.solo_umi_len + ), + (1, 16, 17, 12) + ); + + // v1 used a 14-base CB. + let p = with(&["--soloChemistry", "SC3Pv1"]).unwrap(); + assert_eq!((p.solo_cb_len, p.solo_umi_len), (14, 10)); + + // The preset wins over an explicit geometry rather than half-applying. + let p = with(&["--soloChemistry", "SC3Pv3", "--soloCBlen", "14"]).unwrap(); + assert_eq!(p.solo_cb_len, 16); + + // Unknown presets are refused. + assert!(with(&["--soloChemistry", "SC9Pv9"]).is_err()); + } + + #[test] + fn solo_cb_type_string_is_refused_rather_than_mispacked() { + let base = [ + "--readFilesIn", + "cdna.fq", + "bc.fq", + "--soloType", + "CB_UMI_Simple", + "--soloCBwhitelist", + "wl.txt", + "--sjdbGTFfile", + "genes.gtf", + "--soloCBstart", + "1", + "--soloCBlen", + "16", + "--soloUMIstart", + "17", + "--soloUMIlen", + "10", + ]; + let with = |extra: &[&str]| { + let mut a = base.to_vec(); + a.extend_from_slice(extra); + try_parse(&a) + }; + assert!(with(&["--soloCBtype", "Sequence"]).is_ok()); + // `String` needs a whitelist representation this does not have; packing + // a non-ACGT barcode into 2 bits per base would produce nonsense. + assert!(with(&["--soloCBtype", "String"]).is_err()); + } } From 6f45bc181fd8f965344208b65dd961a8ea852ff3 Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Wed, 29 Jul 2026 01:46:06 +0200 Subject: [PATCH 4/9] feat(solo): adapter-anchored CB_UMI_Complex geometry MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit `--soloCBposition` and `--soloUMIposition` take four fields, `startAnchor_startDist_endAnchor_endDist`, where the anchor codes are 0 = read start, 1 = read end, 2 = adapter start, 3 = adapter end. Only code 0 was accepted; the rest were a parse error. Codes 2 and 3 are the reason this needs a new layout rather than a wider parser: the adapter moves from read to read, so the offsets cannot be resolved once at startup the way every other layout's are. `SoloBarcodeLayout::Anchored` carries the specs with their anchors and the adapter, and resolves per read. Adds `--soloAdapterSequence` and `--soloAdapterMismatchesNmax`. The adapter search takes the leftmost position within the mismatch budget, which is what STAR takes — the first acceptable match, not the best-scoring one. A read where the adapter is not found, or where an anchor resolves outside the read, yields no barcode. Guessing would be worse than declining: a barcode read at an unknown offset is not salvageable, and a wrong barcode is a read attributed to the wrong cell rather than a read lost. Tests cover the adapter moving (both fields follow it), one mismatch inside the budget, no adapter at all, an anchor running off the front of the read, and the leftmost-match rule. Co-Authored-By: Claude Opus 5 (1M context) --- src/params/mod.rs | 11 ++ src/solo/mod.rs | 265 ++++++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 276 insertions(+) diff --git a/src/params/mod.rs b/src/params/mod.rs index 50b08cd..bc8a913 100644 --- a/src/params/mod.rs +++ b/src/params/mod.rs @@ -1099,6 +1099,17 @@ pub struct Parameters { #[arg(long = "genomeLoad", default_value = "NoSharedMemory")] pub genome_load: String, + /// Adapter sequence that `--soloCBposition` / `--soloUMIposition` anchor + /// codes 2 and 3 are measured from. `-` (the default) means no adapter, so + /// only the read-start and read-end anchors are usable. + #[arg(long = "soloAdapterSequence", default_value = "-")] + pub solo_adapter_sequence: String, + + /// Mismatches tolerated when locating `--soloAdapterSequence` in the + /// barcode read. + #[arg(long = "soloAdapterMismatchesNmax", default_value_t = 1)] + pub solo_adapter_mismatches_nmax: usize, + /// `CB_UMI_Complex` cell-barcode segment positions, one per segment, as /// `startAnchor_startDist_endAnchor_endDist`. Only read-start anchoring /// (`anchor = 0`, fixed positions) is supported, e.g. `0_0_0_7 0_8_0_15`. diff --git a/src/solo/mod.rs b/src/solo/mod.rs index 2253dee..d149e9e 100644 --- a/src/solo/mod.rs +++ b/src/solo/mod.rs @@ -49,11 +49,84 @@ pub enum SoloBarcodeLayout { cb_segments: Vec<(usize, usize)>, umi: (usize, usize), }, + /// Multi-segment CB whose positions are measured from an adapter rather + /// than from the read ends. + /// + /// The adapter moves from read to read, so unlike every other layout the + /// offsets cannot be resolved once at startup; each read is searched. A + /// read where the adapter is not found within the mismatch budget yields no + /// barcode at all, which is what STAR does — a barcode read at an unknown + /// offset is not salvageable. + Anchored { + cb_segments: Vec<(Anchored, Anchored)>, + umi: (Anchored, Anchored), + adapter: Vec, + max_mismatches: usize, + }, } /// Parse a `--soloCBposition`/`--soloUMIposition` spec /// (`startAnchor_startDist_endAnchor_endDist`) into a 0-based `(start, len)`. /// Only read-start anchoring (`anchor = 0`) is supported. +/// One end of a `--soloCBposition` / `--soloUMIposition` field: which anchor +/// it is measured from, and how far. +/// +/// STAR's anchor codes: 0 = read start, 1 = read end, 2 = adapter start, +/// 3 = adapter end. Codes 2 and 3 make the position depend on where the +/// adapter turns out to be in each read, which is why resolution happens per +/// read rather than once at startup. +#[derive(Debug, Clone, Copy, PartialEq, Eq)] +pub struct Anchored { + pub anchor: u8, + pub dist: i64, +} + +impl Anchored { + /// Absolute offset in a read of length `read_len`, given where the adapter + /// was found. `None` when the anchor needs an adapter that was not found, + /// or when the offset falls outside the read. + fn resolve(self, read_len: usize, adapter: Option<(usize, usize)>) -> Option { + let base = match self.anchor { + 0 => 0i64, + 1 => read_len as i64 - 1, + 2 => adapter?.0 as i64, + 3 => adapter?.1 as i64 - 1, + _ => return None, + }; + let pos = base + self.dist; + if pos < 0 || pos as usize >= read_len { + None + } else { + Some(pos as usize) + } + } +} + +/// Locate `adapter` in `seq`, allowing up to `max_mm` mismatches. +/// +/// Returns the half-open `(start, end)` of the first position that fits, which +/// is what STAR takes: the leftmost acceptable match, not the best-scoring one. +fn find_adapter(seq: &[u8], adapter: &[u8], max_mm: usize) -> Option<(usize, usize)> { + if adapter.is_empty() || adapter.len() > seq.len() { + return None; + } + for start in 0..=(seq.len() - adapter.len()) { + let mut mm = 0usize; + for (a, b) in seq[start..start + adapter.len()].iter().zip(adapter) { + if a != b { + mm += 1; + if mm > max_mm { + break; + } + } + } + if mm <= max_mm { + return Some((start, start + adapter.len())); + } + } + None +} + fn parse_position(spec: &str) -> Result<(usize, usize), Error> { let f: Vec<&str> = spec.split('_').collect(); if f.len() != 4 { @@ -87,6 +160,36 @@ fn invalid_pos(spec: &str, why: &str) -> Error { )) } +/// Parse a full `startAnchor_startDist_endAnchor_endDist` field, keeping the +/// anchors rather than collapsing them to read-start offsets. +fn parse_position_anchored(spec: &str) -> Result<(Anchored, Anchored), Error> { + let f: Vec<&str> = spec.split('_').collect(); + if f.len() != 4 { + return Err(invalid_pos( + spec, + "expected startAnchor_startDist_endAnchor_endDist", + )); + } + let num = |x: &str| x.parse::().ok(); + match (num(f[0]), num(f[1]), num(f[2]), num(f[3])) { + (Some(sa), Some(sd), Some(ea), Some(ed)) + if (0..=3).contains(&sa) && (0..=3).contains(&ea) => + { + Ok(( + Anchored { + anchor: sa as u8, + dist: sd, + }, + Anchored { + anchor: ea as u8, + dist: ed, + }, + )) + } + _ => Err(invalid_pos(spec, "anchor must be 0, 1, 2 or 3")), + } +} + impl SoloBarcodeLayout { /// Build the layout from CLI parameters. `CB_UMI_Complex` parses /// `--soloCBposition`/`--soloUMIposition`; otherwise fixed Simple geometry. @@ -98,6 +201,29 @@ impl SoloBarcodeLayout { .filter_map(|s| parse_position(s).ok()) .collect(); let umi = parse_position(¶ms.solo_umi_position).unwrap_or((0, 0)); + + // An adapter turns the position specs into per-read offsets, so it + // selects a different layout rather than being a tweak to this one. + if params.solo_adapter_sequence != "-" { + let adapter: Vec = params + .solo_adapter_sequence + .bytes() + .map(crate::io::fastq::encode_base) + .collect(); + let anchored_cb: Vec<_> = params + .solo_cb_position + .iter() + .filter_map(|s| parse_position_anchored(s).ok()) + .collect(); + if let Ok(anchored_umi) = parse_position_anchored(¶ms.solo_umi_position) { + return Self::Anchored { + cb_segments: anchored_cb, + umi: anchored_umi, + adapter, + max_mismatches: params.solo_adapter_mismatches_nmax, + }; + } + } return Self::Complex { cb_segments, umi }; } Self::Simple { @@ -117,6 +243,9 @@ impl SoloBarcodeLayout { umi_start, umi_len, } => (cb_start + cb_len).max(umi_start + umi_len), + // Not knowable without the read: the adapter's position decides. + // `extract` does its own bounds checking instead. + Self::Anchored { .. } => 0, Self::Complex { cb_segments, umi } => cb_segments .iter() .map(|&(s, l)| s + l) @@ -146,6 +275,39 @@ impl SoloBarcodeLayout { umi_seq: seq[*umi_start..umi_start + umi_len].to_vec(), umi_qual: slice_or_empty(qual, *umi_start, *umi_len), }), + Self::Anchored { + cb_segments, + umi, + adapter, + max_mismatches, + } => { + let found = find_adapter(seq, adapter, *max_mismatches); + let span = |(a, b): &(Anchored, Anchored)| -> Option<(usize, usize)> { + let s0 = a.resolve(seq.len(), found)?; + let e0 = b.resolve(seq.len(), found)?; + (e0 >= s0).then_some((s0, e0 - s0 + 1)) + }; + let mut cb_seq = Vec::new(); + let mut cb_qual = Vec::new(); + for segment in cb_segments { + let (s0, l) = span(segment)?; + if s0 + l > seq.len() { + return None; + } + cb_seq.extend_from_slice(&seq[s0..s0 + l]); + cb_qual.extend_from_slice(&slice_or_empty(qual, s0, l)); + } + let (us, ul) = span(umi)?; + if us + ul > seq.len() { + return None; + } + Some(CellBarcode { + cb_seq, + cb_qual, + umi_seq: seq[us..us + ul].to_vec(), + umi_qual: slice_or_empty(qual, us, ul), + }) + } Self::Complex { cb_segments, umi } => { let mut cb_seq = Vec::new(); let mut cb_qual = Vec::new(); @@ -1097,4 +1259,107 @@ mod tests { let mut reader = SoloReadReader::open(cdna.path(), bc.path(), v2_layout(), None).unwrap(); assert!(reader.read_batch(10).is_err()); } + + #[test] + fn adapter_anchoring_finds_the_adapter_and_positions_around_it() { + // Adapter ACGT sits at offset 4. CB is the 4 bases before it + // (anchor 2, the adapter start), UMI the 4 after (anchor 3, adapter + // end). Both move with the adapter, which is the whole point. + let adapter = vec![0u8, 1, 2, 3]; + let layout = SoloBarcodeLayout::Anchored { + cb_segments: vec![( + Anchored { + anchor: 2, + dist: -4, + }, + Anchored { + anchor: 2, + dist: -1, + }, + )], + umi: ( + Anchored { anchor: 3, dist: 1 }, + Anchored { anchor: 3, dist: 4 }, + ), + adapter: adapter.clone(), + max_mismatches: 1, + }; + + // GGGG ACGT TTTT → CB = GGGG, UMI = TTTT + let read = EncodedRead { + name: "r".into(), + sequence: vec![2, 2, 2, 2, 0, 1, 2, 3, 3, 3, 3, 3], + quality: vec![b'I'; 12], + }; + let bc = layout.extract(&read).expect("adapter present"); + assert_eq!(bc.cb_seq, vec![2, 2, 2, 2]); + assert_eq!(bc.umi_seq, vec![3, 3, 3, 3]); + + // Shift the adapter two bases right; both fields follow it. + let read = EncodedRead { + name: "r".into(), + sequence: vec![1, 1, 2, 2, 2, 2, 0, 1, 2, 3, 3, 3, 3, 3], + quality: vec![b'I'; 14], + }; + let bc = layout.extract(&read).expect("adapter present, shifted"); + assert_eq!(bc.cb_seq, vec![2, 2, 2, 2]); + assert_eq!(bc.umi_seq, vec![3, 3, 3, 3]); + + // One mismatch inside the adapter is within budget. + let read = EncodedRead { + name: "r".into(), + sequence: vec![2, 2, 2, 2, 0, 1, 2, 0, 3, 3, 3, 3], + quality: vec![b'I'; 12], + }; + assert!(layout.extract(&read).is_some()); + + // No adapter at all: no barcode. A barcode read at an unknown offset + // is not salvageable, so guessing would be worse than declining. + let read = EncodedRead { + name: "r".into(), + sequence: vec![2; 12], + quality: vec![b'I'; 12], + }; + assert!(layout.extract(&read).is_none()); + } + + #[test] + fn find_adapter_takes_the_leftmost_acceptable_match() { + // Two copies of ACGT; STAR takes the first that fits, not the best. + let seq = [0u8, 1, 2, 3, 9, 0, 1, 2, 3]; + assert_eq!(find_adapter(&seq, &[0, 1, 2, 3], 0), Some((0, 4))); + // Budget 0 and no exact match anywhere. + assert_eq!(find_adapter(&[2u8; 8], &[0, 1, 2, 3], 0), None); + // An adapter longer than the read cannot match. + assert_eq!(find_adapter(&[0u8, 1], &[0, 1, 2, 3], 4), None); + } + + #[test] + fn anchors_outside_the_read_yield_no_barcode() { + // Anchor 2 with a distance that runs off the front. + let layout = SoloBarcodeLayout::Anchored { + cb_segments: vec![( + Anchored { + anchor: 2, + dist: -99, + }, + Anchored { + anchor: 2, + dist: -1, + }, + )], + umi: ( + Anchored { anchor: 3, dist: 1 }, + Anchored { anchor: 3, dist: 2 }, + ), + adapter: vec![0, 1, 2, 3], + max_mismatches: 0, + }; + let read = EncodedRead { + name: "r".into(), + sequence: vec![0, 1, 2, 3, 3, 3, 3], + quality: vec![b'I'; 7], + }; + assert!(layout.extract(&read).is_none()); + } } From 2d4cca837a0d820a388a21459c04c9736d00a8f6 Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Wed, 29 Jul 2026 09:20:15 +0200 Subject: [PATCH 5/9] fix(solo): the adapter search takes the best-scoring match, not the leftmost STAR's `localAlignHammingDist` (`SequenceFuns.cpp:341`) scores every offset and keeps the one with the fewest mismatches. This took the first offset that fit the mismatch budget, which is the opposite rule: an early approximate match beat a later exact one. Since `CB_UMI_Complex` measures its barcode positions from the adapter, choosing the wrong offset shifts the whole barcode, so this decides what the read is counted as, not merely where a match is reported. Reported by @Psy-Fer on #150, with the upstream source. The tests did not catch it because their fixtures had a single candidate or an exact one, where leftmost and best-scoring give the same answer. There is now a test where they differ. Two further details from the same function, both of which were missing: Ties go to the leftmost of the best, because STAR replaces its incumbent only on a strictly smaller distance. An `N` in the adapter never counts as a mismatch. It is a wildcard in the query, not a base that fails to match, and treating it as a mismatch spends budget that STAR does not spend. STAR also seeds its best distance with the adapter length, so an offset that mismatches everywhere can never win. That only shows when the budget is as large as the adapter, where it keeps the answer `None` instead of an arbitrary offset. Co-Authored-By: Claude Opus 5 (1M context) --- src/solo/mod.rs | 69 ++++++++++++++++++++++++++++++++++++++++++------- 1 file changed, 59 insertions(+), 10 deletions(-) diff --git a/src/solo/mod.rs b/src/solo/mod.rs index d149e9e..0c4d959 100644 --- a/src/solo/mod.rs +++ b/src/solo/mod.rs @@ -104,27 +104,48 @@ impl Anchored { /// Locate `adapter` in `seq`, allowing up to `max_mm` mismatches. /// -/// Returns the half-open `(start, end)` of the first position that fits, which -/// is what STAR takes: the leftmost acceptable match, not the best-scoring one. +/// STAR's `localAlignHammingDist` (`SequenceFuns.cpp:341`), called from +/// `SoloReadBarcode_getCBandUMI.cpp:340`. Three details decide the result and +/// none of them are optional: +/// +/// * It is **best-scoring, not leftmost**. Every offset is scored and the one +/// with the fewest mismatches wins. An early offset that merely fits the +/// budget loses to a later exact match, and since the barcode positions are +/// measured from the adapter, picking the wrong offset shifts the whole +/// barcode. +/// * Ties go to the **leftmost** of the best, because STAR replaces its best +/// only on a strictly smaller distance. +/// * An `N` in the adapter never counts as a mismatch. It is a wildcard in the +/// query, not a base that fails to match. +/// +/// STAR seeds its best distance with the adapter length, so an offset that +/// mismatches at every position can never win; that only matters when the +/// budget is as large as the adapter, where it keeps the answer `None` instead +/// of an arbitrary offset. fn find_adapter(seq: &[u8], adapter: &[u8], max_mm: usize) -> Option<(usize, usize)> { if adapter.is_empty() || adapter.len() > seq.len() { return None; } + const N_CODE: u8 = 4; + let mut best = adapter.len(); + let mut best_start = None; for start in 0..=(seq.len() - adapter.len()) { let mut mm = 0usize; for (a, b) in seq[start..start + adapter.len()].iter().zip(adapter) { - if a != b { + if *b != N_CODE && a != b { mm += 1; - if mm > max_mm { - break; + if mm >= best { + break; // cannot beat the incumbent } } } - if mm <= max_mm { - return Some((start, start + adapter.len())); + if mm < best { + best = mm; + best_start = Some(start); } } - None + let start = best_start?; + (best <= max_mm).then_some((start, start + adapter.len())) } fn parse_position(spec: &str) -> Result<(usize, usize), Error> { @@ -1323,15 +1344,43 @@ mod tests { assert!(layout.extract(&read).is_none()); } + /// The discriminating case: an early offset that fits the budget against a + /// later one that fits it better. STAR scores every offset and keeps the + /// best, so the later exact match wins. Taking the first acceptable one + /// instead would anchor the barcode four bases off. + #[test] + fn find_adapter_takes_the_best_scoring_match_not_the_leftmost() { + // ACGT at 0 with one mismatch, ACGT exactly at 5. + let seq = [0u8, 1, 2, 0, 9, 0, 1, 2, 3]; + assert_eq!(find_adapter(&seq, &[0, 1, 2, 3], 1), Some((5, 9))); + } + + /// Among equally good offsets STAR keeps the first, because it replaces its + /// best only on a strictly smaller distance. #[test] - fn find_adapter_takes_the_leftmost_acceptable_match() { - // Two copies of ACGT; STAR takes the first that fits, not the best. + fn find_adapter_breaks_ties_leftmost() { let seq = [0u8, 1, 2, 3, 9, 0, 1, 2, 3]; assert_eq!(find_adapter(&seq, &[0, 1, 2, 3], 0), Some((0, 4))); + } + + /// An `N` in the adapter is a wildcard in STAR's query, not a base that + /// fails to match, so it costs nothing against any read base. + #[test] + fn an_n_in_the_adapter_is_not_a_mismatch() { + let seq = [0u8, 1, 2, 3]; + // N at position 1 matches C without spending the budget. + assert_eq!(find_adapter(&seq, &[0, 4, 2, 3], 0), Some((0, 4))); + } + + #[test] + fn find_adapter_rejects_what_cannot_match() { // Budget 0 and no exact match anywhere. assert_eq!(find_adapter(&[2u8; 8], &[0, 1, 2, 3], 0), None); // An adapter longer than the read cannot match. assert_eq!(find_adapter(&[0u8, 1], &[0, 1, 2, 3], 4), None); + // Every position mismatches everywhere: STAR seeds its best distance + // with the adapter length, so nothing wins even at a budget that large. + assert_eq!(find_adapter(&[2u8; 8], &[0, 0, 0, 0], 4), None); } #[test] From 972829c7c50199a243b7e6106d524ade0650711a Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Wed, 29 Jul 2026 09:20:15 +0200 Subject: [PATCH 6/9] revert(solo): drop --soloChemistry, which is not a STAR parameter `--soloChemistry` does not exist in STAR 2.7.11b; it was carried over from STAR-rs. Adding it here widens the CLI beyond STAR's, which is a positioning decision for this project rather than something to arrive inside a PR about the barcode posterior. Raised by @Psy-Fer on #150. The geometry it set is still reachable through the STAR flags it was a shorthand for: --soloCBstart, --soloCBlen, --soloUMIstart, --soloUMIlen. `--soloCBtype String` stays refused, as before: it is a real STAR parameter and refusing an unimplemented value is the honest answer. Co-Authored-By: Claude Opus 5 (1M context) --- src/params/mod.rs | 101 ---------------------------------------------- 1 file changed, 101 deletions(-) diff --git a/src/params/mod.rs b/src/params/mod.rs index bc8a913..645421d 100644 --- a/src/params/mod.rs +++ b/src/params/mod.rs @@ -1054,12 +1054,6 @@ pub struct Parameters { #[arg(long = "soloCBwhitelist", num_args = 1.., default_values_t = vec!["None".to_string()])] pub solo_cb_whitelist: Vec, - /// 10x chemistry preset. Sets the CB/UMI geometry so it does not have to be - /// spelled out with `--soloCBstart` and friends. `-` (the default) leaves - /// the geometry to those flags. - #[arg(long = "soloChemistry", default_value = "-")] - pub solo_chemistry: String, - /// How cell barcodes are represented: `Sequence` (2-bit packed ACGT, the /// default). #[arg(long = "soloCBtype", default_value = "Sequence")] @@ -1581,46 +1575,6 @@ impl Parameters { )); } - // --soloChemistry presets. STAR applies these before validating the - // geometry, so a preset and an explicit --soloCBstart cannot disagree: - // the preset wins and says so, rather than half-applying. - if params.solo_chemistry != "-" { - let geometry = match params.solo_chemistry.as_str() { - // (cb_start, cb_len, umi_start, umi_len), 1-based as STAR takes them. - "CR_2" | "CR_3" | "CR_3.1" | "CR_4" | "SC3Pv1" => Some((1, 14, 15, 10)), - "SC3Pv2" => Some((1, 16, 17, 10)), - "SC3Pv3" | "SC3Pv4" | "SC5P" => Some((1, 16, 17, 12)), - _ => None, - }; - let Some((cb_start, cb_len, umi_start, umi_len)) = geometry else { - return Err(command.error( - ErrorKind::InvalidValue, - format!( - "unknown --soloChemistry '{}'; expected SC3Pv1, SC3Pv2, SC3Pv3, \ - SC3Pv4, SC5P, or -", - params.solo_chemistry - ), - )); - }; - for (flag, given) in [ - ("--soloCBstart", params.solo_cb_start), - ("--soloCBlen", params.solo_cb_len), - ("--soloUMIstart", params.solo_umi_start), - ("--soloUMIlen", params.solo_umi_len), - ] { - if given != 0 { - log::warn!( - "--soloChemistry {} overrides the geometry given by {flag}", - params.solo_chemistry - ); - } - } - params.solo_cb_start = cb_start; - params.solo_cb_len = cb_len; - params.solo_umi_start = umi_start; - params.solo_umi_len = umi_len; - } - // --soloCBtype: only 2-bit packed ACGT barcodes are supported. `String` // needs a different whitelist representation entirely, so refuse rather // than silently pack a non-ACGT barcode into nonsense. @@ -2720,61 +2674,6 @@ mod tests { assert!(!p.out_sam_attributes.contains(SamAttributes::XS)); } - #[test] - fn solo_chemistry_presets_set_the_geometry() { - let base = [ - "--readFilesIn", - "cdna.fq", - "bc.fq", - "--soloType", - "CB_UMI_Simple", - "--soloCBwhitelist", - "wl.txt", - "--sjdbGTFfile", - "genes.gtf", - ]; - let with = |extra: &[&str]| { - let mut a = base.to_vec(); - a.extend_from_slice(extra); - try_parse(&a) - }; - - // v2: 16-base CB, 10-base UMI. - let p = with(&["--soloChemistry", "SC3Pv2"]).unwrap(); - assert_eq!( - ( - p.solo_cb_start, - p.solo_cb_len, - p.solo_umi_start, - p.solo_umi_len - ), - (1, 16, 17, 10) - ); - - // v3 lengthened the UMI to 12. - let p = with(&["--soloChemistry", "SC3Pv3"]).unwrap(); - assert_eq!( - ( - p.solo_cb_start, - p.solo_cb_len, - p.solo_umi_start, - p.solo_umi_len - ), - (1, 16, 17, 12) - ); - - // v1 used a 14-base CB. - let p = with(&["--soloChemistry", "SC3Pv1"]).unwrap(); - assert_eq!((p.solo_cb_len, p.solo_umi_len), (14, 10)); - - // The preset wins over an explicit geometry rather than half-applying. - let p = with(&["--soloChemistry", "SC3Pv3", "--soloCBlen", "14"]).unwrap(); - assert_eq!(p.solo_cb_len, 16); - - // Unknown presets are refused. - assert!(with(&["--soloChemistry", "SC9Pv9"]).is_err()); - } - #[test] fn solo_cb_type_string_is_refused_rather_than_mispacked() { let base = [ From 9f049f832b889ec097b88fdaf3c0170c661c9301 Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Thu, 30 Jul 2026 14:38:12 +0200 Subject: [PATCH 7/9] docs(divergence): record that all four homopolymer UMIs are rejected rustar-aligner has always rejected poly-A/C/G/T UMIs in every geometry. STAR rejects only poly-A under CB_UMI_Complex, because the four packed constants are built in the SoloReadBarcode constructor from pSolo.umiL while that path does not learn its UMI length until it has extracted one (SoloReadBarcode.cpp:16-21 against getCBandUMI.cpp:353-354), leaving all four constants zero. Poly-A packs to zero and is still caught; the other three are not. CONTRIBUTING requires a deliberate divergence to be recorded with the STAR source it departs from, and this one was not. Verified against the 2.7.11b source rather than asserted. Co-Authored-By: Claude Opus 5 (1M context) --- DIVERGENCE.md | 31 +++++++++++++++++++++++++++++++ 1 file changed, 31 insertions(+) diff --git a/DIVERGENCE.md b/DIVERGENCE.md index bd957ed..424993f 100644 --- a/DIVERGENCE.md +++ b/DIVERGENCE.md @@ -30,6 +30,37 @@ This is the reason faithfulness is reported **tie-adjusted**. On the 10k yeast b **Source.** `src/rng.rs`, `src/align/read_align.rs` (`per_read_seed`, `shuffle_tied_prefix`), `src/params/mod.rs` (`MultimapperOrder`). STAR: `ReadAlign_multMapSelect.cpp`, `ReadAlignChunk` RNG seeding. +### 1.2 Homopolymer UMIs are rejected whichever base repeats + +**What STAR does.** STAR means to reject a UMI made of one repeated base, and under `CB_UMI_Simple` it does. Under `CB_UMI_Complex` it rejects only poly-A, letting poly-C, poly-G and poly-T through. + +The mechanism is an initialisation-order defect, not a rule. The four packed constants are built in the `SoloReadBarcode` constructor (`SoloReadBarcode.cpp:16-21`): + +```cpp +for (uint32 jj=0;jj<4;jj++) { + homoPolymer[jj]=0; + for (uint32 ii=0; ii Date: Thu, 30 Jul 2026 14:42:40 +0200 Subject: [PATCH 8/9] test(solo): name the cbMinP threshold and oneExact guard as their own tests Both rules were asserted inside resolve_multi_prefers_higher_prior, whose name describes only the ordering. The threshold decides between correcting a barcode and dropping the read, and oneExact is not implied by it (a sole candidate holds the whole posterior and is still refused), so each gets a test that says what it checks. Co-Authored-By: Claude Opus 5 (1M context) --- src/solo/count.rs | 41 +++++++++++++++++++++++++++++++++++++++++ 1 file changed, 41 insertions(+) diff --git a/src/solo/count.rs b/src/solo/count.rs index cea1ad1..4176747 100644 --- a/src/solo/count.rs +++ b/src/solo/count.rs @@ -2154,4 +2154,45 @@ mod tests { // is refused even when it dominates the posterior. assert_eq!(resolve_multi_cb(&cands[..1], &[0, 0], 0.0), None); } + + /// Two candidates of equal quality, so the exact-count priors decide the + /// posterior, and STAR accepts the winner only once its share reaches + /// `CB_MIN_P` (0.975, `SoloReadBarcode_getCBandUMI.cpp`). Named apart from + /// the ordering test above because this is the threshold rather than the + /// ordering: it decides between correcting a barcode and dropping the read. + #[test] + fn cb_posterior_below_min_p_is_rejected() { + use crate::solo::whitelist::CbCandidate; + let cands = vec![ + CbCandidate { + wl_index: 0, + mismatch_pos: 1, + mismatch_qual: b'I', + }, + CbCandidate { + wl_index: 1, + mismatch_pos: 2, + mismatch_qual: b'I', + }, + ]; + // 10/13 = 0.769: a clear winner on ordering, still short of the bar. + assert_eq!(resolve_multi_cb(&cands, &[10, 3], 0.0), None); + // 1000/1003 = 0.997 clears it. + assert_eq!(resolve_multi_cb(&cands, &[1000, 3], 0.0), Some(0)); + } + + /// STAR's `oneExact`: a correction needs at least one candidate that was + /// seen exactly in the whitelist counts. A sole candidate with no exact + /// support holds the whole posterior and is still refused, so this guard is + /// not implied by the threshold. + #[test] + fn a_cb_never_seen_exactly_is_refused_without_pseudocounts() { + use crate::solo::whitelist::CbCandidate; + let one = vec![CbCandidate { + wl_index: 0, + mismatch_pos: 1, + mismatch_qual: b'I', + }]; + assert_eq!(resolve_multi_cb(&one, &[0, 0], 0.0), None); + } } From eb5956abd0fd23880cdee915115bceb636e1b4d3 Mon Sep 17 00:00:00 2001 From: Benjamin Demaille Date: Thu, 30 Jul 2026 20:12:18 +0200 Subject: [PATCH 9/9] docs(changelog): record the solo barcode flags and the adapter fix Co-Authored-By: Claude Opus 5 (1M context) --- CHANGELOG.md | 14 ++++++++++++++ 1 file changed, 14 insertions(+) diff --git a/CHANGELOG.md b/CHANGELOG.md index b6e4893..e0367b4 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -13,6 +13,13 @@ Sections commonly used: Features, Bug fixes, Other changes. ### Features +- **STARsolo barcode resolution**: `--soloCBtype`, + `--soloAdapterSequence` / `--soloAdapterMismatchesNmax` for + adapter-anchored `CB_UMI_Complex` geometry, and + `--soloOutFormatFeaturesGeneField3`. Barcode correction now applies + STAR's `cbMinP` posterior threshold (0.975) and its `oneExact` + guard, so a low-confidence correction is dropped rather than accepted. + - **STARsolo single-cell quantification (`--soloType`)** — the 10x Chromium / plate-based count-matrix pipeline, ported from STAR and verified against real STARsolo (#90). @@ -99,6 +106,13 @@ Sections commonly used: Features, Bug fixes, Other changes. ### Bug fixes +- The STARsolo adapter search took the leftmost match inside the + mismatch budget; STAR's `localAlignHammingDist` takes the + best-scoring one. Since `CB_UMI_Complex` measures barcode positions + from the adapter, the wrong offset shifted the whole barcode. Ties now + go to the leftmost of the best, and an `N` in the adapter is a + wildcard rather than a mismatch. + - **STARsolo `Gene` assignment now requires exon concordance**, matching STARsolo: a read counts toward a gene only when every aligned block lies within the gene's exons, rather than merely overlapping one. This