DuckDBのspatial拡張には、PostGISと同じ名前の ST_AsMVT・ST_AsMVTGeom・ST_TileEnvelope があり、データベースサーバーを立てずにMapbox Vector Tile(MVT)を1枚ずつ生成できます。MVTは地図の図形と属性をタイル単位のバイナリ(Protocol Buffers)に詰めた形式で、MapLibre GL JSなどのブラウザ側ライブラリが描画します。
この記事では、DuckDB 1.5.6(2026年9月28日公開)で実際にタイルを生成し、デコードして中身を確かめた結果をもとに、SQLの組み立て方、空のタイルが返る原因、1タイルあたりの応答時間を左右する索引の条件、Flaskでの配信、PostGISとの違いを順に説明します。ベクタタイルとラスタタイルの違いや配信方式の選び方はWeb GISとは?タイル配信とベクタタイル・ライブラリ選定で扱っています。
まとめ|DuckDBでMVTを生成するときの要点
ST_AsMVTとST_AsMVTGeomはDuckDB 1.4.0(2025年9月16日)から使えます。spatial拡張は自動で読み込まれないため、LOAD spatialが必要です。- 座標は先にWebメルカトル(EPSG:3857)へ変換します。
ST_Transform(geom, 'EPSG:4326', 'EPSG:3857')をalways_xy := trueなしで書くと日本のように経度が90度を超える地域の座標はinfになり、タイル生成のSQLはエラーにならずに中身が空のタイルを返します。 - 100万点のテーブルでは、R-tree索引を張って
WHERE ST_Intersects(列, ST_TileEnvelope(...))と列を直接書いたときだけ索引が使われ、z14の1タイルが234.6msから2.1msになりました。今回の1.5.6の検証では、サブクエリで列に別名を付けたSQLで索引が使われませんでした。 - 低ズームは1枚に点が集中します。z7では28,656点・429,063バイトになり、格子に集約して202,255バイトまで減らしました。
- Flaskで配信するときは、リクエストごとに
cursor()で接続を分けます。1つの接続を共有した実験では、16並列200リクエストのうち51件が500エラーになりました。
DuckDBで使えるMVT関連関数と対応バージョン
MVTの生成に使う関数は3つです。ST_TileEnvelope は2025年6月10日、ST_AsMVT と ST_AsMVTGeom は2025年9月10日にduckdb-spatialのリポジトリへマージされ、DuckDB 1.4.0から利用できます。1.4系は1年間のコミュニティサポートが付く長期サポート(LTS)版で、1.4.5でも同じ3関数が使えることを確認しました。
| 関数 | 種類 | 役割 | 既定値 |
|---|---|---|---|
| ST_TileEnvelope(z, x, y) | スカラー | タイルの範囲をEPSG:3857の矩形で返す | 引数は3つのみ |
| ST_AsMVTGeom(geom, bounds, extent, buffer, clip_geom) | スカラー | 図形をタイル内の座標へ変換し切り抜く | 4096 / 256 / true |
| ST_AsMVT(row, layer_name, extent, geom列名, ID列名) | 集約 | 行の集合を1枚のタイル(BLOB)にする | ‘layer’ / 4096 |
DuckDB 1.5.0(2026年3月9日)では、GEOMETRY型がspatial拡張からDuckDB本体へ移り、GEOMETRY('OGC:CRS84') のように座標参照系(CRS)を型に持てるようになりました。関数の使い方は1.4系と変わりませんが、後述の座標変換では読み込んだ形式によって結果が変わります。DuckDB本体の特徴や版の状況はDuckDBとは?特徴・用途とSQLite・PostgreSQLとの違い、インストール手順はDuckDBの使い方を参照してください。
1枚のタイルを生成するSQLの組み立て方
タイルを1枚作るSQLは、範囲を決める、図形をタイル座標へ移す、行をまとめてバイナリにする、の3段で書きます。次の例は、GeoJSONの駅データからz14・x14552・y6451(東京駅周辺)のタイルを作ります。
INSTALL spatial;
LOAD spatial;
-- EPSG:3857 に変換した列を持つテーブルを作る
CREATE TABLE stations AS
SELECT OGC_FID AS id,
name,
ST_Transform(geom, 'EPSG:4326', 'EPSG:3857', always_xy := true) AS geom_3857
FROM ST_Read('stations.geojson');
-- z=14, x=14552, y=6451 のタイルを1枚作る
SELECT ST_AsMVT(
{'geom': ST_AsMVTGeom(geom_3857, ST_Extent(ST_TileEnvelope(14, 14552, 6451)), 4096, 256, true),
'id': id,
'name': name},
'stations', 4096, 'geom', 'id') AS tile
FROM stations
WHERE ST_Intersects(geom_3857, ST_TileEnvelope(14, 14552, 6451));
東京駅周辺に1,000点を置いたGeoJSONでこのSQLを実行し、Pythonの mapbox-vector-tile でデコードすると、タイル内の456点が12,551バイトに収まり、各点は [2534, 2938] のように0〜4096の整数座標で入っていました。
ST_TileEnvelopeとST_AsMVTGeomによるタイル座標変換
ST_AsMVTGeom の第2引数はBOX_2D型です。ST_TileEnvelope が返すのはGEOMETRY型のため、そのまま渡すと Binder Error: No function matches the given name and argument types 'ST_AsMVTGeom(GEOMETRY, GEOMETRY)' で止まります。ST_Extent() で包むのが正しい書き方で、::BOX_2D へのキャストは Unimplemented type for cast になります。
第4引数のbufferは、タイルの外側にどれだけ図形を残すかの幅です。半径3kmの円を同じタイルに入れると、clip_geom がtrueなら座標は-256〜4352の範囲に切り抜かれて40バイト、falseなら-1127〜8921まで残って156バイトになりました。線や面を隣のタイルと継ぎ目なく描くにはある程度のbufferが要り、点だけのレイヤーでも、円の半径やラベルの大きさによって必要なbufferは変わります。タイル境界で欠けないよう、抽出範囲もbufferに対応して広げて確認します。
ST_AsMVTに渡すSTRUCTと属性の型の制約
ST_AsMVT の第1引数は、図形1列と属性の列を並べたSTRUCTです。属性に使える型はVARCHAR・FLOAT・DOUBLE・INTEGER・BIGINT・BOOLEANの6つに限られ、それ以外は実行時にエラーになります。
- DATE・DECIMAL(10,2)・SMALLINT・INTEGER[] は
ST_AsMVT: property column "p" has unsupported type "SMALLINT"の形で拒否されました。SMALLINTも通らないため、::INTEGERへのキャストが要ります。日付はstrftimeで文字列にします。 - 値がNULLの属性はエラーにならず、その地物のプロパティから省かれます。
- 図形の列を2つ入れると
only one geometry column is allowed in the input rowで止まります。 - 第5引数に指定した列(INTEGERかBIGINT)は地物IDになり、属性には残りません。第2引数を省いたときのレイヤー名は
layerです。
2つのレイヤーを1枚のタイルに入れたいときは、レイヤーごとに ST_AsMVT を呼び、BLOBを || で連結します。実際に連結した65バイトのタイルをデコードすると、レイヤーaとbの2つが読めました。
COPYによるMVTファイルの書き出し
配信サーバーを置かずに静的ファイルとして並べる場合は、COPY ... TO にBLOB形式を指定します。MVTの仕様は拡張子を .mvt とすることを推奨しています。
COPY (
SELECT ST_AsMVT({'geom': ST_AsMVTGeom(geom_3857, ST_Extent(ST_TileEnvelope(10, 909, 403))),
'name': name}, 'stations')
FROM stations
WHERE ST_Intersects(geom_3857, ST_TileEnvelope(10, 909, 403))
) TO 'tiles/10/909/403.mvt' (FORMAT BLOB);
ディレクトリは自動では作られないため、tiles/10/909/ を先に作っておきます。全ズームを書き出すとファイル数はズームが1上がるごとに最大4倍になるので、静的に並べる用途は対象範囲とズームを絞った場合に限ります。
座標系の指定ミスで空のタイルが返る原因
エラーが出ないまま地図に何も表示されない場合は、座標変換、タイルの抽出範囲、配信の成否、描画側のレイヤー名を順に確認します。この章では座標変換によって空のタイルが生成される例を扱います。
EPSG:4326は、定義上「緯度、経度」の順に座標を持つ座標系です。一方、GeoJSONや多くのデータは「経度、緯度」の順で持っています。ST_Transform(ST_Point(139.767, 35.681), 'EPSG:4326', 'EPSG:3857') を実行すると、139.767を緯度として扱うため結果は POINT (inf inf) になりました。always_xy := true を付けると POINT (15558791.27 4256815.78) と正しい値が返ります。経度が90度以下の地域では inf にならず、ST_Point(10, 20) が POINT (2226389.82 1118889.97)(正しくは POINT (1113194.91 2273030.93))のように経度と緯度が入れ替わった有限の値になるため、値が有限であることだけでは正常と判断できません。DuckDBの関数説明もこの引数について、EPSG:4326との変換で経度・緯度の順に解釈させるためのものと書いています。
座標が inf になった点は ST_Intersects でどのタイルとも交わらないため、ST_AsMVT には0行が渡ります。このとき返るのはNULLではなく、レイヤー名だけを持つ20バイト程度のBLOBです。上の例では、変換を忘れた場合が21バイト、always_xy を忘れた場合が23バイトでした。NULLを判定してもこの失敗は検出できないため、開発中はタイルのバイト数か地物の数を確かめてください。
ST_ReadとGeoParquetで変換結果が変わる理由
DuckDB 1.5.6では、読み込んだ形式によってGEOMETRY型に付くCRSが違い、引数2つの ST_Transform(geom, 'EPSG:3857') の結果も変わりました。
| 読み込み元 | 列の型 | ST_Transform(geom, ‘EPSG:3857’) | always_xy := true |
|---|---|---|---|
| ST_Read(GeoJSON) | GEOMETRY(‘EPSG:4326’) | POINT (inf inf) | 正しい値 |
| GeoParquet(DuckDBで書き出し) | GEOMETRY(‘OGC:CRS84’) | 正しい値 | 正しい値 |
| 1.4.5でGeoParquetを読む | GEOMETRY(CRSなし) | 引数3つが必要 | 正しい値 |
OGC:CRS84は経度・緯度の順で定義された座標系なので、軸の入れ替えが起きません。読み込み元が混在しても結果を揃えるには、always_xy := true を常に付けるのが確実です。
ST_Read はGDALを使ってGeoJSON・Shapefile・GeoPackage・FlatGeobufなどを読む関数で、公式の説明に「GDALはシングルスレッド」とあり並列化されません。GeoParquetはGDALを通らず、DuckDB本体のParquet読み込みで直接GEOMETRY列として読めます。繰り返し読む大きなデータではGeoParquetへの事前変換を候補にし、変換コストと対象クエリの読み込み時間を比較して採用します。座標系の考え方そのものはGISデータとは?形式・座標系・取り込み設計、Parquetの仕組みはParquetとはで解説しています。
1タイルあたりの応答時間:R-tree索引が効く条件
日本付近の範囲(東経129.5〜145.5度、北緯31〜45度)にランダムな100万点を置き、1タイルの生成時間を計りました(Intel Core i9-9880H・16GBメモリ、各5〜7回の中央値)。
| 条件 | z14(1点) | z10(452点) | 実行計画 |
|---|---|---|---|
| リクエストごとに座標変換 | 382ms | – | 全件走査 |
| 3857の列を事前に作成 | 236ms | 242ms | 全件走査 |
| R-tree索引+サブクエリで別名 | 218ms | 225ms | SEQ_SCAN |
| R-tree索引+列を直接WHEREに書く | 2.1ms | 29.6ms | RTREE_INDEX_SCAN |
| 索引なし・ST_Intersects_Extent | 67.9ms | – | 全件走査 |
差を生んだのはSQLの書き方です。FROM (SELECT g3857 AS g FROM p) WHERE ST_Intersects(g, ...) のように別名を付けたサブクエリを挟むと、索引があっても実行計画はSEQ_SCANのままでした。FROM p WHERE ST_Intersects(g3857, ST_TileEnvelope(...)) と書き直すと、EXPLAIN に RTREE_INDEX_SCAN が出て、z14は約110分の1の時間になりました。
-- 変換済みの列に R-tree 索引を張る
CREATE INDEX stations_rtree ON stations USING RTREE (geom_3857);
-- 索引が使われているかを確かめる
EXPLAIN
SELECT count(*) FROM stations
WHERE ST_Intersects(geom_3857, ST_TileEnvelope(14, 14552, 6451));
公式ドキュメントは、索引が使われる条件を「WHERE句でST_Intersectsなどの空間述語を使い、引数の一方がクエリ計画時に値の決まる定数であること」としています。z・x・yを $z のようなプレースホルダで渡した場合も、EXPLAIN には RTREE_INDEX_SCAN が出て、z14の1タイルは4.4msでした。後述のFlaskの例はこの書き方です。また、R-tree索引は一度メモリに載ると索引を削除するまで解放されず、その分は memory_limit に数えられます。データが大きい場合はこのメモリ量も見積もりに入れてください。
低ズームでタイルが大きくなる問題と間引き
ズームを下げるほど1枚のタイルが広い範囲を覆うため、点の数とバイト数が急増します。同じ100万点で、関東から中部にかかるz7のタイル(113, 50、東経137.8〜140.6度・北緯34.3〜36.6度)は28,656点・429,063バイトになり、gzipで圧縮しても136,799バイトありました。ズーム14の1タイルが36バイトだったのと比べると、同じSQLのままでは低ズームだけが極端に重くなります。
低ズームでは個々の点を送らず、タイル座標上の格子にまとめて件数を属性にする方法が使えます。次のSQLは4096×4096のタイル座標を32座標単位四方の格子に区切り、格子ごとに1点と件数 n を出します。
WITH t AS (
SELECT ST_AsMVTGeom(geom_3857, ST_Extent(ST_TileEnvelope(7, 113, 50)), 4096, 64, true) AS g
FROM stations
WHERE ST_Intersects(geom_3857, ST_TileEnvelope(7, 113, 50))
)
SELECT ST_AsMVT({'geom': ST_Point(gx * 32 + 16, gy * 32 + 16), 'n': n}, 'stations')
FROM (
SELECT floor(ST_X(g) / 32)::INTEGER AS gx,
floor(ST_Y(g) / 32)::INTEGER AS gy,
count(*)::INTEGER AS n
FROM t
GROUP BY ALL
);
この方法でz7のタイルは13,508点・202,255バイトになりました。ST_AsMVTGeom を通した後の座標はすでにタイル座標なので、集約した点はそのまま ST_AsMVT に渡せます。格子を64座標単位四方に広げると4,102点・61,491バイトまで減りました。どのズームから集約に切り替えるかは、描画側で n を円の大きさに使うかどうかと合わせて決めます。
Flaskでタイルを配信しMapLibreで表示する実装
配信前に、書き込み可能な接続で tiles.duckdb を開き、前述のSQLで stations テーブルとR-tree索引を作成して接続を閉じます。次のコードを app.py に保存し、duckdbとFlaskを導入した環境で flask --app app run を実行すると、z・x・yに対応するMVTを返せます。MVTの仕様は、配信時のMIMEタイプを application/vnd.mapbox-vector-tile とすることを推奨しています。
import duckdb
from flask import Flask, Response
app = Flask(__name__)
db = duckdb.connect("tiles.duckdb", read_only=True)
db.execute("LOAD spatial")
SQL = """
SELECT ST_AsMVT({'geom': ST_AsMVTGeom(geom_3857, ST_Extent(ST_TileEnvelope($z, $x, $y)), 4096, 64, true),
'name': name}, 'stations')
FROM stations
WHERE ST_Intersects(geom_3857, ST_TileEnvelope($z, $x, $y))
"""
@app.get("/tiles/<int:z>/<int:x>/<int:y>.mvt")
def tile(z, x, y):
if not (0 <= z <= 22 and 0 <= x < 2 ** z and 0 <= y < 2 ** z):
return Response(status=404)
cur = db.cursor() # リクエストごとに接続を分ける
try:
blob = cur.execute(SQL, {"z": z, "x": x, "y": y}).fetchone()[0]
finally:
cur.close()
return Response(blob, mimetype="application/vnd.mapbox-vector-tile",
headers={"Access-Control-Allow-Origin": "*",
"Cache-Control": "public, max-age=3600"})
Flaskの開発サーバーは複数スレッドでリクエストを処理します。cur = db.cursor() の代わりに db を直接使う版で16並列・200リクエストを送ったところ、149件は成功し、51件が TypeError: 'NoneType' object is not subscriptable の500エラーになりました。別スレッドの実行が結果を上書きしたためで、cursor() でスレッドごとに接続を分けた版は200件すべて成功しています。
ブラウザ側では、MapLibre GL JSのCSSと幅・高さを指定した id="map" の要素を用意します。この例では地図ページをタイルAPIと同じオリジンで配信し、vectorソースにタイルのURLを渡し、レイヤーの source-layer に ST_AsMVT で付けたレイヤー名を指定します。
import * as maplibregl from './maplibre-gl.mjs';
const map = new maplibregl.Map({
container: 'map',
center: [139.767, 35.681],
zoom: 10,
style: {
version: 8,
sources: {
stations: { type: 'vector', tiles: [location.origin + '/tiles/{z}/{x}/{y}.mvt'], minzoom: 0, maxzoom: 14 }
},
layers: [
{ id: 'stations', type: 'circle', source: 'stations', 'source-layer': 'stations',
paint: { 'circle-radius': 2, 'circle-color': '#d33' } }
]
}
});
MapLibre GL JS 6.11.2をChrome 154(ヘッドレス)で開き、100万点のテーブルを配信したサーバーにz10でつなぐと、6枚のタイルが要求され、querySourceFeatures で2,768地物が読み込まれました。MapLibre 6系はESM形式のみで配布されているため、npmパッケージ maplibre-gl の dist にある .mjs ファイル群(ワーカー用を含む)とCSSをページと同じ場所に置き、<script type="module"> の中で上のimportを書きました。source-layer の名前が ST_AsMVT の第2引数と一致していないと、タイルは取得されても何も描かれません。
PostGISのST_AsMVTとの違い
関数名と引数の並びはPostGISに合わせてありますが、PostGISのSQLをそのまま流すと動かない箇所があります。
| 項目 | DuckDB spatial | PostGIS |
|---|---|---|
| ST_AsMVTGeomの範囲指定 | BOX_2D(ST_Extentで包む) | box2d(ST_TileEnvelopeを直接渡せる) |
| ST_TileEnvelopeの引数 | z, x, y のみ | bounds指定可・marginは3.1.0〜 |
| レイヤー名の既定値 | layer | default |
| 属性に使える型 | 6型に限定 | 任意の列(JSONBで可変属性も可) |
| 地物IDの型 | INTEGER・BIGINT | smallint・integer・bigint |
| 範囲の絞り込み | ST_Intersects+R-tree索引 | && 演算子+GiST索引 |
ST_AsMVTGeom の既定値(extent 4096・buffer 256・clip_geom true)は両者で同じで、DuckDBで引数を省いた結果と明示した結果は一致しました。PostGISの公式例は ST_TileEnvelope の margin で絞り込み範囲をbuffer分だけ広げていますが、DuckDBには margin がないため、ラベルを隣のタイルにまたがって描きたい場合は ST_Buffer などで範囲を自分で広げます。
DuckDBで動的に生成するか、事前に生成するかの判断
DuckDBでのタイル生成が向くのは、データが数百万件以下で、更新が日次程度のバッチに限られ、配信先が社内ツールや分析用の地図である場合です。サーバー1プロセスとDuckDBファイル1つで動き、PostGISを用意する必要がありません。上の実測どおり、索引が効く高ズームなら1枚数msで返せます。
反対に、DuckDBを配信サーバーにすべきでないのは、地図を見ている間にもデータが書き換わる場合です。タイルサーバーが read_only=True でファイルを開いている間、別プロセスが同じファイルを読み書きで開くと IO Error: Could not set lock on file で失敗しました(読み取り専用どうしなら2つ目も開けます)。更新のたびにサーバーを止めるか、新しいファイルを書いてから差し替える運用になります。GeoParquetを直接読む構成ならこのロックは起きませんが、R-tree索引が使えなくなります。
- 更新が頻繁で、書き込みと配信を同時に続ける必要がある場合は、PostGISで動的に生成する構成を選びます。
- 不特定多数に公開し、データがほぼ変わらない場合は、tippecanoeなどでMBTilesやPMTilesに事前生成し、CDNから静的に配る方が運用は軽く済みます。
- ブラウザだけで完結させたい場合は、DuckDB-WASMで同じSQLを実行し、MapLibreに渡す方法もあります(DuckDB WASMとは)。
配信方式の比較と、PostGISやPMTilesを含む構成の選び方はWeb GISとは?タイル配信とベクタタイル・ライブラリ選定で詳しく扱っています。
よくある質問
ST_AsMVTとは何をする関数ですか?
行の集合を1枚のMapbox Vector Tile(MVT)のバイナリにまとめる集約関数です。PostGISでは2.4.0から、DuckDBではspatial拡張を読み込んだ1.4.0以降で使えます。図形は事前に ST_AsMVTGeom でタイル内の座標(0〜4096)へ変換しておきます。
DuckDBのどのバージョンからST_AsMVTを使えますか?
1.4.0(2025年9月16日)からです。LTSの1.4.5と最新の1.5.6で動作を確認しました。PythonではPyPIのduckdb 1.5.0以降がPython 3.10以上を求めるため、Python 3.9の環境では1.4系が入ります。どちらでも INSTALL spatial; LOAD spatial; が必要です。
DuckDBのST_Transformでinfが返るのはなぜですか?
EPSG:4326を「緯度、経度」の順で解釈し、経度の値を緯度として変換しているためです。ST_Transform(geom, 'EPSG:4326', 'EPSG:3857', always_xy := true) と指定すると経度・緯度の順で扱われ、正しい値になります。検証した1.5.6の東京駅周辺のGeoJSONでも、ST_Read で読み込んだ列を always_xy なしで変換すると inf になりました。
GeoParquetはGDALを通さずに読み込めますか?
読み込めます。DuckDBはGeoParquetの図形列をParquetの読み込み機能でGEOMETRY型として直接読みます。ST_Read はGDAL経由でシングルスレッドで動くため、Shapefile・GeoPackage・GeoJSONなどを繰り返し読む場合はGeoParquetへの事前変換を検討し、対象データで読み込み時間を比較してください。
ST_AsMVTの結果が空のタイルになるのはなぜですか?
多くは座標変換の誤りで、図形がタイルの範囲と交わらず0行が集約されています。このとき結果はNULLではなく、レイヤー名だけを持つ20バイト前後のBLOBです。ST_AsText で変換後の座標が inf や度単位の小さな値になっていないか、ST_AsMVT の第2引数とMapLibreの source-layer が一致しているかを確かめてください。