scMultiOme

Streamlined, HPC-ready workflow for joint scATAC + scRNA (10x Multiome) analysis using Seurat + Signac.

Given a set of Cell Ranger ARC outputs, the pipeline performs per-sample QC, independent RNA and ATAC dimensionality reduction, weighted-nearest-neighbor (WNN) multimodal integration, motif annotation, and automated cell-type annotation via SingleR reference mapping — all dispatched as a SLURM array so every sample runs in parallel.

---

Workflow

flowchart TD
    A[scMultiOme-meta.s] -->|sbatch --array=1-N| B[run-scMultiOme.s]
    B -->|per sample i| C[scATAC-WNN.R]
    C --> D[QC + filtering]
    D --> E[RNA: SCTransform -> PCA -> UMAP]
    D --> F[ATAC: TF-IDF -> SVD -> UMAP]
    E --> G[WNN integration -> clustering]
    F --> G
    G --> H[Motif annotation - JASPAR2020]
    G --> I[getCellDex-annos.R: SingleR labels]
    H --> J[(Outputs: QC / DimReduction / SingleR-Annotations)]
    I --> J

Loading

Each line of your input manifest becomes one element of a SLURM array, so a 12-sample run launches 12 independent jobs.

---

Repository layout

Path	Description
scMultiOme-meta.s	Top-level launcher. Reads the config and submits the array job. This is the only script you run by hand.
run-scMultiOme.s	Per-task SLURM worker. Resolves the sample/output paths for its array index and calls the R pipeline.
scATAC-WNN.R	Core analysis: QC, dim. reduction, WNN integration, motif annotation.
getCellDex-annos.R	SingleR annotation helper (getCellDex()), sourced by scATAC-WNN.R.
configFile.txt	Template config: species, input/output manifests, HPC resources.
blacklists/	ENCODE blacklist regions (mm10 BED + cached GRanges; hg38 v2 BED).

---

Requirements

• An HPC cluster with SLURM and either conda/mamba (recommended) or a site R module.
• All R/Bioconductor dependencies are bundled in environment.yml — see Quick start below. This removes any need for a site-specific R module.

Full package list (installed automatically by environment.yml)

• R packages: Seurat, Signac, hdf5r, dplyr, tidyverse, ggplot2, glue, patchwork
• Bioconductor: EnsDb.Mmusculus.v79, EnsDb.Hsapiens.v86, BSgenome.Mmusculus.UCSC.mm10, BSgenome.Hsapiens.UCSC.hg38, celldex, SingleR, SingleCellExperiment, JASPAR2020, TFBSTools, motifmatchr

Both mouse (mm10) and human (hg38) are supported; the species is set in the config and controls which genome, annotations, blacklist, and SingleR references are used.

---

Inputs

You provide two plain-text manifests, matched line-by-line (line i of one corresponds to line i of the other):

• inputs.txt — one Cell Ranger ARC output directory per sample, one per line.
• outs.txt — the destination directory for that sample's results, one per line.

# inputs.txt                          # outs.txt
/path/cellranger/sampleA/outs         /path/results/sampleA
/path/cellranger/sampleB/outs         /path/results/sampleB

Each Cell Ranger directory listed in inputs.txt must contain:

• filtered_feature_bc_matrix.h5 (RNA + ATAC counts)
• per_barcode_metrics.csv
• atac_fragments.tsv.gz (with its .tbi index)

---

Configuration

Populate configFile.txt before running:

Variable	Meaning
species	mm10 or hg38 (aliases like mouse/human also accepted).
conda_env	Conda env to activate for R (default scmultiome). Blank = use module load instead.
inputs	Path to the input manifest (default inputs.txt).
outs	Path to the output manifest (default outs.txt).
batch_size	Auto-derived from the line count of inputs — sets the array size.
module1_time	Walltime per sample (SLURM --time, e.g. 0-3 = 3 h).
module1_mem	Memory per sample (e.g. 100GB).
module1_Path	Worker script (default run-scMultiOme.s).
qc_*	Per-cell QC thresholds (see below). Edit to retune filtering.

---

Quick start

1. Create the environment (once)

The bundled environment.yml installs R and every dependency, so you don't need a site module load r/...:

conda config --set channel_priority strict     # recommended for bioconda envs
mamba env create -f environment.yml             # or: conda env create -f environment.yml
conda activate scmultiome

mamba is strongly recommended here — solving this many Bioconductor packages with classic conda is slow. The env is named scmultiome; the worker activates it automatically via the conda_env variable in configFile.txt. (To keep using an HPC R module instead, leave conda_env blank — the worker then falls back to module load.)

2. Configure and launch

From the repository directory (so blacklists/ and the helper script resolve correctly):

# edit configFile.txt (species, conda_env, QC thresholds, resources)
# create inputs.txt + outs.txt (one matched line per sample)
sbatch scMultiOme-meta.s

That's it — scMultiOme-meta.s submits one array job covering every sample.

---

Outputs

For each sample, results land under its outs.txt directory in three subfolders:

QC/

• QC-prefilter.pdf, ATAC-QC-panel.pdf, QC-panel-adj-*.pdf, postQC-panel-*.pdf
• TSS-plot.pdf, nucleosome-histogram-byGroup.pdf
• qc-summary.csv (cells pre/post-QC and percent lost)
• postQC-preNorm.rds

DimReduction/

• RNA-UMAP.pdf, ATAC-UMAP.pdf, WNN-UMAP.pdf
• scRNA-ATAC.rds (object with WNN graph, clusters, and motif info)

SingleR-Annotations/

• Mouse: ImmGen broad + fine labels, marker tables, marker heatmaps, annotated UMAPs, broad-fine-annos-obj.rds.
• Human: the same, computed against both HumanPrimaryCellAtlas and DatabaseImmuneCellExpression references (each in its own subdirectory).

QC thresholds

Cells are retained when they satisfy all of the following. These are set in configFile.txt (the values below are the shipped defaults); edit any line there to retune. If a line is removed or blank, scATAC-WNN.R falls back to the same default.

configFile.txt variable	Metric	Default
qc_min_atac_frags / qc_max_atac_frags	atac_peak_region_fragments	3,000 – 20,000
qc_min_pct_reads_in_peaks	pct_reads_in_peaks	> 15%
qc_max_blacklist_ratio	blacklist_ratio	< 0.1
qc_max_nucleosome_signal	nucleosome_signal	< 4
qc_min_tss_enrichment	TSS.enrichment	> 2
qc_min_nCount_RNA / qc_max_nCount_RNA	nCount_RNA	1,000 – 25,000
qc_max_percent_mt	percent.mt	< 20%

---

Notes

• All dependencies install from environment.yml (conda-forge + bioconda). Per the Signac maintainers, when using conda you should install R packages through conda rather than install.packages() inside R.
• The SingleR reference datasets (ImmGen, HumanPrimaryCellAtlas, DICE) are fetched by celldex from ExperimentHub on first use, so compute nodes need network access (or a pre-populated ExperimentHub cache) for the annotation step.
• Blacklist regions are from the Boyle Lab ENCODE blacklist. The mm10 GRanges object is cached in blacklists/ to avoid re-parsing.
• The RNA/ATAC integration follows the Seurat WNN multiome vignette and the Signac PBMC vignette.

References

• Hao et al. Integrated analysis of multimodal single-cell data. Cell (2021).
• Stuart et al. Single-cell chromatin state analysis with Signac. Nat. Methods (2021).
• Aran et al. Reference-based analysis of lung single-cell sequencing (SingleR). Nat. Immunol. (2019).