混合ガンマ分布モデルをデータにフィットさせる
Rのmixtoolsパッケージでは,混合ガンマ分布モデルもデータへフィットさせることができる。
#λ= .25,.75,α - 4, 10,β= 1, .2が正解 set.seed(0) dat<-c(rgamma(100,4,1),rgamma(300,10,.2)) hist(dat,breaks=20, main="",col="lightblue") model<-gammamixEM(dat) model$gamma.pars model$lambda
まあまあ,て感じ。これ結構混ざっちゃうと,うまくいかないかもね。
混合ガンマ分布の確率密度関数は,
dmixgamma<-function(x,lambda,a1,a2,b1,b2){ y<-lambda*dgamma(x,a1,b1)+(1-lambda)*dgamma(x,a2,b2) y }
これで絵を描くとこう。
hist(dat,breaks=20, main="",col="lightblue",freq=F,ylim=c(0,.06),xlim=c(0,100)) x<-seq(0,100,.01) lines(x,dmixgamma(x,0.25,3.05,11.54,1/1.23,1/4.37),lwd=2,col=2)
…これ,単変量で,しかも前もって2つの混合だってわかっているんだったら,EMアルゴリズムでなくても,この程度の単純さなら普通の(?)L-BFGSとかでいけそうだよね。
そしたら結局はfitdist,mle,bbmle,optimとかでも同じことできそう。
指数正規合成分布とガンマ分布の混合とか,そういうのも自由にできるはず。おもったより簡単ぽい。