C++实现深度神经网络解决多分类逻辑回归问题

如何用C++实现一个深度神经网络(DNN, Deep Neural Network)处理多分类逻辑回归问题(MLR, Multiclass Logistic Regression)且不使用第三方库

唔……

今天就来尝试一下: 使用$\text{sigmoid}$激活函数+$\text{Softmax}$逻辑回归+交叉熵损失函数+反向传播实现一个DNN。整个过程只使用C++标准库而不使用任何第三方库

建议食用本文前请掌握神经网络的基本知识,或者先忽略在”约定和初始化”部分所看到的陌生函数和名词,在后文中会有相应解释。

P.S.: 由于学习Go语言的缘故,我开始左大括号不换行了QwQ。

问题

给定$k$维空间的$K$个点集(在几何距离上相近)。每次查询一个新点,问这个点最有可能属于哪一个点集。

$k$在一定范围内增大不会影响我们的解法,因此这里我们先从简处理$k = 2$的情况,即$2$维平面的点集。


方案

如何表示一个划分方案?

用一个$K \times 1$向量表示,第$i$个数字表示这个点属于点集$i$的概率。

对于初始点集作为给出的标准答案来说一个属于$i$点集的点属于点集$i$的概率当然是$1$。


损失函数

如何评价一个划分方案的好坏?

在监督学习(提供标准答案)的情况下,常用的有两类评价函数,一种是平方差损失函数($MSE$),公式如下:

这个$1/2$是为了求导方便。$MSE$是一种很容易想到而且有效性显而易见的损失函数。不过另一种相对更好的评价函数是交叉熵损失函数($CE$),公式如下:

注意这里$\sum_{i=1}^Ky_i = 1$。$Loss$越小的划分方案我们认为越优,直观理解就是对于$y_i$比较大的项,追求$\ln(\frac{1}{\hat y_i})$较小即$\frac{1}{\hat y_i}$较小即$\hat y_i$较大,这与$y_i$是一致的。


代码实现:

1
2
inline db MSELoss(){db S = 0; for(int j = 0; j < K; S += (Out[j]-Ans[j])*(Out[j]-Ans[j]), j++); return S/2/K;}
inline db CELoss(){db S = 0; for(int j = 0; j < K; S += CE(Ans[j],Out[j]), j++); return -S;}


约定和初始化

为了处理好细节问题,现在这里对神经网络的结构做一下说明和实现。

一开始我试图用一个二维数组来表示神经网络,每个神经元对应一个坐标$(i,j)$。这种方法因为内存连续性,常数比较优秀。但是不能很好的处理每层节点数不同的情况,无法自定义空间复杂度已经最大层数和最大节点数。

还有一种常见的方法使用矩阵表示所有的转移。这样的实现简洁不易出错,代码更易理解。但是灵活性较差,往往需要使用第三方库(例如OpenCV)。而手动实现矩阵则凭空增加了代码复杂度。

于是我选择用vectorstruct组合来完成这一任务。这种做法是最直观的(相对于矩阵方法更贴近于神经网络的物理意义而不是数学意义)。


定义一个神经网络有$N+1$层,分别记为第$0,1,2,…,N$层。

其中第$0$层为输入层,只有$k = 2$个节点,输入形式为$x/\text{RANGE},y/\text{RANGE}$,其中$x$为坐标范围,这样保证输入的每个数都与$1$在一个数量级上,防止$\text{sigmoid}$函数将输入”抹平”。

所有$0$到$N-1$层的激活函数都是$\text{sigmoid}$,但第$N$层的激活函数为$\text{Softmax}$,作为输出层,有恰好$K$个节点。

第$i$层有$L[i].n+1$个节点,其中$0$节点代表偏置,因此$0$节点的输出总是$1$,偏置由边权调整。以后用$(i,j)$表示第$i$层第$j$个节点。

每个节点有$h$表示未经激活函数处理的输入,$x$表示经过激活函数处理的输出,$d$表示反向传播的误差。以上内容按层记录,另外每层$L[i].w[j][k]$表示这样一条边: $(i-1,k) \xrightarrow{L[i].w[j][k]} (i,j)$。

$\eta$是学习率,代码中为etaInOut以及Ans是对外的接口,为了增强通用性,这些下标从$0$开始,注意这意味着接口和神经网络有$1$的下标偏差。


因此基本结构的代码实现如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
#define db double
#define range(x) (x).begin(),(x).end()
#define sigmoid(x) (1./(1.+exp(-(x))))
#define CE(a,b) ((a)*log(max(b,1e-12))+(1-(a))*log(max(1-(b),1e-12)))
typedef vector<db> Arr; mt19937 _rd(1061109589); auto real_rd = std::bind(std::normal_distribution<db>(0.,.01),mt19937(998244353));
struct NeuralNetwork{
struct Layer{
Arr h,x,d; vector<Arr> w; int n;
inline Layer(int _n, int c){
n = _n; x = h = d = Arr(n+1); x[0] = 1; w.clear();
for(int i = 0, j; i <= n; w[i][0] = 0.5, i++)
for(w.push_back(Arr(c+1)), j = 1; j <= c; w[i][j] = real_rd(), j++);
}
}; vector<Layer> L; Arr Out, Ans, In; int N, K; db eta;
inline NeuralNetwork(vector<int> Num, db _eta){L.clear(); eta = _eta; N = -1; K = 0; for(int&j:Num) AddLayer(j);}
}

这里的real_rd是一个正态分布$\mu = 0, \sigma = 0.01$,即$N(0,0.01)$。

所有的边权被随机赋值为一个接近于$0$而很小的数,而偏置常数项被设定为中值$0.5$。

$\eta$设置为$0.01 \sim 0.1$。

我也不知道为啥这样玄学调参

注意$w$的初始化要求比较高: 显然不能全部初始化为$0$,这样整个网络永远都是$0$;也不能全部初始化为一个数,这样同一层所有节点是完全对称的地位,结果也会完全一样,节点数目就形同虚设了;也不能全部随机初始化,因为$\sum wx$很大会导致$\text{sigmoid}$函数值全都趋向于$1$。


神经元和前向传播

如何实现一个神经元?

之前提到过$(i-1,k) \xrightarrow{L[i].w[j][k]} (i,j)$,至于为什么使用$L[i].w[j][k]$而不是$L[i-1].w[k][j]$,是为了后面代码方便。接下来简记为$w_{i,j,k}$,其他同理。

具体地,神经元$(i,j)$受到上一层神经元的影响是每个神经元的输出($x_{i-1,k}$)乘以边权:

可以看到,因为采用$w_{i,j,k}$而不是$w_{i-1,k,j}$,这里可以记为向量点积(inner_product)。$w_{i,j,0}$表示的是偏置常数项。实现的时候,偏置常数项的输出$x_{i,0}$永远是$1$。

然后每个神经元将其受到的影响经过激活函数处理后输出:

其中$\text{sigmoid}(x) = \frac{1}{1+e^{-x}}$,也可以简记为$\sigma(x)$。注意到$\sigma$函数值域是$(0,1)$,只能实现二分类(判断靠近$0$还是靠近$1$),一个它的拓展是$\text{Softmax}$函数,是一个序列到等长序列的映射:

事实上就是经过$e^x$放大以后求出每个数字占整体的比例,来得到一个概率分布。


代码实现:

运行神经网络:

1
2
3
4
5
6
inline void Run(){
for(int i = (copy(range(In),++L[0].x.begin()),1), j; i <= N; i++)
for(j = 1; j <= L[i].n; L[i].x[j]=sigmoid(L[i].h[j]=inner_product(range(L[i-1].x),L[i].w[j].begin(),0.)), j++);
Out.resize(K); for(int j = 0; j < K; Out[j] = exp(L[N].h[j+1]), j++); db S = accumulate(range(Out),0.); for(db&j:Out) j /= S;
}
inline int Judge(){return (int)(max_element(range(Out))-Out.begin());}

Judge()就是输出概率分布中最大值作为神经网络最后的分类,只在检测正确率的时候用到。这里的前向传播和最后的$\text{Softmax}$函数合并了,也可以单独提取出$\text{Softmax}$函数:

1
inline void Softmax(){Out.resize(K); for(int j = 0; j < K; Out[j] = exp(L[N].h[j+1]), j++); db S = accumulate(range(Out),0.); for(db&j:Out) j /= S;}


梯度下降反向传播

系好安全带,开始飙车了!

现在每个$w$都是一开始设定好的,这个神经网络根本没有学习能力。如何让神经网络具有学习能力?

关键就是通过大量数据调整所有$w_{i,j,k}$的值,使得最后的$Loss$最小。

一个很常见的思路就是求导: 每次贪心地将每个$w_{i,j,k}$尝试调大或者调小。至于具体调大还是调小,当然就是向着改变以后$Loss$会降低的方向调整。$Loss$降低的方向就是$Loss$关于$w_{i,j,k}$的偏导数正向的反方向。即每次调整:

这种方法就是梯度下降算法。$\eta$是学习率(就是下降的步长)。

先来尝试推导最后一层节点的偏导数。首先回忆之前的交叉熵损失函数以及$\text{Softmax}$函数:

就可以得到:

最后一个等式是因为$Ans$是一个概率分布,$\sum_{i=1}^K Ans_i = 1$。于是就有:

利用链导法则推广到一般情况: 如果记$d_{i,j} = - \eta \frac{\partial Loss}{\partial h_{i,j}}$,那么就有$\Delta w_{i,j,k} = d_{i,j}x_{i-1,k}$,而:

只要能求出后半部分偏导数就做完了:

然后这里就需要用到$\sigma(x)$的美妙性质了:

因此:

大成功!!!


代码实现:

1
2
3
4
5
6
7
inline void Adjust(){
for(int j = 1; j <= K; L[N].d[j] = -eta*(Out[j-1]-Ans[j-1]), j++);
for(int i = N-1, j, k; i; i--) for(fill(range(L[i].d),0.), j = 1; j <= L[i+1].n; j++)
for(k = 0; k <= L[i].n; L[i].d[k] += L[i+1].d[j]*L[i+1].w[j][k]*L[i].x[k]*(1-L[i].x[k]), k++);
for(int i = N, j, k; i; i--) for(j = 1; j <= L[i].n; j++)
for(k = 0; k <= L[i-1].n; L[i].w[j][k] += L[i].d[j]*L[i-1].x[k], k++);
}

注意一下前向传播和反向传播复杂度都是$O(Nn^2)$的。而反向传播由于链式法则求导运算复杂得多,并且最重要的是反向传播不满足内存连续性。因此反向传播这个部分可以说是整个神经网络的复杂度瓶颈: 一个微小的改动可能造成巨大的常数差异。例如:

1
L[i].d[k] += L[i+1].d[j]*L[i+1].w[j][k]*L[i].x[k]*(1-L[i].x[k])

如果按照之前的公式直接写为:

1
L[i].d[k] += L[i+1].d[j]*L[i+1].w[j][k]*sigmoid(L[i].h[k])*(1-sigmoid(L[i].h[k]))

由于多调用了两次sigmoid函数,而exp函数又十分缓慢,仅仅这一点差别可以让运行时间相差数倍至数十倍(亲测血的教训)。


神经网络

于是我们就完成了,现在放出刚刚完成的完整的深度神经网络DNN代码:

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
// NN.h
#include <bits/stdc++.h>
using namespace std;
#define db double
#define range(x) (x).begin(),(x).end()
#define sigmoid(x) (1./(1.+exp(-(x))))
#define CE(a,b) ((a)*log(max(b,1e-12))+(1-(a))*log(max(1-(b),1e-12)))
typedef vector<db> Arr; mt19937 _rd(1061109589); auto real_rd = std::bind(std::normal_distribution<db>(0.,.01),mt19937(998244353));
struct NeuralNetwork{
struct Layer{
Arr h,x,d; vector<Arr> w; int n;
inline Layer(int _n, int c){
n = _n; x = h = d = Arr(n+1); x[0] = 1; w.clear();
for(int i = 0, j; i <= n; w[i][0] = 0.5, i++)
for(w.push_back(Arr(c+1)), j = 1; j <= c; w[i][j] = real_rd(), j++);
}
}; vector<Layer> L; Arr Out, Ans, In; int N, K; db eta; inline void AddLayer(int n){++N; L.emplace_back(n,K); K = n;}
inline NeuralNetwork(vector<int> Num, db _eta){L.clear(); eta = _eta; N = -1; K = 0; for(int&j:Num) AddLayer(j);}
inline void Run(){
for(int i = (copy(range(In),++L[0].x.begin()),1), j; i <= N; i++)
for(j = 1; j <= L[i].n; L[i].x[j]=sigmoid(L[i].h[j]=inner_product(range(L[i-1].x),L[i].w[j].begin(),0.)), j++);
Out.resize(K); for(int j = 0; j < K; Out[j] = exp(L[N].h[j+1]), j++); db S = accumulate(range(Out),0.); for(db&j:Out) j /= S;
}
inline int Judge(){return (int)(max_element(range(Out))-Out.begin());}
inline db CELoss(){db S = 0; for(int j = 0; j < K; S += CE(Ans[j],Out[j]), j++); return -S;}
inline void Adjust(){
for(int j = 1; j <= K; L[N].d[j] = -eta*(Out[j-1]-Ans[j-1]), j++);
for(int i = N-1, j, k; i; i--) for(fill(range(L[i].d),0.), j = 1; j <= L[i+1].n; j++)
for(k = 0; k <= L[i].n; L[i].d[k] += L[i+1].d[j]*L[i+1].w[j][k]*L[i].x[k]*(1-L[i].x[k]), k++);
for(int i = N, j, k; i; i--) for(j = 1; j <= L[i].n; j++)
for(k = 0; k <= L[i-1].n; L[i].w[j][k] += L[i].d[j]*L[i-1].x[k], k++);
}
};
#undef db
#undef range

瞧!一个深度神经网络的核心代码居然只有30行左右!


训练

如何训练神经网络?

核心问题是如何造数据。这里自己YY了一下,采用这种造数据的方式: 随机一个点集数目$K$,对于每个点集,随机一个矩形区域,在里面生成大量点。然后反复选择其中一部分作为初始给出,另一部分作为检测数据。

下面是一个实例:

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
#include <bits/stdc++.h>
#include "NN.h"
using namespace std;
#define TESTS 5
#define RANGE 100000
#define SET_RANGE 100
#define NUM_PER_SET_MIN 4000
#define NUM_PER_SET_MAX 5000
#define db double
#define rint register int
inline int rf(){int r;int s=0,c;for(;!isdigit(c=getchar());s=c);for(r=c^48;isdigit(c=getchar());(r*=10)+=c^48);return s^45?r:-r;}
struct P{int x,y,i;}; vector<P> Set[32], Homework, Exam; mt19937 rd(time(0)); int K = rd()%30+1, T = TESTS, _; db Loss, Score;
int main(){
for(rint i = K, j, LX, LY; i; i--)
for(LX = rd()%RANGE+1, LY = rd()%RANGE+1, j = rd()%(NUM_PER_SET_MAX-NUM_PER_SET_MIN)+NUM_PER_SET_MIN; j--; )
Set[i].push_back(P{rd()%SET_RANGE+LX,rd()%SET_RANGE+LY,i-1});
for(rint i, j; T--; printf("%.6lf %.6lf\n",Loss/Exam.size(),Score/Exam.size())){
for(vector<P>().swap(Homework), vector<P>().swap(Exam), i = 1; i <= K; i++)
for(shuffle(Set[i].begin(),Set[i].end(),rd), j = 0; j < Set[i].size(); j++)
(j<NUM_PER_SET_MIN-10?Homework:Exam).push_back(Set[i][j]);
NeuralNetwork NN({2,50,K},.1);
// NeuralNetwork NN({2,100,K},.1);
// NeuralNetwork NN({2,200,K},.1);
// NeuralNetwork NN({2,30,30,K},.1);
// NeuralNetwork NN({2,50,50,K},.1);
// NeuralNetwork NN({2,100,100,K},.1);
Loss = Score = 0;
shuffle(Homework.begin(),Homework.end(),rd); _ = 0;
for(auto p:Homework){
NN.In = {(db)p.x/RANGE,(db)p.y/RANGE};
NN.Ans = Arr(K); NN.Ans[p.i] = 1;
NN.Run(); NN.Adjust();
/*
printf("Training #%d:\n",++_);
for(j = 0; j < K; printf("%4.4lf ",NN.Out[j]), j++); puts("");
for(j = 0; j < K; printf("%4.4lf ",NN.Ans[j]), j++); puts("");
printf("Output = %d, Answer = %d\n\n",NN.Judge(),p.i);
//*/
}
shuffle(Exam.begin(),Exam.end(),rd); _ = 0;
for(auto p:Exam){
NN.In = {(db)p.x/RANGE,(db)p.y/RANGE};
NN.Ans = Arr(K); NN.Ans[p.i] = 1;
NN.Run(); Loss += NN.CELoss(); Score += NN.Judge()==p.i;
/*
printf("Testing #%d:\n",++_);
for(j = 0; j < K; printf("%4.4lf ",NN.Out[j]), j++); puts("");
for(j = 0; j < K; printf("%4.4lf ",NN.Ans[j]), j++); puts("");
printf("Output = %d, Answer = %d\n\n",NN.Judge(),p.i);
//*/
}
} return 0;
}

这里先离题说一个小Trick:

1
2
//*
//*/

可以方便地实现多行注释,只要删除第一个/就可以转变为注释状态,加上第一个/就转变为不注释状态。


回到正题,当提供的数据在$4000 * K$左右时,

令人出乎意料的是$NN(\lbrace2,50,K\rbrace,.1)$作为最小的神经网络竟然表现最好,基本可以保持正确率在$90 \%$以上,且绝大多数情况下都是$100 \%$正确。

$NN(\lbrace2,200,K\rbrace,.1)$表现也很好,基本可以保持正确率在$90 \%$以上,且多数情况下是$100 \%$正确。

$NN(\lbrace2,30,30,K\rbrace,.1)$基本可以保持正确率在$70 \%$以上。

$NN(\lbrace2,100,100,K\rbrace,.1)$基本可以保持正确率在$85 \%$以上,经常可以达到$95 \%$甚至$100 \%$。

可能是由于训练数据集大小、神经网络对于这一类问题的适应性、随机的不确定性等原因,单隐层神经网络表现得比多隐层神经网络要好,而且隐层节点数适度减少并不一定会导致正确率下滑。


扫描二维码即可在手机上查看这篇文章,或者转发二维码来分享这篇文章:


文章作者: Magolor
文章链接: https://magolor.cn/2019/08/02/2019-08-02-blog-01/
版权声明: 本博客所有文章除特别声明外,均采用 CC BY-NC-SA 4.0 许可协议。转载请注明来自 Magolor
扫描二维码在手机上查看或转发二维码以分享Magolor的博客