mle2公式调用的自定义密度函数定义错误 [英] Error with custom density function definition for mle2 formula call
问题描述
我想定义自己的密度函数,以用于R
的bbmle
包对mle2
的公式调用中.模型的参数是估计值,但是我无法在返回的mle2
对象上应用residuals
或predict
之类的函数.
I want to define my own density function to be used in the formula call to mle2
from R
's bbmle
package. The parameters of the model are estimated but I cannot apply functions like residuals
or predict
on the returned mle2
object.
这是一个示例,其中我为简单的泊松模型定义了一个函数.
This is an example in which I define a function for a simple Poisson model.
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
上调用残差时,出现错误,提示:
In myfit
, the lambda is estimated correctly but when I call residuals on myfit
, I get an error which says:
Error in myf(9.77598906811668) :
argument "lambda" is missing, with no default
另一方面,如果我使用R
的内置dpois
函数简单地按如下所示拟合模型,则将计算残差:
On the other hand, if I simply fit the model as follows using R
's built-in dpois
function, the residuals are computed:
myfit <- mle2(hpoisson ~ dpois(lambda), start = list(lambda=9), data=data.frame(hpoisson))
residuals(myfit)
任何人都可以告诉我myf
的函数定义中我做错了什么吗?
Could anyone please tell me what I am doing wrong in the function definition of myf
?
谢谢
推荐答案
文档中没有很清楚地说明,但是使用自定义密度函数有一些先决条件:
It's not very clearly explained in the documents, but there are a few prerequisites for using custom density functions:
- 函数名称必须以
d
开头,必须具有第一个自变量x
,并且必须具有命名自变量log
. (log
参数必须做一些明智的事情:尤其是mle2
将使用log=TRUE
调用该函数,并且该函数最好返回对数似然!)通常,尽管它不是必需的,直接计算对数似然然后在log=FALSE
处取幂,而不是计算可能性并在log=TRUE
处记录它,在数字上更明智(在某些情况下,例如零膨胀模型,真的不可行).例如,将我的dmyf()
定义与OP代码中的myf()
定义进行比较... - 为了使用诸如
predict
之类的其他方法,您必须定义一个名称以s
开头的附加函数.它返回指定参数的弯矩列表,摘要统计信息等-参见下面的示例,该示例是从bbmle::spois
复制而来的.
- the function's name must start with
d
, must have first argumentx
, and must have a named argumentlog
. (Thelog
argument has to do something sensible: in particular,mle2
will call the function withlog=TRUE
, and the function had better return the log-likelihood!) In general, although it's not required, it's more numerically sensible to compute the log-likelihood directly and then exponentiate iflog=FALSE
, rather than computing the likelihood and logging it iflog=TRUE
(there are cases, such as zero-inflated models, where this isn't really feasible). For example, compare mydmyf()
definition with themyf()
definition in the OP's code ... - in order to use additional methods such as
predict
you have to define an additional function whose name starts withs
; it returns a list of moments, summary statistics, etc. for a specified parameter -- see example below, which is copied frombbmle::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)
这篇关于mle2公式调用的自定义密度函数定义错误的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!