RustでDeepLearning入門

# RustでDeepLearning 

Rustを使い慣れることを目標に[「ゼロから作るDeep Learning」](https://www.oreilly.co.jp/books/9784873117584/)を参照しながら開発しています.
今回,行列Wを設定すれば推論し、ロスを計算する所まで実装したので、その内容を説明したいと思います.

## ライブラリの使い方
「ゼロから作るDeep Learning」の1章では基本的なライブラリの使い方について説明しています. 
ここではrustの行列演算ライブラリ `ndarray` で遊びます.
ゼロから作るDeep Learningに合わせるならmatplotlib相当のものを探すべきですが、
どれもまだまだ枯れていないように見えたことと,グラフ描画はしないので,一つに絞りました.

実際ndarrayの使い方について説明していきます.

### データの生成
[ソース](https://github.com/rust-ndarray/ndarray/blob/master/src/free_functions.rs#L222)を見ると以下のメソッドを使って,
二次元のndarrayが生成できます.

```rust
pub fn arr2<A: Clone, V: FixedInitializer<Elem = A>>(xs: &[V]) -> Array2<A>
where
    V: Clone,
{
    Array2::from(xs.to_vec())
}
```

rustの文法はなれないと難しいので、少し解説します.
このメソッドの引数の型 `&[V]`,復帰値の型は `Array2<A>`です.
`&` は参照で, Vは`FixedInitializer<Elem = A`>となります.
`Elem` はndarrayのTypeで, `FixedInitializer`はそれを含むTraitです.
`FixedInitializer`はサイズ固定の場合の初期化処理をします.
`Elem`であれば特に問題がないので、intでもfloatでも特に気にせず同じ文法で扱えるようです.

実際に生成してみます

```rust
use ndarray::arr2;
let a = arr2(&[[1, 2, 3],
               [4, 5, 6]]);
println!("{:?}", a);
```

`arr2` という文字からもわかるように二次元のものです.
少し見た限り,データの生成の場合は以下がよく使うものに見えました.

- `array!`
- `Array::from_vec`
- `arr2`
- `arr1`

最初は次元を明記した方がわかりやすいと思うので,arr2を使います.

### 値の取得
二次元配列の場合は
`data[[i, j]]`でアクセスできます.それ含め3通りメソッドが提供されています.
また、範囲外だと`None`になるようです.

```rust
use ndarray::arr2;

let a = arr2(&[[1., 2.],
               [3., 4.]]);

assert!(
    a.get([0, 1]) == Some(&2.) &&
    a.get((0, 2)) == None &&
    a[(0, 1)] == 2. &&
    a[[0, 1]] == 2.
);
```

### 値のループ
ループについて説明します.

基本は以下でいけるようです.

```rust
    let a = arr3(&[[[ 50,  1,  2],
                [ 3,  4,  5]],
               [[ 6,  7,  8],
                [ 9, 10, 11]]]);
    println!("{:?}", a.iter().next().unwrap());
```

`iter`, `iter_mut` :配列の内部値の参照を取得します.実際には`Iter`オブジェクトが返ります.
正直`Iter`オブジェクトの挙動が掴みきれませんが, .next()で次の値をOption型で得られます.

単位でループする場合は以下でできます.

```rust
    let mut a = Array::zeros((10, 10));
    for mut row in a.genrows_mut() {
        println!("{:}", row);
        row.fill(1.);
        println!("{:}", row);
    }
    println!("{:}", a);
```
`genrows`,  `genrows_mut` は列をループします
ちなみに`row.fill(1.)`で値を書き換えるので, `a.cols()`はmutableでないと動きません.


複数のndarrayをループさせるときには以下のようにします.

```rust
    use ndarray::Zip;
    let mut a = Array::zeros((10, 10));
    for mut row in a.genrows_mut() {
        row.fill(1.);
    }
    let mut b = Array::zeros(a.rows());

    Zip::from(a.genrows())
        .and(&mut b)
        .apply(|a_row, b_elt| {
            *b_elt = a_row[a.cols() - 1] + a_row[0];
        });
    println!("b {:?}", b);

    let mut c = Array::zeros(10);
    Zip::from(&mut c)
        .and(&b)
        .apply(|c_elt, &b_elt| {
            *c_elt = b_elt + 1.;
        });
    println!("c {:?}", c);

    use ndarray::Array2;
    type M = Array2<f64>;
    let mut a = M::zeros((12, 8));
    let b = M::from_elem(a.dim(), 1.);
    let c = M::from_elem(a.dim(), 2.);
    let d = M::from_elem(a.dim(), 3.);


    Zip::from(&mut a)
        .and(&b)
        .and(&c)
        .and(&d)
        .apply(|w, &x, &y, &z| {
            *w += x + y * z;
        });
    println!("{:?}", a);
```

zipでは複数の配列をループさせられます.

- 個数を複数にするには`and`で加えればよいです.
- またzipのloopは`Array2`に対して実施しても,全要素に対してループすることになります.
- `apply`をするとArrayView型になります.
- `from_elem`は一つ目の引数は次元で,二つ目の引数の値で埋めます.

### スライス

```rust
let a = arr2(&[[1., 2.],
            [3., 4.]]);

prinltlin("{:?}", a.slice(s![.., 0..1])
```
- `s!`でindex相当を取り,
- `slice`で実際に指定した範囲のデータをアクセスをします.

### 行列の和、差、積、要素積
ここから実際に行列だと思って演算をします.

```rust
    let a = arr2(&[[1., 2.],
                [3., 4.]]);
    let b = arr2(&[[5., 6.],
                [3., 4.]]);
    let c = &a + &b;
    let d = &c - &b;
    let e = &d.dot(&b);
    let f = &d * e; // 要素積
    let g = &f / e; 
    println!("matrix f {}", f);
    println!("matrix g {}", g);
```



## パーセプトロン
二章では二層、三層のシンプルなパーセプトロンについて説明します.
ここで二層のパーセプトロンを具体的には$x_1, x_2$という入力に対して,$w_1, w_2, \theta$というパラメータを定め,
出力を以下で定めたものです.

```math
y = \left\{
\begin{array}{ll}
0 & w_1x_1 + w_2 x_2 < \theta \\
1 & w_1x_2 + w_2x_2 \ge \theta
\end{array}
\right.
```

二層のパーセプトロンで

- `AND`
- `OR`
- `NAND`
が実現でき,逆に
- `XOR`
が二層のパーセプトロンが実現できないことを示します.

まず`AND`ですが
$w_1 = w_2 = 0.5$とし,$\theta = 0.7$とすればよいです.
実際計算してみると$(x_1, x_2)$が(0,0),(1,0),(0,1)のときは高々0.5しか行かず,(1,1)のときは1となります.
$\theta$の条件を考えるとこれが`AND`になることがわかります.

また`OR`は$w_1 = w_2 = 0.5$とし,$\theta = 0.3$とすればよいです.
`NAND`は`AND`を反転させたものなので,
$w_1 = w_2 = -0.5$とし,$\theta = - 0.7$とすればよいです.
とすればよいです.

一方で`XOR`は実現できません.なぜなら
$(1, 0),(0, 1)$でともに
$w_1x_2 + w_2x_2 \ge \theta$を満たしている.
つまり,$w_1 \ge \theta, w_2 \ge \theta$の時,($\theta > 0$の時)$w_1 + w_2 \ge \theta$となります.
逆に$\theta \le 0$の時,$0 \ge \theta$となります.
なので
$(1, 0),(0, 1)$でともに出力が1なら$(1,1),(0, 0)$少なくとも一方が1になることがわかります.
なので、これを実現することは不可能です.

ただし、これを三層のパーセプトロンにすれば実現できます.

```math
\begin{array}{ll}
s_1 =  NAND(x_1, x_2) \\
s_2 =  OR(x_1, x_2) \\
y = AND(s_1, s_2)
\end{array}
```

とすればNANDを実現することができます,

これは三層のパーセプトロンが二層では表せない関数を作れるということで真に表現力が豊かであることがわかります.
(二層のパーセプトロン$F$に対し、$AND(F(x),F(x))=F(x)$となるので、二層のものも全て表すことができます.)

では実際にこれらをrustで実装してみましょう.

```rust
fn over(x: f64, theta: f64) -> f64 {
    return if x >= theta { 1.0 } else { 0.0 }
}

pub fn and(x: &Array1<f64>) -> f64{
    let weights = arr1(&[0.5, 0.5]);
    let theta = 0.7;
    println!("{:?}", weights);
    return over(x.dot(&weights), theta)
}

pub fn or(x: &Array1<f64>) -> f64{
    let weights = arr1(&[0.5, 0.5]);
    let theta = 0.3;
    println!("{:?}", weights);
    return over(x.dot(&weights), theta)
}

pub fn nand(x: &Array1<f64>) -> f64{
    let weights = arr1(&[-0.5, -0.5]);
    let theta = -0.7;
    println!("{:?}", weights);
    return over(x.dot(&weights), theta)
}

pub fn xor(x: &Array1<f64>) -> f64{
    let s1 = or(&x);
    let s2 = nand(&x);
    let s = arr1(&[s1, s2]);
    return and(&s);
}
```


## ニューラルネットワーク
活性化関数と多層パーセプトロンを実装し,実際に一度推論をさせてみます.

### 活性化関数

ニューラルネットワークは入力に対して,行列をかけ,活性化関数をかけを繰り返します.
具体的にある活性化関数としては

- 階段関数
- Sigmoid
- ReLu
- Softmax
等があります.

それらについて実際に見つつ,実装します.

ここで階段関数とは
$f:\mathbb{R} \to \{0, 1\}$で

```math
x \mapsto
\left \{ \begin{array}{ll}
1  &  x \ge \theta \\
0  &  otherwise \\
\end{array}
\right .
```
となる関数のことです.

rustでの実装です。

```rust
pub fn step(x: f64)-> f64 {
    let res = if x > 0.0 {1.0} else {0.0};
    println!("res {:}", res);
    return res
}
```

Sigmoidは
$$
x \mapsto \frac{1}{1 + \exp(-x)}
$$
となる関数のことです.

```rust
pub fn sigmoid(x: f64)-> f64 {
    let res = 1.0 / f64::exp(- x);
    println!("res {:}", res);
    return res
}
```

階段関数とSigmoidの入力が __一次元__ なことに注意してください.

逆に
ReLuは
$$
(x_i)^n_{i=1} \mapsto (\max \{0, x_i\})^n_{i=1}
$$

Softmaxは
$$
(x_i)^n_{i=1} \mapsto \left( \frac{\exp(x_i)}{\sum^n_{j=1} \exp(x_j)} \right)^n_{i=1}
$$

となり多次元の入力です.


実際に実装してみます.
softmaxは`exp`は値が大きくなりすぎる場合があるため,最も値の大きい$x$で分母、分子を割っていることに注意してください.

```rust
pub fn relu(x: &Array1<f64>) -> Array1<f64> {
    return  x.mapv(relu_dim1);
}

fn relu_dim1(x: f64) -> f64 {
    let res = if x > 0.0 {x} else {0.0};
    return res
}

pub fn softmax(x: &Array1<f64>) -> Array1<f64> {
    let y = x - x.to_vec().iter().fold(0.0/0.0, |m, v| v.max(m));
    let exp_y = exp(&y);
    let inv = 1. / &exp_y.sum();
    return exp_y * inv;
}
```

ニューラルネットワークを実装します。
ニューラルネットワークは上で定義した活性化関数と行列の掛け算を組み合わせたて、値を出力します。

実装するとこうなります.

```rust
pub fn forward(){
    let x = arr1(&[1.0, 0.5]);
    let mut weights = HashMap::new();
    let mut biases = HashMap::new();
    weights.insert("W1",
                   arr2(&[[0.1, 0.3, 0.5], [0.2, 0.4, 0.6]]) );
    biases.insert("b1",
                   arr1(&[0.1, 0.2, 0.3]));
    weights.insert("W2",
                   arr2(&[[0.1, 0.4], [0.2, 0.5], [0.3, 0.6]]) );
    biases.insert("b2",
                   arr1(&[0.1, 0.2]));
    weights.insert("W3",
                   arr2(&[[0.1, 0.3], [0.2, 0.4]]) );
    biases.insert("b3",
                   arr1(&[0.1, 0.2]));
    let a1 = x.dot(&weights["W1"]) + &biases["b1"];
    let z1 = relu(&a1);
    let a2 = z1.dot(&weights["W2"]) + &biases["b2"];
    let z2 = relu(&a2);
    let a3 = z2.dot(&weights["W3"]) + &biases["b3"];
    let z3 = softmax(&a3);
    println!("{:?}", z3);
}
```


## ニューラルネットワークの学習
4章では学習させる所まで実装させていますが、ここでは
ロス関数を定義し,それに基づき精度を評価することと
ミニバッチ単位でロスを計算する方法だけ実装します.

これで推論時の結果評価含め最低限できます.
学習のためには他にもままだ必要ですが,最低限使い方がわかったため,ここで終わりにします.

### ロス関数
機械学習ではロス関数を最小化するようにアルゴリズムをチューニングします.
分類だからといって、あってるあっていないだけでなく、どのぐらいずれていたかを評価できるような手法を使います.
そうした評価指標とは異なるもので最適化する最大の目的は勾配の計算のしやすさです.
機械学習を含め反復的な最適化の基本的な手法は勾配をもとにパラメータを変更していくので、
適切な勾配を設定することが重要であったりします.

よく使われる損失関数として2つ,二乗誤差とクロスエントロピーを紹介します.

- $(t_1, \ldots, t_n)$が推論結果
- $(y_1, \ldots, y_n)$が正解だとします.

二乗誤差はそのまま$\sum (x_i - t_i)^2$と定めたものです.
rustで実装するとこうなります.

```rust
pub fn _squared_error(t: &ArrayView1<f64>, y: &ArrayView1<f64>) -> f64{
     (t-y).dot(&(t-y))
}
```

### クロスエントロピーの説明
クロスエントロピーはあまり丁寧に説明されているものが少ないので,長いですが,僕の理解を説明しておきます.

### クロスエントロピーの定義

クロスエントロピーの定義を最初に書くと
$$- \sum y_i log t_i$$
となります.

ただし分類タスクの場合は$y_i$が1 or 0であることに注意すると,正解が$i$の場合は
クロスエントロピーは$- \log t_i$となります.

基本的に、確率分布同士の組に対し,そのずれを測るロスです.
なので推論結果の合計も1,正解の合計も1の時に使います.


これは情報量$-\log t$を$y$に関する確率で期待値を取ったものになります.

#### 情報量
情報量自体の定義も説明します.
情報量はイメージで考えると

- 珍しいものに対しては大きく,普通のものには小さく
- 二つ情報があった場合は情報量は足し合わされる

ものです。

これを数学的な条件と書くと

- 単調減少(確率が高い場合はよくあることなので情報量は小さい)
- 連続
- 独立な確率同士の積が発生する確率はもとの情報量の和になる.

という条件がありえます.

実はこの条件を満たすという制約だけから定数倍を除くと
$- \log p$しかありません.

そのため,情報量は$ - \log p$となります.

情報量的に見ると,クロスエントロピーは予測の平均的な情報量を判定する手法になります.

#### KL divergenceとの関係


分類タスクの場合はクロスエントロピーと$KL divergence$と一致します(log0が現れるので極限で考える必要があります).
それは、$q$を真の分布、$p$を予測の分布とすると,

```math
KL divergence(q,p) = \sum_{i=1}^n q(i) \log \frac{q(i)}{p(i)}
```
で,今の場合
$\sum_{i=1}^n q(i)\log q(i) = 0$となることから,わかります.

KL divergenceは確率分布同士の距離を測る自然な手段であることからも,クロスエントロピーを計算するのが自然なのが,わかっていただけるかと思います.

KL divergenceは二つの確率分布が一致した時に最小になるので,クロスエントロピーも予測の確率分布が真の分布に一致した場合が最小になります.

さらに,$t_i$が$y_i$に近ければ近いほどよい値になります.

__注意(ロスの使い分け)__

ロスの使い分けは絶対的なものはないのですが,分類タスクで二乗誤差を取らない理由の一つは
全ての成分で差が1以下で少し変動ても大きい損失にならないからです.

では、実際に二乗誤差とクロスエントロピーを実装しましょう.



```rust
pub fn squared_error(ts: &ArrayView2<f64>, ys: &ArrayView2<f64>) -> f64 {
    let mut diff = 0.;
    Zip::from(ts.genrows()).and(ys.genrows()).apply(|t, y| {
        diff += _squared_error(&t, &y);
    });
    return diff;
}

pub fn cross_entropy(ts: &ArrayView2<f64>, ys: &ArrayView2<f64>) -> f64 {
    let mut diff = 0.;
    Zip::from(ts.genrows()).and(ys.genrows()).apply(|t, y| {
        diff += _cross_entropy(&t, &y);
    });
    return diff;
}

pub fn _squared_error(t: &ArrayView1<f64>, y: &ArrayView1<f64>) -> f64{
     (t-y).dot(&(t-y))
}

pub fn _cross_entropy(t: &ArrayView1<f64>, y: &ArrayView1<f64>) -> f64 {
    let mut diff = 0.;
    Zip::from(t).and(y).apply(|ti, yi| {
        diff -= yi * ti.log2();
    });
    return diff;
}
```


これで推論し、その損失を評価することができるようになりました.
今回の記事はここまでとなりますが余力があれば実用的なdeep learningライブラリも作れたらと思います.