r 中具有零膨胀泊松分布的混合 GAM

Question

你好。我目前遇到使用零膨胀数据拟合模型的问题。您可以在此处查看响应变量的分布：

随机效应与变量“site”相关。
数据框在这里：

df <- structure(list(Am06 = c(0, 0, 0, 129, 0, 0, 0, 0, 66, 0), 
                     site = c("BEARavn", "LUISav", "PONTav", "TONKam", "GRANav", "VILOav", "VILOam", "PONTam", "LUCEav", "LUCEav"), 
                     smax = c(19.0284, 14.601, 22.537, 23.236273, 16.392, 21.724, 14.263, 20.507, 20.615, 19.11825), 
                     up_365 = c(0, 4.66666666666667, 2.83333333333333, 365, 11.4583333333333, 0.583333333333333, 0, 3.75, 365, 365), 
                     do_60 = c(0.208333333333333, 0, 57.1666666666667, 0, 36, 0, 0, 0, 0, 0)), 
                     row.names = c(234L, 1469L, 2385L, 2716L, 1075L, 2894L, 2849L, 2341L, 1347L, 1338L), 
                     class = "data.frame")

这是我的代码，带有

mgcv::gam

（没有随机效应）以及以下错误：

mod1 <- mgcv::gam(list(as.integer(log(1+df[,1])) ~ s(smax)+ s(up_365)+ 
                  s(do_365)+ s(site, bs = "re"), ~ 1),
                  data = df, family = mgcv::ziplss())

名称错误（dat）<- object$term : 'names' attribute [1] must be the same length as the vector [0]

gamm4V

也没有随机效应和错误

mod2 <- gamm4V(as.integer(log(1+df[,1])) ~ s(smax) + s(up_365) + 
               s(do_365), data = df, family=ziP ,random=~(1|site))

as(value, fieldClass, strict = FALSE) 中的错误：内部问题 在 as() 中：“extended.family” is(object, "family") 为 TRUE，但是 元数据断言“是”关系为 FALSE

然后我尝试使用具有随机效果的列表系统（如跨栏），但我有另一种类型的错误，我不明白。这是我想要拟合的最重要的模型：

mod3 <- gam(list(as.integer(log(1+df[,1])) ~ s(smax) + s(up_365) + 
                 s(do_365) + s(site, bs = "re"),
                 ~ s(smax) + s(up_365) + s(do_365)),
            data = df, family = ziplss, method = "REML")

名称错误（dat）<- object$term : 'names' attribute [1] must be the same length as the vector [0]

唯一有效的（具有随机效果）似乎是这个：

mod4 <- gam(as.integer(log(1+df[,1])) ~ s(smax) + s(up_365) +
            s(do_365), data = df, random=~(1|site), family=ziP)

您能告诉我如何避免这些错误并至少适合这些模型之一吗？

谢谢你

P.

Answer 1

您只能将

glmer()

与 gamm4 一起使用，而

ziP()

或

ziplss()

不在其中。

你必须将你的角色随机效应变量转换为因子；您不能将它们保留为字符，否则您会收到

gam()

看到的错误（第一个和第三个错误）。

random

不是

gam()

的论证；你的最后一个例子没有抛出错误是因为

...

位于

gam()

的函数定义中，并且它默默地吃了

random

参数并且

gam

从未使用过它（出于明显的原因）。

我不明白你为什么要转换响应变量；无需在泊松模型中记录变换。

你想要的是：

df <- transform(df, site = factor(site))
mod1 <- mgcv::gam(list(Am06 ~ s(smax) + s(up_365)+ s(do_365) +
                              s(site, bs = "re"),
                            ~ 1),
                  data = df, family = mgcv::ziplss(),
                  method = "REML")

但我怀疑这就是你想要的模型；恒定的零通货膨胀很少有用，因此请考虑使用模型的均值（比率）部分中的部分或全部协变量来对二项式障碍（列表中有

~ 1

的部分）进行建模。

请注意，

ziplss()

适合零膨胀泊松模型作为障碍模型；观察/未观察的二项式障碍和只能生成 >0 整数的截断泊松部分。请参阅

?ziplss

。

最后，仅仅因为数据中有很多零并不意味着它是零膨胀的。零通胀不是你根据数据来检查的，而是根据模型来检查的。因此，我建议您退后一步，将泊松模型和负二项式模型拟合到数据中，然后评估这些模型的零通胀情况。假设 df` 按照上面的代码块进行转换...

mod_p <- gam(Am06 ~ s(smax) + s(up_365)+ s(do_365) +
                    s(site, bs = "re"),
             data = df, family = poisson(), method = "REML")


mod_nb <- gam(Am06 ~ s(smax) + s(up_365)+ s(do_365) +
                     s(site, bs = "re"),
              data = df, family = nb(), method = "REML")

然后检查这些模型的零通胀，也许可以通过根图（我的

gratia

包中的 rootogram() 中有一个 GAM 版本），或者使用 DHARMa 包及其测试。

仅当这些诊断显示出低估零计数的问题时，才继续使用零改变或膨胀模型。