以下是我的代码：
library(fitdistrplus)
s <- c(11, 4, 2, 9, 3, 1, 2, 2, 3, 2, 2, 5, 8,3, 15, 3, 9, 22, 0, 4, 10, 1, 9, 10, 11,
2, 8, 2, 6, 0, 15, 0 , 2, 11, 0, 6, 3, 5, 0, 7, 6, 0, 7, 1, 0, 6, 4, 1, 3, 5,
2, 6, 0, 10, 6, 4, 1, 17, 0, 1, 0, 6, 6, 1, 5, 4, 8, 0, 1, 1, 5, 15, 14, 8, 1,
3, 2, 9, 4, 4, 1, 2, 18, 0, 0, 10, 5, 0, 5, 0, 1, 2, 0, 5, 1, 1, 2, 3, 7)
o <- fitdist(s, "gamma", method = "mle")
summary(o)
plot(o)
错误说：
fitdist(s, "gamma", method = "mle") 中的错误：函数 mle 无法估计参数，错误代码为 100
Gamma 分布不允许零值（对于 0 的响应，可能性将评估为零，并且对数似然将是无限的），除非形状参数正好是 1.0（即指数分布 - 见下文）。 ..这是一个统计/数学问题，而不是编程问题。您将不得不找到一些明智的做法来处理零值。根据您的应用程序的意义，您可以（例如）
fitdist(s[s>0], ...))fitdist(replace(s,which(s==0),0.1),...)其中哪一个（如果有的话）最好取决于您的应用程序。
@Sandipan Dey 的第一个答案（在数据集中保留零）似乎是有道理的，但实际上它被困在等于 1 的形状参数上。
o <- fitdist(s, "exp", method = "mle")
给出与@Sandipan 的代码相同的答案（除了它估计速率=0.2161572，尺度参数的倒数=4.626262，它是为 Gamma 分布估计的——这只是参数化的变化）。如果您选择拟合指数而不是伽玛，那很好 - 但您应该故意这样做，而不是偶然......
为了说明包含零的拟合可能无法按预期工作，我将构建自己的负对数似然函数并显示每种情况的似然面。
mfun <- function(sh,sc,dd=s) {
-sum(dgamma(dd,shape=sh,scale=sc,log=TRUE))
}
library(emdbook) ## for curve3d() helper function
包含零的表面：
cc1 <- curve3d(mfun(x,y),
## set up "shape" limits" so we evaluate
## exactly shape=1.000 ...
xlim=c(0.55,3.55),
n=c(41,41),
ylim=c(2,5),
sys3d="none")
png("gammazero1.png")
with(cc1,image(x,y,z))
dev.off()
在这种情况下，表面仅定义为 shape=1（即指数分布）；白色区域代表无限对数似然。并不是 shape=1 是最合适的，而是它是唯一合适的...
零排除表面：
cc2 <- curve3d(mfun(x,y,dd=s[s>0]),
## set up "shape" limits" so we evaluate
## exactly shape=1.000 ...
xlim=c(0.55,3.55),
n=c(41,41),
ylim=c(2,5),
sys3d="none")
png("gammazero2.png")
with(cc2,image(x,y,z))
with(cc2,contour(x,y,z,add=TRUE))
abline(v=1.0,lwd=2,lty=2)
dev.off()
mle只需提供要使用 using计算的 gamma 分布参数（比例、形状）的初始值以及参数optim的下限，它应该可以工作。
o <- fitdist(s, "gamma", lower=c(0,0), start=list(scale=1,shape=1))
summary(o)
#Fitting of the distribution ' gamma ' by maximum likelihood
#Parameters :
# estimate Std. Error
#scale 4.626262 NA
#shape 1.000000 NA
#Loglikelihood: -250.6432 AIC: 505.2864 BIC: 510.4766
根据@Ben Bolker 的评论，我们可能希望首先排除零点：
o <- fitdist(s[s!=0], "gamma", method = "mle", lower=c(0,0), start=list(scale=1,shape=1))
summary(o)
#Fitting of the distribution ' gamma ' by maximum likelihood
#Parameters :
# estimate Std. Error
#scale 3.401208 NA
#shape 1.622378 NA
#Loglikelihood: -219.6761 AIC: 443.3523 BIC: 448.19
