ポアソン分布の机上統計実験


 最近、統計の本とか読むようになり、改めて感じますのが、自分に数学的センスが皆無だということ。それとそのセンスって努力云々だけじぢゃ埋められないものもあるということ。
でも、見たりすればなんとな~く理解することはできるので、どうやって自分に視覚的に理解させるのかは、重要なポイントだったりします。

 というわけで、kuboさん的統計実験をやってみますた。お題は最近よく使うポアソン分布の過分散。それをばらつかせたら、どうばらつくのかな~という実験。スクリプトは↓な感じです。

test.x <- seq(0,990,by=10)
test.x <- rep(test.x,10)
id <- 1:1000
id <- as.character(id)
katamuki <- NULL
seppenn <- NULL

for(i in 1:100){

#test.y <- exp(0.005*test.x+1)
test.y <- exp(0.005*test.x+rnorm(length(test.x),1,0.3))
random.y <- rpois(length(test.x),lambda=test.y)
par(new=T)
plot(test.x,random.y,xlim=c(0,1000),ylim=c(0,500),xlab=”標高”,ylab=”個体数密度”,type=”n”)
res.lm <- glm(random.y~test.x,family=poisson)
res.coef <- coef(res.lm)
res.coef <- as.vector(res.coef)
katamuki <- append(katamuki,res.coef[2])
seppenn <- append(seppenn,res.coef[1])

par(new=T)
line <- function(x) {exp(res.coef[2]*x+res.coef[1])}
x <- seq(0,1000,length=100)
plot(x,line(x),col=”red”,xlim=c(0,1000),ylim=c(0,500),type=”n”,xlab=””,ylab=””)
lines(x,line(x),lwd=0.001,col=”blue”)

}

 む~、確率分布で表現できるバラツキって、こんなもんか。0付近でばらつかせるには、回帰式自体をいぢらんといかんのだろうか。そしたら、やっぱり階層ベイズでということになるのかも。
 まだまだ勉強不足っす。