Report10 (データ解析)

情報工学実験 II
データ解析 直線の当てはめ
055706E 上里盛真
1
1
線形回帰・最小二乗法
目的変数と説明変数の間錦を当てはめ、目的変数が説明変数によってどれくらい説明できるかを定量的に分析する
ことを、回帰分析といい、一次式を用いる回帰分析を特に線形回帰と言う。
最小二乗法とは、回帰分析で用いられる推計方法の一つで、目的変数の測定値と説明変数の測定値、および回帰式
を用いて求めた目的変数の推定値の差の二乗の平均が最小になるように探索する方法である。誤差の和は少ない方が
好ましいが、誤差には正負があるのでこれを二乗し解決している。
2
カイ二乗関数
回帰分析では式の信頼性を保つために信頼区間 (推定誤差範囲) というものが使われる。そのためカイ二乗関数と呼
∑N
i 2
ばれる式が用いられている。線形回帰におけるカイ二乗関数の式は x2 (aj ) = i=1 ( ∆
σi ) となる。
3
線形回帰当てはめプログラム
線形回帰なので回帰式は
Y = a1 + a2 x
(1)
これを最小二乗法にかけるのでカイ二乗関数が最小になる a1 , a2 を求める
χ2 (aj ) =
N
∑
1
2
2 (a1 + a2 xi − yi )
σ
i=1 i
(2)
最小を求めるには上式を微分し値が 0 になる点を探索する
N
∑
∂χ2
1
=2
2 (a1 + a2 xi − yi ) = 0
∂a1
σ
i=1 i
(3)
N
∑
∂χ2
1
=2
2 (a1 + a2 xi − yi )xi = 0
∂a2
σ
i=1 i
(4)
つまり下式を求める
これを整理して
∑
∑
a1 S + a2
x−
y=0
∑
∑
∑
a1
x + a2
x2 −
xy = 0
∑ ∑ 2 ∑ ∑
y x − x xy
∑
∑
S x2 − ( x)2
∑
∑ ∑
S xy − x y
∑
∑
a2 =
S x2 − ( x)2
a1 =
さらに誤差範囲は
√
σa1 =
(6)
(7)
(8)
∑
S
∑
√
σa2 =
(5)
S
∑
この辺は講義資料に詳しく載っている。
2
x2
∑
− ( x)2
(9)
S
∑
x2 − ( x)2
(10)
x2
3.1
lin-reg.cpp
一部抜粋
for( i=1; i<=nData; i++ ) {
sigmaTerm = 1/(sigma(i)*sigma(i));
/*前準備
式 (5.10) の計算を行う
*/
s = s + sigmaTerm ;
sx = sx + x(i)*sigmaTerm ;
sy = sy + y(i)*sigmaTerm ;
sxy = sxy + ((x(i)*y(i))*sigmaTerm) ;
sxx = sxx + ((x(i)*x(i))*sigmaTerm) ;
}
//* 切片 a_fit(1) と 傾き a_fit(2) の計算
a_fit(1) = ((sy*sxx) - (sx*sxy)) / ((s*sxx) - (sx*sx));
a_fit(2) = ((s*sxy) - (sy*sx)) / ((s*sxx) - (sx*sx));
//* 切片と傾きの誤差範囲を計算
sig_a(1) = sqrt((sxx)/((s*sxx) - (sx*sx)));
sig_a(2) = sqrt((s)/((s*sxx) - (sx*sx)));
//* 推定値とカイ二乗関数の計算
chisqr = 0.0;
sigmaTerm = 0.0;
for( i=1; i<=nData; i++ ) {
sigmaTerm = sigma(i)*sigma(i);
//推定値
yy(i) = a_fit(1) + a_fit(2)*x(i);
//カイ2乗関数
chisqr = chisqr + ((yy(i)-y(i)) * (yy(i)-y(i))) / sigmaTerm;
}
3.2
実行結果
きちんと回帰式を計算できていることが分かる。
3
図 1: y = 6x2 + 4x + 2 , 信頼区間 0.1
図 2: y = 6x2 + 4x + 2 , 信頼区間 10
図 3: y = 8x2 + 9x + 10 , 信頼区間 0.1
4
4
参考資料
• wikipedia
http://ja.wikipedia.org/wiki/メインページ
• 講義資料ページ
http://jura.nc.ie.u-ryukyu.ac.jp/~tom/jikken/J2-txt-2004/J2-txt-2004.html
5