最小二乘法曲线拟合及其MATLAB代码（使用和不使用曲线拟合工具箱）

赭红色的锆石

已于 2024-07-27 01:02:19 修改

阅读量835

点赞数 15

文章标签：最小二乘法 matlab 算法数学建模线性回归

于 2024-07-24 22:33:54 首次发布

本文链接：https://blog.csdn.net/Zircon_Xray/article/details/140669657

版权

来源：西电数模2024暑期双创课程作业

题目：设由实验测得的数据如下

$x_i$	-3	-2	-1	0	1	2	3
$y_i$	4	2	3	0	-1	-2	-5

试求一条二次曲线（权取 $\omega_i=1$ ），对它们进行最小二乘拟合。

一、最小二乘法

最小二乘法的基本思想

在科学实验和统计分析中，常常需要从一组测试数据中确定变量 x 与 y 之间的近似函数关系 $y=S\left ( x \right )$ ，即根据给定的n个点 $\left ( x_i,y_i \right )$ 确定一条曲线，使得 $\rho_i=S\left ( x_i \right )-y_i$ 总体上尽可能小，这种构造近似函数的方法称为曲线拟合， $S\left ( x \right )$ 称为拟合函数。

这种方法相比插值法的好处主要有三点：

①若用插值法，通过这几个已知点所求得的插值多项式必定是高次插值多项式，而高次插值是数值不稳定的。

②由于数据本身存在的误差，利用插值所得到的插值多项式必保留了所有测量误差，所得结果与实际问题误差较大。

③对数据拟合问题，一般地讲，并不要求所得到的近似解析表达式通过所有已知点，而只要求尽可能通过它们附近，这样可抵消原数据组中测量误差。

那么让 $\sum_{i=1}^{n}\rho_i$ 最小可以吗——偏差有正有负，在求和时可能互相抵消。

那么让 $\sum_{i=1}^{n}\left | \rho_i \right |$ ——可以，但有绝对值符号，不便于分析。

最终，我们选择 $\sum_{i=1}^{n}\rho_i^2$ 最小作为我们的标准。

拟合函数 $S\left ( x \right )$

一般而言，是由m个线性无关函数 $\varphi_1\left ( x \right )$ ， $\varphi_2\left ( x \right )$ ，...， $\varphi_m\left ( x \right )$ （基函数）线性组合而成，即

$S\left ( x \right )=a_1\varphi_1\left ( x \right )+a_2\varphi_2\left ( x \right )+\cdots a_m\varphi_m\left ( x \right )\quad\left ( m<n-1 \right )$

常用的基函数有 $1,x,x^2,\cdots,x^m$ ， $\sin x,\sin2x,\cdots ,\sin mx$ ， $e^{\lambda_1}x,e^{\lambda_2}x,\cdots,e^{\lambda_m}x$

二、用最小二乘法求解矛盾方程组

矛盾方程组

求解线性方程组时，通常要求未知数的个数与方程式个数相等。若方程式的个数多于未知数的个数，方程无解，称为矛盾（超定）方程组。即：

$\sum\limits_{j = 1}^m {{a_{ij}}{x_j} = {b_i}} \mathop {}\nolimits_{} (i = 1,{\rm{ }}2, \cdots ,{\rm{ }}n;{\rm{ }}m < n)$

求解过程

最小二乘法的基本思想：矛盾方程组无精确解，只能寻求某种意义下的近似解，这种近似解并非对精确解之近似，而是寻求各未知数的一组值，使方程组中各式能近似相等。令

${R_i} = \sum\limits_{j = 1}^m {{a_{ij}}{x_j}} - {b_i}\mathop {}\nolimits_{} \mathop {}\nolimits_{} (i = 1,{\rm{ }}2, \cdots ,{\rm{ }}n)$

按最小二乘法原则，求各个方程式误差平方和

$Q = \sum\limits_{i = 1}^n {{R_i}^2} = \sum\limits_{i = 1}^n {{{[\sum\limits_{j = 1}^m {{a_{ij}}{x_j}} - {b_i}]}^2}}$

若 ${x_j}{\rm{ }}(j = 1,2, \cdots ,m)$ 取值，使误差平方和Q达到最小，则称这组值是矛盾方程组的最优近似解。而Q可以看做m个 $x_j$ 自变量的二次函数，且连续。故存在一组数 ${x_1},{\rm{ }}{x_2},{\rm{ }}...,{\rm{ }}{x_m}$ ，使得Q达到最小（极值问题），即

${\raise0.7ex\hbox{${\partial Q}$} \!\mathord{\left/ {\vphantom {{\partial Q} {\partial {x_k}}}}\right.\kern-\nulldelimiterspace} \!\lower0.7ex\hbox{${\partial {x_k}}$}} = 0\mathop {}\nolimits_{} \mathop {}\nolimits_{} (k = 1,{\rm{ }}2, \cdots ,{\rm{ }}m)$

可以求解出极值条件为：

$\sum\limits_{j = 1}^m {(\sum\limits_{i = 1}^n {{a_{ij}}{a_{ik}}} ){x_j}} = \sum\limits_{i = 1}^n {{a_{ik}}{b_i}}\quad (k = 1,{\rm{ }}2, \cdots ,{\rm{ }}m)$

这是m个未知量，m个方程组的线性方程组——对应于原矛盾方程组的正规方程组。显然，这个方程组的解也是原矛盾方程组的最优近似解。

记

$A = \left[ \begin{array}{l} {a_{11}}{\rm{ }}{a_{12}} \cdots {\rm{ }}{a_{1m}}\\ {a_{21}}{\rm{ }}{a_{22}} \cdots {\rm{ }}{a_{2m}}\\ \vdots {\rm{ }} \vdots {\rm{ }} \vdots \\ {a_{n1}}{\rm{ }}{a_{n2}} \cdots {\rm{ }}{a_{nm}} \end{array} \right]$ ${\bf{x}} = {({x_1},{\rm{ }}{x_2}, \cdots ,{x_m})^T}$ ${\bf{b}} = {({b_1},{\rm{ }}{b_2}, \cdots ,{b_n})^T}$

则正规方程组为

${A^T}A{\bf{x}} = {A^T}{\bf{b}}$

解之，得到最优近似解。

三、用多项式作最小二乘曲线拟合

取基函数为 ${\varphi _0}(x) = 1,{\rm{ }}{\varphi _1}(x) = x,{\rm{ }}{\varphi _2}(x) = {x^2},{\rm{ }} \cdots ,{\rm{ }}{\varphi _m}(x) = {x^m}$ ，则拟合多项式为

$P(x) = {a_0} + {a_1}x + {a_2}{x^2} + \cdots + {a_m}{x^m}\mathop {}\nolimits_{} \mathop {}\nolimits_{} \mathop {}\nolimits_{} (m < n - 1)$

对照上面的内容我们可以得到矛盾方程 $A{\bf{\alpha }} = {\bf{y}}$ ，其中

${\bf{\alpha }} = {\left( {{a_0},{a_1}, \cdots ,{a_m}} \right)^T}$ ， ${\bf{y}} = {\left( {{y_1},{y_2}, \cdots ,{y_n}} \right)^T}$ ， $A = \left[ {\begin{array}{*{20}{c}} 1&{{x_1}}&{{x_1}^2}& \cdots &{{x_1}^m}\\ 1&{{x_2}}&{{x_2}^2}& \cdots &{{x_2}^m}\\ {}& \cdots & \cdots & \cdots &{}\\ 1&{{x_n}}&{{x_n}^2}& \cdots &{{x_n}^m} \end{array}} \right]$

对应的正规方程组为： ${A^T}A{\bf{\alpha }} = {A^T}{\bf{y}}$ ，方程的解 $\alpha$ 即为多项式系数。

正规方程的建立

观察到

$\begin{array}{c} {A^T}A = \left[ {\begin{array}{*{20}{c}} 1&1&1& \cdots &1\\ {{x_1}}&{{x_2}}&{{x_3}}& \cdots &{{x_n}}\\ {{x_1}^2}&{{x_2}^2}&{{x_3}^2}& \cdots &{{x_n}^2}\\ \vdots &{}& \vdots &{}& \vdots \\ {{x_1}^m}&{{x_2}^m}&{{x_3}^m}& \cdots &{{x_n}^m} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} 1&{{x_1}}&{{x_1}^2}& \cdots &{{x_1}^m}\\ 1&{{x_2}}&{{x_2}^2}& \cdots &{{x_2}^m}\\ 1&{{x_3}}&{{x_3}^2}& \cdots &{{x_3}^m}\\ \vdots &{}& \vdots &{}& \vdots \\ 1&{{x_n}}&{{x_n}^2}& \cdots &{{x_n}^m} \end{array}} \right]\\ = \left[ {\begin{array}{*{20}{c}} n&{\sum\limits_{i = 1}^n {{x_i}} }&{\sum\limits_{i = 1}^n {{x_i}^2} }& \cdots &{\sum\limits_{i = 1}^n {{x_i}^m} }\\ {\sum\limits_{i = 1}^n {{x_i}} }&{\sum\limits_{i = 1}^n {{x_i}^2} }&{\sum\limits_{i = 1}^n {{x_i}^3} }& \cdots &{\sum\limits_{i = 1}^n {{x_i}^{m + 1}} }\\ \cdots &{}& \cdots &{}& \cdots \\ {\sum\limits_{i = 1}^n {{x_i}^m} }&{\sum\limits_{i = 1}^n {{x_i}^{m + 1}} }&{\sum\limits_{i = 1}^n {{x_i}^{m + 2}} }& \cdots &{\sum\limits_{i = 1}^n {{x_i}^{2m}} } \end{array}} \right]\\ \end{array}$

其是对称正定的，因此，求解正规方程只需要计算 $n,\sum\limits_{i = 1}^n {{x_i}} ,{\rm{ }}\sum\limits_{i = 1}^n {{x_i}^2} , \cdots ,{\rm{ }}\sum\limits_{i = 1}^n {{x_i}^m} ,{\rm{ }}\sum\limits_{i = 1}^n {{x_i}^{m + 1}} , \cdots ,{\rm{ }}\sum\limits_{i = 1}^n {{x_i}^{2m}}$ 和 ${A^T}{\bf{y}} = \left[ \begin{array}{l} \sum\limits_{i = 1}^n {{y_i}} \\ \sum\limits_{i = 1}^n {{y_i}{x_i}} \\ \sum\limits_{i = 1}^n {{y_i}{x_i}^2} \\ {\rm{ }} \vdots \\ \sum\limits_{i = 1}^n {{y_i}{x_i}^m} \end{array} \right]$

因此总的步骤为：

1.计算正规方程组的系数矩阵和常数项各元素

$\begin{array}{l} \sum\limits_{i = 1}^n {{x_i}^0 = } n,{\rm{ }}\sum\limits_{i = 1}^n {{x_i}} ,{\rm{ }}\sum\limits_{i = 1}^n {{x_i}^2} ,{\rm{ }} \cdots ,{\rm{ }}\sum\limits_{i = 1}^n {{x_i}^{2m}} \\ \sum\limits_{i = 1}^n {{y_i}} ,{\rm{ }}\sum\limits_{i = 1}^n {{y_i}{x_i}} ,{\rm{ }}\sum\limits_{i = 1}^n {{y_i}{x_i}^2} ,{\rm{ }} \cdots ,{\rm{ }}\sum\limits_{i = 1}^n {{y_i}{x_i}^m} \end{array}$

2. 利用迭代法或高斯列主元法求正规方程组的解

${a_0},{a_1}, \cdots ,{a_m}$