零膨胀泊松模型不适合

Question

这是数据和设置:

library(fitdistrplus)
library(gamlss)

finalVector <- dput(finalVector)
c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 2, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 1, 2, 1, 1, 1, 1, 1, 
1, 1, 2, 1, 1, 2, 1, 1, 1, 1, 2, 1, 3, 2, 1, 1, 1, 1, 1, 1, 2, 
2, 1, 4, 2, 3, 1, 2, 3, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 
2, 1, 2, 2, 1, 1, 4, 1, 2, 2, 1, 1, 1, 1, 1, 2, 1, 1, 2, 2, 1, 
2, 1, 1, 4, 2, 2, 1, 1, 2, 1, 1, 1, 1, 1, 1)

countFitPoisson <- fitdist(finalVector,"pois",  method = "mle", lower = 0)
countFitZeroPoisson <- fitdist(finalVector, 'ZIP', start = list( ##mu  = mean of poisson, sigma = prob(x = 0))
                                                           mu = as.numeric(countFitPoisson$estimate), 
                                                            sigma = as.numeric(as.numeric(countFitPoisson$estimate))
                                                          ), method = "mle", lower= 0) 

第一个呼叫成功运行.决赛说它没有估计,我不知道为什么.谢谢!

编辑:

假设我正确地完成了代码(不确定),那么我唯一能想到的是没有足够的零来适应模型？

Answer 1

您的数据实际上并非零膨胀,因此拟合模型不会带来改进.fitdistr我没有使用这种方法,而是在glm下面使用和扩展回归模型.但是,所有回归(子)模型只使用常数(或截距)而没有任何实际的回归量.为了可视化,我使用countregR-Forge提供的包(它包含pscl计数数据回归的后继函数)来使用根图.

首先,让我们看看Poisson拟合:

(mp <- glm(finalVector ~ 1, family = poisson))
## Call:  glm(formula = finalVector ~ 1, family = poisson)
## 
## Coefficients:
## (Intercept)  
##      -0.284  
## 
## Degrees of Freedom: 181 Total (i.e. Null);  181 Residual
## Null Deviance:      200.2 
## Residual Deviance: 200.2        AIC: 418.3

这对应于拟合平均值exp(-0.284),即约0.753.如果您比较观察到的和适合的频率,这非常适合数据:

library("countreg")
rootogram(mp)

这表明计数0,1,2的拟合基本上是完美的,并且对于3,4,5只有很小的偏差,但这些频率无论如何都非常低.因此,从这一点来看,似乎不需要扩展模型.

但是为了正式地将模型与其他模型进行比较,可以考虑使用零膨胀泊松(如您所尝试的)障碍泊松或负二项式.零膨胀模型产生:

(mzip <- zeroinfl(finalVector ~ 1 | 1, dist = "poisson"))
## Call:
## zeroinfl(formula = finalVector ~ 1 | 1, dist = "poisson")
## 
## Count model coefficients (poisson with log link):
## (Intercept)  
##     -0.2839  
## 
## Zero-inflation model coefficients (binomial with logit link):
## (Intercept)  
##      -9.151  

因此,计数平均值基本上与之前相同,并且零膨胀的概率基本上为零(plogis(-9.151)约为0.01%).

障碍模型的工作方式类似,但可以使用0-vs-greater的零删失泊松模型和正数的截断泊松模型.然后它也嵌套在泊松模型中,因此可以轻松地进行Wald测试:

mhp <- hurdle(finalVector ~ 1 | 1, dist = "poisson", zero.dist = "poisson")
hurdletest(mhp)
## Wald test for hurdle models
## 
## Restrictions:
## count_((Intercept) - zero_(Intercept) = 0
## 
## Model 1: restricted model
## Model 2: finalVector ~ 1 | 1
## 
##   Res.Df Df Chisq Pr(>Chisq)
## 1    181                    
## 2    180  1 0.036     0.8495

这也清楚地表明没有多余的零和简单的泊松模型就足够了.

作为最终检查,我们还可以考虑负二项模型:

(mnb <- glm.nb(finalVector ~ 1))
## Call:  glm.nb(formula = finalVector ~ 1, init.theta = 125.8922776, link = log)
## 
## Coefficients:
## (Intercept)  
##      -0.284  
## 
## Degrees of Freedom: 181 Total (i.e. Null);  181 Residual
## Null Deviance:      199.1 
## Residual Deviance: 199.1        AIC: 420.3

这又具有几乎相同的均值和巨大的θ参数,其足够接近无穷大(=泊松).因此,总体而言泊松模型就足够了,并且不需要考虑任何扩展.可能性几乎没有变化,附加参数(零通胀,零障碍,θ-离散)不会产生任何改善:

AIC(mp, mzip, mhp, mnb)
##      df      AIC
## mp    1 418.2993
## mzip  2 420.2996
## mhp   2 420.2631
## mnb   2 420.2959