2014-05-13

Quantile Regression

Wikiはこちら
日本語記事で扱っているのはこちら
観察データ $Y=(y_1,...,y_n);y_1 \le y_2 \le ... \le y_n$ があったとする
クオンタイルの回帰(直)線を引く前に、クオンタイル値の推定について考える
今、 $p$ -クオンタイル値を $z(q)$ としてうまいこと求めたい
$z(q)$ が $p$ -クオンタイル値であるとき、観測値が $z(q)$ より小さければ、その差の大きさを $1-p$ で重みづけして足し合わせ、 $z(q)$ より大きければ、その差の大きさを $p$ で重みづけして足し合わせ、その合算を取ることにする。これの最小化をすることで $z(q)$ を求めましょう、という方法をとる
これは、 $z(q)$ の値が『よい値』であると判断するときに、 $|y_i-z(q)|^2$ を用いて定める関数の最小化をするのではなくて、 $|y_i-z(q)|$ を用いて定める関数の最小化をします、という話
実際、を使った関数はどういうものかと言うと次のようになるのだが、
- $(1-p)\sum_{y_i < w} (w-y_i) + p\sum_{y_i > w} (y_i-w)$
どうして、それがよいかを説明しているのが、WikiのページのQuantilesの項(こちら)
- 上の式を満足するようなwというのは、クオンタイル値の定義式 $F_Y(q_{\tau}) = \tau$ に落ち着く、と書いてある
Rでやってみよう。
- 少数の値が得られているときに、上記の関数がどういう挙動をするかを、0,0.01,0.02,...,1の分位点を与える値を、yの最小値から最大値を100等分した値に設定したときに、上記関数の値がどうなるかをしらみつぶしに計算する
その値が0,0.01,...,1に沿って「谷底」になっているところを結ぶのが、上記手法での「クオンタイル値」の求め方

# クオンタイル回帰のためのパッケージ
library(quantreg)
# ちょっとしたデータを作る
y <- c(1,3,4,6,7,8,10)
# クオンタイル0,0.01,...,1
ps <- seq(from=0,to=1,length=101)
# yの値の最小値から最大値までの等分点
qs <- seq(from=min(y),to=max(y),length=101)

# i行目列の値は、ps[i]に相当するクオンタイル値をqs[j]にしたときの『関数』の値
v <- matrix(0,length(ps),length(qs))
for(i in 1:length(ps)){
	for(j in 1:length(qs)){
		s <- which(y < qs[j])
		t <- which(y > qs[j])
		v[i,j] <- sum((ps[i]-1)*(y[s]-qs[j])) + sum(ps[i]*(y[t]-qs[j]))
	}
}
# 赤の強い谷を左下方右上にたどれば、それが、この手法のクオンタイル値のコース
image(ps,qs,v)
# そこに、いわゆる、標本クオンタイル関数の結果をプロットしてみる
q.1 <- quantile(y,p)
q.2 <- kuantile(y,p,type=7) # type=7がデフォルト

points(p,q.1,col=3,pch=20,type="b")
points(p,q.2,col=4,pch=20,type="b")

クオンタイル値の推定値を求める背景の確認が終わったので、それを用いた回帰(直)線を引きたい
横軸xの値に応じて、縦軸yの値が広がりをもって観察されているときに、xの値が与えられたときに、そのxの値に対してyの値のクオンタイル値が知りたいとする。
そのときにxの値に対して各クオンタイル値が回帰(直)線になるようにして見ましょう、という話
引くべき線は直線、その評価は $|y_i-(a x_i +b)|$ のように「推定値からの距離(二乗ではなく)」を基準にします(がそれ以外はいわゆる回帰)

n <- 1000
x <- runif(n )
y <- x
for(i in 1:n){
y[i] <- y[i] + rnorm(1,0,x[i]^2)
}

plot(x,y)

ps <- c(0,0.05,0.1,0.25,0.5,0.75,0.9,0.95,1)

for(i in 1:length(ps)){
	# ps[i]クオンタイルの線形回帰線を引く
	abline(rq(y~x,ps[i]), col=i+1)
}
# 線形回帰の枠組みだから、切片なしの回帰直線もlm()と同様のやり方でできる
for(i in 1:length(ps)){
	# ps[i]クオンタイルの線形回帰線を引く、ただし切片は0
	abline(rq(y~x-1,ps[i]), col=i)
}