誤差と残差の違いから、回帰係数の検定まで、流れを振り返ってみます。
残差は誤差に冪等行列をかけたもので表され、残差平方和は誤差の2次形式で表されます。
冪等行列が理解できたときは感動しました。
(参考資料) 京都産業大学 岩城 秀樹 様 『統計学のレビュー』 http://www.cc.kyoto-su.ac.jp/~iwaki/Stat0.pdf
慶應大学 麻生 良文 様 『重回帰分析』 http://fs1.law.keio.ac.jp/~aso/ecnm/pp/reg2.pdf
統計検定1級過去問 2014人文2
重回帰分析
データはn 個、説明変数は2個とする。(説明変数がp個の場合も同様に考える)
Y=β1 1+β2X2+β3X3+u を成分で表現すると、( u∼N(0,I3) は誤差項 )
(Y1⋮Yn)=(1x21x31⋮⋮⋮1x2nx3n) (β1β2β3)+(u1⋮un) 行列表現では Y=Xβ+u
X=(1X2X3)=(1x21x31⋮⋮⋮1x2nx3n) rank(X)=3 とする。 (full rank)
ただし、1=(1⋮1),X2=(x21⋮x2n),X3=(x31⋮x3n),Y=(Y1⋮Yn)
X′X は3次正方行列で正則とする。 rank(X′X)=3
[参考]
rank(A) = min(m, n) のとき、A は最大階数(full rank)を持つといいフルランク行列という。
回帰係数の推定
Q=∑iu2i=∑i{Yi−(β1⋅1+β2x2i+β3x3i)}2
平方和は行列の積で表すことができる。
Q=(Y−Xβ)′(Y−Xβ)=(Y′−β′X′)(Y−Xβ)=Y′Y−Y′Xβ−β′X′Y+β′X′Xβ
Y′Xβ, β′X′Y はスカラーだから、Y′Xβ=β′X′Y
Q=(Y−Xβ)′(Y−Xβ)=Y′Y−2β′X′Y+β′X′Xβ
行列の微分では、 ∂∂xxTAx=(A+AT)x Aが対称行列の場合は、∂∂xxTAx=2Ax
X′X は対称行列なので、∂Qdβ=−2X′Y+2X′Xβ
ここで、 ∂Qdβ=0 とすると、 X′Xβ=X′Y (正規方程式)
Q を最小にするβ を ˆβ とすると、 ˆβ=(X′X)−1X′Y を得る。
残差と誤差
β=ˆβ のとき、 e=Y−Xˆβ を残差という。
誤差と残差の違いは、
誤差は、真の β による値と実測値の差、 u=Y−Xβ
残差は、推定値 ˆβ による値と実測値の差、e=Y−Xˆβ
残差の性質としては、 ˆY=Xˆβ とすると、
X′(Y−Xˆβ)=0 より、 X′e=0 ˆY′e=ˆβ′X′e=0
ˆY′e=(ˆY,e) は内積を意味するから、 ˆY と残差 e は無相関。
同様に、 X と残差 e は無相関。よって、
残差 e は、 1=(1⋮1),X2=(x21⋮x2n),X3=(x31⋮x3n) とも無相関となる。
即ち、1′e=0 , X′2e=0 , X′3e=0
[注意]
誤差の平方和 u′u はn個の正規分布の平方和であるから、 u′uσ2∼χ2(n)
残差の平方和はどのような分布に従うかが、次に考えるべき問題になる。
(全変動)=(回帰変動)+(残差変動)
予測値 ˆY=Xˆβ
実測値 Y=Xˆβ+e ∴Y=ˆY+e
残差 e=Y−Xˆβ
残差平方和 e′e=(Y−Xˆβ)′(Y−Xˆβ)
全変動
(Y−¯Y⋅1)′(Y−¯Y⋅1)=(ˆY+e−¯Y⋅1)′(ˆY+e−¯Y⋅1)=(ˆY−¯Y⋅1+e)′(ˆY−¯Y⋅1+e)=(ˆY−¯Y⋅1)′(ˆY−¯Y⋅1)+e′e+e′(ˆY−¯Y⋅1)+(ˆY−¯Y⋅1)′e
ここで、 X′Xˆβ=X′Y より、 X′(Y−Xˆβ)=0 ∴X′e=0 また、1′e=0 だから、
(ˆY−¯Y⋅1)′e=ˆY′e−¯Y⋅1′e=ˆY′e−¯Y⋅0=(Xˆβ)′e=ˆβ′X′e=0
同様に、 e′(ˆY−¯Y⋅1)=0 よって、
(Y−¯Y⋅1)′(Y−¯Y⋅1)=(ˆY−¯Y⋅1)′(ˆY−¯Y⋅1)+e′e
これは、(全変動)=(回帰変動)+(残差変動) を意味する。
冪等行列
I=In (n次正方行列)とする。
e=Y−Xˆβ=Y−X(X′X)−1X′Y=[I−X(X′X)−1X′]Y
[I−X(X′X)−1X′]Y=[I−X(X′X)−1X′](Xβ+u)=[I−X(X′X)−1X′](Xβ)+[I−X(X′X)−1X′]u=Xβ−Xβ+[I−X(X′X)−1X′]u=[I−X(X′X)−1X′]u
∴e=[I−X(X′X)−1X′]u
と表される。
ここに登場する行列は興味深い性質を持つ。
・ X(X′X)−1X′ は対称な冪等行列で階数は3である。
・ I−X(X′X)−1X′ は対称な冪等行列で階数はn-3である。
冪等行列では、階数はトレースに等しい。
また、tr(AB)=tr(BA) 、tr(A+B)=tr(A)+tr(B) という性質を用いると、
tr(X(X′X)−1⋅X′)=tr(X′X(X′X)−1)=tr(I3)=3
tr(I−X(X′X)−1⋅X′)=n−3 これが階数になる。
残差平方和の分布
M=I−X(X′X)−1X′ とおけば、M′=M 、M2=M 対称冪等行列となる。(rank(M)=n−3 )
残差は、 e=Mu を満たすので、
e′e=u′M′Mu=u′MMu=u′M2u=u′Mu
よって、残差平方和 SSE は、
SSE=e′e=u′Mu と表される。 (残差は誤差の二次形式)
理論上、 ui∼N(0,σ2) で、各 ui どうしは独立である。
z=1σu とおくと、z∼N(0,In) z はn次元標準正規分布に従うことになり、
SSEσ2=e′eσ2=z′Mz と表される。
Mは対称行列だから直交行列Pで対角化でき、Mは冪等行列でもあるので固有値は 0 か 1。
固有値1に対する固有ベクトルを左から並べてPを構成すれば、P−1MP=(In−3000)
e′eσ2=z′Mz=z′P(P−1MP)P−1z において v=P−1z とすると、
e′eσ2=v′(In−3000)v
ここで、 z∼N(0,In) であり、 E(vv′)=P−1E(zz′)P=P−1InP=In だから、 v∼N(0,In)
よって、
e′eσ2=v21+v22+⋯+v2n−3 vi∼N(0,1) は互いに独立
独立な (n−3) 個の標準正規分布の平方和で表されたので、
SSEσ2∼χ2(n−3) が成立する。
[参考]
対称な冪等行列の対角化については、
甲南大学 村澤 康友 様 『対称行列の固有値分解」
https://ymurasawa.web.fc2.com/ue-ln10.pdf を参照
回帰係数の分布
回帰係数の推定値は、ˆβ=(X′X)−1X′Y
ˆβ−β=(X′X)−1X′Y−(X′X)−1 (X′X)β=(X′X)−1X′(Y−Xβ)=(X′X)−1X′u
ˆβ=β+(X′X)−1X′u より、 ˆβ は u1,u2,⋯,un の線型結合。
よって、 ˆβ はn次元正規分布に従う。
Cov(ˆβ)=E((ˆβ−β) (ˆβ−β)′)=(X′X)−1X′E(uu′)X(X′X)−1=(X′X)−1X′σ2IX(X′X)−1=σ2(X′X)−1
ˆβ∼N(β,σ2(X′X)−1) ˆβ−β∼N(0 , σ2(X′X)−1)
ˆβ の成分 ˆβi については、
(X′X)−1 のij 成分を qij とすると、
σ2(X′X)−1 のii 成分は σ2qii だから、 ˆβi∼N(βi,σ2qii)
よって、 Z=ˆβi−βiσ√qii は、 N(0,1) に従う。
回帰係数のt検定
σ2 が未知のときは、OLS推定量 s2=1n−3e′e を用いて t 検定量を構成する。
定理
Z,Y が独立で Z∼N(0,1) 、 Y∼χ2(m) のとき、 T=Z√Y/m は、自由度mのt分布に従う。
Zについては、上述のものを用い、
Z=ˆβi−βiσ√qii∼N(0,1)
Yについては、カイ二乗分布に従う統計量として、
Y=SSEσ2∼χ2(n−3)
このとき、 Y/m に対応する部分は、
Y/m=SSEσ2(n−3)=s2σ2
よって、 T=ˆβi−βiσ√qii÷sσ=ˆβi−βis√qii∼t(n−3)
Z,Y の独立性をいうには、 ˆβ−β と e′e の独立をいえばよいが難しい。
u の二次形式 e′e=u′Mu
uの一次変換 ˆβ−β=(X′X)−1X′u
慶應大学 長倉先生の資料p5「正規分布に従う確率変数の一次関数と二次形式の独立」によれば
補題 x∼N(μ,σ2In) 、Bは(q,n)行列、Aは(n,n)対称行列とするとき、
BA=0 ならば、Bx と x′Ax は独立であるという。
今、B=(X′X)−1X′ 、A=I−X(X′X)−1X′ であり、
BA=(X′X)−1X′(I−X(X′X)−1X′)=(X′X)−1X′−(X′X)−1X′X(X′X)−1X′=0
よって、 ˆβ−β と e′e の独立はいえる。補題の証明については、ここでは立ち入らない。
[参考]
京都産業大学 岩城秀樹 先生 「統計学のレビュー」のp10に
『 ˆβ と e は独立 』であることの説明がある。
慶應大学 麻生良文 先生の解説がもっとも詳しいが、独立の証明は省かれている。
やはり、線形代数の知識が必要ということである。
http://fs1.law.keio.ac.jp/~aso/ecnm/pp/reg2.pdf
回帰係数のF検定
慶應大学 長倉大輔 様
http://user.keio.ac.jp/~nagakura/zemi/F_test.pdf による。
Z=ˆβi−βiσ√qii∼N(0,1) 、 Y=SSEσ2∼χ2(n−3) のとき、
t検定量を T=Z√Y/(n−3) とすることで、 ˆβi を検定することができた。
この検定量を拡張し、いくつかの係数を同時に検定できるようにしよう。
T2=Z2Y/(n−3) とおくと、Z2∼χ2(1) 、Y=SSEσ2∼χ2(n−3)
独立なカイ二乗分布の商は、F分布に従うので、
Z2∼χ2(1) 、Y=SSEσ2∼χ2(n−3) のとき、
F検定量 F=Z2/1Y/(n−3)∼F(1,n−3) 自由度 (1,n−3) のF分布に従う。
つまり、(F検定量)=(t検定量の2乗) であり、これを用いて回帰係数 ˆβi を検定することができる。
これを一般の ˆβ の検定に拡張できないだろうか?
分子が係数の差の成分の平方であるから、係数の差の二次形式を考える。
長倉先生の「t 統計量とF統計量の分布」ノート版 p4 に的確な解説がある。
http://user.keio.ac.jp/~nagakura/zemi/t_and_F_test.pdf
『正規分布の分散共分散行列の逆行列による2次形式の分布』
xはq行1列、x∼N(0,Σ) 、rank(Σ)=q とするとき、 x′Σ−1x∼χ2(q)
これを用いると、
ˆβ−β∼N(0 ,σ2(X′X)−1) 、rank(σ2(X′X)−1)=3 、ˆβ−β=(X′X)−1X′u だから、
W=(ˆβ−β)′σ−2(X′X)(ˆβ−β)=u′X(X′X)−1(X′X)(X′X)−1X′u=u′X(X′X)−1X′u
M=X(X′X)−1X′ は対称な冪等行列で、階数は3(冪等行列の項を参照)
よって、uの3個の平方和で表されるので、 W∼χ2(3) 自由度3のカイ二乗分布に従う。
もし、帰無仮説 H0:β=0 が真とすると、ˆβ∼N(0 ,σ2(X′X)−1) だから、
W=ˆβ′σ−2(X′X)ˆβ の値は小さくなる。
よって、 F=W/3Y/(n−3)∼F(3,n−3) における F 検定量の値は小さくなる。
反対に、F の値が大きければ、帰無仮説が偽である可能性が高い。
このように、 F の値によって 回帰係数 ˆβ の検定を行うことができる。
[注意]
F 値は、分散分析の場合とは意味合いが異なる。
分散分析では、 F 値は全変動の中での回帰変動の割合であるが、
回帰分析では、カイ二乗分布に従う独立な確率変数の商の分布を考えている。
F統計量と決定係数の関係については「統計学のレビュー」(岩城先生)p12 を参照
富山大学 唐渡広志 先生のサイトにも、F検定の解説がある。
http://www3.u-toyama.ac.jp/kkarato/2017/econometrics/handout/Econometrics-2017-14-1117.pdf
[更新日]