はじめに
機械学習で登場する確率分布について色々な角度から理解したいシリーズです。
この記事では、1次元混合ガウス分布(一変量混合正規分布)の定義を確認します。
【前の内容】
【他の記事一覧】
【この記事の内容】
1次元混合ガウス分布の定義式
1次元混合ガウス分布(Gaussian Mixture Distribution)・一変量混合正規分布(Normal Mixture Distribution)の定義を確認します。
定義式
混合ガウス分布は、複数のガウス分布を線形結合した分布です。ガウス分布については「1次元ガウス分布の定義式の確認 - からっぽのしょこ」を参照してください。
1次元混合ガウス分布は、パラメータ$\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol{\sigma}$を用いて、次の式で定義されます。
ここで、$K$はクラスタ数、$\pi_k$はクラスタ$k$となる確率、$\mu_k$はクラスタ$k$の平均パラメータ、$\sigma_k$は標準偏差パラメータ、$\sigma_k^2$は分散パラメータです。
各クラスとなる確率をまとめたベクトル$\boldsymbol{\pi} = (\pi_1, \pi_2, \cdots, \pi_K)$を混合比率(混合パラメータ)と言います。各成分は0から1の値で、全ての成分の和が1になる必要があります。この条件を数式で表すと、$0 \leq \pi_k \leq 1$または$\pi_k \in (0, 1)$、$\sum_{k=1}^K \pi_k = 1$です。
全てのクラスタの平均パラメータをまとめて$\boldsymbol{\mu} = \{\mu_1, \mu_2, \cdots, \mu_K\}$、標準偏差パラメータをまとめて$\boldsymbol{\sigma} = \{\sigma_1, \sigma_2, \cdots, \sigma_K\}$で表します。$\mu_k$は実数、$\sigma_k, \sigma_k^2$は正の実数を満たす必要があります。
$x$は実数をとります。
混合ガウス分布は、混合比率$\boldsymbol{\pi}$を用いて$K$個のガウス分布の加重平均で定義されています。
1次元混合ガウス分布の確率変数$x$は、割り当てられたクラスタを潜在変数$\mathbf{s}$で表すことで、次の式でも定義できます。
1 of K表現(one-hotベクトル)の潜在変数$\mathbf{s} = (s_1, s_2, \cdots, s_K)$を導入します。割り当てられたクラスタに対応する成分が1で、それ以外の成分が0である$K$次元ベクトル$\mathbf{s}$によってクラスタを表します。各成分が0か1の値をとり、(値が1なのは1つだけなので)全ての成分の和が1になる条件を数式で表すと、$s_k \in \{0, 1\}$、$\sum_{k=1}^K s_k = 1$です。
$\mathbf{s}$は、パラメータ$\boldsymbol{\pi}$を持つカテゴリ分布の確率変数として定義できます。確率$\pi_k$で$s_k = 1$(クラスタが$k$)になります。カテゴリ分布については「カテゴリ分布の定義式 - からっぽのしょこ」を参照してください。
0乗は1なので($a^0 = 1$なので)、$s_k\ (k = 1, \dots, K)$を指数として用いて総乗$\prod_k$の形にすることで、割り当てられたクラスタの分布のみを取り出すことができます。
クラスタが$k$のとき($s_k = 1$のとき)
となります。
つまり、クラスタが$k$のとき($s_k = 1$のとき)、$x$はクラスタ$k$の分布に従います。
混合分布は複数の分布を用いて定義されていますが、確率変数は1つのクラスタに従うように設計されています。
ここまでは、混合ガウス分布の定義を数式から確認しました。次は、グラフから確認します。
グラフの確認
「K個のガウス分布」と「クラスタ数がKの混合ガウス分布」の関係をグラフで確認します。混合ガウス分布のグラフ作成については「分布の作図」を参照してください。
利用するパッケージを読み込みます。
# 利用パッケージ library(tidyverse)
この記事では、基本的にパッケージ名::関数名()の記法を使うので、パッケージを読み込む必要はありません。ただし、作図コードがごちゃごちゃしないようにパッケージ名を省略しているため、ggplot2は読み込む必要があります。
magrittrパッケージのパイプ演算子%>%ではなく、ベースパイプ(ネイティブパイプ)演算子|>を使っています。%>%に置き換えても処理できます。
混合ガウス分布のパラメータ$\boldsymbol{\pi}, \boldsymbol{\mu}, \boldsymbol{\sigma}$を指定します。
# クラスタ数を指定 K <- 3 # 混合比率を指定 pi_k <- c(0.45, 0.25, 0.3) # K個の平均パラメータを指定 mu_k <- c(-4, 0, 2.6) # K個の標準偏差パラメータを指定 sigma_k <- c(1, 1.2, 0.8)
クラスタ数Kと各パラメータの値を指定します。
・K個のガウス分布の計算コード(クリックで展開)
設定したパラメータに応じて、混合ガウス分布の確率変数$x$がとり得る値を作成します。
# 確率変数の値を作成 x_vals <- seq( from = min(mu_k - sigma_k*3), to = max(mu_k + sigma_k*3), length.out = 401 ) head(x_vals)
## [1] -7.00 -6.97 -6.94 -6.91 -6.88 -6.85
$x$の値(x軸の値)をx_valsとします。この例では、クラスタごとに平均を中心に標準偏差の3倍を範囲としたときの最小値から最大値を範囲とします。
クラスタごとに、$x$の点ごとの確率密度を計算して、混合比率を掛けます。
# クラスタごとにガウス分布を計算 dens_cluster_df <- tidyr::expand_grid( k = 1:K, # クラスタ番号 x = x_vals # x軸の値 ) |> # クラスタごとに確率変数の値を複製 dplyr::group_by(k) |> # クラスタごとの計算用にグループ化 dplyr::mutate( dens_k = dnorm(x = x_vals, mean = mu_k[unique(k)], sd = sigma_k[unique(k)]), # クラスタごとの確率密度 dens_k_weighted = pi_k[unique(k)] * dens_k, # 重み付け確率密度 k = factor(k, levels = 1:K), # 色分け用に因子化 param_label_cluster = paste0( "list(", "k==", unique(k), ", mu==", mu_k[unique(k)], ", sigma==", sigma_k[unique(k)], ")" ), # パラメータラベル param_label_weighted = paste0( "list(", "k==", unique(k), ", pi==", pi_k[unique(k)], ", mu==", mu_k[unique(k)], ", sigma==", sigma_k[unique(k)], ")" ) # パラメータラベル ) |> dplyr::ungroup() # グループ化を解除 dens_cluster_df
## # A tibble: 1,203 × 6 ## k x dens_k dens_k_weighted param_label_cluster param_label_wei… ## <fct> <dbl> <dbl> <dbl> <chr> <chr> ## 1 1 -7 0.00443 0.00199 list(k==1, mu==-4, sigm… list(k==1, pi==… ## 2 1 -6.97 0.00485 0.00218 list(k==1, mu==-4, sigm… list(k==1, pi==… ## 3 1 -6.94 0.00530 0.00238 list(k==1, mu==-4, sigm… list(k==1, pi==… ## 4 1 -6.91 0.00578 0.00260 list(k==1, mu==-4, sigm… list(k==1, pi==… ## 5 1 -6.88 0.00631 0.00284 list(k==1, mu==-4, sigm… list(k==1, pi==… ## 6 1 -6.85 0.00687 0.00309 list(k==1, mu==-4, sigm… list(k==1, pi==… ## 7 1 -6.82 0.00748 0.00337 list(k==1, mu==-4, sigm… list(k==1, pi==… ## 8 1 -6.79 0.00814 0.00366 list(k==1, mu==-4, sigm… list(k==1, pi==… ## 9 1 -6.76 0.00885 0.00398 list(k==1, mu==-4, sigm… list(k==1, pi==… ## 10 1 -6.73 0.00961 0.00432 list(k==1, mu==-4, sigm… list(k==1, pi==… ## # … with 1,193 more rows
クラスタごとに、確率密度と混合比率によって重み付けした確率密度を計算します。ガウス分布の確率密度は、dnorm()で計算できます。
また、パラメータを表示するためのラベル(expression()の記法の文字列)を作成しておきます。
・K個のガウス分布の作図コード(クリックで展開)
パラメータの値を数式で表示するための文字列ベクトルを作成します。
# パラメータラベルを作成 param_label_cluster <- dens_cluster_df[["param_label_cluster"]] |> # ラベル列を抽出 unique() # 重複を削除 param_label_cluster
## [1] "list(k==1, mu==-4, sigma==1)" "list(k==2, mu==0, sigma==1.2)" ## [3] "list(k==3, mu==2.6, sigma==0.8)"
データフレームからラベル列を取り出して、unique()で重複を削除します。
$K$個のガウス分布のグラフを作成します。
# 混合ガウス分布を作図 ggplot() + geom_line(data = dens_cluster_df, mapping = aes(x = x, y = dens_k, color = k, linetype = "cluster"), size = 1) + # クラスタごとの確率密度 scale_color_hue(breaks = 1:K, labels = parse(text = param_label_cluster), name = "cluster") + # 色の凡例ラベル scale_linetype_manual(breaks = "cluster", values = "solid", labels = "cluster k distribution", name = "density") + # 線の種類 theme(legend.text.align = 0) + # 図の体裁 guides(linetype = guide_legend(override.aes = list(size = 0.5))) + # 凡例の体裁 labs(title = "Gaussian Distribution", subtitle = parse(text = paste0("K==", K)), x = "x", y = "density")
クラスタごとに色分けして描画します。
クラスタごとに1つの山になります。それぞれの面積が1(それぞれ積分すると1)になります。
・K個の重み付けガウス分布の作図コード(クリックで展開)
同様に、パラメータラベルを取り出します。
# パラメータラベルを作成 param_label_weighted <- dens_cluster_df[["param_label_weighted"]] |> # ラベル列を抽出 unique() # 重複を削除 param_label_weighted
## [1] "list(k==1, pi==0.45, mu==-4, sigma==1)" ## [2] "list(k==2, pi==0.25, mu==0, sigma==1.2)" ## [3] "list(k==3, pi==0.3, mu==2.6, sigma==0.8)"
$K$個の重み付けしたガウス分布のグラフを、元の分布と重ねて作成します。
# 混合ガウス分布を作図 ggplot() + geom_line(data = dens_cluster_df, mapping = aes(x = x, y = dens_k_weighted, color = k, linetype = "weighted"), size = 1) + # クラスタごとの割り引き確率密度 geom_line(data = dens_cluster_df, mapping = aes(x = x, y = dens_k, color = k, linetype = "cluster"), size = 1) + # クラスタごとの確率密度 scale_color_hue(breaks = 1:K, labels = parse(text = param_label_weighted), name = "cluster") + # 色の凡例ラベル scale_linetype_manual(breaks = c("cluster", "weighted"), values = c("dotdash", "solid"), labels = c("cluster k distribution", "weighted distribution"), name = "density") + # 線の種類 theme(legend.text.align = 0) + # 図の体裁 guides(linetype = guide_legend(override.aes = list(size = 0.5))) + # 凡例の体裁 labs(title = "Gaussian Mixture Distribution", subtitle = parse(text = paste0("K==", K)), x = "x", y = "density")
クラスタごとに色分けして、元の分布を点入りの破線で示します。
混合比率(各分布の重み)は0から1の値なので、それぞれ山が小さくなります。値を掛けただけなので、横軸方向には変化しません。
混合比率を掛けたことで、各山の面積が混合比率に一致し、全ての山の面積の合計が1になります。
・混合ガウス分布の計算と作図コード(クリックで展開)
混合ガウス分布の確率密度を計算します。
# 混合ガウス分布を計算 dens_mixture_df <- dens_cluster_df |> dplyr::group_by(x) |> # 点ごとの計算用にグループ化 dplyr::summarise(density = sum(dens_k_weighted), .groups = "drop") # K個の分布の加重平均を計算 dens_mixture_df
## # A tibble: 401 × 2 ## x density ## <dbl> <dbl> ## 1 -7 0.00199 ## 2 -6.97 0.00218 ## 3 -6.94 0.00238 ## 4 -6.91 0.00260 ## 5 -6.88 0.00284 ## 6 -6.85 0.00309 ## 7 -6.82 0.00337 ## 8 -6.79 0.00366 ## 9 -6.76 0.00398 ## 10 -6.73 0.00432 ## # … with 391 more rows
x列でグループ化することでx_valsの値($x$の点)ごとに、重み付けした確率密度を合計します。
混合ガウス分布のパラメータラベルを作成します。
# パラメータラベルを作成 param_label_mixture <- paste0( "list(", "K==", K, ", pi==(list(", paste0(pi_k, collapse = ", "), "))", ", mu==(list(", paste0(mu_k, collapse = ", "), "))", ", sigma==(list(", paste0(sigma_k, collapse = ", "), "))", ")" ) param_label_mixture
## [1] "list(K==3, pi==(list(0.45, 0.25, 0.3)), mu==(list(-4, 0, 2.6)), sigma==(list(1, 1.2, 0.8)))"
混合ガウス分布のグラフを、$K$個のガウス分布と重ねて作成します。
# 混合ガウス分布を作図 ggplot() + geom_line(data = dens_mixture_df, mapping = aes(x = x, y = density, linetype = "mixture"), color = "hotpink", size = 1) + # 混合分布の確率密度 geom_line(data = dens_cluster_df, mapping = aes(x = x, y = dens_k_weighted, color = k, linetype = "weighted"), size = 1) + # クラスタごとの割り引き確率密度 geom_line(data = dens_cluster_df, mapping = aes(x = x, y = dens_k, color = k, linetype = "cluster"), size = 1) + # クラスタごとの確率密度 scale_color_hue(breaks = 1:K, labels = parse(text = param_label_weighted), name = "cluster") + # 色の凡例ラベル scale_linetype_manual(breaks = c("cluster", "weighted", "mixture"), values = c("dotdash", "dashed", "solid"), labels = c("cluster k distribution", "weighted distribution", "mixture distribution"), name = "density") + # 線の種類 theme(legend.text.align = 0) + # 図の体裁 guides(linetype = guide_legend(override.aes = list(size = 0.5))) + # 凡例の体裁 labs(title = "Gaussian Mixture Distribution", subtitle = parse(text = param_label_mixture), x = "x", y = "density")
元の分布を点入りの破線、重み付けした分布を破線で示します。
$K$個の重み付けした分布に関して、x軸の値ごとに積み上げた形になります。
・混合ガウス分布の作図コード(クリックで展開)
最後に、混合ガウス分布だけのグラフを確認します。
# 混合ガウス分布を作図 ggplot() + geom_line(data = dens_mixture_df, mapping = aes(x = x, y = density, linetype = "mixture"), color = "hotpink", size = 1) + # 混合分布の確率密度 scale_linetype_manual(breaks = "mixture", values = "solid", labels = "mixture distribution", name = "density") + # 線の種類 guides(linetype = guide_legend(override.aes = list(size = 0.5))) + # 凡例の体裁 labs(title = "Gaussian Mixture Distribution", subtitle = parse(text = param_label_mixture), x = "x", y = "density")
$K$個の山を持つ分布なのが分かります。ただし、クラスタ間において平均が近い場合は、$K$個以下の山になることもあります。
この記事では、1次元混合ガウス分布の定義を確認しました。
参考文献
- 須山敦志『ベイズ推論による機械学習入門』(機械学習スタートアップシリーズ)杉山将監修,講談社,2017年.
おわりに
いよいよ混合分布編に突入です！
現在、これまでにブログ内で扱った分布に関してまとめています。1次元の混合ガウス分布は登場しなかったと思うのですが、多次元の方の可視化が微妙な出来ばえになったので、分かりやすいようにその前にやっておくことにしました。この記事の図と合わせて見れば問題ないと思います。
それと、1年だか2年前に混合分布を扱ったときに、調べてみて混乱した内容だったからでもあります。確率分布なんだから積分すると1になるだろうと思っていたところ、K(以上)になるような図で説明されている記事を見掛けました。K個の分布を覆うように曲線を描いて混合分布であると説明されてい(たように受け取り)ました。
そういった悩んだ記憶を上書きする意味でも書いておきたかったという事情もあります。
【次の内容】
