前言：
本文属于学习笔记性质。为了让自己更深入地理解卷积神经网络，我只用numpy、pandas等几个库手搓了一个识别MNIST数字的CNN。500张图单次训练，准确率70-80%。
注意：
1.代码并非原创，主要参考了下面的文章，我按自己的思路进行了一些改动。
(29条消息) python神经网络案例——CNN卷积神经网络实现mnist手写体识别_python cnn_腾讯数据架构师的博客-CSDN博客
2.可能有一些错误，欢迎批评指正。
3.有些地方非常话痨，还请见谅。
本网络的架构：
输入28*28分辨率的图像，
卷积层1包含8个5*5的卷积核，输出8张24*24的图像，
池化层1进行2*2最大池化，输出8张12*12的图像。
卷积层2包含16组，每组8个5*5的卷积核，输出16张8*8的图像，
池化层2进行2*2最大池化，输出16张4*4的图像。
全连接层将池化层2的结果拉直形成1*256的向量，输出1*10的one-hot编码。
0.所调用的库
import pandas as pd
import numpy as np为了方便观察loss和特征图，可以再加一个：
import matplotlib.pyplot as plt1.激活函数
常见激活函数主要包括relu、sigmoid、tanh、softmax等， 此处使用relu（卷积层）和softmax（最后的全连接层）。
relu函数写的非常冗长，也可以用np.nditer()简化。
dreluy：
输入：反向传播时传到relu激活层的张量。
输出：对该张量进行relu求导操作后的张量。dsoftmaxy类似。
（注：导数既可以写成是x的函数也可以写成是y的函数，例如y = x^2, dy/dx = 2x = 2*sqrt(y)，因为求导操作是反向传播时才需要的，所以这里把导函数写作y的函数）
def relu(ip_array):
op_array=np.zeros((ip_array.shape))
if ip_array.ndim == 2:
for i in range(ip_array.shape[0]):
for j in range(ip_array.shape[1]):
op_array[i][j] = max(0,ip_array[i][j])
elif ip_array.ndim == 3:
for i in range(ip_array.shape[0]):
for j in range(ip_array.shape[1]):
for k in range(ip_array.shape[2]):
op_array[i][j][k] = max(0,ip_array[i][j][k])
elif ip_array.ndim == 4:
for i in range(ip_array.shape[0]):
for j in range(ip_array.shape[1]):
for k in range(ip_array.shape[2]):
for l in range(ip_array.shape[3]):
op_array[i][j][k][l] = max(0,ip_array[i][j][k][l])
return op_array
def dreluy(ip_array):
op_array=np.zeros((ip_array.shape))
if ip_array.ndim == 2:
for i in range(ip_array.shape[0]):
for j in range(ip_array.shape[1]):
if ip_array[i][j] >0:
op_array[i][j] = 1
else:
op_array[i][j] = 0
elif ip_array.ndim == 3:
for i in range(ip_array.shape[0]):
for j in range(ip_array.shape[1]):
for k in range(ip_array.shape[2]):
if ip_array[i][j][k] >0:
op_array[i][j][k] = 1
else:
op_array[i][j][k] = 0
elif ip_array.ndim == 4:
for i in range(ip_array.shape[0]):
for j in range(ip_array.shape[1]):
for k in range(ip_array.shape[2]):
for l in range(ip_array.shape[3]):
if ip_array[i][j][k][l] >0:
op_array[i][j][k][l] = 1
else:
op_array[i][j][k][l] = 0
return op_array
def softmax(x):
#softmax只允许输入一个列向量，输出一个列向量
y=np.zeros(x.shape)
sum1 = np.sum(np.exp(x))
for i in range(len(x)):
y[i]=np.exp(x[i])/sum1
return y
def dsoftmaxy(y):
dy = np.eye(len(y))-y
for i in range(len(dy)):
dy[i] = dy[i]*y[i]
return dy
2.下面是卷积和池化操作需要的一些工具函数：
2.1. 截取待操作区域的函数（进行卷积和池化操作时用）：
vertex_i和vertex_j是待操作区域的顶角坐标，输出的是以此为顶角，尺寸与卷积核或池化层尺寸相同的的一片区域。
def area(ip_array,k,i,j,filter_size,stride):
'''
ip_array:第k层i*j图片
k：深度，如果为1，本函数返回的是顶角为i,j的待卷积/池化的单片区域，
否则是第k层的区域
'''
vertex_i = stride*i
vertex_j = stride*j
if ip_array.ndim == 2:
return ip_array[vertex_i: vertex_i + filter_size,vertex_j: vertex_j + filter_size]
if ip_array.ndim == 3:
return ip_array[k,vertex_i: vertex_i + filter_size,vertex_j: vertex_j + filter_size]2.2. 获取一片区域中最大元素所在的位置（池化层的反向传播用）：
分情况，如果由多个相同的最大元素，则只输出第一个元素的坐标
def array_argmax(ip_array):
coord = np.where(ip_array == np.max(ip_array))
if len(coord[0]) != 1 :
return coord[0][0],coord[0][1]
else:
return coord[0],coord[1]2.3. 补零（卷积层的反向传播用）：
def zp(ip_array,p):
'''
输入：4*4*3的图片，需要补零的层数（例如1层）
输出：补零后的5*5*3的图片，最外圈为零
'''
if p==0:
return ip_array
else:
if ip_array.ndim == 2:
op_array = np.zeros((ip_array.shape[1] + 2*p,ip_array.shape[2] + 2*p))
for k in range(1):
for i in range(p,ip_array.shape[1]-p):
for j in range(p,ip_array.shape[2]-p):
op_array[i][j] = ip_array[i-p][j-p]
return op_array
elif ip_array.ndim == 3:
op_array = np.zeros((ip_array.shape[0],ip_array.shape[1] + 2*p,ip_array.shape[2] + 2*p))
for k in range(ip_array.shape[0]):
for i in range(p,op_array.shape[1]-p):
for j in range(p,op_array.shape[2]-p):
op_array[k][i][j] = ip_array[k][i-p][j-p]
return op_array
2.4. 卷积：
def conv(ip_array,kernel,b,stride=1):
'''
输入：
无深度：4*4的图片，2*2的核，步长
有深度：4*4*2的图片，一个2*2*2的核，步长
输出：一张卷积后的2*2图片，深度被积分掉
'''
if ip_array.ndim == 2:
op_rows = int((ip_array.shape[0] - kernel.shape[1])/stride + 1)
op_cols = int((ip_array.shape[1] - kernel.shape[1])/stride + 1)
op_array = np.zeros((op_rows,op_cols))
for k in range(1):
for i in range(op_rows):
for j in range(op_cols):
op_array[i][j] = (area(ip_array, k, i, j, kernel.shape[1],stride=1)*kernel).sum() + b
if ip_array.ndim == 3:
op_rows = int((ip_array.shape[1] - kernel.shape[2])/stride + 1)
op_cols = int((ip_array.shape[2] - kernel.shape[2])/stride + 1)
op_array = np.zeros((ip_array.shape[0],op_rows,op_cols))
for k in range(ip_array.shape[0]):
for i in range(op_rows):
for j in range(op_cols):
op_array[k][i][j] = (area(ip_array, k, i, j, kernel.shape[2],stride=1)*kernel[k]).sum() + b
op_array = np.sum(op_array,axis=0)
return op_array
2.5.生成卷积核：
生成的卷积核是一个四维张量，一共有num组，depth个size*size的方形卷积核。
num就是当前卷积层一共要用多少组卷积核进行卷积。depth就是原图有几个通道，或者上一个卷积层输出了几张图像。
class Kernel:
'''
输入：核的size、depth和数量
get_W:返回5*5的权重矩阵W，深度3
get_b:返回一个偏置b，深度3
update：得到微分矩阵后更新W和b
'''
def __init__(self, size, depth, nums):
self.nums = nums
self.depth = depth
self.val = 1/np.sqrt(size**2*depth+size**2*nums)
self.W = np.random.uniform(-self.val,self.val,(nums, depth, size, size))
self.b = np.zeros(self.nums)
def get_W(self):
return self.W
def get_b(self):
return self.b
3.以下是网络的主体：
3.1.卷积层
class Conv:
def __init__(self,ip_rows,ip_cols,depth,ker_size,ker_num,act_name="relu",lr=0.005,stride=1):
self.ip_rows = ip_rows
self.ip_cols = ip_cols
self.depth = depth
self.ker_size = ker_size
self.stride = stride
self.actor = act_name
self.ker_num = ker_num
self.kernel_groups = Kernel(self.ker_size,self.depth,self.ker_num)
self.op_rows = int((self.ip_rows - ker_size)/stride + 1)
self.op_cols = int((self.ip_cols - ker_size)/stride + 1)
self.op_array = np.zeros((self.ker_num,self.op_rows,self.op_cols))
self.lr=lr
def act(self,x):
if self.actor == "relu":
return relu(x)
elif self.actor == "sigmoid":
return sigm(x)
def dact(self,y):
if self.actor == "relu":
return dreluy(y)
elif self.actor == "sigmoid":
return dsigmy(y)
def forward(self,ip_array):
for l in range(self.ker_num):
self.op_array[l] = conv(ip_array,self.kernel_groups.get_W()[l],self.kernel_groups.get_b()[l],self.stride)
self.op_array = self.act(self.op_array)
return self.op_array
def backward(self,ip_array,grad_array,first_layer="N"):
#grad_array是上一层传到本层的梯度，一共是ker_num张
self.delta_array = np.zeros(grad_array.shape)
for num in range(self.ker_num):
self.delta_array[num] = self.dact(self.forward(ip_array))[num]*grad_array[num]
self.dW = np.zeros((self.ker_num,self.depth,self.ker_size,self.ker_size))
self.db = np.zeros(self.ker_num)
for num in range(self.ker_num):
self.db[num] = np.sum(self.delta_array[num])
self.kernel_groups.get_b()[num] += -self.lr*self.db[num]
for dep in range(self.depth):
if ip_array.ndim == 2:
self.dW[num][dep] = conv(ip_array,self.delta_array[num],0,self.stride)
if ip_array.ndim == 3:
self.dW[num][dep] = conv(ip_array[dep],self.delta_array[num],0,self.stride)
self.kernel_groups.get_W()[num][dep] += -self.lr*self.dW[num][dep]
self.padded_grad = zp(self.delta_array,self.ker_size-1) # ker_num个
self.next_grad = np.zeros(ip_array.shape) # depth个
if first_layer == "N":
for depth in range(ip_array.shape[0]):
for num1 in range(self.ker_num):
self.next_grad[depth] += conv(self.padded_grad[num1],np.rot90(self.kernel_groups.get_W()[num1][depth],2),0)
return self.next_grad
else:
pass这里主要说一下反向传播。
假设某卷积层的输入(ip_array)为2通道4*4的图，用3组2通道2*2的卷积核进行卷积，则输出(op_array)是3张3*3的图。则反向传播时，传到本层的梯度(grad_array)应该是
接下来要求的是：
根据正向传播过程：
可以推导出反向传播的公式（这里的星号有的表示点积，有的表示卷积，懒得区分了）：
对b_j，公式1输出的是一张3*3的array，而b_j是一个数，如何对b进行更新呢？
实际上，正向传播时，op_j中的每个元素都是ip1与W_1j、ip_2与W_2j的分别卷积并加和，然后再加上b_j得到的。因此，最终L的对b_j的导数是上面的公式输出的array中所有元素的加和。
如果不好理解，不妨看下面的例子：
y1 = W1∗x1+b, y2 = W2∗x2 + b, z = f(y1) + g(y2)，
因此
dz/db = dz/dy1*dy1/db + dz/dy2*dy2/db = f'(y1) + g'(y2)
对W_ij，从以下例子可以发现，A在op_array的四个元素中系数依次是1，2，4，5，也就是说op_array每个元素对A的导数应该分别是1，2，4，5；即A的导数是op1’*1+op2’*2+op3’*4+op4’*5，这正是梯度矩阵与ip_i卷积的结果。因此，公式2的最后一个偏微分是第i张ip_array，即ip_i，最后一个星号对应的是卷积操作。
类似地，对ip_i，由于op_j是ip_i和W_ij卷积得到，因此公式3的最后一个偏微分似乎应该是W_ij，最后一个星号是卷积操作。但实际操作比较复杂：首先需要把
小结，因为卷积操作比点乘复杂，而且卷积核参数过多（组数*每组通道数*行*列），所以卷积层的反向传播有点麻烦。我的心得是，只要记住传播到某层的梯度肯定与你要更新的东西具有相同的形状，就不会被绕晕或者忘记某个操作。
3.2.池化层
class Maxpooling:
'''
输入：输入图片的行数、列数、深度，池化尺寸和步长。
forward功能：输入4*4*k的待处理图片，输出2*2*k的池化结果
backward功能（上池化）：以k=1为例。
正向传播时，输入4*4的图片A，输出2*2的结果。
反向传播时，先创建一个与A同形状的零矩阵D。导数矩阵del(2,2)中的格点[m][n]对应D中的一个2*2区域：
D[m*stride : m*stride + pool_size][n*stride : n*stride + pool_size]，
D在该区域中仅有一个非零元素x，其位置就是A在该区域中最大元素所处的位置，其数值就是del[m][n]
返回D。
'''
def __init__(self,ip_rows,ip_cols,depth,pool_size,stride=2):
self.ip_rows = ip_rows
self.ip_cols = ip_cols
self.depth = depth
self.pool_size = pool_size
self.stride = stride
self.op_rows = int((self.ip_rows - pool_size)/stride + 1)
self.op_cols = int((self.ip_cols - pool_size)/stride + 1)
def forward(self,ip_array):
self.op_array = np.zeros((self.depth,self.op_rows,self.op_cols))
for k in range(self.depth): #池化在卷积后面，因此输入的图片一定是多层的
for i in range(self.op_rows):
for j in range(self.op_cols):
self.op_array[k][i][j] = np.max(area(ip_array[k],k,i,j,self.pool_size,self.stride))
return self.op_array
def backward(self,ip_array,sens_array):
self.back_array = np.zeros((ip_array.shape))
for k in range(self.depth):
for m in range(self.op_rows):
for n in range(self.op_cols):
i,j = array_argmax(area(ip_array[k],k,m,n,self.pool_size,self.stride))
self.back_array[k][m*self.stride+i,n*self.stride+j] = sens_array[k][m][n].copy()
return self.back_array
池化层的反向传播就容易多了，已经在注释里说明了。
3.3.全连接层
class NN:
'''
单层全连接BP神经网络。
初始化：输入维度、神经元个数（即输出维度）、激活函数为的名字（sigmoid，relu或softmax）
foward(X):对input进行一次正向传播
act/dact:根据输入的激活函数名字指定激活函数
backward(grad):设本层输出为z2，则应输入∂L/∂z2，该函数据此更新W2和b2，并给出∂L/∂z1
'''
def __init__(self,ipsize,opsize,act_name,lr=0.002):
self.ipsize = ipsize
self.opsize = opsize
self.actor = act_name
self.val = 1/np.sqrt(self.ipsize + self.opsize)
self.W = np.random.uniform(-self.val,self.val,(self.opsize,self.ipsize))
self.b = np.zeros(self.opsize)
self.z = np.zeros(self.opsize)
self.lr = lr
def act(self,x):
if self.actor == "relu":
return relu(x)
elif self.actor == "sigmoid":
return sigm(x)
elif self.actor == "softmax":
return softmax(x)
def dact(self,y):
if self.actor == "relu":
return dreluy(y)
elif self.actor == "sigmoid":
return dsigmy(y)
elif self.actor == "softmax":
return dsoftmaxy(y)
def forward(self,X):
self.X = X
self.h = np.dot(self.W,self.X) + self.b
self.y = self.act(self.h)
return self.y
def backward(self,grad):
if self.actor == "softmax":
self.del2 = np.dot(self.dact(self.y),grad)
else:
self.del2 = self.dact(self.y)*grad
self.del1 = np.dot(self.W.T,self.del2)
self.W += -self.lr*np.outer(self.del2,self.X.T)
self.b += -self.lr*self.del2
return self.del1
class CNN:
def __init__(self):
self.conv1 = Conv(28,28,1,5,8,lr=0.01)
self.pool1 = Maxpooling(24,24,8,2)
self.conv2 = Conv(12,12,8,5,16,lr=0.005)
self.pool2 = Maxpooling(8,8,16,2)
self.fc1 = NN(256,10,"softmax",lr=0.002)
def forward(self,one_pic):
self.conved1 = self.conv1.forward(one_pic)
self.pooled1 = self.pool1.forward(self.conved1)
self.conved2 = self.conv2.forward(self.pooled1)
self.pooled2 = self.pool2.forward(self.conved2)
self.res = self.fc1.forward(self.pooled2.flatten().reshape(-1))
return self.res
def backward(self,one_pic,Y):
self.forward(one_pic)
self.grad0 = -Y/self.res
self.gradfc1 = self.fc1.backward(self.grad0)
self.gradfc1 = self.gradfc1.reshape(16,4,4)
self.gradp2 = self.pool2.backward(self.conved2,self.gradfc1)
self.gradc2 = self.conv2.backward(self.pooled1,self.gradp2)
self.gradp1 = self.pool1.backward(self.conved1,self.gradc2)
self.conv1.backward(one_pic,self.gradp1,first_layer="Y")4.整合
损失函数为交叉熵
class CNN:
def __init__(self):
self.conv1 = Conv(28,28,1,5,8,lr=0.01)
self.pool1 = Maxpooling(24,24,8,2)
self.conv2 = Conv(12,12,8,5,16,lr=0.005)
self.pool2 = Maxpooling(8,8,16,2)
self.fc1 = NN(256,10,"softmax",lr=0.002)
def forward(self,one_pic):
self.conved1 = self.conv1.forward(one_pic)
self.pooled1 = self.pool1.forward(self.conved1)
self.conved2 = self.conv2.forward(self.pooled1)
self.pooled2 = self.pool2.forward(self.conved2)
self.res = self.fc1.forward(self.pooled2.flatten().reshape(-1))
return self.res
def backward(self,one_pic,Y):
self.forward(one_pic)
self.grad0 = -Y/self.res
self.gradfc1 = self.fc1.backward(self.grad0)
self.gradfc1 = self.gradfc1.reshape(16,4,4)
self.gradp2 = self.pool2.backward(self.conved2,self.gradfc1)
self.gradc2 = self.conv2.backward(self.pooled1,self.gradp2)
self.gradp1 = self.pool1.backward(self.conved1,self.gradc2)
self.conv1.backward(one_pic,self.gradp1,first_layer="Y")5.训练测试部分
注：MINST_DATA和MINST_TRAIN是把MINST的数据库转化为了CSV格式。该格式中，第一列为实际数字，第2-785列为展开后的28*28像素。
定义一个用来看loss的函数：
def plot_curve(data):
fig=plt.figure()
plt.plot(range(len(data)),data,color="red")
plt.legend(["loss"],loc="upper right")
plt.xlabel("step")
plt.ylabel("loss")
plt.show()if __name__ == '__main__':
#数据预处理，将图片进行Z-score标准化，标签进行one-hot转码
num1 = 600
iter = 1
CNN1 = CNN()
data_csv=pd.read_csv("C://mnist_train.csv")
X=data_csv.iloc[:num1,1:].values/255
for i in range(X.shape[0]):
miu=np.mean(X[i])
sigma=np.std(X[i])
for j in range(X.shape[1]):
X[i][j]=(X[i][j]-miu)/sigma.copy()
X=X.reshape(num1,28,28)
Y=data_csv.iloc[:num1,0].values
Y_oh=np.zeros((len(Y),10))
for i in range(len(Y)):
Y_oh[i][Y[i]]=1
#随便看一张图
plt.matshow(X[0])
loss=[]
for epoch in range(iter):
y=np.zeros(len(Y))
right=0
for i in range(X.shape[0]):
CNN1.backward(X[i],Y_oh[i])
if np.argmax(CNN1.forward(X[i])) == Y[i]:
right += 1
if i%10==0:
print(np.argmax(CNN1.forward(X[i])),Y[i])
loss.append(np.sum(-Y_oh[i]*np.log(CNN1.forward(X[i]))))
if i%100==0:
plot_curve(loss)
print("epoch no.%s"%epoch,"accurancy:",right/len(Y))
#看一下训练完的卷积核长啥样
for i in range(4):
plt.matshow(CNN1.conv1.kernel_groups.get_W()[i])
for i in range(2):
for j in range(2):
plt.matshow(CNN1.conv2.kernel_groups.get_W()[i][j])
#开始测试
print("\ntest begins:\n")
num2 = 100
test_csv=pd.read_csv("C://mnist_train.csv")
X1=test_csv.iloc[:num2,1:].values/255
for i in range(X1.shape[0]):
miu=np.mean(X1[i])
sigma=np.std(X1[i])
for j in range(X1.shape[1]):
X1[i][j]=(X1[i][j]-miu)/sigma.copy()
X1=X1.reshape(num2,28,28)
Y1=test_csv.iloc[:num2,0].values
right1 = 0
for i in range(X1.shape[0]):
if np.argmax(CNN1.forward(X1[i])) == Y1[i]:
right1 += 1
if i%4==0:
print(np.argmax(CNN1.forward(X1[i])),Y1[i])
#看一下特征图长啥样
for i in range(8):
plt.matshow(CNN1.conv1.forward(X1[i]))
print("test accurancy:",right1/len(Y1))
plt.show()
训练结果：
打印结果：
epoch no.0 accurancy: 0.7816666666666666
test begins:
test accurancy: 0.8
第一个卷积层的前四个卷积核：
后两个核看起来比较明显，是检测斜线特征的，看一下对应的特征图：
文章出处登录后可见！
已经登录？立即刷新