Code to estimate the trajectory of the barycenter of a global land-surface property - for now anything related to vegetation phenology.
If you use this repository, please cite the Green Wave paper (Miguel D. Mahecha, Guido Kraemer, Martin Reinhardt, David Montero, Fabian Gans, Ana Bastos, Hannes Feilhauer, Ida Flik, Chaonan Ji, Teja Kattenborn, Mirco Migliavacca, Milena Mönks, Johannes Quaas, Sebastian Sippel, Sophia Walther, Sebastian Wieneke, Christian Wirth and Gustau Camps-Valls, Accelerated North-East Shift of the Global Green Wave Trajectory, CITATION).
@article{Mahecha2025GreenWave,
title = {Accelerated North-East Shift of the Global Green Wave Trajectory},
author = {Mahecha, Miguel D. and Kraemer, Guido and Reinhardt, Martin and
Montero, David and Gans, Fabian and Bastos, Ana and Feilhauer, Hannes and
Flik, Ida and Ji, Chaonan and Kattenborn, Teja and Migliavacca, Mirco and
Mönks, Milena and Quaas, Johannes and Sippel, Sebastian and Walther, Sophia and
Wieneke, Sebastian and Wirth Christian and Camps-Valls, Gustau},
journal = {Proceedings of the National Academy of Sciences (PNAS)},
year = {2026},
doi = {10.1073/pnas.2515835123}
}A global “green wave” seasonally traverses Earth’s surface reflecting large-scale vegetation phenology, driven primarily by solar radiation and modulated by climate variability and ecosystem processes.
We derive a single global trajectory of the green wave by tracking the centroid of terrestrial greenness (or any other biophysical variable) over time. The methods provides a complementary perspective to conventional LSP products by focusing on the coherent large-scale displacement of the biosphere’s phenological center—a kind of speedometer and compass for the global ecosystem.
The resulting trajectory, expressed in interpretable units such as kilometers and days, offers a unified way to quantify global phenology and compare dynamics across Earth observation products and Earth System Models.
- Code to compute green-wave trajectory from analysis-ready data cubes containing relevant variables e.g. a vegetation index time series, its interpolation to daily values, and functions to compute the derived metrics such as mean seasonal trajectories, derivatives, and moments like the viridistice and equiviridis.
- Scripts to reproduce the analytic figures of the manuscript (not the conceptual illustrations)
- Scripts demonstrate how to compute statistics across data products or CMIP model ensembles.
- Computed trajectories for all data sets presented in the paper, dervied mean seasonal trajectories, and statistics
We cannot redistribute the original input datasets (e.g., GIMMS, GLASS, MODIS, Sentinel-2 or CMIP model outputs) because we are not the data providers. Downloading the data from the original data repositories should be done by the users. This repository assumes that an analysis-ready data cube containing the required vegetation index products is already available.
A full description of the data preprocessing workflow is given in the paper. To keep this repository concise, those details are not duplicated here.
Install Julia 1.11 from here, clone the GreenWaveTrajectory repository
git clone https://github.com/mafla/GreenWaveTrajectory.git
cd GreenWaveTrajectoryThis repository contains a Julia environment (via Project.toml and Manifest.toml) plus scripts you can run directly.
From the root of the repository:
julia> using Pkg
julia> Pkg.activate(".")This installs all packages listed in the environment:
julia> Pkg.instantiate()Speed up the first run of the scripts by precompiling all dependencies:
julia> Pkg.precompile()The data needs to be postprocessed locally, to compute the statistics and figures shown in the paper. You can run the postprocessing script:
cd GreenWave/run
julia postprocessing.jlTo compute the green wave trajectory from your own data, you can use the script GreenWave/runscript.jl. It is a template script that you can modify to point to your own data files and set parameters as needed.
cd GreenWave
julia runscript.jlThe figures in the paper can be reproduced by running the scripts in the respective figure folders. The scripts assume that the postprocessing step has been completed and the necessary data files are available in the postprocessing folder. Each figure comes with an own julia environment to ensure reproducibility. You need to instantiate the environment in each figure folder before running the scripts.
cd figures/Figure_02
julia makefigure.jlECEF: Earth-Centered, Earth-Fixed coordinates LLA: Latitude, Longitude, Altitude
Workflow:
- Create Mask for LandSurface and Latitude impact
- Get coords for each pixel in ECEF
- scale the coords to dataarray
flowchart TD
A[YAXArray]-- WeightedCenterTimeSeries -->C[DataFrame X, Y, Z];
C-- smooth -->D[DataFrame X, Y, Z, Xsmooth, Ysmooth, Zsmooth, Latsmooth, Lonsmooth]
D-- interpolate -->E[DataFrame X, Y, Z, dXdt, dYdt, dZdt, speed, Lat, Lon]
E-- getMSCinterpolated -->F[DataFrame <br> X <br> Y<br> Z<br> dXdt<br> dYdt<br> dZdt<br> speed<br> Lat<br> Lon]
E-- getStats --> G[DataFrame <br> Year<br> maxDateZ<br> maxValueZ<br> minDateZ<br> minValueZ<br> maxDateLat<br> maxValueLat<br> minDateLat<br> maxUpSpeedDate<br> maxUpSpeed<br> maxUpSpeedZ<br> maxUpSpeedLat<br> maxDownSpeedDate<br> maxDownSpeed<br> maxDownSpeedZ<br> maxDownSpeedLat]
F-- getStatsMSC --> H[DataFrame <br> maxDateZ<br> maxValueZ<br> minDateZ<br> minValueZ<br> maxDateLat<br> maxValueLat<br> minDateLat<br> maxUpSpeedDate<br> maxUpSpeed<br> maxUpSpeedZ<br> maxUpSpeedLat<br> maxDownSpeedDate<br> maxDownSpeed<br> maxDownSpeedZ<br> maxDownSpeedLat]
flowchart TD
A[X, Y, Z]-- SRDSmoothing -->B[XSmooth, YSmooth, ZSmooth]
B-- Spline interpolation -->C[X, Y, Z, dXdt, dYdt, dZdt, speed, lat, lon, geospeed]
B-- groupby doy, mean -->D[Xmsc, Ymsc, Zmsc]
D-- spline interpolation -->E[Xmsc, Ymsc, Zmsc, dXdtmsc, dYdtmsc, dZdtmsc, speedmsc, latmsc, lonmsc, geospeedmsc]
- add description about the
GREEN_WAVE_PATH - describe in detail how the repo is organized
- add where to set the data directories for the data products
- change logging to @info "$(now()): MESSAGE"