机器学习第五周-梯度下降

最新推荐文章于 2023-08-14 15:24:42 发布

笨鸟的GPS

最新推荐文章于 2023-08-14 15:24:42 发布

阅读量321

点赞数

分类专栏：机器学习

原文链接：https://mp.weixin.qq.com/s/OUslRwKGpS29gncsiyAPyg

版权

机器学习专栏收录该内容

11 篇文章 0 订阅

订阅专栏

学习内容

梯度下降
手动实现
线性回归下的梯度下降
批量和随机梯度下降

一、梯度下降

机器学习的“三板斧”
第一步：定义一个函数几何
第二步：判断函数的好坏
第三步：选择最好的函数
我们的目标是让损失函数最小化，梯度下降是目前机器学习、深度学习解决最优化问题的算法中，最核心、应用最广的方法

梯度下降(Grandient Descent,GD),不是一个机器学习算法，而是一种基于搜索的最优化方法。其作用是用来对原始模型的损失函数进行优化，以便寻找到最优的参数，使得损失函数的值最小。

更聪明的算法，从损失值出发，去更新参数，且要大幅降低计算次数。梯度下降算法作为一个聪明很多的算法，抓住了参数与损失值之间的导数，也就是能够计算梯度（gradient）,通过导数告诉我们此时此刻某参数应该朝什么方向，以怎样的速度运动，能安全高效降低损失值，朝最小损失值靠近

多元函数的导数（derivative）就是梯度（gradient），分别对每个变量进行微分，然后用逗号分割开，梯度就是用括号抱起来，梯度其实是一个向量，比如我们所说的损失函数L的梯度为
在这里插入图片描述

梯度指向误差值增加最快的方向，导数为0（梯度为0向量）的点，就是优化问题的解，每次搜索的步长为某个特定的数值，直到梯度与0向量非常接近为止。

场景假设：选择一个方向下山，那么这个方向怎么选？每次该怎么走？选方向在算法中是以随机方式给出的，总结起来就是随机选择一个方向，然后每次迈步都选择最陡的方向，直到这个方向能达到的最低点

数学推导
在这里插入图片描述

关于参数
在学术上，我们称之为“学习率”(learning rate)，是模型训练时的一个很重要的超参数，能直接影响算法的正确性和效率
首先，学习率不能太大。因此从数学角度上来说，一阶泰勒公式只是一个近似的公式，只有在学习率很小，也就是很小时才成立。并且从直观上来说，如果学习率太大，那么有可能会“迈过”最低点，从而发生“摇摆”的现象（不收敛），无法得到最低点
其次，学习率又不能太小。如果太小，会导致每次迭代时，参数几乎不变化，收敛学习速度变慢，使得算法的效率降低，需要很长时间才能达到最低点

致命问题
从理论上，它只能保证达到局部最低点，而非全局最低点。解决方案：首先随机产生多个初始参数集，分别对每个初始参数集使用梯度下降法，直到函数值收敛于某个值，从这些值中找出最小值。对于梯度下降来说，初始点的位置，也是一个超参数

二、手动实现

1.Scipy库中的derivative方法
scipy.misc.derivative(func, x0, dx=1.0, n=1, args=(), order=3)[source]

def f(x):
    return x**3 + x**2
derivative(f, 1.0, dx=1e-6)

输出：
4.9999999999217337

2.Sympy库中的diff方法
sympy是符号化运算库，能够实现表达式的求导。所谓符号化，是将数学公式以直观符号的形式输出

from sympy import *

# 符号化变量
x = Symbol('x')

func = 1/(1+x**2)

print("x:", type(x))
print(func)
print(diff(func, x))
print(diff(func, x).subs(x, 3))
print(diff(func, x).subs(x, 3).evalf())

输出结果：
x: <class ‘sympy.core.symbol.Symbol’>
1/(x**2 + 1)
-2*x/(x **2 + 1)**2
-3/50
-0.0600000000000000

https://www.cnblogs.com/zyg123/p/10535121.html

3.模拟实现梯度下降

import numpy as np
import matplotlib.pyplot as plt
from scipy.misc import derivative

def lossFunction(x):
    return (x-2.5)**2-1

# 在-1到6的范围内构建140个点
plot_x = np.linspace(-1,6,141)
# plot_y 是对应的损失函数值
plot_y = lossFunction(plot_x)

plt.plot(plot_x,plot_y)
plt.show()

在这里插入图片描述
已知梯度下降的本质是多元函数的导数，这了定义一个求导的方法，使用的是scipy库中的derivative方法。

"""
算法：计算损失函数J在当前点的对应导数
输入：当前数据点theta
输出：点在损失函数上的导数
"""
def dLF(theta):
    return derivative(lossFunction, theta, dx=1e-6)

在这里插入图片描述

theta = 0.0
eta = 0.1
epsilon = 1e-6
while True:
    # 每一轮循环后，要求当前这个点的梯度是多少
    gradient = dLF(theta)
    last_theta = theta
    # 移动点，沿梯度的反方向移动步长eta
    theta = theta - eta * gradient
    # 判断theta是否达到最小值
    # 因为梯度在不断下降，因此新theta的损失函数在不断减小
    # 看差值是否达到了要求
    if(abs(lossFunction(theta) - lossFunction(last_theta)) < epsilon):
        break
print(theta)
print(lossFunction(theta))

输出：
2.498732349398569
-0.9999983930619527

下面创建存放theta点位置的列表，绘图

def gradient_descent(initial_theta, eta, epsilon=1e-6):
    theta = initial_theta
    theta_history.append(theta)
    while True:
        # 每一轮循环后，要求当前这个点的梯度是多少
        gradient = dLF(theta)
        last_theta = theta
        # 移动点，沿梯度的反方向移动步长eta
        theta = theta - eta * gradient
        theta_history.append(theta)
        # 判断theta是否达到损失函数最小值的位置
        if(abs(lossFunction(theta) - lossFunction(last_theta)) < epsilon):
            break

def plot_theta_history():
    plt.plot(plot_x,plot_y)
    plt.plot(np.array(theta_history), lossFunction(np.array(theta_history)), color='red', marker='o')
    plt.show()

eta=0.1
theta_history = []
gradient_descent(0., eta)
plot_theta_history()
print("梯度下降查找次数：",len(theta_history))

在这里插入图片描述
我们发现，在刚开始时移动比较大，这是因为学习率是一定的，再乘上梯度本身数值大（比较陡），后来梯度数值小（平缓）所以移动的比较小。且经历了34次查找。

eta=0.01
theta_history = []
gradient_descent(0., eta)
plot_theta_history()
print("梯度下降查找次数：",len(theta_history))

在这里插入图片描述

eta=0.9
theta_history = []
gradient_descent(0., eta)
plot_theta_history()
print("梯度下降查找次数：",len(theta_history))

在这里插入图片描述可见在一定范围内将学习率调大，还是会逐渐收敛的。但是我们要注意，如果学习率调的过大，一步迈到“损失函数值增加”的点上去了，在错误的道路上越走越远（如下图所示），就会导致不收敛，会报OverflowError的异常。

为避免报错，源代码修改

def lossFunction(x):
    try:
        return (x-2.5)**2-1
    except:
        return float('inf') #正无穷

def gradient_descent(initial_theta, eta, n_iters, epsilon=1e-6):
    theta = initial_theta
    theta_history.append(theta)
    i_iters = 0
    while i_iters < n_iters:
        gradient = dLF(theta)
        last_theta = theta
        theta = theta - eta * gradient
        theta_history.append(theta)
        if(abs(lossFunction(theta) - lossFunction(last_theta)) < epsilon):
            break
        i_iters += 1

三、线性回归下的梯度下降

在这里插入图片描述

上式梯度大小实际上是和样本数量m相关，m越大，累加之和越大，这是不合理的；为了使与m无关，除以m

梯度下降代码

def fit_gd(self, X_train, y_train, eta=0.01, n_iters=1e4):
    """根据训练数据集X_train, y_train, 使用梯度下降法训练Linear Regression模型"""
    assert X_train.shape[0] == y_train.shape[0], \
        "the size of X_train must be equal to the size of y_train"

    def J(theta, X_b, y):
        try:
            return np.sum((y - X_b.dot(theta)) ** 2) / len(y)
        except:
            return float('inf')
        
    def dJ(theta, X_b, y):
        return X_b.T.dot(X_b.dot(theta) - y) * 2. / len(y)

    def gradient_descent(X_b, y, initial_theta, eta, n_iters=1e4, epsilon=1e-8):

        theta = initial_theta
        cur_iter = 0

        while cur_iter < n_iters:
            gradient = dJ(theta, X_b, y)
            last_theta = theta
            theta = theta - eta * gradient
            if (abs(J(theta, X_b, y) - J(last_theta, X_b, y)) < epsilon):
                break

            cur_iter += 1

        return theta

    X_b = np.hstack([np.ones((len(X_train), 1)), X_train])
    initial_theta = np.zeros(X_b.shape[1])
    self._theta = gradient_descent(X_b, y_train, initial_theta, eta, n_iters)

    self.intercept_ = self._theta[0]
    self.coef_ = self._theta[1:]

    return self

使用真实数据

import numpy as np
from sklearn import datasets

boston = datasets.load_boston()
X = boston.data
y = boston.target

X = X[y < 50.0]
y = y[y < 50.0]

首先使用线性回归中的标准方程解进行计算，得到相应系数与截距的最优值，然后计算回归评分（R方）

from myAlgorithm.LinearRegression import LinearRegression
from myAlgorithm.model_selection import train_test_split

X_train, X_test, y_train, y_test = train_test_split(X, y, seed=666)

lin_reg1 = LinearRegression()
%time lin_reg1.fit_normal(X_train, y_train)
lin_reg1.score(X_test, y_test)

输出：
CPU times: user 25.7 ms, sys: 25.5 ms, total: 51.2 ms
Wall time: 116 ms
0.8129802602658466

然后创建一个lin_reg2 ,使用梯度下降法得到最优参数：

lin_reg2 = LinearRegression()
lin_reg2.fit_gd(X_train, y_train, )
lin_reg2.coef_
"""
输出：
array([nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan,nan])
"""

执行过程会出现overflow的报错，并且得到的theta为空
这是因为，在一个真实的数据中，查看前3行数据，就可以观察到，每一个特征所对应的规模是不一样的，有一些很小，有一些很大，使用默认的学习率可能过大，导致不收敛

传递一个很小的学习率

lin_reg2.fit_gd(X_train, y_train, eta=0.000001)
lin_reg2.score(X_test, y_test)
"""
输出：
0.27556634853389195
"""

我们发现得到的R方评分很差，这可能是因为学习率设置的过小，步长太小没有到最优值。为了验证这个假设，可以再进行一次验证：

lin_reg2.fit_gd(X_train, y_train, eta=0.000001, n_iters=1e6)
lin_reg2.score(X_test, y_test)
"""
输出：
0.75418523539807636
"""

因为真实数据整体不在一个规模上，需要在梯度下降之前进行数据归一化操作
在这里插入图片描述

数据归一化

from sklearn.preprocessing import StandardScaler

standardScaler = StandardScaler()
standardScaler.fit(X_train)
X_train_std = standardScaler.transform(X_train)
lin_reg3 = LinearRegression()
lin_reg3.fit_gd(X_train_std, y_train)
X_test_std = standardScaler.transform(X_test)
lin_reg2.score(X_test, y_test)
"""
输出：
0.8129802602658466
"""

这和“正规方程解”得到的score一样，这就说明我们找到了损失函数的最小值。

四、速度更快的随机梯度下降法

每次更新参数时是需要计算所有样本的，通过对整个数据集的所有样本的计算来求解梯度的方向。这种计算方法被称为：批量梯度下降法BGD(Batch Gradient Descent)。但是这种方法在数据量很大时需要计算很久

随机梯度下降法SGD（stochastic gradient descent），随机梯度下降是每次迭代使用一个样本来对参数进行更新。虽然不是每次迭代得到的损失函数都向着全局最优方向，但是大的整体的方向是向全局最优解的，最终的结果往往是在全局最优解附近

在这里插入图片描述

随机下降与学习率的取值：其过程就是每次随机取出一个i,得到一个向量，沿着这个随机产生的向量的方向进行搜索，不停的迭代，得到的损失函数的最小值

随机梯度下降法的搜索过程如下图，如果是批量搜索，那么每次都是沿着一个方向前进的，但是随机梯度下降法由于不能保证随机选择的方向是损失函数减小的方向，更不能保证一定是减小速度最快的方向，所以搜索路径就会呈现下图的态势。随机梯度下降有着不可预知性。
在这里插入图片描述
但实验结论告诉我们，通过随机梯度下降法，依然能够达到最小值的附件（用精度换速度）

随机梯度下降的过程，学习率取值很重要，设置学习率是逐渐递减的，设计一个函数，使学习率随着下降循环次数的增加而减小。

我们想到最简单的表示方法是一个倒数的形式，不过这样也会有问题，如果循环次数比较小的时候，学习率下降的太快了，比如循环次数为1变为2，则学习率减少50%。因此我们可以将分子变为常数，并在分母上增加一个常数项来来缓解初始情况下，学习率变化太大的情况，且更灵活。

批量梯度下降代码：

import numpy as np
import matplotlib.pyplot as plt

m = 100000
x = np.random.normal(size=m)
X = x.reshape(-1,1)
y = 4. * x + 3. + np.random.normal(0, 3, size=m)

def J(theta, X_b, y):
    try:
        return np.sum((y-X_b.dot(theta)) ** 2) / len(y)
    except:
        return float('inf')

def dJ(theta, X_b, y):
    return X_b.T.dot(X_b.dot(theta) - y) * 2. / len(y)

def gradient_descent(X_b, y, initial_theta, eta, n_iters=1e4, epsilon=1e-6):
    theta = initial_theta
    cur_iter = 0
    while(cur_iter < n_iters):
        gradient = dJ(theta, X_b, y)
        last_theta = theta
        theta = theta - eta * gradient
        if (abs(J(theta, X_b, y) - J(last_theta, X_b, y)) < epsilon):
            break
        cur_iter += 1
    return theta

%%time
X_b = np.hstack([np.ones((len(X),1)),X])
initial_theta = np.zeros(X_b.shape[1])
eta = 0.01
theta = gradient_descent(X_b, y, initial_theta, eta)
theta
"""
CPU times: user 2.3 s, sys: 493 ms, total: 2.79 s
Wall time: 1.73 s
array([3.00328612, 4.00425536])
"""

随机梯度下降代码：

# 传递的不是整个矩阵X_b，而是其中一行X_b_i;传递y其中的一个数值y_i
def dJ_sgd(theta, X_b_i, y_i):
    return X_b_i.T.dot(X_b_i.dot(theta) - y_i) * 2.

def sgd(X_b, y, initial_theta, n_iters):
    t0 = 5
    t1 = 50
    # 
    def learning_rate(cur_iter):
        return t0 / (cur_iter + t1)
    theta = initial_theta

    for cur_iter in range(n_iters):
        # 随机找到一个样本（得到其索引）
        rand_i = np.random.randint(len(X_b))
        gradient = dJ_sgd(theta, X_b[rand_i], y[rand_i])
        theta = theta - learning_rate(cur_iter) * gradient
    return theta

%%time
X_b = np.hstack([np.ones((len(X),1)),X])
initial_theta = np.zeros(X_b.shape[1])
theta = sgd(X_b, y, initial_theta, n_iters=len(X_b)//3)
print(theta)

"""
输出：
[2.9287233  4.05019715]
CPU times: user 296 ms, sys: 6.49 ms, total: 302 ms
Wall time: 300 ms
"""

通过简单的例子，我们可以看出，在批量梯度下降的过程中消耗时间1.73 s，而在随机梯度下降中，只使用300 ms。并且随机梯度下降中只考虑的三分一的样本量，且得到的结果一定达到了局部最小值的范围内。

随机梯度的样本量代码改进：
在之前，只是简单的只考虑的三分一的样本量，实际上我们已经考虑所有的样本量n次，然后在每次考虑样本量时采用随机的方式。我们首先要了解：表示对所有样本的循环次数，所有样本的个数为m。因此我们在循环时，应该使用双重循环，即外层循环为每次循环所有样本，内层循环为在所有样本中进行随机选择，次数为样本数量。那么要注意：在计算学习率时，次数就变为

for i_iter in range(n_iters):
                # 将原本的数据随机打乱，然后再按顺序取值就相当于随机取值
                indexes = np.random.permutation(m)
                X_b_new = X_b[indexes,:]
                y_new = y[indexes]
                for i in range(m):
                    gradient = dJ_sgd(theta, X_b_new[i], y_new[i])
                    theta = theta - learning_rate(i_iter * m + i) * gradient

Sklearn中的SGD

from sklearn.linear_model import SGDRegressor

sgd_reg = SGDRegressor()    # 默认n_iter=5
%time sgd_reg.fit(X_train_std, y_train)
sgd_reg.score(X_test_std, y_test)

sgd_reg = SGDRegressor(n_iter=100)
%time sgd_reg.fit(X_train_std, y_train)
sgd_reg.score(X_test_std, y_test)

sklearn中的算法实现，其实和我们的实现有很大的区别，进行了很多的优化。

总结：
批量梯度下降法BGD(Batch Gradient Descent)。
优点：全局最优解；易于并行实现；
缺点：当样本数据很多时，计算量开销大，计算速度慢。
针对于上述缺点，其实有一种更好的方法：随机梯度下降法SGD（stochastic gradient descent），随机梯度下降是每次迭代使用一个样本来对参数进行更新。
优点：计算速度快；
缺点：收敛性能不好

##############################################################################################################
@ 2019.12.08 木居居士的机器学习小组第五周打卡
安利公益监督学习组织 - 【公众号】数据科学家联盟
https://mp.weixin.qq.com/s/1WWmbLZucz9vIp-4tKKQ5Q
感谢木东大佬、饼干大佬、南头大佬、星空妹砸、Desitiny、DD的无私付出，抱拳ing~

笨鸟的GPS

关注

0
点赞
踩
0

收藏

觉得还不错? 一键收藏
0
评论
机器学习第五周-梯度下降

学习内容梯度下降手动实现线性回归下的梯度下降批量和随机梯度下降一、梯度下降机器学习的“三板斧”第一步：定义一个函数几何第二步：判断函数的好坏第三步：选择最好的函数我们的目标是让损失函数最小化，梯度下降是目前机器学习、深度学习解决最优化问题的算法中，最核心、应用最广的方法梯度下降(Grandient Descent,GD),不是一个机器学习算法，而是一种基于搜索的最优化方法...
复制链接

扫一扫