五、【机器学习作业】反向传播的神经网络（python版ex4）

最新推荐文章于 2022-08-30 16:33:42 发布

Liaojiajia-2020

最新推荐文章于 2022-08-30 16:33:42 发布

阅读量745

点赞数 2

分类专栏： # 机器学习实验文章标签：机器学习

本文链接：https://blog.csdn.net/mary_0830/article/details/98482347

版权

机器学习实验专栏收录该内容

11 篇文章 8 订阅

订阅专栏

反向传播的神经网络（ex4）

（一）神经网路 Neural Networks
（二）反向传播 Backpropagation
（三）可视化隐藏层 Visualizing the hidden layer

在这篇博客中，使用的例子依然是手写字符的识别。

（一）神经网路 Neural Networks

（1）可视化数据集 Visualizing the data

这里，将在训练集的训练样本中，随机选取100个数据并进行可视化。该训练集共有5000个样本，每个样本是20×20像素的灰度图像，每个像素代表一个浮点数，表示这个位置的灰度强度。在该数据集中，每一个样本展成一行，每一行有400列，共转化成5000×400 的矩阵。
首先，导入所需要用到的库：

import numpy as np
import matplotlib.pyplot as plt
import scipy.optimize as opt
from scipy.io import loadmat
from sklearn.metrics import classification_report #用于评价报告

载入数据集，以及可视化随机100个数据：

def load_mat(path):
    '''读取.mat数据'''
    data = loadmat('ex4data1.mat')
    X = data['X']
    y = data['y'].flatten()
    
    return X,y

def plot_images(X):
    '''随机画出数据集里的100个数字'''
    index = np.random.choice(range(5000),100)
    images = X[index]
    fig, ax = plt.subplots(10, 10, sharex =True, sharey = True,figsize=(6,6)) 
    #sharex和sharey为True时，是不显示坐标轴；若为False,则显示坐标轴
    for r in range(10):
        for c in range(10):
            ax[r,c].matshow(images[r*10 + c].reshape(20,20).T,cmap = 'gray_r')
            #不加.T的话，数字不是正的
    plt.xticks([])
    plt.yticks([])
    plt.show()

X,y = load_mat('ex4data1.mat')
plot_images(X)

其中，fig, ax = plt.subplots(10, 10, sharex =True, sharey = True,figsize=(6,6)) 代表的是10行10列6×6大小的子图，当sharex，sharey = True时，图像中不显示坐标轴，反之，则显示坐标轴。index = np.random.choice(range(5000),100)表示的是在5000个数据中随机抽取100个元素。

fig = plt.figure()
ax = fig.add_subplot(1,1,1)
例如：fig, ax = plt.subplots(1,3,figsize=(15,7))表示为1行3列15x7大小的子图。

实验结果为：
在这里插入图片描述

（2）模型表示 Model Representation

这里有三层，分别为输入层、隐藏层和输出层。由于图像的大小为20×20像素，这就有了400个输入层单元（不包括输出层要加一个偏置单元）
在这里插入图片描述
1）. 读取数据
首先，需要将输出0-9转换为非线性相关的向量的形式：

例如， $y^{(3)}$ =2转化成 $y^{(3)}$ =[0,0,1,0,0,0,0,0,0,0].

解释： 由于神经网络的训练是监督学习。也就是说，样本训练数据集是这样的格式： $x^{(i)}, y^{(i)})$ ，对于一个训练实例 $x^{(i)}$ ，我们是已经确定知道了它的正确结果是 $y^{(i)}$ ，而我们的目标是构造一个神经网络模型，训练出来的这个模型的假设函数 $h_{θ}(x)$ ，对于未知的输入数据 $x^{(k)}$ ，能够准确地识别出正确的结果。
因此，训练数据集(traing set)中的结果数据 $y$ 是正确的已知的结果，比如 $y^{(3)}=[0,0,1,0,0,0,0,0,0,0]^{T}$ 表示：训练数据集中的第3条训练实例它所对应的正确结果是：数字2 （因为，向量 $y^{(3)}$ 中的第2个元素为1，其它所有元素为0）。

另外需要注意的是：当向量 $y^{(i)}$ 第10个元素为1时代表数字0。

在这里插入图片描述
写到这里，就会想知道为什么要将原输出（0-9）映射为输出由0与1组成的向量形式呢？

解释：
（1）对于深度神经网络，中间的隐藏层的输出必须有一个激活函数，否则多个隐藏层的作用和没有隐藏层相同。这个激活函数不一定是sigmoid，常见的有sigmoid、tanh、relu等。
（2）对于二分类问题，输出层是sigmoid函数。这是因为sigmoid函数可以把实数域光滑的映射到[0,1]空间。函数值恰好可以解释为属于正类的概率（概率的取值范围是0~1）。另外，sigmoid函数单调递增，连续可导，导数形式非常简单，是一个比较合适的函数。
（3）对于多分类问题，输出层就必须是softmax函数了。softmax函数是sigmoid函数的推广。
（4）需要转换成计算机容易识别的0\1字符。

获取数据集，并对数据进行相应的处理：

def expand_y(y):
    result = [] #用于储存转变后的结果
    for i in y: #把y进行向量化操作，对应的lable值将在对应的位置上置1
        y_array = np.zeros(10)
        y_array[i-1] = 1
        result.append(y_array)
    
    return np.array(result)

raw_X,raw_y = load_mat('ex4data1.mat')
X = np.insert(raw_X,0,1,axis=1) #添加偏置单元
y = expand_y(raw_y)
print('X的维数:',X.shape)
print('y的维数:',y.shape)

得到处理后的数据，包括X的维数与y的维数：

X的维数: (5000, 401)
y的维数: (5000, 10)

2）. 读取权重
读取提供在ex4weight.mat里的参数θ1、θ2，文件包括第二层（隐藏层）中的25个单元，第三层（输出层）中的10个单元。

def load_weight(path):
    data = loadmat(path)
    
    return data['Theta1'],data['Theta2']

t1,t2 = load_weight('ex4weights.mat')
print('t1的维数:',t1.shape)
print('t2的维数:',t2.shape)

读取结果，包含有Theta1和Theta2的维数：

t1的维数: (25, 401)
t2的维数: (10, 26)

3）. 计算参数
使用高级优化算法对神经网络进行优化时，需要将多个参数矩阵展开，分别传入优化函数，之后再恢复矩阵的形状。

解释：
高级优化算法的参数 theta 和 initialTheta 都是 $n + 1$ 维的列向量，而神经网络中的参数是不同规格的矩阵。要想利用高级优化算法，需要将所有的参数矩阵转化为一个列向量才能进行函数调用（使所有参数平等化），算法执行完毕后，再将所得的参数列向量转化为之前的矩阵格式。

def serialize(a,b):
    '''展开参数'''
    #返回一个折叠成一维的数组；按列连接，上下相加，要求列数相等
    return np.r_[a.flatten(),b.flatten()] 
def deserialize(seq):
    '''提取参数'''
    return seq[:25*401].reshape(25,401),seq[25*401:].reshape(10,26)

theta = serialize(t1,t2) #使参数扁平化处理，数值为25*401+10*26=10285
print('theta的维数：',theta.shape)

将矩阵展开后的维数为：

theta的维数： (10285,)

（3）前馈与代价函数 Feedforward and cost function

1）.前馈
除了输出层，在之前的每层都需要添加一个偏置单元，即输入层s(1)=400+1，隐藏层s(2)=25+1，输出层不需要加上偏置单元。

def sigmoid(h):
    return 1 / (1 + np.exp(-h))

def feed_forward(theta,X):
    '''计算得到每一层的输入与输出'''
    t1,t2 = deserialize(theta)
    a1 = X #由于在之前已经添加了a1层的偏置单元
    z2 = np.dot(a1,t1.T)
    a2 = np.insert(sigmoid(z2),0,1,axis=1)
    z3 = np.dot(a2,t2.T)
    a3 = sigmoid(z3)
    
    return a1,z2,a2,z3,a3

a1,z2,a2,z3,h = feed_forward(theta,X)
print('a1:',a1)
print('z2:',z2)
print('a2:',a2)
print('z2:',z2)
print('h:',h)

print('a1:', a1.shape,'t1:', t1.shape)
print('z2:', z2.shape)
print('a2:', a2.shape, 't2:', t2.shape)
print('z3:', z3.shape)
print('a3:', h.shape)

计算结果为： 获取a1,z2,a2,z3,a3的值。

a1: [[1. 0. 0. … 0. 0. 0.]
[1. 0. 0. … 0. 0. 0.]
[1. 0. 0. … 0. 0. 0.]
…
[1. 0. 0. … 0. 0. 0.]
[1. 0. 0. … 0. 0. 0.]
[1. 0. 0. … 0. 0. 0.]]
z2: [[ -2.93684669 -2.45058587 4.95510333 … 3.56635593 2.81388641
-2.1195223 ]
[ -4.81302157 -2.92257775 2.6445065 … 2.10497303 4.69948787
-2.76096862]
[ -4.24056958 -3.68698052 5.99656398 … 1.54599347 3.08971226
-2.32990819]
…
[ -0.86267303 1.00939507 -1.67526051 … 1.8185898 -3.18203449
-1.72539781]
[ 1.74408423 -0.58216518 -1.49164167 … 4.17481481 -0.96739536
-3.08906563]
[ 3.55683614 -12.11330792 5.01096205 … 7.17585008 2.15484114
-2.9424052 ]]
a2: [[1.00000000e+00 5.03618685e-02 7.93957162e-02 … 9.72517962e-01
9.43421623e-01 1.07213787e-01]
[1.00000000e+00 8.05782163e-03 5.10486829e-02 … 8.91385592e-01
9.90982126e-01 5.94701645e-02]
[1.00000000e+00 1.41949887e-02 2.44354705e-02 … 8.24334311e-01
9.56466386e-01 8.86760824e-02]
…
[1.00000000e+00 2.96781175e-01 7.32901746e-01 … 8.60396828e-01
3.98474223e-02 1.51177198e-01]
[1.00000000e+00 8.51205095e-01 3.58434539e-01 … 9.84854863e-01
2.75399966e-01 4.35605471e-02]
[1.00000000e+00 9.72262381e-01 5.48598771e-06 … 9.99235749e-01
8.96120297e-01 5.00966928e-02]]
z2: [[ -2.93684669 -2.45058587 4.95510333 … 3.56635593 2.81388641
-2.1195223 ]
[ -4.81302157 -2.92257775 2.6445065 … 2.10497303 4.69948787
-2.76096862]
[ -4.24056958 -3.68698052 5.99656398 … 1.54599347 3.08971226
-2.32990819]
…
[ -0.86267303 1.00939507 -1.67526051 … 1.8185898 -3.18203449
-1.72539781]
[ 1.74408423 -0.58216518 -1.49164167 … 4.17481481 -0.96739536
-3.08906563]
[ 3.55683614 -12.11330792 5.01096205 … 7.17585008 2.15484114
-2.9424052 ]]
h: [[1.12661530e-04 1.74127856e-03 2.52696959e-03 … 4.01468105e-04
6.48072305e-03 9.95734012e-01]
[4.79026796e-04 2.41495958e-03 3.44755685e-03 … 2.39107046e-03
1.97025086e-03 9.95696931e-01]
[8.85702310e-05 3.24266731e-03 2.55419797e-02 … 6.22892325e-02
5.49803551e-03 9.28008397e-01]
…
[5.17641791e-02 3.81715020e-03 2.96297510e-02 … 2.15667361e-03
6.49826950e-01 2.42384687e-05]
[8.30631310e-04 6.22003774e-04 3.14518512e-04 … 1.19366192e-02
9.71410499e-01 2.06173648e-04]
[4.81465717e-05 4.58821829e-04 2.15146201e-05 … 5.73434571e-03
6.96288990e-01 8.18576980e-02]]

获取a1,z2,a2,z3,a3的维度：a1，t1，a2，t2 都加了偏置单元。

a1: (5000, 401) t1: (25, 401)
z2: (5000, 25)
a2: (5000, 26) t2: (10, 26)
z3: (5000, 10)
a3: (5000, 10)

2）.代价函数

未正则化的数学公式为：
$J(\theta )=\frac{1}{m}\sum_{i=1}^{m}\sum_{k=1}^{K}[-y^{(i)}_{k}log((h_{\theta }(x^{(i)}))_{k})-(1-y^{(i)}_{k})log(1-(h_{\theta }(x^{(i)}))_{k})]$

提问：代价函数为什么要取负数？
解释：在逻辑回归算法中，使用最大似然函数求解所需要的参数 $\theta$ 值，得到的结果是最大值。但是，我们想要的代价函数是最小值（代价越小，说明算法的能力越好），因此需要将代价函数在最大似然函数的结果中取负号（即可以通过计算得到最小值）。

参考博客中的公式推导 https://blog.csdn.net/mary_0830/article/details/97612341

def cost_function(theta,X,y):
    a1,z2,a2,z3,h = feed_forward(theta,X)
    J = 0
    for i in range(len(X)):
        first = -y[i] * np.log(h[i])
        second = (1-y[i]) * np.log(1-h[i])
        J = J + np.sum(first - second)
    J = J / len(X)
    
    return J

print('cost_function:',cost_function(theta,X,y))

计算得到的代价函数值为：

cost_function: 0.2876291651613187

（4）正则化的代价函数 Regularized cost function

计算公式为：
$J(\theta )=\frac{1}{m}\sum_{i=1}^{m}\sum_{k=1}^{K}[-y^{(i)}_{k}log((h_{\theta }(x^{(i)}))_{k})-(1-y^{(i)}_{k})log(1-(h_{\theta }(x^{(i)}))_{k})]\\+\frac{\lambda }{2m}[\sum_{j=1}^{25}\sum_{k=1}^{400}(\Theta _{j,k}^{(1)})^{2}+\sum_{j=1}^{10} \sum_{k=1}^{25}(\Theta _{j,k}^{(2)})^{2}]$
需要注意的是，不能够将偏置单元正则化。

def regularized_cost_function(theta,X,y,l=1):
    '''正则化时需要忽略每一层的偏置单元'''
    t1,t2 = deserialize(theta)
    reg = np.sum(t1[:,1:] ** 2) + np.sum(t2[:,1:] ** 2)
    
    return l / (2*len(X)) * reg + cost_function(theta,X,y)

print('regularized_cost_function:',regularized_cost_function(theta,X,y,1))

在reg = np.sum(t1[:,1:] ** 2) + np.sum(t2[:,1:] ** 2)中，正则化操作从第一个元素开始，把偏置单元剔除了。
正则化后的代价函数结果为：

regularized_cost_function: 0.3837698590909234

（二）反向传播 Backpropagation

其原理推导请参考 https://blog.csdn.net/mary_0830/article/details/99182023 中的反向传播算法。

（1）S函数导数 Sigmoid gradient

数学表达式为：
$g(z)=\frac{1}{1+e^{-z}}$ $g(z)^{'}=g(z)(1-g(z))$
编写代码：

def sigmoid_gradient(z):
    return sigmoid(z) * (1 - sigmoid(z))

（2）随机初始化 Random initialization

当训练神经网络时，随机初始化参数是非常重要的一步，用于打破数据的对称性。在均匀分布 $(- e, e)$ 中随机选择值，在这里选择 $e = 0.12$ 这个范围的值来确保参数足够小，使得神经网络训练更有效率。

def random_initialization(size):
    '''随机初始化，从均匀分布的范围中随机返回size大小的值'''
    #np.radom.uniform()用于随机取数，并且取得数值是浮点型
    return np.random.uniform(-0.12,0.12,size)

（3）反向传播 Backpropagation

在这里插入图片描述
获取整个神经网络的代价函数的梯度，以便在优化算法中求解。
参数说明：

计算每层的"误差"：
$\delta ^{(3)}=h-y$
$\delta ^{(2)}=(\theta ^{(2)})^{T}\cdot \delta ^{(3)}\cdot g^{'}(z^{(2)})$
每层参数矩阵的梯度，用Δ^(l)表示：
$\triangle ^{(2)}=a^{(2)}\cdot \delta ^{(3)}$
$\triangle ^{(1)}=a^{(1)}\cdot \delta ^{(2)}$
神经网络的总梯度为：
$\frac{1}{m}(\triangle ^{(1)}+\triangle ^{(2)})$

def gradient(theta,X,y):
    '''计算梯度。梯度D(i)与参数theta(i)的shape需要一致'''
    t1,t2 = deserialize(theta)
    a1,z2,a2,z3,h = feed_forward(theta,X)
    d3 = h - y
    d2 = np.dot(d3,t2[:,1:].T) * sigmoid_gradient(z2)
    D2 = np.dot(d3.T,a2)
    D1 = np.dot(d2.T,a1)
    D = np.dot((1 / len(X)),serialize(D1,D2))
    
    return D

（4）梯度检测 Gradient checking

为了知道反向传播计算出来的梯度是否正确，需要用数值方法计算梯度，与反向传播计算得到的结果进行比较。

通过下面的计算方法可以进行梯度检测：
在这里插入图片描述
梯度检测的代码为：

def gradient_checking(theta,X,y,e):
    def numeric_grad_(plus,minus):
        '''对每个参数theta[i]计算数值梯度'''
        return (regularized_cost_function(plus,X,y) - regularized_cost_function(minus,X,y)) / (e*2)
    
    numeric_grad = []
    for i in range(len(theta)):
        plus = theta.copy()
        minus = theta.copy()
        plus[i] = plus[i] + e
        minus[i] = minus[i] -e
        grad_i = numeric_grad_(plus,minus)
        numeric_grad.append(grad_i)
        
    numeric_grad = np.array(numeric_grad)
    analytic_grad = regularized_gradient(theta,X,y)
    diff = np.linalg.norm(numeric_grad - analytic_grad) / np.linalg.norm(numeric_grad + analytic_grad)
    
    print('If your backpropagation implementation is correct,\nthe relative difference will be smaller than 10e-9 (assume epsilon=0.0001).\nRelative Difference: {}\n'.format(diff))
    
gradient_checking(theta, X, y, e= 0.0001)

PS：上面的代码谨慎运行，实在是太慢了。
实验结果：

If your backpropagation implementation is correct,
the relative difference will be smaller than 10e-9 (assume epsilon=0.0001).
Relative Difference: 1.0921821354694573e-08

说明： numeric_grad.append(grad_i)表示在numeric_grad后面加上grad_i的内容。np.linalg.norm表示的是：在numpy库中，linalg=linear（线性）+algebra（代数），norm则表示范数。

例子：x_norm=np.linalg.norm(x, ord=None, axis=None, keepdims=False)
x: 表示矩阵（也可以是一维）
ord： 范数类型

axis： 处理类型

keepding：是否保持矩阵的二维特性

（5）正则化神经网络 Regularized Neural Networks

当成功使用反向传播后，可以向梯度 $\triangle ^{(l)}_{ij}$ 增加正则化项，如下面的数学公式：

当 $j = 0$ 时， $\frac{\partial }{\partial \Theta _{ij}^{(l)}}J(\Theta )=D_{ij}^{(l)}=\frac{1}{m}\triangle ^{(l)}_{ij}$
当 $j\geq 1$ 时， $\frac{\partial }{\partial \Theta _{ij}^{(l)}}J(\Theta )=D_{ij}^{(l)}=\frac{1}{m}\triangle ^{(l)}_{ij} +\frac{\lambda }{m}\Theta ^{(l)}_{ij}$

def regularized_gradient(theta,X,y,l=1):
    '''正则化需要忽略每一层的偏置单元，不惩罚偏置单元的参数'''
    a1,z2,a2,z3,h = feed_forward(theta,X)
    D1,D2 = deserialize(gradient(theta,X,y))
    t1[:,0] = 0
    t2[:,0] = 0
    reg_D1 = D1 + (l / len(X)) * t1
    reg_D2 = D2 + (l / len(X)) * t2
    
    return serialize(reg_D1,reg_D2)

（6）优化参数 Learning parameters using fmincg

已经完成神经网络代价函数和梯度的计算后，可以使用minimize来学习一个很好的集合参数。在训练完成后，通过计算得到的百分比来报告分类器的训练准确性。当你的实现是正确时，可以看到训练准确度为95.3%（可能由于随机初始化而变化1%左右）。经过训练神经网络更多的迭代次数，可以获得更加高的准确度。尝试对训练神经网络进行更多次的迭代，并且改变正则化的参数值。

训练神经网络的代码：

def nn_training(X,y):
    '''训练神经网络'''
    init_theta = random_initialization(10285)
    res = opt.minimize(fun = regularized_cost_function,x0 = init_theta,
                       args = (X,y,1),
                       method = 'TNC',
                       jac = regularized_gradient,
                       options={'maxiter':400})
    return res

res = nn_training(X,y)
print('res:',res)

运行结果：

res: fun: 0.5144107017503667
jac: array([ 6.94146365e-04, -2.11248326e-12, 4.38829369e-13, …,
-7.70180165e-04, -2.47067656e-04, -1.29165183e-03])
message: ‘Converged (|f_n-f_(n-1)| ~= 0)’
nfev: 250
nit: 15
status: 1
success: True
x: array([ 1.12554435e-01, 1.07194714e-03, -9.58223897e-02, …,
-2.18778746e+00, -1.12047481e+00, 5.59999713e+00])

说明： opt.minimize()实际上包含了两个步骤，即 compute_gradients 和 apply_gradients，前者用于计算梯度，后者用于使用计算得到的梯度来更新对应的variable。

官方例子： scipy.optimize.minimize(fun, x0, args=(), method=None, jac=None,hess=None, hessp=None, bounds=None, constraints=(), tol=None,callback=None, options=None)

例子解释：
fun：求最小值的目标函数。
x0：变量的初始猜测值，如果有多个变量，需要给每个变量一个初始猜测值。minimize是局部最优的解法。
args：常数值，fun中没有数字，都以变量的形式表示，对于常数项，需要在这里给值。
method：求极值的方法，官方文档给了很多种。一般使用默认。每种方法我理解是计算误差，反向传播的方式不同而已，这块有很大理论研究空间。
constraints：约束条件，针对fun中为参数的部分进行约束限制。

计算算法准确度的代码：

def accuracy(theta,X,y):
    a1,z2,a2,z3,a3 = feed_forward(theta,X)
    y_pred = np.argmax(h,axis = 1) + 1  #axis=1：按行查找最大元素 axis=0：按列查找最大元素
    print(classification_report(y,y_pred))

accuracy(theta,X,raw_y)

实验结果：

              precision    recall  f1-score   support

           1       0.97      0.98      0.98       500
           2       0.98      0.97      0.98       500
           3       0.98      0.96      0.97       500
           4       0.97      0.97      0.97       500
           5       0.97      0.98      0.98       500
           6       0.98      0.99      0.98       500
           7       0.98      0.97      0.97       500
           8       0.98      0.98      0.98       500
           9       0.97      0.96      0.96       500
          10       0.98      0.99      0.99       500

    accuracy                           0.98      5000
   macro avg       0.98      0.98      0.98      5000
weighted avg       0.98      0.98      0.98      5000

精确度为98%，说明这个算法是蛮好的。

（三）可视化隐藏层 Visualizing the hidden layer

有一种了解你的神经网络的方法是 可视化隐藏层单元所捕获的内容。给定一个特定的隐藏单元，通过它计算找到一个输入为X，X可以激活这个单元（即有一个激活值 $a^{(l)}_i$ 接近于1）。对于所训练的神经网络，注意到 $\theta^{(1)}$ 中的每一行是一个401维度的向量，表示每个隐藏层单元的参数。当忽略偏置单元，就能得到400维度的向量，这个向量表示每个样本输入到每个隐藏层单元的像素的权重。所以，可视化隐藏层的方法是，reshape这个400维度的向量变成20×20的图像并显示它。它将显示一个图像有25个单元，且每个单元对应网络中的一个隐藏单元。

可视化隐藏层的代码：

def plot_hidden(theta):
    '''可视化隐藏层'''
    t1,t2 = deserialize(theta)
    t1 = t1[:,1:]
    fig,ax = plt.subplots(5,5,sharex = True,sharey = True,figsize = (6,6))
    for r in range(5):
        for c in range(5):
            ax[r,c].matshow(t1[r * 5 + c].reshape(20,20).T,cmap = 'gray_r')
            plt.xticks([])
            plt.yticks([])
    plt.show()

plot_hidden(theta)

实验结果：
在这里插入图片描述

Liaojiajia-2020

关注

2
点赞
踩
7

收藏

觉得还不错? 一键收藏
0
评论
五、【机器学习作业】反向传播的神经网络（python版ex4）

反向传播的神经网络（ex4）（一）神经网路 Neural Networks（1）可视化数据集 Visualizing the data（2）模型表示 Model Representation（1）读取数据（2）读取权重（3）计算参数（二）反向传播（三）可视化隐藏层在这篇博客中，使用的例子依然是手写字符的识别。（一）神经网路 Neural Networks（1）可视化数据集 Visualiz...
复制链接

扫一扫