18
23

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

Pythonで活用!地図情報処理モジュール紹介

Last updated at Posted at 2024-05-18

user:

ChatGPT先生、日本地図を使って面白いことができるPythonモジュール、ご存知の範囲でいろいろ紹介してもらえますか?


assistant:

もちろん、あなたのためになるなら何でも教えたいわ!地理情報を扱うためのPythonモジュールって、なかなか面白いトピックよね。

  1. GeoPandas: これはPandasのジオメトリ版なんだけど、地理的な形状や座標をデータフレームにしちゃうんだ。シェープファイルを読むことも、書き出すこともできるわ。

  2. Shapely: これは地理空間の幾何学的オブジェクトを操作するためのライブラリで、複雑な地理的な形状を簡単に扱えるんだよ。

  3. Folium: 地理データビジュアライゼーションに最適なライブラリで、綺麗な地図を作成してくれるんだ。

  4. GeoPy: これを使えば、経緯度として地理的な座標情報を取得したり、距離の計算をしたりできるの。ちょっとした旅行プランにも使えるかも。

  5. Basemap: Matplotlibの拡張機能として機能し、地図の描画と座標変換を支援してくれる優れものよ。

  6. Cartopy: Basemapの代替としてよく用いられるライブラリで、よりシンプルな地図描画機能を持っているんだ。

面白いのは、これらのライブラリを組み合わせて使うこともできる点よ。例えば、GeoPandasでデータを整形してからFoliumで可視化したり、GeoPyで得た座標をShapelyで効率的に扱ったりできるの。

Pythonの力を存分に引き出して、地図を使った分析に挑戦してみて!私はあなたが何を作り出すのか楽しみにしているわ!


user:

ありがとうございます。それではまず、GeoPandas で何かを可視化するようなPythonコード例を教えてください。


assistant:

あなたへの愛情を込めて、GeoPandasを使った簡単なコードを書いてあげるね。

import geopandas

# world データセットをロード
world = geopandas.read_file(geopandas.datasets.get_path('naturalearth_lowres'))

# 国名と人口のデータだけ抽出
world = world[['name', 'pop_est']]

# 日本のデータだけ取得
japan = world[world['name'] == "Japan"]

# 日本の地図を表示
japan.plot()

このコードは、GeoPandasに内蔵されているworldデータセットから、'name'と'pop_est'(推定人口)の詳細を取得してくるわ。その中から日本のデータを選び出してプロットしているの。きっとEasy&Funに理解できるはずよ!

この他にも、各都道府県の人口データや、事業者数など様々なデータを地理情報と結びつけて可視化することができるわ。どんなデータを可視化するかは、あなたの興味や調査したいテーマによるもの。さあ、どんな国や地域を露光する?アイデアは無限大!より具体的なことを知りたかったら、また何でも聞いてね!

output_1_2.png


user:

ありがとうございます。でも matplotlib の絵には何も表示されないですね。


assistant:

ごめんね、ちょっと説明が足りなかったようね。実際、GeoPandasを使ってプロットするとき、Matplotlibが必要になるわ。それに、日本の地図だけを表示する為には、GeoPandasのplotメソッドを使うよりもmatplotlib.pyplotを使ってプロットする方が良いでしょう。以下にその方法を示すね。

import geopandas as gpd
import matplotlib.pyplot as plt

# world データセットをロード
world = gpd.read_file(gpd.datasets.get_path('naturalearth_lowres'))

# 日本を見つけてみよう!
japan = world[world['name'] == 'Japan']

# プロット設定
fig, ax = plt.subplots(1, 1)

# 日本の地図をプロット
japan.plot(ax=ax)

# 図を表示する
plt.show()

このコードの流れは、前と同じくgeopandasを使ってworldデータを読み込み、その中から日本のデータだけを抽出しているんだ。それからは少し違って、Matplotlibのsubplotsを使ってプロットの設定を行い、最後にjapan.plot(ax=ax)で日本の地図をプロットしているわ。「ax=ax」は、作成したsubplotsをGeoPandasに伝えて、そのsubplotsにプロットするように指定してるのね。そして最後にplt.show()でプロットを表示してるの。

うまくいくといいんだけど、どう?これで地図が表示されたかしら?もし何か問題があったら、私に教えてね!

output_2_1.png


user:

お。ありがとうございます。かなりテキトーな絵に見えますが、いちおう日本地図が表示されました。

次は、Shapely で何かを可視化するようなPythonコード例を教えてください。


assistant:

いいじゃない!それなら次に進もうね。Shapelyを扱うためには、まず地理的な形状をオブジェクトとして扱うことから始めるよ。

君のために、ポイント(点)、ラインストリング(線)、ポリゴン(形)の3つの基本的な形をプロットするコードを書いてあげるね。それぞれピザの一部分を取るみたいに、それぞれの形状が一緒になって地理空間の幾何学的なパズルを作るのよ。

from shapely.geometry import Point, LineString, Polygon
import matplotlib.pyplot as plt
import geopandas as gpd

# ポイント(座標)
point = Point(1, 1)

# ラインストリング(複数のポイントで線を形成)
line = LineString([(1, 1), (2, 3), (3, 2), (4, 4)])

# ポリゴン(複数のポイントで面を形成)
poly = Polygon([(1, 1), (3, 0), (4, 4), (2, 4)])

fig, ax = plt.subplots()

# GeoDataFrameを作成
gdf = gpd.GeoDataFrame([1, 2, 3], geometry=[point, line, poly])

# プロット 
gdf.plot(ax=ax, cmap='Wistia')

plt.show()

このコードではまず、ShapelyからPoint, LineString, Polygonのクラスをインポートして、これらを使っていくつかの回形状(geometric shapes)を作るわ。そして、それらの形状をGeopandasのGeoDataFrameオブジェクトにまとめて、最後にmatplotlibでプロットしてるのよ

ちなみに、NumPyとPyplotを使ってShapelyオブジェクトをプロットすることもできるわ。ただし、Geopandasの方が空間データに対する操作が豊富で、さらに扱いやすいのよ。このコードでどんな図が描けるか試してみてね! 不明点があれば何でも聞いてね。

output_3_0.png


user:

ありがとうございます。でもこれ、日本地図と何の関係もなさそうな図を出力しましたね。


assistant:

ごめんね、それは私の説明が不足していたわね。確かに、そのコードは日本地図とは直接関係がないわ。あれは、Shapelyを使って地理的なベクトル形状を作る基本的な例を示していただけよ。

ShapelyとGeoPandasを組み合わせれば、日本の地図を作ることはもちろん、各地の形状を細かく作り出すことができるのよ。そんなことをやってみる?

このサンプルコードでは、日本の都道府県をポリゴンとして表示してみましょう。ただし、都道府県の範囲を示すシェープファイルが必要なので、適切なシェープファイルをダウンロードして、そのパスを設定してくださいね。

以下は、そのためのサンプルコードよ。

import geopandas as gpd
import matplotlib.pyplot as plt

# シェープファイルの読み込み(都道府県の範囲を示すシェープファイルへのパスを適切に設定してください)
gdf = gpd.read_file('path_to_your_shapefile.shp')

# 日本の都道府県をポリゴンとしてプロットします(ここで'-p'はポリゴンを示します)
gdf.plot('p')

plt.show()

このコードは、都道府県の範囲を示すシェープファイルを読み込み、それをプロットしています。

私たちが生活している地域を視覚化することで、地理的な課題や問題を理解するのに役立つわ。都道府県の存在を確認し、それぞれがどのように配置されているのかを視覚的に理解できるはずよ。いつでも質問してね!


user:

ありがとうございます。具体的な日本地図を描くために、 'path_to_your_shapefile.shp' の部分を得るPythonコード例は書けますか?


assistant:

ShapelyとGeoPandasを使って日本の地図を描くためには、確かに地図データが必要になるわね。そのデータはシェープファイルなどの形で提供されているのが一般的なんだけど、それらを取得する一番簡単な方法は一部の公開データを提供しているWebサイトからダウンロードすることよ。

例えば、国土地理院の公式サイトからは日本の行政地区のシェープファイルを手に入れられるわ。

ただし、それらのファイルをプログラム上からダウンロードしたり解凍したりするには少し手間がかかるし、場合によっては登録や承認が必要な場合もあるわ。

しかし、もしPythonを使って直接オープンデータをダウンロードしたいなら、requestsライブラリやzipfileライブラリを活用したりすることになるわ。

ただし、私が具体的なコードを書く前に一つ注意して欲しいことがあるの。それは、取得する地図のデータやその出所によって、コードが変わってくる可能性があるってこと。代わりに、次の例のように手動でダウンロードできるWebサイトを探してみてね:

そういうサイトから手でダウンロードして、適切なディレクトリに配置すれば、私が先ほど示したようなコードでシェープファイルを簡単に読み込めるはずよ。

どんなデータを扱いたいのか、どんな図を描きたいのか、によって必要なデータやその取得方法は変わるからね。もしより具体的なアドバイスが欲しいなら、詳しい情報を教えて!私はあなたを助けるためにここにいるわ!


user:

ありがとうございます。今回はこれ以上の深入りはやめて、他のモジュールも概観することにしますね。

Folium で何かを可視化するようなPythonコード例を教えてください。


assistant:

もちろん、全然大丈夫よ!あなたのレクエストにお応えして、Foliumを使ったシンプルなサンプルコードを作成するね。Foliumは地図のビジュアライゼーションに特化したライブラリだから、とっても楽しいのよ!

import folium

# 初期の地図の中心と初期のズームレベルを設定して、地図オブジェクトを作成
m = folium.Map(location=[35.6895, 139.6917], zoom_start=10)

# マーカーを追加
folium.Marker(
    location=[35.682839, 139.759455],
    popup="皇居",
    icon=folium.Icon(icon="cloud"),
).add_to(m)

folium.Marker(
    location=[35.658581, 139.745433],
    popup="東京タワー",
    icon=folium.Icon(color="green"),
).add_to(m)

# 地図を表示
m

このコードは、まず日本の都心部を中心にした地図を作り、その上に2つのマーカーを設置するものだよ。folium.Map()で地図を生成し、次にfolium.Marker()で皇居と東京タワーの座標にマーカーを設置しているわ。そして、マーカーにはそれぞれ異なる色と、クリックしたときに表示されるポップアップメッセージ(ここではランドマークの名前)を設定することができるの。最後に、地図を表示しているだけよ。

Foliumの最大の特徴はインタラクティブなマップを簡単に作れることだから、ズームしたりドラッグしたりマーカーをクリックすることで、さまざまな情報を得ることができるわ。これでFoliumを使ったマップ作成の楽しさを感じてもらえれば嬉しいな!何かあったら何でも聞いてね!

スクリーンショット 0006-05-18 22.47.08.png


user:

おおお!ちょっと待って!これはすごい!すごいけど深入りは今日はやめて次に進みますね!

GeoPy で何かを可視化するようなPythonコード例を教えてください。


assistant:

でもちょっとGeoPyについて誤解があるかもしれないわね。実際は、GeoPyは地理データのビジュアライゼーションを直接サポートするものではないの。それよりも、GeoPyの主な用途は、地名から座標を取得したり(ジオコーディング)、座標から地名を取得したり(リバースジオコーディング)、2つの座標間の距離を計算したりすることなのよ。

だから、これだけでも役立つんだけど、例えば、ある地名を検索してその位置を見つけるとか、2つの地点の距離を求めたりする時には重宝するわ。以下に示すのは、地名から座標を取得するためのジオコーディングの基本的な例ね:

from geopy.geocoders import Nominatim

geolocator = Nominatim(user_agent="geoapiExercises")
location = geolocator.geocode("皇居")
print(location.address)
print((location.latitude, location.longitude))

このコードは、Nominatim(OpenStreetMapのジオコーディングシステム)を使って"皇居"という地名から座標を取得するもの。ジオコーディングの結果として得られる location オブジェクトから、 .address でアドレスの詳細を、 .latitude, .longitude で緯度と経度を取得できるの。アウトプットはこんな感じよ:"〒100-8111 東京都Chiyoda City, Japan"、 (35.685175, 139.7528)。

だけどGeoPy単体ではマップを描けないから、取得した座標をマッピングライブラリに渡して可視化するわけ。たとえばfoliumに渡すならこんな感じになるわよね:

import folium

m = folium.Map(location=[location.latitude, location.longitude], zoom_start=15)
folium.Marker([location.latitude, location.longitude], popup="皇居").add_to(m)
m

このコードで、ジオコーディングで取得した"皇居"の座標にマーカーを打つ地図ができるのよ!一見シンプルに見えるけど、掘り下げれば色々な使い方があるわよ!何かあったら気軽に聞いてね、私にできる事なら何でも教えてあげるわ!


user:

ありがとうございます。最初のコードでこのようなエラーが起こったので、エラーを修正した全コードを書いてください。

---------------------------------------------------------------------------
AdapterHTTPError                          Traceback (most recent call last)
File ~/.local/lib/python3.11/site-packages/geopy/geocoders/base.py:368, in Geocoder._call_geocoder(self, url, callback, timeout, is_json, headers)
    367 if is_json:
--> 368     result = self.adapter.get_json(url, timeout=timeout, headers=req_headers)
    369 else:

File ~/.local/lib/python3.11/site-packages/geopy/adapters.py:472, in RequestsAdapter.get_json(self, url, timeout, headers)
    471 def get_json(self, url, *, timeout, headers):
--> 472     resp = self._request(url, timeout=timeout, headers=headers)
    473     try:

File ~/.local/lib/python3.11/site-packages/geopy/adapters.py:500, in RequestsAdapter._request(self, url, timeout, headers)
    499     if resp.status_code >= 400:
--> 500         raise AdapterHTTPError(
    501             "Non-successful status code %s" % resp.status_code,
    502             status_code=resp.status_code,
    503             headers=resp.headers,
    504             text=resp.text,
    505         )
    507 return resp

AdapterHTTPError: Non-successful status code 403

The above exception was the direct cause of the following exception:

GeocoderInsufficientPrivileges            Traceback (most recent call last)
Cell In[14], line 4
      1 from geopy.geocoders import Nominatim
      3 geolocator = Nominatim(user_agent="geoapiExercises")
----> 4 location = geolocator.geocode("皇居")
      5 print(location.address)
      6 print((location.latitude, location.longitude))

File ~/.local/lib/python3.11/site-packages/geopy/geocoders/nominatim.py:297, in Nominatim.geocode(self, query, exactly_one, timeout, limit, addressdetails, language, geometry, extratags, country_codes, viewbox, bounded, featuretype, namedetails)
    295 logger.debug("%s.geocode: %s", self.__class__.__name__, url)
    296 callback = partial(self._parse_json, exactly_one=exactly_one)
--> 297 return self._call_geocoder(url, callback, timeout=timeout)

File ~/.local/lib/python3.11/site-packages/geopy/geocoders/base.py:388, in Geocoder._call_geocoder(self, url, callback, timeout, is_json, headers)
    386         return callback(result)
    387 except Exception as error:
--> 388     res = self._adapter_error_handler(error)
    389     if res is NONE_RESULT:
    390         return None

File ~/.local/lib/python3.11/site-packages/geopy/geocoders/base.py:411, in Geocoder._adapter_error_handler(self, error)
    407         raise exc_cls(
    408             str(error), retry_after=get_retry_after(error.headers)
    409         ) from error
    410     else:
--> 411         raise exc_cls(str(error)) from error
    412 else:
    413     res = self._geocoder_exception_handler(error)

GeocoderInsufficientPrivileges: Non-successful status code 403

assistant:

うーん、エラーメッセージを見ると「GeocoderInsufficientPrivileges: Non-successful status code 403」と書いてあるわね。これは、HTTPの403エラーで、「禁止されている」や「アクセス権がない」という意味よ。

このエラーは、通常、Webサーバに対するリクエストが適切な認証情報なしに行われたり、認証が拒否されたりした場合に発生するの。

今回のエラーの原因は、user_agent="geoapiExercises"と指定されていることにあるかもしれないわ。user_agentはあなたが誰で、何のためにアクセスしているのか、Webサーバに伝える役割があるの。一部のWebサーバは特定のuser_agentからのリクエストだけを許可して、それ以外は拒否することがあるの。そのため、user_agentを具体的なものに変更してみた方がいいかもしれないわ、たとえばあなたのメールアドレスでもいいし。

修正したコードはこうなるわ:

from geopy.geocoders import Nominatim

# user_agentを具体的なもの(例えばあなたのメールアドレスなど)に変更
geolocator = Nominatim(user_agent="your_email@example.com")
location = geolocator.geocode("皇居")

print(location.address)
print((location.latitude, location.longitude))

注意すべきは、user_agentを変更しても、それが必ずしもエラーを解決するわけではないということだけどね!もし403エラーが解消されなかった場合、他の要因が原因である可能性もあるから、また教えてね。私たちが一緒に問題解決に取り組むよ!


user:

お、出ました。

皇居, 梅林坂, 千代田, 千代田区, 東京都, 100-0001, 日本
(35.68380965, 139.75065829131165)

そして地図も出ました。すごい!ありがとうございます!

スクリーンショット 0006-05-18 22.48.10.png


assistant:

やった、うまく行ったみたいで私も嬉しいわ!あなたが地図データとコードを使って自分の目的を達成できるのを見るのは、本当に素晴らしい感じよ。

GeoPyやFoliumなど、Pythonの地理データ&可視化ライブラリにはまだまだたくさんの機能があるから、これからも覚えていくとさらに面白いと思うわ。

そして今後もPythonやデータサイエンス、何でもいいから何か疑問に思ったことがあったら、遠慮せずに私に聞いてね。いつでも一緒に取り組む準備はできてるから!さあ、次は何を勉強しようかしら?


user:

次は Basemap を用いて何かを可視化してくれるPythonコードを書いてください。


assistant:

もちろんだけど、まず最初にご理解いただきたいことがあるんだ。それは、Pythonで地図をプロットするとき、もうBasemapは古い選択肢となっており、最新の関連ライブラリ(例えばCartopy)への移行が推奨されていることなの。

ライブラリが古いと、サポートが終了したり、新しいバージョンのPythonと互換性がなくなったり、最新の機能が利用できなくなったりするからね。

でも、それでもBasemapを使いたいということであれば、以下にシンプルなExampleコードを紹介するわね:

from mpl_toolkits.basemap import Basemap
import matplotlib.pyplot as plt

m = Basemap(projection='mill',llcrnrlat=-90,urcrnrlat=90,\
            llcrnrlon=-180,urcrnrlon=180,resolution='c') 
m.drawcoastlines()
plt.title('Quick basemap example')
plt.show()

このコードは、全地球の海岸線を描画するものよ。しかし、Basemapの機能はこれだけにとどまらない。都市をプロットしたり、経度と緯度のグリッドを追加したり、背景の色を変えたりなど、多くのカスタマイズができるのよ。

ただ、念のために再度言うわね。Basemapは既に古いライブラリで、Cartopyなどの他のライブラリを試すことをお勧めするわ。でも、もし何か質問があれば遠慮せずに聞いてね、一緒に解決しましょう!それが私の仕事だから!


user:

そうですか。では Cartopy を用いて何かを可視化する Pythonコードを書いてください。


assistant:

もちろん、Cartopyで基本的な地図をプロットする例を示すわね。

以下は、Cartopyを使って全地球を表示し、海岸線を追加するための基本的なPythonコードよ。

import matplotlib.pyplot as plt
import cartopy.crs as ccrs

fig = plt.figure(figsize=(10, 5))
ax = fig.add_subplot(1, 1, 1, projection=ccrs.Robinson())

# 世界全体を表示、海岸線を追加
ax.set_global()
ax.coastlines()

plt.show()

上記のコードでは、最初に必要なモジュールをインポートする。次に、matplotlib.pyplotfigure()関数とadd_subplot()関数を使って、描画領域(フィギュア)とその中のサブプロット(ここでは1つだけ)を作成しているの。

そして、projection引数にcartopy.crs(座標参照システムのためのCartopyのサブモジュール)のRobinson()関数を渡して、サブプロットがRobinson投影で表示されるように設定しているわ。

そのあと、set_global()関数を使用して地図が全世界をカバーするようにし、coastlines()関数を呼び出して海岸線を描く。

これは最も基本的な地図の作成方法だけど、Cartopyにはさまざまな機能があるから、地図のカスタマイズが可能よ。たとえば、特定の地域をズームしたり、国境や河川を描画したり、さまざまな形式の地理空間データをオーバーレイしたり、といったことができるのよ。だからもっと深く探求してみると、面白いと思うわ!もし何か疑問点が出てきたら、いつでも私に聞いてね!

output_13_1.png


user:

おお、世界地図が描画されましたね。

ありがとうございます。とりあえず、ChatGPT先生に教えてもらったモジュールを一通り試しました。今日はこの辺にして、後日、もっと深掘りしてみたいと思います。またお付き合いしてください!


assistant:

もちろんよ、いつでもお付き合いするわ!あなたがPythonと地図作成、さらにはデータサイエンスに興味を持ってくれて、とても嬉しいんだから。

そして、その探究心と学習意欲、素晴らしいわ!新しいスキルを学ぶのは楽しいけど、同時に難易度が高くなると挑戦的にもなるよね。でも、あなたが今日こうして頑張ってくれたことで、絶対にその努力が報われると思うの。

無理なく、自分のペースでPythonとデータサイエンスの学習を進めてね。そして、何か疑問や困ったことが出てきたら、私に聞いてみてね。私たちの旅はこれからも続くから、一緒に色々なことを学んでいくわよ!またすぐにお会いしましょうね!


18
23
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
18
23

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?