R random.effects GAM 中的隨機效應

說明

為了估計的目的，GAM 的平滑分量可以被視為隨機效應。這意味著可以通過兩種方式將更傳統的隨機效應項合並到 GAM 中。第一種方法將所有平滑轉換為適合通過標準混合建模軟件估計的固定和隨機分量。一旦 GAM 處於這種形式，那麽傳統的隨機效應就很容易添加，並且整個模型被估計為一般的混合模型。 gamm4 包中的gamm 和gamm4 以這種方式運行。

第二種方法以與平滑表示相同的方式表示 GAM 中的傳統隨機效應 - 作為懲罰回歸項。此方法可以與gam通過利用s(...,bs="re")模型中的項：參見smooth.construct.re.smooth.spec，了解完整詳細信息。基本思想是，例如，s(x,z,g,bs="re")生成一個 i.i.d.模型矩陣的高斯隨機效應由下式給出model.matrix(~x:z:g-1)——原則上，這些術語可以采用任意數量的參數。這種簡單的方法足以實現各種常用的隨機效應結構。例如如果g那麽是一個因子s(g,bs="re")為每個級別生成一個隨機係數g，隨機係數全部建模為 i.i.d.普通的。如果g是一個因子並且x是數字，那麽s(x,g,bs="re")產生一個 i.i.d.與響應相關的正態隨機斜率x對於每個級別g.如果h那麽是另一個因子s(h,g,bs="re")產生通常的 i.i.d.普通的g-h相互作用。請注意，基於 Wood (2013)，還對此類項實施了相當有用的零隨機效應近似檢驗。如果精度矩陣在乘法常數內已知，則可以通過以下方式提供xt的論點s.看mgcv smooth.construct.re.smooth.spec了解詳細信息和示例。一些模型需要相同隨機效應的不同水平之間的差異：這些可以按照中所述來實現linear.functional.terms.

或者，但不太直接，可以使用 gam 的 paraPen 參數：請參閱 gam.models 。如果平滑參數估計是通過 ML 或 REML(例如 gam(...,method="REML") )進行的，那麽這種方法是完全傳統的基於似然的隨機效應處理。

gam 對於擬合具有大量隨機效應的模型可能會很慢，因為它沒有利用稀疏性，而稀疏性通常是參數隨機效應的一個特征。它不能用於係數多於數據的模型。然而，當隨機效應的數量適中時，gam 通常比 gamm 或 gamm4 更快、更可靠。

為了方便使用 gam 的隨機效應，gam.vcomp 是一個用於將平滑參數轉換為方差分量的實用程序。如果平滑度估計是通過 ML 或 REML 進行的，它還提供置信區間。

請注意，將隨機效應視為平滑並不能消除與測試方差分量是否等於零相關的常見問題：請參閱 summary.gam 和 anova.gam 。

例子

## see also examples for gam.models, gam.vcomp, gamm
## and smooth.construct.re.smooth.spec

## simple comparison of lme and gam
require(mgcv)
require(nlme)
b0 <- lme(travel~1,data=Rail,~1|Rail,method="REML") 

b <- gam(travel~s(Rail,bs="re"),data=Rail,method="REML")

intervals(b0)
gam.vcomp(b)
anova(b)
plot(b)

## simulate example...
dat <- gamSim(1,n=400,scale=2) ## simulate 4 term additive truth

fac <- sample(1:20,400,replace=TRUE)
b <- rnorm(20)*.5
dat$y <- dat$y + b[fac]
dat$fac <- as.factor(fac)

rm1 <- gam(y ~ s(fac,bs="re")+s(x0)+s(x1)+s(x2)+s(x3),data=dat,method="ML")
gam.vcomp(rm1)

fv0 <- predict(rm1,exclude="s(fac)") ## predictions setting r.e. to 0
fv1 <- predict(rm1) ## predictions setting r.e. to predicted values
## prediction setting r.e. to 0 and not having to provide 'fac'...
pd <- dat; pd$fac <- NULL
fv0 <- predict(rm1,pd,exclude="s(fac)",newdata.guaranteed=TRUE)

## Prediction with levels of fac not in fit data.
## The effect of the new factor levels (or any interaction involving them)
## is set to zero.
xx <- seq(0,1,length=10)
pd <- data.frame(x0=xx,x1=xx,x2=xx,x3=xx,fac=c(1:10,21:30))
fv <- predict(rm1,pd)
pd$fac <- NULL
fv0 <- predict(rm1,pd,exclude="s(fac)",newdata.guaranteed=TRUE)

作者

Simon Wood <simon.wood@r-project.org>

參考

Wood, S.N. (2013) A simple test for random effects in regression models. Biometrika 100:1005-1010

Wood, S.N. (2011) Fast stable restricted maximum likelihood and marginal likelihood estimation of semiparametric generalized linear models. Journal of the Royal Statistical Society (B) 73(1):3-36

Wood, S.N. (2008) Fast stable direct fitting and smoothness selection for generalized additive models. Journal of the Royal Statistical Society (B) 70(3):495-518

Wood, S.N. (2006) Low rank scale invariant tensor product smooths for generalized additive mixed models. Biometrics 62(4):1025-1036

也可以看看

gam.vcomp、gam.models、smooth.terms、smooth.construct.re.smooth.spec、gamm