はじめに
はじめまして！データサイエンティストの渡邊です。株式会社GEOTRA(@GEOTRA )にてGIS（Geographic Information System：地理情報システム）データを用いた分析を行っています。
本記事ではgeopandasで座標とポリゴン間の距離を算出する方法を紹介します。
背景
- 「あるエリアから◯km以内の座標を抽出する」といった分析を行いたいときがある
- 単にそのエリアの代表点をとってその点との距離を算出するのもいいが、geopandasを用いれば座標とポリゴン間の距離をすぐに算出できる
やること
-
geopandas.GeoSeries.distanceを用いる- 返り値の単位はメートル
-
POINTとPOLYGON間の距離はもちろん、POINTとMULTIPOLYGON間の距離も算出できる
サンプルコード
import geopandas as gpd # 0.10.2
from shapely.geometry import Point, Polygon, MultiPolygon # 1.8.0
from keplergl import KeplerGl #0.3.2
import random
# POINTの例 皇居周辺のランダムな座標
random_points = [Point(random.uniform(139.7, 139.8), random.uniform(35.65, 35.75)) for _ in range(50000)]
# MULTIPOLYGONの例（ここでは2つのポリゴンを含むMULTIPOLYGONを作成）
polygon1 = [[139.74762183584198,35.69430755236971],[139.75036385428467,35.69512122576469],[139.75231490586677,35.693493870673066],[139.755162386557,35.69220909319844],[139.75521511768082,35.69109560264534],[139.75595335341572,35.6903675427258],[139.75964453208658,35.68972513139477],[139.75964453208658,35.689211198605214],[139.76033003669582,35.68908271489071],[139.7613846591735,35.687283921153686],[139.76059369231666,35.6856135763648],[139.76006638107694,35.68574206566698],[139.759802725458,35.68377187362315],[139.75800986724536,35.68420018037952],[139.75753528713128,35.68300091566584],[139.75579516004245,35.683557719381916],[139.75579516004245,35.68300091566584],[139.75621700903469,35.68282959066395],[139.75611154678708,35.68222995025875],[139.75495146206185,35.6822727818656],[139.7543186885745,35.68017400608707],[139.75532057992837,35.679660011766856],[139.75484599981428,35.678846180655036],[139.75358045284148,35.679060347542176],[139.75025839203707,35.67816084275526],[139.7479909537104,35.68064516463877],[139.7469890623565,35.68064516463877],[139.74614536437403,35.68248693955575],[139.74509074189632,35.68304374685881],[139.74630355774534,35.68672714344338],[139.74603990212643,35.68844029321275],[139.75046931653227,35.68972513139477],[139.7479909537104,35.69066733290892],[139.74762183584198,35.691481043442074],[139.74836007157688,35.692380398057296],[139.74762183584198,35.69430755236971]]
polygon2 = [[139.7494178796471,35.7171301381481],[139.745660622306,35.720104370005345],[139.75026326254869,35.723459780355356],[139.75439624562256,35.720180630900565],[139.7494178796471,35.7171301381481]]
multipolygon = MultiPolygon([Polygon(polygon1), Polygon(polygon2)])
# POINTのジオデータフレーム
gdf_point = gpd.GeoDataFrame(geometry=random_points, crs='epsg:4326')
gdf_point['geometry'] = gdf_point['geometry'].to_crs('epsg:6677') # 平面直角座標への変換
# MULTIPOLYGONのジオデータフレーム
gdf_polygon = gpd.GeoDataFrame(geometry=[multipolygon], crs='epsg:4326')
gdf_polygon['geometry'] = gdf_polygon['geometry'].to_crs('epsg:6677') # 平面直角座標への変換
# 距離(m)の算出
gdf_point['distance'] = gdf_point.distance(gdf_polygon.unary_union)
# keplerglで可視化
kepler_map = KeplerGl(height=800)
kepler_map.add_data(data=gdf_polygon, name='Polygon')
kepler_map.add_data(data=gdf_point, name='Point')
kepler_map
二つのポリゴンからの距離がしっかり算出されているのが確認できる
注意点
座標参照系（CRS）の設定は対象のエリアに合わせて設定する必要がある(東京都はepsg:6677)
# POINTのジオデータフレーム
gdf_point = gpd.GeoDataFrame(geometry=random_points, crs='epsg:4326')
gdf_point['geometry'] = gdf_point['geometry'].to_crs('epsg:6677') # 平面直角座標への変換
# MULTIPOLYGONのジオデータフレーム
gdf_polygon = gpd.GeoDataFrame(geometry=[multipolygon], crs='epsg:4326')
gdf_polygon['geometry'] = gdf_polygon['geometry'].to_crs('epsg:6677') # 平面直角座標への変換
おわりに
最後までお読みくださりありがとうございました。
GEOTRAのQiitaでは、引き続き地理空間情報に関する情報を中心に技術系のテクニックを発信していきます。今後も皆さんのお役に立てるコンテンツを配信できればと思っておりますので、皆様のいいね・フォローをお待ちしています！
また、こちらの記事を読んでGEOTRAに興味を持ってくださった方は、是非Wantedlyから会社紹介記事もご覧ください！GEOTRAは、学生インターンを含め、現在一緒に働く仲間を募集しております。興味をお持ちの方は、Wantedly上の募集ページやストーリー、各種サイトをご確認ください！
企業HP： https://www.geotra.jp/
note： https://note.com/2022geotra/
Wantedly： https://www.wantedly.com/stories/s/geotra_introduction
LinkedIn： https://www.linkedin.com/company/%E6%A0%AA%E5%BC%8F%E4%BC%9A%E7%A4%BEgeotra/