第４回MIKU勉強会

線形回帰関数について

#lm()はＲの線形回帰関数
lm(formula, data, subset, weights, na.action,
   method = "qr", model = TRUE, x = FALSE, y = FALSE, qr = TRUE,
   singular.ok = TRUE, contrasts = NULL, offset, ...)

つまり、lm()関数を使うと言うのは、weight = a * height +b の式に乗せる、ということ。

# glm() はＲの一般化線形回帰関数
glm(formula　●~▲　●と▲の関数, family = gaussian, data, weights, subset,
    na.action, start = NULL, etastart, mustart, offset,
    control = list(...), model = TRUE, method = "glm.fit",
    x = FALSE, y = TRUE, contrasts = NULL, ...)

つまり、glm()関数を使うと言うのは、y=\sum_{i=0}^n a_i x^iにあてはめること。 

data(women)
x<-women$height
y<-women$weight

plot(x,y)

# Function to estimate polynomial coefficients and coordinates estimated
# 多項式近似の係数と、推定座標を返す関数
glmPolynomial<-function(y,x,dg){# y:dependent values (従属値列),x:descriptive values (説明変数),dg:degrees of polynomials (多項式の次数)
# Xs: Matrix consisted of x^1,x^2,...x^{dg}
# Xs: xの1,2,...,dg 乗の値でできた行列
Xs<-matrix(0,length(x),dg) 
for(i in 1:dg){
Xs[,i]<-x^i
}
# glm() is generalized linear regression function in R
# glm() はＲの一般化線形回帰関数
fit<-glm(y~Xs[,1:dg]) 
# beta: coefficients estimated
# beta: 推定係数ベクトル
beta<-fit$coefficients
# Xs2 has an additional column for x^0
# Xs2 には、x^0のための列を追加
Xs2<-cbind(rep(1,length(x)),Xs)
# predY is a matrix of values for individual degrees of x
# predY は xの各次数の項
predY<-t(t(Xs2)*beta)
# apply(): Sum up all columns corresponding to 0 to dg degrees
# apply()関数で0:dg次数のすべてを足し合わせる
sumupPredY<-apply(predY,1,sum)
list(beta=beta,x=x,predY=sumupPredY)
}

x<-women$height
y<-women$weight
xlim<-c(min(x),max(x))
ylim<-c(min(y),max(y))
plot(x,y,xlim=xlim,ylim=ylim,col="red")

dgs<-1:5
for(i in dgs){
outglm<-glmPolynomial(y,x,i)
par(new=TRUE)
# gray() changes contrast 
# gray() 関数は線の濃淡を変える
#plot(outglm$x,outglm$predY,xlim=xlim,ylim=ylim,col=gray(i/(max(dgs)*1.5)),type="l")
plot(outglm$x,outglm$predY,xlim=xlim,ylim=ylim,col=i,type="l")
}