1
1

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

@rino_yume(rino yume)

GeoParquetを使って大きなGISデータを扱う――GeoPandasとDuckDBで試す

1
Posted at

GeoParquetを使って大きなGISデータを扱う――GeoPandasとDuckDBで試す

GISデータの受け渡しでは、GeoPackageやShapefileを使うことが多いと思います。
編集や目視確認ならそれで困りませんが、建物ポリゴン、人口メッシュ、観測点などが増えてくると、別の選択肢も欲しくなります。

GeoParquetは、そうした大きめのベクトルデータを分析したり、クラウド上へ置いたりするときに使いやすい形式です。この記事では、仕組みを整理したうえで、GeoPandasとDuckDBを使った基本的な読み書きを試します。

この記事の内容は2026年8月3日時点で確認しています。GeoParquet 2.0が現行仕様ですが、サンプルコードはツール間の互換性を優先し、GeoParquet 1.1＋WKB＋bbox coveringで作成します。

GeoParquetとは

GeoParquet は、Apache Parquetへ地理空間データの扱い方を加えた形式です。
Point、LineString、Polygonなどのジオメトリに加えて、CRS、主ジオメトリ列、地物の種類、空間範囲などをメタデータとして持たせます。

Parquet自体は列指向のファイル形式です。CSVのように1行ずつ読むのではなく、必要な列だけを選んで読みやすい構造になっています。

GeoParquetファイル
├─ Row Group 1
│  ├─ id列
│  ├─ category列
│  ├─ value列
│  └─ geometry列
├─ Row Group 2
│  ├─ id列
│  ├─ category列
│  ├─ value列
│  └─ geometry列
└─ Footer
   ├─ スキーマ
   ├─ 列統計
   ├─ Row Groupの位置
   └─ GeoParquetメタデータ

Parquetでは、行をまとめた単位をRow Group、各Row Group内の列データをColumn Chunk、さらに細かな圧縮単位をPageと呼びます。

たとえば100列ある建物データから、建物ID、高さ、用途だけを集計する場合、geometryを含む残りの列を読まずに処理できます。DuckDBでは、この列選択が自動的にParquet Readerへ渡されます。

SELECT
    building_id,
    height,
    usage
FROM read_parquet('buildings.parquet');

GeoParquetのバージョン

2026年8月時点の現行仕様はGeoParquet 2.0です。

GeoParquet 2.0では、Apache Parquet 2.11で追加されたGEOMETRYとGEOGRAPHYの論理型を使います。ジオメトリはWKBで格納され、CRSや空間統計をParquet側の型情報として扱えるようになりました。

一方、GeoPandasの安定版APIでは、to_parquet()のschema_versionとしてGeoParquet 1.1までが案内されています。既存環境を含めて確実に読み書きしたい場合は、しばらくGeoParquet 1.1を使う場面も残ります。

項目 GeoParquet 1.1 GeoParquet 2.0
geometry WKBまたはGeoArrow ParquetのGEOMETRY / GEOGRAPHY論理型
CRS geoメタデータ内のPROJJSON 論理型のCRS＋GeoParquetメタデータ
空間絞り込み 行単位のbbox coveringを利用 geometry列の空間統計を利用
互換性 現行ツールで扱いやすい 対応が進んでいる途中
この記事のサンプル 使用する 概要のみ説明

新規システムでGeoParquet 2.0を採用する場合は、GeoPandas、PyArrow、GDAL、DuckDB、QGISなど、実際に使うツールの対応バージョンをそろえてから決めるのが安全です。

GeoPackageやGeoJSONとの違い

GeoParquetは便利ですが、すべてのGISファイルを置き換えるものではありません。

形式 向いている用途 補足
GeoParquet 大量データの分析、配布、クラウド保存 列選択、圧縮、分割読み込みに向く
GeoPackage QGISでの編集、複数レイヤ管理 SQLiteベースで更新しやすい
GeoJSON Web API、小規模データ交換 ブラウザから扱いやすいが大容量には重い
FlatGeobuf BBOX指定でのベクトル配信 空間インデックスを持つ単一レイヤ形式
Shapefile 古いGISソフトとの受け渡し 複数ファイル構成や文字数などの制約がある

実務では、次のように役割を分けると扱いやすくなります。

編集用       : GeoPackage
分析・保管用 : GeoParquet
Web配信用    : MVT / GeoJSON
検索・集計   : DuckDB / PostGIS

処理の流れ

ChatGPT Image 2026年8月3日 13_16_07.png

サンプルの構成

この記事で使うサンプル一式は、次の構成です。
ZIPファイルをダウンロードして、解凍してください。

geoparquet_sample_r001/
├─ README.md
├─ requirements.txt
├─ data/
│  ├─ input/
│  └─ output/
├─ src/
│  ├─ create_sample.py
│  ├─ convert_to_geoparquet.py
│  ├─ inspect_geoparquet.py
│  ├─ read_bbox.py
│  └─ query_duckdb.py
├─ sql/
│  └─ example_queries.sql
└─ tests/
   └─ test_workflow.py

外部データを用意しなくても確認できるように、静岡市周辺を想定した架空のPointデータを生成します。実在の施設や観測所を示すものではありません。

環境を作る

Python 3.13を使用します。

Windows PowerShell

py -3.13 -m venv .venv
.\.venv\Scripts\Activate.ps1
python -m pip install --upgrade pip
pip install -r requirements.txt

WSL / Ubuntu

python3.13 -m venv .venv
source .venv/bin/activate
python -m pip install --upgrade pip
pip install -r requirements.txt

requirements.txtは次の内容です。

geopandas>=1.1,<2.0
pyarrow>=18
shapely>=2.1
pandas>=2.2
numpy>=2.0
tzdata>=2025.2
duckdb>=1.5
pytest>=8.3

サンプルGeoParquetを作る

実行します。

python src/create_sample.py

出力先は次のファイルです。

data/output/sample_points.parquet

保存部分は次のようになっています。

from pathlib import Path

import geopandas as gpd


def write_geoparquet(
    gdf: gpd.GeoDataFrame,
    output_path: Path,
    row_group_size: int = 50_000,
    compression_level: int = 15,
) -> None:
    """GeoDataFrameをGeoParquet 1.1.0として保存します。

    CRSが未設定のまま保存すると、読み手が座標の意味を正しく判断できません。
    そのため、書き出す前に必ずCRSを確認します。
    """
    if gdf.crs is None:
        raise ValueError(
            "CRSが未設定です。元データの仕様を確認してください。"
        )

    output_path.parent.mkdir(parents=True, exist_ok=True)

    # WKBは現在のGISツールで最も互換性を確保しやすい符号化です。
    # write_covering_bbox=Trueを指定すると、各地物の外接矩形がbbox列へ
    # 書き込まれ、GeoPandasのbbox読み込みやDuckDBの一次絞り込みに使えます。
    gdf.to_parquet(
        output_path,
        index=False,
        compression="zstd",
        compression_level=compression_level,
        geometry_encoding="WKB",
        write_covering_bbox=True,
        schema_version="1.1.0",
        row_group_size=row_group_size,
    )

ここでは次の設定にしています。

  • schema_version="1.1.0"
  • geometry_encoding="WKB"
  • write_covering_bbox=True
  • compression="zstd"
  • row_group_size=50_000

配布向けGeoParquetの公式ガイドでは、ZSTD圧縮、空間的な並べ替え、5万～15万行程度のRow Groupが推奨されています。ただし、最適値は地物の複雑さと検索方法で変わるため、実データで計測して決めます。

なぜ空間順に並べ替えるのか

GeoParquetは、GeoPackageのR-treeやFlatGeobufのPacked Hilbert R-treeと同じ空間インデックスを持つわけではありません。

GeoParquet 1.1ではbbox列の統計、GeoParquet 2.0ではgeometry列の空間統計を使い、対象範囲と交差しないRow Groupなどを読み飛ばします。

全国の地物がランダムな順番で並んでいると、1つのRow Groupの範囲が全国へ広がり、読み飛ばしの効果が小さくなります。近い地物を近い行へ並べておくことが重要です。

サンプルでは分かりやすさを優先し、単純なグリッド番号で並べています。

import numpy as np


# 0.01度ごとの簡易グリッド番号を作ります。
# 実運用ではHilbert曲線、GeoHash、メッシュコードなども検討します。
gdf["_sort_y"] = np.floor(
    (gdf["latitude"] - 34.90) / 0.01
).astype("int32")
gdf["_sort_x"] = np.floor(
    (gdf["longitude"] - 138.30) / 0.01
).astype("int32")

gdf = (
    gdf.sort_values(["_sort_y", "_sort_x", "feature_id"], kind="stable")
    .drop(columns=["_sort_x", "_sort_y"])
    .reset_index(drop=True)
)

メタデータを確認する

PyArrowを使うと、ParquetのFooterとGeoParquetのgeoメタデータを確認できます。

python src/inspect_geoparquet.py data/output/sample_points.parquet

中心部分は次のコードです。

import json
from pathlib import Path

import pyarrow.parquet as pq


def read_geo_metadata(path: Path) -> dict | None:
    """ParquetのKey-Value Metadataからgeoメタデータを取得します。"""
    metadata = pq.read_metadata(path)
    key_value_metadata = metadata.metadata or {}
    geo_raw = key_value_metadata.get(b"geo")

    if geo_raw is None:
        return None

    return json.loads(geo_raw.decode("utf-8"))

最低限、次の項目を確認しておきます。

  • GeoParquetの仕様バージョン
  • 主ジオメトリ列
  • geometryの符号化
  • CRS
  • geometry type
  • bbox coveringの有無
  • Row Group数
  • 圧縮方式

GeoPandasでBBOX指定して読む

write_covering_bbox=Trueで保存したGeoParquet 1.1は、GeoPandasのread_parquet()へbboxを渡して部分読み込みできます。

python src/read_bbox.py \
  data/output/sample_points.parquet \
  --bbox 138.37 34.97 138.43 35.03

コードは次の形です。

import geopandas as gpd


# bboxはxmin, ymin, xmax, ymaxの順です。
# 指定する座標はGeoParquet側のCRSへ合わせます。
gdf = gpd.read_parquet(
    "data/output/sample_points.parquet",
    columns=[
        "feature_id",
        "category",
        "status",
        "value",
        "geometry",
    ],
    bbox=(138.37, 34.97, 138.43, 35.03),
)

columnsを指定する場合でも、GeoDataFrameとして読み込むにはgeometry列が必要です。属性だけを読みたい場合は、pandas.read_parquet()やDuckDBを使います。

DuckDBで属性集計する

DuckDBはParquetを直接検索できます。いったんDuckDBファイルへ取り込まなくても、SQLからそのまま集計できます。

SELECT
    category,
    COUNT(*) AS record_count,
    ROUND(AVG(value), 2) AS average_value
FROM read_parquet('data/output/sample_points.parquet')
GROUP BY category
ORDER BY record_count DESC;

このSQLではgeometryを参照していないため、DuckDBは必要な属性列だけを読みます。

Pythonから実行する場合は次のとおりです。

from pathlib import Path

import duckdb


path = Path("data/output/sample_points.parquet").resolve()
path_sql = path.as_posix().replace("'", "''")

connection = duckdb.connect(database=":memory:")
try:
    sql = f"""
        SELECT
            category,
            COUNT(*) AS record_count,
            ROUND(AVG(value), 2) AS average_value
        FROM read_parquet('{path_sql}')
        GROUP BY category
        ORDER BY record_count DESC
    """
    result = connection.execute(sql).fetchdf()
    print(result)
finally:
    connection.close()

DuckDBでBBOX検索する

GeoParquet 1.1のbboxは、次のようなSTRUCT列です。

bbox
├─ xmin
├─ ymin
├─ xmax
└─ ymax

DuckDBでは、通常の数値条件として候補を絞れます。

SELECT
    feature_id,
    category,
    status,
    value,
    longitude,
    latitude
FROM read_parquet('data/output/sample_points.parquet')
WHERE bbox.xmax >= 138.37
  AND bbox.xmin <= 138.43
  AND bbox.ymax >= 34.97
  AND bbox.ymin <= 35.03
ORDER BY feature_id;

この判定は外接矩形同士の交差です。Pointではそのまま最終結果として使えますが、LineStringやPolygonでは候補が多めに残る場合があります。

厳密に判定する場合は、DuckDBのspatial拡張を読み込み、ST_Intersects()を追加します。

INSTALL spatial;
LOAD spatial;

SELECT
    feature_id,
    category,
    value
FROM read_parquet('data/output/sample_points.parquet')
WHERE bbox.xmax >= 138.37
  AND bbox.xmin <= 138.43
  AND bbox.ymax >= 34.97
  AND bbox.ymin <= 35.03
  AND ST_Intersects(
        ST_GeomFromWKB(geometry),
        ST_MakeEnvelope(138.37, 34.97, 138.43, 35.03)
      );

DuckDBのバージョンやGeoParquetの形式によって、geometry列が最初からGEOMETRYとして認識される場合があります。その場合はST_GeomFromWKB(geometry)ではなく、geometryをそのまま空間関数へ渡します。添付サンプルでは列型を確認して式を切り替えています。

GeoPackageをGeoParquetへ変換する

添付サンプルのconvert_to_geoparquet.pyを使います。

python src/convert_to_geoparquet.py \
  data/input/buildings.gpkg \
  data/output/buildings.parquet \
  --layer buildings \
  --spatial-sort

平面直角座標系へ変換して保存する例です。

python src/convert_to_geoparquet.py \
  data/input/buildings.gpkg \
  data/output/buildings_epsg6676.parquet \
  --layer buildings \
  --target-crs EPSG:6676 \
  --make-valid \
  --spatial-sort

CRSが未設定の場合、スクリプトは変換を止めます。座標値だけを見てEPSGコードを推測すると、位置ずれを見落としやすいためです。

--make-validを付けると不正ジオメトリを修復します。ただし、修復後にGeometryCollectionへ変わることがあります。変換前後で地物数、geometry type、空間範囲を確認してください。

GDALで変換する

GDALのParquetドライバを使う方法もあります。

ogr2ogr \
  -f Parquet \
  data/output/buildings.parquet \
  data/input/buildings.gpkg \
  buildings \
  -lco COMPRESSION=ZSTD \
  -lco WRITE_COVERING_BBOX=YES \
  -lco SORT_BY_BBOX=YES \
  -lco ROW_GROUP_SIZE=100000

SORT_BY_BBOX=YESは空間的な並べ替えに便利ですが、大きなデータでは一時領域も必要になります。処理時間と一時ディスク容量を確認して使います。

実務で確認しておきたい点

CRSを曖昧にしない

GeoParquet 1.1では、crs項目が省略された場合、仕様上はOGC:CRS84として扱われます。CRS不明を意味するわけではありません。

元データのCRSが分からないときに、安易にset_crs("EPSG:4326")を実行しない方が安全です。set_crs()は座標値を変換せず、CRS情報だけを設定します。座標変換にはto_crs()を使います。

安定したIDを持たせる

GeoParquetは一括分析や配布に向きます。別ファイルの属性結果を結合したり、更新差分を追ったりする場合は、行番号ではなく安定した一意IDを持たせます。

building_id
mesh_id
station_id
scenario_id

geometryを毎回読まない

属性集計だけなら、geometry列は不要です。列指向形式の利点を生かすため、必要な列だけをSELECTします。

SELECT
    municipality_code,
    damage_rank,
    COUNT(*) AS building_count
FROM read_parquet('buildings.parquet')
GROUP BY municipality_code, damage_rank;

1ファイルを大きくしすぎない

配布ガイドでは、2GBを超える大きなデータは空間分割も検討するよう案内されています。

たとえば、地域メッシュや都道府県コードで分けます。

data/
├─ prefecture_code=22/
│  ├─ mesh=5238/part-000.parquet
│  └─ mesh=5239/part-000.parquet
└─ prefecture_code=23/
   └─ mesh=5237/part-000.parquet

DuckDBは複数のParquetをまとめて読めます。

SELECT *
FROM read_parquet('data/**/*.parquet');

Hive形式でフォルダーを分けると、条件に合わないファイル自体を読み飛ばせます。

頻繁な編集には使わない

GeoParquetは、地物を1件ずつ追加・更新する用途には向きません。編集途中のデータはGeoPackageやPostGISで管理し、分析・配布段階でGeoParquetへ書き出す方が扱いやすくなります。

どのような案件で使いやすいか

GeoParquetは、次のようなデータと相性が良いです。

  • PLATEAUの建物データ
  • 国勢調査やe-Statの地域メッシュ
  • 河川、水路、道路などの全国ベクトル
  • AMeDASや水位観測所の時系列属性付き地点
  • 浸水想定、被害判定、危険度評価の計算結果
  • DuckDB、Spark、DataFrame系ツールで集計する中間データ

OpenLayersへGeoParquetをそのまま渡すより、DuckDBで表示範囲を絞り、FastAPIからMVTやGeoJSONとして返す構成の方が、ブラウザ側を軽くしやすくなります。

GeoParquet
    ↓
DuckDBでBBOX検索・属性集計
    ↓
FastAPI
    ├─ MVT
    ├─ GeoJSON
    └─ 集計JSON
    ↓
OpenLayers

テスト

添付サンプルには、次の内容を確認するテストを入れています。

  • GeoParquet 1.1として保存できること
  • geoメタデータが存在すること
  • WKBとbbox coveringが設定されていること
  • 指定したRow Group数になること
  • GeoPandasでBBOX読み込みできること
  • DuckDBで属性集計できること
pytest -q

まとめ

GeoParquetは、GIS編集用というより、大量のベクトルデータを分析・保管・配布するための形式です。

最初に試す設定としては、次の組み合わせが扱いやすいと思います。

GeoParquet 1.1
WKB
ZSTD
bbox coveringあり
空間順に並べ替え
Row Groupは5万～10万行から確認
CRSと一意IDを明示

GeoParquet 2.0へ移行できる環境では、ParquetネイティブのGEOMETRY / GEOGRAPHY型と空間統計を利用できます。一方で、既存ツールを含めた互換性が必要なら、1.1を選ぶ理由もまだあります。

形式だけを先に決めるのではなく、QGISなどのアプリケーションで編集するのか、DuckDBで集計するのか、Web地図へ配信するのかを整理してから使い分けるのが良さそうです。

今後、いろいろな使い方について記事を書いていく予定です。

参考リンク

1
1
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
1
1

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?