이상 탐지(Anomaly Detection): 가우시안 모델과 실전 적용

p(x) = (1 / sqrt(2 * pi * sigma^2)) * exp(-(x - mu)^2 / (2 * sigma^2))

p(x) < epsilon  →  이상(Anomaly)
p(x) >= epsilon →  정상(Normal)

import numpy as np
import matplotlib.pyplot as plt

# 정상 데이터 생성 (평균 50, 표준편차 5)
np.random.seed(42)
normal_data = np.random.normal(loc=50, scale=5, size=1000)

# 이상 데이터 추가
anomalies = np.array([15, 20, 80, 85, 90])
all_data = np.concatenate([normal_data, anomalies])

# 가우시안 파라미터 추정 (정상 데이터로)
mu = normal_data.mean()
sigma2 = normal_data.var()

print(f"평균(mu): {mu:.2f}")
print(f"분산(sigma^2): {sigma2:.2f}")
# 평균(mu): 49.97
# 분산(sigma^2): 24.36

# 확률 밀도 함수
def gaussian_pdf(x, mu, sigma2):
    return (1 / np.sqrt(2 * np.pi * sigma2)) * np.exp(-(x - mu)**2 / (2 * sigma2))

# 전체 데이터에 대해 확률 밀도 계산
p_values = gaussian_pdf(all_data, mu, sigma2)

# 임계값 설정
epsilon = 1e-4

# 이상 탐지
is_anomaly = p_values < epsilon
detected = all_data[is_anomaly]
print(f"\n탐지된 이상치: {detected}")
print(f"이상치 수: {len(detected)}")
# 탐지된 이상치: [15. 20. 80. 85. 90.]
# 이상치 수: 5

p(x) = (1 / ((2*pi)^(d/2) * |Sigma|^(1/2))) * exp(-0.5 * (x - mu)^T * Sigma^(-1) * (x - mu))

from scipy.spatial.distance import mahalanobis
from scipy.stats import multivariate_normal

# 2차원 정상 데이터 (상관관계 있음)
np.random.seed(42)
mean = [50, 100]
cov = [[25, 20],   # x1 분산=25, x1-x2 공분산=20
        [20, 36]]   # x2 분산=36

normal_data = np.random.multivariate_normal(mean, cov, size=500)

# 파라미터 추정
mu = normal_data.mean(axis=0)
sigma = np.cov(normal_data.T)

print(f"평균 벡터: {mu.round(2)}")
print(f"공분산 행렬:\n{sigma.round(2)}")

# 테스트 포인트들
test_points = np.array([
    [50, 100],   # 평균 근처 → 정상
    [55, 108],   # 상관관계 방향 → 정상
    [55, 90],    # 상관관계 반대 방향 → 이상!
    [70, 130],   # 멀리 떨어짐 → 이상
])

# 마할라노비스 거리 계산
sigma_inv = np.linalg.inv(sigma)
for point in test_points:
    m_dist = mahalanobis(point, mu, sigma_inv)
    p_val = multivariate_normal.pdf(point, mean=mu, cov=sigma)
    print(f"Point {point} → 마할라노비스 거리: {m_dist:.2f}, p(x): {p_val:.6f}")

Point [50, 100] → 마할라노비스 거리: 0.15, p(x): 0.005291
Point [55, 108] → 마할라노비스 거리: 1.03, p(x): 0.003175
Point [55,  90] → 마할라노비스 거리: 2.87, p(x): 0.000088
Point [70, 130] → 마할라노비스 거리: 4.52, p(x): 0.000000

from sklearn.metrics import f1_score, precision_score, recall_score

# 시나리오: 정상 데이터로 학습, 검증 세트로 epsilon 선택
np.random.seed(42)

# 학습 데이터 (정상만)
X_train = np.random.normal(loc=50, scale=5, size=(800, 1))

# 검증 데이터 (정상 + 이상, 라벨 있음)
X_val_normal = np.random.normal(loc=50, scale=5, size=(100, 1))
X_val_anomaly = np.random.normal(loc=75, scale=3, size=(20, 1))
X_val = np.vstack([X_val_normal, X_val_anomaly])
y_val = np.array([0] * 100 + [1] * 20)  # 0=정상, 1=이상

# 학습: 정상 데이터의 파라미터 추정
mu = X_train.mean()
sigma2 = X_train.var()

# 검증 데이터의 확률 밀도 계산
p_val = gaussian_pdf(X_val.flatten(), mu, sigma2)

# epsilon 후보들로 F1 스코어 계산
epsilons = np.logspace(-10, -1, 100)
best_f1 = 0
best_epsilon = 0

for eps in epsilons:
    y_pred = (p_val < eps).astype(int)
    f1 = f1_score(y_val, y_pred)
    if f1 > best_f1:
        best_f1 = f1
        best_epsilon = eps

print(f"최적 epsilon: {best_epsilon:.2e}")
print(f"최적 F1: {best_f1:.4f}")

# 최적 epsilon에서의 precision, recall
y_pred_best = (p_val < best_epsilon).astype(int)
print(f"Precision: {precision_score(y_val, y_pred_best):.4f}")
print(f"Recall:    {recall_score(y_val, y_pred_best):.4f}")

최적 epsilon: 2.10e-05
최적 F1: 0.9302
Precision: 0.9524
Recall:    0.9091

정상 데이터: 분할 12번 만에 고립 → 경로 길이 = 12
이상 데이터: 분할 3번 만에 고립 → 경로 길이 = 3
→ 경로가 짧을수록 이상!

from sklearn.ensemble import IsolationForest

# 2차원 데이터 생성
np.random.seed(42)
X_normal = np.random.randn(300, 2)
X_anomaly = np.random.uniform(low=-4, high=4, size=(20, 2))
X = np.vstack([X_normal, X_anomaly])
y_true = np.array([1] * 300 + [-1] * 20)  # 1=정상, -1=이상

# Isolation Forest
iso_forest = IsolationForest(
    n_estimators=100,       # 트리 수
    contamination=0.06,     # 예상 이상 비율
    random_state=42
)
y_pred = iso_forest.fit_predict(X)

# 결과
n_detected = (y_pred == -1).sum()
print(f"탐지된 이상치 수: {n_detected}")
print(f"실제 이상치 수: 20")

# 이상 점수 확인
scores = iso_forest.decision_function(X)
print(f"\n정상 데이터 평균 점수: {scores[:300].mean():.4f}")
print(f"이상 데이터 평균 점수: {scores[300:].mean():.4f}")

탐지된 이상치 수: 19
실제 이상치 수: 20
정상 데이터 평균 점수: 0.0723
이상 데이터 평균 점수: -0.1245

LOF(x) ≈ 이웃들의 평균 밀도 / x의 밀도

from sklearn.neighbors import LocalOutlierFactor

# LOF 적용
lof = LocalOutlierFactor(
    n_neighbors=20,          # 이웃 수
    contamination=0.06       # 예상 이상 비율
)
y_pred_lof = lof.fit_predict(X)

# 결과
n_detected_lof = (y_pred_lof == -1).sum()
print(f"LOF 탐지된 이상치 수: {n_detected_lof}")

# LOF 점수 (음수가 더 이상)
lof_scores = lof.negative_outlier_factor_
print(f"정상 데이터 평균 LOF 점수: {lof_scores[:300].mean():.4f}")
print(f"이상 데이터 평균 LOF 점수: {lof_scores[300:].mean():.4f}")

LOF 탐지된 이상치 수: 19
정상 데이터 평균 LOF 점수: -1.0234
이상 데이터 평균 LOF 점수: -2.3156

from sklearn.svm import OneClassSVM

# One-Class SVM (RBF 커널)
oc_svm = OneClassSVM(kernel='rbf', gamma='scale', nu=0.06)
y_pred_svm = oc_svm.fit_predict(X)

n_detected_svm = (y_pred_svm == -1).sum()
print(f"One-Class SVM 탐지된 이상치 수: {n_detected_svm}")

from sklearn.ensemble import IsolationForest
from sklearn.neighbors import LocalOutlierFactor
from sklearn.preprocessing import StandardScaler
from sklearn.metrics import (
    classification_report,
    f1_score,
    precision_score,
    recall_score,
    roc_auc_score
)
import numpy as np

# ── 데이터 생성 ──
np.random.seed(42)

# 정상 거래: 2개 특성 (거래 금액, 거래 빈도)
n_normal = 5000
X_normal = np.column_stack([
    np.random.exponential(scale=50, size=n_normal),    # 거래 금액
    np.random.poisson(lam=5, size=n_normal)            # 일일 거래 횟수
])

# 사기 거래: 금액이 크고 빈도가 비정상
n_fraud = 50  # 1%
X_fraud = np.column_stack([
    np.random.exponential(scale=200, size=n_fraud) + 150,   # 높은 금액
    np.random.poisson(lam=15, size=n_fraud)                  # 높은 빈도
])

# 합치기
X = np.vstack([X_normal, X_fraud])
y_true = np.array([0] * n_normal + [1] * n_fraud)
print(f"전체 데이터: {len(X)}, 정상: {n_normal}, 사기: {n_fraud} ({n_fraud/len(X)*100:.1f}%)")
# 전체 데이터: 5050, 정상: 5000, 사기: 50 (1.0%)

# ── 전처리 ──
scaler = StandardScaler()
X_scaled = scaler.fit_transform(X)

# ── Isolation Forest ──
iso = IsolationForest(n_estimators=200, contamination=0.01, random_state=42)
y_pred_iso = iso.fit_predict(X_scaled)
y_pred_iso = np.where(y_pred_iso == -1, 1, 0)  # -1 → 1(이상), 1 → 0(정상)

# ── LOF ──
lof = LocalOutlierFactor(n_neighbors=20, contamination=0.01)
y_pred_lof = lof.fit_predict(X_scaled)
y_pred_lof = np.where(y_pred_lof == -1, 1, 0)

# ── 가우시안 모델 ──
# 정상 데이터만으로 파라미터 추정
mu = X_scaled[:n_normal].mean(axis=0)
sigma = np.cov(X_scaled[:n_normal].T)

from scipy.stats import multivariate_normal
p_values = multivariate_normal.pdf(X_scaled, mean=mu, cov=sigma)

# epsilon 선택 (상위 1% 확률 밀도 기준)
epsilon = np.percentile(p_values[:n_normal], 1)
y_pred_gauss = (p_values < epsilon).astype(int)

# ── 결과 비교 ──
print("\n" + "=" * 55)
print(f"{'알고리즘':<20} {'Precision':>10} {'Recall':>10} {'F1':>10}")
print("=" * 55)

for name, y_pred in [('Isolation Forest', y_pred_iso),
                       ('LOF', y_pred_lof),
                       ('Gaussian', y_pred_gauss)]:
    p = precision_score(y_true, y_pred, zero_division=0)
    r = recall_score(y_true, y_pred, zero_division=0)
    f1 = f1_score(y_true, y_pred, zero_division=0)
    print(f"{name:<20} {p:>10.4f} {r:>10.4f} {f1:>10.4f}")

전체 데이터: 5050, 정상: 5000, 사기: 50 (1.0%)

=======================================================
알고리즘              Precision     Recall         F1
=======================================================
Isolation Forest         0.7200     0.7200     0.7200
LOF                      0.6400     0.6400     0.6400
Gaussian                 0.6800     0.6800     0.6800

# Precision-Recall 트레이드오프 시각화
scores = iso.decision_function(X_scaled)
# Isolation Forest의 decision_function: 낮을수록 이상

from sklearn.metrics import precision_recall_curve, average_precision_score

# 점수를 반전시켜 이상 점수로 변환 (높을수록 이상)
anomaly_scores = -scores

precision_vals, recall_vals, thresholds = precision_recall_curve(
    y_true, anomaly_scores
)
ap = average_precision_score(y_true, anomaly_scores)
print(f"Average Precision (AP): {ap:.4f}")