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.

The derived half of a metadata record against the authored half. two halves, two mechanismsderived at write timeโœ“extent, CRS, geometry typesโœ“feature count and schemaโœ“null counts per columnโœ“content checksumauthored, and gatedโœ“title and feature definitionโœ“content dateโœ“lineage narrativeโœ“licence, attribution, contactNever type what a program can compute; always gate what it cannot spatialworkflow.io
Anything a program can compute should never be typed; anything it cannot must be gated.

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.

One authored file and a written dataset producing three metadata artefacts. one source, three outputsauthored YAMLreviewed like coderead the writtenfilenot the framemerge + gateraises if incompletesidecar + tags +STACgeneratedReading the written file is what catches what the driver changed spatialworkflow.io
One authored file in, three artefacts out, none of them typed twice.

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.

The checks a metadata publish gate performs. the publish gateevery required authored field presentlicence is a known identifierthe dataset has a CRSthe dataset is not emptybbox is west, south, east, northpystac passed a reversed onewarn instead of raisea warning shipsFive assertions, and one anti-pattern spatialworkflow.io
The bbox order check is here because the STAC schema check passed a reversed box.

Edge cases or notes

  • Multi-layer files need a record per layer. Plus one for the container.
  • total_bounds on 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=str or 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.

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.