GIS Vector File Formats Compared: Shapefile, GeoPackage, GeoJSON, Parquet
On this page
Problem statement
The same 338,611 land parcels β HM Land Registry's INSPIRE polygons for Birmingham β written five ways, with the size on disk and the time to write and read each back:
print(compare_formats(parcels).to_string(index=False))
format MB write_s read_s columns_unchanged
geoparquet 45.62 0.9 0.5 7/7
flatgeobuf 177.80 2.9 1.3 7/7
shapefile 179.76 3.8 1.0 4/7
geopackage 193.97 2.8 1.5 7/7
geojson 389.75 12.7 5.5 7/7
compare_formats() is defined in Measure it on your own data below; the run's load average is printed with the data. GeoParquet is a quarter of the size of the next smallest and the fastest to write and read. GeoJSON is twice the size of anything else and the slowest. The shapefile comes out smaller than the GeoPackage here β and three of its seven columns do not survive unchanged: two names are cut to ten characters and a timestamp comes back as text.
Size and speed are the visible differences. The ones that cost more are what each format can carry besides geometry β field names, types, layers, an index, a CRS β because those decide what comes back unchanged.
Quick answer
Pick by the constraint that actually binds:
- GeoPackage β the sensible default for interchange: one file, many layers, UTF-8, real types
- GeoParquet β the default inside a pipeline: smallest, fastest, preserves pandas dtypes
- GeoJSON β for the web and APIs, in WGS 84, expect it to be large
- FlatGeobuf β streaming and cloud-native reads with a built-in spatial index
- Shapefile β only when a recipient requires it; expect renamed columns and converted types
parcels.to_file("out/parcels.gpkg", layer="parcels", driver="GPKG") # interchange
parcels.to_parquet("out/parcels.parquet", compression="zstd", write_covering_bbox=True) # pipeline
parcels.to_crs(4326).to_file("out/parcels.geojson", driver="GeoJSON",
RFC7946="YES", COORDINATE_PRECISION=6) # web
parcels.to_file("out/parcels.fgb", driver="FlatGeobuf") # streaming
print(sorted(f.name for f in Path("out").iterdir()))
['compare', 'parcels.fgb', 'parcels.geojson', 'parcels.gpkg', 'parcels.parquet']
The one habit worth adopting: keep a full-fidelity master in GeoPackage or GeoParquet, and derive every other product from it rather than converting between products.
The five formats at a glance
Each format, and what it costs you
Shapefile: what it costs you
It is over thirty years old, it is everywhere, and every one of its limits is a real-world bug. Round-tripping the parcels:
back = gpd.read_file("out/compare/shapefile/d.shp")
print("renamed:", sorted(set(parcels.columns) - set(back.columns)), "->", sorted(set(back.columns) - set(parcels.columns)))
print({c: f"{parcels[c].dtype} -> {back[c].dtype}" for c in ["VALIDFROM", "BEGINLIFESPANVERSION"] if c in back}
or {c: str(back[c].dtype) for c in back.columns if c.startswith(("VALID", "BEGINL"))})
renamed: ['BEGINLIFESPANVERSION', 'NATIONALCADASTRALREFERENCE'] -> ['BEGINLIFES', 'NATIONALCA']
{'VALIDFROM': 'datetime64[ms, UTC] -> str'}
The list of limits worth memorising, each checked with GDAL 3.12:
- field names: 10 characters β
NATIONALCADASTRALREFERENCEbecomesNATIONALCA; of two long names sharing a start,population_2025becomespopulationandpopulation_2026becomespopulati_1 - text fields: 254 characters β longer values are cut, with one warning per layer
- no datetime field β GDAL writes datetimes as a text field, with a warning that "DateTime" was requested; they read back as strings
- booleans β stored as a one-character integer field; GDAL reads them back as booleans, other software sees 1 and 0
- no 64-bit integer precision in practice β large integers are written as doubles
- one geometry type per file β polygons and lines need separate files
- 2 GB limit per
.shpand per.dbf - encoding is a sidecar β a missing
.cpgmeans GDAL reads the text as Latin-1 - several files that must travel together
Use it when a recipient's software requires it, rename the fields yourself before writing, and keep the real version elsewhere.
GeoPackage: the sensible default
import pyogrio
wards.to_file("out/atlas.gpkg", layer="wards", driver="GPKG")
parcels.to_file("out/atlas.gpkg", layer="parcels", driver="GPKG")
print(pyogrio.list_layers("out/atlas.gpkg"))
print({c: str(t) for c, t in gpd.read_file("out/atlas.gpkg", layer="parcels", rows=1).dtypes.items()})
[['wards' 'Polygon']
['parcels' 'Polygon']]
{'gml_id': 'str', 'INSPIREID': 'int32', 'LABEL': 'int32', 'NATIONALCADASTRALREFERENCE': 'int32', 'VALIDFROM': 'datetime64[ms, UTC]', 'BEGINLIFESPANVERSION': 'datetime64[ms, UTC]', 'geometry': 'geometry'}
GeoPackage is an SQLite database with an OGC-standardised schema. That gives you a single file, UTF-8 throughout, no practical limit on field-name or text length, real date and datetime types (the timezone survives), multiple layers of different geometry types, an R-tree spatial index, and room for rasters and metadata tables. QGIS, ArcGIS, GDAL and anything else built in the last decade reads it.
Its weaknesses are SQLite's: one writer at a time (a second writer waits or fails with "database is locked"), and locking that is unreliable on some network shares.
GeoJSON: for the web, with rules
RFC 7946 says GeoJSON coordinates are WGS 84 longitude and latitude. GDAL does not enforce that by default β it will write projected coordinates and label them with the pre-standard crs member:
wards.to_file("out/wards_bng.geojson", driver="GeoJSON")
print(open("out/wards_bng.geojson").read(200).split('"features"')[0])
wards.to_crs(4326).to_file("out/wards.geojson", driver="GeoJSON", RFC7946="YES", COORDINATE_PRECISION=6)
print(gpd.read_file("out/wards.geojson").crs, f"{os.path.getsize('out/wards_bng.geojson') / 1e3:,.0f} kB ->",
f"{os.path.getsize('out/wards.geojson') / 1e3:,.0f} kB")
{
"type": "FeatureCollection",
"name": "wards_bng",
"crs": { "type": "name", "properties": { "name": "urn:ogc:def:crs:EPSG::27700" } },
EPSG:4326 167 kB -> 99 kB
A web client following the standard ignores the crs member and puts the metre coordinates somewhere near 0Β°, 0Β°. Write RFC7946="YES" from WGS 84. GeoJSON is plain JSON, so every language and every browser reads it without a GIS library; the costs are size β text with repeated keys, twice the size of any other format for the parcels β and the WGS 84 rule. For line-delimited streaming, GeoJSONSeq (.geojsonl) writes one feature per line.
FlatGeobuf: streaming and cloud-native
FlatGeobuf is a flat binary sequence of features with an optional packed Hilbert R-tree at the front. A reader can fetch the index, work out which byte ranges it needs, and request only those β over HTTP, from a static file:
import json
import subprocess
import sys
url = "/vsicurl/https://flatgeobuf.org/test/data/UScounties.fgb"
probe = f'import geopandas as gpd; print(len(gpd.read_file("{url}", bbox=(-105.3, 39.6, -104.6, 40.1))))'
run = subprocess.run([sys.executable, "-c", probe], capture_output=True, text=True,
env=dict(os.environ, CPL_VSIL_SHOW_NETWORK_STATS="YES"))
output = run.stdout + run.stderr
stats = json.loads(output[output.index("{"):])["methods"]["GET"]
print(output.split()[0], "counties;", stats["count"], "range requests,", f"{stats['downloaded_bytes']:,} bytes of a 14,100,008-byte file")
7 counties; 7 range requests, 196,608 bytes of a 14,100,008-byte file
That is the best way to serve large vector data without a tile server or database. FlatGeobuf is write-once: there is no editing a feature in place.
GeoParquet: the pipeline format
slim = gpd.read_parquet("out/parcels.parquet", columns=["INSPIREID", "geometry"])
print(slim.shape, {c: str(t) for c, t in gpd.read_parquet("out/parcels.parquet").dtypes.items() if c != "geometry"})
(338611, 2) {'gml_id': 'str', 'INSPIREID': 'int32', 'LABEL': 'int32', 'NATIONALCADASTRALREFERENCE': 'int32', 'VALIDFROM': 'datetime64[ms, UTC]', 'BEGINLIFESPANVERSION': 'datetime64[ms, UTC]'}
Parquet stores columns rather than rows, compresses each column with a codec suited to its type, and records the schema β including pandas extension dtypes and timezone-aware timestamps β in the file. GeoParquet adds a standard metadata block describing the geometry column and its CRS. Here it was the smallest file, the fastest read, and it returned every column with its original type.
The limitation is ecosystem: QGIS reads GeoParquet only when its GDAL was built with Parquet support, and the GDAL inside the pyogrio wheels has no Parquet driver at all β GeoPandas reads Parquet through pyarrow instead. It is a pipeline and analytics format, not a delivery format.
Measure it on your own data
The comparison at the top of the page is compare_formats(), defined in the setup block below: it writes a GeoDataFrame in all five formats, reads each back, and reports size, time and how many columns kept their name and type. Run it on your own data, because the ratios depend on it. Birmingham's 69 electoral wards β few features, detailed outlines, two short attributes β come out differently:
print(compare_formats(wards, out=Path("out/compare_wards")).to_string(index=False))
format MB write_s read_s columns_unchanged
geoparquet 0.05 0.0 0.0 3/3
shapefile 0.07 0.0 0.0 3/3
flatgeobuf 0.07 0.0 0.0 3/3
geojson 0.16 0.0 0.0 3/3
geopackage 0.17 0.0 0.0 3/3
Numbers from your own data beat any table in an article. For 69 features the GeoPackage is the largest file, because an SQLite database with its metadata tables and spatial index costs close to 100 kB before any data; for 338,611 parcels that overhead disappears into the total. GeoParquet is the smallest either way.
Birmingham's parcels in five formats
The data used on this page
- Parcels β HM Land Registry's INSPIRE Index Polygons for Birmingham, 338,611 freehold parcels in British National Grid, with an ID, two copies of it and two timestamps. Downloaded as GML and saved as a GeoPackage once. The Land Registry republishes the file monthly, so the count changes slightly.
- Wards β Birmingham's electoral wards (December 2024, generalised) from the ONS boundary service.
The block also defines compare_formats().
import io
import os
import time
import warnings
import zipfile
from pathlib import Path
import geopandas as gpd
import pandas as pd
import requests
os.makedirs("data/raw", exist_ok=True)
Path("out").mkdir(exist_ok=True)
if not os.path.exists("data/raw/parcels.gpkg"):
reply = requests.Session().get("https://use-land-property-data.service.gov.uk/datasets/"
"inspire/download/Birmingham.zip", timeout=600)
reply.raise_for_status()
zipfile.ZipFile(io.BytesIO(reply.content)).extract("Land_Registry_Cadastral_Parcels.gml")
gpd.read_file("Land_Registry_Cadastral_Parcels.gml").to_file("data/raw/parcels.gpkg", layer="parcels")
os.remove("Land_Registry_Cadastral_Parcels.gml")
if not os.path.exists("data/raw/wards.gpkg"):
reply = requests.get("https://services1.arcgis.com/ESMARspQHYMw9BZ9/arcgis/rest/services/"
"Wards_December_2024_Boundaries_UK_BGC/FeatureServer/0/query",
params={"where": "LAD24NM='Birmingham'", "outFields": "WD24CD,WD24NM",
"outSR": 27700, "f": "geojson"}, timeout=120)
reply.raise_for_status()
gpd.read_file(io.BytesIO(reply.content)).to_file("data/raw/wards.gpkg")
parcels = gpd.read_file("data/raw/parcels.gpkg")
wards = gpd.read_file("data/raw/wards.gpkg")
print(len(parcels), "parcels,", len(wards), "wards; load average", os.getloadavg()[0])
def compare_formats(gdf, out=Path("out/compare")):
writers = {
"shapefile": lambda p: gdf.to_file(p / "d.shp", driver="ESRI Shapefile"),
"geopackage": lambda p: gdf.to_file(p / "d.gpkg", driver="GPKG"),
"geojson": lambda p: gdf.to_crs(4326).to_file(p / "d.geojson", driver="GeoJSON"),
"flatgeobuf": lambda p: gdf.to_file(p / "d.fgb", driver="FlatGeobuf"),
"geoparquet": lambda p: gdf.to_parquet(p / "d.parquet", compression="zstd"),
}
readers = {
"shapefile": lambda p: gpd.read_file(p / "d.shp"),
"geopackage": lambda p: gpd.read_file(p / "d.gpkg"),
"geojson": lambda p: gpd.read_file(p / "d.geojson"),
"flatgeobuf": lambda p: gpd.read_file(p / "d.fgb"),
"geoparquet": lambda p: gpd.read_parquet(p / "d.parquet"),
}
rows = []
for name, write in writers.items():
folder = out / name
folder.mkdir(parents=True, exist_ok=True)
with warnings.catch_warnings():
warnings.simplefilter("ignore") # the shapefile's truncation warnings
start = time.perf_counter(); write(folder); write_s = time.perf_counter() - start
size = sum(f.stat().st_size for f in folder.iterdir())
start = time.perf_counter(); back = readers[name](folder); read_s = time.perf_counter() - start
same = [c for c in gdf.columns if c in back.columns and str(back[c].dtype) == str(gdf[c].dtype)]
rows.append({"format": name, "MB": round(size / 1e6, 2), "write_s": round(write_s, 1),
"read_s": round(read_s, 1), "columns_unchanged": f"{len(same)}/{len(gdf.columns)}"})
return pd.DataFrame(rows).sort_values("MB")
338611 parcels, 69 wards; load average 5.8740234375
Example 1: what survives a round trip
from shapely.geometry import Point
sample = gpd.GeoDataFrame({
"population_2025": [1, 2], "population_2026": [3, 4],
"text": ["a" * 300, "MΓΌnchen"],
"flag": pd.array([True, False], dtype="boolean"),
"when": pd.to_datetime(["2026-08-11 14:30", "2026-08-12 09:15"], utc=True),
"count": pd.array([1, None], dtype="Int64"),
}, geometry=[Point(407000, 287000), Point(408000, 288000)], crs=27700)
for name, path, write, read in [
("shapefile", "out/t.shp", lambda p: sample.to_file(p), gpd.read_file),
("geopackage", "out/t.gpkg", lambda p: sample.to_file(p, driver="GPKG"), gpd.read_file),
("geojson", "out/t.geojson", lambda p: sample.to_crs(4326).to_file(p, driver="GeoJSON"), gpd.read_file),
("geoparquet", "out/t.parquet", lambda p: sample.to_parquet(p), gpd.read_parquet),
]:
with warnings.catch_warnings():
warnings.simplefilter("ignore")
write(path)
back = read(path)
changed = {c: str(back[c].dtype) for c in back.columns if c != "geometry"
and (c not in sample or str(back[c].dtype) != str(sample[c].dtype))}
print(f"{name:10} text {back['text'].str.len().max()} chars | changed: {changed}")
shapefile text 254 chars | changed: {'population': 'int64', 'populati_1': 'int64', 'flag': 'bool', 'when': 'str', 'count': 'float64'}
geopackage text 300 chars | changed: {'flag': 'bool', 'when': 'datetime64[ms, UTC]', 'count': 'float64'}
geojson text 300 chars | changed: {'population_2025': 'int32', 'population_2026': 'int32', 'flag': 'bool', 'when': 'datetime64[ms, UTC]', 'count': 'float64'}
geoparquet text 300 chars | changed: {}
One script shows the whole comparison: the shapefile renames both population columns, cuts the text to 254 characters and turns the timestamp into text; the GeoPackage and GeoJSON keep names and timestamps but turn the nullable Int64 and boolean into plain NumPy types, and the GeoJSON even narrows the 64-bit integers to 32-bit because the values fit; GeoParquet changes nothing.
Example 2: a multi-layer GeoPackage as a deliverable
import pyogrio
dest = Path("out/delivery.gpkg")
dest.unlink(missing_ok=True)
wards.to_file(dest, layer="wards", driver="GPKG")
parcels.clip(wards[wards["WD24NM"] == "Ladywood"]).to_file(dest, layer="parcels_ladywood", driver="GPKG")
for name, geom_type in pyogrio.list_layers(dest):
info = pyogrio.read_info(dest, layer=name)
print(f"{name:18} {geom_type:14} {info['features']:>7,} features {info['crs']}")
wards Polygon 69 features EPSG:27700
parcels_ladywood MultiPolygon 2,622 features EPSG:27700
Example 3: partial reads
Reading only a 5 km square of Birmingham from each format:
BBOX = (405000, 285000, 410000, 290000)
reads = {
"geopackage": lambda: gpd.read_file("out/compare/geopackage/d.gpkg", bbox=BBOX),
"flatgeobuf": lambda: gpd.read_file("out/compare/flatgeobuf/d.fgb", bbox=BBOX),
"shapefile": lambda: gpd.read_file("out/compare/shapefile/d.shp", bbox=BBOX),
"geojson": lambda: gpd.read_file("out/compare/geojson/d.geojson", bbox=tuple(gpd.GeoSeries.from_xy(
[BBOX[0], BBOX[2]], [BBOX[1], BBOX[3]], crs=27700).to_crs(4326).total_bounds)),
"geoparquet + bbox column": lambda: gpd.read_parquet("out/parcels.parquet", bbox=BBOX),
}
for name, read in reads.items():
start = time.perf_counter(); part = read()
print(f"{name:25} {len(part):>7,} parcels in {time.perf_counter() - start:.2f} s")
try:
gpd.read_parquet("out/compare/geoparquet/d.parquet", bbox=BBOX)
except ValueError as err:
print("geoparquet, no bbox column:", str(err)[:75], "...")
geopackage 27,543 parcels in 0.11 s
flatgeobuf 27,543 parcels in 0.10 s
shapefile 27,543 parcels in 0.16 s
geojson 27,570 parcels in 4.79 s
geoparquet + bbox column 27,547 parcels in 0.19 s
geoparquet, no bbox column: Specifying 'bbox' not supported for this Parquet file (it should either hav ...
Every format can answer a bounding-box read, but not equally. GeoPackage and FlatGeobuf use a spatial index; a shapefile skips records using the bounding box stored in each record header; GeoJSON has to parse the whole file first, which is why it took about forty times as long. GeoParquet can only filter by bounding box if it was written with write_covering_bbox=True, and it filters on each feature's bounding box, so it can return a few features the others exclude. (The GeoJSON read uses a box in degrees, so its count differs slightly too.)
Example 4: choose the format from the audience
AUDIENCE = {
"analyst_python": {"format": "parquet", "crs": None, "note": "keeps dtypes exactly"},
"gis_desktop": {"format": "gpkg", "crs": 27700, "note": "opens everywhere"},
"web_map": {"format": "geojson", "crs": 4326, "note": "RFC7946=YES, simplify first"},
"large_web_layer": {"format": "fgb", "crs": 4326, "note": "bbox reads over HTTP"},
"legacy_client": {"format": "shp", "crs": 27700, "note": "rename fields first"},
}
for who, spec in AUDIENCE.items():
print(f"{who:16} -> .{spec['format']:8} {spec['note']}")
analyst_python -> .parquet keeps dtypes exactly
gis_desktop -> .gpkg opens everywhere
web_map -> .geojson RFC7946=YES, simplify first
large_web_layer -> .fgb bbox reads over HTTP
legacy_client -> .shp rename fields first
Same geometry, different containers
Every format here stores the same Simple Features geometry, usually as WKB inside; the measurements differ because the containers do. The shapefile's attribute table is the dBase format of the 1980s, with ten-character names and no timestamp type, which is why NATIONALCADASTRALREFERENCE came back as NATIONALCA and VALIDFROM as text: 4 of 7 columns unchanged. A GeoPackage is an SQLite database with a standard schema, so names, types and an R-tree index survive, at the price of being the largest binary file here.
The index decides partial reads. Reading the Ladywood ward's bounding box took a tenth of a second or so from GeoPackage, FlatGeobuf and the shapefile, which skips records by the box in each record header, and over four seconds from GeoJSON, which must parse every feature first (Example 3). FlatGeobuf goes further: with its index at the front of the file, a client fetched seven counties in 7 HTTP range requests, 196,608 bytes of a 14.1 MB file, which is what makes a static file behave like a feature service. GeoParquet's advantage is columnar storage, the subject of GeoParquet and columnar storage.
So there is no best format, only a best pairing of format and consumer β and a pipeline should keep one high-fidelity master and derive deliveries from it, because a chain of conversions through lossy formats is where names, types and precision quietly go.
Truncation, concurrency and other format traps
- Shapefile field names collide after truncation:
population_2025becomespopulationandpopulation_2026becomespopulati_1β GDAL cuts the second one shorter and appends a counter. Rename deliberately. - GeoPackage and concurrency: SQLite allows one writer at a time; others wait or fail with "database is locked". Do not write to one GeoPackage from parallel processes.
- GeoJSON must be WGS 84: RFC 7946 removed the
crsmember, but GDAL still writes one for a projected layer unless you passRFC7946="YES". Reproject before writing. .geojsonl/ GeoJSONSeq: One feature per line, streamable, far friendlier for very large exports than a single JSON document.- Parquet needs GDAL built with Parquet support: true of recent QGIS installers, not of the GDAL inside the pyogrio wheels. Check the recipient's setup before delivering it.
- FlatGeobuf is append-friendly, not edit-friendly: It is designed for write-once, read-many.
- KML, GML, MapInfo TAB: Still common in some sectors. GDAL reads them all β convert on ingest rather than working in them.
Internal links
- How to Read and Write GeoPackage Files in Python
- How to Batch Convert Shapefiles to GeoPackage in Python
- How to Reduce GIS File Size in Python Without Wrecking the Data
- Why Are My Shapefile Column Names Truncated? How to Fix It
- GeoPandas to_file() Fails on a Column Type: How to Fix It
- WKT, WKB and GeoJSON: How Geometry Is Actually Stored
- Metadata disappears when you convert the file β measured losses per format
FAQ
Which vector format should I use by default?
GeoPackage for anything you hand to another person, GeoParquet for anything that stays inside a Python pipeline. Both are lossless for attributes, which shapefile is not.
Is the shapefile really that bad?
It works, and it is universal. The problems are specific and predictable: 10-character field names, 254-character text, timestamps written as text, one geometry type per file, and several files that must stay together.
Why is my GeoJSON so large?
It is text: coordinates as decimal strings and key names repeated per feature. For the Birmingham parcels it was twice the size of any other format. Reduce COORDINATE_PRECISION (6 decimal places of a degree is about 0.1 m), simplify the geometry, and rely on gzip in transport β GeoJSON compresses very well.
When should I use FlatGeobuf?
When clients need to read part of a large dataset directly from static hosting. Its front-loaded spatial index lets a bbox query fetch only the relevant byte ranges, even over HTTP β seven range requests and under 200 kB of a 14 MB file for the Denver-area counties above.
Can QGIS open GeoParquet?
Versions whose GDAL was built with the Parquet driver can, which includes the current official installers. Check before delivering; GeoPackage remains the safe interchange choice.
How do I keep timezone-aware timestamps?
Use GeoParquet, which preserves the pandas dtype exactly. GeoPackage stores a datetime but may normalise the offset; shapefile drops the time entirely.
Does converting between formats lose data?
It can. Every conversion is limited by the target format, so a chain through shapefile permanently loses long field names and text. Keep a master in a high-fidelity format and derive from it.