Tutorials › Statistical Methods for Single-Cell RNA Sequencing › Worked Analysis of a Spatial Tissue Section

Spatial Transcriptomics · Tutorial 655 of 700

Worked Analysis of a Spatial Tissue Section

Follow a hypothetical spot-count section from coordinate checks to deconvolution and spatial testing, and learn how to report proportion uncertainty and sensitivity to analysis choices.

Advanced 12 min read

What You'll Learn

  • Check spatial coordinates, spot layout, and tissue membership before building a neighborhood graph.
  • Fit and interpret a constrained two-cell-type mixture for spot-level marker counts.
  • Report bootstrap uncertainty for deconvolved proportions and distinguish it from other sources of uncertainty.
  • Calculate Moran’s I and describe a permutation test for spatial patterning.
  • Compare results across neighborhood graphs and reference profiles without treating sensitivity as proof of correctness.

From a Tissue Section to a Defensible Analysis

A spatial analysis is a chain of linked decisions. Coordinates determine which spots count as neighbors; a reference profile helps determine each spot’s estimated cell-type composition; and the resulting proportion map is tested against a spatial null model. An error or strong assumption at an early step can change the final interpretation. The aim here is to make those dependencies visible in a compact, hypothetical analysis.

The example uses six spots arranged in two rows and three columns. Two marker genes distinguish cell types A and B. We check the coordinates, estimate cell-type proportions from the marker profiles, summarize uncertainty, and test whether the estimated proportion of A has positive spatial autocorrelation. The test is an analysis of this map, not evidence that one cell type causes another to occur nearby.

Definition: A spot-count tissue section is a set of spatially located measurements, each typically aggregating RNA from multiple cells. A deconvolved proportion is an estimated fraction of the cells or cellular contribution in a spot assigned to a specified cell type, under a chosen reference and model. A spatial neighborhood graph specifies which pairs of spots are treated as neighbors and, through its weights, how strongly those pairs contribute to a spatial statistic.

Check Coordinates Before Modeling

Suppose the six spot centers have coordinates \((x,y)\), in arbitrary grid units, and all six pass the tissue mask. The intended layout is:

Spot\(x\)\(y\)Estimated proportion of A
1110.8
2210.7
3310.2
4120.8
5220.7
6320.2

Before fitting a model, verify that coordinates are finite, unique, and expressed in a consistent scale and orientation; that count rows match coordinate records; and that the tissue mask excludes off-tissue locations. Plotting the spots is a basic but important diagnostic: it can reveal swapped coordinate columns, duplicated locations, unexpected gaps, or an incorrect image orientation. Coordinate checks do not establish that spots are independent or that a chosen graph is biologically appropriate. They establish that the spatial structure being analyzed is the intended one.

Worked Example: Defining Neighbors From the Grid

For a simple rook graph, treat two spots as neighbors if their grid positions share an edge: they differ by one unit in exactly one coordinate. Spot 1 is adjacent to spots 2 and 4; spot 2 is adjacent to spots 1, 3, and 5; spot 3 is adjacent to spots 2 and 6. The bottom row has the corresponding horizontal edges.

There are four horizontal edges, two in each row, and three vertical edges, one in each column. Thus there are seven undirected edges. If the adjacency matrix \(W=(w_{ij})\) records each undirected edge in both directions, then \(w_{ij}=w_{ji}=1\) for neighbors, zero otherwise, and \(S_0=\sum_i\sum_j w_{ij}=14\). This count is a useful check: an accidental missing or extra edge changes both the graph and the spatial statistic.

Deconvolution and Proportion Uncertainty

For a transparent worked fit, suppose the two reference profiles, on the same normalized expression scale as the spot measurements, are A: \((8,2)\) and B: \((2,8)\) for marker genes 1 and 2. Let \(p_i\) denote the proportion assigned to A in spot \(i\); the proportion assigned to B is \(1-p_i\). The linear mixture predicts marker expression \((2+6p_i,\;8-6p_i)\). Estimate \(p_i\) by minimizing squared residuals subject to \(0\leq p_i\leq1\).

Model and assumptions: For a spot with observed two-marker vector \(y_i\), the mixture model is \(y_i=p_i r_A+(1-p_i)r_B+\varepsilon_i\), where \(r_A\) and \(r_B\) are the reference profiles and \(0\leq p_i\leq1\). The fitting example uses ordinary least squares, which targets squared error and is most directly justified when errors are approximately mean-zero with comparable variance across the two markers. The profiles must be on a compatible scale, the selected markers must be informative for the cell types, and the reference should represent the cells contributing to the spot. These are modeling assumptions, not consequences of the optimization.

Worked Example: Fitting a Spot With Two Markers

Suppose spot 5 has observed marker vector \(y_5=(6.2,3.8)\). Its residual sum of squares at a candidate \(p\) is

$$ Q(p)=[6.2-(2+6p)]^2+[3.8-(8-6p)]^2=2(4.2-6p)^2. $$

The constrained minimum occurs at \(p=4.2/6=0.7\), which is within \([0,1]\). The fitted profile is \((2+6(0.7),8-6(0.7))=(6.2,3.8)\), so the residuals are \((0,0)\). The estimate is therefore \(\hat p_{5,A}=0.7\), with estimated B proportion \(0.3\). This exact fit is a consequence of the deliberately simple two-marker example; real data generally have residual error, and an exact fit with few markers is not evidence that the reference or model is correct.

A proportion estimate alone does not communicate its uncertainty. One practical approach is a bootstrap that resamples cells within each reference cell type, rebuilds the reference profiles, refits the spots under the same constraints, and takes quantiles of the resulting estimates. This captures variability induced by reference sampling under the chosen procedure. It does not automatically capture uncertainty from an incorrect cell-type list, unmodeled tissue structure, or a wrong neighborhood graph.

Worked Example: Reporting Bootstrap Uncertainty

Suppose 1,000 reference-resampling bootstrap fits for spot 5 give a median A proportion of 0.69 and a 2.5th and 97.5th percentile of 0.61 and 0.77. Report the estimate as 0.70 with a 95% bootstrap interval of 0.61 to 0.77, and state how the interval was obtained. The percentile interval describes the central range of estimates produced by this bootstrap procedure; it is not a 95% probability statement that the fixed, unknown proportion lies in this particular interval.

Inspect the bootstrap distribution as well as its endpoints. A pile-up at zero or one indicates that the constraint is frequently active, making a symmetric standard-error summary potentially misleading. A wide or multimodal distribution can indicate weak marker separation or instability in the reference profiles. If several spots are resampled jointly, preserve the intended dependence structure; resampling each spot’s measurements as though they were independent can understate uncertainty when they share reference estimates.

Testing a Spatial Pattern in the Estimated Map

A common summary of spatial autocorrelation is Moran’s \(I\). For a map \(z_i\), let \(\bar z\) be its mean, let \(W=(w_{ij})\) be a specified spatial weights matrix, and let \(S_0=\sum_i\sum_j w_{ij}\). Moran’s \(I\) compares the weighted cross-products of neighboring deviations with the total squared deviation. Positive values indicate that similar values tend to be neighbors relative to the permutation reference; the statistic alone is not a calibrated significance test.

Definition: For \(n\) spots and a weights matrix \(W\), Moran’s statistic is
$$ I=\frac{n}{S_0}\frac{\sum_i\sum_j w_{ij}(z_i-\bar z)(z_j-\bar z)}{\sum_i(z_i-\bar z)^2}. $$
When \(W\) is symmetric and stores both directions of each undirected edge, the numerator counts each undirected edge product twice. State the graph and weights convention when reporting \(I\).

Worked Example: Moran’s I for the Rook Graph

Use the A proportions \(0.8,0.7,0.2,0.8,0.7,0.2\). Their mean is \(\bar z=3.4/6=0.5667\). The centered values, repeated in each row, are \(7/30,4/30,-11/30\), and their squared sum across all six spots is

$$ \sum_i(z_i-\bar z)^2=2\frac{7^2+4^2+(-11)^2}{900}=\frac{372}{900}=0.4133. $$

For the seven undirected edges, the four horizontal-edge products sum to \(2[(7/30)(4/30)+(4/30)(-11/30)]=-32/900=-0.03556\). The three vertical-edge products sum to \((7^2+4^2+(-11)^2)/900=186/900=0.2067\). Thus the undirected edge-product sum is \(154/900=0.1711\), and the directed numerator in Moran’s formula is \(2(0.1711)=0.3422\). With \(n=6\) and \(S_0=14\),

$$ I=\frac{6}{14}\frac{0.3422}{0.4133}\approx0.355. $$

This positive statistic describes the observed arrangement relative to the chosen graph, but does not by itself give a p-value. To test against spatial randomness, hold the graph fixed and permute the observed spot values among locations. Under this null, values are exchangeable over the tested locations; if tissue regions or covariates make that implausible, an unrestricted permutation test may be miscalibrated. A permutation test also does not account for uncertainty in the deconvolved values unless that uncertainty is incorporated into the testing procedure.

Sensitivity, Diagnostics, and Alternatives

A complete analysis asks whether its conclusions are robust to defensible alternatives. For example, a queen graph also connects diagonal grid neighbors, whereas the rook graph does not. A distance-based graph instead connects spots within a stated radius, and a \(k\)-nearest-neighbor graph fixes the number of neighbors but may create uneven physical distances. Each graph defines a different spatial question. Recompute the statistic and its calibrated test under plausible choices; do not select a graph only because it produces a small p-value.

Reference sensitivity is equally important. In the hypothetical fit above, suppose an alternative reference gives A: \((9,1)\) and B: \((1,9)\). For spot 5’s observed vector \((6.2,3.8)\), the predicted profile is \((1+8p,9-8p)\), and least squares gives \(\hat p=(6.2-1)/8=0.65\). This differs from 0.70 under the first reference. Across the six spots, the corresponding estimates transform from \(0.8,0.7,0.2\) in each row to \(0.725,0.650,0.275\) in each row. The pattern remains ordered in this particular example, but both the estimated magnitudes and the uncertainty should be reported as reference-dependent.

Useful diagnostics include plotting fitted and observed marker profiles, checking residuals by marker and spot, inspecting estimates at the zero and one boundaries, and mapping bootstrap interval widths. For spatial tests, verify the graph visually, check the number and distribution of neighbors, and compare results across prespecified plausible graphs. A sudden change in conclusions after a small radius change, or strong differences between reasonable references, signals fragility rather than a reason to hide the less favorable result.

Moran’s \(I\) is not the only option. A local statistic can identify where clustering appears, but introduces multiple-testing and interpretation issues. A spatial regression can adjust for measured covariates, but requires a specified model and its own assumptions. Permutation tests are flexible for a global statistic when exchangeability is credible; they do not repair confounding or poor deconvolution. The appropriate alternative depends on whether the target is a global pattern, local regions, or an association adjusted for covariates.

Key takeaway: Treat the final spatial result as conditional on coordinate processing, the deconvolution reference, the uncertainty procedure, and the graph. Report those choices, the estimated proportions with uncertainty, the spatial statistic and its calibration, and sensitivity to reasonable alternatives.

Check Your Understanding

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

  1. Why should a coordinate plot be checked before choosing a neighborhood graph?
  2. Under the two-profile model in this tutorial, what A proportion minimizes squared error for an observed marker vector \((5,5)\)?
  3. What does a bootstrap interval from resampling reference cells describe, and what does it not necessarily capture?
  4. For the rook graph above, why is the directed numerator twice the sum over the seven undirected edges?
  5. Why can an unrestricted permutation test be questionable when tissue regions differ systematically?

Solutions. (1) A plot can reveal coordinate swaps, duplicates, gaps, or incorrect orientation that would produce the wrong neighbor relationships. (2) The model prediction is \((2+6p,8-6p)\); the two residuals are \(3-6p\) and \(-3+6p\), so \(Q(p)=2(3-6p)^2\), minimized at \(p=0.5\), within the constraint. (3) It summarizes variation induced by the specified reference-resampling and refitting procedure; it does not automatically cover reference mismatch, omitted cell types, or graph uncertainty. (4) A symmetric weights matrix represents every undirected edge as two directed entries, \((i,j)\) and \((j,i)\), each contributing the same product. (5) If values are not exchangeable across regions, unrestricted permutations break the observed region structure, so the resulting null distribution may not represent the intended null.

Clinical Biostats