Courseraの機械学習コースの9週目です。今回は異常検知とレコメンデーションについてです。プログラム演習が今回までなので、実装系の細かい内容に関しては今回が最後になるっぽいです。
記事「Coursera機械学習入門オンライン講座虎の巻(文系社会人にオススメ)」を最初に読むことをお勧めします。

サマリ

実績時間：約9時間
目安時間：約6時間(正確には6時間3分)
※コースに各章の時間目安が書かれていて、それを足しただけです。

一覧

No.	章	目安時間	内容
1	Anomaly Detection	108分	異常検知
1.1	Density Estimation	29分	確率密度による推測
1.2	Building an Anomaly Detection System	32分	異常検知システムの構築
1.3	Multivariate Gaussian Distribution	27分	多変量ガウス分布
1.4	Review	20分	異常検知のクイズ
2	Recommender Systems	255分	レコメンデーション
2.1	Predicting Movie Ratings	21分	映画評価値の予測
2.2	Collaborative Filtering	18分	協調フィルタリング
2.3	Low Rank Matrix Factorization	16分	低ランク行列分解
2.4	Review	200分	レコメンデーションのクイズと異常検知も含んだプログラム演習

1. Anomaly Detection

1.1. Density Estimation

1.1.1. Problem Motivation

Video:7分
異常検知(Anomaly Detection)の使い所についてです。
異常検知は教師なし学習でありつつ、教師あり学習とも似ている点があります。
飛行機のエンジンを例として、以下の2つの特徴量を持ちます。

$x_1$：エンジンから生成される熱量
$x_2$：振動

2つの特徴量から正常か異常かを学習させ、新しいエンジンの特徴量$x_{test}$に対して異常かを判断させてい場合が異常検知のユースケースです。

正常である確率を$p(x_{test})$の関数を使って計算し、しきい値$\epsilon$によって正常・異常を識別します。

$p(x_{test}) < \epsilon \quad :異常$
$p(x_{test}) \geq \epsilon \quad :正常$

1.1.2. Gaussian Distribution

Video:10分
ガウス分布(Gaussian Distribution)について。
ガウス分布(正規分布)は以下の式です。

p(x) = \frac{1}{\sqrt{2\pi}\sigma}\exp\left(-\frac{(x-\mu)^2}{2\sigma^2}\right)

平均($\mu$)と分散($\sigma^2$)の値をパラメータにとり、分布が変わります。

統計学で学ぶガウス分布/正規分布と基本同じです。私は「スバラシク実力がつくと評判の統計学キャンパス・ゼミ」や「まずはこの一冊から意味がわかる統計学」で勉強しました。
「基本同じ」と言ったのは、分散については不偏分散を使わないことが多いそうです。以下の式の$n-1$の部分に$n$を使います。これは、データ数が多ければどちらでも大差ないからだそうです。

{S^2 = \frac{1}{n-1} \sum_{i=1}^n (X_i - \overline{X})^2}

以下の記事参照

1.1.3. Algorithm

Video:12分
異常検知のアルゴリズム。
以下の計算をして異常かどうかを判断します。

ステップ1 平均値算出

平均値$\mu_j$を求める

{\mu_j = \frac{1}{m} \sum^m_{i=1}x_j^{(i)}}

ステップ2 分散算出

分散$\sigma_j^2$を求める

{\sigma_j^2 = \frac{1}{m} \sum^m_{i=1} (x_j^{(i)} - \mu_j)^2}

ステップ3 確率計算としきい値による異常判断

ガウス分布による確率の積$p(x_{test})$を求め、しきい値$\epsilon$と比較し異常判断

{p(x_{test}) = \Pi^{n}_{j=1} p(x_j; \mu_j, \sigma^2_j) = \Pi^n_{j=1} \frac{1}{\sqrt{2\pi} \sigma_j} exp(-\frac{(x_j - \mu_j)^2}{2\sigma_j^2})} \\
if \quad p(x) < \epsilon \quad then \quad 異常

1.2. Building an Anomaly Detection System

1.2.1. Developing and Evaluating an Anomaly Detection System

Video:13分
異常検知システムの開発と評価。
以下の内容を学習。

訓練：交差検証：テストを6:2:2の割合で分割
正解ラベル$y$で、正常の場合$y=0$、異常の場合$y=1$とします
訓練データは正解ラベルを含めない(教師あり学習ではないし、そもそもほとんどが正常データのはずだから)
交差検証、テストデータは正解ラベルを含めてF値等で評価(不均衡データなので正答率では評価しない)
交差検証で評価時にしきい値$\epsilon$の値を変えながらF値を見て最適なしきい値を決めます

※不均衡データに対して正答率を使わないことは6週目に学習しました。

1.2.2. Anomaly Detection vs. Supervised Learning

Video:7分
異常検知と教師あり学習のガイドライン。
異常検知でも評価時に正解ラベル$y$を使いましたが、だったら教師あり学習でもいいのではないか、という疑問が出ます。
下のスライドがどちらを使うべきかのガイドラインです。
大雑把には異常のパターンを学習できれば教師あり学習です。「異常のパターンを学習」するためには、一定量のデータ数があり、パターン化できる異常である必要があります。

以下はそれぞれの具体例です。

1.2.3. Choosing What Features to Use

Video:12分
どの特徴量を使うかの選択。
特徴量がガウス分布(正規分布)に従っていない場合があります。その場合は$log$関数やルートでガウス分布となるように変換をします。実務的にはヒストグラムを確認しながら探索的に変換をします。

異常のデータに対する確率密度$p(x)$が大きかった場合にどうするかです。その場合は、新しい特徴量を作り、再度モデルを作って確認するのが定石のようです。

1.3. Multivariate Gaussian Distribution (Optional)

1.3.1. Multivariate Gaussian Distribution

Video:13分
多変量ガウス分布(Multivariate Gaussian Distribution)について。
前の章で習った通常のガウス分布では対応できない場合があります。それは下図に左のような分布を持ったデータです。特徴量$x_1$と$x_2$が関連していて中心から同心円状に分布しないような場合です。

このような場合に多変量ガウス分布を使います。モデル$p(x_1),p(x_2)$と独立して確率を定義せずに$p(x)$とひとつにまとめます。そして、パラメータに$\mu$と$\Sigma$(共分散行列)を使います。これにより、分布を下図のように変化させます。

1.3.2. Anomaly Detection using the Multivariate Gaussian Distribution

Video:14分
多変量ガウス分布(Multivariate Gaussian Distribution)を使った異常検知。
以下の計算をして異常かどうかを判断します。

ステップ1 平均値算出

平均値$\mu_j$を求めます(オリジナルと同じ)

{\mu_j = \frac{1}{m} \sum^m_{i=1}x_j^{(i)}}

ステップ2 分散算出

共分散行列$\Sigma$を求める

{\Sigma = \frac{1}{m} \sum^m_{i=1} (x^{(i)} - \mu) (x^{(i)} - \mu)^T}

ステップ3 確率計算としきい値による異常判断

ガウス分布による確率の積$p(x)$を求め、しきい値$\epsilon$と比較し異常判断

{p(x) = \frac{1}{(2\pi)^{\frac{n}{2}} | \Sigma |^{\frac{1}{2}}} exp(-\frac{1}{2}(x - \mu)^T \Sigma^{-1}(x - \mu))} \\
if \quad p(x) < \epsilon \quad then \quad 異常

オリジナルのガウス分布と多変量ガウス分布をどのように使い分けるかのガイドラインです。

オリジナルの特徴
- 新たな特徴量をマニュアルで作らなければいけない場合があります
- 計算コストが小さい
多変量ガウス分布の特徴
- 特徴量間の関係性を自動で検知してくれます
- 計算コストが高い
- データ数が特徴量の種類数の10倍程度必要

1.4. Review

1.4.1. Lecture Slides

Reading:10分
異常検知のスライドPDF。

1.4.2. Quiz: Anomaly Detection

Quiz:5問, Reading:10分
簡単な内容です。
でも勘違いをして1問ミスしてしまい、80%でした。

2. Recommender Systems

2.1. Predicting Movie Ratings

2.1.1. Problem Formulation

Video:7分
レコメンダシステム(Recommender System)とはなにか。
Amazonなどで出てくるおすすめ商品を表示するやつです。
アカデミアの世界ではたいして話題にならない領域らしいですが、ビジネスの世界では非常に重視されるようです。
例えば下表のような人ごとの5段階での映画評価のマトリックスがあったとします。

	山田	田中	佐藤
映画1	3	?	5
映画2	4	2	?

以下の表記をします。

表記	意味	表での例
$n_u$	ユーザ数	3
$n_m$	評価対象映画数	2
$r^{(i,j)}$	ユーザ$j$が映画$i$を評価済の場合1	$r^{(1,3)}=1, r^{(2,1)}=0$
$y^{(i,j)}$	ユーザ$j$による映画$i$の評価	$y^{(1,3)}=5, y^{(2,1)}=undefined$

2.1.2. Content Based Recommendations

Video:14分
コンテントベースのレコメンデーション(Content Based Recommendations)についてです。
映画などの内容に着目して、レコメンデーションをする手法です。
例えば下表の2種類の映画があったとして、特徴量($x^{(i)}$)を設定します。下表の例だと映画1はアクションムービーで、映画2はラブロマンスものという意味です。

	特徴1(アクション)	特徴2(ロマンス)
映画1	1	0
映画2	0.1	0.9

そして、映画の特徴量を($x^{(i)}$)、学習すべきパラメータを$\theta^{(j)}$とし、予測するユーザ評価は以下の式で表します。
※$\theta^{(j)}$の$j$はユーザ$j$を表しています。例えば1人目が田中太郎さんだとしたら$\theta^1$は田中太郎さんの評価を予測するパラメータです。

{y^{(i,j)} = (\theta^{(j)})^{T}} x^{(i)}

以下の式で最適化をします。右の$\lambda$が入っている部分は正則化項です。

でもって、最急降下法でパラメータを最適化します。

コンテントベースのレコメンデーション手法の弱点は、すべての評価対象が特徴を持つことができない点です。全映画にロマンスやアクションなどの特徴を定義しておくことは現実的ではありません。

2.2. Collaborative Filtering

2.2.1. Collaborative Filtering

Video:10分
協調フィルタリングです。
コンテントベースの場合は、映画の特徴量$x^{(i)}$からユーザごとの評価に対するパラメータ$\theta^{(i)}$を学習しました。協調フィルタリングでは、ユーザごとの評価に対するパラメータ$\theta^{(i)}$から映画の特徴量$x^{(i)}$の学習が追加され$\theta \to x \to \theta \to x \to ...$と交互に学習していきます。
目的関数は以下の式です。

2.2.2. Collaborative Filtering Algorithm

Video:8分
協調フィルタリングのアルゴリズム。
映画の特徴量$x^{(i)}$とユーザごとの評価に対するパラメータ$\theta^{(i)}$を同時に学習する方法です。以下の手順です。

Step 1. パラメータ初期化

パラメータ$x^{(1)}, ..., x^{(n_m)}, \theta^{(1)}, ..., \theta^{(n_u)}$を乱数で初期化します。5週目に習ったSymmetry Breakingです。

Step 2. 学習

以下の目的関数で両パラメータを同時に学習します。

{J = \frac{1}{2} \sum_{(i,j):r(i,j)=1} ((\theta^{(j)})^{\mathrm{T}}x^{(i)} - y^{(i,j)})^2 + \frac{\lambda}{2}\sum^{n_m}_{i=1} \sum^{n}_{k=1}(x_k^{(i)})^2
+ \frac{\lambda}{2}\sum^{n_u}_{i=1} \sum^{n}_{k=1}(\theta_k^{(j)})^2}

最急降下法によるパラメータ更新。

{x_k^{(i)} := x_k^{(i)} - \alpha (\sum_{j:r(i,j)=1} ((\theta^{(j)})^{\mathrm{T}} x^{(i)} - y^{(i,j)})\theta^{(j)}_k + \lambda x^{(i)}_k)
\\
\theta_k^{(i)} := \theta_k^{(j)} - \alpha (\sum_{i:r(i,j)=1} ((\theta^{(j)})^{\mathrm{T}} x^{(i)} - y^{(i,j)})x^{(i)}_k + \lambda \theta^{(j)}_k)
}

Step 3. 評価値算出

学習したパラメータを使って各ユーザの評価値を算出します。

{y^{(i,j)} = (\theta^{(j)})^{T}} x^{(i)}

2.3. Low Rank Matrix Factorization

2.3.1. Vectorization: Low Rank Matrix Factorization

Video:8分
低ランク行列分解。

映画評価表に対して$Y$という行列を定義します。

そして、評価値の予測するには、以下の行列で表しました。

\displaystyle \begin{bmatrix} (x^{(1)})^T(\theta^{(1)}) & \ldots & (x^{(1)})^T(\theta^{(n_u)})\\ \vdots & \ddots & \vdots \\ (x^{(n_m)})^T(\theta^{(1)}) & \ldots & (x^{(n_m)})^T(\theta^{(n_u)})\end{bmatrix}

ごちゃごちゃしているので、$X$と$\Theta$を使ってすっきりと表します。

X = \begin{bmatrix} - & (x^{(1)})^T & - \\ & \vdots & \\ - & (x^{(n_m)} & - \end{bmatrix},\ \Theta = \begin{bmatrix} - & (\theta^{(1)})^T & - \\ & \vdots & \\ - & (\theta^{(n_u)} & - \end{bmatrix}\\
Y = X\Theta^T

こうしてランクが低くなるので、低ランク行列分解と言います。

映画を$x^{(i)} \in \mathbb{R}^n$の特徴量で表すと、映画同士の類似性を計算できます。$||x^{(i)}-x^{(j)}||$で計算し、値(距離)が小さければ類似しているということになります。

2.3.2. Implementational Detail: Mean Normalization

Video:8分
平均標準化(mean Normalization)を使って実装を良くする話。
もし、Eveが新ユーザ登録して、映画未評価だったと仮定します。平均標準化をせずに、?の値に0を使うと評価予測値が0になります。

しかし、平均標準化をすることで未評価に0を使っても平均値が出力されます(平均を0とするので)。

2週目に学習したFeature Scalingと考え方は同じです。ただ、評価の場合は値の最大値と最小値が決まっているのが少し違います。

2.4. Review

2.4.1. Lecture Slides

Reading:10分
レコメンデーションのスライドPDF。

2.4.2. Quiz: Recommender Systems

Reading:10分
数式系が難しく時間がかかってしまいました。3回目で80%のギリギリ合格。
あとで理解を深めてからもう一度挑戦したい・・・

2.4.3. Anomaly Detection and Recommender Systems

Programming Assignment:180分
異常検知とレコメンデーションのプログラム演習です。
GitHubを見てやったので1時間程度で終了。忙しかったこともあり、そんなに真面目にやっていないです。

感想

今回の内容は、仕事の忙しさなどにより2週間以上かけてしまいました。そのため、思い出すのに時間がかかったりと非効率的なコストをかけてしまっています。
エッセンスだけ理解して、実際に使う時になったらきちんと理解しようと考え、短時間で済ませています。なので数式などは、あまり理解していません(理解する努力も大してしていません)。

関連記事

私より詳しくまとめてくれています。

他の週のリンクです。