Model Validation Demo¶
This notebook demonstrates TerraFlow's model validation workflow, including:
- Spatial block cross-validation with buffer zones (Roberts et al. 2017)
- Cohen's kappa against a reference dataset
- Moran's I on score residuals
- Validation block in report.json
In [1]:
Copied!
import json
import tempfile
import textwrap
from pathlib import Path
import numpy as np
import pandas as pd
import rasterio
from rasterio.transform import from_origin
from terraflow.pipeline import run_pipeline
from terraflow.validation import run_validation
# ── Synthetic inputs ─────────────────────────────────────────────────────────
tmp = Path(tempfile.mkdtemp())
# Synthetic 5×5 raster in EPSG:4326
raster_path = tmp / "raster.tif"
arr = np.arange(25, dtype="float32").reshape(5, 5)
with rasterio.open(
raster_path, "w", driver="GTiff",
height=5, width=5, count=1, dtype="float32",
crs="EPSG:4326", transform=from_origin(-100.0, 40.0, 0.01, 0.01),
) as ds:
ds.write(arr, 1)
# Synthetic climate CSV
climate_path = tmp / "climate.csv"
pd.DataFrame({
"lat": [40.0, 40.01, 40.02],
"lon": [-100.0, -99.99, -99.98],
"mean_temp": [18.0, 19.0, 20.0],
"total_rain": [100.0, 120.0, 140.0],
}).to_csv(climate_path, index=False)
# Synthetic reference labels (lat/lon within ROI)
ref_path = tmp / "reference.csv"
pd.DataFrame({
"lat": [39.97, 39.98, 39.99, 40.00, 40.01],
"lon": [-100.03, -100.01, -99.99, -99.97, -99.95],
"label": ["low", "medium", "high", "high", "medium"],
}).to_csv(ref_path, index=False)
# Config with validation section
config_path = tmp / "config.yml"
config_path.write_text(textwrap.dedent(f"""\
raster_path: {raster_path}
climate_csv: {climate_path}
output_dir: {tmp}/outputs
roi:
type: bbox
xmin: -101.0
ymin: 39.0
xmax: -99.0
ymax: 41.0
model_params:
v_min: 0.0
v_max: 25.0
t_min: 0.0
t_max: 40.0
r_min: 0.0
r_max: 300.0
w_v: 0.4
w_t: 0.3
w_r: 0.3
max_cells: 50
validation:
n_blocks_side: 2
buffer_deg: 0.01
reference_csv: {ref_path}
"""))
# Run pipeline to produce features.parquet, then validate
df = run_pipeline(config_path)
print(f"Pipeline complete: {len(df)} cells, run_id={df['run_id'].iloc[0]}")
report_path = run_validation(config_path)
print(f"Updated report: {report_path}")
import json
import tempfile
import textwrap
from pathlib import Path
import numpy as np
import pandas as pd
import rasterio
from rasterio.transform import from_origin
from terraflow.pipeline import run_pipeline
from terraflow.validation import run_validation
# ── Synthetic inputs ─────────────────────────────────────────────────────────
tmp = Path(tempfile.mkdtemp())
# Synthetic 5×5 raster in EPSG:4326
raster_path = tmp / "raster.tif"
arr = np.arange(25, dtype="float32").reshape(5, 5)
with rasterio.open(
raster_path, "w", driver="GTiff",
height=5, width=5, count=1, dtype="float32",
crs="EPSG:4326", transform=from_origin(-100.0, 40.0, 0.01, 0.01),
) as ds:
ds.write(arr, 1)
# Synthetic climate CSV
climate_path = tmp / "climate.csv"
pd.DataFrame({
"lat": [40.0, 40.01, 40.02],
"lon": [-100.0, -99.99, -99.98],
"mean_temp": [18.0, 19.0, 20.0],
"total_rain": [100.0, 120.0, 140.0],
}).to_csv(climate_path, index=False)
# Synthetic reference labels (lat/lon within ROI)
ref_path = tmp / "reference.csv"
pd.DataFrame({
"lat": [39.97, 39.98, 39.99, 40.00, 40.01],
"lon": [-100.03, -100.01, -99.99, -99.97, -99.95],
"label": ["low", "medium", "high", "high", "medium"],
}).to_csv(ref_path, index=False)
# Config with validation section
config_path = tmp / "config.yml"
config_path.write_text(textwrap.dedent(f"""\
raster_path: {raster_path}
climate_csv: {climate_path}
output_dir: {tmp}/outputs
roi:
type: bbox
xmin: -101.0
ymin: 39.0
xmax: -99.0
ymax: 41.0
model_params:
v_min: 0.0
v_max: 25.0
t_min: 0.0
t_max: 40.0
r_min: 0.0
r_max: 300.0
w_v: 0.4
w_t: 0.3
w_r: 0.3
max_cells: 50
validation:
n_blocks_side: 2
buffer_deg: 0.01
reference_csv: {ref_path}
"""))
# Run pipeline to produce features.parquet, then validate
df = run_pipeline(config_path)
print(f"Pipeline complete: {len(df)} cells, run_id={df['run_id'].iloc[0]}")
report_path = run_validation(config_path)
print(f"Updated report: {report_path}")
INFO:terraflow:Loaded config from /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/config.yml
INFO:terraflow:Computed run fingerprint LE4QHd2w_bFCsKWgP4GjlhcPTeRZp7dEJs7PPOdRgMY (config=9c91f32d9396a3cf1a02349bc1fb4378b6f1c5fda2a24ee53e732cb390d234e3, inputs=2)
INFO:terraflow:DataCatalog built: raster=raster.tif (CRS EPSG:4326, shape (5, 5)), climate=climate.csv (3 rows, vars=['mean_temp', 'total_rain'])
INFO:terraflow:Loaded raster from /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/raster.tif
INFO:terraflow:Loaded climate CSV from /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/climate.csv with 3 rows
INFO:terraflow:Climate variables: ['mean_temp', 'total_rain']
INFO:terraflow:Climate CSV validated successfully: 3 valid records
INFO:terraflow:Loaded raster: /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/raster.tif (CRS: EPSG:4326)
INFO:terraflow:Loaded climate data: /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/climate.csv
INFO:terraflow:Clipped raster to ROI
INFO:terraflow.climate:ClimateInterpolator initialised: strategy='spatial', interpolation_method='linear', records=3, variables=['mean_temp', 'total_rain']
INFO:terraflow:Initialized climate interpolator with strategy='spatial', method='linear'
INFO:terraflow:Sampled 25 cells from 25 valid cells in ROI
INFO:terraflow.climate:Linear interpolation failed (QH6154 Qhull precision error: Initial simplex is flat (facet 2 is coplanar with the interior point)
While executing: | qhull d Qt Qc Qbb Qz Q12
Options selected for Qhull 2019.1.r 2019/06/21:
run-id 1123646479 delaunay Qtriangulate Qcoplanar-keep Qbbound-last
Qz-infinity-point Q12-allow-wide _pre-merge _zero-centrum Qinterior-keep
Pgood _max-width 0.02 Error-roundoff 1.4e-13 _one-merge 9.7e-13
Visible-distance 2.8e-13 U-max-coplanar 2.8e-13 Width-outside 5.5e-13
_wide-facet 1.7e-12 _maxoutside 1.1e-12
The input to qhull appears to be less than 3 dimensional, or a
computation has overflowed.
Qhull could not construct a clearly convex simplex from points:
- p3(v4): 40 -1e+02 1e+02
- p1(v3): 40 -1e+02 0.1
- p2(v2): 40 -1e+02 0
- p0(v1): 40 -1e+02 0.21
The center point is coplanar with a facet, or a vertex is coplanar
with a neighboring facet. The maximum round off error for
computing distances is 1.4e-13. The center point, facets and distances
to the center point are as follows:
center point 40.01 -99.99 25.08
facet p1 p2 p0 distance= -1.1e-08
facet p3 p2 p0 distance= -6.3e-16
facet p3 p1 p0 distance= 2.5e-15
facet p3 p1 p2 distance= -1.2e-14
These points either have a maximum or minimum x-coordinate, or
they maximize the determinant for k coordinates. Trial points
are first selected from points that maximize a coordinate.
The min and max coordinates for each dimension are:
0: 40 40.02 difference= 0.02
1: -100 -2.225e-308 difference= 100
2: 0 100 difference= 100
If the input should be full dimensional, you have several options that
may determine an initial simplex:
- use 'QJ' to joggle the input and make it full dimensional
- use 'QbB' to scale the points to the unit cube
- use 'QR0' to randomly rotate the input for different maximum points
- use 'Qs' to search all points for the initial simplex
- use 'En' to specify a maximum roundoff error less than 1.4e-13.
- trace execution with 'T3' to see the determinant for each point.
If the input is lower dimensional:
- use 'QJ' to joggle the input and make it full dimensional
- use 'Qbk:0Bk:0' to delete coordinate k from the input. You should
pick the coordinate with the least range. The hull will have the
correct topology.
- determine the flat containing the points, rotate the points
into a coordinate plane, and delete the other coordinates.
- add one or more points to make the input full dimensional.
), falling back to nearest-neighbor
INFO:terraflow.climate:Linear interpolation failed (QH6154 Qhull precision error: Initial simplex is flat (facet 2 is coplanar with the interior point)
While executing: | qhull d Qt Qc Qbb Qz Q12
Options selected for Qhull 2019.1.r 2019/06/21:
run-id 1123646479 delaunay Qtriangulate Qcoplanar-keep Qbbound-last
Qz-infinity-point Q12-allow-wide _pre-merge _zero-centrum Qinterior-keep
Pgood _max-width 0.02 Error-roundoff 1.4e-13 _one-merge 9.7e-13
Visible-distance 2.8e-13 U-max-coplanar 2.8e-13 Width-outside 5.5e-13
_wide-facet 1.7e-12 _maxoutside 1.1e-12
The input to qhull appears to be less than 3 dimensional, or a
computation has overflowed.
Qhull could not construct a clearly convex simplex from points:
- p3(v4): 40 -1e+02 1e+02
- p1(v3): 40 -1e+02 0.1
- p2(v2): 40 -1e+02 0
- p0(v1): 40 -1e+02 0.21
The center point is coplanar with a facet, or a vertex is coplanar
with a neighboring facet. The maximum round off error for
computing distances is 1.4e-13. The center point, facets and distances
to the center point are as follows:
center point 40.01 -99.99 25.08
facet p1 p2 p0 distance= -1.1e-08
facet p3 p2 p0 distance= -6.3e-16
facet p3 p1 p0 distance= 2.5e-15
facet p3 p1 p2 distance= -1.2e-14
These points either have a maximum or minimum x-coordinate, or
they maximize the determinant for k coordinates. Trial points
are first selected from points that maximize a coordinate.
The min and max coordinates for each dimension are:
0: 40 40.02 difference= 0.02
1: -100 -2.225e-308 difference= 100
2: 0 100 difference= 100
If the input should be full dimensional, you have several options that
may determine an initial simplex:
- use 'QJ' to joggle the input and make it full dimensional
- use 'QbB' to scale the points to the unit cube
- use 'QR0' to randomly rotate the input for different maximum points
- use 'Qs' to search all points for the initial simplex
- use 'En' to specify a maximum roundoff error less than 1.4e-13.
- trace execution with 'T3' to see the determinant for each point.
If the input is lower dimensional:
- use 'QJ' to joggle the input and make it full dimensional
- use 'Qbk:0Bk:0' to delete coordinate k from the input. You should
pick the coordinate with the least range. The hull will have the
correct topology.
- determine the flat containing the points, rotate the points
into a coordinate plane, and delete the other coordinates.
- add one or more points to make the input full dimensional.
), falling back to nearest-neighbor
INFO:terraflow:Interpolated climate for 25 cells using strategy='spatial'
INFO:terraflow:Closed raster dataset
INFO:terraflow:Artifacts written to /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/outputs/runs/LE4QHd2w_bFCsKWgP4GjlhcPTeRZp7dEJs7PPOdRgMY (fingerprint=LE4QHd2w_bFCsKWgP4GjlhcPTeRZp7dEJs7PPOdRgMY, cells=25, total=0.17s)
INFO:terraflow:Running validation on /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/outputs/runs/LE4QHd2w_bFCsKWgP4GjlhcPTeRZp7dEJs7PPOdRgMY
/Users/chandhini/akhil/TerraFlow/terraflow/validation.py:313: UserWarning: Reference points are up to 1.4 degrees distance from nearest cell — possible extent mismatch. Ensure reference CSV uses WGS84 degrees and covers the same region as the pipeline output. kappa = _compute_kappa(df, reference_df) INFO:terraflow:Cohen's kappa computed from 5 reference points: 0.0000
INFO:terraflow:Validation block written to /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/outputs/runs/LE4QHd2w_bFCsKWgP4GjlhcPTeRZp7dEJs7PPOdRgMY/report.json
Pipeline complete: 25 cells, run_id=LE4QHd2w_bFCsKWgP4GjlhcPTeRZp7dEJs7PPOdRgMY Updated report: /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/outputs/runs/LE4QHd2w_bFCsKWgP4GjlhcPTeRZp7dEJs7PPOdRgMY/report.json
In [2]:
Copied!
import json
with open(report_path) as f:
report = json.load(f)
val = report.get("validation", {})
for key, value in val.items():
print(f"{key}: {value}")
import json
with open(report_path) as f:
report = json.load(f)
val = report.get("validation", {})
for key, value in val.items():
print(f"{key}: {value}")
method: spatial_block_cv citation: Roberts et al. 2017, Ecography n_blocks_side: 2 buffer_deg: 0.01 n_folds: 4 mean_fold_accuracy: 0.8472222222222222 cohen_kappa: 0.0 morans_i_residuals: -0.03518103862225815 kriging_loocv_rmse: None reference_dataset: /var/folders/h_/_bjpnb894kgf18frynwf6bww0000gn/T/tmp51f29au5/reference.csv n_reference_points: 5 note: model has no free parameters; fold accuracy reflects spatial label consistency, not fit generalization
Interpreting Results¶
Mean fold accuracy: Measures how consistently the suitability labels hold across spatial blocks. Because TerraFlow's model has no free parameters (scores are deterministic given config), fold accuracy reflects spatial label distribution consistency rather than model generalization.
Cohen's kappa: Measures agreement between TerraFlow's classification and the reference dataset, corrected for chance. Values: <0 = less than chance, 0 = chance, 0.2-0.4 = fair, 0.4-0.6 = moderate, 0.6-0.8 = substantial, 0.8-1.0 = near-perfect.
Moran's I: Measures spatial autocorrelation in score residuals. Positive values indicate spatial clustering of similar residuals; values near 0 indicate spatial randomness.
Kriging LOOCV RMSE: Per-variable leave-one-out cross-validation error for kriging interpolation, computed during the pipeline run (not by validation).
In [3]:
Copied!
loocv = report.get("kriging_loocv")
if loocv:
print("Kriging LOOCV RMSE:")
for var, rmse in loocv.items():
print(f" {var}: {rmse:.6f}")
else:
print("No kriging LOOCV data (kriging was not used in this run)")
loocv = report.get("kriging_loocv")
if loocv:
print("Kriging LOOCV RMSE:")
for var, rmse in loocv.items():
print(f" {var}: {rmse:.6f}")
else:
print("No kriging LOOCV data (kriging was not used in this run)")
No kriging LOOCV data (kriging was not used in this run)
References¶
- Roberts, D.R. et al. (2017). Cross-validation strategies for data with temporal, spatial, hierarchical, or phylogenetic structure. Ecography, 40(8), 913-929.
- Cliff, A.D. & Ord, J.K. (1981). Spatial Processes: Models and Applications. Pion.