カテゴリカル分布完全ガイド — 多項分布との関係・共役性・LDAへの応用

サイコロを1回だけ振って「3の目が出た」と記録する。じゃんけんで「グー・チョキ・パー」のどれを出したかを1試行で観測する。MNISTの画像1枚に「0〜9」のどれかのラベルを付ける。スパムフィルタが1通のメールを「スパム / 普通 / 重要」のどれかに分類する。

これらに共通する構造は何でしょうか。「K個の選択肢の中から、ちょうど1つが選ばれる」 という確率現象です。コインの表裏という2択を扱うベルヌーイ分布を、より一般のK択に拡張したもの — それがカテゴリカル分布 (categorical distribution、カテゴリ分布) です。

カテゴリカル分布が活きる応用先は驚くほど広いです。たとえば、

  • ニューラル分類器のsoftmax出力 — 画像認識・自然言語処理の最終層は、ロジットをsoftmaxで正規化したカテゴリカル分布です。クロスエントロピー損失は、その分布の負対数尤度に他なりません。
  • トピックモデル (LDA) — 文書に対する潜在トピック割当も、各トピックからの単語生成も、すべてカテゴリカル分布から引かれます。
  • 強化学習の方策 — 離散行動空間のpolicy $\pi(a \mid s)$ もカテゴリカル分布で、Categorical方策勾配の出発点になります。
  • 多クラスナイーブベイズ・ロジスティック回帰 — 観測ラベルの生成モデルとして必ず登場します。

つまり、機械学習・統計モデリングを少しでも掘ろうとすると、カテゴリカル分布は避けて通れない離散分布の基本部品なのです。

本記事の内容

  • カテゴリカル分布の直感的イメージ (K面ダイス) と定式化 (one-hot 表現 / インデックス表現)
  • 期待値・分散・共分散・エントロピーの導出
  • ベルヌーイ・二項・多項分布との関係を「家系図」で整理
  • ディリクレ分布との共役性とベイズ更新の導出
  • 最尤推定 (count/N) とラプラス平滑化を含むベイズ推定
  • Python (NumPy/SciPy/PyTorch) でのサンプリングと推定の実装
  • LDA・softmaxクロスエントロピー・多クラス分類への応用
  • 概念図・関係マップを含む10枚の図解

alt

前提知識

この記事を読む前に、以下の記事を読んでおくと理解が深まります。

カテゴリカル分布の直感 — K面ダイスを1回振る

最初に「定義式ドン置き」を避けたいので、まずは身近な道具で直感を作りましょう。手元にK面のダイスがあるとします。$K=6$ なら普通のサイコロ、$K=20$ ならボードゲームの20面ダイス、$K=2$ なら表裏のあるコイン (ベルヌーイ分布の世界に戻る)。

このダイスは必ずしも公平とは限らず、k番目の面が出る確率を $\pi_k$ とします。たとえば $K=3$ のダイスで $\pi_1 = 0.5, \pi_2 = 0.3, \pi_3 = 0.2$ なら、3つの面のうち1番目の面が出やすく、3番目の面は出にくい、という偏ったダイスです。当然ですが、確率は全部足すと1になります。

$$ \sum_{k=1}^{K} \pi_k = 1, \qquad \pi_k \geq 0 \;(\text{すべての } k) $$

このダイスを1回だけ振って、出た面の番号 (または「どれが出たか」のラベル) を記録する。これがカテゴリカル分布の試行です。1回だけというのが重要で、n回振って各面が何回出たかを記録すると、それはすぐ後で見る多項分布になります。

ここで読者の頭にすでに2つの疑問が浮かんでいるはずです。

  1. 「出た面」を数学的にどう書けば、後の式変形が楽になるのか?
  2. 同じ「K択を1回」を表すのに、ベルヌーイ分布をそのまま拡張する書き方と、確率質量関数で書く書き方の2つがある気がする。どちらが本質か?

この2つの疑問を、次のセクションで定式化に変換していきます。

alt

カテゴリカル分布の数学的定義

インデックス表現 — もっとも素朴な書き方

最も自然なのは、「出たカテゴリ番号」をスカラー値 $z \in \{1, 2, \dots, K\}$ で表す方法です。確率質量関数は単純に、

$$ p(z = k \mid \boldsymbol{\pi}) = \pi_k, \qquad k = 1, 2, \dots, K $$

これでも形式的には何も問題ありません。「$z=3$ になる確率は $\pi_3$」と言うだけです。直感は完璧で、機械学習のコードでもラベル変数は通常この形 (整数のクラスID) で扱われます。

しかし、後で期待値・対数尤度・勾配を計算しようとすると、$z$ がスカラーだと和を取りにくい場面が出てきます。たとえば「N個のサンプルの対数尤度の和」を1つの式にまとめたいとき、$z$ がスカラーだと「$z_i = k$ のとき $\pi_k$ を選ぶ」という条件分岐が式に紛れ込み、解析が煩雑になります。

one-hot 表現 — K次元ベクトルとして書く

これを回避するために、機械学習・統計学で広く使われているのが one-hot 表現 (1 of K representation) です。出たカテゴリが $k$ なら、$K$ 次元ベクトル

$$ \boldsymbol{s} = (s_1, s_2, \dots, s_K)^\top, \qquad s_j = \begin{cases} 1 & (j = k) \\ 0 & (j \neq k) \end{cases} $$

を使って表します。たとえば $K=6$ で「3の目」が出たなら、

$$ \boldsymbol{s} = (0, 0, 1, 0, 0, 0)^\top $$

このベクトルの各要素は次の性質を満たします。

$$ s_k \in \{0, 1\}, \qquad \sum_{k=1}^{K} s_k = 1 $$

1つの成分だけが1で、残りはすべて0。「K個のうち1つを選ぶ」という事実が、ベクトルの構造に組み込まれているのがミソです。one-hot 表現のもとで、カテゴリカル分布の確率質量関数は次のように書けます。

$$ \mathrm{Cat}(\boldsymbol{s} \mid \boldsymbol{\pi}) = \prod_{k=1}^{K} \pi_k^{\,s_k} $$

最初は奇妙な書き方に見えますが、よく見るとこれは「$s_k = 1$ となる1つの $k$ だけ $\pi_k^1 = \pi_k$ が残り、それ以外は $\pi_k^0 = 1$ で消える」という賢い書き方です。つまり結果は単に $\pi_k$ になります — インデックス表現の $p(z=k) = \pi_k$ と完全に一致します。

なぜ one-hot 表現が便利か

one-hot 表現の真価は、対数尤度を取ったときに現れます。

$$ \log \mathrm{Cat}(\boldsymbol{s} \mid \boldsymbol{\pi}) = \sum_{k=1}^{K} s_k \log \pi_k $$

これがそのまま、ニューラル分類器の損失関数として有名なクロスエントロピー (cross-entropy) の形になります。深層学習で CrossEntropyLoss を使うとき、内部では実はこの式 (の符号反転) が計算されているわけです。ベクトルの内積1発で書けるので、行列演算と相性がよく、GPU上での実装も非常に効率的になります。

alt

ここまでで「K択を1回」を表す確率分布の数学的な書き方が2通り (インデックス / one-hot) あって、どちらも同じものを表していること、そして one-hot 表現が解析・実装の両面で扱いやすいことが見えました。次に、この分布の基本統計量 (期待値・分散) を求めて、その性質を定量的に押さえます。

期待値・分散・共分散・エントロピー

期待値: $\mathbb{E}[\boldsymbol{s}] = \boldsymbol{\pi}$

one-hot ベクトル $\boldsymbol{s}$ の期待値を成分ごとに計算します。

$$ \mathbb{E}[s_k] = \sum_{\boldsymbol{s}} s_k \cdot \mathrm{Cat}(\boldsymbol{s} \mid \boldsymbol{\pi}) $$

$\boldsymbol{s}$ の取りうる値は $\boldsymbol{e}_1, \boldsymbol{e}_2, \dots, \boldsymbol{e}_K$ ($K$ 個の単位ベクトル) です。$\boldsymbol{s} = \boldsymbol{e}_j$ のとき $s_k = \delta_{jk}$ (クロネッカーのデルタ) で、$\mathrm{Cat}(\boldsymbol{e}_j \mid \boldsymbol{\pi}) = \pi_j$。よって、

$$ \mathbb{E}[s_k] = \sum_{j=1}^{K} \delta_{jk}\, \pi_j = \pi_k $$

ベクトルでまとめると、

$$ \mathbb{E}[\boldsymbol{s}] = \boldsymbol{\pi} = (\pi_1, \pi_2, \dots, \pi_K)^\top $$

直感的にも明らかで、「平均的にどのカテゴリが出るか」と聞かれたら「$\pi_k$ の確率で k 番目」と答えるしかないからです。

分散と共分散: $\mathrm{Var}[s_k] = \pi_k (1 – \pi_k)$、$\mathrm{Cov}[s_k, s_l] = -\pi_k \pi_l$

$s_k$ は 0/1 の値しか取らないので、$s_k^2 = s_k$ という単純な性質があります。これを使うと、

$$ \mathbb{E}[s_k^2] = \mathbb{E}[s_k] = \pi_k $$

したがって分散は、

$$ \mathrm{Var}[s_k] = \mathbb{E}[s_k^2] – \mathbb{E}[s_k]^2 = \pi_k – \pi_k^2 = \pi_k (1 – \pi_k) $$

これはベルヌーイ分布の分散 $\pi(1-\pi)$ と完全に同じ形です — 当然で、各成分 $s_k$ が「k 番目が出るか出ないか」というベルヌーイ変数になっているからです。

次に、異なる成分間の共分散を計算します。$k \neq l$ のとき、$s_k$ と $s_l$ は同時に1にはなれない (どちらか1つしか出ないので) ので、$s_k s_l = 0$ が常に成り立ちます。したがって、

$$ \mathbb{E}[s_k s_l] = 0 $$

共分散の定義に代入すると、

$$ \mathrm{Cov}[s_k, s_l] = \mathbb{E}[s_k s_l] – \mathbb{E}[s_k]\mathbb{E}[s_l] = 0 – \pi_k \pi_l = -\pi_k \pi_l $$

共分散がになるのは、「k が出れば l は出ない」という排他関係を反映しています。one-hot ベクトルの成分は互いに足し合わせると常に1 (制約) なので、ある成分が大きくなれば他の成分は必然的に小さくなる — その負の相関が共分散に現れているわけです。

エントロピー: $H[\boldsymbol{\pi}] = -\sum_k \pi_k \log \pi_k$

カテゴリカル分布のエントロピーは、定義どおり、

$$ H[\boldsymbol{\pi}] = -\sum_{k=1}^{K} \pi_k \log \pi_k $$

これは情報理論でシャノンエントロピーとして現れる最も基本的な量です。「結果が出るまで何ビット (またはナット) の不確実性があるか」を測ります。

最大値は、すべての $\pi_k$ が等しい (一様分布 $\pi_k = 1/K$) のときで、

$$ H_{\max} = -\sum_{k=1}^{K} \frac{1}{K} \log \frac{1}{K} = \log K $$

逆に、ある $\pi_k$ が1に近づくとエントロピーは0に近づきます。完全に確定した結果には情報量がない、というシャノンの直感をそのまま表しています。

ここまでで、one-hot 表現を使うとカテゴリカル分布の統計量がきれいな閉形式で書けることがわかりました。次に、この分布が他の離散分布とどのように関連しているのかを整理し、「ベルヌーイ・二項・多項・カテゴリカル」の4つの位置関係を1枚の家系図にまとめます。

ベルヌーイ・二項・多項分布との関係

4つの分布を「2軸」で整理する

確率の入門書では、ベルヌーイ・二項・カテゴリカル・多項という4つの分布が次々と紹介され、それぞれ似て非なる定義が並んで初学者を混乱させがちです。ここで一度整理しておくと、4つは 「カテゴリ数 K」「試行回数 n」 の2軸で見事に分類できます。

試行回数 \ カテゴリ数 $K = 2$ (2値) $K \geq 2$ (多値)
$n = 1$ (1回) ベルヌーイ分布 $\mathrm{Bern}(\pi)$ カテゴリカル分布 $\mathrm{Cat}(\boldsymbol{\pi})$
$n \geq 1$ ($n$ 回) 二項分布 $\mathrm{Bin}(n, \pi)$ 多項分布 $\mathrm{Mult}(n, \boldsymbol{\pi})$

つまり、

  • ベルヌーイ → 二項: 同じ試行を $n$ 回繰り返して、「成功回数」を数える
  • ベルヌーイ → カテゴリカル: K択に拡張する ($K=2$ なら一致)
  • カテゴリカル → 多項: K択を $n$ 回繰り返して、「各カテゴリの出現回数」を数える
  • 二項 → 多項: 2択を多値化する

alt

限定形としての関係を式で確認

$K=2$ のカテゴリカル分布は、ベルヌーイ分布とまったく同じになります。$\boldsymbol{\pi} = (\pi, 1-\pi)$ とし、one-hot 表現 $\boldsymbol{s} = (s_1, s_2)$ で書くと、

$$ \mathrm{Cat}(\boldsymbol{s} \mid \boldsymbol{\pi}) = \pi^{s_1} (1-\pi)^{s_2} $$

$s_1 + s_2 = 1$ なので $s_2 = 1 – s_1$ を代入し、$s_1$ を新たに $x$ と書き直せば、

$$ \pi^{x} (1-\pi)^{1-x} $$

これはまさにベルヌーイ分布の確率質量関数です。

一方、$n$ 個の独立同分布なカテゴリカル変数 $\boldsymbol{s}^{(1)}, \dots, \boldsymbol{s}^{(n)}$ について、その $\boldsymbol{m} = \sum_{i=1}^{n} \boldsymbol{s}^{(i)}$ を考えます。$m_k$ は「n回中、k 番目のカテゴリが出た回数」になり、$\sum_k m_k = n$ を満たします。$\boldsymbol{m}$ の従う分布が多項分布で、

$$ \mathrm{Mult}(\boldsymbol{m} \mid n, \boldsymbol{\pi}) = \frac{n!}{m_1! m_2! \cdots m_K!} \prod_{k=1}^{K} \pi_k^{m_k} $$

カテゴリカル分布は、多項分布で $n=1$ とおいた特別な場合でもあります。$n=1$ なら $m_k$ は0か1で、上の式の多項係数 $n!/(m_1! \cdots m_K!) = 1!/(\cdots) = 1$ (どれか1つだけ1で残りは0) となり、$\prod_k \pi_k^{m_k}$ がそのまま残ります。これがカテゴリカルの定義式そのものです。

機械学習での使い分け

機械学習の文脈では、これら4分布の使い分けは次のように現れます。

  • 2値分類 (シグモイド出力) → ベルヌーイ分布 + バイナリクロスエントロピー
  • 多クラス分類 (softmax出力、1サンプル) → カテゴリカル分布 + クロスエントロピー
  • マルチラベル分類 (n個の独立な2値) → n個のベルヌーイ
  • 多クラス、複数試行集計 (例: 文書のbag-of-words) → 多項分布

カテゴリカル分布は、特に1つの観測ラベルの生成モデルを書くときに必須です。ニューラル分類器は何百万枚の画像を学習しますが、各画像のラベルを生成しているのは「1サンプル分のカテゴリカル分布」なのです。

これで、カテゴリカル分布が他の離散分布とどう関係するかが整理できました。次は、このカテゴリカル分布をベイズ流に扱う上で必須となる共役事前分布 — ディリクレ分布を導入します。

ディリクレ分布との共役性

なぜ共役性が嬉しいのか

ベイズ推定の枠組みでは、パラメータ $\boldsymbol{\pi}$ にも分布を考えます。観測データ $\mathcal{D}$ を見る前の信念を表す事前分布 $p(\boldsymbol{\pi})$ と、データから定まる尤度 $p(\mathcal{D} \mid \boldsymbol{\pi})$ をベイズの定理で組み合わせると、観測後の信念である事後分布 $p(\boldsymbol{\pi} \mid \mathcal{D})$ が得られます。

$$ p(\boldsymbol{\pi} \mid \mathcal{D}) = \frac{p(\mathcal{D} \mid \boldsymbol{\pi})\, p(\boldsymbol{\pi})}{p(\mathcal{D})} $$

ここで、事前分布と事後分布が同じ族に属するような事前分布のことを共役事前分布 (conjugate prior) と呼びます。共役性があれば、事後分布が解析的に書けるので、計算が劇的に楽になります — MCMCも変分推論もいらず、「パラメータの足し算」で更新が終わります。

カテゴリカル分布 (および多項分布) の共役事前分布が ディリクレ分布 (Dirichlet distribution) です。

ディリクレ分布の定義

$K$ 次元の確率ベクトル $\boldsymbol{\pi} = (\pi_1, \dots, \pi_K)$ で $\pi_k \geq 0$, $\sum_k \pi_k = 1$ を満たすもの全体 — つまり K-1 次元の単体 (simplex) — の上の分布として、ディリクレ分布は次で定義されます。

$$ \mathrm{Dir}(\boldsymbol{\pi} \mid \boldsymbol{\alpha}) = \frac{1}{B(\boldsymbol{\alpha})} \prod_{k=1}^{K} \pi_k^{\,\alpha_k – 1} $$

ここで $\boldsymbol{\alpha} = (\alpha_1, \dots, \alpha_K)$ はすべて正の集中度パラメータ (concentration parameter) で、$B(\boldsymbol{\alpha})$ は正規化定数 (多変量ベータ関数):

$$ B(\boldsymbol{\alpha}) = \frac{\prod_{k=1}^{K} \Gamma(\alpha_k)}{\Gamma\!\left(\sum_{k=1}^{K} \alpha_k\right)} $$

特に $K=2$ のディリクレ分布は、ちょうどベータ分布になります — ベータ分布の素直な多変量拡張だと思って差し支えありません。

共役性の証明 — パラメータの足し算で更新が完結する

$N$ 個のカテゴリカル観測 $\boldsymbol{s}^{(1)}, \dots, \boldsymbol{s}^{(N)}$ (i.i.d.) があるとします。各カテゴリ $k$ の出現回数を $N_k = \sum_{i=1}^{N} s_k^{(i)}$ と書きます ($\sum_k N_k = N$)。尤度は、

$$ p(\mathcal{D} \mid \boldsymbol{\pi}) = \prod_{i=1}^{N} \prod_{k=1}^{K} \pi_k^{s_k^{(i)}} = \prod_{k=1}^{K} \pi_k^{N_k} $$

事前分布をディリクレ $\mathrm{Dir}(\boldsymbol{\pi} \mid \boldsymbol{\alpha})$ とすると、事後分布は (正規化定数を無視して)、

$$ p(\boldsymbol{\pi} \mid \mathcal{D}) \propto p(\mathcal{D} \mid \boldsymbol{\pi})\, p(\boldsymbol{\pi}) = \left(\prod_{k=1}^{K} \pi_k^{N_k}\right) \left(\prod_{k=1}^{K} \pi_k^{\alpha_k – 1}\right) $$

積を1つにまとめると、

$$ p(\boldsymbol{\pi} \mid \mathcal{D}) \propto \prod_{k=1}^{K} \pi_k^{\,\alpha_k + N_k – 1} $$

この形は、ちょうどディリクレ分布の核そのものです。新しい集中度パラメータを $\alpha_k’ = \alpha_k + N_k$ と書けば、

$$ p(\boldsymbol{\pi} \mid \mathcal{D}) = \mathrm{Dir}(\boldsymbol{\pi} \mid \boldsymbol{\alpha} + \boldsymbol{N}) $$

となります。事後分布もディリクレ — 確かに共役です。更新ルールは驚くほどシンプルで、

観測ごとに、対応するカテゴリの集中度パラメータを1ずつ増やすだけ

$\boldsymbol{\alpha}$ の各成分は「観測前に各カテゴリを何回見たことにするか」という疑似カウント (pseudo-count) と解釈できます。たとえば $\boldsymbol{\alpha} = (1, 1, \dots, 1)$ なら「各カテゴリを1回ずつ見たことにする」一様事前。$\boldsymbol{\alpha} = (10, 10, \dots, 10)$ なら「各カテゴリを10回ずつ」と、より強い事前信念を表します。

alt

事後平均と事後モード

事後分布 $\mathrm{Dir}(\boldsymbol{\alpha} + \boldsymbol{N})$ の平均値は、ディリクレ分布の公式から、

$$ \mathbb{E}[\pi_k \mid \mathcal{D}] = \frac{\alpha_k + N_k}{\sum_{j=1}^{K} (\alpha_j + N_j)} = \frac{\alpha_k + N_k}{\alpha_0 + N} $$

ここで $\alpha_0 = \sum_k \alpha_k$ と書きました。これはスムージング (smoothing) されたカウント比率です。データの素朴な比率 $N_k / N$ に、疑似カウント $\alpha_k$ が「足し込まれて」いる形になっています。

事後モード (MAP推定値) は、$\alpha_k > 1$ のもとで、

$$ \hat{\pi}_k^{\mathrm{MAP}} = \frac{\alpha_k + N_k – 1}{\alpha_0 + N – K} $$

事前 $\boldsymbol{\alpha} = (1, 1, \dots, 1)$ (一様事前) のときはちょうど MAP = MLE になりますが、$\boldsymbol{\alpha} = (2, 2, \dots, 2)$ のときはラプラス平滑化 (Laplace smoothing) と呼ばれ、$N_k = 0$ のカテゴリにも非ゼロの確率を割り当てる効果があります — ナイーブベイズ分類器の常套手段ですね。

ここまでで、カテゴリカル分布をベイズ推定する道具が揃いました。次は、ベイズではなく最尤推定で $\boldsymbol{\pi}$ を点推定する方法を、ラグランジュ未定乗数法で導出します。

最尤推定 — count/N の正当化

最尤推定の問題設定

$N$ 個の i.i.d. 観測 $\boldsymbol{s}^{(1)}, \dots, \boldsymbol{s}^{(N)}$ が与えられたとき、$\boldsymbol{\pi}$ の最尤推定量を求めます。対数尤度は、

$$ \ell(\boldsymbol{\pi}) = \log \prod_{i=1}^{N} \mathrm{Cat}(\boldsymbol{s}^{(i)} \mid \boldsymbol{\pi}) = \sum_{i=1}^{N} \sum_{k=1}^{K} s_k^{(i)} \log \pi_k = \sum_{k=1}^{K} N_k \log \pi_k $$

ここで $N_k$ は前と同じく「カテゴリ $k$ の観測回数」です。これを $\boldsymbol{\pi}$ について最大化したい。ただし制約 $\sum_k \pi_k = 1$ があります。

ラグランジュ未定乗数法による導出

制約付き最大化なので、ラグランジュ関数を組みます。

$$ \mathcal{L}(\boldsymbol{\pi}, \lambda) = \sum_{k=1}^{K} N_k \log \pi_k – \lambda \left( \sum_{k=1}^{K} \pi_k – 1 \right) $$

$\pi_k$ で偏微分してゼロとおきます (各 k について)。

$$ \frac{\partial \mathcal{L}}{\partial \pi_k} = \frac{N_k}{\pi_k} – \lambda = 0 \quad \Longrightarrow \quad \pi_k = \frac{N_k}{\lambda} $$

ここで制約 $\sum_k \pi_k = 1$ を使うと、$\sum_k N_k / \lambda = 1$ より $\lambda = \sum_k N_k = N$。よって、

$$ \hat{\pi}_k^{\mathrm{MLE}} = \frac{N_k}{N} $$

直感どおりの結果です — 「k 番目が出た回数 / 全試行回数」 がそのまま $\pi_k$ の最尤推定値。ベルヌーイ分布での $\hat{\pi} = m/N$ の自然な多次元拡張になっています。

MLEの問題点 — ゼロ確率の罠

最尤推定はシンプルで美しいのですが、実用上の弱点があります。観測されなかったカテゴリの確率がゼロになることです。たとえば、新製品のレビューに「☆1, ☆2, ☆3, ☆4, ☆5」の5カテゴリがあるとして、まだ100件のレビューしか集まっていない状況で、「☆1」がたまたま0件だったとしましょう。MLEは $\hat{\pi}_{\star 1} = 0$ と答えます。これは「未来も絶対に ☆1 のレビューは来ない」と言っているに等しく、明らかに過信です。

この問題を防ぐ最もシンプルな方法が、先ほど見たベイズ推定 (ディリクレ事前) です。$\boldsymbol{\alpha} = (1, 1, \dots, 1)$ のラプラス平滑化を入れれば、観測ゼロのカテゴリにも (微小だが) 確率が割り振られ、ゼロ確率の罠を避けられます。

alt

ここまでで、頻度論 (MLE) とベイズの2つの推定アプローチを見ました。理論はそろったので、ここからは実装に進みます。

Python実装 — サンプリングと推定

NumPy/SciPy でカテゴリカル分布を扱う

カテゴリカル分布は、NumPyでは np.random.choice、SciPyでは scipy.stats.multinomial (n=1 のとき) として扱えます。最小例から始めましょう。

import numpy as np
from scipy import stats

# パラメータ: K=4 のカテゴリカル分布
pi = np.array([0.5, 0.25, 0.15, 0.10])
K = len(pi)

# (A) インデックス表現でサンプリング
rng = np.random.default_rng(0)
z_samples = rng.choice(K, size=1000, p=pi)  # 値域 0..K-1
print("インデックス表現の最初10サンプル:", z_samples[:10])

# (B) one-hot 表現でサンプリング (1 試行の多項分布)
s_samples = rng.multinomial(n=1, pvals=pi, size=1000)  # 形状 (1000, K)
print("one-hot 表現の最初3サンプル:")
print(s_samples[:3])

# 同等性の確認: one-hot を argmax するとインデックスに戻る
z_from_onehot = s_samples.argmax(axis=1)
print("two methods agree?", np.array_equal(z_from_onehot[:100],
                                            rng.choice(K, size=100, p=pi)) == False)
# (上は別 seed で進んでいるので False で正常)

このコードでは、同じパラメータ $\boldsymbol{\pi}$ から (A) インデックス表現(B) one-hot 表現の2通りでサンプリングしています。np.random.choicenp.random.multinomial(n=1) は本質的に同じ分布を別の形式で出力しているだけで、argmax で相互変換できます。実装で使いやすい方を選べばよいでしょう — 損失関数の計算には one-hot、ラベル管理には整数インデックスが便利です。

経験度数と最尤推定

サンプルから最尤推定値 $\hat{\pi}_k = N_k / N$ を計算します。

import numpy as np
import matplotlib.pyplot as plt

rng = np.random.default_rng(42)
true_pi = np.array([0.5, 0.25, 0.15, 0.10])
K = len(true_pi)

# 試行回数を増やしながら推定値の収束を観察
Ns = [10, 50, 200, 1000, 5000, 20000]
all_samples = rng.choice(K, size=max(Ns), p=true_pi)

print(f"{'N':>6} | " + " | ".join([f"pi_{k+1}" for k in range(K)]) + " | RMSE")
print("-" * 50)
for N in Ns:
    counts = np.bincount(all_samples[:N], minlength=K)
    pi_hat = counts / N
    rmse = np.sqrt(np.mean((pi_hat - true_pi) ** 2))
    print(f"{N:>6} | " + " | ".join([f"{p:.3f}" for p in pi_hat]) + f" | {rmse:.4f}")

出力例 (シードを固定しているので決定的に再現できます) を見ると、$N$ が小さいうちは推定値が真値から大きくぶれていますが、$N=20000$ ではほぼ真値 $(0.5, 0.25, 0.15, 0.10)$ に一致しています。RMSE は $1/\sqrt{N}$ で減少していて、これは i.i.d. 観測における最尤推定の漸近理論と整合します。

ベイズ推定 — ディリクレ事後

次に、ディリクレ事前を入れたベイズ推定を実装します。前述のとおり、事後の集中度パラメータは「事前 + カウント」の足し算で済みます。

import numpy as np
from scipy import stats

rng = np.random.default_rng(7)
true_pi = np.array([0.5, 0.25, 0.15, 0.10])
K = len(true_pi)

# 事前: 一様 (Laplace事前)
alpha_prior = np.ones(K)

# 観測データ
N = 100
samples = rng.choice(K, size=N, p=true_pi)
counts = np.bincount(samples, minlength=K)

# 事後パラメータ = 事前 + カウント
alpha_post = alpha_prior + counts
print("事前 alpha:", alpha_prior)
print("カウント N_k:", counts)
print("事後 alpha:", alpha_post)

# 事後平均 (ベイズ点推定の標準的な選び方)
post_mean = alpha_post / alpha_post.sum()
print(f"事後平均 (Bayes): {post_mean}")
print(f"最尤推定  (MLE)  : {counts / N}")
print(f"真値              : {true_pi}")

# 事後分布から1万サンプル引いてヒストグラム
pi_samples = rng.dirichlet(alpha_post, size=10000)  # 形状 (10000, K)
print("\n事後サンプルから計算した 95% 信用区間:")
for k in range(K):
    lo, hi = np.quantile(pi_samples[:, k], [0.025, 0.975])
    print(f"  pi_{k+1}: [{lo:.3f}, {hi:.3f}] (真値 {true_pi[k]:.3f})")

このコードのポイントは2つあります。1つ目、事後パラメータの計算が alpha_prior + counts の1行で完結していること — これがディリクレ-カテゴリカル共役の威力です。2つ目、np.random.dirichlet で事後分布から直接サンプリングしてしまえば、任意の汎関数の事後分布 (例: 信用区間や最大値の確率) もモンテカルロで簡単に計算できます。出力を見ると、真値はどのカテゴリでも 95% 信用区間にきちんと収まっており、ベイズ推定が点推定よりも豊富な情報を提供してくれることがわかります。

alt

PyTorch でのカテゴリカル分布

深層学習のフレームワークでは、カテゴリカル分布は明示的にクラスとして提供されています。PyTorch なら torch.distributions.Categorical です。

import torch
import torch.nn.functional as F

# ロジット (生のスコア、softmax前) から作る
logits = torch.tensor([2.0, 1.0, 0.1, -1.0])
dist = torch.distributions.Categorical(logits=logits)

# 確率
print("確率 pi:", dist.probs)               # softmax(logits)
print("対数確率 log pi:", dist.logits)       # 正規化された log probs

# サンプリング
samples = dist.sample((1000,))
print("サンプル形状:", samples.shape)
print("経験度数:", torch.bincount(samples, minlength=4) / 1000.0)

# ある観測の対数尤度 (これがクロスエントロピー)
y = torch.tensor(0)  # 真クラス
log_p_y = dist.log_prob(y)
print(f"log P(y={y.item()}) = {log_p_y.item():.4f}")

# PyTorch の cross_entropy と一致するか確認
ce = F.cross_entropy(logits.unsqueeze(0), y.unsqueeze(0))
print(f"F.cross_entropy = {ce.item():.4f}  (-log_p_y = {-log_p_y.item():.4f})")

F.cross_entropy-dist.log_prob(y) の値がぴったり一致します。これは「分類器の損失関数は、出力カテゴリカル分布の負対数尤度そのもの」という事実を、コードで直接確認している形です。深層学習で CrossEntropyLoss を呼ぶたびに、裏では K 次元のカテゴリカル分布の最尤推定問題を解いているわけです。

ここまでで理論と実装が一通りそろいました。最後に、カテゴリカル分布が実際の応用 — トピックモデルや深層分類器 — でどう活躍しているかを見て、本記事を締めます。

応用 — LDA・softmax・クロスエントロピー

LDA (Latent Dirichlet Allocation) におけるカテゴリカル分布

トピックモデルの代表格であるLDAは、文書集合から潜在トピックを抽出する確率モデルです。生成過程は次のように記述されます。

  1. 各文書 $d$ について、トピック分布 $\boldsymbol{\theta}_d \sim \mathrm{Dir}(\boldsymbol{\alpha})$
  2. 各トピック $k$ について、単語分布 $\boldsymbol{\phi}_k \sim \mathrm{Dir}(\boldsymbol{\beta})$
  3. 文書 $d$ の各単語位置 $n$ について、
  4. トピック割当 $z_{d,n} \sim \mathrm{Cat}(\boldsymbol{\theta}_d)$
  5. 単語 $w_{d,n} \sim \mathrm{Cat}(\boldsymbol{\phi}_{z_{d,n}})$

LDAの中で2回カテゴリカル分布が出てきます。トピック割当 $z$ の生成と、単語 $w$ の生成です。そして、それぞれのカテゴリカル分布のパラメータがディリクレ分布から生成されます — 共役性が効いてGibbsサンプリングの形が美しく書けるのが、LDAが成功した大きな理由の1つです。

「カテゴリカル + ディリクレ」の組み合わせは、LDAだけでなく隠れマルコフモデル (HMM) の出力分布混合カテゴリカルモデルLatent class analysis など、潜在変数モデル全般の基礎になっています。

alt

softmax クロスエントロピー損失との関連

ニューラル分類器の最終層は、ロジットベクトル $\boldsymbol{z} \in \mathbb{R}^K$ を softmax で正規化してカテゴリカル分布のパラメータ $\boldsymbol{\pi}$ に変換します。

$$ \pi_k = \mathrm{softmax}(\boldsymbol{z})_k = \frac{e^{z_k}}{\sum_{j=1}^{K} e^{z_j}} $$

訓練時の損失は、真ラベル $y$ (one-hot $\boldsymbol{t}$) とのクロスエントロピー:

$$ \mathcal{L}(\boldsymbol{z}, y) = -\sum_{k=1}^{K} t_k \log \pi_k = -\log \pi_y $$

これは前のセクションで導いたカテゴリカル分布の対数尤度 (の符号反転) そのものです。つまり、深層分類器の訓練は

「観測ラベル $y$ がカテゴリカル分布 $\mathrm{Cat}(\boldsymbol{\pi}(\boldsymbol{z}))$ から生成されたとして、$\boldsymbol{z}$ について最尤推定を行う」

という枠組みに他なりません。ImageNetの100万枚を分類するモデルも、本質はカテゴリカル分布の最尤推定なのです。

さらに、勾配計算で softmaxlog をまとめた log_softmax を使うのは、数値安定性の理由 (大きなロジットに対するオーバーフロー回避) ですが、これも「カテゴリカル分布の対数尤度をいかに正確に計算するか」という1点に集約されます。

alt

強化学習の方策とサンプリング

離散行動空間における強化学習の方策 $\pi(a \mid s)$ は、状態 $s$ で行動 $a$ を選ぶ確率 — まさにカテゴリカル分布です。REINFORCEActor-Critic のような方策勾配法は、この方策カテゴリカル分布から行動をサンプリングし、報酬で重み付けして対数尤度の勾配を更新します。サンプリングは PyTorch なら Categorical(logits=...).sample() の1行、ログ確率は .log_prob(action) で取れる — カテゴリカル分布の理解がそのまま実装力に直結します。

K と N のスケーリング

最後に、カテゴリカル分布の難しさがカテゴリ数 K に依存することに触れておきましょう。

  • エントロピーは最大で $\log K$ — K が増えると不確実性が指数的でなく対数的に増える
  • 最尤推定のRMSEはおおよそ $O(\sqrt{K/N})$ — K が増えると同じ精度を出すのに比例して多くのサンプルが必要
  • 語彙が数万に及ぶ言語モデルでは、巨大なKに対応するために hierarchical softmax や noise-contrastive estimation などの近似手法が使われる

つまり、カテゴリカル分布は「離散分布の単純な拡張」でありながら、$K$ が大きいときの取り扱いは現代の機械学習における重要な研究テーマの1つでもあります。

alt

まとめ

本記事では、カテゴリカル分布について「直感 → 定義 → 統計量 → 関連分布 → 共役性 → 推定 → 実装 → 応用」というロードマップで丁寧に辿ってきました。要点を整理しておきます。

  • K面ダイスを1回振るときの分布。インデックス表現とone-hot表現の2通りで書け、後者の方が解析・実装で扱いやすい
  • one-hot表現の確率質量関数 $\mathrm{Cat}(\boldsymbol{s} \mid \boldsymbol{\pi}) = \prod_k \pi_k^{s_k}$。対数尤度がそのままクロスエントロピーの形になる
  • 期待値 $\mathbb{E}[\boldsymbol{s}] = \boldsymbol{\pi}$、分散 $\mathrm{Var}[s_k] = \pi_k(1-\pi_k)$、共分散 $\mathrm{Cov}[s_k, s_l] = -\pi_k \pi_l$ (排他性ゆえ負)、エントロピー $H = -\sum_k \pi_k \log \pi_k$ で最大値 $\log K$
  • ベルヌーイ ($K=2$) ・二項 (ベルヌーイの $n$ 回和)・多項 (カテゴリカルの $n$ 回和) と「カテゴリ数・試行回数」の2軸で整理できる
  • 共役事前分布はディリクレ。事後は $\mathrm{Dir}(\boldsymbol{\alpha} + \boldsymbol{N})$ — パラメータの足し算で更新が完結
  • 最尤推定は $\hat{\pi}_k = N_k / N$。ゼロ確率の罠を防ぐにはラプラス平滑化 (一様ディリクレ事前) が有効
  • 応用先: LDA・HMM・softmax分類器・強化学習方策など、離散ラベルが出てくる場面ほぼすべて

カテゴリカル分布は「離散分布の基本構成要素」です。多項分布もLDAもsoftmax分類器も、結局はカテゴリカル分布を組み立て直したものに過ぎません。ここを掴んでおけば、後続の確率モデルがすべて「カテゴリカルとガウスのレゴブロックの組合せ」に見えてきます。

次のステップとして、以下の記事も参考にしてください。