在 Gamma GLMM 中处理空间自相关

Question

问题概述

我正在使用一个数据集来检查植物生长与 各种因素，包括 14 个不同的海冰范围 (SeaIce) 采样点。为了分析数据，我使用伽玛广义 具有对数链接函数的线性混合模型 (GLMM)。然而，我是 关心模型残差和 I 中的空间自相关 我不确定解决此问题的最佳方法。

为了解决这个问题，我正在利用

DHARMa

包装连同

lme4

对于 模型拟合。

library(DHARMa)
library(lme4)
library(MASS)
library(gstat)
library(dplyr)
library(sf)
library(sp)

您可以在此 GitHub 中找到该模型的数据集 存储库。

FinalDataset <- read.csv("https://raw.githubusercontent.com/derek-corcoran-barrios/SeaIceQuestion/master/FinalDataset.csv")

当前模型结构

响应变量

• Growth：代表每个个体每年的增长增量 植物。

预测因素

• SeaIce.s： 每个地点每年的海冰范围，是 利息（按比例计算）。
• 年龄： 每个测量年份的每株植物的年龄，包括 由于生物关系而与 SeaIce.s 进行交互 植物年龄和生长之间。

随机效应

• 地点：代表 14 个采样点的因子，具有不同的 他们之间的距离。我包括一个随机斜率

  (SeaIce.s | Site)

  考虑每个增长与海冰关系的不同斜率 网站。
• 年份： 该研究涵盖 1983 年至 2015 年，并包括 作为随机效应

  (1 | year)

  。
• 个体： 代表测量的每株植物的因子。

型号

fullmod <- glmer(Growth ~ SeaIce.s * age + (SeaIce.s | Site) + (1 | year) + (1 | Individual), data = FinalDataset, family = Gamma(link = "log"), control = glmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e5)))

主要关注点

我主要关心的是如何正确考虑空间自相关 在这个模型内。我发现使用这个的时候有问题 dHARMA 中的代码。

res2_null <- simulateResiduals(fullmod)
res3_null <- recalculateResiduals(res2_null, group = FinalDataset$Site)
locs_null <- FinalDataset %>% group_by(Site) %>% summarise(across(c(Latitude, Longitude), mean))
testSpatialAutocorrelation(res3_null, x = locs_null$Longitude, y = locs_null$Latitude)

## 
##  DHARMa Moran's I test for distance-based autocorrelation
## 
## data:  res3_null
## observed = 0.343448, expected = -0.076923, sd = 0.150983, p-value =
## 0.005365
## alternative hypothesis: Distance-based autocorrelation

从这个结果可以看出，存在很高的空间自相关性。

我尝试过的一些事情：

尝试的解决方案1

我尝试过的解决方案之一这个 问题 将纬度和经度作为固定效应包括在内。然而，当我测试时 使用 testSpatialAutocorrelation() 中的空间自相关 DHARMa 包，我仍然得到一个显着的莫兰 I 模型 包括纬度和经度。

LonLatmod <- glmer(Growth ~ SeaIce.s * age + Latitude + Longitude + (SeaIce.s | Site) + (1 | year) + (1 | Individual), data = FinalDataset, family = Gamma(link = "log"), control = glmerControl(optimizer = "bobyqa", optCtrl = list(maxfun = 2e5)))

res <-  simulateResiduals(LonLatmod)
res2 <- recalculateResiduals(res, group = FinalDataset$Site)
locs_latlon <- FinalDataset %>% group_by(Site) %>% summarise(across(c(Latitude, Longitude), mean))
testSpatialAutocorrelation(res2, x =  locs_latlon$Longitude, y = locs_latlon$Latitude)

## 
##  DHARMa Moran's I test for distance-based autocorrelation
## 
## data:  res2
## observed = 0.382948, expected = -0.076923, sd = 0.150649, p-value =
## 0.002269
## alternative hypothesis: Distance-based autocorrelation

尝试的解决方案2

我还尝试在 MASS 包中使用 glmmPQL() 。但我还没有去过 由于各种错误，即使没有，也能够使其工作 添加关联结构。下面的例子：

formula_glmmPQL <- as.formula("Growth ~ SeaIce.s * age")
model_glmmPQL <- glmmPQL(formula_glmmPQL,
                         random = list(~ SeaIce.s|Site, 
                                      ~ 1|year, 
                                      ~1| Individual),
                         data = FinalDataset,
                         na.action=na.omit,
                         family = Gamma(link = "log"))

#Error in pdFactor.pdLogChol(X[[i]], ...) :
#NA/NaN/Inf in foreign function call (arg 3)

尝试的解决方案3

我也尝试过在 gstat 包中使用 variogram() ，但我 不确定如何解释变异函数的形状以及要改变什么 我的模型基于此。

FinalDatasetSF <- st_as_sf(FinalDataset, coords = c("Longitude", "Latitude"), crs = st_crs(4326))

null_mod <- variogram(log(Growth) ~ 1, FinalDatasetSF)

Abn_fit_null <- fit.variogram(null_mod, model = vgm(1, "Sph", 
                                              700, 1))
plot(null_mod, model=Abn_fit_null)

我还探索了 nlme 包中的 corAR1() 函数，但它 似乎与 glmer() 模型不兼容。

编辑

我添加了 @SarahS 的解决方案，有或没有纬度和 经度作为固定效应

require(glmmTMB)
# Set up the necessary variables
FinalDataset$pos <- numFactor(FinalDataset$Latitude, FinalDataset$Longitude)
FinalDataset$group <- factor(rep(1, nrow(FinalDataset)))
# Fit model
TestA <- glmmTMB(Growth ~ SeaIce.s * age + Latitude + Longitude + (SeaIce.s | Site) + (1 | year) + (1 | Individual) + 1 + exp(pos + 0 | group), data = FinalDataset, family = Gamma(link = "log"))


TestB <- glmmTMB(Growth ~ SeaIce.s * age + Latitude + Longitude + (SeaIce.s | Site) + (1 | year) + (1 | Individual) + 1 + exp(pos + 0 | group), data = FinalDataset, family = Gamma(link = "log"))

然而，这些都不能解决空间自相关问题，如下所示

res <-  simulateResiduals(TestA)
res2 <- recalculateResiduals(res, group = FinalDataset$Site)

locs_latlon <- FinalDataset %>% group_by(Site) %>% summarise(across(c(Latitude, Longitude), mean))

testSpatialAutocorrelation(res2, x =  locs_latlon$Longitude, y = locs_latlon$Latitude)

## 
##  DHARMa Moran's I test for distance-based autocorrelation
## 
## data:  res2
## observed = 0.409445, expected = -0.076923, sd = 0.151346, p-value =
## 0.001311
## alternative hypothesis: Distance-based autocorrelation

还有这里

res <-  simulateResiduals(TestB)
res2 <- recalculateResiduals(res, group = FinalDataset$Site)

testSpatialAutocorrelation(res2, x =  locs_latlon$Longitude, y = locs_latlon$Latitude)

## 
##  DHARMa Moran's I test for distance-based autocorrelation
## 
## data:  res2
## observed = 0.409445, expected = -0.076923, sd = 0.151346, p-value =
## 0.001311
## alternative hypothesis: Distance-based autocorrelation

征求建议

我将非常感谢任何有关如何的建议或见解 有效解决我当前建模中的空间自相关问题 接近。

Answer 1

您尝试过使用glmmTMB中的空间自相关调整吗？

我相信代码是

install.packages("glmmTMB"); require(glmmTMB)
# Set up the necessary variables
FinalDataset$pos <- numFactor(FinalDataset$Latitude, FinalDataset$Longitude)
FinalDataset$group <- factor(rep(1, nrow(FinalDataset)))
# Fit model
glmmTMB(Growth ~ SeaIce.s * age + Latitude + Longitude + (SeaIce.s | Site) + (1 | year) + (1 | Individual) + 1 + exp(pos + 0 | group), data = FinalDataset, family = Gamma(link = "log"))