Skip to content

多元线性回归

1、基本概念

  线性回归是机器学习中有监督机器学习下的一种算法。 回归问题主要关注的是因变量(需要预测的值,可以是一个也可以是多个)和一个或多个数值型的自变量(预测变量)之间的关系。

  需要预测的值:即目标变量,target,y,连续值预测变量。

  影响目标变量的因素:X1...Xn,可以是连续值也可以是离散值。

  因变量和自变量之间的关系:即模型,model,是我们要求解的。

1.1、连续值

1.2、离散值

1.3、简单线性回归

  前面提到过,算法说白了就是公式,简单线性回归属于一个算法,它所对应的公式。

  y=wx+b

  这个公式中,y 是目标变量即未来要预测的值,x 是影响 y 的因素,w,b 是公式上的参数即要求的模型。其实 b 就是咱们的截距,w 就是斜率嘛! 所以很明显如果模型求出来了,未来影响 y 值的未知数就是一个 x 值,也可以说影响 y 值 的因素只有一个,所以这是就叫简单线性回归的原因。

  同时可以发现从 x 到 y 的计算,x 只是一次方,所以这是算法叫线性回归的原因。 其实,大家上小学时就已经会解这种一元一次方程了。为什么那个时候不叫人工智能算法呢?因为人工智能算法要求的是最优解!

1.4、最优解

  Actual value:真实值,一般使用 y 表示。

  Predicted value:预测值,是把已知的 x 带入到公式里面和出来的参数 w,b 计算得到的,一般使用 y^ 表示。

   Error:误差,预测值和真实值的差距,一般使用 ε 表示。

  最优解:尽可能的找到一个模型使得整体的误差最小,整体的误差通常叫做损失 Loss。

  Loss:整体的误差,Loss 通过损失函数 Loss function 计算得到。

1.5、多元线性回归

  现实生活中,往往影响结果 y 的因素不止一个,这时 x 就从一个变成了 n 个,x1...xn 同时简单线性回归的公式也就不在适用了。多元线性回归公式如下:

  y^=w1x1+w2x2++wnxn+b

  b是截距,也可以使用w0来表示

  y^=w1x1+w2x2++wnxn+w0

  y^=w1x1+w2x2++wnxn+w01

  使用向量来表示,X表示所有的变量,是一维向量;W表示所有的系数(包含w0),是一维向量,根据向量乘法规律,可以这么写:

  y^=WTX

2、正规方程

2.1、最小二乘法

最小二乘法(Least Squares Method)是一种常见的数学优化技术,用于拟合一组数据点的最佳线性回归方程。它的目标是找到一条直线,使得该直线与数据点之间的误差平方和最小。

具体而言,最小二乘法的目标是最小化每个数据点与拟合直线之间的垂直距离的平方和,即最小化所有数据点到拟合直线的残差平方和。残差是指每个数据点的实际值与拟合直线在该点处的预测值之间的差异

最小二乘法可用于回归分析,即对一组数据进行拟合,并预测未来数据点的值。此外,它还可以用于求解线性方程组,如线性回归模型的系数或参数估计等。

  最小二乘法可以将误差方程转化为有确定解的代数方程组(其方程式数目正好等于未知数的个数),从而可求解出这些未知参数。这个有确定解的代数方程组称为最小二乘法估计的正规方程。公式如下:

θ=(XTX)1XTy 或者 W=(XTX)1XTy ,其中的Wθ 即使方程的解!

公式是如何推导的?

最小二乘法公式如下:

J(θ)=12i=0n(hθ(xi)yi)2

使用矩阵表示:

之所以要使用转置T,是因为,矩阵运算规律是:矩阵A的一行乘以矩阵B的一列!

2.2、多元一次方程举例

1、二元一次方程

2、三元一次方程

3、八元一次方程

Python
# 上面八元一次方程对应的X数据
X = np.array([[  0 ,14 , 8 ,  0 ,  5,  -2,   9,  -3],
 [ -4 , 10 ,  6 ,  4 ,-14 , -2 ,-14  , 8],
 [ -1 , -6  , 5 ,-12 ,  3 , -3 ,  2 , -2],
 [  5 , -2  , 3 , 10  , 5 , 11 ,  4  ,-8],
 [-15 ,-15  ,-8 ,-15 ,  7 , -4, -12 ,  2],
 [ 11 ,-10 , -2 ,  4  , 3 , -9 , -6 ,  7],
 [-14 ,  0 ,  4 , -3  , 5 , 10 , 13 ,  7],
 [ -3 , -7 , -2 , -8  , 0 , -6 , -5 , -9]])
# 对应的y
y = np.array([ 339 ,-114  , 30 , 126, -395 , -87 , 422, -309])
display(X,y)

2.3、矩阵转置公式与求导公式

转置公式如下:

  • (mA)T=mAT,其中m是常数
  • (A+B)T=AT+BT
  • (AB)T=BTAT
  • (AT)T=A

假设我们有一个矩阵 A,它是一个 3×2 的矩阵,其中包含了 6 个元素。我们可以将其表示为:

A=(123456)

现在,如果我们想要转置矩阵 A,我们可以将 A 中的行和列互换。这样就得到了一个新的矩阵 AT,它的行数和 A 的列数相等,列数和 A 的行数相等。在这种情况下,AT 是一个 2×3 的矩阵。我们可以使用下面的公式来计算 AT

Ai,jT=Aj,i

也就是说,新矩阵中的第 i 行和第 j 列的元素等于原矩阵中的第 j 行和第 i 列的元素。

下面是具体的演示:

AT=(135246)

在上面的例子中,我们首先写出了矩阵 A,然后使用公式 Ai,jT=Aj,i 将其转置为矩阵 AT。转置矩阵的行和列与原矩阵相反,其中原矩阵的第 i 行变为了新矩阵的第 i 列,原矩阵的第 j 列变为了新矩阵的第 j 行。

求导公式如下:

  • XTX=I$$
  • AXTX=A
  • AXX=AT
  • XAX=AT
  • XTAXX=(A+AT)X; A不是对称矩阵
  • XTAXX=2AX; A是对称矩阵

2.4、推导正规方程 θ 的解

  1. 矩阵乘法公式展开
  1. 进行求导(注意X、y是已知量,θ 是未知数):
  • yTy 是常量

根据2.3、矩阵转置公式与求导公式可知

  1. 根据上面求导公式进行运算:
  • 矩阵运算分配律
  1. 令导数J(θ)=0
  1. 矩阵没有除法,使用逆矩阵进行转化:

到此为止,公式推导出来了~

2.5、凸函数判定

判定损失函数是凸函数的好处在于我们可能很肯定的知道我们求得的极值即最优解,一定是全局最优解。

如果是非凸函数,那就不一定可以获取全局最优解~

来一个更加立体的效果图:

判定凸函数的方式: 判定凸函数的方式非常多,其中一个方法是看黑塞矩阵是否是半正定的。

黑塞矩阵(hessian matrix)是由目标函数在点 X 处的二阶偏导数组成的对称矩阵。

对于我们的式子来说就是在导函数的基础上再次对θ来求偏导,结果就是 XTX。所谓正定就是 XTX 的特征值全为正数,半正定就是 XTX 的特征值大于等于 0, 就是半正定。

这里我们对 J(θ) 损失函数求二阶导数的黑塞矩阵是 XTX ,得到的一定是半正定的,自己和自己做点乘嘛!

X 为任意大小的实矩阵时,XTX 一定是半正定的。

要证明 XTX 是半正定的,可以考虑对于任意非零向量 v,证明 vT(XTX)v 非负。注意到 vT(XTX)v=(Xv)T(Xv)=||Xv||20,其中 ||Xv|| 表示向量 Xv 的范数,等于 (Xv)T(Xv)。因此,XTX 是半正定的。

这里不用数学推导证明这一点。在机器学习中往往损失函数都是凸函数,到深度学习中损失函数往往是非凸函数,即找到的解未必是全局最优,只要模型堪用就好!机器学习特点是:不强调模型 100% 正确,只要是有价值的,堪用的,就Okay!

3、线性回归实战

3.1、使用正规方程进行求解

3.1.1、简单线性回归

一元一次方程,在机器学习中一元表示一个特征,b表示截距,y表示目标值。

Python
import numpy as np
import matplotlib.pyplot as plt
# 转化成矩阵
X = np.linspace(0,10,num = 30).reshape(-1,1)
# 斜率和截距,随机生成
w = np.random.randint(1,5,size = 1)
b = np.random.randint(1,10,size = 1)
# 根据一元一次方程计算目标值y,并加上“噪声”,数据有上下波动~
y = X * w + b + np.random.randn(30,1)
plt.scatter(X,y)
# 重新构造X,b截距,相当于系数w0,前面统一乘以1
X = np.concatenate([X,np.full(shape = (30,1),fill_value= 1)],axis = 1)
# 正规方程求解
θ = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(y).round(2)
print('一元一次方程真实的斜率和截距是:',w, b)
print('通过正规方程求解的斜率和截距是:',θ)
# 根据求解的斜率和截距绘制线性回归线型图
plt.plot(X[:,0],X.dot(θ),color = 'green')

效果如下(random.randn是随机生成正太分布数据,所以每次执行图形会有所不同):

3.1.2、多元线性回归

二元一次方程,x1x2 相当于两个特征,b是方程截距

Python
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d.axes3d import Axes3D # 绘制三维图像
# 转化成矩阵
x1 = np.random.randint(-150,150,size = (300,1))
x2 = np.random.randint(0,300,size = (300,1))
# 斜率和截距,随机生成
w = np.random.randint(1,5,size = 2)
b = np.random.randint(1,10,size = 1)
# 根据二元一次方程计算目标值y,并加上“噪声”,数据有上下波动~
y = x1 * w[0] + x2 * w[1] + b + np.random.randn(300,1)
fig = plt.figure(figsize=(9,6))
ax = Axes3D(fig)
ax.scatter(x1,x2,y) # 三维散点图
ax.view_init(elev=10, azim=-20) # 调整视角
# 重新构造X,将x1、x2以及截距b,相当于系数w0,前面统一乘以1进行数据合并
X = np.concatenate([x1,x2,np.full(shape = (300,1),fill_value=1)],axis = 1)
w = np.concatenate([w,b])
# 正规方程求解
θ = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(y).round(2)
print('二元一次方程真实的斜率和截距是:',w)
print('通过正规方程求解的斜率和截距是:',θ.reshape(-1))
# # 根据求解的斜率和截距绘制线性回归线型图
x = np.linspace(-150,150,100)
y = np.linspace(0,300,100)
z = x * θ[0] + y * θ[1] + θ[2]
ax.plot(x,y,z ,color = 'red')

效果如下:

3.2、机器学习库scikit-learn

3.2.1、scikit-learn简介

3.2.2、scikit-learn实现简单线性回归

L2范数是指向量中每个元素的平方和的平方根。在数学中,L2范数也称为欧几里得范数,通常用 ||x||2 表示,其中x是一个向量。L2范数可用于衡量向量的大小或长度。

L2范数在机器学习和数据科学中经常用于正则化和优化问题中,例如岭回归和支持向量机等算法。在这些情况下,通过对损失函数添加L2正则化项,可以避免过拟合并提高模型的泛化能力。同时,L2范数也被用来作为相似性度量的度量标准,如余弦相似性。

||x||2=(|xi|2)

其中,||x||2 表示 x 的 L2 范数,x_i 表示 x 中的第 i 个元素,∑ 表示对所有元素求和,|x_i| 表示 x_i 的绝对值,sqrt 表示平方根。

Python
from sklearn.linear_model import LinearRegression
import numpy as np
import matplotlib.pyplot as plt
# 转化成矩阵
X = np.linspace(0,10,num = 30).reshape(-1,1)
# 斜率和截距,随机生成
w = np.random.randint(1,5,size = 1)
b = np.random.randint(1,10,size = 1)
# 根据一元一次方程计算目标值y,并加上“噪声”,数据有上下波动~
y = X * w + b + np.random.randn(30,1)
plt.scatter(X,y)
# 使用scikit-learn中的线性回归求解
model = LinearRegression()
model.fit(X,y)
w_ = model.coef_
b_ = model.intercept_
print('一元一次方程真实的斜率和截距是:',w, b)
print('通过scikit-learn求解的斜率和截距是:',w_,b_)
plt.plot(X,X.dot(w_) + b_,color = 'green')

3.2.3、scikit-learn实现多元线性回归
Python
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d.axes3d import Axes3D
# 转化成矩阵
x1 = np.random.randint(-150,150,size = (300,1))
x2 = np.random.randint(0,300,size = (300,1))
# 斜率和截距,随机生成
w = np.random.randint(1,5,size = 2)
b = np.random.randint(1,10,size = 1)
# 根据二元一次方程计算目标值y,并加上“噪声”,数据有上下波动~
y = x1 * w[0] + x2 * w[1] + b + np.random.randn(300,1)
fig = plt.figure(figsize=(9,6))
ax = Axes3D(fig)
ax.scatter(x1,x2,y) # 三维散点图
ax.view_init(elev=10, azim=-20) # 调整视角
# 重新构造X,将x1、x2以及截距b,相当于系数w0,前面统一乘以1进行数据合并
X = np.concatenate([x1,x2],axis = 1)
# 使用scikit-learn中的线性回归求解
model = LinearRegression()
model.fit(X,y)
w_ = model.coef_.reshape(-1)
b_ = model.intercept_
print('二元一次方程真实的斜率和截距是:',w,b)
print('通过scikit-learn求解的斜率和截距是:',w_,b_)
# # 根据求解的斜率和截距绘制线性回归线型图
x = np.linspace(-150,150,100)
y = np.linspace(0,300,100)
z = x * w_[0] + y * w_[1] + b_
ax.plot(x,y,z ,color = 'green')

4、线性回归算法推导

高斯与最小二乘的故事

4.1、深入理解回归

  回归简单来说就是“回归平均值”(regression to the mean)。但是这里的 mean 并不是把历史数据直接当成未来的预测值,而是会把期望值当作预测值。 追根溯源回归这个词是一个叫高尔顿的人发明的,他通过大量观察数据发现:父亲比较高,儿子也比较高;父亲比较矮,那么儿子也比较矮!正所谓“龙生龙凤生凤老鼠的儿子会打洞”!但是会存在一定偏差~

  父亲是 1.98,儿子肯定很高,但有可能不会达到1.98   父亲是 1.69,儿子肯定不高,但是有可能比 1.69 高

  大自然让我们回归到一定的区间之内,这就是大自然神奇的力量。

  高尔顿是谁?达尔文的表弟,这下可以相信他说的十有八九是对的了吧!

  人类社会很多事情都被大自然这种神奇的力量只配置:身高、体重、智商、相貌……

  这种神秘的力量就叫正态分布。大数学家高斯,深入研究了正态分布,最终推导出了线性回归的原理:最小二乘法

  接下来,我们跟着高斯的足迹继续向下走~

4.2、误差分析

  误差 εi 等于第 i 个样本实际的值 yi 减去预测的值 y^ ,公式可以表达为如下:

  εi=|yiy^|

  εi=|yiWTxi|

  假定所有的样本的误差都是独立的,有上下的震荡,震荡认为是随机变量,足够多的随机变量叠加之后形成的分布,它服从的就是正态分布,因为它是正常状态下的分布,也就是高斯分布均值是某一个值,方差是某一个值。 方差我们先不管,均值我们总有办法让它去等于零 0 的,因为我们这里是有截距b, 所有误差我们就可以认为是独立分布的,1<=i<=n,服从均值为 0,方差为某定值的高斯分布。线性回归算法中我们假设误差符合均值为0,方差为定值的正态分布!!!

4.3、最大似然估计

  最大似然估计(maximum likelihood estimation, MLE)一种重要而普遍的求估计量的方法。最大似然估计明确地使用概率模型,其目标是寻找能够以较高概率产生观察数据的系统发生树。最大似然估计是一类完全基于统计的系统发生树重建方法的代表。

  是不是,有点看不懂,太学术了,我们举例说明~

  假如有一个盒子,里面有黑白两种颜色的球,数目多少不知,两种颜色的比例也不知。我们想知道盒子中白球和黑球的比例,但我们不能把盒子中的球全部拿出来数。现在我们可以每次任意从已经摇匀的盒子中拿一个球出来,记录球的颜色,然后把拿出来的球再放回盒子中。这个过程可以重复,我们可以用记录的球的颜色来估计盒子中黑白球的比例。假如在前面的一百次重复记录中,有七十次是白球,请问盒子中白球所占的比例最有可能是多少?

请告诉我你的答案!

很多小伙伴,甚至不用算,凭感觉,就能给出答案:70%

下面是详细推导过程:

  • 最大似然估计,计算

  • 白球概率是p,黑球是1-p(罐子中非黑即白)

  • 取一个请问是白球的概率是多少?

    • p
  • 取两个球,两个球都是白色,概率是多少?

    • p2
  • 取5个球都是白色,概率是多少?

    • p5
  • 取10个球,9个是白色,一个是黑色,概率是多少呢?

    取出 10 个球,其中 9 个是白色,1 个是黑色的概率表示如下:

    P(x=k)=C109p9(1p)

    其中,C109表示从 10 个元素中选择 9 个元素的组合数目。

    所以,取出 10 个球,其中 9 个是白色,1 个是黑色的概率为:

    • C101=C109 这个两个排列组合公式是相等的~
    • C109p9(1p)=C101p9(1p)
  • 罐子取100个球,70次是白球,30次是黑球,概率是多少?

  • P=C10030p70(1p)30
  • 最大似然估计,什么时候P最大呢,最大才符合真实观测,符合实际情况?

    C10030是常量,可以去掉

    p > 0,1- p > 0,所以上面概率想要求最大值,那么求导数即可!

  • P=70p69(1p)30+p7030(1p)29(1)

    令导数为0:

  • 0=70p69(1p)30+p7030(1p)29(1)

    公式化简:

  • 0=70(1p)p30
  • 0=70100p
  • p = 70%

4.4、高斯分布-概率密度函数

最常见的连续概率分布是正态分布,也叫高斯分布,而这正是我们所需要的,其概率密度函数如下:

公式如下:

随着参数μ和σ变化,概率分布也产生变化。 下面重要的步骤来了,我们要把一组数据误差出现的总似然,也就是一组数据之所以对应误差出现的整体可能性表达出来了,因为数据的误差我们假设服从一个高斯分布,并且通过截距项来平移整体分布的位置从而使得 μ=0,所以样本的误差我们可以表达其概率密度函数的值如下:

简化如下:

4.5、误差总似然

和前面黑球白球问题类似,也是一个累乘问题~

  根据前面公式εi=|yiWTxi|可以推导出来如下公式:

公式中的未知变量就是WT,即方程的系数,系数包含截距~如果,把上面当成一个方程,就是概率P关于W的方程!其余符号,都是常量!

现在问题,就变换成了,求最大似然问题了!不过,等等~

累乘的最大似然,求解是非常麻烦的!

接下来,我们通过,求对数累乘问题,转变为累加问题(加法问题,无论多复杂,都难不倒我了!)

4.6、最小二乘法MSE

PW=i=0n12πσe(yiWTxi)22σ2

根据对数,单调性,对上面公式求自然底数e的对数,效果不变~

loge(PW)=loge(i=0n12πσe(yiWTxi)22σ2)

接下来 log 函数继续为你带来惊喜,数学上连乘是个大麻烦,即使交给计算机去求解它也得哭出声来。惊喜是:

  • loga(XY)=logaX+logaY
  • logaXY=logaXlogaY
  • logaXn=nlogaX
  • loga(X1X2Xn)=logaX1+logaX2++logaXn
  • logxxn=n(nR)
  • loga1X=logaX
  • logaNyx=yxlogaN

loge(PW)=loge(i=0n12πσe(yiWTxi)22σ2)

      累乘问题变成累加问题~

乘风破浪,继续推导--->

      

      

上面公式是最大似然求对数后的变形,其中πσ都是常量,而(yiWTxi)2肯定大于!上面求最大值问题,即可转变为如下求最小值问题:

L(W)=12i=0n(y(i)WTx(i))2 L代表Loss,表示损失函数,损失函数越小,那么上面最大似然就越大~

有的书本上公式,也可以这样写,用J(θ)表示一个意思,θ角色就是W:

J(θ)=12i=1n(y(i)θTx(i))2=12i=1n(θTx(i)y(i))2

进一步提取:

其中:

  y^=hθ(X)=Xθ 表示全部数据,是矩阵,X表示多个数据,进行矩阵乘法时,放在前面

  y^i=hθ(x(i))=θTx(i) 表示第i个数据,是向量,所以进行乘法时,其中一方需要转置

因为最大似然公式中有个负号,所以最大总似然变成了最小化负号后面的部分。 到这里,我们就已经推导出来了 MSE 损失函数J(θ),从公式我们也可以看出来 MSE 名字的来 历,mean squared error,上式也叫做最小二乘法!

4.7、归纳总结升华

  这种最小二乘法估计,其实我们就可以认为,假定了误差服从正太分布,认为样本误差的出现是随机的,独立的,使用最大似然估计思想,利用损失函数最小化 MSE 就能求出最优解!所以反过来说,如果我们的数据误差不是互相独立的,或者不是随机出现的,那么就不适合去假设为正太分布,就不能去用正太分布的概率密度函数带入到总似然的函数中,故而就不能用 MSE 作为损失函数去求解最优解了!所以,最小二乘法不是万能的~

Released under the MIT License.