衛星から夜の地球を観測したデータは「夜間光データ」と呼ばれています。

以前、夜間光データと土地価格に相関があるかを調査した記事(https://sorabatake.jp/27232/)を公開しました。結果は、夜間光と土地価格に相関はありそうだが、それぞれの変化率には相関がなさそうというものでした。

せっかくならもっとたくさんのデータとの相関を調査してみようということで、e-Stat(政府統計の総合窓口)で公開されているありとあらゆるデータとの相関を調査してみることにしました。

400種類以上の統計データと夜間光との相関係数を算出し、正の相関と負の相関のある統計データをそれぞれランキングで発表したいと思います！　

扱うデータについて

前回の記事同様、夜間光データはGoogle Earth Engine(GEE)で公開されているVIIRS(※1)を用います。

※1 アメリカの衛星であるSuomi-NPPに搭載されているVisible/Infrared Imager and Radiometer Suite(VIIRS)

夜間光と相関を調査する統計データはe-Statから取得します。
e-Statは日本の統計が閲覧できる政府統計のポータルサイトです。(https://www.e-stat.go.jp/)各省庁が公表する統計データがまとまっており、検索したり地図上に表示できたりします。

今回は以下のページから取得できる、都道府県別の統計データを対象に相関を調査します。
https://www.e-stat.go.jp/regional-statistics/ssdsview/prefectures

A:人口・世帯、B:自然環境、C:経済基盤、D:行政基盤、E:教育、F:労働、G:文化・スポーツ、H:居住、I:健康・医療、J:福祉・社会保障、K:安全、と幅広いデータを取得することができます。

ここから2019年のデータを可能な限り取得し、2019年の夜間光データの平均値と相関を調査していきます。

事前の準備として、上記のページから統計データのcsvをダウンロードしておきます。
25個ずつしかデータが選択できないため、下の写真のように複数のcsvをダウンロードし、フォルダにまとめておきます。
※e-StatのデータはAPI経由でもダウンロード可能です。
https://www.e-stat.go.jp/developer

夜間光テーブルの作成

それではまずは都道府県ごとに夜間光の平均値を求めてみます。

まずは必要なライブラリをインポートします。

import numpy as np
import pandas as pd
import warnings
import json
import glob
from scipy import stats
warnings.filterwarnings('ignore')

# reminder that if you are installing libraries in a Google Colab instance you will be prompted to restart your kernal

try:
   import geemap, ee
   import seaborn as sns
   import matplotlib.pyplot as plt
   import geopandas as gpd
   from matplotlib import cm

except ModuleNotFoundError:
   if 'google.colab' in str(get_ipython()):
       print("package not found, installing w/ pip in Google Colab...")
       !pip install geemap seaborn matplotlib
       !pip install geopandas
   else:
       print("package not found, installing w/ conda...")
       !conda install mamba -c conda-forge -y
       !mamba install geemap -c conda-forge -y
       !conda install seaborn matplotlib -y
   import geemap, ee
   import seaborn as sns
   import matplotlib.pyplot as plt
   import geopandas as gpd
   from matplotlib import cm

次に、GEEの認証とモジュールのインポートを行います。
GEE認証を実行すると認証に必要なURLが表示されるので、URLへアクセスしてGoogleアカウントを指定し、表示された認証コードをGoogle Colabのボックスへコピペしてください。

try:
        ee.Initialize()
except Exception as e:
        ee.Authenticate()
        ee.Initialize()

都道府県ごとのポリゴンデータを取得します。
GADM(https://gadm.org/download_country_v3.html)からファイルをダウンロードして、 GeoPandasで読み込みます。ここで後にgeometryを取得しやすくするために、読み込んだGeoPandasのテーブルをjson形式に変換したものも作成しておきます。

prefectures_shp_path = 'shpファイルのパスに変更してください'
jpnShp = gpd.read_file(prefectures_shp_path, crs='EPSG:4326')
df = json.loads(jpnShp.to_json())
jpnShp.plot()

日本のポリゴンデータを取得することができました。

取得したポリゴンデータから、都道府県名のリストを取得します。

pres = jpnShp['NL_NAME_1'].unique()
print(pres)

無事都道府県のリストを取得できています。

今回扱う夜間光データであるVIIRSの月間コンポジット(月間で集約してノイズが除去されたもの)のImageCollectionを取得します。今回は2019年度の統計データとの相関係数を調査するため、夜間光データの期間も2019年4月1日から2020年3月31日までを指定します。

viirs = ee.ImageCollection("NOAA/VIIRS/DNB/MONTHLY_V1/VCMSLCFG").filterDate("2019-04-01", "2020-03-31").select('avg_rad').median()

都道府県ごとのgeometryを取得する関数を定義します。
引数で都道府県名を与えることで、先ほど作成したjson形式のデータからGEEのgeometryを作成し、返します。

def get_prefectures_geo(prefecture):
 features = [f for f in df['features'] if f['properties']['NL_NAME_1']==prefecture]
 return ee.Geometry.MultiPolygon(features[0]['geometry']['coordinates'])

夜間光データのテーブルを作成します。
ループで都道府県ごとの夜間光データの平均値を取得し、夜間光データテーブル(df_nl)に格納します。

cols_nl = ['pres', '2019-NL']
df_nl = pd.DataFrame(index=[], columns=cols_nl)

for i in range(len(pres)):
 aoi = get_prefectures_geo(pres[i])
 mean = ee.Number(viirs.reduceRegion(reducer=ee.Reducer.mean(), geometry=aoi, scale=450).get('avg_rad')).getInfo()
 df_nl = df_nl.append(pd.Series([pres[i], mean],index=df_nl.columns), ignore_index=True)

print(df_nl)

夜間光テーブル(df_nl)が作成できました。

統計データのテーブルの作成

先ほどe-Statからダウンロードしておいたデータのテーブルを作成します。csvを順々にダウンロードしていき、統計データのカラムのみを抽出してpandasのテーブル(df)のカラムに加えていきます。

folder_path = 'csvをまとめておいたフォルダパスに変更してください'

path = sorted(glob.glob(folder_path + "/*.csv"))

for i in range(len(path)):
   if i == 0:
       try:
           df =  pd.read_csv(path[i], encoding="CP932").drop(columns=["調査年","/項目"])
       except UnicodeDecodeError:
           df =  pd.read_csv(path[i], encoding="UTF8").drop(columns=["調査年","/項目"])
      
   else:
       try:
           df_app = pd.read_csv(path[i], encoding="CP932").drop(columns=["調査年","/項目"])
       except UnicodeDecodeError:
           df_app = pd.read_csv(path[i], encoding="UTF8").drop(columns=["調査年","/項目"])

       df = pd.merge(df, df_app, on='地域', how='inner')

このように都道府県ごとの統計データをもったテーブルが作成されます。

相関調査・ランキング発表

先ほど作成した夜間光テーブルと統計データテーブルを結合します。

# 表の結合
df_merge = pd.merge(df_nl, df, left_on='pres', right_on='地域', how='inner').drop(columns=['地域'])

夜間光と統計データの相関係数を算出し、相関係数が大きい順にソートして、csvで出力します。

# 相関マトリクス
df_merge.corr()

# 相関マトリクスを配列格納
df_rslt = df_merge.corr()

# 夜間光と相関係数が大きい順にソート
df_rslt_sort = df_rslt.sort_values('2019-NL', ascending=False)

# 夜間光以外の列を削除
df_rslt_sort_value = df_rslt_sort['2019-NL']

# CSV出力
df_rslt_sort_value.to_csv('出力するファイルパスに変更してください')

出力したCSVから相関係数が正のデータと相関係数が負のデータのランキングをそれぞれ発表します。

相関係数は一般的に0.7以上であれば強い正の相関、-0.7以下であれば強い負の相関と言われているため、その範囲に該当するデータ内でランキングをお見せしたいと思います。

※今回は因果関係を調査したわけではなく、あくまで夜間光と相関がありそうなデータを探索してみただけのため、教育や医療分野など夜間光と関係性がなさそうなデータもあります。

まずは表1に正の相関のあるデータのランキングを示します。相関係数が0.7以上のデータは57種類ありました。1,2位は人口密度、3,4,5,6位は医療関係、7,8,9位は教育関係になっています。

人口が多ければ人口密度も高くなり、医療機関も多く必要になり、学校もたくさん必要になります。10位以降を見てもわかる通り、消防署、郵便局、公園など人口が多い場所に必要な施設の数との強い相関が見られます。

夜間光の平均値が高い場所は、栄えている都道府県であり人口が多いことから、これらのデータとの強い相関が見られていると思われます。

人口と直接関係がなさそうなデータとしては、26位の民営賃貸住宅の家賃、32位の国内銀行預金残高、39位のきまって支給する現金給与月額等の個人のお金関係や、46位の財政力指数と言った行政レベルの経済指標などが目立ちます。夜間光の平均値が高い場所は、個人の所得や貯金額、行政レベルの経済指標が高い傾向にあることがわかります。

表１　夜間光データと正の相関のあるデータTOP57

次に表2に負の相関のあるデータのランキングを示します。相関係数が-0.7以下のデータは4種類ありました。

1位の従属人口指数は、生産年齢人口（15～64歳人口）に対する年少人口（15歳未満人口）と老年人口（65歳以上人口）の割合です。つまりこの従属人口指数が夜間光と負の相関があるということは、夜暗い都道府県ほど生産年齢人口に対して年少・老年人口が多い（高齢化が進んでいる）ということを示唆しています。

2位の人口10万人あたりの民生委員(児童委員)数とはそれぞれ以下のような方々です。

・民生委員
　厚生労働大臣から委嘱され、それぞれの地域において、常に住民の立場に立って相談に応じ、必要な援助を行い、社会福祉の増進に努める方々であり、「児童委員」も兼ねている。

・児童委員
　地域の子どもたちが元気に安心して暮らせるように、子どもたちを見守り、子育ての不安や妊娠中の心配ごとなどの相談・支援等を行う。

夜間光との関係は不明ですが、夜明るいほど委員数が少ないという傾向がありました。

3位は地方交付税割合です。地方交付税とは、「地方公共団体間の財源の不均衡を調整し、どの地域に住む国民にも一定の行政サービスを提供できるよう財源を保証するため、国から地方公共団体に交付される資金」のことです。

つまり地方交付税割合が高いということは、その都道府県の財源は国からの交付金の割合が高いということになり、都道府県自体の財政が低い水準であるということです。

よって夜間光と負の相関があるということは、夜暗い地域ほど地方交付税割合が高い、すなわち財政が低い水準であるということです。夜間光は経済指標と相関があると言われており、その通りの結果になっていることがわかります。

4位は着工新設持ち家比率です。夜暗い地域ほど新しく戸建てを建てた人の割合が高いということです。地方は家を建てる土地があり、都心はマンション等の集合住宅が多いことから、これもイメージ通りの結果になっているのではないでしょうか。

表２　夜間光データと負の相関のあるデータ

まとめ

夜間光とe-Statで公開されているありとあらゆる統計データについて相関を調査してみました。

今回は都道府県レベルを対象に相関調査を行いましたが、市区町村レベルで調査を行えばより詳細な調査を行うことができます。また、対象の統計データをe-Statに絞りましたが、その他政府以外が公開しているような統計データと比較してみると、また新たな発見がありそうです。

最近はe-Statのように気軽に統計データをダウンロードすることができるようになりました。今回は夜間光データを扱いましたが、その他の衛星データを使って統計データとの関係を調査してみても良さそうです。みなさんも衛星データと相関のある面白いデータを探ってみてはいかがでしょうか。

コラム：擬似相関に注意

ここまで夜間光と相関のある統計データをみてきましたが、医療関係や教育関係など夜の明るさと因果関係があるのか疑わしい項目が多く見られます。それらの項目を見ると、人口との相関がありそうです。そこで、それぞれの統計データの人口との相関係数も算出し、夜間光との相関係数との関係を確認してみます。

# 総人口と相関が強い順にソート
df_rslt_sort_population = df_rslt.sort_values('#A011000_総人口【万人】', ascending=False)

# 総人口以外の列を削除
df_rslt_sort_population = df_rslt_sort_population['#A011000_総人口【万人】']

# CSV出力
df_rslt_sort_population.to_csv('出力するファイルパス')

出力したcsvから夜間光との相関係数、人口との相関係数の散布図を作成しました。
夜間光との相関係数が高いものは、人口との相関係数も高くなっていることがわかります。つまり、人口との見かけ上の相関(※)となっていると思われます。

2つのデータに相関関係があっても、因果関係があるとは限らないということには注意しておきましょう。