KAN保姆级解析(2)--KAN源码中的BSpline插值

引言

  KAN源码中KANlayer模块涉及到了许多coef2curve和coef2curve操作，这里实际就是在利用Bspline曲线进行插值或者说拟合，以下介绍相关的BSpline插值知识。后一篇中我们将结合相关代码详细解释Bspline是如何融入到KAN的计算过程中的。

Bspline曲线

 对于Bspline曲线的构造我们只需要给定一系列的控制点$p_0,...,p_n$（control points）以及节点向量 U = { u 0 , . . , u m + 1 } U=\{u_0,..,u_{m+1}\} U={u0​,..,um+1​}（knot vector），以及次数$k$即可进行计算。
 首先，BSpline曲线的基础计算公式是：
$$C(u)=\sum _{i=0}^n{p}_i{N}_{i,k}(u)$$
其中$p_i$是人为选取或者求取的控制点，$N_{i,k}(t)$是$k$次（$k+1$阶）基函数（这里，$k$次的意思是$N_{i,k}(t)$中包含了$t$的$k$次项$t^k$，这里我们假定$k$从0开始取值）。
 在这一公式的基础上我们只需要确定控制点以及基函数，就可以计算给定$t$时的曲线点位置。后者可以通过Cox-de Boor递归公式来确定：

$$\begin{array}{rl}& N_{i,0}\left(u\right)=\{\begin{array}{cc}1& \mathrm{if}{u}_i\le u<u_{i+1}\mathrm{and}{u}_i<u_{i+1}\\ & \\ 0& \mathrm{otherwise}\end{array}\\ & N_{i,j}(u)=\frac{u-u_i}{u_{i+j}-u_i}{N}_{i,j-1}(u)+\frac{u_{i+j+1}-u}{u_{i+j+1}-u_{i+1}}{N}_{i+1,j-1}(u),\end{array}$$

节点向量

这里需要注意的是，节点向量的设置是比较有讲究的，假定需要$m+1$个节点，通常有以下几个方案

1. 在$[0,1]$范围内进行划分（通常就是等距）；这样形成了open Bspline曲线，不经过控制点折线的第一段和最后一段；
2. 前$k+1$个节点设置为0，$u_0=u_1=,..,=u_k=0$，中间进行采样（通常就是等距），后$k+1$个节点设置为1，$u_{m-k-1}=,..=u_m=1$;这样形成了 clamped B-spline曲线，它会穿过第一个和最后一个控制点，并且和第一段最后一段折线相切；
3. 重复多个节点。这样形成closed B-spline曲线.

次数、控制点数目和节点数目的关系¹

 次数$k$，控制点数目$n+1$和节点数目$m+1$需要满足$m=n+k+1$，也就是说如果给出了5个控制点，要求计算次数为3的B样条曲线，我们需要自行设置长度为9的节点向量。这是因为由Cox-de Boor递归公式我们可以知道计算 $N_{i,k}(u)$，我们需要计算 $N_{i,k-1}(u)$ 和 $N_{i+1,p-1}(u)$，可以画出如下计算流图：

从而得到$m=n+k+1$。
通过推导我们也可以进一步得到以下两个性质：

1. 基函数 $N_{i,k}(u)$ 在区间 $[u_i,u_{i+k+1})$ 上非零。换句话说，$N_{i,k}(u)$ 在 $k+1$ 个结节区间 $[u_i,u_{i+1}),[u_{i+1},u_{i+2}),\dots ,[u_{i+k},u_{i+k+1})$ 上非零。
2. 在任意结节区间 $[u_i,u_{i+1})$ 上，最多有 $k+1$ 个度数为 $k$ 的基函数是非零的，即：$N_{i-k,k}(u),N_{i-k+1,k}(u),N_{i-k+2,k}(u),\dots ,N_{i-1,k}(u)$ 和 $N_{i,k}(u)$。

Bspline插值

 KAN中首先使用到的技术就是BSpline插值，整体而言，使用Bspline进行曲线插值的步骤可以总结为如下：

步骤1：数据准备

准备需要拟合的数据点集合 D = { D 0 , D 1 , . . , D n } \mathbb{D}=\{\mathbf{D_0},\mathbf{D_1},..,\mathbf{D_n}\} D={D0​,D1​,..,Dn​}，通常是二维或三维空间中的点。

步骤2：节点向量的选择

 选择节点向量，这决定了B样条基函数的形状和数量。节点向量可以是等距的或非等距的，通常情况下选择非等距节点能够更灵活地控制曲线的形状，简单情况下，我们进行等距采样。
 而在节点向量结构方面，我们会选用clamped结构（为什么这么选，其实是没懂的，感觉插值下反正都会经过第一个和最后一个控制点，使用open好像也没差）,有 $m+1$ 个结节 $u_0=u_1=\dots =u_p=0$, $u_{p+1}$, $u_{p+2}$, $\dots$, $u_{m-p-1}$, $u_{m-p}=u_{m-p+1}=\dots =u_m=1$。

步骤3：构建插值矩阵

 根据数据点和节点向量构建插值矩阵。插值矩阵的每一行对应一个数据点，每一列对应一个基函数在该数据点处的取值。插值矩阵每个元素的求取是通过Codx算法得到的，流程如下：

算法描述

输入：n、k、m、u 和 m+1 个闭合结节 { u 0 , . . . , u m } \{ u_0, ..., u_m \} {u0​,...,um​}
输出：系数 $N_{0,k}(u)$、$N_{1,k}(u)$、…、$N_{n,k}(u)$ 存储在 $N[0]$、$N[1]$、…、$N[n]$ 中

算法步骤²

1. 将 $N[\mathrm{0..}n]$ 初始化为 0； // 初始化

2. 如果 $u=u_0$，则 // 排除特殊情况

   • $N[0]=1.0$;
   • 返回

3. 否则如果 $u=u_m$，则

   • $N[n]=1.0$;
   • 返回

4. 结束如果

5. // 现在 $u$ 在 $u_0$ 和 $u_m$ 之间

6. 令 $u$ 在结节区间 $[u_k,u_{k+1})$ 中；

7. $N[k]:=1.0$；// 度数 0 的系数

8. 对于 $d$ 从 1 到 $k$ 循环执行：// 度数 $d$ 从 1 到 $k$

   • 开始
   • $N[k-d]=\frac{u_{k+1}-u}{u_{k+1}-u_{k-d+1}}\ast N[(k-d)+1]$; // 右侧（西南角）项
   • 对于 $i$ 从 $k-d+1$ 到 $k-1$ 循环执行：// 计算内部项
     • $N[i]:=\frac{u-u_i}{u_{i+d}-u_i}\ast N[i]+\frac{u_{i+d+1}-u}{u_{i+d+1}-u_{i+1}}\ast N[i+1]$;
   • $N[k]=\frac{u-u_k}{u_{k+d}-u_{i+k}}\ast N[k]$; // 左上角项
   • 结束

9. // 数组 $N[\mathrm{0..}n]$ 中存储了系数。

对每前$n+1$节点都进行这样的操作（后面的节点会因为计算项数不足被舍掉，所以只会有n+1个这样的求系数操作），我们就能得到插值矩阵：
[
M_k = \begin{bmatrix}
N_{0,k}(u_0) & N_{1,k}(u_0) & \cdots & N_{n,k}(u_0) \
N_{0,k}(u_1) & N_{1,k}(u_1) & \cdots & N_{n,k}(u_1) \
\vdots & \vdots & \ddots & \vdots \
N_{0,k}(u_n) & N_{1,k}(u_n) & \cdots & N_{n,k}(u_n)
\end{bmatrix}
]
其中，$u_0,u_1,\dots ,u_n$ 是插值点，$N_{i,k}(u)$ 是$k$次B样条基函数在该节点处的取值。

步骤4：控制点的计算

使用最小二乘法，求解线性方程组$D=NP$其中 $N$是插值矩阵，$P$是控制点的向量，$D$是数据点坐标向量。

步骤5：B样条曲线的计算

对于给出的任意待插值坐标序列$\mathbf{u}$，使用Cox-de Boor递归公式得到在$n$个基函数在该位置上的取值，使用插值点和计算得到的控制点$p$进行加权求和，计算code构建B样条曲线。

总结

 使用Bspline进行插值实际只有两个关键步骤：

1. 构造节点向量
2. 计算控制点
   然后基于Cox-de Boor递归公式就可以计算出插值曲线在任意位置的值了。

参考

---

1. B-spline Basis Functions: Definition ↩︎

2. Curve Global Interpolation ↩︎