Data Visualizationとグラフ描画
Spiegelhalter (2019)は、data visualization (dataviz)の概念の導入において、いくつかのグラフ描画を説明している。以下において、ビンの中のゼリービーンの個数の推測データ、異性パートナーの人数調査、病院ごとの手術成功率、明白な関係があるが無相関の例、人口増加の折れ線グラフ、人口動態の散布図の各グラフ描画のPythonスクリプト例を試みた。
Pythonスクリプトの作成においては、Spiegelhalter(2019), p. xviに記載されているウェブサイトからダウンロードしたファイルを参考・利用した。
ビンの中のゼリービーンの個数を推測させることをSpiegelhalterらが行ったYou-Tubeでの調査では、915名による回答が得られ、推定個数値は、219個から31,337個であった。推定個数の分布を表すグラフを描画するスクリプトをリスト１のように作成した。実行すると図１のグラフが表示される。データファイル02-1-bean-data-full-x.csvは、Spiegelhalter(2019), p. xviに記載されているウェブサイトからダウンロードしたファイルに含まれていたものである。
図１
左側に素データのstrip-chart、箱図、ヒストグラムが表示され、右側に素データの対数変換値のグラフが表示されている。
リスト１ 個数の推定数の分布
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sb
import numpy as np
# ダウンロードしたファイルの読み込み
d_csv =
pd.read_csv('02-1-bean-data-full-x.csv', header=None)
data = d_csv.values
x_ck = np.array([200, 500, 2000, 5000,
20000])
fig, ax = plt.subplots(nrows=3,
ncols=2, sharex='col', figsize=(12,7))
fig.suptitle('Guess at number of beans
in jar', fontsize=20)
sb.stripplot(data, orient='h',
edgecolor='b', linewidth=1, alpha=0.5,
ax=ax[0][0])
sb.stripplot(np.log10(data),
orient='h', edgecolor='b',
linewidth=1, alpha=0.5, ax=ax[0][1])
ax[1][0].boxplot(data, vert=False)
ax[1][1].boxplot(np.log10(data),
vert=False)
ax[2][0].hist(data, bins=20)
ax[2][0].set_xlabel('Raw data',
fontsize=16)
ax[2][1].hist(np.log10(data), bins=20)
ax[2][1].set_xlabel('Logarithmic
scale', fontsize=16)
ax[2][1].set_xticks(np.log10(x_ck),
[f'{v:d}' for v in x_ck])
plt.tight_layout()
plt.show()
性生活についての調査（the National Sexual Attitudes and Lifestyle Survey: Natsal）の３回目のデータ（Natsal-3）における男女別の異性パートナーの数を棒グラフで表示するスクリプトをリスト２のように作成した。リスト２のスクリプトを実行すると、図２のグラフが表示される。データファイル02-4-sexual-partners-counts-x.csvは、Spiegelhalter(2019), p. xviに記載されているウェブサイトからダウンロードしたファイルに含まれていたものである。
図２
切りの良い数、10、20、30、40の比率が不連続に多いことを、Spiegelhalter(2019)は指摘している。
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sb
import numpy as np
# ダウンロードしたファイルの読み込み
d_csv =
pd.read_csv('02-4-sexual-partners-counts-x.csv')
data = d_csv.values
NumPartners = d_csv['NumPartners']
MenPercnt = d_csv['MenPercent']
WomenPercnt = d_csv['WomenPercent']
plt.xlim(0,50)
plt.bar(NumPartners, MenPercnt,
width=0.4, align='edge',
color='g')
plt.bar(NumPartners, WomenPercnt,
width=-0.4, align='edge',
color='r')
plt.ylabel('Percentage', fontsize=14)
plt.xlabel('Reported numberof lifetime
opposite-sex partners',
fontsize=14)
plt.text(25, 8, 'Men 35-44', c='g',
fontsize=16)
plt.text(25, 10, 'Women 35-44', c='r',
fontsize=16)
plt.tight_layout()
plt.show()
Spiegelhalter (2019), p.57には、病院ごとの手術数と30日生存確率の散布図が掲載されている。このグラフを描くスクリプトをリスト３のように作成してみた。リスト３のスクリプトを実行すると図３のグラフが表示される。データファイル02-5-child-heart-surgery-1991-x.csvと02-4-sexual-partners-counts-x.csvは、Spiegelhalter(2019), p. xviに記載されているウェブサイトからダウンロードしたファイルに含まれていたものである。
図３
上のグラフが1991―1995のデータであり、下のグラフは2012―2015のものである。
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sb
import numpy as np
# ダウンロードしたファイルの読み込み
#
Survival in under-1s, 1991-1995
#
d_csv1991 =
pd.read_csv('02-5-child-heart-surgery-1991-x.csv')
data1991 = d_csv1991.values
Hospitals = d_csv1991['Hospital']
Operations = d_csv1991['Operations']
ThirtyDay =
d_csv1991['ThirtyDaySurvival']
plt.figure(figsize=(7,7))
plt.subplot(211)
for Hosp, Ope, Thirty in
zip(Hospitals, Operations, ThirtyDay):
plt.scatter(Ope, Thirty,
label=Hosp)
plt.legend(bbox_to_anchor=(1.0, 1.0),
fontsize=8)
plt.title('Survival in under-1s,
1991-1995')
plt.xlabel('Number of operations')
plt.ylabel('% 30-day survival')
plt.tight_layout()
#
#
ダウンロードしたファイルの読み込み
#
Survival for all children, 2012-2015
#
d_csv2012 =
pd.read_csv('02-5-child-heart-surgery-2012-x.csv')
print(d_csv2012)
data2012 = d_csv2012.values
Hospitals = d_csv2012['Hospital']
Operations = d_csv2012['Operations']
ThirtyDay = d_csv2012['ThirtyDaySurvival']
plt.subplot(212)
for Hosp, Ope, Thirty in
zip(Hospitals, Operations, ThirtyDay):
plt.scatter(Ope, Thirty,
label=Hosp)
plt.legend(bbox_to_anchor=(1.0, 1.0),
fontsize=8)
plt.title('Survival for all children,
2012-2015')
plt.xlabel('Number of operations')
plt.ylabel('% 30-day survival')
plt.tight_layout()
plt.show()
Spiegelhalter (2019)、p.59に、関係が見られるパターンであるが無相関である散布図が掲載されている。２変量の間に明瞭な逆U字型の関係があるが、相関係数が０であるパターンを描くスクリプトをリスト４のように作成してみた。リスト４のスクリプトを実行すると、図４のグラフが表示される。
図４
変数Yは、変数Xの逆U字型の関数として関係があるが、相関係数は０である。
import numpy as np
import matplotlib.pyplot as plt
x = np.linspace(1,100,100)
y = (x - np.mean(x))**2 * (-0.01)
y = 100 + y
# y = 100 - 0.01
* (x - x_mean)**2
r = np.corrcoef(x,y)
plt.title(f'r = {r[0,1]:.3f}',
fontsize=24)
plt.plot(x, y, 'o', lw=0)
plt.xlabel('X', fontsize=16)
plt.ylabel('Y', fontsize=16)
plt.show()
人口増加を描く折れ線を素データに対するものと対数変換値に対するものを描くスクリプトを、リスト５およびリスト６に示すように作成してみた。
リスト５のスクリプトを実行すると図５のグラフが表示される。
図５
データ02-7-pop-past-future-2.7-necessities-mod.csvは、、Spiegelhalter(2019), p. xviに記載されているウェブサイトからダウンロードしたファイルに含まれていたファイル02-7-pop-past-future-2.7-necessities.csvにおいて、列変数名NorthernAmericaとLatinAmericaandCarribeanをN AmericanとS Americanに変更したものである。
リスト６のスクリプトを実行すると、図６のグラフが表示される。
図６
リスト６におけるデータファイル02-7-pop-past-future million in 51.csvは、Spiegelhalter(2019), p. xviに記載されているウェブサイトからダウンロードしたファイルに含まれていたものである。
図６のグラフは、対数変換された値に対するグラフである。色付きの折れ線は、地域別の人口を表し、灰色の折れ線は国別の人口を表している。
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sb
import numpy as np
#
# ダウンロードしたファイルの読み込み
# 地域別データ
#
d_csv_std = pd.read_csv('02-7-pop-past-future-2.7-necessities-mod.csv')
d_keys = d_csv_std.keys()
data = d_csv_std.values
names = d_keys[1:8]
data = data.T[:8]
for i in range(7):
plt.plot(data[0], data[i+1],
label=names[i])
plt.xlabel('Year')
plt.ylabel('Millions')
plt.legend(bbox_to_anchor=(1.0, -0.1),
ncols=4)
plt.tight_layout()
plt.show()
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sb
import numpy as np
#
#
ダウンロードしたデータの読み込み
#
国別データ
#
d_csv_all =
pd.read_csv('02-7-pop-past-future million in 51.csv')
data_all = d_csv_all.values
data_all = data_all.T
year_all = data_all[0]
data_all =
np.log10(data_all[1:]/1000)
# 対数変換
for i in range(len(data_all)):
plt.plot(year_all,
data_all[i], c='k', alpha=0.1)
#
# ダウンロードしたファイルの読み込み
# 地域別データ
#
d_csv_std =
pd.read_csv('02-7-pop-past-future-2.7-necessities-mod.csv')
d_keys = d_csv_std.keys()
data = d_csv_std.values
names = d_keys[1:8]
data = data.T[:8]
year = data[0]
data = np.log10(data[1:8])
# 対数変換
for i in range(7):
plt.plot(year, data[i],
lw=5, label=names[i])
plt.xlabel('Year')
plt.ylabel('Millions')
#
#
対数変換の目盛り用
#
pops = [1, 5, 10, 50, 100, 500, 1000,
5000, 10000]
log_pops = np.log10(pops)
plt.yticks(log_pops, pops)
plt.legend(bbox_to_anchor=(1.0, -0.1),
ncols=4)
plt.tight_layout()
plt.show()
人口動態の散布図を描くスクリプトをリスト７のように作成してみた。入力データファイル02-8-pop-past-future-x.xlsxは、Spiegelhalter(2019), p. xviに記載されているウェブサイトからダウンロードしたファイルに含まれていた02-8-pop-past-future-x.csvをExcelのファイルに変換したものである。ダウンロードしたCSVファイルをpandas.read_csvで読み込もうとすると文字コードエラーが出たが、Excelファイルに変換するとpandas.read_excelで読み込めた。CSVファイルからExcelファイルへの変換は、CSVファイルをExcelで読み込んだ後、Excelファイルとして保存して行った。Excelファイルとしての保存は、ファイルの保存名のファイル拡張子として.xlsxを選ぶとExcelファイルとして保存される。
リスト７のスクリプトを実行すると図７のグラフが表示される。
図７
縦軸は、1951年と2015年の人口比であり、横軸が地域名である。横軸の値は、円が重ならないようにランダムなジッターが加えられている。
円の大きさは、1951年の人口に応じて設定されている。人口の多い国は大きい円で表されていて目立つようになっている。
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sb
import numpy as np
#
# データファイルの読み込み
#
with
pd.ExcelFile('02-8-pop-past-future-x.xlsx') as f_xlsx:
data_xlsx =
pd.read_excel(f_xlsx)
v_data_xlsx = data_xlsx.values
# ndarray型に変換
Country = v_data_xlsx.T[1]
# 国（地域）名の取り出し
Country[62] = "Cote
d'Ivoire"
# 国名データの修正
Pops = data_xlsx[1951].values /
1000 # 1951年の人口
# 2015年の人口の1951年の人口に対する比
Increase = data_xlsx[2015].values /
data_xlsx[1951].values
#
#
Indexの値により国を地域に分類するときの参照に用いる配列
#
idx_lst = np.arange(0, 273+1)
idx_africa =
np.concatenate((idx_lst[15:34+1],idx_lst[36:44+1],
idx_lst[46:58+1],idx_lst[60:76+1]))
idx_asia =
np.concatenate((idx_lst[79:86+1], idx_lst[89:93+1],
idx_lst[95:103+1], idx_lst[105:115+1],
idx_lst[117:134+1]))
idx_europe =
np.concatenate((idx_lst[137:146+1], idx_lst[148:160+1],
idx_lst[162:177+1], idx_lst[179:187+1]))
idx_northam = idx_lst[241:245+1]
idx_latinam =
np.concatenate((idx_lst[190:215+1], idx_lst[217:224+1],
idx_lst[226:239+1]))
idx_oceania =
np.concatenate((idx_lst[248:249+1], idx_lst[251:255+1],
idx_lst[257:263+1], idx_lst[265:273+1]))
Index = data_xlsx['Index'].values
Continent = np.zeros(273) # 初期値0
for i in range(273):
if
np.count_nonzero(idx_africa - Index[i]) < len(idx_africa):
Continent[i] = 1 # Africa
elif
np.count_nonzero(idx_asia - Index[i]) < len(idx_asia):
Continent[i] = 2 # Asia
elif
np.count_nonzero(idx_europe - Index[i]) < len(idx_europe):
Continent[i] = 3 # Europe
elif
np.count_nonzero(idx_latinam - Index[i]) < len(idx_latinam):
Continent[i] = 4 # Latin America
elif
np.count_nonzero(idx_northam - Index[i]) < len(idx_northam):
Continent[i] = 5 # North America
elif
np.count_nonzero(idx_oceania - Index[i]) < len(idx_oceania):
Continent[i] = 6 # Oceania
idx_Classified = Continent > 0
#
Continentに値が割り振られたか
idx_Bigger = Pops > 1
# Popsの値が１より大きいか
#
#
上記２つの条件が満たされているか
#
idx_Keepers =
np.logical_and(idx_Classified, idx_Bigger)
#
#
上記2つの条件を満たすデータを選ぶ
#
Kpt_Continent = Continent[idx_Keepers]
Kpt_Increase = Increase[idx_Keepers]
Kpt_Pops = Pops[idx_Keepers]
Kpt_Country = Country[idx_Keepers]
#
# 散布図における円の重なりを防ぐためのランダムジッター
#
Kpt_Continent +=
np.random.uniform(-0.3, 0.3, len(Kpt_Continent))
v_min = np.log10(np.min(Kpt_Pops))
v_max = np.log10(np.max(Kpt_Pops))
def s_marker(v):
# vの値に対応する円の大きさ
v = np.log10(v)
return 10 + 10 * (v - v_min)
/ (v_max - v_min)
#
#
散布図において名前を表示する国のリスト
#
Targets = {"Cote d'Ivoire",
"Saudi Arabia", "Japan", "Georgia",
"Albania",
"Papua New Guinea"}
ContinentNames = ["Africa",
"Asia", "Europe", "Latin\nAmerica",
"North\nAmerica", "Oceania"]
plt.figure(figsize=(9, 7))
for i in range(len(Kpt_Continent)):
plt.plot(Kpt_Continent[i],
Kpt_Increase[i], 'o', c='k',
markersize=s_marker(Kpt_Pops[i]), fillstyle='none')
if Kpt_Country[i] in
Targets:
plt.text(Kpt_Continent[i], Kpt_Increase[i], Kpt_Country[i],
ha='center',
va='bottom', fontsize=12)
plt.plot([], 'o', c='k',
markersize=s_marker(5), fillstyle='none', label='\n5')
plt.plot([], 'o', c='k',
markersize=s_marker(10), fillstyle='none', label='\n10')
plt.plot([], 'o', c='k', markersize=s_marker(50),
fillstyle='none', label='\n50')
plt.plot([], 'o', c='k',
markersize=s_marker(100), fillstyle='none', label='\n100')
plt.legend(bbox_to_anchor=(1.15, 1),
title='Population\n(millions)')
plt.xticks(np.arange(1, 6.1),
ContinentNames, fontsize=16)
plt.yticks(fontsize=14)
plt.ylabel("Relative increase in
population, 1951-2015", fontsize=16)
plt.tight_layout()
plt.show()
岡本安晴（2019） いまさら聞けないPythonでデータ分析．丸善出版
Spiegelhalter, D. (2019). The art of statistics: Learning from data. Pelican