ポアソン回帰の最尤推定と対数尤度計算
ポアソン回帰は、カウントデータ(非負整数値)を扱う回帰手法で、目的変数がポアソン分布に従うと仮定します。
ポアソン回帰モデルの構成
ポアソン回帰モデル:
$$Y_i \\sim \\text{Poisson}(\\mu_i)$$
$$\\log(\\mu_i) = \\beta_0 + \\beta_1 x_i$$
ここで:
- $Y_i$:カウント型の目的変数
- $\\mu_i = E[Y_i]$:期待値(平均パラメータ)
- $\\log(\\mu_i)$:対数リンク関数
Step 1: ポアソン分布の確率質量関数
ポアソン分布の確率質量関数:
$$P(Y_i = y_i | \\mu_i) = \\frac{\\mu_i^{y_i} e^{-\\mu_i}}{y_i!}$$
対数尤度への寄与:
$$\\log P(Y_i = y_i | \\mu_i) = y_i \\log \\mu_i - \\mu_i - \\log(y_i!)$$
Step 2: 期待値の計算
与えられたパラメータ:$\\beta_0 = 0$、$\\beta_1 = 1$
モデル:$\\log(\\mu_i) = 0 + 1 \\cdot x_i = x_i$
したがって:$\\mu_i = \\exp(x_i)$
各データ点の期待値:
- $\\mu_1 = \\exp(0) = 1$
- $\\mu_2 = \\exp(1) = e \\approx 2.718$
- $\\mu_3 = \\exp(2) = e^2 \\approx 7.389$
Step 3: 対数尤度の計算
全体の対数尤度:
$$\\log L = \\sum_{i=1}^3 [y_i \\log \\mu_i - \\mu_i - \\log(y_i!)]$$
データ点1: $(x_1, y_1) = (0, 1)$、$\\mu_1 = 1$
$$\\log L_1 = 1 \\cdot \\log(1) - 1 - \\log(1!) = 0 - 1 - 0 = -1$$
データ点2: $(x_2, y_2) = (1, 3)$、$\\mu_2 = e \\approx 2.718$
$$\\log L_2 = 3 \\cdot \\log(2.718) - 2.718 - \\log(3!)$$
$$= 3 \\times 1 - 2.718 - \\log(6) = 3 - 2.718 - 1.792 = -1.51$$
データ点3: $(x_3, y_3) = (2, 9)$、$\\mu_3 = e^2 \\approx 7.389$
$$\\log L_3 = 9 \\cdot \\log(7.389) - 7.389 - \\log(9!)$$
$$= 9 \\times 2 - 7.389 - \\log(362880) = 18 - 7.389 - 12.801 = -2.19$$
Step 4: 全体の対数尤度
$$\\log L = \\log L_1 + \\log L_2 + \\log L_3$$
$$= -1 + (-1.51) + (-2.19) = -4.70$$
ポアソン回帰の対数尤度関数
一般形:
$$\\log L(\\boldsymbol{\\beta}) = \\sum_{i=1}^n \\left[ y_i (\\mathbf{x}_i^T \\boldsymbol{\\beta}) - \\exp(\\mathbf{x}_i^T \\boldsymbol{\\beta}) - \\log(y_i!) \\right]$$
特徴:
- 線形性:$y_i$ に対して線形
- 指数項:$\\exp(\\mathbf{x}_i^T \\boldsymbol{\\beta})$ による非線形性
- 最適化:通常、反復重み付き最小二乗法(IRLS)で解く
Step 5: 最尤推定の実装手順
1. スコア関数(勾配)
$$\\frac{\\partial \\log L}{\\partial \\beta_j} = \\sum_{i=1}^n (y_i - \\mu_i) x_{ij}$$
2. ヘッセ行列(2次微分)
$$\\frac{\\partial^2 \\log L}{\\partial \\beta_j \\partial \\beta_k} = -\\sum_{i=1}^n \\mu_i x_{ij} x_{ik}$$
3. ニュートン・ラフソン法による更新
$$\\boldsymbol{\\beta}^{(t+1)} = \\boldsymbol{\\beta}^{(t)} - \\mathbf{H}^{-1} \\mathbf{s}$$
Step 6: モデル診断
ポアソン回帰の診断手法
- 過分散の検定:分散が平均より大きいかチェック
- ピアソン残差:$r_i^P = \\frac{y_i - \\hat{\\mu}_i}{\\sqrt{\\hat{\\mu}_i}}$
- 偏差残差:$r_i^D = \\text{sign}(y_i - \\hat{\\mu}_i) \\sqrt{2[y_i \\log(y_i/\\hat{\\mu}_i) - (y_i - \\hat{\\mu}_i)]}$
- 情報量基準:AIC, BIC によるモデル選択