How to write a metadata record for a dataset in Python
On this page
Problem statement
Writing metadata by hand fails for a reason that has nothing to do with discipline: most of the fields are facts about the data, so a human copying them in is a human introducing errors. Extents drift after a clip, feature counts go stale after a filter, and the CRS in the record is the one the dataset had two versions ago.
The fix is to split the record in two. The derivable half โ extent, CRS, schema, counts, checksums, geometry types โ is computed at write time from the file being written. The authored half โ title, what one feature is, lineage, licence, contact โ is a small YAML file a person maintains, and a gate refuses to publish without it.
This guide builds both halves and the gate, for a real dataset: the points of interest that OpenStreetMap contributors have mapped as areas in Rutland โ schools, car parks, pubs, playgrounds โ from Geofabrik's extract of 1 January 2026.
Quick answer
import geopandas as gpd, hashlib, pathlib, datetime, json, yaml
def write_metadata(data_path, authored_path, out_path=None, layer=None):
gdf = gpd.read_file(data_path, layer=layer)
authored = yaml.safe_load(pathlib.Path(authored_path).read_text())
REQUIRED = ["title", "feature_definition", "content_date", "lineage", "licence", "contact"]
missing = [k for k in REQUIRED if not authored.get(k)]
if missing:
raise ValueError(f"{authored_path}: missing {missing}")
b = gdf.total_bounds
record = {
**authored,
"derived": {
"file": pathlib.Path(data_path).name,
"layer": layer,
"features": int(len(gdf)),
"geometry_types": sorted(gdf.geom_type.dropna().unique().tolist()),
"crs": gdf.crs.to_string() if gdf.crs is not None else None,
"bbox": [round(float(v), 6) for v in b],
"schema": {c: str(t) for c, t in gdf.dtypes.items() if c != gdf.geometry.name},
"null_counts": {c: int(gdf[c].isna().sum())
for c in gdf.columns if c != gdf.geometry.name},
"generated": datetime.datetime.now(datetime.timezone.utc).isoformat(timespec="seconds"),
},
}
out = pathlib.Path(out_path or (str(data_path) + ".meta.json"))
out.write_text(json.dumps(record, indent=2, default=str))
return record
record = write_metadata("data/pois.gpkg", "data/pois.meta.yml")
print({k: v for k, v in record["derived"].items() if k not in ("schema", "null_counts", "generated")})
print("null counts:", record["derived"]["null_counts"])
{'file': 'pois.gpkg', 'layer': None, 'features': 538, 'geometry_types': ['MultiPolygon'], 'crs': 'EPSG:4326', 'bbox': [-0.797341, 52.533354, -0.431659, 52.752769]}
null counts: {'osm_id': 0, 'code': 0, 'fclass': 0, 'name': 345}
The derived block is regenerated on every write, so it cannot disagree with the data. The authored block is reviewed in a pull request like code.
Step-by-step solution
1. Write the authored half first, as YAML next to the data
%%writefile data/pois.meta.yml
title: Points of interest mapped as areas, Rutland
feature_definition: >
One polygon per OpenStreetMap object that Geofabrik exports to its
pois_a layer (amenities, shops, leisure and tourism features drawn as
outlines rather than points), within the Rutland extract.
content_date: 2026-01-01
lineage:
- OpenStreetMap, Geofabrik extract rutland-260101-free.shp.zip
- layer gis_osm_pois_a_free_1, read with GeoPandas and written unchanged to GeoPackage
licence: ODbL-1.0
attribution: ยฉ OpenStreetMap contributors
contact: gis@example.org
accuracy:
statement: Not assessed. OpenStreetMap publishes no positional accuracy for individual features.
YAML rather than JSON because a human edits it, and multi-line strings are readable. The accuracy field says what is known โ nothing, for OSM โ rather than leaving a blank or inventing a figure; a statement that accuracy was not assessed is itself metadata.
2. Derive the rest at write time
Compute the derived block from the file you just wrote, not from the in-memory frame โ that is what catches a writer that dropped a column or promoted a geometry type:
print("in memory: ", pois.geom_type.value_counts().to_dict())
for path in ("data/pois.gpkg", "data/pois.shp"):
back = gpd.read_file(path)
print(f"{path:17}", back.geom_type.value_counts().to_dict(), "| columns:", [c for c in back.columns if c != "geometry"])
in memory: {'Polygon': 536, 'MultiPolygon': 2}
data/pois.gpkg {'MultiPolygon': 538} | columns: ['osm_id', 'code', 'fclass', 'name']
data/pois.shp {'Polygon': 536, 'MultiPolygon': 2} | columns: ['osm_id', 'code', 'fclass', 'name']
The frame holds 536 Polygons and 2 MultiPolygons. The GeoPackage holds 538 MultiPolygons โ GDAL declared the layer MultiPolygon because the input mixed the two, and promoted every Polygon to fit. A record derived from the frame would describe geometry nobody received. (A layer of Polygons alone is written as Polygons; it is the mix that triggers the promotion.)
3. Include the schema and the null counts
The schema is what breaks consumers when it changes, and the null counts tell a user whether a column is worth reading. Both are one line to compute and are the fields people thank you for.
4. Add a content checksum, not a byte checksum
A byte hash of a GeoPackage changes every time you write it. Hash the canonicalised data instead โ see How to checksum spatial datasets so you can prove they match.
5. Embed the record where the format allows
A sidecar is easy to separate from the data. Put at least the title, licence, version and content date inside the file as well, so a copy on somebody's laptop still knows what it is. How to store metadata inside a GeoPackage and How to embed provenance tags in a GeoTIFF cover the two formats that make this easy.
6. Emit the standard shapes from the same record
Do not author STAC or ISO directly; generate them. Metadata standards compared has the mapping.
7. Gate the publish step
The check is three lines and it is the only thing that keeps records complete once the deadline is close.
Code examples
The data used on this page
Points of interest mapped as areas, from Geofabrik's OpenStreetMap extract for Rutland of 1 January 2026 (ยฉ OpenStreetMap contributors, ODbL), and the Copernicus GLO-30 DEM over Rutland for the raster example (ยฉ DLR e.V. 2010โ2014 and ยฉ Airbus Defence and Space GmbH 2014โ2018, provided under COPERNICUS by the European Union and ESA).
import datetime
import hashlib
import json
import pathlib
from pathlib import Path
import geopandas as gpd
import pandas as pd
import rasterio
import requests
import yaml
from rasterio.windows import from_bounds
Path("data").mkdir(exist_ok=True)
if not Path("data/rutland.shp.zip").exists():
Path("data/rutland.shp.zip").write_bytes(requests.get(
"https://download.geofabrik.de/europe/united-kingdom/england/rutland-260101-free.shp.zip", timeout=300).content)
pois = gpd.read_file("/vsizip/data/rutland.shp.zip/gis_osm_pois_a_free_1.shp")
pois.to_file("data/pois.gpkg")
pois.to_file("data/pois.shp")
if not Path("data/dem.tif").exists():
url = ("/vsicurl/https://copernicus-dem-30m.s3.amazonaws.com/Copernicus_DSM_COG_10_N52_00_W001_00_DEM/"
"Copernicus_DSM_COG_10_N52_00_W001_00_DEM.tif")
with rasterio.open(url) as src:
window = from_bounds(-0.83, 52.52, -0.42, 52.76, src.transform)
profile = src.profile | {"driver": "GTiff", "width": int(window.width), "height": int(window.height),
"transform": src.window_transform(window)}
with rasterio.open("data/dem.tif", "w", **profile) as dst:
dst.write(src.read(1, window=window), 1)
print(f"{len(pois)} points of interest mapped as areas; files:", sorted(p.name for p in Path("data").iterdir()))
538 points of interest mapped as areas; files: ['dem.tif', 'pois.cpg', 'pois.dbf', 'pois.gpkg', 'pois.meta.yml', 'pois.prj', 'pois.shp', 'pois.shx', 'rutland.shp.zip']
Example 1 โ derive a raster's block too
import numpy as np
import rasterio
def derived_raster(path):
with rasterio.open(path) as src:
stats = []
for i in range(1, src.count + 1):
band = src.read(i, masked=True)
stats.append({
"band": i,
"dtype": src.dtypes[i - 1],
"nodata": src.nodatavals[i - 1],
"valid_fraction": float(band.count() / band.size),
"min": None if band.count() == 0 else round(float(band.min()), 2),
"max": None if band.count() == 0 else round(float(band.max()), 2),
"mean": None if band.count() == 0 else round(float(band.mean()), 2),
})
return {
"size": [src.width, src.height],
"bands": src.count,
"crs": src.crs.to_string() if src.crs else None,
"pixel_size": [abs(src.transform.a), abs(src.transform.e)],
"bbox": [round(v, 6) for v in src.bounds],
"band_stats": stats,
"tags": src.tags(),
}
print(json.dumps(derived_raster("data/dem.tif"), indent=1))
{
"size": [
984,
863
],
"bands": 1,
"crs": "EPSG:4326",
"pixel_size": [
0.00041666666666666664,
0.0002777777777777778
],
"bbox": [
-0.83,
52.520278,
-0.42,
52.76
],
"band_stats": [
{
"band": 1,
"dtype": "float32",
"nodata": null,
"valid_fraction": 1.0,
"min": 9.42,
"max": 217.53,
"mean": 89.23
}
],
"tags": {
"AREA_OR_POINT": "Area"
}
}
valid_fraction is the raster equivalent of a null count and answers "is this tile mostly nodata?" without opening it in a viewer. Read it together with nodata: here it is 1.0 because the Copernicus DEM declares no nodata value at all, so nothing can be counted as missing โ which is itself worth recording.
Example 2 โ a schema block that can be diffed between versions
def schema_block(gdf):
out = {}
for col, dtype in gdf.dtypes.items():
if col == gdf.geometry.name:
continue
s = gdf[col]
entry = {"type": str(dtype), "nulls": int(s.isna().sum())}
if pd.api.types.is_string_dtype(s) or s.dtype == object:
uniques = s.dropna().unique()
if len(uniques) <= 25:
entry["values"] = sorted(map(str, uniques))
else:
entry["distinct"] = int(len(uniques))
entry["max_length"] = int(s.dropna().astype(str).str.len().max() or 0)
elif pd.api.types.is_numeric_dtype(s):
entry["min"], entry["max"] = float(s.min()), float(s.max())
out[col] = entry
return out
for column, entry in schema_block(gpd.read_file("data/pois.gpkg")).items():
print(f"{column:7}", {k: (v[:6] + ["โฆ"] if isinstance(v, list) and len(v) > 6 else v) for k, v in entry.items()})
osm_id {'type': 'str', 'nulls': 0, 'distinct': 538, 'max_length': 10}
code {'type': 'int32', 'nulls': 0, 'min': 2001.0, 'max': 2964.0}
fclass {'type': 'str', 'nulls': 0, 'distinct': 55, 'max_length': 16}
name {'type': 'str', 'nulls': 345, 'distinct': 190, 'max_length': 48}
Recording the value list for low-cardinality columns turns a code list into something a consumer can validate against, and makes a new category visible in a diff. None of these columns is low-cardinality enough โ fclass has 55 values โ so the block records the count instead, and the fact that 345 of the 538 features have no name.
Example 3 โ the publish gate
REQUIRED = ["title", "feature_definition", "content_date", "lineage", "licence", "contact"]
KNOWN_LICENCES = {"CC0-1.0", "CC-BY-4.0", "CC-BY-SA-4.0", "ODbL-1.0",
"ODC-BY-1.0", "OGL-UK-3.0", "proprietary"}
def gate(record):
problems = []
problems += [f"missing {k}" for k in REQUIRED if not record.get(k)]
if record.get("licence") not in KNOWN_LICENCES:
problems.append(f"licence {record.get('licence')!r} is not a known identifier")
d = record.get("derived", {})
if not d.get("crs"):
problems.append("dataset has no CRS")
if d.get("features", 0) == 0:
problems.append("dataset is empty")
bbox = d.get("bbox") or [0, 0, 0, 0]
if bbox[0] > bbox[2] or bbox[1] > bbox[3]:
problems.append(f"bbox is not in west, south, east, north order: {bbox}")
if problems:
raise ValueError("metadata gate failed:\n " + "\n ".join(problems))
return True
print("as written:", gate(record))
broken = {**record, "licence": "OSM licence", "contact": "",
"derived": {**record["derived"], "bbox": record["derived"]["bbox"][2:] + record["derived"]["bbox"][:2]}}
try:
gate(broken)
except ValueError as exc:
print(exc)
as written: True
metadata gate failed:
missing contact
licence 'OSM licence' is not a known identifier
bbox is not in west, south, east, north order: [-0.431659, 52.752769, -0.797341, 52.533354]
The bbox order check is worth having explicitly: a reversed bounding box passed pystac's validation against the STAC schema when this site tested it, in the STAC validation guide.
Explanation
Why the derived block must come from the written file
The file is what ships, and writers change things: a shapefile truncates column names to ten characters, and a GeoPackage written from a layer that mixes Polygons and MultiPolygons promotes every one of them to MultiPolygon โ 538 of 538 here. A record derived from the in-memory frame documents a dataset that nobody received.
Why the authored fields need a gate rather than a template
A template with empty fields ships with empty fields. The only mechanism that reliably produces complete records is one that stops the publish, because that is the only moment at which somebody has both the knowledge and the motivation to fill them in.
Why null counts and value lists earn their place
They are the two things a consumer would otherwise compute themselves, and they are the two that change silently between versions. A column that was 2% null and is now 60% null is a broken upstream join, and the metadata diff is where it shows up first.
Why YAML for the authored half
The authored half is a small document with multi-line prose in it, maintained by a person and reviewed in a pull request. JSON is unpleasant to write multi-line strings in; YAML is not. The generated record is JSON because programs read it.
Edge cases or notes
- Multi-layer files need a record per layer. Plus one for the container.
total_boundson an empty frame returns[nan, nan, nan, nan]. Guard it, or the gate's bbox check compares NaNs and passes.- Round the bbox deliberately. Six decimal places is ample and keeps the record diffable.
- Datetime objects are not JSON. Pass
default=stror convert explicitly. - Do not store credentials in lineage. Record the endpoint, not the token.
- Keep the authored file in version control even when the data is not.
- Regenerate on every release. A record that is older than the data is worse than none.
- A record is a deliverable. Ship it in the zip.
Internal links
- Spatial metadata explained: what a dataset must tell you โ which fields and why
- How to read the metadata already inside your spatial files โ what you can recover from existing files
- How to store metadata inside a GeoPackage โ embedding the record
- How to embed provenance tags in a GeoTIFF โ the raster equivalent
- How to checksum spatial datasets so you can prove they match โ the content hash field
- Metadata standards compared: ISO 19115, STAC and Frictionless โ emitting the standard shapes
- How to validate a GeoDataFrame against a schema before analysis โ the consumer side of the schema block
- Metadata extents and dates do not match the data โ what the derived block prevents
FAQ
How do I generate metadata for a spatial dataset in Python?
Compute the derivable fields โ extent, CRS, schema, counts, geometry types, checksum โ from the written file, and merge them with a small authored YAML file holding title, feature definition, lineage, licence and contact.
Which fields should a person write?
Title, what one feature represents, the content date, lineage, licence, attribution and a contact. Everything else can be derived.
Should the record be JSON or YAML?
Author in YAML because a human maintains it; generate JSON because programs consume it.
Why derive from the file rather than the DataFrame?
Because writers change things โ shapefiles truncate column names, GeoPackages promote geometry types โ and the record must describe what was actually shipped.
How do I stop records from going stale?
Regenerate the derived block on every write, and gate the publish step so an incomplete authored block raises rather than warns.
What checksum should go in the record?
A content hash of the canonicalised data. A byte hash of a GeoPackage or GeoJSON changes between identical writes.