Modeling the spatial spread of emerging infectious diseases in highly heterogeneous environments

John M. Drake

Odum School of Ecology, University of Georgia

Ecological Society of America Annual Meeting

Salt Lake City, UT

July 29, 2026

West Nile virus, 1999–2005

Too coarse to explain, or too costly to fit

Heterogeneity is the rule — climate, land use, hosts, and vectors all change from county to county. Three ways we usually model spread through it:

Diffusion — traveling waves on a homogeneous plane; the heterogeneity is washed out.
Gravity models — fit well, but phenomenological: one exponent conflates production, movement, and establishment.
Coupled metapopulation ODEs — full dynamics at every location, but overparameterized — fittable only with rich panel data (Li, Ionides, King, Pascual & Ning 2024).

But we rarely have panel data. What we usually do have is first-report data — the first WNV case in a county, the first mussel in a lake. We want a model we can fit to that.

Why not just write down the full mechanistic model? Each option costs something. Diffusion — reaction–diffusion waves — assumes a homogeneous plane, and the landscape isn’t homogeneous. The gravity model fits well but is phenomenological: one distance exponent conflates how infection is produced, how it moves, and where it establishes. The opposite extreme is a full network of coupled ODEs — the mechanistic ideal, and King and colleagues showed it can be fit for early-pandemic COVID — but it needs dense case time series at every node. For most invasions, West Nile included, we don’t have that. What we have is first-report data: the week a county first reports the virus, the year a lake first reports a mussel. So we need a mechanistic model we can fit to first-report data. That’s what I’ll build.

What we infer is which landscape features drive invasion — and how

The spatial units are handed to us by the data — counties for WNV, lakes for mussels. Spatial resolution isn’t something we estimate.

What is unknown: which features of the landscape shape invasion, and through what functional form.

So we want a framework where each candidate mechanism can be added parametrically — every feature enters as estimable terms we can test.

A framework for spread in heterogeneous landscapes

Spread decomposes into production, movement, and observation

Every invasion event decomposes into a source that produces infection, a redistribution process that carries it across space, and a destination that receives it.

Each process enters at its own biological scale

For an uninfected county j at week t, the force of infection sums over all currently infected sources i:

\lambda_{ij,t} = \underbrace{K\vphantom{S_{i,t}}}_{\text{scale}}\;\underbrace{D_{j,t}}_{\text{destination}}\;\underbrace{R_{ij,t}}_{\text{redistribution}}\;\underbrace{S_{i,t}}_{\text{source}} \qquad \Lambda_{j,t} = \sum_{i \,\in\, \text{infected}} \lambda_{ij,t}

Each component is parameterized independently
Covariates enter at the appropriate mechanism
The classical gravity model is a special case

A likelihood turns the framework into inference

Each week, every at-risk county experiences infection (or invasion) as a Poisson process with rate \Lambda_{j,t}:

-\log \mathcal{L} = \sum_{t}\left[\; -\!\!\sum_{j \,\in\, \text{new}(t)} \log\!\left(1 - e^{-\Lambda_{j,t}}\right) \;+\!\! \sum_{j \,\text{ remains at risk at } t} \Lambda_{j,t}\;\right]

Fit directly to surveillance data by maximum likelihood → parameter estimates, standard errors, and formal model comparison (AIC).

West Nile virus in North America

Macroscopic patterns

Cumulative U.S. counties reporting West Nile virus, 1999–2009 (inset: year of first report): fast continental spread, seasonal steps, and an early long-distance jump to Florida and Georgia in 2001.

Baseline model overpredicts the early part of the invasion and underpredicts the late part of the invasion

A homogeneous model — no landscape heterogeneity, no seasonality, one constant rate — spreads the invasion evenly across the decade, overpredicting the early years and underpredicting the sharp increase in 2002.

Hypothesis 1: Heterogeneity in Culex habitat

Modeled habitat suitability for the three predominant Culex vectors (Gorris et al. 2021), averaged to a single surface and aggregated by county — highest across the upper Midwest and California’s Central Valley.

Space and Culex habitat: a big gain, but the timing is still wrong

\lambda_{ij,t}=K\,(\textit{Culex}_i)^{e^{s_1}}\,(\textit{Culex}_j)^{e^{d_1}}\,\textit{dist}_{ij}^{\,-e^{r_1}},\qquad K=10^{\alpha}

A power-law distance kernel plus Culex suitability at the source drops AIC from 36,431 to 27,206 — a large gain — yet the simulated invasion still runs ahead of the real one and climbs too smoothly.

Culex habitat at the destination adds nothing

Also letting the destination depend on Culex suitability leaves the fit unchanged — AIC 27,208, essentially identical — so receptivity at the destination carries no signal.

Hypothesis 2: Temperature-dependent transmission

Temperature at the source fixes the timing

\lambda_{ij,t}=K\,[\,\textit{Culex}_i\,\mathcal{R}_{\text{rel}}(T_{i,t})\,]^{e^{s_1}}\,(\textit{Culex}_j)^{e^{d_1}}\,\textit{dist}_{ij}^{\,-e^{r_1}}

Scaling the source by the Culex thermal-performance term drops AIC to 24,134 and finally reproduces the seasonal steps — for no extra parameter.

Hypothesis 3: Migratory bird movement gives the best model

\lambda_{ij,t}=K\,[\,\textit{Culex}_i\,\mathcal{R}_{\text{rel}}(T_{i,t})\,]^{e^{s_1}}\,(\textit{Culex}_j)^{e^{d_1}}\,(\,\textit{dist}_{ij}^{\,-e^{r_1}}+e^{r_2}\,\textit{bird}_{ij,t}\,)

Adding weekly movement of 27 WNV-competent bird species to the redistribution kernel gives the best fit — AIC 23,921.

Temperature still leaves one thing out: that early jump to Florida and Georgia, the virus in the Southeast in 2001, years ahead of the main wave. Distance alone can’t produce a jump that far. We know how it got there — migrating birds — and they aren’t in the model yet. So we add them. We take the twenty-seven North American bird species the CDC has detected West Nile in, from its ArboNET dead-bird surveillance. For each, a BirdFlow model, fit to eBird’s weekly abundance maps, predicts how the species moves week to week; we weight it by population size from Rosenberg and colleagues’ 2019 estimates. Summed, that’s a weekly county-to-county flow of competent birds, added to the kernel alongside distance. Now a county can be seeded by a nearby source or a migratory link to a distant one. It’s the best model — AIC twenty-three thousand nine hundred — and the bird term supplies the long-range jump that reached Florida and Georgia early. One caveat: the timing of the main wave still isn’t quite right — the invasion holds back through 2002 and breaks out in 2003, and the simulation runs a little ahead. But the sequence is complete: space, habitat, temperature, birds.

The model reproduces the wave — and shows us what’s missing

Difference between the best model’s simulated first-report time and the observed, county by county: mostly pale (the east-to-west wave is recovered), but early across the interior (blue) and late at the northeastern origin (red).

Space sets the shape; temperature sets the timing

Model	What it adds	k	NLL	AIC	Delta AIC
Homogeneous null	no space, no heterogeneity	1	18215	36431	12510
Culex source	+ distance kernel + source suitability	3	13600	27206	3284
Culex destination	+ destination suitability	4	13600	27208	3286
Thermal source	+ temperature at source	4	12063	24134	212
Thermal + bird kernel	+ bird migration	5	11956	23921	0

Space + habitat — a distance kernel and source suitability are the huge first gain over the null
Thermal source — temperature at the source fixes the timing, for no extra parameter
Thermal + bird kernel — adding bird migration gives the best model

A general template for emerging diseases in complex landscapes

The framework

Decompose spread into source · redistribution · destination · observation
Fit directly to surveillance by likelihood
Let model comparison find the minimal structure

Where to intervene

Source — transmission control
Redistribution — movement & vectors
Destination — susceptibility

Heterogeneity is the rule.

Step back from West Nile — the framework is the contribution. Any invasion in a heterogeneous landscape can be written this way: a source, a redistribution, a destination, and an observation model, fit to the first-report data we typically have, with model comparison deciding how much structure the system demands. Because the mechanisms are separate, it maps onto management: vector control acts on the source, movement and vector dispersal on redistribution, host susceptibility on the destination. The same four terms can be written for the next emerging pathogen and fit to its first reports, and model comparison identifies which part of the landscape carries the signal — and where an intervention would act. That’s the contribution: not a model of West Nile, but a way to ask, for any invasion in a varying landscape, how much of that variation the spread used.

Acknowledgements

Collaborators — Eric Marty, JP Schmidt, Suzanne O’Regan

Data & methods

Dispersal / S–R–D framework — Jongejans et al. 2015, Theoretical Ecology 8:207
Culex suitability — Gorris et al. 2021, Parasites & Vectors 14:547
Mosquito thermal biology — Mordecai et al. 2019, Ecology Letters 22:1690
Bird abundances — Rosenberg et al. 2019, Science 366:120
Bird movement (BirdFlow) — Fuentes et al. 2023, Methods in Ecology and Evolution 14:923
Metapopulation inference — Li et al. 2024, J. R. Soc. Interface 21:20240217

Backup

Fitted parameters (Thermal + bird kernel)

Parameter	Symbol	Estimate	Interpretation
Scale	\alpha	-4.24	K = 10^{\alpha}
Source exponent	e^{s_1}	0.37	sub-linear in suitability × temperature
Distance exponent	e^{r_1}	1.63	power-law decay, heavy tail
Bird term	e^{r_2}	\approx 1.3\times10^{3}	migration flux added to kernel
Destination exponent	e^{d_1}	0.55	weak destination effect

The observation layer: right-censoring and detection

We fit to infected counties; every not-yet-infected county contributes a survival (zero-arrival) term each week.
The likelihood cleanly separates process from observation, so reporting can be modeled rather than assumed perfect.
Extensions: explicit detection probability and subsampling calibration of the non-infected pool (as in the invasive-species version of this framework).

The model sequence

Baseline gravity distance kernel, suitability fixed at source · Fitted source free source suitability · Fitted destination free destination
Thermal source temperature seasonality at source (best “for free”) · Thermal both ends temperature at destination (rejected)
Thermal + bird kernel bird migration in the kernel (best overall)
Each model adds one mechanism to the one before it — a nested sequence, so AIC differences are honest.