Engineering Skills

製品開発エンジニアがデータ解析のノウハウを垂れ流します

てこ比

今回はてこ比について。各観測値が推定値に与える影響、および全体の平均からどの程度ずれているかを示す指標で、0 から 1 までの値をとります。

モデル式(最小二乗法)

まずはモデル式の多項式近似パラメータ推定について、最小二乗法で求めます。モデルは下記で、以下データの数を[math] \displaystyle n [/math]、近似パラメータの数を[math] \displaystyle p [/math]とします。

[math] \displaystyle y ={\beta}_0 + \sum_{i=1}^{n}{{\beta}_i x_i} [/math]

これを行列式で表します。

[math] \mathbf{y}=\left \{ \begin{matrix} y_1 \\ y_2 \\ \vdots \\ y_n \end{matrix} \right \} [/math]

[math] \mathbf{X}=\left \{ \begin{matrix} 1 & x_{11} & x_{12} & \cdots & x_{1p} \\ 1 & x_{21} & x_{22} & \cdots & x_{2p} \\ \vdots & \vdots & \vdots & \ddots & \vdots \\ 1 & x_{n1} & x_{n2} & \cdots & x_{np} \\ \end{matrix} \right \} [/math]

[math] \mathbf{\beta}=\left \{ \begin{matrix} \beta_1 \\ \beta_2 \\ \vdots \\ \beta_n \end{matrix} \right \} [/math]

[math] \mathbf{\varepsilon}=\left \{ \begin{matrix} \varepsilon_1 \\ \varepsilon_2 \\ \vdots \\ \varepsilon_n \end{matrix} \right \} [/math]

と置くと、近似モデルの行列式表示は下記のようになります。

[math] \displaystyle \mathbf{y}= \mathbf{X} \beta + \varepsilon [/math]

近似モデルからの推定値を[math] \displaystyle \mathbf{\hat{y}} [/math]と置くと、

[math] \displaystyle \mathbf{\hat{y}}= \mathbf{X} \beta [/math]

最小二乗法から[math] \displaystyle \mathbf{\beta} [/math]の推定値[math] \displaystyle \mathbf{b} [/math]が下記です。

[math] \displaystyle \mathbf{b}= (\mathbf{X}^{\mathit{T}} \mathbf{X})^{-1} \mathbf{X}^{\mathit{T}} \mathbf{y} [/math]

ハット行列

ここまでの結果から[math] \displaystyle \mathbf{\hat{y}} [/math][math] \displaystyle \mathbf{y} [/math]の関係は下記のようになります。

[math] \displaystyle \mathbf{\hat{y}}= \mathbf{X} (\mathbf{X}^{\mathit{T}} \mathbf{X})^{-1} \mathbf{X}^{\mathit{T}} \mathbf{y} = \mathbf{H}\mathbf{y} [/math]

[math] \displaystyle \mathbf{H} [/math] は応答変数 [math] \displaystyle \mathbf{y}[/math]から 回帰推定値 [math] \displaystyle \mathbf{\hat{y}}[/math] を生成する行列なので、ハット行列(hat matrix)と 呼ばれています。ハット行列は "射影行列" とも呼ばれます。

ハット行列は冪等(べきとう)行列でもあります。冪等行列とは、自分自身との積が自分自身に一致する行列のことです。。

[math] \displaystyle \mathbf{H}^{2} = \mathbf{X} (\mathbf{X}^{\mathit{T}} \mathbf{X})^{-1} \mathbf{X}^{\mathit{T}} \mathbf{X} (\mathbf{X}^{\mathit{T}} \mathbf{X})^{-1} \mathbf{X}^{\mathit{T}} = \mathbf{X} \mathbf{I} (\mathbf{X}^{\mathit{T}} \mathbf{X})^{-1} \mathbf{X}^{\mathit{T}} = \mathbf{X} (\mathbf{X}^{\mathit{T}} \mathbf{X})^{-1} \mathbf{X}^{\mathit{T}} = \mathbf{H} [/math]

つまり [math] \displaystyle \mathbf{H}^{2}=\mathbf{H} [/math]

てこ比(leverage)

ハット行列の対角成分 [math] \displaystyle h_{ii} [/math] が大きいと [math] \displaystyle i [/math]番目のデータは推定値に大きな 影響を及ぼします。これを「てこ比(leverage)」と呼びます。ここまでの結果からハット行列を用いて下記のように書けます。

[math] \displaystyle \mathbf{\hat{y}} = \mathbf{H}\mathbf{y} [/math]

上記から [math] \displaystyle \hat{ y_{i} } [/math]

[math] \displaystyle \hat{ y_{i} } = h_{i1} y_{1} + h_{i2} y_{2} + \cdots + h_{ij} y_{j} + h_{in} y_{n} [/math]

以上より、てこ比[math] \displaystyle h_{ii} [/math]は下記のようになり

[math] \displaystyle h_{ii} = \frac {\partial \hat{ y_{i} } } {\partial y_{i} } [/math]

ハット行列の対角成分になります。

各観測点毎に定義されるてこ比の平均値は回帰モデルのパラメータ数を[math] \displaystyle p [/math]、観測データ数を[math] \displaystyle n [/math]とすると、[math] \displaystyle \frac{(p + 1)}{n} [/math] なので、例えば[math] \displaystyle \frac{2(p + 1)}{n} [/math] より大きいデータには注意を払う必要があります。このデータが回帰係数の推定に大きな影響を与えているので、外れ値などの可能性があるからです。

てこ比の性質は、

・説明変数のデータを変えずに目的変数[math] \displaystyle \hat{ y_{i} } [/math]の値を1だけ変えたときの予測値[math] \displaystyle y_{i} [/math]の変化量

[math] \displaystyle 0 \leq h_{ii} \leq 1 [/math]

・合計は[math] \displaystyle (p + 1) [/math]

[math] \displaystyle h_{ii} = \frac{1}{n} + \frac{D_i}{n-1} [/math] (Dはマハラノビスの距離)

標準化誤差

標準化誤差の前に、回帰の残差ベクトル[math] \displaystyle \mathbf{\hat{e}} [/math]について考えてみると、次のようになります。

[math] \displaystyle \mathbf{\hat{e}} = \mathbf{y} - \mathbf{\hat{y}} = ( \mathbf{I} - \mathbf{H} ) \mathbf{y} [/math]

標準化残差は、残差([math] \displaystyle e_i [/math] )を誤差変動の平均平方の平方根(標準偏差の推定値)で割ったものになります。[math] \displaystyle e_i [/math]の分散を考えてみると、

[math] \displaystyle \mathbf{Var[\hat{e}]} = ( \mathbf{I} - \mathbf{H} ) {\mathbf{\sigma}}^2 [/math]

上記より

[math] \displaystyle \mathbf{\hat{e_i}} = N(0, {{\sigma}^2}(1-h_{ii}) ) [/math]
[math] \displaystyle \frac {\mathbf{\hat{e_{i}}}} {\sigma \sqrt{1-h_{ii}}} = N(0, 1 ) [/math]

残差分散[math] \displaystyle {\sigma}^2 [/math]の推定値を [math] \displaystyle {s_e}^2 [/math] とすると,

[math] \displaystyle \mathbf{r_{i}} = \frac {\mathbf{\hat{e_{i}}}} {{s_e} \sqrt{1-h_{ii}}} = \frac {\mathbf{\hat{e_{i}}}} {\sqrt{MSE(1-h_{ii})}} \sim t(n-p-1) [/math]

一つの基準として、標準化残差の絶対値が2より大きい場合は標準化残差が大きいと見なされます。

回帰診断

以上を踏まえた回帰分析の診断方法の一つに、てこ比ー標準化誤差をプロットすることがあります。

Y=Xにノイズを乗せた外れ値がなさそうな場合の一例は下記です。赤線は、横軸でてこ比が[math] \displaystyle \frac{2(p + 1)}{n} [/math] に縦線を、縦軸で標準化残差の絶対値が2の場所に横線を引いています。

f:id:OceanOne:20210625235401p:plainf:id:OceanOne:20210625235815p:plain

次に、明らかな外れ値を加えた一例を下記に示します。先ほど定義した赤線箇所から、外れ値らしきデータは逸脱していることが分かります。

f:id:OceanOne:20210625235504p:plainf:id:OceanOne:20210625235522p:plain

まとめ

回帰分析を診断するための、てこ比と標準化誤差の利用法/解釈について記載しました。今回も、こちらのツールに実装しています。