Raw Data Processing in Single-Cell Sequencing

1. Key Steps in Raw Data Processing

2. Starting Point

• Raw input: lane-demultiplexed FASTQ files (derived from base calling and quality estimation).
• Issues like base-calling errors and index hopping can affect data quality but can be mitigated computationally or experimentally (e.g., dual indexing).

3. Raw Data Quality Control (QC)

• Use tools like FastQC to assess sequencing read quality (quality scores, base content, N content).
• QC is essential even if downstream tools include built-in checks.
• For single-cell data, transcript-derived reads (typically read 2 in 10x Chromium data) are biologically meaningful for QC metrics.
• Technical reads (barcodes, UMIs) may not show typical biological patterns but metrics like N base fraction remain relevant.
• Multiple FastQC reports can be aggregated with MultiQC.

4. Alignment and Mapping

• Purpose: determine the origin loci of each sequenced fragment (genomic or transcriptomic).
• Raw sequences include:
  • Cell Barcodes (CB): unique cell identifiers.
  • Unique Molecular Identifiers (UMIs): distinguish molecules to correct amplification bias.
  • Raw cDNA sequences: actual read sequences.
• Accurate mapping is critical; errors lead to inaccurate count matrices.
• scRNA-seq datasets are large (hundreds of millions to billions of reads), making mapping computationally demanding.
• Many RNA-seq aligners are protocol-agnostic and do not handle scRNA-seq specifics (CBs, UMIs).
• Specialized tools for scRNA-seq mapping and processing include:

• These tools differ in internal methods, file handling (BAM vs. in-memory), and computational trade-offs.

5. Types of Mapping Algorithms

Mapping approaches can be classified by:

1. Type of mapping performed
2. Reference sequence used

a) Main Mapping Types

b) Alignment-based Methods

• Use heuristics to find candidate loci.
• Score nucleotide-level alignments using dynamic programming.
• Variants include:
  • Global alignment: entire query and reference aligned.
  • Local alignment: aligns subsequences.
  • Semi-global (fitting) alignment: most of the query aligns to a substring of the reference.
• Soft-clipping allows partial alignment by reducing penalties at read ends.
• These variations modify dynamic programming rules but keep overall complexity similar.

Accurate raw data processing, especially mapping and UMI resolution, is essential to generate reliable count matrices for downstream single-cell RNA-seq analyses.

1. Raw Data Quality Control

a) Alignment Heuristics and Optimizations

• Banded alignment: Limits dynamic programming computation to a diagonal band, ignoring low-score alignments [Chao et al., 1992].
• X-drop and Z-drop heuristics: Early pruning of unpromising alignments to save computation [Zhang et al., 2000; Li, 2018].
• Wavefront alignment: Achieves optimal alignments with reduced time and memory, especially for high-scoring matches [Marco-Sola et al., 2020].
• Vectorization and parallelism: Use of instruction-level parallelism and difference encoding to speed up dynamic programming [Farrar, 2007; Suzuki and Kasahara, 2018].

b) CIGAR String

• Encodes alignment operations as a concise string in SAM/BAM files.
• Example:
  • 3M2D4M = 3 matches/mismatches, 2 deletions, 4 matches/mismatches.
  • Extended form 3=2D2=2X distinguishes matches (=) and mismatches (X).
• Used to represent the backtrace of the alignment.

c) Alignment-Based vs Lightweight-Mapping Methods

• Spliced alignment: Handles reads spanning splice junctions.
• Contiguous alignment: Aligns reads continuously without splicing.

d) Reference Sequence Choices for Mapping

• Full genome mapping requires splice-aware aligners (e.g., STAR, STARsolo).
• Transcriptome mapping uses contiguous alignment, faster but misses non-transcript reads.
• Augmented transcriptome includes additional sequences to capture unspliced or intronic reads, improving assignment and reducing artifacts.

Mapping choice impacts sensitivity, computational cost, and biological interpretation, especially for single-nucleus data and RNA velocity analyses.

e) Summary of Key Points

• Heuristics like banded alignment and X-drop optimize alignment speed without sacrificing accuracy.
• CIGAR strings encode detailed alignment operations essential for downstream analysis.
• Alignment-based methods provide quality scores; lightweight methods trade accuracy for speed.
• Reference choice (genome, transcriptome, augmented transcriptome) affects mapping completeness and computational resources.
• Augmented transcriptomes are crucial for capturing reads outside annotated transcripts, especially in single-nucleus data.

1. Barcoded Beads and Cell Barcodes (CBs)

• Barcoded beads tag RNA from individual cells with unique oligonucleotide cell barcodes (CBs).
• These CBs enable in silico demultiplexing of sequencing reads into individual cellular bins.
• Beads contain high-diversity DNA barcodes for parallel barcoding of molecular content.

2. Alignment Orientation Filtering

• Some chemistries require reads to align in a specific orientation to the reference transcript (e.g., 10x Genomics 3' Chromium expects reads on the forward strand).
• Reads mapped in the reverse-complement orientation may be filtered out based on user settings.
• This "stranded" protocol must be documented when applicable.

3. Errors in Barcoding

Key sources of discrepancy between observed CBs and actual cells:

4. Strategies for Cell Barcode Identification and Correction

1. Correction against a Known Whitelist

   • CBs are matched against a known permit list of barcodes.
   • Barcodes not in the list are corrected by finding the closest match (using Hamming or edit distance).
   • Ambiguities arise when corrupted barcodes are close to multiple whitelist entries.
   • Tie-breaking can use base quality scores or restrict correction to barcodes observed above a frequency threshold.
   • Reads with ambiguous corrections are usually discarded.

2. Knee or Elbow-Based Methods

   • Construct a cumulative frequency plot of barcodes ranked by read/UMI counts.
   • Identify an inflection point ("knee" or "elbow") separating high-quality barcodes from likely erroneous ones.
   • Barcodes above the knee form a data-driven permit list for correction.
   • Flexible for chemistries with or without external whitelist.
   • May be overly conservative or fail if no clear knee exists.

3. Filtering Based on Expected Cell Count

   • User provides an estimate of expected cell number.
   • Barcodes are sorted by frequency; those near the expected count frequency are considered valid.
   • Remaining barcodes are corrected against this valid set.
   • Useful when knee detection is unreliable.

4. Filtering Based on Forced Number of Valid Cells

   • User manually sets the number of valid barcodes.
   • Barcodes above this threshold are accepted; others corrected accordingly.
   • Simple but assumes the threshold accurately reflects real cell count.

5. UMI Resolution and Deduplication

• After CB correction, reads are assigned to corrected CBs or discarded.
• To quantify gene abundance per cell, UMI deduplication is performed to remove PCR amplification bias.
• The goal is to identify reads and UMIs derived from each original pre-PCR molecule.
• The deduplicated UMI counts per gene per cell serve as the raw expression estimates for downstream analysis.

> UMI deduplication is essential to accurately estimate molecule counts by removing PCR duplicates within each cell barcode.

1. Cell Barcode Correction

UMI (Unique Molecular Identifier) resolution is crucial in single-cell RNA sequencing (scRNA-seq) to accurately count original molecules by correcting errors and ambiguities in UMI tags.

2. Key Concepts in UMI Resolution

• Reads and UMIs:

  • Reads mapped to a gene locus are called reads of that reference.
  • UMIs tagging these reads are UMIs of that reference.
  • Reads tagged by the same UMI form the reads of that UMI.
  • Each read has one UMI but can map to multiple genes (multimapping).

• Ideal UMI scenario:

  • Each UMI uniquely tags one pre-PCR molecule.
  • Reads of a UMI map to a single gene.
  • Molecule count = number of distinct UMIs per gene.

3. Practical Challenges Necessitating UMI Correction

4. Evidence and Impact

• UMI errors are common and correlate with PCR amplification level.
• Multimapping occurs at non-trivial rates, especially in gene families with similar sequences.
• Discarding multimapping reads biases gene abundance estimates.

5. Approaches to UMI Resolution

• Discarding multi-gene equivalence classes (ECs):

  • Count molecules by tallying ECs uniquely assigned to genes.
  • Simple but loses information from multimapping reads.

• Statistical inference (e.g., EM algorithm):

  • Probabilistically assign multi-gene UMIs to genes.
  • Uses gene-unique ECs as evidence to apportion ambiguous UMIs.
  • Maximizes likelihood of observed data given gene counts.
  • Applied independently per cell to build the count matrix.

6. Quantification Summary

• The count vector per cell is generated after UMI resolution.
• Statistical models improve accuracy but are limited by sparse per-cell data.
• UMI resolution is essential for reliable gene abundance estimation in scRNA-seq.

To retain: UMI errors and multimapping complicate molecule counting; statistical methods like EM improve UMI assignment and gene quantification, crucial for accurate single-cell data analysis.

1. UMI Resolution

UMI (Unique Molecular Identifier) resolution is a critical step in single-cell RNA-seq data processing to accurately count unique RNA molecules and avoid inflation due to PCR duplicates.

a) Key Concepts

• UMI deduplication aims to collapse reads sharing the same cell barcode (CB), gene, and UMI sequence into a single count, representing one original molecule.
• Errors in UMI sequences caused by sequencing or PCR can lead to overestimation of molecule counts if not corrected.
• UMI resolution methods use error correction algorithms to group similar UMIs likely originating from the same molecule.

b) Common UMI Resolution Strategies

c) Important Considerations

• UMI length and complexity affect the probability of collisions (different molecules sharing the same UMI).
• Sequencing depth and error rates influence the choice of UMI resolution method.
• Over-aggressive collapsing can underestimate molecule counts; under-collapsing inflates counts.
• Some tools integrate UMI resolution with cell barcode correction and gene assignment for joint optimization.

To retain accurate molecule counts, UMI resolution must balance error correction with avoiding over-collapsing distinct molecules.

d) Summary Workflow in UMI Resolution

1. Group reads by cell barcode and gene.
2. Construct UMI similarity graph based on sequence differences (usually 1 mismatch).
3. Apply error correction rules (e.g., directional adjacency) to merge UMIs.
4. Count unique UMIs post-correction as proxy for original molecules.

This step is essential before generating the final gene count matrix used for downstream analyses such as clustering, differential expression, and trajectory inference.

1. Overview of Preprocessing Tools for Single-Cell RNA-seq

• Common tools for raw data processing in single-cell and single-nucleus RNA-seq are generally robust and accurate.
• The choice of preprocessing tool often has less impact on downstream analyses (e.g., clustering) than other analysis steps.
• Lightweight mapping restricted to spliced transcriptomes can increase spurious mapping and false gene expression.
• Tool choice depends on:
  • Task requirements
  • Computational resources

• alevin-fry: Memory-frugal, small index size, suitable for single-nucleus RNA-seq including unspliced transcripts.
• Among alignment-based pipelines, STARsolo is preferred over Cell Ranger for speed, memory efficiency, and similar results.

2. Real-World Example: Processing with alevin-fry

a) Data and Reference Preparation

• Input data: FASTQ files from single-cell sequencing.
• Reference: Genome sequences (FASTA) + gene annotations (GTF).
• Example uses chromosome 5 of human genome (GRCh38, GENCODE v32/Ensembl 98).
• Reads subset from a human brain tumor dataset (10x Genomics).

b) Software Setup

• Use conda to create environment and install:
  • simpleaf (wrapper for alevin-fry pipeline)
  • alevin-fry
  • salmon
  • pyroe
• For Apple silicon (M1/M2), use Rosetta2 translation layer to ensure compatibility.

c) Commands Summary

3. Simplified Raw Data Processing Pipeline with simpleaf

• simpleaf encapsulates the alevin-fry pipeline into two main commands:

• Output structure after indexing:
  • ref/ folder: splici reference files
  • index/ folder: salmon index built on splici reference

4. Environment Variable and Configuration

• Set environment variable for configuration storage:

mkdir alevin_fry_home
export ALEVINF_RY_HOME='alevin_fry_home'

• Use simpleaf set-paths to detect and configure paths to required tools.

Key point:
The choice of preprocessing tool should balance accuracy, speed, memory usage, and specific analysis goals; lightweight mapping methods like alevin-fry are efficient for standard analyses, while alignment-based methods are preferred for detailed genomic mapping.

1. Count Data Representation

Count matrix files generated by alevin-fry:

• The number of rows in quants_mat_cols.txt = 3 × number of annotated genes (S, U, A variants).
• Each row in quants_mat_rows.txt corresponds to a cell barcode after correction/filtering.

2. Loading count data into analysis tools

• Use pyroe.load_fry() (Python) or fishpond::loadFry() (R) to load the count matrix as an AnnData object.
• Default: the X layer sums spliced (S) and ambiguous (A) counts per gene.
• Alternative: load separate layers for unspliced (U), spliced (S), and ambiguous (A) counts for more detailed analysis.

Including intronic (unspliced) counts can increase sensitivity and improve downstream analyses.

3. Typical pipeline steps for count data generation

• Output quantification folder contains count matrices and metadata.
• simpleaf wraps these steps into a simplified command-line interface.

4. Key points on count data representation

• Count matrix is sparse: only non-zero entries stored in quants_mat.mtx.
• Gene identifiers include splice status suffixes: enables separate quantification of spliced, unspliced, and ambiguous reads.
• Cell barcodes are corrected and filtered to remove errors and low-quality cells.
• Quantification outputs can be flexibly loaded to include or exclude intronic counts depending on analysis goals.

The count matrix structure supports advanced single-cell RNA-seq analyses that leverage splicing information for improved biological insight.

1. Key Tools for Single-Cell RNA-Seq Raw Data Processing

Several tools are widely used for processing raw single-cell RNA-seq data, each with specific strengths:

2. Alevin-fry Pipeline Workflow

Alevin-fry is a modular pipeline for scRNA-seq data processing with key steps:

1. Cell Barcode Correction
   Command example:

   alevin-fry generate-permit-list -u CB_permit_list -d expected_orientation
   

   • Corrects cell barcodes by generating a permit list.
   • Handles strand specificity (e.g., -d fw for forward strand).

2. Filter Mapping Information
   Command example:

   alevin-fry collate -i gpl_out_dir -r alevin_map_dir -t num_threads
   

   • Collates mapping data for downstream processing.

3. UMI Resolution and Quantification
   Command example:

   alevin-fry quant -r resolution -m txp_to_gene_mapping -i gpl_out_dir
   

   • Resolves UMIs and quantifies gene expression.
   • Uses transcript-to-gene mapping files (e.g., 3col.tsv in splici_ref folder).

3. Workflow Wrappers and Resources

• Quantaf: Nextflow workflow for large-scale processing of single-cell and single-nucleus datasets using alevin-fry. Provides preprocessed public datasets.
• Simpleaf: Wrapper for alevin-fry enabling the entire pipeline execution (reference generation to quantification) with minimal commands.
• Pyroe (Python) and Roe (R): Helper packages for processing alevin-fry quantification outputs.
• Galaxy Project Tutorials: Accessible tutorials for scRNA-seq raw data processing and quality control.

4. Quality Control and Doublet Detection Tools

• dropkick: Automated QC and cell identification for droplet-based scRNA-seq.
• miQC: Probabilistic framework for adaptive quality control.
• scds and Solo: Computational annotation and identification of doublets in scRNA-seq data.

5. Important References for Methods and Tools

• Alignment and quantification methods: STAR, Minimap2, kallisto|bustools, alevin-fry.
• Doublet detection: Solo, scds, DoubletDecon.
• Quality control: dropkick, miQC.
• Data processing pipelines: nf-core scrnaseq, Quantaf, Simpleaf.

To remember: Efficient and accurate scRNA-seq raw data processing relies on modular pipelines like alevin-fry combined with barcode correction, mapping filtering, and UMI resolution, supported by wrappers and QC tools for streamlined workflows.

1. Practical Implementation with alevin-fry

alevin-fry is a modular, efficient, and memory-constant tool designed for preprocessing droplet-based single-cell RNA-seq data. It focuses on accurate gene abundance estimation by integrating alignment and quantification steps.

a) Key Features of alevin-fry

• Modularity: Separates indexing, mapping, and quantification into distinct steps, allowing flexibility and customization.
• Efficiency: Uses lightweight mapping algorithms (e.g., selective alignment) to speed up processing without sacrificing accuracy.
• Constant Memory Usage: Designed to handle large datasets with fixed memory requirements, enabling scalability.

b) Workflow Overview

1. Indexing: Build a transcriptome index from reference sequences to enable fast mapping.
2. Mapping: Map reads to the transcriptome using selective alignment or pseudoalignment, capturing potential multi-mapping reads.
3. Quantification: Resolve ambiguous mappings and estimate gene-level abundances, accounting for unique molecular identifiers (UMIs) and cell barcodes.
4. Post-processing: Filter out empty droplets and low-quality cells, correct errors in barcodes and UMIs, and generate count matrices for downstream analysis.

c) Advantages Over Other Tools

d) Practical Tips

• Use selective alignment mode in alevin-fry for improved accuracy in mapping reads with splicing or sequencing errors.
• Perform barcode and UMI correction to reduce noise and improve quantification reliability.
• Apply empty droplet filtering to distinguish true cells from background noise.
• Generate gene count matrices compatible with downstream tools like Seurat or Scanpy.

To retain: alevin-fry enables fast, memory-efficient, and accurate preprocessing of droplet-based single-cell RNA-seq data by modularizing indexing, mapping, and quantification, with integrated error correction and filtering steps.