Banner

3.1 Foundations of Spatial Regression¶

Part: 3 — Spatial Regression (Chapter 3.1 of 7)

Dataset: diabetes_atlantic.shp (666 counties, Atlantic seaboard)

Focus: Why OLS fails on spatial data, substantive vs. nuisance dependence, the model taxonomy, and data preparation

Library: spatialstats.regression, spatialstats.weights (PySpatialStats)

Table of Contents¶

  1. Overview

  2. Python Setup

  3. Why OLS Fails on Spatial Data

  4. Two Mechanisms: Substantive vs. Nuisance Dependence

  5. A Taxonomy of Spatial Regression Models

  6. Data Preparation

  7. Summary and Conclusions

  8. Resources

1. Overview¶

Part 2 established that Dia_pct is spatially autocorrelated across counties. This chapter asks the natural next question: what does that mean for a regression model that explains diabetes prevalence from covariates like obesity and physical inactivity? The short answer — plain OLS is no longer safe, in one of two distinct ways — sets up everything that follows in Part 3.

Step What we do
Diagnosis Show concretely why spatially dependent residuals break OLS's standard-error and, in one case, coefficient guarantees
Two mechanisms Distinguish substantive dependence (real spillovers between counties) from nuisance dependence (an unmodelled spatially structured omitted variable)
Taxonomy Place OLS, SLX, SAR (SLM), SEM, SARAR, SDM, and SDEM on one map, and note which this library implements
Data Prepare the regression dataset used throughout Part 3: Dia_pct on six health/access covariates plus an urban/rural indicator, over 666 Atlantic-seaboard counties

By the end, Chapter 3.2 fits the OLS baseline on this exact dataset and starts accumulating the diagnostic evidence that Chapters 3.3–3.5 use to choose a model.

2. Python Setup¶

In [1]:
import sys, warnings
from pathlib import Path

import numpy as np
import pandas as pd
import geopandas as gpd
import matplotlib.pyplot as plt

warnings.filterwarnings("ignore", category=FutureWarning)
warnings.filterwarnings("ignore", category=DeprecationWarning)


def find_repo_root(start: Path = Path.cwd()) -> Path:
    for p in [start, *start.parents]:
        if (p / "pyproject.toml").exists() and (p / "data").is_dir():
            return p
    raise FileNotFoundError("Run this notebook from inside the PySpatialStats repository.")


ROOT = find_repo_root()
DATA = ROOT / "data"
try:
    import spatialstats
except ImportError:
    sys.path.insert(0, str(ROOT))
    import spatialstats

from spatialstats.weights import Queen
from spatialstats.esda import Moran
from spatialstats.regression import OLS

plt.rcParams.update({
    "figure.dpi": 100, "savefig.dpi": 100, "font.size": 10, "axes.titlesize": 11,
    "axes.titleweight": "bold", "axes.spines.top": False, "axes.spines.right": False,
    "axes.grid": True, "grid.alpha": 0.25, "legend.frameon": False,
})
C = {"blue": "#2b6cb0", "red": "#c0392b", "green": "#2f855a", "orange": "#dd6b20",
     "grey": "#718096", "purple": "#6b46c1", "teal": "#2c7a7b"}

print(f"spatialstats : version {spatialstats.__version__}")
print(f"data folder  : {DATA.relative_to(ROOT)}/")

spatialstats : version 0.2.0
data folder  : data/

3. Why OLS Fails on Spatial Data¶

OLS assumes each observation's error term is drawn independently. Part 2's whole apparatus — Moran's I, LISA — exists because that assumption routinely fails for areal data: neighbouring counties share climate, policy, infrastructure, and demographics, so their unexplained variation is correlated too. Two distinct failure modes follow, depending on where the dependence actually lives:

  • If the dependence is in the errors (an omitted spatially structured variable, or measurement error that clusters geographically), OLS coefficients remain unbiased, but the reported standard errors are wrong — usually too small, because the effective sample size is less than $n$ independent observations would give. Confidence intervals are too narrow and p-values overstate significance.
  • If the dependence reflects a genuine spillover — one county's outcome directly depends on its neighbours' outcomes, not just on shared covariates — then the true model has $y$ on both sides implicitly (through a term like $\rho W y$), and omitting that term biases every coefficient, not just the standard errors, because $Wy$ is correlated with the included regressors whenever they are themselves spatially patterned.

Chapter 3.2 puts numbers on the first failure mode directly, and Chapter 3.3's Lagrange Multiplier tests are built specifically to tell the two apart.

4. Two Mechanisms: Substantive vs. Nuisance Dependence¶

Substantive dependence Nuisance dependence
Interpretation A real causal channel: outcomes in one place influence outcomes in neighbouring places An artefact of the model: some spatially structured driver is missing from $X$, so its influence shows up in the residuals instead
Formal structure $y = \rho W y + X\beta + \varepsilon$ — the spatial lag (SAR/SLM) model $y = X\beta + u,\ u = \lambda W u + \varepsilon$ — the spatial error (SEM) model
What $\beta$ means Not a marginal effect — a change in $x_k$ anywhere propagates through $\rho$ to every other county (Chapter 3.4's impacts()) Retains its ordinary OLS interpretation; $\lambda \ne 0$ says "something spatially structured is missing", not "your $\beta$'s are wrong"
Correct response Model the spillover explicitly (spatial lag or Durbin model) Add the missing covariate if you can identify it, or model the nuisance directly (spatial error model)

This is Part 2, Chapter 2.1's dependence-vs-heterogeneity distinction, made operational: a regression's LM tests (Chapter 3.3) are a formal, data-driven way to ask which of the two — if either — is present, rather than guessing from theory alone.

5. A Taxonomy of Spatial Regression Models¶

$$ \text{OLS} \;\longrightarrow\; \text{SLX} \;\longrightarrow\; \text{SAR (SLM)} \;\longrightarrow\; \text{SEM} \;\longrightarrow\; \text{SARAR} \;\longrightarrow\; \text{SDM} \;\longrightarrow\; \text{SDEM} $$

  • OLS (Ordinary Least Squares): The basic, familiar regression model. It assumes each area's outcome is independent of its neighbors, except for shared covariates. No special spatial effects are included. Use this as a baseline.

  • SLX (Spatial Lag of X): Like OLS, but also includes "neighbor averages" of each predictor as extra inputs. For example, it doesn't just use a county’s own obesity rate, but also its neighbors’ average obesity rate. Helps capture the effect of neighboring conditions.

  • SAR / SLM (Spatial Lag Model): Adds the outcome values of neighboring areas ("Wy") directly into the regression. Here, a change in one county can ripple through to others—outcomes actually influence each other, beyond the covariates.

  • SEM (Spatial Error Model): Instead of outcomes influencing each other, this model assumes that some unmeasured, spatially-patterned factor is affecting the residuals (“leftover errors”). Neighboring counties’ errors are correlated, but outcomes themselves don’t cause each other.

  • SARAR: Combines SAR/SLM and SEM in one model, meaning both outcome spillovers and spatially-patterned residuals are present. Not implemented here.

  • SDM (Spatial Durbin Model): Combines SAR/SLM and SLX—outcomes and predictors can both spill over to neighbors. Good when both direct and indirect neighborhood effects may matter.

  • SDEM: Combines SEM and SLX—neighboring effects come through the predictors and the residuals, but not the outcome itself.

Model Structure Available in this library
OLS $y = X\beta + \varepsilon$ ✅ OLS
SLX $y = X\beta + WX\theta + \varepsilon$ (lagged covariates only, still OLS-estimable) ✅ OLS(lag_x=True) / SLX
SAR / SLM $y = \rho Wy + X\beta + \varepsilon$ ✅ SpatialLag
SEM $y = X\beta + u,\ u = \lambda Wu + \varepsilon$ ✅ SpatialError
SARAR $y = \rho Wy + X\beta + u,\ u = \lambda Wu + \varepsilon$ (both mechanisms at once) ❌ not implemented
SDM $y = \rho Wy + X\beta + WX\theta + \varepsilon$ (SAR + lagged covariates) ✅ SpatialDurbin (SpatialLag(slx=True))
SDEM $y = X\beta + WX\theta + u,\ u = \lambda Wu + \varepsilon$ (SEM + lagged covariates) ✅ SpatialError(slx=True)

Every model but SARAR is available; Chapter 3.4 works through SLM, SEM, and the Durbin variants (SDM, SDEM) in detail, and Chapter 3.5 shows the "common factor" relationship that connects SDM back to SEM.

6. Data Preparation¶

The running example for the rest of Part 3: diagnosed diabetes prevalence (Dia_pct) explained by obesity, physical inactivity, exercise access, primary-care access, food-environment quality, social vulnerability, and an urban/rural indicator, across 666 counties along the Atlantic seaboard — a larger, more geographically extensive dataset than Part 2's 217-county Northeast subset, chosen here so the regression has more power and a genuinely non-trivial spatial footprint.

In [2]:
counties = gpd.read_file(DATA / "diabetes_atlantic.shp")
print(f"diabetes_atlantic.shp: {counties.shape[0]} counties across {counties['STATE'].nunique()} states/DC")
print(counties["STATE"].value_counts())

diabetes_atlantic.shp: 666 counties across 11 states/DC
STATE
Georgia                 159
Virginia                128
North Carolina          100
Pennsylvania             67
New York                 62
West Virginia            55
South Carolina           46
Maryland                 24
New Jersey               21
Delaware                  3
District of Columbia      1
Name: count, dtype: int64

UrbRural is a string column ("Urban"/"Rural") and must be coded numerically before it can enter a regression.

In [3]:
covariates = ["Obesity", "PhysInact", "Acc_Exer", "PCP_rate", "FoodEnvIdx", "SVI"]
X = counties[covariates].copy()
X["Urban"] = (counties["UrbRural"] == "Urban").astype(int)
y = counties["Dia_pct"]

print(f"y: Dia_pct, range [{y.min():.2f}, {y.max():.2f}], mean {y.mean():.2f}")
X.describe().T.assign(missing=X.isna().sum())

y: Dia_pct, range [5.86, 16.96], mean 9.14
Out[3]:
count mean std min 25% 50% 75% max missing
Obesity 666.0 27.720706 4.792474 16.8000 23.5600 27.96000 31.260000 41.1800 0
PhysInact 666.0 21.540345 3.784371 13.2800 18.7450 21.17000 23.830000 35.1600 0
Acc_Exer 666.0 64.620270 22.826274 0.0000 50.0000 66.40000 82.350000 100.0000 0
PCP_rate 666.0 55.706156 38.313176 5.0000 32.8000 47.60000 70.550000 536.6000 0
FoodEnvIdx 666.0 7.541141 0.976839 2.6000 7.0000 7.60000 8.200000 10.0000 0
SVI 666.0 0.565418 0.275070 0.0013 0.3521 0.56605 0.812575 0.9987 0
Urban 666.0 0.714715 0.451890 0.0000 0.0000 1.00000 1.000000 1.0000 0

Spatial weights are needed from the first model onward (Chapter 3.2's residual diagnostics already require w), so we build them once here: Queen contiguity, row-standardised — the same default used throughout Part 2.

In [4]:
w = Queen(counties)
print(w.summary())

Spatial Weights (W)
===================
  Observations        666
  Transform           row
  Min neighbors       1
  Mean neighbors      5.423
  Max neighbors       11
  Islands             0
  Components          1
  Nonzero weights     3612

No islands on this larger, more contiguous dataset — unlike Part 2's Northeast subset (which had Nantucket as an isolated county), every Atlantic county has at least one Queen neighbour, so none of Part 2's island-handling caveats apply here. This matters for Part 3 specifically: every model below needs $W$ to represent a fully row-stochastic, invertible system (via $(I - \rho W)^{-1}$ for the lag model), which an island would complicate.

In [5]:
global_moran = Moran(y, w, permutations=999, seed=1)
print(f"Moran's I of raw Dia_pct: {global_moran.I:.4f}  (p_sim={global_moran['p_sim']:.4f})")

Moran's I of raw Dia_pct: 0.3707  (p_sim=0.0010)

Diabetes prevalence is itself strongly autocorrelated, exactly as in Part 2 — but that alone does not tell us why a regression of Dia_pct on these covariates would need a spatial model: if the covariates themselves fully explain the spatial pattern, plain OLS residuals could still be independent. Chapter 3.2 checks this directly, on the residuals rather than the raw outcome.

7. Summary and Conclusions¶

This chapter set up the conceptual and data foundation for spatial regression.

What we did¶

  1. Explained the two distinct ways OLS can fail on spatial data: nuisance dependence (wrong standard errors, unbiased coefficients) and substantive dependence (biased coefficients from an omitted spatial-lag term).
  2. Mapped that distinction onto the spatial error (SEM) and spatial lag (SAR/SLM) models respectively.
  3. Placed every model in this library's taxonomy — OLS, SLX, SAR, SEM, SDM, SDEM — noting that SARAR alone is not implemented.
  4. Prepared the Part 3 regression dataset: Dia_pct on six covariates plus a coded Urban indicator, across 666 Atlantic-seaboard counties, with a row-standardised, island-free Queen W.
  5. Confirmed the outcome variable is itself strongly spatially autocorrelated (Moran's I = 0.371, $p_\text{sim}=0.001$) — a stronger raw signal than Part 2's 217-county Northeast subset (I = 0.240) — before any covariate has even been considered.

Practical takeaways¶

  • Do not assume which mechanism (substantive vs. nuisance) is at work from theory alone — Chapter 3.3's LM tests exist because the two are formally distinguishable from the same OLS residuals.
  • A string categorical column silently breaks a regression if not coded numerically first; always check X.dtypes before fitting.
  • Confirming zero islands in $W$ up front avoids a class of failures (a non-invertible $(I-\rho W)$) that would otherwise surface confusingly deep inside Chapter 3.4's maximum-likelihood fitting.

Next steps¶

Chapter 3.2 (OLS Baseline and Diagnostics) fits the OLS model on exactly this data and runs the classical and spatial diagnostics — condition number, normality, heteroskedasticity, and residual Moran's I — that Chapter 3.3 builds on.

8. Resources¶

Textbooks¶

Reference Notes
Anselin, L. (1988). Spatial Econometrics: Methods and Models. Kluwer Academic. The foundational text for this entire Part; origin of the substantive/nuisance distinction
LeSage, J. & Pace, R. K. (2009). Introduction to Spatial Econometrics. CRC Press. Modern standard reference; origin of the impacts framework used in Chapter 3.4
Elhorst, J. P. (2014). Spatial Econometrics: From Cross-Sectional Data to Spatial Panels. Springer. Clear taxonomy of the SLX/SAR/SEM/SDM/SDEM model family

Key papers¶

Paper Contribution
Anselin, L. (1988). Lagrange multiplier test diagnostics for spatial dependence and spatial heterogeneity. Geographical Analysis, 20, 1–17. Formal test statistics used in Chapter 3.3
LeSage, J. P. & Pace, R. K. (2009). Chapter 2, spatial model taxonomy. In Introduction to Spatial Econometrics. The OLS→SLX→SAR→SEM→SARAR→SDM→SDEM taxonomy reproduced in Section 5

In this library¶

Object Role
spatialstats.regression.OLS Baseline model; lag_x=True gives SLX
spatialstats.regression.SpatialLag, SpatialError SAR/SLM and SEM, by maximum likelihood
spatialstats.regression.SpatialDurbin SpatialLag(slx=True), i.e. SDM
spatialstats.weights.Queen Spatial weights, reused directly from Part 2
Chapter 3.2 lm_tests, moran_residuals, classical diagnostics
Chapter 3.4 Full treatment of SLM, SEM, SDM, and impacts()