glm起始值不接受log-link

And*_*ndi 5 r offset glm

我想运行带有日志链接和偏移量的高斯GLM.出现以下问题:

y <- c(1,1,0,0)
t <- c(5,3,2,4)
Run Code Online (Sandbox Code Playgroud)

没问题:

exp(coef(glm(y~1 +  offset(log(t)), family=poisson)))
Run Code Online (Sandbox Code Playgroud)

family=gaussian,需要指定起始值,它在这里工作:

exp(coef(glm(y~1, family=gaussian(link=log), start=0)))
Run Code Online (Sandbox Code Playgroud)

但在这里不起作用:

exp(coef(glm(y~1 +  offset(log(t)), family=gaussian(link=log), start=0)))
Run Code Online (Sandbox Code Playgroud)

eval中的错误(expr,envir,enclos):找不到有效的起始值:请指定一些"

有没有人看到什么是错的(希望只是在我的编码中)?

Jam*_*mes 10

我看起来好像start没有被识别出来offset.您试图在y值中取0的日志-Inf.glm在没有得到一些帮助的情况下寻找解决方案时,显然无法解决这个问题start.您的y值中的小扰动将允许解决方案.

exp(coef(glm(I(y+.Machine$double.eps)~1 + offset(log(t)), family=gaussian(link=log))))
(Intercept) 
  0.1481481
Run Code Online (Sandbox Code Playgroud)


Ben*_*ker 10

以下是一些考古学的结果,它解释了正在发生的事情,在glm功能的深处:

调试(使用debug("glm"))和单步调试该函数表明它在以下调用时失败:

if (length(offset) && attr(mt, "intercept") > 0L) {
  fit$null.deviance <- eval(call(if (is.function(method)) "method" else method, 
    x = X[, "(Intercept)", drop = FALSE], y = Y, weights = weights, 
    offset = offset, family = family, control = control, 
    intercept = TRUE))$deviance
}
Run Code Online (Sandbox Code Playgroud)

这是尝试计算模型的零偏差.仅评估是否存在截距项和偏移项(我不确定原因;可能是前一次调用计算的默认空偏差glm在这种情况下是错误的并且必须重新计算?).它调用glm.fit(默认值method),但没有起始值,因为对于仅拦截模型,这些通常是不必要的.

现在调试里面glm.fit看看会发生什么:我们得到

  if (is.null(etastart) && is.null(start) && is.null(mustart) && 
    ((family$link == "inverse" && any(y == 0)) || (family$link == 
        "log" && any(y <= 0))))
    stop("cannot find valid starting values: please specify some")
Run Code Online (Sandbox Code Playgroud)

我们看到因为没有传递起始值,因为使用了日志链接,并且因为某些gaussian()值等于零,所以拟合失败.因此,如果(并且仅在?)指定了偏移量和截距,使用了日志链接,并且响应中存在零值,则应该发生这种情况.

如果你y; 添加线

    start = start[1], etastart = etastart[1], mustart = mustart[1],
Run Code Online (Sandbox Code Playgroud)

符合零偏差的呼叫(即上面所示的呼叫); 而且dump("glm",file="glmtemp.R"),似乎工作正常......我认为这应该是一个合理的通用解决方案.如果有人想在R开发列表中提出这个问题,请随意.

  • 挖得好!:) 我认为有一个更好的解决办法:只需添加“mustart = fit$fitted.values”。我已在 https://bugs.r-project.org/bugzilla/show_bug.cgi?id=16877 提交了补丁。 (2认同)