PCA(Principal Component Analysis) 이상탐지 - 실습, 코드구현

 

https://soochan8.tistory.com/12

PCA(Principal Component Analysis) 이상탐지 - 개념, 수식

 

 

PCA는 파이썬의 scikit-learn 라이브러리에 잘 구현되어 있다. 보통은 그냥 딸깍 해서 쓰면 되지만,

 

처음 배울 때 확실히 이해했는지 확인해보고 싶다면 numpy정도만 사용해 직접 구현을 해 보시는 것을 추천한다.

 

이해했다고 생각했지만 구현하면서 알게 된 것도 많다.

 

 

 

이론편에서 있는 step과 최대한 비슷하게 구성했다.

 

정말 간단한 3차원의 데이터 10개를 임의로 생성한 후, 표준화하고 공분산행렬을 계산한다.

import numpy as np
np.random.seed(0)

# 1. 샘플 데이터 준비 (3차원 데이터, 10개의 샘플)
X = np.array([[2.5, 2.4, 3.2], 
              [0.5, 0.7, 1.0], 
              [2.2, 2.9, 3.1], 
              [1.9, 2.2, 2.9], 
              [3.1, 3.0, 3.5], 
              [2.3, 2.7, 3.0], 
              [2.0, 1.6, 2.5], 
              [1.0, 1.1, 1.2], 
              [1.5, 1.6, 2.0], 
              [1.1, 0.9, 1.5]])

# 2. 표준화
X_mean = np.mean(X, axis=0)
X_std = np.std(X, axis=0)
X_standardized = (X - X_mean) / X_std

# 3. 공분산 행렬 계산
cov_matrix = np.cov(X_standardized.T)  #변수들을 행으로 배치해야 하기 떄문에 전치시킴

print("표준화된 데이터 행렬\n", X_standardized)
print("\n 공분산 행렬 \n", cov_matrix)

 

<출력>

표준화된 데이터 행렬
 [[ 0.92627881  0.61016865  0.9461558 ]
 [-1.7585873  -1.506743   -1.62365008]
 [ 0.52354889  1.23278973  0.82934645]
 [ 0.12081898  0.36112022  0.59572773]
 [ 1.73173864  1.35731394  1.29658388]
 [ 0.6577922   0.9837413   0.71253709]
 [ 0.25506228 -0.38602507  0.12849029]
 [-1.08737078 -1.00864614 -1.39003137]
 [-0.41615425 -0.38602507 -0.4555565 ]
 [-0.95312747 -1.25769457 -1.03960329]]

 공분산 행렬 
 [[1.11111111 1.0288103  1.07518529]
 [1.0288103  1.11111111 1.06198995]
 [1.07518529 1.06198995 1.11111111]]

 

 

고유값분해를 하면 고유벡터와 고유값이 나온다.

 

이 중 고유값이 큰 순서대로 정렬해 첫번째 주성분 축(고유벡터)을 찾고,

여기에 표준화된 데이터를 projection한다. 데이터를 벡터에 프로젝션하는건 내적연산 (np.dot)으로 가능하다.

 

10x3 데이터와 3x1 벡터 (고유값이 가장 큰 첫번째 고유벡터)  내적해 10x1의 데이터를 만드는 것이다.

pc를 몇 개를 쓸 건지에 따라서 고유벡터 수를 선택할 수 있다 

# 4. 고유값과 고유벡터 계산
eigenvalues, eigenvectors = np.linalg.eig(cov_matrix)

# 5. 고유값 정렬 및 첫 번째 주성분 선택
sorted_indices = np.argsort(eigenvalues)[::-1]  # 고유값 내림차순 정렬
eigenvectors = eigenvectors[:, sorted_indices]  # 정렬된 고유벡터

k = 1
W = eigenvectors[:, :k]  # 첫 번째 주성분 축 (3x1)

# 6. 첫번째 eigenvector에 데이터를 projection
Z = np.dot(X_standardized, W) 

print("eigenvalues: ",eigenvalues, "\n eigenvectors \n" ,eigenvectors)
print("\n PC Score: \n", Z) # 각 샘플이 주성분 축에서 얼마나 떨어져있는가?

 

<출력>

eigenvalues:  [3.22188767 0.08335334 0.02809232] -> 3개 주성분의 고유값

 eigenvectors -> 각 주성분에 대응되는 고유벡터. 이 축에 데이터를 projection한다. (주성분공간으로 투영) 
 [[-0.57616336 -0.64223112 -0.50554423]
 [-0.5737167   0.75832621 -0.30950041]
 [-0.58213824 -0.11171637  0.8053785 ]]

 PC Score: -> 위 코드에서 Z에 해당하는 부분. 10개 데이터를 PC1에 투영시켰을 때 값
 [[-1.43454533e+00]
 [ 2.82286599e+00]
 [-1.49171602e+00]
 [-6.23588059e-01]
 [-2.53126909e+00]
 [-1.35817966e+00]
 [-2.87627888e-04]
 [ 2.01437075e+00]
 [ 7.26438714e-01]
 [ 1.87591033e+00]]

 

PC 1에 투영시키는 것의 의미는

아래 그림처럼 우리 샘플 데이터는 3차원 공간에 있으니, 빨간색 축이 첫번째 고유벡터(PC1)이라고 가정한다면

저 축 위에다가 점들을 projection시키는 것이다. 위 결과의 PC Score가 그 값이고..

https://d-craftshop.tistory.com/5

2차원의 경우에는 이렇다. 3차원도 동일하게 PC를 하나만 선택했을 때는 선에,

두 개를 선택했다면 2차원 평면에 주성분을 나타낼 수 있다.

https://ddongwon.tistory.com/114

 

 

 

여기가 핵심이다.

 

우리는 주성분 공간에 투영된 1x10 데이터 Z를 가지고 있다.

이걸 원래 차원인 3차원으로 복원시키려면? 아까 내적했던 고유벡터를 전치시킨 다음에 또 곱해주면 된다. 나는 이걸 이해하는데 꽤 오래 걸렸다..

 

# 7. 복원 
X_reconstructed = np.dot(Z, W.T)  # Z에 eigenvector의 전치를 내적 -> 다시 3차원데이터로 복원
X_reconstructed = X_reconstructed * X_std + X_mean  # 정규화 과정의 역

print("원본 데이터 \n" , X_standardized)
print(" \n 재구축된 데이터 \n" , X_reconstructed)

 

<출력>

원본 데이터 
 [[ 0.92627881  0.61016865  0.9461558 ]
 [-1.7585873  -1.506743   -1.62365008]
 [ 0.52354889  1.23278973  0.82934645]
 [ 0.12081898  0.36112022  0.59572773]
 [ 1.73173864  1.35731394  1.29658388]
 [ 0.6577922   0.9837413   0.71253709]
 [ 0.25506228 -0.38602507  0.12849029]
 [-1.08737078 -1.00864614 -1.39003137]
 [-0.41615425 -0.38602507 -0.4555565 ]
 [-0.95312747 -1.25769457 -1.03960329]]
 
 재구축된 데이터 
 [[2.42569734 2.57093379 3.10492875]
 [0.59844464 0.60942937 0.9831793 ]
 [2.45023462 2.59727387 3.13342069]
 [2.07763986 2.19730387 2.70077514]
 [2.8964039  3.07622406 3.65149869]
 [2.3929217  2.53575006 3.06687069]
 [1.81012345 1.91013252 2.39014334]
 [0.94544544 0.98192497 1.38610592]
 [1.49821732 1.57531006 2.02796758]
 [1.00487173 1.04571743 1.4551099 ]]

 

이건 임의로 생성한 데이터라 의미는 없을테지만

이런식으로 처음과 같은 3x10형태로 데이터가 복원된 것을 확인할 수 있다.

 

이후 원래 데이터와의 차이(L2 norm)를 이용해 재구축오차를 계산한다.

# 8. 샘플별 재구축 오차 계산
sample_errors = np.linalg.norm(X - X_reconstructed, axis=1)
print("샘플별 재구축 오차:\n", sample_errors)

 

<출력>

샘플별 재구축 오차:
 [0.20923142 0.1348233  0.39417942 0.26693395 0.26497798 0.20021032
 0.37987336 0.22705358 0.03734912 0.17971666]

 

이것을 이상치 스코어로 사용하는 것이다..

 

여하튼 PCA는  이런 식으로 구현할 수 있다. 

 

 

 

 

---

 

 

PCA기반 이상탐지 - 간단한 예제

조금 더 높은 차원의 데이터를 생성해 이상탐지하는 과정을 테스트해보았다.

아래 코드에서는 사이킷런의 pca를 이용했다.

n_samples = 300
n_features = 10 

# 밀도가 다른 두 개의 정상데이터 군집을 생성
X1, _ = make_blobs(n_samples=150, centers=1, cluster_std=0.50, random_state=0, center_box=(-5.0, 5.0), n_features=n_features)  # 밀도가 높은 군집
X2, _ = make_blobs(n_samples=150, centers=1, cluster_std=1.5, random_state=1, center_box=(5.0, 15.0), n_features=n_features)  # 밀도가 낮은 군집
X = np.vstack((X1, X2))

# 고르게 분포된 이상치 추가
outliers = np.random.uniform(low=-10, high=20, size=(20, n_features)) 
X = np.vstack((X, outliers))

df = pd.DataFrame(X, columns=[f'Feature {i+1}' for i in range(n_features)])

 

<출력>

 

feature 1,2,3만 시각화해보면 이런 모습이다. 세 개의 차원만 확인한 것이긴 하지만, 분포로부터 벗어난 이상치들이 생성된 것을 대략적으로 확인할 수 있다.

Scree plot을 보면 너무 간단해서인지 첫 번쨰 주성분이 전체 분산의 70% 이상을 설명하는 모습이다. 

 

 

from pyod.models.pca import PCA  # PCA 기반 이상치 탐지 모델

# PCA 기반 이상치 탐지 모델
pca_model = pyodPCA(contamination=0.06)
pca_model.fit(X)

# 이상치 예측 (1: 이상치, 0: 정상)
y_pred = pca_model.labels_

pca = SKPCA(n_components=2)
X_pca = pca.fit_transform(X)

# 재구축 오차 계산 (10D -> 2D -> 10D 복원)
X_reconstructed = pca.inverse_transform(X_pca)
reconstruction_error = np.mean(np.square(X - X_reconstructed), axis=1)


plt.figure(figsize=(8, 6))

# 정상 데이터와 이상치 다른 색으로 표시
plt.scatter(X_pca[:, 0], X_pca[:, 1], c=y_pred, cmap='coolwarm', s=10, label='Data Points')

# 이상치에만 재구축 오차 텍스트 표시
for i, txt in enumerate(np.round(reconstruction_error, 2)):
    if y_pred[i] == 1:  # 이상치인 경우에만
        plt.annotate(txt, (X_pca[i, 0], X_pca[i, 1]), fontsize=8, color='black')

plt.title('PCA-based Outlier Detection with Reconstruction Error Score')
plt.xlabel('Principal Component 1')
plt.ylabel('Principal Component 2')
plt.show()

 

주성분 1,2를 축으로 2차원그래프를 그리고, 이상치로 판단한 데이터들의 이상치스코어를 표기했다.

이때 재구축오차의 상위 6%값을 임계값으로 지정하여 이를 이상치로 간주하였다.

 

임의로 생성한 데이터기에, contamination 값을 바꿔 가며  PC1, PC2의 그래프를 보고 임계값을 정했지만, 실제 데이터에 적용한다면 불량/정상 라벨링이 되어 있을 것이기에 재현율, 정밀도 등을 비교해가며 임계값을 조정하면 될 것이다.