논문 ‘Kernel Principal Component Analysis and its Applications in Face Recognition and Active Shape Models’ 리뷰
커널PCA 개념과 코드 실습 중심으로
본 페이지는 Kernel PCA를 다룬 논문 ‘Kernel Principal Component Analysis and its Applications in Face Recognition and Active Shape Models’에 대한 리뷰 입니다. 고려대 일반대학원 IME654(2020년 2학기) 과제 일환으로 작성됐습니다. 일부 이미지는 온라인 상 공유된 자료들을 참고했으며 소스는 모두 기재돼 있습니다.
(논문 링크: https://arxiv.org/abs/1207.3538)
이 논문은 뉴욕 Rensselaer Polytechnic Institute 박사과정 중이던 Quan Wang님이 2012년 발표했습니다. 저자는 현재 구글 엔지니어로 활동 중이시네요.
Kernel PCA(Schölkopf et al., 1999)를 처음 제시한 논문은 아니지만, Kernel PCA를 설명할 때 상당히 자주 인용되고 있습니다. 이유는 정확히 모르겠지만 PCA, Kernel PCA가 명쾌하게 설명돼 있어 좋았습니다.
논문은 개략적으로 1) PCA와 2) PCA 개념을 활용한 Active Shape Model, 3) Kernel PCA, 4) Kernel PCA를 기반으로 한 Active Shape Model 순으로 서술합니다. PCA는 차원축소 중 feature extraction 방법론 중 하나이고, Active shape model은 이미지 내 특정 shape을 찾는(object-fitting) 통계적 모델의 일종입니다. 이번 과제의 핵심은 Kernel PCA 부분이기 때문에, 논문의 서술 순서를 따라 차례로 요약한 뒤 커널 PCA의 코드 구현까지 시행해보겠습니다.
1. PCA(주성분분석)와 Active shape model
앞에서 언급했듯이 PCA는 고차원 데이터를 저차원 데이터로 변환하는 차원축소 방법론 중 하나입니다. 변환하는 과정에서 원 데이터의 ‘분산’을 최대한 보존하는 변수를 찾는 게 PCA의 목적입니다.
PCA가 데이터를 변환할 때는 사영(project)하는 방법을 씁니다. 그래픽하게는 (그림 1.1)처럼 나타낼 수 있습니다. 세모는 각각 (x, y)좌표를 갖는 2차원에서의 객체이고, 직선은 1차원의 축(주성분)입니다. 세모를 축에 꽂아 내리면 직선 위에서 각각 실수값을 갖게 되는데 이게 새로운 변수값이 됩니다. 새로운 변수값이 원 데이터의 분포 또는 변동을 잘 나타내야 하기 때문에 PCA가 하고자 하는 사영은 (그림 1.1)의 왼쪽보다 오른쪽에 가깝습니다.
사영하는 과정을 식으로 나타내면 아래와 같습니다. D차원의 데이터를 M개 주성분 축에 사영해 M차원의 데이터로 변환하는 과정입니다.
이제 ‘무슨 축에 사영하느냐’(무슨 A를 쓸 것인가)는 문제만 남습니다. 목적은 원 데이터의 분산을 최대한 보존하는 것이기 때문에, 이는 ‘Y의 분산을 최대화하는 A를 찾자’는 최적화 문제가 됩니다. 유도과정을 생략하고 결론만 적자면 “ 원 데이터 공분산 행렬의 고유벡터들을 A(=u)로 쓰면 Y의 분산이 최대화되더라 ” 가 됩니다. 이를 수식으로 표현하면 아래와 같습니다.
보통 PCA 설명은 여기까지 인데 논문에서는 한가지 포인트가 더 나옵니다. 사영(차원축소)한 결과물인 Y로부터 원 데이터 X를 ‘복원’하는 부분인데요. 정확히는 주성분을 다시 내적해주는 방식으로 원 데이터를 근사(approximate)하게 됩니다. 이는 Active shape model에서 중요하게 다뤄집니다.
Active shape model은 PCA의 논리를 차용해 이미지에서 특정 shape를 찾는 모델입니다. shape는 검출하고자 하는 부분의 경계에 점들을 찍어 나타내며, 이 점들의 좌표를 정확히 예측하는 게 목적입니다.
학습 이미지에서 shape을 이루는 점들은 특정 분포를 따른다고 볼 수 있습니다. 예컨대 얼굴 shape를 이루는 점 50개의 좌표를 한 벡터로 놓으면, 학습 이미지 1000장에 해당하는 벡터 1000개(각 100차원)는 어떤 분포를 이룰 것입니다. 이 분포를 학습시킨 뒤 test 이미지의 shape을 분포로부터 크게 벗어나지 않게 한다면 실제 shape에 가까운 점들을 찍을 수 있겠죠. 이때 1000개 벡터에 PCA를 적용하면 변동이 가장 큰 좌표 구역(ex. 귀의 위치)이 어디인지 나타낸 분포의 축(주성분)을 알 수 있을 것입니다.
테스트 이미지의 점들을 학습시킨 분포에 맞게 조정하는 과정은 PCA에서 원 데이터를 복원하는 과정과 같습니다. 학습한 벡터들의 평균 벡터와 주성분(공분산 행렬의 고유벡터)을 통해 원 벡터를 복원할 수 있습니다. (그림 1.3)은 왼쪽부터 (1) 학습 이미지의 평균 shape벡터 (2) 주성분 (3) 파라미터를 최적화해 test 이미지에 가장 가깝게 원 shape 벡터를 복원한 모습을 나타냅니다. 이 과정을 수식으로 표현하면 아래와 같습니다.
2. Kernel PCA
다시 PCA로 돌아가 커널PCA로 발전시켜보겠습니다. PCA로 차원축소를 해도 linear한 분류가 어려운 데이터들이 있습니다. 이때 커널 PCA를 사용하면 데이터를 고차원으로 옮겨 선형 분류가 가능한 형태로 변환하는 것과 같은 효과를 얻을 수 있기 때문에 용이합니다.
우선 데이터를 기존 D차원에서 더 고차원인 M차원으로 변환한 ɸ(x)를 상정합니다. ɸ(x)에 PCA를 적용하면 아래와 같이 식이 전개됩니다.
마지막 식을 ɸ(x)에 대해 풀면 연산이 굉장히 복잡해지지만, 커널트릭은 이를 ‘커널 함수’라고 불리는 k(x, x) 로 대체해 매우 간단한 연산으로 치환시킵니다. 중간 유도과정은 생략하고 우리가 찾고자 하는 주성분 v를 구하는 식으로 정리하면 아래와 같습니다.
즉, 커널함수를 통해 ɸ(x)를 명시적으로 구하지 않고서도 데이터를 고차원으로 옮긴 후 PCA를 적용한 값을 구할 수 있다는 게 커널PCA의 핵심입니다.
커널함수로는 대표적으로 2가지를 사용합니다. 다항커널함수(Polynomial kernel function)과 가우시안 또는 RBF커널함수입니다. 이는 아래 코드 구현을 통해 소개할 예정입니다.
논문은 결과적으로 커널 PCA를 Active shape model에 적용할 경우 모델 성능이 높아진다는 것으로 마무리됩니다. PCA를 적용할 때 shape 벡터의 분포를 선형적으로만 분석할 수 있다면, 커널 PCA로는 더욱 정교한 분석이 가능하기 때문에 예측 성능이 높아진 것으로 이해됩니다.
그럼 마지막으로 논문 일부를 코드로 구현해보겠습니다. 논문에선 3가지 실험이 나오며, 여기선 이중 첫번째 실험인 가상의 3차원 구(球) 데이터에 linear PCA와 커널 PCA를 적용한 결과를 살펴보겠습니다.
# 라이브러리 불러오기
import numpy as npimport pandas as pdfrom mpl_toolkits.mplot3d import Axes3Dimport matplotlib.pyplot as plt# 구 형태의 데이터셋 생성
def create_sphere(cx,cy,cz, r, n): phi = np.random.uniform(0, 2*np.pi, size=n) theta = np.random.uniform(0, np.pi, size=n) r_xy = r*np.sin(theta) x = cx + np.cos(phi) * r_xy y = cy + np.sin(phi) * r_xy z = cz + r * np.cos(theta) return np.stack([x,y,z])
#(0,0,0)이 중심이고 반지름이 100인 구 위의 1000개의 데이터(Class1)sphere1 = create_sphere(0, 0, 0, 100, 1000) #(0,0,0)이 중심이고 반지름이 30인 구 위의 200개 데이터(Class2)
sphere2 = create_sphere(0, 0, 0, 30, 200)print(np.shape(sphere1))print(np.shape(sphere2))
fig = plt.figure(figsize=(8,8))ax = fig.add_subplot(111, projection='3d')ax.scatter(sphere1[0,:],sphere1[1,:],sphere1[2,:], c='b', s=1, label='Class1')ax.scatter(sphere2[0,:],sphere2[1,:],sphere2[2,:], c='r', s=1, label='Class2')ax.set_xlabel("X")ax.set_ylabel('Y')ax.set_zlabel('Z')plt.legend()plt.show()
# 데이터 결합 및 정제
X1_1 = sphere1[0,:]X2_1 = sphere1[1,:]X3_1 = sphere1[2,:]Y_1 = 100X1_2 = sphere2[0,:]X2_2 = sphere2[1,:]X3_2 = sphere2[2,:]Y_2 = 30X1 = np.concatenate((X1_1,X1_2), axis=0)X2 = np.concatenate((X2_1,X2_2), axis=0)X3 = np.concatenate((X3_1,X3_2), axis=0)X = np.vstack((X1, X2, X3))y = np.zeros(np.shape(X)[1])y[:1000] = Y_1y[1000:] = Y_2data = np.vstack((X, y))data = pd.DataFrame(data).Tdata.columns=["X1", "X2", "X3", "y"]X = data[["X1","X2","X3"]]y = data["y"]
# Standard(Linear) PCA #중요 내용이 아니기 때문에 패키지로 처리from sklearn.decomposition import PCAlinear_pca = PCA(n_components=2)linear_pca.fit(X)pca_X = linear_pca.transform(X)x1 = pca_X[:,0]x2 = pca_X[:,1]plt.figure(figsize=(6,6))plt.scatter(x1[y==100], x2[y==100], c='b', alpha=0.5, s=1, label='Class1')plt.scatter(x1[y==30], x2[y==30], c='r', alpha=0.5, s=1, label='Class2')plt.legend()plt.show()
# Kernel PCA with Polynomial kerneldef poly_kpca(X, degree, n_components, coef0=1): K = (np.dot(X, X.T) + coef0) ** degree N = K.shape[0] one_n = np.ones((N, N)) / N K = K - one_n.dot(K) - K.dot(one_n) + one_n.dot(K).dot(one_n) eigvals, eigvecs = eigh(K) pc = np.column_stack((eigvecs[:, -i] for i in range(1, n_components + 1))) return pcpc = poly_kpca(X, 5, 2) #degree=5로 설정x1 = pc[:,0]x2 = pc[:,1]plt.figure(figsize=(6,6))plt.scatter(x1[y==100], x2[y==100], c='b', alpha=0.5, label='Class1')plt.scatter(x1[y==30], x2[y==30], c='r', alpha=0.5, label='Class2')plt.legend()plt.show()
#분류 성능 확인 - Accuracy score : 0.8292from sklearn.model_selection import train_test_splitfrom sklearn.linear_model import SGDClassifierX_train, X_test, y_train, y_test = train_test_split(pc, y, test_size=0.2, shuffle=True)model = SGDClassifier(loss="hinge", penalty="l2", max_iter=100)result = model.fit(X_train, y_train)y_pred = result.predict(X_test)from sklearn.metrics import accuracy_scoreprint("Accuracy score : {0:0.4f}".format(accuracy_score(y_test, y_pred)))
# Kernel PCA with Gaussian kernel
from scipy.spatial.distance import pdist, squareformfrom scipy.linalg import eighdef rbf_kpca(X, gamma, n_components): distance = pdist(X, 'sqeuclidean') #샘플 간 유클리디안 거리^2 pairwise_distance = squareform(distance) #거리를 정방대칭행렬로 변환 K = np.exp(-pairwise_distance / (2 * (gamma**2))) #커널 행렬 계산 N = K.shape[0] one_n = np.ones((N, N)) / N K = K - one_n.dot(K) - K.dot(one_n) + one_n.dot(K).dot(one_n) #커널 행렬을 중앙에 맞춤 eigvals, eigvecs = eigh(K) #커널행렬의 고유값, 고유벡터를 오름차순으로 반환 pc = np.column_stack((eigvecs[:, -i] for i in range(1, n_components + 1))) return pcpc = rbf_kpca(X, gamma=20, n_components=2) #gamma=20으로 설정x1 = pc[:,0]x2 = pc[:,1]plt.figure(figsize=(6,6))plt.scatter(x1[y==100], x2[y==100], c='b', alpha=0.5, label='Class1')plt.scatter(x1[y==30], x2[y==30], c='r', alpha=0.5, label='Class2')plt.legend()plt.show()
# 분류 성능 확인 - Accuracy score : 0.9958
from sklearn.model_selection import train_test_splitfrom sklearn.linear_model import SGDClassifierX_train, X_test, y_train, y_test = train_test_split(pc, y, test_size=0.2, shuffle=True)model = SGDClassifier(loss="hinge", penalty="l2", max_iter=100)result = model.fit(X_train, y_train)y_pred = result.predict(X_test)from sklearn.metrics import accuracy_scoreprint("Accuracy score : {0:0.4f}".format(accuracy_score(y_test, y_pred)))
# gamma(hyperparameter) test#gamma=1일 때 0.8625에서 시작해 gamma=21에서 1.0에 수렴acc = []for i in range(1, 50): pc = rbf_kpca(X, gamma=i, n_components=2) X_train, X_test, y_train, y_test = train_test_split(pc, y, test_size=0.2, shuffle=True) model = SGDClassifier(loss="hinge", penalty="l2", max_iter=100) result = model.fit(X_train, y_train)y_pred = result.predict(X_test)accuracy = accuracy_score(y_test, y_pred)acc.append(accuracy)print(f'gamma={i} accuracy score={accuracy}')
