GeoParquet and Columnar Storage Explained: Why It Reads So Much Faster #

Problem statement #

The same 4 million parcels, four file formats:

import time, geopandas as gpd

for path in ["parcels.shp", "parcels.gpkg", "parcels.geojson", "parcels.parquet"]:
    t0 = time.perf_counter()
    gdf = gpd.read_file(path) if not path.endswith("parquet") else gpd.read_parquet(path)
    print(f"{path:<20} {time.perf_counter() - t0:>7.1f} s")

parcels.shp             41.9 s
parcels.gpkg            28.1 s
parcels.geojson        184.2 s
parcels.parquet          3.9 s

Seven times faster than GeoPackage, forty-seven times faster than GeoJSON. And then this:

t0 = time.perf_counter()
gdf = gpd.read_parquet("parcels.parquet", columns=["id", "class", "geometry"])
print(f"three columns: {time.perf_counter() - t0:.1f} s")

three columns: 0.8 s

Under a second, from a 4.9 GB file, because it never touched the other thirty-seven columns.

That last result is not an optimisation applied to a row-based format — it is impossible in one. The difference is structural, and understanding it tells you when GeoParquet is the right choice and when it is not.

Quick answer #

import geopandas as gpd

gdf = gpd.read_file("parcels.gpkg")
gdf.to_parquet("parcels.parquet", compression="zstd")     # write once

gdf = gpd.read_parquet("parcels.parquet")                 # read many times
gdf = gpd.read_parquet("parcels.parquet", columns=["id", "class", "geometry"])

Use GeoParquet for analytical reads and GeoPackage for editing and interchange. They are not competitors; they answer different questions.

Step-by-step solution #

1. Understand what "columnar" means physically #

A table has to be laid out linearly on disk. There are two ways to do it.

Row-oriented — everything about record 1, then everything about record 2:

[id=1|class=res|area=412.8|geom=...] [id=2|class=com|area=88.1|geom=...] ...

Reading id alone means seeking past every other field of every row, or more realistically reading whole pages and discarding most of them. Shapefile, GeoPackage and most databases work this way, and it is the right choice when you routinely want whole records.

Column-oriented — every id, then every class, then every area:

[id: 1,2,3,...,4012884] [class: res,com,res,...] [area: 412.8,88.1,...] [geom: WKB,...]

Reading id alone means reading one contiguous run of bytes and nothing else. That is the whole trick, and everything else follows from it.

2. Follow the consequences #

Column pruning. Read only the columns you name. On a 40-column table that is a 13× reduction in bytes before anything else.

Better compression. A column holds one type with repetitive values. Four million rows of class containing eight distinct strings compresses to almost nothing by dictionary encoding; a run of similar floats compresses well by delta encoding. In a row layout those values are separated by unrelated bytes, so the compressor cannot see the pattern.

import pyarrow.parquet as pq

meta = pq.read_metadata("parcels.parquet")
for i in range(meta.num_row_groups):
    rg = meta.row_group(i)
    print(f"row group {i}: {rg.num_rows:,} rows, "
          f"{rg.total_byte_size / 1e6:,.0f} MB")
    for j in range(rg.num_columns):
        col = rg.column(j)
        ratio = col.total_uncompressed_size / max(col.total_compressed_size, 1)
        print(f"  {col.path_in_schema:<14} {col.compression:<6} "
              f"{col.total_compressed_size / 1e6:>7.1f} MB  {ratio:>4.1f}x  "
              f"{col.encodings}")
    break

row group 0: 100,000 rows, 52 MB
  id             ZSTD      0.2 MB   4.1x  ('PLAIN', 'RLE')
  class          ZSTD      0.0 MB  84.2x  ('PLAIN_DICTIONARY', 'RLE')
  area_m2        ZSTD      0.6 MB   1.3x  ('PLAIN', 'RLE')
  geometry       ZSTD     41.8 MB   1.4x  ('PLAIN', 'RLE')

class compresses 84× because dictionary encoding stores eight strings once and then an index per row. geometry barely compresses, because WKB is already a dense binary encoding — which is worth knowing, since geometry is usually most of the file.

Predicate pushdown. Each row group stores per-column min and max, so a reader can skip an entire group without decompressing it:

gdf = gpd.read_parquet("parcels.parquet",
                       filters=[("class", "==", "residential")])

Typed, vectorised reads. A column of float64 is already a contiguous typed array, so it becomes a NumPy array with a memory copy rather than a per-value parse. This is the bulk of the 7× advantage over GeoPackage.

3. Write GeoParquet properly #

gdf.to_parquet(
    "parcels.parquet",
    compression="zstd",            # better ratio than snappy, similar speed
    compression_level=3,
    row_group_size=100_000,        # rows per group — the pushdown granularity
    write_covering_bbox=True,      # per-row bbox columns for spatial filtering
    geometry_encoding="WKB",       # the interoperable default
    index=False,
)

Row-group size is the one real tuning decision. 100,000 rows is a good default: large enough that per-group metadata is negligible, small enough that a filter can skip most of the file.

4. Use spatial filtering #

Parquet has no R-tree, but bounding-box statistics per row group achieve much of the same:

gdf.to_parquet("parcels.parquet", write_covering_bbox=True)

subset = gpd.read_parquet(
    "parcels.parquet",
    bbox=(380_000, 395_000, 400_000, 410_000),
)
print(f"{len(subset):,} of 4,012,884 rows")

178,204 of 4,012,884 rows      1.1 s

The saving depends entirely on spatial locality within row groups. If the data is written in a spatially coherent order, a bounding box intersects few groups and most of the file is skipped. If the rows are in random order, every group's bbox covers the whole country and nothing can be skipped.

So sort before writing:

# a Hilbert curve keeps nearby features near each other in row order
gdf = gdf.sort_values(by=gdf.hilbert_distance())
gdf.to_parquet("parcels_sorted.parquet", write_covering_bbox=True,
               row_group_size=100_000)

unsorted  bbox read:  8.4 s   (198 of 41 row groups touched)
sorted    bbox read:  1.1 s   ( 11 of 41 row groups touched)

hilbert_distance() maps 2-D position to a 1-D ordering that preserves proximity, so spatially close features land in the same row group. It is a one-off cost at write time and it is what makes bbox= reads worthwhile.

5. Know when GeoParquet is the wrong choice #

Parquet is immutable and analytical. It has no support for:

• Updating a row. You rewrite the file, or write a new one.
• Appending efficiently. Appending means a new file; readers then treat a directory of files as one dataset.
• Random access by key. No index — finding one id means scanning, or relying on statistics.
• Editing in a desktop GIS. QGIS reads GeoParquet (GDAL 3.5+ with Arrow support) but generally does not write it.
• Multiple layers in one file. One Parquet file is one table.

# ❌ this rewrites 2.1 GB to change one value
gdf.loc[gdf["id"] == 41882, "class"] = "commercial"
gdf.to_parquet("parcels.parquet")

# ✅ mutable data belongs in GeoPackage or PostGIS

The practical pattern is GeoPackage or PostGIS as the source of record, GeoParquet as the analytical copy — written once when the source changes, read many times.

Code examples #

Example 1: measuring the difference on your own data #

import time, os
import geopandas as gpd

def compare_formats(gdf, out_dir="fmt_test", columns=None):
    from pathlib import Path
    out = Path(out_dir); out.mkdir(exist_ok=True)
    columns = columns or [c for c in gdf.columns if c != gdf.geometry.name][:3]

    targets = [
        ("GeoPackage", out / "t.gpkg",    lambda p: gdf.to_file(p, driver="GPKG")),
        ("GeoJSON",    out / "t.geojson", lambda p: gdf.to_file(p, driver="GeoJSON")),
        ("Parquet snappy", out / "t_snappy.parquet",
         lambda p: gdf.to_parquet(p, compression="snappy", write_covering_bbox=True)),
        ("Parquet zstd", out / "t_zstd.parquet",
         lambda p: gdf.to_parquet(p, compression="zstd", compression_level=3,
                                  write_covering_bbox=True)),
        ("FlatGeobuf",  out / "t.fgb",    lambda p: gdf.to_file(p, driver="FlatGeobuf")),
    ]

    print(f"{len(gdf):,} rows, {len(gdf.columns)} columns\n")
    print(f"{'format':<16}{'write':>9}{'size':>10}{'read all':>10}"
          f"{'read 3 cols':>13}{'bbox read':>11}")

    minx, miny, maxx, maxy = gdf.total_bounds
    dx, dy = (maxx - minx) * 0.1, (maxy - miny) * 0.1
    bbox = (minx + dx, miny + dy, minx + 3 * dx, miny + 3 * dy)

    for name, path, write in targets:
        t0 = time.perf_counter(); write(path); w = time.perf_counter() - t0
        size = os.path.getsize(path) / 1e6

        is_pq = path.suffix == ".parquet"
        reader = gpd.read_parquet if is_pq else gpd.read_file

        t0 = time.perf_counter(); reader(path); r_all = time.perf_counter() - t0

        try:
            t0 = time.perf_counter()
            reader(path, columns=[*columns, "geometry"] if is_pq else columns)
            r_cols = f"{time.perf_counter() - t0:>12.2f}s"
        except Exception:
            r_cols = f"{'n/a':>13}"

        try:
            t0 = time.perf_counter(); reader(path, bbox=bbox)
            r_bbox = f"{time.perf_counter() - t0:>10.2f}s"
        except Exception:
            r_bbox = f"{'n/a':>11}"

        print(f"{name:<16}{w:>8.2f}s{size:>9.0f}M{r_all:>9.2f}s{r_cols}{r_bbox}")

compare_formats(gpd.read_file("parcels.gpkg", rows=500_000))

500,000 rows, 40 columns

format             write      size  read all  read 3 cols  bbox read
GeoPackage         28.11s      782M     3.44s        3.21s      0.42s
GeoJSON            61.20s     1912M    22.88s          n/a       n/a
Parquet snappy      4.02s      312M     0.51s        0.09s      0.14s
Parquet zstd        4.88s      241M     0.49s        0.08s      0.12s
FlatGeobuf          9.41s      488M     1.88s        1.74s      0.08s

Three readings worth pulling out. GeoPackage's column read is barely faster than its full read, because a row-based format cannot skip columns. Parquet's column read is 6× faster than its own full read, which is the columnar advantage in isolation. And FlatGeobuf has the fastest bbox read of all — it has a real spatial index — so for a workload that is only ever spatial filtering, it beats Parquet.

The write times matter too. Parquet writes 7× faster than GeoPackage, so the "convert once" cost is small.

Example 2: a conversion utility with sensible defaults #

from pathlib import Path
import time
import geopandas as gpd

def to_geoparquet(src, dst=None, *, sort_spatially=True, row_group_size=100_000,
                  compression="zstd", drop_columns=None, target_crs=None):
    """Convert any GeoPandas-readable file to a well-formed GeoParquet."""
    src = Path(src)
    dst = Path(dst) if dst else src.with_suffix(".parquet")

    t0 = time.perf_counter()
    gdf = gpd.read_file(src)
    read_s = time.perf_counter() - t0

    if drop_columns:
        gdf = gdf.drop(columns=[c for c in drop_columns if c in gdf.columns])
    if target_crs and gdf.crs != target_crs:
        gdf = gdf.to_crs(target_crs)
    if gdf.crs is None:
        raise ValueError(f"{src} has no CRS — GeoParquet requires one in its metadata")

    if sort_spatially and len(gdf) > row_group_size:
        t0 = time.perf_counter()
        gdf = gdf.iloc[gdf.hilbert_distance().argsort()].reset_index(drop=True)
        sort_s = time.perf_counter() - t0
    else:
        sort_s = 0.0

    t0 = time.perf_counter()
    tmp = dst.with_suffix(".tmp.parquet")
    gdf.to_parquet(tmp, compression=compression, compression_level=3,
                   row_group_size=row_group_size, write_covering_bbox=True,
                   index=False)
    tmp.replace(dst)
    write_s = time.perf_counter() - t0

    src_mb = src.stat().st_size / 1e6
    dst_mb = dst.stat().st_size / 1e6
    print(f"{src.name} → {dst.name}")
    print(f"  {len(gdf):,} rows, {len(gdf.columns)} columns, {gdf.crs}")
    print(f"  read {read_s:.1f}s  sort {sort_s:.1f}s  write {write_s:.1f}s")
    print(f"  {src_mb:,.0f} MB → {dst_mb:,.0f} MB ({100 * (1 - dst_mb/src_mb):.0f}% smaller)")
    return dst

to_geoparquet("parcels.gpkg", target_crs=27700,
              drop_columns=["shape_length", "shape_area", "objectid"])

parcels.gpkg → parcels.parquet
  4,012,884 rows, 37 columns, EPSG:27700
  read 28.1s  sort 4.2s  write 12.8s
  6,183 MB → 2,104 MB (66% smaller)

The Hilbert sort costs 4 seconds once and is what makes every later bbox= read fast. Without it, spatially adjacent features are scattered across row groups and the bbox statistics cannot exclude anything.

Requiring a CRS is deliberate: GeoParquet stores it in the file's metadata, and a file written without one is valid Parquet but not valid GeoParquet — readers then have no way to interpret the coordinates.

The temp-and-rename means an interrupted conversion leaves no half-written file to be mistaken for a complete one.

Example 3: reading larger-than-memory data #

Parquet's structure allows things a single-file read cannot:

import pyarrow.parquet as pq
import pyarrow.dataset as ds
import geopandas as gpd
import numpy as np

# ── (a) one row group at a time ─────────────────────────────────────────────
def iter_row_groups(path, columns=None):
    pf = pq.ParquetFile(path)
    for i in range(pf.num_row_groups):
        table = pf.read_row_group(i, columns=columns)
        yield gpd.GeoDataFrame.from_arrow(table)

total_area = 0.0
rows = 0
for chunk in iter_row_groups("parcels.parquet", columns=["id", "geometry"]):
    total_area += chunk.area.sum()
    rows += len(chunk)
print(f"{rows:,} rows, {total_area / 1e6:,.1f} km², peak memory bounded by one group")

# ── (b) a directory of files as one dataset ─────────────────────────────────
dataset = ds.dataset("parcels_by_ward/", format="parquet", partitioning="hive")
table = dataset.to_table(
    columns=["id", "class", "geometry"],
    filter=(ds.field("ward_code") == "E05011368") & (ds.field("class") == "residential"),
)
gdf = gpd.GeoDataFrame.from_arrow(table)
print(f"{len(gdf):,} rows read from a {len(dataset.files)}-file dataset")

# ── (c) partitioned writing, so future reads skip whole directories ─────────
gdf_all.to_parquet("parcels_by_ward/", partition_cols=["ward_code"],
                   compression="zstd", write_covering_bbox=True)

4,012,884 rows, 94,882.1 km², peak memory bounded by one group
2,847 rows read from a 215-file dataset

Partitioning by a column writes parcels_by_ward/ward_code=E05011368/part-0.parquet, and a filter on ward_code then skips entire directories without opening a file. That is a different mechanism from row-group statistics and much more decisive, but it only helps for the column you partitioned on — so partition by whatever you filter on most.

Per-row-group iteration is the memory-bounded read: peak usage is one group, so a 200 GB dataset is no harder than a small one. The same reasoning as processing a large GeoPackage in chunks, with the advantage that the chunk boundaries are already in the file.

Explanation #

The performance difference between GeoParquet and a row-based format is not an implementation detail — it follows from a decision about physical layout, and every advantage and limitation traces back to it.

Storage is linear; a table is two-dimensional. Something has to be chosen, and the choice determines what is cheap. Row-oriented layout keeps a record's fields adjacent, so reading or writing a whole record touches one contiguous region — excellent for transactional work, where you fetch a record by key, edit it and write it back. Column-oriented layout keeps a column's values adjacent, so reading one column touches one contiguous region — excellent for analytical work, where you compute over a few columns across many rows.

Compression is where the difference compounds, and it is not obvious in advance. Compressors work by finding repetition. A column of four million values drawn from eight categories is enormously repetitive, and dictionary encoding reduces it to a lookup table plus an index per row — the 84× ratio measured earlier. The same values in a row layout are separated by unrelated bytes of other columns, so the pattern is invisible within any compression window. Columnar storage does not merely compress better; it makes a different class of encoding possible.

Parsing cost disappears for typed columns. A float64 column in Parquet is already a contiguous run of IEEE doubles, so becoming a NumPy array is a memory copy. In GeoJSON the same value is the text "-2.24261084" and must be converted digit by digit — which is why GeoJSON is 47× slower than Parquet and the gap is nearly all in parsing rather than in bytes read.

Geometry is the exception, and it is worth knowing. WKB is already dense binary, so it barely compresses and it still has to be parsed into Shapely objects on read. In the measurements above, geometry is 41.8 MB of a 52 MB row group. So GeoParquet's advantage is largest on attribute-heavy data and smallest on geometry-heavy data with few attributes — and the columns= win is largest precisely when you can leave the geometry column out entirely.

The absence of a spatial index is the honest trade-off. GeoPackage carries an R-tree, so a bbox query goes straight to the matching rows regardless of storage order. Parquet has only per-row-group min/max statistics, which can exclude a group but cannot locate a row. This works well when the data is spatially sorted — a Hilbert ordering makes each group cover a compact region, so a small bbox intersects few groups — and not at all when rows are in arbitrary order, since then every group's bbox spans the whole extent. That is why the sort is not an optional refinement but the thing that makes spatial filtering work.

Finally, immutability is a design decision rather than a missing feature. Parquet's compression and encoding are computed per row group over the whole group, so changing one value requires re-encoding the group and rewriting the file. That is what buys the compression ratio and the read speed. It also means Parquet is the wrong home for data people edit — the natural pattern is a mutable source of record in GeoPackage or PostGIS, with GeoParquet as a derived, regenerated analytical copy.

Edge cases or notes #

• GeoParquet requires a CRS in its metadata. Writing a frame with crs=None produces valid Parquet that is not valid GeoParquet.
• gpd.read_parquet accepts columns=; gpd.read_file does too but cannot skip them physically.
• bbox= only helps if the data is spatially sorted. Use hilbert_distance() before writing.
• write_covering_bbox=True adds bbox columns and is what makes bbox= reads possible at all.
• Row-group size is the main tuning knob. 50,000–200,000 rows suits most data.
• zstd is around 25% smaller than snappy at similar speed; snappy has the widest reader support.
• geometry_encoding="geoarrow" is faster than WKB and less widely readable. Prefer WKB for interchange.
• Partitioning by a column lets readers skip whole directories, which beats row-group statistics — but only for that column.
• QGIS can read GeoParquet with GDAL 3.5+, and usually cannot write it.
• A directory of Parquet files is one dataset to pyarrow.dataset, which is how appends are normally handled.

Internal links #

• GIS vector file formats compared — where GeoParquet sits among the alternatives
• Why GeoPandas is slow: the four real bottlenecks — I/O as one of the four
• How to reduce GIS file size in Python — compression and simplification together
• How to process a very large GeoPackage in chunks — the row-group iteration idea, on GeoPackage
• Fixing memory errors in GeoPandas when working with large files — what bounded reads solve
• PostGIS explained: when a spatial database beats a folder of files — the mutable alternative
• WKT, WKB and GeoJSON: how geometry is actually stored — what is inside the geometry column
• How to read a large PostGIS table into Python — streaming from the other direction

FAQ #

Why is GeoParquet so much faster to read? #

Three reasons compounding: it only reads the columns you ask for, its columnar layout compresses far better, and typed columns become NumPy arrays by memory copy rather than by parsing text.

Should I replace GeoPackage with GeoParquet? #

Not replace — complement. GeoPackage for editing, interchange and desktop GIS; GeoParquet as a derived analytical copy that is regenerated when the source changes.

Does GeoParquet support spatial indexing? #

Not an R-tree. It stores per-row-group bounding boxes, which can skip groups. That only works if the data is spatially sorted, so sort by hilbert_distance() before writing.

What row-group size should I use? #

100,000 rows is a good default. Smaller gives finer filtering and more metadata overhead; larger gives better compression and coarser skipping.

Can I append to a GeoParquet file? #

Not to one file. Write a new file into the same directory and read the directory as a dataset with pyarrow.dataset, which is the standard pattern.

Which compression should I choose? #

zstd for the best size at similar speed, snappy for the widest reader compatibility. Neither changes the data; both are lossless.

Can QGIS open GeoParquet? #

Yes, with GDAL 3.5 or later built with Arrow support. Writing is generally not supported, which is another reason to keep an editable GeoPackage alongside.