机器学习第四章线性回归算法进阶4.2梯度下降法求解多变量线性回归（《大话Python机器学习》学习笔记）

第四章线性回归算法进阶

4.2梯度下降法求解多变量线性回归

  梯度下降法是对最小二乘法进行优化求解回归的一种算法，它采用了迭代的形式来寻找成本函数J(θ)的最小值。其中J(θ)：

4.2.1梯度下降的含义

定义：
  来自于数学中的微积分，通过对多元函数参数求偏导数，把求得的各参数的偏导数以向量的形式写出来就是梯度。
几何意义：
  函数变化增加最快的地方
  梯度上升：
    沿着梯度向量的方向更容易找到函数的最大值。
  梯度下降：
    沿着梯度向量相反的方向，梯度减小最快，更容易找到函数的最小值。
梯度下降基本含义：
  是一种最优化算法，通过采用迭代的方法求解目标函数得到最优解；
  是在成本函数的基础上，利用梯度迭代求出局部最优解。

4.2.2梯度下降相关概念

• 学习率（Learning Rate）：也称为步长，它决定了在梯度下降迭代的过程中，每一步沿梯度负方向前进的长度。
• 特征（Feature）：样本中输入部分
• 假设函数（Hypothesis Function）：在监督学习中，为了拟合输入样本而使用的假设函数，记作h_θ(x)
• 损失/成本函数（Loss Function）：为了评估模型拟合的好坏，用损失函数来度量拟合的程度。线性回归中，损失函数通常为样本输出和假设函数的差取平方。
    在假设函数中，α决定了逼近最低点的速率；若α过大，可能会出现在最低点附近反复震荡的情况，无法收敛；若α过小，而导致逼近的速率太慢，需要迭代多次才能逼近最低点。因此，当梯度下降法选择α时，可以大致以3的倍数再以10的倍数的规律来选出一组α，从而找到使得J(θ)快速下降收敛的α。
  
  

4.2.3梯度下降法求解线性回归算法

求解思路：
  对θ取一随机初始值，可以是全零的向量，然后不断迭代改变θ的值使其成本函数根据梯度下降的方向减小，知道收敛求出某θ值使得J(θ)最小或是局部最小。
  更新规则：

  J(θ)：对θ的偏导决定了梯度下降的方向，将J(θ)带入更新规则中得到：

  根据多变量线性回归的表达式h_θ(x)=θ₀+θ₁x₁+θ₂x₂+θ₃x₃+…+θ_nx_n，在对θ₀，θ₁，θ₂，θ₃，…θ_n取一个随机初始值的基础上，得到初始预测值h_θ(x)，然后计算更新规则式（上图），得到新的θ_j，在此基础上计算成本函数或者损失函数J(θ)的值。
  接下来再根据更新规则得到新的θ_j，如此循环，最终得到成本函数或者损失函数J(θ)的值达到最小，且处于收敛状态时，即为梯度下降法求解线性回归算法得到的最优解。
批量梯度下降法
  对于上图算式，由于每一次迭代都需要遍历所有的训练数据，如果训练数据庞大，那么每次迭代运算的复杂度比较高，使得收敛速度很慢。
随机梯度下降法
  为了解决批量梯度下降法的缺点，随机梯度下降法在更新参数时，不必遍历全部训练数据，只要随机选取其中一个训练数据进行参数更新，收敛速度快，减少复杂度。

1. 由于迭代所使用的训练数据不是全部数据，因此每一步下降的过程中，成本函数J(θ)值的优化不是按照规律逐步减少的，反而会有波动以及大幅跳跃的现象。
2. 随机梯度下降法每一次更新规则之后，J(θ)值的下降是极其不规则，而且有向上的波动。
     好处：对于类似的盆地区域（很多局部极小值点）该波动的特点可能会使得优化的方向从当前的局部极小值点跳到另一个更好的局部极小值点，从而使得非凸函数最终收敛于一个较好的局部极值点，甚至全局极值点。
     而在批量梯度下降法中，存在多个局部最优解的情况下，有时得到的最优解不一定是全局最优解，而是局部最优解，对于这种情况，可能需要随机初始化θ，得到多个最优值，再从中筛选。

4.2.4梯度下降法的Python实现：影厅观影人数的拟合（二）

1. 导入数据

df = pd.read_csv('D:/PythonProject/machine/data/3_film.csv')
print(df)

filmnum filmsize ratio quality
0 45 106 17 6
1 44 99 15 18
2 61 149 27 10
3 41 97 27 16
4 54 148 30 8
… … … … …
121 49 115 23 2
122 37 78 20 4
123 64 165 48 8
124 41 84 16 3
125 40 102 24 1
[126 rows x 4 columns]

2. 插入1列全为1的数组

# 在df第1列和第2列之家插入一列全是1的数组
df.insert(1,'Ones',1)
print(df)

filmnum Ones filmsize ratio quality
0 45 1 106 17 6
1 44 1 99 15 18
2 61 1 149 27 10
3 41 1 97 27 16
4 54 1 148 30 8
… … … … … …
121 49 1 115 23 2
122 37 1 78 20 4
123 64 1 165 48 8
124 41 1 84 16 3
125 40 1 102 24 1
[126 rows x 5 columns]

3. 选取特征向量与响应向量，并划分数据

# 计算df的列数
cols = df.shape[1]
# 取数据df的第2列之后的数据作为X变量
X = df.iloc[:,1:cols]
# 取数据df的第一列为y变量
y = df.iloc[:,0:1]
# 把X，y转换成数组形式，便于计算
X = np.array(X.values)
y = np.array(y.values)
# 以25%的数据构建测试样本，剩余作为训练样本
X_train,X_test,y_train,y_test = train_test_split(X,y,test_size=0.25,random_state=1)
print(X_train.shape,X_test.shape,y_train.shape,y_test.shape)

(94, 4) (32, 4) (94, 1) (32, 1)

4. 构建成本函数

def computeCost(X,y,theta):
    inner = np.power(((X * theta.T) - y),2)
    return np.sum(inner) / (2 * len(X))

5. 构建梯度下降法求解函数

"""
    X、y：输入变量
    theta：参数
    alpha：学习率
    iters：梯度下降迭代次数
"""
def gradientDescent(X,y,theta,alpha,iters):
    # 构建零值矩阵
    temp = np.matrix(np.zeros(theta.shape))
    # 计算需要求解参数的个数
    parameters = int(theta.ravel().shape[1])
    # 构建iters个0的数组
    cost = np.zeros(iters)

    for i in range(iters):
        # 计算hθ(x)-y
        error = (X * theta.T) - y
        # 对于theta中每一个元素依次计算：
        for j in range(parameters):
            # 计算两矩阵相乘(hθ(x)-y)x
            term = np.multiply(error,X[:,j])
            # 更新法则
            temp[0,j] = theta[0,j] - ((alpha / len(X)) * np.sum(term))

        theta = temp
        # 基于求出来的theta求解成本函数
        cost[i] = computeCost(X,y,theta)
    return theta,cost

6. 设定相关参数的初始值，并代入gradientDescent()函数中求解

# 设定学习率
alpha = 0.000001
# 设定迭代次数
iters = 100
theta = np.matrix(np.array([0,0,0,0]))
# 采用gradientDescent()函数来优化求解，带入初始值并求解
g,cost = gradientDescent(X,y,theta,alpha,iters)
print(g)

[[0.00354631 0.40050803 0.08201549 0.0328456 ]]

7. 对预测集X_test进行预测，并对结果进行评价

# 求出预测集y_test的预测值
y_hat = X_test * g.T
# 设置图片尺寸
plt.figure(figsize=(10,6))
# 创建t变量
t = np.arange(len(X_test))
# 绘制y_test曲线
plt.plot(t,y_test,'r',linewidth=2,label='y_test')
# 绘制y_hat曲线
plt.plot(t,y_hat,'g',linewidth=2,label='y_hat')
# 设置图例
plt.legend()
plt.show()

8. 对预测结果进行评价

# 拟合优度R2的输出方法二
print('r2_2={}'.format(r2_score(y_test,y_hat)))
# 计算MAE
print('MAE={}'.format(metrics.mean_absolute_error(y_test,y_hat)))
# 计算MSE
print('MSE={}'.format(metrics.mean_squared_error(y_test,y_hat)))
# 计算RMSE
print('RMSE={}'.format(np.sqrt(metrics.mean_squared_error(y_test,y_hat))))

r2_2=0.8072533281432706
MAE=5.262018943705062
MSE=52.2456418234778
RMSE=7.22811467974034