ディープラーニング (Deep Learning:CNN )

はじめに

Deep Learning(ディープラーニング)という言葉は、もはや殆どの人が耳にしたことがあると思います。ディープラーニングは何かを調べていくと、今いち何か掴めないという人が多いのではないでしょうか?今回はディープラーニングの1種であるCNN( Convolution Neural Networks)を中心に画像処理を例として ディープラーニングとはなにか、どう学習させていくかを説明いたします。

前段は特徴量作成ステージ。学習により畳み込みフィルターが作られる。後段はClassification。普通のNNとなる。Fully connectedであり、最後にはSOFTMAX関数を用いる。

Deep learningの正体

ディープラーニングはDeep Neural networks(DNN)とも呼ばれており、その正体はニューラルネットワーク(NN)そのものです。ディープラーニングという言葉が流行り始めた際、中身がニューラルネットワークという事を知って、本当に驚きました。というのも、NNはとても古く伝統ある学習アルゴリズムでしたが、ランダムフォレストやSVMの台頭により、誰も利用しない、古い学習アルゴリズムとして扱われていたからです。厳密にはNNとDNNの違いは、特徴量の作成を学習化にとりいれた事でした。今まで画像処理の分野では技術者は頑張って特徴量を作り、その特徴量をNNなどの学習機で学習させていましたが、DNNでは、画像を突っ込めば前段で特徴量を勝手に作ってくれます。(厳密には特徴量を作成するフィルターを作成する。そのフィルターをかませることで特徴量マップを作成する)これは、とても斬新でした。

Forward process

さて、ディープラーニングの一種であるCNNを例にどのようなに学習しているかを説明します。CNNがわかると、DNNの構造がわかり、例えば最近の技術であるDCGANの理解も進みます。

CNNはFowardプロセスとして前段と後段の2つに分けることができます。まずは前段の畳込みとプーリングについて説明します。

前段( 畳込みとプーリング )

1.畳み込み処理

3x3の場合は、上記のように畳み込み処理を行う。畳み込み処理をして出てきた値に対して、 ReLU( max( 0, x ) ) を実施している例が多い。フィルターのサイズについては、経験的に決める。

CNNでやっていることは実は単純です。CNNでは、入力画像にフィルターによる畳込み演算をし、それを特徴量マップとして次の層に渡します。この特徴マップづくりが特徴量を自動でつくる、と言われている部分でもあります。CNNは多層であるということを言いましたが、次の層でも、同様に畳込み演算を実施して、更に特徴量マップを作成します。そして次の層へ…ということを繰り返し実施するのです。何層にするかは自分で決めればいいのですが、大体2層くらいで落ち着かせる事が多くあります。

CNNの肝は特徴量マップ作りですが、これは言い換えればフィルター(畳み込みの重み)を学習により自動で求めることにあります。繰り返しフィルターの重みをを更新していく処理により精度を高めていきます。フィルターの更新にはバックプロパゲーション(誤差伝播)により更新していくのですが、こちらは後述します。

なお、畳込み処理ではReLU処理をする場合がほとんどです。ReLUとはmax( 0, x )というただのMAX処理で、畳こみ値がマイナスであれば0にする操作となります。

2.プーリング処理

Min poolingの例。画素の一番低いピクセル(黒に近いもの)を選択する。Max Poolingであれば白に一番近いものが選ばれれる。Poolingにより特徴マップの次元は削減される。Max Pooling, Min Poolingの他にAverage Pooling等がある。

CNNではPoolingという作業(間引き)を行っています。Poolingは特徴量をマップを作った直後(あるいはReLU処理後)に実施します

最も使われるPooling操作は畳み込み演算後の特徴マップから、ウィンドウの中で最大値だけをとるMaxPoolingです。Poolingを行う最大の目的は次元の削減にあります。またPooling操作は必須ではありません。事実、DCGANなど最近のDNNの一種では行っていません。ただ、伝統的な?CNNでは利用されています。また、至るところでも解説されているので覚えておくといいでしょう。

前段における注意事項

Note 1:入力画像がカラー画像の場合

入力画像がカラー画像のとき、つまりR,G,Bチャンネルが3つあるときですが、フィルターが3次元となります。例えばウィンドウサイズが3x3であれば、RGBを含め3x3x3のフィルターを用意します。このフィルターを適用して、特徴量マップの1つの要素ができます。

Note 2:何枚の特徴量マップを作ればいいのか

各層において、何枚の特徴量マップを作ればいいのかは、自由に決めます。64,128など、霧のいい数字を用いることが多いです。

NOTE 3: 何層作ればいいのか

2層の例が圧倒的に多いですが、こちらも経験的に決めているようです。

後段処理(Classification)

後段処理では、伝統的なニューラルネットワークを行います。例えば、最後に14x14の特徴量MAPが3つできたとします(14x14x3)。これを1次元にすると588次元です(Flattering処理)。これを用いて、普段のNNで使われるテクニックを使い、指定のクラスの次元まで落とし込むのです。SOFTMAXがよく使われています。この方法については、従来のNNそのものであるので割愛します。

Backward process

誤差伝搬はどのようにやるのでしょうか。後段処理についてはNNで今までと変わりません。前段の誤差伝搬処理について解説していきます。

誤差伝搬を行いますが、今までとの逆の操作をしていくことになります。前段処理は次のような流れでした。

画像 –> 特徴量マップ –> ReLU –> MaxPooling

本当は数回特徴量作成とPoolingをやりますが、ここでは簡単に一層で作成したと仮定し説明します。さて、逆操作では以下になります。

d Max Pooling –> d ReLU–>d Convolution(特徴マップ逆操作 )

d Max Pooling

Forward processでのMax Poolingは最大値だけを抽出し、次元を縮小する操作でした。Backward Processでは、次元縮小に利用したフィルターを用いて、後段から来た誤差をもとの大きさにもどしてやります。

d11, d12, d21, d22は後段のNNからの逆伝搬で算出された誤差。

Backword processでは、このMAX Filterを使うので、前段で作ったMAX Filterは保持しておくことが必要です。

d ReLU

ReLUはMAX関数。0より小さければ0,それ以外は1にする関数となる

この操作を図示すると次のとおりです。単純に0より上であれば1にする操作となります。

例えば4x4の特徴量Mapにd ReLUを実施すると次のようになる。

d Convolution

特徴量MAPの作成に利用したフィルターを使い、復元していきます。

はみ出した領域は、0でパディングしておく(つまり、はみ出した部分の重みは考慮しない)

180度転置したものが伝搬で使われるのかについてはForward And Backpropagation in Convolutional Neural Network[1]に詳しく述べられていますので、参照ください。参考サイトは英語ですが、わかりやすい数式で述べられているので、興味ある人は読めるはずです。

参考文献

[1] medium.com/@2017csm1006/forward-and-backpropagation-in-convolutional-neural-network-4dfa96d7b37e

https://towardsdatascience.com/a-comprehensive-guide-to-convolutional-neural-networks-the-eli5-way-3bd2b1164a53

おすすめ記事

Matrix Factorization(MF)

はじめに

Matrix Factorization(MF)について解説していきます。

現在世にある推薦システムは殆どMFをベースに作られているといっても過言ではありません。Simon Funkという人がNetflix Prizeに応募するために開発し、とても話題になりました。SVDを用いて推薦システムを構築したのですが、あるトリックを用いたSVDであったため、中のエンジンはFunk-SVDなんて言われたりしています。他にもNMFを利用したりSVD++など改良が進んでいますが、メインのやり方は本稿の説明でわかるかと思います。

結局、やることは従来のSVDとおなじになる。ただし、求め方で工夫を行う

何故SVDを改良する必要があったのかと言えば、SVDは疎でない(スッカスカでない)行列の解析には素晴らしく良い結果をもたらしていました。しかしながら、商品や映画の推薦となると、全てのユーザーが全てのアイテムを推薦している訳はありません。実際に映画のレーティングの行列を作ると95%のがNaN(未評価)であったということからも、推薦のもとになる行列はとても疎な行列です。この状況を改善しようとした試みがFunkSVDだったというわけです。

概要

Matrix Factorizationの流れを説明します。次のようなレーティングマトリックスを2つ(Item & User )に分けるものです。他のサイトと差別化するために、本ブログでは 一番わかり易い 表ベースの解説をします。今、以下のようなレーティングマトリックスがあったとします。

Movie 1Movie 2Movie 3Movie 4
User 15NaN83
User 277NaN4
User 3NaN677
User 410NaNNaNNaN

上記の行列を2つの行列q, pに分ける事をMFでは行います。

\hat{r}_{ui} = q_i^Tp_u

\hat{r}は推定されたレーティングという意味し、実際のレーティングはrとします。

p_uはユーザー、q_iはアイテムを表す行列です。今回の目的は、ユーザー特徴量行列とアイテム特徴量行列をつくり、それらを掛けると元の行列にする事です。

これは実際のレーティングrと推定レーティング\hat{r}の2乗誤差を最小にするというな行列を分解する,という意味でもあります。式で表せば次のようになります。

minimum( p,q ) \sum_{u,i\in K}(r_{ui}-q_i^Tp_u)^2

p, qを作る

ユーザ特徴量行列p、アイテム特徴量行列qは今は未知なので適当なランダム値で埋めてやります。

p=
Factor 1Factor 2Factor 3
User 11.40.4-1.2
User 21.31.5-0.3
User 30.80.10.6
User 4-0.70.51.1
q=
Movie 1Movie 2Movie 3Movie 4
Factor 10.60.51.3-0.5
Factor 2-1.1-0.40.71.3
Factor 310.30.5-0.1

さて、説明のために、User1のMovie3に注目します。

実際の評価rはレーティングマトリックスをみてみると”8″でした。

ユーザ特徴量行列p(青)ではUser 1は[1.4, 0.4, -1.2]

アイテム特徴量行列q(緑)のMovie3は[1.3, 0.7, 0.5]

です。推定された評価\hat{r}は[1.4 * 1.3 + 0.4 * 0.7 + 0.5 * -1.2 ] = 1.6となります。なお、算出しなくてもいいのですが2乗誤差はは(8 - 1.6)^2=40.96となりました。

更新処理

誤差を元にp, qの値を更新します。p,qの行列から参照した要素は次のとおりでした。

From p1.40.4-1.2
From q1.30.70.5

新しい特徴値は次の計算式を用いて求めます。

p_i = p_i + 2 \alpha ( r - \hat{r} ) q_i

q_i = q_i + 2 \alpha ( r - \hat{r} ) p_i

ここでαは学習率と呼ばれ、0.1くらいにするケースが多いです。これに習ってpの最初の要素を計算してみると、次のようになります。(p_0は1.4, q_0は1.3を指します)。何故この更新式になるかは[2]を参照ください。

1.4 + 2 \times 0.1 \times ( 8 - 1.6 ) \times 1.3 = 3.064

出た値を用いてp_0の値を更新します。

p3.0640.4-1.2

同様にp, qの残りに対しても更新処理を行います。なお、qの更新式ではpの値も使いますが、既に更新された値を使います。

上記例ではp_0の値は3.064になりましたので、q_0を求める際には、新しい数値を使います。あとは繰り返し処理を繰り返し行い、 すべてのp, q行列を更新していきます。以上となります。

これを学習データを用いて適当に取ってきた値で特徴量行列を更新していきます。

ちなみに本当の評価rがNaN値であれば重みの更新処理をしなくても構いません。(というかしません。)。これがNaN値を考慮したSVDであると言える理由です。

重み更新は、ニューラルネットワークなどではおなじみなのですがStochastic Gradient Descent (SGD):確率的勾配降下法と呼ばれています。ランダムに(確率的に)適当な値を引っ張ってきて、重みを更新している(勾配をもとに、2乗誤差を低くするように降下させている)からです。差分がゼロ(微分がゼロ)になれば、最適解に到達したと言えます。

なぜ特異値分解:SVDと呼ばれるのか

今までの操作のどこが特異値分解なのでしょうか。

特異値分解とは、ある行列をエラーを最小とするように2つに分割する操作でした。今回分割する際に行ったことは2乗誤差が最小になるように分割していく操作です。SVDはもともと最小二乗近似を得る方法とも言えます。 実際に開発者のFunkも自身のブログで次のように述べています。

Singular value decomposition is just a mathematical trick for finding those two smaller matrices which minimize the resulting approximation error–specifically the mean squared error (rather convenient!).

さて、大まかな説明はこれで終わりなのですが、更にいろいろな工夫をすることで精度が高くなります。SVDのかわりにNMFを使うことだったり、今回はマトリックスを分割する際のエラーの計算方法として単純にpqに分けるようにしましたが、実際に精度を出すためには正則化したり、バイアス項を突っ込む必要があります。

さらなる改良

正則化

結局重みは次の式を最小化するということでした。

minimum( p,q ) \sum_{u,i\in K}(r_{ui}-q_i^Tp_u)^2

ここで上記式に正則化をします。(過学習を防ぎ、一般的に使えるようにするためにペナルティ項をおく)

minimum( p,q ) \sum_{u,i\in K}(r_{ui}-q_i^Tp_u)^2 + \lambda (\|q_i\|^2 +\|p_u\|^2)

最初の項目に\lambda (\|q_i\|^2 +\|p_u\|^2)をおいただけです。この意味ですが、1つの映画を酷評したユーザーがいて、その映画以外は評価していなかったとします。すると、このユーザーから他の映画へのすべての評価が非常に低くなります。この正則化をおけば、 例えば q_iに大きな値を与えることを抑制することになります。

バイアス項

推測されるレーティングは次の式でもとまるものでした。

\hat{r} = q_i^Tp_u

ここに3つのバイアスを加えます。

  • \mu 全アイテムの平均レーティング
  • b_i アイテムiの平均レーティングから\mu を引いたもの
  • b_u はユーザーuの平均レーティングから\mu を引いたもの

最終的に、レーティングの式は

\hat{r} = q_i^Tp_u+\mu + b_i + b_u

となります。そして、こちらにも正則化項をいれます。

minimum( p,q ) \sum_{u,i\in K}(r_{ui}-q_i^Tp_u)^2 + \lambda (\|q_i\|^2 +\|p_u\|^2+b_i^2+b_u^2)

これが最終的なSVDの式となります。長くなりましたが、これでMFの説明は終わりとなります。

最後に

その他、NMFやSVD++など、エンジンを変えて更新式を実行したい場合は、次のURLに詳細に書かれています。[1]

もしNMFによる実装方法などを記事化してほしい方がいれば、コメントをいただけますと幸いです。

索引

関連記事

Word2Vec

はじめに

Word2Vecは近年で最も注目された技術の一つです。主にNLP( Natural Language Processing:自然言語処理 )の分野の技術なのですが、広告や商品ページなどの推薦エンジン技術として利用されています。実際に楽天などでは本技術を使っているようです[1]。
コンセプトはとてもシンプルです。今回はWord2Vecがどういう技術かをpythonをもとにnumpyを使ったソースコードの実装例をもとに説明します。(gensimなどのライブラリは使いません!)。実際内部でどのような処理が行われているか参考にしてください。

どのようなモデルか

このような構成をしています。すべて本ブログで解説しますので、読み終わる頃にはよくわかります。

Word2Vecは何か

単語をベクトル化する、これをWord2Vecが行います。単語のベクトル化ができることで、四則演算が使えてしまうのです。。

例えば、仮に単語の「神様」:[5, 1, 3]というベクトルで表した場合、「神様」と似ているベクトルを探すこと(例えばコサイン積をとって1に近いものを探すこと)ができます。他にも「神様」:[5,1,3]へ新しく、例えば「海」:[1,2,3]というベクトルを足し込んだ結果を探すと、「リヴァイアサン」:[6,3,6]がヒットするかもしれません。外積をとって直行する言葉を探すなんてこともできそうですね。

その昔、類似度を探す技術と言えば共起行列を用いて探る方法が一般的でした。文章から単語を抽出してともに出現する単語から共起行列を作成してSVDなどで解析し、特徴量空間にマッピングして類似度を算出するのです。ただ、この方法では行列のサイズが大きくなると計算がおいつかなくなる問題がありました。2013年にGoogleの研究者がWord2Vecを提案し、業界がざわつきました。

単語をベクトル化することは「word embedding」と呼ばれます。

Word2Vecの実装に向けて

Word2Vecのアルゴリズムは大きく2つあり、違いは次のとおりです。

  • Continuous Bag-Of-Words (CBOW)
    • Skip-gramよりも数倍学習が早く、頻出単語に関しては若干の高い精度を誇っている。
  • Continuous Skip-gram (SG)
    • 小さい学習データから珍しい単語やフレーズを上手に表現できる。

2つともに近くの単語( context words )からの確率分布を用いた手法であります。近傍の単語から新しい単語の類推をしていきます。今回はSG( Skip-gram )を説明します。

学習方法

学習データは文章を使います。今回は次のフレーズを元にWord2Vecを説明します。

「Although the world is full of suffering, it is full of the overcoming of it.」
(世界は苦しいことでいっぱいだけれども、それに打ち勝つことでもあふれている。) : マザーテレサ

上記英語フレーズを学習データとします。実践では長い文章を用意します。

さて、Skip-gram では学習にニューラルネットワークを使います。教師データとなるトレーニングデータは単語の近傍の単語(context wordsとよばれる)です。具体的に行きましょう。例えば、先のフレーズの「world」であれば、[Although, the, is, full]がcontext wordsとして学習セットになります。

近傍の単語数はWindow Sizeともよばれ、上記例でのWindow Sizeは5となる

フレーズ内の、すべての単語に対し同様に学習セットを用意します。たとえば、Althoughの近傍の単語context wordsは[the, world]となります。後述しますがovercomingの近傍ワードは[ of,the, of, it ]ですが、重複は除かれるルールになるので [ of,the, it ] となります。

Vocabulary & One-hot Encoding

Vocabulary

具体的な学習データの作り方ですがVocabulary(ボキャブラリ)を構築し、one-hotエンコーディングを実施する必要があります。

「Vocabulary」とは重複のない単語群のデータです。まずは文章を空白区切りします。

[‘Although’, ‘the’, ‘world’, ‘is’, ‘full’, ‘of’, ‘suffering,’, ‘it’, ‘is’, ‘full’, ‘of’, ‘the’, ‘overcoming’, ‘of’, ‘it’]

そして、重複を除きます。

[‘Although’, ‘the’, ‘world’, ‘is’, ‘full’, ‘of’, ‘suffering,’, ‘it’, ‘overcoming’]

‘it’, ‘is’ , ‘full’ , ‘of’ , ‘the’の重複が取り除かれ1つになり、9次元のデータがとれました。これで「Vocabulary」 が完成しました。

One-hot encoding( One hot vector )

続いてone-hotエンコーディングを行います。図を見るとわかりますが、単純に単語があるところを1にし、それ以外を0にするだけのベクトルを作るという意味です。

One-hot vectorとも呼ばれています。さて、以上で下準備は終わりました。

最初にworldのcontext wordsは「Although, the, is, full」と説明しました。このデータはone-hotベクトルで表現すると [list([0,0,1,0,0,0,0,0,0]), list( [ 1,0,0,0,0,0,0,0,0 ], [ 0,1,0,0,0,0,0,0,0 ], [ 0,0,0,1,0,0,0,0,0 ], [ 0,0,0,0,1,0,0,0,0 ] ) ]というペアになります。

Training

学習の流れですが、学習の全体像は次のとおりとなります。

今、例としてHidden層の次元数は10としている。最初の重み行列W1のサイズはしたがって9×10となる。W2も10×9となる

Word2Vecでは、ニューラルネットによる学習において、上記のうち、2つの重みマトリックス(W1, W2)を更新していきます。(その方法は次の章で記載します。)

学習終了後ですが、使われるものは重み行列W1です。実はこれが特徴量行列となっています。例えば、”although”[ 1,0,0,0,0,0,0,0,0]があったとして、これと学習後の重み行列W1をかけ合わせて出てきたものが”although”の特徴量ベクトルです。実際には掛け合わせると言っても one-hotベクトル内で 1が立っている行をW1から抜き出すという操作になります。

さて、では具体的にどのようにWord2Vecでは学習されているか解説します。

Trainig – Forward pass

まず、一番最初のプロセスとして行列W1,W2はランダムな0-1の範囲で初期化をしておきます。そしてニューラルネットのforward関数としては次の関数を使うことになります。

def forward_pass(onehot_input, w1, w2):
    h = np.dot(w1.T, onehot_input)
    u = np.dot(w2.T, h)
    y = self.softmax(u)
    return y, h, u

onehot_inputは先の例で言うところの [0,0,1,0,0,0,0,0,0] (入力単語 worldのone hot表現)です。重み行列W1とonehot_inputでドット積を実行してHiddenベクトルh(隠れ層)が求められます。続いて、重み行列W2とHiddenベクトルhでドット積を実行し、出力ベクトルuを算出しsoftmax関数をかけて予測値ベクトル yを算出します。

結局のところforward_pass関数におけるリターンバリューとしては[予測値ベクトル y、Hiddenベクトルh、出力層ベクトルu]の3つを得ます。 softmax関数により必ず値が0〜1の範囲になります。

def softmax(x):
    e_x = np.exp(x - np.max(x))
    return e_x / e_x.sum(axis=0)

Training- Error

forward_pass関数の結果から予測ベクトルyが返ってきました。このベクトルの値はsoftmax関数を噛ましているので必ず0-1の範囲に収まっています。さて今、次の予測ベクトルyが返ってきたとします。

y = [0.888, 0.801, 0.766, 0.552, 0.659, 0.371, 0.517, 0.087, 0.212]

ここで、出力すべきであった 「although, the, is, full」のone-hotベクトル[ 1,0,0,0,0,0,0,0,0 ], [ 0,1,0,0,0,0,0,0,0 ], [ 0,0,0,1,0,0,0,0,0 ], [ 0,0,0,0,1,0,0,0,0 ] それぞれと 予測ベクトルyの差分を取ります。期待された結果からどれだけ間違えたかを差分として求めるのです。具体的に示します。

予測の確率であるyとone-hot表現された単語との差分を、オレンジ色の枠線内が表している。

こうして出来た、それぞれの単語で求められたエラーベクトルをすべて足し込み1つのエラーベクトルをeを作ります。このエラーベクトルeはBack Propagation処理で利用し、重み行列を更新していきます。

Training- BackPropagation

ニューラルネットワークでおなじみのバックプロパゲーション(誤差逆伝播)を行います。目的はエラーベクトルeを利用して重み行列(W1, W2)を更新することにあります。

具体的な流れは次のとおりです。

  1. Hidden ベクトルhとErrorベクトルeで外積を取り、重み行列w2に対する差分行列dW2を作成します。
  2. 重み行列W2とErrorベクトルeで内積を取り、ベクトルmを作ります。入力ベクトルonehot_inputとベクトルmで外積を取り、重み行列w1に対する差分行列dW1を作成します。

それぞれの重み行列W1, W2から、 それぞれ差分行列dW1 *学習率 ,dW2 *学習率を引き、新しい重み行列を得ます。学習率は1%くらいに設定しておきます。

learning_ratio = 0.01
def backprop( w1, w2, e, h, x ):
    dW2 = np.outer(h, e)
    dW1 = np.outer(x, np.dot(w2, e.T))
    # 重みの更新
    w1 = w1 - (learning_ratio * dW1)
    w2 = w2 - (learning_ratio * dW2)

以上で終わりです。これをぶん回していくことで学習していきます。学習後、W1の重み行列が特徴量行列となっています。

Trainig-Loss

損失関数を定義し、学習進み具合をチェックできます。損失関数は実際に学習作業と関係ありませんが、ぶん回す回数(エポック数)を決める上で、発散したりしていないか、きちんと収束しているかを確認する上で重要となります。

損失関数は2つの部分で構成されています。1つめは出力層のすべての要素の合計、の負の値です(softmaxの前)。 2 つめは、近傍の数であるcontext_wordsの数と、出力層のすべての要素(指数関数の後)の合計のログの乗算です。

loss += -np.sum(u) + len(context) * np.log(np.sum(np.exp(u)))

類似語の探索

特徴量ベクトルは、最初に説明したとおりone-hotベクトルと学習後の重み行列W1とのドット積によりもとまります。特徴量ベクトルを取得したら、後は、四則演算となります。例えば類似ワードは内積値が1日買い物となります。ベクトル同士の四則演算をして、得た結果かから色々と考えてみると面白いかと思います。

参考文献

https://towardsdatascience.com/an-implementation-guide-to-word2vec-using-numpy-and-google-sheets-13445eebd281

おすすめ記事

特異値分解 (SVD: Single value decomposition)

特異値分解とは

特異値分解 (SVD)はPCAとほとんど同じと思ってください。SVDはPCAが適用できないような行列に対して、データ削減や特徴抽出が可能な方法となります。

PCAを行う行列は正方行列であり、ランクがフルである、つまり逆行列を持つ正則行列である必要がありました。SVDでは疑似逆行列を利用するので正方行列ではない行列やランク落ちしている行列からでも、計算処理が可能となります。

MAT_A = \begin{pmatrix} a & b \\ c & d \\ \end{pmatrix}

正則行列は、逆行列を持つ正方行列を指します。正方行列とは上記のようにN x N行列のことで、横と縦の数が一致した行列を指します。

データは常に正則行列であるとは限りません。PCAは正則行列である必要がある制約条件がありました。SVDは疑似逆行列を利用することで、PCAに似た計算処理を可能にするのです。

SVDはどんな時に用いられるのか

SVDが必要な背景としては、例えば推薦システムなどで用いられます。推薦システムは例えばユーザに映画を紹介したい、あるいはユーザに商品を紹介したい等のシステムです。この際に、ユーザー数と映画数が同じでないと計算上不都合が生じてしまうのです。実際には要素数が同じなどありえません。そこで、疑似行列を使った処理(SVD)を用いる、ということです。

どのようにしてSVDを用いて推薦システムを作るのか解説していきます。

SVDによる推薦システム

今ユーザー3名が5つの映画にそれぞれ次のような評価をつけたとします。

review_user1 = [1,4,5,4,5]
review_user2 = [1,3,1,4,4]
review_user3 = [1,2,5,4,1]
mat = numpy.array( [review_user1, review_user2, review_user3 ] )

上記のデータをSVDします。numpyのライブラリを使うことで一発で求まることができます。

def learn( mat ):
    u,s,v=numpy.linalg.svd( mat )
    return [u, s, v]

U, S, Vが求まりました。この関係ですが、元の行列matはU S V^tで求めるように分解をしています。Sは 対角行列 となっており、対角成分に特異値が並んでいます。対角成分は左上が最大で右下が最小となるようになっています。なんだか面白い行列が出てきました。

さて、このUSVを求めたら類似度はどのように求まるかというと、つぎのようにすれば求まります。eには現在自分の映画5個における評価を入れておきます。e= numpy.array([1,5,4,4,5])

def getSimilar( usv, e, dim =2):
    u, s, v = usv
    vv = numpy.array( numpy.matrix(v).T )    
    U = numpy.matrix( u[ :, 0:dim ] )
    V = numpy.matrix( vv[ :, 0:dim ] )
    S = numpy.matrix( numpy.diag( s[:dim] ) )
    E = numpy.matrix( e )
    pos = E * U * S.I
    result = []
    for i, user in enumerate(V):
        val = float(user *pos.T)
        val /= (numpy.linalg.norm( user ) * numpy.linalg.norm( pos ))
        result.append( [val+1, i ])
    return sorted( result, reverse=True )

getSimilar関数はコサイン積をとっています。Dimは次元数という意味で、PCAで言うところの寄与数の高いところ順に内積を取るイメージを取ってください。これにより自分と似ているユーザーが返ってきます。

最後に

実は今紹介したSVDによる推薦システムは、データにSVDを適用し、コサイン積を求めたものでした。推薦システムとしては1つ大きな問題が残っております。それは評価していないアイテム、例えば見ていない映画の評価を0とするとした評価がSVD適用時に考慮されてしまっているということです。見ていないものが、つまらなかった(評価0)は少し辺です。そして、スッカスカの巨大行列をSVDするっていうのは何か違和感があります。そこで評価していないものはスキップするというようなアイディアを取り入れた方法を使う必要があります。これが matrix factorization(MF)という手法になります。Netflixの推薦エンジンのコンペ
で圧倒的なスコアを叩き出したエンジンがまさにMFを使っていました。MFについて、次回ご説明します。

参考文献

https://ja.wikipedia.org/wiki/%E7%89%B9%E7%95%B0%E5%80%A4%E5%88%86%E8%A7%A3

おすすめ記事

ガボールフィルター(Gabor Filter)

ガボールフィルターは、以下の式で表される。

gaussian(σ, z) * exp ( i * z )

ここで、gaussian項はガウス分布(正規分布)の関数、exp項のiは虚数を意味し、zは(x, y):位置情報 を表している。

上記のガウス分布は、2次元の場合

gaussian (σ, x, y ) = ( 1 / 4 * π * σ ) * exp ( -1 * (x2 + y2) / (4 * σ2) )

として表される。

上記ガウシアンの式をgnuplotで描画した例は以下の通り

Gabolフィルターのカーネル。gauss(x) = 1 / (2 * sqrt(pi)) * exp(- x * x / 4 )

σ = 1(標準偏差) として、gauss(x) = 1 / (2 * sqrt(pi)) * exp(- x * x / 4 ) をプロット(シグマを1としたため、式が簡略化)

ガボール gaussian(σ, z) * exp ( i * z )の分解

式には、exp(i * z)の虚数が入っている。このパートを、実部と虚部に分ける。(exp(ix) は cosx + i sinxと分解できる)

画像処理でガボールフィルターを用いて処理を施す場合は、この二つのパートに分離して作業をすることが多い。

ガボールフィルターには、ガウシアン項とExp項があり、Exp項には虚数があり、Cos項 と Sin項に分離できる。(exp(iθ) = cosθ + isinθ)(下記で説明する回転のsin, cosとは違うので混同しないように。)

それぞれの項目と、ガウス関数を掛け合わせたのは以下の通り。

ところで、本やWebでは、ガボールカーネルの関数は、上記で示した gaussian(σ, z) * exp ( i * z ) ではなく、

gaussian(σ, z) * (exp ( i * z ) – exp ( σ2 ) )

と表される。このexp ( σ2 )を引いている理由は、

gaussian(σ, z) * exp ( i * z ) を積分した結果が0にならない(直流成分が0とならない)ためである。(これはCos項:実数部が原因している。)

exp ( σ2 )を引いた画像は以下の通り。

Cos項とガウシアン関数の畳み込みカーネル:実部実数部カーネル(cos * gaussian) 3 – dimention

Sin項とガウシアン関数を畳み込みカーネル:虚部虚数部カーネル(sin * gaussian) 3 – dimention
splot cos(x) * gauss(x,y)
splot sin(x) * gauss(x,y)
gauss(x,y) = (1 / (4 * pi)) * exp(-1 * (x ** 2 + y ** 2) / (4))
splot cos(x) * gauss(x,y)
splot sin(x) * gauss(x,y)

として描画。(ガウス関数はσを1とおいた。そのため式が簡略化されている。)注意:実際に、直流をゼロにするためには、Cos項は、

      gauss(x,y) * ( cos(x) –  exp(-1) )

とすべき。

尚、カーネルにスケーリング(a倍)を施したい場合は、

 1 / sqrt(a) * gaussian( σ, z / a ) * (exp ( i * z / a ) – exp ( σ2 ) ) 

一般的に書くと、aj 倍したい場合、

 1 / pow(a, 1 / j) * gaussian( σ, z / aj ) * (exp ( i * z / aj ) – exp ( σ2 ) ) 

となる。

なお、回転を施したい場合は、zに回転行列をzにかければよい。

z = (x, y)とすると、

  [x’]     [ cosθ    sinθ ]  [x]

  [y’]     [ -sinθ   cosθ ]  [y]

z’ = (x’, y’)

参考文献

https://ja.wikipedia.org/wiki/%E3%82%AC%E3%83%9C%E3%83%BC%E3%83%AB%E3%83%95%E3%82%A3%E3%83%AB%E3%82%BF

おすすめ記事

正準相関分析(Canonical Correlation Analysis:CCA)

はじめに

1つの目的変数を、複数の変数から説明するのが重回帰分析でした。正準相関分析は、複数の目的変数に対して、複数の変数がどう説明するのか見る方法です。

重回帰分析の際には、1つの変数を説明するために複数の変数を使っていましたが、正準相関分析は複数の変数を説明するために複数の変数を使います。重回帰分析の一般化、としてみることが出来ます。

数学的な意味の正準相関分析

正準相関分析とは、例えば2つの次元に貼っているベクトル( Va, Ua ):左上, (Vb, Ub ):左下があるとして、同じ空間(右図)に射影し、そのベクトルの相関を最大(コサイン積を最大)にするという見方もできる。上図のWa, 3次元から2次元への、Wbは2次元から2次元へのそれぞれ射影行列となる。索引:https://gregorygundersen.com/blog/2018/07/17/cca/

数学正準相関分析は、「複数の変数を2つの合成変数にまとめ、その合成変数間の相関係数を最大にする」という手法です。

つまり数に重み係数をつけて足し合わせたものを変量(正準変量、合成変量とも呼ばれる)として、その変量同士の相関関係である正準相関係数を最大にするような重み相関係数を求めます。

変量は回帰分析においてのy=axにおいてのyの部分です。回帰分析の場合、yの値が決まっていますが、正準相関分析ではyの値は決まっておらず、解析結果によって求めます。

正準相関分析は共通の空間(正準空間)に写像する手法です。

複数の情報源に共通して含まれる情報だけを抽出することによって、情報を統合するための次元圧縮法とも考えることが出来ます。

正準相関分析の具体例

理論だと少しむずかしいので具体的な例で説明します。「数学点数日本史点数握力遠投距離 」 のデータが5人分あったとします。

people1=[90, 40, 50,30]
people2=[80, 50, 55,10]
people3=[30, 80, 40,5 ]
people4=[50, 50, 35,50]
people5=[70, 80, 30,40]

同じ人たちの「走り幅跳び距離、英語点数」のデータもあります。

people1=[1.5, 40]
people2=[3.1, 30]
people3=[1.7, 50]
people4=[0.5, 80]
people5=[1.3, 90]

走り幅跳び距離と英語点数が、どう影響を受けるかということを考えることができます。

2つグループの相関関係を無理矢理最大化するような式をそれぞれ作ってやる、という事が正準相関分析となります。先に説明したとおり、解析結果後にそれぞれのグループで合成変量がでてきます。y=aXでいうところのyの部分です。このyの値をどう捉えるかはその観測者次第となります。

分析結果については、PCAのように寄与率と、最も相関があるように出来た重みが出てきます。N次元のデータであれば、N個の軸での解析が可能です。ただ、2つのグループの次元が違う場合がほとんどです。その場合は、最も長い次元にあわされますが、経験的に5−6個の軸から見れば大体解析できます。以下にCCAのコードを記載します。(察しの通り、正方行列ではないので特異値分解を中でおこなっています)

from sklearn.cross_decomposition import CCA
    cca = CCA(n_components=5)
    cca.fit(X1, X2) 

n_componentsは、寄与率の高い上から何個までの重み群と寄与率取ってくるのか、指定するものです。解析した異次元数を表しています。

CCAが中で何をしているかは、数学的な説明が多くなる為詳細を省きます。正準相関分析に関して解き方などの理論が気になる方は、君山 由良先生の多変量回帰分析・正準相関分析・多変量分散分析―多変量間の相関と因果関係の因子 (統計解説書シリーズ A-16) をおすすめします。非常に丁寧に解説されております。

参考文献

おすすめサイト

回帰分析(重回帰分析)

回帰分析の概念はとても簡単です。データにフィッティングする直線を見つけるという事だけです。例えば、身長と体重の相関関係を知りたいとします。直感的にですが、身長が高くになるに連れて、体重も増えそうですね。

体重 = a * 身長

となるような直線を推定する、これが回帰分析です。ここでよく見ると、y = a * xという式と同じです。( yが体重、xが身長)。つまり、回帰分析は傾きaを求めるための分析方法ということです。

さて、今の回帰分析は単回帰と呼ばれています。もっと複雑に体重を求めたいかもしれません。例えば説明変数「性別、宗教、胸囲、頭の多さ、身長」から、「体重」を求める、というものです。このように複数の説明変数から、体重を求めたい場合(説明変数の傾きを求めたい場合)は重回帰分析と呼ばれています。式は次のような形となります。

y = a_1x_1 + a_2x_2 + a_3x_3 + ...

説明の通り複数の説明変数がある場合が、重回帰分析と呼ばれています。重回帰分析でどうやって係数群をもとめるかは最小二乗法のコラムをご覧ください。

少し話はそれますが先程上げた「性別、宗教、胸囲、頭の多さ、身長」から「体重」を求めるという行為は実はやってはいけません。コンプライアンス的にという意味ではなく、重回帰分析のルールとして次のがあるためです。重回帰分析には次のルールがあります。

  • 説明変数は独立であること。
  • 説明変数は連続値であること。

身長と、胸囲は相関関係があり独立ではありませんね。身長が高くなると胸囲があがるという相関関係がありそうです。また、性別は0か1ですので連続値ではありません。宗教もイスラム、キリスト、仏教、神など連続ではないです。

なお、ルールを無視しても精度は出ます。ただ適切ではないです。論文などの公式なところに使えば突っ込まれるだけで、社内のデータなどでは使っても問題はありません。それでは、多変量解析のデータを解析する上でどういう方法を使えばいいかについて、迷った場合は以下を参照ください。

従属変数従属変数が量的独立変数が量的使用できる多変量解析の種類
ありYESYES重回帰分析
ありYESNO(質的)分散分析、多重分類分析、数量化理論Ⅰ類
ありNOYES重回帰分析、判別分析、ロジスティック回帰分析
ありNONO(質的)対数線形モデル、数量化理論Ⅱ類
なしYES因子分析、主成分分析
なしNO(質的)クラスター分析、多次元尺度法、数量化理論Ⅲ類

従属変数とは、y=a_1x_1 a_2x_2でいうところの、yが存在するかどうかを表します。質的なデータとは0,1で区別されるようなデータ、量的とは連続なデータか、0〜1の範囲でくまなく分布する連続するデータか、という意味を示します。

例えば、宗教(質的)と性別(質的)から癌になる確率(量的)を求めたい場合は数量化理論Ⅰ類を用いることとなります。

さて、重回帰分析では、一個の従属変数と複数個の独立変数の線形合成変数の相関が最大となるような重みを求めることを目的としています。(従属変数y、独立変数x、重みa)説明変数が一つなら単回帰分析。2つ以上なら重回帰分析となり、単回帰分析とは直線近似を指し、重回帰分析はそれの一般化で超平面を指します。

今、目的変数yに対し、説明変数xがどれだけの重みを持っているのかということを求めるのが重回帰、単回帰分析であると説明しましたが、改めて式にして表せば、

(単回帰分析) Y=aX

(重回帰分析) Y=aX // 太字になっているということで、ベクトルを表している。

終わりに

目的変数がひとつであることが回帰分析であると説明しました。重回帰分析は単回帰分析の一般化です。重回帰分析を更に一般化すると、正準相関分析となります。正準相関分析については、また別途解説いたします。

参考文献

https://ja.wikipedia.org/wiki/%E5%9B%9E%E5%B8%B0%E5%88%86%E6%9E%90

おすすめの記事

ロジスティック回帰分析

ロジスティック回帰分析とは、モデルがlog( \frac{1} {1 - p} ) = a_1 x_1 + a_2 x_2 + a_3 x_3 ...となるような{a_1,a_2,a_3}ベクトルを探すモデルです。(上記pはベクトル、x1..xnは今既知の値です。)

一旦モデルを作成さえすれば、あとは実測値を入れるだけで良いのがポイントで、このあたりは回帰分析と同じとなります。それでは一体何が違うのでしょうか。

ロジスティック回帰分析では、2値のデータの予測をしたいときに利用します。(男か女か、イスラム教かそれ以外か、日本人か外人かなど)。回帰分析では、値の推定に正規分布を利用しているため、正規分布に従わない、たとえば0か1かしか取らないものの推定が出来ないのです。少し専門的な言い方をすると、量的変数の解析には回帰分析、質的変数の解析(ただし2値に限る)にはロジスティック回帰分析を使います。量的変数は連続値のデータで、質的変数は離散的なデータと捉えて問題ありません。

出力が2値ですので、面白い使い方として、たとえばこれをしたら売れるか(1)、売れないか(0)みたいな出力結果を得られます。ですので、例えばマーケッターの方があるターゲットに売れるか売れないかの2値をとる計算などにロジスティック回帰分析が使えます。

以下のような関数を使います。

P = exp( Z ) / ( 1 + exp( Z ) ) =  1 / ( 1 + exp( Z ) )

ここで、上記の中のZは次のように表します。

Z = B_0 + B_1*x_1+ B_2*x_2+ B_3*x_3 …

補足ですがPの式はよくみるとシグモイド関数です。このBの係数(偏回帰係数
とよばれる)を求めるために式を変形します。

ln( P / (1-P )) = Z

この式をみるとわかるとおり、あとは、最小二乗法を用いて偏回帰係数を求めるだけで、回帰分析と全く一緒となります。これにより、ロジスティック回帰分析が出来ます。なお、最尤法(コンピューターをぶん回してフィッティング関数を見つける方法)を利用する方法でも解は求まるのですが、ほとんど最小2乗法を利用した計算法で事が足ります。

参考文献

おすすめの記事

最尤法

最尤法は、ある観察された離散的なデータを、連続的な関数にフィッティング(推定)する事を目的とした方法です。

今観測されたデータを正規分布や二項分布、ガンマ、ポアソンなど何でもいいのですが、連続関数にうまく当てはめてやる。うまく当てはまったということは、その際にその分布のパラメータの推定を行えていることを意味します。

今自分で適当に仮定した確率モデル(例えば正規分布など)があった場合、観測された実際の分布との適合度を示すものを「尤度」と呼びます。

尤度は高ければ高いほど、その仮定した確率モデルのフィッティングがうまく言っている事を表す。(尤度は、ある理論分布(正規分布など)を仮定したとき、現実に得られた実データと 同じ結果がその分布から生じる確率)

仮定する確率モデルのことを「尤度関数」と呼ぶ。(つまり観測結果をフィッティングした関数のこと)

尤度関数は単純に正規分布の場合や、複数の分布の混合の形で表される場合もあります。

最尤推定、最尤法とは、尤度を最大化するような尤度関数を決めてやる(パラメータを推定する)ことです。

参考文献

https://ja.wikipedia.org/wiki/%E6%9C%80%E5%B0%A4%E6%8E%A8%E5%AE%9A

おすすめの記事

FPマイニング(頻出パターンマイニング )

FPマイニングとはFrequency Pattern Mining の略で、物事の関連を見出すことを目的とします。たとえば、「おむつを買った人とはビールを買う」などそうした関連を見出すことを目的とします。どういうところで活躍しているのでしょうか?たとえば広告の技術などで活躍します。例えば、AというサイトをみてBというサイトを見たい人はCという商品が欲しい。それではCの広告を出そう!という事になります。こうした関連性を見出したいときはFPマイニングが利用されます。それでは説明に入っていきます。

アソシエーションルール

「アソシエーションルール」という言葉は、日本語で、連関ルール、関連ルール、相関ルール(⇐ 全部同じ意味を表す)と呼ばれています。記号で表すと以下のようになります。

{X}→{Y}

XであればYであるという意味で。よく例題で用いられるのはビールと紙おむつです。

  {ビール} →{紙おむつ}

上記の意味は「ビールを買った人はおむつを買う」という意味です。

FPマイニングはそのまんま実装すれば非常に簡単ですが問題があります。データ量が多くなると、とたんに難しくなるということです。これは、全ての組み合わせを探索する必要があるため、組み合わせ爆発が発生することによります。(要素数をXとすると2Xパターンが存在)

そのため工夫としてAprioriアルゴリズムやFP-Growth といった手法がとられており、実際にはFP-Treeを構築するFPGrowthがというアルゴリズムが用いられています。

アソシエーションルールの評価方法について

アソシエーションルールの評価方法には全部で3つの指標があります。

  1. 支持度( support )
  2. 確信度( confidence )
  3. リフト( lift )

結局は確信度が大事なわけなのですが、確信度は支持度(Support, sup)で算出されるため、Supを効率的に算出することを目標としています。

まずは、上記1,2,3それぞれの式がどのように表されるか説明いたします。

1. ルールX ⇒ Y の支持度

supp(X ⇒ Y) = σ(X ∪ Y)/ M

2. ルールX ⇒ Y の確信度

conf (X ⇒ Y) = σ(X ∪ Y)/ σ(X) もしくは、

conf (X ⇒ Y) = supp(X ⇒ Y) / supp(X)

3. ルールX ⇒ Y のリフト

lift(X ⇒ Y) = conf (X ⇒ Y) / supp(Y)

ここで、上記のσの意味は全体の集合から、要素Xを含むのが幾つあったのかを表すために用いられています。

例えば、集合が3つあったとして、それぞれの要素X, Y, Z, Wは次のようになっているとします。

{X, Y}, {X,Z}, {X, Y, Z, W} 

ここで、σ( X )は 3となります。そして、XとYをふくんでいる集合の数はσ( X ∪ Y ) = 2 となります。なお、集合の数は3つなので、式中にあるMは3となります。
それぞれの「支持」「確信」「リフト」に関して説明は以下のとおりです。

  • 支持度については全体の集合からそのルールが幾つあったのかという確率である。
  • 確信度については、ルールXがあったときに Yも発生する確率である。(XであればYも買う確率)
  • リフトについては 確信度を支持度で割っている。支持度がひくければ低いほど値が高く(データの中から出現度合いが低いとリフトが上がる)、確信も高いと高くなる。これは、確率P r ( X∪Y ) /Pr(X)Pr(Y)の近似値としても知られている。ちなみ、これは相関である。つまりLIFTは相関を表している
    • Liftは同じ確信があるときに、どちらの商品が人気あがあるのかがよくわかる。例えば、0.3/0.1, 0.3 / 0.05だと後者の方が大きい。除数は支持率(全体の数である)

3つの指標において、最終的に関連が高いものは次の二つの事が言えます。

  • 確信度conf が高い
  • 支持度support が高い

FPマイニングはこうした数値を出して、確信度が高い!支持度が低かった!など分析していくことです。

前述のとおり、全探索すれば関連性は見いだせるのですが全てでるので良いのですが、実際には不可能です。そのため、以下のような簡略化するアルゴリズムが用いられます。

Support(頻度)さえでれば、Confidenceが求まることになります。Confidence( ビール →おむつ ) = support( {ビール、おむつ} ) / support( ビール )

ですので、どのように効率的にSupportをだすかが鍵となります。

Apriori アルゴリズム

全データ数を見て、その関連を見つけるのは大変です。(組み合わせ爆発問題)。そのためAprioriアルゴリズムによって効率的に上記のアソシエーションルールを満たすものを探し出します。以下はその様子を図示した内容です。

集合から、要素を取り出し、その数をカウントしてC1を作成します。しきい値(この場合は1とした)より低いものを無視して、L1を作ります。L1の要素のすべての組み合わせからC2を作ります。ここで、C2のセットで1つ減らした要素( {1 2}であれば{1}と{2}が) L1に入っていない場合、C2の要素としてカウントは出来ません。C2作成に当たりその例はありませんでした。ただC3でその事例が発生します。

C2からスキャンして、先程と同様にしきい値より下のものを削除しL2をつくります。以下同様にC3をつくります。先に説明したとおり、C3の部分集合がL2に含まれていない場合は、その部分集合は入れてはいけません。そのため、{2 3 5}のみが候補として残ります。

さて、Supportさえ出てしまえば、Confidenceを算出することが出来ます。Aprioriアルゴリズムの問題点はメモリです。100個のアイテムから候補を作成するのに2の100乗の候補が必要になり、メモリを圧迫するうえ、データのスキャンにも時間が必要となります。

参考文献

おすすめの記事