mgcv:如何识别 gam 和 gamm 模型中的确切节点值?
mgcv: How to identify exact knot values in a gam and gamm model?
我正在使用 gams 来拟合资源选择函数,以确定候鸟的能量开发阈值。我的模型如下所示:
m4 <-gam(used~ti(propwells_buff_res1500, bs = "cr", k = 5) +
ti(year, bs = "cr", k = 5) +
ti(propwells_buff_res1500, year, bs = "cr", k = 5),
family = binomial(link = "cloglog"), data=mov, gamma=1.4, method="ML")
used 是动物使用的位置,propwells_buff_res1500 是随机生成的 "available" 点(由 1500 米半径的圆圈缓冲),它们内部的能量发展量不同。我已将结限制为 5,但是,我希望能够提取精确的结值,因为据我所知,结值代表阈值……也就是动物使用下降的表面干扰百分比。
我希望这是有道理的。如果没有,我只想知道如何获得结值。从 plot(m4) 中,我可以看到非线性线的斜率开始发生变化的位置,但了解确切的值将非常有帮助。
到目前为止,我已经尝试过:
smooth <- m4$smooth[[3]]
smooth$knots
##this knot option isn't available to me,
##I saw it in an old post from 2016, figured out that XP should replace knots
smooth$XP
##and all this returns is list()
非常感谢任何帮助,谢谢。
为了得到结,您可以提取边际平滑项的 xp
分量(注意它是小写的 xp
因为有一个 XP
在顶层平滑,这是另外一回事)。
这是一个例子
library('mgcv')
## simulate some data
set.seed(1729)
df <- gamSim(2) # this is a bivariate example
## fit the model
mod <- gam(y ~ ti(x, bs = 'cr', k = 5) +
ti(z, bs = 'cr', k = 5) +
ti(x, z, bs = rep('cr', 2), k = 5),
data = df$data, method = 'REML')
## extract the 3rd smooth
sm <- mod[['smooth']][[3]]
边缘碱基在sm$margin
,这只是两个平滑对象的列表:
r$> str(sm$margin, max = 1)
List of 2
$ :List of 21
..- attr(*, "class")= chr [1:2] "cr.smooth" "mgcv.smooth"
..- attr(*, "qrc")=List of 4
.. ..- attr(*, "class")= chr "qr"
..- attr(*, "nCons")= int 1
$ :List of 21
..- attr(*, "class")= chr [1:2] "cr.smooth" "mgcv.smooth"
..- attr(*, "qrc")=List of 4
.. ..- attr(*, "class")= chr "qr"
..- attr(*, "nCons")= int 1
其中每个都有一个 xp
组件:
sm_x <- sm$margin[[1]]
sm_z <- sm$margin[[2]]
因此 x
的边缘 CRS 的节点是:
r$> sm_x$xp
0% 25% 50% 75% 100%
0.0005697084 0.2477067126 0.4704501621 0.7121602102 0.9960833385
并且 z
是
r$> sm_z$xp
0% 25% 50% 75% 100%
0.007381999 0.244705125 0.488819070 0.717802322 0.991505836
为什么要这些值?它们位于观察到的协变量值的五分位数处:
r$> with(df$data, quantile(x, probs = seq(0, 1, length = 5)))
0% 25% 50% 75% 100%
0.0005697084 0.2477067126 0.4704501621 0.7121602102 0.9960833385
r$> with(df$data, quantile(z, probs = seq(0, 1, length = 5)))
0% 25% 50% 75% 100%
0.007381999 0.244705125 0.488819070 0.717802322 0.991505836
这就是 mgcv 为 CRS 基础放置结的方式。可以使用 place.knots()
:
恢复确切位置
r$> with(df$data, place.knots(x, 5))
[1] 0.0005697084 0.2477067126 0.4704501621 0.7121602102 0.9960833385
r$> with(df$data, place.knots(z, 5))
[1] 0.007381999 0.244705125 0.488819070 0.717802322 0.991505836
但是从边缘平滑对象拉结更安全,因为用户总是可以通过 knots
参数指定结给 gam()
。
我正在使用 gams 来拟合资源选择函数,以确定候鸟的能量开发阈值。我的模型如下所示:
m4 <-gam(used~ti(propwells_buff_res1500, bs = "cr", k = 5) +
ti(year, bs = "cr", k = 5) +
ti(propwells_buff_res1500, year, bs = "cr", k = 5),
family = binomial(link = "cloglog"), data=mov, gamma=1.4, method="ML")
used 是动物使用的位置,propwells_buff_res1500 是随机生成的 "available" 点(由 1500 米半径的圆圈缓冲),它们内部的能量发展量不同。我已将结限制为 5,但是,我希望能够提取精确的结值,因为据我所知,结值代表阈值……也就是动物使用下降的表面干扰百分比。
我希望这是有道理的。如果没有,我只想知道如何获得结值。从 plot(m4) 中,我可以看到非线性线的斜率开始发生变化的位置,但了解确切的值将非常有帮助。
到目前为止,我已经尝试过:
smooth <- m4$smooth[[3]]
smooth$knots
##this knot option isn't available to me,
##I saw it in an old post from 2016, figured out that XP should replace knots
smooth$XP
##and all this returns is list()
非常感谢任何帮助,谢谢。
为了得到结,您可以提取边际平滑项的 xp
分量(注意它是小写的 xp
因为有一个 XP
在顶层平滑,这是另外一回事)。
这是一个例子
library('mgcv')
## simulate some data
set.seed(1729)
df <- gamSim(2) # this is a bivariate example
## fit the model
mod <- gam(y ~ ti(x, bs = 'cr', k = 5) +
ti(z, bs = 'cr', k = 5) +
ti(x, z, bs = rep('cr', 2), k = 5),
data = df$data, method = 'REML')
## extract the 3rd smooth
sm <- mod[['smooth']][[3]]
边缘碱基在sm$margin
,这只是两个平滑对象的列表:
r$> str(sm$margin, max = 1)
List of 2
$ :List of 21
..- attr(*, "class")= chr [1:2] "cr.smooth" "mgcv.smooth"
..- attr(*, "qrc")=List of 4
.. ..- attr(*, "class")= chr "qr"
..- attr(*, "nCons")= int 1
$ :List of 21
..- attr(*, "class")= chr [1:2] "cr.smooth" "mgcv.smooth"
..- attr(*, "qrc")=List of 4
.. ..- attr(*, "class")= chr "qr"
..- attr(*, "nCons")= int 1
其中每个都有一个 xp
组件:
sm_x <- sm$margin[[1]]
sm_z <- sm$margin[[2]]
因此 x
的边缘 CRS 的节点是:
r$> sm_x$xp
0% 25% 50% 75% 100%
0.0005697084 0.2477067126 0.4704501621 0.7121602102 0.9960833385
并且 z
是
r$> sm_z$xp
0% 25% 50% 75% 100%
0.007381999 0.244705125 0.488819070 0.717802322 0.991505836
为什么要这些值?它们位于观察到的协变量值的五分位数处:
r$> with(df$data, quantile(x, probs = seq(0, 1, length = 5)))
0% 25% 50% 75% 100%
0.0005697084 0.2477067126 0.4704501621 0.7121602102 0.9960833385
r$> with(df$data, quantile(z, probs = seq(0, 1, length = 5)))
0% 25% 50% 75% 100%
0.007381999 0.244705125 0.488819070 0.717802322 0.991505836
这就是 mgcv 为 CRS 基础放置结的方式。可以使用 place.knots()
:
r$> with(df$data, place.knots(x, 5))
[1] 0.0005697084 0.2477067126 0.4704501621 0.7121602102 0.9960833385
r$> with(df$data, place.knots(z, 5))
[1] 0.007381999 0.244705125 0.488819070 0.717802322 0.991505836
但是从边缘平滑对象拉结更安全,因为用户总是可以通过 knots
参数指定结给 gam()
。