Linear Regression

Figure2.1 in The Elements of Statistical Learning を書いてみました。 前回同様Input SpaceはX = {\bf R}^{2}であることに注意!

library(ElemStatLearn)

x = mixture.example$x
g = mixture.example$y
xnew = mixture.example$xnew

X = cbind(matrix(1, dim(x)[1], 1), x)
Xnew = cbind(matrix(1, dim(xnew)[1], 1), xnew)
beta = solve(t(X) %*% X) %*% t(X) %*% g

yhat = Xnew %*% beta 
px1 = mixture.example$px1
px2 = mixture.example$px2
contour(px1, px2, matrix(yhat, length(px1), length(px2)), levels=0.5, axes=FALSE, xlab="", ylab="", labels="", main="Linear Regression of 0/1 Response")
points(x, col=ifelse(g==1, "coral", "cornflowerblue"))
points(expand.grid(px1, px2), pch=".", cex=1.2, col=ifelse(yhat>0.5, "coral", "cornflowerblue"))
box()

file:///Users/hiroyukitanaka/Dropbox/myR/plots/esl2.1.png

広告

k-Nearest Neighbor

The Elements of Statistical LearningのFigure2.3を書いてみました。 参考URL

ElemStatLearnっていうライブラリにデータがあるみたいです。感激しました。 上のサンプルではclassライブラリ?のknnに丸投げしてますが、勉強のために自分で書いてみました。

rm(list=ls())
library(ElemStatLearn)
euclidean_distance = function(x, y){
  return(sqrt(t(x-y) %*% (x-y)))
}
N = function(x, k, train){
  return(sort(apply(train, 1, function(di) (sqrt(t(di-x) %*% (di-x)))), index.return=TRUE)$ix[1:k])
}
knn = function(x, k, train, g){
  sum(g[N(x, k, train)]) / k
}
k = 15
x = mixture.example$x
g = mixture.example$y
xnew = mixture.example$xnew
px1 = mixture.example$px1
px2 = mixture.example$px2
px = as.matrix(expand.grid(px1, px2), length(px1), length(px2))
#よくわからなければ↓を実行してみると良い
#----------------------------------------
#sample = matrix(px[1,], ncol=1) #着目点
#points(t(sample)) #着目点をプロット(黒点)
#points(x[N(sample, k, x), ], col="red")#着目点のk個のNeighbor
#knn(sample, k, x, g)
#----------------------------------------
#predlabels = ifelse(apply(px, 1, knn, k=k, train=x, g=g)>0.5, 1, 0)
prob = apply(px, 1, knn, k=k, train=x, g=g)
#par(mar=rep(2:4)) 
contour(px1, px2, matrix(prob, length(px1), length(px2)), levels = 0.5, labels="", xlab="", ylab="", main="15-nearnest neighbor", axes=FALSE)
points(x, col=ifelse(g==1, "coral", "cornflowerblue"))
points(expand.grid(px1, px2), pch=".", cex=1.2, col=ifelse(prob > 0.5, "coral", "cornflowerblue"))
box()

結果: file:///Users/hiroyukitanaka/Dropbox/myR/plots/esl2.3.png

注意したいのは、これinput spaceが2次元だってことですね。 僕はinput spaceが1次元(x)でyを予測する、って考えていたので非常に混乱してました。 この問題はclassificationなのでトレーニングデータにおける、3次元目は0,1で与えられてます(これがg)。

つまりG = {0, 1}(当たり前だが…)。 そしてx \in R^{2}に対してY(x)を計算し、境界を0.5としてclassificationするというわけですね。 次はregressionしてみます。input spaceが2次元だと言うことがわかったので、すんなりいけそうです。(多分 あー研究が進まないからこういうおもしろいことをしてしまう( ´ゝ`)