从感知机到神经网络：Python实现与测试

本文为吴恩达神经网络课程打的学习笔记，包含了本人自己的完整推导过程，并且为便于理解，将代码进行了部分重构。

简介

神经网络算法可以算是一种仿生学，其基本数据处理单元——感知机——模仿的是生物神经系统内的神经元，它能够接受来自多个源的信号输入，然后将信号转化为便于传播的信号在进行输出(在生物体内表现为电信号)。

感知机

感知机结构：

{z=wT∗x+bhw(x)=step(z)

其中z为神经元将接收到的信号进行整合，step()函数将整合后的信号进行转化并输出。

早期常用到的step函数有：

heaviside(x)={01x<0x≥0

sgn(x)=⎧⎩⎨−101x<0x=0x>0

多层感知机(Multi-Layer Perceptron)

多层感知机由输入层、输出层和隐含层组成：

实际上在机器学习发展的早期，感知机模型只能用来解决一些鸡毛蒜皮的问题，主要原因就是受step()函数的限制，因为早期的step()函数只能产生离散的输出，并且掩盖了原信号(即无法从感知机的输出来推断原信号)。下图为感知机模型解决异或问题的方案：

反向传播算法的出现

在很长一段时间里人们都没有找到训练MLP的方法，直到1986年D. E. Rumelhart等人发明了反向传播算法，MLP才具有学习能力。

反向传播算法具体过程是这样的：

1. 对于一个训练用例，网络将其作为输出并计算输出，即预测值
2. 测量网络的输出误差
3. 计算在隐含层的最后一层中哪些单元对误差产生了贡献
4. 再从那些产生误差的单元往回追溯，找出隐含层中所有导致误差的单元，直到回溯到输入层为止
5. 调节那些单元的权重分布

为了使这个算法正常工作，作者还对MLP架构做出了一个至关重要的改进：将step函数替换为逻辑函数σ(x)=11+e(−x)，这个函数不仅是全局可导的，而且具有连续输出。

目前常用的激活函数有三种：

⎧⎩⎨⎪⎪σ(x)=11+e(−x)tanh(x)=2σ(2x)−1ReLU(x)=max(0,x)

下面是它们的图象与导数：

单感知机模型的实现

外部库

import numpy as npimport matplotlib.pyplot as pltimport h5pyimport scipyfrom PIL import Imagefrom scipy import ndimage%matplotlib inlinenp.random.seed(1)

引入数据

所用数据为h5格式的关于猫的图片数据，每张图片像素为64*64，RGB通道。其中猫的图片标记为y=1，不是猫的图片标记为y=0。

def load_dataset():    train_dataset = h5py.File('datasets/train_catvnoncat.h5', "r")    train_set_x_orig = np.array(train_dataset["train_set_x"][:]) # your train set features    train_set_y_orig = np.array(train_dataset["train_set_y"][:]) # your train set labels    test_dataset = h5py.File('datasets/test_catvnoncat.h5', "r")    test_set_x_orig = np.array(test_dataset["test_set_x"][:]) # your test set features    test_set_y_orig = np.array(test_dataset["test_set_y"][:]) # your test set labels    classes = np.array(test_dataset["list_classes"][:]) # the list of classes    train_set_y_orig = train_set_y_orig.reshape((1, train_set_y_orig.shape[0]))    test_set_y_orig = test_set_y_orig.reshape((1, test_set_y_orig.shape[0]))    return train_set_x_orig, train_set_y_orig, test_set_x_orig, test_set_y_orig, classes

train_set_x_orig, train_set_y, test_set_x_orig, test_set_y, classes = load_dataset()print("X_train shape:\t{}\nY_train shape:\t{}\nX_test shape:\t{}\nY_test shape:\t{}\n".format      (train_set_x_orig.shape,train_set_y.shape,test_set_x_orig.shape,test_set_y.shape))

输出为：
X_train shape: (209, 64, 64, 3)
Y_train shape: (1, 209)
X_test shape: (50, 64, 64, 3)
Y_test shape: (1, 50)

index = 30        #更改index值以查看不同样本plt.imshow(train_set_x_orig[index])print ("y = " + str(train_set_y[0,index]))

输出为：

m_train=train_set_x_orig.shape[0]m_test=test_set_x_orig.shape[0]num_px=train_set_x_orig.shape[1]print("m_train:{}\tm_test:{}\tnum_px:{}".format(m_train,m_test,num_px))

输出为：m_train:209 m_test:50 num_px:64

数据的矩阵化

在多样本输入条件下，对于每一个样本

xi=⎡⎣⎢⎢⎢feature1feature2...featuren⎤⎦⎥⎥⎥

，都需进行计算wT∗x+b，最简单的方法就是使用循环让模型针对每一个样本进行训练，这样做的时间复杂度为O(样本数∗每个样本的训练次数)。

由于每个样本的计算之间没有依赖关系，可以将数据矩阵化实现并行计算：
X=[x1x2…xm]
，其中每个样本xi均为列向量。

感知机的最终输出
A=[σ(z1)σ(z2)…σ(zm)]，其中zi表示样本转换后的数字信号。

train_set_x_flatten = train_set_x_orig.reshape(m_train,num_px*num_px*3).Ttest_set_x_flatten = test_set_x_orig.reshape(m_test,num_px*num_px*3).Ttrain_set_x = train_set_x_flatten/255test_set_x = test_set_x_flatten/255print("X_train_flatten shape:\t{}\nX_test_flatten shape:\t{}\n".format(train_set_x.shape,test_set_x.shape))

输出为：
X_train_flatten shape: (12288, 209)
X_test_flatten shape: (12288, 50)

实现激活函数

σ(x)=11+e(−x)

def sigmoid(z):    s = 1/(1+np.exp(-z))    return s

参数矩阵初始化

随即初始化权重矩阵w，而偏置常量b设为1。

#参数：输入层的节点数dim,同样本特征数#返回：随机初始化的参数字典parametersdef initialize_parameters(dim):    parameters={        "w":np.random.randn(dim,1)*0.01,        "b":1    }    return parameters

前向传播

前向传播就是针对输入给出相应的输出并计算误差的过程，此处使用交叉熵L(a,y)=−ylna−(1−y)ln(1−a)来评估误差：

⎧⎩⎨⎪⎪⎪⎪⎪⎪Z=WT∗X+bA=σ(Z)Loss=1m∑i=1mL(ai,y)

#参数：包含W,b的参数字典parameters，样本数据X，样本标签Y#返回：损失值cost，包含A的缓存cachedef forward_propagation(parameters, X, Y):    w=parameters['w']    b=parameters['b']    m = X.shape[1]    Z=np.dot(w.T,X)+b    A = sigmoid(Z)    cache=A    cost = ((-1/m)*(Y*np.log(A)+(1-Y)*np.log(1-A))).sum(axis=1)    cost = np.squeeze(cost)    assert(cost.shape == ())    return cost,cache

反向传播

反向传播算法实际上是一个优化问题，它需要找到一组能最小化Loss的解(W,b)。最简单的方法就是求导，直接计算满足一阶导为零且二阶导大于或小于零的参数即可，实际上在线性模型中就有这么一种方法称为正规方程。

另外一种逼近解的方法叫梯度下降法，对于某参数θ的计算过程如下：

θnext_step=θcur_step−α∂Loss∂θ

其中α称为学习率。

具体对于参数w跟b的计算为：

∂L∂A=1m∗(−YA+1−Y1−A)

∂L∂Z===∂L∂A∗∂A∂Z1m∗(−yA+1−Y1−A)∗A∗(1−A)1m(A−Y)

∂L∂W==∂L∂Z∗∂Z∂W1m(A−Y)∗XT

∂L∂b==∂L∂Z∗∂Z∂b1m(A−Y)

我们只需要∂L∂W=1mX∗(A−Y)T和∂L∂b=1m(A−Y)即可。为了减少计算量，可将前向传播时计算过的A缓存起来。

#参数：样本数据X，样本标签Y，包含A的缓存cache#返回：包含dw,db的梯度字典gradsdef backward_propagation(X,Y,cache):    m=X.shape[1]    A=cache    grads = {"dw": np.dot(X,(A-Y).T)/m,             "db": ((A-Y)/m).sum(axis=1)}    return grads

参数的迭代优化

对于给定的学习率与迭代次数，对参数W，b进行优化。optimize函数返回优化好的参数与记录下的损失值cost以供可视化。

θnext_step=θcur_step−α∂Loss∂θ

def optimize(parameters, X, Y, num_iterations, learning_rate):    costs = []    w=parameters['w']    b=parameters['b']    for i in range(num_iterations):        cost,cache=forward_propagation(parameters,X,Y)        grads=backward_propagation(X,Y,cache)        dw = grads["dw"]        db = grads["db"]        w = w-learning_rate*dw        b = b-learning_rate*db        parameters={            'w':w,            'b':b        }        #每一百次迭代记录一次cost         if i % 100 == 0:            costs.append(cost)    params = {"w": w,              "b": b}    return params, costs

预测函数

预测函数以给定的参数(W,b)与样本X为参数，会返回给定参数下的模型对样本X的预测值。

def predict(parameters, X):    w = parameters["w"]    b = parameters["b"]    m = X.shape[1]    Y_prediction = np.zeros((1,m))    Z=np.dot(w.T,X)+b    A = sigmoid(Z)    for i in range(A.shape[1]):        Y_prediction[0, i] = 1 if A[0,i]>0.5 else 0    #注意取矩阵元素与取数组元素的区别    return Y_prediction

评估模型

可视化模型在训练过程中的学习曲线。

def show_info(parameters,X_train,Y_train,X_test,Y_test,costs,learning_rate):    Y_prediction_test = predict(parameters, X_test)    Y_prediction_train = predict(parameters, X_train)    print("train accuracy: "  + str(np.sum((Y_prediction_train == Y_train)/X_train.shape[1])))    print("test accuracy: "  + str(np.sum((Y_prediction_test == Y_test)/X_test.shape[1])))    plt.plot(np.squeeze(costs))    plt.ylabel('cost')    plt.xlabel('iterations (per tens)')    plt.title("Learning rate =" + str(learning_rate))    plt.show()

整合模型

def model(X_train, Y_train, X_test, Y_test, num_iterations = 2000, learning_rate = 0.005):    #1.参数初始化    parameters = initialize_parameters(X_train.shape[0])    #2.迭代优化参数    parameters, costs = optimize(parameters, X_train, Y_train, num_iterations, learning_rate)    show_info(parameters,X_train,Y_train,X_test,Y_test,costs,learning_rate)    return parameters

训练模型

parameters = model(train_set_x, train_set_y, test_set_x, test_set_y, num_iterations = 2000, learning_rate = 0.005)

输出为：

多层网络实现

网络模型

此处实现一个L层的神经网络，最后一层的激活函数使用σ(x)，其余层均使用ReLU(x)。整个模型框图如下图所示：

参数矩阵维度分析

对于L层的网络模型，易得每一层(第i层)都有这么几个参数(i>0)：

• 神经元数，n[i]
• 权重参数，W[i]
• 偏置参数，b
• 线性输出，Z[i]=W[i]∗A[i−1]+b（注意，此处不再使用WT∗X）
• 激活函数，
  g[i](x)={σ(x)ReLU(x)i=Li<L
• 输出，A[i]=g(Z[i])

下面简单推导一下每一层的参数矩阵的维度。

令每个样本均为列向量，则易得

dim(A[i])=(n[i],m)

其中m为样本数量。

然后从输入层开始，输入层的矩阵维度为dim(A[0])=(n[0],m)，而

dim(W[1]∗A[0])=dim(A[1])=(n[1],m)

所以dim(W[1])=(n[1],n[0])；继续沿层数推下去会发现一个规律，得出

dim(W[i])=(n[i],n[i−1])

参数矩阵初始化

#三层模型layers_dims = [12288, 7,1]         #输入层，隐含层，输出层

def initialize_parameters(layer_dims):    parameters = {}    L = len(layer_dims)            # number of layers in the network    for i in range(1, L):        parameters['W' + str(i)] = np.random.randn(layer_dims[i],layer_dims[i-1])*0.01        parameters['b' + str(i)] = np.ones((layer_dims[i],1))        assert(parameters['W' + str(i)].shape == (layer_dims[i], layer_dims[i-1]))        assert(parameters['b' + str(i)].shape == (layer_dims[i], 1))    return parameters

激活函数

def sigmoid(Z):    A = 1/(1+np.exp(-Z))    return Adef relu(Z):    A = np.maximum(0,Z)    assert(A.shape == Z.shape)    return A

每一层神经元的运算

完成每一层所有神经元的运算，支持两种激活函数。

#参数：上一层的输出A_prev，当前层的参数W,b，激活方式activation#返回：当前层的输出A，缓存A_prev,W,Zdef perceptron_forward(A_prev, W, b, activation):    Z=np.dot(W,A_prev)    if activation == "sigmoid":        A = sigmoid(Z)    elif activation == "relu":        A = relu(Z)    assert (A.shape == (W.shape[0], A_prev.shape[1]))    cache=(A_prev,W,Z)    return A, cache

前向传播(L层的神经元运算)

def forward_propagation(X,Y, parameters):    caches = []    A = X        #第零层的输出就为X    L = len(parameters)//2                  # number of layers in the neural network    for l in range(1, L):        A_prev = A        A, cache = perceptron_forward(A_prev,parameters["W"+str(l)],parameters["b"+str(l)],"relu")        caches.append(cache)    AL, cache = perceptron_forward(A,parameters["W"+str(L)],parameters["b"+str(L)],"sigmoid")    caches.append(cache)    assert(AL.shape == (1,X.shape[1]))    m = Y.shape[1]    cost = -np.sum(np.multiply(Y,np.log(AL))+np.multiply(1-Y,np.log(1-AL)))/m    return AL, caches,cost

反向传播推导

首先从最后一层开始，

dA[L]==∂L∂A[L]1m∗(−YA[L]+1−Y1−A[L])

dZ[L]===∂L∂Z[L]dA[L]∗∂A[L]∂Z[L]dA[L]∗A[L]∗(1−A[L])

dW[L]===∂L∂W[L]dZ[L]∗∂Z[L]∂W[L]dZ[L]∗A[L−1].T

db[L]=1mdZ[L]

于是对于使用\sigma(x)的最后一层，得到了dW[L]=dZ[L]∗A[L−1].T和db[L]=1mdZ[L]，继续往前推导。

dA[L−1]===∂L∂A[L−1]dZ[L]∗∂Z[L]∂A[L−1]W[L].T∗dZ[L]

dZ[L−1]===∂L∂Z[L−1]dA[L−1]∗∂A[L−1]∂Z[L−1](0,dA[L−1])

dW[L−1]===∂L∂W[L−1]dZ[L−1]∗∂Z[L−1]∂W[L−1]dZ[L−1]∗A[L−2].T

db[L−1]=1mdZ[L−1]

可得对于使用ReLU(x)的L-1层，每层的梯度为dW[i]=dZ[i]∗A[i−1].T，db[i]=1mdZ[i]。

以上结论可变换为：

⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪dA[L]dZ[L]dW[i]db[i]dA[i−1]dZ[i−1]=1m∗(−YA[L]+1−Y1−A[L])=dA[L]∗A[L]∗(1−A[L])=dZ[i]∗A[i−1].T=1mdZ[i]=W[i].T∗dZ[i]=(0,dA[i−1])(1)(2)(3)(4)(5)(6)

在多层神经网络的反向传播过程中，为快速求出第i层的dW[i]与db[i]，需要缓存前向传播过程中的A[i−1]、W[i]与Z[i]，还有反向传播过程中的dA[i]。

每一层神经元的运算

#参数：当前层的dA，包含A_prev,W,Z的cache，激活方式activation#返回：上一层的梯度dA_prev，当前层的梯度dW,dbdef perceptron_backward(dA, cache, activation):    A_prev,W,Z=cache    if activation == "relu":        dZ = np.array(dA, copy=True)        dZ[Z <= 0] = 0    elif activation == "sigmoid":        s = 1/(1+np.exp(-Z))        dZ = dA * s * (1-s)    assert (dZ.shape == Z.shape)    m = A_prev.shape[1]    dW = np.dot(dZ,A_prev.T)    db = np.sum(dZ,axis=1,keepdims=True)/m    dA_prev = np.dot(W.T,dZ)    return dA_prev, dW, db

反向传播

def backward_propagation(AL, Y, caches):    grads = {}    L = len(caches) # the number of layers    m = AL.shape[1]    dAL = - (np.divide(Y, AL) - np.divide(1 - Y, 1 - AL))/m    current_cache = caches[L-1]        #取L层的缓存数据    grads["dA" + str(L-1)], grads["dW" + str(L)], grads["db" + str(L)] = perceptron_backward(dAL,current_cache,"sigmoid")    for l in reversed(range(L - 1)):        current_cache = caches[l]        dA_prev, dW, db = perceptron_backward(grads["dA"+str(l+1)],current_cache,"relu")        grads["dA" + str(l)] = dA_prev        grads["dW" + str(l + 1)] = dW        grads["db" + str(l + 1)] = db    return grads

迭代优化

def update_parameters(parameters, grads, learning_rate):    L = len(parameters) // 2 # number of layers in the neural network    for l in range(L):        parameters["W" + str(l+1)] = parameters["W"+str(l+1)]-learning_rate*grads["dW"+str(l+1)]        parameters["b" + str(l+1)] = parameters["b"+str(l+1)]-learning_rate*grads["db"+str(l+1)]    return parameters

预测函数

若\sigma(x)的输出>0.5则判为正例，否则判为反例。

def predict(X, y, parameters):    m = X.shape[1]    n = len(parameters) // 2 # number of layers in the neural network    p = np.zeros((1,m))    probas, _ , _ = forward_propagation(X,y, parameters)    for i in range(0, probas.shape[1]):        if probas[0,i] > 0.5:            p[0,i]=1 if probas[0,i]>0.5 else 0    print("Accuracy: "  + str(np.sum((p == y)/m)))    return p

整合模型

def L_layer_model(X, Y, layers_dims, learning_rate = 0.005, num_iterations = 2000, print_cost=False):    costs = []    parameters = initialize_parameters(layers_dims)    for i in range(0, num_iterations):        AL, caches,cost = forward_propagation(X,Y, parameters)        grads = backward_propagation(AL, Y, caches)        parameters = update_parameters(parameters, grads, learning_rate)        if i % 100 == 0:            costs.append(cost)    # plot the cost    plt.plot(np.squeeze(costs))    plt.ylabel('cost')    plt.xlabel('iterations (per tens)')    plt.title("Learning rate =" + str(learning_rate))    plt.show()    return parameters

训练模型

parameters = L_layer_model(train_set_x, train_set_y, layers_dims, learning_rate=0.005,num_iterations = 2000, print_cost = True)

输出为：

评估模型

predict(train_set_x,train_set_y,parameters)predict(test_set_x,test_set_y,parameters)

输出为：
Accuracy: 0.99043062201
Accuracy: 0.76

具有隐含层的神经网络比单感知机要提高了六个百分点的准确率，不过两个模型都存在明显的过拟合现象。