transform(GeoDataFrame) uses geopandas' per-geometry affine_transform (Python loop): 41 s for 5.5 M polygons vs 1.8 s with shapely.transform

[!NOTE]
This whole message is AI-generated. The issue was automatically discovered and reported by an AI agent (Claude) during an autonomous bug hunt on the spatialdata code base. It has not been verified or triaged by a human yet; the needs: triage label is set so that a maintainer can confirm it. The reproduction script below was executed by the agent in an isolated environment (see Environment) and its output is pasted verbatim.

Summary

geopandas delegates affine_transform to shapely.affinity.affine_transform per geometry in Python. On the 5.48 M Visium HD bins: spatialdata.transform 41.3 s, geopandas affine_transform 42.7 s, vectorised shapely.transform 1.8 s with identical coordinates. The synthetic script below (2 M squares) reproduces the ratio.

Severity (agent's assessment): medium/high for real data — on the hot path of get_extent(exact=True), transform_to_coordinate_system, aggregate (when transformations differ), transform_to_data_extent and plotting

Where: src/spatialdata/_core/operations/transform.py, GeoDataFrame overload (data.geometry.affine_transform(shapely_notation))

Expected behaviour

Vectorised transformation (seconds, not minutes).

Reproduction

Save as repro.py and run uv run repro.py (the PEP 723 header pins spatialdata to the commit the bug was found on; replace the URL fragment with @main to test the current main branch).

# /// script
# requires-python = ">=3.12"
# dependencies = [
#     "spatialdata @ git+https://github.com/scverse/spatialdata.git@ccf1ea048d054b6624214bf618008a9f9ae223e0",
# ]
# ///
"""transform(GeoDataFrame) uses geopandas' per-geometry affine_transform (Python loop); shapely.transform is ~20x faster."""
import time
import warnings
import numpy as np
import shapely
import geopandas as gpd
from spatialdata import transform
from spatialdata.models import ShapesModel
from spatialdata.transformations import Affine, get_transformation

warnings.simplefilter("ignore")
n = 2_000_000
rng = np.random.default_rng(0)
xs, ys = rng.uniform(0, 5000, n), rng.uniform(0, 5000, n)
m = np.array([[0.9, 0.1, 5.0], [-0.1, 0.9, 7.0], [0.0, 0.0, 1.0]])
shapes = ShapesModel.parse(gpd.GeoDataFrame({"geometry": shapely.box(xs, ys, xs + 2, ys + 2)}), transformations={"global": Affine(m, input_axes=("x", "y"), output_axes=("x", "y"))})
t0 = time.time(); out = transform(shapes, to_coordinate_system="global"); t_sd = time.time() - t0
A, b = m[:2, :2], m[:2, 2]
t0 = time.time(); vec = shapely.transform(shapes.geometry.values, lambda c: c @ A.T + b); t_vec = time.time() - t0
same = np.allclose(shapely.get_coordinates(vec[:10000]), shapely.get_coordinates(out.geometry.values[:10000]))
print(f"{n:,} polygons: spatialdata.transform {t_sd:.1f}s | shapely.transform (vectorised) {t_vec:.1f}s | identical coordinates: {same}")
bug = same and t_sd > 5 * t_vec
print("VERDICT:", "BUG REPRODUCED (large avoidable cost)" if bug else "NOT REPRODUCED")

Observed output

2,000,000 polygons: spatialdata.transform 16.6s | shapely.transform (vectorised) 0.7s | identical coordinates: True
VERDICT: BUG REPRODUCED (large avoidable cost)

Possible fix direction (unverified)

A, b = matrix[:-1, :-1], matrix[:-1, -1]
transformed = gpd.GeoSeries(shapely.transform(data.geometry.values, lambda c: c @ A.T + b), index=data.index, crs=data.crs)

(shapely.transform handles Polygon/MultiPolygon/holes/Points and include_z for 3D geometries.)

Environment

uv run repro.py with the PEP 723 metadata in the script (fresh, isolated environment; spatialdata built from main @ ccf1ea0 (2026-08-28); Python 3.13, latest releases of the dependencies at run time: pandas 3.0, anndata 0.13, zarr 3.3, dask 2026.8, numpy 2.5, geopandas 1.1, shapely 2.1). macOS (arm64). Also reproduced in a second environment with pandas 2.3.3 / anndata 0.12.11 / numpy 2.4.4 / zarr 3.2.1.

Automatically generated; discovered by an AI agent (Claude) and not yet reviewed by a human.

Start with src/spatialdata/_core/operations/transform.py and inspect the GeoDataFrame overload that calls data.geometry.affine_transform(shapely_notation). Run the provided repro.py with uv to establish the baseline, then verify the transformation remains coordinate-equivalent while avoiding the reported per-geometry cost. Done means the benchmark improves substantially without changing the output coordinates.