はじめに

こんばんは。
MYJLab Advent Calender 2024の2日目、本日の担当はM1の yoshiki495 です。

本記事では、「ゼロから作るDeep Learning ❺」の要点をまとめてみようかと思います。
統計学および機械学習において基礎的かつ重要な概念から具体的な生成モデルの実装方法までが記載されていますので、本書片手に要点を確認していただけると幸いです。

1: 正規分布

1.1 確率の基礎

現実世界では不確実性が常に存在します。
この不確実性を定量的に扱うために「確率」が用いられます。

1.1.1 確率変数と確率分布

確率変数は、結果が確率的に決まる変数です。
例えば、サイコロの目は確率変数としてモデル化できます。
サイコロの各目が出る確率は均等であり、この場合の確率分布は一様分布と呼ばれます。

表 1-1 サイコロの目の確率分布

サイコロの目 x	1	2	3	4	5	6
確率 $p(x)$	$1/6$	$1/6$	$1/6$	$1/6$	$1/6$	$1/6$

確率分布は、すべての可能な値に対してその発生確率を示します。
確率分布が成立するためには、各確率が $0$ 以上 $1$ 以下であり、全ての確率の和が $1$ である必要があります。

確率分布の条件:

$0 \leq p(x_k) \leq 1 \quad (k=1,2,\dots,N)$
$\sum_{k=1}^{N} p(x_k) = 1$

1.1.2 確率分布の種類

確率分布は大きく「離散型確率分布」と「連続型確率分布」に分けられます。

離散型確率分布: 確率変数が離散値を取る場合（例: サイコロの目）
連続型確率分布: 確率変数が連続値を取る場合（例: 身長、気温）

連続型確率分布では、特定の点における確率ではなく、区間における確率を計算します。
これは「確率密度」として表現され、積分によって求められます。

1.1.3 期待値と分散

確率分布の中心的な指標として「期待値（平均）」と「分散」があります。

期待値 $E[X]$: 確率変数の平均値
- 離散型: $E[X] = \sum_{k=1}^{N} x_k p(x_k)$
- 連続型: $E[X] = \int_{-\infty}^{\infty} x p(x) dx$
分散 $\text{Var}[X]$: 確率変数のばらつきの指標
- 離散型: $\text{Var}[X] = \sum_{k=1}^{N} (x_k - \mu)^2 p(x_k)$
- 連続型: $\text{Var}[X] = \int_{-\infty}^{\infty} (x - \mu)^2 p(x) dx$

分散の平方根を取ったものを「標準偏差」と呼び、データの散らばり具合を直感的に理解することができる値となります。

1.2 正規分布

正規分布は統計学および機械学習において最も重要な連続型確率分布です。
これは「ガウス分布」とも呼ばれ、自然現象や測定誤差など多くの現象がこの分布に従います。

1.2.1 正規分布の確率密度関数

正規分布の確率密度関数は以下の式で定義されます。

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

ここで、

$\mu$ : 平均
$\sigma$ : 標準偏差

図 1-5 平均0、標準偏差1の正規分布

正規分布は平均を中心に左右対称の釣り鐘型の曲線を描きます。

1.2.2 正規分布のコード

以下は、正規分布の確率密度関数をPythonで実装した例です。

import numpy as np
import matplotlib.pyplot as plt

def normal(x, mu=0, sigma=1):
    y = 1 / (np.sqrt(2 * np.pi) * sigma) * np.exp(-(x - mu)**2 / (2 * sigma**2))
    return y

x = np.linspace(-10, 10, 1000)
y0 = normal(x, mu=-3)
y1 = normal(x, mu=0)
y2 = normal(x, mu=5)

plt.plot(x, y0, label='$\mu=-3$')
plt.plot(x, y1, label='$\mu=0$')
plt.plot(x, y2, label='$\mu=5$')
plt.legend()
plt.xlabel('x')
plt.ylabel('y')
plt.show()

図1-6 平均が-3, 0, 5の正規分布（標準偏差はすべて1）

平均を変えることで、分布のピーク位置が移動します。

1.3 中心極限定理

正規分布が多くの現象で見られる理由は、「中心極限定理」によるものです。
中心極限定理は、任意の分布から独立にサンプリングした大数の平均が正規分布に近づくという定理です。

1.3.1 中心極限定理とは

中心極限定理によれば、任意の確率分布 $p(x)$ から独立に生成された $N$ 個のサンプル ${x^{(1)}, x^{(2)}, \dots, x^{(N)}}$ のサンプル平均 $\bar{x}$ は、サンプルサイズ $N$ が大きくなるにつれて、平均 $\mu$、分散 $\sigma^2 / N$ の正規分布に近づきます。

1.4 サンプル和の確率分布

中心極限定理を理解するために、サンプル平均ではなく「サンプル和」について考えます。サンプル和 $s$ は、以下の式で表されます。

$s=x^{(1)} + x^{(2)} + \dots + x^{(N)}$

1.4.1 サンプル和の期待値と分散

サンプル和の期待値と分散は以下の通りです。

$E[s]=N\mu$

$\text{Var}[s] = N\sigma^2$

1.4.2 コードで確かめる

以下のコードは、一様分布からサンプル和の分布を実験的に確認します。

import numpy as np
import matplotlib.pyplot as plt

x_sums = []
N = 5

for _ in range(10000):
    xs = np.random.rand(N)  # 一様分布からの乱数
    s = np.sum(xs)         # 和を求める
    x_sums.append(s)

# 正規分布
mu = N / 2
sigma = np.sqrt(N / 12)
x_norm = np.linspace(0, N, 1000)
y_norm = (1 / (np.sqrt{2 * np.pi} * sigma)) * np.exp(-(x_norm - mu)**2 / (2 * sigma**2))

# グラフの描画
plt.hist(x_sums, bins='auto', density=True, alpha=0.6, color='g', label='Sample Sum')
plt.plot(x_norm, y_norm, color='r', label='Normal Distribution')
plt.title(f'N={N}')
plt.xlabel('Sum of Samples')
plt.ylabel('Probability Density')
plt.legend()
plt.show()

図1-12 サンプル和のヒストグラム（サンプルサイズは5）

この結果、サンプル和の分布が正規分布に近づくことが確認できます。

1.4.3 一様分布の平均と分散

一様分布 $U(0,1)$ の平均と分散は以下の通りです。

$\mu = \frac{1}{2}, \quad \sigma^2 = \frac{1}{12}$

これらは積分によって計算され、サンプル和の理論的期待値と分散に一致します。

2: 最尤推定

2.1 生成モデルの概要

最尤推定（Maximum Likelihood Estimation, MLE）は、観測データに最も適合するモデルのパラメータを推定する手法です。

2.1.1 生成モデルとは

生成モデルの目的はデータの生成過程をモデル化し、観測データの確率分布を推定することです。
具体的には、観測データを最もよく説明する確率分布を見つけます。

2.1.2 母集団とサンプル

統計学では、全体のデータを「母集団」と呼び、その中から抽出した一部を「サンプル」と呼びます。
母集団の分布を推定するために、サンプルデータを用います。

2.2 実データを使った生成モデルの実装

ここでは、身長データを用いて生成モデルを実装し、最尤推定によって正規分布のパラメータを推定します。

2.2.1 身長データセットの読み込み

以下のコードで、身長データを読み込み、ヒストグラムをプロットして分布を確認します。

import os
import numpy as np
import matplotlib.pyplot as plt

file = os.path.join("step02", "height.txt")
data = np.loadtxt(file)
plt.hist(data, bins="auto", density=True)
plt.xlabel("Height (cm)")
plt.ylabel("Probability Density")
plt.show()

2.2.2 正規分布による生成モデル

身長データが正規分布に従うと仮定し、そのパラメータ（平均 $\mu$ と標準偏差 $\sigma$）を推定します。

import numpy as np
import matplotlib.pyplot as plt

# 正規分布の定義
def normal(x, mu=0, sigma=1):
    y = 1 / (np.sqrt{2 * np.pi} * sigma) * np.exp(-(x - mu)**2 / (2 * sigma**2))
    return y

# サンプルデータ
x = np.linspace(-10, 10, 1000)
y0 = normal(x, mu=-3)
y1 = normal(x, mu=0)
y2 = normal(x, mu=5)

plt.plot(x, y0, label='$\mu=-3$')
plt.plot(x, y1, label='$\mu=0$')
plt.plot(x, y2, label='$\mu=5$')
plt.legend()
plt.xlabel('x')
plt.ylabel('y')
plt.show()

図2-6 サンプルのデータ分布と正規分布の比較

実データのヒストグラムと推定した正規分布が重なっており、生成モデルの精度を確認できます。

2.3 最尤推定の理論

最尤推定は、観測データの尤度を最大化するパラメータを求める手法です。
以下にその理論的背景と具体的な計算方法を示します。

尤度関数 (Likelihood Function)

観測データ $D = {x^{(1)}, x^{(2)}, \dots, x^{(N)}}$ が与えられたとき、そのデータが特定のパラメータ $\theta$ によって生成される確率を尤度と呼びます。

$L(\theta) = p(D | \theta)$

対数尤度 (Log-Likelihood)

計算を簡便にするため、尤度の対数を取ったものを対数尤度と呼びます。

$\log L(\theta) = \sum_{i=1}^{N} \log p(x^{(i)} | \theta)$

微分を使って最大値を探す

関数の最大化には微分を用います。以下はその具体例です。

# 関数の定義
def f(x):
    return -2 * x**2 + 3 * x + 4

# 微分
def df_dx(x):
    return -4 * x + 3

# 最大値を求める
x_max = 3 / 4
y_max = f(x_max)
print(f"最大値はx={x_max}のとき、y={y_max}")

2.3.1 正規分布の最大尤推定

正規分布の密度関数を用いて、最尤推定を行います。
観測データ $D = {x^{(1)}, x^{(2)}, \dots, x^{(N)}}$ に対して、以下の式を最大化します。

$\log L(\mu, \sigma) = -\frac{N}{2} \log (2\pi\sigma^2) - \frac{1}{2\sigma^2} \sum_{i=1}^{N} (x^{(i)} - \mu)^2$

パラメータ推定式

平均 $\mu$: $\mu = \frac{1}{N} \sum_{n=1}^{N} x^{(n)}$
標準偏差 $\sigma$: $\sigma = \sqrt{ \frac{1}{N} \sum_{n=1}^{N} (x^{(n)} - \mu)^2 }$

2.4 生成モデルの用途

生成モデルは、観測データの特徴を捉えて新しいデータを生成することが可能です。
以下のコード例では、身長データから推定した正規分布を基に新たなデータを生成し、元のデータと比較しています。

import os
import numpy as np
import matplotlib.pyplot as plt

# データの読み込み
path = os.path.join(os.path.dirname(__file__), 'height.txt')
xs = np.loadtxt(path)

# パラメータの推定
mu = np.mean(xs)
sigma = np.std(xs)

# 新しいサンプルの生成
samples = np.random.normal(mu, sigma, 10000)

# ヒストグラムのプロット
plt.hist(xs, bins='auto', density=True, alpha=0.7, label='Original Data')
plt.hist(samples, bins='auto', density=True, alpha=0.7, label='Generated Data')
plt.xlabel('Height (cm)')
plt.ylabel('Probability Density')
plt.legend()
plt.show()

図2-8: 生成モデルによるデータ生成結果

生成されたデータの分布が元のデータと一致していることが確認でき、生成モデルの有効性が示されます。

3: 多次元正規分布

3.1 NumPyと多次元配列

NumPy多次元配列を扱うPythonライブラリとして、多く利用されています。

3.1.1 多次元配列とベクトルの基礎

多次元配列（テンソル）は、複数の値を1つの構造で扱うためのデータ形式です。
次元や軸といった概念があり、スカラー（0次元）、ベクトル（1次元）、行列（2次元）といった形で扱われます。

スカラー: 単一の値
ベクトル: 値が1つの方向に並ぶデータ構造
行列: 値が2つの方向（行と列）に並ぶデータ構造

これらをNumPyで扱う方法は以下となっています。

import numpy as np
# ベクトル
x = np.array([1, 2, 3])
print(f"ベクトル: {x}, 次元: {x.ndim}, 形状: {x.shape}")

# 行列
W = np.array([[1, 2], [3, 4]])
print(f"行列: {W}, 次元: {W.ndim}, 形状: {W.shape}")

3.1.2 要素ごとの演算

NumPyでは、多次元配列間での要素ごとの演算（アダマール積など）が簡単に行えます。

A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])
print(A + B)  # 要素ごとの和
print(A * B)  # 要素ごとの積

3.1.3 ベクトルの内積と行列積

ベクトルや行列の積も簡単に計算できます。

# ベクトルの内積
a = np.array([1, 2])
b = np.array([3, 4])
print(np.dot(a, b))

# 行列の積
A = np.array([[1, 2], [3, 4]])
B = np.array([[5, 6], [7, 8]])
print(np.dot(A, B))

3.2 多次元正規分布

3.2.1 多次元正規分布の理論

多次元正規分布は、以下の確率密度関数で定義されます：

$N(x; \mu, \Sigma) = \frac{1}{(2\pi)^{D/2} |\Sigma|^{1/2}} \exp\left(-\frac{1}{2}(x - \mu)^T \Sigma^{-1} (x - \mu) \right)$

$x$: $D$次元ベクトル
$\mu$: 平均ベクトル
$\Sigma$: 共分散行列

3.2.2 多次元正規分布の実装

def multivariate_normal(x, mu, cov):
    D = len(x)
    det = np.linalg.det(cov)
    inv = np.linalg.inv(cov)
    norm_const = 1 / np.sqrt((2 * np.pi)**D * det)
    return norm_const * np.exp(-0.5 * (x - mu).T @ inv @ (x - mu))

3.3 多次元正規分布の可視化

3.3.1 3Dプロットと等高線

2次元正規分布の例を、3Dグラフと等高線で描画します。

コード例

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

# 平均と共分散
mu = np.array([0, 0])
cov = np.array([[1, 0.5], [0.5, 1]])

# グリッドを作成
x = np.linspace(-3, 3, 100)
y = np.linspace(-3, 3, 100)
X, Y = np.meshgrid(x, y)

# Z値を計算
Z = np.array([[multivariate_normal(np.array([xi, yi]), mu, cov) for xi, yi in zip(x_row, y_row)] for x_row, y_row in zip(X, Y)])

# グラフを描画
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(X, Y, Z, cmap="viridis")
plt.show()

3.4 多次元正規分布の最尤推定

3.4.1 最尤推定の式

多次元正規分布の最尤推定は、以下の式を用いて行います：

平均: $\hat{\mu} = \frac{1}{N} \sum_{i=1}^N x_i$
共分散行列: $\hat{\Sigma} = \frac{1}{N} \sum_{i=1}^N (x_i - \hat{\mu})(x_i - \hat{\mu})^T$

コード例

python
コードをコピーする
python
コードをコピーする
def estimate_parameters(data):
    mu = np.mean(data, axis=0)
    cov = np.cov(data, rowvar=False)
    return mu, cov

3.4.2 データを使った実装

サンプルデータから推定されたパラメータを使い、多次元正規分布をプロットします。

ステップ 4: 混合ガウスモデル (GMM)

これまで学んだ正規分布を拡張し、複数の正規分布を組み合わせた「混合ガウスモデル (GMM)」について学びます。GMMは、多峰性分布をモデル化するために使われ、機械学習やデータ生成で重要な役割を果たします。本章では以下を解説します。

身の回りにある多峰性分布の例
混合ガウスモデルの理論
GMMを用いたデータ生成
パラメータ推定とその課題

4.1 身の回りにある多峰性分布

4.1.1 男女混合の身長分布

男女別に見ると、身長は正規分布に従う傾向がありますが、性別を区別せずに集計すると、異なる平均値を持つ2つの正規分布が混ざった分布（多峰性分布）が得られます。

4.1.2 アリの体長分布

アリの体長は、役割（女王アリ、働きアリなど）によって異なるため、1つの分布ではモデル化できない、多峰性分布の例となります。

4.2 混合ガウスモデルの理論

4.2.1 混合ガウスモデル (GMM) の概要

GMMは、複数の正規分布を加重平均で組み合わせたモデルです。
データが以下の手順で生成されると仮定します：

複数の正規分布の中から1つを選ぶ
選ばれた正規分布からデータを生成する

GMMの確率密度関数は次の式で表されます：

$p(x) = \sum_{k=1}^{K} \phi_k \mathcal{N}(x; \mu_k, \Sigma_k)$

ここで、

$K$: 正規分布の数
$\phi_k$: 第 $k$ の正規分布の重み ($\sum_{k=1}^K \phi_k = 1$)

$\sum_{k=1}^K \phi_k = 1$
$\mathcal{N}(x; \mu_k, \Sigma_k)$: 第 $k$ の正規分布

4.3 GMMを用いたデータ生成

4.3.1 データ生成の手順

$\phi_k$ に基づいてどの正規分布を使うか決定
選ばれた正規分布からサンプルを生成

import numpy as np
import matplotlib.pyplot as plt

# GMMのパラメータ
mus = np.array([[2.0, 4.0], [6.0, 8.0]])
covs = np.array([[[1.0, 0.5], [0.5, 1.0]], [[1.5, -0.5], [-0.5, 1.5]]])
phis = np.array([0.4, 0.6])

def sample_gmm(mus, covs, phis, n_samples=500):
    samples = []
    for _ in range(n_samples):
        k = np.random.choice(len(phis), p=phis)  # 分布選択
        sample = np.random.multivariate_normal(mus[k], covs[k])  # サンプル生成
        samples.append(sample)
    return np.array(samples)

# サンプル生成
samples = sample_gmm(mus, covs, phis)
plt.scatter(samples[:, 0], samples[:, 1], alpha=0.5)
plt.title("Generated Data from GMM")
plt.show()

4.4 パラメータ推定と課題

4.4.1 GMMのパラメータ推定

GMMのパラメータ推定は以下の3つを求めることを指します：

$\phi_k$: 各分布の重み
$\mu_k$: 各分布の平均
$\Sigma_k$: 各分布の共分散行列

観測データ $D = {x^{(1)}, x^{(2)}, \dots, x^{(N)}}$ が与えられた場合、GMMの尤度は次のように定義されます：

$p(D; \theta) = \prod_{n=1}^N \sum_{k=1}^K \phi_k \mathcal{N}(x^{(n)}; \mu_k, \Sigma_k)$

ここで、$\theta = (\phi, \mu, \Sigma)$ はGMMのパラメータをまとめたものです。

4.4.2 最尤推定の課題

尤度の対数を最大化する問題：

$\log p(D; \theta) = \sum_{n=1}^N \log \left( \sum_{k=1}^K \phi_k \mathcal{N}(x^{(n)}; \mu_k, \Sigma_k) \right)$

この式は「log-sum」形を持つため解析的に解けません。

4.4.3 EMアルゴリズム

GMMのパラメータ推定には、EMアルゴリズム（Expectation-Maximization Algorithm）が用いられます。

Eステップ: 各データ点がどの分布に属するかの確率を計算
Mステップ: その確率に基づいてパラメータを更新

<第5章以降は執筆中です>

おわりに

本記事は途中式や実装などをかなり省略し、要点のみを押さえた内容となっていますので、省略された部分が気になる方はご自身の理解のためにも本書を読むことをお勧めします。
本書自体は復習などが頻繁に行われており、統計学などの基礎知識が少ない私でも読み進めることが可能な内容でした。

今後は他のシリーズにも挑戦してみようと思います。

ではこれにて。