Go to Qiita Advent Calendar 2024 Top
LoginSignup
5
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?

More than 1 year has passed since last update.

RDBMS-GIS(MySQL,PostgreSQLなど)Advent Calendar 2022

Day 19
@myonyomu

PostGISの空間検索速度とインデックスの話、あるいは&&のおまじないの現在について

Last updated at Posted at 2022-12-18

この記事はRDBMS-GIS(MySQL,PostgreSQLなど) Advent Calendar 2022の19日目のエントリーです。

カレンダー上の他の人の記事にMySQLが多くて「カレンダーの種類間違えたんじゃないか」とビビり散らかしてましたが、PostGISの文字を指差し確認したので大丈夫だと判断して出します。

先に言いたいことについてまとめ

PostGISでST_Intersects等の空間演算で、Bounding Boxでの検索を重ねて利用することでオプティマイザにインデックスの利用を明示することができるというテクニックがありました。

-- ST_Intersectsだけでも使ってインデックスを利用して検索してくれそうだが無理。

SELECT
  geometry
FROM
  geometry_table
WHERE
  -- Boundary box同士を比較。不要だが、Indexを利用するためには入れなければいけない。 
  geometry_query && geometry 
  -- 実際に入れたい検索条件。あるジオメトリに交差しているジオメトリを出したい。
  AND ST_Intersects(ST_GeomFromText('POLYGON((0 0, 90 0, 90 90, 0 90, 0 0))', 4326), geometry);

3年ほど前にもブログで検証してます。

最近お仕事でPostgreSQLのバージョンアップ作業をやってきて、
「このおまじない的に入れられてきた&&オペレーターは、PostGISのメジャーバージョンが上がった現代でも同じように利用しても良いのか？Release Note見る限りIndex周りや空間演算系にちょいちょい機能追加入ってない？でもメジャーバージョン上がってるし一つ一つ機能見ていくには大量にありすぎるしなあ……」1
と感じる機会があったため、「もしいまだにおまじないが必要ならば知らない人にこのテクニックを周知できるし、必要なければクエリの修正点が見えるためヨシ！」の心持ちでこの記事を書き始めました。

なお結論としては、一応おまじないを抜きにできるけど、リファレンス先に読んでおけば確認する必要なかったかもしれないと、リファレンス見て適切な型をオプティマイザに認識させないと遅くなることもあるです。

前置き：RDBにおける地理空間の検索

GISはその名の通り、地理情報をあつかうシステムのことです。
地理情報を扱う、ということは、RDBに店舗・顧客の地点位置情報や、都道府県・市町村の行政区画などの凹多角形を格納する必要性が出てくることもあります。
また、情報を格納するということは、当然のようにシステムで検索する必要性も出てきます。

GISにおいては、検索を行うためのパラメータに特定カラムの値の一致を見るだけでなく、空間的な関係性を用いるというパターンも存在します。
SQLにおける表現では、前者はWHERE point_name = '皇居'のようにWHERE句で検索パラメータを利用することが可能ですが、空間検索においては、「地点Aから半径2km以内に存在する地点すべて」というようなパラメータを文字列や数値を利用して表現することは難しいです。

image.png

この課題を解決するために、PostGISでは、PostgreSQLにおける空間検索を可能にする様々な拡張機能が提供されています。
例えば、そもそもRDBのカラムに地理情報を格納するためのGeometry型の定義や、上記の「地点Aから半径2km以内に存在する地点すべて」という検索を行うことができるST_Intersectsという関数などがその機能に当たります。

地理空間の検索におけるインデックス

RDBの検索ではインデックスを用いることで検索を高速化することが可能ですが、地理空間情報においても同様にインデックスを用いた検索の高速化を図ることができます。

ただし、PostGISで利用されているインデックスの種類は、B-treeではなくGiSTと呼ばれるインデックスで2、ありとあらゆる説明を簡略化して乱暴に言ってしまうと、「図形や空間情報を扱うときに向いているよ」「&&,<<のような演算子の挙動を制御して、インデックスを利用するようにするよ」というものです。
ただし、ここで気を付けなければいけないのは、GiSTで利用される演算子は空間図形そのもの同士を比較しているのではなく、空間図形の「バウンディングボックス」と呼ばれる外接する四角形を利用して判定を行っている、という点です。
例として、&&オペレータは二つの空間図形のバウンディングボックス同士を比較し、バウンディングボックス同士が重なっていればtrueを返す、という挙動となっています。

image.png

このようなバウンディングボックス同士の演算においてはオプティマイザが優先的にインデックスを利用するのですが、ST_IntersectsのようなGiSTで提供されている演算子以外の空間演算関数を利用するとうまくインデックスが利用されなくなるため、先にバウンディングボックスでの演算で篩にかけてから実際に行いたい空間演算を行う、というパターンが定着していました。

……と、いうのがPostgreSQL11系、PostGIS2系までの話。
前置きが長くなりましたが、今回の記事では「そんなことしなくても、今提供されているバージョンだったら空間演算関数だけでうまくIndexを使ってくれるんじゃないの？」というのを検証します。

検証環境

サーバー建てるのも面倒なので、PostGISが提供しているDocker Imageを利用します。

  • postgis/postgis:15-3.3
    • PostgreSQL: 15.1
    • PostGIS: 3.3.2

Dockerfile

FROM postgis/postgis:15-3.3

docker-compose.yaml

version: '3'
services:
  db:
    container_name: postgis-test
    build: .
    ports:
      - 5433:5432
    environment:
      POSTGRES_PASSWORD: postgres
    volumes:
      - ./data:/work
    working_dir: /work

workディレクトリ内にsqlを用意しておき、事前にテーブルを作ってデータを入れておきます。
CSVインポートを利用したいため、WKTからGeometryに変換して入れるためのカラムも作成しています。
また、単なる検証のため、データベースやユーザーの設定を簡略化してpostgresユーザーのまま実行してます。

CREATE TABLE geometry_test_table (
    mesh_id INTEGER PRIMARY KEY,
    population_count NUMERIC,
    geometry_polygon_wkt TEXT,
    geometry_polygon GEOMETRY(POLYGON, 4326)
);

CREATE INDEX
    geometry_polygon_idx
ON
    geometry_test_table
USING
    GIST(geometry_polygon);

今回入れるデータとしては、検索対象の数が稼げれば何でもよいので、国土交通省のページでダウンロード可能な500ｍメッシュのデータの北海道のものを利用しました。3
image.png
北海道のメッシュのみを入れると、おおよそ45000件。この中から検索することで検証してみます。

検証

ポリゴン同士の交差判定(ST_Intersects)

釧路のあたりを適当な範囲を取って検索

EXPLAIN ANALYZE
SELECT
    mesh_id
FROM
    geometry_test_table
WHERE
    ST_Intersects(ST_GeomFromText(
            'POLYGON((143.74278419330417 43.494825159846755,145.5575882704526 43.494825159846755,145.5575882704526 42.90695306110539,143.74278419330417 42.90695306110539,143.74278419330417 43.494825159846755))'
            , 4326
        ),
        geometry_polygon
    );

                                                                                                                           QUERY PLAN

-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
 Gather  (cost=1104.56..56831.69 rows=3649 width=4) (actual time=2.118..170.826 rows=3625 loops=1)
   Workers Planned: 1
   Workers Launched: 1
   ->  Parallel Bitmap Heap Scan on geometry_test_table  (cost=104.56..55466.79 rows=2146 width=4) (actual time=0.576..2.279 rows=1812 loops=2)
         Filter: st_intersects('0103000020E6100000010000000500000045475BE3C4F76140D0634B6E56BF454045475BC3D7316240D0634B6E56BF454045475BC3D7316240353AB4091774454045475BE3C4F76140353AB4091774454045475BE3C4F76140D0634B6E56BF4540'::geometry, geometry_polygon)
         Heap Blocks: exact=154
         ->  Bitmap Index Scan on geometry_polygon_idx  (cost=0.00..103.65 rows=3649 width=0) (actual time=0.825..0.826 rows=3625 loops=1)
               Index Cond: (geometry_polygon && '0103000020E6100000010000000500000045475BE3C4F76140D0634B6E56BF454045475BC3D7316240D0634B6E56BF454045475BC3D7316240353AB4091774454045475BE3C4F76140353AB4091774454045475BE3C4F76140D0634B6E56BF4540'::geometry)
 Planning Time: 0.744 ms
 Execution Time: 171.492 ms

Bitmap Heap Scanを行っています。

&&オペレータで先に篩をかけてからのポリゴン同士の交差判定

おまじないつきの方法です。

EXPLAIN ANALYZE
SELECT
    mesh_id
FROM
    geometry_test_table
WHERE
    geometry_polygon && ST_GeomFromText(
            'POLYGON((143.74278419330417 43.494825159846755,145.5575882704526 43.494825159846755,145.5575882704526 42.90695306110539,143.74278419330417 42.90695306110539,143.74278419330417 43.494825159846755))'
            , 4326
    )
    AND ST_Intersects(ST_GeomFromText(
            'POLYGON((143.74278419330417 43.494825159846755,145.5575882704526 43.494825159846755,145.5575882704526 42.90695306110539,143.74278419330417 42.90695306110539,143.74278419330417 43.494825159846755))'
            , 4326
        ),
        geometry_polygon
    );

                                                                                     QUERY PLAN
-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
 Bitmap Heap Scan on geometry_test_table  (cost=11.28..8095.12 rows=293 width=4) (actual time=0.755..4.524 rows=3625 loops=1)
   Recheck Cond: (geometry_polygon && '0103000020E6100000010000000500000045475BE3C4F76140D0634B6E56BF454045475BC3D7316240D0634B6E56BF454045475BC3D7316240353AB4091774454045475BE3C4F76140353AB4091774454045475BE3C4F76140D0634B6E56BF4540'::geometry)
   Filter: st_intersects('0103000020E6100000010000000500000045475BE3C4F76140D0634B6E56BF454045475BC3D7316240D0634B6E56BF454045475BC3D7316240353AB4091774454045475BE3C4F76140353AB4091774454045475BE3C4F76140D0634B6E56BF4540'::geometry, geometry_polygon)
   Heap Blocks: exact=154
   ->  Bitmap Index Scan on geometry_polygon_idx  (cost=0.00..11.21 rows=293 width=0) (actual time=0.663..0.679 rows=3625 loops=1)
         Index Cond: ((geometry_polygon && '0103000020E6100000010000000500000045475BE3C4F76140D0634B6E56BF454045475BC3D7316240D0634B6E56BF454045475BC3D7316240353AB4091774454045475BE3C4F76140353AB4091774454045475BE3C4F76140D0634B6E56BF4540'::geometry) AND (geometry_polygon && '0103000020E6100000010000000500000045475BE3C4F76140D0634B6E56BF454045475BC3D7316240D0634B6E56BF454045475BC3D7316240353AB4091774454045475BE3C4F76140353AB4091774454045475BE3C4F76140D0634B6E56BF4540'::geometry))
 Planning Time: 0.509 ms
 Execution Time: 5.323 ms

ST_Intersects内にST_GeomFromTextを利用した時と同様にBitmap Heap Scanを行っていますが、作成されたBitmapにRecheckを行い、条件を絞っています。

WITH句付きでポリゴン同士の交差判定

上記とやっていることは同じですが、検索に利用するパラメータをWITH句に切り出し、スカラー値にしています。

EXPLAIN ANALYZE
WITH geometry_query AS (
    SELECT
        ST_GeomFromText(
            'POLYGON((143.74278419330417 43.494825159846755,145.5575882704526 43.494825159846755,145.5575882704526 42.90695306110539,143.74278419330417 42.90695306110539,143.74278419330417 43.494825159846755))'
            , 4326
        )
)

SELECT
    mesh_id
FROM
    geometry_test_table
WHERE
    ST_Intersects((SELECT * FROM geometry_query), geometry_polygon);

                                                                  QUERY PLAN
-----------------------------------------------------------------------------------------------------------------------------------------------
 Index Scan using geometry_polygon_idx on geometry_test_table  (cost=0.29..149.38 rows=5 width=4) (actual time=0.144..4.326 rows=3625 loops=1)
   Index Cond: (geometry_polygon && $0)
   Filter: st_intersects($0, geometry_polygon)
   InitPlan 1 (returns $0)
     ->  Result  (cost=0.00..0.01 rows=1 width=32) (actual time=0.001..0.001 rows=1 loops=1)
 Planning Time: 0.114 ms
 Execution Time: 4.889 ms

Indexの状況が最初に確認され、&&オペレータによるインデックスの探索が行われているように見えます。
Index Only Scanに近そうな実行計画です。
実際の速度としては、クエリ内で関数を呼び出すよりもずっと早いです。

WITH句と&&オペレータおまじないの同時使用

最後に、全部盛りも試してみます。

EXPLAIN ANALYZE
WITH geometry_query AS (
    SELECT
        ST_GeomFromText(
            'POLYGON((143.74278419330417 43.494825159846755,145.5575882704526 43.494825159846755,145.5575882704526 42.90695306110539,143.74278419330417 42.90695306110539,143.74278419330417 43.494825159846755))'
            , 4326
        )
    )

SELECT
    mesh_id
FROM
    geometry_test_table
WHERE
    (SELECT * FROM geometry_query) && geometry_polygon
    AND ST_Intersects((SELECT * FROM geometry_query), geometry_polygon);

                                                                  QUERY PLAN
----------------------------------------------------------------------------------------------------------------------------------------------
 Index Scan using geometry_polygon_idx on geometry_test_table  (cost=0.33..33.35 rows=1 width=4) (actual time=0.291..4.841 rows=3625 loops=1)
   Index Cond: ((geometry_polygon && $1) AND (geometry_polygon && $2))
   Filter: st_intersects($2, geometry_polygon)
   CTE geometry_query
     ->  Result  (cost=0.00..0.01 rows=1 width=32) (actual time=0.038..0.038 rows=1 loops=1)
   InitPlan 2 (returns $1)
     ->  CTE Scan on geometry_query  (cost=0.00..0.02 rows=1 width=32) (actual time=0.040..0.041 rows=1 loops=1)
   InitPlan 3 (returns $2)
     ->  CTE Scan on geometry_query geometry_query_1  (cost=0.00..0.02 rows=1 width=32) (actual time=0.000..0.001 rows=1 loops=1)
 Planning Time: 0.520 ms
 Execution Time: 5.618 ms

こちらも高速ではありますが、同一のインデックススキャンをCTEを元に二回行っているため、若干の無駄があります。

検証からの結論

さて、ここまでで4件のパターンを見て速度を比較してきました。

手法 Planning Time Execution Time
追加なし 0.744 ms 171.492 ms
&&オペレータ 0.509 ms 5.323 ms
WITH句利用 0.114 ms 4.889 ms
&&とWITH句 0.520 ms 5.618 ms

&&オペレータの利用でも速度は出ていますが、それよりもWITH句を利用してスカラー値を関数のパラメータにした方が効率的なスキャンが行われています。
また、With句の場合は&&オペレータを利用せずにST_Intersectsを利用していても、オプティマイザは適切に演算子を選んでインデックススキャンを行っているように見えます。

Index Cond: (geometry_polygon && $0)

というわけで、

  • ST_Intersectsを代表とする空間演算関数のみを利用していても、適切にGiSTで提供されている演算子を利用してインデックススキャンが行われる。
  • ただし、空間演算関数に渡す引数には、オプティマイザが認識できる正しい型が必要になる

というのが結論です。

今回の場合は、WITH句を利用したスカラー値はPolygonであることが確定していました。しかし、そもそもST_GeomFromTextは、ポリゴンのみならず、LineやPointといったGeometryで表現可能な型すべてに変換可能です。
PostgreSQL14系のリファレンスによると、GiSTで提供される演算子はインデックス可能な演算子の型が決まっています。そのため、実際に変換後にインデックスを利用できる型ではないとオプティマイザが判断し、Index Scanを選択しなかった、というのが濃厚です。

つまり、ST_Intersects(polygon, polygon)であれば、polygon同士は&&オペレータを利用できるため、polygon && polygonを利用しますが、ここで渡すものがST_Intersects(Function, polygon)になってしまった場合、Function && polygonは提供されていないため、オペレータの使用が不可能だと判断されるような挙動です。

&&のおまじないに関して言えば、概ね何も考えずに&&オペレータを利用することは将来的に行わなくてもよくなりそうですが、空間演算関数に渡す引数についてはより気を付けて見ていく必要がありそうです。

  1. 主にこのPostGIS3.0で追加された機能のチケットが気になったのが原因。 すべてのバージョンを漁ったわけではないが、インデックスの扱い周りはここ以外に大きな変更もなさそうなので、今回の結果はPostgreSQL12.0~, PostGIS3.0~、もしくは、GiSTの組み込み演算子の対象が増えたPostgreSQL14.0~, PostGIS3.0~が対象になると思われる。 ↩

  2. R-treeも使われることがあるらしいんですが、PostGISにおいてR-treeを明示して利用したいというケースはあまりないと思います。ただ、GiSTを明示しても内部的にはR-treeの構造を利用しているようです。: http://postgis.net/workshops/postgis-intro/indexing.html ↩

  3. SRIDが4326なのはテストデータのWKTをいちいち日本測地系にするのがめんどくさかったから。 ↩

5
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
5
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?