Edward Turk

Master's research · Thesis in progress

Mapping soil moisture from vegetation phenology

Can satellite-observed vegetation predict soil moisture across an entire watershed?

My master's research predicts soil moisture from space. The output is a continuous 10 m map of the upper East River valley near Crested Butte, Colorado, one surface per year from 2018 to 2025. I built it for terrain that is difficult, and in places impossible, to instrument directly.

The question came from the field, and it started with methane. On the Slate River floodplain we ran collaborative flux campaigns in 2024 and 2025 with groups from Stanford and Cal State East Bay, taking roughly sixty chamber measurements across five sites laid out along a distance-from-river transect. The pattern was clear: methane flux rose with soil moisture, and the wettest sites were the active ones. That's a satisfying result over a thousand square metres of floodplain. It's also a trap. The processes that produce methane happen in anoxic microsites millimetres across. Moisture gradients running for kilometres switch those microsites on and off. You can't put a chamber on every hillside of a 6,000 km² watershed.

So I stopped asking about methane and asked a mapping question instead. If you can't measure the subsurface everywhere, can you infer it? The whole project rests on vegetation as a sensor: roots connect subsurface moisture to the canopy, plants integrate the water available to them over a whole season, and a satellite can watch that integration happen. Greenness is not moisture. But the timing and shape of a season of greenness carries information about the moisture underneath it, and unlike a soil probe, a satellite covers every pixel of the valley, every year, for free.

Type
Master's research
Status
Thesis in progress
Domain
Water resources
Tech
Sentinel-2 · R · Python · GEE
Two-panel map: mean predicted 10 m soil moisture across the study area, and the year-to-year coefficient of variation of those predictions.
The output: continuous 10 m soil moisture for ground nobody can walk, and a companion surface showing where the prediction is steady from year to year and where it swings. (a) Mean predicted volumetric water content, 2018–2025, from the 2018-trained model. (b) Inter-annual coefficient of variation of the same predictions. Because the topographic and vegetation-community inputs are held fixed across years, panel (b) is a sensitivity surface for the year-varying satellite term. It is not a map of model accuracy and not observed field variability.

How it works

I turn ten years of satellite imagery into a small set of vegetation-timing measurements, then relate those to field soil-moisture readings with a model built to be tested on sites it has never seen. That last clause is the design constraint that shaped everything else: a map is only worth making if the model behind it works somewhere nobody sampled.

A word on each half. On the satellite side I fit a third-order Fourier harmonic regression to every pixel's NDVI record, one fit per year, which turns a ragged, cloud-gapped series of scene-dates into a smooth seasonal curve you can actually measure. From that curve come the phenology metrics: how green a pixel was in June, how green it got at its peak, how green it stayed in September, and how far it swung between them. The field half is 437 soil-moisture observations across 17 sites. The topographic baseline the vegetation signal has to beat is HAND (height above nearest drainage), computed from a LiDAR DEM. That's the standard physically-grounded answer to "where does water collect?"

The two halves meet in a pooled linear mixed-effects model with a random intercept per site, fit in R with lme4. I validate it leave-one-site-out: hold a whole site back, fit on the rest, predict the held-out site from the fixed effects only. No site intercept, because a new location on a map doesn't come with one. That's a deliberately unflattering test, and it's the only one that resembles what the model has to do once I deploy it as a raster.

  1. ~419 Sentinel-2 scene-dates (2016–2025)
  2. per-pixel 3rd-order Fourier harmonic regression
  3. phenology metrics: June mean, Peak NDVI, Sep mean, seasonal amplitude
  4. 437 soil-moisture observations · 17 sites
  5. pooled linear mixed-effects model (site random intercepts)
  6. leave-one-site-out validation
  7. annual 10 m VWC maps
Model specification and supporting analyses

Specification. The deployed model is VWC ~ wetness index + vegetation community + Peak NDVI + (1 | site), fit with lme4 in R. Vegetation community enters as a between-site fixed effect: sites in different plant communities sit at different moisture levels, and telling the model which community a location belongs to keeps that offset off the spectral term.

Metric selection. I entered candidate spectral predictors one at a time and scored them under leave-one-site-out cross-validation rather than by in-sample fit. Peak NDVI won on every out-of-site measure.

The PCMA residual. The pipeline also produces a phenology-conditioned moisture anomaly (PCMA): a per-year regression of late-season NDVI on early-season phenology and vegetation class, whose per-pixel residual I read as a moisture anomaly. It's a genuinely useful diagnostic product. The residual field is spatially coherent (variogram range ~274 m, roughly 65% of its variance spatially structured) and stable year to year. But in the head-to-head it lost. Peak NDVI out-predicted the PCMA residual at held-out sites on every metric, and adding the residual on top of Peak NDVI bought nothing. The simpler metric carries the transferable signal. I'm reporting that as a result, not hiding it as a failure.

The snowmelt bridge. A separate analysis at SNOTEL station 380 supplies the mechanistic why. Deep soil moisture (20 in) tracks snowmelt duration (r = +0.909, p = 0.012) and 7-day maximum melt rate (r = +0.831, p = 0.040), but not peak snow-water equivalent. Melt dynamics recharge the deep soil, not snowpack size, and deep moisture is what feeds the early-season productivity the phenology metrics later detect. Six water years. Exploratory, and labelled that way in the manuscript.

Full pipeline, notebooks and figure crosswalk: see the repository.

Three panels: a map of mean phenology-conditioned moisture anomaly residuals, a map of how consistently each pixel's residual keeps the same sign across years, and a semivariogram.
The PCMA residual field is a coherent, repeatable pattern rather than noise. A real anomaly surface, even though it lost the head-to-head to Peak NDVI. (a) Mean PCMA residual, 2018–2025 (blue wetter than expected, red drier). (b) Modal-sign agreement across the eight annual residuals. (c) Empirical semivariogram of the mean residual. A diagnostic of the residual field's internal structure, not independent validation against field moisture.
Four panels relating snowmelt duration, melt rate and peak snow-water equivalent to soil moisture at two depths at a SNOTEL station.
It is how fast and how long the snow melts, not how much fell, that recharges the deep soil the vegetation later draws on. SNOTEL station 380: deep (20 in) soil moisture against snowmelt duration and 7-day maximum melt rate, and against peak snow-water equivalent. Shallow (2 in) moisture shown for contrast. Six water years, exploratory, small-n.

What it found

Bar chart of the percentage of model skill lost when each predictor is dropped, shown separately for within-site fit and for leave-one-site-out prediction.
The inversion, in one chart: Peak NDVI contributes least when the model is scored on sites it was fit on, and most when it is scored on sites it was not. Drop-one (LOVO) contribution for Peak NDVI, topographic wetness and vegetation community, in-sample versus leave-one-site-out, 2018: the percentage of full-model skill lost when each predictor is removed. Predictors are correlated, so each value is the skill a term adds on top of the other two in this model, not an independent share of variance.

In the maps

The argument runs as a sequence of figures: where the data came from, why vegetation should carry a moisture signal at all, what the satellite actually measures, where the spectral term helps, whether it survives being tested on unseen sites, and what falls out the far end.

Study-area map showing soil-moisture sampling locations and model predictor surfaces across the East River watershed.
The sampling network sits in the upper East River valley near Crested Butte, Colorado, a snowmelt-driven system where the wet ground and the walkable ground are not the same ground. Field sampling locations over a 10 m hillshade, with the East River and Slate River. Insets over the Pump House cluster show the two model predictor surfaces: annual peak NDVI from Sentinel-2, and the HAND-based topographic wetness index. Any counts printed on this figure are figure-level and differ from the modelled subset described in the method.
A gloved hand beside a small pit dug into stony meadow ground, with water standing in the bottom of the hole and excavated cobbles piled at its edge.
Soil core sampling in the field. Water stands in the open hole, and the material dug out of it is coarse and stony.
Concept diagram: plant roots connect subsurface soil moisture to the canopy, so satellite-observed greenness carries information about the moisture below.
Roots link what is below the surface to what a satellite can see above it. That is the whole premise in one picture. Conceptual diagram of "vegetation as a sensor." Illustrative, not a data product or a measured result.
A season of NDVI observations for one location with a smooth harmonic curve fitted through them.
A ragged, cloud-gapped year of satellite observations becomes a smooth seasonal curve you can take measurements off. Third-order Fourier harmonic regression fitted to a per-pixel Sentinel-2 NDVI time series, 2016–2025. Illustrative example for one location. Not a watershed-wide signal and not a moisture measurement.
Scatter plot comparing per-site R-squared from a topographic-wetness-only model against a model that also includes a spectral metric, with a one-to-one reference line.
Adding the best site-specific spectral metric raises adjusted per-site R² where the topographic baseline is weakest, and which metric helps is site-specific. Per-site adjusted R², wetness-only (x-axis) versus wetness plus the best site-specific spectral metric (y-axis), with a 1:1 line. This is within-site fit, not transfer to new sites, and per-site values are exploratory at small n.
Observed versus predicted soil moisture for sites held out of model fitting, plotted against a one-to-one line.
Every point here is a site the model had never seen when it made the prediction. Leave-one-site-out observed versus predicted volumetric water content, 2018, from the pooled mixed model; held-out sites are predicted from the fixed effects only. One validation year shown.
Median predicted soil moisture along elevation and slope gradients, comparing the full model with a topography-only baseline.
Spectral metrics correct topographic underestimation on hillslopes. The lift lands on steeper, higher terrain, exactly where terrain alone has least to say. Median predicted VWC with interquartile ribbons along elevation and slope, full model (Peak NDVI + topographic wetness + community) versus a topography-and-community null. Descriptive full-model versus topographic-baseline comparison, not uniform field validation, and the highest slope and elevation bins are small-n.

A companion strand: imaging the microsites directly

Running alongside the mapping work is a strand that attacks the same question from the opposite end of the scale. Using synchrotron µ-X-ray fluorescence at Stanford's SLAC, we imaged Fe(II)/Fe(III) redox chemistry inside intact soil cores (the millimetre-scale anoxic microsites themselves) and correlated that chemistry with lab methane-flux incubations on the same material.

This is scientific context for why moisture is the variable worth mapping. It is not validation of the maps and it is not the deliverable. It's here because it's the other half of how I think about the problem: I've collected the field and lab data that GIS layers are made of (flux chambers, soil cores, porewater chemistry), not just downloaded them.

Synchrotron micro-X-ray fluorescence images of an intact soil core showing the distribution of Fe(II) and Fe(III).
The microsites that make methane, imaged directly. This is the scale the satellite will never see. Synchrotron µ-XRF mapping of Fe(II)/Fe(III) redox chemistry in intact soil cores at Stanford's SLAC. Companion laboratory science and scientific context. It does not validate the mapping products.
A white flux chamber sitting on grass at the edge of a willow thicket, with a ribbed hose, a marker flag, a clipboard holding a data sheet and a pair of work gloves beside it. Conifers and a pointed peak stand in the background under a blue sky.
Methane flux gas sampling in the field. A survey chamber seated on meadow ground, with a paper data sheet beside it.

Honest limits & status

Honest limits

  • A first-author manuscript is in preparation: not submitted, not under review, not published. It is intended for MDPI Remote Sensing.
  • Prediction skill degrades above roughly 24% volumetric water content: the wettest ground is systematically under-predicted, a ceiling a linear model cannot clear with these predictors.
  • Calibrated on one watershed. Leave-one-site-out is cross-validation inside a single sampling network, not independent regional validation. Transfer beyond this basin is untested.

Outputs

Interested in this work?

I'm an early-career environmental data scientist in the SF Bay Area.

Get in touch See all work