数值计算方法 Chapter3. 曲线拟合的最小二乘法

393 阅读2分钟

本文已参与「新人创作礼」活动,一起开启掘金创作之路。

1. 线性拟合和二次拟合函数

最小二乘法本质上就是求一个事先定义一个函数,然后使用已知的采样点结果拟合函数的参数,使得所有采样点的均方误差最小。

即:

φ(x)=argmin iφ(xi)yi2\varphi(x) = argmin\ \sum_{i}|\varphi(x_i) - y_i|^2

1. 线性拟合

我们假定拟合曲线为:

φ(x)=ax+b\varphi(x) = ax + b

则有拟合误差为:

Q=i(axi+byi)2Q = \sum_{i}(ax_i+b-y_i)^2

要使得拟合误差QQ最小,则我们有QQ对于参数a,ba,b的偏导均为0,因此,我们即有:

{Qa=i2xi(axi+byi)=0Qb=i2(axi+byi)=0\left\{ \begin{aligned} \frac{\partial Q}{\partial a} &= \sum_i 2x_i(ax_i + b - y_i) = 0 \\ \frac{\partial Q}{\partial b} &=\sum_i 2(ax_i + b - y_i) = 0 \end{aligned} \right.

可以解得:

{a=nixiyi(ixi)(iyi)nixi2(ixi)2b=1niyi1n(ixi)a\left\{ \begin{aligned} a &= \frac{n \cdot \sum_{i} x_iy_i - (\sum_i x_i)(\sum_i y_i)}{n \cdot \sum_i{x_i^2} - (\sum_i x_i)^2} \\ b &= \frac{1}{n} \sum_i y_i - \frac{1}{n} (\sum_i x_i) \cdot a \end{aligned} \right.

2. 二次拟合函数

类似的,我们可以得到二次拟合函数的最小二乘法的结果。

定义拟合函数为:

φ(x)=ax2+bx+c\varphi(x) = ax^2 + bx + c

则有拟合误差:

Q=i(axi2+bxi+c)2Q = \sum_i(ax_i^2 + bx_i + c)^2

同样,我们可以解得,在极值条件下,有:

{Qa=i2xi2(axi2+bxi+cyi)=0Qb=i2xi(axi2+bxi+cyi)=0Qc=i2(axi2+bxi+cyi)=0\left\{ \begin{aligned} \frac{\partial Q}{\partial a} &= \sum_i 2x_i^2(ax_i^2 + bx_i + c - y_i) = 0 \\ \frac{\partial Q}{\partial b} &=\sum_i 2x_i(ax_i^2 + bx_i + c - y_i) = 0 \\ \frac{\partial Q}{\partial c} &=\sum_i 2(ax_i^2 + bx_i + c - y_i) = 0 \end{aligned} \right.

亦即:

(ixi4ixi3ixi2ixi3ixi2ixiixi2ixin)(abc)=(ixi2yiixiyiiyi)\begin{pmatrix} \sum_i x_i^4 & \sum_i x_i^3 & \sum_i x_i^2 \\ \sum_i x_i^3 & \sum_i x_i^2 & \sum_i x_i \\ \sum_i x_i^2 & \sum_i x_i & n \end{pmatrix} \begin{pmatrix} a \\ b \\ c \end{pmatrix} = \begin{pmatrix} \sum_i x_i^2y_i \\ \sum_i x_iy_i \\ \sum_i y_i \end{pmatrix}

可解得:

(abc)=(ixi4ixi3ixi2ixi3ixi2ixiixi2ixin)1(ixi2yiixiyiiyi)\begin{pmatrix} a \\ b \\ c \end{pmatrix} = \begin{pmatrix} \sum_i x_i^4 & \sum_i x_i^3 & \sum_i x_i^2 \\ \sum_i x_i^3 & \sum_i x_i^2 & \sum_i x_i \\ \sum_i x_i^2 & \sum_i x_i & n \end{pmatrix}^{-1} \cdot \begin{pmatrix} \sum_i x_i^2y_i \\ \sum_i x_iy_i \\ \sum_i y_i \end{pmatrix}

类似的也可以推至n阶的情况,这里就不一一展开了。

3. aebxa\cdot e^{bx}型函数

形如aebxa\cdot e^{bx}的函数直接用最小二乘法倒是没法直接求解,不过可以通过一定的函数变换转换成nn阶函数形式,然后我们就可以仿照上述方式进行求解了。

比如这里的aebxa\cdot e^{bx}类型的函数,我们做函数变换,令z=ln(y)z= ln(y),则有:

z=bx+ln(a)z = bx + ln(a)

即可使用最小二乘法进行求解。

类似的,还可以有函数y=1ax+by = \frac{1}{ax+b},我们令z=1yz = \frac{1}{y},即可有:

z=ax+bz = ax + b

2. 解矛盾方程组

书中这一小节事实上就是给前面最小二乘法的内容提供一些理论上的支持,没有啥更多的内容,因此,我们仅在这里摘录书中的定理如下:

定理3.1
(1) AAmmnn列的矩阵,bb为列向量,ATAX=ATbA^TAX = A^Tb称为矛盾方程组AX=bAX = b的法方程,法方程恒有解。
(2) XXminAXb22min||AX-b||_2^2的解,当且仅当XX满足ATAX=ATbA^TAX = A^Tb,即XX是法方程的解。