DeepLearning:五、BP神经网络in python

下面这段代码来自IMPLEMENTING A NEURAL NETWORK FROM SCRATCH IN PYTHON – AN INTRODUCTION中的代码，将上一章节中的BP神经网络实现出来。代码并没有考虑效率，但是代码很容易理解整个BP过程，这里我将我的理解记录下来并分享出来，加深印象。
目录

• BP神经网路回顾
• BP神经网络in python
  
  • 数据
  • 神经网路结构
  • 神经网络模型
  • 预测
  • 结果可视化
• 这里使用的所有代码

BP神经网络回顾

我们先看一下单个神经元的结构：

我们在看一下一个这里需要构建的神经网络的结构：

再回顾一下BP误差传播算法：

用数学的形式表示就是：

z1a1z2a2=xW1+b1=tanh(z1)=a1W2+b2=y^=softmax(z2)

$$\begin{aligned} z_1 & = xW_1 + b_1 \\ a_1 & = \tanh(z_1) \\ z_2 & = a_1W_2 + b_2 \\ a_2 & = \hat{y} = \mathrm{softmax}(z_2) \end{aligned}$$
其中$x$表示输入向量，$z_i$是$i$层的输入，$a_i$是应用激活函数的输出。$W_1,b_1,W_2,b_2$是这个网络的模型参数。应用BP公式求导（这是是人工求导的，在theano中是可以自动求导）得：
$$\begin{aligned} & \delta_3 = \hat{y} - y \\ & \delta_2 = (1 - \tanh^2 z_1) \circ \delta_3W_2^T \\ & \frac{\partial{L}}{\partial{W_2}} = a_1^T \delta_3 \\ & \frac{\partial{L}}{\partial{b_2}} = \delta_3\\ & \frac{\partial{L}}{\partial{W_1}} = x^T \delta2\\ & \frac{\partial{L}}{\partial{b_1}} = \delta2 \\ \end{aligned}$$

BP算法代码in python

数据

def generate_data():
    np.random.seed(0)
    X, y = datasets.make_moons(200,noise=0.20)
    return X, y

这里面使用的是scikit-learn中的函数，产生200个数据，下图是这个图的可视化结果：

神经网络结构

class Config:
    nn_input_dim = 2  # 输入的维度
    nn_output_dim = 2  # 输出维度
    # 梯度下降参数
    epsilon = 0.01  # 学习率
    reg_lambda = 0.01  # 正则化长度

如网络图所示，这里的网络有2维输入，2维输出，学习步长是0.01（如果这里不懂，可以看一下梯度下降法的原理），正则化长度。

神经网络模型

# 这个function是学习神经网络的参数以及建立模型
# - nn_hdim: 隐藏层的节点数
# - num_passes: 梯度下降法使用的样本数量
def build_model(X, y, nn_hdim, num_passes=20000):
    # Initialize the parameters to random values. We need to learn these.
    num_examples = len(X)
    np.random.seed(0)
    W1 = np.random.randn(Config.nn_input_dim, nn_hdim) / np.sqrt(Config.nn_input_dim)
    b1 = np.zeros((1, nn_hdim))
    W2 = np.random.randn(nn_hdim, Config.nn_output_dim) / np.sqrt(nn_hdim)
    b2 = np.zeros((1, Config.nn_output_dim))

    # 最后返回的模型，主要就是每一层的参数向量
    model = {}

    # 梯度下降法
    for i in range(0, num_passes):

        # 正向传播过程
        z1 = X.dot(W1) + b1
        a1 = np.tanh(z1)
        z2 = a1.dot(W2) + b2
        exp_scores = np.exp(z2)
        probs = exp_scores / np.sum(exp_scores, axis=1, keepdims=True)

        # 误差反向传播过程
        delta3 = probs
        delta3[range(num_examples), y] -= 1
        dW2 = (a1.T).dot(delta3)
        db2 = np.sum(delta3, axis=0, keepdims=True)
        delta2 = delta3.dot(W2.T) * (1 - np.power(a1, 2))
        dW1 = np.dot(X.T, delta2)
        db1 = np.sum(delta2, axis=0)

        # 添加正则项 (b1 and b2 不需要做正则化)
        dW2 += Config.reg_lambda * W2
        dW1 += Config.reg_lambda * W1

        # 梯度下降参数更新
        W1 += -Config.epsilon * dW1
        b1 += -Config.epsilon * db1
        W2 += -Config.epsilon * dW2
        b2 += -Config.epsilon * db2

        # 更新模型参数
        model = {'W1': W1, 'b1': b1, 'W2': W2, 'b2': b2}

    return model

前向传播过程数学表达：

z1a1z2a2=xW1+b1=tanh(z1)=a1W2+b2=y^=softmax(z2)

$$\begin{aligned} z_1 & = xW_1 + b_1 \\ a_1 & = \tanh(z_1) \\ z_2 & = a_1W_2 + b_2 \\ a_2 & = \hat{y} = \mathrm{softmax}(z_2) \end{aligned}$$
误差反向传播数学表达：
$$\begin{aligned} & \delta_3 = \hat{y} - y \\ & \delta_2 = (1 - \tanh^2 z_1) \circ \delta_3W_2^T \\ & \frac{\partial{L}}{\partial{W_2}} = a_1^T \delta_3 \\ & \frac{\partial{L}}{\partial{b_2}} = \delta_3\\ & \frac{\partial{L}}{\partial{W_1}} = x^T \delta2\\ & \frac{\partial{L}}{\partial{b_1}} = \delta2 \\ \end{aligned}$$
梯度下降更新参数：

这里是整个神经网络的核心的地方，代码的解释我个人觉得使用数学公式表示的更加明了清晰。这段代码中还有将loss值打印出来的函数，我将这段代码删除了，方便代码更加紧凑，更加容易理解。

预测

def predict(model, x):
    W1, b1, W2, b2 = model['W1'], model['b1'], model['W2'], model['b2']
    # Forward propagation
    z1 = x.dot(W1) + b1
    a1 = np.tanh(z1)
    z2 = a1.dot(W2) + b2
    exp_scores = np.exp(z2)
    probs = exp_scores / np.sum(exp_scores, axis=1, keepdims=True)
    return np.argmax(probs, axis=1)

这里就是一个前向传播过程：

z1a1z2a2=xW1+b1=tanh(z1)=a1W2+b2=y^=softmax(z2)

$$\begin{aligned} z_1 & = xW_1 + b_1 \\ a_1 & = \tanh(z_1) \\ z_2 & = a_1W_2 + b_2 \\ a_2 & = \hat{y} = \mathrm{softmax}(z_2) \end{aligned}$$

结果可视化

def visualize(X, y, model):
    # plt.scatter(X[:, 0], X[:, 1], s=40, c=y, cmap=plt.cm.Spectral)
    # plt.show()
    plot_decision_boundary(lambda x:predict(model,x), X, y)
    plt.title("Logistic Regression")


def plot_decision_boundary(pred_func, X, y):
    # 设置最小最大值并填充
    x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
    y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5
    h = 0.01
    # 生成数据网格
    xx, yy = np.meshgrid(np.arange(x_min, x_max, h), np.arange(y_min, y_max, h))
    # 预测整个数据网格上的数据
    Z = pred_func(np.c_[xx.ravel(), yy.ravel()])
    Z = Z.reshape(xx.shape)
    # 绘制数据点以及边界
    plt.contourf(xx, yy, Z, cmap=plt.cm.Spectral)
    plt.scatter(X[:, 0], X[:, 1], c=y, cmap=plt.cm.Spectral)
    plt.show()

将数据点以及边界绘制出来。

运行

def main():
    X, y = generate_data()
    model = build_model(X, y, 3, print_loss=True)
    visualize(X, y, model)


if __name__ == "__main__":
    main()

运行结果

这里使用到的所有代码

下面这是这里使用到的所有代码，为了紧凑，删除了一些代码，如果你想看源代码看这里。想看原文看这里

import numpy as np
from sklearn import datasets, linear_model
import matplotlib.pyplot as plt


class Config:
    nn_input_dim = 2 
    nn_output_dim = 2  
    epsilon = 0.01  
    reg_lambda = 0.01  


def generate_data():
    np.random.seed(0)
    X, y = datasets.make_moons(200, noise=0.20)
    return X, y


def visualize(X, y, model):
    plot_decision_boundary(lambda x:predict(model,x), X, y)
    plt.title("Logistic Regression")


def plot_decision_boundary(pred_func, X, y):
    x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
    y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5
    h = 0.01
    xx, yy = np.meshgrid(np.arange(x_min, x_max, h), np.arange(y_min, y_max, h))
    Z = pred_func(np.c_[xx.ravel(), yy.ravel()])
    Z = Z.reshape(xx.shape)
    plt.contourf(xx, yy, Z, cmap=plt.cm.Spectral)
    plt.scatter(X[:, 0], X[:, 1], c=y, cmap=plt.cm.Spectral)
    plt.show()

def predict(model, x):
    W1, b1, W2, b2 = model['W1'], model['b1'], model['W2'], model['b2']
    z1 = x.dot(W1) + b1
    a1 = np.tanh(z1)
    z2 = a1.dot(W2) + b2
    exp_scores = np.exp(z2)
    probs = exp_scores / np.sum(exp_scores, axis=1, keepdims=True)
    return np.argmax(probs, axis=1)



def build_model(X, y, nn_hdim, num_passes=20000, print_loss=False):
    num_examples = len(X)
    np.random.seed(0)
    W1 = np.random.randn(Config.nn_input_dim, nn_hdim) / np.sqrt(Config.nn_input_dim)
    b1 = np.zeros((1, nn_hdim))
    W2 = np.random.randn(nn_hdim, Config.nn_output_dim) / np.sqrt(nn_hdim)
    b2 = np.zeros((1, Config.nn_output_dim))

    model = {}

    for i in range(0, num_passes):

        z1 = X.dot(W1) + b1
        a1 = np.tanh(z1)
        z2 = a1.dot(W2) + b2
        exp_scores = np.exp(z2)
        probs = exp_scores / np.sum(exp_scores, axis=1, keepdims=True)

        delta3 = probs
        delta3[range(num_examples), y] -= 1
        dW2 = (a1.T).dot(delta3)
        db2 = np.sum(delta3, axis=0, keepdims=True)
        delta2 = delta3.dot(W2.T) * (1 - np.power(a1, 2))
        dW1 = np.dot(X.T, delta2)
        db1 = np.sum(delta2, axis=0)

        dW2 += Config.reg_lambda * W2
        dW1 += Config.reg_lambda * W1

        W1 += -Config.epsilon * dW1
        b1 += -Config.epsilon * db1
        W2 += -Config.epsilon * dW2
        b2 += -Config.epsilon * db2

        model = {'W1': W1, 'b1': b1, 'W2': W2, 'b2': b2}

    return model


def main():
    X, y = generate_data()
    model = build_model(X, y, 3)
    visualize(X, y, model)


if __name__ == "__main__":
    main()