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

Question

我想知道如何使用 MLE 的方法将任何分布拟合到给定的数据集。作为一个特定的例子，当 https://rivista-statistica.unibo.it/article/viewFile/6836/7039 中描述的广义林德利分布应用于数据时，任何人都可以建议一个工作代码，该代码将为 $\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？那么如何使用它来拟合直方图上的分布呢？

Answer 1

我将根据 Parari et al. 中的一个问题（Proschan 的空调数据（Table 7.2，结果为Table7.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))
}

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

Data fitting by the method of maximum likelihood for a new distribution

r

histogram

mle