Speeding Up GeoPandas I/O with pyogrio and Arrow

A gpd.read_file() call that takes two minutes on a multi-million-feature layer is usually not limited by your disk, your CPU or the driver — it is limited by how many times the Python interpreter is asked to touch a single feature. This guide is for anyone whose pipeline spends most of its wall clock in the read step: it sits under GeoPandas DataFrames Explained in Mastering Core Geospatial Python Libraries, and covers the pyogrio engine that GeoPandas 1.x defaults to, the use_arrow=True transfer path, the read-time pushdowns that skip work entirely, and a benchmark harness honest enough to publish.

Why This Approach / What Goes Wrong

Both I/O engines call the same GDAL/OGR C library, open the same file with the same driver and produce the same GeoDataFrame. The difference is entirely in the hand-off. Fiona exposes OGR as a feature iterator: every feature is materialised as a Python dict of properties plus a GeoJSON-ish geometry mapping, and GeoPandas then builds the geometry column by constructing one Shapely object per row. On a 2.4-million-feature layer that is 2.4 million interpreter round-trips before any analysis starts, and the per-feature overhead — dict allocation, attribute boxing, reference counting — dominates the actual decoding.

pyogrio removes the loop rather than optimising it. It asks OGR for the whole layer in bulk from C, fills typed NumPy arrays column by column, and converts the geometry column with a single vectorised shapely.from_wkb call over a buffer of WKB blobs. Python sees a handful of finished arrays, not N features. That is why the speed-up is roughly a constant factor across file sizes instead of something that only shows up on small files: you have deleted a per-row cost, not amortised it.

Where the per-feature Python loop disappears in the GeoPandas I/O stack Two stacks sit on the same GDAL and OGR driver layer. On the left, the fiona engine returns a feature iterator, which becomes one Python dictionary and one Shapely object per feature, and only then a GeoDataFrame assembled from those dictionaries. On the right, pyogrio performs a bulk vectorised read that produces columnar batches converted to NumPy or Arrow in a handful of C calls, and the GeoDataFrame is simply wrapped around those ready-made arrays. The driver work at the bottom is identical in both stacks, so the per-feature interpreter round-trips are the entire difference. The same driver, two hand-offs read the stacks bottom-up: only the middle tier differs engine="fiona" engine="pyogrio" GeoDataFrame assembled out of N Python dicts GeoDataFrame wrapped around ready-made arrays one dict + one Shapely object per feature N interpreter round-trips · the bottleneck columnar batches → NumPy / Arrow a handful of C calls, whatever N is fiona.open() lazy feature iterator pyogrio.read_dataframe() bulk vectorised read GDAL / OGR drivers identical C library, identical file parsing, identical bytes off disk Nothing about the file changed — only how many times Python was asked to look at a feature.
The driver tier is shared; the per-feature Python tier exists only on the Fiona path, and deleting it is where the factor of five to ten comes from.

use_arrow=True shortens the hand-off one step further. GDAL 3.6 added an Arrow array stream interface, so instead of pyogrio calling OGR's per-feature field accessors in C and packing NumPy arrays itself, the driver emits Arrow record batches — typed, columnar, dictionary-encoded for repeated strings — and pyarrow hands those buffers to pandas with far less copying. The win is largest on wide attribute tables and on string-heavy layers; on a narrow layer whose cost is dominated by parsing complex polygons into Shapely, it is modest. It is a transfer-format change, not a magic flag.

Granularity of the hand-off in three read paths Three lanes decode the same one million features. The Fiona lane shows a long row of small cells, one Python dictionary per feature, meaning one million interpreter round-trips. The default pyogrio lane shows the same long row of small cells but they are C-level field accesses, so there is no Python cost per feature. The pyogrio lane with use_arrow set to true shows only a few wide bars, because GDAL hands over whole columnar record batches. The file, the driver and the resulting GeoDataFrame are identical in all three lanes; only the granularity of the hand-off changes. One million features, three hand-off granularities each cell is one unit of work crossing the boundary into Python Fiona engine one dict per feature … × 1,000,000 1,000,000 interpreter round-trips pyogrio, default one C field access per feature … × 1,000,000 stays inside C · zero Python per feature pyogrio + Arrow one record batch per chunk … × ~16 batches columnar buffers handed over whole Same file, same driver, same resulting GeoDataFrame — only the size of each unit crossing the boundary changes.
Fiona pays per feature in Python, default pyogrio pays per feature in C, and the Arrow path pays per batch — which is why the third step is a smaller jump than the first.

The largest win, though, is not a faster read: it is a read that never happens. columns=, bbox=, mask=, where= and rows= are all evaluated by GDAL before anything reaches Python, so unwanted fields are never decoded and unwanted features are never touched. On an indexed format the bounding-box filter also lets the driver skip whole pages of the file. Choosing the container matters as much as choosing the engine — see GeoParquet vs Shapefile for storage — and if the file is larger than RAM the question stops being speed and becomes memory, which is the subject of GeoPandas vs Fiona for large files.

Prerequisites

conda install -c conda-forge "geopandas=1.0.*" "pyogrio=0.9.*" "pyarrow=16.*" \
    "shapely=2.0.*" "fiona=1.9.*"

Install the whole set from conda-forge rather than mixing wheels: pyogrio, Fiona and pyproj each bind GDAL and PROJ, and two different GDAL builds in one environment produce import-time ABI errors or, worse, a silently different driver set on each engine — which quietly invalidates any benchmark you run.

Step-by-Step Implementation

1. Confirm which engine is actually running.

Never benchmark or tune before you know what you have. GeoPandas 1.x selects pyogrio when it is importable and falls back to Fiona otherwise, so a stale environment can silently give you the slow path.

import geopandas as gpd
import pyogrio

print(gpd.options.io_engine)              # None -> auto-select (pyogrio if installed)
print(pyogrio.__version__)                # 0.9.0
print(pyogrio.__gdal_version_string__)    # 3.9.2  (>= 3.6 required for Arrow)

gpd.options.io_engine = "pyogrio"         # pin it explicitly for reproducible runs

2. Read the layer's metadata without reading the layer.

pyogrio.read_info() opens the dataset, reads the header and closes it. You get the feature count, field names, dtypes and CRS for the price of a header read — enough to plan the pushdowns in the next steps.

info = pyogrio.read_info("buildings.fgb")
print(info["features"], "features")        # 2412885 features
print(info["crs"])                         # EPSG:25832
print(list(info["fields"])[:5])            # ['osm_id', 'use', 'height_m', 'levels', 'source']

3. Turn on the Arrow transfer path.

use_arrow=True is accepted by gpd.read_file() and forwarded to pyogrio. Set it per call, or once for the whole process with the PYOGRIO_USE_ARROW=1 environment variable when you cannot edit every call site.

buildings = gpd.read_file("buildings.fgb", use_arrow=True)
print(type(buildings), buildings.shape)    # <class 'geopandas.geodataframe.GeoDataFrame'> (2412885, 62)

4. Read only the columns you need.

A 62-column building layer read for a height histogram wastes 60 columns of decoding, allocation and memory. columns= is a pyogrio-only argument that selects fields inside the driver. When you want attributes with no geometry at all, drop to the pyogrio API directly and pass read_geometry=False — it returns a plain pandas DataFrame and skips WKB parsing entirely.

# Geometry plus two attributes, selected inside GDAL
heights = gpd.read_file(
    "buildings.fgb",
    columns=["use", "height_m"],
    use_arrow=True,
)

# Attribute-only: no geometry column, no WKB parsing, returns a pandas DataFrame
attrs = pyogrio.read_dataframe("buildings.fgb", columns=["use", "height_m"],
                               read_geometry=False, use_arrow=True)

5. Push the spatial filter down — in the file's own CRS.

bbox= hands a rectangle to GDAL, which uses the format's spatial index to skip non-matching features before decoding. The trap is projection: with the pyogrio engine the bounding box must already be in the dataset's CRS — unlike the Fiona engine, pyogrio does not reproject a GeoSeries for you. Convert your area of interest explicitly against the CRS reported in step 2. Bounding-box tuples are always (minx, miny, maxx, maxy) in x/y order regardless of what axis order the EPSG authority declares, because pyproj transformations behind to_crs() are built with always_xy=True; passing latitude first silently matches nothing.

from shapely.geometry import box

# AOI drawn in lon/lat, reprojected to the file's CRS (UTM 32N here) before use
aoi = gpd.GeoSeries([box(11.52, 48.10, 11.62, 48.17)], crs="EPSG:4326")
aoi = aoi.to_crs(info["crs"])
minx, miny, maxx, maxy = aoi.total_bounds

core = gpd.read_file(
    "buildings.fgb",
    bbox=(minx, miny, maxx, maxy),      # already in EPSG:25832
    columns=["use", "height_m"],
    use_arrow=True,
)
print(len(core), "of", info["features"], "features decoded")

Use mask=<geometry> instead of bbox= when the area of interest is a real polygon rather than a rectangle; the two are mutually exclusive and passing both raises. For attribute cuts, where= sends an OGR SQL predicate into the driver, and rows=slice(0, 5000) gives you a cheap sample for schema checks.

Three read-time pushdowns narrowing a large layer before it reaches Python A shrinking stack of bars. The full layer is twelve point four gigabytes, forty eight million features and sixty one columns. Applying a bounding box lets the spatial index skip whole pages and leaves one point nine gigabytes and six point two million features. Adding an OGR SQL where clause on the use field leaves four hundred megabytes and one point one million features. Selecting two columns leaves forty eight megabytes. Each filter is evaluated inside GDAL, so the discarded bytes are never decoded into Python objects at all. Every filter you push down is decoding that never happens whole layer · 12.4 GB · 48 M features · 61 columns bbox=(minx, miny, maxx, maxy) spatial index skips whole pages of the file 1.9 GB · 6.2 M features where="use = 'residential'" OGR SQL evaluated in the driver, not in pandas 0.4 GB · 1.1 M features columns=["use", "height_m"] 59 fields never leave the driver 48 MB what the GeoDataFrame actually holds Bar widths are proportional to bytes decoded — the fastest read is the one that touches the fewest of them.
Engine choice changes the cost per feature; pushdown changes how many features and fields there are to pay for at all.

6. Write with the same engine.

to_file() routes through pyogrio too, and the write path benefits for the same reason: the whole column is handed to OGR in bulk instead of one dst.write(feature) call per row. Writes are typically two to four times faster than the Fiona path, and the gap is smaller than on reads because the driver's own encoding work — building a spatial index, compressing, fsyncing — is unchanged.

# Bulk write through pyogrio; FlatGeobuf and GeoPackage both carry a spatial index
core.to_file("munich_residential.fgb", driver="FlatGeobuf")

# Or the pyogrio API directly, for explicit layer/driver control
pyogrio.write_dataframe(core, "munich.gpkg", layer="residential", driver="GPKG")

# GeoParquet does not go through OGR at all — pandas writes it via pyarrow
core.to_parquet("munich_residential.parquet")

7. Time it, do not guess.

Measure with a warm-up run, several repeats and a median — a single cold-cache timing mostly measures your disk. Compare identical work: same file, same columns, same process.

import statistics
import time

import geopandas as gpd

def timed(fn, repeats=5):
    """Median wall-clock seconds over `repeats` runs, after one warm-up."""
    fn()                                            # warm the OS page cache
    samples = []
    for _ in range(repeats):
        t0 = time.perf_counter()
        fn()
        samples.append(time.perf_counter() - t0)
    return statistics.median(samples)

path = "buildings.fgb"
cases = {
    "fiona":                  lambda: gpd.read_file(path, engine="fiona"),
    "pyogrio":                lambda: gpd.read_file(path, engine="pyogrio"),
    "pyogrio + arrow":        lambda: gpd.read_file(path, engine="pyogrio", use_arrow=True),
    "arrow + 2 columns":      lambda: gpd.read_file(path, engine="pyogrio", use_arrow=True,
                                                    columns=["use", "height_m"]),
}

baseline = None
for name, fn in cases.items():
    elapsed = timed(fn)
    baseline = baseline or elapsed
    print(f"{name:<20} {elapsed:7.2f} s   {baseline / elapsed:5.1f}x")
# fiona                 142.10 s     1.0x
# pyogrio                23.14 s     6.1x
# pyogrio + arrow        15.82 s     9.0x
# arrow + 2 columns       5.41 s    26.3x

Four habits keep numbers like these defensible. Warm the cache before timing, or you are ranking your storage rather than the engines. Report a median of several runs, because the first read of a session also pays for GDAL driver registration and PROJ database initialisation. Time each case in the same interpreter on the same file, since a format conversion between runs changes the experiment. And resist extrapolating across shapes of data: a wide, string-heavy attribute table and a narrow layer of hundred-vertex polygons sit at opposite ends of the Arrow speed-up range, so a factor measured on one says little about the other. If a result looks implausible, print gpd.options.io_engine and pyogrio.__gdal_version_string__ inside the timed function — a silent fall back to Fiona, or an environment with two GDAL builds, explains most surprising benchmarks.

Verification

A benchmark is only meaningful if every case returned the same data. Assert parity first, then report the timings — and state the machine, the GDAL version and the cache state alongside any number you publish.

import geopandas as gpd
from geopandas.testing import assert_geodataframe_equal

ref = gpd.read_file("buildings.fgb", engine="fiona")
fast = gpd.read_file("buildings.fgb", engine="pyogrio", use_arrow=True)

assert len(ref) == len(fast), "engines disagree on feature count"
assert ref.crs == fast.crs, "CRS lost or altered on one path"
assert set(ref.columns) == set(fast.columns), "column sets differ"

# Same geometries, same attributes; dtypes may differ (Arrow gives nullable types)
assert_geodataframe_equal(ref, fast[ref.columns], check_dtype=False)
print(f"parity OK: {len(ref):,} features, {ref.crs}")
# parity OK: 2,412,885 features, EPSG:25832
Measured read times for one layer under four read strategies A horizontal bar chart of median read times for the same two point four million feature FlatGeobuf. The Fiona engine takes one hundred and forty two seconds and is the baseline. The pyogrio engine takes twenty three point one seconds, about six times faster. Adding use_arrow brings it to fifteen point eight seconds, about nine times faster. Restricting the read to two columns brings it to five point four seconds, about twenty six times faster. All runs are medians of five warm-cache repeats on one machine. Reading the same 2.4 M-feature FlatGeobuf four ways median of 5 warm-cache runs · GDAL 3.9 · local NVMe · your numbers will differ engine="fiona" 142.1 s · baseline engine="pyogrio" 23.1 s · 6.1× faster + use_arrow=True 15.8 s · 9.0× faster + columns=[2 of 62] 5.4 s · 26× faster 0 s 50 s 100 s 150 s The engine swap is the big step; Arrow is a real but smaller one; not reading columns you never use beats both.
Ranked honestly, the ordering is engine, then pushdown, then Arrow — and the last two compound.

Edge Cases & Debugging

Frequently Asked Questions

Do I have to change my code to use pyogrio? No. From GeoPandas 1.0 the engine is selected automatically and pyogrio wins whenever it is importable, so an unmodified gpd.read_file() is already on the fast path. The changes worth making are the deliberate ones: pinning gpd.options.io_engine = "pyogrio" so a benchmark cannot silently fall back, adding use_arrow=True, and passing columns= and bbox=.

Is use_arrow=True always faster? Not always, and rarely dramatically on its own. It removes per-field conversion between GDAL and Python, so it shines on wide tables with many string columns and does little on a two-column layer of complex polygons where WKB-to-Shapely conversion dominates. Treat it as a measured optimisation on your own data, not a default assumption — the harness in step 7 answers the question in a minute.

Does pyogrio also make GeoParquet faster? GeoParquet does not go through pyogrio at all. gdf.to_parquet() and gpd.read_parquet() are implemented on pyarrow directly, which is why they are usually the fastest option available and support their own column and row-group filtering. Compare the containers in GeoParquet vs Shapefile for storage before assuming the engine is your bottleneck.

Can pyogrio read a file bigger than RAM? No — it still materialises the full result set, so a whole-layer read of an oversized file fails the same way it always did. Pushdowns help when the result fits even though the file does not; when neither fits, you need a streaming or chunked strategy instead, which is exactly the trade-off laid out in GeoPandas vs Fiona for large files.

Related