Python3入门机器学习(六)- 梯度下降法

1. 梯度下降法简介

1-1

以下是定义了一个损失函数以后,参数theta对应的损失函数J的值对应的示例图,我们需要找到使得损失函数值J取得最小值对应的theta(这里是二维平面,也就是我们的参数只有一个)

在直线方程中,导数代表斜率 在曲线方程中,导数代表切线斜率 导数代表theta单位变化时,J相应的变化

1-2

1-3

η太小,会减慢收敛学习速度

1-4

η太大,甚至导致不收敛

1-5

其他注意事项

  • 并不是所有函数都有唯一的极值点

1-6 解决方案:

  • 多次运行,随机化初始点
  • 梯度下降法的初始点也是一个超参数

1-7


2 梯度下降法模拟

2.1 实现

import numpy as np
import matplotlib.pyplot as plt
# 简单模拟一个损失函数
plot_x = np.linspace(-1,6,141)
plot_y = (plot_x-2.5)**2-1
# 绘制我们模拟的损失函数
plt.plot(plot_x,plot_y)

2.1-1

def dJ(theta):
    """损失函数的导数"""
    return 2*(theta-2.5)

def J(theta):
    """损失函数"""
    try:
        return (theta-2.5)**2-1
    except:
        return float('inf')

def gradient_descent(initial_theta,eta,n_iters = 1e4,epsilon=1e-8):
    """
    梯度下降法封装
    initial_theta:初始化的theta值
    eta:学习率η
    n_iters: 最大循环次数
    epsilon: 精度
    """
    theta = initial_theta
    # theta_history 保存theta的变化值
    theta_history.append(initial_theta)
    i_iters = 0
    
    while i_iters<n_iters:
        """
        如果theta两次变化之间的损失函数值的变化小于我们定义的精度
        则可以说明我们已经找到了最低的损失函数值和对应的theta
        
        如果循环次数超过了我们设置的循环次数,
        则说明可能由于η设置的过大导致无止境的循环
        """
        gradient = dJ(theta)
        last_theta = theta
        theta = theta - eta * gradient
        theta_history.append(theta)
        
        if (abs(J(theta)-J(last_theta)) < epsilon):
            break
        i_iters += 1

def plot_theta_history():
    plt.plot(plot_x,J(plot_x))
    plt.plot(np.array(theta_history),J(np.array(theta_history)),color='r',marker='+')
    print("the size of theta_history is %d"%len(theta_history))

2.2 使用不同η学习率测试并观察我们的梯度下降法的结果

eta = 0.1
theta_history = []
gradient_descent(0.,eta)
plot_theta_history()

2.2-1

eta = 0.01
theta_history = []
gradient_descent(0.,eta)
plot_theta_history()

2.2-2

eta = 0.001
theta_history = []
gradient_descent(0.,eta)
plot_theta_history()

2.2-3

eta = 0.8
theta_history = []
gradient_descent(0.,eta)
plot_theta_history()

2.2-4

可以发现,只要η不超过一个限度,我们编写的函数都可以在有限次数之后找到最优解,并且η越小,学习的次数越多 下面来看一下,如果eta取较大值1.1,会出现什么情况

eta = 1.1
theta_history = []
gradient_descent(0.,eta)
# 数据量太大会报错
# plot_theta_history()
print(len(theta_history))
# 输出10001
theta_history[-1]
# 输出 nan(not a number)

可以看出当我们的eta取1.1,函数会循环直至终止,这是由于,我们的η设置过大,导致每次循环过后,损失函数j的值都向大的方向变化

eta = 1.1
theta_history = []
gradient_descent(0.,eta,n_iters = 10)
plot_theta_history()

2.2-5


3 多元线性回归中的梯度下降法

3-1

一个三维空间中的梯度下降法(x,y为系数,z为损失函数)

3-2

推导过程

3-3

3-4

上面推导出的式子的大小是和样本数有关的,m越大,结果越大,这是不合理的,我们希望和m无关

3-5


4 线性回归中的梯度下降法的实现

4.1 比较笨的方法实现

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(X_b)
            except:
                return float('inf')

        def dJ(theta, X_b, y):
            res = np.empty(len(theta))
            res[0] = np.sum(X_b.dot(theta) - y)

            for i in range(1, len(theta)):
                res[i] = np.sum((X_b.dot(theta) - y).dot(X_b[:, i]))

            return res * 2 / len(X_b)

        def gradient_descent(X_b, y, initial_theta, eta, n_iters=n_iters, epsilon=1e-8):
            """
            梯度下降法封装
            X_b: X特征矩阵
            y: 结果向量
            initial_theta:初始化的theta值
            eta:学习率η
            n_iters: 最大循环次数
            epsilon: 精度
            """
            theta = initial_theta
            i_iters = 0

            while i_iters < n_iters:
                """
                如果theta两次变化之间的损失函数值的变化小于我们定义的精度
                则可以说明我们已经找到了最低的损失函数值和对应的theta
                
                如果循环次数超过了我们设置的循环次数,
                则说明可能由于η设置的过大导致无止境的循环
                """
                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

                i_iters += 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)

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

        return self

4.2 测试我们的算法

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(666)
x = 2*np.random.random(size = 100)
# 定义截距为4 斜率为3
y = x * 3. + 4. + np.random.normal(size=100)
X = x.reshape(-1,1)
plt.scatter(x,y)

4.2-1

from machine_learning.LinearRegression import LinearRegression
lin_reg = LinearRegression()
lin_reg.fit_gd(X,y)

lin_reg.interception_
# 4.021457858204859
lin_reg.coef_
# array([3.00706277])

4.3 向量化

4.3-1

4.3-2

4.3-3

修改之前的求导函数

def dJ(theta, X_b, y):
            # res = np.empty(len(theta))
            # res[0] = np.sum(X_b.dot(theta) - y)
            #
            # for i in range(1, len(theta)):
            #     res[i] = np.sum((X_b.dot(theta) - y).dot(X_b[:, i]))
            #
            # return res * 2 / len(X_b)
            return X_b.T.dot(X_b.dot(theta) - y) * 2. / len(X_b)

使用真实的数据测试

4.3-4

使用真实的数据,调整eta和iters,要么由于eta太小导致无法得出真实的结果,要么由于eta太大导致训练时间加长,这是由于数据的规模在不同的特征上不同,所以我们需要对数据进行归一化

4.4 数据归一化

4.4-1

4.4-2

4.4-3

如果样本数非常多,那么即使使用梯度下降法也会导致速度比较慢,因为在梯度下降法中,每一个样本都要参与运算。这时候需要采用随机梯度下降法,我们将在下一小节进行介绍


5 随机梯度下降法

5.1 随机梯度下降法介绍

5.1

批量梯度下降法带来的一个问题是η的值需要设置的比较小,在样本数比较多的时候导致不是速度特别慢,这时候观察随机梯度下降法损失函数的求导公式,可以发现,我们对每一个Xb都做了求和操作,又在最外面除以了m,那么可以考虑将求和和除以m的两个运算约掉,采用每次使用一个随机的Xb

5.2

5.3

由于我们使用的事随机梯度下降法,所以导致我们的最终结果不会像批量梯度下降法一样准确的朝着一个方向运算,而是曲线行下降,这时候我们就希望,越到下面,η值相应减小,事运算次数变多,从而精确计算结果

5-4

这里使用了模拟退火的思想

5-5

5.2 随机梯度下降法实现

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(t):
        return t0 / (t + 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
        eta = 0.01
        X_b = np.hstack([np.ones((len(X), 1)), X])
        initial_theta = np.zeros(X_b.shape[1])
        # 随机的检查了3分之一个样本总量的样本
        _theta = sgd(X_b, y, initial_theta, n_iters=len(X_b)//3)
        # 输出 CPU times: user 318 ms, sys: 5.22 ms, total: 323 ms
        # Wall time: 337 ms
        # _theta: array([3.03182269, 3.93118623])

5.3 随机梯度下降法的封装和测试

    def fit_sgd(self, X_train, y_train, n_iters=5, t0=5, t1=50):
        """
        根据训练数据集X_train, y_train, 使用随机梯度下降法训练Linear Regression模型
        :param X_train:
        :param y_train:
        :param n_iters: 在随机梯度下降法中,n_iters代表所有的样本会被看几圈
        :param t0:
        :param t1:
        :return:
        """
        assert X_train.shape[0] == y_train.shape[0], \
            "the size of X_train must be equal to the size of y_train"
        assert n_iters >= 1

        def dJ_sgd(theta, X_b_i, y_i):
            """
            去X_b,y 中的随机一个元素进行导数公式的计算
            :param theta:
            :param X_b_i:
            :param y_i:
            :return:
            """
            return X_b_i * (X_b_i.dot(theta) - y_i) * 2.

        def sgd(X_b, y, initial_theta, n_iters, t0=5, t1=50):

            def learning_rate(t):
                """
                计算学习率,t1 为了减慢变化速度,t0为了增加随机性
                :param t: 第t次循环
                :return:
                """
                return t0 / (t + t1)

            theta = initial_theta
            m = len(X_b)

            for cur_iter in range(n_iters):
                # 对X_b进行一个乱序的排序
                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(cur_iter * m + i) * gradient

            return theta

        X_b = np.hstack([np.ones((len(X_train), 1)), X_train])
        initial_theta = np.random.randn(X_b.shape[1])
        self._theta = sgd(X_b, y_train, initial_theta, n_iters, t0, t1)

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

        return self

模拟数据进行测试

5.3-1

真实数据波士顿房价进行测试

5.3-2

5.4 使用Sklearn中的 随机梯度下降法

5.4-1

需要注意的是sklearn中的梯度下降法比我们自己的算法要复杂的多,性能和计算准确度上都比我们的要好,我们的算法只是用来演示过程,具体生产上的使用还是应该使用Sklearn提供的


6 梯度下降法 的调试

6.1 梯度下降法调试的原理

可能我们计算出梯度下降法的公式,并使用python编程实现,预测的过程中程序并没有报错,但是可能我们需要求的梯度的结果是错误的,这个时候需要怎么样去调试发现错误呢。

首先以二维坐标平面为例,一个点(O)的导数就是曲线在这个点的切线的斜率,在这个点两侧各取一个点(AB),那么AB两点对应的直线的斜率应该大体等于O的切线的斜率,并且这A和B的距离越近,那么两条直线的斜率就越接近 事实上,这也正是导数的定义,当函数y=f(x)的自变量x在一点x0上产生一个增量Δx时,函数输出值的增量Δy与自变量增量Δx的比值在Δx趋于0时的极限a如果存在,a即为在x0处的导数,记作f'(x0)或df(x0)/dx

6-1

扩展到多维维度则如下

6-2

梯度下降法调试的实现

np.random.seed(666)
X = np.random.normal(size=(1000,10))
X_b = np.hstack([np.ones((len(X),1)),X])

# 真实的θ值
true_theta = np.arange(1,12,dtype=float)
# np.random.normal(size=1000) 添加噪音
y = X_b.dot(true_theta) + np.random.normal(size=1000)
# 实现J(θ)函数
def J(theta,X_b,y):
    try:
        return np.sum((y-X_b.dot(theta))**2)/len(X_b)
    except:
        return float('inf')
    
# 实现数学推导出的dJ(θ)
def d_J_main(theta,X_b,y):
    return X_b.T.dot(X_b.dot(theta) - y) * 2. /len(X_b)

# 实现debug模式的dJ(θ)
def d_J_debug(theta,X_b,y,epsilon=0.01):
    res = np.empty(len(theta))
    for i in range(len(theta)):
        theta_1 = theta.copy()
        theta_1[i] += epsilon
        theta_2 = theta.copy()
        theta_2[i] -= epsilon
        res[i] = (J(theta_1,X_b,y)-J(theta_2,X_b,y))/(2*epsilon)
    return res
# 批量梯度下降法,d_J为求导函数,作为一个参数传入,用于切换求导策略
def gradient_descent(d_J,X_b, y, initial_theta, eta, n_iters=1e4, epsilon=1e-8):
            theta = initial_theta
            i_iters = 0

            while i_iters < n_iters:
                gradient = d_J(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

                i_iters += 1

            return theta
X_b = np.hstack([np.ones((len(X), 1)), X])
initial_theta = np.zeros(X_b.shape[1])
eta = 0.01
# 使用d_J_debug调试模式求出theta
%time theta = gradient_descent(d_J_debug,X_b, y, initial_theta, eta)
print(theta)
# 使用数学解求出theta
%time theta = gradient_descent(d_J_main,X_b, y, initial_theta, eta)
print(theta)
# 输出结果
CPU times: user 531 ms, sys: 214 ms, total: 745 ms
Wall time: 613 ms
[ 0.94575233  1.98082712  3.06882065  3.94835863  4.97139932  5.9859077
  7.01077392  7.99250414  8.99151383  9.97525811 10.99758484]
CPU times: user 67 ms, sys: 27.6 ms, total: 94.6 ms
Wall time: 76.7 ms
[ 0.94575233  1.98082712  3.06882065  3.94835863  4.97139932  5.9859077
  7.01077392  7.99250414  8.99151383  9.97525811 10.99758484]

由此可以看出,我们的d_J_debug和d_J_main的结果是相近的,所以我们的d_J_main的数学推导是没问题的。 我们可以在真正的机器学习之前,先使用d_J_debug这种调试方式来验证一下我们的d_J_main的结果是否正确,然后再进行机器学习。

d_J_debug是通用的,可以放在任何求导的debug过程中,所以可以作为我们机器学习的工具箱来使用


7.梯度下降法的总结

7.1 小批量

  • 批量梯度下降法
  • 随机梯度下降法 下面来看下二者的对比

维度

批量梯度下降法

随机梯度下降法

计算方式

每次对所有的样本看一遍才可以计算出梯度

每一次只需观察一个样本

速度

稳定性

高,一定可以先向损失函数下降的方式前进

低,每一次的方式不确定,甚至向反方向前进

综合二者的优缺点,有一种新的梯度下降法

7-1

小批量梯度下降法:即,我们每一次不看全部样本那么多,也不是只看一次样本那么少,每次只看k个样本 对于小批量梯度下降法,由多了一个超参数

def fit_lit_sgd(self, X_train, y_train, n_iters=5, t0=5, t1=50,k=10):
        """
        根据训练数据集X_train, y_train, 使用随机梯度下降法训练Linear Regression模型
        :param X_train:
        :param y_train:
        :param n_iters: 在随机梯度下降法中,n_iters代表所有的样本会被看几圈
        :param t0:
        :param t1:
        :param k: 小批量随机下降法的超参数k
        :return:
        """
        assert X_train.shape[0] == y_train.shape[0], \
            "the size of X_train must be equal to the size of y_train"
        assert n_iters >= 1

        def dJ_sgd(theta, X_b_k, y_k):
            """
            去X_b,y 中的随机选择k个元素进行导数公式的计算
            :param theta:
            :param X_b_i:
            :param y_i:
            :return:
            """
            return np.sum((X_b_k * (X_b_k.dot(theta) - y_k) ))* 2/len(X_b_k).

        def sgd(X_b, y, initial_theta, n_iters, t0=5, t1=50):

            def learning_rate(t):
                """
                计算学习率,t1 为了减慢变化速度,t0为了增加随机性
                :param t: 第t次循环
                :return:
                """
                return t0 / (t + t1)

            theta = initial_theta
            m = len(X_b)

            for cur_iter in range(n_iters):
                # 每次看k个元素
                i =0
                while i < m:
                    X_b_new = X_b[i:i+k]
                    y_new = y[i:i+k]
                    gradient = dJ_sgd(theta, X_b_new, y_new)
                    theta = theta - learning_rate(cur_iter * m + i+k) * gradient
                    i = i+k

            return theta

        X_b = np.hstack([np.ones((len(X_train), 1)), X_train])
        initial_theta = np.random.randn(X_b.shape[1])
        self._theta = sgd(X_b, y_train, initial_theta, n_iters, t0, t1)

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

        return self

7.2 随机

随机梯度下降法的优点

7.2-1

7.3 梯度上升法

7.3-1

7.3-2

7.3-3

本文参与腾讯云自媒体分享计划,欢迎正在阅读的你也加入,一起分享。

发表于

我来说两句

0 条评论
登录 后参与评论

相关文章

来自专栏目标检测和深度学习

李飞飞等提出新的迭代视觉推理框架,在ADE上实现8.4 %的绝对提升

译者 | 梁红丽 张蔚敏 编辑 | 明 明 出品 | AI科技大本营 【AI科技大本营导读】近日,李飞飞等提出了一种新的迭代视觉推理框架。该框架超越了目前缺乏推...

396110
来自专栏机器之心

教程 | 如何估算深度神经网络的最优学习率

38450
来自专栏AI科技大本营的专栏

李飞飞等提出新的迭代视觉推理框架,在ADE上实现8.4 %的绝对提升

译者 | 梁红丽 张蔚敏 编辑 | 明 明 【AI科技大本营导读】近日,李飞飞等提出了一种新的迭代视觉推理框架。该框架超越了目前缺乏推理能力的识别系统。该框架包...

36370
来自专栏ATYUN订阅号

了解学习速率以及它如何提高深度学习的表现

学习速率是深度学习中的一个重要的超参数,如何调整学习速率是训练出好模型的关键要素之一。这篇文章将着重说明以下几点: 什么是学习速率? 它的意义是什么? 如何系统...

39750
来自专栏专知

【干货】计算机视觉实战系列06——用Python做图像处理

【导读】专知成员Hui上一次为大家介绍主成分分析(PCA)、以及其在图像上的应用,这一次为大家详细讲解SciPy库的使用以及图像高斯模糊实战。 【干货】计算机视...

436140
来自专栏机器之心

盘点 | 对比图像分类五大方法:KNN、SVM、BPNN、CNN和迁移学习

选自Medium 机器之心编译 参与:蒋思源、黄小天、吴攀 图像分类是人工智能领域的基本研究主题之一,研究者也已经开发了大量用于图像分类的算法。近日,Shiyu...

1.2K80
来自专栏智能算法

结合Scikit-learn介绍几种常用的特征选择方法(下)

5 两种顶层特征选择算法 之所以叫做顶层,是因为他们都是建立在基于模型的特征选择方法基础之上的,例如回归和SVM,在不同的子集上建立模型,然后汇总最...

69740
来自专栏和蔼的张星的图像处理专栏

DSST详解

有一段时间没有看tracking了,前面一个月老师没有找,我也没有看文章,主要去看c++和cs231n去了。上周一老师找了我一次,于是赶紧把tracking又拾...

24120
来自专栏机器学习算法工程师

fine-gained image classification

我们在路边看到萌犬可爱至极,然后却不知道这个是哪种狗;看见路边的一个野花却不知道叫什么名字,吃着一种瓜,却不知道是甜瓜还是香瓜傻傻分不清……

11920
来自专栏智能算法

深度学习三人行(第7期)----深度学习之避免过拟合(正则化)

今天我们一起学习下深度学习中如何避免过拟合,我们多多交流,共同进步。本期主要内容如下:

21640

扫码关注云+社区

领取腾讯云代金券