如何基于GAM拟合结果的估计公式手动计算预测值？

背景与示例代码 #

以下是使用R语言mgcv包拟合广义可加模型（GAM）并获取预测值的代码：

library('mgcv')

a0 = 1
a1 = 5

f2 = function(x) x^3

ysim1 = function(n = 500) {
  set.seed(10)
  x1 = runif(n)
  x2 = runif(n)
  e = rnorm(n)
  f = f2(x1)
  y = a0 + a1*x1 + f + e
  data.frame(y = y, x1, x2, f2 = f)
}

df1 = ysim1()
head(df1)

m = gam(y ~ x1 + s(x2), data = df1)
pred1 = predict(m)[1]

运行后得到第一个样本的预测值：

> pred1
       1 
3.846924

对比线性回归场景，我们可以直接提取系数结合自变量计算预测值：

ysim2 = function(n = 500) {
 set.seed(10)
 x = runif(n)
 e = rnorm(n)
 y = a0 + a1*x + e
 data.frame(y, x)
}

df2 = ysim2()
head(df2)

m2 = lm(y~x, df2)

pred2 = predict(m2)[1]

预测值：

> pred2
     1 
3.5459 

手动计算方式：coef(m2)[[1]] + coef(m2)[[2]]*df2$x[1]

GAM模型的手动计算方法 #

GAM包含平滑项，无法像线性回归那样直接用简单线性组合计算，但可以通过**线性预测器矩阵（lpmatrix）**实现手动计算：

1. 提取模型系数：

gam_coef <- coef(m)

2. 获取第一个样本对应的线性预测器矩阵：
   这个矩阵包含线性项x1和平滑项s(x2)的基函数取值，通过predict函数的type = "lpmatrix"参数生成：

lp_matrix <- predict(m, newdata = df1[1, ], type = "lpmatrix")

3. 计算预测值：
   将矩阵与系数做矩阵乘法，得到最终预测值：

manual_pred <- as.vector(lp_matrix %*% gam_coef)

执行后结果与pred1完全一致：

> manual_pred
[1] 3.846924

原理说明 #

GAM的预测本质是计算线性预测器，即模型系数与设计矩阵的乘积。其中平滑项s(x2)由mgcv自动转换为一组基函数，设计矩阵中会包含这些基函数的列。通过lpmatrix获取对应样本的设计矩阵后，即可像线性回归一样完成手动计算。

内容的提问来源于stack exchange，提问作者M. Er