NPP-HGSTFormer

NPP-HGSTFormer (Net Primary Productivity Heterogeneous Graph Spatiotemporal Transformer) reconstructs annual NPP at 30 m spatial resolution from 500 m observations. The model combines multi-source remote sensing and meteorological data in a heterogeneous graph, learns spatial and cross-resolution relationships with a Heterogeneous Graph Transformer (HGT), and models interannual dynamics with a temporal Transformer.

The current implementation was developed for Wangcheng District, Changsha, China, using annual data from 2016 to 2024. The framework is designed for research on fine-scale ecosystem productivity, carbon-cycle monitoring, and geospatial super-resolution.

Highlights

• Reconstructs 30 m NPP from 500 m MODIS NPP observations.
• Integrates 30 m NDVI, 500 m land use/land cover (LULC), and ERA5 meteorological variables.
• Represents multi-source, multi-resolution data as annual heterogeneous graphs.
• Uses HGT layers to learn spatial, cross-scale, and land-cover relationships.
• Uses a temporal Transformer to learn dependencies among annual graph representations.
• Trains without 30 m NPP labels through an aggregation-consistency loss.

Model overview

flowchart LR
    A["NDVI time series<br/>30 m"] --> D["Annual heterogeneous graph"]
    B["NPP + ERA5<br/>500 m"] --> D
    C["LULC<br/>500 m"] --> D
    D --> E["HGT spatial encoder"]
    E --> F["Annual NDVI-node embeddings"]
    F --> G["Temporal Transformer"]
    G --> H["MLP prediction head"]
    H --> I["Predicted NPP<br/>30 m"]
    I --> J["Mean aggregation<br/>to 500 m"]
    B --> K["Observed NPP<br/>500 m"]
    J --> L["Aggregation-consistency loss"]
    K --> L

Loading

For each year, the code constructs three node types and three active edge relations:

Node type	Resolution	Feature dimension	Features
ndvi	30 m	4	Temporal mean, standard deviation, maximum, and minimum of NDVI
npp	500 m	5	NPP, 2 m air temperature (t2m), soil moisture (swvl1), surface solar radiation (ssrd), and total precipitation (tp)
lulc	500 m	18	One-hot-encoded IGBP land-cover class

Edge relation	Mapping	Purpose
("ndvi", "spatial", "ndvi")	Four-neighbor pixel connections	Preserve local spatial continuity at 30 m
("ndvi", "pixel_to", "npp")	Many fine pixels to one coarse cell	Connect 30 m predictors to 500 m supervision
("lulc", "cover_to", "npp")	One coarse LULC cell to one NPP cell	Add land-cover constraints

ERA5 variables are concatenated with the NPP-node features in the active graph-building workflow. The repository also contains standalone ERA5-node and ERA5-to-NPP edge utilities, but they are not called by the current graph builder.

Repository structure

NPP-HGSTFormer/
|-- readme.md
`-- SuperResolution/
    `-- model/
        |-- nodes/
        |   |-- npp_nodes.py
        |   |-- ndvi_nodes.py
        |   |-- lulc_nodes.py
        |   `-- era5_nodes.py
        |-- edges/
        |   |-- ndvi_ndvi_edges.py
        |   |-- npp_ndvi_edges.py
        |   |-- npp_lulc_edges.py
        |   `-- npp_era5_edges.py
        |-- pyg/
        |   `-- pyg_build.py
        |-- models/
        |   `-- model.py
        `-- train/
            `-- trainer.py

Data requirements

The data are not distributed with this repository. Prepare annual files for 2016-2024 using the following directory structure and filenames:

data/
|-- NPP/
|   `-- MOD17A3HGF.061_Npp_500m_doy{year}001000000_aid0001_sf.tif
|-- NDVI/
|   `-- HLSS30.020_30m_aid0001_49N_{year}_NDVI.nc
|-- lulc/
|   `-- MCD12Q1.061_LC_Type1_doy{year}001000000_aid0001.tif
`-- ERA5/
    `-- ERA5_{year}.nc

Expected variables and array shapes:

• The NDVI NetCDF file must contain a variable named NDVI with shape [time, height, width]. The graph builder converts it to four annual statistical features.
• The ERA5 NetCDF file must contain t2m, swvl1, ssrd, and tp as two-dimensional annual grids.
• The NPP and LULC rasters must have identical grid dimensions and aligned cell ordering.
• The NDVI height and width must be integer multiples of the corresponding NPP dimensions.
• All annual files must use consistent spatial extents, node ordering, and array shapes.

The source products used in the associated study are:

• MODIS MOD17A3HGF annual NPP, 500 m.
• Harmonized Landsat and Sentinel-2 S30 (HLSS30), used to derive 30 m NDVI.
• MODIS MCD12Q1 annual land cover, 500 m.
• ERA5 meteorological variables.

NASA products can be located through NASA Earthdata Search, and ERA5 data are available from the Copernicus Climate Data Store.

Installation

Create and activate a virtual environment first.

python -m venv .venv

On Windows:

.venv\Scripts\Activate.ps1

On macOS or Linux:

source .venv/bin/activate

Install PyTorch using the command appropriate for your operating system and CUDA environment from the official PyTorch installation page. Then install PyTorch Geometric and the remaining dependencies:

python -m pip install torch_geometric
python -m pip install numpy scipy rasterio xarray netCDF4

If PyTorch Geometric requires platform-specific wheels in your environment, follow its official installation guide.

Configuration

The current research code uses explicit filesystem paths. Before running it, update:

1. DATA_DIR and OUT_DIR in SuperResolution/model/pyg/pyg_build.py.
2. GRAPH_DIR in SuperResolution/model/train/trainer.py.

The default temporal split is:

• Training: 2016-2022
• Testing: 2023-2024

Usage

Run commands from the SuperResolution directory so that the model namespace can be resolved.

1. Build annual heterogeneous graphs

cd SuperResolution
python -m model.pyg.pyg_build

This writes one serialized PyTorch Geometric graph per year (2016.pt through 2024.pt) to the configured OUT_DIR.

2. Train and evaluate the model

python -m model.train.trainer

The trainer automatically uses a CUDA device when available. It saves the lowest-training-loss checkpoint as best_model.pth in the current working directory and prints coarse-scale RMSE, MAE, and R2 metrics for the training and testing sequences.

Default implementation settings

Component	Setting
Hidden dimension	64
HGT layers	2
HGT attention heads	4
Temporal Transformer layers	2
Temporal attention heads	4
Transformer dropout	0.1
Optimizer	Adam
Learning rate	1e-4
Epochs	50
Objective	Mean-squared aggregation-consistency loss

The prediction head uses the final temporal embedding of each NDVI node to produce one 30 m NPP value. During training and evaluation, the predictions are averaged within their associated 500 m NPP cells and compared with the observed coarse-resolution NPP.

Important implementation notes

• This repository is a research implementation and currently has no command-line configuration layer.
• The training script saves the best checkpoint but does not implement early stopping.
• Evaluation metrics are calculated after aggregating the fine-resolution prediction to the coarse NPP grid.
• The current evaluate function evaluates the final graph in the supplied temporal sequence.
• The graph builder uses explicit Python loops for edge generation; large study areas may require substantial memory and preprocessing time.
• Reproducing manuscript results requires the same preprocessing, spatial alignment, masks, and input data used in the study.

Citation

If you use this repository, please cite the associated manuscript:

Wang, L., Guo, Y., Ren, H., and Yin, L. A Dual-Transformer Network for Spatiotemporal Super-Resolution Reconstruction of Ecosystem Net Primary Production (NPP) Based on Heterogeneous Graph: A Case Study of the Wangcheng District, Changsha, Central China.

Publication details and a DOI will be added when available.

Acknowledgements

This work was supported by the Hunan Provincial Natural Resource Science and Technology Planning Program (Grant No. HBZ20240144), the Open Project of the Technology Innovation Center for Natural Ecosystem Carbon Sink, Ministry of Natural Resources (Grant No. CS2023D04), and the Hunan Provincial Natural Science Foundation (Grant No. 2023JJ60188).