如何分析纵向计数数据:考虑 GLMM 中的时间自相关?

机器算法验证 r 混合模式 自相关 错误 面板数据
2022-01-30 06:15:26

您好统计大师和 R 编程向导,

我有兴趣将动物捕获建模为环境条件和一年中的某一天的函数。作为另一项研究的一部分,我统计了三年内大约 160 天的捕获次数。在这些日子里,我每天都有温度、降雨量、风速、相对湿度等。因为数据是从相同的 5 个地块重复收集的,所以我使用地块作为随机效应。

我的理解是 nlme 可以很容易地解释残差中的时间自相关,但不能处理像 lme4 这样的非高斯链接函数(它不能处理自相关?)。目前,我认为在日志(计数)上使用 R 中的 nlme 包可能会起作用。所以我现在的解决方案是运行类似的东西:

m1 <- lme(lcount ~ AirT + I(AirT^2) + RainAmt24 + I(RainAmt24^2) + RHpct + windspeed + 
      sin(2*pi/360*DOY) + cos(2*pi/360*DOY), random = ~1|plot, correlation =
      corARMA(p = 1, q = 1, form = ~DOY|plot), data = Data)

其中 DOY = 一年中的某一天。最终模型中可能会有更多的交互,但这是我的总体思路。我也可以尝试用类似的东西进一步模拟方差结构

weights = v1Pow

我不确定是否有更好的方法来处理泊松混合模型回归或任何东西?我刚刚在 Kedem 和 Fokianos 的“时间序列分析的回归模型”的第 4 章中找到了数学讨论。目前它有点超出我的能力,尤其是在应用程序中(用 R 编码)。我还在 Zuur 等人中看到了 MCMC 解决方案。BUGS 语言(使用 winBUGS 或 JAG)的混合效果模型书(第 23 章)。那是我最好的选择吗?R中是否有一个简单的MCMC包可以处理这个问题?我不太熟悉 GAMM 或 GEE 技术,但如果人们认为它们会提供更好的洞察力,我愿意探索这些可能性。我的主要目标是创建一个模型来预测给定环境条件下的动物捕获量。其次,我想解释一下动物在活动方面的反应。

任何关于进行(哲学)的最佳方式、如何在 R 或 BUGS 中编码的想法都将不胜感激。我对 R 和 BUGS(winBUGS)相当陌生,但正在学习。这也是我第一次尝试解决时间自相关问题。

谢谢,丹

1个回答

记录转换您的响应是一种选择,尽管并不理想。通常首选 GLM 框架。如果您不熟悉 GLM,请先查看它们,然后再查看混合模型扩展。对于计数数据,泊松或负二项分布假设可能是合适的。如果方差高于表示过度分散的平均值(https://en.wikipedia.org/wiki/Overdispersion),则表示负二项式。参数估计的解释对两者是等价的。

在我的经验中,R 中存在几个选项,其中 lme4 最常被引用。

#glmer
library(lme4) 
glmer(count ~ AirT + I(AirT^2) + RainAmt24 + I(RainAmt24^2) + RHpct + windspeed + sin(2*pi/360*DOY) + cos(2*pi/360*DOY) + (1|plot), family=poisson, data = Data) 
# use glmer.nb with identical syntax but no family for negative binomial.

# glmmADMB with negative binomial
install.packages("glmmADMB", repos=c("http://glmmadmb.r-forge.r-project.org/repos", getOption("repos")),type="source") 
require(glmmADMB)
glmmadmb(count ~ AirT + I(AirT^2) + RainAmt24 + I(RainAmt24^2) + RHpct + windspeed + sin(2*pi/360*DOY) + cos(2*pi/360*DOY) + (1|plot), 
           family="nbinom", zeroInflation=FALSE, data=Data)

# glmmPQL, requires an estimate for theta which can be obtained from a 
# glm model in which the correlation structure is ignored.
library(MASS)
glmmPQL(count ~ AirT + I(AirT^2) + RainAmt24 + I(RainAmt24^2) + RHpct + windspeed + sin(2*pi/360*DOY) + cos(2*pi/360*DOY) , random = list(~1 | plot), data = Data,family = negative.binomial(theta = 4.22, link = log))

这些链接也可能有帮助:

https://udrive.oit.umass.edu/xythoswfs/webui/_xy-11096203_1-t_yOxYgf1s http://www.cell.com/trends/ecology-evolution/pdf/S0169-5347(09)00019-6.pdf

两者都是由 lme4 的作者 Ben Bolker 编写的。

我没有测试这些例子,但它们应该让你知道从哪里开始。如果您希望验证它们的实施,请提供数据。