Spatial Random Effects in Count Regression

Why Add a Spatial Random Effect?

A spot-count model can include measured predictors such as tissue region, cell-type proportion, or a treatment indicator, as well as an offset for sequencing depth. Even after these terms are included, nearby spots may have similar counts for reasons the measured covariates do not capture. A spatial random effect represents this remaining structured variation. It is not a substitute for measured covariates: it describes spatial pattern left over after conditioning on them.

The model below is for one gene across $n$ spots. The count at spot $i$ is $Y_i$, its library-size offset is $L_i$, and its measured covariate vector is $x_i$. The random effect $u_i$ is a latent, spot-specific deviation in log mean expression. The model combines a negative-binomial observation distribution with a joint distribution for the vector of spatial effects.

Definition: A spatial random-effects negative-binomial regression is a hierarchical model in which, conditional on spatial effects $u=(u_1,\ldots,u_n)$, the spot counts are independent and $$ Y_i\mid u_i,x_i \sim \operatorname{NB}(\mu_i,\theta), \qquad \log(\mu_i)=\log(L_i)+x_i^\top\beta+u_i. $$ Here $\mu_i=E(Y_i\mid u_i,x_i)$, $\theta>0$ is the negative-binomial size parameter, and $\beta$ is the vector of regression coefficients. In the mean-size parameterization used here, $\operatorname{Var}(Y_i\mid u_i,x_i)=\mu_i+\mu_i^2/\theta$. The spatial effects have a joint mean-zero distribution, commonly $u\sim N(0,\tau^2R)$, where $R$ is a specified spatial correlation matrix and $\tau^2$ controls the effect variance.

The offset $L_i$ is a known positive exposure, often a library-size measure. Its coefficient is fixed at one: doubling $L_i$, while holding the other terms fixed, doubles the conditional mean. The covariates $x_i$ include an intercept, and the random effects have mean zero under their distribution. The particular spatial correlation structure is a modeling choice; this tutorial uses $R$ generically rather than prescribing a specific structure.

Assumptions and Interpretation

The model assumes that the negative-binomial distribution is appropriate for counts conditional on the latent effects and covariates. Its conditional variance allows more variability than a Poisson distribution, with $\theta$ controlling the amount: as $\theta$ grows, the variance approaches $\mu_i$. Conditional independence means that after $u$, $x$, and the offset are specified, the remaining count fluctuations at different spots are independent. Dependence between spots is represented through the joint distribution of $u$.

The spatial effects are assumed to have mean zero and a covariance or correlation structure that reflects the intended notion of spatial similarity. The covariance matrix must be positive semidefinite. The assumed structure should be chosen in light of the spot layout and scientific question; an incorrect structure can leave residual dependence or impose an unrealistic pattern. A Gaussian distribution for $u$ is a distributional assumption about latent log-mean deviations, not a claim that the observed counts are Gaussian.

With $u_i$ held fixed, $\exp(\beta_j)$ is the multiplicative change in the conditional expected count for a one-unit increase in covariate $j$, holding the other covariates, the offset, and the spatial effect fixed. Likewise, $\exp(u_i)$ multiplies the mean relative to the mean at $u_i=0$. It is a residual spatial deviation, not a direct measure of a particular biological mechanism.

A subtle distinction is that these are conditional interpretations. If $u_i$ is Gaussian with variance $\tau^2$, and $m_i=L_i\exp(x_i^\top\beta)$, then the mean after averaging over the random effect is $$ E(Y_i\mid x_i,L_i)=E_u(m_i e^{u_i})=m_i e^{\tau^2/2}. $$ Thus $m_i$ is the mean at the reference value $u_i=0$, not the mean averaged over the random-effect distribution. If the variance of $u_i$ is the same across spots, this averaging changes the intercept-scale mean but not the multiplicative covariate ratio.

Worked Example: Combining a Covariate, Offset, and Spatial Effect

Consider a hypothetical model with $\log(\mu_i)=\log(L_i)+\beta_0+\beta_1x_i+u_i$. Let $\exp(\beta_0)=0.0002$, $\exp(\beta_1)=1.5$, and consider a spot with $L_i=20{,}000$, $x_i=1$, and $u_i=0.1$. Compute its conditional mean count.

Substituting these values gives $$ \mu_i=20{,}000(0.0002)(1.5)e^{0.1} =6(1.10517) =6.631. $$

The model therefore gives a conditional expected count of approximately 6.63 for this hypothetical spot. At the same library size and spatial effect, changing $x_i$ from 0 to 1 multiplies the conditional mean by 1.5. The positive spatial effect multiplies it by $e^{0.1}$, approximately 1.105, relative to the value at $u_i=0$. Neither multiplier is a causal effect without additional design and assumptions.

Likelihood and the Role of the Spatial Distribution

Given a value of $u$, the negative-binomial probability mass function for a count $y_i$ is $$ P(Y_i=y_i\mid u_i,x_i)= \frac{\Gamma(y_i+\theta)}{\Gamma(\theta)y_i!} \left(\frac{\theta}{\theta+\mu_i}\right)^\theta \left(\frac{\mu_i}{\theta+\mu_i}\right)^{y_i}. $$ The conditional likelihood is the product of these probabilities across spots, because conditional independence is assumed. The random effects cannot simply be treated as unrelated extra coefficients: their joint distribution supplies the spatial structure and regularizes their values toward the mean-zero distribution.

For likelihood-based inference on $\beta,\theta,\tau^2$, and any parameters in $R$, the random effects are integrated out. In schematic form, the marginal likelihood is $$ P(y\mid x,L)=\int \left\{\prod_{i=1}^{n}P(Y_i=y_i\mid u_i,x_i,L_i,\theta)\right\} f(u;\tau^2R)\,du. $$ The integral averages the conditional count likelihood over the assumed distribution of spatial effects. In a Bayesian analysis, the same hierarchy is used to form a joint posterior for effects and model parameters. In either approach, uncertainty in the latent effects and model parameters matters; plugging in one estimated spatial surface and proceeding as if it were known can understate uncertainty.

A short derivation shows why the model creates dependence in the counts. Write $m_i=L_i\exp(x_i^\top\beta)$, so that $\mu_i=m_i e^{u_i}$. For jointly Gaussian effects with covariance matrix $C=\tau^2R$, the normal moment-generating identity gives $E(e^{u_i})=e^{C_{ii}/2}$ and $E(e^{u_i+u_j})=e^{(C_{ii}+C_{jj}+2C_{ij})/2}$. By conditional independence and the law of total covariance, $$ \operatorname{Cov}(Y_i,Y_j\mid x,L) =\operatorname{Cov}(\mu_i,\mu_j) =m_i m_j e^{(C_{ii}+C_{jj})/2}\left(e^{C_{ij}}-1\right), \qquad i\ne j. $$ Therefore, positive latent covariance induces positive marginal count covariance, even though the counts are independent conditional on $u$. The negative-binomial size parameter, by contrast, governs extra variation within a spot conditional on its latent mean. These two mechanisms are related in the resulting marginal variability but are not interchangeable.

Worked Example: Translating a Spatial Effect Into a Mean Ratio

Suppose two hypothetical neighboring spots have the same offset and measured covariates, and their fitted spatial effects are $u_1=0.2$ and $u_2=-0.2$. Let the common mean at $u=0$ be 10 counts. Find their conditional means and compare them.

For spot 1, $\mu_1=10e^{0.2}=12.214$, approximately 12.2 counts. For spot 2, $\mu_2=10e^{-0.2}=8.187$, approximately 8.19 counts. Their ratio is $$ \frac{\mu_1}{\mu_2}=\frac{10e^{0.2}}{10e^{-0.2}}=e^{0.4}=1.492. $$

Under the fitted conditional model, spot 1 has an expected count about 1.49 times that of spot 2, holding their measured predictors and offsets equal. This difference is attributed to their estimated residual spatial effects. It does not establish why the spots differ, and it may be sensitive to the assumed spatial structure and measured covariates.

Worked Count Likelihood

A small calculation illustrates how the negative-binomial part evaluates a proposed conditional mean. It does not estimate a spatial surface by itself: a full spatial model also accounts for the distribution of $u$ and integrates or otherwise makes inference over the latent effects.

Worked Example: Evaluating a Conditional Negative-Binomial Likelihood

Suppose two hypothetical spots have observed counts $y_1=0$ and $y_2=2$. At $u_1=u_2=0$, take both conditional means to be 2 and the common size parameter to be $\theta=2$. Compute the conditional log likelihood.

For $y_1=0$, the combinatorial factor is 1 and $$ \log P(Y_1=0\mid \mu_1=2,\theta=2) =2\log\left(\frac{2}{2+2}\right) =-2\log 2 =-1.386294. $$ For $y_2=2$, the gamma-function factor is $\Gamma(4)/\{\Gamma(2)2!\}=3$, so $$ \log P(Y_2=2\mid \mu_2=2,\theta=2) =\log 3+2\log\left(\frac{2}{4}\right)+2\log\left(\frac{2}{4}\right) =\log 3-4\log 2 =-1.673977. $$

Conditional independence makes the total log likelihood the sum, $-1.386294-1.673977=-3.060271$, or approximately $-3.060$. This value evaluates the conditional count model at the specified means and effects. It is not, on its own, evidence that $u=0$ is the best spatial fit: assessing a proposed latent vector must also account for its probability under the spatial-effect distribution.

Diagnostics, Alternatives, and Limitations

A spatial random effect is useful only if its assumptions and fit are plausible. Examine fitted counts and residuals against library size, measured covariates, and spatial location. Residual spatial pattern after fitting suggests that the chosen spatial distribution or covariate set may not capture the dependence. Conversely, an estimated spatial effect that closely tracks a measured predictor can indicate that the two are difficult to distinguish. The effect should be interpreted as residual structure conditional on the model, not as an automatic correction for confounding.

Useful checks include comparing observed and fitted count distributions, checking whether variation beyond the fitted negative-binomial model remains, and inspecting whether residuals retain spatial autocorrelation. A handful of extreme counts, many zeros, or systematic regional differences can expose lack of fit, but no single diagnostic identifies its cause. The earlier material on Moran’s I and empirical variograms supplies tools for describing spatial pattern; such summaries can help diagnose residual structure, but they do not validate the count distribution or establish that the random-effect covariance is correct.

The main alternatives answer somewhat different modeling needs. A negative-binomial regression without a spatial effect is simpler and may suffice when conditional residuals show no meaningful spatial dependence, but it treats spots as independent after covariates. Adding spatially varying fixed predictors can be more interpretable when the relevant spatial features are measured, but cannot represent unmeasured residual dependence by itself. A spatial random effect accommodates residual dependence, at the cost of relying on a distributional and covariance assumption and potentially absorbing variation that would otherwise be attributed to covariates. A Poisson model is a more restrictive count alternative because its conditional variance equals its mean; it is not a replacement for spatial structure.

The model also has limits. A mean-zero latent effect is not identifiable as a separate overall intercept without a centering convention or a mean-zero distributional constraint. The spatial effect can be confounded with broad-scale covariates, especially when those covariates vary smoothly across tissue. Offset errors can distort all estimated rates, and a single covariance structure may not represent boundaries, disconnected tissue regions, or multiple spatial scales. Finally, the latent effect does not explain the biological process that generated the residual pattern; it accounts for that pattern statistically.

Common mistake: “The spatial random effect is positive at this spot, so the gene is more highly expressed there because of the measured cell type.” This conclusion assigns a residual term to a specific cause that the model has not identified. The correction is to say that the fitted conditional mean is higher at that spot because of its positive residual spatial effect, after accounting for the included covariates and offset. A causal or cell-type-specific interpretation needs evidence and modeling beyond this residual term.

Specify the count model.
Choose the negative-binomial mean-size parameterization, list measured predictors, and define the positive library-size offset.

Specify the latent effect.
State its mean, covariance structure, and scale, and ensure the covariance matrix is valid for the spots being modeled.

Fit and interpret conditionally.
Distinguish covariate effects at fixed spatial effect from means averaged over its distribution.

Check residual behavior.
Assess count fit, remaining overdispersion, spatial residual pattern, and potential overlap between covariates and the estimated spatial effect.

Check Your Understanding

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

If $u_i$ increases by 0.3 while all other terms remain fixed, by what factor does the conditional mean change?
Why does a spatial random effect induce marginal dependence when counts are conditionally independent?
In the model, what does the library-size offset assume about the effect of doubling library size?
How should a positive fitted spatial effect be interpreted when measured covariates and the offset are held fixed?

Solutions. (1) The mean is multiplied by $e^{0.3}=1.350$, approximately, because the effect enters the log mean additively. (2) The counts share correlated latent effects; averaging over those effects induces covariance between their conditional means, and hence between the counts. (3) It assumes the conditional expected count doubles, since the coefficient of $\log(L_i)$ is fixed at one. (4) It indicates a higher conditional expected count than at $u_i=0$, given the included terms; it does not identify a biological cause or establish causation.