mle2 公式调用的自定义密度函数定义错误

Question

我想定义我自己的密度函数，用于从 R 的 bbmle 包调用 mle2 的公式。估计了模型的参数，但我无法在返回的 mle2 对象上应用 residuals 或 predict 等函数。

这是我为简单泊松模型定义函数的示例。

library(bbmle)

set.seed(1)
hpoisson <- rpois(1000, 10)

myf <- function(x, lambda, log = FALSE) {
  pmf <- (lambda^x)*exp(-lambda)/factorial(x)
  if (log)
    log(pmf)
  else
    pmf
}

myfit <- mle2(hpoisson ~ myf(lambda), start = list(lambda=9), data=data.frame(hpoisson))
residuals(myfit)

在 myfit 中，lambda 被正确估计，但是当我在 myfit 上调用残差时，我收到一条错误消息：

Error in myf(9.77598906811668) : 
  argument "lambda" is missing, with no default

另一方面，如果我简单地使用 R 的内置 dpois 函数如下拟合模型，则计算残差：

myfit <- mle2(hpoisson ~ dpois(lambda), start = list(lambda=9), data=data.frame(hpoisson))
    residuals(myfit)

谁能告诉我 myf 的函数定义中我做错了什么？

谢谢

Answer 1

文档中解释的不是很清楚，但是使用自定义密度函数有几个先决条件：

函数名称 必须以 d 开头，必须有第一个参数 x，并且必须有一个命名参数 log。（log 参数必须做一些明智的事情：特别是，mle2 会用 log=TRUE 调用函数，函数最好 return 对数似然！）一般来说，虽然这不是必需的，但直接计算对数似然然后在 log=FALSE 时求幂比在 [=16] 时计算似然并记录它在数值上更明智=]（在某些情况下，例如零膨胀模型，这是不可行的）。例如，将我的 dmyf() 定义与 OP 代码中的 myf() 定义进行比较 ...
为了使用额外的方法，例如predict，你必须定义一个名称以s开头的额外函数；它 return 是指定参数的时刻列表、汇总统计信息等 -- 请参见下面的示例，该示例是从 bbmle::spois.

library("bbmle")
set.seed(1)
hpoisson <- rpois(1000, 10)

dmyf <- function(x, lambda, log = FALSE) {
    logpmf <- x*log(lambda)-lambda-lfactorial(x)
    if (log) return(logpmf)  else return(exp(logpmf))
}
smyf <- function(lambda) {
    list(title = "modified Poisson",
         lambda = lambda, mean = lambda,
         median = qpois(0.5, lambda),
         mode = NA, variance = lambda, sd = sqrt(lambda))
}
myfit <- mle2(hpoisson ~ dmyf(lambda),
              start = list(lambda=9), data=data.frame(hpoisson))
residuals(myfit)

Answer 2

不是真正的答案，但需要更多帮助：

我用它来尝试制作一个 "custom" beta 二项式函数来模仿 bbmle 小插图第一位中的那个。

set.seed(1001)
x1 <- rbetabinom(n=1000, prob=0.1, size=50, theta=10)
dmybetabinom <- function(x, N, theta, p, log=FALSE) {
    (choose(N,x)*beta(N-x+theta*(1-p),x+theta*p))/beta(theta*(1-p),theta*p)
}

该函数的工作方式类似于 dbetabinom:

dbetabinom(0:9,size=9,theta=4, prob=0.5) [1] 0.04545455 0.08181818 0.10909091 0.12727273 0.13636364 0.13636364 0.12727273 0.10909091 [9] 0.08181818 0.04545455 dmybetabinom(0:9,N=9,theta=4, p=0.5) [1] 0.04545455 0.08181818 0.10909091 0.12727273 0.13636364 0.13636364 0.12727273 0.10909091 [9] 0.08181818 0.04545455

但是当我尝试在其上使用 mle2 函数时，我遇到了这个错误：

m0fa <- mle2(x1~dmybetabinom( N=50, theta, p), start=list(p=0.2, theta=9), data=data.frame(x1) )`

Error in optim(par = c(0.2, 9), fn = function (p)  :   non-finite finite-difference value [1] `

mle2 公式调用的自定义密度函数定义错误

Error with custom density function definition for mle2 formula call

r

poisson

mle