識別モデルと生成モデル

はじめに

分類問題に適用される機械学習の手法は、以下の3つに大別できる（下図参照）。

識別関数を用いるもの
識別モデルを用いるもの
生成モデルを用いるもの

識別関数は1つの点を推定するが、識別モデルと生成モデルは確率を推定する。各手法の代表的なアルゴリズムも上図に示した。

今回の記事では、2分類問題を取り上げ、識別モデルと生成モデルの2つのアプローチでこれを解き、それぞれの手法の特長を解説する。

ソースコード

今回のソースコードはここにあるsample.ipynbである。

識別モデル

いま、観測データ $\{X,Y\}$ が与えられているとする。ここで

(1) $\begin{eqnarray*} X&=&(x_1,\cdots,x_N), \;\;x_n \in \mathcal{R}^1\\ Y&=&(y_1,\cdots,y_N), \;\;\;y_n \in \{0,1\} \end{eqnarray*}$

である。パラメータ $\theta=(\alpha,\beta)$ を導入し、同時確率分布 $p(X,Y,\theta)$ を考える。ベイズの定理から次式を得る。

(2) $\begin{equation*} p(\theta|X,Y)=\frac{p(Y|X,\theta)p(\theta)}{p(Y|X)} \end{equation*}$

最初に尤度 $p(Y|X,\theta)$ をBernoulli分布を用いてモデル化する。

(3) $\begin{eqnarray*} p(Y|X,\theta)&=&\prod_{n=1}^N p(y_n|x_n,\alpha,\beta) \\ p(y_n|x_n,\alpha,\beta) &=& {\rm Bern}(y_n|\nu_n) \\ &=& \nu_n^{y_n}(1-\nu_n)^{(1-y_n)} \\ \nu_n&=&f_s(z) \\ z&=&\alpha+\beta x_n \\ f_s(z)&=&\frac{1}{1+\exp{(-z)}} \end{eqnarray*}$

次に事前分布 $p(\theta)$ をモデル化する。

(4) $\begin{eqnarray*} p(\theta)&=&p(\alpha)p(\beta) \\ p(\alpha)&=&\mathcal{N}(\alpha|0,\sigma^2_0) \\ p(\beta)&=&\mathcal{N}(\beta|0,\sigma^2_0) \end{eqnarray*}$

ここで、 $\mathcal{N}(z|\mu,\sigma^2)$ は平均 $\mu$ 、標準偏差 $\sigma$ の正規分布を表す。先に与えた関数 $f_s(z)$ はシグモイド関数である（下図参照）。

この関数の値は $y_n$ が1になる確率に相当するので、0と1を分離する境界線は $f_s(z)=0.5$ となる $z$ の値、 $z=0$ から求めることができる。すなわち

(5) $\begin{equation*} \alpha + \beta x_n=0 \end{equation*}$

故に

(6) $\begin{equation*} x_n=-\frac{\alpha}{\beta}\equiv b_d \end{equation*}$

が境界線となる。以上の準備のあと事後確率、 $p(\alpha|X,Y)$ と $p(\beta|X,Y)$ を求めることになる。 $\sigma_0$ はあらかじめ与える定数である。

先に見た尤度 $p(Y|X,\theta)$ は、 $X$ を与えて、それが属するクラス $Y$ を識別するモデルとなっている。これが識別モデルと呼ばれる所以である。

生成モデル

同時確率分布 $p(X,Y,\theta)$ に戻る。この分布はベイズの定理より以下のように変形することもできる。

(7) $\begin{equation*} p(\theta|X,Y)=\frac{p(X|Y,\theta)p(\theta)}{p(X|Y)} \end{equation*}$

最初に尤度 $p(X|Y,\theta)$ をモデル化する。

(8) $\begin{equation*} p(X|Y,\theta)=\prod_{n=1}^N p(x_n|y_n,\theta) \end{equation*}$

ここで、 $y_n$ が0のときの分布と1のときの分布に分けて考える。

(9) $\begin{eqnarray*} p(x_n|y_n=0,\theta)&=&\mathcal{N}(x_n|\mu_0,\sigma_c^2) \\ p(x_n|y_n=1,\theta)&=&\mathcal{N}(x_n|\mu_1,\sigma_c^2) \end{eqnarray*}$

どちらも正規分布とする。その標準偏差は共通の値とし、平均だけを異なるものとする。次に、事前分布 $p(\theta)$ をモデル化する。

(10) $\begin{eqnarray*} p(\theta)&=&p(\mu_i)p(\sigma_c),\;\;i=0,1 \\ p(\mu_i)&=&\mathcal{N}(\mu_i|0,\sigma^2_3) \\ p(\sigma_c)&=&\mathcal{N}(\sigma_c|0,\sigma^2_4),\;\;\sigma_c \geq 0 \end{eqnarray*}$

$\sigma_3,\sigma_4$ はあらかじめ与える定数である。2つの $x_n$ の分布を正規分布とし、その標準偏差を同じものとしたので、2つの分布の境界線は次式で与えられる。

(11) $\begin{equation*} b_d=\frac{\mu_0+\mu_1}{2} \end{equation*}$

以上の準備のあと事後確率、 $p(\mu_0|X,Y)$ 、 $p(\mu_1|X,Y)$ 、 $p(\sigma_c|X,Y)$ を求めることになる。

ここでは尤度として2つの $X$ の分布を考えた。すなわち、各クラス $Y$ に属するサンプル $X$ を生成するモデルとなっている。これが生成モデルと呼ばれる所以である。

PyMC3による実装

ここまでの計算を、データセット「iris」を用いてPyMC3により行う。このデータセットにはアヤメの３品種

setosa
versicolor
virgnica

が50個ずつ集められており、4つの特徴量

がく片の長さ：sepal length
がく片の幅：sepal width
花びらの長さ：petal length
花びらの幅：petal width

の値が格納されている。データの先頭の様子は以下の通りである。

sepal_lengthについてのバイオリン図を次に示す。縦軸はsepal_lengthの値、横軸は3つの品種を表す。

バイオリン図とは、ヒストグラムを縦軸に沿って描画し、それを左右に展開したものである。今回は4つの特徴量の中のsepal lengthを $x$ 、setosaとversicolorの2品種を $y$ として用いる。

最初に識別モデルを考える。コードは以下の通り。

import pymc3 as pm

with pm.Model() as model_0:
    alpha = pm.Normal('alpha', mu=0, sd=10)
    beta = pm.Normal('beta', mu=0, sd=10)
    z = alpha + pm.math.dot(x, beta)
    fs = pm.Deterministic('fs', 1 / (1 + pm.math.exp(-z)))
    bd = pm.Deterministic('bd', -alpha / beta)
    yl = pm.Bernoulli('yl', p=fs, observed=y)
    trace = pm.sample(5000)

import pymc3 as pm

with pm.Model() as model_0:

alpha = pm.Normal('alpha', mu=0, sd=10)

beta = pm.Normal('beta', mu=0, sd=10)

z = alpha + pm.math.dot(x, beta)

fs = pm.Deterministic('fs', 1 / (1 + pm.math.exp(-z)))

bd = pm.Deterministic('bd', -alpha / beta)

yl = pm.Bernoulli('yl', p=fs, observed=y)

trace = pm.sample(5000)

4行目： $p(\alpha)$ を定義する。 $\sigma_0=10$ とした。
5行目： $p(\beta)$ を定義する。
6行目： $z$ を定義する。
7行目： $f_s(z)$ を定義する。
8行目： $b_d$ を定義する。
9行目：尤度 $p(Y|X,\alpha,\beta)$ を定義する。
10行目：MCMCを行う。

収束具合を見るため、次のコードでGelman-Rubinテストを行う。

pm.gelman_rubin(trace)

1	pm.gelman_rubin(trace)

計算される値は全て1.1未満となることを確認できる（結果は略）。1.1未満であれば収束したとみなして良い。次に、事後確率から算出される要約統計量を次のコードで計算する。

varnames = ['alpha', 'beta', 'bd']
pm.summary(trace, varnames)

1 2	varnames = ['alpha', 'beta', 'bd'] pm.summary(trace, varnames)

出力は以下の通り。

「mean」欄の「bd」の値5.42が、0と1を分ける境界線である。下図は、境界線の平均値とシグモイド関数の平均値を描画したものである。

黒丸はsepal lengthの観測値（下側がsetosa）、濃い赤のラインが境界線 $b_d$ の平均値、薄い赤の領域が $b_d$ の95%HPDである。また、濃い青のラインがシグモイド関数 $f_s$ の平均値、薄い青の領域が $f_s$ の95%HPDである。ある量が95%HPDの領域内にあるとは、95%の確率でその領域内に存在することを意味する。HPDはHighest Posterior Density（最高事後密度）の略である。

次に、生成モデルの場合を考える。コードは以下の通り。

with pm.Model() as model:
    mus = pm.Normal('mus', mu=0, sd=10, shape=2)
    sigma = pm.HalfNormal('sigma', sd=5)
    setosa = pm.Normal('setosa', mu=mus[0], sd=sigma, observed=x[:50])
    versicolor = pm.Normal('versicolor', mu=mus[1], sd=sigma, observed=x[50:])
    bd = pm.Deterministic('bd', (mus[0] + mus[1]) / 2)
    trace = pm.sample(5000)

with pm.Model() as model:

mus = pm.Normal('mus', mu=0, sd=10, shape=2)

sigma = pm.HalfNormal('sigma', sd=5)

setosa = pm.Normal('setosa', mu=mus[0], sd=sigma, observed=x[:50])

versicolor = pm.Normal('versicolor', mu=mus[1], sd=sigma, observed=x[50:])

bd = pm.Deterministic('bd', (mus[0] + mus[1]) / 2)

trace = pm.sample(5000)

2行目： $p(\mu_i)$ を定義する。 $\sigma_3=10$ とした。
3行目： $p(\sigma_c)$ を定義する。 $\sigma_4=5$ とした。
4,5行目：2つの分布 $p(x_n|y_n=0,\theta)$ と $p(x_n|y_n=1,\theta)$ を定義する。
6行目： $b_d$ を定義する。
7行目：MCMCを行う。

先と同様に、Gelman-Rubinテストは合格する（詳細は略）。要約統計量は以下の通り。

境界線の平均値は5.47程度になることが分かる。これは先に求めた識別モデルでの結果、5.42に近い値である。境界線の平均値と95%HPDを図示したものが下図である。

識別モデルでの結果に比べ、95%HPDが少し狭くなっている。すなわち、生成モデルの方が確からしさの高い値になる。

識別モデルと生成モデルの比較

識別モデルと生成モデルの違いは、それぞれの尤度を比較すると明確になる。すなわち、前者の尤度は $p(Y|X,\theta)$ であり、これは $X$ が観測されたときクラス $Y$ が実現する確率を表している。識別を行うと言う目的に直接アプローチする手法である。一方、後者では、その尤度 $p(X|Y,\theta)$ を見ても分かる通り、識別に直接アプローチせず、 $X$ の分布を最初に求めるている。この分布を通して未知パラメータの事後分布を求めることで本来の目的である識別を行っている。生成モデルの場合、 $X$ の分布が得られるので、擬似データの生成や、外れ値検知などにも応用することができる。
今回の例では、生成モデルの方が確度の高い境界線を得ることができた。その理由は、setosaとversicolorのサンプルが正規分布で良く近似できたためであると思われる。サンプルの分布が正規分布に従わない場合は、識別モデルの方が良いアプローチとなる。
今回生成モデルで用いた境界線 $b_d$ はとても簡単な式 $(\mu_0+\mu_1)/2$ で表すことできた。その理由は、2つの正規分布の標準偏差を同じものとしたためである。サンプルの分布を正規分布とし、その標準偏差を等しくする解析手法を線形判別分析と呼ぶ。

Kumada Seiya

仕事であろうとなかろうと勉強し続ける、その結果”中身”を知ったエンジニアになれる

AnacondaのNumPyとpipによるNumPyの速度の違い

Amazon SageMakerのAPIを使った賢いWebアプリの作り方

識別モデルと生成モデル

はじめに

ソースコード

識別モデル

生成モデル

PyMC3による実装

識別モデルと生成モデルの比較

Kumada Seiya

最近の記事

レイアウト解析

ChatGPTを用いた画像異常検知

クォータニオンと回転行列の関係

OpenCVで深層学習

OpenCVのcv::Matに関する小ネタ集

テンソルネットワークの入り口（その2）

線形回帰をPythonで数式から逃げずに実装してみた

Dropout

PandasをDaskで簡単に並列処理して高速化しよう

Unsupervised Representation Learning by Predicting Image Rotationsの紹介

主成分分析とは何なのか、とにかく全力でわかりやすく解説する

k近傍法による時系列データの異常検知

Google Vision APIでOCR

アーカイブ

カテゴリー