Compute chromVAR deviations from an explicit sequence reference

Compute GC-corrected motif deviations and store their Z-scores in a new assay without changing source counts.

Usage

RunChromVAR(
  object,
  genome,
  genome.id,
  assay = NULL,
  new.assay = "chromvar",
  background.peaks = NULL,
  niterations = 50L,
  w = 0.1,
  bs = 50L,
  seed = 1L
)

Arguments

object: A Seurat object with aligned peak counts and Signac motifs. Motif IDs must survive Seurat assay creation unchanged (for example, underscore-containing IDs are rejected rather than silently renamed).
genome: A loaded BSgenome, DNAStringSet or indexed Rsamtools FaFile. Character package names are not accepted and no reference is downloaded.
genome.id: Explicit non-empty reference label, identical to the source assay and any known label on the loaded reference. Unlabeled reference provenance is the caller's declaration, not an assembly inference. Known chromosome lengths must agree and source peaks must fit the reference.
assay: Source chromatin assay; defaults to the active assay.
new.assay: Unused name for the chromVAR Z-score assay.
background.peaks: Optional integer matrix with one row per source peak, in exactly the same order and with matching row names. Values are 1-based source-peak indices. With NULL, use GC and accessibility matching. For explicit backgrounds, the caller is responsible for appropriate matching; this function validates the current peak-index universe only.
niterations, w, bs: Parameters passed unchanged to getBackgroundPeaks.
seed: Seed for background sampling; caller RNG is restored.

Value

The Seurat object with a new assay. Its data layer contains Z-scores; Misc(new.assay, "scop_chromvar") contains deviations, GC bias, actual background peaks, ordered source identifiers and parameters. This receipt is a calculation-time snapshot: later subsetting or RenameCells operations do not synchronize its nested experiments or identifiers. Source metadata and default assay are unchanged; no count QC is created from scores.

Details

This calls chromVAR directly and retains the deviations object and actual background indices in assay Misc, while the assay data layer contains Z-scores. No peaks or cells are silently dropped, and non-finite scores are not replaced with zero. Every source peak and cell must have positive total counts. At least two source-assay cells are required. Z-scores quantify motif-associated accessibility deviations, not TF expression, occupancy, or causality. Explicit stored sparse zeros are removed locally before motif indexing. The underlying computeDeviations uses the caller's BiocParallel configuration; this function does not register a global parallel backend.

Examples

(function() {
  old <- BiocParallel::bpparam()
  on.exit(BiocParallel::register(old, default = TRUE))
  BiocParallel::register(BiocParallel::SerialParam(), default = TRUE)
  n <- 200L
  counts <- withr::with_seed(321, Matrix::Matrix(matrix(
    1 + stats::rpois(n * 4, 4), n, 4,
    dimnames = list(paste0("p", seq_len(n)), paste0("cell", 1:4))
  ), sparse = TRUE))
  ranges <- GenomicRanges::makeGRangesFromDataFrame(data.frame(
    seqnames = "chr1", start = (seq_len(n) - 1L) * 100L + 1L, end = seq_len(n) * 100L
  ))
  atac <- ReadATACData(counts, ranges = ranges, genome = "toy", verbose = FALSE)
  gc <- 10L + seq_len(n) %% 81L
  reference <- Biostrings::DNAStringSet(c(chr1 = paste0(vapply(
    gc,
    function(k) paste0(strrep("G", k), strrep("A", 100L - k)), ""
  ), collapse = "")))
  # Toy annotations illustrate the interface, not a biological TF claim.
  incidence <- Matrix::Matrix(cbind(Grich = gc > 50, Arich = gc < 30), sparse = TRUE)
  rownames(incidence) <- rownames(atac)
  pwm <- stats::setNames(rep(list(matrix(1, 4, 3,
    dimnames = list(c("A", "C", "G", "T"), NULL)
  )), 2), colnames(incidence))
  Signac::Motifs(atac[["peaks"]]) <- Signac::CreateMotifObject(data = incidence, pwm = pwm)
  scored <- RunChromVAR(atac, reference, genome.id = "toy", niterations = 10, w = 0.2, bs = 10)
  head(SeuratObject::LayerData(scored[["chromvar"]], layer = "data"))
})()
#>            cell1      cell2     cell3       cell4
#> Grich -0.5791645 -0.3256361 0.4755358 -0.04483244
#> Arich  0.0420846 -0.3166115 0.3583383 -0.30367846