GeoParquet vs Shapefile for Storage
The Shapefile has been the default vector format for thirty years, and it shows: 10-character field names, a 2 GB ceiling, a CRS in a separate file that goes missing, and multiple sidecar files per dataset. GeoParquet fixes all of it. This guide compares the two for storage and interchange and shows how to migrate with GeoPandas. It is for anyone choosing a working format for a Python pipeline. It sits under Cloud-Native Geospatial Formats in Geospatial Data Ingestion & Processing Workflows.
Why This Approach / What Goes Wrong
The Shapefile is not one file but a family — .shp (geometry), .shx (index), .dbf (attributes), and an optional .prj (CRS) — and its limits cause real data loss at every layer:
- Field names truncate to ten characters. The DBF attribute table, inherited from dBASE III, caps column names at ten bytes, so
population_densitysilently becomespopulatioand collides withpopulation. There is no way to recover the original names automatically once written. - The CRS lives in a droppable sidecar. The
.prjfile is optional and frequently absent after a copy or an email, leaving data with no coordinate system — the single most common ingestion bug, addressed in Coordinate Reference System Transformations and diagnosed at library level in Coordinate Systems with PyProj. - A hard 2 GB ceiling. Both
.shpand.dbfuse 32-bit byte offsets, so either file fails past 2 GB — a limit modern building-footprint or parcel datasets hit routinely. - Lossy attribute types. DBF has no native boolean, no timezone-aware datetime, no int64, and a fixed numeric field width, so dtypes are coerced — often to strings — on write. DBF also has no true
NULL: a missing number is written as blanks and read back as0orNaNdepending on the driver, so "no data" and "zero" become indistinguishable in exactly the columns where the difference matters. - Geometry types are promoted on write. A Shapefile declares one geometry type for the whole file and stores every polygon as a multipart record, so a
Polygoncolumn comes back asMultiPolygonafter a round trip and a layer holding both lines and points cannot be written at all. If your data legitimately mixes types, see handling mixed geometry types in a GeoDataFrame — the fix is a format that tolerates them, not a cast. - Ring orientation flips. The Shapefile specification wants clockwise exterior rings; GeoJSON's RFC 7946 wants counter-clockwise. Converting between them repeatedly with tools that "helpfully" normalise orientation produces diffs in files whose geometry never changed, which is maddening when the files are under version control.
GeoParquet stores everything in one compressed columnar file: full-length field names, native Arrow dtypes, and the CRS embedded in the file schema as a PROJJSON block. Because it is Apache Parquet with a geo metadata key, features are grouped into row groups whose per-column statistics include a bounding box, so a reader can skip whole groups that miss a query window — the basis for querying GeoParquet with DuckDB Spatial without downloading the file. The only reason to write a Shapefile today is a downstream tool that accepts nothing else.
The decision is rarely GeoParquet or Shapefile in the abstract — it is which format owns which role. Use GeoParquet as the working and interchange format between pipeline stages, and treat Shapefile purely as an export target for a legacy consumer.
| Concern | Shapefile | GeoParquet |
|---|---|---|
| Files per dataset | 3-4 sidecars | 1 |
| Field-name length | 10 chars | unlimited |
| Size ceiling | 2 GB per .shp/.dbf |
none (chunked) |
| CRS storage | separate .prj, droppable |
embedded PROJJSON |
| Attribute types | DBF-coerced | native Arrow dtypes |
| Compression | none | ZSTD / Snappy |
| Partial reads | full-file scan | row-group bbox skip |
Which GeoParquet version you write matters
"GeoParquet" is a metadata specification layered on Parquet, and the version stamped into the file changes what readers can do with it. Version 1.0 fixed the geo metadata key, WKB geometry encoding, and per-column CRS as PROJJSON. Version 1.1 added two things that show up directly in Python: GeoArrow encodings, which store coordinates as native nested Arrow lists instead of opaque WKB blobs, and a covering entry that names a bounding-box column so a reader knows which columns to test for spatial pruning.
That second one is the practical difference. Row groups always carry min/max statistics, but a reader has to be told which columns describe the envelope before it can use them as a spatial filter. A GeoParquet 1.0 file has no covering declaration, so gpd.read_parquet(..., bbox=...) has nothing to prune with and quietly falls back to reading everything and filtering in memory — correct results, none of the savings. Write the covering column and the same call reads a fraction of the file.
Where Shapefile still wins, and where GeoPackage fits
Two honest concessions. Shapefile is still the only vector format that literally every desktop GIS, survey instrument, and municipal permitting portal accepts without argument, so it remains a legitimate export target. And GeoParquet is not universally readable yet: support depends on the GDAL build, not just the GDAL version — the Parquet driver is optional and many Linux distribution packages ship without it.
# Does this GDAL actually have the Parquet driver compiled in?
ogrinfo --formats | grep -i parquet # (Parquet) Geoparquet -raw- (rw+v)
If that prints nothing, GDAL was built without Arrow support and QGIS on that machine cannot open a .parquet regardless of its version number. When you need a single file that everything reads, has no 2 GB ceiling, keeps the CRS internally, supports multiple layers and mixed geometry types, and stores real NULLs, the answer is GeoPackage — not Shapefile. What it does not give you is columnar column pruning or row-group skipping, so it is an interchange format rather than an analytical one.
Prerequisites
geopandas>=1.0— forwrite_covering_bbox,bbox=filtering on read, and GeoArrow encoding;0.14reads and writes GeoParquet 1.0 but has none of thosepyarrow>=15— the GeoParquet engine that backsto_parquet/read_parquet
conda install -c conda-forge "geopandas=1.0.*" "pyarrow=15.*"
Install from conda-forge, not pip, so GDAL, PROJ, and the bindings stay ABI-compatible — the same environment discipline used across the Shapefile & GeoJSON Parsing stage. If you are pinned to geopandas 0.14 for other reasons, everything below still works except the covering-bbox write and the bbox= read; substitute a DuckDB query for windowed reads in the meantime.
Step-by-Step Implementation
1. Read a legacy Shapefile and inspect the damage.
import geopandas as gpd
legacy = gpd.read_file("census_legacy.shp")
print(legacy.columns.tolist()) # ['populatio', 'median_inc', 'geometry'] — truncated
print(legacy.crs) # None if the .prj was missing
print(legacy.dtypes["median_inc"]) # object — DBF coerced the numeric column to string
2. Repair names, dtypes, and CRS, then write GeoParquet. Assert the known source CRS explicitly rather than guessing — here EPSG:25832 (ETRS89 / UTM 32N), a metric projected CRS suited to area work, not Web Mercator.
legacy = legacy.rename(columns={"populatio": "population", "median_inc": "median_income"})
legacy["median_income"] = legacy["median_income"].astype("float64") # undo DBF string coercion
if legacy.crs is None:
legacy = legacy.set_crs(epsg=25832) # set_crs asserts; it does NOT reproject
legacy.to_parquet(
"census.parquet",
compression="zstd", # smaller than the snappy default at similar decode speed
) # full names, real dtypes, and CRS now embedded in one file
3. Sort spatially before writing so row-group skipping actually helps. Row-group bbox pruning only pays off when nearby features sit near each other in the file; write in insertion order and every group's bbox spans the whole extent.
census = gpd.read_parquet("census.parquet")
# Hilbert-curve ordering groups spatially-adjacent rows into the same row group
census = census.sort_values(
by="geometry",
key=lambda geom: geom.hilbert_distance(total_bounds=census.total_bounds),
)
census.to_parquet("census.parquet", compression="zstd", row_group_size=50_000)
4. Write the covering bounding-box column and size the row groups deliberately. Sorting alone is not enough — the file has to declare where its envelope columns are, and the row group has to be small enough that skipping one saves real bytes but large enough that compression and metadata overhead stay reasonable.
census.to_parquet(
"census.parquet",
compression="zstd",
row_group_size=50_000, # ~64-128 MB per group for a polygon layer
write_covering_bbox=True, # geopandas>=1.0 — adds the bbox struct + covering metadata
schema_version="1.1.0", # emit GeoParquet 1.1 metadata
)
# Now a window read prunes instead of scanning
turin = gpd.read_parquet("census.parquet", bbox=(7.60, 45.00, 7.80, 45.10))
print(len(turin), "features intersect the window")
Row-group size is a genuine trade-off rather than a tuning ritual. Very small groups (a few thousand rows) give fine-grained skipping but multiply footer metadata, and the footer is read in full on every open — a file with tens of thousands of row groups can spend more time parsing statistics than reading data. Very large groups compress better and keep the footer tiny, but a single overlapping feature forces the whole group to decode. For polygon layers, aim for groups in the 64–128 MB range; for dense point data, larger groups are usually fine because the per-row cost is so low.
5. Compare on-disk size and round-trip fidelity.
import os
shp_bytes = sum(
os.path.getsize(f"census_legacy{ext}") for ext in (".shp", ".shx", ".dbf", ".prj")
)
pq_bytes = os.path.getsize("census.parquet")
print(f"Shapefile set: {shp_bytes/1e6:.1f} MB GeoParquet: {pq_bytes/1e6:.1f} MB")
# Shapefile set: 88.4 MB GeoParquet: 19.7 MB
6. Measure the read, not just the file size. Size is the headline, but the reason to migrate is the read path. Time the three access patterns you actually use — full load, column subset, and spatial window — because they diverge sharply.
import time
import geopandas as gpd
def timed(label, fn):
t0 = time.perf_counter()
out = fn()
print(f"{label:<28} {time.perf_counter() - t0:6.2f}s {len(out):>8,} rows")
timed("Shapefile, full load", lambda: gpd.read_file("census_legacy.shp"))
timed("GeoParquet, full load", lambda: gpd.read_parquet("census.parquet"))
timed("GeoParquet, 2 columns", lambda: gpd.read_parquet(
"census.parquet", columns=["population", "geometry"]))
timed("GeoParquet, bbox window", lambda: gpd.read_parquet(
"census.parquet", bbox=(7.60, 45.00, 7.80, 45.10)))
# Shapefile, full load 31.80s 1,200,000 rows
# GeoParquet, full load 4.10s 1,200,000 rows
# GeoParquet, 2 columns 1.60s 1,200,000 rows
# GeoParquet, bbox window 0.30s 14,802 rows
The full-load gap is mostly decompression versus DBF parsing and is roughly constant. The interesting numbers are the last two: column pruning and row-group skipping are structural wins a Shapefile cannot offer at any file size, because a .shp/.dbf pair has no internal index and no column independence — every read is a full sequential scan. That is also why the gap widens with dataset size rather than staying proportional.
7. Export Shapefile only when a legacy consumer demands it. Accept the truncation knowingly, and keep GeoParquet as the source of truth.
gpd.read_parquet("census.parquet").to_file("for_legacy_tool.shp") # names re-truncated to 10 chars
For a consumer that needs one portable file but not specifically a Shapefile, prefer GeoPackage — to_file("census.gpkg", driver="GPKG") — which keeps full field names, real dtypes, the CRS, and mixed geometry types, and has no size ceiling.
Verification
Confirm GeoParquet preserved the full schema, dtypes, and CRS that Shapefile lost — and that the file carries row-group statistics for skipping.
import geopandas as gpd
import pyarrow.parquet as pq
restored = gpd.read_parquet("census.parquet")
assert "population" in restored.columns and "median_income" in restored.columns
assert restored["median_income"].dtype == "float64" # not coerced back to string
assert restored.crs.to_epsg() == 25832
assert len(restored) == len(legacy)
print("CRS preserved:", restored.crs.to_epsg()) # CRS preserved: 25832
meta = pq.ParquetFile("census.parquet").metadata
print("Row groups:", meta.num_row_groups) # Row groups: 3 — enables bbox skipping
Then confirm the file actually declares a covering column, because this is the check that separates "GeoParquet that skips" from "GeoParquet that scans":
import json
import pyarrow.parquet as pq
geo = json.loads(pq.ParquetFile("census.parquet").schema_arrow.metadata[b"geo"])
print(geo["version"]) # 1.1.0
print(geo["columns"]["geometry"]["covering"]["bbox"])
# {'xmin': ['bbox', 'xmin'], 'ymin': ['bbox', 'ymin'],
# 'xmax': ['bbox', 'xmax'], 'ymax': ['bbox', 'ymax']}
print(geo["columns"]["geometry"]["crs"]["id"]) # {'authority': 'EPSG', 'code': 25832}
A KeyError on covering means the file was written without write_covering_bbox=True (or by a 1.0-era writer), and every bbox= read against it will silently degrade to a full scan.
Edge Cases & Debugging
- CRS is
Noneafter reading a Shapefile. The.prjwas missing;set_crsto the known source EPSG before anything else. - Mangled column names. Shapefile truncation; rename explicitly — there is no way to recover the original names automatically.
to_parquetfails.pyarrownot installed or too old.- Other tools can't read GeoParquet. Older GIS software predates it; export a Shapefile or GeoPackage for those, keep GeoParquet internally.
- Categorical/datetime columns differ after round trip. Parquet preserves dtypes Shapefile coerces to strings — usually an improvement, but check downstream assumptions.
- Multi-layer needs. Shapefile is one layer per file; if you need many layers in one file, use GeoPackage, not Shapefile.
pyarrowversion skew between writer and reader. A file written by a newerpyarrowcan fail on an older reader; pin the same major version across the whole pipeline.- Row-group skipping still reads everything. The features were never spatially sorted (step 3) or the writer emitted no bbox statistics — Hilbert-sort before writing so nearby features share a group.
PolygonbecameMultiPolygonafter a Shapefile round trip. Expected: the format stores every polygon as multipart. Call.explode(index_parts=False)if single-part geometry is required downstream.- Numeric columns come back with
0where the source had blanks. DBF cannot storeNULL; the zeros are fabricated. Re-derive nulls from a companion flag column, or migrate before the information is lost. - QGIS or ArcGIS refuses to open the
.parquet. The bundled GDAL was compiled without the Arrow/Parquet driver. Check withogrinfo --formats | grep -i parqueton that machine; if it is missing, export GeoPackage for that consumer rather than chasing the version number. - Footer parsing dominates the read time. Too many tiny row groups — each carries per-column statistics that are read in full on open. Raise
row_group_sizeuntil the footer is a small fraction of the file. to_parquetrejectswrite_covering_bboxas an unknown argument. You are ongeopandas<1.0; upgrade, or generate the bbox columns manually and accept that only DuckDB, notread_parquet, will use them.- Diffs appear in files nobody edited. Ring-orientation normalisation on a Shapefile round trip. Keep the canonical copy in GeoParquet and regenerate exports rather than storing both under version control.
Frequently Asked Questions
Is GeoParquet a replacement for GeoPackage as well as for Shapefile? No — they solve different problems and it is worth keeping both. GeoParquet is columnar and optimised for analytical reads: column pruning, row-group skipping, cheap remote access. GeoPackage is a SQLite database optimised for editing, multiple layers in one file, and near-universal desktop support. Use GeoParquet as the analytical and interchange store between pipeline stages; use GeoPackage when a human is going to open the file in a desktop GIS and edit it.
Should I write WKB or GeoArrow geometry encoding? WKB unless you have measured a reason not to. GeoArrow avoids a serialize/deserialize step and can be meaningfully faster for point-heavy data, but it is a GeoParquet 1.1 feature and older readers — including some cloud query engines — will reject or misread the file. WKB is the interoperable default and the decode cost is rarely the bottleneck compared to network transfer. Revisit if you are handing billions of points to an Arrow-native consumer.
How do I store a dataset too large for one Parquet file?
Write a partitioned dataset — a directory of Parquet files under Hive-style region=…/year=… subdirectories — rather than one enormous file. Engines prune whole partitions from the directory names before they look at any footer, which is a cheaper filter than row-group statistics, and each part stays independently readable and rewritable. Partition on the column you filter by most often, and keep parts in the hundreds of megabytes; thousands of tiny parts reintroduce the small-file problem the format was meant to solve.
Does GeoParquet preserve geometry precision exactly? Yes. WKB stores IEEE 754 doubles and Parquet compresses them losslessly, so coordinates round-trip bit-for-bit. Shapefile also stores doubles, so precision is not where it loses — the losses are in field names, dtypes, nulls, and the CRS. If you want smaller files by reducing precision, do it deliberately with a simplification or coordinate-rounding step before writing, not by hoping the format does it for you.
Can I append to a GeoParquet file? Not to a single file — Parquet's footer is written last and describes the whole file, so appending means rewriting. This is the practical argument for partitioned datasets: appending becomes writing a new part file into the directory, which is atomic and cheap. If you need row-level updates rather than appends, you have outgrown a file format and want a database, which is the trade-off examined in DuckDB Spatial vs PostGIS for Analytics.
What should I do with the Shapefiles I already have? Convert once, verify with the assertions above, and keep the originals archived rather than in the working path. The conversion is where you get to fix what the format broke — restoring field names, re-typing columns, asserting the CRS — and doing it as a one-off scripted migration means the repairs are recorded and repeatable. Leaving Shapefiles in the pipeline means paying the truncation and CRS-loss tax on every future read, and the bulk cleanup patterns for that are in automating Shapefile cleanup with Python.