1. 概述
线性回归(Linear Regression)是机器学习和统计学中最基础和最广泛应用的模型,是一种对自变量和因变量之间关系进行建模的回归分析1。最小二乘法以误差的平方和作为目标函数,通过求解方程组来获取模型的参数,成为线性回归的普遍方法。
2. 原理步骤
例如有测量数据集
{
(
x
i
,
z
i
)
∣
i
=
1
,
…
,
n
}
\{(x_i,z_i)|i=1,\dots,n\}
{(xi,zi)∣i=1,…,n},观察发现符合抛物线规律
z
=
θ
3
x
2
+
θ
2
x
+
θ
1
\begin{align} z=\theta_3x^2+\theta_2x+\theta_1 \end{align}
z=θ3x2+θ2x+θ1
需要确定系数
θ
3
,
θ
2
,
θ
1
\theta_3,\theta_2,\theta_1
θ3,θ2,θ1,你可以选三组数据代入公式求解,但那并不代表所有数据的分布,需要使用线性回归找到最合适的解。那就把所有点代入抛物线表达式得到方程组:
θ
3
x
1
2
+
θ
2
x
1
+
θ
1
=
z
1
θ
3
x
2
2
+
θ
2
x
2
+
θ
1
=
z
2
⋮
θ
3
x
n
2
+
θ
2
x
n
+
θ
1
=
z
n
\begin{aligned} \theta_3x_{1}^{2}+\theta_2x_1+\theta_1&=z_1 \\ \\ \theta_3x_{2}^{2}+\theta_2x_2+\theta_1&=z_2 \\ \vdots& \\ \theta_3x_{n}^{2}+\theta_2x_{n}+\theta_1&=z_{n} \end{aligned}
θ3x12+θ2x1+θ1θ3x22+θ2x2+θ1⋮θ3xn2+θ2xn+θ1=z1=z2=zn
一般
n
n
n 远大于待求系数的个数,方程组无精确解,写成矩阵的形式,习惯上将高次的系数放后边:
[
1
x
1
x
1
2
1
x
2
x
2
2
⋮
⋮
⋮
1
x
n
x
n
2
]
[
θ
1
θ
2
θ
3
]
≠
[
z
1
z
2
⋮
z
n
]
\begin{align} \begin{bmatrix} 1 & x_1 & x_{1}^{2} \\ 1 & x_2 & x_{2}^{2} \\ \vdots & \vdots & \vdots \\ 1 & x_{n} & x_{n}^{2} \\ \end{bmatrix} \begin{bmatrix} \theta_1\\ \theta_2\\ \theta_3\\ \end{bmatrix} \ne \begin{bmatrix} z_1 \\ z_2\\ \vdots\\ z_{n} \end{bmatrix} \end{align}
11⋮1x1x2⋮xnx12x22⋮xn2
θ1θ2θ3
=
z1z2⋮zn
定义式(2)左边为估计值
Y
Y
Y
Y
=
[
y
1
y
2
⋮
y
n
]
=
[
1
x
1
x
1
2
1
x
2
x
2
2
⋮
⋮
⋮
1
x
n
x
n
2
]
[
θ
1
θ
2
θ
3
]
=
A
⋅
Θ
Y= \begin{bmatrix} y_1 \\ y_2 \\ \vdots \\ y_n \end{bmatrix}= \begin{bmatrix} 1 & x_1 & x_{1}^{2} \\ 1 & x_2 & x_{2}^{2} \\ \vdots & \vdots & \vdots \\ 1 & x_{n} & x_{n}^{2} \\ \end{bmatrix} \begin{bmatrix} \theta_1\\ \theta_2\\ \theta_3\\ \end{bmatrix}=A \sdot \Theta
Y=
y1y2⋮yn
=
11⋮1x1x2⋮xnx12x22⋮xn2
θ1θ2θ3
=A⋅Θ
A
A
A 为测量自变量组成的矩阵,由回归对象的数学模型式(1)确定,可以记为将
n
n
n个式(1)按行依次堆满矩阵得到;
Θ
\Theta
Θ 是待求系数向量。现在目标变成了找到一组
Θ
\Theta
Θ,使
∣
∣
Z
−
Y
∣
∣
2
2
||Z-Y||_{2}^{2}
∣∣Z−Y∣∣22 (向量元素的平方和)最小 ,即:
m
i
n
J
=
(
Z
−
A
Θ
)
T
(
Z
−
A
Θ
)
\begin{align} \begin{matrix} min & J=(Z-A\Theta)^T(Z-A\Theta) \end{matrix} \end{align}
minJ=(Z−AΘ)T(Z−AΘ)
可以使用之前介绍的牛顿拉夫逊法来求解,也可以用其它梯度下降法,然而这种形式的目标函数可以求出解的表达式,使用更方便,式(3)展开为:
J
=
(
Z
T
−
Θ
T
A
T
)
(
Z
−
A
Θ
)
=
Z
T
Z
−
Z
T
A
Θ
−
Θ
T
A
T
Z
+
Θ
T
A
T
A
Θ
=
Z
T
Z
−
2
Z
T
A
Θ
+
Θ
T
A
T
A
Θ
∵
(
Θ
T
A
T
Z
)
T
=
Z
T
(
Θ
T
A
T
)
T
=
Z
T
A
Θ
→
Θ
T
A
T
Z
=
Z
T
A
Θ
\begin{align} J&=(Z^T-\Theta^TA^T)(Z-A\Theta) \nonumber\\ \nonumber\\ &=Z^TZ-Z^TA\Theta-\Theta^TA^TZ+\Theta^TA^TA\Theta \nonumber\\ \nonumber\\ &=Z^TZ-2Z^TA\Theta+\Theta^TA^TA\Theta \\ \nonumber\\ &\because (\Theta^TA^TZ)^T=Z^T(\Theta^TA^T)^T=Z^TA\Theta \to \Theta^TA^TZ=Z^TA\Theta \nonumber \end{align}
J=(ZT−ΘTAT)(Z−AΘ)=ZTZ−ZTAΘ−ΘTATZ+ΘTATAΘ=ZTZ−2ZTAΘ+ΘTATAΘ∵(ΘTATZ)T=ZT(ΘTAT)T=ZTAΘ→ΘTATZ=ZTAΘ
式(4)对
Θ
\Theta
Θ 求偏导并令其为零:
∂
J
∂
Θ
=
0
−
2
(
Z
T
A
)
T
+
(
Θ
T
A
T
A
)
ˊ
Θ
+
Θ
T
A
T
A
(
Θ
)
ˊ
=
−
2
(
Z
T
A
)
T
+
A
T
A
Θ
+
Θ
T
A
T
A
=
−
2
A
T
Z
+
2
A
T
A
Θ
=
0
Θ
=
(
A
T
A
)
−
1
A
T
Z
\begin{align} \frac{\partial J}{\partial \Theta}&=0-2(Z^TA)^T+(\Theta^TA^TA\acute{)} \Theta +\Theta^TA^TA(\Theta \acute{)} \nonumber\\ \nonumber\\ &=-2(Z^T A)^T+A^T A\Theta+\Theta^T A^T A \nonumber\\ \nonumber\\ &=-2A^T Z+2A^T A\Theta \nonumber\\ \nonumber\\ &=0 \\ \nonumber\\ \Theta&=(A^TA)^{-1}A^TZ \end{align}
∂Θ∂JΘ=0−2(ZTA)T+(ΘTATA)ˊΘ+ΘTATA(Θ)ˊ=−2(ZTA)T+ATAΘ+ΘTATA=−2ATZ+2ATAΘ=0=(ATA)−1ATZ
式(6)就是经常见到的最小二乘法求解模型系数公式,只要获取足量的数据,矩阵
A
A
A与系统自变量有关,向量
Z
Z
Z与系统的应变量有关,直接一步就可以得到模型的参数了。
3. 实例测试
测试最小二乘法对函数 z = x 2 − 2 x + 1 z=x^2-2x+1 z=x2−2x+1 的拟合能力。
![]() 41组数据 |
![]() 401组数据 |
分别取41和401组数据时 [ a , b , c ] [a,b,c] [a,b,c] 分别为 [ 0.9817 , − 1.8774 , 1.0976 ] [0.9817,-1.8774,1.0976] [0.9817,−1.8774,1.0976] 和 [ 1.0128 , − 1.9979 , 1.0160 ] [1.0128,-1.9979,1.0160] [1.0128,−1.9979,1.0160]。
4. 总结
推导了最小二乘法的求解公式,在线性回归中,重要的是回归的数学模型如何建立。值得一提的一点,式(1)抛物线含有平方项,属于非线性多项式模型,但是在求解未知系数时,其实是在求解线性方程组(5),所以尽管模型是非线性,过程还是属于线性回归。若模型采用多层神经网络,则式(5)是非线性方程组,宜采用梯度下降法直接求式(3)极值,这时候就叫深度学习了!
5. 详细代码
TestLSM.m
%
% Auther: sanfan66
%
clc
clear all
close all
X=[-2:0.01:2]';
Z0=TestFun(X);
n=length(X);%数据组数
rng(0,'v5uniform')%随机数种子
Z=Z0+randn(size(Z0));
nTheta=3;%3个待求系数
A=ones(n,nTheta);
A(:,2)=X;
A(:,3)=X.^2;
theta=inv(A'*A)*A'*Z;
Y=A*theta;
plot(X,Z0,'b',X,Z,'O',X,Y,'r')
legend('理论','测量','最小二乘')
xlim([-2 2]);
ylim([-4 12]);
function Z=TestFun(X)
Z=X.*X-2*X+1;
end
邱锡鹏. 神经网络与深度学习. 2020:机械工业出版社. ↩︎