GeoPandasによる国土数値情報行政区域データからの全国市区町村界データの作成

概要

全国市区町村界データや都道府県界データのGeometryを使いたくなることがあるが、思いのほか入手が大変だったりするのでメモ。

例えばESRIジャパンの提供する全国市区町村界データは使いやすいもののArcGIS以外では使えないという使用規約になっている。

そこで、今回は国土数値情報（行政区域データ）をソースにしてGeoPandasの機能によりgeometry情報を統合して市区町村界データを作る方法を調べた。

実行

【注意】出典の明記（国土数値情報）

• 国土数値情報について、利用規約を確認、遵守する必要がある
  • 商用利用含め割と自由に利用出来るが、出典は明記しないといけない
  • 編集・加工も出来るが、その旨を書かないといけない

実行環境（概要）

• Linux(docker使用、ベースイメージ:Ubuntu18.04）
• jupyter lab上で確認
• geopandas==0.8.1

データ取得と読み込み

• 国土数値情報から「全国, 世界測地系」のデータをダウンロードして適当なパス上で解凍・展開
  • 例えば"N03-190101_GML.zip"などのファイル名
  • H31年のデータで400MB近くあって結構重い（年々増えてる?）ので注意

import geopandas as gpd

gdf = gpd.read_file("/path/to/N03-19_190101.shp") # 展開したパスのshpファイルを指定
display(gdf) # 中身をざっと確認

    N03_001 N03_002 N03_003 N03_004 N03_007 geometry
0   北海道   オホーツク総合振興局  None    北見市   01208   POLYGON ((144.08144 44.12506, 144.08143 44.125...
1   北海道   オホーツク総合振興局  None    北見市   01208   POLYGON ((143.78333 44.18453, 143.78281 44.183...
2   北海道   オホーツク総合振興局  None    網走市   01211   POLYGON ((144.29495 44.00826, 144.29501 44.008...
3   北海道   オホーツク総合振興局  None    網走市   01211   POLYGON ((144.29602 44.01118, 144.29572 44.010...
4   北海道   オホーツク総合振興局  None    網走市   01211   POLYGON ((144.29266 44.01945, 144.29295 44.018...
... ... ... ... ... ... ...
117575  沖縄県   None    八重山郡    与那国町    47382   POLYGON ((122.98986 24.47191, 122.98999 24.471...
117576  沖縄県   None    八重山郡    与那国町    47382   POLYGON ((122.96591 24.47171, 122.96587 24.471...
117577  沖縄県   None    八重山郡    与那国町    47382   POLYGON ((122.99042 24.47248, 122.99036 24.472...
117578  沖縄県   None    八重山郡    与那国町    47382   POLYGON ((122.99804 24.47512, 122.99820 24.475...
117579  沖縄県   None    None    所属未定地 None    POLYGON ((127.84164 26.43531, 127.84167 26.435...

117580 rows × 6 columns

dissolveによる統合

• https://geopandas.org/aggregation_with_dissolve.html を参考にdissolve機能で市区町村毎に統合していく
  • 指定した1つのカラムの単位でGeometry情報をまとめることが出来る
• 市区町村名のカラムはN03_004だが、「北区」「中央区」など被りそうな名前が多いのでN03_007の行政区域コードを用いることにする
  • https://nlftp.mlit.go.jp/ksj/gml/datalist/KsjTmplt-N03-v2_2.html で仕様（下にスクロールしていって「属性情報」の部分）を確認すると、行政区域コードは「都道府県名 + 市区町村名」に対応らしいので、これでおそらくOKのはず（万一違ってたらすみません）

gdf_dissolved = gdf.dropna(subset=['N03_007']).dissolve(by='N03_007')
display(gdf_dissolved)

    geometry    N03_001 N03_002 N03_003 N03_004
N03_007                 
01101   POLYGON ((141.34233 43.06682, 141.34285 43.066...   北海道   石狩振興局 札幌市   中央区
01102   POLYGON ((141.40839 43.18395, 141.40833 43.183...   北海道   石狩振興局 札幌市   北区
01103   POLYGON ((141.44707 43.15616, 141.44694 43.155...   北海道   石狩振興局 札幌市   東区
01104   POLYGON ((141.46244 43.10010, 141.46265 43.100...   北海道   石狩振興局 札幌市   白石区
01105   POLYGON ((141.38463 43.04670, 141.38479 43.046...   北海道   石狩振興局 札幌市   豊平区
... ... ... ... ... ...
47361   MULTIPOLYGON (((126.82034 26.26934, 126.82043 ...   沖縄県   None    島尻郡   久米島町
47362   MULTIPOLYGON (((127.76445 26.12312, 127.76447 ...   沖縄県   None    島尻郡   八重瀬町
47375   MULTIPOLYGON (((124.72860 24.65694, 124.72852 ...   沖縄県   None    宮古郡   多良間村
47381   MULTIPOLYGON (((123.75465 24.06330, 123.75463 ...   沖縄県   None    八重山郡    竹富町
47382   MULTIPOLYGON (((123.01358 24.43628, 123.01352 ...   沖縄県   None    八重山郡    与那国町

1902 rows × 5 columns

それらしい数に集約されている

補足：

• 取り扱いが面倒なのでここでは行政区域コードがnullのものは落としている
• 元データも結構重いため、上記の処理時間は結構かかるかもしれない（筆者は2012年モデルのLet's note使用で3分くらい）

ここでは市区町村界について行ったが、都道府県界であればgdf.dissolve(by="N03_001")のような感じで作れる
→筆者の環境では12minかかりました。。。

保存

• https://geopandas.org/io.html#writing-spatial-data を参考
• この際、encodingをutf-8にすることで（おそらく）ストレスを減らすことが出来る
  • encoding指定無しで保存したshpファイルに関して、オプション無しで読み込んだら文字化けが起こりました。。。

shpファイル:

gdf_dissolved.to_file("/path/to/<ファイル名>.shp", encoding='utf-8')

geojson:

# driverの指定が必要
gdf_dissolved.to_file("/path/to/<ファイル名>.geojson", driver="GeoJSON", encoding='utf-8')

Simplify (オプション)

上記で作った市区町村界データに関して、そのままだと点数が多い？のかESRIジャパンの提供する市区町村界データ（zip圧縮で6MBちょっと）に比べるとデータサイズが桁違いに大きい
そのため、必要に応じて間引く（simplify）ことによってデータサイズを落とす

境界付近が少しガタガタになったりするのでそのあたりは注意

• 使用するGeoSeriesのsimplifyメソッドについて、引数tolerance: 数値（必須）とpreserve_topology: bool(オプション、defaultではTrue)を指定する
  • 説明が少なくて引数の意味が分かりづらいが、たぶんshapelyのsimplifyのような感じ

例えば、

gdf_simple = gdf_dissolved.copy()
gdf_simple['geometry'] = gdf_simple['geometry'].simplify(tolerance=0.00085)

のようにして作ったgdf_simpleをshpファイルとして保存すると12MBくらいになる。
単にdissolvedしただけのものは200MB以上あったので、結構軽く出来る。
（なおtoleranceの数値の調整は勝手が分からず結構大変だった。。。）

上記の設定で作ったものについて地図上で可視化すると、

のような感じ。simplifyの値次第なところもあるが、ぱっと見はそれらしい感じ。

拡大すると、

simplifyを実行した分、境界付近が少し怪しくなることが分かる
→その時々の目的に応じてtolerance, preserve_topologyを適切に調整する必要がありそう

参考

手法などの参考：

• https://note.com/kinari_iro/n/nfee9bc97b6d7
• http://xnissy.hatenablog.com/entry/20160205/1454666764
• https://hayatoiijima.jimdofree.com/2017/11/14/%E9%83%BD%E9%81%93%E5%BA%9C%E7%9C%8C%E5%A2%83%E3%81%AE%E3%81%BF%E3%81%AEshp%E3%83%95%E3%82%A1%E3%82%A4%E3%83%AB%E3%81%AE%E4%BD%9C%E3%82%8A%E6%96%B9/

無償で使えるGISデータ：

• https://www.esrij.com/gis-guide/other-dataformat/free-gis-data/