声音传播的公式最小二乘法公式

2020年9月25日发(作者：祝汝佐)
最小二乘法公式
历史
编辑

1801年，意大利天文学家朱赛普 ·皮亚齐发现了第一颗小行星谷神星。经过40天的跟
踪观测后，由于谷神星运行至太阳背后，使得皮亚 齐失去了谷神星的位置。随后全世界的科
学家利用皮亚齐的观测数据开始寻找谷神星，但是根据大多数人 计算的结果来寻找谷神星
都没有结果。时年24岁的高斯也计算了谷神星的轨道。奥地利天文学家海因里 希·奥尔伯斯
根据高斯计算出来的轨道重新发现了谷神星。
高斯使用的最小二乘法的方法发表于1809年他的著作《天体运动论》中。
法国科学家勒让德于1806年独立发明“最小二乘法”，但因不为世人所知而默默无闻。


二乘法(2张)

勒让德曾与高斯为谁最早创立最小二乘法原理发生争执。
1829年，高斯提供了最小二乘法 的优化效果强于其他方法的证明，因此被称为高斯-
马尔可夫定理。（来自于wikipedia） [1]

适用领域
编辑

代数，数学学科
线性最小二乘的基本公式
编辑

考虑超定方程组（超定指未知数小于方程个数）：

其中m代表有m个等式，n代表有 n 个未知数

，m>n ；将其进行向量化后为：




，

，

显然该方程组一般而言没有解，所以为了选取最合适的

让该等式尽量成立，引入残差平方和函数S

（在统计学中，残差平方和函数可以看成n倍的均方误差MSE）
当

时，

取最小值，记作：

通过对

进行微分求最值，可以得到：

如果矩阵

非奇异则

有唯一解
[1]
：

原理
编辑

在我们研究两个变量（x,y）之间的相互关系时，通常可以得到一系列成对的数据
（x1,y1.x 2,y2... xm,ym）；将这些数据描绘在x -y直角坐标系中，若发现这些点在一条直
线附近，可以令这条直线方程如（式1-1）。

（式1-1）
其中：a0、a1 是任意实数
为建立这直线方程就要确定a0和 a1，应用《最小二乘法原理》，将实测值Yi与利用
计算值Yj（Yj=a0+a1Xi）（式1-1 ）的离差（Yi-Yj）的平方和

最小为“优化判据”。
令：φ =

（式1-2)
把（式1-1）代入（式1-2）中得：
φ =

（式1-3)
当

最小时，可用函数 φ 对a0、a1求偏导数，令这两个偏导数等于零。
∑2(a0 + a1*Xi - Yi）=0（式1-4)
∑2Xi（a0 +a1*Xi - Yi）=0（式1-5)
亦即：
na0 + （∑Xi ) a1 = ∑Yi （式1-6)
（∑Xi ) a0 + （∑Xi^2 ) a1 = ∑（Xi*Yi) （式1-7)
得到的两个关于a0、 a1为未知数的两个方程组，解这两个方程组得出：
a0 = （∑Yi) n - a1（∑Xi) n （式1-8)
a1 = [n∑(Xi Yi) - （∑Xi ∑Yi)] (n∑Xi^2 -∑Xi∑Xi)（式1-9)
这时把a0、a1代入（式1-1）中， 此时的(式1-1）就是我们回归的一元线性方程即：
数学模型。
在回归过程中，回归的关联式不可能全部通过每个回归数据点（x1,y1. x2,y2...xm, ym），
为了判断关联式的好坏，可借助相关系数“R”，统计量“F”，剩余标准偏差“S”进行判断 ；“R”
越趋近于 1 越好；“F”的绝对值越大越好；“S”越趋近于 0 越好。
R = [∑XiYi - m （∑Xi m）（∑Yi m)] SQR{[∑Xi2 - m （∑Xi m)2][∑Yi2 - m （∑Yi
m)2]} （式1-10) *
在（式1-1 0）中，m为样本容量，即实验次数；Xi、Yi分别为任意一组实验数据X、
Y的数值。
[1]

方法
编辑

以最简单的一元线性模型来解释最小二乘法。 什么是一元线性模型呢？监督学习中，
如果预测的变量是离散的，我们称其为分类（如决策树，支持向量 机等），如果预测的变量
是连续的，我们称其为回归。回归分析中，如果只包括一个自变量和一个因变量 ，且二者的
关系可用一条直线近似表示，这种回归分析称为一元线性回归分析。如果回归分析中包括两< br>个或两个以上的自变量，且因变量和自变量之间是线性关系，则称为多元线性回归分析。对
于二维 空间线性是一条直线；对于三维空间线性是一个平面，对于多维空间线性是一个超平
面。
[3]

对于一元线性回归模型, 假设从总体中获取了n组观察值（X1，Y1），（X2，Y2）， …，
（Xn，Yn）。对于平面中 的这n个点，可以使用无数条曲线来拟合。要求样本回归函数尽
可能好地拟合这组值。综合起来看，这条 直线处于样本数据的中心位置最合理。 选择最佳
拟合曲线的标准可以确定为：使总的拟合误差（即总残 差）达到最小。有以下三个标准可以
选择：
[3]

（1）用“残差和最小 ”确定直线位置是一个途径。但很快发现计算“残差和”存在相互抵
消的问题。
（2）用“残差绝对值和最小”确定直线位置也是一个途径。但绝对值的计算比较麻烦。
（3 ）最小二乘法的原则是以“残差平方和最小”确定直线位置。用最小二乘法除了计算
比较方便外，得到的 估计量还具有优良特性。这种方法对异常值非常敏感。
[3]

最常用的是普通最小二乘法（ Ordinary Least Square，OLS）：所选择的回 归模型应
该使所有观察值的残差平方和达到最小。（Q为残差平方和）- 即采用平方损失函数。
[3]

样本回归模型：

其中e
i
为样本（X
i,
Y
i
）的误差。
[3]

平方损失函数：

则通过Q最小确定这条直线，即确定β0和β1 ，把它们看作是Q的函数，就变成了一
个求极值的问题，可以通过求导数得到。求Q对两个待估参数的偏 导数：
[3]


根据数学知识我们知道，函数的极值点为偏导为0的点。
[3]

解得：

这就是最小二乘法的解法，就是求得平方损失函数的极值点。
[3]

公式
编辑


拟合
编辑

对给定数据点集合

，在取定的函数类

中，求

，使误差的平方和

最小，

。从几何意义上讲，就是寻求与给定点集

的距离平方和为最小的曲线y=p(x) 。函数p(x)称为拟合函数或最小二乘解，求拟合函
数p(x)的方法称为曲线拟合的最小二乘法。< br> [1]

最小二乘法的矩阵形式
最小二乘法的矩阵形式为：

其中

为

的矩阵，

为

的列向量，

为

的列向量。如果

（方程的个数大于未知量的个数），这个方程系统称为矛盾方程组（Over Determined
System），如果

（方程的个数小于未知量的个数），这个系统就是Under Determined System。
正常来看，这个方程是没有解的，但在数值计算领域，我们通常是计算

，解出其中的

。比较直观的做法是求解

，但通常比较低效。其中一种常见的解法是对

进行QR分解（

），其中

是

正交矩阵（Orthonormal Matrix），

是

上三角矩阵（Upper Triangular Matrix），则有

用MATLAB命令
1
x=R(Qb)

可解得

。
[1]

最小二乘法的Matlab实现
① 一次函数线性拟合使用polyfit（x,y,1）
②多项式函数线性拟合使用 polyfit（x,y,n），n为次数
拟合曲线
x=[0.5,1.0,1.5,2.0,2.5,3.0]，
y=[1.75,2.45,3.81,4.80,7.00,8.60]。
解：MATLAB程序如下：
1
2
3
4
5
6
x=[0.5,1.0,1.5,2.0,2.5,3.0];

y=[1.75,2.45,3.81,4.80,7.00,8.60];

p=polyfit(x,y,2)

x1=0.5:0.5:3.0;

y1=polyval(p,x1);

plot(x,y,'*r',x1,y1,'-b')

计算结果为：
1
p =0.5614 0.8287 1.1560

即所得多项式为y=0.5614x^2+0.8287x+1.15560
③非线性函数使用
1
2
lsqcurvefit(fun,x0,x,y)

a=nlinfit(x,y,fun,b0)

最小二乘法在交通运输学中的运用
交通发生预测的目的是建立分区产生的交通量与分区土地利用、社会经济特征等变量
之间的定量 关系，推算规划年各分区所产生的交通量。因为一次出行有两个端点，所以我们
要分别分析一个区生成的 交通和吸引的交通。交通发生预测通常有两种方法：回归分析法和
聚类分析法。
[1]

回归分析法是根据对因变量与一个或多个自变量的统计分析，建立因变量和自变量的< br>关系，最简单的情况就是一元回归分析，一般式为：Y=α+βX式中Y是因变量，X是自变量，
α和β是回归系数。若用上述公式预测小区的交通生成，则以下标 i 标记所有变量；如果
用它研究分区交通吸引，则以下标 j 标记所有变量。而运用公式的过程中需要利 用最小二
乘法来求解，上述公式中的回归系数根据最小二乘法可得：
其中，式中的X拔是规划 年的自变量值，Y拔是规划年分区交通生成（或吸引）预测
值。
[1]

课题
编辑

从前面的学习中, 我们知道最小二乘法可以用来处理一组数据, 可以从一组测定的数
据中寻求变量之间的依赖关系, 这种函数关系称为经验公式. 本课题将介绍最小二乘法的
精确定义及如何寻求点与点之间近似成线性关系时的经验公式. 假定实验测得变量之间的
n个数据, 则在 平面上, 可以得到 n个点 , 这种图形称为“散点图”, 从图中可以粗略看出
这些点大致散落在某直线近旁, 我们认为 与 之间近似为一线性函数, 下面介绍求解步骤.
考虑函数 y=a+bx, 其中a 和 b是待定常数. 如果离散点完全的在一直线上,可以认为
变量之间的关系为一元函数 . 但一般说来, 这些点不可能在同一直线上. 但是它只能用直
线来描述时, 计算值与实际值会产生偏差. 当然要求偏差越小越好, 但由于偏差可正可负,
因此不能认为总偏差 时, 拟合函数很好地反映了变量之间的关系,但是因为此时每个偏差
的绝对值可能很大. 为了改进这一缺陷, 就考虑用平均值来代替 . 但是由于绝对值不易作
解析运算, 因此, 进一步用残差平方和函数来度量总偏差. 偏差的平方和最小可以保证每个
偏差都不会很大. 于是问题归结为确定拟合函数中的常数和使残差平方和函数最小. 通过
这种方法确定系数的方法称为最小二乘法.
由极值原理得 , 即
解此联立方程得
(*)
问题 I 为研究某一化学反应过程中, 温度 ℃)对产品得率 (%)的影响, 测得数据如下:
温度 ℃)
100 110 120 130 140 150 160 170 180 190
得率 (%)
45 51 54 61 66 70 74 78 85 89
(1) 利用“ListPlot”函数, 绘出数据 的散点图(采用格式: ListPlot[{ , , …, },
Prolog->AbsolutePointSize[3]] );
(2) 利用“Line”函数, 将散点连接起来, 注意观察有何特征? (采用格式:
Show[Graphics[Line[{ , , …, }]] , Axes->True ])
(3) 根据公式(*), 利用“Apply”函数及集合的有关运算编写一个小的程序, 求经验公
式
(程序编写思路为: 任意给定两个集合A (此处表示温度)、B(此处表示得率), 由公式(*)
可定义两个二元函数(集合A和B为其变量)分别表示 和 . 集合A元素求和: Apply[Plus,A]
表示将加法施加到集合A上, 即各元素相加, 例如Apply[Plus,{1,2,3}]=6;Length[A]表示集合
A 元素的个数, 即为n; A.B表示两集合元素相乘相加;A*B表示集合A与B元素对应相乘得
到的新的集合.)
(4) 在同一张图中显示直线 及散点图;
(5) 估计温度为200时产品得率.
然而, 不少实际问题的观测数据 , , …, 的散点图明显地不能用线性关系来描叙, 但确
实散落在某一曲线近旁, 这时可以根据散点图的轮廓和实际经验, 选一条曲线来近似表达
与 的相互关系.
问题 II 下表是美国旧轿车价格的调查资料, 今以 表示轿车的使用年数, (美元)表示相
应的平均价格, 求 与 之间的关系.
使用年数
1 2 3 4 5 6 7 8 9 10
平均价格
2651 1943 1494 1087 765 538 484290 226 204
(1) 利用“ListPlot”函数绘出数据 的散点图, 注意观察有何特征?
(2) 令 , 绘出数据 的散点图, 注意观察有何特征?
(3) 利用“Line”函数, 将散点 连接起来, 说明有何特征?
(4) 利用最小二乘法, 求 与 之间的关系;
(5) 求 与 之间的关系;
(6) 在同一张图中显示散点图及 关于 的图形.
思考与练习
编辑

1. 假设一组数据 : , , …, 变量之间近似成线性关系, 试利用集合的有关运算, 编写一
简单程序: 对于任意给定的数据集合 , 通过求解极值原理所包含的方程组, 不需要给出 、
计算的表达式, 立即得到 、 的值, 并就本课题 I (3)进行实验.
注: 利用Transpose函数可以得到数据A的第一个分量的集合, 命令格式为:
先求A的转置, 然后取第一行元素, 即为数据A的第一个分量集合, 例如
(A即为矩阵 )
= (数据A的第一个分量集合)
= (数据A的第二个分量集合)
B-C表示集合B与C对应元素相减所得的集合, 如 = .
2. 最小二乘法在数学上称为曲线拟合, 请使用拟合函数“Fit”重新计算 与 的值, 并与
先前的结果作一比较.
注: Fit函数使用格式:
设变量为x, 对数据A进行线性拟合, 如对题1中的A拟合函数为:
实例
编辑

数据编号
实验次数w
x
y
1
2
0.1
1.1
2
1
0.2
1.9
3
1
0.3
3.1
4
1
0.4
3.9
要拟合得到形如y = a + b x 的函数，求解函数中系数的方程组为


其中，

为权重，对应每个实验点的实验次数，4个实验点只有第一个点重复做了 一次且得到
相同结果（如果结果不同则另算一个实验点），其它都没有重复实验，因此总次数为5次。






解得

故拟合方程为