前言
最小二乘法是一种常用的数学方法,用于解决数据拟合问题。在数据拟合问题中,我们需要找到一条曲线或者一个平面来拟合给定的数据。最小二乘法的主要思想是寻找一个拟合函数,使得该函数与所有观测值之间的差的平方和最小。最小二乘法的背景可以追溯到18世纪欧洲的天文学家和数学家吕格朗日,他首次利用最小二乘法求解天体轨迹问题。此后,最小二乘法被广泛应用于数据拟合、统计分析、信号处理、图像处理等领域,并成为现代数学和应用数学中不可或缺的工具。
最小二乘原理是指在数学上寻找误差平方和最小的解决方案,即将所有观测值与拟合值之间的差的平方和最小化。这个原理在统计学中被广泛应用,用于确定回归模型中的系数。
最小二乘估计是根据最小二乘原理,通过具体的计算方法来求解回归模型中的系数,是一种常用的参数估计方法。最小二乘估计可以用于线性回归、多项式回归、非线性回归等各种回归模型的参数估计。它可以通过求解矩阵方程或者最小化目标函数的方法来得到系数的估计值。因此,最小二乘估计是通过具体的计算方法来实现最小二乘原理的应用。
本文由最小二乘原理引入到最小二乘估计
一、最小二乘原理
1.一维例子
假设测量某个物体的长度,
n
n
n个人测得的长度分别为
x
1
,
x
2
,
.
.
.
x
n
x_1,x_2,...x_n
x1,x2,...xn,然后通过这
n
n
n个数据得到该物体最接近真实值的长度,可以按照统计学原理,定义误差平方和:
f
(
x
)
=
(
x
1
−
x
)
2
+
(
x
2
−
x
)
2
+
.
.
.
+
(
x
n
−
x
)
2
.
f(x)=(x_1-x)^2+(x_2-x)^2+...+(x_n-x)^2.
f(x)=(x1−x)2+(x2−x)2+...+(xn−x)2.为了使其误差平方和最小,及求
f
(
x
)
f(x)
f(x)的最小值,令
f
(
x
)
f(x)
f(x)的导数为零,得到
f
′
(
x
)
=
−
2
(
x
1
−
x
)
−
2
(
x
2
−
x
)
−
.
.
.
−
2
(
x
n
−
x
)
=
0.
f'(x)=-2(x_1-x)-2(x_2-x)-...-2(x_n-x)=0.
f′(x)=−2(x1−x)−2(x2−x)−...−2(xn−x)=0.求解可以得到
x
=
x
1
+
x
2
+
.
.
.
+
x
n
n
=
x
ˉ
.
x=\frac{x_1+x_2+...+x_n}{n}=\bar{x}.
x=nx1+x2+...+xn=xˉ.由于
f
′
′
(
x
)
=
2
n
>
0
f''(x)=2n>0
f′′(x)=2n>0,可知
f
(
x
ˉ
)
f(\bar{x})
f(xˉ)是误差函数
f
(
x
)
f(x)
f(x)的最小值。这就是最小二乘的基本原理。
2.二维例子
现有
n
n
n个点
(
x
1
,
y
1
)
,
(
x
2
,
y
2
)
,
.
.
.
,
(
x
n
,
y
n
)
(x_1,y_1),(x_2,y_2),...,(x_n,y_n)
(x1,y1),(x2,y2),...,(xn,yn)很接近某条直线。假设这条待定的直线为
y
=
k
x
+
b
.
y=kx+b.
y=kx+b.
点
(
x
1
,
y
1
)
,
(
x
2
,
y
2
)
,
.
.
.
,
(
x
n
,
y
n
)
(x_1,y_1),(x_2,y_2),...,(x_n,y_n)
(x1,y1),(x2,y2),...,(xn,yn)不可能都满足
y
=
k
x
+
b
y=kx+b
y=kx+b,因此,目标是估计参数
a
,
b
a,b
a,b,使得每个点到该直线的距离和为最小值,由于点到直线的距离公式计算复杂,在这里采用点到直线的纵坐标差
ε
\varepsilon
ε,方便计算。
ε
j
=
y
j
−
k
x
j
−
b
\varepsilon_j =y_j-kx_j-b
εj=yj−kxj−b
目标函数为
J
1
(
k
,
b
)
=
∑
j
=
1
n
ε
j
2
=
(
y
1
−
k
x
1
−
b
)
2
+
(
y
2
−
k
x
2
−
b
)
2
+
.
.
.
+
(
y
n
−
k
x
n
−
b
)
2
J_1(k,b)=\sum_{j=1}^{n}\varepsilon _j^2=(y_1-kx_1-b)^2+(y_2-kx_2-b)^2+...+(y_n-kx_n-b)^2
J1(k,b)=j=1∑nεj2=(y1−kx1−b)2+(y2−kx2−b)2+...+(yn−kxn−b)2
现在要使目标函数
J
1
(
k
,
b
)
J_1(k,b)
J1(k,b)最小,这就构成了最小二乘拟合问题。
分别对
k
k
k和
b
b
b求偏导,并令其为零得到
∂
J
1
(
k
,
b
)
∂
k
=
−
2
x
1
(
y
1
−
k
x
1
−
b
)
−
2
x
2
(
y
2
−
k
x
2
−
b
)
−
.
.
.
−
2
x
n
(
y
n
−
k
x
n
−
b
)
=
0
,
∂
J
1
(
k
,
b
)
∂
b
=
−
2
(
y
1
−
k
x
1
−
b
)
−
2
(
y
2
−
k
x
2
−
b
)
−
.
.
.
−
2
(
y
n
−
k
x
n
−
b
)
=
0.
\begin{aligned}\frac{\partial J_1(k,b)}{\partial k}&=-2x_1(y_1-kx_1-b)-2x_2(y_2-kx_2-b)-...-2x_n(y_n-kx_n-b)=0,\\ \frac{\partial J_1(k,b)}{\partial b}&=-2(y_1-kx_1-b)-2(y_2-kx_2-b)-...-2(y_n-kx_n-b)=0.\end{aligned}
∂k∂J1(k,b)∂b∂J1(k,b)=−2x1(y1−kx1−b)−2x2(y2−kx2−b)−...−2xn(yn−kxn−b)=0,=−2(y1−kx1−b)−2(y2−kx2−b)−...−2(yn−kxn−b)=0.或
(
∑
j
=
1
n
x
j
2
)
k
+
(
∑
j
=
1
n
x
j
)
b
=
∑
j
=
1
n
x
j
y
j
,
(
∑
j
=
1
n
x
j
)
k
+
n
b
=
∑
j
=
1
n
y
j
.
\begin{aligned}(\sum_{j=1}^{n}x_j^2)k+(\sum_{j=1}^{n}x_j)b&=\sum_{j=1}^{n}x_jy_j ,\\ (\sum_{j=1}^{n}x_j)k+nb&=\sum_{j=1}^{n}y_j .\end{aligned}
(j=1∑nxj2)k+(j=1∑nxj)b(j=1∑nxj)k+nb=j=1∑nxjyj,=j=1∑nyj.因此求得的参数拟合值为
[
k
b
]
=
[
∑
j
=
1
n
x
j
2
∑
j
=
1
n
x
j
∑
j
=
1
n
x
j
n
]
−
1
[
∑
j
=
1
n
x
j
y
j
∑
j
=
1
n
y
j
]
.
\begin{bmatrix} k \\ b \end{bmatrix}=\begin{bmatrix} \sum_{j=1}^{n}x_j^2 & \sum_{j=1}^{n}x_j\\ \sum_{j=1}^{n}x_j & n \end{bmatrix}^{-1}\begin{bmatrix} \sum_{j=1}^{n}x_jy_j \\ \sum_{j=1}^{n}y_j \end{bmatrix}.
[kb]=[∑j=1nxj2∑j=1nxj∑j=1nxjn]−1[∑j=1nxjyj∑j=1nyj].
3.最小二乘估计
为了处理高维(大于两个未知参数)线性拟合问题,定义由观测数据构成的信息向量
φ
j
\pmb{\varphi_j}
φj和拟合模型的参数向量
θ
\pmb{\theta}
θ如下,
φ
j
=
[
x
j
1
]
,
θ
=
[
k
b
]
\pmb{\varphi_j}=\begin{bmatrix} x_j \\ 1\end{bmatrix},\pmb{\theta}=\begin{bmatrix} k\\ b\end{bmatrix}
φj=[xj1],θ=[kb]可以得到
y
j
=
k
x
j
+
b
+
ε
j
=
[
x
j
1
]
[
k
b
]
+
ε
j
=
φ
j
T
θ
+
ε
j
.
y_j=kx_j+b+\varepsilon_j=\begin{bmatrix} x_j & 1\end{bmatrix}\begin{bmatrix} k \\ b\end{bmatrix}+\varepsilon_j=\pmb{\varphi_j}^T\pmb{\theta}+\varepsilon_j.
yj=kxj+b+εj=[xj1][kb]+εj=φjTθ+εj.上式称为线性回归模型,在系统辨识中称为辨识模型,也是最小二乘格式。偏差
ε
i
\varepsilon_i
εi称为噪声或随机干扰。
下面求解向量形式的目标函数,将
J
1
(
k
,
b
)
J_1(k,b)
J1(k,b)写为向量形式
J
2
(
θ
)
J_2(\pmb \theta)
J2(θ),可得到
J
2
(
θ
)
=
∑
j
=
1
n
ε
j
2
=
∑
j
=
1
n
(
y
j
−
φ
j
T
θ
)
2
J_2(\pmb \theta)=\sum_{j=1}^{n}\varepsilon _j^2=\sum_{j=1}^{n}(y_j-\pmb{\varphi_j}^T\pmb{\theta})^2
J2(θ)=j=1∑nεj2=j=1∑n(yj−φjTθ)2
对
θ
\pmb{\theta}
θ求偏导并令其为零可得
∂
J
2
(
θ
)
∂
θ
=
−
2
∑
j
=
1
n
φ
j
(
y
j
−
φ
j
T
θ
)
=
0
\frac{\partial J_2(\pmb \theta)}{\partial \pmb \theta}=-2\sum_{j=1}^{n}\pmb{\varphi_j}(y_j-\pmb{\varphi_j}^T\pmb{\theta})=0
∂θ∂J2(θ)=−2j=1∑nφj(yj−φjTθ)=0
则
∑
j
=
1
n
φ
j
y
j
−
(
∑
j
=
1
n
φ
j
φ
j
T
)
θ
=
0
\sum_{j=1}^{n}\pmb{\varphi_j}y_j-(\sum_{j=1}^{n}\pmb{\varphi_j}\pmb{\varphi_j}^T)\pmb{\theta}=0
j=1∑nφjyj−(j=1∑nφjφjT)θ=0
故
θ
\pmb{\theta}
θ的最小二乘估计为
θ
^
=
(
∑
j
=
1
n
φ
j
φ
j
T
)
−
1
∑
j
=
1
n
φ
j
y
j
.
\hat{\pmb{\theta}}=(\sum_{j=1}^{n}\pmb{\varphi_j}\pmb{\varphi_j}^T)^{-1}\sum_{j=1}^{n}\pmb{\varphi_j}y_j.
θ^=(j=1∑nφjφjT)−1j=1∑nφjyj.
现在定义
Y
n
=
[
y
1
y
2
.
.
.
y
n
]
∈
R
n
,
H
n
=
[
φ
1
T
φ
2
T
.
.
.
φ
n
T
]
∈
R
2
n
,
ε
i
=
[
ε
1
ε
2
.
.
.
ε
n
]
∈
R
n
.
\pmb{Y_n}=\begin{bmatrix} y_1\\ y_2\\...\\y_n\end{bmatrix}\in R^n,\pmb{H_n}=\begin{bmatrix} \pmb{\varphi_1}^T\\ \pmb{\varphi_2}^T\\...\\\pmb{\varphi_n}^T\end{bmatrix}\in R^{2n},\pmb{\varepsilon_i}=\begin{bmatrix} \varepsilon_1\\ \varepsilon_2\\...\\\varepsilon_n\end{bmatrix}\in R^n.
Yn=
y1y2...yn
∈Rn,Hn=
φ1Tφ2T...φnT
∈R2n,εi=
ε1ε2...εn
∈Rn.由于
H
n
T
H
n
=
[
φ
1
φ
2
.
.
.
φ
n
]
[
φ
1
T
φ
2
T
.
.
.
φ
n
T
]
=
∑
j
=
1
n
φ
j
φ
j
T
,
H
n
T
Y
n
=
[
φ
1
φ
2
.
.
.
φ
n
]
[
y
1
y
2
.
.
.
y
n
]
=
∑
j
=
1
n
φ
j
y
j
,
\pmb{H_n}^T\pmb{H_n}=\begin{bmatrix} \pmb{\varphi_1}& \pmb{\varphi_2}&...&\pmb{\varphi_n}\end{bmatrix}\begin{bmatrix} \pmb{\varphi_1}^T\\ \pmb{\varphi_2}^T\\...\\\pmb{\varphi_n}^T\end{bmatrix}=\sum_{j=1}^{n}\pmb{\varphi_j}\pmb{\varphi_j}^T,\\ \pmb{H_n}^T\pmb{Y_n}=\begin{bmatrix} \pmb{\varphi_1}& \pmb{\varphi_2}&...&\pmb{\varphi_n}\end{bmatrix}\begin{bmatrix} y_1\\ y_2\\...\\y_n\end{bmatrix}=\sum_{j=1}^{n}\pmb{\varphi_j}y_j,
HnTHn=[φ1φ2...φn]
φ1Tφ2T...φnT
=j=1∑nφjφjT,HnTYn=[φ1φ2...φn]
y1y2...yn
=j=1∑nφjyj,所以最小二乘估计最终可以写为
θ
^
=
(
∑
j
=
1
n
φ
j
φ
j
T
)
−
1
∑
j
=
1
n
φ
j
y
j
=
(
H
n
T
H
n
)
−
1
H
n
T
Y
n
=
θ
^
(
n
)
\hat{\pmb{\theta}}=(\sum_{j=1}^{n}\pmb{\varphi_j}\pmb{\varphi_j}^T)^{-1}\sum_{j=1}^{n}\pmb{\varphi_j}y_j=(\pmb{H_n}^T\pmb{H_n})^{-1}\pmb{H_n}^T\pmb{Y_n}=\hat{\pmb{\theta}}(n)
θ^=(j=1∑nφjφjT)−1j=1∑nφjyj=(HnTHn)−1HnTYn=θ^(n)
参考文献
丁锋.系统辨识新论[M].北京:科学出版社,2013:199.