机器学习神经网络---反向传播Python(吴恩达课后作业)

任务

依旧是手写数字识别，通过神经网络的反向传播来训练参数，得到最优化模型

一、引入库和数据集

import numpy as np
import scipy.io as io
data=io.loadmat('C:/Users/Administrator/Desktop/新建文件夹/ex4data1.mat')
X=data['X']
y=data['y'].ravel()  # 将y一维化
data2=io.loadmat('C:/Users/Administrator/Desktop/新建文件夹/ex4weights.mat')

X.shape    #(5000, 400)
y.shape    #(5000,)

二、获取输入变量和输出变量以及权重

1.输入变量

X=np.insert(X,0,1,axis=1)
X.shape    # (5000, 401)

2.输出变量
因为输出结果有十种类型：0,1,2,3,4,5,6,7,8,9,所以需要将y拓展。

def y_expand(y):
    result=[]
    for i in y:
        y_array=np.zeros(10)
        y_array[i-1]=1
        result.append(y_array)
    return np.array(result)
y=y_expand(y)
y.shape   #  (5000, 10)

3.权重

theta1=data2['Theta1']
theta2=data2['Theta2']
theta1.shape  #(25, 401)
theta2.shape   #(10, 26)

三、展开参数

1.把所有参数展开成一个长向量的形式

def serialize(a,b):
    return np.r_[a.ravel(),b.ravel()]
'''def serialize(a,b):
    return np.append(a.ravel(),b.ravel())'''

theta=serialize(theta1,theta2)
theta.shape  #(10285,)

2.拆分

def deserialize(theta):
    theta1=theta[:25*401].reshape(25,401)
    theta2=theta[25*401:].reshape(10,26)
    return theta1,theta2

三、前向传播

# X.shape(5000,401)
# t1,t2=deserialize(theta)
# t1.shape  (5,401)
# t2.shape (10,26)
def sigmoid(z):
    return 1/(1+np.exp(-z))
    
def feed_forward(theta,X):
    t1,t2=deserialize(theta)
    a1=X
    z2=a1@t1.T
    a2=sigmoid(z2)
    a2=np.insert(a2,0,1,axis=1)
    z3=a2@t2.T
    a3=sigmoid(z3)
    return a1,z2,a2,z3,a3   

四、代价函数

1.代价函数

def cost(theta,X,y):
    a1,z2,a2,z3,h=feed_forward(theta,X)
    J=np.sum(y*np.log(h)+(1-y)*np.log(1-h))/(-len(X))
    return J   # 向量化形式
'''
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+=first+second
    J=J/(-len(X))
    return J
  '''

2.正则化后的代价函数

def costReg(theta,X,y,lam=1):
    t1,t2=deserialize(theta)
    first=np.sum(np.power(t1[:,1:],2))   #忽略每层的偏置项
    second=np.sum(np.power(t2[:,1:],2))  #忽略每层的偏置项
    reg=lam*(first+second)/(2*len(X))
    return cost(theta,X,y)+reg

五、反向传播(梯度矩阵)

1.激活函数sigmoid的导数

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

2.梯度矩阵

def gradient(theta,X,y):
    t1,t2=deserialize(theta)
    a1,z2,a2,z3,a3=feed_forward(theta,X)
    #每一层的误差
    d3=a3-y  #(5000,10)
    d2=d3@t2[:,1:] *sigmoid_gradient(z2)    #(5000,25)
    #梯度
    D2=d3.T@a2/len(X)    #(10,26)
    D1=d2.T@a1/(len(X))    #(25,401)
    D=serialize(D1,D2)
    return D  

3.正则化后的梯度矩阵

def gradientReg(theta,X,y,lam=1):
    D=gradient(theta,X,y)
    D1,D2=deserialize(D)
    t1,t2=deserialize(theta)
    #不惩罚偏置单元的参数
    D1[:,1:]=D1[:,1:]+t1[:,1:]*lam/(len(X)) 
    D2[:,1:]=D2[:,1:]+t2[:,1:]*lam/(len(X))
    D_Reg=serialize(D1,D2)
    return D_Reg

'''
def regularized_gradient(theta, X, y, l=1):
    """不惩罚偏置单元的参数"""
    t1,t2=deserialize(theta)
    D1, D2 = deserialize(gradient(theta, X, y))
    t1[:,0] = 0   
    t2[:,0] = 0#这样做的话已经改变了theta的值，所以调用regularized_gradient(theta,X,y,1)后，theta[:25*401].reshape(25,401)[:,0]的值变为0.
这样用高级优化算法时，计算出错！
    reg_D1 = D1 + (l / len(X)) * t1
    reg_D2 = D2 + (l / len(X)) * t2
    return serialize(reg_D1, reg_D2)'''    
#改进为如下：（使用副本）
'''
def regularized_gradient(theta, X, y, l=1):
    """不惩罚偏置单元的参数"""
    t1,t2=deserialize(theta)
    D1, D2 = deserialize(gradient(theta, X, y))
    t11=t1.copy()
    t22=t2.copy()
    t11[:,0] = 0   
    t22[:,0] = 0
    reg_D1 = D1+ (l / len(X)) * t11
    reg_D2 = D2 + (l / len(X)) * t22
    return serialize(reg_D1, reg_D2)'''

六、梯度检测

将数值计算后的梯度与反向传播正则化后的梯度进行比较

def gradient_checking(theta,X,y,e):
    def a_gradient_estimation(plus,minus):
        return (costReg(plus,X,y)-costReg(minus,X,y))/(2*e)   #对每个参数计算数值梯度
    gradient_estimation=[]
    for i in range(len(theta)):
        plus=theta.copy()   # 要使用副本，否则会改变原来的theta
        minus=theta.copy()
        plus[i]=plus[i]+e
        minus[i]=minus[i]-e
        gradient_i=a_gradient_estimation(plus,minus)
        gradient_estimation.append(gradient_i)
   
    gradient_estimation=np.array(gradient_estimation)
    analytic_gradient=gradientReg(theta,X,y)
    diff=np.linalg.norm(gradient_estimation-analytic_gradient)/np.linalg.norm(gradient_estimation+analytic_gradient)  #np.linalg.norm()是指求范数，默认是二范数
   '''diff应该是一个很小很小的值，也就是说gradient_estimation和analytic_gradient很接近，如果diff较大，则说明出错了'''
    return diff

调用该函数，结果如下：

说明反向传播的执行过程是正确的

七、参数theta随机初始化

目的：打破数据的对称性

init_theta = np.random.uniform(-0.12,0.12,10285)
#从服从的均匀分布范围中随机选取size大小的值
'''
def random_init(size):
    return np.random.uniform(-0.12,0.12,size)

init_theta=random_init(10285)  #25*401+10*26=10285
'''

八、利用优化函数优化参数

import scipy.optimize as opt
res = opt.minimize(fun =costReg,x0=init_theta,args=(X,y,1),method='TNC',
                   jac=gradientReg)                               

注：
（1）nfev是指function
evaluation的数目，即目标优化函数被调用的次数，这里为2442次，success=True表明优化过程是顺利的。
（2）不建议在opt.minimize()中加入options={‘maxiter’:400}，原因是限制了调用次数。如下：

res = opt.minimize(fun =costReg,x0=init_theta,args=(X,y,1),method='TNC',
                   jac=gradientReg,options={'maxiter':400})

可以发现，success=False,nfev=400,message告诉我们 ‘Max. number of function evaluations reached’。

九、利用优化后的参数重新进行预测，并计算准确率

raw_y=data['y']
_,_,_,_,h=feed_forward(res.x,X)  #res.x表示优化后的参数
y_pred=np.argmax(h,axis=1)+1  #这里可看上一节
y_pred=y_pred.reshape(raw_y.shape) #要执行这一步的原因看下面
accuracy=np.mean(y_pred==raw_y)
accuracy   #0.9962

注：

从图片可以看出他们的维度不一样，如果去掉这一步，accuracy将等于0.1

'''
raw_y=data['y']
_,_,_,_,h=feed_forward(res.x,X)  #res.x表示优化后的参数
y_pred=np.argmax(h,axis=1)+1  #这里可看上一节在这里插入代码片
accuracy=np.mean(y_pred==raw_y)
accuracy   # 0.1
'''

十、补充

神经网络的的过程：
1.引入数据集和权重
2.参数展开
3.前向传播
4.定义正则化代价函数
5.反向传播（求梯度矩阵）
6.梯度检测
7.高级优化算法，得到优化后的权重
8.预测和求精度