通过最大似然法拟合新分布的数据答案

【问题标题】：Data fitting by the method of maximum likelihood for a new distribution通过最大似然法拟合新分布的数据
【发布时间】：2019-02-26 06:01:07
【问题描述】：

我想知道如何使用 MLE 方法将任何分布拟合到给定的数据集。作为一个特定的例子，当https://rivista-statistica.unibo.it/article/viewFile/6836/7039 中描述的广义 Lindley 分布应用于数据时，任何人都可以建议一个工作代码，该代码将为 $\theta$ 和 $\beta$ 给出正确的 MLE 结果：5.1、6.3、第 10.8、12.1、18.5、19.7、22.2、23、30.6、37.3、46.3、53.9、59.8、66.2 156？然后如何使用它来拟合直方图上的分布？

【问题讨论】：

标签： r histogram mle

【解决方案1】：

我将根据Parari et al. 中的一个问题（来自 Proschan 的空调数据（表 7.2，表 7.5 中的结果）回答您现在已删除的一个非常相似的问题

一般来说，您需要知道合理的起始值才能进行通用的最大似然估计（这无疑是一个先有鸡还是先有蛋的问题，但有很多解决方案 [使用矩量法，根据问题选择合理的值，观察图表来选择值等]。在这种情况下，由于作者给了你答案，你可以使用这些参数来选择正确数量级的起始值。更一般地说，您需要了解一些有关合理数量级的知识才能开始使用。

请注意，一般最大似然估计存在很多繁琐的数值问题：例如见第 7 章here ...

数据 (x) 和似然函数 (dBEPL) 定义如下。我正在定义密度函数并使用mle2()的公式接口...

dd  <- data.frame(x)
library(bbmle)
## parameters from table
ovec <- list(alpha=.7945,beta=0.1509,omega=6.7278,
              a=0.2035,b=0.2303)
## starting vals -- right order of magnitude
svec <- list(alpha=0.8,beta=0.2,omega=10,a=0.2,b=0.2)
m1 <- mle2( x ~ dBEPL(alpha,beta,omega,a,b),
      data=dd,
      start=svec,
      control=list(parscale=unlist(svec)))
coef(m1)  ## NOT the same as reported, but close
par(las=1,bty="l")
hist(x,col="gray",freq=FALSE,breaks=40,ylim=c(0.,0.014))
with(as.list(coef(m1)),curve(dBEPL(x,alpha,beta,omega,a,b,log=FALSE),add=TRUE,
                             col=2,lwd=2,n=201))

x <- scan(textConnection("
194 413 90 74 55 23 97 50 359 50 130 487 57 102 
15 14 10 57 320 261 51 44 9 254 493 33 18 209 41 
58 60 48 56 87 11 102 12 5 14 14 29 37 186 29 104 
7 4 72 270 283 7 61 100 61 502 220 120 141 22 603 35 
98 54 100 11 181 65 49 12 239 14 18 39 3 12 5 32 9 438 
43 134 184 20 386 182 71 80 188 230 152 5 36 79 59 33 
246 1 79 3 27 201 84 27 156 21 16 88 130 14 118 44 15 42 
106 46 230 26 59 153 104 20 206 566 34 29 26 35 5 82 31 
118 326 12 54 36 34 18 25 120 31 22 18 216 139 67 310 
346 210 57 76 14 111 97 62 39 30 7 44 11 63 23 22 23 14 
18 13 34 16 18 130 90 163 208 124 70 16 101 52 
208 95 62 11 191 14 71"))


dBEPL <- function(x,alpha,beta,omega,a,b,log=TRUE) {
    r <- log(alpha*beta^2*omega/(beta(a,b)*(beta+1)))+
           log(1+x^alpha)+(alpha-1)*log( x)-beta* x^alpha+(omega*a-1) *
           log(1-(1+beta* x^alpha/(beta+1))*exp(-beta* x^alpha))+
           (b-1)*log(1-(1-(1+beta* x^alpha/(beta+1))*exp(-beta* x^alpha))^omega)
    if (log) return(r) else return(exp(r))
}

【讨论】：

这是一个很好的答案，我真的很感激！因此，据我了解，m.l.e.s 由 svec 中列出的起始值决定，因此使用起始值 alpha=0.8,beta=0.2,omega=10,a=0.2,b= 时获得的 m.l.e.s 的值0.2 可能与 ovec 中给出的不同？如果没有给我们 ovec 中的条目：alpha=.7945、beta=0.1509、omega=6.7278、a=0.2035 和 b=0.2303，您能否说明如何使用矩量法获得近似值？跨度>
一般来说，人们希望不同的起始值会收敛到几乎相同的答案。我猜这个函数参数化过大/不稳定/在某些区域具有非常平坦的似然面。至于时刻的方法：这样做需要我花更长的时间来弄清楚所有参数的含义。
谢谢。这是一个有趣的案例：pdfs.semanticscholar.org/7713/…。当将上述方法应用于 LE 分布时，我们获得了 NaN 警告，并且（使用估计量的已知值 9.499·10^(−4) 和 2·10^(−6)），由 bbmle 产生的 mles 完全不同表中的那些（尽管参数为正，但 lambda 为负）。这是为什么呢？
很难说。如果您将 完全可重现 示例作为新问题发布，我（或其他人）可能会回答...
我已将问题移至stackoverflow.com/questions/55032275/nan-errors-with-bbmle