0

我有一些实际数据，恐怕有些令人讨厌。

它本质上是一个正负二项分布（没有任何零计数）。但是，有些异常值似乎会导致一些错误的计算发生（可能是下溢或 NaN？）前 8 个左右的条目是合理的，但我猜最后几个会导致拟合出现一些问题。

这是数据：

> df
   counts  t
1    1968  1
2     217  2
3      55  3
4      26  4
5      11  5
6       5  6
7       8  7
8       3  8
9       1 10
10      1 11
11      1 12
12      1 13
13      1 15
14      1 18
15      1 26
16      1 59

该命令运行一段时间，然后吐出错误消息

> vglm(counts ~ t, data=df, family = posnegbinomial)
Error in if (take.half.step) { : missing value where TRUE/FALSE needed

但是，如果我重新运行这个切断异常值，我会得到一个 posnegbinomial 的解决方案

> vglm(counts ~ t, data=df[1:9,], family = posnegbinomial)
Call:
vglm(formula = counts ~ t, family = posnegbinomial, data = df[1:9,])

Coefficients:
(Intercept):1 (Intercept):2             t 
    7.7487404     0.7983811    -0.9427189 

Degrees of Freedom: 18 Total; 15 Residual
Log-likelihood: -36.21064

如果我尝试家庭 pospoisson（正泊松：没有零值），我会收到类似的错误“参数不可解释为逻辑”。

我确实注意到 Stackoverflow 中有许多类似的问题，关于需要 TRUE/FALSE 的缺失值，但对于其他 R 包。这向我表明，也许包编写者需要更好地预测计算可能会失败。

4

1 回答 1

1

我认为您的近端问题是，您的极值的负二项式的预测均值非常接近于零，以至于它们下溢到零，这是包作者没有预料到/保护的。（关于非线性优化/拟合要实现的一件事是，总是有可能通过给它极端数据来破坏拟合方法......）

我无法让它在 中工作VGAM，但我会提供一些其他建议。

plot(log(counts)~t,data=dd)

并观察数据以获得参数值的初始估计（至少对于平均模型）：

m0 <- lm(log(counts)~t,data=subset(dd,t<10))

我以为我可以vglm()通过设置起始值来开始工作，但这实际上并没有成功，即使我从其他平台获得了相当好的值（见下文）。

glmmADMB

该glmmADMB软件包可以通过以下方式处理正 NB family="truncnbinom"：

library(glmmADMB)
m1 <- glmmadmb(counts~t, data=dd, family="truncnbinom")

（有一些警告信息......）

bbmle::mle2()

这需要更多的工作：它在标准模型上失败了，但是如果我在预测的平均值上设置一个下限就可以了......

library(VGAM)  ## for dposnegbin
library(bbmle)
m2 <- mle2(counts~dposnegbin(size=exp(logk),
                         munb=pmax(exp(logeta),1e-7)),
           parameters=list(logeta~t),
           data=dd,
           start=list(logk=0,logeta=0))

再次警告信息。

比较glmmADMB, mle2, 简单截断lm拟合...

cc <- cbind(coef(m2),
  c(log(m1$alpha),coef(m1)),
  c(NA,coef(m0)))
dimnames(cc) <- list(c("log_k","log_int","slope"),

                 c("mle2","glmmADMB","lm"))

##               mle2   glmmADMB         lm
## log_k    0.8094678  0.8094625         NA
## log_int  7.7670604  7.7670637  7.1747551
## slope   -0.9491796 -0.9491778 -0.8328487

原则上这也是可能的，但它遇到了与...glmmTMB相同的问题。vglm()

于 2016-02-14T23:49:57.333 回答