snRNA-seq Computation & QC
- What does the computational pipeline do between a folder of FASTQ files and a count matrix you can cluster?
- How do FastQC and Cell Ranger metrics each contribute to deciding whether a library is trustworthy?
- Which pool failed QC, and why — and how would we know if we hadn’t been told?
- After QC, what do normalization, batch correction, clustering, annotation, and pseudobulk DE actually do?
- Describe the snRNA-seq computational pipeline from raw reads to differential expression.
- Read CellRanger
qc_report.htmland FastQC outputs and interpret per-pool metrics across the four NS&B 11-KT pools. - Justify, with concrete metrics, why Pool 4 is excluded and Pool 3 (EO) is the workshop’s focus.
- Distinguish bulk-style sample failure (fixable upstream) from analysis decisions (handled downstream in Seurat).
- Recognize the four major analytical steps you will run hands-on: normalization (SCTransform), dimensionality reduction (PCA + Harmony), clustering (Louvain), and pseudobulk differential expression (DESeq2).
Lecture Video
Overview
This pre-class session has two halves. The first is conceptual: the video above walks you through the snRNA-seq computational pipeline, from FASTQ reads through pseudobulk DE, so you can recognize each step before you run it. The second is hands-on: a QC interpretation worksheet that you fill in by clicking through the actual Cell Ranger and FastQC reports for all four pools on the workshop server. By the time you walk into Friday’s session, you should be able to defend — in plain language — which pools to include in analysis and why.
Friday opens with a 15-minute review of your worksheet answers. Come ready to point at specific numbers in the reports, not just say “Pool 4 looked bad.”
The Pipeline at a Glance
Single-nucleus RNA-seq turns physical molecules into a tractable analytical object — a gene-by-nucleus count matrix — through a chain of computational steps. Each step has a name, a tool, and a decision point:
| Step | Tool | Decision you control |
|---|---|---|
| Alignment + UMI counting | Cell Ranger | Pre-mRNA reference; cell-calling thresholds |
| Quality control | FastQC + Cell Ranger metrics | Which samples to keep, which to exclude |
| Normalization | SCTransform | Number of variable features |
| Dimensionality reduction | PCA | Number of principal components |
| Batch correction | Harmony | Which metadata defines “batch” |
| Clustering | Louvain (Seurat) | Resolution parameter |
| Cell-type annotation | Curated markers + module scores | Which marker lists, which thresholds |
| Differential expression | DESeq2 (pseudobulk) | Design matrix, contrast direction |
Each step is short to explain but easy to fumble. The video covers the conceptual content; Friday and Saturday’s hands-on episodes walk through the code.
QC Is the Gate
The most important fact about the pipeline: Cell Ranger does not know whether the molecules in your tube were the molecules you wanted. It aligns reads to the reference and calls cells from the barcode-rank distribution. If those reads were ribosomal RNA, you get a count matrix of ribosomal RNA. There is no error flag.
This is why QC isn’t downstream cleanup — it’s epistemic. You are asking: did the experiment work at all? Two tools together answer the question:
- FastQC reads the FASTQ files before alignment. It catches base-call quality, GC content shifts, adapter contamination, and overrepresented sequences. A sample that fails FastQC has a problem in what was in the tube.
- Cell Ranger reports metrics after alignment: confident mapping rate, intronic fraction, cells called, median genes per cell. A sample that passes FastQC but fails Cell Ranger has a problem in capture or composition, not sequencing quality.
The healthy ranges for snRNA-seq:
| Metric | Healthy range | What low means |
|---|---|---|
| Q30 RNA read (R2) | > 85% | Base-call accuracy problem |
| Confident genome mapping | > 60% | Contamination, wrong reference, or degradation |
| Intronic reads | 25–45% | Reflects nuclear pre-mRNA capture |
| Median genes / cell | > 500 | Real nuclei vs. empty droplets |
Q30 measures sequencer performance — how accurately the base calls were made. You can have 93% Q30 on a library that is 90% ribosomal RNA. Those rRNA reads were sequenced accurately; they’re still rRNA. Always read Q30 together with mapping rate and intronic fraction.
QC Worksheet — Four Pools
Open the actual reports for each pool in side tabs, walk through them, and fill in the table below.
The reports live on the workshop server inside your RStudio session — they are not published to this public site. To open each one, find the file in the Files pane on the right side of RStudio (navigate into data/qc/...), click it, and choose View in Web Browser. That opens the rendered HTML in a new tab served through your authenticated RStudio session.
| Pool | Tissue | Fish | CellRanger | FastQC (flow cell A / B) |
|---|---|---|---|---|
| 1 | Hindbrain | BB48, BB50 | data/qc/pool1_hindbrain/outs/qc_report.html |
data/qc/pool1_hindbrain/fastqc_fcA/ / fastqc_fcB/ |
| 2 | Midbrain | BB48, BB49INJ, BB50, BB50INJ | data/qc/pool2_midbrain/outs/qc_report.html |
data/qc/pool2_midbrain/fastqc_fcA/ / fastqc_fcB/ |
| 3 | Electric organ (EO) | BB48, BB49INJ, BB50, BB50INJ | data/qc/pool3_eo/outs/qc_report.html |
data/qc/pool3_eo/fastqc_fcA/ / fastqc_fcB/ |
| 4 | Skin / Knollenorgan | BB48, BB49INJ, BB50, BB50INJ | data/qc/pool4_skin/outs/qc_report.html |
data/qc/pool4_skin/fastqc_fcA/ / fastqc_fcB/ |
Per-pool comparison table
Read each pool’s qc_report.html and FastQC output. Fill in the cells below. Some are filled in for you as examples — you’ll fill the blanks.
| Metric | Pool 1 (Hindbrain) | Pool 2 (Midbrain) | Pool 3 (EO) | Pool 4 (Skin) |
|---|---|---|---|---|
| Confident genome mapping (%) | 74 | 76 | _____ | _____ |
| Intronic reads (%) | 36 | 38 | _____ | _____ |
| Cells called | 4,761 | 5,234 | _____ | _____ |
| Median genes / cell | 1,500–1,700 | 1,600–2,200 | _____ | _____ |
| FastQC: per-base sequence content | _____ | _____ | _____ | _____ |
| FastQC: GC content | _____ | _____ | _____ | _____ |
| FastQC: overrepresented sequences | _____ | _____ | _____ | _____ |
Short-answer prompts
Based on the metrics you filled in:
- Which pool would you exclude from analysis, and what’s the single strongest piece of evidence?
- Is the exclusion decision a “lower-quality data” call or a “this is not gene-expression data” call? Defend your wording.
- What would you say in a Methods section of a paper to justify the exclusion?
Pool 4 is excluded. The strongest single piece of evidence is the 1% intronic fraction — snRNA-seq from real nuclei runs 25–45% intronic, and 1% means almost nothing was nuclear pre-mRNA. Combined with the 11% confident mapping and the 77,427 phantom cells called, the data is not low-quality gene expression — it is not gene-expression data at all.
A Methods sentence: “Pool 4 (skin/knollenorgan) was excluded due to rRNA contamination arising from low nuclear yield during tissue preparation; Kraken2 analysis ruled out exogenous contamination.”
The intronic-reads percentage is one of the most diagnostic Cell Ranger metrics for snRNA-seq specifically. Why does the healthy range run 25–45% for snRNA-seq vs. <5% for scRNA-seq?
snRNA-seq captures nuclear RNA, which includes large amounts of pre-mRNA — transcripts still being processed and still carrying their introns. scRNA-seq captures cytoplasmic RNA, which is fully spliced (introns removed). A high intronic fraction in snRNA-seq is the validation stamp for nuclear capture; in scRNA-seq it would be a contamination flag. Pool 4’s 1% intronic rate means whatever was captured wasn’t nuclear.
Click into Pool 2’s qc_report.html and find the barcode-rank plot. Then do the same for Pool 4. Compare them visually.
- What does the “knee” in Pool 2’s plot represent biologically?
- Pool 4 has no knee. What does that mean about what’s in the droplets?
- Why does the absence of a knee lead to 77,427 spurious cell calls?
- The knee is the gap between barcodes that contained a real nucleus (high UMI count) and empty droplets (low UMI count). Cells contain 10–100× more mRNA than empty droplets, creating the gap.
- No knee means high-UMI and low-UMI barcodes form a continuous ramp — there is no clean signal distinguishing nucleus-containing droplets from empty ones. In Pool 4, free rRNA molecules loaded into essentially every droplet, smearing the distribution.
- Cell Ranger’s barcode-rank algorithm anchors on the knee to call cells. Without one, it picks an inflection point anywhere on the ramp and calls everything above it — 77,427 barcodes with 25–73 genes each, none of them real cells.
Forensics box
When FastQC flagged Pool 4 with ~7.5% overrepresented “No Hit” sequences, the natural reaction is to assume contamination from another organism — bacterial skin microbiome being the obvious candidate.
Three independent lines of evidence ruled this out:
- BLAST. The “No Hit” consensus sequence returned 100% match to B. brachyistius 18S/28S rRNA. The unmapped reads are the fish’s own ribosomal RNA — not foreign.
- Chemistry. The 10x Genomics poly(dT) capture chemistry requires a poly(A) tail. Bacterial mRNA isn’t polyadenylated. Even if skin microbiome bacteria were present, their transcripts couldn’t be captured.
- Kraken2. A microbial classifier ran against the unmapped reads from both Pool 2 (midbrain, clean) and Pool 4 (skin, failed). In both pools, ~81–90% of unmapped reads were unclassified — identical fractions. The classified hits were reagent background (Homo sapiens handling contaminant, trace Klebsiella) and present equally in both pools. No skin-specific microbial enrichment.
The mechanism: skin nuclei are rare and fragile. Mechanical lysis protocols designed for soft brain tissue don’t efficiently release nuclei from densely-packed keratinocytes. Low yield → few real nuclei → the abundant, stable rRNA pool dominated what was captured. Result: 89% rRNA library, no knee, 77,427 phantom cells. This is a sample-preparation problem, not a bioinformatics problem.
Beyond QC: What Friday and Saturday Will Cover
Once you’ve decided which pools are real, the rest of the pipeline asks: what cell types are in this tissue, and which genes change with treatment?
Normalization removes the technical component of nucleus-to-nucleus variation. SCTransform fits a regularized negative binomial model to each gene and returns variance-stabilized residuals.
PCA compresses ~3,000 highly variable genes into 20–30 principal components — enough to capture meaningful variance, few enough to cluster on.
Harmony corrects GEM-well batch effects in the PCA embedding. With four samples on two capture days, batch is real; Harmony aligns samples within each cell type without flattening cell-type structure.
Clustering uses Louvain community detection on a kNN graph in the corrected embedding. The resolution parameter controls how finely the graph is partitioned — a judgment call, not an answer.
Cell-type annotation matches clusters to curated marker gene lists (electrocyte: scn4aa, scn4ab, kcna7a; muscle: pax7b, myomesins; Schwann: col28a1a; endothelial: cdh5; motor neuron: mnx1, chata). Module scores quantify how well each cluster matches each gene set.
Pseudobulk + DESeq2 runs the central statistical test: aggregate raw counts per (cell-type × sample), then test treatment as a fixed effect. This treats fish, not nucleus, as the unit of replication — statistically correct for snRNA-seq.
We have 2 fish per condition × ~1,400 nuclei per fish = ~5,600 nuclei. A graduate student in your lab proposes testing 11-KT vs. vehicle with a Wilcoxon rank-sum test on the ~2,800 individual 11-KT vs. ~2,800 vehicle nuclei. Predict the consequence before reading the solution.
Nuclei from the same fish are not independent — they share fish-level variation in diet, microbiome, stress, circadian state, and countless other factors unrelated to 11-KT. Treating each nucleus as a replicate inflates the degrees of freedom ~1,400-fold, deflating p-values dramatically. The result: many genes look strongly “significant” but the effect rides on one unusual fish. This is pseudoreplication.
Pseudobulk aggregates raw counts to the (cell-type × fish) level, restoring n=2 per group. DESeq2 then estimates dispersion and tests with honest degrees of freedom. The padj values will be much larger than you’d see in the pseudoreplicated test — that’s the correct statistical behavior, not a failure of the method.
Reference documents (10X Genomics)
For deeper reading on the upstream chemistry, run protocols, and yield planning that feed the Cell Ranger output you read above:
- 10x Chromium GEM-X Universal 3′ v4 — 4-plex OCM User Guide — the exact chemistry used for the four EO samples.
- Chromium Nuclei Isolation Kit user guide (CG000505 Rev B) — alternate kit-based nuclei prep (we use the open-protocol workflow described in Thursday’s episode).
- Demonstrated Protocol: Nuclei isolation for snRNA-seq (CG000124 Rev G) — 10X’s reference nuclei-isolation protocol for context.
- 10X Sequencing Handbook — sequencing planning and per-sample read budgets.
- Flow Cell Calculator (CG000604 Rev B) — Excel calculator we used to spec the four-pool flow cell loading.
Keypoints
- Cell Ranger converts FASTQs into a count matrix but does not validate whether the underlying biology was sound — that’s what QC is for.
- FastQC catches problems in what was sequenced; Cell Ranger catches problems in what was captured. Both must look right.
- The intronic fraction (25–45% healthy) is the single most diagnostic snRNA-seq metric — it distinguishes nuclear capture from contamination.
- Pool 4 was excluded because intronic 1%, mapping 11%, and 77,427 phantom cells indicate the library is host ribosomal RNA, not gene expression.
- Pool 3 (EO) is the workshop’s analytical focus — the library is real, CellRanger’s cell calls look reasonable, and the 11-KT × vehicle 2v2 design is intact.
- The downstream pipeline (SCTransform → PCA → Harmony → Louvain → curated annotation → pseudobulk DESeq2) is what Friday and Saturday’s hands-on sessions cover.
- Nuclei from the same fish are not independent replicates — always aggregate to the biological replicate level before testing for treatment effects.