What a Test for Spatial Variation Answers
A gene may vary across spots because its mean expression changes with location, because neighboring spots have correlated expression after measured covariates are accounted for, or both. A gene-wise test asks whether the observed spatial pattern is inconsistent with a specified null model. It does not, by itself, say that the pattern is large, biologically important, or caused by location.
Let \(y_{gi}\) denote the expression response for gene \(g\) at spot \(i\), for \(i=1,\ldots,n\). The response might be a normalized expression value or a transformed value suitable for a Gaussian model. For count data, the response distribution and mean-variance relationship must instead be represented by an appropriate count likelihood. A test is meaningful only relative to its model: the null hypothesis, spatial alternative, covariates, and error assumptions all matter.
For a mean-trend test, a typical null is that all coefficients for spatial predictors are zero after adjustment for the nonspatial covariates. For an autocorrelation test, the null is instead that the residual values have no spatial association under the test’s specified reference model. State which question is being tested before interpreting a p-value.
Testing Spatial Terms in Regression
Suppose the response is modeled using a nonspatial design matrix \(X\) and a matrix \(Z\) of spatial predictors, such as coordinate terms or a prespecified spatial basis. A Gaussian linear model for one gene is
Here \(\boldsymbol{\beta}_g\) represents nonspatial effects and \(\boldsymbol{\gamma}_g\) represents the spatial terms. For a joint test of \(q\) spatial coefficients, the hypotheses are \(H_0:\boldsymbol{\gamma}_g=\mathbf{0}\) versus \(H_a:\boldsymbol{\gamma}_g\ne\mathbf{0}\). This asks whether the spatial terms improve the model after accounting for \(X\).
The classical partial \(F\) test requires a full-rank design, normally distributed errors with mean zero, a common error variance, and independent errors across spots conditional on the predictors. If \(RSS_0\) and \(RSS_1\) are the residual sums of squares for the null and full models, respectively, and the full model has \(p_1\) fitted coefficients, then
Under \(H_0\) and the stated assumptions, this statistic has an \(F\) distribution with \(q\) and \(n-p_1\) degrees of freedom. The numerator measures the reduction in residual variation per added spatial coefficient; the denominator estimates residual variance in the full model. The test is an omnibus test: it can detect a nonzero combination of spatial terms without identifying a single coefficient as the whole pattern.
Worked Example: Testing a Linear Spatial Trend
Suppose a hypothetical section has 10 spots. For one gene, log-normalized expression is regressed on a centered spatial coordinate, with an intercept. Assume for this example that the Gaussian linear-model conditions hold, including independent, equal-variance errors. The estimated slope is \(0.60\) expression units per coordinate unit, with standard error \(0.20\). Test for a linear trend at \(\alpha=0.05\).
The hypotheses are \(H_0:\beta_{\text{space}}=0\) and \(H_a:\beta_{\text{space}}\ne0\). The coefficient statistic is
There are \(10-2=8\) residual degrees of freedom, because the model estimates an intercept and a slope. The two-sided p-value is \(P(|t_8|\ge 3.00)\), approximately \(0.0171\). Thus, under the model, the observed slope or a more extreme one would be unusual if the true linear slope were zero. At the 0.05 level, this is evidence against \(H_0\).
The estimated trend is \(0.60\) log-expression units per coordinate unit. That is an effect estimate, not a p-value. A 95% confidence interval is approximately \(0.60\pm 2.306(0.20)\), or \((0.139,1.061)\) log-expression units per coordinate unit. The interval describes uncertainty in the slope under the model; it does not establish that the trend is biologically consequential.
Spatial regression can also test several basis terms jointly. The basis and its complexity should be specified in a way that can represent the patterns of interest without fitting arbitrary noise. If the residuals remain spatially correlated, the ordinary partial \(F\) reference distribution may not be valid: treating correlated spots as independent can make standard errors too small and evidence appear stronger than it is.
Worked Example: An Omnibus Test of Two Spatial Terms
Consider hypothetical expression measurements at \(n=12\) spots. A null model contains an intercept and a batch indicator; a full model adds two spatial basis terms. Suppose \(RSS_0=48\) and \(RSS_1=30\). The full model has \(p_1=4\) coefficients: an intercept, a batch coefficient, and two spatial coefficients. Assume the Gaussian independent, equal-variance error conditions are reasonable.
There are \(q=2\) added terms and \(12-4=8\) residual degrees of freedom. The partial \(F\) statistic is
Under \(H_0\), compare this with \(F_{2,8}\). The upper-tail p-value is approximately \(0.153\), so the result is not significant at \(\alpha=0.05\). The full model reduces the residual sum of squares by \(18/48=0.375\), or 37.5%, but that reduction alone is not sufficient evidence against the null under this test. Nor does the nonsignificant result prove that the spatial effect is zero; the data may be too limited or noisy to estimate it precisely.
Testing Residual Spatial Autocorrelation
A spatial trend in the mean and spatial autocorrelation in residuals are different targets. Moran’s \(I\) summarizes whether values at neighboring spots tend to be similar or dissimilar. For a gene, it is often applied to residuals after fitting the nonspatial covariates relevant to the question. Applying it to unadjusted expression can detect spatial structure caused by those covariates rather than residual spatial dependence.
Let \(e_i\) be the residual for spot \(i\), and let \(W=(w_{ij})\) be a specified spatial weight matrix with no self-weights (\(w_{ii}=0\)). Write \(S_0=\sum_i\sum_j w_{ij}\), and let \(\bar e\) be the mean residual. Moran’s statistic is
With nonnegative weights, a relatively large positive value indicates that neighboring residuals tend to be alike; a negative value indicates a tendency toward local dissimilarity. The value is not an effect size on the original expression scale. Its interpretation depends on the weights and the null reference distribution. A test requires assumptions appropriate to that reference distribution, including how values are treated under the null; the statistic alone is not a p-value.
Worked Example: Calculating Moran’s I for One Gene
For illustration, take four hypothetical spots arranged in a path, with equal weights of 1 between adjacent spots and 0 otherwise. Suppose the residuals for a gene are \((2,1,-1,-2)\). They sum to zero, so \(\bar e=0\). The weight matrix is symmetric and has six directed nonzero entries, so \(S_0=6\).
The weighted cross-product counts each undirected edge twice. The three edge products are \(2(1)=2\), \(1(-1)=-1\), and \((-1)(-2)=2\). Therefore the numerator sum is \(2(2-1+2)=6\). The sum of squared residuals is \(2^2+1^2+(-1)^2+(-2)^2=10\), giving
The positive value indicates that residuals on adjacent spots tend to be alike for this arrangement. It does not establish statistical significance: four spots provide little information, and a p-value requires a suitable null distribution for this statistic and analysis. Nor does \(I=0.400\) mean expression changes by 0.400 units per spot; Moran’s \(I\) is a dimensionless measure of spatial association.
Evidence, Effect Size, and Many Genes
A p-value measures how incompatible the data are with a specified null model, conditional on its assumptions. It is not the probability that the null hypothesis is true, and it does not measure the size of the spatial pattern. With many genes, even a modest evidence threshold can produce false discoveries. A gene-wise testing workflow should therefore report the test statistic, effect estimate where available, uncertainty, raw p-value, and a multiple-testing summary.
For example, the Benjamini–Hochberg procedure controls the expected false discovery rate under its applicable dependence conditions. Given \(m\) p-values sorted as \(p_{(1)}\le\cdots\le p_{(m)}\), find the largest rank \(k\) satisfying \(p_{(k)}\le k\alpha/m\). Reject the hypotheses at ranks \(1,\ldots,k\). This is a decision rule across a family of tests; it does not convert the p-value into an effect measure.
Worked Example: Applying a False-Discovery-Rate Threshold
Suppose six hypothetical genes have sorted p-values \(0.001, 0.008, 0.021, 0.040, 0.200,\) and \(0.600\). At target false discovery rate \(\alpha=0.05\), the Benjamini–Hochberg rank thresholds are \(0.05r/6\), for ranks \(r=1,\ldots,6\): \(0.00833, 0.0167, 0.0250, 0.0333, 0.0417,\) and \(0.0500\).
At rank 1, \(0.001\le0.00833\); at rank 2, \(0.008\le0.0167\); and at rank 3, \(0.021\le0.0250\). At rank 4, \(0.040>0.0333\), and the later ranks also fail their thresholds. Thus the largest passing rank is \(k=3\), and the first three genes are rejected by this procedure. This does not mean that each rejected gene has a 5% probability of being a false discovery. The procedure controls an expected proportion over the tested family under its conditions.
Assumptions, Diagnostics, and Method Choice
Before trusting a gene-wise result, check that the model’s assumptions and the biological question align. For Gaussian regression, inspect residual-versus-fitted patterns for unequal variance or mean misspecification, residual distributions for severe departures from normality when sample sizes are small, and residual maps or residual autocorrelation for unmodeled spatial structure. Check that the spatial design terms are not redundant with covariates and that the fitted model is identifiable. For autocorrelation tests, document the graph or weights, examine whether a few high-degree spots dominate the statistic, and ensure the null reference applies to the residuals being tested.
A Gaussian regression test of coordinate terms is suited to a specified change in the conditional mean. Moran’s \(I\) is a global summary of spatial association and can detect patterns not well described by a simple linear trend, but it does not identify the shape or location of a pattern. A spatial covariance or random-effect model can represent residual dependence directly, but inference for a variance component at zero may have a nonstandard null distribution because zero lies on the boundary of the allowable parameter space. Use a reference distribution justified for the model rather than automatically applying an ordinary coefficient test.
For count data, a Gaussian model on raw counts is generally inappropriate because counts are nonnegative and their variability commonly depends on the mean. A count regression can test nested spatial terms through a likelihood-based statistic when the likelihood and regularity conditions justify its reference distribution. A spatial random-effect test needs particular care at a zero variance null. Whatever the response model, a test of spatial association is not a causal test: location may be associated with measured or unmeasured tissue features.
A second mistake is to interpret a nonsignificant result as proof of no spatial variation. A wide confidence interval, limited number of spots, inappropriate spatial scale, or misspecified model can obscure a real pattern. A test addresses only the alternatives its model can detect; a linear coordinate term, for instance, may miss a localized or nonmonotone pattern.
Decide whether the question concerns a spatial mean pattern or residual spatial autocorrelation, and state the gene-wise null and alternative.
Specify the response scale, covariates, spatial terms or weights, and assumptions that justify the test statistic’s null distribution.
Inspect residuals and report interpretable effect estimates and uncertainty where the method provides them.
Summarize the family of tests with an appropriate adjustment, and interpret evidence separately from effect size.
Check Your Understanding
Use the model assumptions and distinctions in this tutorial to answer the questions. Brief solutions follow.
- A regression includes an intercept and three spatial basis terms. State the null hypothesis for an omnibus test of spatial mean variation.
- A gene has Moran’s \(I=0.35\). What can be concluded from this value alone, and what cannot?
- Why can an ordinary independent-error \(F\) test be unreliable when residuals are spatially correlated?
- In the six-gene example, why are the first three genes selected but not the fourth?
- A gene has a very small p-value but a slope close to zero in expression units per coordinate unit. Explain how both results can be true.
Solutions. (1) The null is that all three spatial coefficients equal zero; the alternative is that at least one is nonzero. (2) The positive value indicates neighboring values tend to be alike under the chosen weights, but it supplies neither a significance conclusion nor an expression-scale effect size without an appropriate null reference and additional interpretation. (3) Correlation violates the independent-error assumption and can make the estimated standard error and reference distribution inappropriate, often overstating evidence. (4) Rank 3 is the largest rank with \(p_{(r)}\le 0.05r/6\); rank 4 fails because \(0.040>0.0333\). (5) A precisely estimated small slope can be far enough from zero relative to its standard error to yield a small p-value. Statistical evidence and effect magnitude answer different questions.