在mgcv::gam 中有一种方法可以做到这一点(你的Q2),通过predict.gam 方法和type = "lpmatrix"。
?predict.gam 甚至有一个例子,我在下面复制:
library(mgcv)
n <- 200
sig <- 2
dat <- gamSim(1,n=n,scale=sig)
b <- gam(y ~ s(x0) + s(I(x1^2)) + s(x2) + offset(x3), data = dat)
newd <- data.frame(x0=(0:30)/30, x1=(0:30)/30, x2=(0:30)/30, x3=(0:30)/30)
Xp <- predict(b, newd, type="lpmatrix")
##################################################################
## The following shows how to use use an "lpmatrix" as a lookup
## table for approximate prediction. The idea is to create
## approximate prediction matrix rows by appropriate linear
## interpolation of an existing prediction matrix. The additivity
## of a GAM makes this possible.
## There is no reason to ever do this in R, but the following
## code provides a useful template for predicting from a fitted
## gam *outside* R: all that is needed is the coefficient vector
## and the prediction matrix. Use larger `Xp'/ smaller `dx' and/or
## higher order interpolation for higher accuracy.
###################################################################
xn <- c(.341,.122,.476,.981) ## want prediction at these values
x0 <- 1 ## intercept column
dx <- 1/30 ## covariate spacing in `newd'
for (j in 0:2) { ## loop through smooth terms
cols <- 1+j*9 +1:9 ## relevant cols of Xp
i <- floor(xn[j+1]*30) ## find relevant rows of Xp
w1 <- (xn[j+1]-i*dx)/dx ## interpolation weights
## find approx. predict matrix row portion, by interpolation
x0 <- c(x0,Xp[i+2,cols]*w1 + Xp[i+1,cols]*(1-w1))
}
dim(x0)<-c(1,28)
fv <- x0%*%coef(b) + xn[4];fv ## evaluate and add offset
se <- sqrt(x0%*%b$Vp%*%t(x0));se ## get standard error
## compare to normal prediction
predict(b,newdata=data.frame(x0=xn[1],x1=xn[2],
x2=xn[3],x3=xn[4]),se=TRUE)
这会贯穿整个过程,甚至是在 R 或 GAM 模型之外完成的预测步骤。您将不得不稍微修改示例以执行您想要的操作,因为示例评估模型中的所有项,并且除了样条之外您还有另外两个项 - 本质上您执行相同的操作,但仅针对样条项,即涉及为样条查找Xp 矩阵的相关列和行。然后你还应该注意样条曲线居中,所以你可能也可能不想撤消它。
对于您的 Q1,为示例中的 xn 向量/矩阵选择适当的值。这些对应于模型中nth 项的值。因此,将您想要固定的值设置为某个平均值,然后改变与样条关联的值。
如果您在 R 中执行所有这些操作,则仅在您拥有数据的样条协变量值处评估样条曲线会更容易进入另一个模型。您可以通过创建一个数据框来进行预测,然后使用
predict(mod, newdata = newdat, type = "terms")
其中mod 是拟合的 GAM 模型(通过 mgcv::gam),newdat 是包含模型中每个变量的列的数据框(包括参数项;设置您不想要的项变化到某个恒定的平均值[比如数据集中变量的平均值]或某个水平(如果是一个因素)。 type = "terms" 部分将为newdat 中的每一行返回一个矩阵,其中包含模型中每个项的拟合值的“贡献”,包括样条项。只需取该矩阵中与样条线对应的列 - 再次居中。
也许我误解了您的 Q1。如果要控制结,请参阅mgcv::gam 的 knots 参数。默认情况下,mgcv::gam 在数据的极端处放置一个结,然后剩余的“结”均匀分布在该区间上。 mgcv::gam 没有找到结 - 它会为您放置它们,您可以通过 knots 参数控制它们放置的位置。