「ベイズ統計モデリングによるデータ分析入門」をPythonで写経 ～ Vol.3 データの可視化

書籍の著者　馬場真哉 先生

この記事は、書籍「RとStanではじめるベイズ統計モデリングによるデータ分析入門」第２部第３章「ggplot2によるデータの可視化」の Python 写経活動記録です。

ベイズ統計モデリングでよく使われる可視化手法を学びます！
書籍が採用する「R の可視化ライブラリ ggplot2」に相当する Python の可視化ライブラリを「seaborn」に定めて、動かしていきます！

では書籍を開いてベイズ統計モデリングの旅に出かけましょう🚀

はじめに

このブログシリーズは書籍「RとStanではじめるベイズ統計モデリングによるデータ分析入門」（講談社、「テキスト」と呼びます）の Python 写経です。

テキストの紹介と引用表記はリンク先の記事に掲載しています。

準備

■ 記事の範囲
この記事はテキスト第２部第３章の以下の節を取り扱います。

■ 利用データ
テキスト・サポートサイトのデータファイルを引用しています。
Jupyter Notebook ファイルと同一フォルダ内の「data」フォルダにデータファイルを格納しています。

■ ライブラリのインポート
Jupyter Notebook 形式でコードを記述します。
この記事で用いるライブラリをインポートします。

# インポート

# 数値計算
import pandas as pd

# 可視化
import matplotlib.pyplot as plt
import seaborn as sns
sns.set_theme()                          # ggplot風のスタイル
plt.rcParams['font.family'] = 'Meiryo'   # または import japanize_matplotlib 等

【コード補足】
ggplot 風の見た目にする目的で「sns.set_theme() 」を実行します。

第３章 seabornによるデータの可視化

データの読み込み

データの csv ファイルを pandas のデータフレーム形式で変数 fish に読み込みます。
read_csv 関数を利用します。

# p.103 データの読み込み
fish = pd.read_csv('./data/2-2-1-fish.csv')
print('fish.shape: ', fish.shape)
fish.head(3)

【実行結果】
length（体長）列に 400 個のデータが含まれています。

🔷🔷🔷

seaborn ライブラリの２つのインターフェース

seaborn は Python の可視化ライブラリです。
現在、２つのインターフェース（書き方）があります。

① Function インターフェース
作図するグラフの種類ごとに関数化されています。
例：ヒストグラム：histplot 関数

② Objects インターフェース
ggplot2 のオブジェクトの「足し合わせ」に似た書き方をします。
公式ドキュメントには「オブジェクトインターフェースは現在実験段階であり、未完成です」と書かれています（バージョン 0.13.2）。

この記事は ① Function インターフェース をご紹介し、最後にまとめて ② Objects インターフェース の実装例を添えます。

🔷🔷🔷

ヒストグラムとカーネル密度推定

カーネル密度推定の曲線は、MCMC サンプルを用いたパラメータの事後分布の描画などで頻出します！

ヒストグラムを描画します。
histplot 関数を利用します。
図 2.3.1 の左側に相当します。

# p.104 図2.3.1 ヒストグラム
sns.histplot(data=fish, x='length')
plt.title('ヒストグラム', loc='left');

【実行結果】

カーネル密度推定の結果を描画します。
kdeplot 関数を利用します。
図 2.3.1 の右側に相当します。

# p.104 図2.3.1 カーネル密度推定
sns.kdeplot(data=fish, x='length')
plt.title('カーネル密度推定', loc='left');

【実行結果】

ヒストグラムとカーネル密度推定の曲線を重ね合わせします。
histplot 関数の引数 kde=True を指定します。

# p.105 図2.3.2 グラフの重ね合わせ
sns.histplot(data=fish, x='length', kde=True, stat='density')
plt.title('グラフの重ね合わせ', loc='left');

【実行結果】

ヒストグラムとカーネル密度推定の曲線グラフを横に並べます。

# p.105 図2.3.1 グラフの一覧表示

# 描画領域の設定
fig, ax = plt.subplots(1, 2, figsize=(8, 4), tight_layout=True)
# ヒストグラムの描画
sns.histplot(data=fish, x='length', ax=ax[0])
ax[0].set_title('ヒストグラム', loc='left');
# カーネル密度推定の描画
sns.kdeplot(data=fish, x='length', ax=ax[1])
ax[1].set_title('カーネル密度推定', loc='left');

【実行結果】

【コード補足】
可視化ライブラリ matplotlib の subplots 関数を利用します。
axes と呼ばれる図（figure）内のサブプロットを表のように定義します。
今回は行 1、列 2 を次のように指定します。

fig, ax = plt.subplots(1, 2, ...)

axes は配列 ax に設定されました。
１つ目のヒストグラムは ax[0] を指定して描画します。
２つ目のカーネル密度推定の曲線は ax[1] を指定して描画します。

🔷🔷🔷

箱ひげ図とバイオリンプロット

箱ひげ図やバイオリンプロットは、MCMC サンプルを用いたパラメータの事後分布を「グループ別」で比較したいときなどで利用します！

新しいデータを seaborn ライブラリから pandas のデータフレーム形式で変数 iris に読み込みます。

# データの読み込み seabornのdataset利用

# seabornのdatasetよりデータを取得
iris = sns.load_dataset('iris')
# 結果の表示
print('iris.shape: ', iris.shape)
iris.head(3)

【実行結果】
変数（列）が５つ、標本サイズが 150 の「あやめデータセット」です。

箱ひげ図とバイオリンプロットを横に並べて描画します。
箱ひげ図は boxplot 関数、バイオリンプロットは violinplot 関数を利用します。
図 2.3.3 に相当します。

# p.107 図2.3.3 箱ひげ図とバイオリンプロット

# 描画領域の設定
fig, ax = plt.subplots(1, 2, figsize=(8, 4), tight_layout=True)
# 箱ひげ図の描画
sns.boxplot(data=iris, x='species', y='petal_length', ax=ax[0])
ax[0].set_title('箱ひげ図', loc='left')
# バイオリンプロットの描画
sns.violinplot(data=iris, x='species', y='petal_length', ax=ax[1])
ax[1].set_title('バイオリンプロット', loc='left');

【実行結果】

🔷🔷🔷

散布図

あやめの種類別に色分けをした散布図を描画します。
scatterplot 関数を利用します。
引数 hue に色分けしたい変数 'species'  を与えます。
図 2.3.4 に相当します。

# p.107 図2.3.4 色分けした散布図
sns.scatterplot(data=iris, x='petal_width', y='petal_length', hue='species');

【実行結果】

🔷🔷🔷

折れ線グラフ

折れ線グラフは時系列などの「推移」が大事なケースにおいて、観測値とベイズ統計モデルによる予測値を比較するときなどで活躍します！

新しいデータの csv ファイルを pandas の データフレーム nile に読み込みます。
このファイルは WEB サイトから直接取得できます。

# p.108 ナイル川流量データの読み込み WebサイトのCSVファイルを直接読み込む

# ナイル川流量データのURL
url = 'https://vincentarelbundock.github.io/Rdatasets/csv/datasets/Nile.csv'
# CSVファイルの読み込み
nile = pd.read_csv(
    url,                     # csvファイルの置き場
    usecols=[1, 2],          # 2～3列目を利用
    names=['year', 'nile'],  # 列名の変更
    index_col=0,             # 利用する列のうち最初の列をインデックスにする
    skiprows=1               # １行目を読み込まない
)
# 結果の表示
print('nile.shape: ', nile.shape)
nile.head(3)

【実行結果】
1871～1970年のナイル川の流量 value です。

折れ線グラフを描画します。
lineplot 関数を利用します。
図 2.3.5 に相当します。

# p.109 図2.3.5 折れ線グラフ
plt.figure(figsize=(10, 4))
sns.lineplot(data=nile, x=nile.index, y='nile');

【実行結果】

🔷🔷🔷

Objects インターフェースの例

ここまで描いたグラフを seaborn の Objects インターフェースで再度描いてみましょう。

📊 追加インポート
seaborn の Objects インターフェースのライブラリを「so」の名称で読み込みます。
あわせて、matplotlib ライブラリもインポートします。

# 追加インポート
import seaborn.objects as so   # Objects インターフェース
import matplotlib as mpl       # matplotlib ライブラリ

📊 ヒストグラム
ggplot2 の「+」による足し合わせのように「.」（ドット）でグラフ描画の要素＝オブジェクトをつなぎ合わせる構文です。

# p.104 図2.3.1 ヒストグラム
(
    so.Plot(data=fish, x='length')     # データ
    .add(so.Bars(), so.Hist())         # レイヤー（マーク：Bars, 計算：Hist）
    .label(title='ヒストグラム')        # ラベル（タイトル）
    .theme({'font.family': 'Meiryo'})  # テーマ（フォント）
)

【実行結果】

📊 カーネル密度推定の描画

# p.104 図2.3.1 カーネル密度推定
(
    so.Plot(data=fish, x='length')              # データ
    .add(so.Line(), so.KDE())                   # レイヤー（マーク：Line, 計算：KDE）
    .label(title='カーネル密度推定', y='density')  # ラベル（タイトル、y軸ラベル）
    .theme({'font.family': 'Meiryo'})           # テーマ（フォント）
)

【実行結果】

📊 ヒストグラムとカーネル密度推定の曲線を重ね合わせ

# p.105 図2.3.1 グラフの重ね合わせ
(
    so.Plot(data=fish, x='length')
    .add(so.Bars(alpha=0.5), so.Hist('density'))   # ヒストグラム
    .add(so.Line(), so.KDE())                      # KDEプロット
    .label(title='グラフの重ね合わせ', y='density')
    .theme({'font.family': 'Meiryo'})
)

【実行結果】

📊 ヒストグラムとカーネル密度推定の曲線の横並び描画
matplotlib の subfigures を用いて２つの描画領域 sf1, sf2 を設定します。

# p.105 図2.3.1 グラフの一覧表示

# figureの設定
f = mpl.figure.Figure(figsize=(8, 4))
# subfigureの設定
sf1, sf2 = f.subfigures(1, 2)
# ヒストグラムの描画
(
    so.Plot(data=fish, x='length')
    .add(so.Bars(), so.Hist())
    .label(title='ヒストグラム')
    .theme({'font.family': 'Meiryo'})
    .on(sf1) # subfigure
    .plot()
)
# カーネル密度推定の描画
(
    so.Plot(data=fish, x='length')
    .add(so.Line(), so.KDE())
    .label(title='カーネル密度推定', y='density')
    .theme({'font.family': 'Meiryo'})
    .on(sf2) # subfigure
    .plot()
)

【実行結果】

📊 箱ひげ図とバイオリンプロット
Objects インターフェースには、箱ひげ図とバイオリンプロットを描画するオブジェクトは無さそうです。

📊 散布図

# p.107 図2.3.4 色分けした散布図
(
    so.Plot(data=iris, x='petal_width', y='petal_length')
    .add(so.Dot(), color='species')
    .label(title='散布図')
    .theme({'font.family': 'Meiryo'})
)

【実行結果】

📊 折れ線グラフ

# p.109 図2.3.5 折れ線グラフ
(
    so.Plot(data=nile, x=nile.index, y='nile')
    .add(so.Line())
    .label(title='折れ線グラフ')
    .theme({'font.family': 'Meiryo'})
    .layout(size=(10, 4))
)

【実行結果】

今回の記事は以上です。
楽しかったですね！

シリーズの記事

次の記事（PyMC版）

次の記事（Stan版）

前の記事

目次

ブログの紹介

note で８つのシリーズ記事を書いています。
ぜひ覗いていってくださいね！

１．のんびり統計

統計検定２級の問題集を手がかりにして、確率・統計をざっくり掘り下げるブログです。
雑談感覚で大丈夫です。ぜひ覗いていってくださいね。
統計検定２級公式問題集CBT対応版に対応しています。
Python、EXCELのサンプルコードの配布もあります。

２．統計・データ分析とつながる

シリーズ「統計・データ分析とつながる」は、統計・データ分析との「つながり」を発掘して、コラム風に仕立てたブログシリーズです。
生成 AI の力を借りながら、統計・データ分析の入り口をイメージして、自由気ままに書きました。
たとえば…
・日常生活と統計のつながり
・統計検定２級からその先へのつながり
気楽にお読みいただけたら嬉しいです🍀

３．実験！たのしいベイズモデリング１＆２をPyMC Ver.5で

書籍「たのしいベイズモデリング」・「たのしいベイズモデリング２」の心理学研究に用いられたベイズモデルを PyMC Ver.5で描いて分析します。
この書籍をはじめ、多くのベイズモデルは Ｒ言語＋Stanで書かれています。
PyMCの可能性を探り出し、手軽にベイズモデリングを実践できるように努めます。
身近なテーマ、イメージしやすいテーマですので、ぜひぜひPyMCで動かして、一緒に楽しみましょう！

４．実験！岩波データサイエンス1のベイズモデリングをPyMC Ver.5で

書籍「実験！岩波データサイエンスvol.1」の４人のベイジアンによるベイズモデルを PyMC Ver.5で描いて分析します。
この書籍はベイズプログラミングのイロハをざっくりと学ぶことができる良書です。
楽しくPyMCモデルを動かして、ベイズと仲良しになれた気がします。
みなさんもぜひぜひPyMCで動かして、一緒に遊んで学びましょう！

５．楽しい写経 ベイズ・Python等

ベイズ、Python、その他の「書籍の写経活動」の成果をブログにします。
主にPythonへの翻訳に取り組んでいます。
写経に取り組むお仲間さんのサンプルコードになれば幸いです🍀

６．RとStanではじめる心理学のための時系列分析入門 を PythonとPyMC Ver.5 で

書籍「RとStanではじめる心理学のための時系列分析入門」の時系列分析をPythonとPyMC Ver.5 で実践します。
この書籍には時系列分析のテーマが盛りだくさん！
時系列分析の懐の深さを実感いたしました。
大好きなPythonで楽しく時系列分析を学びます。

７．データサイエンスっぽいことを綴る

統計、データ分析、AI、機械学習、Pythonのコラムを不定期に綴っています。
統計・データサイエンス書籍にまつわる記事が多いです。
「統計」「Python」「数学とPython」「R」のシリーズが生まれています。

８．Python機械学習プログラミング実践記

書籍「Python機械学習プログラミング PyTorch & scikit-learn編」を学んだときのさまざまな思いを記事にしました。
この書籍は、scikit-learnとPyTorchの教科書です。
よかったらぜひ、お試しくださいませ。

最後までお読みいただきまして、ありがとうございました。