この記事は朝日新聞社 Advent Calendar 2023の19日目の記事です。🎉

昨日は　@yamano1212 の「音声合成を使った音声認識のためのデータ作成とその結果」でした。

メディア研究開発センターの新妻です。

皆さんは特に理由もなく、とりあえずデータを地図にマッピングしたくなりませんか？
今日は景気ウォッチャーというデータを地図にマッピングしてみたくなったので、やっちゃいました。

そもそも景気ウォッチャーとは

景気ウォッチャーとは、内閣府が毎月実施している「街角の景況感を判断するためのアンケート調査」を公開したデータで、下記のページで公開されています。

この調査結果は街角景気とも言われ、様々な業種の経営者や労働者にインタビューをして、３ヶ月前と比較した景気の現状や先行きについて5段階評価をしてもらったというものです。

そのため、このデータは客観的な指標ではなく日本で働く皆さんの「実感」に基づいた景気が反映されていると考えられます。

今回使うデータ

今回は公開されたばかりの令和5年11月調査の統計表データを使います。
テキストデータである判断理由集などもございますが、今回は地図にマップにプロットをするため季節調整済みの「全国の分野・業種別、地域別DIの推移」というデータを使います。

この"DI"とは、Diffusion Indexの略です。
これは、質問に対する各選択肢の回答数の百分比に対して、点数をかけて指標化したものです。
このデータでは、5段階評価のうち最も良い回答を1、最も悪い回答を0とし、その間を0.25刻みで点数を割り当て、その回答の構成比率をかけた指標となっています。

単純に解釈するだけであれば、50であれば景気は横ばい、それ以上なら上向きになっており、それ以下なら下向きだと回答されていると考えてください。

DIに関する詳細は下記ページに投げます。

また、景気ウォッチャーの指標については、このページがわかりやすいと思います。

データ準備

まず、下記ページより次のリンクをクリックしてExcelデータをダウンロードします。

今回はダウンロードしたExcelデータのうち、季節調整済みの「景気の現状判断（方向性）ＤＩ」を利用しました。

使用ライブラリ

今回の実装には下記のライブラリを使っています。

geopandas
pandas
xlrd
pyogrio
matplotlib
japanize-matplotlib

景気ウォッチャーをpandasで読み込む

Excel形式ですが、A列からはじまっていなかったり、カラム名が1行に揃っていなかったりと、ちょっとした整形が施されているため、正しくデータを読み込むためには少しコツがいります。
下記のように、読み込んでから一部の行や列を整形するとDataFrameに合わせた形になりました。

df = pd.read_excel('{watcher5.xlsのパス}', sheet_name='地域別（現状)', skiprows=2, usecols='B:R')

# 列名が別々の行になっているので、1行にまとめる
new_columns = df[0:3].bfill().iloc[0].values
new_columns[0] = '年'
new_columns[1] = '月'
df.columns = new_columns
df.drop([0, 1, 2], inplace=True)

# 年が1月の行にしか埋まっていないので、2月以降を1月の年で埋める
df['年'].ffill(inplace=True)

# 年月をint型に変換
df['年'] = df['年'].str.replace('年', '').astype(int)
df['月'] = df['月'].astype(int)

# 必要な列のみ抽出
df.drop(columns=['全国', '関東', '東京都'], inplace=True)

# 年月をindexにする
df['年月'] = ((df['年'].astype(str) + df['月'].astype(str).str.zfill(2))).astype(int)
df.set_index('年月', inplace=True)

下記のようなDataFrameがになっていれば、正しく整形できていそうです。

日本地図の境界データを読み込んで地域単位にまとめる

次にgeopandasを使って日本の境界データを読み込みましょう。
境界データは、国土交通省が公開している国土数値情報ダウンロードサイトの「行政区域データ」を使うことにしましょう。

そのため、上記ページより次のダウンロードボタンを探してダウンロードしましょう。
（ちなみに2023年のデータで427MBあります）

ダウンロードができたら圧縮ファイルを解凍して、その中身を整形します。
初期のデータは市区町村単位の境界データになっているため、景気ウォッチャーの区分に合わせて地方の単位にまで境界データを結合していきます。

境界データの結合には、geopandasのdissolve関数を使います。

実際の結合コードはこちらです。（実行完了までそこそこ時間がかかります。）

import geopandas as gpd

# 景気ウォッチャーの地方区分
# https://www5.cao.go.jp/keizai3/watcher/watcher_mokuteki.html
REGION2PREFS = {
    '北海道': ['北海道'],
    '東北': ['青森県', '岩手県', '宮城県', '秋田県', '山形県', '福島県'],
    '南関東': ['東京都', '神奈川県', '千葉県', '埼玉県'],
    '北関東': ['茨城県', '栃木県', '群馬県'],
    '甲信越': ['山梨県', '長野県', '新潟県'],
    '東海': ['愛知県', '静岡県', '岐阜県', '三重県'],
    '北陸': ['富山県', '石川県', '福井県'],
    '近畿': ['大阪府', '京都府', '兵庫県', '奈良県', '滋賀県', '和歌山県'],
    '中国': ['岡山県', '広島県', '鳥取県', '島根県', '山口県'],
    '四国': ['香川県', '徳島県', '愛媛県', '高知県'],
    '九州': ['福岡県', '佐賀県', '長崎県', '大分県', '熊本県', '宮崎県', '鹿児島県'],
    '沖縄': ['沖縄県'],
}

PREF2REGION = {k:v for v, ks in REGION2PREFS.items() for k in ks}

gdf = gpd.read_file('{解凍した中身の入ってるフォルダのパス}', encoding='shift-jis')

# 都道府県でまとめる
gdf_pref = gdf.dissolve(by='N03_001', as_index=False)

# 都道府県を元に地方区分を割り当て
gdf_pref['region'] = gdf_pref['N03_001'].map(PREF2REGION)

# 地方でまとめる
gdf_region = gdf_pref.dissolve(by='region', as_index=False)

# 不要なカラムを落とす
gdf_region = gdf_region.drop(columns=['N03_001', 'N03_002', 'N03_003', 'N03_004', 'N03_007'])

gdf_region.to_file('{境界データの出力先のパス}', encoding='utf-8', index=False, driver='GeoJSON')

境界データと景気ウォッチャーを結合してプロットする

今回は２０２３年１１月に集計分のデータのみを使うので、まずは一覧にしてみましょう。

df.loc[202311, gdf['region'].tolist()].reset_index(drop=False).sort_values(by=202311, ascending=False)

境界データをgeopandasで読みこんで、このデータと結合させていきましょう。

# 今回のgeopandasはpyogrioをエンジンにしないとうまく読み込めない
gpd.options.io_engine = "pyogrio"

gdf = gpd.read_file(open('{前工程で作った境界データのパス}', encoding='utf-8'))

# 東京を中心とした日本測地系に変換
# これをしないとプロットした地図が歪みます
gdf = gdf.to_crs(epsg=4301)

# 景気ウォッチャーのDataFrameから2023年11月分の指標だけ取り出す
df_202310 = df.loc[202311, gdf['region'].tolist()].reset_index(drop=False)

# 境界データと景気ウォッチャーの2023年11月分データを結合
merged_gdf = gdf.merge(df_202311, left_on='region', right_on='index', how='left')

これで、プロットをするためのデータが完成しました。
あとは、日本地図にデータをプロットするだけです。

# マップにプロット
fig, ax = plt.subplots(figsize=(12, 7))
merged_gdf.plot(column=202311, legend=True, cmap='OrRd', ax=ax, legend_kwds={'loc': 'lower right'})

plt.axis('off')
ax.set_aspect('equal', 'datalim')

上記のコードで、2023年11月の地域別景気の現状判断DIをプロットしたものがこちらです。

日本地図は国土数値情報（国土交通省 https://nlftp.mlit.go.jp/ksj/gml/datalist/KsjTmplt-N03-v3_1.html ）より2023年12月11日に取得

（あとから思ったんですけど、50を白色としてそれより大きい数字を暖色、それより小さい数字を寒色にするほうがわかりやすかったですね。）

結果を読み取ってみる

基本的には３ヶ月前からの比較で回答をしてもらう指標であるため、これ単体で判断をすることは難しいと思いますが、この地図と表を読む限りでは、景気が先月よりも上向いていると感じている人が多くなっているのは沖縄だけで、九州は微増、一方で北陸は微減、その他の地方ではあまり変わらずという感じだと言えそうです。

さらにここから、（pysalなどをいれて）空間統計量を計算したり、別ページで公開されている景気判断理由をテキスト解析したりすることで、この結果をより細かく掘り下げることができそうです。
また、現時点では2023年11月の現状判断のみをプロットしていますが、時系列も考慮することでもっと立体的な分析ができるかもしれません。
今回の記事はここまでですが、また別の機会に上記の分析などを続きとしてやりたいなと思います。

明日は @kwbta のEMNLP参加記だそうです！