[科学计算]曲线拟合——最小二乘法

🚀宁可十年不将军，不可一日不拱卒!

问题：已知 $x_1,x_2,\cdots ,x_m;y_1,y_2,\cdots ,y_m.$ 求一个简单易算的近似函数 $P(x)\approx f(x)$ 。 但是此时 $m$ 很大，而且 $y_i\ne f(x_i)$ ，因此这个时候数值插值就无可奈何了，没有必要再令 $P(x_i)=y_i.i=1,2,\cdots ,m.$ 了，否则会造成过拟合现象，而偏离真实 $f(x)$ 。我们的目标仅仅是 $P(x_i)-y_i$ 总体上尽可能的小。

常见的作法有：

• 使 $\underset{0\le i\le m}{\max }|P(x_i)-y_i|$ 最小，这种方法太复杂，究竟为什么太复杂呢？需要进一步考虑。
• 使 $\sum _{i=1}^m|P(x_i)-y_i|$ 最小，不可导，求解困难。
• 使 $\sum _{i=1}^m|P(x_i)-y_i{|}^2$ 最小，此时对参数的估计是极大似然估计。详细参照另一篇文章。

最小二乘法拟合多项式

确定多项式的格式为 $P(x)=a_0+a_{1}x+a_2{x}^2+\cdots +a_n{x}^n$ ，由此多项式产生的一组观测值为 $(x_i,y_i)(i=1,2,\cdots ,m)$ ，注意观测值是有一定的误差的。我们需要通过一系列的操作，学习到参数 $a_0,a_1,\cdots ,a_n$ 。经过分析，我们选择了使 $\psi (a_0,a_1,\cdots ,a_n)=\sum _{i=1}^m|P(x_i)-y_i{|}^2$ 达到最小值作为评估标准，$\psi$ 最小时的参数 $a_0,a_1,\cdots ,a_n$ 就是我们待求的。有一个大的前提是 $n<<m$ 。
$$\begin{array}{r}\psi (a_0,a_1,\cdots ,a_n)=\sum _{i=1}^m|a_0+a_1{x}_i+a_2{x}_i^2+\cdots +a_n{x}_i^n-y_i{|}^2\end{array}$$
从上式可知， $\psi$ 是一个关于 $a_0,a_1,\cdots ,a_n$ 的多元函数。根据费马定理，由于 $\psi$ 可导，因此在极值点处的必要条件是 $\frac{\partial \psi }{\partial a_i}=0,(i=1,2,\cdots ,n)$ 。这个时候大家是不是看出来选择评估函数时可导的重要性了！

接下来其实就是一个对公式进行整理简化的过程。下式中 $k=1,2,\cdots ,n$ 。
$$\begin{array}{rl}\frac{\partial \psi }{\partial a_k}& =\frac{\partial }{\partial a_k}\sum _{i=1}^m|a_0+a_1{x}_i+a_2{x}_i^2+\cdots +a_n{x}_i^n-y_i{|}^2\\ & =\sum _{i=1}^{m}2(a_0+a_1{x}_i+a_2{x}_i^2+\cdots +a_n{x}_i^n-y_i{|}^2{x}_i^k\\ & =2\sum _{i=1}^m(\sum _{j=0}^n{a}_j{x}_i^j-y_i)x_i^k\\ & =2[\sum _{i=1}^m\sum _{j=0}^n{a}_j{x}_i^{j+k}-\sum _{i=1}^m{y}_i{x}_i^k]\\ & =2[\sum _{j=0}^n{a}_j\sum _{i=1}^m{x}_i^{j+k}-\sum _{i=1}^m{y}_i{x}_i^k]\end{array}$$
很关键的一步替换来了：记 $b_k=\sum _{i=1}^m{x}_i^k,c_k=\sum _{i=1}^m{y}_i{x}_i^k$ , 因此
$$\begin{array}{rl}\frac{\partial \psi }{\partial a_k}& =2[\sum _{j=0}^n{a}_j{b}_{j+k}-c_k]\end{array}$$
令上式为 $0$ ，可得出
$$\begin{array}{r}\sum _{j=0}^n{a}_j{b}_{j+k}=c_k\end{array}$$
将 $n+1$ 个如上的等式全列出来构成一个线性方程组（称为法方程组）
$$\{\begin{array}{rl}\sum _{j=0}^n{a}_j{b}_{j+0}& =c_0\\ \sum _{j=0}^n{a}_j{b}_{j+1}& =c_1\\ ⋮& \\ \sum _{j=0}^n{a}_j{b}_{j+n}& =c_n\end{array}\Rightarrow \left[\begin{array}{cccc}{b}_0& b_1& \cdots & b_n\\ b_1& b_2& \cdots & b_{n+1}\\ ⋮& ⋮& \ddots & ⋮\\ b_n& b_{n+1}& \cdots & b_{2n}\end{array}\right]\left[\begin{array}{c}{a}_0\\ a_1\\ ⋮\\ a_n\end{array}\right]=\left[\begin{array}{c}{c}_0\\ c_1\\ ⋮\\ c_n\end{array}\right]$$
那么这个法方程组是否存在唯一解呢？它的解仅是极值点，那么这个极值点是否又是最小值点呢？这两个问题是有严格的定理证明的。

$定理：最小二乘法拟合多项式存在唯一$

$定理：最小二乘拟合多项式的法方程组的解确是\psi 的最小点。$

最小二乘法的一般形式

现实中，需要被拟合的函数形式各样，有的一些确实可以通过变换转换为多项式的拟合，但大多数情况下这种变换并不容易，因此，为了不失一般性，我们接下来讨论最小二乘法的一般形式。

若已知被拟合函数的形式为 $S(x)=a_0{\psi }_0(x)+a_1{\psi }_1(x)+\cdots +a_n{\psi }_n(x)(n<<m)$ 。这里 ${\psi }_i(x)(i=1,2,\cdots ,n)$ 已知，$a_i(i=1,2,\cdots ,n)$ 参数未知。

用最小二乘法拟合曲线，就是在众多的 $S(x)$ 中，求得一函数 $S^{\ast }(x)$ 使得 $||\delta |{|}_2^2$ 最小。
$$\begin{array}{r}||\delta |{|}_2^2=\sum _{i=0}^m{w}_i{\delta }_i^2=\sum _{i=0}^{m}w(x_i)(S(x_i)-y_i{)}^2\end{array}$$
其中 $w(x_i)>0$ 是不同的点 $(x_i,y_i)$ 的比重，$y_i=f(x_i)$ 是观测值, 其中 $f(x)=S(x)+\epsilon (x)$。

上述的问题可以转化为一个多元函数 $I(a_0,a_1,\cdots ,a_n)$ 求极小值 $(a_0^{\ast },a_1^{\ast },\cdots ,a_n^{\ast })$ 的问题。
$$\begin{array}{r}I(a_0,a_1,\cdots ,a_n)=\sum _{i=0}^{m}w(x_i)[\sum _{j=0}^n{a}_j{\psi }_j(x_i)-f(x_i){]}^2\end{array}$$
由求多元函数极值的必要条件，有 $k=0,1,\cdots ,n$
$$\begin{array}{rl}\frac{\partial I}{\partial a_k}& =2\sum _{i=0}^{m}w(x_i)[\sum _{j=0}^n{a}_j{\psi }_j(x_i)-f(x_i)]{\psi }_k(x_i)\\ & =2[\sum _{i=0}^{m}w(x_i)\sum _{j=0}^n{a}_j{\psi }_j(x_i){\psi }_k(x_i)-\sum _{i=0}^{m}w(x_i)f(x_i){\psi }_k(x_i)]\\ & =2[\sum _{j=0}^n{a}_j\sum _{i=0}^{m}w(x_i){\psi }_j(x_i){\psi }_k(x_i)-\sum _{i=0}^{m}w(x_i)f(x_i){\psi }_k(x_i)]\end{array}$$
若记 $({\psi }_j,{\psi }_k)=\sum _{i=0}^{m}w(x_i){\psi }_j(x_i){\psi }_k(x_i);(f,{\psi }_k)=\sum _{i=0}^{m}w(x_i)f(x_i){\psi }_k(x_i)$ ，由 $\frac{\partial I}{\partial a_k}=0$ 得
$$\begin{array}{r}\sum _{j=0}^n{a}_j({\psi }_j,{\psi }_k)=(f,{\psi }_k)\end{array}$$
从而用矩阵的形式可以得到一个法方程组
$$\begin{array}{r}\left[\begin{array}{cccc}({\psi }_0,{\psi }_0)& ({\psi }_1,{\psi }_0)& \cdots & ({\psi }_n,{\psi }_0)\\ ({\psi }_0,{\psi }_1)& ({\psi }_1,{\psi }_1)& \cdots & ({\psi }_n,{\psi }_1)\\ ⋮& ⋮& \ddots & ⋮\\ ({\psi }_0,{\psi }_n)& ({\psi }_1,{\psi }_n)& \cdots & ({\psi }_n,{\psi }_n)\end{array}\right]\left[\begin{array}{c}{a}_0\\ a_1\\ ⋮\\ a_n\end{array}\right]=\left[\begin{array}{c}(f,{\psi }_0)\\ (f,{\psi }_1)\\ ⋮\\ (f,{\psi }_n)\end{array}\right]\end{array}$$

思考与疑问

1. 拟合与插值的区别是什么？
2. 拟合需要提前知道被拟合函数的形式，也就是说需要有一个模型假设，那么这个作模型假设有什么依据吗？