Is the growing season shifting — are plants greening up earlier or staying green longer?

What you can answer

Whether greenup is starting earlier over the last N years (MCD12Q2 Greenup day-of-year trend, pixel-by-pixel or AOI-mean).
Whether senescence/dormancy is happening later (MCD12Q2 Senescence, Dormancy day-of-year trends → growing-season length = Dormancy − Greenup).
Whether peak greenness is shifting in time or magnitude (MCD12Q2 MidGreenup/Peak dates and EVI_Amplitude/EVI_Area per cycle).
How seasonal productivity is changing (MOD15A2H LAI/FPAR integrated over the season; MOD13Q1 NDVI/EVI seasonal integral as a greenness-productivity proxy).
Whether the trend is consistent with warming springs — you can correlate phenology shift against an independent temperature record you bring in.

What you can NOT answer with these datasets alone

Direct attribution to CO₂, warming, or management — phenology shift is observed, not explained; needs climate/management ancillary data.
Net carbon uptake (NEP/GPP) — LAI/FPAR and NDVI are greenness/structure proxies, not flux. Use a flux product (e.g. MOD17 GPP or eddy-covariance towers) for actual carbon.
Crop-specific phenology — MODIS 250–500m mixes fields; MCD12Q2 reports a generic land-surface phenology, not “corn vs soybean” emergence.
Sub-seasonal timing finer than the cadence — MOD13Q1 is 16-day, MOD15A2H is 8-day; the green-up date has an uncertainty floor near the composite period.
Whether two green-up cycles per year are crop double-cropping or noise — MCD12Q2 reports up to 2 cycles but does not label them.

Code template (Python, cloud-direct, ~30 lines)

import earthaccess
import xarray as xr
import numpy as np

earthaccess.login(strategy="netrc")

# 1. Define AOI + time window
aoi = (-96.0, 39.0, -87.0, 44.0)  # W, S, E, N — US Corn Belt (Iowa/Illinois)
years = (2001, 2023)

# 2. Pull MCD12Q2 annual phenology metrics (greenup / dormancy day-of-year)
pheno = earthaccess.search_data(
    short_name="MCD12Q2",
    bounding_box=aoi,
    temporal=(f"{years[0]}-01-01", f"{years[1]}-12-31"),
)
# ...open HDF-EOS, read SDS "Greenup" and "Dormancy" (band 0 = first cycle).
# Values are days since 1970-01-01; convert to day-of-year per year.
# season_length = Dormancy - Greenup

# 3. Pull MOD15A2H LAI/FPAR for seasonal productivity context
lai = earthaccess.search_data(
    short_name="MOD15A2H",
    bounding_box=aoi,
    temporal=(f"{years[0]}-01-01", f"{years[1]}-12-31"),
)
# ...open, apply QC (FparLai_QC), scale Lai_500m (scale 0.1), integrate over season

# 4. Optional continuity past MODIS with VIIRS LAI/FPAR
viirs_lai = earthaccess.search_data(
    short_name="VNP15A2H",
    bounding_box=aoi,
    temporal=("2012-01-01", f"{years[1]}-12-31"),
)

# 5. Fit a linear trend of AOI-mean Greenup DOY vs year → days/decade earlier
#    Repeat for Dormancy and season_length; plot all three.

Expected output

Time-series: AOI-mean green-up day-of-year per year, with a fitted trend line (days/decade earlier or later).
Time-series: AOI-mean growing-season length (Dormancy − Greenup), per year.
Map: per-pixel green-up trend (days/decade), showing where the shift is strongest.
Companion plot: seasonal-integrated LAI per year, to show whether an earlier/longer season also means more cumulative greenness.

Caveats

MCD12Q2 phenology dates are stored as days since 1970-01-01 (not day-of-year) — you must convert, and watch for the fill value (32767).
MCD12Q2 reports up to two cycles per pixel per year; for single-season temperate crops use cycle 0, but double-cropped or irrigated pixels may populate cycle 1.
LAI/FPAR saturate at high canopy density (LAI ≳ 5–6), so peak-summer productivity differences in dense corn may be compressed.
MODIS Terra (MOD*) has known sensor degradation late in the record; cross-check trends against the VIIRS continuity record before claiming a multi-decadal signal.
A 23-year record is short for climate attribution — report the trend with its confidence interval, not as a settled fact.

Cross-DAAC composition

This is a single-DAAC workflow: all four products (MCD12Q2, MOD15A2H, MOD13Q1, VNP15A2H) live at LP DAAC under one Earthdata Login. The complexity here is temporal stitching (MODIS → VIIRS continuity) and HDF-EOS subdataset handling, not multi-DAAC auth. See recipes/r01-three-daac-composition.mdx for the general cross-DAAC pattern when you add a climate record from another archive.