Pythonで世界地図を描く

#Google Colabで世界地図を書いてみる。
Google Colab等のJupyter notebookを使うと、世界地図を簡単に描ける。本コードは、そのままで表示するだけなら、Google Colabにログインして、コピペで動くと思う。ただし、世界銀行 (The World Bank) のデータを入れ替える場合、ファイル名等の編集が必要になる。
なお、ベースにした資料は、参考資料の1番である。

##Pythonスクリプトを動かす前の環境設定
以下のスクリプトを動かし、必要なデータ(シェープファイルおよび、各国のデータ(今回は緑化率))を持ってくる。

```python
!wget http://www.naturalearthdata.com/http//www.naturalearthdata.com/download/10m/cultural/ne_10m_admin_0_countries.zip
!unzip ne_10m_admin_0_countries.zip
!wget http://api.worldbank.org/v2/en/indicator/AG.LND.FRST.ZS?downloadformat=csv -O AG_LND_FRST_ZS.zip
!unzip AG_LND_FRST_ZS.zip
!apt install proj-bin libproj-dev libgeos-dev
!pip install https://github.com/matplotlib/basemap/archive/v1.1.0.tar.gz
!pip install geonamescache
```
##プログラム本体
実際に実行するプログラムは以下の通りである。なお、ここで、CSVファイルに別データを入れれば、表示が異なるデータを世界地図にマッピングできる。

```python
#
# The original source code is follows. 
# http://ramiro.org/notebook/basemap-choropleth/
#
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from geonamescache import GeonamesCache
from matplotlib.patches import Polygon
from matplotlib.collections import PatchCollection
from mpl_toolkits.basemap import Basemap

filename = 'API_AG.LND.FRST.ZS_DS2_en_csv_v2.csv'
year = '2015'
shapefile = 'ne_10m_admin_0_countries'
cols = ['Country Name', 'Country Code', year]
title = 'Forest area as percentage of land area in {}'.format(year)
imgfile = 'imgtest.png'
descripton = '''
Forest area is land under natural or planted stands of trees of at least 5 meters in situ, whether productive or not, and excludes tree stands in agricultural production systems (for example, in fruit plantations
and agroforestry systems) and trees in urban parks and gardens. Countries without data are shown in grey. Data: World Bank - worldbank.org • Original Code: Ramiro Gómez - ramiro.org'''.strip()

gc = GeonamesCache()
iso3_codes = list(gc.get_dataset_by_key(gc.get_countries(), 'iso3').keys())
df = pd.read_csv(filename, skiprows=4, usecols=cols)
df.set_index('Country Code', inplace=True)
df = df.ix[iso3_codes].dropna()
values = df[year]
num_colors = 9
cm = plt.get_cmap('Greens')
scheme = [cm(i / num_colors) for i in range(num_colors)]
bins = np.linspace(values.min(), values.max(), num_colors)
df['bin'] = np.digitize(values, bins) - 1
df.sort_values('bin', ascending=False).head(10)

fig = plt.figure(figsize=(22, 12))
ax = fig.add_subplot(111, facecolor='w', frame_on=False)
fig.suptitle('Forest area as percentage of land area in {}'.format(year), fontsize=30, y=.95)

m = Basemap(lon_0=0, projection='robin')
m.drawmapboundary(color='w')

m.readshapefile(shapefile, 'units', color='#444444', linewidth=.2)
for info, shape in zip(m.units_info, m.units):
    iso3 = info['ADM0_A3']
    if iso3 not in df.index:
        color = '#dddddd'
    else:
        color = scheme[df.loc[iso3]['bin']]
    patches = [Polygon(np.array(shape), True)]
    pc = PatchCollection(patches)
    pc.set_facecolor(color)
    ax.add_collection(pc)

# Cover up Antarctica so legend can be placed over it.
ax.axhspan(0, 1000 * 1800, facecolor='w', edgecolor='w', zorder=2)

# Draw color legend.
ax_legend = fig.add_axes([0.35, 0.14, 0.3, 0.03], zorder=3)
cmap = mpl.colors.ListedColormap(scheme)
cb = mpl.colorbar.ColorbarBase(ax_legend, cmap=cmap, ticks=bins, boundaries=bins, orientation='horizontal')
cb.ax.set_xticklabels([str(round(i, 1)) for i in bins])

# Set the map footer.
plt.annotate(descripton, xy=(-.8, -3.2), size=14, xycoords='axes fraction')

plt.savefig(imgfile, bbox_inches='tight', pad_inches=.2)
```

# 再生可能エネルギー利用率で見た世界地図
今度も世界銀行のデータのうち、再生可能エネルギー利用率のグラフを書いてみる。
##環境設定
以下のスクリプトを動かし、必要なデータ(シェープファイルおよび、各国のデータ(今回は再生可能エネルギー))を持ってくる。

```python
!wget http://www.naturalearthdata.com/http//www.naturalearthdata.com/download/10m/cultural/ne_10m_admin_0_countries.zip
!unzip ne_10m_admin_0_countries.zip
!wget http://api.worldbank.org/v2/en/indicator/EG.FEC.RNEW.ZS?downloadformat=csv -O EG_FEC_RNEW_ZS.zip
!unzip EG_FEC_RNEW_ZS.zip
!apt install proj-bin libproj-dev libgeos-dev
!pip install https://github.com/matplotlib/basemap/archive/v1.1.0.tar.gz
!pip install geonamescache
```
##プログラム本体
再生可能エネルギー利用率なので、こちらも緑のままとした。なお、年を変えると、該当データも変わる(csvデータ参照)。

```python
#
# The original source code is follows. 
# http://ramiro.org/notebook/basemap-choropleth/
#
import matplotlib as mpl
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
from geonamescache import GeonamesCache
from matplotlib.patches import Polygon
from matplotlib.collections import PatchCollection
from mpl_toolkits.basemap import Basemap

filename = 'API_EG.FEC.RNEW.ZS_DS2_en_csv_v2.csv'
year = '2015'
shapefile = 'ne_10m_admin_0_countries'
cols = ['Country Name', 'Country Code', year]
title = 'Renewable energy percentage of land area in {}'.format(year)
imgfile = 'imgtest.png'
descripton = '''
Renewable energy percentage world map by country.
Countries without data are shown in grey. Data: World Bank - worldbank.org • Original Code: Ramiro Gómez - ramiro.org'''.strip()

gc = GeonamesCache()
iso3_codes = list(gc.get_dataset_by_key(gc.get_countries(), 'iso3').keys())
df = pd.read_csv(filename, skiprows=4, usecols=cols)
df.set_index('Country Code', inplace=True)
df = df.ix[iso3_codes].dropna()
values = df[year]
num_colors = 9
cm = plt.get_cmap('Greens')
scheme = [cm(i / num_colors) for i in range(num_colors)]
bins = np.linspace(values.min(), values.max(), num_colors)
df['bin'] = np.digitize(values, bins) - 1
df.sort_values('bin', ascending=False).head(10)

fig = plt.figure(figsize=(22, 12))
ax = fig.add_subplot(111, facecolor='w', frame_on=False)
fig.suptitle(title, fontsize=30, y=.95)

m = Basemap(lon_0=0, projection='robin')
m.drawmapboundary(color='w')

m.readshapefile(shapefile, 'units', color='#444444', linewidth=.2)
for info, shape in zip(m.units_info, m.units):
    iso3 = info['ADM0_A3']
    if iso3 not in df.index:
        color = '#dddddd'
    else:
        color = scheme[df.loc[iso3]['bin']]
    patches = [Polygon(np.array(shape), True)]
    pc = PatchCollection(patches)
    pc.set_facecolor(color)
    ax.add_collection(pc)

# Cover up Antarctica so legend can be placed over it.
ax.axhspan(0, 1000 * 1800, facecolor='w', edgecolor='w', zorder=2)

# Draw color legend.
ax_legend = fig.add_axes([0.35, 0.14, 0.3, 0.03], zorder=3)
cmap = mpl.colors.ListedColormap(scheme)
cb = mpl.colorbar.ColorbarBase(ax_legend, cmap=cmap, ticks=bins, boundaries=bins, orientation='horizontal')
cb.ax.set_xticklabels([str(round(i, 1)) for i in bins])

# Set the map footer.
plt.annotate(descripton, xy=(-.8, -3.2), size=14, xycoords='axes fraction')

plt.savefig(imgfile, bbox_inches='tight', pad_inches=.2)
```
##出力例
緑の色が濃いところが、再生エネルギーの利用率が高いところである。
![renewable_energy_2015.png](https://qiita-image-store.s3.amazonaws.com/0/153471/dca48b21-7826-e91e-89c9-c3584bfbb569.png)

#参考資料
1.[Creating a Choropleth Map of the World in Python using Basemap](http://ramiro.org/notebook/basemap-choropleth/)
2.[The World Bank Indicators](https://data.worldbank.org/indicator/)