使用splines软件包中的bs函数指定b样条拟合的自由度

Question

我正在使用splines软件包的bs函数来创建一条b样条平滑曲线，用于图形目的。(至少有一个报告指出，Excel使用三阶b样条曲线来绘制平滑的折线图，我希望能够复制这些曲线。)我无法理解bs函数所需的参数。根据bs文档改编的代表性代码如下：

require(splines)
require(ggplot2)
n <- 10
x <- 1:10
y <- rnorm(n)
d <- data.frame(x=x, y=y)
summary(fm1 <- lm(y ~ bs(x, degree=3)), data=d)
x.spline <- seq(1, 10, length.out=n*10)
spline.data <- data.frame(x=x.spline, y=predict(fm1, data.frame(x=x.spline)))
ggplot(d, aes(x,y)) + geom_point + geom_line(aes(x,y), data=spline.data)

bs文档中的示例代码在对bs的调用中指定了df=5，而没有指定degree。我不知道我有多少自由度。我只知道我想要一个三阶b样条。我已经尝试指定不同的df值，而不是degree，或者除了degree之外，我得到了截然不同的结果。这就是为什么我怀疑df规范是这里的问题所在。在这种情况下，我如何计算df？

帮助文件建议df =长度(结)+度数。如果我将内点视为节点，这将为本例提供df=11，这将生成错误消息和无意义的样条拟合。

提前谢谢你。

我显然不清楚我的意图。我正在尝试这样做：How can I use spline() with ggplot?，但是使用b样条。

Answer 1

你不应该尝试去适应每一个点。我们的目标是找到一个可接受的拟合摘要，但这取决于有限数量的节点。在缺省值3以上增加多项式的次数没有太大的价值。只有10分，你肯定不想要df=11。试试df=5，结果应该是相当平坦的。rms/Hnisc软件包的作者Frank Harrell更喜欢限制三次样条，因为在极端情况下的预测是线性的，因此比其他多项式基出现的更少野性。

我更正了几个拼写错误，并添加了一个knots参数以使您的代码正常工作：

require(splines)
require(ggplot2); set.seed(trunc(100000*pi))

n <- 10
x <- 1:10
y <- rnorm(n)
d <- data.frame(x=x, y=y)
summary(fm1 <- lm(y ~ bs(x, degree=3, knots=2)), data=d)
x.spline <- seq(1, 10, length.out=n*10)
spline.data <- data.frame(x=x.spline, y=predict(fm1, data.frame(x=x.spline)))
ggplot(d, aes(x,y)) + geom_point() + geom_line(aes(x,y), data=spline.data)

我从改变随机种子的练习中走出来，认为Frank Harrell知道他在说什么。在使用他的软件包时，我在极端情况下不会得到相同的行为。

Answer 2

我做了更多的工作，并提出了以下建议。首先，道歉。我要找的是一条平滑的样条曲线，而不是回归样条曲线。我没有合适的词汇量来表达这个问题。虽然bs()的帮助文件中的示例似乎提供了这一点，但该函数没有为我的样本数据提供相同的行为。在stats包中还有另一个函数smooth.spline，它提供了我需要的东西。

set.seed(tunc(100000*pi))
n <- 10
x <- 1:n
xx <- seq(1, n, length.out=200)
y <- rnorm(n)
d <- data.frame(x=x, y=y)
spl <- smooth.spline(x,y, spar=0.1)
spline.data <- data.frame(y=predict(spl,xx))
ggplot(d,aes(x,y)) + geom_point() + geom_line(aes(x,y), spline.data)
spl2 <- smooth.spline(x, y, control=
            list(trace=TRUE, tol=1e-6, spar=0.1, low=-1.5, high=0.3))
spline.data2 <- data.frame(predit(spl2,xx))
ggplot(d,aes(x,y)) + geom_point() + geom_line(aes(x,y), spline.data2)

对smooth.spline的两个调用代表了两种方法。第一个方法手动指定平滑参数，第二个方法迭代到最优解。我发现我必须对优化进行适当的约束，才能得到我想要的解决方案。

结果应与Excel折线图使用的b样条曲线相匹配。我有一些合作者认为Excel图形是标准，我至少需要与之相匹配的性能。