线性回归及梯度下降
我们第一个学习的算法是线性回归,通过本篇可以看到这个算法的概况,更重要的是将会了解监督学习的完整流程。
模型函数
我们先从一个例子开始,假如我们拿到了这样一张表格,上面包含了纽约若干程序员的年薪:
那么如果你的工作经验是3.5年,年薪应该是多少呢?
显而易见,只有工作经验和年薪是不同,那我们就把这两项抽取出来,首先通过可视化简单的看看一下这两个的关系。
我们将要用来描述这个回归问题的标记如下:
- 𝑚 代表训练集中实例的数量
- 𝑥 代表特征/输入变量
- 𝑦 代表目标变量/输出变量
- (𝑥,𝑦) 代表训练集中的实例
-
代表第 𝑖 个观察实例
h 代表学习算法的解决方案或函数也称为假设(hypothesis)
这是一个监督学习算法,我们把数据集喂给我们的学习算法,学习算法会输出一个函数,通常表示为小写 h 表示。h 代表 hypothesis(假设),表示 h 根据输入的 𝑥 值来得出 𝑦 值,𝑦 值对应年薪。 因此,h 是一个从 𝑥 到 𝑦 的函数映射。
因而,要解决年薪预测问题,我们实际上是要将已知的数据喂给我们的学习算法,进而学习得到一个假设 h,那么, 对于我们的房价预测问题,我们该如何表达 h?
要回答这个问题之前,我们先来看一下数据集的大概样子,是否可以看出什么?
x = range(11)
y = [103100, 104900, 106800, 108700, 110400,
112300, 114200, 116100, 117800, 119700, 121600]
fig = plt.figure(figsize=(10,6))
ax = fig.add_subplot(111)
ax.scatter(x, y)
plt.show()
从上面简单的可视化,很明显工作经验(x)越长,年薪(y)越高,并且看起来是符合一个线性关系。一种可能的表达方式(模型函数)是:
因为只含有一个输入变量(特征),因此这样的问题叫作单变量线性回归问题。
代价函数(Cost Fuction)
既然已经有了模型函数,我们现在要做的便是为我们的模型选择合适的参数 θ_0 和 θ_1,我们选择的参数决定了我们得到的直线相对于我们的训练集的准确程度,模型所预测的值与训练集中实际值之间的差距(下图中蓝线所指)就是建模误差(modeling error),误差线就是下图中的蓝线
为了方便我们把建模误差的平方和定义如下, 也就是我们说的代价函数:
image
我们的目标便是选择出可以使得建模误差的平方和能够最小的模型参数。 也就是使得代价函数最小
梯度下降
梯度下降是一个用来求函数最小值的算法,我们将使用梯度下降算法来求出代价函数 𝐽(𝜃0, 𝜃1) 的最小值。
梯度下降的思想是:开始时我们随机选择一个参数的组合(𝜃0, 𝜃1, . . . . . . , 𝜃𝑛),计算代价函数,然后我们寻找下一个能让代价函数值下降最多的参数组合。我们持续这么做直到一个局部最小值(local minimum),因为我们并没有尝试完所有的参数组合,所以不能确 定我们得到的局部最小值是否便是全局最小值(global minimum),选择不同的初始参数组合,可能会找到不同的局部最小值。
批量梯度下降(batch gradient descent)算法的公式为:
其中 α 是学习率(learning rate),它决定了我们沿着能让代价函数下降程度最大的方向向下迈出的步子有多大,在批量梯度下降中,我们每一次都同时让所有的参数减去学习速率乘以代价函数的导数。
接下来我们就用梯度下降的方法来拟合例子中的数据。
x = np.array([[0,106000],
[1,107200],
[2,108400],
[3,109600],
[4,110800],
[5,112000],
[6,113200],
[7,114400],
[8,115600],
[9, 116800],
[10, 118000]]
)
m, n = np.shape(x)
x_data = np.ones((m,n))
x_data[:,:-1] = x[:,:-1]
y_data = x[:,-1]
m, n = np.shape(x_data)
theta = np.ones(n)
cost_arr = []
def batchGradientDescent(maxiter,x,y,theta,alpha):
xTrains = x.transpose()
for i in range(0,maxiter):
hypothesis = np.dot(x,theta)
loss = (hypothesis-y)
gradient = np.dot(xTrains,loss)/m
theta = theta - alpha * gradient
cost = 1.0/2*m*np.sum(np.square(np.dot(x,np.transpose(theta))-y))
cost_arr.append(cost)
if cost < 0.0001:
break
return theta
result = batchGradientDescent(10000,x_data,y_data,theta,0.01)
print ("final: ",result)
newy = np.dot(x_data,result)
fig = plt.figure(figsize=(16, 8))
ax1 = fig.add_subplot(121)
ax1.plot(range(len(cost_arr)), cost_arr, 'ko')
ax2 = fig.add_subplot(122)
ax2.plot(x[:,0],newy, 'k--')
ax2.plot(x[:,0],x[:,1], 'ro')
plt.show()
final = [ 1200.00034431 105999.99760913]
这里的两个值就是我们模型函数中的𝜃0, 𝜃1
梯度下降的曲线和最后拟合的直线如下图:
此次我们使用超参 α = 0.01,一共进行了6500多次的尝试最终收敛到0,为了让收敛的速度加快,我们可以调整超参,比如 α = 0.03,就只进行了2000多次就收敛到0。
α = 0.03
多变量线性回归
之前我们研究的都是单变量的回归模型,当我们增加更多的特征,例如年龄、职位等就构成了一个含有多个变量的模型,模型中的特征为
增加更多特征后,我们引入一些新的注释:
𝑛 代表特征的数量
𝑥(𝑖)代表第 𝑖 个训练实例,是特征矩阵中的第𝑖行,是一个向量(vector)。
比方说,上图的
这个公式中有 𝑛 + 1个参数和 𝑛 个变量,为了使得公式能够简化一些,引入 𝑥0 = 1,则公式转化为:
此时模型中的参数是一个𝑛 + 1维的向量,任何一个训练实例也都是𝑛 + 1维的向量,特征矩阵 𝑋 的维度是 𝑚 ∗ (𝑛 + 1)。 因此公式可以简化为:
其中上标 𝑇 代表矩阵转置。
正规方程
最小二乘法可以将误差方程转化为有确定解的代数方程组(其方程式数目正好等于未知数的个数),从而可求解出这些未知参数。这个有确定解的代数方程组称为最小二乘法估计的正规方程。
前面我们讲过,Cost Fuction 是:
其中:
将向量表达形式转为矩阵表达形式,则有:
其中𝑋为𝑚行𝑛列的矩阵(𝑚为样本个数,𝑛为特征个数),𝜃为𝑛行1列的矩阵,𝑦为𝑚行1列的矩阵,对𝐽(𝜃)进行如下变换:
接下来对𝐽(𝜃)偏导,需要用到以下几个矩阵的求导法则:
所以有:
令
则有:
Python中对于矩阵的各种操作可以通过 Numpy 库的一些方法来实现,非常方便。但在这个代码实现中需要注意:X矩阵不能为奇异矩阵,否则是无法求解矩阵的逆的。下面是最小二乘法的代码实现部分。
def standRegres(xArr,yArr):
"""
函数说明:计算回归系数theta
Parameters:
xArr - x数据集
yArr - y数据集
Returns:
theta - 回归系数
"""
xMat = np.mat(xArr); yMat = np.mat(yArr).T
#根据文中推导的公示计算回归系数
xTx = xMat.T * xMat
if np.linalg.det(xTx) == 0.0:
print("矩阵为奇异矩阵,不能求逆")
return
theta = xTx.I * (xMat.T*yMat)
return theta
最小二乘法 vs 梯度下降法
二者都对损失函数的回归系数进行了求偏导,并且所得到的推导结果是相同的,那么究竟哪里不同呢?
最小二乘法通过使推导结果等于0,从而直接求得极值,而梯度下降则是将推导结果带入迭代公式中,一步一步地得到最终结果。简单地说,最小二乘法是一步到位的,而梯度下降是一步步进行的。