Transformer(トランスフォーマー)

今回はNLP界でも話題沸騰中のトランスフォーマーの仕組みについて解説していこうと思います。Transformerの前に、Attentionの技術の理解が必要で、Attentionがどんな物か簡単に(人のYoutubeビデオで)紹介し、そのあとTransformerの技術説明をしていこうと思います。TransformerはDeep Learning 研究の中でも非常に大きなブレイクスルーを起こした技術の一つです。Transformerの派生版である GoogleのBERTやGPT-3など聞いたことがあると思いますが、翻訳のレベルが数段上がりました。

本ブログではできるだけ数式を使わず、簡単に説明したいと思っているのですが、所々に難しいところがあると思います。疑問点など有りましたらご連絡をいただければ幸いです。

Attentionとは

Transformerを学ぶ前にTransformerの目玉の一つでもあるAttentionという技術はどういう技術か紹介します。

「Attention」は注目という意味となります。つまるところAttentionとは、入力の注目する部分を見つける!ということになります。物を学習する時に、画像処理で言えば例えば犬だけに注目して処理したもの(犬を全景処理で切り出したもの)と、画像全体で処理したものを使っての学習では、効率や結果が変わりそうです。この犬だけに注目したい!というのがAttentionという技術となります。この技術は文章にも適用ができ、例えば翻訳の部分では翻訳する際にどこにまず注目した方が良いのか、や、文章の単語同士の関わりの指標みたいなのを抽出することが出来るようになるのです。

言葉で説明してもわかりにくいかもしれません。この技術をYoutubeで15分で紹介している方が居ました。素晴らしく解りやすいので、是非ご覧になってください。

Attentionは対象物のみに集中する、という意味で逆に言えばそれ以外を無視するという雑な言い方も出来ます。実際に我々の脳でも同じことをやっています。例えば集合写真で人数を数えてくれって言われたら、頭の数を数えてほかは無視するようにすると思います。このように、人間と同じようなことをするという方法がAttentionという方法です。

Seqence to sequence with the Attention

さて、Seq2Seqモデルで翻訳を前回の記事では実装しました。Seq2Seqは素晴らしい技術ですが長い文章の翻訳が弱いという弱点があります。それはエンコーダーからデコーダーへ渡されるものがContextベクトルだけであるため、例えば翻訳の最初の部分の表現が弱くなってしまう可能性があるという問題点があったためです。Seq2Seqのモデルはわかりやすく表すと下図のような仕組みとなっていました。なおSeq2Seqに関しては前回の記事を参照ください

Seq2Seqの流れ。コンテキストがバケツリレーでエンコーダーからデコーダーへと渡されている。図は参考文献[1]より

Attentionはこうした問題に対応するべく産まれました[ 2 Bahdanau et al, 2015 ]。Seq2Seqで実際に、どのようなことをAttentionを実現しているかと言えば、雑に言うと、seq2seqのエンコードの部分で利用したHiddenベクトルを、まとめてデコーダー部分にわたしている、というのが特徴になります。上のSeq2Seqのモデルと比較してみてください。

Seq2Seq with Attention方式。Encoderで利用したHiddenベクトルもわたし、そのベクトル群から注目するワードは何かという処理をDecoder内で行う(Attention処理)

Seqence to seqence でAttentionの実現には、デコーダー内でまとめて渡されたHidden ベクトルの内どれが重要であるかを単語単位で調べていき、重み付けをしてSoftmaxし、加重平均しています。詳細な流れですが長くなるので説明しません。というのも、その方法よりもTransformerの方が精度があがっているため、いまさらSeq2SeqのAttentionの実装方式にフォーカスするよりは、Transformerを学んだほうが有用であるためです。(なにより、Transformerの説明の記事なので)

Seq2Seq with Attentionの実現方法については、リクエストが有れば別の機会にご紹介できればと思います。なお、Seqence to sequence のAttentionの実装方法も奥が深く、Global, local attentionなどの手法も提案されています。

Transformer

本題のトランスフォーマーです。トランスフォーマーが産まれたのは、雑に言えばSeq2Seq+Attentionの並列化をして高速に学習したかったからです。

並列化と聞いてまず思い浮かぶのがCNNなどコンボリューションで、並列処理が得意なGPUに向けて行列演算に落とし込めば並列化が見込めます。単純にはSeq2Seqはバケツリレー形式で情報を運んでいる(前の入力に対しての処理をしている)ので並列化が叶いません。

そこで、コンテキストを送ってAttentionを計算していたところから、Encoder,Decoder各層内でself-attention層をつくり、行列演算を実現することで並列化を目指します。

Attentionという名前がついていますが、今までのSeq2Seqで扱っていた方式と全く別なので注意してください。Seq2Seqでは全体でAttentionを実現しているイメージでしたが、Transformerでは層という単位に落とし込んで実現しています。下記に簡単にTransformerの全体図を示します。

トランスフォーマー全体図。例では左のエンコーダーは2段。右のデコーダー1段構成。実際の論文ではEncoder6段、Decoder6段の構成となる。
右側DecoderのOutputでは最後にSoftmax処理が入る。Encoder,Decoderで、Self-attention層があり今までのSeq2Seqとは様相が異なる。

上図がTransformerの全体図になります。Encoder部とDecoder部の2部構成はSeq2Seqと同じです。Encoder,Decoderとも多段で接続(オリジナルの論文では6段)で構成しています。

Transformerでの肝はずばりSelf-attention層となっています。Self-attention層は同じセンテンス内の単語群が、翻訳したい単語にどれだけ影響しているかというものを検出します。

EncoderとDecoderともSelf-attetionやNormalizationなど殆ど同じ層を使いますが、唯一Decoderにしかないものがあり、Encoder-decoder attention層です。後ほど説明しますが、Encoder側からバイパスしてもらう特徴量(最後の段のEncoderから2つの特徴行列(Key, Valueを生成する行列))を用いいてSelf attentionを行います。これにより入力のどのパーツに注目するべきかを考えることができます。なお、全てのDecoder群のEncoder-decoder attention層では、最後のEncoderの特徴(key, valueの生成行列 )を使うことになります。

Transformer 説明の流れ

全体の流れがサラッとわかったところでそれぞれの層がどのような処理をしているかを記載していきます。それぞれの層の説明、特にSelf-attentionはなかなか濃ゆいと思うのでじっくり理解していってください。

Encoder、Decoderそれぞれでは

  1. 事前処理Embedding
  2. 事前処理Positional Encoding
  3. self-attention層( + Encoder-Decoder attention層( decoderのみ ) )
  4. Normalization層
  5. Feed forward層

というような処理を順次行っています。それぞれの処理について説明していきます。

Embedding処理

Embeddingは文字の埋め込みを行います。Word2Vecの様な処理です。具体的には例えば512次元のベクトルに単語を変換します。例えば I am an studentであれば、

I = [ 0.1, 0.2, 0.03, 0.78 …. 0.04 ], am = [ 0.41, 06., 0.033, 0.378 …. 0.3 ], an = [ 0.41, 0.52, 0.073, 0.5678 …. 0.34 ], student = [ 0.81, 0.62, 0.023, 0.478 …. 0.54 ]

のように一つ一つの単語に対してベクトルに変換をします。pytrochではEmbeddingは関数一つで実現でき、上記例では512次元のベクトルに変換するのにembedding = nn.Embedding(4, 512)としておき、embeddingを呼ぶことで実現できます。Embeddingの方法やアルゴリズムについては[10]を参照してみてください。要するにここで文字列をベクトル化してやって学習するようにするのです。

Positional Encoding(PE)処理

Positional Encoding (PE)は順序の定義を行います。今回のTransformerでは行列演算により並列化し従来よりも高速に計算をするという目的がありました。並列に行列計算を行うため、文章の単語に対して順序(もとにいた位置)を考慮させる必要があります。Seq2Seqのようにバケツリレーではないので、バラバラにしてももとの位置が解るような情報量を与える必要があります。

順序でいえば、例えば0,1,2,3とindex値を割り振ればいいと思うかもしれませんが、文字数が大きくなると、値も大きなってしまいます。そもそも-1〜1の範囲で正規化して学習するため整数のindex値は相性が悪そうです。じゃあ正規化をすれば良いのではないかと思うのですが、今度は別の問題が発生します。例えば、5文字のセンテンスがあったとして0.8は四番目の文字を意味します。所が10文字の0.8は8番目になります。つまり、同じ値に場所が違う問題が発生します。

では、文字列の次元を固定(例えば512次元)して、割り振ればいいと思うかもしれませんが、ほとんどの文章は20単語くらいでしょう。そうすると、ものすごい無駄な空間(20次元以降の情報は使われることは殆ど無い)となります。

さて、こうした問題を解決したのが余弦、正弦を使った順序の定義です。数式は次のとおり表されます。

上記式を例えば128次元(d_model)と仮定すると、以下のような画像となります。[11]

positional encoding
各行が文字の位置のベクトルとなる。埋め込みされたベクトルと足し合わせる

Embeddingしたinputを、このPositional Encodeing の値と足し込めたものをインプットとして利用します。

x = Embedding( input ) + PositionalEncoding( input )

Self-attention層

肝のSelf attention層です。センテンスにある各単語同士の注目度を算出していきます。注目度の算出は雑に言うと内積で求めていきます。

Self-Attentionでは、Embedding+PEされたベクトルから、Query, Key, Valueという3つのベクトルを作成します。この3つのベクトルを作るためには3つの行列が必要なのですが、Self-attention層ではその3つの行列を学習で作ることが目標です。

注目度である内積値はQuery, Key, Valueのベクトルから求めます。3つのベクトルの次元数は全て同じで、入力のEmbedding+PEされた次元数より少なくなります。例えば、オリジナルの論文では512次元のEmbedding+PEベクトルがあるとすると、Query, Key, Valueベクトルは64次元まで減っています。

さてさて、ここでいうQuery, Key, Valueとはなんでしょうか?

わたしもそうでしたが、Query, Key, Valueって何で3つに分解するんだろう、それぞれの持つ意味は何だろうかと最初は考えては調べていき、最後は思考停止に陥っていました。今は取り敢えず深く考えずに3つに分解する、それを以ってなんとかSelf-attentionを実現する、という考えで読み勧めてください。

上で説明した話の概念図を下に示します。

今ここで, You love me( お前は俺を愛してる)の3ワードのセンテンスを例としましょう。説明のため、各ワードは16次元でembedding+PE されているとします。まずは、16次元のワードを3次元のベクトルに落とし込みます。(3ワードの”3”とは無関係です。短縮次元はいくつでもいいのです。オリジナルの論文では8で割って512次元を64次元にしていました)。それぞれの生成行列を定義します。初期値はランダムです。

Query, Key, Valueの生成行列を学習することをSelf-attention層では目的とする。各ワードにかけ合わせ、それぞれの単語に対するquery, key, value を作る。16次元を3次元に変換する。

上記処理を行ったあとは、ワードに対するそれぞれのquery, key, valueベクトルを用いて内積値の計算を実施します。まず、Youを例に取り、You, love, meのKeyベクトルとValueベクトルを使って結果zを取り出します。計算処理の内容については下記のとおりです。

self-attention層計算
Youのqueryベクトルに対して、You, Love, MeそれぞれのKeyベクトルと内積値を取る。現在次元数は16次元から3次元(繰り返すが単語数の3とは関係ないので注意)になったので、そのsqrt(現在の次元数3)で除算し、それとValueベクトルとかけ合わせ足し合わせる。なおSoftmaxの合計は必ず1(0.67+0.22+0.11=1)であるため、加重平均ともとれる。これで算出された値をz1ベクトルとする。これを、Love(q2), Me(q3)に対しても行いz2, z3も求める

これが、Self-attention層で行われている処理です。

なおこの処理を複数回行うのがMulti-head attentionという処理になります。今z1だけが出てきましたが、z2, z3 …と出力していきます。オリジナルの論文では8つのマルチヘッドに分割していました。最初に512次元から64次元に落とすために8で除算していたので、そのため8なのかもしれません。なお、8つのマルチヘッドで分割していた場合、それに伴って生成行列も8個ずつ出来ます。

出力されたz1 .. znを連結(concat)して、もとの入力次元に戻す(linear)作業をすることで、最終的な処理は終わりになります。例えば、今回16次元の入力から3次元になりました。8つのヘッドを使ってConcatするとz1 .. z8までの合計で24次元になります。それを入力次元の16次元に戻す処理、つまり射影マトリックスをかける処理(nn.Linear)をすることで処理が終わります。

ノーマライゼーション層(+ residual の説明 )

標準化、正規化、規格化など、呼ばれるノーマライゼーションにより安定した学習が出来るようになります。中で実装している内容ですが、値を平均から引き、標準偏差で割るような操作をしてやります。平均が0で、分散が1になります。まさにノーマライゼーションですね。

言葉だとわかりにくいので、どういう処理をしているのか、ある人[6]がソースコードを書いていました。参考にしてみてください

class Norm(nn.Module):
    def __init__(self, d_model, eps = 1e-6):
        super().__init__()
    
        self.size = d_model
        # create two learnable parameters to calibrate normalisation
        self.alpha = nn.Parameter(torch.ones(self.size))
        self.bias = nn.Parameter(torch.zeros(self.size))
        self.eps = eps
    def forward(self, x):
        norm = self.alpha * (x - x.mean(dim=-1, keepdim=True)) \
        / (x.std(dim=-1, keepdim=True) + self.eps) + self.bias
        return norm

residual(残差処理) の説明

残差処理といいますが、Encoderを構築する際に、Multi header attention で出てきた出力値に、元々の入力Inputを足し合わせます。そしてNormalization処理をします。ソースコード風で簡単に表すと

y = Normlization( MulitHeadderAttention( input ) + input )

これを残差処理といいます。下記のオリジナルの論文を見ると、Add&Norm内で入力値を持ち込んでいるので残差処理をしているのがわかります。

オリジナルの論文の画像。Add&Normの処理が残差処理として利用される。入力のinputと、attentionで出てきた出力を足しこみ平均化(標準化)処理を行っている

Feed Forward 層

Feed Fowardでは、2層の全結合層からなるニューラルネットワークを実装しています。また、Dropoutも定義することで過学習を抑制します。Dropoutに関しては、次の記事が明るいです。[14]

数式ではつぎのようなことを行っています。

FFN(x)=max(0,xW1+b1)W2+b2

図では次のとおりです。いえば、ただのニューラルネットワークです。

よく見られる2層のニューラルネットワーク。中間層のサイズを大きくして表現力を上げる作業をここでは行う。尚、Dropoutも実施する。

ソースコードだと次のような単純な処理をしているだけです。

class FeedForward(nn.Module):
    def __init__(self, d_model, d_ff=2048, dropout = 0.1):
        super().__init__() 
        # We set d_ff as a default to 2048
        self.linear_1 = nn.Linear(d_model, d_ff)
        self.dropout = nn.Dropout(dropout)
        self.linear_2 = nn.Linear(d_ff, d_model)
    def forward(self, x):
        x = self.dropout(F.relu(self.linear_1(x)))
        x = self.linear_2(x)
        return x

ここまでで、Encoderで使われている全てのパーツの説明がおわりました!実はEncoderもDecoderもほとんど同じパーツを使います。唯一,Encoder-decoder attention層だけがちがうのですが、その解説を行います。

Encoder-DecoderAttention

改めてTransformerの全体図をみてみると、次のようになっています。

Transformer 全体図 [8]より。左のEncoder、右のDecoderとも多段で縦に接続されているが、説明上それぞれ一段構成になっている。

ブロックを見てみると殆どの解説が終わっている部分ばかりです。唯一Encoder-DecoderAttentionだけ説明をしていませんでした。ただ、Encoder-Decoder Attentionは次の画像がとてもわかり易いかと思います。

encoder で使われたValue, Keyを渡す。 デコーダーではquery を使う。

やっていることはSelf-attentionと一緒です。ただし、Encoder-Decoder層では、最後のEoncoderパートで使われたValue, Keyベクトルを生成するマトリックスを渡し、Self-attention層で計算した方法と全く同じロジックで計算するのです。これだけとなります。

すべてのパーツの説明が終わりました。後は実装してみましょう!

実装

実装については、ブログの公式Githubにアップロード予定です。

最後に言い訳…

2021年頃にTransformerの記事が書きたいと思って、Twitterで書くと公言していました。しかし、コロナ渦でてんやわんやだったり、そもそもTransformerのボリュームが有りすぎて、画像等の用意を考えると時間的になかなか書けませんでした。そうこうしている間に2022年では素晴らしい記事が出てきてしまい、書くのは意味ないと思ってしまい、全く手が進みませんでした。しかし最近改めて日本語のサイトを見ると、日本語では全体的に網羅されていない、数学的すぎる、簡潔すぎる、論文の解説に集中している、というような本ブログのようなまとめ解説は有用なのかもしれない、と思い、思い切って書いてみました。皆さんの勉強の参考になってくれれば幸いです。

参考文献

[1] https://jalammar.github.io/visualizing-neural-machine-translation-mechanics-of-seq2seq-models-with-attention/

[2] https://towardsdatascience.com/transformers-explained-visually-part-3-multi-head-attention-deep-dive-1c1ff1024853

[3] https://ai.googleblog.com/2017/08/transformer-novel-neural-network.html

[4] https://jalammar.github.io/illustrated-transformer/

[5] https://mchromiak.github.io/articles/2017/Sep/12/Transformer-Attention-is-all-you-need/#.XIWlzBNKjOR

[6]https://towardsdatascience.com/how-to-code-the-transformer-in-pytorch-24db27c8f9ec

[7]https://machinelearningmastery.com/a-gentle-introduction-to-positional-encoding-in-transformer-models-part-1/

[8]https://towardsdatascience.com/transformers-explained-visually-part-2-how-it-works-step-by-step-b49fa4a64f34

[9]https://towardsdatascience.com/how-to-code-the-transformer-in-pytorch-24db27c8f9ec

[10]https://pytorch.org/tutorials/beginner/nlp/word_embeddings_tutorial.html

[11]https://kazemnejad.com/blog/transformer_architecture_positional_encoding/

[12]https://medium.com/swlh/elegant-intuitions-behind-positional-encodings-dc48b4a4a5d1

[13]https://data-science-blog.com/blog/2021/04/07/multi-head-attention-mechanism/

[14]https://medium.com/axinc/dropout%E3%81%AB%E3%82%88%E3%82%8B%E9%81%8E%E5%AD%A6%E7%BF%92%E3%81%AE%E6%8A%91%E5%88%B6-be5b9bba7e89

[15] https://pytorch.org/docs/stable/generated/torch.nn.LayerNorm.html

[16]https://kikaben.com/transformers-encoder-decoder/

Variational AutoEncoder( VAE )

近年、ディープラーニング業界、はたまた画像処理業界ではGANとVAEの2つの技術で話題が持ちきりとなりました。GANについては前回解説しました。今回はそんなVAEのアルゴリズムについて解説をしていきます。GAN同様に最後にPyTorchによる実装例を紹介していきます。

VAEとは

VAEは情報を圧縮して圧縮した情報から元の情報に戻す、というような仕組みをもった、AE(Auto encoder)と言われるものの一種です。AEはただ単にデータの圧縮と再構築をするだけでしたが、VAEは確率モデルの理論を導入しています。VAEは確率分布を使ったモデルということは、未知のデータを確率的に作成できることになります。VAEはGenerative model(生成モデル)と言われています。トレーニングデータからProbability Density Function (PDF)を推定するモデルであるためす。 

GANと違って何ができるのか

VAEの一つの特徴は次元削減です。次元削減といえばPCASVDなどを思い浮かべるかもしれません。VAEは同様にLatent space(潜在空間)に次元圧縮し、またそこから復元するということをしています。

PCAのように圧縮した次元において直行している必要はありません。この辺も普通の次元圧縮とのアプローチとは異なっています。

赤はAEの軸。青はPCAの軸。PCAは軸同士は直行しているが、AEに関しては必ずしも直行している必要はない。

AEとVAEは違いがあまり無く感じますが、前述の通り確率を用いたところに違いがあります。AEは点(Single point)としてデータを潜在空間にマッピングしていましたが、VAEでは潜在空間に落とし込むときにガウス分布に従って落とし込ませます。デコーダーはその点を拾ってデコードするため、見方を変えれば、確率的にサンプリングしてzを拾っている、ということになるのです。このzをデコードして元の入力に近づけるようにするのです。

さて、VAEは確率モデリングと言われます。確率モデリングとは、例えば何かのデータxが分布していたとして、その分布を確率で表現するというモデルです。あるデータxが正規分布に従っていたとしたら、データxの確率分布は正規分布となり、P(X)=正規分布( μ, σ )と言う感じで表現します。P(X)と書くと思わずPは関数と勘違いするかもしれませんが、確率を表すものなので注意が必要です。

数式

VAEはデータを確率モデル化をすることを目標とします。データをX、その確率分布をP(X)と定義すれば、P(X)を見つけること、すなわちP(X)の最大化がVAEの目的です。そしてP(X)に関しては次のように表現することができます。[1]

P(X) = \int P(X \vert z) P(z) dz = \int P(X,z)

上記の式はXが入力画像と考え、その画像を表現する潜在ベクトルzを少なくとも1つ見つけたいという意味があります。

また、事前分布P(z)から\{z_i\}_{i=1}^nをサンプルした際に、P(X)を次のように近似できます。

P(X) \approx \frac{1}{n}\sum_{i=1}^n P(x|z_i)

P(X)を求めればいいのですが、一筋縄には行きません。それはXが高次元であることもそうですが、その確率を求めるのにはすべてのzの空間を舐め回すような非常にたくさんのサンプリングが必要なだけでなく、組み合わせをチェックするような処理が必要となり、総当りでP(X)を求めることは現実的ではありません。

そのため、P(X)を求めるには別のアプローチを考えます。もし事後確率のP(z|X)がもとめられれば、p(x|X)=\int P(x|z)p(z|X) dz より未知の画像xを作り出すような確率分布を求めることができそうです。そうなると 事後確率のP(z|X)を求めるだけとなり、簡単に思えますが事後確率はベイズの定理により次の式になります。

P(z|X) = P(X|z)\cfrac{P(z)}{ P(X)}

よく見るとあの厄介なP(X)が分母にあります。やはり求めることができません。

ここで、諦めずに別のアプローチで P(z|X) を求めていくようにしていきます。

P(z|X) を求める

P(z|X)を求めるにあたり少し工夫します。ここでVAEの名前の由来ともなっていますが、 Variational Inference(VI) という手法を使います。P(z|X)を推定するためにQ(z|X)という分布を考えます。Q(z|X)はP(z|X)の近似で、ガウス分布などの簡単な分布関数の組み合わせによりP(z|X)を近似します。

関数の組み合わせにより、近似していく。青の点線へ緑の分布関数で近似していく様子。 https://towardsdatascience.com/bayesian-inference-problem-mcmc-and-variational-inference-25a8aa9bce29

どれだけ似ているのか、という指標については類似尺度としてKLダイバージェンスを使います。同じ分布であれば0に、異なれば値が大きくなっていくような関数です。KLダイバージェンス次のように書くことが出来ます。

KL(P||Q) = E_{x \sim P(x)} log \cfrac{P(x)}{Q(x)}=\int_{-\infty}^{\infty}P(x)\ln \cfrac{P(x)}{Q(x)}dx

ここから式が多く出てくるので、見やすくするために単純にP(z|X )をP, Q( z|X)をQと表示します。Qの近似は次のようにかくことができます。

\begin{aligned} D_{KL}[Q\Vert P] &= \sum_z Q \log (\cfrac{Q}{P}) \\ &=E [ \log (\cfrac{Q}{P}) ] \\ &= E[\log Q - \log P] \end{aligned}

最後はlogの変換公式により、割り算を引き算にしただけです。さてP である P(z|X ) は P(z|X) = P(X|z)P(z) / P(X) のように表現できたので、先の式に当てはめてみましょう。

\begin{aligned} D_{KL}[Q\Vert P] &= E[\log Q - (\log (P(X \vert z) \cfrac{P(z)}{P(X)})] \\ &= E[\log Q- \log P(X \vert z) - \log P(z) + \log P(X)] \end{aligned}

ここでzに関する期待値でくくっている項に注目するとP(X)があります。P(X)はzに関係ありません。そのため括弧の外に出すことができます。

\begin{aligned} D_{KL}[Q\Vert P] &= E[\log Q- \log P(X \vert z) - \log P(z)] + \log P(X) \\ D_{KL}[Q\Vert P] - \log P(X) &= E[\log Q- \log P(X \vert z) - \log P(z)] \end{aligned}

ここからがトリッキーなのですが、右辺に注目するともう一つのKLダイバージェンスを見つけることができます。先の式を両辺にマイナス倍します。

\begin{aligned} \log P(X) - D_{KL}[Q\Vert P] &= E[-\log Q + \log P(X \vert z) + \log P(z)] \\ &= E[\log P(X \vert z) -( \log Q - \log P(z))] \\ &= E[\log P(X \vert z)] -E[( \log Q - \log P(z))] \\ &= E[\log P(X \vert z)]- D_{KL}[Q\Vert P(z) ] \end{aligned}

なんとPとQを近づけようとした結果、本来の目的であったP(X)に関する式がでてきてしまいました。PとQを縮約して記載していたので、正しく展開してかいてみます。

\log P(X) - D_{KL}[Q(z \vert X) \Vert P(z \vert X)] = E[\log P(X \vert z)] - D_{KL}[Q(z \vert X) \Vert P(z)]

本当の目的はP(X)の最大化でした。結局ですが、これを最終的なVAEの目的の関数として定義することになります。

最終的に得た式は非常に興味深い構造になっています。

  1. Q(z|X) はデータXを受取り、潜在空間にzを投影
  2. zは潜在変数
  3. P(X|z)は潜在変数zからデータ生成

つまりQ(z|X)はエンコーダ、zは潜在変数(エンコードされたデータ),そしてP(X|z)はデコーダです。まさにオートエンコーダーそのものとなりました。

算出された式を観察して、左辺と右辺がありますが、左辺にあるlog(p)最大化することが目的でした。そのためには右辺を最大化していくこと、つまり E[\log P(X \vert z)] を大きくして D_{KL}[Q(z \vert X) \Vert P(z)]を小さくしていくことで、左辺は大きくなっていきます。

そのため、今度は目的を右辺の最大化に絞っていきます。

右辺を最大化する。

右辺には以下の2つの式があります。

  1. E[\log P(X \vert z)]
  2. D_{KL}[Q(z \vert X) \Vert P(z)]

右辺=数式1-数式2、ですので数式1を最大化、数式2は最小化していけば、右辺は大きくなり目的が達成されます。

まず 数式 1についてですが、よく見ると潜在変数zを受取りXを生成する教師つきの学習そのものです。ですので、学習によりなんとかなりそうです。

さて、厄介なのが 数式 2です。ここで一つの仮定を起きます。P(z)は正規分布N(0, 1)と仮定するのです。そして、Xからzを生成する分布もパラメータ\mu(X), \sigma(X)付きの正規分布となります。 平均と分散はXを中心としたという意味です。そして、KLダイバージェンスは次のように表されます。

D_{KL}[N(\mu(X), \Sigma(X)) \Vert N(0, 1)] = \frac{1}{2} \, \left( \textrm{tr}(\Sigma(X)) + \mu(X)^T\mu(X) - k - \log \, \det(\Sigma(X)) \right)

kはガウシアンの次元数、traceは対角要素の和を表します。そして、detは対角要素の積\det \left({\mathbf A}\right) = \prod_{i \mathop = 1}^n a_{ii}です。

導出に関しては、ここでは重要でないため割愛しますが最終的には次のようになります。(導出に関して興味ある方は[2]を参照してください。)

D_{KL}[N(\mu(X), \Sigma(X)) \Vert N(0, 1)] = \frac{1}{2} \sum_k \left( \exp(\Sigma(X)) + \mu^2(X) - 1 - \Sigma(X) \right)

この項目を実装するときにはロス関数の中で利用します。実際との差分を計算するためですね。

全ての要素の説明ができました。あとはやるだけです。

実装

エンコーダーとデコーダーの実装は次のとおりです。

エンコーダー

class Encoder( nn.Module ):
    def __init__( self ):
        super().__init__()
	self.common = nn.Sequential(
            nn.Linear( 784, 400 ),
            nn.ReLU(),
            )
	self.model1 = nn.Sequential(
            self.common,
            nn.Linear( 400, 20 )
            )
        self.model2 = nn.Sequential(
            self.common,
            nn.Linear( 400, 20 )
            )
    def forward( self, img ):
	img_flat = img.view( img.size( 0 ), -1 )
        return self.model1( img_flat ), self.model2( img_flat )

デコーダー

class Decoder( nn.Module ):
    def __init__( self ):
	super().__init__()
	self.model = nn.Sequential(
            nn.Linear( 20, 400 ),
            nn.ReLU(),
            nn.Linear( 400, 784 ),
            nn.Sigmoid(),
            )
    def forward( self, z ):
        return self.model( z )

エンコーダーとデコーダーを用いいたVAEを次のように実装していきます。

VAE

class VAE( nn.Module ):
    def __init__( self ):
        super().__init__()
        self.encoder = Encoder()
	self.decoder = Decoder()

    def _reparameterization_trick( self, mu, logvar ):
        std = torch.exp( 0.5 * logvar )
        eps = torch.randn_like( std )
        return mu + eps * std

    def forward( self, _input ):
        mu, sigma = self.encoder( _input )
	z         = self._reparameterization_trick( mu, sigma )
        return self.decoder( z ), mu, sigma

説明ではzの取得は正規分布でのサンプリングを仮定しました。ところが実際にサンプリングするとバックプロパゲーションが出来ないため学習が出来ません。そこでreparameterization trickというのを用いいます。zを次の式で近似します。

z = \mu(X) + \Sigma^{\frac{1}{2}}(X) \, \epsilon

where, \epsilon \sim N(0, 1)

左図はZを表現するために本来のサンプリングで実装した図。左ではバックプロパゲーションができない。そのため、足し算と掛け算に依る表現に正規分布に従うノイズの足しこみを行い誤差伝搬を可能にする。[3]

最後に損失関数とVAEを使ったコードは次のとおりです。

# Kingma and Welling. Auto-Encoding Variational Bayes. ICLR, 2014                                                                                                                                                                 
# 入力画像をどのくらい正確に復元できたか?                                                                                                                                                                                        
def VAE_LOSS( recon_x, x, mu, logvar ):
    # 数式では対数尤度の最大化だが交差エントロピーlossの最小化と等価                                                                                                                                                              
    BCE = F.binary_cross_entropy(recon_x, x.view(-1, 784), size_average=False)
    # 潜在空間zに対する正則化項. # P(z|x) が N(0, I)に近くなる(KL-distanceが小さくなる)ようにする                                                                                                                               
    KLD = -0.5 * torch.sum(1 + logvar - mu.pow(2) - logvar.exp())
    return BCE + KLD

def main():
    epoch_size = 50
    vae = VAE()
    vae.cuda()
    Tensor = torch.cuda.FloatTensor
    dataloader=get_dataloader()
    optimizer = torch.optim.Adam( vae.parameters(), lr=1e-3 )

    for epoch in range( epoch_size ):
        for i, ( imgs, _ ) in enumerate(dataloader):
            optimizer.zero_grad()
            real_images          = Variable( imgs.type( Tensor ) )
            gen_imgs, mu, logvar = vae( real_images )
            loss                 = VAE_LOSS( gen_imgs, real_images, mu, logvar ).cuda()
            loss.backward()
            optimizer.step()

VAE_LOSSでは得た画像が目的とした式と同じになるような差分を計算しています。

長くなりましたがこれがVAEの全貌です。今回のソースコードはGithubにあげてあります。

https://github.com/octopt/techblog/blob/master/vae/main.py

参考文献

  1. https://en.wikipedia.org/wiki/Law_of_total_probability
  2. https://wiseodd.github.io/techblog/2016/12/10/variational-autoencoder/
  3. https://towardsdatascience.com/understanding-variational-autoencoders-vaes-f70510919f73

おすすめの記事

最小二乗法による重回帰分析

重回帰分析の値を最小二乗法でえる方法を示す。早速だが、次の式に適用することで目的の係数群が手に入る。

\mathbf{w} = (\mathbf{X}^{T}\mathbf{X})^{-1} \mathbf{X}^{T}\mathbf{y}

上記x, y, wは次の通りの関係を持っている。

y_i=a_{0i} + x_{0i}w_{0i} + x_{1i} w_{1i} + ... x_{ni} w_{ni}+ f

aは切片、fは残差と言われる。ここでyはm個あると仮定する。回帰分析のところで説明したとおり、切片aと重みwベクトルを求めるのが回帰分析の目的となる。最初の\mathbf{X}については次のようになる。

\mathbf{X} = \begin{bmatrix} \mathbf{x}_0 & \mathbf{x}_1 & … &\mathbf{x}_n \end{bmatrix}

\mathbf{x}_0 = \begin{bmatrix} x_{00} \\ x_{01} \\ ... \\ x_{0m} \end{bmatrix}

となっている。

上記ではやや理解がしづらいと思うので、実際に具体例をしめす。なお、他のサイトでは切片を求めない場合もある。当サイトでは、今回は切片を考慮に入れて説明しよう。

LSMを用いた重回帰分析の仕方(計算方法具体例)

具体的に説明のためにあるサイト[2]から取得したデータを用いる

営業所売上額広告費販売員数
A85006
B95008
C1370010
D1140013
E1480011
F17120013
G?130014

上記で営業所Gの売上高を予想したいというのが今回の目的である。(そのために重みベクトルと切片を求める。)上記の表で例えば営業所Aの売上高は 8 = a_0 + 500 * w_0 + 6 * w_1 ということになり、最小二乗法を使うことで未知の切片のa_0と重みである\mathbf{w}が求まる。

まずは上記表を元にa0, x0, x1のベクトル、並びにyのベクトルを定義する。

import numpy
def main():
    y = numpy.array( [[8, 9, 13, 11, 14, 17]]).T
    a0 = numpy.array( [[ 1, 1, 1, 1, 1, 1]] ).T # 切片が不要な場合はここは不要
    x1 = numpy.array( [[ 500, 500, 700, 400, 800, 1200]] ).T
    x2 = numpy.array( [[ 6, 8, 10, 13, 11, 13]] ).T

a0に関しては全て1にする。切片であるために1以外の数は適切ではない。そして\mathbf{X}は次のように定義される。

X = numpy.hstack( [a0, x1, x2] ) # 切片が不要であれば最初のa0をなくす。
# [[   1  500    6]
# [   1  500    8]
# [   1  700   10]
# [   1  400   13]
# [   1  800   11]
# [   1 1200   13]]

肝心のXは上記のように出力される。なお、切片が不要な場合は最初の「1」がなくなった行列となる。さて、定義できたのであとは最初に示した式通りに書いていく。

XTX = X.T.dot(X) # X^T.X                                                                                                                                                                                                      
invXTX = numpy.linalg.inv( XTX )
result = ( invXTX.dot( X.T ) ).dot( y )
print( result )
# [[1.14814815]
#  [0.00786008]
#  [0.53909465]]

出力値として、a_0=1.14814815, x_0=0.00786008, x_1=0.53909465となった。

この値を元にGの営業所の売上予想は 1.14814815 + 1300 \times 0.00786008 + 14 \times 0.53909465=18.9135772となる。

参考文献

  1. https://k-san.link/linear-regression-01/
  2. https://istat.co.jp/ta_commentary/multiple
  3. https://datachemeng.com/wp-content/uploads/ordinaryleastsquares.pdf

おすすめのサイト

数量化理論Ⅳ類(MDS)

数量化理論Ⅳ類はMDS,多次元尺度構成法とも呼ばれます。対象物の類似度を測る手法です。

あるアイテムが似ている、似ていないというのを平面上にプロットすることができます。平面上だけでなくても、直線、或いは3次元上でプロットできますが、多くは視覚的にわかりやすい2次元にプロットする例が多いです。なお、「似ている程度」「似 ていない程度」の量のことを数量化理論Ⅳ類では親近性と呼びます。

MDSの説明ですが、以下のドキュメントが非常にわかりやすく書かれています。導出が気になった方は、是非次のドキュメントに目を通してください。

http://www.page.sannet.ne.jp/yo-skmt/MultiEx/No13-15MU-Text02.pdf

女性ファッション雑誌をMDSで可視化

MDSの利用シーンとしての具体例を考えます。例えば女性専門雑誌がいくつかあるとします(VIVI, an anなど)。女性向けファッション雑誌の内容は男性からするとどれも一緒に見えますが、女性から見ればファッションの種類によって雑誌がかわります。では具体的にどの雑誌とどの雑誌が似ているのか、という事を男性が知りたいときにはどうすればいいのでしょうか。(男性は例え本の中身を見てもさっぱりですよね?)

まさにこんなときがMDSの出番と言えるでしょう。MDSを使えば次のようにマッピングできます。

上記画像は https://u-site.jp/research/methods/mds/ より。

雑誌の中身を見て、これとこれが似ているなど判断しなくてもMDSによりわかりやすく可視化することができます。

その他にも、例えばM1グランプリの採点から、採点者の傾向や、お笑い芸人のタイプなどを可視化出来ます。採点表はネットに公開されているので、興味があったらMDSで実験してみたら面白いかもしれません。

MDSの計算方法

対象物間に「距離」という概念を持ち込む必要があります。対象物間どうしの「距離」の概念を何らかの方法で定義できれば、 固有値分解を実施することでMDSが完了します。距離とは何か、ということを疑問に思った方は距離の公理を参照ください。大抵の場合は馴染み深い距離(ユークリッド距離)を利用するケースが多いです。

ちなみに、どのような行列に対しての固有値分解(あるいは特異値分解)をすることで分析手法の名前が変わってきます[2]

  1. 相関行列の固有値分解⇒主成分分析or因子分析 (数量化理論III類)
  2. 分散比行列の固有値分解 ⇒ 判別分析 (数量化理論II類)
  3. 頻度行列の特異値分解 ⇒コレスポンデンス分析(数量化理論III類)
  4. 距離行列の固有値分解 ⇒ MDS (数量化理論IV類)

MDSは対象間同士の距離行列さえ求めてしまえば、あとはやる(主成分分解する)だけで可視化できるというのはある意味興味深いです。

距離行列の作成

実際の例で計算をして見ます。今、A,B,Cの三人がいたとして、それぞれの互いにどう思っているかを評価してもらいます。(親近性の調査)

AはBはC は
AのことNaN33
B のこと 4 NaN 3
C のこと 32 NaN

上記のように、互いに評価をしています。なおこの行列を親近性行列と呼びます。

MDSとは、例えば3人を親密な人は距離が近く、そうでない人同士が遠くなるように位置を求めることです。2次元でも、3次元でも考えられますが、計算できる次元数は{対象者−1}次元となります。なお、親近性行列の要素は e_{11}, e_{12}, e_{13}, e_{21} ... e_{32}, e_{33}のように、 e_{ij}として表すことにします。

1次元の直線でA,B,Cの親近性を表現すること考えてみます。A, B, Cのそれぞれの位置をx1, x2, x3とします。ここで、位置x1, x2, x3の制約条件を次のとおりに定めます。

  • x1, x2, x3の平均は0
  • x1, x2, x3の分散は1

そして関数Qというものを定義します。関数Qは次のような定義です

m=特定の2人の親近性\times 特定の 2人の距離の2乗

つまり次のようになります。

m=e_{ij}\times (x_{i} - x_{j})^2

全員分のmの総和にマイナスしたものをQとして定義します。

Q=-\sum{ m_k }

Qの最大化を行うと、親近性の大きい個体の距離は短く、親近性の小さい個体の距離がなくなるようになります。なお、距離の定義は今回は線分で見るのでユークリッド距離を使いますが対象により変わります。要は距離の公理の「距離」が定義できるものであれば何でも良いです。

さてこれらを用いて計算ですが、細かい導出が坂元保秀教授により大変素晴らしくまとめられています。[1]

数式が多くなるので、結論を書きますと親近性行列から行列Hを定義します。行列Hとはh_{12} = e_{12} + e_{21}という形で、対角線同士を足し合わせる行列となります。対角線の要素[ latex ] e_{11}, e_{22}, e_{33} [/latex]については、e_{11} = -(e_{12}+e_{13}) \times 2 のように、自分意外の行を全て足しこみ、二倍してマイナスにした数字です。以上の操作を図にすると次のとおりです。

行列Hの作り方。親近性行列と、その天地を足す。加えて、親近性行列の行をすべて足しこんだものを2倍してマイナスをかけたものを対角の要素とする

行列Hがでたあとは、固有値分解をすることでx1,x2,x3が算出されます。いま、Pythonを用いて行列Hを固有値分解してみましょう。

import numpy as np
import math

a = [-12, 7, 6]
b = [7, -14, 5]
c = [ 6, 5, -5]
mat = np.array( [ a, b, c ] )
e_value, e_mat = np.linalg.eig( mat )
x1, x2, x3 = e_mat.T[ 0 ]
print x1, x2, x3

e_valueは固有値、e_matには固有ベクトルが入っています。e_matに格納されている値がまさにx1, x2, x3となります。それぞれの値は[-0.521333010475 -0.44553760869 -0.727810505167]となっています。

ここで先程の制約条件を思い出してください。平均が0、分散が1ということで定義しました。この条件が合っていることが確認できるかと思います。(厳密には要素数が3ですので、x1, x2, x3にsqrt(3)をそれぞれ掛けた値で分散が1になるのですが、スケーリングするということだけですので気にしなくて構いません。)

もし直線でなく平面でプロットしたい場合は固有ベクトルの1行目と2行目を使います。今一行目でx1, x2, x3としましたが、2行目をy1, y2, y3としてプロットすればよいだけです。3次元の際は、3行目をzとして利用します。

参考文献

[1] http://www.page.sannet.ne.jp/yo-skmt/MultiEx/No13-15MU-Text02.pdf

[2]https://www.nikkei-r.co.jp/glossary/id=1605

おすすめ記事

単純ベイズ分類器(Naive Bayes Classifier )

はじめに

ベイズ定理を利用した分類機を、ナイーブベイズ分類器とよびます。スパムフィルターとして使われており、聞いたことがある人も少なくないと思います。ベイズ分類器は実装が簡単な上に、出力値が確率であるため扱いやすく、ポピュラーな学習モデルです。 当サイトでは長い理論説明はせず、必要最低限の数式をで説明していこうと思います。

具体的なベイズ定理

男性のタバコを吸う人が肺癌になる確率 を求めます。タバコを吸ったら肺がんになるという確率は普通はわかり得ません。ただこの確率が、世の中の肺がん率、喫煙率、肺癌患者の喫煙率が解れば、ベイズ推定で求めることができるのです。

ベイズ式は以下のように表されます。

Pr( A | B ) = Pr( B | A ) * Pr( A ) / Pr( B )

それぞれは次のとおり定義します。

  1. Pr( A ) = 肺癌率
  2. Pr( B ) = 喫煙率
  3. Pr( A | B ) = 喫煙者の人が肺癌になる確率
  4. Pr( B | A ) = 肺癌の人が喫煙者である確率

なお、サイトに寄っては次のようにも書かれたりしています。

P(A|B) =\frac{P(A \cap B )}{P(B)} = \frac{P(B|A)P(A)}{P(B)}

罹患率等のデータの取得

実際に喫煙者の人が肺癌になる確率を求めます。喫煙率、肺がん率、肺癌患者の喫煙率を入手する必要があり、ネットで検索してみました。

男性の肺癌患者の30%はタバコによるガンだそうです[1]。そこで、式4. 喫煙者の人が肺癌になる確率は次のように定義できます。

4. Pr( B | A )=0.3 (喫煙者の人が肺癌になる確率)

続いて、肺癌率は男性の場合は7.4%[2], 男性の喫煙率は29.4%[3]なので, Pr(A), Pr(B)はそれぞれ次のように定義できます。

1. Pr(A) = 0.074(男性肺癌率)

2. Pr( B ) = 0.294 ( 男性喫煙率)

必要な確率は全て知り得たので、最終的なタバコを吸うと肺癌になる確率を求めます。最初の数式にそれぞれの確率を当てはめます。

Pr( A | B ) = 0.3 * 0.074 / 0.294 = 0.075

上記の通り、約7.5%の確率でタバコを吸うと肺癌になることが判明しました。通常、タバコを吸って肺癌になる確率は解るはずはありませんが、ベイズの定理を用いて求めることができるようになるのです。しかしながら、喫煙は8%弱肺癌に寄与するというのはなんとも面白い結果となりました。

ベイズ分類器

ベイズの定理を見たところで、実際にをどうやって分類器として力を発揮できるかを考えます。

ベイズ分類器では、左側のベイズ理論を右側のように仮定していきます。

P(Y|X) = \frac{P(X|Y) \times P(Y)}{P(X)}

上記式は次のようにベイズ分類機では意味すると当てはめています。

Posterior = \frac{Likelifood \times prior }{ evidence }

なお、Likefoodは尤度、Priorは事前分布、Evidenceは周辺尤度と呼ばれています。

尤度であるP(X|Y) 、それから事前分布であるP(Y)は学習データから入手します。今次のような学習データがあったとします。次の学習データセットはエジンバラ大学の教材[4]を参考にしました。以下の表は頻度表といわれます。

天気気温湿度風速ゴルフ
晴れ暑い高いなしNO
晴れ暑い高いありNO
曇り暑い高いなしYES
雨天普通高いなしYES
雨天寒い普通なしYES
雨天寒い普通ありNO
曇り寒い普通ありYES
晴れ普通高いなしNO
晴れ寒い普通なしYES
雨天普通普通なしYES
晴れ普通普通ありYES
曇り普通高いありYES
曇り暑い普通なしYES
雨天普通高いありNO

ここで、X=( 天気、気温、湿度、風速), そしてY=( ゴルフ )となります。YのYESはゴルフした、Noはゴルフをしなかったという意味です

P( X | Y ) と P( Y )を上記の表から当てはめます。

P(Y=YES)=\frac{ YES 数}{ YES 数+NO数} = \frac{9}{14}

P(X=晴れ|Y=Play) = \frac{ YES 数 \cap 晴れ数}{ YES 数} = \frac{2}{9}

P(X|Y)のXは1つだけで説明しましたが、Xが複数条件の時も知りたくなるかもしれません例えば外が晴れ&寒い時にはゴルフするのかどうかなどです。一件わかりますが、寒い&晴れているなどのXの条件数が多くなるほど計算は複雑になります。状態数が多くなるためです。

そこでXの要素はすべて独立しているという仮定を置きます。そうすることで次のように書き換えることが可能となります。

P(X=晴れ,X=寒い| Y ) = P(X=晴れ|Y) P(X=寒い|Y)

独立という仮定をおくことで、分割することができるようになるのです。こうすることでパラメータ数は非常にが少なくなり(具体的には指数関数が線形になり)高次元のデータでも扱うことができるようになります。

データが独立という仮定を置かれている所以がナイーブ(雑)と言われている理由です。

迷惑メールフィルター実装例

実際にベイズ推定を用いた迷惑メールフィルターがどのように作られるのか説明します。

迷惑メール 30通中を受信したとします。迷惑メールの中で次の単語「バイアグラ」、「限定品」、「機械学習」という単語が含まれた表、頻度表を作成します。次のようになりました。

バイアグラ限定品機械学習
回数7170

迷惑メールが30通ですので、ここから尤度表を作成することができます。単に母数で除算するだけです。

バイアグラ限定品機械学習
YES7/3017/300/30
NO23/3013/3030/30

尤度表から、例えば「バイアグラ」が入り、「限定品」が入っていない場合の迷惑メールの確率 P( 迷惑メール | バイアグラ ∧~限定品)を求めることができます。

\frac{P( バイアグラ \cap\lnot限定品 | 迷惑メール) P(迷惑メール)}{ P( バイアグラ\cap \lnot 限定品 )}

さて、ここで∧の条件は独立である場合分解できるということを説明しました。そのため、上記式は次のように分解できます。

\frac{P( バイアグラ| 迷惑メール) P( \lnot 限定品 | 迷惑メール) P(迷惑メール)}{ P( バイアグラ ) P( \lnot 限定品 ) }

続いて迷惑メールでない 正常なメール100通のうち、先程の単語についての尤度表を作成してみます。

バイアグラ限定品機械学習
YES1/1001/1007/100
NO99/10029/10023/100

ゆう土俵から、正常なメールの確率は次のとおりとなります。

\frac{P( バイアグラ| 正常メール) P( \lnot 限定品 | 正常メール) P( 正常 メール)}{ P( バイアグラ ) P( \lnot 限定品 ) }

これで準備が整いました。メールに「バイアグラ」& NOT「限定品」 の条件 で迷惑メールと正常メールのどちらに分類されるかの確率を求めます。今、分母が共通していますので、分子の部分だけに着目して計算をします。(何故分母は無視したのかは後述します)

迷惑メールの尤度= P( バイアグラ| 迷惑メール) P( \lnot 限定品 | 迷惑メール) P(迷惑メール) = 7/30 * 13/30 * 30/130= 0.023

正常メールの尤度= P( バイアグラ| 正常 メール) P( \lnot 限定品 | 正常 メール) P( 正常 メール) = 1/100 * 29/100 * 100/130 = 0.0022

実はこれだけで、この種のメッセージが迷惑メールである確率をもとめることができます。次式で求まります。

0.023/(0.023 + 0.0022 )=0.912

よく見るとわかりますが、算出された合計の尤度で、該当尤度を除算するだけです。結果として、バイアグラを含み、限定品を含まないメールは、おおよそ91.2%が迷惑メールであるということがわかりました。このように尤度をすべて足しこんだ値を、該当の尤度で割ることで、それぞれのメールの確率がわかります。

分母を省略した理由ですが、分母が共通しています。確率を計算するときに打ち消し合って消えるので無視できるのです。

参考文献

  1. https://www.sankei.com/life/news/180104/lif1801040007-n1.html
  2. https://www.haigan.gr.jp/guideline/2018/1/1/180101000000.html
  3. https://ganjoho.jp/reg_stat/statistics/stat/smoking.html
  4. http://www.inf.ed.ac.uk/teaching/courses/inf2b/learnSlides/inf2b12-learnlec06.pdf

次に読むのおすすめ記事

正準相関分析(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乗の候補が必要になり、メモリを圧迫するうえ、データのスキャンにも時間が必要となります。

参考文献

おすすめの記事