C言語で０からニューラルネットワークを自作してみた

2023/06/23追記：
頂いたコメントを元にコードを修正しました。主に直した点は下記の通りです。
ーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーー
1.静的なメモリを用いることでMNIST全てを用いて学習が行えるようにしました。
（通常のメモリだとオーバーフローしてしまいますが、静的にメモリを取ることで大容量のデータを扱えるようになりました）
2.変数の初期化を明示しました。
3.正解率の求め方を改めました（間違った計算をしていました）
その他軽微なミスを修正しました。
ーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーーー




こんにちにゃんです。
水色桜（みずいろさくら）です。
今回はC言語を用いてニューラルネットワークを自作して、MNISTを分類してみようと思います。
今回作成したモデルの精度を以下に示します。
10回の学習で55.0％の精度を達成できています。
```
1回目計算終了
error : 0.183634
accuracy : 37.899647
2回目計算終了
error : 0.163894
accuracy : 42.598858
3回目計算終了
error : 0.157277
accuracy : 46.738163
4回目計算終了
error : 0.156364
accuracy : 48.517864
5回目計算終了
error : 0.153185
accuracy : 51.317394
6回目計算終了
error : 0.149610
accuracy : 52.837139
7回目計算終了
error : 0.145736
accuracy : 53.437038
8回目計算終了
error : 0.143767
accuracy : 54.056934
9回目計算終了
error : 0.143105
accuracy : 54.506859
10回目計算終了
error : 0.142516
accuracy : 55.116756
```
記事中で何か不明な点、間違いなどありましたら、コメントいただけると嬉しいです。


# C言語とは
1972年にAT&Tベル研究所のデニス・リッチーらが主体となって開発した汎用プログラミング言語。コンパイラ言語であるため、ほかのプログラミング言語と比べて高速に動作します。（Pythonのように機械学習のライブラリが充実していないため、機械学習を行うためにはフルスクラッチする必要があります。そのため、C言語は機械学習にはほとんど使われません。）

# [MNIST](https://github.com/pjreddie/mnist-csv-png)とは
[MNIST(Mixed National Institute of Standards and Technology database)](https://github.com/pjreddie/mnist-csv-png)とは、手書き数字画像60,000枚と、テスト画像10,000枚を集めた、画像データセットです。さらに、手書きの数字「0〜9」に正解ラベルが与えられるデータセットでもあり、画像分類問題で人気の高いデータセットです。

# ニューラルネットワークとは

人間の脳のしくみ（ニューロン間のあらゆる相互接続）から着想を得たもので、脳機能の特性のいくつかをコンピュータ上で表現するために作られた数学モデルです。（[Udemyメディアから引用](https://udemy.benesse.co.jp/data-science/ai/neural-network.html)）
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2847349/f1212ea0-7f28-68e2-557a-ae1aae032711.png)
（[Udemyメディアから引用](https://udemy.benesse.co.jp/data-science/ai/neural-network.html)）

ニューラルネットワークは入力層、中間層、出力層の３つから構成されます。本記事で作成するニューラルネットワークのモデルは次のような感じです。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2847349/cd4328ec-fb94-a8a0-c9f7-5007c67dcecc.png)

ニューラルネットワークにおいては、前の層からの出力に重みとバイアスを加え、それをアクティベーション関数と呼ばれる非線形関数にかけます。これにより、複雑な分布も表現することが可能になります。ニューラルネットワークは学習を行うことで出力層で人間が望む結果（正しい答え、正解）が出るようにします。学習によって中間層の重みとバイアスを最適化していきます。

# 数学的な背景

ニューラルネットワークを構成するうえで、中間層と出力層のアクティベーション関数を選ぶ必要があります。本記事では中間層のアクティベーション関数はLeaky ReLU、出力層のアクティベーション関数はsoftmax関数を用います。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2847349/8011635f-3f11-6fb3-7a43-5b270897e6e2.png)
Leaky ReLU関数は先行研究において優れた成果を出している関数であり、勾配消失（勾配が小さくなり、ニューラルネットワークの学習が進まなくなってしまう問題）が起きにくいという性質を有しています。
softmax関数は出力をすべて足し合わせると１になるように設計された関数であり、分類問題に用いられている関数です。入力に対して０～１の任意の値を返すため、擬似的に確率として扱うことができます。
交差エントロピー誤差はsoftmax関数と一緒に用いると逆伝播が簡単な形に書けるように設計された関数です。一見複雑そうに見えますが、微分すると簡単な形になります。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2847349/d1798e18-17e2-cb5d-530e-7c74a1b9f27d.png)

誤差逆伝播法はニューラルネットワークにおいて勾配を効率的に求めるための手法です。微分の連鎖率を用いて図のように逆側から勾配を求めていきます。
誤差逆伝播法の導出方法は下の画像の通りです。交差エントロピー誤差を出力で微分した値と、出力を入力で微分した値を掛け合わせると、結果としてとっても簡単な形に書き表すことができます。出力から教師データを引いたものになるので、実装もしやすくなっています（そうなるようにsoftmax関数と交差エントロピー誤差は設計されています。）
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/2847349/b0948ed4-c25a-661a-d57e-91126394d224.png)

# 実装
今回実装したコード全体を以下（neural_network.c）に示します。お忙しい方はコードをコピペして使ってみてください（C言語をコンパイルできる環境が必要です。）。

<details><summary>neural_network.c</summary>

```C:neural_network.c
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>

#define max_num_train 50000
#define max_num_test 10000
#define M1 32
#define M2 32 /*隠れ層の厚さ*/
#define N 100  /*学習の繰り返し回数*/
#define rate 0.03f  /*学習率*/

float activ(float x){
    return 1.0/(1.0+exp(-x));
}

float activ_deri(float x){
    return activ(x)*(1.0-activ(x));
}

int main(void)
{
/*ファイルの読み込み*/
FILE *fp;
fp=fopen("mnist_train.txt","r");
if(fp==NULL){
printf("file open error\n");
}
else{
printf("file opened\n");
}

char line[2000]={};

/*データを配列に格納*/
/*計算量が大きいのでmax_numの数だけ使用する*/
static float num[max_num_train][785];
const char *list;
for(int j=0;j<max_num_train;j++){
    fgets(line, 2000, fp);
    list=strtok(line,",");
    num[j][0]=atof(list+0);
    for(int i=1;i<785;i++){
    num[j][i]=atof((list+i))/255.0;
    }
    
}

fclose(fp);

/*ファイルの読み込み*/
FILE *fp_test;
fp_test=fopen("mnist_test.txt","r");
if(fp_test==NULL){
printf("file open error\n");
}
else{
printf("file opened\n");
}

char line_test[2000]={};
const char *list_test;
/*データを配列に格納*/
/*計算量が大きいのでmax_numの数だけ使用する*/
static float num_test[max_num_test][785];
for(int j=0;j<max_num_test;j++){
    fgets(line_test, 2000, fp_test);
    list_test=strtok(line_test,",");
    num_test[j][0]=atof(list_test+0);
    for(int i=1;i<785;i++){
    num_test[j][i]=atof((list_test+i))/255.0;
    }
}

fclose(fp_test);

static float W1[784][M1]={};  /*一層目の重み*/
static float W2[M1][M2]={};  /*二層目の重み*/
static float W3[M2][10]={};  /*三層目の重み*/

/*重みをランダムな数で初期化する*/
for(int i=0;i<784;i++){
    for(int j=0;j<M1;j++){
        W1[i][j]=((float)rand()/RAND_MAX*2.0f-1.0f)/10.0f;
    }
}
for(int i=0;i<M1;i++){
    for(int j=0;j<M2;j++){
        W2[i][j]=((float)rand()/RAND_MAX*2.0f-1.0f)/10.0f;
    }
}
for(int i=0;i<M2;i++){
    for(int j=0;j<10;j++){
        W3[i][j]=((float)rand()/RAND_MAX*2.0f-1.0f)/10.0f;
    }
}

/*正解データをone-hotベクトル化する*/
static float ans[10][max_num_train]={};
for(int j=0;j<max_num_train;j++){
    if(fabsf(num[j][0]-0)<1){
        ans[0][j]=1.0f;
    }
    if(fabsf(num[j][0]-1)<1){
        ans[1][j]=1.0f;
    }
    if(fabsf(num[j][0]-2)<1){
        ans[2][j]=1.0f;
    }
    if(fabsf(num[j][0]-3)<1){
        ans[3][j]=1.0f;
    }
    if(fabsf(num[j][0]-4)<1){
        ans[4][j]=1.0f;
    }
    if(fabsf(num[j][0]-5)<1){
        ans[5][j]=1.0f;
    }
    if(fabsf(num[j][0]-6)<1){
        ans[6][j]=1.0f;
    }
    if(fabsf(num[j][0]-7)<1){
        ans[7][j]=1.0f;
    }
    if(fabsf(num[j][0]-8)<1){
        ans[8][j]=1.0f;
    }
    if(fabsf(num[j][0]-9)<1){
        ans[9][j]=1.0f;
    }
}

static float ans_test[10][max_num_test]={};
for(int j=0;j<max_num_test;j++){
    if(fabsf(num_test[j][0]-0)<1){
        ans_test[0][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-1)<1){
        ans_test[1][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-2)<1){
        ans_test[2][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-3)<1){
        ans_test[3][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-4)<1){
        ans_test[4][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-5)<1){
        ans_test[5][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-6)<1){
        ans_test[6][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-7)<1){
        ans_test[7][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-8)<1){
        ans_test[8][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-9)<1){
        ans_test[9][j]=1.0f;
    }
    
}

/*ここから順伝播と逆伝播を実装していく*/
for(int o=0;o<N;o++){
static float z1[M1][max_num_train]={};  /*一層目の出力を活性化関数に入力したもの*/
static float z2[M2][max_num_train]={};  /*二層目の出力を活性化関数に入力したもの*/
static float z3[10][max_num_train]={};  /*三層目の出力を活性化関数に入力したもの*/
       
/*ここからは順伝播の実装*/
for(int n=0;n<max_num_train;n++){
    float b1[M1]={};  /*一層目のバイアス*/
    float b2[M2]={};  /*二層目のバイアス*/
    float b3[10]={};  /*三層目のバイアス*/
    static float a1[M1][max_num_train]={};  /*一層目の出力*/
    static float a2[M2][max_num_train]={};  /*二層目の出力*/
    static float a3[10][max_num_train]={};  /*三層目の出力*/

    for(int l=0;l<784;l++){
        for(int m=0;m<M1;m++){
            a1[m][n]=a1[m][n]+num[n][l+1]*W1[l][m]+b1[m];  /*一層目の計算*/
        }
    }
    for(int i=0;i<M1;i++){
        z1[i][n]=activ(a1[i][n]);  /*活性化関数sigmoidを適用する*/
        }
    
    for(int l=0;l<M1;l++){
        for(int m=0;m<M2;m++){
            a2[m][n]=a2[m][n]+z1[l][n]*W2[l][m]+b2[m];  /*二層目の計算*/
        }
    }
    for(int i=0;i<M2;i++){
        z2[i][n]=activ(a2[i][n]);  /*活性化関数sigmoidを適用する*/
        }

    for(int m=0;m<10;m++){
        for(int l=0;l<M2;l++){
            a3[m][n]=a3[m][n]+z2[l][n]*W3[l][m]+b3[m];  /*三層目の計算*/
        }
    }
    
    /*ここからSOFTMAX関数の実装*/
    float u=0.00001f;
    float t=-RAND_MAX;
    for(int i=0;i<10;i++){
        if(t<a3[i][n]){
            t=a3[i][n];
        }
    }
    for(int i=0;i<10;i++){
        u += exp(a3[i][n]-t);  /*数値がオーバーフローしないように最大値tのexpで割る*/
    } 
    for(int i=0;i<10;i++){
        z3[i][n] = exp(a3[i][n]-t)/(u+0.00001f);  /*数値がオーバーフローしないように最大値tのexpで割る*/
    }
    /*ここから逆伝播の実装*/
    float delta3[10]={};
    for(int j=0;j<10;j++){
        delta3[j]=(z3[j][n]-ans[j][n])*activ_deri(z3[j][n]);
    }

    float delta2[M1]={};
    for (int i = 0; i < M2; i++) {
        double error = 0.0f;
        for (int j = 0; j < 10; j++) {
            error += W3[i][j] * delta3[j];
        }
        delta2[i] = error * activ_deri(z2[i][n]);
    }
    float delta1[M1]={};
    for (int i = 0; i < M1; i++) {
        double error = 0.0;
        for (int j = 0; j < M2; j++) {
            error += W2[i][j] * delta2[j];
        }
        delta1[i] = error * activ_deri(z1[i][n]);
    }
    
    // 重みの更新
    for (int i = 0; i < 10; i++) {
        for (int j = 0; j < M2; j++) {
            W3[j][i] -= rate * delta3[i] * z2[j][n];
        }
        b3[i] -= rate*delta3[i];
    }
    
    for (int i = 0; i < M2; i++) {
        for (int j = 0; j < M1; j++) {
            W2[j][i] -= rate * delta2[i] * z1[j][n];
        }
        b2[i] -= rate*delta2[i];
    }
    
    for (int i = 0; i < M1; i++) {
        for (int j = 0; j < 784; j++) {
            W1[j][i] -= rate * delta1[i] * num[n][j+1];
        }
        b1[i] -= rate*delta1[i];
    }
}
static float z3_test[10][max_num_test]={};  /*三層目の出力をSOFTMAX関数に入力したもの*/
static float z2_test[M2][max_num_test]={};  /*二層目の出力を活性化関数に入力したもの*/
static float z1_test[M1][max_num_test]={};  /*一層目の出力を活性化関数に入力したもの*/
float b1[M1]={};  /*一層目のバイアス*/
float b2[M2]={};  /*二層目のバイアス*/
float b3[10]={};  /*三層目のバイアス*/
    
/*ここからは更新したパラメータを用いて順伝播を計算し、誤差と正解率を求める*/
for(int n=0;n<max_num_test;n++){
    static float a1_test[M1][max_num_test]={};  /*一層目の出力*/
    static float a2_test[M2][max_num_test]={};  /*二層目の出力*/
    static float a3_test[10][max_num_test]={};  /*三層目の出力*/
    
    for(int m=0;m<M1;m++){
        for(int l=0;l<784;l++){
            a1_test[m][n]=a1_test[m][n]+num_test[n][l+1]*W1[l][m]+b1[m];
        }
    }
    
    for(int i=0;i<M1;i++){
        z1_test[i][n]=activ(a1_test[i][n]);
    }

    for(int m=0;m<M1;m++){
        for(int l=0;l<M2;l++){
            a2_test[l][n]=a2_test[l][n]+z1_test[m][n]*W2[m][l]+b2[l];
        }
    }
    
    for(int i=0;i<M2;i++){
        z2_test[i][n]=activ(a2_test[i][n]);
    }


    for(int m=0;m<10;m++){
        for(int l=0;l<M2;l++){
            a3_test[m][n]=a3_test[m][n]+z2_test[l][n]*W3[l][m]+b3[m];
        }
    }


    float u=0.0001f;
    float t=-RAND_MAX;
    for(int i=0;i<10;i++){
        if(t<a3_test[i][n]){
            t=a3_test[i][n];
        }
    }
    for(int i=0;i<10;i++){
        u += exp(a3_test[i][n]-t);
    } 
    for(int i=0;i<10;i++){
        z3_test[i][n] = exp(a3_test[i][n]-t)/(u+0.000001f);
    }
}

/*交差エントロピー誤差を求める*/
float score=0.0f;
for(int i=0;i<max_num_test;i++){
    for(int j=0;j<10;j++){
    score=score-ans_test[j][i]*log(z3_test[j][i]+0.000001f)/max_num_test/10;
    }
}

/*正解率を求める*/
float acc=0.0f;
for(int j=0;j<max_num_test;j++){
        int zindex=0;
        for(int i=1;i<10;i++){
            if(z3_test[zindex][j]<z3_test[i][j]){
                zindex=i;
            }
        }
        if(fabsf(zindex-num_test[j][0])<1){
            acc=acc+100.0f/max_num_test;
        }
}


/*結果の表示*/
printf("%d回目計算終了\n",o+1);
printf("error : %f\n",score);
printf("accuracy : %f\n",acc);

}
return 0;
}




```

</details>

まず、ファイルのロードを行い、データを配列に格納します。

```C:file_load.c
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>

#define max_num_train 50000
#define max_num_test 10000
#define M1 32
#define M2 32 /*隠れ層の厚さ*/
#define N 100  /*学習の繰り返し回数*/
#define rate 0.03f  /*学習率*/

float activ(float x){
    return 1.0/(1.0+exp(-x));
}

float activ_deri(float x){
    return activ(x)*(1.0-activ(x));
}

int main(void)
{
/*ファイルの読み込み*/
FILE *fp;
fp=fopen("mnist_train.txt","r");
if(fp==NULL){
printf("file open error\n");
}
else{
printf("file opened\n");
}

char line[2000]={};

/*データを配列に格納*/
/*計算量が大きいのでmax_numの数だけ使用する*/
static float num[max_num_train][785];
const char *list;
for(int j=0;j<max_num_train;j++){
    fgets(line, 2000, fp);
    list=strtok(line,",");
    num[j][0]=atof(list+0);
    for(int i=1;i<785;i++){
    num[j][i]=atof((list+i))/255.0;
    }
    
}

fclose(fp);

/*ファイルの読み込み*/
FILE *fp_test;
fp_test=fopen("mnist_test.txt","r");
if(fp_test==NULL){
printf("file open error\n");
}
else{
printf("file opened\n");
}

char line_test[2000]={};
const char *list_test;
/*データを配列に格納*/
/*計算量が大きいのでmax_numの数だけ使用する*/
static float num_test[max_num_test][785];
for(int j=0;j<max_num_test;j++){
    fgets(line_test, 2000, fp_test);
    list_test=strtok(line_test,",");
    num_test[j][0]=atof(list_test+0);
    for(int i=1;i<785;i++){
    num_test[j][i]=atof((list_test+i))/255.0;
    }
}

fclose(fp_test);

```

次に、各種変数を定義し、各層の重みを乱数で初期化します。
正解データ（０～９）は後の計算で用いれるようにone-hotベクトル化します。
（one-hotベクトルとは、一つの要素が１でそれ以外の要素が０であるようなベクトルのこと。）
```C:rand.c
static float W1[784][M1]={};  /*一層目の重み*/
static float W2[M1][M2]={};  /*二層目の重み*/
static float W3[M2][10]={};  /*三層目の重み*/

/*重みをランダムな数で初期化する*/
for(int i=0;i<784;i++){
    for(int j=0;j<M1;j++){
        W1[i][j]=((float)rand()/RAND_MAX*2.0f-1.0f)/10.0f;
    }
}
for(int i=0;i<M1;i++){
    for(int j=0;j<M2;j++){
        W2[i][j]=((float)rand()/RAND_MAX*2.0f-1.0f)/10.0f;
    }
}
for(int i=0;i<M2;i++){
    for(int j=0;j<10;j++){
        W3[i][j]=((float)rand()/RAND_MAX*2.0f-1.0f)/10.0f;
    }
}

/*正解データをone-hotベクトル化する*/
static float ans[10][max_num_train]={};
for(int j=0;j<max_num_train;j++){
    if(fabsf(num[j][0]-0)<1){
        ans[0][j]=1.0f;
    }
    if(fabsf(num[j][0]-1)<1){
        ans[1][j]=1.0f;
    }
    if(fabsf(num[j][0]-2)<1){
        ans[2][j]=1.0f;
    }
    if(fabsf(num[j][0]-3)<1){
        ans[3][j]=1.0f;
    }
    if(fabsf(num[j][0]-4)<1){
        ans[4][j]=1.0f;
    }
    if(fabsf(num[j][0]-5)<1){
        ans[5][j]=1.0f;
    }
    if(fabsf(num[j][0]-6)<1){
        ans[6][j]=1.0f;
    }
    if(fabsf(num[j][0]-7)<1){
        ans[7][j]=1.0f;
    }
    if(fabsf(num[j][0]-8)<1){
        ans[8][j]=1.0f;
    }
    if(fabsf(num[j][0]-9)<1){
        ans[9][j]=1.0f;
    }
}

static float ans_test[10][max_num_test]={};
for(int j=0;j<max_num_test;j++){
    if(fabsf(num_test[j][0]-0)<1){
        ans_test[0][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-1)<1){
        ans_test[1][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-2)<1){
        ans_test[2][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-3)<1){
        ans_test[3][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-4)<1){
        ans_test[4][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-5)<1){
        ans_test[5][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-6)<1){
        ans_test[6][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-7)<1){
        ans_test[7][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-8)<1){
        ans_test[8][j]=1.0f;
    }
    if(fabsf(num_test[j][0]-9)<1){
        ans_test[9][j]=1.0f;
    }
    
}
```

次に、順伝播とSOFTMAX関数を実装します。
```C:SOFTMAX.c
/*ここから順伝播と逆伝播を実装していく*/
for(int o=0;o<N;o++){
static float z1[M1][max_num_train]={};  /*一層目の出力を活性化関数に入力したもの*/
static float z2[M2][max_num_train]={};  /*二層目の出力を活性化関数に入力したもの*/
static float z3[10][max_num_train]={};  /*三層目の出力を活性化関数に入力したもの*/
       
/*ここからは順伝播の実装*/
for(int n=0;n<max_num_train;n++){
    float b1[M1]={};  /*一層目のバイアス*/
    float b2[M2]={};  /*二層目のバイアス*/
    float b3[10]={};  /*三層目のバイアス*/
    static float a1[M1][max_num_train]={};  /*一層目の出力*/
    static float a2[M2][max_num_train]={};  /*二層目の出力*/
    static float a3[10][max_num_train]={};  /*三層目の出力*/

    for(int l=0;l<784;l++){
        for(int m=0;m<M1;m++){
            a1[m][n]=a1[m][n]+num[n][l+1]*W1[l][m]+b1[m];  /*一層目の計算*/
        }
    }
    for(int i=0;i<M1;i++){
        z1[i][n]=activ(a1[i][n]);  /*活性化関数sigmoidを適用する*/
        }
    
    for(int l=0;l<M1;l++){
        for(int m=0;m<M2;m++){
            a2[m][n]=a2[m][n]+z1[l][n]*W2[l][m]+b2[m];  /*二層目の計算*/
        }
    }
    for(int i=0;i<M2;i++){
        z2[i][n]=activ(a2[i][n]);  /*活性化関数sigmoidを適用する*/
        }

    for(int m=0;m<10;m++){
        for(int l=0;l<M2;l++){
            a3[m][n]=a3[m][n]+z2[l][n]*W3[l][m]+b3[m];  /*三層目の計算*/
        }
    }
    
    /*ここからSOFTMAX関数の実装*/
    float u=0.00001f;
    float t=-RAND_MAX;
    for(int i=0;i<10;i++){
        if(t<a3[i][n]){
            t=a3[i][n];
        }
    }
    for(int i=0;i<10;i++){
        u += exp(a3[i][n]-t);  /*数値がオーバーフローしないように最大値tのexpで割る*/
    } 
    for(int i=0;i<10;i++){
        z3[i][n] = exp(a3[i][n]-t)/(u+0.00001f);  /*数値がオーバーフローしないように最大値tのexpで割る*/
    }
```

次に、誤差逆伝播法を用いて重みの更新を行っていきます。
```C:back.c
/*ここから逆伝播の実装*/
    float delta3[10]={};
    for(int j=0;j<10;j++){
        delta3[j]=(z3[j][n]-ans[j][n])*activ_deri(z3[j][n]);
    }

    float delta2[M1]={};
    for (int i = 0; i < M2; i++) {
        double error = 0.0f;
        for (int j = 0; j < 10; j++) {
            error += W3[i][j] * delta3[j];
        }
        delta2[i] = error * activ_deri(z2[i][n]);
    }
    float delta1[M1]={};
    for (int i = 0; i < M1; i++) {
        double error = 0.0;
        for (int j = 0; j < M2; j++) {
            error += W2[i][j] * delta2[j];
        }
        delta1[i] = error * activ_deri(z1[i][n]);
    }
    
    // 重みの更新
    for (int i = 0; i < 10; i++) {
        for (int j = 0; j < M2; j++) {
            W3[j][i] -= rate * delta3[i] * z2[j][n];
        }
        b3[i] -= rate*delta3[i];
    }
    
    for (int i = 0; i < M2; i++) {
        for (int j = 0; j < M1; j++) {
            W2[j][i] -= rate * delta2[i] * z1[j][n];
        }
        b2[i] -= rate*delta2[i];
    }
    
    for (int i = 0; i < M1; i++) {
        for (int j = 0; j < 784; j++) {
            W1[j][i] -= rate * delta1[i] * num[n][j+1];
        }
        b1[i] -= rate*delta1[i];
    }
}
```


最後に交差エントロピー誤差と正解率を求めます。
```C:error.c
static float z3_test[10][max_num_test]={};  /*三層目の出力をSOFTMAX関数に入力したもの*/
static float z2_test[M2][max_num_test]={};  /*二層目の出力を活性化関数に入力したもの*/
static float z1_test[M1][max_num_test]={};  /*一層目の出力を活性化関数に入力したもの*/
float b1[M1]={};  /*一層目のバイアス*/
float b2[M2]={};  /*二層目のバイアス*/
float b3[10]={};  /*三層目のバイアス*/
    
/*ここからは更新したパラメータを用いて順伝播を計算し、誤差と正解率を求める*/
for(int n=0;n<max_num_test;n++){
    static float a1_test[M1][max_num_test]={};  /*一層目の出力*/
    static float a2_test[M2][max_num_test]={};  /*二層目の出力*/
    static float a3_test[10][max_num_test]={};  /*三層目の出力*/
    
    for(int m=0;m<M1;m++){
        for(int l=0;l<784;l++){
            a1_test[m][n]=a1_test[m][n]+num_test[n][l+1]*W1[l][m]+b1[m];
        }
    }
    
    for(int i=0;i<M1;i++){
        z1_test[i][n]=activ(a1_test[i][n]);
    }

    for(int m=0;m<M1;m++){
        for(int l=0;l<M2;l++){
            a2_test[l][n]=a2_test[l][n]+z1_test[m][n]*W2[m][l]+b2[l];
        }
    }
    
    for(int i=0;i<M2;i++){
        z2_test[i][n]=activ(a2_test[i][n]);
    }


    for(int m=0;m<10;m++){
        for(int l=0;l<M2;l++){
            a3_test[m][n]=a3_test[m][n]+z2_test[l][n]*W3[l][m]+b3[m];
        }
    }


    float u=0.0001f;
    float t=-RAND_MAX;
    for(int i=0;i<10;i++){
        if(t<a3_test[i][n]){
            t=a3_test[i][n];
        }
    }
    for(int i=0;i<10;i++){
        u += exp(a3_test[i][n]-t);
    } 
    for(int i=0;i<10;i++){
        z3_test[i][n] = exp(a3_test[i][n]-t)/(u+0.000001f);
    }
}

/*交差エントロピー誤差を求める*/
float score=0.0f;
for(int i=0;i<max_num_test;i++){
    for(int j=0;j<10;j++){
    score=score-ans_test[j][i]*log(z3_test[j][i]+0.000001f)/max_num_test/10;
    }
}

/*正解率を求める*/
float acc=0.0f;
for(int j=0;j<max_num_test;j++){
        int zindex=0;
        for(int i=1;i<10;i++){
            if(z3_test[zindex][j]<z3_test[i][j]){
                zindex=i;
            }
        }
        if(fabsf(zindex-num_test[j][0])<1){
            acc=acc+100.0f/max_num_test;
        }
}


/*結果の表示*/
printf("%d回目計算終了\n",o+1);
printf("error : %f\n",score);
printf("accuracy : %f\n",acc);

}
return 0;
}




```

# 終わりに
今回はC言語を用いてニューラルネットワークを作成する方法について書いてきました。
０からニューラルネットワークを作成してみたいという方の参考になれば嬉しいです。
では、ばいにゃん～。