【python】随机产生20个点用梯度下降法线性回归拟合---BGD+SGD+MBGD

最新推荐文章于 2024-07-27 12:20:46 发布

空空7

最新推荐文章于 2024-07-27 12:20:46 发布

阅读量1k

点赞数 2

分类专栏： Python 文章标签： python 机器学习逻辑回归

本文链接：https://blog.csdn.net/qq_45726331/article/details/115423288

版权

Python 专栏收录该内容

9 篇文章 3 订阅

订阅专栏

怎么拟合一道直线呢？先把直线方程设出来，h为预测函数

现在需要求解最佳的 $\theta_0$ 和 $\theta_1$ 得到最佳的直线
什么是最佳直线，就是用预测函数 $h$ 算出的预测值与实际值的差之和最小（只考虑y值）， $J$ 为代价函数

m是数据集中数据点的个数，也就是样本数；
½是一个常量，这样是为了在求梯度的时候，二次方乘下来的2就和这里的½抵消了，自然就没有多余的常数系数，方便后续的计算，同时对结果不会有影响。
那么现在就是求上述代价函数 $J$ （二元函数）的最小值，怎么求呢？用梯度下降法求，求出最佳的 $\theta_0$ 和 $\theta_1$ 使得代价函数 $J$ 取得最小值。按照如下的迭代公式求解（注意 $\theta$ 的上标表示先后两次迭代）

$\alpha$ 为学习率or步长，定义每一步的长度，不能太短也不能太长
为了简化计算，我们将 $J$ 和 $\nabla J$ 用矩阵形式表示，如下

如何随机产生点

# 20个点
m = 20
# 构造矩阵
x0 = np.ones((m,1))
x1 = np.arange(1, m+1).reshape(m, 1)
x = np.hstack((x0, x1)) # 得到x的值
y = np.random.randint(0, 30, size=(m, 1))
y.sort(axis=0) # 得到随机产生y的值，排序是为了使其类似直线分布

效果如图：
在这里插入图片描述
代价函数

# 求损失值函数
def loss_function(p, x, y):
    temp = np.dot(x, p) - y
    result = 1 / (2*m) * np.dot(np.transpose(temp), temp)
    return result

计算梯度函数

# 求梯度函数
def grad_function(p, x, y):
    temp = np.dot(x, p) - y
    result = 1 / m * np.dot(np.transpose(x),temp)
    return result

梯度迭代过程

# 存储迭代次数
count = 0 
# 设置p的初始值
p = np.array([1, 1]).reshape(2, 1)
# 计算初始梯度
grad = grad_function(p, x, y)
# 计算初始损失值
loss = loss_function(p, x, y)
# 存储损失值
lloss = []
loss = loss_function(p,x,y).flatten().tolist() # 先扁平化得到一维矩阵再转换为列表
lloss.extend(loss) # 这里用extend不用append
count += 1
# 设定精度
e = 1e-5
# 开始循环迭代
while np.all(np.absolute(grad) > e):
    p = p - algha * grad # 迭代公式
    # 更新grad和loss
    grad = grad_function(p, x, y)
    loss = loss_function(p, x, y)
    loss = loss_function(p,x,y).flatten().tolist()
    lloss.extend(loss)
    count += 1

完整代码如下：

import numpy as np 
import matplotlib.pyplot as plt
from numpy.core.numeric import ones 
# 显示中文
plt.rcParams['font.sans-serif']=['SimHei']
plt.rcParams['axes.unicode_minus'] = False
# 20个点
m = 20
# 构造矩阵
x0 = np.ones((m,1))
x1 = np.arange(1, m+1).reshape(m, 1)
x = np.hstack((x0, x1)) # 得到x的值
y = np.random.randint(0, 30, size=(m, 1))
y.sort(axis=0) # 得到随机产生y的值，排序是为了使其类似直线分布
# 设定学习率
algha = 0.01
# 求损失值函数
def loss_function(p, x, y):
    temp = np.dot(x, p) - y
    result = 1 / (2*m) * np.dot(np.transpose(temp), temp)
    return result
# 求梯度函数
def grad_function(p, x, y):
    temp = np.dot(x, p) - y
    result = 1 / m * np.dot(np.transpose(x),temp)
    return result

# 存储迭代次数
count = 0 
# 设置p的初始值
p = np.array([1, 1]).reshape(2, 1)
# 计算初始梯度
grad = grad_function(p, x, y)
# 计算初始损失值
loss = loss_function(p, x, y)
# 存储损失值
lloss = []
loss = loss_function(p,x,y).flatten().tolist() # 先扁平化得到一维矩阵再转换为列表
lloss.extend(loss) # 这里用extend不用append
count += 1
# 设定精度
e = 1e-5
# 开始循环迭代
while np.all(np.absolute(grad) > e):
    p = p - algha * grad # 迭代公式
    # 更新grad和loss
    grad = grad_function(p, x, y)
    loss = loss_function(p, x, y)
    loss = loss_function(p,x,y).flatten().tolist()
    lloss.extend(loss)
    count += 1
# 绘图
fig,ax = plt.subplots(1,2) # 一行两列的子图
fig.suptitle('基于梯度下降算法的线性回归拟合')
fig.dpi = 150
# 得到最佳的p
[p0, p1] = p
print('p:', p0,p1)
# 绘制拟合图像
ax[0].set_title('拟合图像')
ax[0].scatter(x1,y)
ax[0].plot(x1, p0 + p1*x1, color='r')
ax[0].grid(True) # 绘制网格线
# 绘制损失值图像
ax[1].set_title('损失值图像')
li = list(range(count))
print('迭代次数：', count)
ax[1].plot(li,lloss)
ax[1].grid(True)
plt.show()

效果如下：
在这里插入图片描述
上述代码中numpy库函数调用的一些解释
Python numpy.transpose 详解
 08_Numpy初始化生成相同元素值的ndarray数组
 np.dot()函数的用法详解

3. 其他改进形式梯度下降法（SGD+MBGD）

用关键段代码进行一个比较
BGD的非矩阵形式
每次使用全量的样本数据进行计算，得到m次梯度的和，用这个和乘以步长进行迭代
count+1是计算一次m个梯度的和，用和乘以步长更新一次 $\theta$
![在这里插入图片描述](https://img-blog.csdnimg.cn/20210404205201136.png

代码关键段：（完整附在后面）

while count < loop_max:
    count += 1
    # 在标准梯度下降中，权值更新的每一步对多个样例求和，需要更多的计算
    sum_m = np.zeros(2)
    for i in range(m):
        diff = (np.dot(p, x[i]) - y[i]) * x[i] # 一维向量进行内积，就是对应元素相乘再相加，二维及以上就是矩阵乘法运算
        sum_m = sum_m + diff     # 当alpha取值过大时,sum_m会在迭代过程中会溢出

    p = p - alpha * sum_m       # 注意步长alpha的取值,过大会导致振荡
    #w = w - 0.005 * sum_m      # alpha取0.005时产生振荡,需要将alpha调小

    # 判断是否已收敛，前后两次w的值的差小于某个阈值
    if np.linalg.norm(p - error) < epsilon:
        break
    else:
        error = p

SGD
每次从训练集中随机选择一个样本计算得到梯度，用该梯度乘以步长进行迭代
count+1是计算一个数据的梯度，用梯度乘以步长更新一次 $\theta$ ，这个过程进行m次，即总共更新m次 $\theta$
在这里插入图片描述

while count < loop_max:
    count += 1
    # 遍历训练数据集，不断更新p
    for i in range(m):
        # 训练集代入,计算误差值
        diff = (np.dot(p, x[i]) - y[i]) * x[i]  
        # 采用随机梯度下降算法,更新一次权值只使用一个训练数据
        p = p - alpha * diff
    # 终止条件：前后两次计算出的权向量的绝对误差充分小
    if np.linalg.norm(p - error) < epsilon:     
        break
    else:
        error = p

MBGD
综合了BGD与SGD，在每次更新速度与更新次数中间取得一个平衡，其每次更新从训练集中随机选择b,b<m个样本计算梯度和，梯度和乘以步长进行迭代。
count+1是计算b个数据的梯度和，用1/b乘以梯度和乘以步长更新 $\theta$ 。这个过程进行m/b次，即更新m/b次 $\theta$
在这里插入图片描述

while count < loop_max:
    count += 1
    for i in range(1,m,minibatch_size):
        sum_m = np.zeros(2)
        for k in range(i-1,i+minibatch_size-1,1):
            diff = (np.dot(p, x[k]) - y[k]) *x[k]
            sum_m = sum_m + diff  #当alpha取值过大时,sum_m会在迭代过程中会溢出
        p = p- alpha * (1.0/minibatch_size) * sum_m #注意步长alpha的取值,过大会导致振荡

    # 判断是否已收敛
    if np.linalg.norm(p- error) < epsilon:
        break
    else:
        error = p

完整BGD代码：

# 非矩阵形式实现
import numpy as np
from scipy import stats
import matplotlib.pyplot as plt

# 构造训练数据
x1 = np.arange(0., 10., 0.2)
m = len(x1)                                      # 训练数据点数目
x0 = np.full(m, 1.0)
x = np.vstack([x0, x1]).T               # 将偏置b作为权向量的第一个分量
y = 2 * x1 + 5 + np.random.randn(m)

# 两种终止条件
loop_max = 100000   # 最大迭代次数(防止死循环)
epsilon = 1e-5

# 初始化权值
np.random.seed(0)
p = np.random.randn(2)

alpha = 0.001      # 步长(注意取值过大会导致振荡,过小收敛速度变慢)
diff = 0.
error = np.zeros(2)
count = 0          # 循环次数

while count < loop_max:
    count += 1
    # 在标准梯度下降中，权值更新的每一步对多个样例求和，需要更多的计算
    sum_m = np.zeros(2)
    for i in range(m):
        diff = (np.dot(p, x[i]) - y[i]) * x[i] # 一维向量进行内积，就是对应元素相乘再相加，二维及以上就是矩阵乘法运算
        sum_m = sum_m + diff     # 当alpha取值过大时,sum_m会在迭代过程中会溢出

    p = p - alpha * sum_m       # 注意步长alpha的取值,过大会导致振荡
    #w = w - 0.005 * sum_m      # alpha取0.005时产生振荡,需要将alpha调小

    # 判断是否已收敛，前后两次w的值的差小于某个阈值
    if np.linalg.norm(p - error) < epsilon:
        break
    else:
        error = p
print ('loop count = %d' % count,  '\tw:[%f, %f]' % (p[0], p[1]))

plt.scatter(x1, y)
plt.plot(x1, p[1] * x1 + p[0], 'r')
plt.show()

完整SGD代码：

import numpy as np
import matplotlib.pyplot as plt

# 构造训练数据(另一种构造方法)
# x1 = np.arange(0.,10.,0.2)
# m = len(x1)  # 训练数据点数目
# x0 = np.full(m, 1.0)
# x = np.vstack([x0, x1]).T  # 将偏置b作为权向量的第一个分量
# y = 2 *x + 5 +np.random.randn(m) # 基于函数y=2x+5

m = 20
x0 = np.full(m, 1.0)
x1 = np.arange(1, m+1)  
x = np.vstack([x0, x1]).T  
y = np.random.randint(0, 30,size=20)
y.sort(axis=0)

# 两种终止条件
loop_max = 10000   # 最大迭代次数(防止死循环)
epsilon = 1e-5

# 设置p的初始值
np.random.seed(0)
p = np.random.randn(2)
# 设定学习率(注意取值过大会导致振荡,过小收敛速度变慢)
alpha = 0.001   
# 精度   
error = np.zeros(2)
# 循环次数
count = 0         
# 开始迭代
while count < loop_max:
    count += 1
    # 遍历训练数据集，不断更新p
    for i in range(m):
        # 训练集代入,计算误差值
        diff = (np.dot(p, x[i]) - y[i]) * x[i]  
        # 采用随机梯度下降算法,更新一次权值只使用一个训练数据
        p = p - alpha * diff
    # 终止条件：前后两次计算出的权向量的绝对误差充分小
    if np.linalg.norm(p - error) < epsilon:     
        break
    else:
        error = p
[p0, p1] = p
print('p:', p0,p1)
print('迭代次数：', count)
# 绘图
plt.scatter(x1, y)
plt.plot(x1,  p0 + p1 * x1, 'r')
plt.show()

完整MBGD代码

import numpy as np
import matplotlib.pyplot as plt

# 构造训练数据
# x1 = np.arange(0.,10.,0.2)
# m = len(x1)  # 训练数据点数目
# x0 = np.full(m, 1.0)
# x = np.vstack([x0, x1]).T  # 将偏置b作为权向量的第一个分量
# y = 2 * x1 + 5 +np.random.randn(m)

m = 20
x0 = np.full(m, 1)
x1 = np.arange(1, m+1)  
x = np.vstack([x0, x1]).T  
y = np.random.randint(0, 30,size=20)
y.sort(axis=0)

# 两种终止条件
loop_max = 10000  #最大迭代次数(防止死循环)
epsilon = 1e-5

# 设置p的初始值
np.random.seed(0)
p = np.random.randn(2)
#步长(注意取值过大会导致振荡即不收敛,过小收敛速度变慢)
alpha = 0.001  
diff = 0.
error = np.zeros(2)
# 循环次数
count = 0  
# 每次更新的样本数
minibatch_size = 5  
while count < loop_max:
    count += 1
    # 在minibatch梯度下降中，权值更新的每一步对多个样例求和，需要更多的计算
    for i in range(1,m,minibatch_size):
        sum_m = np.zeros(2)
        for k in range(i-1,i+minibatch_size-1,1):
            diff = (np.dot(p, x[k]) - y[k]) *x[k]
            sum_m = sum_m + diff  #当alpha取值过大时,sum_m会在迭代过程中会溢出
        p = p- alpha * (1.0/minibatch_size) * sum_m #注意步长alpha的取值,过大会导致振荡

    # 判断是否已收敛
    if np.linalg.norm(p- error) < epsilon:
        break
    else:
        error = p

[p0, p1] = p
print('p:', p0,p1)
print('迭代次数：', count)
# 绘制
plt.scatter(x1, y)
plt.plot(x1, p[1]* x1 + p[0],'r')
plt.show()

参考梯度下降算法及python实现（学习笔记）

空空7

关注

2
点赞
踩
11

收藏

觉得还不错? 一键收藏
2
评论
【python】随机产生20个点用梯度下降法线性回归拟合---BGD+SGD+MBGD

随机产生20个点假设一条直线，符合预期，现需求解最佳的p0和p1（待定系数）什么是最佳直线，就是用拟合直线算出的预测值与实际值的差之和最小（只考虑y）—目标函数那么现在就是求上述二元函数的最小值，怎么求呢，用梯度下降法求，得到最佳的p0和p1使得目标函数值最小梯度下降法的迭代公式用矩阵的形式表示...
复制链接

扫一扫