# -*- coding: utf-8 -*- ''' Created on 2018年1月26日 @author: Jason.F @summary: 多層感知器實現，加梯度檢驗 訓練集：http://yann.lecun.com/exdb/mnist/ train-images-idx3-ubyte: training set images train-labels-idx1-ubyte: training set labels t10k-images-idx3-ubyte: test set images t10k-labels-idx1-ubyte: test set labels ''' import pandas as pd import numpy as np import time import sys import os import struct from scipy.special import expit import matplotlib.pyplot as pltclass MLPGradientCheck(object):def __init__(self,n_output,n_features,n_hidden=30,l1=0.0,l2=0.0,epochs=500,eta=0.001,alpha=0.0,decrease_const=0.0,shuffle=True,minibatches=1,random_state=None):np.random.seed(random_state)self.n_output=n_output #輸出層數量self.n_features=n_features #輸入層數量self.n_hidden=n_hidden #隱層數量self.w1,self.w2 = self._initial_weights() #初始化權值系數self.l1=l1 #l1正則化系數self.l2=l2 #l2正則化系數self.epochs=epochs #迭代次數self.eta=eta #學習速率self.alpha=alpha #動量學習進度的參數，在上一輪迭代基礎上增加一個因子，用于加快權重更新的學習self.decrease_const=decrease_const #用于降低自適應學習速率n的常熟d,隨著迭代次數的增加而遞減顆更好地收斂self.shuffle=shuffle self.minibatches=minibatches#每批次訓練樣本數def _ecode_labels(self,y,k):onehot = np.zeros((k,y.shape[0]))for idx ,val in enumerate(y):onehot[val,idx]=1.0return onehotdef _initial_weights(self):w1 = np.random.uniform(-1.0,1.0,size=self.n_hidden*(self.n_features+1))w1 = w1.reshape(self.n_hidden,self.n_features+1)w2 = np.random.uniform(-1.0,1.0,size=self.n_output*(self.n_hidden+1))w2 = w2.reshape(self.n_output,self.n_hidden+1)return w1,w2def _sigmoid(self,z):#expit is equivalent to 1.0/(1.0+np.exp(-z))return expit(z)def _sigmoid_gradient(self,z):sg =self._sigmoid(z)return sg*(1-sg)def _add_bias_unit(self,X,how='column'):if how == 'column':#列X_new = np.ones((X.shape[0],X.shape[1]+1))X_new[:,1:]=Xelif how =='row':#行X_new = np.ones((X.shape[0]+1,X.shape[1]))X_new[1:,:]=Xelse:raise AttributeError('`how` must be `column` or `row`')return X_newdef _feedforwrd(self,X,w1,w2):a1=self._add_bias_unit(X, how='column')z2=w1.dot(a1.T)a2=self._sigmoid(z2)a2=self._add_bias_unit(a2, how='row')z3=w2.dot(a2)a3=self._sigmoid(z3)return a1,z2,a2,z3,a3def _L2_reg(self,lambda_,w1,w2):return (lambda_/2.0)*(np.sum(w1[:,1:]**2)+np.sum(w2[:,1:]**2))def _L1_reg(self,lambda_,w1,w2):return (lambda_/2.0)*(np.abs(w1[:,1:]).sum()+np.abs(w2[:,1:]).sum())def _get_cost(self,y_enc,output,w1,w2):term1 = -y_enc *(np.log(output))term2 = (1-y_enc) * np.log(1-output)cost = np.sum(term1-term2)L1_term =self._L1_reg(self.l1, w1, w2)L2_term =self._L2_reg(self.l2, w1, w2)cost =cost + L1_term +L2_termreturn costdef _get_gradient(self,a1,a2,a3,z2,y_enc,w1,w2):#backpropagationsigma3 = a3-y_encz2 = self._add_bias_unit(z2, how='row')sigma2 = w2.T.dot(sigma3) * self._sigmoid_gradient(z2)sigma2 = sigma2[1:,:]grad1 = sigma2.dot(a1)grad2 = sigma3.dot(a2.T)#regularizegrad1[:,1:] += (w1[:,1:] * (self.l1+self.l2))grad2[:,1:] += (w2[:,1:] * (self.l1+self.l2))return grad1,grad2def predict(self,X):a1,z2,a2,z3,a3 = self._feedforwrd(X, self.w1, self.w2)y_pred = np.argmax(z3,axis=0)return y_preddef fit(self,X,y,print_progress=False):self.cost_=[]X_data,y_data =X.copy(),y.copy()y_enc = self._ecode_labels(y, self.n_output)delta_w1_prev =np.zeros(self.w1.shape)delta_w2_prev =np.zeros(self.w2.shape)for i in range(self.epochs):#adaptive learning rateself.eta /= (1+self.decrease_const*i)if print_progress:sys.stderr.write('\rEpoch:%d/%d'%(i+1,self.epochs))sys.stderr.flush()if self.shuffle:idx = np.random.permutation(y_data.shape[0])X_data,y_data = X_data[idx],y_data[idx]mini = np.array_split(range(y_data.shape[0]),self.minibatches)for idx in mini:#feedbacka1,z2,a2,z3,a3 = self._feedforwrd(X[idx], self.w1, self.w2)cost = self._get_cost(y_enc=y_enc[:,idx], output=a3, w1=self.w1, w2=self.w2)self.cost_.append(cost)#compute gradient via backpropagationgrad1,grad2 = self._get_gradient(a1=a1,a2=a2,a3=a3,z2=z2,y_enc=y_enc[:,idx],w1=self.w1,w2=self.w2)##start gradient checking()grad_diff = self._gradient_checking(X=X[idx], y_enc=y_enc[:,idx], w1=self.w1,w2=self.w2, epsilon=1e-5, grad1=grad1, grad2=grad2)if grad_diff<=1e-7:print ('Ok:%s'%grad_diff)elif grad_diff <=1e-4:print ('Warning:%s'%grad_diff)else:print ('Problem:%s'%grad_diff)##end gradient checking#update weights[alpha*delta_w_prev] for momentum learningdelta_w1,delta_w2 = self.eta *grad1,self.eta*grad2self.w1 -= (delta_w1 +(self.alpha * delta_w1_prev))self.w2 -= (delta_w2 +(self.alpha * delta_w2_prev))delta_w1_prev,delta_w2_prev=delta_w1,delta_w2return self def _gradient_checking(self,X,y_enc,w1,w2,epsilon,grad1,grad2):'''Apply gradient checking (for debugging only)Returns:relative_error:float,Relative error between the numerically approximated gradients and the backpropagated gradients.'''num_grad1 = np.zeros(np.shape(w1))epsilon_ary1 = np.zeros(np.shape(w1))for i in range(w1.shape[0]):for j in range(w1.shape[1]):epsilon_ary1[i,j]=epsilona1,z2,a2,z3,a3 =self._feedforwrd(X, w1- epsilon_ary1, w2)cost1 = self._get_cost(y_enc, a3, w1- epsilon_ary1, w2)a1,z2,a2,z3,a3 =self._feedforwrd(X, w1+ epsilon_ary1, w2)cost2 = self._get_cost(y_enc, a3, w1+ epsilon_ary1, w2)num_grad1[i,j] =(cost2-cost1)/(2*epsilon)epsilon_ary1[i,j]=0num_grad2 = np.zeros(np.shape(w2))epsilon_ary2 = np.zeros(np.shape(w2))for i in range(w2.shape[0]):for j in range(w2.shape[1]):epsilon_ary2[i,j]=epsilona1,z2,a2,z3,a3 =self._feedforwrd(X, w1, w2- epsilon_ary2)cost1 = self._get_cost(y_enc, a3, w1, w2- epsilon_ary2)a1,z2,a2,z3,a3 =self._feedforwrd(X, w1, w2+ epsilon_ary2)cost2 = self._get_cost(y_enc, a3, w1, w2+ epsilon_ary2)num_grad2[i,j] =(cost2-cost1)/(2*epsilon)epsilon_ary2[i,j]=0num_grad =np.hstack((num_grad1.flatten(),num_grad2.flatten()))#數值梯度grad = np.hstack((grad1.flatten(),grad2.flatten()))#解析梯度norm1 = np.linalg.norm(num_grad - grad)norm2 = np.linalg.norm(num_grad)norm3 = np.linalg.norm(grad)relative_error=norm1/(norm2+norm3)return relative_errordef load_mnist(path,kind='train'):#load mnist data from pathlabels_path = os.path.join(path,'%s-labels.idx1-ubyte'%kind)images_path = os.path.join(path,'%s-images.idx3-ubyte'%kind)with open(labels_path,'rb') as lbpath:magic,n =struct.unpack('>II',lbpath.read(8))labels = np.fromfile(lbpath,dtype = np.uint8)with open(images_path,'rb') as imgpath:magic,num,rows,cols =struct.unpack('>IIII',imgpath.read(16))images = np.fromfile(imgpath,dtype = np.uint8).reshape(len(labels),784)#28X28像素return images,labels if __name__ == "__main__": start = time.clock() #導入數據集homedir = os.getcwd()#獲取當前文件的路徑X_train,y_train = load_mnist(homedir+'\\mnist', kind='train')print ('Rows:%d,columns:%d'%(X_train.shape[0],X_train.shape[1]))X_test,y_test = load_mnist(homedir+'\\mnist', kind='t10k')print ('Rows:%d,columns:%d'%(X_test.shape[0],X_test.shape[1]))#將特征矩陣的784像素向量還原成18X28圖像fig,ax = plt.subplots(nrows=2,ncols=5,sharex=True,sharey=True)ax=ax.flatten()for i in range(10):img = X_train[y_train==i][0].reshape(28,28)ax[i].imshow(img,cmap='Greys',interpolation='nearest')ax[0].set_xticks([])ax[0].set_yticks([])plt.tight_layout()plt.show()#繪制相同數字的多個示例fig,ax = plt.subplots(nrows=5,ncols=5,sharex=True,sharey=True)ax=ax.flatten()for i in range(25):img = X_train[y_train==7][i].reshape(28,28)ax[i].imshow(img,cmap='Greys',interpolation='nearest')ax[0].set_xticks([])ax[0].set_yticks([])plt.tight_layout()plt.show()'''#將數據存儲為CSV格式np.savetxt('train_img.csv',X_train,fmt='%i',delimiter=',')#指定存儲數據類型為整型，分隔符為,np.savetxt('train_labels.csv',y_train,fmt='%i',delimiter=',')np.savetxt('test_img.csv',X_test,fmt='%i',delimiter=',')np.savetxt('test_labels.csv',y_test,fmt='%i',delimiter=',')#從csv加載數據X_train = np.genfromtxt('train_img.csv',dtype=int,delimiter=',')y_train = np.genfromtxt('train_labels.csv',dtype=int,delimiter=',')X_test = np.genfromtxt('test_img.csv',dtype=int,delimiter=',')y_test = np.genfromtxt('test_labels.csv',dtype=int,delimiter=',')'''#創建感知器模型nn_check = MLPGradientCheck(n_output=10,n_features=X_train.shape[1],n_hidden=50,l2=0.0,l1=0.0,epochs=10,eta=0.001,alpha=0.0,decrease_const=0.0,shuffle=True,minibatches=1,random_state=1)nn_check.fit(X_train,y_train,print_progress=True)print (nn_check.fit(X_train[:5],y_train[:5],print_progress=False))#梯度檢驗結果，計算成本大end = time.clock() print('finish all in %s' % str(end - start))

結果：

執行太長時間了，12個小時就跑了2個循環。

總結

以上是生活随笔為你收集整理的【Python-ML】神经网络-多层感知器增加梯度检验的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。