MATLAB Smoothing Spline 拟合

参考

The Elements of Statistical Learning (chapter 5.4)
MATLAB - Smoothing Splines
MATLAB - fit

1. 基础

Smoothing Spline 可以用于离散数据的函数拟合。考虑下面的问题：在所有存在二阶连续导数的函数中寻找拟合函数$f(x)$，可以使下面式子的值最小，$RSS$可以理解为惩罚系数。
R S S ( f , λ ) = ∑ i = 1 N { y i − f ( x i ) } 2 + λ ∫ { f ′ ′ ( t ) } d t RSS(f,\lambda)=\sum_{i=1}^{N}\{y_i-f(x_i)\}^2+\lambda \int\{f^{''}(t)\}dt RSS(f,λ)=i=1∑N​{yi​−f(xi​)}2+λ∫{f′′(t)}dt
式中$\lambda$是一个常数，称为平滑系数（smoothing parameter）。式子中的前一部分用来衡量拟合曲线与原数据的近似程度，后半部分称为曲率惩罚。$\lambda$的选取原则可以参考下面两种极限情况。

• $\lambda =0$：$f$可以是任意函数。
• $\lambda =\infty$： 简单的最小二乘直线拟合，因为此时任意非零二阶导数产生的惩罚系数都是都是$\infty$而不能被接受，要使RSS最小，$f(x)=0$。

这里我们可以思考二阶导数的含义，即函数在某一点斜率的变化率，二阶导数越大，函数曲线的走向就越曲折，因此当$\lambda =0$，拟合函数$f(x)$可以无限弯曲，经过所有样本点，使误差的平方和为0，从而使RSS最小。而当$\lambda =\infty$，拟合函数$f(x)$无限平滑，以至于成为一条直线。我们总能在$\lambda \in [0,\infty )$中找到符合要求的拟合函数$f(x)$。

2. MATLAB 应用

MATLAB中所使用的Smoothing Splines代价函数与上面的定义非常类似。
$$p\sum _i{w}_i(y_i-s(x_i){)}^2+(1-p)\int (\frac{d^{2}x}{d{x}^2}{)}^{2}dx$$
其中权重$w_i$如果没有定义，默认都是1。
$p$的定义域是0到1，即$p\in [0,1]$当$p$的取值越小，拟合结果越平滑，而当$p$的取值越大，拟合结果就越倾向于经过所有的点。具体原因已经在上文中解释，不再赘述。

2.1 拟合

p=0.01

x = [0.3; 0.4; 0.5; 0.6; 0.7; 0.8; 0.9; 1.0];
y = [0; 143; 298; 482; 762; 1290; 2880; 4016];

f = fit(x, y, 'smoothingspline', 'SmoothingParam', 0.01);
plot(f, x, y);

---

p=1

x = [0.3; 0.4; 0.5; 0.6; 0.7; 0.8; 0.9; 1.0];
y = [0; 143; 298; 482; 762; 1290; 2880; 4016];

f = fit(x, y, 'smoothingspline', 'SmoothingParam', 0.01);
plot(f, x, y);

2.2 函数调用

对拟合生成的函数 f ,跟正常函数一样输入参数即可。

x = [0.3; 0.4; 0.5; 0.6; 0.7; 0.8; 0.9; 1.0];
y = [0; 143; 298; 482; 762; 1290; 2880; 4016];

f = fit(x, y, 'smoothingspline', 'SmoothingParam', 1);
f(0.6)

ans =

   482

函数返回结果与原数据一致（取决于拟合的程度）。