インクを水に一滴落とすと、インクは徐々に広がり、最終的に水全体に均一に拡散します。この過程は不可逆に見えます — 拡散したインクを元の一滴に戻すことは、通常できません。しかし、もし拡散の各ステップを「逆再生」できたらどうでしょうか。均一に広がったインクから、元の美しい一滴の形を復元できるはずです。
この「拡散過程の逆転」というアイデアが、拡散モデル（Diffusion Model） の核心です。データに段階的にノイズを加えて純粋なノイズにする前方過程と、ノイズからデータを段階的に復元する逆過程を学習することで、高品質なデータ生成を実現します。
2020年にホーらが提案したDDPM（Denoising Diffusion Probabilistic Model） は、拡散モデルの理論を洗練させ、GANに匹敵する画像生成品質を達成しました。現在ではStable DiffusionやDALL-E 2などの画像生成AIの基盤技術として、産業応用でも広く使われています。
拡散モデルを理解すると、以下のような応用が開けます。
- 画像生成: テキストからの画像生成（Text-to-Image）の基盤技術
- 画像編集: 画像の一部を修復するインペインティング、解像度向上
- 音声合成: WaveGradなど、音声波形の高品質生成
- 3Dモデル生成: Point-Eなど、3Dオブジェクトの生成
本記事の内容
- 拡散モデルの基本的な発想と直感
- 前方拡散過程の数学的定式化
- 逆拡散過程と変分下界（ELBO）の導出
- ノイズ予測としての損失関数の簡略化
- PyTorchによる2次元データ生成の実装
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
- 正規分布 — ガウス分布の性質
- KLダイバージェンス — 確率分布間の距離
- VAEの実装チュートリアル — 変分推論とELBOの基礎
- マルコフ連鎖 — 状態遷移の確率過程
拡散モデルとは — ノイズの追加と除去
基本的な発想
拡散モデルの発想は驚くほどシンプルです。
前方過程（Forward process）: きれいなデータ $\bm{x}_0$ に少しずつガウスノイズを加えていき、$T$ ステップ後には純粋なガウスノイズ $\bm{x}_T \sim \mathcal{N}(\bm{0}, \bm{I})$ にします。この過程は固定された（学習しない）マルコフ連鎖です。
逆過程（Reverse process）: 純粋なノイズ $\bm{x}_T$ から始めて、少しずつノイズを除去していき、最終的にきれいなデータ $\bm{x}_0$ を復元します。この過程はニューラルネットワークでパラメータ化され、学習によって獲得します。
日常的な例えで考えましょう。精密なジグソーパズルの完成図（データ）があります。ピースを1つずつランダムな位置に動かしていくと（前方過程）、最終的にはピースがバラバラになります（ノイズ）。逆過程は、バラバラのピースから元の完成図を復元する手順を学習することに相当します。
重要なのは、各ステップで加えるノイズの量が小さければ、逆過程の各ステップも「ほぼガウス的」になるということです。これにより、逆過程をガウス分布でパラメータ化できます。
VAEとの関係
拡散モデルは、実はVAEの一般化として理解できます。VAEは1ステップでデータを潜在空間にエンコードしますが、拡散モデルは $T$ ステップの階層的な潜在変数 $\bm{x}_1, \bm{x}_2, \ldots, \bm{x}_T$ を持つ階層VAEと見なせます。
ただし、通常のVAEとの大きな違いが2つあります。
- エンコーダ（前方過程）が固定されている: VAEではエンコーダも学習しますが、拡散モデルでは前方過程は事前に決められたノイズスケジュールに従う固定された過程です
- 潜在変数がデータと同じ次元: VAEは通常、データよりも低い次元の潜在空間を使いますが、拡散モデルの潜在変数 $\bm{x}_t$ はデータ $\bm{x}_0$ と同じ次元を持ちます
拡散モデルの直感的な理解を得たところで、次に前方過程の数学的定式化に進みましょう。
前方拡散過程
ノイズスケジュール
前方過程は、データ $\bm{x}_0 \sim q(\bm{x}_0)$ から始めて、各ステップ $t$ でガウスノイズを加えるマルコフ連鎖です。
$$ \begin{equation} q(\bm{x}_t | \bm{x}_{t-1}) = \mathcal{N}(\bm{x}_t; \sqrt{1 – \beta_t} \, \bm{x}_{t-1}, \beta_t \bm{I}) \end{equation} $$
ここで $\beta_t \in (0, 1)$ は各ステップのノイズの強さを制御するノイズスケジュールです。$\beta_t$ が小さいほど、各ステップで加えるノイズが少なく、遷移が緩やかになります。
各遷移の式を言葉で説明すると、$\bm{x}_t$ は $\bm{x}_{t-1}$ を $\sqrt{1-\beta_t}$ 倍に縮小し（信号を少し弱め）、分散 $\beta_t$ のノイズを加えたものです。$\sqrt{1-\beta_t}$ で縮小する理由は、ステップを重ねても全体の分散が発散しないようにするためです。
DDPMでは $\beta_1 = 10^{-4}$ から $\beta_T = 0.02$ まで線形に増加させる線形スケジュールが使われます。$T = 1000$ ステップが典型的な設定です。
任意のステップへの直接サンプリング
前方過程の大きな利点は、$\bm{x}_0$ から任意のステップ $t$ の $\bm{x}_t$ を中間ステップを経由せずに直接サンプリングできることです。
$\alpha_t = 1 – \beta_t$ および $\bar{\alpha}_t = \prod_{s=1}^{t} \alpha_s$ とおくと
$$ \begin{equation} q(\bm{x}_t | \bm{x}_0) = \mathcal{N}(\bm{x}_t; \sqrt{\bar{\alpha}_t} \, \bm{x}_0, (1 – \bar{\alpha}_t) \bm{I}) \end{equation} $$
すなわち
$$ \bm{x}_t = \sqrt{\bar{\alpha}_t} \, \bm{x}_0 + \sqrt{1 – \bar{\alpha}_t} \, \bm{\epsilon}, \quad \bm{\epsilon} \sim \mathcal{N}(\bm{0}, \bm{I}) $$
この式の導出を確認しましょう。$\bm{x}_1 = \sqrt{\alpha_1} \, \bm{x}_0 + \sqrt{\beta_1} \, \bm{\epsilon}_1$ から始めます。
$\bm{x}_2$ は $\bm{x}_1$ に基づいて
$$ \bm{x}_2 = \sqrt{\alpha_2} \, \bm{x}_1 + \sqrt{\beta_2} \, \bm{\epsilon}_2 $$
$\bm{x}_1$ を代入すると
$$ \bm{x}_2 = \sqrt{\alpha_2} (\sqrt{\alpha_1} \, \bm{x}_0 + \sqrt{\beta_1} \, \bm{\epsilon}_1) + \sqrt{\beta_2} \, \bm{\epsilon}_2 $$
$$ = \sqrt{\alpha_1 \alpha_2} \, \bm{x}_0 + \sqrt{\alpha_2 \beta_1} \, \bm{\epsilon}_1 + \sqrt{\beta_2} \, \bm{\epsilon}_2 $$
ここで、独立なガウス確率変数の和の分散は加法的なので、$\sqrt{\alpha_2 \beta_1} \, \bm{\epsilon}_1 + \sqrt{\beta_2} \, \bm{\epsilon}_2$ は $\mathcal{N}(\bm{0}, (\alpha_2 \beta_1 + \beta_2) \bm{I})$ にしたがいます。$\alpha_2 \beta_1 + \beta_2 = \alpha_2(1-\alpha_1) + (1-\alpha_2) = 1 – \alpha_1\alpha_2 = 1 – \bar{\alpha}_2$ を確認すると、確かに $q(\bm{x}_2|\bm{x}_0) = \mathcal{N}(\sqrt{\bar{\alpha}_2}\,\bm{x}_0, (1-\bar{\alpha}_2)\bm{I})$ が得られます。帰納法で一般の $t$ に拡張できます。
この直接サンプリングは、学習時に非常に重要です。ランダムな $t$ を選んで即座に $\bm{x}_t$ を計算できるため、効率的な学習が可能になります。
$t \to T$ で $\bar{\alpha}_T \to 0$ となるようにスケジュールを設計すれば、$\bm{x}_T \approx \mathcal{N}(\bm{0}, \bm{I})$ となり、前方過程はデータを純粋なノイズに変換できます。
前方過程の数学を理解したところで、次にこの過程を「巻き戻す」逆過程の定式化に進みましょう。
逆拡散過程
逆過程の形
逆過程は、ノイズ $\bm{x}_T \sim \mathcal{N}(\bm{0}, \bm{I})$ から始めて、段階的にノイズを除去してデータ $\bm{x}_0$ を復元する過程です。
ホーらの重要な観察は、$\beta_t$ が十分に小さい場合、$q(\bm{x}_{t-1} | \bm{x}_t)$ もまたガウス分布に近くなるということです。したがって、逆過程もガウス遷移でパラメータ化できます。
$$ \begin{equation} p_\theta(\bm{x}_{t-1} | \bm{x}_t) = \mathcal{N}(\bm{x}_{t-1}; \bm{\mu}_\theta(\bm{x}_t, t), \sigma_t^2 \bm{I}) \end{equation} $$
ここで $\bm{\mu}_\theta(\bm{x}_t, t)$ はニューラルネットワークで予測する平均、$\sigma_t^2$ は固定された分散です。
逆過程の全体は
$$ p_\theta(\bm{x}_{0:T}) = p(\bm{x}_T) \prod_{t=1}^{T} p_\theta(\bm{x}_{t-1} | \bm{x}_t) $$
条件付き逆過程の解析解
$\bm{x}_0$ が既知の場合、条件付き逆過程 $q(\bm{x}_{t-1} | \bm{x}_t, \bm{x}_0)$ はベイズの定理から解析的に計算できます。
$$ q(\bm{x}_{t-1} | \bm{x}_t, \bm{x}_0) = \frac{q(\bm{x}_t | \bm{x}_{t-1}, \bm{x}_0) \, q(\bm{x}_{t-1} | \bm{x}_0)}{q(\bm{x}_t | \bm{x}_0)} $$
マルコフ性から $q(\bm{x}_t | \bm{x}_{t-1}, \bm{x}_0) = q(\bm{x}_t | \bm{x}_{t-1})$ です。3つのガウス分布の比を計算すると、結果もまたガウス分布になります。
$$ q(\bm{x}_{t-1} | \bm{x}_t, \bm{x}_0) = \mathcal{N}(\bm{x}_{t-1}; \tilde{\bm{\mu}}_t(\bm{x}_t, \bm{x}_0), \tilde{\beta}_t \bm{I}) $$
ここで
$$ \tilde{\bm{\mu}}_t(\bm{x}_t, \bm{x}_0) = \frac{\sqrt{\bar{\alpha}_{t-1}} \beta_t}{1 – \bar{\alpha}_t} \bm{x}_0 + \frac{\sqrt{\alpha_t}(1 – \bar{\alpha}_{t-1})}{1 – \bar{\alpha}_t} \bm{x}_t $$
$$ \begin{equation} \tilde{\beta}_t = \frac{1 – \bar{\alpha}_{t-1}}{1 – \bar{\alpha}_t} \beta_t \end{equation} $$
この平均 $\tilde{\bm{\mu}}_t$ は $\bm{x}_0$ と $\bm{x}_t$ の加重平均です。$\bm{x}_0 = (\bm{x}_t – \sqrt{1-\bar{\alpha}_t}\,\bm{\epsilon}) / \sqrt{\bar{\alpha}_t}$ を代入すると
$$ \begin{equation} \tilde{\bm{\mu}}_t(\bm{x}_t, \bm{\epsilon}) = \frac{1}{\sqrt{\alpha_t}}\left(\bm{x}_t – \frac{\beta_t}{\sqrt{1 – \bar{\alpha}_t}} \bm{\epsilon}\right) \end{equation} $$
この式は、平均 $\tilde{\bm{\mu}}_t$ が $\bm{x}_t$ から「ノイズ成分 $\bm{\epsilon}$ を引き去る」操作であることを示しています。ニューラルネットワークの目標は、$\bm{x}_t$ と $t$ から、前方過程で加えられたノイズ $\bm{\epsilon}$ を予測することになります。
変分下界と損失関数の導出
ELBOの導出
データの対数尤度 $\ln p_\theta(\bm{x}_0)$ に対するELBO（変分下界）を導出します。VAEと同様に
$$ \ln p_\theta(\bm{x}_0) \geq \mathbb{E}_q\left[\ln \frac{p_\theta(\bm{x}_{0:T})}{q(\bm{x}_{1:T}|\bm{x}_0)}\right] = \mathcal{L} $$
このELBOは以下の3種類の項に分解できます。
$$ \mathcal{L} = \underbrace{\mathbb{E}_q[\ln p_\theta(\bm{x}_0 | \bm{x}_1)]}_{L_0} – \underbrace{D_{\text{KL}}(q(\bm{x}_T|\bm{x}_0) \| p(\bm{x}_T))}_{L_T} – \sum_{t=2}^{T} \underbrace{\mathbb{E}_q\left[D_{\text{KL}}(q(\bm{x}_{t-1}|\bm{x}_t, \bm{x}_0) \| p_\theta(\bm{x}_{t-1}|\bm{x}_t))\right]}_{L_{t-1}} $$
各項の意味を説明します。
$L_0$（再構成項）: $\bm{x}_1$ から $\bm{x}_0$ を復元する精度を表します。
$L_T$（事前分布との一致）: 前方過程の最終ステップが事前分布 $p(\bm{x}_T) = \mathcal{N}(\bm{0}, \bm{I})$ と一致するかを表します。ノイズスケジュールが適切であれば、この項はほぼゼロであり、定数として扱えます。
$L_{t-1}$（デノイジングマッチング項）: $t = 2, \ldots, T$ の各ステップで、学習する逆過程 $p_\theta(\bm{x}_{t-1}|\bm{x}_t)$ が条件付き逆過程 $q(\bm{x}_{t-1}|\bm{x}_t, \bm{x}_0)$ にどれだけ近いかを表します。この項が学習の核心です。
損失関数の簡略化
$L_{t-1}$ はガウス分布同士のKLダイバージェンスであり、分散が一致する場合、平均の差の二乗に比例します。
$$ L_{t-1} = \frac{1}{2\tilde{\beta}_t} \|\tilde{\bm{\mu}}_t(\bm{x}_t, \bm{x}_0) – \bm{\mu}_\theta(\bm{x}_t, t)\|^2 + \text{const.} $$
$\bm{\mu}_\theta$ を $\bm{\epsilon}$ の予測としてパラメータ化します。
$$ \bm{\mu}_\theta(\bm{x}_t, t) = \frac{1}{\sqrt{\alpha_t}}\left(\bm{x}_t – \frac{\beta_t}{\sqrt{1 – \bar{\alpha}_t}} \bm{\epsilon}_\theta(\bm{x}_t, t)\right) $$
$\tilde{\bm{\mu}}_t$ の表式と比較すると、$L_{t-1}$ は $\bm{\epsilon}_\theta$ とtrue noise $\bm{\epsilon}$ の差に帰着します。
$$ L_{t-1} = \frac{\beta_t^2}{2\tilde{\beta}_t \alpha_t (1-\bar{\alpha}_t)} \|\bm{\epsilon} – \bm{\epsilon}_\theta(\bm{x}_t, t)\|^2 + \text{const.} $$
DDPMでは、時刻依存の係数を無視した簡略化された損失関数が使われます。
$$ \begin{equation} L_{\text{simple}} = \mathbb{E}_{t, \bm{x}_0, \bm{\epsilon}}\left[\|\bm{\epsilon} – \bm{\epsilon}_\theta(\sqrt{\bar{\alpha}_t} \, \bm{x}_0 + \sqrt{1-\bar{\alpha}_t} \, \bm{\epsilon}, t)\|^2\right] \end{equation} $$
この損失関数は驚くほどシンプルです。「ノイズが加えられた画像からノイズを予測する」デノイジングタスクに帰着しているのです。係数を無視することで各タイムステップの重みが均一になり、実験的にはこの簡略化が生成品質を向上させることがわかっています。
損失関数の理論的導出を理解したところで、次に学習とサンプリングのアルゴリズムを具体的に見ていきましょう。
学習とサンプリングのアルゴリズム
学習アルゴリズム
DDPMの学習は以下の手順で行います。
- データ $\bm{x}_0$ をミニバッチからサンプリング
- タイムステップ $t \sim \text{Uniform}(\{1, \ldots, T\})$ をランダムに選択
- ノイズ $\bm{\epsilon} \sim \mathcal{N}(\bm{0}, \bm{I})$ をサンプリング
- ノイズ付きデータ $\bm{x}_t = \sqrt{\bar{\alpha}_t} \, \bm{x}_0 + \sqrt{1-\bar{\alpha}_t} \, \bm{\epsilon}$ を計算
- ノイズ予測ネットワーク $\bm{\epsilon}_\theta(\bm{x}_t, t)$ の出力と $\bm{\epsilon}$ のMSEを計算
- 勾配降下法でパラメータ $\theta$ を更新
各ステップで異なる $t$ が選ばれるため、ネットワークは全てのノイズレベルでのデノイジングを同時に学習します。
サンプリングアルゴリズム
学習済みモデルからサンプルを生成する手順は以下のとおりです。
- $\bm{x}_T \sim \mathcal{N}(\bm{0}, \bm{I})$ をサンプリング
- $t = T, T-1, \ldots, 1$ について以下を繰り返す: – $\bm{z} \sim \mathcal{N}(\bm{0}, \bm{I})$ をサンプリング（$t > 1$ のとき）、$\bm{z} = \bm{0}$（$t = 1$ のとき） – $\bm{x}_{t-1} = \frac{1}{\sqrt{\alpha_t}}\left(\bm{x}_t – \frac{\beta_t}{\sqrt{1-\bar{\alpha}_t}} \bm{\epsilon}_\theta(\bm{x}_t, t)\right) + \sigma_t \bm{z}$
- $\bm{x}_0$ を返す
このサンプリングは $T$ 回の逐次的なステップを必要とするため、GANの1回のforward passと比較すると非常に遅いです。DDPMでは $T = 1000$ が典型的であり、1枚の画像生成に1000回のニューラルネットワーク推論が必要です。この遅さを改善する方法（DDIMなど）が後続の研究で提案されています。
PyTorchによる実装 — 2次元データの拡散モデル
2次元での実装
拡散モデルの動作を理解するために、2次元のガウス混合分布を対象に実装します。
import torch
import torch.nn as nn
import torch.optim as optim
import numpy as np
import matplotlib.pyplot as plt
torch.manual_seed(42)
np.random.seed(42)
# --- 対象データ: 2次元スイスロール ---
def sample_swiss_roll(n_samples):
"""スイスロール型のデータを生成"""
t = 1.5 * np.pi * (1 + 2 * np.random.rand(n_samples))
x = t * np.cos(t)
y = t * np.sin(t)
data = np.stack([x, y], axis=1)
data = data / 10.0 # スケーリング
return torch.tensor(data, dtype=torch.float32)
# --- ノイズスケジュール ---
T = 300 # 拡散ステップ数
beta_start = 1e-4
beta_end = 0.02
betas = torch.linspace(beta_start, beta_end, T)
alphas = 1.0 - betas
alpha_bars = torch.cumprod(alphas, dim=0)
def q_sample(x_0, t, noise=None):
"""前方過程: x_0からx_tを直接サンプリング"""
if noise is None:
noise = torch.randn_like(x_0)
sqrt_alpha_bar = torch.sqrt(alpha_bars[t]).unsqueeze(-1)
sqrt_one_minus_alpha_bar = torch.sqrt(1 - alpha_bars[t]).unsqueeze(-1)
return sqrt_alpha_bar * x_0 + sqrt_one_minus_alpha_bar * noise
ここでは2次元のスイスロールデータを対象としています。スイスロールは渦巻き状の非線形な構造を持つため、モデルが複雑な分布を学習できるかのテストに適しています。ノイズスケジュールは $\beta_1 = 10^{-4}$ から $\beta_T = 0.02$ まで線形に増加する設定です。alpha_bars は $\bar{\alpha}_t = \prod_{s=1}^t \alpha_s$ を cumprod（累積積）で計算しています。
# --- ノイズ予測ネットワーク ---
class NoisePredictor(nn.Module):
def __init__(self, data_dim=2, hidden_dim=256, time_emb_dim=64):
super().__init__()
# 時刻の埋め込み
self.time_emb = nn.Sequential(
nn.Linear(1, time_emb_dim),
nn.SiLU(),
nn.Linear(time_emb_dim, time_emb_dim),
)
# メインネットワーク
self.net = nn.Sequential(
nn.Linear(data_dim + time_emb_dim, hidden_dim),
nn.SiLU(),
nn.Linear(hidden_dim, hidden_dim),
nn.SiLU(),
nn.Linear(hidden_dim, hidden_dim),
nn.SiLU(),
nn.Linear(hidden_dim, data_dim),
)
def forward(self, x, t):
# t を正規化して埋め込み
t_normalized = t.float().unsqueeze(-1) / T
t_emb = self.time_emb(t_normalized)
# x と時刻埋め込みを結合
x_input = torch.cat([x, t_emb], dim=-1)
return self.net(x_input)
ノイズ予測ネットワークは、ノイズが加えられたデータ $\bm{x}_t$ と時刻 $t$ を入力として受け取り、加えられたノイズ $\bm{\epsilon}$ を予測します。時刻 $t$ は専用の埋め込み層で処理されます。これは、異なるノイズレベル（異なる $t$）で異なるデノイジング戦略が必要であるためです。ノイズがわずかな場合（$t$ が小さい）と、ほぼ純粋なノイズの場合（$t$ が大きい）では、ネットワークの振る舞いが変わるべきです。
活性化関数にはSiLU（Sigmoid Linear Unit, $x \cdot \sigma(x)$）を使用しています。SiLUはReLUよりも滑らかであり、拡散モデルの学習で広く使われています。
# --- 学習 ---
model = NoisePredictor()
optimizer = optim.Adam(model.parameters(), lr=3e-4)
n_epochs = 5000
batch_size = 1024
losses = []
for epoch in range(n_epochs):
x_0 = sample_swiss_roll(batch_size)
# ランダムな時刻を選択
t = torch.randint(0, T, (batch_size,))
# ノイズの生成と前方過程
noise = torch.randn_like(x_0)
x_t = q_sample(x_0, t, noise)
# ノイズの予測
noise_pred = model(x_t, t)
# 損失（簡略化されたELBO）
loss = nn.functional.mse_loss(noise_pred, noise)
optimizer.zero_grad()
loss.backward()
optimizer.step()
losses.append(loss.item())
if (epoch + 1) % 1000 == 0:
print(f'Epoch [{epoch+1}/{n_epochs}] Loss: {loss.item():.6f}')
学習ループは理論で導出した簡略化された損失関数をそのまま実装しています。各イテレーションで、ランダムなデータ $\bm{x}_0$、ランダムな時刻 $t$、ランダムなノイズ $\bm{\epsilon}$ を生成し、ノイズ予測のMSEを計算します。
# --- サンプリング ---
@torch.no_grad()
def p_sample(model, x_t, t_idx):
"""逆過程の1ステップ"""
t = torch.full((x_t.shape[0],), t_idx, dtype=torch.long)
beta_t = betas[t_idx]
alpha_t = alphas[t_idx]
alpha_bar_t = alpha_bars[t_idx]
# ノイズ予測
eps_pred = model(x_t, t)
# 平均の計算
mu = (1.0 / torch.sqrt(alpha_t)) * (
x_t - (beta_t / torch.sqrt(1 - alpha_bar_t)) * eps_pred
)
if t_idx > 0:
noise = torch.randn_like(x_t)
sigma = torch.sqrt(beta_t)
return mu + sigma * noise
else:
return mu
@torch.no_grad()
def sample(model, n_samples):
"""逆拡散過程でサンプルを生成"""
x = torch.randn(n_samples, 2)
trajectory = [x.numpy().copy()]
for t in reversed(range(T)):
x = p_sample(model, x, t)
if t % 30 == 0:
trajectory.append(x.numpy().copy())
return x.numpy(), trajectory
サンプリング関数は逆拡散過程を実装しています。$\bm{x}_T \sim \mathcal{N}(\bm{0}, \bm{I})$ から始めて、$t = T-1, T-2, \ldots, 0$ の順にデノイジングステップを繰り返します。trajectory に途中の状態を記録することで、サンプリング過程の可視化が可能です。
# --- 結果の可視化 ---
generated, trajectory = sample(model, 2000)
fig, axes = plt.subplots(2, 4, figsize=(18, 9))
# 上段: 拡散過程の可視化
real_data = sample_swiss_roll(2000).numpy()
noise_levels = [0, T//4, T//2, 3*T//4]
for idx, t_val in enumerate(noise_levels):
ax = axes[0, idx]
t_tensor = torch.full((2000,), t_val, dtype=torch.long)
x_0_tensor = torch.tensor(real_data, dtype=torch.float32)
x_t = q_sample(x_0_tensor, t_tensor).numpy()
ax.scatter(x_t[:, 0], x_t[:, 1], s=2, alpha=0.5, c='steelblue')
ax.set_title(f'Forward: t = {t_val}', fontsize=12)
ax.set_xlim(-3, 3)
ax.set_ylim(-3, 3)
ax.set_aspect('equal')
ax.grid(True, alpha=0.3)
# 下段: 逆拡散過程の可視化
n_traj = len(trajectory)
traj_indices = np.linspace(0, n_traj-1, 4, dtype=int)
for idx, traj_idx in enumerate(traj_indices):
ax = axes[1, idx]
traj_data = trajectory[traj_idx]
ax.scatter(traj_data[:, 0], traj_data[:, 1], s=2, alpha=0.5,
c='coral')
step = T - traj_idx * 30 if traj_idx < n_traj-1 else 0
ax.set_title(f'Reverse: step {traj_idx}/{n_traj-1}', fontsize=12)
ax.set_xlim(-3, 3)
ax.set_ylim(-3, 3)
ax.set_aspect('equal')
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig('diffusion_process.png', dpi=150, bbox_inches='tight')
plt.show()
この可視化から、拡散モデルの前方過程と逆過程の両方を確認できます。
-
上段（前方拡散過程）: 左端のスイスロール構造のデータが、$t$ の増加とともに徐々にノイズ化されていく様子です。$t = 0$ ではスイスロールの渦巻き構造が明確ですが、$t = T/4$ で構造が崩れ始め、$t = 3T/4$ ではほぼガウスノイズに近づいています。これは理論どおり、$\bar{\alpha}_t \to 0$ でデータ成分が消失する過程です
-
下段（逆拡散過程）: ガウスノイズから始めて、学習済みモデルで段階的にデノイジングした結果です。初期のランダムな点群が、ステップを経るごとに構造を獲得し、最終的にスイスロールの形を復元しています。逆過程が前方過程の「巻き戻し」として正しく機能していることが確認できます
# --- 最終結果の比較 ---
fig, axes = plt.subplots(1, 3, figsize=(16, 5))
# (a) 元データ
ax = axes[0]
ax.scatter(real_data[:, 0], real_data[:, 1], s=3, alpha=0.5, c='blue')
ax.set_title('Original Data', fontsize=13)
ax.set_xlim(-3, 3)
ax.set_ylim(-3, 3)
ax.set_aspect('equal')
ax.grid(True, alpha=0.3)
# (b) 生成データ
ax = axes[1]
ax.scatter(generated[:, 0], generated[:, 1], s=3, alpha=0.5, c='red')
ax.set_title('Generated Data', fontsize=13)
ax.set_xlim(-3, 3)
ax.set_ylim(-3, 3)
ax.set_aspect('equal')
ax.grid(True, alpha=0.3)
# (c) 損失の推移
ax = axes[2]
window = 100
smoothed = np.convolve(losses, np.ones(window)/window, mode='valid')
ax.plot(smoothed, 'b-', linewidth=1.5)
ax.set_xlabel('Epoch', fontsize=12)
ax.set_ylabel('MSE Loss', fontsize=12)
ax.set_title('Training Loss', fontsize=13)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.savefig('diffusion_result.png', dpi=150, bbox_inches='tight')
plt.show()
最終結果の比較から、拡散モデルの生成能力を確認できます。
-
左図（元データ）: スイスロール分布の特徴的な渦巻き構造です
-
中央図（生成データ）: 拡散モデルで生成したデータです。スイスロールの渦巻き構造が再現されており、密度の分布も元データに近いことがわかります。GANと異なりモード崩壊の心配がなく、分布全体を均一にカバーしています
-
右図（損失の推移）: MSE損失が学習を通じて安定して減少しています。VAEと同様に、拡散モデルの損失は学習の進行を反映するため、学習状況の監視に有用です。GANのようなゲームベースの不安定性がない点が特徴的です
拡散モデルの理論的考察
なぜ拡散モデルは高品質なのか
拡散モデルがGANやVAEと比較して高品質な生成を実現できる理由は複数あります。
段階的な生成: 1ステップで画像全体を生成するGANやVAEと異なり、拡散モデルは $T$ ステップかけて徐々にデータを形成します。各ステップでの変化量が小さいため、ネットワークが学習すべきタスクが容易です。
全モードのカバー: 拡散モデルの目的関数は対数尤度のELBOであり、データ分布のモードを漏らすことなく学習します。GANのモード崩壊の問題がありません。
安定した学習: 単一の目的関数（ノイズ予測のMSE）を最小化するだけであり、GANのような2つのネットワークの均衡問題がありません。
サンプリング速度の問題
拡散モデルの最大の弱点は、サンプリングに $T$ ステップを要することです。DDPMでは $T = 1000$ であり、GANの約1000倍の計算が必要です。
この問題に対する改善手法として以下が提案されています。
- DDIM（Denoising Diffusion Implicit Models）: 非マルコフな逆過程を使い、ステップ数を大幅に削減（例: 50ステップ）
- Progressive Distillation: 学習済みモデルを蒸留してステップ数を半減させる操作を繰り返す
- Consistency Models: 任意のノイズレベルから直接データを予測するモデル
まとめ
本記事では、拡散モデル（DDPM）の数理を完全に導出し、実装しました。
- 前方過程はデータに段階的にガウスノイズを加え、$\bar{\alpha}_t$ を用いた直接サンプリングが可能
- 逆過程はノイズ予測ネットワーク $\bm{\epsilon}_\theta$ でパラメータ化され、条件付き逆過程とのKLダイバージェンスを最小化する
- 損失関数は $\|\bm{\epsilon} – \bm{\epsilon}_\theta(\bm{x}_t, t)\|^2$ という単純なノイズ予測のMSEに簡略化される
- 安定した学習、モード崩壊のない生成、尤度ベースの目的関数がGANに対する利点
- サンプリング速度の遅さが主な弱点であり、DDIMなどの加速手法で改善される
次のステップとして、以下の記事も参考にしてください。
- スコアマッチングとスコアベース生成モデル — 拡散モデルのスコア関数の視点からの理解
- 画像生成モデルの比較 — GAN, VAE, Diffusionの体系的比較
- VAEの実装チュートリアル — 変分推論の基礎