「データ解析のための統計モデリング入門」をPythonで写経 Vol.23 ～ 11章「空間構造のある階層ベイズモデル」②ベイズ統計モデリング～一次元空間の階層ベイズモデルと欠測データの推定

11章「空間構造のある階層ベイズモデル」

書籍の著者　久保拓弥 先生

書籍「データ解析のための統計モデリング入門」11章「空間構造のある階層ベイズモデル」の Python写経活動記録 です。 

この記事は前回記事に引き続き、一次元空間の空間相関を考慮した階層ベイズモデル に取り組みます。

目的変数データに欠測が含まれている場合、ベイズ統計モデルが欠測を補完してくれます。
ただし、不適切なモデルの場合には欠測の補完も不適切になります。
一次元空間データの欠測補完を２つのベイズ統計モデルで実行し、空間相関を考慮したモデルが適切であることを体感します。

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

はじめに

このブログシリーズは、書籍「データ解析のための統計モデリング入門 一般化線形モデル・階層ベイズモデル・MCMC」（岩波書店、「テキスト」と呼びます）の Python 写経を通じて得た「統計モデリングの楽しさ」をご紹介します。

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

準備

準備

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

■ 利用データ
テキスト・サポートサイトのデータファイルを引用しています。
▶️ サポートサイト

Jupyter Notebook ファイルと同一フォルダ内に「data」フォルダを用意して、data フォルダ配下の章別フォルダにデータファイルを格納しています。

■ ライブラリのインポート
この記事で用いるライブラリをインポートします。

# インポート

# 数値計算
import numpy as np
import pandas as pd

# PyMC
import pymc as pm
import pytensor.tensor as pt
import arviz as az

# Rデータセットの読み込み
import rdata

# 可視化
import matplotlib.pyplot as plt
plt.rcParams['font.family'] = 'Meiryo'  # または import japanize_matplotlib

統計モデリング・サマリー

この記事で扱う統計モデリングの概要です。

■ 統計モデル
一次元空間の空間相関を考慮するベイズ統計モデルと考慮しないモデルを比較します。
空間相関は、ランダム切片「場所差」の事前分布に $${\text{ICAR}}$$ を指定して表現します。

■ モデリング手続き
前回記事の例題データの一部に意図的に欠測を設定します。
空間相関を考慮するモデルは前回記事のモデルと同じです。
空間相関を考慮しないモデルは場所差の事前分布に正規分布を指定します。
今回は欠測の補完に注力します。
いつもの細かなモデリング手続きには言及しないことにします。

欠測のある観測データの作成

◼️ データの読み込み
Y.RData ファイルを rdata ライブラリで変換して、pandas データフレームの data に読み込みます。

# 例題：一次元空間上の個体数分布のデータの読み込み p.242~

# RDataファイルの読み込み
data = rdata.read_rda('./data/ch11/Y.RData', default_encoding='ASCII')
data = pd.DataFrame(data)

# csvファイルの出力
# data.to_csv('./data/ch11/data.csv')

# データフレームの表示
print('data.shape: ', data.shape)
data.head()

【実行結果】
データの個数（標本サイズ）は 50 です。
１本の直線上に等間隔に配置された「50 区画」を観察して記録した「個体数」に関する仮想の観測データです。
区画（あるいは位置）は インデックス $${j}$$ で示されます。

【変数の説明】
区画 $${j}$$ の個体数 y です。

50 の区画で構成される一次元空間のイメージ図です。

🍀🍀🍀

◼️ データの可視化
テキスト p.244 図 11.2 に相当する区画と個体数 y・局所密度 m の関係を可視化します。

# 例題の一次元空間上の架空データの可視化 p.244 図11.2

# 描画設定
bottom, top = -0.8, 23  # y軸の範囲

# 描画領域の設定
plt.figure(figsize=(6.4, 4.8))

# 区画 j と個体数 y の散布図の描画
plt.plot(data.Y, 'o', ms=7, alpha=0.7, label='観測値')

# （通常は分からない）区画 j の局所密度 m の赤点線の描画
plt.plot(data.m, color='tab:red', ls='--', label='平均値')

# 修飾
plt.xlabel('位置 $j$', fontsize=12)
plt.ylabel('個体数 $y_j$', fontsize=12)
plt.ylim(bottom, top)
plt.legend();

【実行結果】
位置は区画のことです。
個体数は各区間でばらついていますが、近くの区画、隣り合う区画では値が似ていることを確認できます。
局所密度（赤点線）は滑らかな曲線であり、やはり、近くの区画、隣り合う区画では値が似ていることを確認できます。

🍀🍀🍀

◼️ 欠測データの作成
例題データの一部を欠測値 NaN に置き換えます。

# 例題の架空データ（一部欠損版） p.253 図11.6

# 欠測値を含むデータの作成
missing = [5, 8, 11, 12, 25, 26, 27, 28, 29]
data_missing = data.copy()
data_missing.loc[missing, 'Y'] = np.nan

# 結果の表示
print('data_missing.shape:', data_missing.shape)
data_missing.head(10)

【実行結果】
データの先頭 10 個を表示しています。
区画インデックス $${j}$$（Python なので 0 始まりです）の 5 と 8 の観測値 Y が欠測値 NaN です。
欠測値に置き換えた区画インデックスは $${[5, 8, 11, 12, 25, 26, 27, 28, 29]}$$ です。

◼️ 欠測データの可視化
欠測データを可視化します。
テキスト p.253 図 11.6 に相当します。

# 欠測データの可視化

# 描画設定
bottom, top = -0.8, 23  # y軸の範囲

# 描画領域の設定
plt.figure(figsize=(6.4, 4.8))

# 観測値の散布図の描画
plt.plot(data_missing.Y, 'o', ms=6, alpha=0.7, label='観測値')

# 欠測値の散布図の描画
plt.plot(missing, data.Y[missing], 'o', ms=7, color='black', label='欠測値')

# 欠測値の背景（グレー）の塗りつぶし描画
for fill_range in [[5, 5], [8, 8], [11, 12], [25, 29]]:
    plt.fill_between(np.arange(fill_range[0] - 0.5, fill_range[1] + 0.6, 0.1),
                     bottom, top, color='gray', alpha=0.2)

# 平均値の描画
plt.plot(data.m, color='tab:green', ls='--', label='平均値')

# 修飾
plt.xlabel('位置 $j$', fontsize=12)
plt.ylabel('個体数 $y_j$', fontsize=12)
plt.ylim(bottom, top)
plt.legend();

【実行結果】
グレイの縦長帯は欠測のある区画です。
黒点が（本来は分からない）欠測区画の観測値です。

２つのベイズ統計モデルで欠損値の補完（予測）を行います！

空間相関を考慮する階層ベイズモデル

ベイズ統計モデルの概要

前回記事の ICAR モデルと同じです。

◼️ 尤度と事前分布

• 区画 $${j}$$ の個体数 $${y_j}$$ はポアソン分布 $${p(y_j \mid \lambda_j)}$$ に従うとします。

• 平均パラメータ $${\lambda_j}$$ は線形予測子と対数リンク関数を用いて $${\log \lambda_j = \beta + r_j}$$ とします。

• $${\beta}$$ は全区画共通の大域的なパラメータであり、平均 0、標準偏差 100 の正規分布（ひらべったい正規分布）に従うとします。

• 場所差 $${r_j}$$ は局所的なパラメータであり、標準偏差パラメータ $${s}$$ の ICAR モデル $${\text{ICAR}(s)}$$ で階層事前分布を表現できるとします。

• $${r_j}$$ の標準偏差パラメータ $${s}$$ は全区画共通の大域的なパラメータであり、区間 $${[0, 10^3]}$$ の一様分布（幅が十分に広い連続一様分布）に従うとします。

◼️ ベイズ統計モデルの数式表現

パラメータの事後分布の推定

PyMC ライブラリでベイズ統計モデリングを実装します。

◼️ 隣接行列 W の作成
隣接行列 W は ICAR のパラメータです。
区画 $${j}$$ の１つ隣の区画に１、その他の区画に０を設定した行列です。

# 隣接行列Wの作成

# 標本サイズの設定
n = len(data_missing)

# 隣接行列の作成
W = np.array([[1 if abs(i - j)==1 else 0 for j in range(n)] for i in range(n)])

# 結果の表示
print('W.shape:', W.shape)
W

【実行結果】
行数・列数は例題データの標本サイズ 50 です。
２行目の区画 1 を見ましょう。
１が立っている区画は両隣に相当する 0 と 2 です。

🍀🍀🍀

◼️ モデルの定義

# モデルの定義 (A)空間相関を考慮しているモデル

# 設定と準備
coords = {'data': data_missing.index}   # 座標ラベルの設定：データ行の識別子

# モデリング
with pm.Model(coords=coords) as model_ma:
    
    # dataの定義: 目的変数=生存種子数Y ※NaNを含むので Data() を使えない
    Y = data_missing.Y

    # 事前分布
    # 切片β: 無情報事前分布 N(0,100)
    beta = pm.Normal('beta', mu=0, sigma=100)
    # 個体差rのICARモデルの標準偏差s: 無情報事前分布 U(0, 100) ※テキストはU(0, 1000)
    s = pm.Uniform('s', lower=0, upper=100)
    # 個体差r: 階層事前分布 ICAR(W)
    r_raw = pm.ICAR('r_raw', W=W, sigma=1, dims='data')  # ※ sigma=1
    r = pm.Deterministic('r', r_raw * s, dims='data')    # ※ r_raw * s

    # 線形予測子: 指数関数でlog(λ)を平均λに変換
    lam = pm.Deterministic('lam', pt.exp(beta + r), dims='data')

    # 尤度関数: 平均λのポアソン分布
    obs = pm.Poisson('obs', mu=lam, observed=Y, dims='data')

【実行結果】なし

【ワーニングの概要】
欠測値を含むことが原因となり、２つのワーニングが表示されます。
コード実行に支障がないので、ワーニングを無視します。

欠測値が含まれる観測値 Y（個体数）のデータ型は float 型です。
観測値 Y は尤度：ポアソン分布に従うので、PyMC内部では numpy を使って int 型（整数型）に変換しますが、float 型の NaN を含んでいるためワーニングを発しています。
このワーニングを無視してもコード実行に支障はありません。

PyMC が欠測値の自動補完を行うことを情報提供しています。

【コードの補足：観測値に欠測値を含む場合】
欠測値を含むデータを pm.Data() で定義できません（エラーになります）。
そこで、尤度のポアソン分布の observed 引数に Data化しない素の観測データを「直接」設定します。

Y = data_missing.Y
obs = pm.Poisson('obs', mu=lam, observed=Y, dims='data')

🍀🍀🍀

◼️ モデルの確認
モデルの数式と有向グラフを可視化します。
数式を表示します。

# モデルの表示
model_ma

【実行結果】

モデルの有向グラフを描画します。

# モデルの可視化
pm.model_to_graphviz(model_ma)

【実行結果】
複雑なモデルです。
$${\texttt{obs\_observed}}$$ は欠測のない観測値（41 個）、$${\texttt{obs\_unobserved}}$$ は欠測のある観測値（9 個）です。
欠測値のある観測値も MCMC サンプルが生成されます。

🍀🍀🍀

◼️ MCMC の実行
MCMC サンプリングを行います。
PyMC 標準のNUTS サンプラーを利用します。
合計 4000 個の MCMC サンプルを得ます。

%%time
# MCMCサンプリング
# thinning:未実施, パラメータの初期値:未設定, サンプリングアルゴリズム:NUTS

with model_ma:
    idata_ma = pm.sample(draws=1000, tune=1000, chains=4, random_seed=42)

【実行結果】
Divergences は０個です。

ワーニングが表示されますが、いったん無視します。

【コードの補足：欠測値がある場合の NUTSサンプラー】
観測値に欠測値を含む場合（自動補完する場合）、PyMC 標準の NUTS サンプラーを使う必要があります。
例えば nutpie などの外部の NUTS サンプラーを指定するとエラーになります。

パラメータ推定値の確認

ざっくり収束の確認などを行います。
観測値が欠損値を含む場合、MCMC 結果の評価指標が悪くなります。
（理由は分かりません…）
今回は MCMC 結果の評価をスルーしたり、しきい値を緩くします。

◼️ $${\widehat{R}}$$ の確認
$${\widehat{R}}$$ が 1.1 以下になっていることを確認します。
前回記事まではしきい値を 1.01 にしていましたが、今回は諸事情で 1.1 に引き上げました…

全パラメータの「$${\widehat{R}}$$ ＞1.1」の個数が０になればOKです。

# r_hat>1.1の確認
# 設定
idata_in = idata_ma      # idata名
threshold = 1.1          # しきい値

# しきい値を超えるR_hatの個数を表示
print((az.rhat(idata_in) > threshold).sum())

【実行結果】
すべてのパラメータの $${\widehat{R}}$$ は1.1 以下です。

◼️ トレースプロットの確認

# トレースプロットの表示
var_names = ['beta', 's', 'lam', 'r']
pm.plot_trace(idata_ma, var_names=var_names,
              backend_kwargs={'tight_layout': True});

【実行結果】
右のトレースプロットは、各チェーンが「毛糸玉」のようにゲジゲジと混ざり合い、ドリフトがないのでＯＫ（収束ＯＫ）としましょう。
（左の４本の chain は若干ずれ気味に見えます…）

🍀🍀🍀

◼️ 事後分布の要約統計量
パラメータの事後分布の推定値を確認します。

# 推論データの要約統計情報の表示
pm.summary(idata_ma, hdi_prob=0.95, var_names=var_names, round_to=3)

【実行結果】
$${\beta, s}$$ の事後平均と 95% HDI は前回記事の「欠測の無いデータの ICAR モデルの結果」と似た値を推定しています。

（参考：前回記事の欠測のないデータの ICAR モデルのsummary()）

空間相関を考慮しない階層ベイズモデル

ベイズ統計モデルの概要

場所差 $${r_j}$$ に空間相関を想定しない、つまり、区画ごとに独立していると仮定したモデルです。

◼️ 尤度と事前分布
空間相関を想定するモデルと相違する場所差 $${r_j}$$ の事前分布を確認します。

• 場所差 $${r_j}$$ は局所的なパラメータであり、平均 0、標準偏差 $${s}$$ の正規分布 $${\text{Normal}(0, s)}$$ を事前分布（階層事前分布）に指定します。

◼️ ベイズ統計モデルの数式表現

パラメータの事後分布の推定

PyMC ライブラリでベイズ統計モデリングを実装します。

◼️ モデルの定義

# モデルの定義 (B)空間相関を考慮していないモデル

# 設定と準備
coords = {'data': data_missing.index}   # 座標ラベルの設定：データ行の識別子

# モデリング
with pm.Model(coords=coords) as model_mb:
    
    # dataの定義: 目的変数=生存種子数Y ※NaNを含むので Data() を使えない
    Y = data_missing.Y

    # 事前分布
    # 切片β: 無情報事前分布 N(0,100)
    beta = pm.Normal('beta', mu=0, sigma=100)
    # 個体差rの正規分布の標準偏差s: 無情報事前分布 U(0, 1000)
    s = pm.Uniform('s', lower=0, upper=1000)
    # 個体差r: 階層事前分布 N(0,s)
    r = pm.Normal('r', mu=0, sigma=s, dims='data')

    # 線形予測子: 指数関数でlog(λ)を平均λに変換
    lam = pm.Deterministic('lam', pt.exp(beta + r), dims='data')

    # 尤度関数: 平均λのポアソン分布
    obs = pm.Poisson('obs', mu=lam, observed=Y, dims='data')

【実行結果】なし
※欠測値に関するワーニングは無視します。

◼️ モデルの確認
モデルの数式と有向グラフを可視化します。
数式を表示します。

# モデルの表示
model_mb

【実行結果】

モデルの有向グラフを描画します。

# モデルの可視化
pm.model_to_graphviz(model_mb)

【実行結果】
複雑なモデルです。
$${\texttt{obs\_observed}}$$ は欠測のない観測値（41 個）、$${\texttt{obs\_unobserved}}$$ は欠測のある観測値（9 個）です。
欠測値のある観測値も MCMC サンプルが生成されます。

🍀🍀🍀

◼️ MCMC の実行
MCMC サンプリングを行います。
PyMC 標準のNUTS サンプラーを利用します。
合計 4000 個の MCMC サンプルを得ます。

%%time
# MCMCサンプリング
# thinning:未実施, パラメータの初期値:未設定, サンプリングアルゴリズム:NUTS

with model_mb:
    idata_mb = pm.sample(draws=1000, tune=1000, chains=4, random_seed=42,
                         target_accept=0.9)

【実行結果】
Divergences は０個です。

ワーニングが表示されますが、いったん無視します。

パラメータ推定値の確認

ざっくり収束の確認などを行います。
観測値が欠損値を含む場合、MCMC 結果の評価指標が悪くなります。
（理由は分かりません…）
今回は MCMC 結果の評価をスルーしたり、しきい値を緩くします。

◼️ $${\widehat{R}}$$ の確認
$${\widehat{R}}$$ が 1.1 以下になっていることを確認します。
全パラメータの「$${\widehat{R}}$$ ＞1.1」の個数が０になればOKです。

# r_hat>1.1の確認
# 設定
idata_in = idata_mb      # idata名
threshold = 1.1          # しきい値

# しきい値を超えるR_hatの個数を表示
print((az.rhat(idata_in) > threshold).sum())

【実行結果】
すべてのパラメータの $${\widehat{R}}$$ は1.1 以下です。

◼️ トレースプロットの確認

# トレースプロットの表示
var_names = ['beta', 's', 'lam', 'r']
pm.plot_trace(idata_mb, var_names=var_names,
              backend_kwargs={'tight_layout': True});

【実行結果】
右のトレースプロットは、各チェーンが「毛糸玉」のようにゲジゲジと混ざり合い、ドリフトがないのでＯＫ（収束ＯＫ）としましょう。
（左の４本の chain は若干ずれ気味に見えます…）

🍀🍀🍀

◼️ 事後分布の要約統計量
パラメータの事後分布の推定値を確認します。

# 推論データの要約統計情報の表示
pm.summary(idata_mb, hdi_prob=0.95, var_names=var_names, round_to=3)

【実行結果】
$${\beta}$$ の事後平均は空間相関を考慮する階層ベイズモデルと似た値を推定しています。
しかし、ばらつきパラメータ $${s}$$ の事後平均は空間相関を考慮する階層ベイズモデルと少々相違している感じです。

（参考：空間相関を考慮する階層ベイズモデルのsummary()）

欠測値補完の結果確認

MCMC サンプルを用いて、平均個体数 $${\lambda}$$ ベースの予測を可視化します。
テキスト p.254 図 11.7 の予測チャートに相当します。
※ $${\lambda}$$ の MCMC サンプルを用いても大丈夫です。

# 空間統計モデルによる欠測データの予測 p.254 図11.7 ※テキストと異なる結果になった

## 設定と準備
# λの予測値の算出
lam_pred_ma = np.exp(betas_ma + rs_ma)
lam_pred_mb = np.exp(betas_mb + rs_mb)
# λの予測値の10%,50%,90%パーセンタイル点の算出
lam_pred_quantiles_ma = np.quantile(lam_pred_ma, q=[0.1, 0.5, 0.9], axis=1)
lam_pred_quantiles_mb = np.quantile(lam_pred_mb, q=[0.1, 0.5, 0.9], axis=1)
# 描画設定
bottom, top = -0.8, 23  # y軸の範囲

## 描画
# 描画領域の設定
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(7.7, 9.3), tight_layout=True)

# (A)空間相関を考慮しているモデル
# 観測値の散布図の描画
ax1.plot(data_missing.Y, 'o', ms=6, alpha=0.7, label='観測値')
# 欠測値の散布図の描画
ax1.plot(missing, data.Y[missing], 'o', ms=7, color='black', label='欠測値')
# 欠測値の背景（グレー）の塗りつぶし描画
for fill_range in [[5, 5], [8, 8], [11, 12], [25, 29]]:
    ax1.fill_between(np.arange(fill_range[0] - 0.5, fill_range[1] + 0.6, 0.1),
                     bottom, top, color='gray', alpha=0.2)
# 観測値のベースになった平均λの折れ線グラフ（緑）の描画
ax1.plot(data.m, color='tab:green', ls='--', label='観測値の平均λ')
# λの予測値（中央値）の折れ線グラフ（赤）の描画
ax1.plot(lam_pred_quantiles_ma[1, :], color='tab:red', label='λの予測値:中央値')
# λの予測値の80%区間の塗りつぶし描画
ax1.fill_between(data.index,
                 lam_pred_quantiles_ma[0, :],
                 lam_pred_quantiles_ma[2, :],
                 color='lightpink', alpha=0.5, label='λの予測値:80%区間')
# 修飾
ax1.set_title('(A) 空間相関を考慮しているモデル', fontsize=12)
ax1.set_xlabel('位置 $j$', fontsize=12)
ax1.set_ylabel('個体数 $y_j$', fontsize=12)
ax1.set_ylim(bottom, top)
ax1.legend(bbox_to_anchor=(1.39, 1));

# (B)空間相関を考慮していないモデル
# 観測値の散布図の描画
ax2.plot(data_missing.Y, 'o', ms=6, alpha=0.7, label='観測値')
# 欠測値の散布図の描画
ax2.plot(missing, data.Y[missing], 'o', ms=7, color='black', label='欠測値')
# 欠測値の背景（グレー）の塗りつぶし描画
for fill_range in [[5, 5], [8, 8], [11, 12], [25, 29]]:
    ax2.fill_between(np.arange(fill_range[0] - 0.5, fill_range[1] + 0.6, 0.1),
                     bottom, top, color='gray', alpha=0.2)
# 観測値のベースになった平均λの折れ線グラフ（緑）の描画
ax2.plot(data.m, color='tab:green', ls='--', label='観測値の平均λ')
# λの予測値（中央値）の折れ線グラフ（赤）の描画
ax2.plot(lam_pred_quantiles_mb[1, :], color='tab:red', label='λの予測値:中央値')
# λの予測値の80%区間の塗りつぶし描画
ax2.fill_between(data.index,
                 lam_pred_quantiles_mb[0, :],
                 lam_pred_quantiles_mb[2, :],
                 color='lightpink', alpha=0.5, label='λの予測値:80%区間')
# 修飾
ax2.set_title('(B) 空間相関を考慮していないモデル', fontsize=12)
ax2.set_xlabel('位置 $j$', fontsize=12)
ax2.set_ylabel('個体数 $y_j$', fontsize=12)
ax2.set_ylim(bottom, top)
ax2.legend(bbox_to_anchor=(1, 1));

【実行結果】

（参考：前回記事の欠測のないデータの ICAR モデルの予測）

【観察】
空間相関を考慮する階層ベイズモデルの方は、前後の区画の値から「滑らかに」予測できている感じがします。
さらには、前回の欠測値のないデータによる ICAR モデルの予測と似た予測ができています。

一方で、空間相関を考慮しない階層ベイズモデルの方は、欠測区画の予測は大きくぶれて（暴れて）いる感じです。
欠測のない区画も予測値は隣の区画と外れた「ジグザク」な予測になっています。

🍀🍀🍀

【テキストの解説】
テキストの解説を引用いたします。

✅️ 空間相関を考慮する階層ベイズモデル

✅️ 空間相関を考慮しない階層ベイズモデル

まとめ

観測値（PyMCのobservedデータ）に欠測がある場合、ベイズ統計モデルが欠測部分を推定して補完してくれることを学びました。

データが生成された現象をベイズ統計モデルが適切に表現できている場合は補完値（予測値）は妥当そうな値になりました。

一方で適切ではないベイズ統計モデルの場合には、補完値（予測値）の妥当性に懸念が生じることも分かりました。

🍀🍀🍀

シリーズ記事 vol.1 に掲載した統計モデリング観を振り返りたいと思います。

💡 統計モデリング観

分析者自身が、統計モデルを選択・作成して、モデルが現象を説明・予測できているか評価する。

この一文に何とか近づけるようにと努力して、ブログ原稿を書いたつもりです。

原稿を書くことはテキストとの格闘です。

テキストを読めば読むほどに疑問が積み上がっていく状況。
ChatGPT のおかげで、サクッと疑問を解決できたこともあれば、ChatGPT のせいで、疑問がますます複雑になって理解不能になったこともありました。

「原稿書き⇒疑問発生⇒疑問解消（or 解決諦め）⇒原稿見直し」のサイクルに長い時間がかかりました。
原稿の文章案とコード案を何度も見返して、修正修正修行の繰り返しの日々だったような気がします。

何が残ったかは分かりませんが、ひとまず、何かが出来上がったことを額面通りに受け止めようと思います。

統計学はなぞの多い学問です…

🍀🍀🍀

最後に。
疑問解消のための ChatGPT との応答を通じて多くの知識を得ました。
この知識を ChatGPT に「教科書」の名でまとめてもらっています。
このブログシリーズで積み上がった「教科書」の一覧です。

今回のブログは以上です。

このブログシリーズは今回が最終回です。

久保先生、楽しく統計モデリングを学ぶことができました。
ありがとうございました。

ブログを読んでくださったみなさん、長丁場で長文のブログを根気よく読んでくださり、ありがとうございました。

シリーズの記事

前の記事

目次

ブログの紹介

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

１．のんびり統計

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

２．実験！たのしいベイズモデリング１＆２を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の教科書です。
よかったらぜひ、お試しくださいませ。

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