从具有随机截距和二次项的 GLMM 生成空间预测(栅格)
Generating a spatial prediction (raster) from a GLMM with a random intercept and a quadratic term
如标题所示,我正在尝试生成描述相对使用概率的预测栅格。为了创建一个可重现的示例,我使用了 maxlike
包中的 MaungaWhau
数据。 MaungaWhau
是一个 list
,它包含两个栅格图层以及一个包含 1000 个位置的图层。所以,有了这些数据和包...
library(maxlike)
library(lme4)
data(MaungaWhau)
我们可以制作两个栅格图层,一个栅格堆栈,以及一个SpatialPoints
object.
elev <- raster(MaungaWhau$elev, xmn=0, xmx=61, ymn=0, ymx=87)
precip <- raster(MaungaWhau$precip, xmn=0, xmx=61, ymn=0, ymx=87)
rs <- stack(elev, precip)
PointDat <- SpatialPoints(MaungaWhau$xy)
然后我制作了一个包含 IndID
的新数据框:每个人的唯一 ID (AAA - DDD); Used
:表示点是否被使用的二进制响应变量(分别为1或0);以及 SpatialPoints
object.
中每个点的 Elev
和 Precip
值
df <- data.frame(IndID = sample(c("AAA", "BBB", "CCC", "DDD"), 1000, replace = T),
Used = sample(c(0,1),1000, replace = T),
Elev = extract(elev, PointDat),
Precip = extract(precip, PointDat))
head(df); tail(df)
> head(df); tail(df)
IndID Used Elev Precip
1 DDD 0 0.3798393 0.6405494
2 DDD 1 0.8830846 1.1174869
3 AAA 0 1.9282864 0.9641432
4 DDD 0 1.5024634 0.4695881
5 BBB 1 1.3089075 -0.1341483
6 BBB 1 0.5733952 0.6246699
然后我构建了一个资源选择模型 (RSF) 并将 IndID 指定为随机效应。另请注意,我为 Elev
.
添加了一个二次项
#Make model
Mod <- glmer(Used ~ Elev + I(Elev^2) + Precip + (1|IndID), family=binomial, data = df, verbos = 1)
summary(Mod)
鉴于已经使用过的和可用的点,我对解释不感兴趣。我的第一个问题是更多的确认。 raster
包小插图指出 "The names in the Raster object should exactly match those expected by the model." 在与 I(Elev^2)
拟合的二次项的实例中,我是否正确预测 Elev
将是 'looking'?似乎是这种情况,因为在下面的 predict
代码中没有与 Elev
相关的错误。
其次,如何处理随机截距项(1|IndID)
?我对边际预测感兴趣,不需要考虑个人。
用下面的代码
#Change names of layers in raster stack to match model
names(rs) <- c("Elev", "Precip")
Pred <- predict(rs, Mod)
我收到错误:
Error in eval(expr, envir, enclos) : object 'IndID' not found
是否可以在不将 IndID
协变量传递给 predict
函数的情况下为 'typical' 个体生成边际预测?换句话说,我想在制作预测表面时忽略 IndID
项和对截距的相关单独调整。
与其依赖 predict
函数,不如通过首先创建固定效应的 beta 对象来手动生成预测
betas <- fixef(Mod)
然后通过将每个栅格乘以各自的 beta 系数来生成预测。
Pred <- betas[1] + (elev * betas[2]) + (elev^2 * betas[3]) + (precip * betas[4])
plot(Pred)
然后很容易添加或删除截距并手动指定 link 函数(例如 logit)。
lme4 (merMod) 对象的 predict
函数默认进行条件预测。
要进行 marginal/unconditional 预测,您需要使用 re.form
参数。您的代码如下所示:
Pred <- predict(rs, Mod, re.form = NA)
如果您还想在响应变量的范围内进行预测,您可以使用 type
参数。有关可用参数的更多详细信息,请参阅帮助页面,?predict.merMod
。
如标题所示,我正在尝试生成描述相对使用概率的预测栅格。为了创建一个可重现的示例,我使用了 maxlike
包中的 MaungaWhau
数据。 MaungaWhau
是一个 list
,它包含两个栅格图层以及一个包含 1000 个位置的图层。所以,有了这些数据和包...
library(maxlike)
library(lme4)
data(MaungaWhau)
我们可以制作两个栅格图层,一个栅格堆栈,以及一个SpatialPoints
object.
elev <- raster(MaungaWhau$elev, xmn=0, xmx=61, ymn=0, ymx=87)
precip <- raster(MaungaWhau$precip, xmn=0, xmx=61, ymn=0, ymx=87)
rs <- stack(elev, precip)
PointDat <- SpatialPoints(MaungaWhau$xy)
然后我制作了一个包含 IndID
的新数据框:每个人的唯一 ID (AAA - DDD); Used
:表示点是否被使用的二进制响应变量(分别为1或0);以及 SpatialPoints
object.
Elev
和 Precip
值
df <- data.frame(IndID = sample(c("AAA", "BBB", "CCC", "DDD"), 1000, replace = T),
Used = sample(c(0,1),1000, replace = T),
Elev = extract(elev, PointDat),
Precip = extract(precip, PointDat))
head(df); tail(df)
> head(df); tail(df)
IndID Used Elev Precip
1 DDD 0 0.3798393 0.6405494
2 DDD 1 0.8830846 1.1174869
3 AAA 0 1.9282864 0.9641432
4 DDD 0 1.5024634 0.4695881
5 BBB 1 1.3089075 -0.1341483
6 BBB 1 0.5733952 0.6246699
然后我构建了一个资源选择模型 (RSF) 并将 IndID 指定为随机效应。另请注意,我为 Elev
.
#Make model
Mod <- glmer(Used ~ Elev + I(Elev^2) + Precip + (1|IndID), family=binomial, data = df, verbos = 1)
summary(Mod)
鉴于已经使用过的和可用的点,我对解释不感兴趣。我的第一个问题是更多的确认。 raster
包小插图指出 "The names in the Raster object should exactly match those expected by the model." 在与 I(Elev^2)
拟合的二次项的实例中,我是否正确预测 Elev
将是 'looking'?似乎是这种情况,因为在下面的 predict
代码中没有与 Elev
相关的错误。
其次,如何处理随机截距项(1|IndID)
?我对边际预测感兴趣,不需要考虑个人。
用下面的代码
#Change names of layers in raster stack to match model
names(rs) <- c("Elev", "Precip")
Pred <- predict(rs, Mod)
我收到错误:
Error in eval(expr, envir, enclos) : object 'IndID' not found
是否可以在不将 IndID
协变量传递给 predict
函数的情况下为 'typical' 个体生成边际预测?换句话说,我想在制作预测表面时忽略 IndID
项和对截距的相关单独调整。
与其依赖 predict
函数,不如通过首先创建固定效应的 beta 对象来手动生成预测
betas <- fixef(Mod)
然后通过将每个栅格乘以各自的 beta 系数来生成预测。
Pred <- betas[1] + (elev * betas[2]) + (elev^2 * betas[3]) + (precip * betas[4])
plot(Pred)
然后很容易添加或删除截距并手动指定 link 函数(例如 logit)。
lme4 (merMod) 对象的 predict
函数默认进行条件预测。
要进行 marginal/unconditional 预测,您需要使用 re.form
参数。您的代码如下所示:
Pred <- predict(rs, Mod, re.form = NA)
如果您还想在响应变量的范围内进行预测,您可以使用 type
参数。有关可用参数的更多详细信息,请参阅帮助页面,?predict.merMod
。