r - 用数据在同一图中绘制 GAM 预测值和 GAMM AR1 时出现问题

Question

我已经用我的数据运行了 GAM，并且我正在将 GAM 的预测值与数据点一起绘制在图表中。不同区域有 15 个相同的图，其中一些存在自相关问题。对于这些，我运行了一个 GAMM AR1 模型，现在我想绘制与其他区域相似的预测值。在下面你会看到两张图，左边的一张是来自 GAM 的预测值和置信区间，以及真实数据。在右侧，您将看到来自 GAMM AR1 的线。

如您所见，GAM 图具有预测值、CI 和带有数据点的“真实数据”x 轴。GAMM AR1 有一条蓝线，“GAMM 值”位于 x 轴上。

我如何绘制来自 GAMM AR1 的预测值，类似于我对 GAM 所做的？请参阅下面的数据和脚本。

数据（数据框“例如”）：

structure(list(Year = c(1970, 1971, 1972, 1973, 1974, 1975, 1976, 
1977, 1978, 1979, 1980, 1981, 1982, 1983, 1984, 1985, 1986, 1987, 
1988, 1989, 1990, 1991, 1992, 1993, 1994, 1995, 1996, 1997, 1998, 
1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, 2009, 
2010, 2011, 2012, 2013, 2014, 2015, 2016), F = c(0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 16, 14, 24, 10, 15, 26, 20, 
15, 19, 13, 18), M = c(0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 
0, 0, 27, 40, 35, 39, 19, 30, 42, 42, 39, 56, 50), U = c(100, 
79, 71, 87, 119, 56, 98, 78, 50, 58, 71, 131, 159, 89, 100, 43, 
28, 89, 108, 95, 110, 131, 114, 45, 49, 56, 52, 51, 69, 81, 85, 
60, 54, 46, 54, 57, 1, 5, 8, 5, 0, 1, 1, 5, 8, 2, 0), Tot = c(100, 
79, 71, 87, 119, 56, 98, 78, 50, 58, 71, 131, 159, 89, 100, 43, 
28, 89, 108, 95, 110, 131, 114, 45, 49, 56, 52, 51, 69, 81, 85, 
60, 54, 46, 54, 57, 44, 59, 67, 54, 34, 57, 63, 62, 66, 71, 68
), ratio = c(NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, 1.6875, 2.85714285714286, 1.45833333333333, 
3.9, 1.26666666666667, 1.15384615384615, 2.1, 2.8, 2.05263157894737, 
4.30769230769231, 2.77777777777778), popsize = c(NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA)), row.names = c(NA, 
-47L), class = c("tbl_df", "tbl", "data.frame"))

GAM 脚本和使用 CI 和数据绘制 pred 值：

m.eg <- gam(Tot~s(Year),family=poisson,data=eg)

YearP=seq(1970,2016,by=1)
meg.pred=predict(m.eg,newdata=data.frame(Year=YearP),type="response",se.fit=T)

plotCI(x=YearP, y=meg.pred$fit,uiw=2*meg.pred$se.fit, type="l",sfrac=0.003,
       ylim=c(40,140),xlim=c(1970,2016),
       col="red",gap=0,lwd=1.6,cex=1.2,las=1, 
       xlab="", ylab="Number of animals")
points(eg$Year,eg$Tot,pch=19,cex=0.9)

GAMM AR1 脚本和线图：

m.eg.ar1 <- gamm(Tot~s(Year),family=poisson,correlation = corAR1(form = ~ Year), data=eg)
plot(eg$Year,predict(m.eg),col="white")
lines(eg$Year,predict(m.eg.ar1$gam),col="blue")

标签： rplotgam

Answer 1

预测值gamm几乎与 for 相同，gam因此您几乎拥有它。唯一的区别是一个gamm对象由一个gam和一个对象组成，lme并且该predict方法仅gam作为参数。所以你必须在predict调用中指定这个。下面的代码适用于您的示例：

m.eg.ar1 <- gamm(Tot~s(Year),family=poisson,correlation = corAR1(form = ~ Year), data=eg)
megar1.pred=predict(m.eg.ar1$gam,newdata=data.frame(Year=YearP),type="response",se.fit=T)
plotCI(x=YearP, y=megar1.pred$fit,uiw=2*megar1.pred$se.fit, type="l",sfrac=0.003,
       ylim=c(40,140),xlim=c(1970,2016),
       col="red",gap=0,lwd=1.6,cex=1.2,las=1, 
       xlab="", ylab="Number of animals")
points(eg$Year,eg$Tot,pch=19,cex=0.9)