R/mgcv:为什么 te() 和 ti() 张量积会产生不同的表面?

机器算法验证 r 广义加法模型 毫克CV
2022-02-28 19:05:03

用于拟合张量积交互mgcv包有两个函数:. 我了解两者之间的基本分工(拟合非线性交互与将此交互分解为主要影响和交互)。我不明白的是为什么并且可能会产生(略微)不同的结果。Rte()ti()te(x1, x2)ti(x1) + ti(x2) + ti(x1, x2)

MWE(改编自?ti):

require(mgcv)
test1 <- function(x,z,sx=0.3,sz=0.4) { 
  x <- x*20
 (pi**sx*sz)*(1.2*exp(-(x-0.2)^2/sx^2-(z-0.3)^2/sz^2)+
             0.8*exp(-(x-0.7)^2/sx^2-(z-0.8)^2/sz^2))
}
n <- 500

x <- runif(n)/20;z <- runif(n);
xs <- seq(0,1,length=30)/20;zs <- seq(0,1,length=30)
pr <- data.frame(x=rep(xs,30),z=rep(zs,rep(30,30)))
truth <- matrix(test1(pr$x,pr$z),30,30)
f <- test1(x,z)
y <- f + rnorm(n)*0.2

par(mfrow = c(2,2))

# Model with te()
b2 <- gam(y~te(x,z))
vis.gam(b2, plot.type = "contour", color = "terrain", main = "tensor product")

# Model with ti(a) + ti(b) + ti(a,b)
b3 <- gam(y~ ti(x) + ti(z) + ti(x,z))
vis.gam(b3, plot.type = "contour", color = "terrain", main = "tensor anova")

# Scatterplot of prediction b2/b3
plot(predict(b2), predict(b3))

在这个例子中差异不是很大,但我只是想知道为什么会有差异。

会话信息:

 > devtools::session_info("mgcv")
 Session info
 -----------------------------------------------------------------------------------
 setting  value                       
 version  R version 3.3.1 (2016-06-21)
 system   x86_64, linux-gnu           
 ui       RStudio (0.99.491)          
 language en_US                       
 collate  en_US.UTF-8                 
 tz       <NA>                        
 date     2016-09-13                  

 Packages      ---------------------------------------------------------------------------------------
 package * version date       source        
 lattice   0.20-33 2015-07-14 CRAN (R 3.2.1)
 Matrix    1.2-6   2016-05-02 CRAN (R 3.3.0)
 mgcv    * 1.8-12  2016-03-03 CRAN (R 3.2.3)
 nlme    * 3.1-128 2016-05-10 CRAN (R 3.3.1)
1个回答

这些表面上是相同的模型,但实际上在拟合时存在一些细微的差异。一个重要的区别是,与模型ti()相比,带有项的模型估计了更多的平滑度参数te()

> b2$sp
te(x,z)1 te(x,z)2 
3.479997 5.884272 
> b3$sp
    ti(x)     ti(z)  ti(x,z)1  ti(x,z)2 
 8.168742 60.456559  2.370604  2.761823

这是因为有更多的惩罚矩阵与这两个模型相关联;ti()模型中,我们每个“项”有一个,而te()模型中只有两个,每个边缘平滑一个。

我看到带有ti()as 的模型用于决定我是否想要y^=β0+s(x,y)或者y^=β0+s(x)+s(y). te()如果我使用术语,我无法比较这些模型,所以我使用ti(). 一旦我确定我是否需要s(x,y)如果我需要,我可以重新调整模型,如果不需要,我可以为每个边际效应te()单独调整s()s(x,y).

method = "ML"请注意,您可以通过使用( 或进行拟合来使模型彼此更接近"REML",但除非所有术语都被完全惩罚,否则您不应该将“固定”效果与它们进行比较"REML",默认情况下它们不是,但会说与select = TRUE)。