神经网络结构、学习与实现

神经网络是人工智能和统计学习等领域研究的热点,其能将样本特征映射至高维空间进行分类/回归。本文梳理了神经网络的结构及其学习算法,并在Python下进行实现及测试。

神经网络的结构与学习

神经元与网络结构

神经网络(Neural Network, NN)是一种模仿生物神经网络(动物的中枢神经系统,特别是大脑)的结构和功能的数学模型或计算模型,用于对函数进行估计或近似 [1]。神经网络的输入是单个样本的特征值向量,输出也是一个向量。神经网络由多个神经元组成。单个神经元的结构如下图:

neuron

单个神经元的输入是多维向量,输出为一个数。具体的,对于每个神经元,其对应的计算过程为:先对多维向量赋权重求和得到$net_j$值,计算为$ net_j = \sum w_i x_i + b_j$ ,其中$w_i$为输入信号的权重参数,$b_j$为偏差量。 再对$net_j$值做非线性变换得到输出值,输出值计算为 $z_j = \sigma (net_j) $,其中$ \sigma(\cdot) $是一个非线性变换函数,其导数通常可以通过原始函数值进行计算。常用的非线性变换函数有:

  • Logistic-Sigmoid函数:$\sigma(x) = \frac{1}{1+e^{-x}}$,其导数为 $\sigma ‘(x) = \sigma(x) (1- \sigma(x))$
  • Tanh-Sigmoid函数:$tanh(x) = \frac{e^x - e^{-x}}{e^x + e^{-x}}$,其导数为 $tanh’(x) = 1 - tanh^2(x)$
  • Softplus函数:$f(x) = log(1+e^x)$,其导数为 $f’(x) = \frac{e^x}{1+e^x} = \frac{1}{1+e^{-x}} = \sigma(x)$
  • ReLu函数:$f(x) = max(0, x)$,其导数在$x>0$时为$1$,在$x<=0$时为$0$

神经网络通常为按层级排列的多层结构组成,例如下图的三层神经网络,三层分别为输入层、隐藏层和输出层。输出层每个神经元对应一个目标值。例如在$c-$类分类问题中,若样本$x$的真实类别为第$k$个类别,则目标输出为 $t_i = 0 \ if \ i != k \ else \ t_i = 1$ 。

NN_structure

神经网络输出

神经网络输出过程会按网络的拓扑顺序算出每个节点值。下面以使用三层神经网络计算不同类别后验概率值为例进行说明。首先,若希望神经网络最后的输出值 $z_k$ 为对应第$k$个类别的概率值,则可对输出层的输出值进行Softmax函数处理,$Softmax(z_k) = \frac{e^{z_k}}{\sum e^{z_j}} = p_k$,其中$z_k$为输出层第$k$个节点的输出值。现在从输入层开始进行对每个节点关联的值进行计算:在已知输入向量$x$及初始化了各层权重的情况下(记第$h$层(有$i$个节点)到$h+1$层(有$j$个节点)之间的权重记为$W^{(h+1)}$,其为一个$j \times i$的矩阵),各层节点值的计算顺序为:$net^{(2)} -> z^{(2)} -> net^{(3)} -> z^{(3)} -> p^{(3)}$,其中$net^{(2)}$表示值向量,即上图中第二层(隐藏层)所有全体$net_j^{(2)}$组成的向量,$z_j = f(net_j)$,$p_k = softmax(z_k)$。将此计算过程向量化,则有:

  • $net^{(2)} = W^{(1)}x$
  • $z^{(2)} = f(net^{(2)})$,
  • $net^{(3)} = W^{(2)}z^{(2)}$
  • $z^{(3)} = f(net^{(3)})$
  • $p^* = exp(z^{(2)})$
  • $p = \frac{p}{sum(p)}$

总的来说上述迭代计算的过程即对应神经网络代表的函数,其输入为$x$,函数中有参数$W^{(h)}$,输出为$p$。只要给定$W^{(h)}$及$x$,以上所有值$net^{(\cdot )}, z^{(\cdot )}, p$ 都是确定的。

NN学习

参数的学习过程是最优化的过程,其过程与感知器的学习相似,即对参数定义损失函数,在搜索空间内找最优的参数值使得损失函数值最小。

定义损失函数与更新规则

损失函数是模型参数$W$的函数,其作用是衡量出在特定参数下,NN最终输出向量$Z$与目标向量 $T$之间的差异。对于回归问题,可以设置单个样本的损失函数为$J(W) = 1/2 \sum (t_k - z_k)^2 = 1/2 (T-Z)^2$;对于分类问题,由于NN的输出值是样本对于每个类别的后验概率$P(c_k | x)$,其目标向量$T$为样本的真实概率 $P(c_k|x) = 1 \ if \ x\in c_k\ ,else P(c_k |x)=0$。衡量两个概率向量差异的一种方式是使用交叉熵函数,对于预测概率$P$和真实概率$T$($\sum_i ti = 1$ ) ,交叉熵计算为$J= -\sum_i t_i log (p_i) $,其中$p_i = softmax(z_i) = \frac{e^{z_i}}{\sum_j e^{z_j}}$。为了方便下面的说明,下面先分别将softmax函数和交叉熵函数和对输入$ a_i$进行求导。

  • Softmax 函数对$a_k$求导有:$\partial p_i / \partial a_k = p_i(1-p_i), i = k $ ,和$\partial p_i / \partial a_k = -p_i p_k, i \neq k $

  • 交叉熵对$a_k$求导有:
    $$
    \frac{\partial J}{\partial z_k} = \sum_u \frac{\partial J}{\partial p_u} \frac{\partial p_u}{\partial z_k } = \sum_u (- \frac{t_u }{p_u} \frac{\partial p_u }{\partial z_k}) = - \frac{t_k }{p_k} \frac{\partial p_k }{\partial z_k} -\sum_o \frac{t_i }{p_i} \frac{\partial p_i }{\partial z_k}
    $$
    其中$o = { i \neq k }$,进而$\frac{\partial J}{\partial z_k} = - \frac{t_k }{p_k} p_k(1-p_k) + \sum_o \frac{t_i }{p_i} p_ip_k = \frac{\partial J}{\partial z_k} = -t_k+ \sum_i t_i p_k = p_k -t_k $。

有了单个样本的函数就可以求权重对样本集的损失值,例如将所有样本的损失值求和或求平均值作为样本集的损失值。在分类问题中,可以设置样本集的损失函数为$J(W) = \sum_k^m J(W, x_k) = \sum_k \sum_i -t_i log (p_i) $,其中$m$为样本集容量。有了损失值就可使用梯度下降法对参数进行迭代求解。下面是权重参数的更新规则(对单个样本而言):

$W(n+1) = W(n) - \eta \frac{\partial J}{\partial W(n)} $,其中$\eta$ 是每次更新的步长,记 $\Delta W(n) = \frac{\partial J}{\partial W(n)}$,则更新规则为$W(n+1) = W(n) - \eta\Delta W(n)$。同理,偏移量$b$的修改规则为$b(n+1) = b(n) - \eta\Delta b(n)$,其中$\Delta b(n) = \frac{\partial J}{\partial b(n)}$

对样本集而言,根据$J(W) = \sum_k^m J(W, x_k) $,这里需要将$\Delta W(n)$修改为$\sum_k^m \Delta W_k(n)$ ,及$\Delta b(n)$修改为$\sum_k^m \Delta b_k(n)$ 。

反向传播算法

有了以上更新规则就可以对权重$W$进行迭代求解,从更新规则来看,需要计算出 $\Delta W$ 。该计算过程可以使用反向传播算法。反向传播算法求导思想一个很好的解释可以见此博客 。总的来说:

  1. 对于复杂网络的求导,导数值是从输入到输出所有路径导数的和,而每条路径的导数可以使用链式法则进行计算;
  2. 反向传播过程即从终点开始$J$,按网络结构的拓扑序的逆序顺序对所有端点求 $\frac{\partial J}{\partial (\cdot )}$ 。

权重学习

有了以上损失函数、更新规则和导数求导的方法,现在可以对神经网络中的权重$W, b$进行求导,并根据导数给出迭代更新方程。以下层数标记的详细规定可见这里

  1. 从隐藏层到输出层

    先假设网络结构为只包含一个隐藏层,即从第二层隐藏层到输出层的权重矩阵为$W^{(2)}$。下面讨论先计算权重 $W^{(2)}$对损失函数的导数,将输出层上的节点标为 $k$ ,隐藏层中的节点标为 $j$,则有$\frac{\partial J}{\partial Wkj} =\frac{\partial J}{\partial net_k} \frac{\partial net_k}{\partial Wkj}$,其中前一项 $\frac{\partial J}{\partial net_k} $对于所有 $\frac{\partial J}{\partial Wk(\cdot) }$都相同,记 $\delta = \frac{\partial J}{\partial net_k}$。对于分类问题,规定输出层的$net_k = a_k$,则有$\delta^{(3)}_k = \frac{\partial J}{\partial a_k} \frac{\partial a_k}{\partial net_k} = (p_k - t_k) \times 1= p_k- t_k $。故有$\frac{\partial J}{\partial Wkj} =\frac{\partial J}{\partial net_k} \frac{\partial net_k}{\partial Wkj} = \delta^{(3)}_k z_j = (p_k-t_k) z_j$。将该计算过程向量化则有$\delta^{(3)} =P - T $, 及$\Delta W^{(2)} = \frac{\partial J}{\partial W^{(2)}} = \delta^{(3)} (z^{(2)})^t $,其中$\Delta W^{(2)}$是一个与$W^{(2)}$有相同形状的矩阵。

  2. 从输入层到隐藏层,隐藏层到隐藏层

    记从隐含层$n-1$到下一隐含层的权重为$W^{(n-1)}$。对于第$n$层隐含层,记其上层(第$n-1$层)的节点为 $i$,该层(第$n$层)节点为$j$,再下层(第$n+1$层)节点为$k$。则$W^{(n-1)}$是一个$n_j \times n_i$的矩阵。现假设该神经网络从$n$层开始所有 $\delta^{(o)}, \Delta W^{(o)}$都已知,现希望求第$n-1$层的$\Delta W^{(n-1)}$ 。将损失函数对矩阵值求导,则有 $\frac{\partial J}{\partial Wji^{(n-1)}} =\frac{\partial J}{\partial net_j^{(n)}} \frac{\partial net_j^{(n)}}{\partial Wji^{(n-1)}}$,同上,记$\delta_j^{(n)} = \frac{\partial J}{\partial net_j^{(n)}}$,则有
    $$
    \delta_j^{(n)} = \frac{\partial J}{\partial net_j^{(n)}} = \frac{\partial J}{\partial z_j^{(n)}} \frac{\partial z_j^{(n)}}{\partial net_j^{(n)}} = \sum_k \frac{\partial J}{\partial net_k^{(n+1)}} \frac{\partial net_k^{(n+1)}}{\partial z_j^{(n)}} \frac{\partial z_j^{(n)}}{\partial net_j^{(n)}} = \sum_k \delta_k^{(n+1)} Wkj^{(n)}f(net_j^{(n)})
    $$
    这样损失函数对该层权重的导数即为有$\frac{\partial J}{\partial Wji^{(n-1)}} =\frac{\partial J}{\partial net_j^{(n)}} \frac{\partial net_j^{(n)}}{\partial Wji^{(n-1)}} = \delta_j^{(n)} z_i^{(n-1)}$。将该层计算过程向量化有$\delta^{(n)} = f’(net^{(n)}) \circ (W^{(n)})^t \delta^{(n+1)}$,其中$\circ$运算为Hadamard变换)。这样,$\Delta W^{(n-1)} = \frac {\partial J}{\partial W^{(n-1)}} = \times \delta^{(n)} (z^{(n-1)})^t$。对于计算从输入层到隐含层的参数导数情况,只需令该式中的$z^{(n-1)} = x$ 即可求出导数。

因此,神经网络的迭代过程可以表示为,先进行正向传播求出每个节点的值,再反向传播求得$\Delta W$进而更新权重。例如,对于前面讨论的三层神经网络,正向传播的计算过程为$net^{(1)} -> z^{(1)} -> net^{(2)} -> z^{(2)} -> p^{(2)}$,反向传播的计算过程为:$\delta^{3} -> \Delta W^{(2)} -> \delta^{(2)} -> \Delta W^{(1)} $。即有:

  • $\delta^{(3)} = P - T $
  • $\Delta W^{(2)} = \delta^{(3)} (z^{(2)})^t $
  • $\delta^{(2)} = f’(net^{(2)}) \circ (W^{(2)})^t \delta^{(3)}$
  • $\Delta W^{(1)} = \delta^{(2)} (x)^t $

同理可推得$\Delta b^{(n)} = \delta^{(n+1)} $。若激活函数为Logistic-Sigmoid记为$f(\cdot)$,由于其导数为$f’(x) = f(x) (1-f(x))$可以表示为函数值的函数,记$df(x) = x(1-x)$,则$ f’(net^{(2)}) = df(z^{(2)}) = z^{(2)} \circ (1- z^{(2)}) $

有了$\Delta W, \Delta b$,就可以根据迭代规则($W(n+1) = W(n) - \eta\Delta W(n)$)对权重进行学习。

NN学习的优化

有了以上多层网络的学习算法即可以开始对网络进行学习。在实际应用时,使用以上结构得到的性能和准确率往往是不够的。一方面由于神经网络能进行高维空间映射,其可以将函数的复杂度不断上升使得其对样本集的学习分类效果达到100%,而这样的得到的复杂分界面通过是过拟合的,即该神经网络应用在测试数据集时得到较差的结果。另一方面,由于神经网络在高维空间中进行梯度下降学习,其学习过程会困于局部最小值或遇到马鞍面的情况,对于这些问题,神经网络都不能很好得解决。下面讨论对NN学习中得优化,其中向量化和加动量能提高NN运行速度,Regulation和权重衰减能缓和过拟合的情况。

  • 向量化:向量化即将函数运行的过程进行向量化运算。上面对NN计算的讨论可以看成是对单个样本$x^{(i)}$运算向量化结果,事实上,对于样本集$X = (x^{(1)}, x^{(2)}, …, x^{(m)})$,其计算过程也可以向量化。首先,对于拥有$d$个特征的单个样本$x^{(1)}$的情况,其前向计算过程为:
    $$
    x^{(1)} [d \times 1] -> net^{(2)} [n2 \times 1]-> z^{(2)} [n2 \times 1]-> net^{(3)} [n3 \times 1]-> z^{(3)} [n3 \times 1] -> p^{(3)} [n_3 \times 1]
    $$
    其中$n2$表示第二层的节点数。若将$x^{(1)}$ 替换为$X \ [d \times N]$,其中$N$为样本集容量。这样,根据前向传播的公式,只需将$b^{(1)} [n2 \times N]$,即可将上述所有矩阵的列向量变成对应独立每个样本得到的值,即可以直接得到 $p^{(3)} [n3 \times N]$,对应$N$个样本下每个类别的概率值。进一步得,根据反向传播公式$\Delta W^{(2)} [n3 \times n2]= \delta^{(3)} (z^{(2)})^t $,由于向量化得到得$ \delta^{(3)} ,z^{(2)} $的维数分别为$[n3 \times N], [n2 \times N]$,反向传播公式的到$\Delta W^{(2)}$等价为$N$个单样本得到的$\Delta W^{(2)}$的叠加!由于$\partial J(W) / \partial w = \sum_k^N \partial J(W, x_k) / \partial w$,计算样本集对于单个权重$w$的导数时需要计算$N$个单样本对权重$w$求导数再就和,故以上向量化的结果可以把这个值直接得出。同理,根据公式$\Delta b^{(2)} = \delta^{(3)} $,若要再向量化中对偏移量也得到类似的结果只需令$\Delta b^{(2)} = sum(\delta^{(3)}, axis=1) $,即$\Delta b_i^{(2)} = \sum_k^N (\delta_i^{(3)})^{(k)}$ 。详细的计算过程可以查看下面的代码实现。
  • Regulation:正规化即加入对参数数值大小的惩罚项,这样可以防止参数过大。例如,对损失函数添加$W$每个数值平方,这样有 $J(new) = J + \lambda/2 ( \sum w_i^2)$,其中$\lambda$ 是调整系数。这样即可在减少输出层与真实输出间差异的同时确保函数不会太过复杂以至于过拟合(过拟合下权重通常有较大的数值),加了规范化项后对$w_i$求导的数值就要加入一项 $\lambda w_i$ 。
  • 加动量:加动量即在更新权重值的时候考虑当前推荐的移动方向(及速度)及历史推荐的移动方向(及速度),该描述体现为将更新公式变成$W(n+1) = W(n) + (1-\alpha) \Delta W(n) + \alpha \Delta W(n-1)$,其中$\alpha$为调节动量大小的参数。加入动量后对于搜索空间中陡峭的区域权重值的变化较快,这样有利于减少函数收敛及在平坦区域上搜索的时间。
  • 权重衰减:权重衰减和正规化规范化都对权重的数值大小进行调整,权重衰减的规则时在更新后,将权重进一步缩减$W(n+1) = (1-\epsilon ) W(n+1)$,其中$\epsilon$为调节参数。

Python 实现

好了,通过上文的讨论,现在终于可以动手写代码实现了。下面实现一个三层神经网络分类器,其输出层使用Softmax函数以获得样本的后验概率估计值。

生成200个呈半月形分布数据,并先使用上次编写的感知器分出一个线性分界面。从图中可以看出感知器得到的线性分解面不能很好解决该数据集的分类问题。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
import numpy as np
from sklearn import datasets
# 生成数据
np.random.seed(0)
X, y = datasets.make_moons(200, noise=0.25)
# 感知器学习,画图
%matplotlib inline
import matplotlib.pyplot as plt
from mlxtend.plotting import plot_decision_regions
ppn = Perceptron(epochs=100, eta=0.01)
newy = np.where(y > 0, 1, -1)
ppn.increment_train(X, newy)
#ppn.batch_train(X, y)
print('Weights: %s, predict rate:%s' % (ppn.w_, ppn.get_predict_rate(X, newy)))
plot_decision_regions(X, newy, clf=ppn)
plt.title('Perceptron')
plt.xlabel('sepal length [cm]')
plt.ylabel('petal length [cm]')
plt.show()

perceptron_data

ANN实现

下面用向量化的操作实现多层神经网络。其中,输出层使用softmax函数,损失函数定义为$J(W) = \sum_k^m J(W, x_k) + \lambda/2 \sum w^2 = \sum_k \sum_i -t_i log (p_i) + \lambda/2 \sum w^2 $。学习算法使用批量下降。其中网络结构及激活函数可以自定义。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
# Artificial Neural Network
# classifier = ANN(eta = step_len, epochs = max_iteration, _lambda = lambda_value)
# classifier.train(X, y)
# make sure X = [nsample X nfeature], y = [ nsample X nclass ]
# y only have two class label {0, 1, 2, etc}, (begin with 0)
# [= 1] position in each row indicate the y_i label, e.g. y_i \in {1}, then y_i = [0, 1, 0, etc.]
class ANN(object):
def __init__(self, eta=0.01, epochs=100, _lambda = 0, \
momentum_alpha = 0.0, decay_r = 0.0):
self.eta = eta # eta should > 0
self.epochs = epochs
self._lambda = _lambda
self.momentum_alpha = momentum_alpha
self.decay_r = decay_r
# can custom to have different activation function for each layer, here assign each layer with
# identical activation function for simplifying the code structure
def init_activation_function(self, f, df):
self.f = f
self.df = df
# network struct, input_n[int], output_n[int], hide_n[n1[int], n2[int], etc.]
# the network will have a first layer that have input_n nodes
# last layer have output_n nodes
# several hide layer: each have n1 nodes, n2 nodes , etc.
# e.g init_network_structure(3, 2, [5, 5])
def init_network_structure(self, input_n, output_n, hide_n=[3], nsample = 1):
np.random.seed(0)
all_n = [input_n] + hide_n + [output_n]
self.all_n = all_n
self.ln = len(all_n)
# init weight and bias, which will use both in train and predict process
# k = 1.0 / nsample
k = 0.01
# weight, init with normal(0, (0.01)^2)
self.w = [np.random.randn(all_n[i+1], all_n[i]) * k for i in range(len(all_n)-1)]
# bias, init with normal(0, (0.01)^2)
self.b = [np.random.randn(all_n[i+1], 1) * k for i in range(len(all_n)-1)]
# init net, d_weight, d_bias, p, etc, which will only use in train process
# forward para.
self.net = [np.zeros((i, nsample)) for i in all_n]
self.z = [np.zeros((i, nsample)) for i in all_n]
# for the predict layer p
self.p = np.zeros((output_n, nsample))
# backpropagation para.
# d_weight & d_bias
self.delta = [np.zeros((i, nsample)) for i in all_n]
self.dw = [np.zeros((all_n[i+1], all_n[i])) for i in range(len(all_n)-1)]
self.db = [np.zeros((all_n[i+1], 1)) for i in range(len(all_n)-1)]
self.dw_pre = [np.zeros((all_n[i+1], all_n[i])) for i in range(len(all_n)-1)]
self.db_pre = [np.zeros((all_n[i+1], 1)) for i in range(len(all_n)-1)]
# i = layer_i
# from second layer to last layer
def forward_propagation(self, i):
# net[n_l X nsample] = [n_l X nsample] + [n_l X 1]
# the self.b will change from to [n_l X 1] to [n_l X nsample] automatically
self.net[i] = self.w[i-1].dot(self.z[i-1]) + self.b[i-1]
self.z[i] = self.f(self.net[i]) if i != self.ln - 1 else self.net[i]
# if is the last output layer
# update output
if i == self.ln - 1:
self.p = np.exp(self.z[i])
sum_p = self.p.sum(0)
# exp(p) so sum_p > 0
self.p = self.p / sum_p
# i = layer_i
# from last layer to first layer
def backpropagation(self, i, T=0):
# if is the last layer
if i == self.ln - 1:
self.delta[i] = self.p - T
else:
# self.df(self.z[i]) == (self.f(self.net[i]))'
self.delta[i] = self.df(self.z[i]) * \
(self.w[i].T).dot(self.delta[i+1])
# dw is the sum of all sample dw_i
# print self.dw[i].shape, self.delta[i+1].shape, self.z[i].T.shape
self.dw[i-1] = self.delta[i].dot(self.z[i-1].T)
# db is the sum of all sample db_i
self.db[i-1] = np.sum(self.delta[i], axis=1, keepdims=True)
# get probability
# X = [n_new_sample X nfeature]
# return p = [nclass X n_new_sample]
def get_p(self, X):
tmp_z = X.T
# from second layer to last layer
for i in range(1, self.ln):
tmp_net = self.w[i-1].dot(tmp_z) + self.b[i-1]
tmp_z = self.f(tmp_net) if i != self.ln - 1 else tmp_net
tmp_p = np.exp(tmp_z)
sum_p = tmp_p.sum(0)
p = tmp_p / sum_p
return p
# transform probability to class label
# return predict = [n_new_sample, ] array
def predict(self, X):
p = self.get_p(X)
## transform p to class label
return np.argmax(p, axis=0)
# calculate sample set loss
def calculate_loss(self, X, y):
p = self.get_p(X)
# Calculating the loss
nsample = X.shape[0]
# sum(J(i)) / m
loss = -1 * np.sum(np.log(p) * y.T) / nsample
# regularization term
loss += self._lambda / 2.0 * sum([np.sum(np.square(_w)) for _w in self.w])
return loss
def calculate_match_rate(self, X, y):
pre_label = self.predict(X)
# Calculating the loss
match_time = sum([y[i][pre_label[i]] for i in range(len(y))])
return 1.0 * match_time / y.shape[0]
def log_shape(self):
print 'layer nodes: ', self.all_n
print 'p :', self.p.shape, self.p
for i in range(self.ln):
print 'layer %d: ' % (i)
print 'net :', self.net[i].shape, self.net[i]
print 'z :', self.z[i].shape, self.z[i]
print 'delta :', self.delta[i].shape, self.delta[i]
if i != self.ln - 1:
print 'w :', self.w[i].shape, self.w[i]
print 'b :', self.b[i].shape, self.b[i]
print 'dw :', self.dw[i].shape, self.dw[i]
print 'db :', self.db[i].shape, self.db[i]
print ' '
# btach train
def train(self, X, y, hide_n = [3], f=None, df=None, \
_lambda = None, eta=None, epochs=None, \
momentum_alpha = 0.0, decay_r = 0.0, print_loss=0):
# init struct first
nsample, nfeature = X.shape
nsample, nclass = y.shape
self.init_network_structure(nfeature, nclass, hide_n, nsample)
# init activation_function
self.init_activation_function(f, df)
self.eta = self.eta if eta is None else eta
self.epochs = self.epochs if epochs is None else epochs
self._lambda = self._lambda if _lambda is None else _lambda
self.momentum_alpha = momentum_alpha
self.decay_r = decay_r
# init input layer, and target nodes
self.z[0] = X.T
target = y.T
for _ in range(self.epochs):
# Forward propagation, for second to last layer
for i in range(1, self.ln):
self.forward_propagation(i)
# Backpropagation, for last-1 to first layer
# update output layer
self.backpropagation(self.ln-1, target)
# for self.ln-2 to 0
for i in range(self.ln-2, 0, -1):
self.backpropagation(i)
# after backpropagation, weight of loss of all sample's is store in dw
# update W and b
for i in range(0, self.ln-1):
update_w = self.eta * (self.dw[i] + self._lambda * self.w[i])
update_d = self.eta * self.db[i]
# self.w[i] -= update_w
# self.b[i] -= update_d
if momentum_alpha > 0 and _ > 0:
self.w[i] -= update_w * (1 - momentum_alpha ) + momentum_alpha * self.dw_pre[i]
self.b[i] -= update_d * (1 - momentum_alpha ) + momentum_alpha * self.db_pre[i]
else:
self.w[i] -= update_w
self.b[i] -= update_d
self.w[i] = self.w[i] * (1.0 - self.decay_r)
#self.b[i] = self.b[i] * (1.0 - self.decay_r)
if momentum_alpha > 0.0:
self.dw_pre[i] = update_w
self.db_pre[i] = update_d
# self.log_shape()
if print_loss and _ % print_loss == 0:
print "Loss after iteration %i: %f; match times:%f" %(_, \
self.calculate_loss(X, y), self.calculate_match_rate(X, y))

测试结果与参数讨论

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
# Logistic-Sigmoid
# f = lambda x: 1.0 / (1 + np.exp(-x)) #apply to net
# df = lambda x: x * (1 - x) #apply to z
# # Tanh-Sigmoid
f = lambda x: np.tanh(x)
df = lambda x: 1 - x*x
# init input data for NN
# 200 X 2
nsample, nfeature = X.shape
nclass = 2
tX = X
ty = np.zeros((nsample, nclass))
for i in range(len(ty)):
ty[i][int(y[i])] = 1
ty = ty
print tX.shape, ty.shape
nn = ANN()
nn_hdim = [3]
nn.train(tX, ty, hide_n = nn_hdim, f=f, df=df, \
eta = 0.001, epochs=20000, _lambda = 0.001, \
momentum_alpha = 0, decay_r = 0, print_loss=1000)
plt.title('Hidden Layer %s, match time: %.1f%%' % (nn_hdim, \
100 * nn.calculate_match_rate(tX, ty)))
plot_decision_regions(X, y, clf=nn)
plt.show()

在不同激活函数、网络结构、参数下进行测试。首先是激活函数,在该问题上tanh函数与logistic函数得到相似的结果。

nn_predict

下面设置$\eta = 0.01, \lambda = 0.001$,迭代次数为30000次在不同的隐藏层节点数下进行测试。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
plt.figure(figsize=(16, 32))
var_layer = [[1], [3], [2], [5], [10], [50]]
var_eta = [0.001] * 6
var_epochs = [30000] * 6
var_f = [f] * 6
var_df = [df] * 6
for i, nn_hdim in enumerate(var_layer):
nn.train(tX, ty, hide_n = var_layer[i], f=f, df=df, \
eta = var_eta[i], epochs=var_epochs[i], _lambda = 0.001, \
momentum_alpha = 0, decay_r = 0, print_loss=10000)
plt.subplot(5, 2, i+1)
if i == 0:
plt.title('hide layer %s, match time: %.1f%%' % (var_layer[i],\
100 * nn.calculate_match_rate(tX, ty)))
else:
plt.title('hide layer %s, match time: %.1f%%' % (var_layer[i],\
100 * nn.calculate_match_rate(tX, ty)))
plot_decision_regions(X, y, clf=nn)
plt.show()

nn_predict

可以看出在隐藏层只有一层的情况下,增加隐藏层的节点数量能提高分类性能,但节点增加到一定数量后如果不能增加迭代次数和调整相应的$\eta,\lambda$性能会饱和。下图为隐藏层节点都为50,在不同迭代次数和步长下学习的结果。可见隐藏节点较多的时候网络倾向于过拟合。

nn_predict

下面是设置不同的隐藏层结构下的结果。图上参数表示 [n1, n2] 表示有两个隐藏层,第一层有n1个节点,第二层有n2个节点。增加隐藏层的层数能提高分类性能,但是会增加函数的复杂度,使得结果倾向于过拟合。

nn_predict

总结

待续

参考

[1] Pattern classification: https://www.amazon.com/Pattern-Classification-Pt-1-Richard-Duda/dp/0471056693

[2] The element of statistical learning: https://www.amazon.com/Elements-Statistical-Learning-Prediction-Statistics/dp/0387848576

[3] Implementing a Neural Network from Scratch in Python: http://www.wildml.com/2015/09/implementing-a-neural-network-from-scratch/

[4] UFLDL: http://ufldl.stanford.edu/wiki/index.php/UFLDL%E6%95%99%E7%A8%8B