下面的有些叙述基于我个人理解, 可能与专业书籍描述不同, 但是最终都是表达同一个意思, 如果有不同意见的小伙伴, 请在评论区留言, 我不胜感激.
参考:
-
周志华-机器学习
-
最小二乘法求解多元回归方程: https://blog.csdn.net/weixin_39445556/article/details/83543945
-
梯度下降算法: https://blog.csdn.net/h_l_dou/article/details/82826653
-
回归方程度量方式: https://blog.csdn.net/qq_15698613/article/details/86478736
前提介绍
常见统计量
均值:
中位数:
将样本数据按数值大小顺序排列, 取中间的变量.
众数:
数据中出现次数最多的数.
方差:
描述数据离散程度.
标准差:
简单线性回归
问题背景
举一个例子说明: 商家卖出产品的数量与做广告数量的线性关系(如下图), 当做的广告越多, 那么卖出的产品就越多; 在这个问题中, 卖出产品这个样本中, 就存在一个属性: 做广告数量, 因此可以用线性方程对这个问题进行模拟, 利用线性方程大概预测出广告与产品之间的数量关系, 但是如何得到最优的线性方程, 下面通过简单线性回归模型求解.
产品数量与广告数量的关系:
线性回归方程:
建立简单线性回归方程
-
简单线性回归模型:
β0: 回归线截距
β1: 回归线斜率
ε: 偏差 -
简单线性回归方程:
E(y) = β0 + β1*x
该方程用于描述在给定x情况下y的期望值. -
估计的简单线性回归方程:
ŷ = b0 + b1*x
描述给定x情况下求ŷ 即y的估计值 -
偏差ε假定:
ε是随机变量, 均值为0
ε的方差对于所有的自变量x都是一样的
ε独立, 满足正态分布 -
求解最优简单线性回归模型的前提:
-
b0, b1求解公式:
-
使用前面的案例求解对应的b0, b1:
分子 = (1-2)(14-20)+(3-2)(24-20)+(2-2)(18-20)+(1-2)(17-20)+(3-2)(27-20)= 6 + 4 + 0 + 3 + 7= 20
分母 = （1-2）^2 + (3-2)^2 + (2-2)^2 + (1-2)^2 + (3-2)^2= 1 + 1 + 0 + 1 + 14
b1 = 20/4 =5
b0 = 20 - 5*2 = 20 - 10 = 10
最终得线性回归方程为: y = 5x +10
Python代码实现简单线性回归
import numpy as np#求b0, b1
def fitSLR(x, y):n = len(x)dinominator = 0numerator = 0for i in range(0, n):numerator += (x[i] - np.mean(x))*(y[i] - np.mean(y))dinominator += (x[i] - np.mean(x))**2b1 = numerator/float(dinominator)b0 = np.mean(y)/float(np.mean(x))return b0, b1def predict(x, b0, b1):return b0 + x*b1x = [1, 3, 2, 1, 3]
y = [14, 24, 18, 17, 27] b0, b1 = fitSLR(x, y)print "intercept:", b0, " slope:", b1x_test = 6y_test = predict(6, b0, b1)print "y_test:", y_test
多元线性回归
问题背景
-
简单线性回归模型只能描述具有一个属性的样本, 当样本具有多个属性的时候就不能求解出线性回归方程; 当简单线性回归模型中, 一个属性是分类型变量, 例如: 汽车运送货物, 汽车具有不同的车型, 车型决定了运送时间, 面对这样的属性的时候, 就不能使用分类型属性作为自变量, 而应该将分类型属性进行转变, 将一个分类型属性转变为多个连续型属性, 将简单线性回归模型转化为多元线性回归, 下面将介绍这种转变方法.
-
与简单线性回归的区别就在于, 它存在存在个自变量(即样本存在多个属性). 下面将使用多元线性回归描述.
多元回归方程
-
多元线性回归模型:
y = β0 + β1x1 + β2x2 + … + βp*xp + ε
其中：β0，β１，β2… βp是参数
ε是误差值 -
多元回归方程:
E(y) = β0 + β1x1 + β2x2 + … + βp*xp -
估计多元线性回归方程:
y_hat = b0＋b１x1 + b2x2 + … + bp*xp
处理分类型数据
正如问题背景中描述, 存在如下案例:
英里数 次数 车型 时间
100 4 1 9.3
50 3 0 4.8
100 4 1 8.9
100 2 2 6.5
50 2 2 4.2
80 2 1 6.2
75 3 1 7.4
65 4 0 6
90 3 0 7.6
车型, 运送次数两个属性决定了运送时间, 并且车型属于分类型属性.
针对上面的问题, 有如下解决方案:
将车型分解转换, 车型包含3个子集{0, 1, 2}, 因此将车型转化为3个属性:(0, 1, 2), 每个属性分别使用0/1表示存在/不存在, 如下:
英里数 次数 0 1 2 时间
100 4 0 1 0 9.3
50 3 1 0 0 4.8
100 4 0 1 0 8.9
100 2 0 0 1 6.5
50 2 0 0 1 4.2
80 2 0 1 0 6.2
75 3 0 1 0 7.4
65 4 1 0 0 6
90 3 1 0 0 7.6
原来存在3个属性{英里数, 运输次数, 车型}, 现在变为{英里数, 运输次数, 0, 1, 2}; 最终每个样本可在属性空间中标记对应位置. 因此最终问题还是求解多元回归方程
处理二值数据
在处理分类型数据中得到的0/1只能表示样本该属性存在或不存在, 并不能在属性空间中呈现一个连续的状态, 如下图:
我们很难找到一个评判标准, 必须根据实际数据来决定判断标准, 并且判断标准很容易受数据影响, 比如: 对于图中的模型使用 h(x)>0.2 表达Yes,
因此我们需要使用Sigmoid函数使曲线平滑化, 将含有离散型数据的函数变为连续性数据的函数.(这一切的一切都是为了简化求解多元回归方程!!!)
Sigmoid方程:
z∈(-∞, +∞), y∈(0, 1)
Sigmoid函数图像:
求解多元回归方程
测试数据为X（x0，x1，x2···xn）,分别代表样本属性
要学习的参数为： Θ（θ0，θ1，θ2，···θn）, 代表每个属性前系数, 也可看作是斜率.
-
线性方程如下:
-
向量表示方式:
-
使用Sigmoid函数处理后的预测函数: (将Z方程带入Sigmoid函数中)
-
使用概率表示:(即样本为正反例的概率)
-
正例(y = 1):
-
反例(y = 0):
-
-
求解最优多元线性回归模型的前提:
与简单线性回归方程一样, 必须保证该模型值最小
i: 样本序号
m: 样本总数
hθ(x): θ属性下, x的预测值
上述模型表达每个样本估计值与真实值的差距, 因此我们需要找到合适的 Θ（θ0，θ1，θ2，···θn）使得该模型最后的值达到最小. -
!!! Cost函数(求解多元线性回归方程的问题核心):
-
描述某一样本不同二值数据下的hθ(x), 即 y = 0 时的 hθ(x) 和 y = 1 时的 hθ(x); 由于前面方程带入Sigmoid方程, 所以h(x)取对数
-
由于Cost函数是分段函数, 因此对Cost求均值(如下操作), 得到J(θ), 这步操作的目的是将Cost分段函数合二为一, 使用一个函数描述所有情况下的样本.
-
-
利用均方误差代价函数可得J(θ)与hθ(x)之间的关系, 如下图, 最终问题就转化为求解J(θ)方程的最小值.
-
寻找合适的Θ（θ0，θ1，θ2，···θn）使得 J(θ) 存在最小值:
针对该问题有如下方式求解(使用均方误差最小化):- 最小二乘函数求导，让导函数为0时的结果，就是最小二乘的解。求导过程如下:
变为矩阵表达式
变为矩阵函数
对函数求偏导, 最终得到问题的解
代码实现:
- 最小二乘函数求导，让导函数为0时的结果，就是最小二乘的解。求导过程如下:
import numpy as np
import matplotlib.pyplot as plt#设置随机种子
seed = np.random.seed(100)#构造一个100行1列到矩阵。矩阵数值生成用rand，得到到数字是0-1到均匀分布到小数。
X = 2 * np.random.rand(100,1) #最终得到到是0-2均匀分布到小数组成到100行1列到矩阵。这一步构建列X1(训练集数据)
#构建y和x的关系。 np.random.randn(100,1)是构建的符合高斯分布（正态分布）的100行一列的随机数。相当于给每个y增加列一个波动值。
y= 4 + 3 * X + np.random.randn(100,1)#将两个矩阵组合成一个矩阵。得到的X_b是100行2列的矩阵。其中第一列全都是1.
X_b = np.c_[np.ones((100,1)),X]#!!!!!!!!!!!!!!!!!问题的核心, X, Y矩阵的转换
theta_best = np.linalg.inv(X_b.T.dot(X_b)).dot(X_b.T).dot(y)
# print(theta_best)
# 生成两个新的数据点,得到的是两个x1的值
X_new = np.array([[0],[2]])# 填充x0的值，两个1
X_new_b = np.c_[(np.ones((2,1))),X_new]print (X_new_b)# 用求得的theata和构建的预测点X_new_b相乘，得到yhat
y_predice = X_new_b.dot(theta_best)
print(y_predice)
# 画出预测函数的图像，r-表示为用红色的线
plt.plot(X_new,y_predice,'r-')
# 画出已知数据X和掺杂了误差的y，用蓝色的点表示
plt.plot(X,y,'b.')
# 建立坐标轴
plt.axis([0,2,0,15,])plt.show()- 梯度下降求解:
梯度下降法则:
α代表学习率, 即梯度下降过程中, 每次迭代的偏向度, 偏向度越大则迭代得越快, 反之越慢
梯度下降案例介绍:
设有如下单变量函数:
函数求导:
初始化(设置初始状态θ以及学习率α):
根据梯度下降迭代计算,得到如下迭代计算过程:
最终得到J(θ)函数最小值, 然后固定θ, 得到问题的解:
变为矩阵表达式
然后对J(θ)求关于θ的偏导
J’(θ) 代入梯度下降法则中, 得到多元回归方程的更新法则:
- 公式分解
代表该属性下, 预测值与实际值的差, 也就是
-------------------------------------------------------------------------------------------------------------------
等价于如下公式
最终通过上述转化, 就可以利用矩阵转变为h(x)的关系式.
代码实现:
import numpy as np
import random
'''
求解非线性回归方程, 找到合适的theta使得 求和(hypothesis - y) 最小
因此利用梯度下降得到最小值,在得到最小值的过程中不断修正theta, 最终可解
'''#theta: 要求的回归方程前系数,向量, theta1,theta2
#alpha: 学习率
#m: 实例个数
#numIterations: 梯度下降迭代次数
def gradientDescent(x, y, theta, alpha, m, numIterations):#求矩阵装置xTrans = x.transpose()for i in range(0, numIterations):#x, theta内积, 即方程所有的估计值hypothesis = np.dot(x, theta)#估计值减去实际值得出丢失部分, 利用丢失的部分求gradient#依次使数据尽量拟合到方程上, loss的值逐渐向最小值处靠近#数学上则是利用求该点曲线斜率,根据斜率来决定需要下降多少,即求gradientloss = hypothesis - y# avg cost per example (the 2 in 2*m doesn't really matter here.# But to be consistent with the gradient, I include it)#最小二乘法求cost得出最小值 //均方误差代价函数cost = np.sum(loss ** 2) / (2 * m)print("Iteration %d | Cost: %f" % (i, cost))# avg gradient per example#np.dot(xTrans, loss)代表在当前样本loss下的总的下降度#要表示当前样本就需要: 总/mgradient = np.dot(xTrans, loss) / m# update#每下降一次就更新一次总的theta(看似是总的theta,实际更新还是从theta0, theta1...一个一个走,每次loss的变化也就改变了对应的theta,loss使得估计数据逐渐靠近真实数据, 从而theta也跟着改变,当迭代次数足够多的时候, theta就得到最接近真实值)theta = theta - alpha * gradientreturn theta#利用随机数得到获得x, y的数据集
def genData(numPoints, bias, variance):#x:numPoints行,2列; y:numPosits行,1列//标签x = np.zeros(shape=(numPoints, 2))y = np.zeros(shape=numPoints)# basically a straight linefor i in range(0, numPoints):x[i][0] = 1x[i][1] = i#bias:偏差, variance:方差y[i] = (i + bias) + random.uniform(0, 1) * variancereturn x, yx, y = genData(100, 25, 10)
#获取x的行数,列数
m, n = np.shape(x)
numIterations= 100000
alpha = 0.0005
#新建一个矩阵
theta = np.ones(n)
theta = gradientDescent(x, y, theta, alpha, m, numIterations)
print(theta)
sklearn求解多元回归方程
from sklearn.linear_model import LinearRegression
import numpy as np
import matplotlib.pyplot as plt
# 解析解求线性回归# 手动构建数据集和y与x的对应关系
x = 2 * np.random.rand(100,1)
y= 4 + 3*x + np.random.randn(100,1)line_reg = LinearRegression()
# 训练数据集,训练完成后，参数会保存在对象line_reg中。
line_reg.fit(x,y)# line_reg.intercept为截距，就是w0，line_reg.coef_为其他参数，coef的全拼为coefficient
print(line_reg.intercept_,line_reg.coef_)x_new = np.array([[0],[2]])
# line_reg.predict(x_new) 为预测结果
print(line_reg.predict(x_new))plt.plot(x_new,line_reg.predict(x_new),'r-')
# 画出已知数据X和掺杂了误差的y，用蓝色的点表示
plt.plot(x,y,'b.')
# 建立坐标轴
plt.axis([0,2,0,15,])plt.show()回归中的相关系数r和决定系数R^2
相关系数
相关系数: 描述两个值线性相关强度的量
取值范围 [-1, 1]: 正向相关: >0, 负向相关：<0, 无相关性：=0
决定系数
决定系数：R平方值
决定系数（coefficient of determination），有的教材上翻译为判定系数，也称为拟合优度。是相关系数的平方。表示可根据自变量的变异来解释因变量的变异部分。
R平方值为0.8，则表示回归关系可以解释80%的变异。
R2=SSR/SST=1-SSE/SST
其中：SST=SSR+SSE，SST (total sum of squares)为总平方和，SSR (regression sum of squares)为回归平方和，SSE (error sum of squares) 为残差平方和。
R平方的局限性：R平方随着自变量的增加会变大，R平方和样本量是有关系的。因此，我们要到R平方进行修正。修正的方法：
Python实现R^2:
import numpy as np
from astropy.units import Ybarn
import mathdef computeCorrelation(X, Y):xBar = np.mean(X)yBar = np.mean(Y)SSR = 0varX = 0varY = 0for i in range(0 , len(X)):diffXXBar = X[i] - xBardiffYYBar = Y[i] - yBarSSR += (diffXXBar * diffYYBar)varX += diffXXBar**2varY += diffYYBar**2SST = math.sqrt(varX * varY)return SSR / SSTtestX = [1, 3, 8, 7, 9]
testY = [10, 12, 24, 21, 34]print computeCorrelation(testX, testY)
总结: 当计算出回归方程, 可以使用相关系数和决定系数来判断当前回归方程的正确率
