Skip to contents

`ms_tally()` is the second link of the msuiter working grammar (`ms_variants() |> ms_tally() |> ...`): it counts the records of an [MsVariants] object against a 2bit reference genome and returns an [MsCatalog] holding one channel table of the requested modality. Reference contexts, DBS-pair exclusion and the per-record skip ledger are owned by the Rust core; this layer only validates arguments, picks the table and assembles the catalog. All violations raise `msuiter_error_*` conditions.

Usage

ms_tally(variants, genome, mode = "SBS96", ...)

Arguments

variants

An [MsVariants] object, as returned by [ms_variants()]. The canonical columns `chrom`, `start`, `ref` and `alt` feed the tally; an optional `sample` column attributes records to samples (columns of the count matrix appear in first-appearance order; when the column is absent, all records are attributed to the single sample id `"sample"`), and an optional `strand` column carries the transcription annotation (`T`/`U`/`B`/`N`) for the strand-aware modes.

genome

Path to an (uncompressed) UCSC 2bit reference genome whose chromosome names match the variant table exactly (no `chr` prefix normalization; unmatched chromosomes are ledgered as `skipped:unknown_chrom` rows, never an error).

mode

Channel table to build: `"SBS96"` (default), `"SBS192"`, `"SBS384"`, `"SBS1536"` or `"DBS78"`. Any other value is an error that lists every legal mode.

...

No passthrough arguments: anything arriving in `...` is rejected as an argument-spelling error (`msuiter_error_input`).

Value

An [MsCatalog] object.

Details

Strand policy (mirrors the Rust `sbs` module; there is no silent degradation): `"SBS192"` requires a `strand` column on the variant table – B/N-strand records have no SBS192 channel, so inventing the annotation would misplace real mutations and is a structured error instead; `"SBS384"` defaults a missing `strand` column to `"N"`; all other modes ignore strand.

Adjacent doublet-base substitutions are counted in DBS78 and excluded from every SBS matrix; indels, invalid alleles, reference mismatches and context-window violations are skipped into the ledger, never silently coerced (SPMG parity). The ledger is switch-independent: it is written for every input record whatever the mode.

The returned catalog's provenance merges the variants' own provenance (caller, matched normal, source, parse statistics, ...) with the tally context: `genome` (build label from the variants), `genome_path` (the 2bit file used), `mode`, `ledger` (one `record TAB destination` line per input record, in input order), `n_skipped` and `n_variants`. On a key collision the tally context wins: the tally-produced keys (`genome`, `genome_path`, `mode`, `ledger`, `n_skipped`, `n_variants`) are written after the merge and overwrite same-named keys coming from the variants' provenance.

See also

[ms_variants()] for the input container, [ms_catalog()] for the low-level constructor.

Examples

# A minimal VCF whose records sit on one short chromosome.
lines <- c(
  "##fileformat=VCFv4.2",
  "#CHROM\tPOS\tID\tREF\tALT\tQUAL\tFILTER\tINFO",
  "chr1\t3\t.\tA\tT\t.\tPASS\t.",
  "chr1\t6\t.\tA\tG\t.\t.\t."
)
path <- tempfile(fileext = ".vcf")
writeLines(lines, path)
v <- ms_variants(path, "GRCh38",
  caller = "synthetic", matched_normal = "none")

# A tiny 2bit reference is built in-memory for the example; point
# `genome` at a real UCSC .2bit file in your own analysis. The chromosome
# reads TCACGACG (2-bit packing, T=0 C=1 A=2 G=3, first base high).
twobit <- tempfile(fileext = ".2bit")
u32 <- function(x) writeBin(as.integer(x), raw(), size = 4L,
  endian = "little")
# record: dnaSize, nBlockCount=0, no blocks, reserved x2, packed dna
rec <- c(u32(8L), u32(0L), raw(0L), raw(0L), u32(0L), u32(0L),
  as.raw(c(0x19, 0xE7)))
# index: nameLen, name, data offset (16-byte header + 9-byte index entry)
idx <- c(as.raw(4L), charToRaw("chr1"), u32(25L))
con <- file(twobit, open = "wb")
writeBin(c(u32(0x1A412743), u32(0L), u32(1L), u32(0L), idx, rec), con)
close(con)

cat1 <- ms_tally(v, twobit, mode = "SBS96")
cat1
#> <MsCatalog> SBS96: 96 channels x 1 samples; skipped 0/2 variants
cat1@provenance$n_skipped # records the ledger, not an error
#> [1] 0
unlink(c(path, twobit))