1

我想用定制的负指数概率密度函数拟合一些字段数据。（动机——我最终希望将现场数据拟合到Bullock Shea 和 Skarpaas 2006表 3 中的许多分布）。

首先，我dnegexp根据这篇文章定义了一个函数： Error with custom density function definition for mle2 formula call

dnegexp <- function(x, mya, myb, log=FALSE){
   logne <- log(mya)-myb*x
   if(log) return(logne) else return(exp(logne))
}

然后我创建了一个rnegexp函数，它使用任何给定的参数从这个分布中生成一个数据集。（请原谅这个功能的不雅 - 这不是这篇文章的主题，但欢迎提出更好的想法......）

rnegexp <- function(n, thisa, thisb){
  rnums <- array(data=NA,dim=n)
  counter=0
  while(counter<=n){
    candidate <- runif(1,0,100)
    if(runif(1,0,1)<dnegexp(candidate,thisa, thisb)) {
      rnums[counter] <- candidate
      counter <- counter+1
    }
  }
  return(rnums)
}

有了这两个函数，我通过创建一个数据集来测试我的工作流程，该数据集遵循具有已知参数的负指数分布：

set.seed(501)
mynegexpvals <- rnegexp(100,0.08, 0.09)

hist(mynegexpvals, freq=FALSE)
mydists <- seq(0,100, by=1)
lines(mydists,dnegexp(mydists, 0.08, 0.09), col="blue", lwd=2)

生成值的直方图与生成 pdf 的行

当我尝试mle2从bbmle包中使用来查找参数时，它给出了无意义的值，即使我给了它确切的生成参数作为起始值：

> library(bbmle)
> mle2(mynegexpvals ~ dnegexp(a,b), start = list(a=0.08, b=0.09), data=data.frame(mynegexpvals))

Call:
mle2(minuslogl = mynegexpvals ~ dnegexp(a, b), start = list(a = 0.08, 
    b = 0.09), data = data.frame(mynegexpvals))

Coefficients:
            a             b 
 2.421577e+12 -6.849330e+12 

Log-likelihood: 6.148807e+15

尝试限制搜索空间会在边界处产生参数估计：

> mle2(mynegexpvals ~ dnegexp(a,b), start = list(a=0.08, b=0.09), data=data.frame(mynegexpvals), method="L-BFGS-B", lower=c(a=0.04, b=0.0001), upper=c(a=1000, b=1000))

Call:
mle2(minuslogl = mynegexpvals ~ dnegexp(a, b), start = list(a = 0.08, 
    b = 0.09), method = "L-BFGS-B", data = data.frame(mynegexpvals), 
    lower = c(a = 0.04, b = 1e-04), upper = c(a = 1000, b = 1000))

Coefficients:
    a     b 
1e+03 1e-04 

Log-likelihood: 690.69 
Warning message:
In mle2(mynegexpvals ~ dnegexp(a, b), start = list(a = 0.08, b = 0.09),  :
  some parameters are on the boundary: variance-covariance calculations based on Hessian may be unreliable

我在这里错过了什么重要的东西吗？我是否定义了我的 dnegexp 函数错误？

4

0 回答 0