The necessity of data-driven optimization of the number of principal components in scRNA-seq analysis - examples on real and simulated data

Technical Summary: The Necessity of Data-Driven Optimization of the Number of Principal Components in scRNA-seq Analysis

Problem Statement
Single-cell RNA-sequencing (scRNA-seq) generates high-dimensional data, necessitating dimensionality reduction techniques like Principal Component Analysis (PCA) for downstream tasks such as clustering, visualization, and trajectory analysis. Canonical pipelines in widely used software packages (e.g., Seurat v5, Scanpy v1.11) typically default to using a fixed number of principal components (PCs), often 50, or suggest a broad range of 10–50. However, the paper argues that this choice lacks rigorous justification. Selecting too few components risks losing local biological diversity, while selecting too many may introduce technical noise and increase computational costs. Furthermore, the optimal number of components is likely dataset- and task-specific, yet current practices rely on static defaults that do not account for evolving data characteristics, protocols, or dataset complexity. While data-driven methods exist to estimate the number of components, they have not been systematically evaluated specifically for scRNA-seq analysis.

Methodology
The authors conducted a comprehensive benchmarking study involving 25 real scRNA-seq datasets (from Tabula Sapiens, Tabula Muris, and other sources) and 37 simulated datasets generated by four different simulators (ScDesign2, Sparsim, Splatter, Symsim). The study evaluated the impact of PCA embedding size on downstream clustering performance and assessed eight data-driven methods for selecting the number of components:

1. Fixed Default: PC50 (50 components).
2. Threshold-based: Scree plot (SCREE), Eigengaps (EG and EG-bound), Explained Variance (EV), and Adjusted Correlation Thresholding (ACT).
3. Hypothesis-based: Kritchman-Nadler (KN).
4. Simulation-based: Jackstraw (JS) and Parallel Analysis (PA).

The evaluation framework focused on clustering quality using the Leiden algorithm, measured by the Adjusted Rand Index (ARI) against ground truth labels (for simulated data) or known cell types (for real data). To ensure fair comparison, the resolution parameter of the Leiden algorithm was optimized to align with the true number of clusters, acknowledging that clustering quality is sensitive to both the number of PCs and the resolution. Additional metrics included neighborhood preservation (trustworthiness score), convergence success rates, and computational time (wall-clock time). Experiments were conducted under varying pre-processing conditions, including standardized vs. non-standardized data and different numbers of highly variable genes (2,000 vs. 5,000).

Key Results

• Impact of Embedding Size: For most datasets, clustering performance exhibits a broad optimal range between 10 and 100 components, where the default of 50 often falls. However, specific datasets (e.g., Koh, Kumar, and certain simulated groups) showed narrow optimal ranges or peaks that did not include 50 components, demonstrating that fixed defaults can lead to suboptimal results for complex or specific data types.
• Performance of Data-Driven Methods:
  • Clustering Accuracy: No single method consistently outperformed all others across all datasets. However, ACT, SCREE, Jackstraw (JS), Parallel Analysis (PA), and the default PC50 generally achieved the best clustering performance.
  • Explained Variance (EV): This method tended to select too many components, introducing noise and resulting in poorer clustering performance on most simulated datasets, though it performed adequately on some real datasets (e.g., Skin, Prostate).
  • Eigengaps (EG): This method generally selected fewer components, which negatively impacted clustering on complex real datasets but performed well on simulated data.
  • Kritchman-Nadler (KN): While theoretically sound, this method suffered from poor convergence rates on simulated datasets, limiting its practical utility.
  • Neighborhood Preservation: EG consistently achieved the highest trustworthiness scores, preserving local neighborhood structures better than other methods, particularly on simulated data.
• Computational Cost: Simulation-based methods (JS, PA) were computationally slower due to repeated eigenvalue calculations but scaled similarly to other methods. The KN method showed potential for "exploding" in computation time due to solving non-linear systems.

Significance and Claims
The paper concludes that while the general recommendation of 10–50 components (including the default 50) performs well for many datasets, it is insufficient as a universal rule given the increasing size, complexity, and diversity of scRNA-seq data. The authors argue that data-driven selection of the number of components is a necessity to avoid suboptimal embeddings that could impair downstream biological conclusions.

The study provides practical recommendations, categorizing methods into tiers based on performance, stability, and convergence. The authors recommend ACT, SCREE, JS, and PA as the most robust approaches. They caution against relying solely on dimensionality-restricted methods (like fixed defaults or EG-bound) which may fail to adapt to shifting data characteristics. Furthermore, the paper highlights the critical interdependence between the number of PCA components and the resolution parameter in clustering algorithms, suggesting that these hyperparameters should be optimized jointly rather than independently.

Ultimately, the work serves as a benchmark for the field, moving beyond the "one-size-fits-all" default approach and advocating for adaptive, data-driven strategies in scRNA-seq analysis pipelines.