読者です 読者をやめる 読者になる 読者になる

R 分布の知識整理(2) t分布1/歪度(わいど)と尖度(せんど)

t分布と正規分布の関係、歪度(わいど)と尖度(せんど)について確認する。説明はここがわかりやすい。t分布とは、

n個のデータx_{1},x_{2},\dots,x_{n}が独立にN(\mu,\sigma^{2})に従うとき、

t=\frac{\overline{x}-\mu}{\sqrt{\frac{\hat{\sigma}^{2}}{n}}

ただし、\overline{x}は標本平均、\hat{\sigma}^{2}は不偏分散

に従う分布である。上記の場合は自由度がn-1

 nの値を増加させるにつれ、正規分布に近づくことを確認する。また、次の式で定義される歪度(わいど)と尖度(せんど)の動きも追ってみる。

歪度は、\frac{\sum^n_{i=1}(x_i-\bar{x})^3}{{s_{xx}}^{\frac{3}{2}}(n-1)}

尖度は、\frac{\sum^n_{i=1}(x_i-\bar{x})^4}{{s_{xx}}^{2}(n-1)}

以下Rのコード。

maxdf <- 100
y <- rep(0,maxdf)
w <- rep(0,maxdf)
s <- rep(0,maxdf)

for(i in 1:maxdf) {
  x <- rt(df=i,10000)
  y[i] <- ks.test(x,"pnorm",mean=mean(x),sd=sd(x))$p.value
  w[i] <- sum({x-mean(x)}^3)/sd(x)^(3/2)/(length(x)-1)
  s[i] <- sum({x-mean(x)}^4)/sd(x)^2/(length(x)-1)
} #文字化け対策で一部()を{}にしています

par(mfrow=c(3,1),)
plot(y,type="l",ylab="p.value of ks.test",xlab="degrees of freedom")
abline(h=0.05,col=2)
abline(v=20,col=2)

plot(w,type="l",ylab="skewness",xlab="degrees of freedom of t-distribution",ylim=c(0,10))
abline(h=0,col=2)
abline(v=20,col=2)

plot(s,type="l",ylab="kurtosis",xlab="degrees of freedom of t-distribution",ylim=c(0,10))
abline(h=3,col=2)
abline(v=20,col=2)
par(mfrow=c(1,1))

結果は以下の通り。

上:自由度が20を超えたあたりから、ks.test(p値が高い程、正規分布に近いと推定出来る)のp値がほぼ完全に0.05より大きくなり、t分布は正規分布に近いと推定が出来る。

中:歪度(わいど)は、正規分布と同様0に収束する

下:尖度(せんど)は、正規分布と同様の3に収束する

f:id:nakhirot:20130918205935p:plain

 【疑問】

以前、時系列の状態空間モデルのノイズ項に正規分布の代わりにt分布を用いるモデルを見たことがある。時系列のモデル構築に用いるデータが20足らずの場合、ノイズはt分布を使用するとrobustになる、ということなのだろうか。後日考えることにする。