パーセンタイル的なあれを数値的に計算する2

前回、回帰もできるし罰則もつけられるよねー

ってことで試した。

tau適当にいじってもそれらしい結果が出てきた。

全く検証してないけど。

多分、定数項には正則化項をつけないのが正しい気がする。

n<-50
a<-3
b<-2.5

x <- runif(n)
y <- a+ b*x + rnorm(n,mean=0,sd=b*x)

dmat <- cbind(rep(1,n),x)

#初期値をOLS解にする
lres<-lm(y~x,data=data.frame(x=x,y=y))
beta<-lres$coefficients

#%点
tau<- 0.9
#正則化項　15くらいまでにしないと発散して解が出ない
lambda<-0.01

#学習率とか計算の回数
learning_rate <- 0.1
itermax<-1000

#1000くらいじゃ収束しない時がある
for(i in 1:itermax){
  #ざんさ
  resd<- y - dmat %*% beta
  #一回微分
  temp<- tau-(resd<0)
  beta_new <- colSums(temp[,1] * (-dmat))/n + lambda * beta
  #目的関数(いらないけど)
  #objF<- sum((tau-(resd<0))*resd)/n
  
  #更新
  beta <- beta - learning_rate*(beta_new)
}

#あってるかどうか見る
plot(x,y,xlim=c(0,1),ylim=c(0,10))
par(new=T)
plot(x,dmat %*% beta,xlim=c(0,1),ylim=c(0,10))