多元线性回归解析解

最新推荐文章于 2024-05-31 14:04:00 发布

ma_no_lo

最新推荐文章于 2024-05-31 14:04:00 发布

阅读量367

点赞数 14

文章标签：线性回归概率论算法 python 机器学习

本文链接：https://blog.csdn.net/ma_no_lo/article/details/138015728

版权

数学推导

基本概念

在整个线性回归的学习中，求得一组最优的权重是我们的目的。

权重W是一维向量，记作W=[w1 w2 .... wn].一般的，权重为n行的列向量。

特征X是一维向量，记作X=[x1 x2 x3 .....xn]，一般的，特征为n行的列向量。

在二维空间中，我们有y = ax+b,现在，我们将其拓展到n维空间，则有y = w0 + w1x1 w2x2 +.....+wnxn，其中，w0称为偏置。整体记作y = （w.T）X + ε

中心极限定理

中心极限定理是概率论中随机变量分布渐进于正态分布，而ε服从正态分布。

最大似然估计（MLE）

给定一个概率分布，已知其概率密度函数或者概率质量函数f，以及分布参数θ，我们可以从这个分布中抽出一个具有n个值的采样，利用f计算其似然函数。记作

L(θ|x1,x2....xn) = f(θ)(x1,x2....xn)

即给定分布参数后出现样本组的概率。

其中，最常见的连续概率分布是正态分布

$f(x|\mu ,\sigma ^{2}) = \frac{1}{\sqrt{2\pi \sigma ^{2}}}*e^{-\frac{(x-\mu )^{2}}{2\sigma ^{2}}}$

所以正态分布的线性回归的最大似然为

$L(\theta ) (\varepsilon1...\varepsilon m ) = f(\varepsilon 1....\varepsilon n|\mu ,\sigma ^{2})$

假设误差服从正态分布，符合中心定理，将μ置零，同时σ为常数。事件相互独立，所以可以写成

$\prod_{i=1}^{m}f(\varepsilon i|\mu ,\sigma ^{2})$ = $\prod_{i=1}^{m}\frac{1}{\sqrt{2\pi \sigma ^{2}}}*e^{-\frac{(\varepsilon i )^{2}}{2\sigma ^{2}}}$

又因为 $\varepsilon i = \left | yi - \hat{y} \right |$ = $\left | yi - W.Txi \right | = \left | yi - \theta .Txi \right |$

带入得到

$L(\theta ) = \prod_{i=1}^{m}\frac{1}{\sqrt{2\pi \sigma ^{2}}}*e^{-\frac{(yi-\theta .Txi )^{2}}{2\sigma ^{2}}}$

我们要求的便是当L(θ)最大时所对应的θ*

$argmaxL(\theta ) =argmax \prod_{i=1}^{m}\frac{1}{\sqrt{2\pi \sigma ^{2}}}*e^{-\frac{(yi-\theta .Txi )^{2}}{2\sigma ^{2}}}$

对数似然函数

$l(\theta ) = lnL(\theta ) = mln\frac{1}{\sqrt{2\pi }\sigma } - \frac{1}{2\sigma ^{2}}\sum_{i=1}^{m}(yi - \theta .Txi)^{2}$

$J(\theta ) = \frac{1}{2}\sum_{i=1}^{m}(h\theta (xi)-yi)^{2}$

J(θ)称为损失函数，当其最小时，对应的θ为所求

解析解求法

因为 $J(\theta ) = \frac{1}{2}\sum_{i=1}^{m}(h\theta (xi)-yi)^{2}$ 等价于一个长度为m的向量乘以他自己，但是向量相乘在矩阵中要为行向量乘列向量，所以要转置其中一个。

$J(\theta ) = \frac{1}{2}(X\theta -y).T(X\theta -y)$

根据转置的性质

得到 $J(\theta ) = \frac{1}{2}(\theta .TX.TX\theta - \theta .TX.Ty - y.TX\theta + y.Ty )$

梯度

方向导数

设二元函数上有一点P(x0,y0)，那么以P点为起的射线

$x = x0 + t\cos \alpha$

$y = y0 + t\cos \beta$

$\alpha +\beta = \frac{\pi }{2}$

$\lim_{x\rightarrow 0}\frac{f(x0+t\cos \alpha,yo+t\cos \beta )-f(x0,y0)}{t} = \frac{\partial f}{\partial l}|(x0,y0)$

f(x,y)在(x0,y0)可微，那么方向导数存在， $\frac{\partial f}{\partial l}|(x0,y0) = \frac{\partial f}{\partial x}(x0,y0)\cos \alpha +\frac{\partial f}{\partial y}(x0,y0)\cos \beta$

那么

$\frac{\partial f}{\partial l}|(x0,y0) = (\frac{\partial f}{\partial x}(x0,y0),\frac{\partial f}{\partial y}(x0,y0))\cdot (cos\alpha ,\cos \beta )$

记 $\bigtriangledown f(x0,y0) = (\frac{\partial f}{\partial x}(x0,y0),\frac{\partial f}{\partial y}(x0,y0))$ ，其中， $\bigtriangledown f(x0,y0)$ 记作梯度。

即梯度是某点下降最快的方向导数。

最优解

最优解一定是一个驻点，但驻点不一定是最优解，这时我们用Hessian矩阵的半正定来判断。

Hessian矩阵

Hessian矩阵是由目标函数在X处的二阶偏导数所构成的对称矩阵。

半正定：特征值全部大于等于0

正定：特征值全部大于0

$J{}''(\theta ) = X.TX$

这相当于把X中的元素平方相加，结果肯定大于等于0。所以目标函数的黑塞矩阵半正定，所以函数是凸函数，驻点即为最优解。

目标函数求导

矩阵求导的性质

$\frac{\partial \theta .TA\theta }{\partial \theta } = 2A\theta$

$\frac{\partial \theta .TA}{\partial \theta } = A$

$\frac{\partial A\theta }{\partial \theta } = A.T$

所以 $J{}'(\theta ) = X.TX\theta - X.Ty$

令其等于0，得到

$\theta = (X.TX)^{-1}\cdot X.Ty$

至此，解析解求解完成。

代码实现

import numpy as np
import matplotlib.pyplot as plt


X1 = 2 * np.random.rand(100,1)#创建一百行一列的列向量
X2 = 3 * np.random.rand(100,1)
#这里的y是真实值，等于y_hat+error
y = 5 + 4*X1 + 3*X2 + np.random.randn(100,1)#error服从正态分布

X_b = np.c_[np.ones((100,1)),X1,X2]
print(X_b)
theta = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)
print(theta)

#使用模型进行预测
X_new = np.array([[0,0],
                 [2,3]])
X_new_b = np.c_[np.ones((2,1)),X_new]
print(X_new)
y_predict  = X_new_b.dot(theta)
print(y_predict)

可视化

fig = plt.figure()
ax3d = fig.add_subplot(projection='3d')
ax3d.plot(X1,X2,y,'b.')
x1 = X_new[:,0]
x2 = X_new[:,1]
Z = y_predict[:,0]
print(x1)
print(x2)
print(Z)
ax3d.plot(x1,x2,Z,color='black')
plt.show()