Tutorials › High-Dimensional Genomic Data Analysis › Choosing Factor Number With Parallel Analysis

Dimension Reduction and Clustering · Tutorial 606 of 795

Choosing Factor Number With Parallel Analysis

You will learn to construct and interpret a dimension-matched parallel-analysis benchmark, apply it to PCA and factor analysis, and recognize when its noise reference is inappropriate.

Intermediate 12 min read

What You'll Learn

  • Define the null reference used in parallel analysis and why it must match the data dimensions and preprocessing.
  • Compare ordered observed eigenvalues with simulated mean or upper-quantile benchmarks.
  • Apply the leading-component stopping rule in concrete PCA and factor-analysis examples.
  • Recognize how centering, scaling, dependence, and batch effects affect the comparison.
  • Distinguish parallel analysis from the Kaiser rule and scree-plot judgment.
  • Diagnose unstable or misleading retention decisions.

Why Compare Eigenvalues With Noise?

A large first principal component can arise even when variables have no meaningful shared structure: finite samples produce uneven eigenvalues by chance. This is especially important in genomic profiles, where the number of measured features can be large relative to the number of samples. Parallel analysis asks whether each observed eigenvalue is larger than would typically arise from noise in a data matrix with the same dimensions and preprocessing.

The method is a retention rule, not a guarantee that retained components represent biological mechanisms. It compares the observed spectrum with a specified noise reference. The conclusion therefore depends on how that reference is generated, and on whether its assumptions are reasonable for the data.

Definition: Parallel analysis compares the ordered eigenvalues from an observed data matrix with corresponding eigenvalues from many simulated noise matrices having the same number of rows and columns and the same relevant preprocessing. A component is retained when its observed eigenvalue exceeds its chosen noise benchmark. Common benchmarks are the simulated mean eigenvalue or a specified upper quantile, often the 95th percentile.

Set Up the Comparison

Let \(X\) be a data matrix with \(n\) samples in rows and \(p\) features in columns. After choosing whether to center and scale the features, calculate the eigenvalues of the resulting covariance or correlation matrix. Write them in descending order as \(\lambda_1\geq\lambda_2\geq\cdots\). For each simulated noise matrix, apply the same preprocessing and calculate its ordered eigenvalues.

The noise generator defines the null reference. One option is to generate independent random values for every sample-feature entry, often from a normal distribution. Another is to permute each feature independently across samples. Independent column-wise permutation preserves each feature’s observed values while breaking its alignment with the other features. It is appropriate only when the row structure permits such permutation; for example, restrictions may be needed for repeated measures or other dependent observations.

If there are \(B\) simulated matrices, let \(\lambda^{(b)}_j\) be eigenvalue \(j\) from simulation \(b\). At each rank \(j\), the simulations provide a reference distribution. A mean-based rule uses the average of those \(B\) values; an upper-quantile rule uses a chosen quantile of them. The choice should be stated because the rules need not select the same number of components.

$$ \bar{\lambda}^{\,\mathrm{noise}}_j = \frac{1}{B}\sum_{b=1}^{B}\lambda^{(b)}_j, \qquad q_{j,\,0.95} = \text{the 95th percentile of } \{\lambda^{(1)}_j,\ldots,\lambda^{(B)}_j\}. $$

A common retention convention compares eigenvalues in order and retains the leading components while each observed eigenvalue exceeds its benchmark. Stop at the first rank that does not exceed the benchmark. This consecutive rule avoids interpreting a later isolated exceedance as evidence for a meaningful extra dimension after an earlier rank has failed. State whether the benchmark is the mean or an upper quantile; “parallel analysis” alone does not specify that choice.

Centering and scaling must be consistent between observed and simulated data. For a correlation-matrix analysis, standardize each nonconstant feature as appropriate in both the observed and simulated matrices. For a covariance-matrix analysis, use the same units and variance treatment in both. A reference generated for a different \(n\), \(p\), or preprocessing scheme is not the same noise benchmark.

Why the Dimensions and Preprocessing Matter

Even independent noise can have unequal sample eigenvalues. The observed spectrum should therefore be compared with the spectrum generated by noise under the same sampling geometry, rather than with a universal cutoff. For centered data, the rows sum to zero, so the centered matrix has rank at most \(n-1\); it also has rank at most \(p\). Consequently, no more than \(\min(n-1,p)\) eigenvalues can be nonzero. This fact matters when \(p\) is much larger than \(n\).

For standardized features, the sample correlation matrix has trace \(p\): its diagonal entries are all 1, and the trace is their sum. Since the trace is also the sum of its eigenvalues, their average across all \(p\) ranks is 1. That does not imply that every eigenvalue should be compared with 1. Sampling variation can make the first few noise eigenvalues larger than 1, particularly in high dimensions. Parallel analysis estimates those rank-specific fluctuations by simulation.

$$ \sum_{j=1}^{p}\lambda_j = \operatorname{tr}(R) = p, \qquad \operatorname{rank}(X_{\mathrm{centered}}) \leq \min(n-1,p). $$
Justification: Under independent-feature noise, the population correlation matrix has no off-diagonal associations. The sample correlation matrix is nevertheless random, so its eigenvalues vary across samples. Simulating matrices with the same \(n\), \(p\), and preprocessing approximates the rank-specific eigenvalue variation expected under that noise model. Comparing \(\lambda_j\) with the simulated reference for the same rank therefore accounts for finite-sample variation that a fixed cutoff of 1 ignores.

Worked PCA Examples

Worked Example: Applying a 95th-Percentile Benchmark

Suppose a hypothetical study has \(n=40\) samples and \(p=5\) standardized features. Its observed correlation-matrix eigenvalues are \(2.35, 1.42, 0.78, 0.30,\) and \(0.15\), which sum to \(5.00\), as required for a correlation matrix. After \(B=2{,}000\) independent-noise simulations, the rank-specific 95th-percentile benchmarks are \(1.55, 1.30, 1.12, 0.98,\) and \(0.85\).

At rank 1, \(2.35>1.55\), so the first observed eigenvalue exceeds its benchmark. At rank 2, \(1.42>1.30\), so the second also exceeds its benchmark. At rank 3, \(0.78<1.12\), so the comparison stops there. The parallel-analysis rule retains two components.

The comparison is appropriate if the simulations used the same dimensions, centering, and standardization as the observed analysis. The result says that the first two eigenvalues exceed the chosen 95th-percentile independent-noise benchmarks; it does not establish that the components are biological or that the 95th percentile is a formal family-wise error guarantee.

Worked Example: A High-Dimensional, Small-Sample Matrix

Consider a hypothetical centered profile matrix with \(n=25\) samples and \(p=100\) features. Its centered rank cannot exceed \(\min(25-1,100)=24\), so at most 24 eigenvalues can be nonzero. Suppose a parallel analysis using matrices of exactly these dimensions and the same preprocessing gives observed leading eigenvalues \(7.20, 4.10, 1.25, 1.10\), with corresponding 95th-percentile noise benchmarks \(5.80, 3.70, 1.30, 1.18\).

The first two comparisons pass: \(7.20>5.80\) and \(4.10>3.70\). The third does not pass, since \(1.25<1.30\). The fourth also does not pass, since \(1.10<1.18\). Under the consecutive rule, retain two components.

The large difference between \(p\) and \(n\) makes dimension matching especially important: a simulation with a different sample size would produce a different noise spectrum. The result is still conditional on the independent-noise reference. If sample dependence or batch effects are present, the simulated null may fail to represent the relevant background structure.

Using Parallel Analysis for Common Factors

PCA describes total variance through eigenvalues of a covariance or correlation matrix. Common-factor analysis instead seeks dimensions that account for shared covariance, separating it conceptually from feature-specific variance. For factor retention, parallel analysis can compare observed eigenvalues of a reduced correlation matrix—one in which the diagonal entries represent estimated communalities—with corresponding eigenvalues from simulated noise data processed using the same factor-analysis procedure. The precise construction varies across implementations, so report which eigenvalues and which simulation procedure were used.

The noise reference should match the target question. Comparing ordinary PCA eigenvalues with a noise reference can help choose principal components, but does not by itself answer how many common factors to retain. For a common-factor analysis, use a factor-oriented parallel-analysis procedure and apply the same communalities or communality-estimation procedure to the observed and simulated data. The rotation used to present a retained factor solution does not determine the parallel-analysis benchmark.

Worked Example: Comparing Common-Factor Eigenvalues

Suppose a hypothetical factor analysis uses a specified reduced-correlation-matrix procedure. The first three observed common-factor eigenvalues are \(1.62, 1.18,\) and \(0.91\). Simulations using the same \(n\), \(p\), preprocessing, and communality-estimation procedure give mean rank-specific benchmarks \(1.35, 1.10,\) and \(0.83\).

Using the mean-benchmark rule, rank 1 passes because \(1.62>1.35\). Rank 2 also passes because \(1.18>1.10\). Rank 3 passes because \(0.91>0.83\). Thus, the leading-three rule retains three common factors, provided the next rank is checked and is the first one that does not exceed its benchmark. These values illustrate the comparison; a complete analysis must calculate enough ranks to identify that stopping point.

This decision is not interchangeable with a PCA decision: the eigenvalues being compared arise from a common-factor construction. It is also conditional on the chosen mean benchmark and on the specified factor-analysis procedure.

Diagnostics, Alternatives, and Limitations

Before interpreting a retention result, check that observed and simulated eigenvalues are computed from the same matrix type, in descending order, after the same centering and scaling. Confirm that all simulations use the intended dimensions, and that the number of nonzero eigenvalues is consistent with the rank bound. If the number of simulations is small, the estimated upper quantiles can be noisy; increasing \(B\) can reduce Monte Carlo variability, though it cannot correct a poorly chosen null model.

The most consequential diagnostic is whether the noise generator represents the structure that should count as background. Independent-column permutations preserve marginal feature values but destroy cross-feature association. They can also destroy row-level dependence. Independent normal simulation additionally assumes a particular marginal distribution. If batches, related samples, or other known dependencies are relevant, a naive independent-noise reference may label those structures as signal. Design a null simulation that preserves appropriate nuisance structure, or adjust for it before analysis when justified. An eigenvalue above a simple independent-noise benchmark is not proof of biological meaning.

Parallel analysis is commonly considered alongside two alternatives. The Kaiser rule retains correlation-matrix components with eigenvalues greater than 1; it is simple, but does not account for rank-specific finite-sample noise. A scree plot displays the observed eigenvalues and asks where their decline visibly changes; it is useful as a descriptive view, but locating an elbow can be subjective. These approaches may disagree because they use different criteria. Report the selected rule and its result rather than presenting the number of retained dimensions as uniquely determined.

ApproachDecision basisMain caution
Parallel analysisObserved rank-specific eigenvalues exceed simulated noise benchmarksDepends on the null generator and benchmark choice
Kaiser ruleCorrelation-matrix eigenvalue exceeds 1Ignores the dimension-specific noise spectrum
Scree plotVisual change in the eigenvalue declineThe apparent elbow can be subjective

Common mistakes include comparing a correlation-matrix spectrum to simulations based on unstandardized data, changing sample or feature dimensions in the null simulations, and treating an exceedance as a test that the component is biologically real. Another error is to declare that rank \(j\) passes when its observed eigenvalue is below the benchmark. Check each inequality explicitly and state the benchmark being used.

Key takeaway: Parallel analysis asks whether each leading eigenvalue is unusually large relative to a specified, dimension-matched noise reference. Its result is useful only when the observed and simulated analyses are comparable and the null model reflects the background structure relevant to the question.

Check Your Understanding

Answer each question, then compare your reasoning with the solutions.

  1. With \(n=25\) centered samples and \(p=100\) features, what is the maximum possible rank of the centered data matrix?
  2. An observed rank-3 eigenvalue is \(1.25\), and its benchmark is \(1.30\). Does it pass?
  3. Why should a noise simulation use the same sample and feature dimensions as the observed data?
  4. Observed eigenvalues are \(1.7, 1.2, 0.1\); benchmarks are \(1.4, 1.3, 0.8\). How many leading components pass the consecutive rule?
  5. What does an eigenvalue passing parallel analysis establish, and what does it not establish?

Solutions. (1) The rank is at most \(\min(n-1,p)=\min(24,100)=24\). (2) No: \(1.25<1.30\), so rank 3 does not exceed its benchmark. (3) Eigenvalue distributions depend on matrix dimensions and preprocessing; a mismatched simulation supplies a different noise reference. (4) Only rank 1 passes: \(1.7>1.4\), while rank 2 fails because \(1.2<1.3\). Rank 3 also fails because \(0.1<0.8\). (5) It establishes that the eigenvalue exceeds the specified independent-noise benchmark. It does not show that the component is biological; batch structure, for example, can also produce an eigenvalue exceedance. A batch-aware null or appropriate adjustment may be needed for that separate question.

Clinical Biostats