目次
序論
10.1 条件独立性
基本概念
コードによる条件独立性の検証
実行結果
10.2 有向グラフモデル
基本概念
10.2.1 例1:単純な有向グラフ(雨→外出→傘)
10.2.2 例2:コンピュータビジョン応用 —— 画像分類における有向グラフ
10.2.3 例3:有向グラフによる推論(結果から原因を推測)
10.2.4 まとめ
10.3 無向グラフモデル
基本概念
10.3.1 例1:単純な無向グラフ(画像ピクセル間依存関係)
10.3.2 例2:コンピュータビジョン応用 —— 画像ノイズ除去(無向グラフ)
10.4 有向グラフモデルと無向グラフモデルの比較
主要な比較(マインドマップ)
視覚的比較コード
10.5 コンピュータビジョンにおけるグラフモデル
主な応用分野
実装例:CRFによる画像セマンティックセグメンテーション(無向グラフの拡張)
10.6 複数未知変数を持つモデルの推論
基本概念
10.6.1 最大事後確率解の求め方
10.6.2 後続分布の周辺分布の計算
10.6.3 周辺分布の最大化
10.6.4 後続分布からのサンプリング
10.7 サンプル生成
10.7.1 有向グラフモデルからのサンプル生成
10.7.2 無向グラフモデルからのサンプル生成
10.8 学習
基本概念
10.8.1 有向グラフモデルの学習
10.8.2 無向グラフモデルの学習
議論
備考
演習問題
まとめ
序論
こんにちは!今日は『コンピュータビジョン:モデル、学習と推論』という書籍の第10章「グラフモデル」について解説します。グラフモデルは不確実性を扱う際の強力なツールであり、複雑な確率関係を直感的に表現し、変数間の依存関係を明確にします。この記事では煩雑な数式を避け、分かりやすい言葉+実行可能なPythonコード+視覚化図表を通じて、グラフモデルの核となる概念と応用を深く理解していただきます。
10.1 条件独立性
基本概念
条件独立性はグラフモデルの中心的な概念です。簡単に言うと、「ある条件Cが与えられたとき、AとBは互いに影響を与えない」という意味です。
たとえば、あなたが外出するかどうか(A)と、スーパーマーケットが割引になるかどうか(B)は通常関連していますが、もし「大雨が降っている(C)」という条件が与えられれば、スーパーマーケットが割引になっても外出しません。このとき、AとBはCの条件下で独立します。
数学的には、条件独立は P(A,B|C) = P(A|C) × P(B|C) と表されます。
条件独立性のコード検証
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import multivariate_normal
# MacシステムのMatplotlibでの日本語表示設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS', 'DejaVu Sans']
plt.rcParams['axes.unicode_minus'] = False
plt.rcParams['font.family'] = 'sans-serif'
plt.rcParams['font.family'] = 'Arial Unicode MS'
plt.rcParams['axes.facecolor'] = 'white'
# 3つの変数:A(外出), B(割引), C(雨)
# データ生成:C=1は雨、C=0は晴れ
np.random.seed(42)
n_samples = 1000
# C(雨)のデータ生成(30%の確率)
C = np.random.binomial(1, 0.3, n_samples)
# Aの生成:晴れの日は80%外出、雨の日は5%外出
A = np.where(C == 0,
np.random.binomial(1, 0.8, n_samples),
np.random.binomial(1, 0.05, n_samples))
# Bの生成:割引確率50%(Aとは無関係)
B = np.random.binomial(1, 0.5, n_samples)
# 条件独立性の検証:P(A=1|B=1,C=1) と P(A=1|C=1) を計算
# 1. C=1かつB=1のときのA=1の確率
mask_c1 = C == 1
mask_c1_b1 = (C == 1) & (B == 1)
p_a1_c1_b1 = np.sum(A[mask_c1_b1] == 1) / len(A[mask_c1_b1])
p_a1_c1 = np.sum(A[mask_c1] == 1) / len(A[mask_c1])
# グラフ表示
fig, ax = plt.subplots(figsize=(8, 5))
categories = ['P(A=1|C=1,B=1)', 'P(A=1|C=1)']
values = [p_a1_c1_b1, p_a1_c1]
bars = ax.bar(categories, values, color=['#1f77b4', '#ff7f0e'])
# 数値ラベル追加
for bar, val in zip(bars, values):
ax.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.01,
f'{val:.3f}', ha='center', va='bottom')
ax.set_ylabel('確率値')
ax.set_title('条件独立性の検証:雨の日において、BとAは独立')
ax.set_ylim(0, 0.1)
plt.show()
# 結果出力
print(f"P(A=1|C=1,B=1) = {p_a1_c1_b1:.3f}")
print(f"P(A=1|C=1) = {p_a1_c1:.3f}")
print("両者はほぼ一致し、条件独立性が確認されました!")
実行結果
このコードは柱状グラフを生成し、雨の日における割引の有無が外出確率に与える影響がほぼないことを視覚的に示します。条件独立性が成立していることが確認できます。
10.2 有向グラフモデル
基本概念
有向グラフモデル(ベイズネットワーク)は変数間の依存関係を矢印付きの辺で表し、矢印の方向は「影響される側」を示します。例えば「雨→外出」は、雨が外出に影響を与えることを意味します。
有向グラフモデルは「家族系譜図」のようなものと考えられます:親ノードが「原因」、子ノードが「結果」で、子ノードはその親ノードのみに依存します。
10.2.1 例1:単純な有向グラフ(雨→外出→傘)
import networkx as nx
import matplotlib.pyplot as plt
# Macフォント設定(重複設定で確実に有効化)
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 有向グラフの作成
G = nx.DiGraph()
# ノードと辺の追加
nodes = ['雨', '外出', '傘']
edges = [('雨', '外出'), ('外出', '傘'), ('雨', '傘')]
G.add_nodes_from(nodes)
G.add_edges_from(edges)
# 有向グラフの描画
fig, ax = plt.subplots(figsize=(8, 5))
pos = nx.spring_layout(G, seed=42) # レイアウト固定
nx.draw(G, pos, with_labels=True, ax=ax,
node_color='#1f77b4', node_size=2000,
font_size=12, font_family='Arial Unicode MS',
arrowstyle='->', arrowsize=20)
ax.set_title('有向グラフモデル例:雨→外出→傘', fontsize=14)
plt.show()
# 結合確率の計算:P(雨,外出,傘) = P(雨) × P(外出|雨) × P(傘|雨,外出)
# 条件確率表の定義
p_rain = 0.3 # P(雨=1)
p_go_out_rain = 0.05 # P(外出=1|雨=1)
p_go_out_no_rain = 0.8 # P(外出=1|雨=0)
p_umbrella_rain_go = 0.95 # P(傘=1|雨=1,外出=1)
p_umbrella_rain_no_go = 0.1 # P(傘=1|雨=1,外出=0)
p_umbrella_no_rain_go = 0.1 # P(傘=1|雨=0,外出=1)
p_umbrella_no_rain_no_go = 0.01 # P(傘=1|雨=0,外出=0)
# 結合確率の計算:雨=1,外出=1,傘=1
joint_prob = p_rain * p_go_out_rain * p_umbrella_rain_go
print(f"P(雨=1,外出=1,傘=1) = {joint_prob:.4f}")
10.2.2 例2:コンピュータビジョン応用 —— 画像分類における有向グラフ
import numpy as np
import matplotlib.pyplot as plt
from sklearn.datasets import load_digits
from sklearn.naive_bayes import GaussianNB
from sklearn.model_selection import train_test_split
from sklearn.metrics import accuracy_score
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 手書き数字データセットの読み込み(有向グラフモデルの古典的な応用:ナイーブベイズ)
digits = load_digits()
X, y = digits.data, digits.target
# データセットの分割
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=42)
# ナイーブベイズ(有向グラフモデルの一種:カテゴリ→特徴、特徴間は条件独立)
gnb = GaussianNB()
gnb.fit(X_train, y_train)
y_pred = gnb.predict(X_test)
# 正確率の計算
accuracy = accuracy_score(y_test, y_pred)
print(f"ナイーブベイズ(有向グラフモデル)分類精度:{accuracy:.4f}")
# 元画像と予測結果の可視化
fig, axes = plt.subplots(2, 5, figsize=(12, 6))
axes = axes.flatten()
# テストサンプルの前10個を表示
for i in range(10):
ax = axes[i]
ax.imshow(X_test[i].reshape(8, 8), cmap='gray')
ax.set_title(f"正解:{y_test[i]}\n予測:{y_pred[i]}")
ax.axis('off')
plt.suptitle('有向グラフモデル(ナイーブベイズ)手書き数字分類結果', fontsize=14)
plt.tight_layout()
plt.show()
10.2.3 例3:有向グラフによる推論(結果から原因を推測)
import numpy as np
import matplotlib.pyplot as plt
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 有向グラフ:風邪→発熱、インフルエンザ→発熱
# 先験確率と条件確率の定義
p_cold = 0.2 # 風邪の先験確率
p_flu = 0.1 # インフルエンザの先験確率
p_fever_cold = 0.9 # 風邪時に発熱する確率
p_fever_flu = 0.95 # インフルエンザ時に発熱する確率
p_fever_no = 0.05 # 風邪もインフルエンザもない時に発熱する確率
# 発熱(結果)が観測された場合、風邪とインフルエンザの確率を推論(原因)
# ベイズの定理を使用:P(風邪|発熱) = P(発熱|風邪)P(風邪) / P(発熱)
# P(発熱) の計算:P(発熱|風邪)P(風邪) + P(発熱|インフルエンザ)P(インフルエンザ) + P(発熱|無)P(無)
p_fever = (p_fever_cold * p_cold) + (p_fever_flu * p_flu) + (p_fever_no * (1 - p_cold - p_flu))
p_cold_fever = (p_fever_cold * p_cold) / p_fever
p_flu_fever = (p_fever_flu * p_flu) / p_fever
# 推論結果の可視化
fig, ax = plt.subplots(figsize=(8, 5))
categories = ['風邪', 'インフルエンザ']
probabilities = [p_cold_fever, p_flu_fever]
bars = ax.bar(categories, probabilities, color=['#2ca02c', '#d62728'])
# 数値ラベル追加
for bar, prob in zip(bars, probabilities):
ax.text(bar.get_x() + bar.get_width()/2, bar.get_height() + 0.01,
f'{prob:.3f}', ha='center', va='bottom')
ax.set_ylabel('事後確率')
ax.set_title('発熱が観測された場合の風邪・インフルエンザ確率推論(有向グラフ推論)')
ax.set_ylim(0, 1)
plt.show()
print(f"P(風邪|発熱) = {p_cold_fever:.3f}")
print(f"P(インフルエンザ|発熱) = {p_flu_fever:.3f}")
10.2.4 まとめ
有向グラフモデルの主な特徴:
- 辺は方向があり、「因果関係」(親→子)を示す;
- 結合確率は各ノードの親ノード条件付き確率の積に分解できる;
- 明確な因果関係がある場面(診断、分類)に適している。
10.3 無向グラフモデル
基本概念
無向グラフモデル(マルコフ随機場)は変数間の「相互影響」を矢印なしの辺で表現し、因果関係はなく、より平等な関係を示します。例えば「画像の隣接ピクセル」は、左のピクセルが明るければ右のピクセルも明るい傾向があり、どちらが原因かは明確ではありません。このような関係には無向グラフが適しています。
10.3.1 例1:単純な無向グラフ(画像ピクセル依存関係)
import numpy as np
import matplotlib.pyplot as plt
import networkx as nx
# Macフォント設定(日本語とマイナス記号の正常表示を保証)
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# ===================== 1. 3x3ピクセルの無向グラフの描画 =====================
G = nx.Graph()
# ピクセルノードの追加(番号0-8、3x3グリッドに対応)
nodes = [f'ピクセル({i},{j})' for i in range(3) for j in range(3)]
G.add_nodes_from(nodes)
# 隣接エッジの追加(上下左右)
edges = []
for i in range(3):
for j in range(3):
# 右のピクセル(同じ行、列+1)
if j < 2:
edges.append((f'ピクセル({i},{j})', f'ピクセル({i},{j+1})'))
# 下のピクセル(同じ列、行+1)
if i < 2:
edges.append((f'ピクセル({i},{j})', f'ピクセル({i+1},{j})'))
G.add_edges_from(edges)
# 無向グラフの描画
fig, ax = plt.subplots(figsize=(10, 8))
pos = nx.spring_layout(G, seed=42) # 固定シードでレイアウトを安定化
nx.draw(G, pos, with_labels=True, ax=ax,
node_color='#ff7f0e', node_size=3000,
font_size=10, font_family='Arial Unicode MS',
edge_color='gray')
ax.set_title('無向グラフモデル例:3x3画像ピクセルの隣接依存関係', fontsize=14)
plt.show()
# ===================== 2. ピクセル依存関係の検証(類似度計算で代替) =====================
np.random.seed(42)
img = np.random.randint(0, 255, (3, 3)) # 3x3ランダムグレースケール画像
print(f"3x3画像のピクセル値:\n{img}")
# 隣接ピクセル間の類似度を計算(相関係数ではなく、単一サンプル向け)
# 類似度式:1 - |ピクセル1 - ピクセル2| / 255 (類似度が1に近いほど依存度が高い)
similarity = []
similarity_details = [] # 詳細情報を格納
for i in range(3):
for j in range(3):
# 右隣ピクセルとの類似度
if j < 2:
pix1 = img[i, j]
pix2 = img[i, j+1]
sim = 1 - abs(pix1 - pix2) / 255.0
similarity.append(sim)
similarity_details.append((f'({i},{j})-({i},{j+1})', pix1, pix2, sim))
# 下隣ピクセルとの類似度
if i < 2:
pix1 = img[i, j]
pix2 = img[i+1, j]
sim = 1 - abs(pix1 - pix2) / 255.0
similarity.append(sim)
similarity_details.append((f'({i},{j})-({i+1},{j})', pix1, pix2, sim))
# 各ピクセル対の類似度を出力
print("\n隣接ピクセル対の類似度(類似度が1に近いほど依存度が高い):")
for pair, p1, p2, sim in similarity_details:
print(f"ピクセル{pair}:{p1} vs {p2} → 類似度 = {sim:.3f}")
# 平均類似度の計算
mean_sim = np.mean(similarity)
print(f"\nすべての隣接ピクセルの平均類似度:{mean_sim:.3f}")
# ===================== 3. 可視化:画像+類似度分布 =====================
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))
# 3x3シミュレーション画像の表示(ピクセル値をラベル)
ax1.imshow(img, cmap='gray', vmin=0, vmax=255)
ax1.set_title('3x3ランダムグレースケール画像', fontsize=12)
for i in range(3):
for j in range(3):
ax1.text(j, i, str(img[i,j]), ha='center', va='center', color='red', fontsize=12)
ax1.set_xticks([])
ax1.set_yticks([])
# 隣接ピクセル類似度分布の表示
ax2.bar(range(len(similarity)), similarity, color='#1f77b4')
ax2.axhline(y=mean_sim, color='red', linestyle='--', label=f'平均類似度:{mean_sim:.3f}')
ax2.set_xlabel('隣接ピクセル対番号')
ax2.set_ylabel('類似度(0~1)')
ax2.set_title('隣接ピクセル対の類似度分布', fontsize=12)
ax2.legend()
ax2.grid(axis='y', alpha=0.3)
ax2.set_ylim(0, 1) # y軸範囲を限定して視覚的に明確に
plt.tight_layout()
plt.show()
10.3.2 例2:コンピュータビジョン応用 —— 画像ノイズ除去(無向グラフ)
import numpy as np
import matplotlib.pyplot as plt
from skimage import data, util, img_as_float
from skimage.restoration import denoise_bilateral
from skimage.metrics import peak_signal_noise_ratio
# Macフォント設定(日本語とマイナス記号の正常表示を保証)
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 1. 画像の読み込みと前処理
img = img_as_float(data.camera()) # 0〜1の浮動小数点に標準化
# キーとなる修正:RandomStateオブジェクトを整数シードに置き換え(新しいバージョンのrngパラメータ要件に準拠)
seed = 42
noisy_img = util.random_noise(img, mode='gaussian', var=0.01, rng=seed)
# 2. MRF(無向グラフ)に基づくバイラテラルフィルタによるノイズ除去
# sigma_color=0.1:色の類似度重み(小さいほどノイズ除去が激しい)
# sigma_spatial=15:空間近傍範囲(大きいほど考慮する隣接ピクセルが多い)
denoised_img = denoise_bilateral(noisy_img, sigma_color=0.1, sigma_spatial=15)
# 3. 可視化比較
fig, axes = plt.subplots(1, 3, figsize=(18, 6))
axes[0].imshow(img, cmap='gray')
axes[0].set_title('元画像', fontsize=14)
axes[0].axis('off')
axes[1].imshow(noisy_img, cmap='gray')
axes[1].set_title('ノイズ付加後の画像', fontsize=14)
axes[1].axis('off')
axes[2].imshow(denoised_img, cmap='gray')
axes[2].set_title('MRF(バイラテラル)によるノイズ除去後の画像', fontsize=14)
axes[2].axis('off')
plt.suptitle('無向グラフモデル(マルコフ随機場)による画像ノイズ除去比較', fontsize=16)
plt.tight_layout()
plt.show()
# 4. PSNRの計算と出力
psnr_noisy = peak_signal_noise_ratio(img, noisy_img, data_range=1.0)
psnr_denoised = peak_signal_noise_ratio(img, denoised_img, data_range=1.0)
print(f"ノイズ画像PSNR:{psnr_noisy:.2f} dB")
print(f"ノイズ除去画像PSNR:{psnr_denoised:.2f} dB")
print(f"ノイズ除去によるPSNR向上:{psnr_denoised - psnr_noisy:.2f} dB")
10.4 有向グラフモデルと無向グラフモデルの比較
主な比較
視覚的比較コード
import numpy as np
import matplotlib.pyplot as plt
import networkx as nx
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 比較図の作成
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(12, 5))
# 有向グラフの描画
G_dir = nx.DiGraph()
G_dir.add_edges_from([('原因', '結果'), ('原因', '中間変数'), ('中間変数', '結果')])
pos_dir = nx.spring_layout(G_dir, seed=42)
nx.draw(G_dir, pos_dir, with_labels=True, ax=ax1,
node_color='#1f77b4', node_size=2000,
font_size=12, font_family='Arial Unicode MS',
arrowstyle='->', arrowsize=20)
ax1.set_title('有向グラフモデル(因果関係)', fontsize=12)
# 無向グラフの描画
G_undir = nx.Graph()
G_undir.add_edges_from([('ピクセルA', 'ピクセルB'), ('ピクセルB', 'ピクセルC'), ('ピクセルA', 'ピクセルC')])
pos_undir = nx.spring_layout(G_undir, seed=42)
nx.draw(G_undir, pos_undir, with_labels=True, ax=ax2,
node_color='#ff7f0e', node_size=2000,
font_size=12, font_family='Arial Unicode MS',
edge_color='gray')
ax2.set_title('無向グラフモデル(相互依存)', fontsize=12)
plt.suptitle('有向グラフ vs 無向グラフ 視覚的比較', fontsize=14)
plt.tight_layout()
plt.show()
10.5 コンピュータビジョンにおけるグラフモデル
主な応用分野
- 画像分類:有向グラフ(ナイーブベイズ、ベイズネットワーク);
- 画像ノイズ除去/修復:無向グラフ(マルコフ随機場);
- 物体検出:混合グラフモデル(因果と依存の組み合わせ);
- セマンティックセグメンテーション:条件付き随機場(CRF、無向グラフの拡張)。
実装例:CRFによる画像セマンティックセグメンテーション(無向グラフの拡張)
import numpy as np
import matplotlib.pyplot as plt
from skimage import data, segmentation, color
from skimage.filters import gaussian
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# サンプル画像の読み込み
img = data.astronaut()
# Lab色空間への変換(CRFに適した形式)
img_lab = color.rgb2lab(img)
# 初期セグメンテーション(スーパーピクセル)
segments = segmentation.slic(img, n_segments=100, compactness=10, sigma=1)
# CRF(無向グラフモデル):スーパーピクセルの隣接依存関係に基づく
# 簡易CRF実装(核心は隣接スーパーピクセルの色類似度)
def crf_simple(img_lab, segments, iterations=5):
n_segments = np.max(segments) + 1
# 各スーパーピクセルの平均色を計算
mean_colors = np.zeros((n_segments, 3))
for i in range(n_segments):
mean_colors[i] = np.mean(img_lab[segments == i], axis=0)
# CRF反復最適化
for _ in range(iterations):
new_segments = segments.copy()
for i in range(img_lab.shape[0]):
for j in range(img_lab.shape[1]):
# 隣接ピクセルのスーパーピクセルを確認
neighbors = []
if i > 0:
neighbors.append(segments[i-1, j])
if i < img_lab.shape[0]-1:
neighbors.append(segments[i+1, j])
if j > 0:
neighbors.append(segments[i, j-1])
if j < img_lab.shape[1]-1:
neighbors.append(segments[i, j+1])
# 色が最も近いスーパーピクセルを選択
if neighbors:
current_color = img_lab[i, j]
distances = [np.linalg.norm(current_color - mean_colors[n]) for n in neighbors]
best_neighbor = neighbors[np.argmin(distances)]
new_segments[i, j] = best_neighbor
segments = new_segments
# 平均色を更新
for i in range(n_segments):
mean_colors[i] = np.mean(img_lab[segments == i], axis=0)
return segments
# CRFセグメンテーションの実行
crf_segments = crf_simple(img_lab, segments)
# 可視化比較
fig, (ax1, ax2, ax3) = plt.subplots(1, 3, figsize=(18, 6))
ax1.imshow(img)
ax1.set_title('元画像', fontsize=12)
ax1.axis('off')
ax2.imshow(segmentation.mark_boundaries(img, segments))
ax2.set_title('初期スーパーピクセルセグメンテーション', fontsize=12)
ax2.axis('off')
ax3.imshow(segmentation.mark_boundaries(img, crf_segments))
ax3.set_title('CRF(無向グラフ)最適化後セグメンテーション', fontsize=12)
ax3.axis('off')
plt.suptitle('CRF(無向グラフモデル)による画像セマンティックセグメンテーション', fontsize=14)
plt.tight_layout()
plt.show()
10.6 複数未知変数を持つモデルの推論
基本概念
グラフモデルに複数の未知変数が含まれる場合、推論の目的は主に以下の4つです:
- 最大事後確率(MAP)の求解:最も可能性が高い変数組み合わせを特定(例:最も可能性のある病気の組み合わせ);
- 周辺分布の求解:単一変数の確率を求める(他の変数は無視);
- 周辺分布の最大化:単一変数の最も可能性が高い値を求める;
- 事後分布からのサンプリング:事後分布に従ったサンプルを生成。
10.6.1 最大事後確率の解の求解
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import minimize
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 事後確率関数の定義:P(x1,x2) = -((x1-2)^2 + (x2-3)^2 + 0.5*x1*x2)
# (負号はminimizeが最小値を求めるため、確率最大値に対応)
def posterior_neg(params):
x1, x2 = params
return (x1-2)**2 + (x2-3)**2 + 0.5*x1*x2
# MAP(最大事後確率)の求解
initial_guess = [0, 0]
result = minimize(posterior_neg, initial_guess, method='L-BFGS-B')
map_x1, map_x2 = result.x
# 事後分布とMAP点の可視化
x1 = np.linspace(-2, 6, 100)
x2 = np.linspace(-1, 7, 100)
X1, X2 = np.meshgrid(x1, x2)
Z = -((X1-2)**2 + (X2-3)**2 + 0.5*X1*X2) # 事後確率
fig, ax = plt.subplots(figsize=(8, 6))
contour = ax.contourf(X1, X2, Z, levels=20, cmap='viridis')
plt.colorbar(contour, ax=ax, label='事後確率')
ax.scatter(map_x1, map_x2, color='red', s=100, label=f'MAP点: ({map_x1:.2f}, {map_x2:.2f})')
ax.set_xlabel('x1')
ax.set_ylabel('x2')
ax.set_title('最大事後確率(MAP)求解例', fontsize=12)
ax.legend()
plt.show()
print(f"MAP解:x1={map_x1:.2f}, x2={map_x2:.2f}")
10.6.2 事後確率分布の周辺分布の求解
import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import quad
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 結合事後確率 P(x1,x2) の定義
def joint_posterior(x1, x2):
return np.exp(-((x1-2)**2 + (x2-3)**2 + 0.5*x1*x2))
# x1の周辺分布 P(x1) = ∫P(x1,x2)dx2
def marginal_x1(x1):
# 積分範囲:x2は-∞から+∞(近似で-10から10)
return quad(lambda x2: joint_posterior(x1, x2), -10, 10)[0]
# x1の周辺分布を計算
x1_vals = np.linspace(-2, 6, 100)
marginal_vals = [marginal_x1(x1) for x1 in x1_vals]
# 周辺分布の可視化
fig, ax = plt.subplots(figsize=(8, 5))
ax.plot(x1_vals, marginal_vals, color='#2ca02c', linewidth=2)
ax.fill_between(x1_vals, marginal_vals, alpha=0.3, color='#2ca02c')
ax.set_xlabel('x1')
ax.set_ylabel('周辺確率 P(x1)')
ax.set_title('事後確率の周辺分布(x1)', fontsize=12)
plt.show()
# 周辺分布のピークを求める
max_idx = np.argmax(marginal_vals)
max_x1 = x1_vals[max_idx]
print(f"x1周辺分布のピーク:{max_x1:.2f}")
10.6.3 周辺分布の最大化
import numpy as np
import matplotlib.pyplot as plt
from scipy.optimize import minimize
from scipy.integrate import quad # キー:quad積分関数のインポート補完
# Macフォント設定(日本語とマイナス記号の正常表示を保証)
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 結合事後確率関数 P(x1,x2) の定義
def joint_posterior(x1, x2):
"""2次元結合事後確率分布の定義"""
return np.exp(-((x1-2)**2 + (x2-3)**2 + 0.5*x1*x2))
# 周辺分布の負関数(最小化することで周辺分布の最大化を実現)
def marginal_neg(x1):
"""
x1周辺分布の負値を計算
引数x1:scalar(最適化器は配列を受け取るため、最初の要素を取り出す)
戻り値:周辺分布の負値(minimizeで最小値を求めることで最大値を求める)
"""
# 修正:最適化器が渡すx1は配列(例:[0])なのでスカラーに変換
x1_scalar = x1[0] if isinstance(x1, (np.ndarray, list)) else x1
# x2で積分してx1の周辺分布を計算
marginal = quad(lambda x2: joint_posterior(x1_scalar, x2), -10, 10)[0]
return -marginal # 負号:minimizeで最小値を求める=周辺分布の最大値を求める
# 周辺分布の最大化(即ち負周辺分布の最小化)
initial_guess = [0] # 初期推定値
result = minimize(marginal_neg, initial_guess, method='L-BFGS-B')
max_marginal_x1 = result.x[0]
# x1の範囲内の周辺分布値を計算(可視化用)
x1_vals = np.linspace(-2, 6, 100)
marginal_vals = []
for x1 in x1_vals:
# 各x1に対してx2で積分して周辺分布P(x1)を計算
marg = quad(lambda x2: joint_posterior(x1, x2), -10, 10)[0]
marginal_vals.append(marg)
marginal_vals = np.array(marginal_vals)
# 結果の可視化
fig, ax = plt.subplots(figsize=(8, 5))
# 周辺分布曲線を描画
ax.plot(x1_vals, marginal_vals, color='#d62728', linewidth=2, label='周辺分布 P(x1)')
# 周辺分布最大化点をマーク
max_marginal_val = -marginal_neg([max_marginal_x1]) # 該当点の周辺分布値を計算
ax.scatter(max_marginal_x1, max_marginal_val, color='red', s=100,
label=f'周辺分布最大化点: x1={max_marginal_x1:.2f}')
# チャートの装飾
ax.set_xlabel('x1', fontsize=11)
ax.set_ylabel('周辺確率 P(x1)', fontsize=11)
ax.set_title('周辺分布最大化例', fontsize=12)
ax.legend(fontsize=10)
ax.grid(axis='y', alpha=0.3)
plt.show()
# 結果出力
print(f"周辺分布最大化のx1値:{max_marginal_x1:.2f}")
print(f"該当点の周辺確率値:{max_marginal_val:.4f}")
10.6.4 事後分布からのサンプリング
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import multivariate_normal
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 事後分布(多次元ガウス分布)の定義
mean = [2, 3]
cov = [[1, -0.25], [-0.25, 1]] # 共分散行列
rv = multivariate_normal(mean, cov)
# サンプリング
np.random.seed(42)
samples = rv.rvs(size=1000)
# サンプリング結果の可視化
fig, ax = plt.subplots(figsize=(8, 6))
# 分布の輪郭を描画
x1 = np.linspace(-1, 5, 100)
x2 = np.linspace(0, 6, 100)
X1, X2 = np.meshgrid(x1, x2)
Z = rv.pdf(np.dstack((X1, X2)))
ax.contour(X1, X2, Z, levels=10, colors='gray', alpha=0.5)
# サンプル点を描画
ax.scatter(samples[:,0], samples[:,1], alpha=0.6, s=20, color='#1f77b4', label='サンプル点')
ax.set_xlabel('x1')
ax.set_ylabel('x2')
ax.set_title('事後分布からのサンプリング例(多次元ガウス)', fontsize=12)
ax.legend()
plt.show()
# サンプルの平均と共分散を検証
print(f"サンプル平均:{np.mean(samples, axis=0)}")
print(f"サンプル共分散:\n{np.cov(samples.T)}")
10.7 サンプル生成
10.7.1 有向グラフモデルからのサンプリング
import numpy as np
import matplotlib.pyplot as plt
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 有向グラフ:天気→活動→気分
# 条件確率表の定義
# P(天気):晴=0.6、曇=0.3、雨=0.1
weather_probs = [0.6, 0.3, 0.1]
weather_labels = ['晴', '曇', '雨']
# P(活動|天気):晴→散歩(0.7)/読書(0.3);曇→散歩(0.2)/読書(0.8);雨→散歩(0.0)/読書(1.0)
activity_probs = [[0.7, 0.3], [0.2, 0.8], [0.0, 1.0]]
activity_labels = ['散歩', '読書']
# P(気分|活動):散歩→嬉しさ(0.9)/悲しみ(0.1);読書→嬉しさ(0.6)/悲しみ(0.4)
mood_probs = [[0.9, 0.1], [0.6, 0.4]]
mood_labels = ['嬉しさ', '悲しみ']
# サンプル生成関数
def sample_directed_graph(n_samples=1000):
samples = []
for _ in range(n_samples):
# 1. 天気をサンプル
weather = np.random.choice(3, p=weather_probs)
# 2. 活動をサンプル(天気に依存)
activity = np.random.choice(2, p=activity_probs[weather])
# 3. 気分をサンプル(活動に依存)
mood = np.random.choice(2, p=mood_probs[activity])
samples.append([weather_labels[weather], activity_labels[activity], mood_labels[mood]])
return np.array(samples)
# サンプル実行
samples = sample_directed_graph(1000)
# サンプル結果の可視化
fig, axes = plt.subplots(1, 3, figsize=(15, 5))
# 天気分布
weather_counts = np.unique(samples[:,0], return_counts=True)
axes[0].bar(weather_counts[0], weather_counts[1], color='#1f77b4')
axes[0].set_title('天気サンプル分布', fontsize=12)
axes[0].set_ylabel('サンプル数')
# 活動分布
activity_counts = np.unique(samples[:,1], return_counts=True)
axes[1].bar(activity_counts[0], activity_counts[1], color='#ff7f0e')
axes[1].set_title('活動サンプル分布', fontsize=12)
# 気分分布
mood_counts = np.unique(samples[:,2], return_counts=True)
axes[2].bar(mood_counts[0], mood_counts[1], color='#2ca02c')
axes[2].set_title('気分サンプル分布', fontsize=12)
plt.suptitle('有向グラフモデルのサンプル結果', fontsize=14)
plt.tight_layout()
plt.show()
# 条件確率の検証出力
happy_walk = np.sum((samples[:,1]=='散歩') & (samples[:,2]=='嬉しさ')) / np.sum(samples[:,1]=='散歩')
print(f"サンプルから得られるP(嬉しさ|散歩) = {happy_walk:.2f}(理論値0.9)")
10.7.2 無向グラフモデルからのサンプリング
import numpy as np
import matplotlib.pyplot as plt
from scipy.stats import multivariate_normal
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# 無向グラフのサンプリング:MCMC(マルコフ連鎖モンテカルロ)
def mcmc_sample_undirected(n_samples=10000, burn_in=1000):
# 無向グラフのエネルギー関数(多次元ガウス)
mean = [0, 0]
cov = [[1, 0.8], [0.8, 1]]
rv = multivariate_normal(mean, cov)
# サンプル点の初期化
current = np.array([0.0, 0.0])
samples = []
for i in range(n_samples + burn_in):
# 提案分布:ガウスランダムウォーク
proposal = current + np.random.normal(0, 0.5, 2)
# 受理確率の計算
p_current = rv.pdf(current)
p_proposal = rv.pdf(proposal)
accept_prob = min(1, p_proposal / p_current)
# 受理または拒否
if np.random.uniform() < accept_prob:
current = proposal
# 燃焼期間をスキップ
if i >= burn_in:
samples.append(current.copy())
return np.array(samples)
# サンプル実行
samples = mcmc_sample_undirected()
# サンプル結果の可視化
fig, ax = plt.subplots(figsize=(8, 6))
# 真の分布を描画
x = np.linspace(-3, 3, 100)
y = np.linspace(-3, 3, 100)
X, Y = np.meshgrid(x, y)
Z = multivariate_normal.pdf(np.dstack((X, Y)), mean=[0,0], cov=[[1,0.8],[0.8,1]])
ax.contour(X, Y, Z, levels=10, colors='gray', alpha=0.5)
# サンプル点を描画
ax.scatter(samples[:,0], samples[:,1], alpha=0.1, s=5, color='#d62728')
ax.set_xlabel('x1')
ax.set_ylabel('x2')
ax.set_title('無向グラフモデルMCMCサンプル結果', fontsize=12)
plt.show()
# サンプル共分散の検証
print(f"サンプル共分散:\n{np.cov(samples.T)}")
print(f"理論共分散:\n[[1, 0.8], [0.8, 1]]")
10.8 学習
基本概念
グラフモデルの「学習」とは、データからモデルのパラメータ(例えば条件確率表、エネルギー関数のパラメータ)を推定することです。
10.8.1 有向グラフモデルの学習
import numpy as np
import matplotlib.pyplot as plt
from sklearn.naive_bayes import GaussianNB
from sklearn.datasets import load_iris
from sklearn.model_selection import train_test_split
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# アイリスデータセットの読み込み
iris = load_iris()
X, y = iris.data, iris.target
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=42)
# 有向グラフモデル(ナイーブベイズ)のパラメータ学習
gnb = GaussianNB()
gnb.fit(X_train, y_train)
# 学習されたパラメータの表示(平均と分散)
print("学習されたクラス平均:")
for i, class_name in enumerate(iris.target_names):
print(f"{class_name}: {gnb.theta_[i]}")
print("\n学習されたクラス分散:")
for i, class_name in enumerate(iris.target_names):
print(f"{class_name}: {gnb.var_[i]}")
# 学習効果の検証
y_pred = gnb.predict(X_test)
accuracy = np.sum(y_pred == y_test) / len(y_test)
print(f"\nモデル精度:{accuracy:.4f}")
# 学習結果の可視化(前2つの特徴量)
fig, ax = plt.subplots(figsize=(8, 6))
colors = ['#1f77b4', '#ff7f0e', '#2ca02c']
for i, color in enumerate(colors):
mask = y_test == i
ax.scatter(X_test[mask, 0], X_test[mask, 1], color=color, label=iris.target_names[i], alpha=0.7)
# 間違いをマーク
mask_wrong = mask & (y_pred != i)
ax.scatter(X_test[mask_wrong, 0], X_test[mask_wrong, 1], color='red', marker='x', s=100)
ax.set_xlabel(iris.feature_names[0])
ax.set_ylabel(iris.feature_names[1])
ax.set_title('有向グラフモデル(ナイーブベイズ)学習結果', fontsize=12)
ax.legend()
plt.show()
10.8.2 無向グラフモデルの学習
import numpy as np
import matplotlib.pyplot as plt
from skimage import data, util, img_as_float
from skimage.restoration import denoise_bilateral
from sklearn.metrics import mean_squared_error
# Macフォント設定
plt.rcParams['font.sans-serif'] = ['Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False
# ===================== 1. 最適化:ノイズ画像数の削減+計算量削減 =====================
# 画像の読み込みと標準化
img = img_as_float(data.camera())
# 最適化1:ノイズ画像を5つから2つへ減らして計算量を大幅削減(コア高速化ポイント)
noisy_imgs = [util.random_noise(img, mode='gaussian', var=0.01 + i*0.005, rng=i)
for i in range(2)] # 5つ→2つ、計算量を2.5倍軽減
# ===================== 2. 最適化:パラメータ範囲の絞り込み+事前処理 =====================
# 最適化2:パラメータ範囲を絞って探索回数を削減(主な傾向だけ保持)
sigma_colors = [0.1, 0.15, 0.2] # 5つ→3つ、計算量をさらに1.7倍軽減
mse_scores = []
# パラメータを走査してMSEを計算(ロジックは変更なし、計算量を削減)
for sigma in sigma_colors:
# 最適化3:バイラテラルフィルタの空間近傍計算量を削減(sigma_spatialを15から10へ、計算速度を約2倍高速化)
denoised_imgs = [denoise_bilateral(noisy_img, sigma_color=sigma, sigma_spatial=10)
for noisy_img in noisy_imgs]
mse = np.mean([mean_squared_error(img, denoised) for denoised in denoised_imgs])
mse_scores.append(mse)
# 最適パラメータの選定
best_sigma = sigma_colors[np.argmin(mse_scores)]
print(f"学習された最適sigma_color値(MRFのbetaに相当):{best_sigma}")
print(f"最小平均二乗誤差(MSE):{np.min(mse_scores):.4f}")
# ===================== 3. 可視化(ロジックは変更なし、効率を向上) =====================
fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 6))
# パラメータ-MSE曲線
ax1.plot(sigma_colors, mse_scores, marker='o', color='#d62728', linewidth=2, label='MSE曲線')
ax1.scatter(best_sigma, np.min(mse_scores), color='red', s=150,
label=f'最適パラメータ={best_sigma}\n最小MSE={np.min(mse_scores):.4f}')
ax1.set_xlabel('sigma_color(MRFの等価betaパラメータ)', fontsize=11)
ax1.set_ylabel('平均二乗誤差(MSE)', fontsize=11)
ax1.set_title('無向グラフモデル(MRF)パラメータ学習', fontsize=13)
ax1.legend(fontsize=10)
ax1.grid(axis='y', alpha=0.3)
# ノイズ除去効果比較(1組のノイズ画像のみ使用して計算量削減)
noisy_img = noisy_imgs[0]
denoised_best = denoise_bilateral(noisy_img, sigma_color=best_sigma, sigma_spatial=10)
denoised_worst = denoise_bilateral(noisy_img, sigma_color=sigma_colors[0], sigma_spatial=10)
compare_img = np.hstack([noisy_img, denoised_worst, denoised_best])
ax2.imshow(compare_img, cmap='gray')
ax2.set_title(f'ノイズ画像 | 最悪パラメータ({sigma_colors[0]}) | 最適パラメータ({best_sigma})', fontsize=12)
ax2.axis('off')
plt.suptitle('無向グラフモデル(MRF)パラメータ学習結果', fontsize=15)
plt.tight_layout()
plt.show()
議論
グラフモデルは確率理論とコンピュータビジョンを結ぶ橋渡しであり、その主な価値は以下の通りです:
- 複雑な確率関係を視覚化し、理解とモデリングの難易度を下げます;
- 条件独立性に基づいて計算を簡略化し、大規模問題を解けるようにします;
- 推論と学習を兼ね備え、データから自動的にモデルを最適化できます。
実際の応用では、有向グラフと無向グラフの選択は変数間の因果関係によります:
- 因果関係がある(診断、分類)→ 有向グラフ;
- 相互依存がある(画像ピクセル、時系列データ)→ 無向グラフ。
備考
- 本記事のすべてのコードはMacシステム上で直接実行可能です。日本語表示が設定されています;
- 複雑な数式は避けており、直感的理解と実践的な効果に焦点を当てています;
- グラフモデルの高度な内容(正確な推論、近似推論)については、原著『コンピュータビジョン:モデル、学習と推論』をご参照ください;
- 依存ライブラリのインストール:
pip install numpy matplotlib scipy scikit-learn scikit-image networkx。
演習問題
- 本文の有向グラフコードを修正し、「咳→発熱→肺炎」の推論プロセスをシミュレートしてください;
- 無向グラフノイズ除去コードのbetaパラメータを調整し、異なるパラメータがノイズ除去効果に与える影響を観察してください;
- 本文のサンプルコードを組み合わせて、簡単な「画像修復」モデルを実装してください(無向グラフで欠損領域を埋めます)。
まとめ
1. 核心概念:条件独立性はグラフモデルの基盤であり、有向グラフは因果関係を、無向グラフは相互依存関係を表現します;
2. 実践的ポイント:有向グラフは分類/推論に適し、無向グラフは画像モデリング/ノイズ除去に適し、CRFは無向グラフにおける視覚的応用の中心です;
3. 重要な能力:グラフモデルの学習はパラメータ推定であり、推論は未知変数の解(MAP、周辺分布、サンプリングなど)を求めます。
この記事がグラフモデルを完全に理解するお手伝いになれば幸いです。質問があればコメント欄でお気軽にどうぞ。