Current best practices in single‐cell RNA‐seq analysis: a tutorial

Key results

1. SuggestedPre-processing: scran pooled size factors for non-full-length data followed by log(x+1) without gene scaling; ComBat when batches share cell compositions; 1,000-5,000 highly variable genes (HVGs); UMAP to visualize, PCA (or diffusion maps for trajectories) to summarize.

   Scran beat other tested normalizations for batch correction (Buttner et al, 2019) and differential expression (Vieth et al, 2019); ComBat handled low-to-medium complexity experiments and beat MNN in simple batch settings (Buttner et al, 2019); Klein et al (2015) saw similar PCA embeddings with 200 to 2,400 HVGs; UMAP favored for speed and scalability (Becht et al, 2018). · Normalization; Batch effects and data integration; Feature selection; Visualization

2. SuggestedDownstream: Louvain modularity clustering on a K-nearest-neighbour graph; automated plus manual marker-based annotation when reference atlases exist; Slingshot for simple and PAGA for complex trajectories; MAST or limma for differential expression on uncorrected data with covariates in a full-rank design.

   Louvain beat other scRNA-seq clustering methods (Duò et al, 2018; Freytag et al, 2018); Saelens et al (2018) ranked Slingshot best for linear to multifurcating paths and recommended PAGA for complex ones; MAST was the best single-cell DE tool in Soneson & Robinson (2018), runs 10-fold to 100-fold faster than weighted bulk tools, and limma-voom adds a further 10-fold speedup. · Clustering; Cluster annotation; Trajectory inference; Differential expression testing

3. SuggestedFive processing stages fall into measured, corrected and reduced layers; measured data suit statistical testing, corrected data visual comparison, and reduced data visualization, clustering and trajectory inference, since correction distorts gene variances that tests rely on.

   Fig EV2: ComBat batch correction and DCA denoising tend to lower gene coefficient of variation in mouse intestinal epithelium (Haber et al, 2017) and mouse embryonic stem cell (Klein et al, 2015) data, so over-corrected genes are more likely called DE; Soneson & Robinson (2018) benchmarked DE methods on raw and normalized inputs only. · Stages of pre‐processed data; Figure legends; Tables

4. SuggestedMarker-gene P-values are inflated because clusters and markers come from the same data, so P-values should not validate clusters, though gene rankings remain usable; the same circularity affects differential testing along pseudotime, and a permutation test accounting for clustering is proposed.

   Clustering random data simulated with splatter (Zappia et al, 2017) still yielded significant marker genes (Appendix); inflated P-values overstate marker numbers, yet, if clusters are meaningful, top-ranked genes stay the best candidates; DE over experimental conditions within clusters is independent of clustering and keeps expected significance. · Cluster annotation; Gene expression dynamics; Differential expression testing

5. ShownApplied to mouse intestinal epithelium data, the workflow illustrates jointly visualized quality control thresholds, ComBat reducing batch separation, Louvain clusters annotated with literature markers, and Slingshot distal and proximal enterocyte trajectories, released as an adaptable R/Python notebook.

   Count-depth cutoff at 1,500 (peak near 1,200 counts) and gene cutoff at 700 (noise peak at approx. 400 genes), with high mitochondrial fractions confined to low-count cells; PCA elbow between PCs 5 and 7 of the first 31; markers Slc12a2 (stem), Arg2 (enterocyte), Tff3 (goblet), Defa24 (Paneth); ComBat made sample separation in UMAP less visible. · Introduction; Quality control; Figure legends