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

草薙の研究ログ

英語教育関係。でも最近は統計(特にR)ネタが中心。

ゼロ過剰ポアソン分布をデータにフィットさせる

ゼロ過剰ポアソン分布では2つのプロセスが考えられている。まず,最初に0の確率がσである二項分布で,値が0でないときのカウントデータがμに従うという。統計モデリングにおいて「ゼロが多いときに使うといい」とはよく聞くものの,厳密にいえば,このような二段階の過程が十分に考えられる事象を対象にするといいと思う。そして,殆どの場合は,先んじている0かそうでないかということにはあまり関心がないときに。

#使うパッケージ
library(fitdistrplus)
library(gamlss)

#ある実測のデータ
d<-c(0,5,6,0,1,8,0,2,0,1,4,2,0,3,8,9,1,1,0,10,1,2,0,2,11,3,0,8,5,11,3,8,8,8,14,7,8,6,5,1,8,7,2,0,6,8,3,3,0,6,2,11,2,0,0,0,0,0,0)

#可視化
hist(d)

#0の割合
length(d[d==0])/length(d)

#まずはポアソンを普通にデータにフィット
fit.pois<-fitdist(d,"pois")
summary(fit.pois)

#次にゼロ過剰ポアソン分布をデータにフィット
fit.ZIP<-fitdist(d,"ZIP",start=1)
summary(fit.ZIP)

#2つのモデルを比較してみる
gofstat(list(fit.pois,fit.ZIP),fitnames=c("Poisson","ZIP"))

#これではσを出してくれない。gamlssを使う
fit.ZIP2<-gamlss(d~1,family=ZIP)
fit.ZIP2


うむ。