Chapter 29: Reproducibility and Experiment Provenance¶
Part VI — Data, Reproducibility & Verification
Learning Objectives¶
By the end of this chapter you will be able to:
- Explain what "provenance" means for a simulation result, beyond just saving the output
- Read every field of an
ExperimentRecordand know what each one is for - Find a
ModelExecutorby name in the registry it joins automatically - Validate an executor automatically with
ExecutorTestHarness, before ever running a real experiment - Use the words reproducibility, replicability and generalizability as the National Academies define them
- Decide when two results must be identical and when they only need to agree within a tolerance — and hash the right thing for each
- Record the environment of a run, and write the criteria that will judge a run before running it
# Standard imports — add chapter-specific imports below
import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
Chapter 22 built one complete ModelExecutor and ran it once. This chapter asks what happens after that one run — six months later, when a reviewer asks "are you sure this is the same input that produced Figure 3," or a collaborator asks "can I rerun exactly what you ran." Chapter 22's checksums were the first piece of the answer, and Chapter 28's run.json did the same for the data preparation; this chapter is the rest of it.
Why Provenance Matters¶
A result without provenance is a claim, not evidence. "I ran a fire model and got this map" is a claim; "I ran forest_fire version a3f9c1e, on input forest.gpkg with SHA-256 eba1aa64..., with initial_fire_density=0.1, and got output with SHA-256 3335784f..." is evidence — every piece of it independently checkable by someone who wasn't in the room when it ran.
Three failure modes provenance specifically guards against, none of them requiring anyone to have acted in bad faith:
- Silent input drift. A shared input file gets edited — a bug fix, a data correction — after an experiment ran against it. Without a checksum, nothing distinguishes "the exact file that produced this result" from "a file that looks similar now."
- Silent code drift. The model's own logic changes between when a result was produced and when someone tries to reproduce it. Without recording which version of the code ran, "rerun the experiment" silently means "run a different experiment that happens to share a name."
- Parameter amnesia. Six months on, even the person who ran the experiment often can't recall the exact parameter values used — which
seed, whichinitial_fire_density— from memory alone.
Reproducible, Replicable, General¶
The words are used loosely, so it is worth fixing them before going on. The National Academies of Sciences, Engineering, and Medicine (2019) define them for computational work:
- Reproducibility: obtaining consistent results using the same input data, computational steps, methods, code and conditions of analysis. This is what provenance serves, and what the rest of this chapter is about.
- Replicability: obtaining consistent results across studies that collect their own data to answer the same question.
- Generalizability: the extent to which results apply in other contexts or populations.
The same report accepts that identical bytes are not always attainable, and that consistent results can then be defined within an accepted range of variation — which is where tolerances come in, later in this chapter. Checking one implementation of a model against another (Chapter 30) is none of the three: it is a test of software, not a replication of a study.
The problem these words name is not hypothetical. In a survey of 1,576 researchers, more than 70% had tried and failed to reproduce another scientist's experiments, and more than half had failed to reproduce their own (Baker, 2016). In computing, a first-hand account describes the code of a former student, in folders with a mysterious numbering scheme and scarcely any documentation, as needing at least a year to run reliably (Barba, 2016). Almost every cause such accounts list — no code, no data, no record of the environment — is something a run can write down at the moment it happens.
ExperimentRecord in Depth¶
Chapter 22 touched three or four of ExperimentRecord's fields in passing. Here's the complete picture — every field exists to close one of the three gaps above:
from dissmodel.executor import ExperimentRecord, DataSource
print("ExperimentRecord fields:")
for name in ExperimentRecord.model_fields:
print(" -", name)
ExperimentRecord fields: - experiment_id - created_at - model_name - model_commit - code_version - resolved_spec - source - input_format - column_map - band_map - parameters - period - output_path - artifacts - metrics - status - logs
Grouped by what they're for:
| Field(s) | Closes which gap |
|---|---|
experiment_id, created_at |
identity — which run is this, and when did it happen |
model_name, model_commit, code_version |
code drift — which model, which exact version of its code |
resolved_spec |
the fully-resolved configuration actually used, defaults and all — not just what the caller explicitly passed |
source (a DataSource: type, uri, collection, version, checksum) |
input drift — exactly which file, and proof via checksum |
input_format, column_map, band_map |
how the raw input was interpreted — which column meant what |
parameters, period |
parameter amnesia — every value the model actually ran with |
output_path, artifacts, metrics |
the result itself, and proof it's the exact bytes save() produced — artifacts maps each file to its checksum, and output_sha256 is a shortcut to its "output" entry |
status, logs |
did it actually complete, and what happened along the way |
model_commit deserves a second look: pinning a commit hash, not just a version string like "1.2.0", means "run against exactly this code" survives even a model package that hasn't cut a formal release yet — the situation of every extension package of Chapter 21 before its first release, and of any model still installed from a branch.
Checksums as Proof, Not Just Metadata¶
A checksum only does its job if it's computed at the moment of use, not trusted from a filename or a timestamp. Chapter 22's load() already did this correctly — load_dataset() hashes the file's actual bytes as it reads them, not the string "forest.gpkg":
from dissmodel.io import load_dataset, save_dataset
from dissmodel.geo import vector_grid
gdf = vector_grid(dimension=(10, 10), resolution=1, attrs={"state": 0})
save_dataset(gdf, "checksum_demo.gpkg")
_, checksum_before = load_dataset("checksum_demo.gpkg")
print("Checksum, unchanged file:", checksum_before[:16], "...")
# Now change the file — flip a single cell's state and re-save
gdf.loc[0, "state"] = 1
save_dataset(gdf, "checksum_demo.gpkg")
_, checksum_after = load_dataset("checksum_demo.gpkg")
print("Checksum, one cell changed:", checksum_after[:16], "...")
print("Same file?", checksum_before == checksum_after)
Checksum, unchanged file: 4537ffb4b08603d1 ... Checksum, one cell changed: b778922a5be7d28b ... Same file? False
One value in one cell, and the checksum changes completely — that's the point. record.source.checksum, recorded once at load() time and never recomputed silently later, is what lets someone six months from now confirm "yes, this really is the file the experiment ran against," rather than a file that merely shares its name.
# A sample input for the executor examples below
sample_gdf = vector_grid(dimension=(15, 15), resolution=1, attrs={"state": 0})
save_dataset(sample_gdf, "forest.gpkg")
'5e1adfbd0270b15fe722d99ef660164fe12f42c3b663c1f25ba8f2ea3c9cd746'
Registering and Discovering Executors¶
A real project accumulates more than one ModelExecutor — the fire model from Chapter 22, whatever Chapter 26's LUCC models eventually add, an ensemble driver from Chapter 32. ExecutorRegistry gives every one of them a name, so calling code can look one up without importing its class directly. Registration is automatic: defining a subclass of ModelExecutor with a name adds it to the registry (__init_subclass__ does it). The explicit register() call below only makes the step visible; it is not needed:
from dissmodel.executor import ModelExecutor, ExecutorRegistry
class ForestFireExecutor(ModelExecutor):
name = "forest_fire"
def load(self, record):
gdf, checksum = load_dataset(record.source.uri)
record.source.checksum = checksum
return gdf
def run(self, data, record):
from dissmodel.core import Environment
from dissmodel_ca.models import FireModel
env = Environment(start_time=0, end_time=record.parameters.get("end_time", 20))
model = FireModel(
gdf=data,
initial_fire_density=record.parameters.get("initial_fire_density", 0.05),
seed=record.parameters.get("seed", 42),
)
model.initialize() # seed the starting fires (see Chapter 22)
env.run()
return data
def save(self, result, record):
uri = record.output_path or "output.gpkg"
record.output_path = uri
record.output_sha256 = save_dataset(result, uri)
record.status = "completed"
return record
ExecutorRegistry.register(ForestFireExecutor)
print("Registered:", ExecutorRegistry.list())
retrieved_cls = ExecutorRegistry.get("forest_fire")
print("Looked up by name:", retrieved_cls)
Registered: ['forest_fire'] Looked up by name: <class '__main__.ForestFireExecutor'>
The platform of Chapter 33 relies on this lookup: a job names its executor as a string (the class key of a registry file), and the worker finds the class here, without importing every executor class directly.
Testing Executors¶
An executor is code, and code that's never been tested is code you don't actually know works. ExecutorTestHarness runs two distinct checks, deliberately kept separate:
from dissmodel.executor.testing import ExecutorTestHarness
harness = ExecutorTestHarness(ForestFireExecutor)
passed = harness.run_contract_tests()
print("\nContract tests passed:", passed)
ExecutorTestHarness — ForestFireExecutor ──────────────────────────────────────────────────── ✅ name attribute exists ✅ name is a non-empty string ✅ name has no whitespace ✅ load() is implemented ✅ run() is implemented ✅ save() is implemented ✅ run() signature is correct ✅ save() signature is correct ──────────────────────────────────────────────────── All 8 checks passed ✅ Contract tests passed: True
Contract tests check the executor's shape, not its behavior: does it declare a name, does it implement load/run/save, do the method signatures match what ModelExecutor expects. None of this runs a single tick of simulation — it's the equivalent of Chapter 5's pytest checking that a function exists and is callable, before ever checking what it computes.
A sample-data cycle goes further: it actually runs load → run → save end to end, against a real ExperimentRecord, and confirms the result looks like a completed experiment should — status == "completed", an output_sha256 actually set:
from dissmodel.executor import DataSource
record = ExperimentRecord(
experiment_id="harness_demo",
model_name="forest_fire",
source=DataSource(type="local", uri="forest.gpkg"),
parameters={"initial_fire_density": 0.1, "seed": 7, "end_time": 5},
)
cycle_passed = harness.run_with_sample_data(record)
print("\nSample-data cycle passed:", cycle_passed)
▶ Running forest_fire... validate()... load()... run()... Running from 0 to 5 (duration: 5)
save()... ✅ Cycle OK — status=completed sha256=d4b68bbebc76... Sample-data cycle passed: True
Neither check looks at the science itself. An earlier draft of this executor forgot model.initialize(), so it simulated a forest with no fire at all, and it still passed both checks with a valid checksum. A contract test proves the executor has the right shape; only a test on the result (for example, assert (result["state"] != 0).any()) proves the model actually did something.
Both checks belong in a project's own test suite — the same pytest-based suite Chapter 5 first introduced — run automatically before any executor is trusted with a real experiment, not eyeballed once by hand and forgotten.
Hash for Identity, Tolerance for Equivalence¶
A checksum answers one question exactly: is this the same file? It answers nothing else. It does not say whether the content is right, and — less obviously — two files can hold the same values and still have different checksums, because the bytes of a file also depend on how it was written. The cell below writes one array of numbers three times as a GeoTIFF, with three compression settings:
import hashlib
import rasterio
from rasterio.transform import from_origin
values = np.random.default_rng(0).random((200, 200)).astype("float32")
def write(path, **options):
with rasterio.open(path, "w", driver="GTiff", height=200, width=200, count=1, dtype="float32",
crs="EPSG:31983", transform=from_origin(0, 2000, 10, 10), **options) as dst:
dst.write(values, 1)
for path, options in [("same_a.tif", {"compress": "deflate", "zlevel": 6}),
("same_b.tif", {"compress": "deflate", "zlevel": 1}),
("same_c.tif", {"compress": "lzw"})]:
write(path, **options)
with rasterio.open(path) as src:
decoded = src.read(1)
file_hash = hashlib.sha256(open(path, "rb").read()).hexdigest()[:12]
value_hash = hashlib.sha256(decoded.tobytes()).hexdigest()[:12]
print(f"{path} file {file_hash}… values {value_hash}…")
same_a.tif file af8d8e8f432c… values 168698eb29e0… same_b.tif file 704f98e19e8a… values 168698eb29e0… same_c.tif file 88ab5aad0cab… values 168698eb29e0…
Three file hashes, one value hash. This is not a contrived case. When the Open LUCC workflow of Chapter 31 was run again on a Mac with an Apple processor, its output GeoTIFF had a different SHA-256 from the one recorded on Linux — and every value in it was equal. The compressed bytes differed in every strip of the file, because each platform's GDAL encodes DEFLATE slightly differently. A check that required identical file hashes everywhere would have rejected a correct run. For identity across platforms, hash the decoded values, not the file.
The second case is the opposite: values that are not identical, from software that is. Floating-point addition is not associative, so the same numbers added in another order can differ in the last digit:
a = (0.1 + 0.2) + 0.3
b = 0.1 + (0.2 + 0.3)
print(repr(a), repr(b), a == b, abs(a - b))
0.6000000000000001 0.6 False 1.1102230246251565e-16
Any change that reorders a sum — another library version, another number of threads, another CPU — can move a result by that much, and a threshold in the model can then turn a last-digit difference into a different decision (Chapter 27's flood rule is an example). TerraME itself shows it: run eight times on the same machine, in the same container, LuccME's first test lab wrote outputs with several different hashes, all within 10⁻¹² of each other — the last digit the file prints (Chapter 30).
So there are two kinds of check, and a record should say which one it uses:
| Question | Check | Example |
|---|---|---|
| Is this the archived file? | Hash (SHA-256 of the file) | a reference output downloaded from a release |
| Is this the same data, on another platform? | Hash of the decoded values | two GeoTIFFs written by different GDAL builds |
| Is this the same result, computed again? | Tolerance (largest absolute difference ≤ ε) | a model run on another machine, a regenerated reference |
A tolerance has to be stated, and it has to be justified: the precision of the file (10⁻¹² for a CSV with twelve decimals), the precision of the arithmetic (about 10⁻⁷ relative for 32-bit floats), or a difference that was measured and explained. A tolerance chosen after seeing the differences is a regression guard — it stops things getting worse — not evidence of equivalence.
Recording the Environment¶
The same code on the same input can still give a different result in a different environment, so a record has to say where it ran. Three things are enough to start: the interpreter, the machine, and the version of every package.
import importlib.metadata as metadata
import json
import platform
import sys
environment = {
"python": sys.version.split()[0],
"platform": platform.platform(),
"machine": platform.machine(),
"packages": {name: metadata.version(name)
for name in ["dissmodel", "disscube", "disslucc", "numpy", "rasterio", "geopandas"]},
}
print(json.dumps(environment, indent=2))
{
"python": "3.13.16",
"platform": "Linux-6.18.44-fc-v139-x86_64-with-glibc2.39",
"machine": "x86_64",
"packages": {
"dissmodel": "0.6.6",
"disscube": "0.6.0",
"disslucc": "0.5.0",
"numpy": "2.5.4",
"rasterio": "1.5.2",
"geopandas": "1.2.0"
}
}
In a real project this is two files written by the run itself rather than a dictionary in a notebook. pip freeze > requirements.lock pins the version of every installed package, not only the ones the project imports directly — the Open LUCC workflow archives its lock file next to the record of the reported run. A short description of the machine (processor, memory, operating system) goes next to it, because time and memory measurements mean nothing without it. One caveat that the Mac rerun made concrete: a lock file pins Python packages, not the compiled libraries inside them (GDAL, GEOS, PROJ), which is where the platform differences came from.
Criteria Before the Run¶
The last piece of provenance is not about the past at all. Before a run whose results will be reported, write down what will count as success — in a file, committed, with the commit recorded. Afterwards, a script reads the file and the run's metrics and says pass or fail. The Open LUCC workflow keeps its criteria in acceptance.toml, frozen at a named commit before the reported runs:
[meta]
status = "frozen"
frozen_commit = "0eef81fddfdbf75b2a0d6db33559831b111fbb85"
[demand] # the allocation meets the demand
max_abs_error_km2 = 1000.0
max_iterations_per_step = 3000
[determinism] # a repeated run gives the same result
max_abs_diff = 0.0
[resolution] # class areas agree across the two grids
max_class_area_diff_fraction = 0.005
Writing criteria first prevents the most common way of fooling oneself — moving the threshold to where the result landed. It does not make the criteria meaningful, and each one should be read for what it leaves open. The cell below applies a demand criterion like the one above to a toy run in which the allocation never moves at all:
import tomllib
criteria = tomllib.loads("""
[demand]
max_abs_error_km2 = 1000.0
""")
# A toy run: the demand of one class grows by 30 km² a year, and the allocation never moves.
demand = np.array([30.0 * year for year in range(1, 15)])
allocated = np.zeros_like(demand)
largest_error = float(np.max(np.abs(allocated - demand)))
ok = largest_error <= criteria["demand"]["max_abs_error_km2"]
print(f"largest error {largest_error:.0f} km² -> {'PASS' if ok else 'FAIL'}")
largest error 420 km² -> PASS
It passes. The criterion is the allocation's own convergence tolerance, so any class whose total change is smaller than the tolerance passes without moving. This is not an invented weakness: the small synthetic case of the Open LUCC workflow (Chapter 31) passes its demand criterion with zero iterations in every year, its largest error growing year after year up to 429 km², under the 1,000 km² tolerance. The criterion shows that the iteration converged within its limit; it does not show that the model did anything. A criterion that is met is a statement about the criterion — reading what it cannot see is part of writing it.
A last story from the ecosystem makes the same point about tests in general. DisSCube 0.5.0 read PRODES maps with a reader that silently dropped the yearly deforestation classes from 2010 onwards. Every synthetic test passed. The bug was found by looking at a map of the BR-163 corridor, where the deforested share grew implausibly little — 18.1% to 20.5% between 2008 and 2024 — and fixed in 0.6.0, which gives 14.5% to 32.5% for the same cells. Records, hashes and criteria tell you what ran; they do not replace looking at the result.
Exercises¶
- Break it on purpose. Remove the
save()method from a copy ofForestFireExecutorand rerunrun_contract_tests(). Which specific check fails, and does the failure message tell you what's missing without you having to guess? - Trace a field to its gap. For
model_commitspecifically, write one sentence describing a real scenario where having it recorded — versus not — would have mattered to someone trying to reproduce a result. - Hash the values. Write a function
value_hash(path)that returns the SHA-256 of every band of a GeoTIFF, decoded, together with its shape and data type. Why should the shape and the type be part of the hash? - Choose the check. For each case, say whether you would compare by file hash, by value hash, or by tolerance, and give the tolerance if one is needed: (a) a reference CSV downloaded from a Zenodo release; (b) the output raster of your model, run again on a colleague's laptop; (c) a regenerated TerraME output that prints twelve decimals.
- A criterion that sees movement. Rewrite the toy demand criterion so that the run of Criteria Before the Run fails. What does your new criterion require, and could a correct model ever fail it?
# Your code here
Summary¶
Key concepts introduced¶
- Provenance as evidence, not a claim: three specific failure modes (input drift, code drift, parameter amnesia) it guards against
- Reproducibility, replicability and generalizability as the National Academies define them, and why checking one implementation against another is none of the three
- Every
ExperimentRecordfield, grouped by which failure mode it closes — includingmodel_commitpinning an exact commit rather than a version string, andartifactsholding the checksum of every output - Checksums computed from actual file bytes at the moment of use, proven to change the instant a single value does
- Hash for identity, tolerance for equivalence: the same values can have different file hashes, and the same software can produce values that differ in the last digit; hash the decoded values across platforms, and state every tolerance and where it comes from
- Recording the environment: interpreter, machine, and a lock file of every package — and what a lock file does not pin
- Acceptance criteria written before the run, and read for what they cannot see
ExecutorRegistryandExecutorTestHarness's two-layer check — and why neither, nor any record, replaces looking at the result
Chapter 30 asks the question a provenance record alone cannot answer: whether the Python components compute what the original TerraME and LuccME software computes, and how to state precisely where they differ.
Further Reading¶
- National Academies of Sciences, Engineering, and Medicine (2019). Reproducibility and Replicability in Science. The National Academies Press. https://doi.org/10.17226/25303
- Baker, M. (2016). 1,500 scientists lift the lid on reproducibility. Nature, 533(7604), 452–454. https://doi.org/10.1038/533452a
- Barba, L. A. (2016). The hard road to reproducibility. Science, 354(6308), 142. https://doi.org/10.1126/science.354.6308.142
- Gundersen, O. E. (2020). The reproducibility crisis is real. AI Magazine, 41(3), 103–106. https://doi.org/10.1609/aimag.v41i3.5318
- Sandve, G. K. et al. (2013). Ten Simple Rules for Reproducible Computational Research. PLOS Computational Biology, 9(10), e1003285
- Wilkinson, M. D. et al. (2016). The FAIR Guiding Principles for scientific data management and stewardship. Scientific Data, 3, 160018. https://doi.org/10.1038/sdata.2016.18
- DisSModel executor API reference: https://dissmodel.github.io/dissmodel/