Python Remind Memo - World Heritage Plotting in Python and ChatGPT

-

国土数値情報の地域資源・観光という項目に「世界文化遺産」、「世界自然遺産」というコンテンツがあったので、図化を試みてみました。(URL:https://nlftp.mlit.go.jp/ksj/gml/datalist/KsjTmplt-A34-v1_2.html)

まずは、GIS情報（Shapeファイル）のダウンロードから始めますが、今回は手動でダウンロード→解凍→GoogleDriveにアップロードするのではなく、GoogleClabの中でこれらの手順を済ませてしまいたいと思いました。
調べてみると、GeoPandasにZIPファイルを直接読み込む機能があるようなので、以下のコードを実行しました。
!pip install geopandas
import geopandas
import pathlib
NOTEBOOK_PATH = pathlib.Path().resolve()
DATA_DIRECTORY = NOTEBOOK_PATH / "data"

print("####", DATA_DIRECTORY)

world_heritages = geopandas.read_file('https://nlftp.mlit.go.jp/ksj/gml/data/A34/A34-17/A34-180316_GML.zip')
world_heritages.head()


ちなみに、'https://nlftp.mlit.go.jp/ksj/gml/data/A34/A34-17/A34-180316_GML.zip'は国土数値情報の該当ページのHTMLコードを解読して、ZIPファイルの保存場所を調べました。
しかし、この結果は次のようなErrorになりました。
ERROR:fiona._env:`/vsizip/vsicurl/https://nlftp.mlit.go.jp/ksj/gml/data/A34/A34-17/A34-180316_GML.zip' not recognized as a supported file format.

ChatGPT先生にこのエラーの理由をお尋ねすると、以下のご教示をいただきました。

このエラーは、geopandas.read_file()メソッドが指定されたURLから直接ZIPファイルを読み込むことができないことを示しています。read_file()メソッドは、URLから直接ファイルを読み込む機能を持っていないため、このエラーが発生しました。

ZIPファイルをダウンロードして解凍するためには、先ほど紹介した方法を使用して手動で処理する必要があります。以下に、改めてZIPファイルをダウンロードして解凍する手順を示します。

python
import urllib.request
import zipfile

zip_url = "https://nlftp.mlit.go.jp/ksj/gml/data/A34/A34-17/A34-180316_GML.zip"

# ZIPファイルをダウンロード
urllib.request.urlretrieve(zip_url, "example.zip")

# ZIPファイルを解凍
with zipfile.ZipFile("example.zip", "r") as zip_ref:
    zip_ref.extractall("destination_folder")

上記のコードの最後の行の"destination_folder"の部分を、ZIPファイルを解凍したいディレクトリのパスに変更してください。これにより、ZIPファイルがダウンロードされ、指定したディレクトリに解凍されます。

もし、解凍後のデータを読み込むためにgeopandasを使用したい場合は、解凍されたファイルのパスを指定してgeopandas.read_file()を呼び出すことができます。

python
import geopandas

shapefile_path = "destination_folder/A34-180316_GML/A34-180316.shp"
world_heritages = geopandas.read_file(shapefile_path)
world_heritages.head()

上記のコードでは、解凍されたデータの中にある.shpファイルを指定してgeopandas.read_file()を使用しています。もちろん、解凍されたデータのファイル形式に応じて、適切なファイルを指定してください。


ここでありがちなエラー。
Error:unsupported operand type(s) for +: 'PosixPath' and 'str'
というのが発生しました。こちらもChatGPT先生にご質問したところ、


おっしゃるエラーは、PosixPath型（またはPath型）のオブジェクトと文字列型のオブジェクトを+演算子で結合しようとした際に発生するものです。Pathオブジェクトは、ファイルシステムのパスを操作するためのクラスであり、文字列とは異なる型です。

このエラーを解決するためには、Pathオブジェクトを文字列に変換してから結合する必要があります。以下のように、str()関数を使用してPathオブジェクトを文字列に変換してください。

python
import pathlib

# DATA_DIRECTORYがPosixPathオブジェクトとして定義されている場合
destination_folder = str(DATA_DIRECTORY) + "/A34-180316_GML"
shapefile_path = str(destination_folder) + "/A34-180316.shp"

または、Pathオブジェクトを/演算子で連結する方法もあります。

python
import pathlib

# DATA_DIRECTORYがPosixPathオブジェクトとして定義されている場合
destination_folder = DATA_DIRECTORY / "A34-180316_GML"
shapefile_path = destination_folder / "A34-180316.shp"

いずれの方法でも、Pathオブジェクトを文字列に変換するか、/演算子を使ってPathオブジェクト同士を結合することで、エラーを回避できます。


だそうです。Pathlibで作成したPathhはPosixPath型という特殊なTypeらしく、PosixPath型（またはPath型）のオブジェクトと文字列型のオブジェクトを+演算子で結合できないとのことで、一度PosixPathをstr()で文字列に変換する作業がいるようです。そんな面倒な手続きを踏まずにPathlibが直接文字列を吐き出せばいいような気がしますが、何かりゆうがあるのでしょうね？
ちょっとChatGPT先生に聞いてみました。


pathlibライブラリは、様々な操作を行うためにパスを表現するための高水準なAPIを提供します。その際に、PosixPath型やWindowsPath型など、環境に応じた具体的なパスオブジェクトの型を返すことがあります。

pathlibの設計哲学の一つは、環境に依存しないパス表現を提供することです。つまり、pathlibはWindowsやUnix系の異なる環境で使用されることを想定しており、それぞれの環境に適したパス表現を返すことができます。例えば、PosixPathはUnix系の環境（Linux、macOSなど）で使われることが多く、WindowsPathはWindows環境で使われることが多いです。

pathlibオブジェクトを文字列に変換するためには、str()関数を使用することができます。これにより、PosixPathやWindowsPathオブジェクトを必要な場所で文字列に変換することができます。

python
import pathlib

path = pathlib.Path("example.txt")
path_str = str(path)

一方で、pathlibはパスを直接操作するために/演算子をオーバーロードしています。これにより、パスオブジェクト同士を結合する際に直感的なコードを書くことができます。

python
import pathlib

path1 = pathlib.Path("dir1")
path2 = pathlib.Path("dir2")
combined_path = path1 / path2

したがって、pathlibライブラリはパスをより柔軟に操作できるために、環境に依存しないパスオブジェクトの型を返す設計となっています。必要に応じて、パスオブジェクトを文字列に変換することで、従来の文字列操作を行うこともできます。


WindowsやUNIXといったOSの違いに左右されないようにしているというのが回答のようですね。納得！！！
本題に戻ります。
A34-180316_GML.zipは解凍してみるとA34a-180316.shp(ポリゴンデータ)、A34b-180316.shp(ポイントデータ)、A34c-180316.shp(ラインデータ)の３つのシェープファイルを持っているようなので、それぞれべつのGeoDataFrameデータとして読み込みます。
import urllib.request
import zipfile

zip_url = "https://nlftp.mlit.go.jp/ksj/gml/data/A34/A34-17/A34-180316_GML.zip"

# ZIPファイルをダウンロード
urllib.request.urlretrieve(zip_url, "world_herirages.zip")

import pathlib
# ZIPファイルを解凍するフォルダー'data'を現在のディレクトリの下に作成
NOTEBOOK_PATH = pathlib.Path().resolve()
DATA_DIRECTORY = NOTEBOOK_PATH / "data"

# ZIPファイルを解凍して、'data'に保存
with zipfile.ZipFile("world_herirages.zip", "r") as zip_ref:
    zip_ref.extractall(DATA_DIRECTORY)

!pip install geopandas
import geopandas
shapefile_path = str(DATA_DIRECTORY) + "/A34-180316_GML/A34a-180316.shp"
world_heritages_poly = geopandas.read_file(shapefile_path)

shapefile_path = str(DATA_DIRECTORY) + "/A34-180316_GML/A34b-180316.shp"
world_heritages_point = geopandas.read_file(shapefile_path)

shapefile_path = str(DATA_DIRECTORY) + "/A34-180316_GML/A34c-180316.shp"
world_heritages_line = geopandas.read_file(shapefile_path)

world_heritages_poly.head()


A34a_001A34a_002A34a_003A34a_004A34a_005A34a_006A34a_007A34a_008A34a_009geometry
011–@—²Ž›’nˆæ‚Ì•§‹³Œš‘¢•¨]“Þ—ÇŒ§15.3–@—²Ž›ú@AúAAúCAúE199312POLYGON ((135.73529 34.61004, 135.73536 34.609...
121–@—²Ž›’nˆæ‚Ì•§‹³Œš‘¢•¨]“Þ—ÇŒ§15.3–@‹NŽ›ú@AúAAúCAúE199312POLYGON ((135.74658 34.62322, 135.74659 34.622...
232•P˜Hé]•ºŒÉŒ§107.0•P˜Héú@AúC199312POLYGON ((134.68559 34.83558, 134.68638 34.836...
343ŒÃ“s‹ž“s‚Ì•¶‰»à]Ž ‰êŒ§A‹ž“s•{1056.0‰ê–Εʗ‹_ŽÐúAAúC199412MULTIPOLYGON (((135.75183 35.06112, 135.75205 ...
453ŒÃ“s‹ž“s‚Ì•¶‰»à]Ž ‰êŒ§A‹ž“s•{1056.0‰ê–ÎŒä‘c_ŽÐúAAúC199412POLYGON ((135.77132 35.03817, 135.77165 35.039...

表示が文字化け！またまたChatGPT先生に質問。
User
world_heritages_point.head()　の表示が文字化けするのですが
ChatGPT

文字化けが発生する場合、GeoPandasのデータフレームで日本語が正しく表示されない可能性があります。文字エンコーディングを指定してデータを読み込むことで、文字化けを解決することができます。

GeoPandasでデータを読み込む際に、encodingパラメータを使用して文字エンコーディングを指定します。一般的に、日本語の文字エンコーディングには'utf-8'または'shift-jis'が使用されます。以下のようにencodingパラメータを指定してデータを読み込んでみてください。

python
world_heritages_point = geopandas.read_file(shapefile_path, encoding='shift-jis')
world_heritages_point.head()

上記の例では、shapefile_pathはデータのパスを指す変数であり、'shift-jis'という文字エンコーディングが使用されています。もし文字化けが解消されない場合、'utf-8'を試してみてください。

python
world_heritages_point = geopandas.read_file(shapefile_path, encoding='utf-8')
world_heritages_point.head()

適切な文字エンコーディングを指定することで、文字化けを解消できるはずです。ただし、データ自体に文字化けが含まれている場合は、それを解決するためには元のデータ自体を修正する必要があります。

ほうほう、ではまず’Shift-jis’でお試し。

A34a_001A34a_002A34a_003A34a_004A34a_005A34a_006A34a_007A34a_008A34a_009geometry
011法隆寺地域の仏教建造物‐奈良県15.3法隆寺@、A、C、E199312POLYGON ((135.73529 34.61004, 135.73536 34.609...
121法隆寺地域の仏教建造物‐奈良県15.3法起寺@、A、C、E199312POLYGON ((135.74658 34.62322, 135.74659 34.622...
232姫路城‐兵庫県107.0姫路城@、C199312POLYGON ((134.68559 34.83558, 134.68638 34.836...
343古都京都の文化財‐滋賀県、京都府1056.0賀茂別雷神社A、C199412MULTIPOLYGON (((135.75183 35.06112, 135.75205 ...
453古都京都の文化財‐滋賀県、京都府1056.0賀茂御祖神社A、C199412POLYGON ((135.77132 35.03817, 135.77165 35.039...

どうやら成功したようです。

これをとりあえずシンプルにplot()してみようと思います。
!pip install japanize_matplotlib
import matplotlib.pyplot as plt
import japanize_matplotlib

fig, ax = plt.subplots()

# ポリゴンデータのプロット
world_heritages_poly.plot(ax=ax, color='blue', label='地域')

# ポイントデータのプロット
world_heritages_point.plot(ax=ax, color='red', label='代表地点', markersize=2)

# ラインデータのプロット
world_heritages_line.plot(ax=ax, color='green', label='線')

# 凡例の作成
proxy_poly = plt.Rectangle((0, 0), 1, 1, fc='blue')
proxy_point = plt.Rectangle((0, 0), 1, 1, fc='red')
proxy_line = plt.Rectangle((0, 0), 1, 1, fc='green')

# 凡例の表示
ax.legend([proxy_poly, proxy_point, proxy_line], ['地域', '代表地点', '線'])

# プロットの表示
plt.show()

なんとなく雰囲気はわかりますが、日本地図を下敷きにしないと、なんのことやらさっぱりわからないですね。これをFoliumの地図上に表示したいと思います。
まずはfoliumの準備。適当な大きさの日本地図を表示します。
import folium

# Foliumの地図オブジェクトを作成
# 地図の中心は東京を指定
map = folium.Map(location=[35.6895, 139.6917], zoom_start=4)

map



インタラクティブな地図はblogerには簡単にアップロードできないので画像のコピペです。
この地図に先ほどの世界文化遺産に関するGeoDataFrameをプロットします。
import folium

# Foliumの地図オブジェクトを作成
# 地図の中心は東京を指定
map = folium.Map(location=[35.6895, 139.6917], zoom_start=4)

# ポリゴンをプロット
folium.GeoJson(world_heritages_poly).add_to(map)

# ポリゴンをプロット
folium.GeoJson(world_heritages_point).add_to(map)

# ラインストリングをプロット
folium.GeoJson(world_heritages_line).add_to(map)

map


これはpointにマーカーを、PolygonとLineは多角形や線（マーカーに隠れて見えません）を描画しています。アップするとそれらも見えます。わかりやすいように多角形(緑)と線(赤)に変更します。
import folium

# Foliumの地図オブジェクトを作成
# 地図の中心は東京を指定
map = folium.Map(location=[35.6895, 139.6917], zoom_start=4)

# ポリゴンをプロット
folium.GeoJson(
    world_heritages_poly,
    name='polygon',
    style_function=lambda feature: {
        'fillColor': 'green',
        'color': 'green',
        'weight': 2,
        'fillOpacity': 0.5
    }
).add_to(map)

# ポリゴンをプロット
folium.GeoJson(world_heritages_point).add_to(map)

# ラインストリングをプロット
folium.GeoJson(
    world_heritages_line,
    name='line',
    style_function=lambda x: {'color': 'red', 'weight': 2}
).add_to(map)
map



東海付近をアップしてみます。
京都・奈良周辺は当然多いとして、案外、富士山周辺にたくさんの世界文化遺産が存在していることに驚きました。でも、これではそれぞれの場所が何なのかがさっぱりわかりません。ぜひ、ポップアップを設定して、マーカーなどをクリックしたら、名称がでるようにしたいですね。
import folium

# Foliumの地図オブジェクトを作成
# 地図の中心は東京を指定
map = folium.Map(location=[35.6895, 139.6917], zoom_start=4)

# ポリゴンをプロット
folium.GeoJson(
    world_heritages_poly,
    name='polygon',
    style_function=lambda feature: {
        'fillColor': 'green',
        'color': 'green',
        'weight': 2,
        'fillOpacity': 0.5
    }
).add_to(map)

# ポイントをプロット
# GeoDataFrameから必要な列を取り出し、ポップアップに設定
for _, row in world_heritages_point.iterrows():
    popup_text = row['A34b_003']  # ポップアップに表示する列の値
    popup = folium.Popup(popup_text, max_width=200)
    folium.Marker([row.geometry.y, row.geometry.x], popup=popup).add_to(map)



# ラインストリングをプロット
folium.GeoJson(
    world_heritages_line,
    name='line',
    style_function=lambda x: {'color': 'red', 'weight': 2}
).add_to(map)
map



「富士山－信仰の対象と芸術の源泉」だそうです。多くのマーカーをクリックしましたが、全部同じポップアップでした。そうか、富士山か…。
しかし、四国に世界文化遺産がないというのは、とても寂しいですね。

これに四国88箇所霊場の位置をプロットして仕上げにしたいと思います。
#Google ColabでGoogle Driveを使用するときのおまじない
from google.colab import drive
drive.mount('/content/drive')

#四国巡礼88箇所の位置情報を記録したShapeファイルを読み取り
gdf_shikoku88 = gpd.read_file('/content/drive/MyDrive/Colab Notebooks/202303_Shikoku88Data/shikoku_88temples.shp')

for index, row in gdf_shikoku88.iterrows():
    popup_text = '第' + str(index) + '番霊場<br>' + row['Name']  # ポップアップに表示する列の値
    popup = folium.Popup(popup_text, max_width=200)
    folium.Marker(
        location=[row.geometry.y, row.geometry.x],
        icon=folium.Icon(icon='cloud', color='red'),  # マーカーの種類を変更
        popup=popup,
    ).add_to(map)

map

水色マーカーは世界文化遺産、赤色マーカーは四国88箇所霊場です！！！
本編はここで完結します。





















コメント

コメントを投稿

このブログの人気の投稿

【完結】ランニング、お食事　２０２２年５月～２０２２年１２月

-
イメージ
２０２２年１２月某日 山口県長門市ラン、南の山並みに風力発電の風車が遠望できます。 ２０２２年１２月某日 山口県長門市黄波戸温泉を目指してラン！トワイライトを通過した夜景です。 ２０２２年１２月２９日 山陽道移動中。徳山西〜徳山東区間にて渋滞１０km。下松ＳＡにてやや遅めの昼食ということで名物牛骨ラーメンと半チャーハンのセット１１５０円を注文でございます。スープはすき焼きに使う牛脂の風味がします。見た目は透明感のあるスープですが、見た目よりもヘビーでした。 ２０２２年１２月２４日 クリスマスイブラン！イルミネーションの中を快走です。 ２０２２年１１月某日 香川大学農学部の収穫祭。おそらく３年ぶり。学食でチキンフライ定食、そして、お大根とお茶、蜜柑、シークワーサーを購入。 ２０２２年１２月某日 時間的には夕ランですが、冬至近辺は夜ランですね。レインボー通りのイルミネーションがキラキラです。 ２０２２年１１月某日 エリエールゴルフクラブ香川にてラウンド中の昼食にカキフライ＆アジフライ。ボリュームあり、美味でございました。冬のラウンドはカキフライに限ります。 ちなみに、ラウンドはグリーンが難しすぎてさんざんでした。ロングホール４ホールのうち、３ホールが３オン４パットなんて…ショットは切れていたんですけどねぇ。月例が開催されていたので、グリーン設定が激ムズ。 ２０２２年１２月初旬 師走の高知の空はひんやりと快晴です。アンパンマン列車が写り込んでいます。 ２０２２年１１月某日 道の駅牟礼源平の里にてオリーブハマチ丼と今季初の牡蠣フライをいただきました。オリーブハマチは真ん中に落とされた存在感たっぷりの黄身とマッチングが間違いない味わいです。牡蠣フライはショーケースに置いていましたが、レンチンで完全復活！どちらも激旨！小鉢はキンピラゴボウです。 １３時半前後にお邪魔し、行列かありましたが、３０分以内でお食事させていただきました。 県外の方にぜひお薦めしたい一品でございます。 ２０２２年１１月某日 福岡県筑紫野市にふるさと納税で、御礼の品は博多もつ鍋太閤さんのもつ鍋です。18000円で4〜5人分のハイボリュームでしかもおいし～い。3人家族の我が家には満腹の量でした。ラーメンも付いていて、翌朝はラーメンの残り、よく夕にはご飯を入れて、三度楽しむことができました。 ２０２２年１１某日
続きを読む

ランニング、グルメ、ドライブ　２０２３年４月〜

-
イメージ
２０２３年６月某日  I rounded at Elleair Golf Club Kagawa. It have renovated one of holes that has been named 'Rapan hole'. It held Rapan Cup Event that is named after 'Rapan Hole'. Rapan hole at Elliair golf club Kagawa is below, that future rabbit shape lown into the bunker. ２０２３年６月某日 神戸フランツさんの南京町店にてマンゴ壷プリンをテイクアウト。マンゴは万能です。めちゃ美味しいです。ここにしかない美味なのか、わたしにはわかりません。 ２０２３年６月某日 高知県仁淀川町　大渡ダム湖畔に佇む茶農家の店あすなろさんにてランチ。奮発してあすなろ御膳１５８０円をオーダー。ワンプレートに土佐ならではのいたどりやシイタケのたたき、茶味噌をトッピングしたおにぎりなどがきれいに盛り付けられていて特に女子が喜ぶ内容。がっつり男子にはポリュームが物足りないかもしれませんが、そこは食後の抹茶ソフトでカバーできます。 お土産は沢渡抹茶大福です。 松山市駅の西端角地に小さなお店を構えるじゃこ天ほか練り物のお店。松山市内で立地良しで、色々な種類の練り物がリーズナブルに購入できる貴重なお店です。宇和島じゃこ天と八幡浜じゃこ天を食べ当てることができないので、口コミするのは僭越ですが、軟骨天は美味いです。 またまた松山市内で「坂の上の雲」ラン！ 今回は花園町にある「子規誕生地跡」を訪問。看板には坂の上の雲の冒頭の一節が紹介されています。 この物語の主人公は、あるいはこの時代の小さな日本ということになるかもしれないが、ともかくわれわれは三人の人物のあとを追わねばならない。そのうちのひとりは、俳人になった。俳句、短歌と行った日本の古い短詩型に新風を入れてその中興の祖になった正岡子規である。 以下解説から抜粋。 近代文学史上大きな足跡を残した正岡子規は、慶応３年（１８６７）年９月１７日、この地に生まれ、東京根岸で３５年の生涯を閉じた。子規の本名は常規（つねのり）であるが、近親者は生涯を通じて「のぼさん」と呼んだ
続きを読む

ランニング、グルメ、クライム　２０２３年７月〜

-
イメージ
２０２３年８月某日 高松市内にあるうなぎ処鰻松さんでうなぎらんちでございます。うなぎは基本的にエクセレント！！！鰻松さんのタレは甘さ控えめのすっきりしたお味だと感じました。焼き加減、肝吸もグッド！また来たいお店です。 ２０２３年１０月某日 インドネシアのお土産コーヒー。小袋に分けられたコーヒー豆です。今日はSUMATORA MANDHELINGです。濃厚な味わいと深い風味を感じます。 ２０２３年９月某日 新日本百名山の一峰。 この日もたくさんの猫ちゃんに出会いました。おやまをぐるりとまわりながら上がるので、勾配もやや緩やかです。 今回、新日本百名山の由来を調べて、比較的高い年齢層を対象にしているということで、回ってみようと思っています。 来来亭　高松松縄店をリピートらんち。濃いめの醤油スープとたっぷり背脂は万人受けすると思います。元気のいいスタッフさん達も好印象。ラーメン・チャーハンセット、美味しくいただきました。 愛媛県西予市城川町の中城本家酒造合同会社さまで購入した「純米吟醸　城川郷」をいただきました。伊豫の地酒、軒先に杉玉が下がる古風な酒蔵。すっきりしていて、深い味わい。たぶん酒蔵でしか買えない酒粕を使ったしっとりケーキもとっても美味。わたしの大好きな酒蔵の1つです。 愛媛県久万高原町の美川地区にある道の駅みかわに併設されたレストランにてらんち！最近みかわ中華なるメニューができて麻婆豆腐や担々麺などセレクションが充実しています。以前、油淋鶏と担々麺をいただいたので今回は麻婆豆腐！みかわ中華の大ファンになっています。 小豆島とりっぷのらんちにBreadShopのすこっぷさんにお邪魔。１５時を過ぎていたので品薄でしたが、美味しそうなカスクートとスパイシーベーコンを購入。 ２０２３年９月某日 愛媛県久万高原町にある高原ゴルフ倶楽部を初らうんど。前半インは実力どおり５０超、後半アウトは５０未満でトータル１００超。うーん、ずーと悩んでいたアプローチ＆パットがこの日はまあまあ良かったのですが、ドライバーに乱れ、アイアンはシャンクが散発し始めて、結局普段どおりのスコアです。 アウトの後半に長く感じる（実際飛ばないわたしにはツーオンは難しい）ホールが多く、難しめ。そのためかアウトのらうんどは３時間かかりました。インはスムーズに回れました。 名物のナイスショット観音様を撮影。お
続きを読む