一、神经网络
在前面的练习中，你实现了神经网络的前馈传播，并利用它预测手写数字的权重。在这个练习中，你将实现反向传播算法来学习神经网络的参数。
1.1 Visualizing the data(可视化数据)
这与您在上一次练习中使用的数据集相同。在ex3data1.mat中有5000个训练示例，其中每个训练示例是20像素乘20像素的数字灰度图像。像素用浮点数表示，表示该位置的灰度强度。20×20像素网格被“展开”成400维矢量。每一次训练数据矩阵X中的一行。
它给出了一个5000×400矩阵X，其中每一行都是手写数字映像的训练示例。训练集的第二部分是包含训练集标签的5000维向量y。“0”数字被标记为“10”，而数字“1”到“9”则被标记为“1”到“9”。’’’
import matplotlib.pyplot as plt
import numpy as np
import scipy.io as scio #引入读取.mat文件的库
import scipy.optimize as opt
#import displayData as dd
plt.ion()
input_layer_size=400 #输入特征的维度 每张数字图片20*20=400
hidden_layer_size=25 #隐含层的个数为25
num_labels=10 #10个标签 注意“0”对应第十个标签 1-9一次对应第1-9个标签
data=scio.loadmat('ex4data1.mat')#读取训练集 包括两部分 输入特征和标签
X=data['X']#提取输入特征 5000*400的矩阵 5000个训练样本 每个样本特征维度为400 一行代表一个训练样本
y=data['y'].flatten() #提取标签 data['y']是一个5000*1的2维数组 利用flatten()将其转换为有5000个元素的一维数组
y_raw=data['y'].flatten()
编写可视化数据集的程序：
def display_data(x):
(m, n) = x.shape #100*400
example_width = np.round(np.sqrt(n)).astype(int) #每个样本显示宽度 round()四舍五入到个位 并转换为int
example_height = (n / example_width).astype(int) #每个样本显示高度 并转换为int
#设置显示格式 100个样本 分10行 10列显示
display_rows = np.floor(np.sqrt(m)).astype(int)#floor(x)，其功能是“向下取整”，下取整是直接去掉小数部分。
display_cols = np.ceil(m / display_rows).astype(int)#ceil(x),其功能是向上取整。
# 待显示的每张图片之间的间隔
pad = 1
# 显示的布局矩阵 初始化值为-1
display_array = - np.ones((pad + display_rows * (example_height + pad),
pad + display_rows * (example_height + pad)))
# Copy each example into a patch on the display array
curr_ex = 0
for j in range(display_rows):
for i in range(display_cols):
if curr_ex > m:
break
# Copy the patch
# Get the max value of the patch
max_val = np.max(np.abs(x[curr_ex]))
display_array[pad + j * (example_height + pad) + np.arange(example_height),
pad + i * (example_width + pad) + np.arange(example_width)[:, np.newaxis]] = \
x[curr_ex].reshape((example_height, example_width)) / max_val
curr_ex += 1
if curr_ex > m:
break
# 显示图片
plt.figure()
plt.imshow(display_array, cmap='gray', extent=[-1, 1, -1, 1])
plt.axis('off')
加载数据在上面的display_data(x)函数上：
m = y.size #训练样本的数量
# 随机抽取100个训练样本 进行可视化
rand_indices = np.random.permutation(range(m)) #获取0-4999 5000个无序随机索引
selected = X[rand_indices[0:100], :] #获取前100个随机索引对应的整条数据的输入特征
display_data(selected) #调用可视化函数 进行可视化
运行结果：
1.2 load weight 读取权重
这里我们提供了已经训练好的参数θ1，θ2，存储在ex4weight.mat文件中。这些参数的维度由神经网络的大小决定，第二层有25个单元，输出层有10个单元(对应10个数字类)。
data=scio.loadmat('ex4weights.mat')
theta1=data['Theta1'] #theta1:(25,401)
theta2=data['Theta2'] #theta2:(10,26)
#当我们使用高级优化方法来优化神经网络时，
#我们需要将多个参数矩阵展开，才能传入优化函数，然后再恢复形状。
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(theta1, theta2) # 扁平化参数，25*401+10*26=10285
1.3 正向传播
1.3.1标签向量化
首先我们要将标签值（1，2，3，4，…，10）转化成非线性相关的向量，向量对应位置（y[i-1]）上的值等于1，例如y[0]=6转化为y[0]=[0,0,0,0,0,1,0,0,0,0]。
def expand_y(y):
result=[]
# 把y中每个类别转化为一个向量，对应的lable值在向量对应位置上置为1
for i in y:
y_array=np.zeros(10)
y_array[i-1]=1
result.append(y_array)
return np.array(result)
X=np.insert(X,0,1,axis=1) # 5000 x 401
y=expand_y(y) #5000 x 10
1.3.2 定义sigmoid()函数
def sigmoid(z):
return 1/(1+np.exp(-z))
1.3.2 前向传播
假设现在有一个三层的神经网络，如图：
def feed_forward(theta,X):
t1, t2 = deserialize(theta)
a1=X
[email protected] # (5000,25)
a2=np.insert(sigmoid(z2),0,1,axis=1) # (5000,26)
[email protected] # (5000,10)
a3=sigmoid(z3)
return a1,z2,a2,z3,a3
a1,z2,a2,z3,h=feed_forward(theta,X)
t1, t2 = deserialize(theta)
1.3.3 Compute Cost (Feedforward)计算成本
神经网络的成本函数(没有正则化)是：
def cost(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(theta,X, y))#0.287629165161
#注意不要将每层的偏置项正则化。
def regularized_cost(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(theta,X,y)
print(regularized_cost(theta,X,y,1))#0.383769859091
1.3.4 Regularized cost function（正规化成本函数）
给出了正则化神经网络的代价函数：
#注意不要将每层的偏置项正则化。
def regularized_cost(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(theta,X,y)
print(regularized_cost(theta,X,y,1))#0.383769859091
二、Backpropagation 反向传播
在本练习中，您将实现反向传播算法来计算神经网络成本函数的梯度。
一旦计算了梯度，就可以通过使用高级优化器(如minimize)最小化成本函数J(Θ)来训练神经网络。
首先实现反向传播算法来计算(非正则化)神经网络参数的梯度。然后验证对于非正则情形的梯度计算是正确的，最后实现正则化神经网络的梯度。
2.1 Sigmoid gradient（函数导数）
Sigmoid 函数为：
对Sigmoid 函数导数得：
def sigmoid_gradient(z):
return sigmoid(z)*(1-sigmoid(z))
2.2 Random initialization 随机初始化
当我们训练神经网络时，随机初始化参数是很重要的，可以打破数据的对称性。一个有效的策略是在均匀分布(−e，e)中随机选择值，我们可以选择 e = 0.12 这个范围的值来确保参数足够小，使得训练更有效率。
def random_init(size):
'''从服从的均匀分布的范围中随机返回size大小的值'''
return np.random.uniform(-0.12, 0.12, size)
2.3 Backpropagation 反向传播
目标：获取整个网络代价函数的梯度。以便在优化算法中求解。
这里面一定要理解正向传播和反向传播的过程，弄清楚各种参数在网络中的维度。
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 (5000, 401) t1 (25, 401)
z2 (5000, 25)
a2 (5000, 26) t2 (10, 26)
z3 (5000, 10)
a3 (5000, 10)'''
先来从后往前计算每层的“误差“。注意到这里的误差用双引号括起来，因为并不是真正的误差。
注意第一层是没有误差的，因为是输入层。
然后来计算每层参数矩阵的梯度，用Δ(l)表示：
最后网络的总梯度为：
到这里反向传播就完成了，接着就可以利用梯度下降法或者更高级的优化算法来训练网络。
def gradient(theta,X,y):
t1,t2=deserialize(theta)
a1,z2,a2,z3,h=feed_forward(theta,X)
d3=h-y #(5000,10)
[email protected][:,1:]*sigmoid_gradient(z2) #(5000,25)
[email protected] # (10,26)
[email protected] #(25,401)
D=(1/len(X))*serialize(D1,D2) #(10285,)
return D
2.4 Gradient checking 梯度检测
在你的神经网络,你是最小化代价函数J(Θ)。执行梯度检查你的参数,你可以想象展开参数Θ(1)Θ(2)成一个长向量θ。通过这样做,你能使用以下梯度检查过程。
def gradient_checking(theta,X,y,e):
def a_numeric_grad(plus,minus):
return (regularized_cost(plus,X,y,l=1)-regularized_cost(minus,X,y,l=1))/(2*e)#对每个参数theta_i#计算数值梯度，即理论梯度。
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=a_numeric_grad(plus,minus)
numeric_grad.append(grad_i)
numeric_grad=np.array(numeric_grad)
analytic_grad=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,0.0001)#这个运行很慢，谨慎运行
运行结果：
If your backpropagation implementation is correct,
the relative difference will be smaller than 10e-9 (assume epsilon=0.0001).
Relative Difference: 0.19495118386913435
注意：np.linalg.norm(求范数)
linalg=linear（线性）+algebra（代数），norm则表示范数。
函数参数：
x_norm=np.linalg.norm(x, ord=None, axis=None, keepdims=False)
①x: 表示矩阵（也可以是一维）
②ord：范数类型
2.5 Regularized Neural Networks 正则化神经网络
在成功地实现了反向传播算法之后，将向梯度添加正则化。
def regularized_gradient(theta,X,y,l=1):
"""不惩罚偏置单元的参数"""
a1,z2,a2,z2,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 #(25,401)
reg_D2=D2+(l/len(X))*t2 #(10,26)
return serialize(reg_D1,reg_D2)
2.6 Learning parameters using fmincg 优化参数
在成功地实现了神经网络代价函数和梯度计算之后，使用minimize来学习一个很好的集合参数。训练后完成后，通过计算得到正确示例的百分比来报告分类器的训练准确性。如果您的实现是正确的，您应该看到报告的训练准确率约为95.3%(这可能由于随机初始化而变化约1%)。通过训练神经网络进行更多的迭代，可以获得更高的训练精度。尝试训练神经网络进行更多的迭代(例如，将MaxIter设置为400)，并改变正则化参数。
def nn_training(X,y):
init_theta=random_init(10285)
res = opt.minimize(fun=regularized_cost,
x0=init_theta,
args=(X, y, 1),
method='TNC',
jac=regularized_gradient,
options={'maxiter': 400})
return res
res = nn_training(X, y)
print(res)
运行结果：
fun: 0.53044032678964603
jac: array([ -4.38866225e-03, -2.11248326e-12, 4.38829369e-13, ...,
-4.84312600e-04, 5.02428977e-04, -3.28169024e-04])
message: 'Converged (|f_n-f_(n-1)| ~= 0)'
nfev: 228
nit: 17
status: 1
success: True
x: array([ 0.85774573, 0.08465805, -0.00913285, ..., -2.81358358,
-2.93387811, -1.68379532])
2.7 accuracy（精确度）
from sklearn.metrics import classification_report
def accuracy(theta,X,y):
a1,z2,a2,z2,h=feed_forward(theta,X)
y_pred=np.argmax(h,axis=1)+1 #y_pred:(5000,) axis=1：按行查找最大元素 axis=0：按列查找最大元素
print(classification_report(y, y_pred))
accuracy(theta,X,y_raw)
运行结果：
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
micro avg 0.98 0.98 0.98 5000
macro avg 0.98 0.98 0.98 5000
weighted avg 0.98 0.98 0.98 5000
注意：classification_report解释
https://blog.csdn.net/genghaihua/article/details/81155200
三、Visualizing the hidden layer 可视化隐藏层
了解你的神经网络正在学习的一种方法是可视化被隐藏的单位捕捉到的表示。给定一个特定的隐藏单元，一种可视化的方法它计算的是找到一个输入x，它将导致它**(也就是说，**值(a(L)i)接近1)。对于您训练的神经网络，请注意Θ(1)的第一行是表示i的参数的401维向量。隐藏单位。如果我们放弃偏置项，我们得到一个400维向量，它表示从每个输入像素到隐藏单元的权重，这是一种将捕获的“表示”可视化的方法。隐藏单元是将这个400维向量重新塑造成20×20的图像并显示它。它将显示一个图像有25个单元，每个单元对应于网络中的一个隐藏单元。
def plot_hidden(theta):
t1,t2=deserialize(theta)
t1=t1[:,1:]
fig,ax_array=plt.subplots(5,5,sharex=True,sharey=True,figsize=(6,6))
for r in range(5):
for c in range(5):
ax_array[r,c].matshow(t1[r*5+c].reshape(20,20),cmap='gray_r')
plt.xticks([])
plt.yticks([])
plt.show()
plot_hidden(theta)
运行结果：