曲线拟合的最小二乘法

本文参考书为马东升著《数值计算方法》

最小二乘法

与插值法的区别：只需“逼近” f(x)f(x)f(x) ，而不用满足插值原则（即经过插值节点）

最小二乘原理

设 y=f(x)y=f(x)y=f(x) 在区间 [a,b][\ a\ ,b\ ][ a ,b ] 上有定义，求 φ(x)\varphi(x)φ(x) ，使得 φ(xi)\varphi(x_i)φ(xi) 与 f(xi)f(x_i)f(xi) 的误差或者距离最小

1-范数意义下的误差

R1=∑i=1n∣φ(xi)−yi∣R_1=\sum\limits_{i=1}^n|\varphi(x_i)-y_i|R1=i=1∑n∣φ(xi)−yi∣
∞\infty∞-范数意义下的误差

R∞=max⁡1≤i≤n∣φ(xi)−yi∣R_{\infty}=\max\limits_{1\le i\le n}|\varphi(x_i)-y_i|R∞=1≤i≤nmax∣φ(xi)−yi∣
2-范数意义下的误差（均方误差）

R22=∑i=1n(φ(xi)−yi)2R_2^2=\sum\limits_{i=1}^n(\varphi(x_i)-y_i)^2R22=i=1∑n(φ(xi)−yi)2

按均方误差达到极小构造拟合曲线的方法称为最小二乘法
正则方程组

系数矩阵为对称矩阵的方程组，称为正则方程组或法方程组

直线拟合

设已知数据点 (xi,yi)(x_i,y_i)(xi,yi) 分布大致为直线，构造拟合直线 y=a+bxy=a+bxy=a+bx ，则需使残差平方和
∑i=1m[yi−(a+bxi)]2\sum\limits_{i=1}^m[\ y_i-(a+bx_i)]^2 i=1∑m[ yi−(a+bxi)]2
为最小

计算时正则方程组为
(∑i=1m1∑i=1mxi∑i=1mxi∑i=1mxi2)(ab)=(∑i=1myi∑i=1mxiyi)\begin{pmatrix}{\sum\limits_{i=1}^m1}&{\sum\limits_{i=1}^mx_i}\\ {\sum\limits_{i=1}^mx_i}&{\sum\limits_{i=1}^mx_i^2}\end{pmatrix}\begin{pmatrix}a\\b\end{pmatrix}=\begin{pmatrix}{\sum\limits_{i=1}^my_i}\\{\sum\limits_{i=1}^mx_iy_i}\end{pmatrix} ⎝⎜⎛i=1∑m1i=1∑mxii=1∑mxii=1∑mxi2⎠⎟⎞(ab)=⎝⎜⎛i=1∑myii=1∑mxiyi⎠⎟⎞
如果有权值 ωi\omega_iωi ，则在每一项乘以 ωi\omega_iωi 即可，如下
(∑i=1mωi∑i=1mωixi∑i=1mωixi∑i=1mωixi2)(ab)=(∑i=1mωiyi∑i=1mωixiyi)\begin{pmatrix} {\sum\limits_{i=1}^m\omega_i}&{\sum\limits_{i=1}^m\omega_ix_i}\\ {\sum\limits_{i=1}^m\omega_ix_i}&{\sum\limits_{i=1}^m\omega_ix_i^2}\end{pmatrix} \begin{pmatrix}a\\b\end{pmatrix}= \begin{pmatrix}{\sum\limits_{i=1}^m\omega_iy_i}\\{\sum\limits_{i=1}^m\omega_ix_iy_i}\end{pmatrix} ⎝⎜⎛i=1∑mωii=1∑mωixii=1∑mωixii=1∑mωixi2⎠⎟⎞(ab)=⎝⎜⎛i=1∑mωiyii=1∑mωixiyi⎠⎟⎞

例：已知数据点

xix_ixi	0	0.2	0.4	0.6	0.8
yiy_iyi	0.9	1.9	2.8	3.3	4.2

求拟合直线 y=a+bxy=a+bxy=a+bx

解：设拟合直线为 y=a+bxy=a+bxy=a+bx ，则有正则方程组
(5221.2)(ab)=(13.16.84)\begin{pmatrix}5 &2\\2 &1.2\end{pmatrix} \begin{pmatrix}a\\b\end{pmatrix}= \begin{pmatrix}13.1\\6.84\end{pmatrix} (5221.2)(ab)=(13.16.84)
解得 a=1.02a=1.02a=1.02 ，b=4b=4b=4

故拟合直线为 y=1.02+4xy=1.02+4xy=1.02+4x

超定方程组的最小二乘解

给定数据序列 (xi,yi)(i=1,2,⋯,n)(x_i,y_i)\ (i=1,2,\cdots,n)(xi,yi) (i=1,2,⋯,n) ，作拟合直线 p(x)p(x)p(x) 过这些点，其中
p(xi)=a0+a1xi=yip(x_i)=a_0+a_1x_i=y_i p(xi)=a0+a1xi=yi
故有 nnn 个方程组却只有两个未知量，故方程组为矛盾方程组（无解）

因而需要寻求在均方误差极小意义下的解

设线性方程组 Ax=bAx=bAx=b ，其中
A=[aij]n×mx=[x1,x2,⋯,xn]Tb=[b1,b2,⋯,bm]TA=[\ a_{ij}\ ]_{n\times m}\ \ x=[\ x_1\ ,x_2\ ,\cdots\ ,x_n\ ]^T\ \ b=[\ b_1\ ,b_2\ ,\cdots\ ,b_m\ ]^T A=[ aij ]n×m x=[ x1 ,x2 ,⋯ ,xn ]T b=[ b1 ,b2 ,⋯ ,bm ]T
若有向量 xxx 使得 ∑i=1m(ai1x1+ai2x2+⋯+aimxm−bi)2\sum\limits_{i=1}^m(a_{i1}x_1+a_{i2}x_2+\cdots+a_{im}x_m-b_i)^2i=1∑m(ai1x1+ai2x2+⋯+aimxm−bi)2 达到最小，则称 xxx 为超定方程组的最小二乘解

法方程

ATAX=ATbA^TAX=A^TbATAX=ATb 称为矛盾方程 AX=bAX=bAX=b 的法方程，当且仅当 XXX 满足 ATAX=ATbA^TAX=A^TbATAX=ATb ，即 XXX 是法方程的解

例：解矛盾方程组
{x1+x2+x3=2x1+3x2−x3=−12x1+5x2+2x3=13x1−x2+5x3=−2\begin{cases}x_1+x_2+x_3=2\\x_1+3x_2-x_3=-1\\2x_1+5x_2+2x_3=1\\3x_1-x_2+5x_3=-2\end{cases} ⎩⎪⎪⎪⎨⎪⎪⎪⎧x1+x2+x3=2x1+3x2−x3=−12x1+5x2+2x3=13x1−x2+5x3=−2
解：法方程为
(1511191136319331)(x1x2x3)=(−36−5)\begin{pmatrix}15 &11&19\\11&36&3\\19&3&31\end{pmatrix}\begin{pmatrix}x_1\\x_2\\x_3\end{pmatrix}=\begin{pmatrix}-3\\6\\-5\end{pmatrix} ⎝⎛1511191136319331⎠⎞⎝⎛x1x2x3⎠⎞=⎝⎛−36−5⎠⎞
解得 x1=1.5917x_1=1.5917x1=1.5917 ，x2=0.5899x_2=0.5899x2=0.5899 ，x3=0.7572x_3=0.7572x3=0.7572

注： ATA^TAT 在前，可记忆为 AX=bAX=bAX=b 左右两边在前面乘以 ATA^TAT

可线性化模型的最小二乘拟合

在实际问题中，变量关系不一定是线性，故可以将拟合函数变为
f(y)=a+bg(x)f(y)=a+bg(x) f(y)=a+bg(x)
即通过变量替换，把非线性问题转化为线性问题

如把双曲线
1y=a+b1x\dfrac 1y=a+b\dfrac 1x y1=a+bx1
转化为
y^=a+bx^\hat y=a+b\hat x y^=a+bx^

多变量的数据拟合

影响 yyy 的是多个变量，设拟合函数为 y=Aαy=A\alphay=Aα ，其中
y=(y1y2⋮ym)A=(1x11x21⋯xn11x12x22⋯xn2⋮⋮⋮⋮⋮1x1mx2m⋯xnm)α=(α0α1⋮αn)y=\begin{pmatrix}y_1\\y_2\\\vdots\\y_m\end{pmatrix}\ \ A=\begin{pmatrix}1 &x_{11}&x_{21}&\cdots&x_{n1}\\1 &x_{12}&x_{22}&\cdots&x_{n2}\\\vdots&\vdots&\vdots&\vdots&\vdots\\1 &x_{1m}&x_{2m}&\cdots&x_{nm}\end{pmatrix}\ \ \alpha=\begin{pmatrix}\alpha_0\\\alpha_1\\\vdots\\ \alpha_n\end{pmatrix} y=⎝⎜⎜⎜⎛y1y2⋮ym⎠⎟⎟⎟⎞ A=⎝⎜⎜⎜⎛11⋮1x11x12⋮x1mx21x22⋮x2m⋯⋯⋮⋯xn1xn2⋮xnm⎠⎟⎟⎟⎞ α=⎝⎜⎜⎜⎛α0α1⋮αn⎠⎟⎟⎟⎞
故有正则方程组
ATAα=ATyA^TA\alpha=A^Ty ATAα=ATy
当 ATAA^TAATA 非奇异时，可求得
α=(ATA)−1ATy\alpha=(A^TA)^{-1}A^Ty α=(ATA)−1ATy

多项式拟合

正则方程组为
(m∑xi∑xi2⋯∑xin∑xi∑xi2∑xi3⋯∑xin+1⋮⋮⋮⋮⋮∑xin∑xin+1∑xin+2⋯∑xi2n)(a0a1⋮an)=(∑yi∑xiyi⋮∑xinyi)\begin{pmatrix}m&\sum x_i&\sum x_i^2&\cdots&\sum x_i^n\\\sum x_i&\sum x_i^2&\sum x_i^3&\cdots&\sum x_i^{n+1}\\\vdots&\vdots&\vdots&\vdots&\vdots\\\sum x_i^n&\sum x_i^{n+1}&\sum x_i^{n+2}&\cdots&\sum x_i^{2n}\end{pmatrix} \begin{pmatrix}a_0\\a_1\\\vdots\\a_n\end{pmatrix} = \begin{pmatrix}\sum y_i\\\sum x_iy_i\\\vdots\\\sum x_i^ny_i\end{pmatrix} ⎝⎜⎜⎜⎛m∑xi⋮∑xin∑xi∑xi2⋮∑xin+1∑xi2∑xi3⋮∑xin+2⋯⋯⋮⋯∑xin∑xin+1⋮∑xi2n⎠⎟⎟⎟⎞⎝⎜⎜⎜⎛a0a1⋮an⎠⎟⎟⎟⎞=⎝⎜⎜⎜⎛∑yi∑xiyi⋮∑xinyi⎠⎟⎟⎟⎞
其中 ∑\sum∑ 是 ∑i=1m\sum\limits_{i=1}^mi=1∑m 的缩写

例：给定数据表

xxx	1	2	3	4	6	7	8
yyy	2	3	6	7	5	3	2

试用最小二乘法求二次拟合多项式

解：设二次拟合多项式
y=a0+a1x+a2x2y=a_0+a_1x+a_2x^2 y=a0+a1x+a2x2
写出正则方程组
(7∑xi∑xi2∑xi∑xi2∑xi3∑xi2∑xi3∑xi4)(a0a1a2)=(∑yi∑xiyi∑xi2yi)\begin{pmatrix}7&\sum x_i&\sum x_i^2\\\sum x_i&\sum x_i^2&\sum x_i^3\\\sum x_i^2&\sum x_i^3&\sum x_i^4\end{pmatrix}\begin{pmatrix}a_0\\a_1\\a_2\end{pmatrix}=\begin{pmatrix}\sum y_i\\\sum x_iy_i\\\sum x_i^2y_i\end{pmatrix} ⎝⎛7∑xi∑xi2∑xi∑xi2∑xi3∑xi2∑xi3∑xi4⎠⎞⎝⎛a0a1a2⎠⎞=⎝⎛∑yi∑xiyi∑xi2yi⎠⎞
解得 a0=−1.3185a_0=-1.3185a0=−1.3185 ，a1=3.4321a_1=3.4321a1=3.4321 ，a2=−0.3864a_2=-0.3864a2=−0.3864

正交多项式及最小二乘拟合

数值计算方法第五章—曲线拟合的最小二乘法相关推荐

工程数学计算方法第五章数值积分
工程数学计算方法第五章数值积分数值积分数值微分微分差商法原理误差分析优化:变步长算法插值法数值积分机械求积公式插值型数值积分代数精度⭐ 节点等距分布复合求积复合求积基 ...
计算机最小二乘法拟合,数值计算_第6章曲线拟合的最小二乘法.doc
数值计算_第6章曲线拟合的最小二乘法.doc 还剩 29页未读, 继续阅读下载文档到电脑,马上远离加班熬夜! 亲,很抱歉,此页已超出免费预览范围啦! 如果喜欢就下载吧,价低环保! 内容要点: 第 ...
曲线拟合最小二乘法对数c语言实现,数值计算_第6章曲线拟合的最小二乘法
<数值计算_第6章曲线拟合的最小二乘法>由会员分享,可在线阅读,更多相关<数值计算_第6章曲线拟合的最小二乘法(23页珍藏版)>请在人人文库网上搜索. 1.第6章曲线拟合的最 ...
科学计算机解多元方程组,南京大学计算机科学与技术系数值计算方法(第3章)3 向量和矩阵的范数病态方程组解线性方程组的迭代法.ppt...
南京大学计算机科学与技术系数值计算方法(第3章)3 向量和矩阵的范数病态方程组解线性方程组的迭代法.ppt 上传人:ct****5 文档编号:43116102 上传时间:2020-01-21 格 ...
数值计算方法第四章—插值法
插值法本文参考书为马东升著<数值计算方法> 代数插值插值已知某些(有限)点的函数值求其余点的函数值代数多项式插值 n+1n+1n+1 个互异节点处满足插值条件 P(xi)=yi,i ...
数值计算方法第三章—线性方程组的数值解法知识点总结
线性方程组的数值解法本文参考书为马东升著<数值计算方法> 高斯消去法顺序高斯消去法通过初等变换消去方程组系数矩阵主对角线以下的元素,而使方程组化为等价的上三角形方程组列主元高斯消去 ...
数值计算方法第六章—数值积分和数值微分
数值积分和数值微分本文参考书为马东升著<数值计算方法> 数值积分概述基本思想用简单的函数代替复杂的 f(x)f(x)f(x) 求得复杂的积分梯形公式 I=∫abf(x)dx≈(b− ...
角速度求积分能得到欧拉角吗_数值计算方法第六章数值积分和数值微分
写在章前:积分与微分的计算,是具有广泛应用的古典问题. 然而,在微积分教材中,只对简单的或特殊的情况,提供了函数的积分或微分的解析表达式,在理论上可以使用Newton-Leibniz公式计算.但对很多 ...
数值计算方法（五）——迭代法求方程根
(一)直接迭代数学描述: 代码实现: /***@name Equation_iteration:方程求根的迭代法*@param1 x0:初始值 **/ double Equation_iterati ...

数值计算方法第五章—曲线拟合的最小二乘法