您可以对样条公式进行逆向工程，而无需进入R代码。 知道这一点就足够了

• 样条是分段多项式函数。

• 次多项式由它们在点处的值确定。$d$$d+1$

• 多项式的系数可以通过线性回归获得。

因此，您只需在每对连续节点（包括数据范围的隐含端点）之间创建间隔的到的幂对预测进行回归。每个这样的结“bin”中的每个样条基元素都有一个单独的公式。例如，在下面的示例中，使用了三个内部结（对于四个结箱）和三次样条（），产生三次多项式，每个具有系数。的幂次比较高$d+1$$x$$x^d$$d=3$$4\times 4=16$$d+1=4$$x$涉及，必须保持系数中的所有精度。正如您可能想象的那样，任何样条基本元素的完整公式都会变得很长！

正如我很久以前提到的，能够将一个程序的输出用作另一个程序的输入（无需人工干预，这会引入不可重现的错误）是一种有用的统计沟通技巧。这个问题提供了一个很好的例子来说明该原则是如何应用的：与其手动复制这个 16 位系数，我们可以拼凑出一种方法，将计算的样条曲线转换为 Excel 可以理解的公式。我们需要做的就是从上述方法中提取样条系数，将它们重新格式化为类似 Excel 的公式，然后将它们复制并粘贴到 Excel 中。$64$RR

这种方法适用于任何统计软件，甚至是源代码不可用的无证专有软件。

这是一个取自问题的示例，但已修改为在三个内部点 ( ) 以及端点处具有结。绘图显示的版本，然后是 Excel 的渲染。在这两种环境中都很少进行自定义（除了指定颜色以大致匹配 Excel 的默认颜色）。$200, 500, 800$$(1, 1000)$RR

（版本中的垂直灰色网格线R显示了内部结的位置。）

这是完整的R代码。这是一个简单的 hack，完全依赖paste函数来完成字符串操作。（更好的方法是创建一个公式模板并使用字符串匹配和替换命令填充它。）

#
# Create and display a spline basis.
#
x <- 1:1000
n <- ns(x, knots=c(200, 500, 800))

colors <- c("Orange", "Gray", "tomato2", "deepskyblue3")
plot(range(x), range(n), type="n", main="R Version",
     xlab="x", ylab="Spline value")
for (k in attr(n, "knots")) abline(v=k, col="Gray", lty=2)
for (j in 1:ncol(n)) {
  lines(x, n[,j], col=colors[j], lwd=2)
}
#
# Export this basis in Excel-readable format.
#
ns.formula <- function(n, ref="A1") {
  ref.p <- paste("I(", ref, sep="")
  knots <- sort(c(attr(n, "Boundary.knots"), attr(n, "knots")))
  d <- attr(n, "degree")
  f <- sapply(2:length(knots), function(i) {
    s.pre <- paste("IF(AND(", knots[i-1], "<=", ref, ", ", ref, "<", knots[i], "), ", 
                   sep="")
    x <- seq(knots[i-1], knots[i], length.out=d+1)
    y <- predict(n, x)
    apply(y, 2, function(z) {
      s.f <- paste("z ~ x+", paste("I(x", 2:d, sep="^", collapse=")+"), ")", sep="")
      f <- as.formula(s.f)
      b.hat <- coef(lm(f))
      s <- paste(c(b.hat[1], 
            sapply(1:d, function(j) paste(b.hat[j+1], "*", ref, "^", j, sep=""))), 
            collapse=" + ")
      paste(s.pre, s, ", 0)", sep="")
    })
  })
  apply(f, 1, function(s) paste(s, collapse=" + "))
}
ns.formula(n) # Each line of this output is one basis formula: paste into Excel

第一个样条输出公式（这里产生的四个）是

"IF(AND(1<=A1, A1<200), -1.26037447288906e-08 + 3.78112341937071e-08*A1^1 + -3.78112341940948e-08*A1^2 + 1.26037447313669e-08*A1^3, 0) + IF(AND(200<=A1, A1<500), 0.278894459758071 + -0.00418337927419299*A1^1 + 2.08792741929417e-05*A1^2 + -2.22580643138594e-08*A1^3, 0) + IF(AND(500<=A1, A1<800), -5.28222778473101 + 0.0291833541927414*A1^1 + -4.58541927409268e-05*A1^2 + 2.22309136420529e-08*A1^3, 0) + IF(AND(800<=A1, A1<1000), 12.500000000002 + -0.0375000000000067*A1^1 + 3.75000000000076e-05*A1^2 + -1.25000000000028e-08*A1^3, 0)"

要在 Excel 中使用此功能，您需要做的就是删除周围的引号并在其前面加上“=”符号。（通过更多的努力，您可以R编写一个文件，当由 Excel 导入时，该文件在所有正确的位置包含这些公式的副本。）将其粘贴到公式框中，然后拖动该单元格直到“A1”引用第一个要计算样条的复制并粘贴（或拖放）该单元格以计算其他单元格的值。我用这些公式填充单元格 B2:E:102，引用值。$x$$x$

rms您可能会发现使用 R包更容易使用三次回归样条的截断幂基础。拟合模型后，您可以使用 中的Function或latex函数检索拟合样条函数的代数表示rms。

您已经执行了以下操作：

> rm(list=ls())
> set.seed(1066)
> x<- 1:1000
> y<- rep(0,1000)
> y[1:500]<- pmax(x[1:500]+(runif(500)-.5)*67*500/pmax(x[1:500],100),0.01)
> y[501:1000]<-500+x[501:1000]^1.05*(runif(500)-.5)/7.5
> df<-as.data.frame(cbind(x,y))
> library(splines)
> spline1 <- glm(y~ns(x,knots=c(500)),data=df,family=Gamma(link="log"))
> 

现在我将向您展示如何以两种不同的方式预测 x=12 的（响应）：首先使用 predict 函数（简单的方法！）

> new.dat=data.frame(x=12)
> predict(spline1,new.dat,type="response")
       1 
68.78721 

第二种方式是直接基于模型矩阵。注意我使用exp的链接功能是日志。

> m=model.matrix( ~ ns(df$x,knots=c(500))) 
> prd=exp(coefficients(spline1) %*% t(m)) 
> prd[12]
[1] 68.78721

请注意，在上面我提取了第 12 个元素，因为它对应于 x=12。如果您想预测训练集之外的 x，那么您可以再次使用 predict 函数。假设我们想要找到 x=1100 的预测响应值

> predict(spline1, newdata=data.frame(x=1100),type="response")
       1 
366.3483