This guide is part of Automated Vector & Raster Cleaning Workflows, the reference for building repeatable Python pipelines that repair, normalize, and standardize spatial data before it reaches an analytical store.

Choosing a Vector Output Format for Spatial Pipelines

The output format is the decision a spatial ETL pipeline lives with longest. Every cleaning stage upstream — geometry repair, CRS normalization, schema harmonization — exists to produce one durable artifact, and the container you write that artifact into determines how cheaply it can be queried, streamed, versioned, and re-read for the next year. For GIS analysts and data engineers, the reflex to write a Shapefile because “everything opens it” quietly imposes a 10-character column-name limit, a 2 GB per-file ceiling, and a four-file bundle that a single mv can corrupt. Meanwhile the format that would have made a downstream DuckDB scan ten times faster, or served a web map directly from object storage, went unchosen.

This guide compares the four formats a modern pipeline actually ships to — GeoParquet, FlatGeobuf, GeoPackage, and Shapefile — as output sinks, not as interchange curiosities. The distinctions that matter are structural: whether storage is columnar or row-oriented, whether a spatial index ships inside the file, whether a client can read a bounding-box slice over an HTTP range request, and whether the coordinate reference system survives the write intact.

Prerequisites & Environment

The examples assume a current GeoPandas stack. GeoParquet write support relies on pyarrow; FlatGeobuf and GeoPackage write support is provided through GDAL’s pyogrio engine, which GeoPandas 1.0+ uses by default.

• Python 3.10+ — modern type hints and compatibility with GeoPandas 1.x.
• Core dependencies — geopandas>=1.0, pyarrow>=14.0 (GeoParquet + ZSTD), pyogrio>=0.7 (GDAL-backed FlatGeobuf/GeoPackage I/O), shapely>=2.0 (vectorized geometry), and optionally duckdb>=0.10 with its spatial extension for in-place GeoParquet analytics.
• GDAL 3.6+ — required for stable FlatGeobuf write support and GeoPackage R-tree indexing through pyogrio.

pip install "geopandas>=1.0" "pyarrow>=14.0" "pyogrio>=0.7" "shapely>=2.0" "duckdb>=0.10"

Confirm the writer engine and Arrow availability before relying on GeoParquet in a job:

import geopandas as gpd
import pyogrio

print(gpd.__version__)           # expect 1.x
print(pyogrio.__gdal_version__)  # expect (3, 6, x) or newer
# GeoParquet requires pyarrow; ImportError here means ZSTD writes will fail
import pyarrow  # noqa: F401

Version / Compatibility Matrix

This is the comparison matrix. Each row is a structural property that determines how the format behaves as a pipeline sink; the rightmost column names the workload each format is built for.

Two rows deserve emphasis for pipeline authors. The Shapefile field-name limit is not cosmetic: the DBF component silently truncates any attribute name to 10 bytes, so population_density and population_dens... collide into the same populati_1-style mangled key, corrupting the schema you spent an upstream stage harmonizing. And GeoParquet’s columnar layout is the single property that unlocks predicate pushdown — a WHERE on one column reads only that column’s chunks — which is why it pairs with the analytical engines covered in PostGIS vs DuckDB for analytical loads.

Step-by-Step Implementation

Each step writes the same cleaned GeoDataFrame to a different sink so the trade-offs are concrete rather than abstract. Assume gdf has already passed geometry repair and CRS normalization and carries a valid gdf.crs.

Step 1 — Write GeoParquet with ZSTD and Preserved CRS

GeoParquet embeds the CRS as PROJJSON inside the file’s geo metadata key, so there is no sidecar to lose. ZSTD gives the best compression ratio for the mixed numeric-and-WKB payload typical of vector data, and it decompresses fast enough that analytical scans stay I/O-bound rather than CPU-bound.

import geopandas as gpd

def write_geoparquet(gdf: gpd.GeoDataFrame, path: str) -> None:
    """Write a cleaned GeoDataFrame to GeoParquet with ZSTD compression.

    CRS is serialized into the GeoParquet `geo` metadata automatically;
    no .prj sidecar is produced or needed.
    """
    if gdf.crs is None:
        raise ValueError("Refusing to write GeoParquet without a CRS set")
    gdf.to_parquet(
        path,
        compression="zstd",       # best ratio for WKB + attribute mix
        compression_level=9,
        geometry_encoding="WKB",  # portable default; "geoarrow" for 1.1 readers
        write_covering_bbox=True, # GeoParquet 1.1 bbox column for spatial filters
        index=False,
    )

write_covering_bbox=True adds the GeoParquet 1.1 bounding-box covering column, letting readers skip row groups that fall outside a query window — the columnar equivalent of a coarse spatial index.

Step 2 — Write FlatGeobuf with a Packed Spatial Index

FlatGeobuf builds a packed Hilbert R-tree over all features and stores it in the file header at write time. That index is what lets a remote client fetch only the features inside a bounding box using HTTP range requests, without a server or a database.

def write_flatgeobuf(gdf: gpd.GeoDataFrame, path: str) -> None:
    """Write a GeoDataFrame to FlatGeobuf with its packed R-tree index.

    The spatial index is written into the header so bbox-filtered range
    reads work directly against the file on object storage.
    """
    if gdf.crs is None:
        raise ValueError("Refusing to write FlatGeobuf without a CRS set")
    gdf.to_file(
        path,
        driver="FlatGeobuf",
        engine="pyogrio",
        SPATIAL_INDEX="YES",  # build the packed Hilbert R-tree in the header
    )

FlatGeobuf enforces a single geometry type per file. If gdf mixes points and polygons, split by geom_type before writing or the driver raises on the first mismatched feature.

Step 3 — Write GeoPackage for Desktop-GIS Interchange

When a human analyst will open the output in QGIS or ArcGIS, GeoPackage is the right sink: one SQLite file, arbitrary column names, an R-tree spatial index, and support for multiple layers and mixed geometry types in the same container.

def write_geopackage(
    gdf: gpd.GeoDataFrame,
    path: str,
    layer: str,
) -> None:
    """Write (or append) a named layer to a single GeoPackage container."""
    if gdf.crs is None:
        raise ValueError("Refusing to write GeoPackage without a CRS set")
    gdf.to_file(
        path,
        driver="GPKG",
        layer=layer,
        engine="pyogrio",
        SPATIAL_INDEX="YES",  # R-tree virtual table for fast desktop queries
    )

Because a GeoPackage holds many layers, a single call can append a roads layer beside an existing parcels layer in the same file — the multi-layer property Shapefile cannot express.

Step 4 — Verify CRS and Schema Round-Trip

The write is not done until the CRS and column names survive a read-back. This check belongs in the pipeline, not in a notebook, because a silent CRS drop is the failure mode that corrupts every downstream spatial join.

def assert_roundtrip(original: gpd.GeoDataFrame, path: str) -> None:
    """Read the written file back and assert CRS and columns are intact."""
    reloaded = gpd.read_file(path) if not path.endswith(".parquet") \
        else gpd.read_parquet(path)
    assert reloaded.crs == original.crs, (
        f"CRS drift: wrote {original.crs}, read {reloaded.crs}"
    )
    lost = set(original.columns) - set(reloaded.columns)
    assert not lost, f"Columns lost on round-trip (likely 10-char truncation): {lost}"

Run assert_roundtrip against a Shapefile sink and it will surface exactly the truncation the format imposes — the diagnostic that motivates the dedicated conversion recipe below.

Step 5 — Encode the Format Decision as a Function

Rather than hard-coding a format in every job, make the sink a policy the pipeline resolves from the consumer’s needs. This keeps the decision explicit and auditable.

from typing import Literal

Consumer = Literal["analytics", "web_stream", "desktop_gis", "legacy_export"]

def resolve_output_writer(consumer: Consumer):
    """Map a downstream consumer to the format best suited to it."""
    return {
        "analytics": write_geoparquet,     # columnar, DuckDB-native
        "web_stream": write_flatgeobuf,    # bbox range reads
        "desktop_gis": write_geopackage,   # single-file, multi-layer
        "legacy_export": write_shapefile,  # only when explicitly required
    }[consumer]

Advanced Patterns & Edge Cases

Migrating a Shapefile Archive to a Columnar Sink

The most common pipeline task is not choosing a format for greenfield data but retiring a directory of legacy Shapefiles into an analytical store. That migration has its own hazards — 10-character truncation, CRS recovery from a possibly-missing .prj, and ZSTD level tuning — covered end to end in converting Shapefiles to GeoParquet with GeoPandas. Treat it as the batch counterpart to Step 1 here.

Reconciling Column Names Before They Reach the Sink

A format that permits arbitrary column names does not absolve the pipeline of harmonizing them. If several sources feed one GeoParquet output, their attribute schemas must agree first — otherwise the columnar file inherits three spellings of the same field. Run the reconciliation described in attribute mapping and schema harmonization upstream of the write, so the format choice preserves a schema that is already clean.

Partitioned GeoParquet for Large Outputs

A single GeoParquet file is fine to tens of millions of features, but beyond that, partition by a spatial or temporal key so readers touch only relevant fragments. Write a directory of part=*.parquet files keyed on a quadkey or date, and DuckDB or pyarrow will treat the directory as one logical dataset with partition pruning — the columnar analogue of tiling.

Performance Optimization

The table below is indicative only — absolute numbers depend on geometry complexity, attribute width, hardware, and compression level. It reflects the relative behavior of one representative dataset (roughly one million polygon features with a dozen attributes) and should be read for ratios between formats, not as a benchmark to quote.

The pattern that matters: GeoParquet wins decisively on size and on selective reads because compression and column pruning compound, while FlatGeobuf wins on bounding-box delivery because its index answers a spatial query without reading the body. Measure on your own data before committing:

import time
from pathlib import Path

def benchmark_write(gdf: gpd.GeoDataFrame, out_dir: str) -> dict[str, dict]:
    """Write the same GeoDataFrame to each format and record size + write time."""
    results: dict[str, dict] = {}
    writers = {
        "data.parquet": lambda p: write_geoparquet(gdf, p),
        "data.fgb": lambda p: write_flatgeobuf(gdf, p),
        "data.gpkg": lambda p: write_geopackage(gdf, p, layer="features"),
    }
    for name, writer in writers.items():
        path = str(Path(out_dir) / name)
        start = time.perf_counter()
        writer(path)
        results[name] = {
            "write_seconds": round(time.perf_counter() - start, 3),
            "bytes": Path(path).stat().st_size,
        }
    return results

Integration into ETL Pipelines

The format decision is a load-stage concern, but it reaches back into the whole pipeline. Write to a temporary path (data.parquet.part) and rename atomically only after assert_roundtrip passes, so an interrupted job never leaves a half-written file that the next run mistakes for complete. Route features that fail the round-trip check — a Shapefile that lost columns, a FlatGeobuf mixed-geometry rejection — to a dead-letter path rather than failing the batch, so one bad layer never blocks the rest.

For analytical targets, the columnar sink is only half the story: the query engine on the other side determines whether column pruning and predicate pushdown are actually exercised. A GeoParquet directory read by DuckDB’s spatial extension gives serverless analytical scans, while the same data loaded into PostGIS gives transactional upserts and a persistent GIST index — the trade-off worked through in the PostGIS versus DuckDB comparison for analytical loads. Keep the sink format and the load target chosen together, because a row-oriented format feeding a columnar engine wastes the engine, and a columnar file re-parsed row-by-row wastes the format.

Troubleshooting Reference

Related

• Converting Shapefiles to GeoParquet with GeoPandas — batch-convert a Shapefile archive while preserving CRS and surviving 10-character column truncation
• Attribute Mapping & Schema Harmonization — reconcile column names and types across sources before they reach a columnar sink
• PostGIS vs DuckDB Spatial for Analytical Loads — pick the load target that consumes the format you write
• Automated Vector & Raster Cleaning Workflows — the full reference for repairing and normalizing spatial data before load