最小二乘拟合方法

本页目录

曲线拟合工具箱在拟合数据时使用最小二乘法。进行拟合需要一个能表示响应数据与预测数据之间关系的模型，模型中有一些需要确定的系数。拟合的目的就是要获得对这些未知系数的估计。[大谦MATLAB，dqmatlab点com]

最小二乘法通过最小化残差的平方和来获得待定系数的估计。第i个数据点的残差定义为测量响应值yi和拟合响应值\({\hat{y}}_{i}\)之间的差值，即

\[{r}_{i}={y}_{i}-{\hat{y}}_{i}\]

残差的平方和用下式给定

\[S=\sum_{i=1}^{n} {{r}_{i}^{2}}=\sum_{i=1}^{n} {({y}_{i}-{\hat{y}}_{i}{)}^{2}}\]

其中，n为参与拟合的数据点的个数，S为误差估计的平方和。支持的最小二乘拟合类型包括：

线性最小二乘

加权线性最小二乘

稳健最小二乘

非线性最小二乘

线性最小二乘

曲线拟合工具箱用线性模型拟合数据时采用的是线性最小二乘方法。线性模型用一个线性方程来定义。为了演示线性最小二乘拟合，假设用1阶多项式模型拟合n个数据点，模型为

y=p1x+p2

为了求解方程中的未知系数p1和p2，构造一个有2个未知参数n个线性方程的系统S。如果n大于未知参数的个数，则等式系统是超定的。

\[S=\sum_{i=1}^{n} {({y}_{i}-({p}_{1}{x}_{i}+{p}_{2}){)}^{2}}\]

因为最小二乘拟合最小化残差的平方和，各系数可以通过对各系数求偏导数并使偏导数等于0来求得。

\[\begin{matrix} \frac{\partial S}{\partial{p}_{1}}=-2\sum_{i=1}^{n} {{x}_{i}({y}_{i}-({p}_{1}{x}_{i}+{p}_{2}))=0} \\ \frac{\partial S}{\partial{p}_{2}}=-2\sum_{i=1}^{n} {({y}_{i}-({p}_{1}{x}_{i}+{p}_{2}))=0} \end{matrix}\]

用b1和b2替换p1和p2，前面的方程变为

\[\begin{matrix} \sum {{x}_{i}({y}_{i}-({b}_{1}{x}_{i}+{b}_{2}))=0} \\ \sum {({y}_{i}-({b}_{1}{x}_{i}+{b}_{2}))=0} \end{matrix}\]

其中，i从1变到n。移项以后，可以写成下面的标准形式：

\[\begin{matrix} {b}_{1}\sum {{x}_{i}^{2}+{b}_{2}\sum {{x}_{i}=\sum {{x}_{i}{y}_{i}}}} \\ {b}_{1}\sum {{x}_{i}+n{b}_{2}=\sum {{y}_{i}}} \end{matrix}\]

求解b1

\[{b}_{1}=\frac{n\sum {{x}_{i}{y}_{i}-\sum {{x}_{i}\sum {{y}_{i}}}}}{n\sum {{x}_{i}^{2}-(\sum {{x}_{i}}{)}^{2}}}\]

用b1的值求解b2

\[{b}_{2}=\frac{1}{n}(\sum {{y}_{i}-{b}_{1}\sum {{x}_{i}}})\]

如上所示。求系数p1和p2只需要简单的几步计算。

加权线性最小二乘

通常假设响应数据的方差不变，如果违背该假设，拟合会更容易受到数据质量的影响。可以用加权最小二乘回归来改进拟合结果，进行这种回归时，在拟合时添加额外的权重因子。加权最小二乘回归最小化下面的误差估计

\[S=\sum_{i=1}^{n} {{w}_{i}({y}_{i}-{\hat{y}}_{i}{)}^{2}}\]

其中，wi为权重。权重确定每个响应值对最终的参数估计有多大影响。高级别的数据点对拟合的影响力比低级别数据点的大。如果权重已知，建议使用本方法。

通过拟合数据和绘残差图，常常可以确定方差是否不变。在图1-1中，数据包含不同质量的点，并假设拟合正确。从残差图(图1-2)中可以发现低质量的数据，它们对应于小预测值生成的散点数据比大预测值的还大的情况。

Document Image
\[\]

图1-1 数据的拟合曲线

Document Image
\[\]

图1-2 拟合残差图

提供的权重应该把响应方差变换到一个常数值。如果知道数据的方差，则权重用下式给定

\[{w}_{i}=1/{\sigma}^{2}\]

如果不知道方差，可以用类似下面的方程求权重的近似值。

\[{w}_{i}={\left( \frac{1}{n}\sum_{i=1}^{n} {({y}_{i}-\bar{y}{)}^{2}} \right)}^{-1}\]

如果数据集包含重复采样数据，本方程会工作得比较好。此时，n为重复采样次数。但是，权重的变化会很大。更好的办法是绘方差的图形并用敏感模型拟合数据。模型的形式不很重要—多项式或幂函数在很多情况下都工作得很好。

稳健最小二乘

最小二乘法的主要缺点是对异常值太敏感。异常值对拟合结果具有很大的影响，因为对残差取平方会将异常值的影响放大。要使异常值的影响最小，可以用稳健最小二乘拟合数据。工具箱提供了下面2个稳健回归方案：

最小绝对残差(LAR)—LAR方案求使残差的绝对差，而不是平方差最小的曲线。所以，极值对拟合的影响更小。

双二次加权—本方案使加权平方和最小，其中每个数据点的权重与该点到拟合线的距离有关。离拟合线更近的点权重更大，离拟合线更远的点权重更小，超出预期范围的点的权重为0。对于大多数情况，双二次加权方案比LAR方案的效果更好，因为它不仅用一般最小二乘法找到了拟合数据的曲线，还使异常值的影响最小。

用双二次加权进行稳健拟合使用的是重复加权的最小二乘迭代算法，按照下面的步骤进行：

用加权最小二乘法拟合模型；

计算调整的残差并进行标准化。调整的残差为

\[{r}_{adj}=\frac{{r}_{i}}{\sqrt{1-{h}_{i}}}\]

ri为一般最小二乘残差，hi为中心化杠杆值。中心化杠杆值通过减小高杠杆值数据点的权重来调整残差。中心化杠杆值很高的数据点对最小二乘拟合具有很大的影响。标准化的调整残差由下式给定。

\[u=\frac{{r}_{adj}}{Ks}\]

K为调和常数，等于4.685，s为MAD/0.6745给定的稳健方差，其中MAD为所有残差中值的绝对差。

用函数u计算稳健性权重。双二次权重由下式给定

\[{w}_{i}=\left\{ \begin{matrix} (1-({u}_{i}{)}^{2}{)}^{2} \\ 0 \end{matrix} \right. \begin{matrix} \left| {u}_{i} \left|\lt1 \\ \left| {u}_{i} \left|\geq1 \end{matrix}\]

注意，如果提供自己的回归权重矢量，最后的权重是稳健性权重和回归权重的乘积。

如果拟合收敛，则完成拟合任务；否则，返回到第1步进行下一次迭代。

非线性最小二乘

曲线拟合工具箱用非线性模型拟合数据时，采用非线性最小二乘公式。非线性模型由非线性方程或线性方程与非线性方程的组合构成。例如幂函数就是非线性的。

按照矩阵形式，非线性模型可以用下面的公式表示：

\[y=f(X,\beta)+\varepsilon\]

其中，y是一个n×1的响应数据矢量；f是β和X的函数；β是m×1的系数矢量；X是n×m的模型设计矩阵；ε是n×1的误差矢量。非线性模型比线性模型更难拟合，因为系数不能用简单的矩阵计算技巧进行计算。一般采用迭代的方法进行计算。

工具箱提供的非线性迭代算法包括：

置信域法—默认算法，如果指定系数约束，必须使用该方法。该方法求解复杂非线性问题的效率比其他几种算法高得多，并且它是普通Levenberg-Marquardt法的改进版本。

Levenberg-Marquardt法—老牌的非线性迭代算法，经过了很多年的考验。如果置信域法得到的结果不理想，并且没有系数约束，试用本法。

Gauss-Newton法—本算法比其他算法快得多，但是它假设数据残差接近0。把该算法放到工具箱中，主要出于教学的目的，解决实际问题时，应该把它作为最后的选择。