キヨシの命題

データサイエンス（統計解析や機械学習など）に関する情報を発信します

NEW ARTICLE
M1チップ搭載のMacbookでRedash Dockerコンテナを起動する
  1. HOME
  2. エンジニアリング
  3. RedShiftのGEOMETRY型を使って緯度経度から都道府県を判定してみる

RedShiftのGEOMETRY型を使って緯度経度から都道府県を判定してみる

RedShiftのGEOMETRY型を使って緯度経度から都道府県を判定してみる

最近TwitterでJTCの内情のツイートばかりをしていて、技術的な発信をおろそかにしていました…。反省です。
ということで、久々のブログ投稿です。

仕事上よく緯度・経度を持つ空間データの分析をするのですが、毎回DWH(AWS RedShift)からPython実行環境にロードし、GeoPandasやShapelyを使って処理していました。

ただ、結局Pythonで処理した結果をまたDWHに戻してBIで集計や可視化をしたり、後工程の作業につなげたりしています。
DWH上で空間データの処理を完結できると途中でPythonを挟む必要がなくなるため無駄な工数やコストを削減できます。

そこで、RedShiftではGEOMETRY型に対して空間関数がサポートされているようなので試してみました。

今回は、RedShiftにGEOMETRY型を含むテーブルを作り、空間関数を使って緯度経度から都道府県を判定してみます。

都道府県ごとのポリゴンデータはこちらで公開されているデータセットを使用します。

japonyol.net
 
47都道府県のポリゴンデータ geojson
https://japonyol.net/editor/article/47-prefectures-geojson.html
市区町村の行政区域データは国土数値情報にあるから、これで事足りる（使用例＝日本全国市区町村の行政区域）。 都道府県のポリゴンデータは案外少なかったので5年前に KML で出していたが、もう要らんだろうと考え削除していた。 それでもやはりデー...

目次

環境

GEOMETRY型を含むテーブルをRedShiftに作成

RedShiftでGEOMETRY型のカラムを持つテーブルを作成する方法はいくつかありますが、AWS公式ドキュメントではCOPYコマンドが最も効率的な方法と紹介されています。

今回は公式ドキュメントの説明に従い、ポリゴンを含むデータセットをS3に格納し、S3からCOPYコマンドでテーブルを作ることにします。

まず、次のPythonコードでこちらから都道府県ごとのポリゴンを保持するgeojsonを取得し、S3にアップロードします。

import geopandas as gpd
from shapely.geometry import Point

path = "https://japonyol.net/editor/article/prefectures.geojson"
prefgeorange = gpd.read_file(path,encoding='SHIFT-JIS')

d = prefgeorange.to_wkb(hex=True)

# upload to S3
bucket_name = "xxxxx"
file_name = "xxxxx"
s3 = boto3.resource("s3")
obj = s3.Object(bucket_name, file_name)
obj.put(Body=data.to_csv(None, index=False, sep="\t"))

ここで、GEOMETRY型にしたいカラムについてはWKB形式もしくはWKT形式でデータを保持しておく必要があります。
GeoPandasでgeojsonを読み込むとShapleyのMultiPolygon型でポリゴンデータがcastされるため、このままファイル出力してS3にアップロードはできません。
今回はWKB形式でMultiPolygon型を文字列に変換するようにしています。

S3にファイルアップロード後、次のCOPYコマンドを実行します。
COPYコマンドを使うために必要なiam roleはこちらを参照してください。

CREATE TABLE sandbox.m_pref_georange(
  pref_name VARCHAR(256) PRIMARY KEY NOT NULL,
  geo_range GEOMETRY
);

copy sandbox.m_pref_georange
from 's3://bucket_name/file_name' 
iam_role 'arn:aws:iam::xxxxxxxxx:role/xxxxxxxxxxxxxxxx'
format as csv
delimiter '\t'
ignoreheader 1
;

この他にCOPYコマンドを経由せず、対象のカラムをWKT型もしくはWKB型に文字列変換してテーブル作った後にそのカラムをGEOMETRY型に変換したテーブルを作成する方法もあります。

RedShiftの空間関数を使ってみる

GEOMETRY型を含むテーブルを作成できたのでいよいよ空間関数を使用してみます。
今回はGeolonia社が公開している全国の町丁目、大字、小字レベルの緯度経度を持つ住所データに対して、空間関数を使って都道府県を判定してみます。
こちらの住所マスタには既に都道府県カラムが含まれているので、都道府県のポリゴンから空間関数で判定した都道府県と一致するかどうかを確かめられます。

まず、住所マスタに都道府県ごとのポリゴンを結合し、付与した都道府県が正しいかを集計します。

select
    case when geolonia.pref = polygon.pref_name then 1 else 0 end as is_correct
    ,count(1) as cnt
from
    (
    select
        都道府県名 as pref
        ,市区町村名 as city
        ,大字町丁目名 as town
        ,緯度 as latitude
        ,経度 as longitude
    from sandbox.m_geolonia
    ) geolonia	
    inner join
        (
        select
            geometry
            ,pref_name
        from scraping.m_pref_georange
        ) polygon
        on
            ST_Within(
                ST_Point(cast(geolonia.longitude as double precision),cast(geolonia.latitude as double precision))
                ,polygon.geometry
            )
group by
	case when geolonia.pref = polygon.pref_name then 1 else 0 end
;

OUTPUT:

ここで、SQLで参照している空間関数について解説します。
他にも空間関数はいろいろあります。詳しくはこちらを参照してください。

関数名説明
ST_Point数値型の緯度、経度をもとにPointというGEOMETRY型に変換します
ST_WithinPointとPolygonを引数にPointが含まれればTrue、含まれなければFalseを返します。これをJOIN条件に入れることで都道府県の判定ができます。

先ほどのSELECT結果を見ると、17件だけポリゴンから付与した都道府県と住所マスタの都道府県が一致していませんでした。
都道府県が一致していなかったデータを確認します。

select
    geolonia.*
    ,polygon.pref_name
from
    (
    select
        都道府県名 as pref
        ,市区町村名 as city
        ,大字町丁目名 as town
        ,緯度 as latitude
        ,経度 as longitude
    from sandbox.m_geolonia
    ) geolonia	
    inner join
        (
        select
            geometry
            ,pref_name
        from scraping.m_pref_georange
        ) polygon
        on
            ST_Within(
                ST_Point(cast(geolonia.longitude as double precision),cast(geolonia.latitude as double precision))
                ,polygon.geometry
            )
where
	geolonia.pref <> polygon.pref_name

OUTPUT:

試しに先頭のレコードについて実際にどの緯度経度のデータで誤ったのかを確認してみます。
「宮城県 刈田郡七ヶ宿町 字杉沢口」は福島県と判定されていますが、Google Mapで住所マスタの緯度経度を検索したところ、福島県と宮城県の県境付近でした。

これは都道府県別のポリゴンの精度によるもので、より緻密なポリゴンを持つマスタを参照すればほぼ完ぺきに緯度経度から都道府県を判定できると思います。
ただ、ポリゴンを緻密にすればするほどデータ量が増えるので処理コストも増えます。
データの扱いやすさとポリゴンの精度はトレードオフの関係にあります。

最後に

RedShiftの空間関数、なかなか便利でした。
ただ、GEOMETRY型を含むテーブルの操作扱いが特殊で、テーブルプレビューのためにselect * from xxx limit 100とかを実行するとエラーが発生します。
空間データを持つカラムならすべてGEOMETRY型にするというのは推奨しませんが、こちらに記載されている空間関数の使用を前提とするような場合であればRedShiftでGEOMETRY型を使うのは多いにアリだと思います。

参考資料