Scilab中的样条系数

Question

我想将分段三次样条拟合到大量数据。我不认为我想要 B 样条，因为我希望样条准确地通过数据点。在 Scilab 中我能看到的唯一方法是使用 splin 和 interp。

但是，我想计算结点之间每个样条曲线的系数（因为我需要获取这些系数并将它们放入不同的软件中）。所有splin 给你的是衍生品。有没有办法获得三次样条系数？或者有没有一种方法可以轻松地从一阶导数生成系数？

Answer 1

是的，可以从您拥有的 y 值和splin 返回的导数值中获取系数。每个区间 [x(i), x(i+1)] 你有 4 个系数要找到，还有 4 个方程：两端的值，两端的导数。最直接的方法是告诉 Scilab 为每个子区间求解这个 4 x 4 系统：这不应该比样条本身的评估花费更多的时间。下面的程序就是这样做的。

x = [0,1,2,3,4,5]   // x values
y = [1,0,1,0,1,0]   // y values
d = splin(x,y)      
n = length(x)-1     // number of subintervals
cfs = zeros(4,n)    // matrix to store coefficients in
for i=1:n
    a = x(i)
    b = x(i+1)
    cfs(:,i) = [1,a,a^2,a^3; 1,b,b^2,b^3; 0,1,2*a,3*a^2; 0,1,2*b,3*b^2] \ [y(i);y(i+1);d(i);d(i+1)]
end

前两个方程1,a,a^2,a^3; 1,b,b^2,b^3 将多项式的值与 y 值相关联；另外两个0,1,2*a,3*a^2; 0,1,2*b,3*b^2 对其导数执行相同的操作。 （公式只是 x 的幂的导数。）

上述脚本的输出：

    1.     1.   - 8.6    13.     13.   
  - 3.4  - 3.4    11.  - 10.6  - 10.6  
    3.1    3.1  - 4.1    3.1     3.1   
  - 0.7  - 0.7    0.5  - 0.3   - 0.3   

每一列都有四个系数：例如，样条的第一部分是1-3.4x+3.1x^2-0.7x^3。由于这是非结样条，splin 命令的默认模式，第二段与第一段相同；并且 last 与倒数第二个相同。

您可以通过绘制片段来检查它是否正常工作：

for i=1:n  
    t = linspace(x(i),x(i+1))
    plot(t,cfs(:,i)'*[ones(t); t; t.^2; t.^3])
end

---

也就是说，用以下形式表示形成样条的多项式会更容易

 p(x) = y(i) + A*(x-x(i)) + B*(x-x(i))*(x-x(i+1)) + C*(x-x(i))^2*(x-x(i+1))

其中系数很容易一一找到，无需求解线性系统：

• A = (y(i+1)-y(i))/(x(i+1)-x(i)) 等于 x(i+1) 处的值
• B = (d(i)-A)/(x(i)-x(i+1))，通过将导数等同于 x(i)
• C = (d(i+1)-A-B*(x(i+1)-x(i)))/(x(i+1)-x(i))^2，通过将导数等同于 x(i+1)

当然，这些系数应该与上述适当的多项式一起使用。这是这个替代版本

for i=1:n
    A = (y(i+1)-y(i))/(x(i+1)-x(i))
    B = (d(i)-A)/(x(i)-x(i+1))
    C = (d(i+1)-A-B*(x(i+1)-x(i)))/(x(i+1)-x(i))^2
    cfs(:,i) = [y(i);A;B;C]
end
// Again, plot for testing
for i=1:n
    t = linspace(x(i),x(i+1))
    plot(t,cfs(:,i)'*[ones(t); t-x(i); (t-x(i)).*(t-x(i+1)); ((t-x(i)).^2).*(t-x(i+1))])
end