※上記の広告は60日以上更新のないWIKIに表示されています。更新することで広告が下部へ移動します。

# GLMM
library(lme4)
glmres7 <- glmer(y.obs ~ year + SST + I(SST^2) + ( 1 | station), data=data2,
                              family=poisson)
# GAM
library(mgcv)
glmres8 <- gam(y.obs ~ year + s(SST) + station, data=data2,
                              family=poisson)
plot(glmres8)

#-------------------------------------------------------

##重要 ndataの作成 説明変数だけのデータ
ndata <- expand.grid(year=unique(data2$year),station=unique(data2$station),SST=mean(data2$SST))
#uniqueは要素を全部、この場合はyearとstationの全ての組み合わせ、SSTは平均を入れている



# glmの場合
p6 <- predict(glmres6,newdata=ndata,type="response")
cbind(p6, ndata)
p6 <- tapply(p6,ndata$year,mean)
plot(p6/mean(p6),type="b",ylim=c(0,1.5))

# GLMM
p7 <- exp(c(1, 1+ as.numeric(coef(glmres7)$station[1,2:6])))
points(p7/mean(p7),type="b",col=2,pch=2)

# GAM
p8 <- predict(glmres8,newdata=ndata,type="response")
p8 <- tapply(p8,ndata$year,mean)
points(p8/mean(p8),type="b",ylim=c(0,1.5),col=3)
|