14

我经常使用 B 样条进行回归。到目前为止，我从来不需要bs详细了解输出：我只需选择我感兴趣的模型，并将其与lm. 但是，我现在需要在外部（非 R）代码中重现 b 样条模型。那么，由 生成的矩阵是什么意思bs呢？例子：

x <- c(0.0, 11.0, 17.9, 49.3, 77.4)
bs(x, df = 3, degree = 1) # generate degree 1 (linear) B-splines with 2 internal knots
#              1         2         3
# [1,] 0.0000000 0.0000000 0.0000000    
# [2,] 0.8270677 0.0000000 0.0000000    
# [3,] 0.8198433 0.1801567 0.0000000    
# [4,] 0.0000000 0.7286085 0.2713915    
# [5,] 0.0000000 0.0000000 1.0000000   
# attr(,"degree")
# [1] 1
# attr(,"knots")
# 33.33333% 66.66667% 
#  13.30000  38.83333 
# attr(,"Boundary.knots")
# [1]  0.0 77.4
# attr(,"intercept")
# [1] FALSE
# attr(,"class")
# [1] "bs"     "basis"  "matrix"

好的，degree正如我在输入中指定的那样，1 也是如此。knots告诉我两个内部结分别位于 x = 13.3000 和 x = 38.8333。看到结处于固定分位数有点惊讶，我希望 R 会为我的数据找到最佳分位数，但这当然会使模型非线性，并且在不知道响应数据的情况下也是不可能的。intercept = FALSE意味着基础中没有包含截距（这是一件好事吗？我一直被教导不要在没有截距的情况下拟合线性模型……好吧，我猜lm无论如何只是添加一个）。

但是，矩阵呢？我真的不明白如何解释它。有三列，我认为这意味着基函数是三个。这是有道理的：如果我有两个内部结K1和K2，我将在左边界结B1和之间有一个样条线K1，另一个样条线在K1和之间K2，最后一个在K2和之间B2，所以...三个基函数，好的。但究竟哪些是基函数？例如，此列是什么意思？

#              1
# [1,] 0.0000000
# [2,] 0.8270677
# [3,] 0.8198433
# [4,] 0.0000000
# [5,] 0.0000000

编辑：这与此问题相似但不完全相同。该问题询问回归系数的解释，但我比这更早一步：我想了解模型矩阵系数的含义。如果我尝试制作与第一个答案中建议的相同的情节，我会得到一个混乱的情节：

b <- bs(x, df = 3, degree = 1)
b1 <- b[, 1]  ## basis 1
b2 <- b[, 2]  ## basis 2
b3 <- b[,3]
par(mfrow = c(1, 3))
plot(x, b1, type = "l", main = "basis 1: b1")
plot(x, b2, type = "l", main = "basis 2: b2")
plot(x, b3, type = "l", main = "basis 3: b3")

在此处输入图像描述

这些不能是 B 样条基函数，因为它们有太多的结（每个函数应该只有一个）。

第二个答案实际上可以让我在 R 之外重建我的模型，所以我想我可以这样做。然而，这个答案也不能准确解释b矩阵的元素是什么：它处理线性回归的系数，我还没有在这里介绍。确实这是我的最终目标，但我也想了解这个中间步骤。

4

2 回答 2

12

矩阵b

#              1         2         3
# [1,] 0.0000000 0.0000000 0.0000000    
# [2,] 0.8270677 0.0000000 0.0000000    
# [3,] 0.8198433 0.1801567 0.0000000    
# [4,] 0.0000000 0.7286085 0.2713915    
# [5,] 0.0000000 0.0000000 1.0000000

实际上只是 的每个点中三个基函数的值的矩阵x，这对我来说应该是显而易见的，因为它与多项式线性模型的解释完全相同。事实上，由于边界结是

bknots <- attr(b,"Boundary.knots")
# [1]  0.0 77.4

和内部结是

iknots <- attr(b,"knots")
# 33.33333% 66.66667% 
#  13.30000  38.83333

那么三个基函数，如此处所示，是：

knots <- c(bknots[1],iknots,bknots[2])
y1 <- c(0,1,0,0)
y2 <- c(0,0,1,0)
y3 <- c(0,0,0,1)
par(mfrow = c(1, 3))
plot(knots, y1, type = "l", main = "basis 1: b1")
plot(knots, y2, type = "l", main = "basis 2: b2")
plot(knots, b3, type = "l", main = "basis 3: b3")

在此处输入图像描述

现在，考虑b[,1]

#              1
# [1,] 0.0000000
# [2,] 0.8270677
# [3,] 0.8198433
# [4,] 0.0000000
# [5,] 0.0000000

这些必须是b1in的值x <- c(0.0, 11.0, 17.9, 49.3, 77.4)。事实上，b1是 0 inknots[1] = 0和 1 in knots[2] = 13.3000，这意味着在x[2](11.0) 中，值必须是11/13.3 = 0.8270677，如预期的那样。同样，由于b1为 0 knots[3] = 38.83333，因此 (17.9) 中的值x[3]必须为(38.83333-13.3)/17.9 = 0.8198433。因为x[4], x[5] > knots[3] = 38.83333,b1那里是 0。其他两列可以给出类似的解释。

于 2016-09-12T12:55:39.360 回答
3

只是对上面@DeltaIV 的出色答案的一个小修正（看起来我无法评论。）

所以在b1他计算的时候b1(x[3])，应该是(38.83333-17.9)/(38.83333-13.3)=0.8198433线性插值。其他一切都很完美。

注意b1应该是这样的

\frac{t}{13.3}I(0<=t<13.3)+\frac{38.83333-t}{38.83333-13.3}I(13.3<=t<38.83333)

于 2019-06-18T21:59:18.727 回答