實驗指導書       下載密碼:u8dl

本篇部落格主要講解，吳恩達機器學習第四周的程式設計作業，作業內容主要是對手寫數字進行識別，是一個十分類問題，要求使用兩種不同的方法實現：一是用之前講過的邏輯迴歸實現手寫數字識別，二是用本週講的神經網路實現手寫數字識別。實驗的原始版本是用Matlab實現的，本篇部落格主要用Python來實現。

目錄

1.實驗包含的檔案

2.資料集

3.利用邏輯迴歸進行手寫數字識別(多分類)

4.邏輯迴歸多分類器的完整專案程式碼

5.利用神經網路進行手寫數字識別

6.神經網路分類器的完整專案程式碼

1.實驗包含的檔案

其中前六個檔案已經寫好了，實驗任務是完成後四個檔案的關鍵程式碼。

2.資料集

ex3data1.mat檔案包含了5000個手寫數字的訓練樣本，每個訓練樣本的輸入是一張20*20的灰度圖，輸入特徵是400維的，輸出是一個標量，取值範圍1-10，1-9依次代表數字1-9，10代表數字0.

3.利用邏輯迴歸進行手寫數字識別(多分類)

• 開啟邏輯迴歸多分類器的主程式ex3.py
• 載入需要的庫

'''系統自帶的庫'''
import matplotlib.pyplot as plt
import numpy as np
import scipy.io as scio   #引入讀取.mat檔案的庫

'''自己編寫的庫'''
import displayData as dd  #引入視覺化資料集的程式
import lrCostFunction as lCF #引入邏輯迴歸多分類代價函式
import oneVsAll as ova    #引入邏輯迴歸訓練程式
import predictOneVsAll as pova  #引入邏輯迴歸預測程式

• 首先載入資料集，並可視化前100個訓練樣本

input_layer_size = 400  # 輸入特徵的維度 每張數字圖片20*20=400
num_labels = 10         # 10個標籤 注意“0”對應第十個標籤   1-9一次對應第1-9個標籤
                       

'''第一部分  載入手寫數字訓練資料集 並可視化部分訓練樣本'''

# 載入訓練集
print('載入並可視化資料 ...')

data = scio.loadmat('ex3data1.mat')  #讀取訓練集 包括兩部分 輸入特徵和標籤
X = data['X']   #提取輸入特徵 5000*400的矩陣  5000個訓練樣本 每個樣本特徵維度為400 一行代表一個訓練樣本
y = data['y'].flatten() #提取標籤 data['y']是一個5000*1的2維陣列 利用flatten()將其轉換為有5000個元素的一維陣列

m = y.size  #訓練樣本的數量

# 隨機抽取100個訓練樣本 進行視覺化
rand_indices = np.random.permutation(range(m)) #獲取0-4999 5000個無序隨機索引
selected = X[rand_indices[0:100], :]  #獲取前100個隨機索引對應的整條資料的輸入特徵

dd.display_data(selected)   #呼叫視覺化函式 進行視覺化

• 檢視視覺化程式displayData.py

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)
    display_cols = np.ceil(m / display_rows).astype(int)

    # 待顯示的每張圖片之間的間隔
    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')

• 檢視視覺化效果

• 編寫邏輯迴歸的代價函式(正則化)並做簡單測試

注意使用向量化程式設計。

'''第2-1部分  編寫邏輯迴歸的代價函式(正則化),做簡單的測試'''

# 邏輯迴歸（正則化）代價函式的測試用例
print('Testing lrCostFunction()')

theta_t = np.array([-2, -1, 1, 2]) #初始化假設函式的引數  假設有4個引數
X_t = np.c_[np.ones(5), np.arange(1, 16).reshape((3, 5)).T/10] #輸入特徵矩陣 5個訓練樣本 每個樣本3個輸入特徵，前面新增一列特徵=1
y_t = np.array([1, 0, 1, 0, 1]) #標籤 做2分類

lmda_t = 3  #正則化懲罰性係數
cost, grad = lCF.lr_cost_function(theta_t, X_t, y_t, lmda_t) #傳入代價函式 

#返回當前的代價函式值和梯度值  與期望值比較 驗證程式的正確性

np.set_printoptions(formatter={'float': '{: 0.6f}'.format})
print('Cost: {:0.7f}'.format(cost))
print('Expected cost: 2.534819')
print('Gradients:\n{}'.format(grad))
print('Expected gradients:\n[ 0.146561 -0.548558 0.724722 1.398003]')

帶有正則化的邏輯迴歸代價函式表示式：

梯度計算：

其中：

在lrCostFunction.py中編寫程式碼計算代價函式和梯度：

#sigmoid函式
def sigmoid(z):
    return 1 / (1 + np.exp(-z))

import numpy as np
from sigmoid import *


#邏輯迴歸的假設函式
def h(X,theta):
    return sigmoid(np.dot(X,theta)) #X: m*(n+1)  theta:(n+1,) 內積返回結果(m*1,)

#計算代價函式
def Compute_cost(theta, X, y, lmd):
    m = y.size  #樣本數

    cost = 0
    myh=h(X,theta)  #得到假設函式值
    
    term1=-y.dot(np.log(myh)) #y:(m*1,)   log(myh):(m*1,)  得到一個數值
    term2=(1-y).dot(np.log(1-myh))#1-y:(m*1,)   log(1-myh):(m*1,) 得到一個數值
    term3=theta[1:].dot(theta[1:])*lmd #正則化項 注意不懲罰theta0 得到一個數值 thrta[1:] (n,)
    cost=(1/m)*(term1-term2)+(1/(2*m))*term3
   
    return cost

#計算梯度值
def Compute_grad(theta,X,y,lmd):
    m = y.size  #樣本數

    grad = np.zeros(theta.shape) #梯度是與引數同維的向量
    myh=h(X,theta)  #得到假設函式值

    reg=(lmd/m)*theta[1:] #reg (n,)
    beta=myh-y    #beta: (m,)
    grad=beta.dot(X)/m #beta:(m,)  X:m*(n+1)  grad:(n+1,)
 
    grad[1:]+=reg
    return grad

def lr_cost_function(theta, X, y, lmd):
   
    cost=Compute_cost(theta, X, y, lmd)
    grad=Compute_grad(theta, X, y, lmd)
    return cost, grad

測試結果：

說明我們計算代價函式和梯度的程式碼正確。

• 訓練邏輯迴歸多分類器

'''第2-2部分  訓練多分類的邏輯迴歸 實現手寫數字識別'''
print('Training One-vs-All Logistic Regression ...')

lmd = 0.1 #正則化懲罰項係數
all_theta = ova.one_vs_all(X, y, num_labels, lmd)  #返回訓練好的引數

• 編寫邏輯迴歸多分類器的訓練程式one_vs_all.py

import scipy.optimize as opt #高階優化函式的包
import lrCostFunction as lCF
from sigmoid import *


#定義一個優化函式 實際上呼叫的是Python內建的高階優化函式 
#可以把它想象成梯度下降法 但是不用手動設定學習率
''' fmin_cg優化函式 第一個引數是計算代價的函式 第二個引數是計算梯度的函式 引數x0傳入初始化的theta值
    args傳入訓練樣本的輸入特徵矩陣X,對應的2分類新標籤y,正則化懲罰項係數lmd
    maxiter設定最大迭代優化次數
'''
def optimizeTheta(theta,X,y,lmd):
    res=opt.fmin_cg(lCF.Compute_cost,fprime=lCF.Compute_grad,x0=theta,\
                    args=(X,y,lmd),maxiter=50,disp=False,full_output=True)
    return res[0],res[1]

def one_vs_all(X, y, num_labels, lmd):
    
    (m, n) = X.shape #m為訓練樣本數  n為原始輸入特徵數

    '''
    邏輯迴歸多分類器的訓練過程：
    用邏輯迴歸做多分類 相當於做多次2分類 每一次把其中一個類別當作正類 其餘全是負類
    手寫數字識別是10分類 需要做十次2分類 
    比如：第一次把數字0當作正類 設定新的標籤為1  數字1-9為負類  設定新的標籤是0 進行2分類
         第一次把數字1當作正類 設定新的標籤為1  數字2-9和0為負類  設定新的標籤是0 進行2分類
         以此類推

    '''
    all_theta = np.zeros((num_labels, n + 1)) #存放十次2分類的 最優化引數
    initial_theta=np.zeros(n+1)   #每一次2分類的初始化引數值
    
    X = np.c_[np.ones(m), X]  #新增一列特徵 值為1

   
    for i in range(num_labels):
        print('Optimizing for handwritten number {}...'.format(i))
        iclass=i if i else 10 #數字0 屬於第十個類別
        logic_Y=np.array([1 if x==iclass else 0 for x in y]) #設定每次2分類的新標籤 
        itheta,imincost=optimizeTheta(initial_theta,X,logic_Y,lmd)
        all_theta[i,:]=itheta
    print('Done')

    return all_theta #返回十次2分類的 最優化引數

• 在訓練集上測試分類器的準確率

'''第3部分  在訓練集上測試 之前訓練的邏輯迴歸多分類器的準確率'''

pred = pova.predict_one_vs_all(all_theta, X) #分類器的預測類別

print('Training set accuracy: {}'.format(np.mean(pred == y)*100))  #與真實類別進行比較 得到準確率

• 編寫測試程式predictOneVsAll.py

def predict_one_vs_all(all_theta, X):
    m = X.shape[0]   #shape[0]返回2維陣列的行數   樣本數
  
    p = np.zeros(m)  #儲存分類器預測的類別標籤

    X = np.c_[np.ones(m), X]  #增加一列1 X:5000*401

    Y=lCF.h(X,all_theta.T)  #all_theta:10*401   Y:5000*10  每一行是每一個樣本屬於10個類別的概率

    p=np.argmax(Y,axis=1) #找出每一行最大概率所在的位置
    
    p[p==0]=10  #如果是數字0的話  他屬於的類別應該是10

    return p

• 分類器在訓練集上的準確率

4.邏輯迴歸多分類器的完整專案程式碼

下載連結     下載密碼:07nt

5.利用神經網路進行手寫數字識別

• 開啟神經網路分類器的主程式ex3nn.py
• 引入必要的包

#系統自帶包
import matplotlib.pyplot as plt
import numpy as np
import scipy.io as scio  #讀取.mat檔案

#手寫的包
import displayData as dd  #視覺化資料
import predict as pd     #神經網路分類器的預測程式

• 載入資料集並可視化

input_layer_size = 400  # 輸入層的單元數  原始輸入特徵數 20*20=400
hidden_layer_size = 25  # 隱藏層 25個神經元
num_labels = 10         # 10個標籤 數字0對應類別10  數字1-9對應類別1-9
                       


'''第1部分 載入資料集並可視化'''


print('Loading and Visualizing Data ...')

data = scio.loadmat('ex3data1.mat') #讀取資料
X = data['X']  #獲取輸入特徵矩陣 5000*400
y = data['y'].flatten()  #獲取5000個樣本的標籤 用flatten()函式 將5000*1的2維陣列 轉換成包含5000個元素的一維陣列
m = y.size  #樣本數 5000

# 隨機選100個樣本 視覺化
rand_indices = np.random.permutation(range(m))
selected = X[rand_indices[0:100], :]

dd.display_data(selected)

• 檢視視覺化程式displayData.py

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)
    display_cols = np.ceil(m / display_rows).astype(int)

    # 待顯示的每張圖片之間的間隔
    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')

• 視覺化效果

• 載入訓練好的神經網路引數

因為本週的課程只講了神經網路的前向傳播，下週講解神經網路的反向傳播演算法，來訓練得到一組最優的引數。所以提供了訓練好的引數，該實驗只需要完成神經網路前向傳播的預測過程，不涉及訓練過程。

'''第2部分 載入訓練好的神經網路引數'''


print('Loading Saved Neural Network Parameters ...')

data = scio.loadmat('ex3weights.mat') #讀取引數資料
#本實驗神經網路結構只有3層 輸入層，隱藏層 輸出層
theta1 = data['Theta1'] #輸入層和隱藏層之間的引數矩陣
theta2 = data['Theta2'] #隱藏層和輸出層之間的引數矩陣

• 實現預測過程

'''第3部分 利用訓練好的引數 完成神經網路的前向傳播 實現預測過程'''
pred = pd.predict(theta1, theta2, X)

print('Training set accuracy: {}'.format(np.mean(pred == y)*100))

• 編寫神經網路的前向傳播程式 predict.py

import numpy as np
from sigmoid import *

def predict(theta1, theta2, x):
    
    #theta1:25*401 輸入層多一個偏置項
    #theta2:10*26  隱藏層多一個偏置項
    m = x.shape[0]  #樣本數
    num_labels = theta2.shape[0]  #類別數

    x=np.c_[np.ones(m),x] #增加一列1   x：5000*401
    p = np.zeros(m)
    z1=x.dot(theta1.T)  #z1:5000*25
    a1=sigmoid(z1)    #a1:5000*25
    a1=np.c_[np.ones(m),a1] #增加一列1 a1:5000*26
    z2=a1.dot(theta2.T) #z2:5000*10
    a2=sigmoid(z2)   #a2:5000*10
    
    p=np.argmax(a2,axis=1) #輸出層的10個單元 第一個對應數字1...第十個對應數字0
 
    p+=1  #最大位置+1 即為預測的標籤
    return p

• 在訓練集上的預測準確率

6.神經網路分類器的完整專案程式碼

下載連結 下載密碼:yi6n