人工神经网络入门教程

1 人工神经网络简介

生物学动机

人工神经网络ANN的研究一定程度上受到了生物学的启发，生物的学习系统由相互连接的神经元（neuron）组成的异常复杂的网格。而人工神经网络由一系列简单的单元相互密集连接构成的，其中每一个单元有一定数量的实值输入，并产生单一的实数值输出。
据估计人类的大脑是由大约(10^{11})次方个神经元相互连接组成的密集网络，平均每个神经元与其他(10^4)个神经元相连。神经元的活性通常被通向其他神经元的连接激活或抑制。

2 神经网络表示

1993年的ALVINN系统是ANN学习的一个典型实例，这个系统使用一个学习到的ANN以正常的速度在高速公路上驾驶汽车。ANN的输入是一个30*32像素的网格，像素的亮度来自一个安装在车辆上的前向摄像机。ANN的输出是车辆行驶的方向。

3 适合神经网络学习的问题

• 实例是用很多属性-值表示的
• 目标函数的输出可能是离散值、实数值或者由若干实数属性或者离散属性组成的向量
• 训练数据可能包含错误
• 可容忍长时间的训练
• 可能需要快速求出目标函数值
• 人类能否理解学到的目标函数不是重要的。

4 感知器

感知器以一个实数值向量作为输入，计算这些输入的线性组合 ，然后如果结果大于某个阈值，就输出1，否则输出-1。

[o(x_1, x_2,...,x_n) = egin{equation} egin{cases} 1 quad if (w_0 + w_1x_1+...+w_nx_n)\ -1 quad otherwise end{cases} end{equation} ]

假设输入空间（特征空间）是(X subseteq R^n)，输出空间是(Y = {-1, +1})。输入(x subseteq X)是实例的特征向量，对应于输入空间的点；输出(y subseteq Y)表示实例的类别。由输入空间到输出空间的如下函数称为感知机：

[f(x)=sign(w cdot x + b) ]

其中w b称为感知机的模型参数，w叫做权值（weight）或者权值向量（weight vector），b叫做偏置（bias）。(wcdot x)表示w和x的内积。sign是符号函数：

[sign(x) = egin{equation}egin{cases} +1, x >= 0 \ -1 , x <= 0 end{cases}end{equation} ]

感知机是一种线性分类模型，它的假设空间是定义在特征空间中的所有线性分类模型。
感知机有如下几何解释：线性方程$$w cdot x + b = 0$$
对应于特征空间的一个超平面S。这个超平面把特征空间划分成两部分，位于两部分的点分别被分为正负两类。

5 多层网络和反向传播算法

5.1 可微阈值单元

• 多个线性单元的连接仍然是产生的是线性函数，而我们更选择选择能够表征非线性函数的网络。
• 感知器单元的不连续阈值使它不可微，所以不适合梯度下降算法。

我们需要的是这样一种单元，它的输出是输入的非线性函数，并且输出是输入的可微函数。sigmoid单元是一种非常类似感知器的单元，并且它基于一个平滑的可微阈值函数。

[o = sigma(w * x) ]

[sigma(y)=frac1{1+e^{-y}} ]

• (sigma)经常被称为Sigmoid函数或者logistic函数。它的输出范围是0到1，随输入单调递增。因为这个函数把非常大的输入值映射到一个小范围的输出，它经常被称为Sigmoid单元的挤压函数。
• sigmoid函数的导数很容易用它的输出表示。

[frac{mathrm{d} sigma(y)}{mathrm{d} y} = sigma(y) * (1 - sigma(y)) ]

5.2 反向传播算法

对于一系列确定的单元互联形成的多层网络，反向传播算法可用来学习这个网络的权值。它采用梯度下降方法试图最小化网络输出值与目标值之间的误差平方。

[E(w) = frac1{2}sum_{din{D}}sum_{kin{outputs}}(t_{kd}-o_{kd})^2 ]

其中outputs是网络输出单元的集合，(t_{kd})和(o_{kd})是与训练样例(d)和第(k)个输出单元相关的输出值。
反向传播算法面临的学习问题是搜索一个巨大的假设空间，这个空间由网络中所有单元的所有可能的权值定义。
在多层网络中，误差曲面可能有多个局部极小值。梯度下降仅能保证收敛到局部极小值，而未必得到全局最小的误差。
包含两层sigmoid单元的前馈网络的反向传播算法：

---

• 创建具有(n_{in})个输入，(n_{hidden})个隐藏单元，(n_{out})个输出单元的网络。
• 初始化所有的网络权值为小的随机值（-0.05-0.05）
• 在遇到终止条件前：
  对于训练样例train_sample中的每一个(x, t)（x网络输入值的向量，t网络输出值的向量）

1. 把实例(x)输入网络，并计算网络中每个单元(u)的输出(o)
2. 对于网络的每一个输出单元(k)，计算它的误差项(delta_k) $$delta_k = o_k(1-o_k)(t_k-o_k)$$
3. 对于网络的每个隐藏单元(h)，计算它的误差(delta_h) $$delta_h = o_h(1-o_h)sum_{kin outputs}w_{kh}delta_k$$
4. 更新每个网络权值(w_{ji}) $$w_{ji}=w_{ji} + Delta w_{ji} $$
   其中(delta_{ji}=eta delta_j x_{ji})

---

梯度下降更新法则依照以下三者来更新每一个权： 学习速率(eta)，该权值涉及到输入(x_{ji})，这个单元输出的误差。
为了直观的理解它，先考虑网络的每一个输出单元k的(delta_k)是怎样计算的。(delta_k)与delta法则中的((t_k-o_k))相似，但乘上了挤压函数的导数(o_k(1-o_k))。
每个隐藏单元h的(delta_h)具有相似的形式。因为训练样例仅对网络的输出提供了目标值，所有缺少直接的目标值来计算隐藏单元的误差值。因此采取下边间接办法计算隐藏单元的误差项：对受隐藏单元h影响的每一个单元的误差(delta_h)进行加权求和，每个误差(delta_k)权值为(w_{kh})，(w_{kh})就是从隐藏单元h到输出单元k的权值。这个权值刻画了隐藏单元h对于输出单元k的误差应“负责”的程度。

5.3 反向传播算法的推导

随机的梯度下降算法迭代处理训练样例，每次处理一个。对于每个训练样例(d)，利用这个样例的误差(E_d)的梯度修改权值。换句话说，对于每一个训练样例(d)，每个权(w_{ji})增加(Delta w_{ji})。

[Delta w_{ji} = -etafrac{partial E_d}{partial w_{ji}} ]

[E_d(w) = frac1{2}sum_{kin{outputs}}(t_k - o_k)^2 ]

• (x_{ji}) = 单元(j)的第(i)个输入
• (w_{ji}) = 与单元(j)的第(i)个输入相关联的权值
• (net_j = sum{w_{ji}x_{ji}})
• (o_j) = 单元(j)计算出的输出
• (t_j) = 单元(j)的目标输出
• (sigma) = sigmoid函数
• (outputs =)网络最后一层的单元的集合
• (Downstream(j))单元的直接输入中包含单元(j) 的输出的单元的集合

权值(w_{ji})仅能通过(net_j)影响网络的其他部分。所以我们用链式规则得到 ：

[frac{partial E_d}{partial w_{ji}} = frac{partial E_d}{partial net_j}frac{partial net_j}{partial w_{ji}} = frac{partial E_d}{partial net_j} x_{ji} ]

剩下的任务就是为$frac{partial E_d}{partial net_j} $导出一个方便的表达式。

1. 输出单元的权值训练法则

[frac{partial E_d}{partial net_j} = frac{partial E_d}{partial o_j} frac{partial o_j}{partial net_j} ]

（1） 考虑公式中的第一项（$frac{partial E_d}{partial o_j} $）

[frac {partial E_d}{partial o_j} = frac {partial}{o_j} frac 1 2 sum_{k in outputs} (t_k - o_k)^2 ]

除了当$k = j (时，所有输出单元k的导数)frac {partial} {partial o_j} (t_k-o_k)^2(为0。 所以我们不必对多个输出单元求和，只需设)k=j$

[ egin{split} frac{partial E_d}{partial o_j} &= frac{partial }{partial o_j}frac12(t_j-o_j)^2 \ &= (t_j-o_j)frac{partial(t_j-o_j)}{partial o_j} \ &= -(t_j-o_j)end{split} ]

（2） 考虑公式中的第二项（(frac{partial o_j}{partial net_j})）
既然(o_j=sigma(net_j))，导数(frac{partial o_j}{partial net_j})就是sigmoid函数的导数。

[frac {partial o_j}{partial net_j} = frac{partialsigma(net_j)}{partial net_j} = o_j(1-o_j) ]

由上可得：

[frac{partial{E_d}}{partial{net_j}}=-(t_j-o_j)o_j(1-o_j) ]

输出单元的随机梯度下降法则：

[Delta w_{ji}=etafrac{partial E_d}{partial w_{ji}}=eta(t_j-o_j)o_j(1-o_j)x_{ji} ]

2. 隐藏单元的权值训练法则
   对于网络中的内部单元或者隐藏单元的情况，推导(w_{ji})必须考虑(w_{ji})间接地影响网络输出，从而影响(E_d)。
   (net_j)只能通过 (Downstream(j))中的单元影响网络输出，再影响(E_d)。

[egin{split}frac{partial E_d}{partial net_j} &=sum_{kin Downstream(j)}frac{partial E_d}{partial net_k}frac{net_k}{net_j} \&=sum_{kin Downstream(j)}-delta_kfrac{net_k}{net_j}\&=sum_{kin Downstream(j)}-delta_kfrac{partial net_k}{partial o_j}frac {partial o_j}{partial net_j}\&=sum_{kin Downstream(j)}-delta_k w_{kj}frac {partial o_j}{partial net_j}\&=sum_{kin Downstream(j)}-delta_k w_{kj}o_j(1-o_j)end{split} ]

用(delta_j)表示(-frac{partial E_d}{partial net_j})，得到：

[delta_j=o_j(1-o_j)sum_kin Donwstream(j)delta_k w_{kj} ]

[Delta w_{ji}=eta delta_j x_{ji} ]

5.4 代码实践

//#include "load_dataset.hpp"

#include <stdint.h>
#include <sys/stat.h>

#include <fstream>
#include <string>

#include <opencv2/opencv.hpp>
using namespace cv;
#include <iostream>
#include <string>
#include <vector>
using namespace std;

#include <time.h>
#include <stdio.h>
#include <stdlib.h>

uint32_t swap_endian(uint32_t val) {
    val = ((val << 8) & 0xFF00FF00) | ((val >> 8) & 0xFF00FF);
    return (val << 16) | (val >> 16);
}

inline cv::Mat str2Mat(const char* pixels, const int row, const int col){
    cv::Mat m(row, col, CV_8UC1);
    for (int y = 0; y < row; y++) {
        memcpy(m.data + y * m.step, pixels + y * col, col);
    }
    return m;
}

void convert_dataset(const string &image_filename, const string &label_filename, std::vector< std::vector<cv::Mat> > &vvMat) {
    vvMat.resize(10);
    // Open files
    std::ifstream image_file(image_filename, std::ios::in | std::ios::binary);
    std::ifstream label_file(label_filename, std::ios::in | std::ios::binary);
    //CHECK(image_file) << "Unable to open file " << image_filename;
    //CHECK(label_file) << "Unable to open file " << label_filename;
    // Read the magic and the meta data
    uint32_t magic;
    uint32_t num_items;
    uint32_t num_labels;
    uint32_t rows;
    uint32_t cols;
    
    image_file.read(reinterpret_cast<char*>(&magic), 4);
    magic = swap_endian(magic);
    //CHECK_EQ(magic, 2051) << "Incorrect image file magic.";
    label_file.read(reinterpret_cast<char*>(&magic), 4);
    magic = swap_endian(magic);
    //CHECK_EQ(magic, 2049) << "Incorrect label file magic.";
    image_file.read(reinterpret_cast<char*>(&num_items), 4);
    num_items = swap_endian(num_items);
    label_file.read(reinterpret_cast<char*>(&num_labels), 4);
    num_labels = swap_endian(num_labels);
    std::cout << "num_items " << num_items << ", num_labels " << num_labels << std::endl;
    //CHECK_EQ(num_items, num_labels);
    image_file.read(reinterpret_cast<char*>(&rows), 4);
    rows = swap_endian(rows);
    image_file.read(reinterpret_cast<char*>(&cols), 4);
    cols = swap_endian(cols);
    std::cout << "rows " << rows << ", cols " << cols << std::endl;
    
    // Storing to db
    char label;
    char* pixels = new char[rows * cols];
    
    for (int item_id = 0; item_id < num_items; ++item_id) {
        image_file.read(pixels, rows * cols);
        label_file.read(&label, 1);
        cv::Mat img = str2Mat(pixels, rows, cols);
        //std::cout << int(label) << std::endl;
        //print_mat(img);
        //string winname = "0";
        //winname[0] += label;
        //imshow(winname, img);
        //waitKey();
        //destroyAllWindows();
        vvMat[label].push_back(img);
        //break;
    }
    
    delete[] pixels;
}



inline double sigmoid(double z) {
    return 1.0 / (1.0 + exp(-z));
}
//sigmoid函数的导数
inline double sigmoid_prime(double z) {
    return sigmoid(z) * (1 - sigmoid(z));
}


const int n_in      = 784;
const int n_hidden  = 30;
const int n_out     = 10;



double w1[n_hidden][n_in];
double w2[n_out][n_hidden];
double b1[n_hidden];
double b2[n_out];
double delta_w1[n_hidden];
double delta_w2[n_out];
double delta_b1[n_hidden];
double delta_b2[n_out];
double eta = 0.01;

void init_wb();
void init_delta();
double evaluate(const vector<vector<Mat>> &vvMatTest, const int num = 1000);
void forward(double* x, double* y, double* h);
void Mat8uc1_darr(const Mat& m, double *arr);

void print_mat(const cv::Mat &m);
void print_array(const double* arr, int size1, int size2 = 0);

void sgd(const vector<vector<Mat>> &vvMatTrain, const int numTrain, const int beginIndex = 0);

int main() {
    string data_path = "/Users/zhangxin/datas/mnist/";
    string images_test = "t10k-images-idx3-ubyte";
    string labels_test = "t10k-labels-idx1-ubyte";
    string images_train = "train-images-idx3-ubyte";
    string labels_train = "train-labels-idx1-ubyte";
    std::vector< std::vector<cv::Mat> > vvMat;
    std::vector< std::vector<cv::Mat> > vvMatTest;
    convert_dataset(data_path + images_test,  data_path + labels_test, vvMatTest);
    convert_dataset(data_path + images_train, data_path + labels_train, vvMat);

    init_wb();
    evaluate(vvMatTest, 2000);
    
    for (int i = 0; i < 1000; i++) {
        printf("%4d --------
", i);
        sgd(vvMat, 100, i * 100);
        evaluate(vvMatTest, 2000);
    }
    //evaluate(vvMatTest, 1000);

    return 0;
    
}

void sgd(const vector<vector<Mat>> &vvMatTrain, const int numTrain, const int beginIndex){
    init_delta();
   for (int j = 0; j < numTrain ; j++) {
        for (int i = 0; i < n_out; i++) {
            double y[n_out] = {0.0};
            double x[n_in] = {0.0};
            double h[n_hidden] = {0.0};
            Mat8uc1_darr(vvMatTrain[i][(beginIndex + j) % vvMatTrain[i].size()], x);
            forward(x, y, h);
            
            double t[n_out] = {0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0};
            t[i] = 1.0;
            
            for (int m = 0; m < n_out; m++) {
                delta_w2[m] = delta_b2[m] = y[m] * (1-y[m]) * (t[m] - y[m]);
                b2[m] += eta * delta_b2[m];
                
                for (int n = 0; n < n_hidden; n++) {
                    w2[m][n] += eta * delta_w2[m] * h[n];
                }
            }
            
            for (int m = 0; m < n_hidden; m++) {
                delta_b1[m] = 0.0;
                for (int k = 0; k < n_out; k++) {
                    delta_b1[m] += w2[k][m] * delta_w2[k];
                }
                delta_w1[m] = delta_b1[m] = delta_b1[m] * h[m] * (1 - h[m]);
                b1[m] += eta * delta_b1[m];

                for (int n = 0; n < n_in; n++) {
                    w1[m][n] += eta * delta_w1[m] * x[n];
                }
                
            }
        }
    }
}


double evaluate(const vector<vector<Mat>> &vvMatTest, const int num){
    int sum = 0;
    int sum_right = 0;
    vector<int> vec_sum(10, 0);
    vector<int> vec_sum_right(10, 0);
    int result[10][10];
    memset(result, 0, sizeof(int) * 100);
    
    for (int i = 0; i < n_out && i < vvMatTest.size(); i++) {
        for (int j = 0; j < num && j < (int)vvMatTest[i].size(); j++) {//
            double x[n_in];
            double y[n_out];
            double h[n_hidden];
            const Mat &img = vvMatTest[i][j];
            
            Mat8uc1_darr(img, x);
//            print_mat(img);
//            print_array(x, img.rows, img.cols);
//            imshow("00", img);
//            waitKey();
//            continue;
            
            forward(x, y, h);
            int predict_label = 0;
            double predict_f = y[0];//y.at<double>(0, 0);
            //printf("%2d %4d :[ 0] %6.3f, ", i, j, predict_f);
            for (int k = 1; k < n_out; k++) {
                double f = y[k];//y.at<double>(0, k);
                //printf("%d %6.3f, ", k, f);
                if (f  > predict_f ) {
                    predict_f = f;
                    predict_label = k;
                }
            }
            //printf(" predict : %d %6.3f
", predict_label, predict_f);
            
            vec_sum[i]++;
            sum++;
            if (predict_label == i) {
                sum_right++;
                vec_sum_right[i]++;
            }
            result[i][predict_label]++;
        }
    }
    std::cout << "sum : " << sum << ", sum_right : " << sum_right << std::endl;
    for (int i = 0; i < n_out; i++) {
        printf(" [%d] [%4d %4d] ", i, vec_sum[i], vec_sum_right[i]);
        for (int j = 0; j < n_out; j++) {
            printf ("%4d ", result[i][j]);
        }
        printf("
");
    }
    return sum_right / double(sum);
}

void forward(double* x, double* y, double* h){
    
    for (int i = 0; i < n_hidden; i++) {
        h[i] = b1[i];
        for (int j = 0; j < n_in; j++) {
            h[i] += x[j] * w1[i][j];
            
        }
        h[i] = sigmoid(h[i]);
    }
    
    for (int i = 0; i < n_out; i++) {
        y[i] = b2[i];
        for (int j = 0; j < n_hidden; j++) {
            y[i] += h[j] * w2[i][j];
        }
        y[i] = sigmoid(y[i]);
    }
    
}

void Mat8uc1_darr(const Mat& m, double *arr) {
    for (int y = 0; y < m.rows; y++) {
        for (int x = 0; x < m.cols; x++) {
            arr[y * m.cols + x] = m.at<uchar>(y,x) / 255.0;
        }
    }
}


void init_delta(){
    memset(delta_w1, 0, sizeof(double) * n_hidden);
    memset(delta_w2, 0, sizeof(double) * n_out);
    memset(delta_b1, 0, sizeof(double) * n_hidden);
    memset(delta_b2, 0, sizeof(double) * n_out);
}


void init_wb() {
    std::srand((unsigned)std::time(0)); // use current time as seed for random generator
    cout << " init w1 : " << endl;
    for (int i = 0; i < n_hidden; i++) {
        cout << " w1[" << i << "]" << endl;
        for (int j = 0; j < n_in; j++) {
            w1[i][j] = std::rand()/double(RAND_MAX) - 0.5;
            w1[i][j] *= 0.1;
            printf("%6.3f ", w1[i][j]);
        }
        printf("
");
    }
    cout << " init w2 : " << endl;
    for (int i = 0; i < n_out; i++) {
        cout << " w1[" << i << "]" << endl;
        for (int j = 0; j < n_hidden; j++) {
            w2[i][j] = std::rand()/double(RAND_MAX) - 0.5;
            w2[i][j] *= 0.1;
            printf("%6.3f ", w2[i][j]);
        }
        printf("
");
    }
    cout << "init b1" << endl;
    for (int i = 0; i < n_hidden; i++) {
        b1[i] = std::rand() / double(RAND_MAX) - 0.5;
        b1[i] *= 0.1;
        printf("%6.3f ", b1[i]);
    }
    printf("
");
    
    cout << "init b2" << endl;
    for (int i = 0; i < n_out; i++) {
        b2[i] = std::rand() / double(RAND_MAX) - 0.5;
        b2[i] *= 0.2;
        printf("%6.3f ", b2[i]);
    }
    printf("
");
    
}



void print_mat(const cv::Mat &m) {
    for (int y = 0; y < m.rows; y++) {
        printf("[%4d] ", y);
        for (int x = 0; x < m.cols; x++) {
            switch(m.type()){
                case CV_8UC1:
                    printf("%3d ", m.data[y * m.step + x]); break;
                case CV_64FC1:
                    printf("%6.3f ", m.at<double>(y,x));    break;
            }
        }
        std::cout << std::endl;
    }
}
void print_array(const double* arr, int size1, int size2) {
    for (int i = 0; i < size1; i++) {
        if( size2 > 1){
            printf("[%4d] ", i);
            for (int j = 0; j < size2; j++) {
                printf("%6.3f ", arr[i * size2 + j]);
            }
            printf("
");
        }else {
            printf("%6.3f ", arr[i]);
        }
    }
    printf("
");
}

在迭代了1000次后得到 准确率94.85%， 输出结果：

999 --------
sum : 10000, sum_right : 9485
 [0] [ 980  965]  965    0    0    1    0    2    6    4    2    0 
 [1] [1135 1114]    0 1114    3    3    0    1    4    2    8    0 
 [2] [1032  961]   10    3  961    4    8    2    8   13   21    2 
 [3] [1010  952]    0    0   15  952    1   15    1   10   14    2 
 [4] [ 982  929]    1    2    3    2  929    1   10    2    6   26 
 [5] [ 892  822]    6    1    3   23    3  822   13    4   10    7 
 [6] [ 958  921]   12    3    4    1    5    8  921    0    4    0 
 [7] [1028  967]    4   10   17    3    8    1    0  967    3   15 
 [8] [ 974  918]    5    3    5   14    6    8    5    8  918    2 
 [9] [1009  936]   10    7    1   14   18    7    1    9    6  936