import numpy as np
import matplotlib.pyplot as plt
from sklearn.preprocessing import StandardScaler
from sklearn.metrics import roc_auc_score, roc_curve

# ===== 1. MVTecAD風の合成データ =====
np.random.seed(42)

def generate_mvtec_like_data(n_normal=800, n_anomaly=80, n_features=20):
    """
    製造検査データを模した合成データ
    正常データ: 低次元多様体上のデータ（実効的な自由度が低い）
    異常データ: 多様体から逸脱したデータ
    """
    # 正常データ: 5次元の潜在因子から20次元に写像
    latent_dim = 5
    W = np.random.randn(latent_dim, n_features) * 0.5
    z_normal = np.random.randn(n_normal, latent_dim)
    X_normal = z_normal @ W + np.random.randn(n_normal, n_features) * 0.1

    # 非線形性を追加
    X_normal += 0.3 * np.sin(X_normal)

    # 異常データの生成（3種類）
    n_each = n_anomaly // 3

    # タイプ1: 特徴量のシフト
    X_anom1 = np.random.randn(n_each, latent_dim) @ W + np.random.randn(n_each, n_features) * 0.1
    X_anom1 += 0.3 * np.sin(X_anom1)
    for i in range(n_each):
        shift_idx = np.random.choice(n_features, 3, replace=False)
        X_anom1[i, shift_idx] += np.random.choice([-1, 1], 3) * np.random.uniform(2, 4, 3)

    # タイプ2: 分布外のデータ
    X_anom2 = np.random.randn(n_each, n_features) * 1.5

    # タイプ3: 相関構造の破壊
    n_rest = n_anomaly - 2 * n_each
    X_anom3 = np.random.randn(n_rest, latent_dim) @ W
    X_anom3 += np.random.randn(n_rest, n_features) * 0.8  # ノイズを増大

    X_anomaly = np.vstack([X_anom1, X_anom2, X_anom3])

    X = np.vstack([X_normal, X_anomaly])
    y = np.array([0] * n_normal + [1] * n_anomaly)

    return X, y

X, y_true = generate_mvtec_like_data()
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)

# 学習/テストの分割
n_train_normal = 600
X_train = X_scaled[y_true == 0][:n_train_normal]
X_test = X_scaled
y_test = y_true

input_dim = X_scaled.shape[1]
print(f"データ形状: {X.shape}, 正常: {np.sum(y_true==0)}, 異常: {np.sum(y_true==1)}")

# ===== 2. 手法1: オートエンコーダ =====
class AnomalyAE:
    """異常検知用オートエンコーダ"""
    def __init__(self, input_dim, hidden_dim=64, latent_dim=8, lr=0.005):
        self.lr = lr
        dims = [input_dim, hidden_dim, hidden_dim // 2, latent_dim,
                hidden_dim // 2, hidden_dim, input_dim]
        self.weights = []
        self.biases = []
        for i in range(len(dims) - 1):
            scale = np.sqrt(2.0 / (dims[i] + dims[i+1]))
            self.weights.append(np.random.randn(dims[i], dims[i+1]) * scale)
            self.biases.append(np.zeros(dims[i+1]))
        self.n_layers = len(self.weights)

    def forward(self, X):
        self.acts = [X]
        self.pres = []
        h = X
        for i in range(self.n_layers):
            z = h @ self.weights[i] + self.biases[i]
            self.pres.append(z)
            if i < self.n_layers - 1:
                h = np.maximum(0, z)
            else:
                h = z
            self.acts.append(h)
        return h

    def backward_and_update(self, X, X_hat, batch_size):
        delta = 2 * (X_hat - X) / batch_size
        for i in range(self.n_layers - 1, -1, -1):
            gW = self.acts[i].T @ delta
            gb = np.sum(delta, axis=0)
            self.weights[i] -= self.lr * gW
            self.biases[i] -= self.lr * gb
            if i > 0:
                delta = (delta @ self.weights[i].T) * (self.pres[i-1] > 0).astype(float)

    def fit(self, X_train, epochs=300, batch_size=64):
        n = X_train.shape[0]
        losses = []
        for epoch in range(epochs):
            idx = np.random.permutation(n)
            eloss = 0
            for s in range(0, n, batch_size):
                e = min(s + batch_size, n)
                Xb = X_train[idx[s:e]]
                Xh = self.forward(Xb)
                loss = np.mean((Xb - Xh) ** 2)
                eloss += loss * (e - s)
                self.backward_and_update(Xb, Xh, e - s)
            losses.append(eloss / n)
        return losses

    def anomaly_score(self, X):
        X_hat = self.forward(X)
        return np.mean((X - X_hat) ** 2, axis=1)

# AEの学習
print("AEの学習...")
ae = AnomalyAE(input_dim, hidden_dim=64, latent_dim=8, lr=0.005)
ae_losses = ae.fit(X_train, epochs=300)
ae_scores = ae.anomaly_score(X_test)
print(f"AE AUC: {roc_auc_score(y_test, ae_scores):.4f}")

# ===== 3. 手法2: Deep SVDD =====
class DeepSVDD:
    """Deep Support Vector Data Description"""
    def __init__(self, input_dim, hidden_dim=64, feature_dim=16, lr=0.005, weight_decay=1e-4):
        self.lr = lr
        self.weight_decay = weight_decay
        self.feature_dim = feature_dim

        # ネットワーク（バイアスなし! — Deep SVDDの重要な制約）
        dims = [input_dim, hidden_dim, hidden_dim // 2, feature_dim]
        self.weights = []
        for i in range(len(dims) - 1):
            scale = np.sqrt(2.0 / (dims[i] + dims[i+1]))
            self.weights.append(np.random.randn(dims[i], dims[i+1]) * scale)
        self.n_layers = len(self.weights)
        self.center = None

    def forward(self, X):
        """順伝播（バイアスなし）"""
        self.acts = [X]
        self.pres = []
        h = X
        for i in range(self.n_layers):
            z = h @ self.weights[i]  # バイアスなし
            self.pres.append(z)
            if i < self.n_layers - 1:
                h = np.maximum(0, z)  # ReLU
            else:
                h = z  # 出力層は線形
            self.acts.append(h)
        return h

    def init_center(self, X_train, batch_size=256):
        """中心 c の初期化: 全データの特徴表現の平均"""
        features = []
        for s in range(0, len(X_train), batch_size):
            e = min(s + batch_size, len(X_train))
            feat = self.forward(X_train[s:e])
            features.append(feat)
        features = np.vstack(features)
        self.center = np.mean(features, axis=0)
        # 中心が原点に近すぎる場合のガード
        self.center[(np.abs(self.center) < 0.01)] = 0.01
        print(f"中心ベクトルのノルム: {np.linalg.norm(self.center):.4f}")

    def fit(self, X_train, epochs=300, batch_size=64):
        """ワンクラスDeep SVDDの学習"""
        self.init_center(X_train)

        n = X_train.shape[0]
        losses = []

        for epoch in range(epochs):
            idx = np.random.permutation(n)
            eloss = 0

            for s in range(0, n, batch_size):
                e = min(s + batch_size, n)
                Xb = X_train[idx[s:e]]
                m = e - s

                # 順伝播
                features = self.forward(Xb)

                # 損失: ||φ(x) - c||² + λ||W||²
                diff = features - self.center
                dist_sq = np.sum(diff ** 2, axis=1)
                loss = np.mean(dist_sq)
                eloss += loss * m

                # 逆伝播
                delta = 2 * diff / m  # d(||φ-c||²)/dφ

                for i in range(self.n_layers - 1, -1, -1):
                    gW = self.acts[i].T @ delta
                    # 重み減衰の追加
                    gW += self.weight_decay * self.weights[i]
                    self.weights[i] -= self.lr * gW
                    if i > 0:
                        delta = (delta @ self.weights[i].T) * (self.pres[i-1] > 0).astype(float)

            losses.append(eloss / n)

        return losses

    def anomaly_score(self, X):
        """異常スコア: ||φ(x) - c||²"""
        features = self.forward(X)
        return np.sum((features - self.center) ** 2, axis=1)

# Deep SVDDの学習
print("\nDeep SVDDの学習...")
svdd = DeepSVDD(input_dim, hidden_dim=64, feature_dim=16, lr=0.005, weight_decay=1e-5)
svdd_losses = svdd.fit(X_train, epochs=300)
svdd_scores = svdd.anomaly_score(X_test)
print(f"Deep SVDD AUC: {roc_auc_score(y_test, svdd_scores):.4f}")

# ===== 4. 手法3: AnoGAN的手法（簡易版） =====
class SimpleGenerator:
    """簡易生成器 G: z → x"""
    def __init__(self, latent_dim, hidden_dim, output_dim, lr=0.001):
        self.lr = lr
        dims = [latent_dim, hidden_dim, hidden_dim, output_dim]
        self.weights = []
        self.biases = []
        for i in range(len(dims) - 1):
            scale = np.sqrt(2.0 / (dims[i] + dims[i+1]))
            self.weights.append(np.random.randn(dims[i], dims[i+1]) * scale)
            self.biases.append(np.zeros(dims[i+1]))
        self.n_layers = len(self.weights)

    def forward(self, z):
        self.acts = [z]
        self.pres = []
        h = z
        for i in range(self.n_layers):
            pre = h @ self.weights[i] + self.biases[i]
            self.pres.append(pre)
            if i < self.n_layers - 1:
                h = np.maximum(0.01 * pre, pre)  # LeakyReLU
            else:
                h = pre  # 線形出力
            self.acts.append(h)
        return h

    def backward(self, grad_output):
        delta = grad_output
        for i in range(self.n_layers - 1, -1, -1):
            gW = self.acts[i].T @ delta / delta.shape[0]
            gb = np.mean(delta, axis=0)
            self.weights[i] -= self.lr * gW
            self.biases[i] -= self.lr * gb
            if i > 0:
                leaky_mask = np.where(self.pres[i-1] > 0, 1.0, 0.01)
                delta = (delta @ self.weights[i].T) * leaky_mask

class SimpleDiscriminator:
    """簡易識別器 D: x → [0, 1]"""
    def __init__(self, input_dim, hidden_dim, lr=0.001):
        self.lr = lr
        dims = [input_dim, hidden_dim, hidden_dim // 2, 1]
        self.weights = []
        self.biases = []
        for i in range(len(dims) - 1):
            scale = np.sqrt(2.0 / (dims[i] + dims[i+1]))
            self.weights.append(np.random.randn(dims[i], dims[i+1]) * scale)
            self.biases.append(np.zeros(dims[i+1]))
        self.n_layers = len(self.weights)

    def forward(self, x):
        self.acts = [x]
        self.pres = []
        h = x
        for i in range(self.n_layers):
            pre = h @ self.weights[i] + self.biases[i]
            self.pres.append(pre)
            if i < self.n_layers - 1:
                h = np.maximum(0.01 * pre, pre)  # LeakyReLU
                self.feature = h  # 中間特徴量を保存（識別損失用）
            else:
                h = 1.0 / (1.0 + np.exp(-np.clip(pre, -500, 500)))  # Sigmoid
            self.acts.append(h)
        return h

    def get_features(self, x):
        """中間層の特徴量を返す"""
        h = x
        for i in range(self.n_layers - 1):
            h = h @ self.weights[i] + self.biases[i]
            h = np.maximum(0.01 * h, h)
        return h

class AnoGANLike:
    """AnoGAN的な異常検知（f-AnoGAN風のエンコーダ付き）"""
    def __init__(self, input_dim, latent_dim=8, hidden_dim=64, lr_g=0.0005,
                 lr_d=0.0005):
        self.latent_dim = latent_dim
        self.G = SimpleGenerator(latent_dim, hidden_dim, input_dim, lr=lr_g)
        self.D = SimpleDiscriminator(input_dim, hidden_dim, lr=lr_d)

        # f-AnoGAN風のエンコーダ
        self.enc_weights = []
        self.enc_biases = []
        enc_dims = [input_dim, hidden_dim, hidden_dim // 2, latent_dim]
        for i in range(len(enc_dims) - 1):
            scale = np.sqrt(2.0 / (enc_dims[i] + enc_dims[i+1]))
            self.enc_weights.append(np.random.randn(enc_dims[i], enc_dims[i+1]) * scale)
            self.enc_biases.append(np.zeros(enc_dims[i+1]))
        self.lr_enc = 0.001

    def encode(self, x):
        """エンコーダ: x → z"""
        h = x
        for i in range(len(self.enc_weights) - 1):
            h = h @ self.enc_weights[i] + self.enc_biases[i]
            h = np.maximum(0, h)
        h = h @ self.enc_weights[-1] + self.enc_biases[-1]
        return h

    def fit_gan(self, X_train, epochs=200, batch_size=64):
        """GANの学習"""
        n = X_train.shape[0]

        for epoch in range(epochs):
            idx = np.random.permutation(n)

            for s in range(0, n, batch_size):
                e = min(s + batch_size, n)
                X_real = X_train[idx[s:e]]
                m = e - s

                # ---- 識別器の更新 ----
                z = np.random.randn(m, self.latent_dim)
                X_fake = self.G.forward(z)

                D_real = self.D.forward(X_real)
                D_fake = self.D.forward(X_fake)

                # 識別器の損失の勾配（Binary Cross Entropy）
                grad_D_real = -(1.0 / (D_real + 1e-8)) / m
                grad_D_fake = (1.0 / (1.0 - D_fake + 1e-8)) / m

                # 簡略化: 識別器の重みを直接更新
                for i in range(self.D.n_layers - 1, -1, -1):
                    if i == self.D.n_layers - 1:
                        # Sigmoid出力層の勾配
                        sig_real = D_real
                        delta_real = grad_D_real * sig_real * (1 - sig_real)
                        sig_fake = D_fake
                        delta_fake = grad_D_fake * sig_fake * (1 - sig_fake)
                    gW = (self.D.acts[i].T @ delta_real + self.D.acts[i].T @ delta_fake) / 2
                    gb = np.mean(delta_real + delta_fake, axis=0) / 2
                    self.D.weights[i] -= self.D.lr * np.clip(gW, -1, 1)
                    self.D.biases[i] -= self.D.lr * np.clip(gb, -1, 1)
                    if i > 0:
                        leaky = np.where(self.D.pres[i-1] > 0, 1.0, 0.01)
                        delta_real = (delta_real @ self.D.weights[i].T) * leaky
                        delta_fake = (delta_fake @ self.D.weights[i].T) * leaky

                # ---- 生成器の更新 ----
                z = np.random.randn(m, self.latent_dim)
                X_fake = self.G.forward(z)
                D_fake = self.D.forward(X_fake)

                # 生成器の勾配: -log(D(G(z)))
                grad_G = -(1.0 / (D_fake + 1e-8)) / m
                sig = D_fake
                delta = grad_G * sig * (1 - sig)

                # D を通じた勾配
                for i in range(self.D.n_layers - 1, -1, -1):
                    if i > 0:
                        leaky = np.where(self.D.pres[i-1] > 0, 1.0, 0.01)
                        delta = (delta @ self.D.weights[i].T) * leaky
                    else:
                        delta = delta @ self.D.weights[i].T

                # Gの逆伝播
                self.G.backward(delta)

    def fit_encoder(self, X_train, epochs=200, batch_size=64):
        """f-AnoGAN風のエンコーダの学習"""
        n = X_train.shape[0]

        for epoch in range(epochs):
            idx = np.random.permutation(n)

            for s in range(0, n, batch_size):
                e = min(s + batch_size, n)
                X_batch = X_train[idx[s:e]]
                m = e - s

                # エンコーダ順伝播
                z_enc = self.encode(X_batch)
                x_recon = self.G.forward(z_enc)

                # 損失: ||x - G(E(x))||²
                recon_loss_grad = 2 * (x_recon - X_batch) / m

                # G を通じた z への勾配（Gの重みは固定）
                delta = recon_loss_grad
                for i in range(self.G.n_layers - 1, -1, -1):
                    if i > 0:
                        leaky = np.where(self.G.pres[i-1] > 0, 1.0, 0.01)
                        delta = (delta @ self.G.weights[i].T) * leaky
                    else:
                        delta = delta @ self.G.weights[i].T

                # エンコーダの逆伝播
                grad_z = delta
                h_list = [X_batch]
                h = X_batch
                for i in range(len(self.enc_weights) - 1):
                    h = h @ self.enc_weights[i] + self.enc_biases[i]
                    h = np.maximum(0, h)
                    h_list.append(h)

                # 最終層
                grad_enc = grad_z
                for i in range(len(self.enc_weights) - 1, -1, -1):
                    gW = h_list[i].T @ grad_enc / m
                    gb = np.mean(grad_enc, axis=0)
                    self.enc_weights[i] -= self.lr_enc * np.clip(gW, -1, 1)
                    self.enc_biases[i] -= self.lr_enc * np.clip(gb, -1, 1)
                    if i > 0:
                        grad_enc = (grad_enc @ self.enc_weights[i].T) * (h_list[i] > 0).astype(float)

    def fit(self, X_train, epochs_gan=200, epochs_enc=200, batch_size=64):
        """GANとエンコーダを順に学習"""
        print("  GAN学習中...")
        self.fit_gan(X_train, epochs=epochs_gan, batch_size=batch_size)
        print("  エンコーダ学習中...")
        self.fit_encoder(X_train, epochs=epochs_enc, batch_size=batch_size)

    def anomaly_score(self, X):
        """異常スコア: ||x - G(E(x))||²"""
        z = self.encode(X)
        x_recon = self.G.forward(z)
        return np.mean((X - x_recon) ** 2, axis=1)

# AnoGAN的手法の学習
print("\nAnoGAN的手法の学習...")
anogan = AnoGANLike(input_dim, latent_dim=8, hidden_dim=64)
anogan.fit(X_train, epochs_gan=200, epochs_enc=200)
anogan_scores = anogan.anomaly_score(X_test)
print(f"AnoGAN-like AUC: {roc_auc_score(y_test, anogan_scores):.4f}")

# ===== 5. 結果の比較と可視化 =====
# AUCスコアの計算
auc_ae = roc_auc_score(y_test, ae_scores)
auc_svdd = roc_auc_score(y_test, svdd_scores)
auc_anogan = roc_auc_score(y_test, anogan_scores)

print("\n" + "=" * 50)
print("AUC-ROC スコアの比較")
print("=" * 50)
print(f"  AE:           {auc_ae:.4f}")
print(f"  Deep SVDD:    {auc_svdd:.4f}")
print(f"  AnoGAN-like:  {auc_anogan:.4f}")

# ROC曲線の比較
fig, axes = plt.subplots(1, 3, figsize=(18, 5))

# (a) ROC曲線
methods = {
    f'AE (AUC={auc_ae:.3f})': ae_scores,
    f'Deep SVDD (AUC={auc_svdd:.3f})': svdd_scores,
    f'AnoGAN-like (AUC={auc_anogan:.3f})': anogan_scores,
}
colors = ['blue', 'red', 'green']
for (name, scores), color in zip(methods.items(), colors):
    fpr, tpr, _ = roc_curve(y_test, scores)
    axes[0].plot(fpr, tpr, color=color, lw=2, label=name)
axes[0].plot([0, 1], [0, 1], 'k--', lw=1, alpha=0.5)
axes[0].set_xlabel('False Positive Rate')
axes[0].set_ylabel('True Positive Rate')
axes[0].set_title('ROC Curve Comparison')
axes[0].legend()
axes[0].grid(True, alpha=0.3)

# (b) スコア分布の比較（箱ひげ図）
all_scores_list = [ae_scores, svdd_scores, anogan_scores]
method_names = ['AE', 'Deep SVDD', 'AnoGAN-like']

normal_scores_list = [s[y_test == 0] for s in all_scores_list]
anomaly_scores_list = [s[y_test == 1] for s in all_scores_list]

positions_normal = [1, 3, 5]
positions_anomaly = [1.6, 3.6, 5.6]

bp1 = axes[1].boxplot([ns / np.percentile(ns, 95) for ns in normal_scores_list],
                       positions=positions_normal, widths=0.4,
                       patch_artist=True, boxprops=dict(facecolor='lightblue'))
bp2 = axes[1].boxplot([ans / np.percentile(ans, 95) for ans in anomaly_scores_list],
                       positions=positions_anomaly, widths=0.4,
                       patch_artist=True, boxprops=dict(facecolor='lightsalmon'))

axes[1].set_xticks([1.3, 3.3, 5.3])
axes[1].set_xticklabels(method_names)
axes[1].set_ylabel('Normalized Anomaly Score')
axes[1].set_title('Score Distribution (Normal vs Anomaly)')
axes[1].legend([bp1['boxes'][0], bp2['boxes'][0]], ['Normal', 'Anomaly'])
axes[1].grid(True, alpha=0.3)

# (c) 学習曲線の比較
axes[2].plot(ae_losses, label='AE', color='blue', alpha=0.7)
axes[2].plot(svdd_losses, label='Deep SVDD', color='red', alpha=0.7)
axes[2].set_xlabel('Epoch')
axes[2].set_ylabel('Training Loss')
axes[2].set_title('Training Loss Curves')
axes[2].legend()
axes[2].grid(True, alpha=0.3)

plt.suptitle('Deep Anomaly Detection: Method Comparison', fontsize=14)
plt.tight_layout()
plt.show()

# ===== 6. Deep SVDDの特徴空間の可視化 =====
from sklearn.decomposition import PCA

# Deep SVDDの特徴表現を取得
features_test = svdd.forward(X_test)

# PCAで2次元に射影
pca = PCA(n_components=2)
features_2d = pca.fit_transform(features_test)
center_2d = pca.transform(svdd.center.reshape(1, -1))[0]

fig, axes = plt.subplots(1, 2, figsize=(14, 5))

# (a) 特徴空間の可視化
scatter_normal = axes[0].scatter(features_2d[y_test == 0, 0], features_2d[y_test == 0, 1],
                                  c='blue', s=10, alpha=0.5, label='Normal')
scatter_anomaly = axes[0].scatter(features_2d[y_test == 1, 0], features_2d[y_test == 1, 1],
                                   c='red', s=30, marker='x', label='Anomaly')
axes[0].scatter(center_2d[0], center_2d[1], c='black', s=200, marker='*',
                label='Center c', zorder=5)

# ハイパースフィアの等高線（2D射影）
normal_dists = np.sum((features_2d[y_test == 0] - center_2d) ** 2, axis=1)
radius = np.percentile(normal_dists, 95)
theta = np.linspace(0, 2 * np.pi, 100)
axes[0].plot(center_2d[0] + np.sqrt(radius) * np.cos(theta),
             center_2d[1] + np.sqrt(radius) * np.sin(theta),
             'k--', lw=2, alpha=0.5, label='95% boundary')

axes[0].set_xlabel('PC1')
axes[0].set_ylabel('PC2')
axes[0].set_title('Deep SVDD Feature Space (PCA projection)')
axes[0].legend()
axes[0].grid(True, alpha=0.3)

# (b) 中心からの距離の分布
dists_normal = np.sqrt(svdd_scores[y_test == 0])
dists_anomaly = np.sqrt(svdd_scores[y_test == 1])

axes[1].hist(dists_normal, bins=40, alpha=0.7, label='Normal', density=True, color='blue')
axes[1].hist(dists_anomaly, bins=20, alpha=0.7, label='Anomaly', density=True, color='red')
axes[1].set_xlabel('Distance from center ||φ(x) - c||')
axes[1].set_ylabel('Density')
axes[1].set_title('Distribution of Distance from Center')
axes[1].legend()
axes[1].grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

# ===== 7. 異常タイプ別の検出性能分析 =====
n_normal = np.sum(y_true == 0)
n_anom_each = (np.sum(y_true == 1)) // 3

# 異常タイプのラベル
anomaly_type = np.zeros(len(y_true), dtype=int)
anomaly_type[n_normal:n_normal + n_anom_each] = 1      # タイプ1: 特徴量シフト
anomaly_type[n_normal + n_anom_each:n_normal + 2*n_anom_each] = 2  # タイプ2: 分布外
anomaly_type[n_normal + 2*n_anom_each:] = 3             # タイプ3: 相関破壊

print("異常タイプ別の平均スコア:")
print(f"{'タイプ':<20} {'AE':>10} {'Deep SVDD':>12} {'AnoGAN-like':>14}")
print("-" * 58)

type_names = ['正常', '特徴量シフト', '分布外', '相関破壊']
for t in range(4):
    mask = anomaly_type == t
    print(f"{type_names[t]:<20} {np.mean(ae_scores[mask]):>10.4f} "
          f"{np.mean(svdd_scores[mask]):>12.4f} {np.mean(anogan_scores[mask]):>14.4f}")

# 各手法の得意・不得意の可視化
fig, axes = plt.subplots(1, 3, figsize=(18, 5))
method_scores = {'AE': ae_scores, 'Deep SVDD': svdd_scores, 'AnoGAN-like': anogan_scores}

for idx, (name, scores) in enumerate(method_scores.items()):
    for t in range(4):
        mask = anomaly_type == t
        axes[idx].hist(scores[mask], bins=30, alpha=0.6,
                       label=type_names[t], density=True)
    axes[idx].set_xlabel('Anomaly Score')
    axes[idx].set_ylabel('Density')
    axes[idx].set_title(f'{name}: Score by Anomaly Type')
    axes[idx].legend(fontsize=8)
    axes[idx].grid(True, alpha=0.3)

plt.tight_layout()
plt.show()

# ===== 8. 最終比較テーブル =====
print("\n" + "=" * 70)
print("深層異常検知手法の比較サマリ")
print("=" * 70)
print(f"{'手法':<16} {'AUC':>8} {'原理':<16} {'推論速度':<10} {'学習容易性':<10}")
print("-" * 70)
results = [
    ('AE', auc_ae, '再構成', '高速', '容易'),
    ('Deep SVDD', auc_svdd, 'ハイパースフィア', '高速', '容易'),
    ('AnoGAN-like', auc_anogan, 'GAN逆写像', '高速(f-)', 'やや困難'),
]
for name, auc, principle, speed, ease in results:
    print(f"{name:<16} {auc:>8.4f} {principle:<16} {speed:<10} {ease:<10}")

print("\n注意: 上記のAUCは合成データに対する結果であり、実データでの性能を")
print("保証するものではありません。データの特性に応じた手法選択が重要です。")