Two Different Questions About Reproducibility
A single-cell analysis can be repeatable in one sense and unstable in another. The same workflow may reproduce its output files exactly when rerun, yet yield different cell embeddings after a new stochastic fit. Conversely, the coordinates in two embeddings may differ while their cell neighborhoods and biological interpretation remain similar. A useful assessment therefore begins by asking what output is being compared, and what kind of variation is relevant to the claim?
Terminology varies across fields. Here, computational repeatability means obtaining the same computational result by rerunning a fully specified analysis on the same inputs. Statistical reproducibility means that the target scientific pattern or conclusion remains sufficiently stable when relevant sources of statistical variation are allowed to change, such as stochastic fitting or the sampled donors. This does not mean that every estimate or cell coordinate must be identical.
These are not universal labels with one standard threshold. A report should define the rerun conditions, the output being checked, and what counts as agreement. For example, a saved count matrix may be required to match exactly, while floating-point summaries may be allowed to differ by a small numerical tolerance. A claim about a cell population may be judged by neighborhood or cluster agreement rather than equality of embedding coordinates.
What Should Match, and What May Vary?
The expected degree of agreement depends on the object. If the input count matrix and its cell and gene identifiers are unchanged, basic input checks such as dimensions and identifiers should match exactly. A deterministic transformation with unchanged code and settings should also produce the same values, apart from documented numerical precision differences. By contrast, stochastic model fitting, sampling-based procedures, or changes in the set of donors can produce legitimate variation.
Before comparing runs, state the target quantity. For an embedding, possible targets include the pairwise distances among cells, local neighborhoods, or downstream group structure. For a donor-level analysis, the target may be an estimated contrast across donors. These are different targets and need not behave alike. In particular, raw embedding axes have no general biological meaning: rotations or reflections can change coordinates without changing the underlying geometry.
Worked Example: Checking an Exact Computational Output
Suppose a hypothetical workflow converts a fixed count matrix into a table of cell and gene identifiers, then computes library sizes by summing counts within each cell. The input contains three cells with total counts 120, 95, and 105. A second run uses the same input and analysis specification and reports totals of 120, 95, and 105 in the same cell order.
For cell \(i\), the library size is \(L_i=\sum_{g=1}^{G} x_{ig}\), where \(x_{ig}\) is the count for gene \(g\) in cell \(i\). The first run gives \(L=(120,95,105)\), and the second gives \(L'=(120,95,105)\). The elementwise differences are \(L-L'=(0,0,0)\), so the summaries agree exactly.
This establishes agreement for these reported totals, not for every part of the workflow. A thorough check would also compare the input dimensions, cell and gene identifiers, and the exact count values. If cell order is permitted to change, compare values after matching by cell identifier rather than comparing rows by position. The conclusion is limited: these checks support computational repeatability of the specified outputs under the tested conditions; they do not show that a biological conclusion is stable to a different sample of donors.
Why Embedding Coordinates Need Special Treatment
An embedding maps cells into a lower-dimensional space to represent some structure in the data. The representation may be useful even when its coordinate axes are not uniquely determined. For example, rotating every point preserves all pairwise distances. Thus, direct subtraction of coordinates from two runs can suggest disagreement even when the geometry is unchanged.
A practical comparison should match the metric to the intended use. Pairwise-distance agreement evaluates global geometry; neighborhood overlap evaluates local structure. If coordinate-by-coordinate comparison is needed, first align the configurations using a suitable transformation, such as a rotation and translation. Alignment should not be used to conceal meaningful changes: if scale, distances, or neighborhoods are part of the scientific claim, those must be assessed directly as well.
Worked Example: Comparing Two Cell Embeddings
Consider two hypothetical two-cell embeddings. Run A places cell 1 at \((1,0)\) and cell 2 at \((0,1)\). Run B places cell 1 at \((0,1)\) and cell 2 at \((-1,0)\). A raw coordinate comparison shows differences, but the second configuration is a rotation of the first.
In run A, the distance between the cells is \(d_A=\sqrt{(1-0)^2+(0-1)^2}=\sqrt{2}\approx1.414\). In run B, it is \(d_B=\sqrt{(0-(-1))^2+(1-0)^2}=\sqrt{2}\approx1.414\). The pairwise distance is exactly preserved.
To see the coordinate change, rotate run A by 90 degrees counterclockwise using \((x,y)\mapsto(-y,x)\). Cell 1 moves from \((1,0)\) to \((0,1)\), and cell 2 moves from \((0,1)\) to \((-1,0)\). The aligned coordinates are precisely those of run B. Therefore, coordinate equality before alignment would be the wrong criterion for this comparison. In a real analysis, two cells are not enough to establish useful embedding stability; the same reasoning would be applied to all cells, alongside neighborhood or distance summaries relevant to the scientific interpretation.
Separate Donor Variation From Run Variation
In single-cell studies, cells from the same donor are not independent biological replicates of that donor. If the scientific claim concerns a condition contrast across people or animals, the donor-level summaries are the units that represent between-donor replication. Repeating a fit on the same donors can reveal algorithmic instability, but it does not create additional donors.
One way to organize the comparison is to index a summary by donor \(d\) and run \(r\). Let \(Y_{dr}\) be a donor-level estimate from run \(r\). Differences across donors within a run describe variation among the sampled donors; differences across runs for the same donor describe sensitivity to the computational fit. These sources should not be pooled as if they were interchangeable.
The identity follows by writing \(Y-E[Y]=(Y-E[Y\mid D])+(E[Y\mid D]-E[Y])\). When the squared expression is averaged, the cross-product has expectation zero because \(E[Y-E[Y\mid D]\mid D]=0\). The two remaining expected squared terms are \(E[\operatorname{Var}(Y\mid D)]\) and \(\operatorname{Var}(E[Y\mid D])\), respectively. The result does not say that these sources can always be estimated precisely from a small number of runs or donors; it clarifies what they represent.
Worked Example: Repeated Donor-Level Contrasts
In a hypothetical study, four donors each have a paired condition contrast in a gene summary, measured in log-normalized expression units. Three stochastic fits produce the following donor-level contrasts:
| Donor | Run 1 | Run 2 | Run 3 |
|---|---|---|---|
| 1 | 0.8 | 0.7 | 0.9 |
| 2 | 0.4 | 0.5 | 0.3 |
| 3 | 0.6 | 0.5 | 0.7 |
| 4 | 0.2 | 0.3 | 0.1 |
The mean contrast in each run is \( (0.8+0.4+0.6+0.2)/4=0.50\), \( (0.7+0.5+0.5+0.3)/4=0.50\), and \( (0.9+0.3+0.7+0.1)/4=0.50\) log-normalized expression units. Thus the overall mean is stable across these three fits. But the donor-level values are not identical: for donor 1 they range from 0.7 to 0.9, and for donor 2 from 0.3 to 0.5. A stable mean alone would miss this variation.
For donor 1, the across-run mean is \(0.8\), and the sample standard deviation is \(s=\sqrt{((0.8-0.8)^2+(0.7-0.8)^2+(0.9-0.8)^2)/(3-1)}=0.10\) log-normalized expression units. The same calculation gives \(0.10\) for each donor. This indicates modest, nonzero run-to-run variation in individual donor estimates for this example.
For comparison, the sample standard deviation across donors in Run 1 is \(s=\sqrt{((0.8-0.5)^2+(0.4-0.5)^2+(0.6-0.5)^2+(0.2-0.5)^2)/(4-1)}\approx0.2582\) log-normalized expression units. The corresponding standard error of the mean across the four donors is \(SE=s/\sqrt{4}\approx0.1291\) log-normalized expression units. This standard error describes uncertainty in the donor-average contrast under the usual sampling interpretation; it is not the standard deviation across computational runs. The three runs do not turn four donors into twelve independent biological replicates.
A Practical Assessment and Its Limits
A reproducibility assessment is most informative when it follows a defined target rather than demanding that every file be identical. First identify the fixed inputs and analysis specification. Then identify which outputs are expected to match exactly, which may differ only within numerical tolerance, and which are expected to vary because of stochastic fitting or sampling. Finally, select a comparison that reflects the scientific claim.
Name the object being assessed: for example, input counts, an embedding's cell neighborhoods, or a donor-level contrast.
State which inputs, analysis settings, and computational conditions are held fixed and which are deliberately varied.
Use exact equality for identifiers or count values when appropriate, numerical tolerances for floating-point summaries, and geometry- or neighborhood-based comparisons for embeddings.
Summarize differences across runs and across donors separately, and state what level of disagreement would change the scientific interpretation.
Useful diagnostics include checking whether cell identifiers are matched correctly, whether a small number of cells account for most embedding changes, and whether donor-level contrasts change sign or magnitude across fits. A high average similarity can conceal unstable subgroups; a stable overall effect can conceal donor-specific changes. Inspecting summaries and relevant individual cases is therefore preferable to reporting one aggregate score without context.
This comparison also has limits. Agreement across repeated fits on the same data does not show that the result will generalize to new donors. Stability of an embedding does not establish that its axes have biological meaning or that an apparent group difference is causal. A reproducibility threshold is purpose-dependent, and there is no single metric that certifies every part of an analysis. The useful conclusion is narrower: under the stated rerun conditions, specified outputs did or did not meet specified agreement criteria.
Common Pitfall
A common mistake is to treat every numerical difference as evidence that an analysis failed, or to treat an identical overall mean as proof that the analysis is stable. The first ignores transformations such as rotation that preserve embedding structure. The second can hide donor-specific changes that cancel in the average. The correction is to define the scientific target first, compare outputs using a suitable measure, and report run-level and donor-level variation separately. Do not interpret computational repeats as independent biological observations.
Check Your Understanding
Answer each question, then compare your reasoning with the solutions.
- Two embeddings have different coordinates but identical pairwise distances. What can be concluded, and what cannot be concluded from this fact alone?
- A workflow is run five times on the same six donors. How many independent donor-level biological replicates are present?
- Across four runs, a donor-level contrast is 0.2, 0.4, 0.6, and 0.8 units. Find its mean and sample standard deviation across runs.
- Why might an unchanged mean contrast across repeated fits fail to demonstrate stability for every donor?
Solutions. (1) The embedding configurations have the same pairwise geometry, so coordinate differences alone do not establish structural disagreement. This does not establish stable neighborhoods in a larger data set, biological validity, or agreement of other workflow outputs. (2) There are six independent donor-level replicates, not thirty; repeated fits do not add donors. (3) The mean is \((0.2+0.4+0.6+0.8)/4=0.50\) units. The sample standard deviation is \(\sqrt{((0.2-0.5)^2+(0.4-0.5)^2+(0.6-0.5)^2+(0.8-0.5)^2)/(4-1)}=\sqrt{0.20/3}\approx0.2582\) units. (4) Donor-specific changes can offset one another in the average, so the mean may stay constant while individual donor estimates move. Compare donor-level results as well as the aggregate.