Spatial Difference-in-Differences Overview

• Spatial Difference-in-Differences is a causal framework that adapts standard DiD techniques to account for spatial spillovers and interference in geographic or network settings.
• It employs various design strategies—such as ring designs, exposure mappings, and continuous functional approaches—to capture treatment effects over space and time.
• Accurate inference requires robust methods like spatial HAC and cluster-robust estimators to adjust for spatially correlated errors and interference bias.

Spatial Difference-in-Differences (Spatial DiD) is a class of causal designs that adapts difference-in-differences to settings in which treatment assignment, exposure, and potential outcomes are indexed by geographic or network proximity. The central difficulty is that treatment is rarely binary and isolated in space: nearby units can be partially exposed, interference violates SUTVA, the boundary between treated and control areas may be unknown or time-varying, and residual shocks can be spatially correlated even when treatment does not spill over (Butts, 2021, Ferman, 2019, Kikuchi, 17 Oct 2025). The literature therefore replaces a single treated indicator with location-specific distances, exposure mappings, boundary-aware local comparisons, or continuous space-time treatment fields, while preserving the core DiD logic of comparing outcome changes under parallel counterfactual evolution (Xu, 2023, Butts, 2021, Butts, 2021, Kikuchi, 13 Oct 2025).

1. Core identification problem

In the canonical two-period setting, the DiD estimator is

$$\hat{\beta}_{DID} = (\bar{Y}_{T}^{post} - \bar{Y}_{T}^{pre}) - (\bar{Y}_{C}^{post} - \bar{Y}_{C}^{pre}),$$

or equivalently a two-way fixed-effects regression,

$$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$

Spatial DiD departs from this benchmark because outcomes may depend on own treatment and on a spatial exposure term, written in one strand as $Y_{it}(D_i,h_i(D))$, where $h_i(D)$ aggregates the effect of nearby treated units (Butts, 2021).

Under random sampling, no anticipation, and parallel counterfactual trends, classical DiD under spillovers satisfies

$$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$

so the standard estimator no longer recovers the total effect on the treated. The bias has two sources: untreated units cease to represent the counterfactual trend because their outcomes are affected by treatment, and treated outcomes combine the effect of own treatment with the effect from the treatment status of close units (Butts, 2021).

Spatially correlated errors are analytically distinct from spillovers. If treatment does not affect untreated units but nearby units share shocks, identification can still hold while inference fails. In that case the error covariance can be written as

$$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$

and cluster-robust variance estimators that assume independence across clusters understate true variance and produce inflated rejection rates (Ferman, 2019). With few treated units, this problem becomes more severe: some single-treated procedures remain asymptotically valid under weak dependence, but with more than one treated unit naive placebo-style resampling can be invalid (Alvarez et al., 2020).

Many modern staggered-adoption DiD estimators are linear combinations of $2\times 2$ DiD contrasts, so the same spatial contamination and spatial inference problems carry over to cohort-specific and event-study settings rather than disappearing with newer estimators (Ferman, 2019).

2. Exposure mappings and causal estimands

A dominant formalization of Spatial DiD summarizes interference through an exposure mapping. In a two-period interference model, neighbor treatments are reduced to an exposure state $G_i = G(i, W_{-i}) \in \mathcal{G}$, and potential outcomes are indexed as $\tilde{y}_{i2}(w_i,g)$. This supports exposure-specific causal targets such as the direct average treatment effect on the treated at exposure level $g$, spillover effects among treated and untreated units, and the misspecification-robust EDATT when the researcher’s exposure map is only an approximation. Identification uses overlap, no anticipation, and a modified parallel trends condition defined within exposure levels rather than within binary treatment groups (Xu, 2023).

A related network- and spatial-dependence framework defines exposure histories through an interference weights matrix. With exposure

$$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$0

the causal estimand is the Average Exposure Effect on the Exposed,

$$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$1

Identification relies on no anticipation, positivity of exposure histories, and conditional parallel trends under a reference exposure history. The proposed doubly robust estimator remains consistent if either the exposure propensity score or the outcome regression trend is correctly specified, and under the paper’s network dependence conditions it is asymptotically normal with variance attaining the semiparametric efficiency bound (Jetsupphasuk et al., 5 Feb 2025).

In staggered adoption with interference, the literature separates own-adoption effects from spillover effects generated by other adopters. One formulation defines the Dynamic Switching Effect, the Control-State Spillover Effect, and the Dynamic Total Effect,

$$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$2

for cohort $$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$3 at event time $$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$4. Identification compares treated cohorts to never-treated units with the same baseline and target-date exposure states, and spillover responses are learned from never-treated units and transported to the exposure distribution faced by treated cohorts (Tagawa, 14 May 2026).

A more continuous representation replaces discrete exposure states by a space-time treatment field. In that approach, the exposure of unit $$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$5 at time $$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$6 is written as

$$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$7

where the kernel $$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$8 may be exponential, logarithmic, or power-law. This embeds interference directly into the regressor of a spatial event-study specification,

$$y_{it} = \alpha_i + \lambda_t + \beta D_{it} + \varepsilon_{it}.$$9

and converts binary treatment timing into dynamic dose variation over space and time (Kikuchi, 17 Oct 2025).

3. Main design families

Spatial DiD is not a single estimator but a family of designs that differ in how they encode spatial exposure and how they construct counterfactual comparisons.

Design family	Main object	Typical maintained restriction
Rings and distance-curve DiD	$Y_{it}(D_i,h_i(D))$0 or $Y_{it}(D_i,h_i(D))$1	Local parallel trends; zero-effect region beyond some distance
Geographic difference-in-discontinuities	$Y_{it}(D_i,h_i(D))$2	Local continuity and time-invariant border discontinuity
Exposure-mapping DiD	DATT, EDATT, spillovers, AEE	Modified or conditional parallel trends within exposure strata
Continuous functional framework	$Y_{it}(D_i,h_i(D))$3, $Y_{it}(D_i,h_i(D))$4, $Y_{it}(D_i,h_i(D))$5	Diffusion/advection scope conditions and exposure support

With geocoded microdata, a common design compares an inner ring of units near a point treatment to an outer ring just further away. In panel form, the first-differences specification is

$Y_{it}(D_i,h_i(D))$6

while a more flexible alternative estimates a treatment effect curve $Y_{it}(D_i,h_i(D))$7 over distance using partitioning-based least squares. The ring design requires local parallel trends and, in its classical form, correct knowledge of the treatment’s spatial reach so that the outer ring is unaffected. The distance-curve approach weakens this by requiring only that treatment effects become zero somewhere within a sufficiently large control radius $Y_{it}(D_i,h_i(D))$8 (Butts, 2021).

At administrative borders, geographic difference-in-discontinuities combines geographic RD with a pre/post difference. The key estimand is

$Y_{it}(D_i,h_i(D))$9

the difference between the post-treatment and pre-treatment border jumps. This design is intended for settings in which sorting and compound treatments make a cross-sectional geographic RD problematic, but these border-aligned discontinuities are stable over time and can therefore be differenced out (Butts, 2021).

A related but distinct local-space design is Spatial First Differences, which differences neighboring units in the cross section:

$h_i(D)$0

It is presented as the spatial analog of first differences in time and as a cross-sectional counterpart to DID. Although it is not itself a panel DiD estimator, it shares the same logic of purging locally smooth unobservables by differencing nearby units rather than the same units over time (Druckenmiller et al., 2018).

4. Boundaries, decay structures, and continuous-functional Spatial DiD

One recent line of work reformulates spatial treatment as a continuous field $h_i(D)$1 defined over space-time:

$h_i(D)$2

Its dynamics are governed by an advection-diffusion equation,

$h_i(D)$3

with diffusion coefficient $h_i(D)$4, advection field $h_i(D)$5, and source term $h_i(D)$6. This framework also defines cumulative exposure

$h_i(D)$7

and spatial gradients $h_i(D)$8, which are policy-relevant quantities absent from standard banded DiD (Kikuchi, 17 Oct 2025).

In this continuous-functional formulation, treatment boundaries are level sets or thresholds of the field. For isotropic exponential decay,

$h_i(D)$9

the boundary corresponding to a relative threshold $$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$0 is

$$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$1

For Gaussian diffusion from an instantaneous point source,

$$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$2

and the sensitivity of the boundary to the diffusion parameter is

$$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$3

The same literature compares exponential, logarithmic, and power-law steady-state decay functions and treats the sign of the decay parameter as a diagnostic of whether a diffusion interpretation is plausible (Kikuchi, 17 Oct 2025).

A closely related physics-based framework derives spatial and temporal treatment-effect boundaries from advection-diffusion-reaction dynamics. In the low-Péclet, low-Reynolds, linear-decay regime, the scalar transport equation reduces to a Helmholtz equation and yields exponential far-field decay. The detectable spatial boundary is

$$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$4

with the approximation

$$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$5

Applicability is tied to explicit scope conditions based on dimensionless numbers: $$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$6, $$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$7, and $$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$8. Negative estimated $$E\big[ \hat{\tau}^{DiD} \big] = \tau_{total} - \tau_{spill}(0),$$9, directional asymmetry, or poor fit to exponential decay are treated as diagnostics of scope failure rather than as mere estimation noise (Kikuchi, 13 Oct 2025).

The empirical healthcare application of the continuous-functional framework uses 32,520 ZIP Code Tabulation Areas, hospital coordinates from HIFLD, CDC PLACES outcomes, and synthetic ZCTA-level demographics. For ACCESS2 at the ZCTA level, the refined exponential fit gives $$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$0 per km with $$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$1, and the 10% detectable boundary is $$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$2 km with 95% CI $$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$3. Model selection strongly favors logarithmic decay over exponential, with $$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$4, while power-law and logarithmic decay are nearly tied with $$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$5. The same application reports 2–13 $$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$6 stronger distance effects for elderly populations and education gradients in which high education reduces distance sensitivity by 5–13 $$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$7 (Kikuchi, 17 Oct 2025).

The air-pollution application analyzes 791 ground-based PM$$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$8 monitors and 189,564 satellite-based NO$$\operatorname{Cov}(\varepsilon_{it},\varepsilon_{jt'}) = \rho\big(d_{ij}, |t-t'|\big),$$9 grid cells in the Western United States over 2019–2021. Within 100 km of coal plants, both pollutants show positive spatial decay: for PM$2\times 2$0, $2\times 2$1 km$2\times 2$2 and $2\times 2$3 km; for NO$2\times 2$4, $2\times 2$5 km$2\times 2$6 and $2\times 2$7 km. Beyond 100 km, estimated decay parameters become negative, which the framework interprets as evidence that urban sources dominate and that diffusion-based Spatial DiD is not appropriate in that regime. The paper states that the framework applies in four of eight analyzed regional splits and rejects in the remainder (Kikuchi, 13 Oct 2025).

The nonparametric distance-curve literature reaches a similar substantive conclusion from a different direction: ring averages can conceal steep local gradients and misstate reach. In the reanalysis of Linden and Rockoff (2008), a two-ring DiD implies an average decline of approximately 7.5% within 0–0.1 miles, whereas the estimated nonparametric $2\times 2$8 shows declines around 20% within a few hundred feet and non-significant effects farther away within the same inner ring; beyond approximately 0.1 miles, $2\times 2$9 is centered at zero (Butts, 2021).

5. Estimation and inference under spatial dependence

Spatial DiD requires separate treatment of identification and inference. Even when parallel trends are credible, cross-sectional dependence in residuals can invalidate standard standard errors. Cluster-robust variance estimators are valid only with many independent clusters, treated and control groups with similar exposure to common shocks, short windows, and weak or short-range spatial correlation. They are invalid under strong cross-cluster dependence, persistent common shocks, long windows, few clusters, or misspecified cluster boundaries (Ferman, 2019).

The dependence of inference distortions on the time horizon is explicit in the theory. When common shocks are more temporally persistent than idiosyncratic errors, size distortions increase with the length of the pre/post window, and analogous distortions grow with event-study lag length. This is why estimators that emphasize short time differences, or many local $G_i = G(i, W_{-i}) \in \mathcal{G}$0 contrasts, can mitigate spatial inference distortions relative to long-horizon TWFE (Ferman, 2019).

Several variance estimators and resampling procedures appear in the literature. Conley spatial HAC estimators are the workhorse when the relevant distance metric and sufficient granularity are available, but they can be unstable or even infeasible when the distance metric has little support. Wild cluster bootstrap can reduce over-rejection with moderate numbers of clusters but does not repair misspecification of the dependence metric. Driscoll–Kraay requires large $G_i = G(i, W_{-i}) \in \mathcal{G}$1 and is often ill-suited to short DiD panels. Randomization or permutation inference is useful when assignment is known at a cluster level and the randomization scheme mirrors the actual assignment structure. Common Correlated Effects can help with spatially correlated shocks under fixed $G_i = G(i, W_{-i}) \in \mathcal{G}$2, but it imposes factor-structure restrictions and can be problematic with heterogeneous treatment effects (Ferman, 2019).

Under explicit interference, inference is built around spatially robust score covariance estimators. One approach uses a SHAC variance estimator in a GMM system for doubly robust DiD with interference, where the score covariance is kernel-weighted over spatial distances and can be conservative for the conditional finite-population variance (Xu, 2023). A related network-dependent framework derives an efficient influence function for the doubly robust estimator and estimates the asymptotic variance with a network HAC-type kernel over graph or geographic distance; dependent-data cross-fitting is proposed to accommodate flexible nuisance estimation under weak dependence (Jetsupphasuk et al., 5 Feb 2025).

Few treated units create an additional complication. When there is a single treated unit, placebo-style methods derived from the empirical distribution of control residuals remain asymptotically valid under weak dependence. With multiple treated units, however, inference must account for treated-treated covariance. Conservative Test 1 approximates a worst-case perfect-correlation scenario among treated units, while Conservative Test 2 uses aggregate-outcome structure $G_i = G(i, W_{-i}) \in \mathcal{G}$3 to achieve less conservative inference when unit sizes $G_i = G(i, W_{-i}) \in \mathcal{G}$4 are available. The same paper recommends Conservative Test 2 as the default for multiple treated units with aggregate outcomes; multiple-testing adjustments based on single-treated valid tests are also available but tend to yield wider confidence sets (Alvarez et al., 2020).

The more recent continuous-functional spatial event-study literature adopts the same broad practice: use spatial HAC or recent robust procedures for long-range spatial dependence, cluster at appropriate geographic levels such as counties or hospital catchment areas, and align bandwidths and cutoffs with estimated boundaries rather than arbitrary rings (Kikuchi, 17 Oct 2025).

6. Diagnostics, applications, and limitations

Diagnostics in Spatial DiD are unusually central because the main threats are often visible in the spatial structure itself. The inference literature recommends Moran’s I, semivariograms, residual correlograms, and distance-binned correlations to detect residual spatial dependence, as well as sensitivity analysis over cluster definitions and Conley bandwidths (Ferman, 2019). Interference-oriented DiD adds placebo pre-trends within exposure strata, overlap checks across $G_i = G(i, W_{-i}) \in \mathcal{G}$5 cells, and robustness to alternative exposure mappings or neighborhood radii (Xu, 2023). The continuous-functional and physics-based approaches add sign-reversal and symmetry diagnostics: positive estimated decay near sources is taken as evidence for diffusion-like attenuation, whereas $G_i = G(i, W_{-i}) \in \mathcal{G}$6, downwind-upwind asymmetry, or poor exponential fit indicate that the maintained propagation mechanism is failing (Kikuchi, 13 Oct 2025, Kikuchi, 17 Oct 2025).

The empirical range of Spatial DiD is broad. In healthcare access, the estimated 37.1 km ACCESS2 boundary implies a concrete treatment zone around hospitals, while the same study argues that targeting elderly and low-education rural populations yields larger returns because their distance sensitivity is much stronger (Kikuchi, 17 Oct 2025). In the Community Health Centers rollout, the staggered-spillover framework estimates negative own-adoption and spillover effects on older-adult mortality, with the Control-State Spillover Effect accounting for about 40% of the Dynamic Total Effect across event times $G_i = G(i, W_{-i}) \in \mathcal{G}$7 to $G_i = G(i, W_{-i}) \in \mathcal{G}$8 (Tagawa, 14 May 2026). In place-based policy evaluation, accounting for spillovers changes conclusions materially: the TVA reanalysis reports a more negative total effect on agriculture and a smaller positive total effect on manufacturing once spillovers by distance are modeled, and in Opportunity Zones negative spillovers just outside treated zones help explain why neighbor comparisons are much larger than eligible-only comparisons (Butts, 2021).

The design family also varies in its external-validity claims. Geographic DiDisc is explicitly local to borders and is designed to difference out time-invariant sorting and compound treatments rather than to estimate broad regional effects (Butts, 2021). Ring and distance-curve designs are local to a radius $G_i = G(i, W_{-i}) \in \mathcal{G}$9 within which local parallel trends are thought plausible (Butts, 2021). Physics-based continuous-functional approaches are broader in spatial ambition, but only under explicit scope conditions such as $\tilde{y}_{i2}(w_i,g)$0, $\tilde{y}_{i2}(w_i,g)$1, and $\tilde{y}_{i2}(w_i,g)$2 (Kikuchi, 13 Oct 2025).

The main limitations recur across frameworks. Cross-sectional diffusion fits are correlational and need panel exposure DiD or instruments to strengthen identification (Kikuchi, 17 Oct 2025). Exposure mappings can be misspecified, in which case interpretation shifts from a causal effect at the true exposure to an EDATT or other estimand defined by the chosen map (Xu, 2023). Unknown dependence metrics make inference fragile, especially with few clusters or limited distance variation (Ferman, 2019). Staggered-adoption spillover designs do not generally identify a global pure direct effect without additional no-interaction structure, so the literature often treats switching, spillover, and total effects as the primary estimands instead (Tagawa, 14 May 2026).

This suggests that Spatial Difference-in-Differences is best understood as a design family rather than a single estimator. Its unifying principle is unchanged from standard DiD: identify causal effects from outcome changes under parallel counterfactual evolution. What changes is the object to which that logic is applied—rings, border neighborhoods, exposure states, network histories, or a continuous treatment field—and the corresponding need to model interference, detect boundaries, and conduct inference in the presence of spatial dependence.