論文一覧に戻る 📚 用語解説(ジャストインタイム型データサイエンス教育)
ランダム効果モデル
Random Effects Model (RE)
個体固有効果を「ランダム変数」として扱うパネルモデル。FE より仮定は強いが効率が良い。
パネル分析REREランダム効果random effects
📍 文脈💡 30秒結論

📍 あなたが今見ているもの

論文中に 「ランダム効果モデル」として登場する用語。

ランダム効果モデル とは:個体固有効果を「ランダム変数」として扱うパネルモデル。FE より仮定は強いが効率が良い。

💡 30秒で分かる結論

🎨 直感で掴む — ランダム効果モデル

SSDSE-B-2026 で 47 都道府県 × 複数年のパネルを考えるとき、 「東京は他と違う」「沖縄は他と違う」というような 個体固有のクセ をどう扱うかが論点になります。 ランダム効果 (Random Effects, RE) モデルは、 この個体クセを 母集団からランダムに引かれた変量 と見なして、 47 個分のダミー変数を立てずに 分散 $\sigma_\mu^2$ を 1 つ推定 することで効率的に処理します。

比喩で言えば、 学校の生徒成績を分析するときに「クラス A」「クラス B」とダミーを大量に作る代わりに、 「クラスごとのバラつき」という 1 つの分散 でまとめてしまう発想です。 ダミーで自由度を消費しないため、 サンプルが少ない場合や時間不変変数(地域の気候など)の係数を推定したいときに重宝します。

📐 定義・数式

$$ y_{it} = \alpha + \boldsymbol{x}_{it}^{\top}\boldsymbol{\beta} + u_i + \varepsilon_{it}, \quad u_i \sim \mathcal{N}(0, \sigma_u^2),\ \varepsilon_{it} \sim \mathcal{N}(0, \sigma_\varepsilon^2) $$

重要な仮定は $\operatorname{Cov}(u_i, \boldsymbol{x}_{it}) = 0$ — 個体固有効果と説明変数が無相関であること(直交性)。 これが満たされれば RE 推定量は BLUE となります。 SSDSE-B-2026 で婚姻率 $x$ と地域固有効果 $u_i$ が相関していそうなら、 RE の前提が崩れます(その場合は FE)。

🔬 数式を言葉で読み解く — 記号 → 意味

記号意味SSDSE-B-2026 解釈
$u_i$個体ランダム効果(平均 0 の確率変数)「都道府県ごとの未観測の地域特性」を 1 つの分散で要約
$\sigma_u^2$個体間分散地域差の大きさ
$\sigma_\varepsilon^2$個体内(時間内)分散同じ都道府県内の年次変動
$\rho = \frac{\sigma_u^2}{\sigma_u^2 + \sigma_\varepsilon^2}$級内相関(ICC)同じ都道府県の観測がどれだけ似るか

推定は GLS(一般化最小二乗)。 OLS との違いは観測値を $1-\hat{\theta}$ だけ平均から引いた「部分組内変換」を行う点で、 FE(完全変換)と OLS(無変換)の中間に位置します。

📖 詳細な解説

この用語は、 統計データ解析・データサイエンスの世界で重要な概念の1つです。 ジャストインタイム型学習では、 必要なときに参照し、 関連概念と合わせて学ぶことで定着を図ります。

基本的な定義

この用語の基本的な意味、 数学的定義、 直感的理解について、 上記の3つの概念マップを通じて、 関連する用語と一緒に把握しましょう。

使い時の判断基準

Python による実装例

▼ コード解説(ランダム効果モデルの基本(PanelOLS RE))
🎯 解説: linearmodels.RandomEffects で個体差 u_i を確率変数として扱う。 個体差が X と無相関なら FE より効率的。
📥 入力例: y = α + βX + u_i + ε_it u_i ~ N(0, σ_u²)
📤 実行例: Between R² + Within R² の合成 全体 R² ≈ 0.95 σ_u² と σ_ε² を分離推定
💬 読み方: RE は「個体差は確率変数で X と独立」と仮定。 仮定が成り立てば SE が小さく検定力が高い。
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
import pandas as pd
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

# SSDSE データの読み込み
df = pd.read_csv('data/raw/SSDSE-B-2026.csv', encoding='cp932')

# 基本統計
df.describe()
df.info()

# 可視化
df.hist(bins=30, figsize=(15, 10))
plt.show()

📖 包括的解説 — この概念を完全マスター

📍 学習の3ステップ

  1. 定義を理解する:この概念は何か? 数式や条件を確認
  2. 具体例を見る:実データ(SSDSE 等)で計算してみる
  3. 応用する:自分のデータに適用、 結果を解釈

🔧 Python実装パターン

▼ コード解説(FE と RE の比較)
🎯 解説: 同じデータで FE と RE を推定して係数を比較。 大きく異なれば仮定が崩れている。
📥 入力例: FE 係数 vs RE 係数 両者の差を計算
📤 実行例: 高齢化率の係数: FE=0.85, RE=0.92 差 0.07 → ハウスマン検定で有意か確認
💬 読み方: 係数差が大きいときは FE が安全。 差が小さくハウスマン検定で棄却されなければ RE で効率を取る。
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
# 基本パターン
import pandas as pd
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt
import seaborn as sns

# データ読み込み
df = pd.read_csv('data/raw/SSDSE-B-2026.csv', encoding='cp932')

# 基本統計量
df.describe()

# 可視化
sns.pairplot(df[['食料費', '教育費', '住居費']])
plt.show()

📚 統計概念マップでの位置

このページの上にある3つの概念マップ(関係マップ、 包含マップ、 ツリーマップ)でこの概念の位置づけが視覚的に分かります。 関連手法を辿って学習を進めましょう。

🎯 SSDSE-B-2026 で挑戦

統計データ活用コンペティションのSSDSE-B-2026データは、 47都道府県の社会経済データ。 この概念を使って以下のような分析ができます:

💡 よく使うコマンド集

機能 Python (pandas) Python (scipy)
要約統計df.describe()stats.describe()
平均df.mean()np.mean()
標準偏差df.std()np.std()
相関df.corr()stats.pearsonr()
t検定stats.ttest_ind()
回帰stats.linregress()
分布フィッティングstats.norm.fit()

🚧 一般的な落とし穴と対策

📊 結果報告の標準フォーマット

🌐 関連分野での応用

🎓 さらに学ぶための文献

🔗 統計用語ネットワーク

この概念は、 他の多くの統計概念と密接に関連しています。 ジャストインタイム型学習では、 必要に応じて関連用語へジャンプしながら全体像を構築します。

主要な関連概念のグループ

グループ 主要概念
記述統計平均、 中央値、 最頻値、 分散、 標準偏差、 共分散、 相関係数
可視化ヒストグラム、 散布図、 箱ひげ図、 ヒートマップ
推測統計標本平均、 標準誤差、 信頼区間、 p値、 有意水準
確率分布正規分布、 t分布、 χ²分布、 F分布、 二項分布
仮説検定t検定、 F検定、 χ²検定、 ノンパラ検定
回帰単回帰、 重回帰、 OLS、 Ridge、 LASSO
分類ロジスティック回帰、 決定木、 SVM、 k-NN
教師なし学習クラスタリング、 PCA、 因子分析
時系列ARIMA、 VAR、 指数平滑法、 自己相関
因果推論DiD、 IV、 傾向スコア、 交絡変数
前処理標準化、 正規化、 欠損値処理、 多重共線性対策
評価R²、 残差、 CV、 RMSE、 効果量

学習順序の推奨

  1. 記述統計(平均、 分散、 標準偏差)
  2. 可視化(ヒストグラム、 散布図)
  3. 確率分布(正規分布)
  4. 推測統計(標準誤差、 信頼区間、 p値)
  5. 仮説検定(t検定、 χ²検定)
  6. 相関と回帰(単回帰、 重回帰)
  7. 多変量解析(PCA、 クラスタリング)
  8. 機械学習(決定木、 RF、 NN)
  9. 時系列・因果推論(応用)

📝 実践練習 — SSDSE-B-2026 で挑戦

初級課題

  1. 東北6県の家計食料費の基本統計量を計算
  2. 食料費のヒストグラムを描く
  3. 食料費と教育費の散布図を描く
  4. 都道府県を「東日本/西日本」に分け、 平均を比較

中級課題

  1. 家計支出 5項目で相関行列を作成、 ヒートマップ可視化
  2. 食料費 → 教育費の単回帰を実行、 残差分析
  3. 家計5項目で PCA を実施、 バイプロット表示
  4. k-means (k=3) で都道府県をクラスタリング、 解釈

上級課題

  1. 地域別の家計パターンに有意差があるか ANOVA で検定
  2. 重回帰で教育費を予測、 多重共線性を VIF で確認
  3. Ridge/LASSO で正則化、 CV で α を最適化
  4. 階層クラスタリングと Ward 法で都道府県を分類、 デンドログラム作成

📚 統計学習の総合ガイド

🎯 学習目標

このページの概念をマスターすることで、 以下のスキルが身につきます:

📊 SSDSE-B-2026 データの構造

このコンペの主要データセット(SSDSE-B-2026)の構造:

🔍 主要な変数群

カテゴリ 変数例
人口総人口、 年齢別人口、 性別人口
人口動態出生数、 死亡数、 合計特殊出生率、 婚姻数
気候気温、 降水量、 降水日数
教育幼小中高校数、 教員数、 生徒数、 大学進学率
経済求職件数、 求人件数、 旅館数
医療病院数、 診療所数、 歯科診療所
家計消費支出、 食料費、 住居費、 教育費等の項目別

💡 ジャストインタイム型学習

このガイドは「必要なときに必要な知識」を提供する設計:

🛠️ Python データサイエンス環境

▼ コード解説(ハウスマン検定(仕様検査))
🎯 解説: ハウスマン検定で「RE の前提が成立するか」を χ² で検定。 棄却されれば FE を採用。
📥 入力例: H0: RE と FE の係数が同じ H1: 異なる(個体効果と X が相関)
📤 実行例: χ² = 12.3, df=3, p=0.006 → H0 棄却、 FE 採用
💬 読み方: p<0.05 で FE 一択。 p≥0.05 でも理論的に X と u_i が相関しそうなら FE が無難。
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
# 必須ライブラリのインストール
pip install pandas numpy scipy statsmodels scikit-learn matplotlib seaborn

# 標準的なインポート
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from scipy import stats
from sklearn.preprocessing import StandardScaler
from sklearn.model_selection import train_test_split
from sklearn.metrics import r2_score, mean_squared_error

# 日本語表示の設定(matplotlib)
plt.rcParams['font.family'] = 'Hiragino Sans'
plt.rcParams['axes.unicode_minus'] = False

# データ読み込み(SSDSE は cp932 エンコーディング)
df = pd.read_csv('data/raw/SSDSE-B-2026.csv', encoding='cp932')
print(df.shape)
print(df.head())
print(df.describe())

🌟 効果的なEDAテンプレート

▼ コード解説(分散成分の推定(σ_u と σ_ε))
🎯 解説: 個体間分散 σ_u² と個体内分散 σ_ε² を分離推定。 σ_u² が大きいほど「県固有要因」が強い。
📥 入力例: df, EntityEffects=True 分散の分解
📤 実行例: σ_u² ≈ 1.5(県固有) σ_ε² ≈ 0.3(年内ばらつき) ICC = σ_u²/(σ_u²+σ_ε²) ≈ 0.83
💬 読み方: ICC(級内相関)が 0.5 を超えれば「同じ県の年同士は強く相関」。 パネル分析が必要なサイン。
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
def quick_eda(df, target=None):
    \"\"\"探索的データ分析の基本テンプレート\"\"\"
    print(f"Shape: {df.shape}")
    print(f"\\nColumn types:\\n{df.dtypes}")
    print(f"\\nMissing values:\\n{df.isnull().sum()}")
    print(f"\\nBasic stats:\\n{df.describe()}")

    # 数値列の可視化
    numeric_cols = df.select_dtypes(include=[np.number]).columns
    df[numeric_cols].hist(bins=20, figsize=(15, 10))
    plt.tight_layout()
    plt.show()

    # 相関ヒートマップ
    if len(numeric_cols) > 1:
        plt.figure(figsize=(12, 10))
        sns.heatmap(df[numeric_cols].corr(), annot=True, fmt='.2f',
                    cmap='RdBu_r', center=0)
        plt.show()

    # ターゲットがあれば散布図行列
    if target and target in df.columns:
        sns.pairplot(df[numeric_cols[:5]], hue=target if df[target].dtype == 'O' else None)
        plt.show()

📈 報告書テンプレート

分析結果を報告する際の標準的な構成:

  1. 背景・目的:なぜこの分析が必要か
  2. データ:出所、 サンプルサイズ、 期間
  3. 方法:使用した統計手法、 仮定
  4. 結果:図表、 統計量、 検定結果
  5. 解釈:結果が何を意味するか
  6. 限界:分析の制約
  7. 結論:要点まとめ、 今後の課題

🗺️ 統計手法選択フローチャート

Q1: 何を知りたい?

Q2: データの種類は?

Q3: サンプルサイズは?

Q4: 仮定は?

📏 効果量の参照表

p値だけでなく効果量も併記するのが現代統計の標準。 主要な指標と Cohen の解釈基準:

統計量 効果量
2群平均差Cohen's d0.20.50.8
相関r0.10.30.5
線形回帰0.020.130.26
ANOVAη² (eta²)0.010.060.14
χ²Cramér's V0.10.30.5
ロジスティックOdds Ratio1.52.54.0

🗺️ 概念マップ — 3つの視点で体系を理解する

ランダム効果モデル がデータサイエンスの体系の中でどこに位置するかを、 3つの異なる視点で可視化します。 同じ情報でも見方を変えると気付きが変わります。

📍 体系階層のパス

🌐 体系階層に未登録

① 🔗 関係マップ — 「他の手法とどう繋がっているか」

中心の概念から放射状に、 前提・兄弟・発展形・応用先などの関係性を矢印で結びます。 横の繋がりを見るのに最適。 ノードをドラッグ、 ホイールでズーム、 クリックで遷移

凡例:現在の用語上位カテゴリ兄弟(並列)前提発展形応用先2階層先

② ⭕ 包含マップ — 「どのカテゴリに含まれているか」

大きな円が小さな円を包含する Circle Packing 図。 「ランダム効果モデル」は緑色でハイライト

📍現在地:統計・データサイエンス

③ 🌳 ツリーマップ — 「面積で見るボリューム比較」

長方形を入れ子に分割した Treemap 図。 各分野の規模感を面積で比較。 「ランダム効果モデル」は緑色でハイライト

🎯 3つのマップの使い分け

マップ 分かること こんな時に見る
🔗 関係マップ手法間の横の関係(前提→発展→応用)「次に何を学べばよい?」 学習順序の判断
⭕ 包含マップ分類体系の入れ子構造(上位⊃下位)「この手法はどんなジャンルに属する?」
🌳 ツリーマップ分野の規模比較(面積=ボリューム)「データサイエンス全体の俯瞰像」

💡 ジャストインタイム学習のヒント:3つの視点を行き来することで、 概念を多角的に理解できます。 包含マップやツリーマップはズーム/ドリルダウンで大分類から細部まで探索できます。

🔖 キーワード索引(拡張)

変量効果モデル周辺の重要語をクイックアクセス:

変量切片モデル 変量傾きモデル ICC(級内相関) 分散成分 縮小推定(shrinkage) BLUP 過小収束 / Singular fit REML vs ML 変量と固定の選択 statsmodels MixedLM pymer4 / lme4 PyMC(ベイズ)

🧮 SSDSE-B 実値計算 — 「地方ブロック」を変量効果にした賃金モデル

47都道府県を 8 地方ブロックに分け、 「ブロック」を変量効果として扱う線形混合モデル。 目的変数は「現金給与総額」、 固定効果は「製造品出荷額(log)」「総人口(log)」。

▼ コード解説(ランダム切片モデル(混合効果))
🎯 解説: statsmodels.MixedLM で線形混合モデル。 県ごとに切片を確率変数として推定。
📥 入力例: groups='都道府県' re_formula='1'
📤 実行例: 切片の分散 ≈ 1.5 各県の切片予測値 BLUP
💬 読み方: ランダム切片は「県ごとに死亡率のベースラインが違う」状況をモデル化。 教育・医療データに頻出。
▼ コード解説(ランダム傾きモデル)
🎯 解説: 切片だけでなく傾きも県ごとに確率変数化。 県によって関係の強さが異なる場合に有効。
📥 入力例: re_formula='1 + 高齢化率'
📤 実行例: 切片と傾きの共分散 傾きの分散 ≈ 0.02
💬 読み方: ランダム傾きは「県ごとに高齢化率の効きが違う」ケースに対応。 多くの社会データでは切片のみで十分。
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
import numpy as np
import pandas as pd
import statsmodels.formula.api as smf

df = pd.read_csv('data/raw/SSDSE-B-2026.csv', encoding='cp932', header=1)
df.columns = [c.strip() for c in df.columns]

def region(code):
    code = int(code)
    if code == 1:                  return '北海道'
    if 2  <= code <= 7:            return '東北'
    if 8  <= code <= 14:           return '関東'
    if 15 <= code <= 23:           return '中部'
    if 24 <= code <= 30:           return '近畿'
    if 31 <= code <= 35:           return '中国'
    if 36 <= code <= 39:           return '四国'
    return '九州沖縄'

df['region']  = df['地域コード'].apply(region)
df['log_mfg'] = np.log(df['製造品出荷額等'])
df['log_pop'] = np.log(df['総人口'])

md  = smf.mixedlm('現金給与総額 ~ log_mfg + log_pop',
                  data=df, groups=df['region'])
res = md.fit(reml=True)
print(res.summary())
icc = res.cov_re.iloc[0,0] / (res.cov_re.iloc[0,0] + res.scale)
print(f'ICC = {icc:.3f}')

典型的な出力例: 固定効果 log_mfg=+0.012, log_pop=+0.18、 ブロック分散 σ_u²≈0.4、 残差分散 σ_e²≈1.2、 ICC≈0.25。 「ブロック間で平均賃金水準に 25% 程度のばらつきがある」と読める。 ICC が 0.05 を下回るなら変量効果を入れるご利益はほぼ無く、 通常 OLS で十分。

BLUP(各ブロックの予測切片偏差)

▼ コード解説(BLUP(最良線形不偏予測量))
🎯 解説: 各県の固有効果 u_i の予測値(BLUP)を取得。 「平均からどれだけ離れているか」を県別に可視化。
📥 入力例: model.random_effects 属性
📤 実行例: 沖縄: u_i = -2.1(平均より低い死亡率) 秋田: u_i = +1.5(高め)
💬 読み方: BLUP は「個体固有のシフト量」。 全国平均からの逸脱を県別に解釈でき、 地域政策の参考になる。
▼ コード解説(年度効果(時間ランダム効果))
🎯 解説: TimeEffects に対応するランダム効果を追加。 年度ごとの全国共通ショックを確率変数化。
📥 入力例: 年度をクラスターとした追加 RE
📤 実行例: 2020 年度のランダム効果 ≈ +0.8(コロナ) 他年度はほぼ 0
💬 読み方: 年度を RE にすると「平均的な年からの逸脱」が分かる。 通常は FE(年度ダミー)で十分。
1
2
3
blup = res.random_effects
for r, s in sorted(blup.items(), key=lambda x: x[1].iloc[0]):
    print(f'{r:<6s}  BLUP切片偏差 = {s.iloc[0]:+.3f}')

⚠️ 変量効果モデルの落とし穴 — 実務で必ず踏む 6 つ

① グループ数が少なすぎて分散成分が推定できない

変量効果に指定したグループの水準数が 5 未満だと、 群間分散 σ_u² の推定が不安定になり、 Singular fit(分散がほぼゼロに収束)が起きる。 SSDSE のように 8 ブロックでも厳しく、 6 ブロック以下なら基本的に変量効果として扱わず、 固定効果(ダミー変数)にした方が良い。 経験則として「変量効果のためには水準数 ≥ 5、 望ましくは ≥ 10」が目安とされる(Gelman & Hill, 2007)。 これを知らずに 3 群で REML を回すと、 信頼区間がほぼ無意味になる。

② 「変量効果 = 完全にランダム」と誤解する

「変量効果」と書かれるからといって、 各群の効果が完全に独立にサンプリングされるわけではない。 実際には「群効果が共通の正規分布 N(0, σ_u²) から抽出された」という階層的仮定を置いており、 群間の縮小(shrinkage)が暗黙的に働く。 そのため極端なサンプル数の群(n=1 や n=2)の推定値は集団平均に強く引き寄せられる。 これを「BLUP の縮小」と呼び、 古典 OLS の群ダミーとは予測値が大きく変わる。 縮小は予測誤差を下げるが、 個別群の効果を強調したい記述目的には不向き。

③ REML と ML の使い分けを意識しない

分散成分推定では REML(制限付き最尤)が ML より不偏に近く既定で推奨されるが、 固定効果の構造を尤度比検定で比較するときは ML に切り替える必要がある。 REML 同士の尤度は固定効果の構造を変えると単純比較できないからである。 statsmodels の MixedLM.fit() は既定 reml=True、 lme4 の lmer も既定 REML=TRUE。 AIC/BIC でモデル選択するときは特に注意する。 一方、 分散成分の有意性検定は REML のままで尤度比検定を実施するのが定石。

④ Singular fit を無視して結論を出す

変量傾きを入れると分散共分散行列が半正定値の境界に張り付き、 「σ_slope ≈ 0」「相関 ±1」となる Singular fit が頻発する。 これは「データから変量傾きを支持する情報がない」シグナルで、 モデルを変量切片のみに簡略化するか、 ベイズ事前分布を使った正則化(rstanarm, brms, PyMC)に切り替えるべき。 警告を無視して結果を読むと SE が過大評価され、 固定効果の検定がおかしくなる。 lme4 では isSingular() で診断可能。

⑤ 自由度・p値の計算法が実装でバラバラ

混合モデルでは「自由度をどう数えるか」が未解決問題で、 R の lme4 はあえて p 値を出さない。 一方、 lmerTest は Satterthwaite 近似、 SAS は Kenward-Roger、 statsmodels は Wald 検定で済ます。 同じデータでも実装によって p 値が 0.04 / 0.06 と境界を跨ぐことがある。 論文では「どの近似を使ったか」を必ず明記する。 ベストプラクティスはブートストラップによる信頼区間、 もしくはベイズ事後信用区間で報告すること。

⑥ 「Fixed vs Random」を経済学とは逆に覚える

統計学の「変量効果モデル」と、 経済学の「固定効果 vs 変量効果」の用語は概念が微妙に違う。 経済学で言う「固定効果」は群ダミー、 「変量効果」は群が誤差項と無相関と仮定したモデルで、 Hausman 検定で選ぶ。 統計学では「変量効果=階層モデルの上位レベル」と呼ぶ。 論文間を読み比べるときに混乱しないよう、 用語の出典分野を確認する習慣をつける。 SSDSE をパネル化して扱う際は経済学流の用法に従う方が論文と整合する。

🐍 Python 実装バリエーション — statsmodels / pymer4 / PyMC

1. statsmodels.formula.api.mixedlm(標準・REML)

▼ コード解説(F 検定(プール OLS vs RE))
🎯 解説: Breusch-Pagan LM 検定で「個体効果があるか」を検定。 p<0.05 なら RE/FE が必要。
📥 入力例: H0: σ_u² = 0(個体効果なし、 プール OLS で OK)
📤 実行例: LM 統計量 = 285.3, p<0.001 → 個体効果あり、 RE/FE 推奨
💬 読み方: プール OLS で済むケースは稀。 BP-LM 検定はパネル分析の最初の診断として推奨。
▼ コード解説(ランダム効果の信頼区間)
🎯 解説: u_i の 95% 信頼区間で「県固有効果が有意か」を判定。
📥 入力例: 各県の BLUP ± 1.96 × SE
📤 実行例: 沖縄: -2.1 ± 0.4(有意に低い) 秋田: +1.5 ± 0.5(有意に高い)
💬 読み方: 0 を含まない CI の県は「他県と統計的に異なる」。 政策ターゲット選定に使える。
1
2
3
4
5
import statsmodels.formula.api as smf
m1 = smf.mixedlm('y ~ x1 + x2', df, groups=df['region']).fit(reml=True)
m2 = smf.mixedlm('y ~ x1 + x2', df, groups=df['region'],
                 re_formula='~x1').fit(reml=True)
print(m1.summary()); print(m2.summary())

2. pymer4 — R の lme4 を Python から呼ぶ

R の lme4 と等価のインターフェースで、 p値が Satterthwaite 近似で出る。 学術論文で lme4 を再現したい場合に便利。

▼ コード解説(予測(新しい県・将来年度))
🎯 解説: predict() で BLUP を考慮した予測。 新しい年度・新しい県(or 既存県の将来)を予測。
📥 入力例: 新規データ X_new
📤 実行例: 県別予測値(BLUP 込み) RE 不要な場合の予測値
💬 読み方: RE モデルでは既存県の予測には BLUP を加算、 新規県の予測には加算しない(平均効果のみ)。
▼ コード解説(RE 推定の効率性)
🎯 解説: RE と FE の SE 比較。 RE の方が SE が小さい(情報量を多く使う)。
📥 入力例: FE の SE と RE の SE
📤 実行例: FE: SE=0.05 RE: SE=0.03 → RE は約 40% 効率が高い
💬 読み方: RE の効率性は「個体間情報も使う」ことから来る。 ただし仮定が崩れれば偏りが入るため、 トレードオフ。
1
2
3
4
from pymer4.models import Lmer
m = Lmer('y ~ x1 + x2 + (1|region)', data=df).fit()
print(m.summary())
print(m.fixef)

3. PyMC でベイズ階層モデル

▼ コード解説(random effects の Python 実装(コードブロック 13))
🎯 解説: random effects を SSDSE-B-2026 都道府県データで実行する Python コード。 47 都道府県 × 112 指標のパネルデータを使い、 公的統計の実値に基づいて手法を可視化する。 教育用ハンズオン教材として scipy/pandas/sklearn の標準ライブラリで完結。
📥 入力例: data/raw/SSDSE-B-2026.csv encoding=shift_jis(または cp932), skiprows=1 47 都道府県 × 5 年(2019-2023)= 235 行 数値特徴量 100+ 列(家計・人口・教育・医療)
📤 実行例: random effects の主要出力を確認 都道府県別・年度別の指標が出力される 例: 秋田・沖縄・東京の差異が顕著に表れる プロット・統計量・モデルパラメータ等の結果を取得
💬 読み方: random effects の結果は SSDSE-B-2026 の地域特性を反映。 都市圏/地方圏/観光圏でパターンが異なり、 47 県の異質性を視覚化・定量化できる。 統計データ活用コンペの実分析に直接転用可能。
▼ コード解説(random effects の Python 実装(コードブロック 14))
🎯 解説: random effects を SSDSE-B-2026 都道府県データで実行する Python コード。 47 都道府県 × 112 指標のパネルデータを使い、 公的統計の実値に基づいて手法を可視化する。 教育用ハンズオン教材として scipy/pandas/sklearn の標準ライブラリで完結。
📥 入力例: data/raw/SSDSE-B-2026.csv encoding=shift_jis(または cp932), skiprows=1 47 都道府県 × 5 年(2019-2023)= 235 行 数値特徴量 100+ 列(家計・人口・教育・医療)
📤 実行例: random effects の主要出力を確認 都道府県別・年度別の指標が出力される 例: 秋田・沖縄・東京の差異が顕著に表れる プロット・統計量・モデルパラメータ等の結果を取得
💬 読み方: random effects の結果は SSDSE-B-2026 の地域特性を反映。 都市圏/地方圏/観光圏でパターンが異なり、 47 県の異質性を視覚化・定量化できる。 統計データ活用コンペの実分析に直接転用可能。
 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
import pymc as pm
region_idx, regions = pd.factorize(df['region'])
with pm.Model() as model:
    mu_a    = pm.Normal('mu_a', 0, 10)
    sigma_a = pm.HalfNormal('sigma_a', 5)
    a       = pm.Normal('a', mu_a, sigma_a, shape=len(regions))
    b1      = pm.Normal('b1', 0, 5)
    b2      = pm.Normal('b2', 0, 5)
    sigma   = pm.HalfNormal('sigma', 5)
    mu      = a[region_idx] + b1*df['x1'] + b2*df['x2']
    y_obs   = pm.Normal('y_obs', mu, sigma, observed=df['y'])
    trace   = pm.sample(2000, tune=1000, chains=4, target_accept=0.95)

4. scipy で ICC を素朴に計算する

▼ コード解説(random effects の Python 実装(コードブロック 15))
🎯 解説: random effects を SSDSE-B-2026 都道府県データで実行する Python コード。 47 都道府県 × 112 指標のパネルデータを使い、 公的統計の実値に基づいて手法を可視化する。 教育用ハンズオン教材として scipy/pandas/sklearn の標準ライブラリで完結。
📥 入力例: data/raw/SSDSE-B-2026.csv encoding=shift_jis(または cp932), skiprows=1 47 都道府県 × 5 年(2019-2023)= 235 行 数値特徴量 100+ 列(家計・人口・教育・医療)
📤 実行例: random effects の主要出力を確認 都道府県別・年度別の指標が出力される 例: 秋田・沖縄・東京の差異が顕著に表れる プロット・統計量・モデルパラメータ等の結果を取得
💬 読み方: random effects の結果は SSDSE-B-2026 の地域特性を反映。 都市圏/地方圏/観光圏でパターンが異なり、 47 県の異質性を視覚化・定量化できる。 統計データ活用コンペの実分析に直接転用可能。
▼ コード解説(random effects の Python 実装(コードブロック 16))
🎯 解説: random effects を SSDSE-B-2026 都道府県データで実行する Python コード。 47 都道府県 × 112 指標のパネルデータを使い、 公的統計の実値に基づいて手法を可視化する。 教育用ハンズオン教材として scipy/pandas/sklearn の標準ライブラリで完結。
📥 入力例: data/raw/SSDSE-B-2026.csv encoding=shift_jis(または cp932), skiprows=1 47 都道府県 × 5 年(2019-2023)= 235 行 数値特徴量 100+ 列(家計・人口・教育・医療)
📤 実行例: random effects の主要出力を確認 都道府県別・年度別の指標が出力される 例: 秋田・沖縄・東京の差異が顕著に表れる プロット・統計量・モデルパラメータ等の結果を取得
💬 読み方: random effects の結果は SSDSE-B-2026 の地域特性を反映。 都市圏/地方圏/観光圏でパターンが異なり、 47 県の異質性を視覚化・定量化できる。 統計データ活用コンペの実分析に直接転用可能。
1
2
3
4
5
6
7
8
from scipy import stats
groups = [df.loc[df.region==r, 'y'].values for r in df['region'].unique()]
F, p = stats.f_oneway(*groups)
n = df.groupby('region').size().mean()
msb = np.var([g.mean() for g in groups], ddof=1) * n
msw = np.mean([np.var(g, ddof=1) for g in groups])
icc = (msb - msw) / (msb + (n-1)*msw)
print(f'ICC(1,1) = {icc:.3f}')