読者です 読者をやめる 読者になる 読者になる

野次馬エンジニア道

野次馬な気持ちでプログラミングをあれこれと綴ります

重回帰分析 - 説明変数の選択

単回帰に続いて今日は重回帰分析

重回帰分析

重回帰モデル

{ \displaystyle
 y_i = \beta_0 + \beta_1x_{i1} + \beta_2x_{i2} + \cdots + \beta_px_{ip} +\varepsilon_i \quad (i=1,\cdots,n) \quad \varepsilon_i \sim N(0, \sigma^2)
}

残差と平方和は、

{ \displaystyle
e_i  = y_i - \hat{y}_i = y_i - (\beta_0 + \beta_1x_{i1} + \beta_2x_{i2} + \cdots + \beta_px_{ip})
}

{ \displaystyle
S_e  =\sum_{i=1}^n e_i^2  = \sum_{i=1}^n \left\{  y_i - (\beta_0 + \beta_1x_{i1} + \beta_2x_{i2} + \cdots + \beta_px_{ip}) \right\} ^2
}

となる。

偏回帰係数の推定

{ S_e} を最小にするような偏回帰係数 { \hat{\beta_0},\hat{\beta_1}, \cdots, \hat{\beta_p}} を求める。

{ \hat{\beta_0},\hat{\beta_1}, \cdots, \hat{\beta_p}} のそれぞれで偏微分して0にすると

{ \displaystyle
 \overline{y}_i = \hat{\beta}_0 + \hat{\beta}_1\overline{x}_{1} + \hat{\beta}_2\overline{x}_{2} + \cdots + \hat{\beta}_p\overline{x}_{p}
}

{ \displaystyle
\hat{\beta}_1S_{11} + \hat{\beta}_2S_{12} + \cdots + \hat{\beta}_pS_{1p} = S_{1y}
}

{ \displaystyle
\hat{\beta}_1S_{21} + \hat{\beta}_2S_{22} + \cdots + \hat{\beta}_pS_{2p} = S_{2y} \\
\vdots \quad \quad \quad \quad \vdots \quad \quad \quad \quad \vdots
}

{ \displaystyle
\hat{\beta}_1S_{p1} + \hat{\beta}_2S_{p2} + \cdots + \hat{\beta}_pS_{pp} = S_{py} \\
}

となる。

ここで説明変数間の共分散 {S_{jk}} と 説明変数と応答変数の共分散を { S_{jy}}

{ \displaystyle
 S_{jk} = S_{kj} = \sum_{i=1}^{n}(x_{ij}-\overline{x}_j)(x_{ik}-\overline{x}_k)
}

{ \displaystyle
S_{jy} = \sum_{i=1}^{n}(x_{ij}-\overline{x}_j)(y_{j}-\overline{y})
}

とすると最小二乗法の解は

 {
\begin{bmatrix}
\hat{\beta}_1 \\
\hat{\beta}_2 \\
\vdots \\
\hat{\beta}_p \\
\end{bmatrix}
 = 
\begin{bmatrix}
S_{11} & S_{12} & \cdots & S_{1p} \\
S_{21} & S_{22} & \cdots & S_{2p} \\
\vdots   & \vdots   & \ddots   & \vdots \\ 
S_{p1} & S_{p2} & \cdots & S_{pp} 
\end{bmatrix} ^{-1}

\begin{bmatrix}
S_{1y} \\
S_{2y} \\
\vdots \\
S_{py} \\
\end{bmatrix}
}

で求まる。

多重共線性

上記の逆行列が求まらない場合、行列式が0に近い場合は説明変数間に相関がある可能性がある(多重共線性)。 その場合は相関関係にある変数を外して解析をやり直す必要がある。

回帰の評価

平方和の分解と決定係数

モデルが応答変数をどの程度説明しているかを決定係数で見る。

決定係数は、y の変動(平方和)を

{ \displaystyle
 \sum_{i=1}^n (y_i-\overline{y})^2 = \sum_{i=1}^n (\hat{y_i}-\overline{y})^2 + \sum_{i=1}^n e_i^2
}

{ \displaystyle
S_{yy} = S_R + S_e
}

データの変動(全平方和) = 予測値の変動(回帰平方和) + 残差変動(残差平方和)

に分解し、{ \sum_{i=1}^n (y_i-\overline{y})^2} で両辺を割って基準化することによって導出されて

決定係数 = 予測値の変動(回帰平方和) / データの変動(全平方和)

{ \displaystyle
R^2 = \frac{S_R}{S_{yy}}
}

で求められる。単回帰の場合さらに決定係数は相関係数の二乗と等しかった。

分散分析表と自由度調整済み決定係数

さらに、データの自由度 { n - 1  = } 予測モデルの自由度 + 残差の自由度 として、各自由度で変動(平方和)を割って 分散(平均平方和)に分解する。

これを表にまとめると

変動要因 変動(平方和) 自由度 分散(平均平方和)
モデル {\sum_{i=1}^n (\hat{y_i}-\overline{y})^2} p { \sum_{i=1}^n (\hat{y_i}-\overline{y})^2/p }
残差 {\sum_{i=1}^n e_i^2 } (n-1) - p { \sum_{i=1}^n e_i^2 /\{(n-1)-p\}}
合計 {\sum_{i=1}^n (y_i-\overline{y})^2} n-1 { \sum_{i=1}^n (y_i-\overline{y})^2/ (n-1) }

このとき、自由度調整済み決定係数は、残差の分散を用いて

{ \displaystyle
R^{*2} = 1 - \frac{ \frac{S_{e}}{n-p-1} }{ \frac{S_{yy}}{n-1} }
}

となる。さらに決定係数と比較して

{ \displaystyle
R^2 > R^{*2}
}

となる。

説明変数の選択

自由度調整済み決定係数は説明変数の組みがモデルをどの程度説明するかを示している。 実際に応答変数に寄与している説明変数だけをモデルに含めるべき。

モデルと残差の変動は、

{ \displaystyle
 \sum_{i=1}^n (\hat{y_i}-\overline{y})^2 \sim \chi^2_(p)
}

{ \displaystyle
 \sum_{i=1}^n e_i^2 / \sigma^2 \sim \chi^2_(n-1-p)
}

となるので、その分散の比は

{ \displaystyle
  F = \frac{ \frac{S_{R}}{p} }{ \frac{ S_{e}}{n-1-p} } \sim F_{p,n-1-p}
}

となる。

決定係数、自由度調整済み決定係数、このF値の基準を元にモデル(説明変数の組み)を選択していく。