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，拉姆達正確地估計，但是當我呼籲myfit殘差，我得到它說的錯誤：

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

在另一方面，如果我只是擬合模型如下使用內置的dpois功能R的，殘差計算：

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，並且該函數最好返回對數似然值！）通常，雖然不是必需的，但它在數值上更合理計算直接對數似然比，然後指數log=FALSE，而不是計算可能性並記錄，如果log=TRUE（有些情況下，如零膨脹模型，這是不可行的）。例如，將我的dmyf()定義與OP代碼中的myf()定義進行比較...
• 爲了使用其他方法（如predict），您必須定義一個名稱以s開頭的附加功能;它將返回一個指定參數的時刻列表，彙總統計信息等 - 請參閱下面的示例，該示例從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

不是一個真正的答案，但需要更多這方面的幫助：

我用這個來嘗試做一個「自定義」β-二項功能模仿一個在第一位的小插曲。

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] `