我经常使用B-splines进行回归。到目前为止,我从来不需要详细了解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意味着在基础中没有包括截取(这是一件好事吗?我一直被教导不要在没有intercept...well的情况下拟合线性模型,lm只是添加了一个)。
但是,矩阵又如何呢?我真的不知道该怎么解释它。有三列,我认为这意味着基函数是三。这是有道理的:如果我有两个内部结点K1和K2,我将在左边界结点B1和K1之间有一条样条曲线,在K1和K2之间有另一条样条曲线,在K2和B2之间有一条样条曲线,so...three基函数,好的。但是到底哪些是基函数呢?例如,此列的含义是什么?
# 1
# [1,] 0.0000000
# [2,] 0.8270677
# [3,] 0.8198433
# [4,] 0.0000000
# [5,] 0.0000000编辑:这与this question相似,但并不完全相同。这个问题问的是回归系数的解释,但我比这更早一步:我想了解模型矩阵系数的含义。如果我尝试按照the first answer中的建议绘制相同的图,我会得到一个混乱的图:
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样条基函数,因为它们有太多的结(每个函数应该只有一个)。
The second answer实际上允许我在R之外重建我的模型,所以我想我可以这样做。然而,这个答案也没有确切地解释b矩阵的元素是什么:它处理的是线性回归的系数,我在这里还没有介绍。这确实是我的最终目标,但我也想了解这个中间步骤。
发布于 2016-09-12 20:55:39
矩阵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 然后,如here所示,三个基函数是:
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这些必须是x <- c(0.0, 11.0, 17.9, 49.3, 77.4)中b1的值。事实上,b1在knots[1] = 0中为0,在knots[2] = 13.3000中为1,这意味着在x[2] (11.0)中,该值必须为11/13.3 = 0.8270677,正如预期的那样。同样,由于knots[3] = 38.83333的b1为0,因此x[3] (17.9)中的值必须为(38.83333-13.3)/17.9 = 0.8198433。从x[4], x[5] > knots[3] = 38.83333开始,b1就是0。对于其他两列也可以给出类似的解释。
发布于 2019-06-19 05:59:19
只是对上面@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)
https://stackoverflow.com/questions/39446842
复制相似问题