最小2乗法の復習
昨日のTokyo.Rで、判別分析とかロバスト推定とかの話を聞いて、久しぶりにその辺りの復習をしようと思い立ちました。GW中に、線形回帰からロジスティック回帰まで復習できればいいなぁと思っています。
さて、線形回帰を行うにあたって、重要な要素である、最小2乗法について復習します。
最小2乗法の理解を深めるために実際に手を動かしてみます。
最小2乗法については、最小二乗法についてなどを見て頂くとして、今回は
これなら分かる最適化数学―基礎原理から計算手法まで | |
金谷 健一 共立出版 2005-09 売り上げランキング : 40962 Amazonで詳しく見る by G-Tools |
まで落とし込んだ後、データが与えられて実際に直線を引くまでを行います。
ここまでの式変形が一つの肝だとは思いますが、今回はそれが分かってるというところからのスターですね。
線形の単回帰分析なので、そこまで難しくはありません。
y=ax+bの直線について最小2乗法適用
参考文献の例題を解いてみます。
鉱石の密度x(g/cm^3)と鉄含有量y(%)の間の関係
x | 2.8 | 2.9 | 3.0 | 3.1 | 3.2 | 3,2 | 3.2 | 3.3 | 3.4 |
---|---|---|---|---|---|---|---|---|---|
y | 30 | 26 | 33 | 31 | 33 | 35 | 37 | 36 | 33 |
このデータについて、y=ax+bの直線を求めてみます。
# XとYの値 Xs <- c(2.8, 2.9, 3.0, 3.1, 3.2, 3.2, 3.2, 3.3, 3.4) Ys <- c(30, 26, 33, 31, 33, 35, 37, 36, 33) # とりあえず点のplot plot(Xs, Ys) #X^2とxyを求める Xs2 <- sum(Xs^2) XY <- sum(Xs * Ys) # aとbについての連立一次方程式を解く est <- solve(matrix(c(Xs2, sum(Xs), sum(Xs), length(Xs)), 2, 2), matrix(c(XY, sum(Ys)), 2, 1)) # 与えられた点と一緒にプロット plot(Xs, Ys) abline(b=est[1,1], a=est[2,1], col="red") # R組み込みのlm関数の結果と比較 test.frame <- data.frame(cbind(Xs, Ys)) test.lm <- lm(Ys~Xs, data=test.frame) abline(test.lm, col="blue")
結果
(当たり前ですが、双方ほぼ同じ結果が得られました)
> est [,1] [1,] 12.067669 [2,] -5.011278 > > > test.lm Call: lm(formula = Ys ~ Xs, data = test.frame) Coefficients: (Intercept) Xs -5.011 12.068