高斯分佈和伽馬分佈的混合

Question

我正在尋找一些腳本/軟件包在R（Python也會這樣做）從高斯和伽瑪分佈的混合物中找出組件分佈參數。到目前爲止，我已經使用R包「mixtools」將 作爲高斯混合數據進行建模，但我認爲它可以通過Gamma加高斯模型更好地建模。

感謝

Answer 1

這裏有一個可能性：

定義效用函數：

rnormgammamix <- function(n,shape,rate,mean,sd,prob) { 
    ifelse(runif(n)<prob, 
      rgamma(n,shape,rate), 
      rnorm(n,mean,sd)) 
} 

（這可能是做多一點點高效...）

dnormgammamix <- function(x,shape,rate,mean,sd,prob,log=FALSE) { 
    r <- prob*dgamma(x,shape,rate)+(1-prob)*dnorm(x,mean,sd) 
    if (log) log(r) else r 
} 

生成假數據：

set.seed(101) 
r <- rnormgammamix(1000,1.5,2,3,2,0.5) 
d <- data.frame(r) 

方法＃1：bbmle包。適合的形狀，比率，對數標準差，對數尺度的概率。

library("bbmle") 
m1 <- mle2(r~dnormgammamix(exp(logshape),exp(lograte),mean,exp(logsd), 
        plogis(logitprob)), 
    data=d, 
    start=list(logshape=0,lograte=0,mean=0,logsd=0,logitprob=0)) 
cc <- coef(m1) 

png("normgam.png") 
par(bty="l",las=1) 
hist(r,breaks=100,col="gray",freq=FALSE) 
rvec <- seq(-2,8,length=101) 
pred <- with(as.list(cc), 
      dnormgammamix(rvec,exp(logshape),exp(lograte),mean, 
          exp(logsd),plogis(logitprob))) 
lines(rvec,pred,col=2,lwd=2) 
true <- dnormgammamix(rvec,1.5,2,3,2,0.5) 
lines(rvec,true,col=4,lwd=2) 
dev.off() 

tcc <- with(as.list(cc), 
      c(shape=exp(logshape), 
       rate=exp(lograte), 
       mean=mean, 
       sd=exp(logsd), 
       prob=plogis(logitprob))) 
cbind(tcc,c(1.5,2,3,2,0.5)) 

擬合是合理的，但參數是相當遙遠 - 我認爲這種模式是不是在這個參數政權（即Gamma和高斯非常強烈的可識別部件可以互換）

library("MASS") 
ff <- fitdistr(r,dnormgammamix, 
    start=list(shape=1,rate=1,mean=0,sd=1,prob=0.5)) 

cbind(tcc,ff$estimate,c(1.5,2,3,2,0.5)) 

fitdistr得到相同的結果mle2，這表明我們 在當地最低。如果我們從真實參數開始，我們將得到 爲合理值並接近真實參數。

ff2 <- fitdistr(r,dnormgammamix, 
    start=list(shape=1.5,rate=2,mean=3,sd=2,prob=0.5)) 
-logLik(ff2) ## 1725.994 
-logLik(ff) ## 1755.458