分からんこと多すぎ

何故か就職できたので技術磨く。一人前の多重線形代数使いを目指しつつ、機械学習、データマイニングをやる

ゼロから始めるDeepLearning_その1_ニューラルネットとは

対象とする人

ディープラーニングすごい! ←聞き飽きた

チュートリアルあるよ! ←ふわっとしすぎて具体的なところが分からん

こういう論文あるよ! ←読めるわけないだろ

そういう人向け。(たぶん学部四年程度向け

ニューラルネット初学者が、書ききるまで怪しいところ満載でも突っ走ります。

ニューラルネット(この記事)
→(AutoEncoder)
→(DenoisingAutoEncoder)
ホップフィールドネットワーク
ボルツマンマシン
Restrictedボルツマンマシン
→(Gaussian Binary - Restricted Boltzmann Machines)
→(DeepBeliefNetwork)
→(DeepNeuralNetworks)
→畳み込みニューラルネット(後日)

までやる。

太線以外は読み飛ばしてOK

本文中では怖い式は使わない。(Appendixに書く)

分からない人がいたら、説明を追加する。

正直算数が覚束ないので、間違ってたらごめんなさい。

ディープラーニングとは

①.ディープラーニングとは、”ニューラルネット(っぽいもの)”を”多層に”積み上げたものによる学習器の総称である。

②.特徴量を自動で作ってくれる。

③.画像認識とかの記録を驚くほど塗り替えた、新進気鋭の技術。

今回は、まずニューラルネットワークの話をする。

ニューラルネットをとりあえず使ってみる

ニューラルネットを使うと、任意の関数が近似できる。

クラス判別の場合は、ある関数が境界面をつくっていると考える。

境界面をつくっているある関数を推定する問題と見れば、任意の関数を近似する問題と同様に解ける。

関数近似

図的にはこんな感じ。

横軸が入力 x 、縦軸が出力 y .

図中の黒線は推定したい関数y=sin(x)

図中の赤点は推定された関数(初期のニューラルネットによる関数)

状況:赤点の横軸上の値 x と黒線上の値 y=sin(x) + noise が与えられている

目的:赤点の縦軸上の値(出力値 y' )が黒線に近づくように未知の関数(今回はsin関数)を求めたい

ニューラルネットの役割:未知の関数そのものを、近似する関数になること

更新を重ねると、以下のように推定した関数が変化していく。

ニューラルネットとは

ニューラルネットとは、多層パーセプトロンのこと。

そもそもパーセプトロンってなに?

  1962年にローゼンブラットが作ったアルゴリズム

  y \left( x \right) = f \left( w^{T} \phi \left( x \right) \right)

  で表される線形識別モデルのこと。

  ここで非線形活性化関数fは

   f \left( a \right) = \begin{cases} +1, \; a \geq 0 \\ -1, \; a < 0  \end{cases}

  というステップ関数である。

  つまり、

  入力データxを、

  関数Φによって特徴量に変換して、

  重みwで足しあわせ、

  ステップ関数で二値に分類するものである。

尤度関数がパラメータの凸関数にならないため、勾配法を何回も解いて解に近づけることになる。

当然過学習もする。(EarlyStoppingという、途中で学習を打ち切る方法で対処したりする)

良いところは、あらゆる関数を近似できることと、
学習し終わったモデル自体は軽量なので、新しいデータを高速に処理できること。

ニューラルネットを図的に書くと、こんな感じになる。

入力変数ベクトルxが与えられた時、線形変換wをいじることで、出力変数ベクトルyを目標ベクトルtに近づけるのが目的となる。

したがって、以下の誤差関数を最小化すれば良い。

一個ずつ説明していく。

ニューラルネットの仕組み

①.まず、d次元のデータベクトルがやってくる

あと、関数の切片とかをいじるために、バイアスパラメータx0=1をくっつける。

②.入力ベクトルが線形変換されて、各活性aに入る

つまり,

a_{j} = \sum_{i=0}^{d} w_{ji}^{(1)} x_{i}

③.各活性がしきい値を超えたら発火する(1の値を取る)

つまり,

z_{j} = h( a_{j} )

非線形関数h(・)はロジスティックシグモイド関数や、tanh関数のようなシグモイド関数を用いる。

例えば,標準シグモイド関数は,こんな関数.

\sigma \left( a \right) = \frac{1}{1 + \exp \left( - a \right)}

非線形変換がないと任意の関数を近似することができないので、何らかの関数が必要とされる。

特にこれらのシグモイド関数が用いられるのは、微分した時の値の計算が簡単なためである。

指数関数族は微分にて最強。

④.隠れ層zが線形変換されて、出力ユニット活性aに入る

ついでに隠れ層にバイアスz0=1をくっつけて、

a_{k} = \sum_{j=0}^{m} w_{kj}^{(2)} z_{j}

⑤.(クラス判別問題の場合)各活性がしきい値を超えたら発火する

y_{k} \left( x,w \right) = \sigma \left( a_{k} \right)

関数近似の場合は、ここでシグモイド関数を用いない)

したがって、①から⑤までの式をくっつけると、

図からも見て取れるように、ニューラルネットワークは、パーセプトロンが二段重なった構造になっている。

パーセプトロンとの違いは、ステップ関数ではなく、連続で非線形シグモイド関数を用いる点である。

連続関数を用いているため、パラメータに関して微分可能であり、高速な学習を可能としている。

また、シグモイド関数非線形領域を用いた場合のみ、ニューラルネットは万能の関数近似器になる。

シグモイド関数の真ん中ら辺は、線形領域なので、ただの線形変換しかしない。
 ただの線形変換によってニューラルネットを作った時、
 モデルの表現力が著しく損なわれることを、ローゼンブラットが証明した)

隠れ層の素子数を、入出力次元数よりも小さくすると、
主成分分析のように特徴量が圧縮される(AutoEncoder)。

逆に、隠れ層の素子数を、入出力次元数よりも大きくすると、
スパースコーディングのように特徴量がスパースになる(SparseAutoEncoder)。

どうやって学習するの?

誤差関数は、次の式で与えられる。

これをwを変化させることで最小化する。

関数の最小値は、微分して勾配が0になる点にあるので、

\nabla E(w) = 0

最急降下法を用いると、wの更新式は以下のようになる。

w_{(n+1)} = w_{(n)} - \eta \nabla E(w_{(n)})

これはすべてのデータを一度に使うバッチ学習というもので、
ニューラルネットの学習にはふさわしくない。
(計算オーダーが、後に紹介する最適化よりもはるかに大きい)

そこで今回は、得られたデータごとに逐次、勾配降下法を行うこととする。
オンライン学習

n番目のデータが来るごとに、誤差E_{n}を計算する。

その和が最終的な誤差関数となる。

 E ( w) = \sum_{n=1}^{N} E_{n} (w)

この誤差関数を最小化していくのが、オンライン勾配降下法である。

オンライン勾配降下法は、確率的勾配降下法(StochasticGradientDescent)と呼ばれている。

SGDの話は、昔ちょっと触れた。

確率的勾配降下法(SGD)の並列化について - 分からんこと多すぎ

w_{(n+1)} = w_{(n)} - \eta \nabla E_{n} (w_{(n)})

要するに、上式のようなパラメータの更新をすれば良い。

Back propagation

上図のような関係を考える。

y^{(l+1)}は、y^{(l)}を重みw^{(l)}で足しあわせて、関数f^(l+1)をかけたものである。

y_{k}^{(l+1)} = f^{(l+1)} ( \sum_{j} w_{kj}^{(l)} y_{j}^{(l)} )

この時、n番目のデータに関する誤差関数E_{n}^{(l)}の重みパラメータw^{(l)}による微分は、以下のようになる。

(l)層での出力ベクトル(シグモイド関数後の値)をy^{(l)}と書くこととする)


(導出は長いので最後に書いた)

f:id:rishida:20140218211520p:plain

ちなみにNは層の総数。

したがってパラメータの更新式は、以下のようになる。

w_{(n+1)}^{(l)} = w_{(n)}^{(l)} - \eta \nabla E_{n} (w_{(n)}^{(l)}) = w_{(n)}^{(l)} - \eta \delta^{(l+1)} {y^{(l)}}^{T}

実装

Rで実装する。

#使い方例
if(FALSE){
  w <- NULL
  for(i in 1:500){
    x <- runif(min = -3, max = 3, n = 100)
    y <- sin(x)
    if( (i < 10) || (490 < i)){
      w <- neural.net(x = matrix(x,nrow = 1), y = matrix(y, nrow=1) , w = w, graph = 'y = sin(x)', noise = 0.1)
    }
    w <- neural.net(x = matrix(x,nrow = 1), y = matrix(y, nrow=1) , w = w, graph = NULL, noise = 0.1)
  }
}

#ニューラルネット関数
#必須項目:
#入力行列x(入力の次元数=行数),出力行列y,重みベクトルのリストw
#オプション:
#中間層素子数nnum,更新重みeta,クラス判別問題か否かのフラグCLASS,
#グラフを書くときのタイトルgraph,ノイズの分散noise
neural.net <- function(x,y,w = NULL, nnum = 5, eta = 0.1, CLASS = FALSE, graph = NULL, noise = 0){
  sigmoid <- function(x){
    return(( 1 / ( 1 + exp( -x))))
  }

  #create weight
  if(is.null(w)){
    w <- list(NULL)
    w[[1]] <- matrix(rnorm(n=nnum * (dim(x)[1]+1),mean=0,sd=(1/sqrt(dim(x)[1]))),nrow=nnum,ncol=(dim(x)[1]+1) )
    w[[2]] <- matrix(rnorm(n=dim(y)[1] * (nnum+1),mean=0,sd=(1/sqrt(dim(x)[1]))),nrow=dim(y)[1],ncol=(nnum+1))
  }

  #mix data
  idx <- sample(c(1:dim(x)[2]),size=dim(x)[2],replace=FALSE)
  x <- matrix(x[,idx],ncol=length(idx))
  y <- matrix(y[,idx],ncol=length(idx))

  #add noise
  y <- y + rnorm(mean=0,sd=noise,n=prod(dim(y)))

  for(i in 1:dim(x)[2]){
    #compute neural.net
    out <- list(NULL)
    out[[1]] <- rbind(1,x[,i])
    out[[2]] <- w[[1]] %*% out[[1]]
    out[[2]] <- rbind(1,sigmoid(out[[2]]))
    out[[3]] <- w[[2]] %*% out[[2]]
    if(CLASS == TRUE){
      out[[3]] <- sigmoid(out[[3]])
    }

    #Back Propagate
    #BP1
    delta <- -1 * (y[,i] - out[[3]]) * (sigmoid(out[[3]])) #* (1 - sigmoid(out[[3]])))
    wdet <- t(out[[2]] %*% t(delta))
    hoge <- y[,i] + (y[,i] - out[[3]])
    #eta <- (matrix(c(w[[2]] * wdet) * c(out[[2]]),nrow=dim(y)[1]) - (wdet * t(out[[2]] %*% out[[3]]))) / (wdet^2)
    w[[2]] <- w[[2]] - eta * wdet
    #BP2
    delta <- ((out[[2]][-1] * (1 - out[[2]][-1])) * as.vector(w[[2]][,-1] * delta))
    w[[1]] <- w[[1]] - eta * t(out[[1]] %*% t(delta))

    #plot
    if(!is.null(graph)){
      if( ((i - 1) %% 10) == 0 ){
        sim.neural(x,y,w,CLASS,graph)
      }
    }
  }

  print(paste('err',out[[3]] - y[,dim(x)[2]]))
  return(w)
}

#推定結果を図示する関数
sim.neural <- function(x,y,w,CLASS = FALSE,graph){
  sigmoid <- function(x){
    return( 1 / ( 1 + exp( -x)))
  }

  #sort
  idx <- order(x)
  x <- matrix(x[,idx],ncol=length(idx))
  y <- matrix(y[,idx],ncol=length(idx))

  #simurate
  out <- list(NULL)
  out[[1]] <- w[[1]] %*% rbind(1,x)
  out[[1]] <- sigmoid(out[[1]])
  out[[2]] <- w[[2]] %*% rbind(1,out[[1]])
  if(CLASS == TRUE){
    out[[2]] <- sigmoid(out[[2]])
  }

  #plot
  plot(x,y,xlim=c(min(x)-0.1,max(x)+0.1),ylim=c(min(y,out[[2]])-0.1,max(y,out[[2]])+0.1),xlab="",ylab="",type='l')
  par(new=TRUE)
  plot(x,out[[2]],col='red',xlim=c(min(x)-0.1,max(x)+0.1),ylim=c(min(y,out[[2]])-0.1,max(y,out[[2]])+0.1),xlab="x",ylab="y",main=graph)
  return(0)
}

飛ばしたBackPropagationの計算。


次回はボルツマンマシンの話。