运行 mle2 函数时出错 (bbmle)
Error when running mle2 function (bbmle)
当 运行 R 中 bbmle
包中的 mle2()
函数时,我收到以下错误:
some parameters are on the boundary: variance-covariance calculations based on Hessian may be unreliable
我想了解这是因为我的数据有问题还是因为正确调用函数有问题。不幸的是,我无法 post 我的真实数据,所以我使用了相同样本大小的类似工作示例。
我使用的自定义 dAction
函数是一个 softmax 函数。优化必须有上限和下限,所以我使用 L-BFGS-B 方法。
library(bbmle)
set.seed(3939)
### Reproducible data
dat1 <- rnorm(30, mean = 3, sd = 1)
dat2 <- rnorm(30, mean = 3, sd = 1)
dat1[c(1:3, 5:14, 19)] <- 0
dat2[c(4, 15:18, 20:22, 24:30)] <- 0
### Data variables
x <- sample(1:12, 30, replace = TRUE)
pe <- dat1
ne <- dat2
### Likelihood
dAction <- function(x, a, b, t, pe, ne, log = FALSE) {
u <- exp(((x - (a * ne) - (b * pe)) / t))
prob <- u / (1 + u)
if(log) return(prob) else return(-sum(log(prob)))
}
### Fit
fit <- mle2(dAction,
start = list(a = 0.1, b = 0.1, t = 0.1),
data = list(x = x, pe = pe, ne = ne),
method = "L-BFGS-B",
lower = c(a = 0.1, b = 0.1, t = 0.1),
upper = c(a = 10, b = 1, t = 10))
Warning message:
In mle2(dAction, start = list(a = 0.1, b = 0.1, t = 0.1), data = list(x = x, :
some parameters are on the boundary: variance-covariance calculations based on Hessian may be unreliable
以下是 summary()
的结果:
summary(fit)
Maximum likelihood estimation
Call:
mle2(minuslogl = dAction, start = list(a = 0.1, b = 0.1, t = 0.1),
method = "L-BFGS-B", data = list(x = x, pe = pe, ne = ne),
lower = c(a = 0.1, b = 0.1, t = 0.1), upper = c(a = 10, b = 1,
t = 10))
Coefficients:
Estimate Std. Error z value Pr(z)
a 0.1 NA NA NA
b 0.1 NA NA NA
t 0.1 NA NA NA
-2 log L: 0.002048047
Warning message:
In sqrt(diag(object@vcov)) : NaNs produced
以及置信区间的结果
confint(fit)
Profiling...
2.5 % 97.5 %
a NA 1.0465358
b NA 0.5258828
t NA 1.1013322
Warning messages:
1: In sqrt(diag(object@vcov)) : NaNs produced
2: In .local(fitted, ...) :
Non-positive-definite Hessian, attempting initial std err estimate from diagonals
我不完全理解你的问题的背景,但是:
这个问题(是否是一个真正的问题在很大程度上取决于我不理解的上述上下文)与您的约束有关。如果我们在没有约束的情况下进行拟合:
### Fit
fit <- mle2(dAction,
start = list(a = 0.1, b = 0.1, t = 0.1),
data = list(x = x, pe = pe, ne = ne))
## method = "L-BFGS-B",
## lower = c(a = 0.1, b = 0.1, t = 0.1),
## upper = c(a = 10, b = 1, t = 10))
我们得到的系数低于您的界限。
coef(fit)
a b t
0.09629301 0.07724332 0.02405173
如果这是正确的,则至少有一个约束将处于活动状态(即,当我们使用下限进行拟合时,至少有一个参数将达到界限 - 事实上,它是所有参数)。当拟合在边界上时,用于计算置信区间(Wald 区间)的最简单机制不起作用。但是,这不会影响您在上面报告的配置文件置信区间估计值。这些是正确的 - 下限报告为 NA
因为置信下限在边界处(如果你愿意,你可以用 0.1 替换它们)。
如果你没想到最优拟合在边界上,那我就不知道怎么回事了,可能是数据问题。
你的对数似然函数没有错,但有点令人困惑,因为你有一个 log
论点,即 returns 负对数似然当 log=FALSE
(默认)和log=TRUE
时的可能性。在我意识到这一点之前,我重写了函数(我还通过尽可能在对数尺度上进行计算,使其在数值上更加稳定)。
dAction <- function(x, a, b, t, pe, ne) {
logu <- (x - (a * ne) - (b * pe)) / t
lprob <- logu - log1p(exp(logu))
return(-sum(lprob))
}
当 运行 R 中 bbmle
包中的 mle2()
函数时,我收到以下错误:
some parameters are on the boundary: variance-covariance calculations based on Hessian may be unreliable
我想了解这是因为我的数据有问题还是因为正确调用函数有问题。不幸的是,我无法 post 我的真实数据,所以我使用了相同样本大小的类似工作示例。
我使用的自定义 dAction
函数是一个 softmax 函数。优化必须有上限和下限,所以我使用 L-BFGS-B 方法。
library(bbmle)
set.seed(3939)
### Reproducible data
dat1 <- rnorm(30, mean = 3, sd = 1)
dat2 <- rnorm(30, mean = 3, sd = 1)
dat1[c(1:3, 5:14, 19)] <- 0
dat2[c(4, 15:18, 20:22, 24:30)] <- 0
### Data variables
x <- sample(1:12, 30, replace = TRUE)
pe <- dat1
ne <- dat2
### Likelihood
dAction <- function(x, a, b, t, pe, ne, log = FALSE) {
u <- exp(((x - (a * ne) - (b * pe)) / t))
prob <- u / (1 + u)
if(log) return(prob) else return(-sum(log(prob)))
}
### Fit
fit <- mle2(dAction,
start = list(a = 0.1, b = 0.1, t = 0.1),
data = list(x = x, pe = pe, ne = ne),
method = "L-BFGS-B",
lower = c(a = 0.1, b = 0.1, t = 0.1),
upper = c(a = 10, b = 1, t = 10))
Warning message:
In mle2(dAction, start = list(a = 0.1, b = 0.1, t = 0.1), data = list(x = x, :
some parameters are on the boundary: variance-covariance calculations based on Hessian may be unreliable
以下是 summary()
的结果:
summary(fit)
Maximum likelihood estimation
Call:
mle2(minuslogl = dAction, start = list(a = 0.1, b = 0.1, t = 0.1),
method = "L-BFGS-B", data = list(x = x, pe = pe, ne = ne),
lower = c(a = 0.1, b = 0.1, t = 0.1), upper = c(a = 10, b = 1,
t = 10))
Coefficients:
Estimate Std. Error z value Pr(z)
a 0.1 NA NA NA
b 0.1 NA NA NA
t 0.1 NA NA NA
-2 log L: 0.002048047
Warning message:
In sqrt(diag(object@vcov)) : NaNs produced
以及置信区间的结果
confint(fit)
Profiling...
2.5 % 97.5 %
a NA 1.0465358
b NA 0.5258828
t NA 1.1013322
Warning messages:
1: In sqrt(diag(object@vcov)) : NaNs produced
2: In .local(fitted, ...) :
Non-positive-definite Hessian, attempting initial std err estimate from diagonals
我不完全理解你的问题的背景,但是:
这个问题(是否是一个真正的问题在很大程度上取决于我不理解的上述上下文)与您的约束有关。如果我们在没有约束的情况下进行拟合:
### Fit
fit <- mle2(dAction,
start = list(a = 0.1, b = 0.1, t = 0.1),
data = list(x = x, pe = pe, ne = ne))
## method = "L-BFGS-B",
## lower = c(a = 0.1, b = 0.1, t = 0.1),
## upper = c(a = 10, b = 1, t = 10))
我们得到的系数低于您的界限。
coef(fit)
a b t
0.09629301 0.07724332 0.02405173
如果这是正确的,则至少有一个约束将处于活动状态(即,当我们使用下限进行拟合时,至少有一个参数将达到界限 - 事实上,它是所有参数)。当拟合在边界上时,用于计算置信区间(Wald 区间)的最简单机制不起作用。但是,这不会影响您在上面报告的配置文件置信区间估计值。这些是正确的 - 下限报告为 NA
因为置信下限在边界处(如果你愿意,你可以用 0.1 替换它们)。
如果你没想到最优拟合在边界上,那我就不知道怎么回事了,可能是数据问题。
你的对数似然函数没有错,但有点令人困惑,因为你有一个 log
论点,即 returns 负对数似然当 log=FALSE
(默认)和log=TRUE
时的可能性。在我意识到这一点之前,我重写了函数(我还通过尽可能在对数尺度上进行计算,使其在数值上更加稳定)。
dAction <- function(x, a, b, t, pe, ne) {
logu <- (x - (a * ne) - (b * pe)) / t
lprob <- logu - log1p(exp(logu))
return(-sum(lprob))
}