Statics / PCA, MCA, FAMD.md

PCA, MCA, FAMD

조회

셋 다 일단 통계분석이긴 한데... 거기까진 감이 좀 오시죠? 근데 저게 다 뭐시여...?

 

- PCA(Principal Component Analysis): 주성분분석
- MCA(Multiple Correspondence Analysis): 다중 대응분석
- FAMD(Factor Analysis of Mixed Data): 혼합 자료의 요인분석 (feat. 안티그래비티 번역)
이렇다.


일단 ADsP(혹은 빅분기나 ADP) 준비하셨던 분들 PCA 뭔지 들어보셨죠? 그죠. 그 차원 축소하는 그거. 나머지 둘은 PCA 칭구칭긔들… 그니까 똑같이 차원 축소 하는 애들이다. 아니 하나만 해도 골치아픈데 저런게 왜 세개나 있는데요? 그건 분석 들어가면서 설명드림.


import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns

from sklearn.decomposition import PCA # 주성분분석
from sklearn.preprocessing import StandardScaler # 을 하려면 스케일러가 필요합니다 
from sklearn.feature_selection import mutual_info_classif # MI (결과는 다른 분 걸 사용함)
from sklearn.preprocessing import LabelEncoder 

import prince # FAMD용

우리 시각화할거+데이터프레임 불러올거 있어서 씨본 판다스 불러야됩니다. 걍 저거 복붙하세요.

 

# 데이터가 쓸만한게 있던가... 
df = pd.read_csv('data/Ibuprofen.csv', sep=';')

# 어차피 분석할거라 대충 채울게요 
df.fillna({'Targets':0.0}, inplace=True)
df.fillna({'Bioactivities':0}, inplace=True)
df.dropna(inplace=True)

df.shape

그래프 설정은 알아서 하시고… 오늘은 그냥 이거 해보자! 라서 켐블에 있는거 대충 불러와서 결측값도 대충 때우고 다 날렸음. 다시 말하는거지만 이 글에서는 아 이런거구나만 익혀두세요.

 

# 수치형, 범주형 칼럼 분리

# 수치형 (정수, 실수 모두 포함)
numeric_columns = df.select_dtypes(include=['number']).columns

# 범주형 (문자열, Object 타입)
# categorical_columns = df.select_dtypes(include=['object']).columns <-이거 하면 이름같은거 다 껴서 안돼요 
# 착한 분석가 여러분들은 직접 확인하시고 픽하십셔
categorical_columns = ['Type', 'Passes Ro3', 'Structure Type']

그럼 이제 하나만 해도 머리 터지는걸 왜 세개나 하는지를 알려드리겠음. 셋 다 차원 축소하는 건 맞는데, PCA는 수치형만 받고 MCA는 범주형만 받는다. FAMD는 잡식성이라 수치형이고 범주형이고 안 가린다. 근데 우리가 데이터를 취급하면서 한쪽만 있을 수가 없잖아요? 그래요. 그겁니다.

 

그리고 저 중에서도 분석에 굳이? 싶은 데이터가 있을 수 있으니 착한 분석가 여러분들은 명세표를 꼭 확인하세요.


PCA

스케일러

# Scaler로 조정
scaler = StandardScaler()
numeric_scale = numeric_columns

transformed_numeric = scaler.fit_transform(df[numeric_scale]) # 변환된 데이터

transformed_numeric = pd.DataFrame(transformed_numeric, columns=numeric_scale) # 데이터프레임으로 변환

아잇 이게머야!!! 스케일러는 또 머예여!!! 이거 안 하면 아마 제 1 주성분에 분자량 들어갈거임. 걔가 제일 크잖아요. 그그 초유의 사태를 방지하기 위해 모든 칼럼의 값을 평균이 0, 표준편차가 1인 범위로 맞추는 게 스탠다드 스케일러다. 스케일러도 종류가 있음.

 

pca = PCA(n_components=0.90)
X_pca = pca.fit_transform(transformed_numeric)

print("성공! 최종 변수 개수:", transformed_numeric.shape[1])
print("압축된 주성분 개수:", X_pca.shape[1])

이거 하면 주성분분석 끝입니다. 저기 n_components에 0.9는 뭔 소리냐면 전체 분산의 90%를 확보하게끔 알아서 하십쇼라는 얘기. 실수 말고 숫자로 넣게 되면 주성분 개수가 됩니다.

 

결과 보기

# 1. 각 주성분의 설명력(분산 비율) 확인
explained_variance = pca.explained_variance_ratio_
cumulative_variance = np.cumsum(explained_variance)

print("--- PCA 결과 성적표 ---")
for i, ratio in enumerate(explained_variance):
    print(f"제 {i+1} 주성분: {ratio:.4f} (누적: {cumulative_variance[i]:.4f})")

# 2. 제1주성분(PC1)에 가장 큰 영향을 준 변수 top 10
# (주성분이 어떤 원래 변수들로 만들어졌는지 확인)
loadings = pd.DataFrame(pca.components_[0], index=numeric_columns, columns=['PC1_Loading'])
print("\n--- PC1을 구성하는 핵심 변수 TOP 10 ---")
print(loadings['PC1_Loading'].abs().sort_values(ascending=False).head(10))
--- PCA 결과 성적표 ---
제 1 주성분: 0.3561 (누적: 0.3561)
제 2 주성분: 0.1519 (누적: 0.5080)
제 3 주성분: 0.1303 (누적: 0.6383)
제 4 주성분: 0.0901 (누적: 0.7285)
제 5 주성분: 0.0732 (누적: 0.8017)
제 6 주성분: 0.0547 (누적: 0.8564)
제 7 주성분: 0.0418 (누적: 0.8982)
제 8 주성분: 0.0288 (누적: 0.9269)

--- PC1을 구성하는 핵심 변수 TOP 10 ---
Heavy Atoms          0.401871
HBA                  0.341076
QED Weighted         0.340439
Molecular Weight     0.315908
AlogP                0.296601
Np Likeness Score    0.290416
#RO5 Violations      0.285224
#Rotatable Bonds     0.254815
HBD                  0.230318
Aromatic Rings       0.213672
Name: PC1_Loading, dtype: float64

자, 이게 결과입니다. 악 내눈! 여기서 눈아프면 밑에꺼 어떻게 볼겨… 진정하십쇼. 원래 주성분분석에서 변수는 제 n 주성분으로 나오는데(이게 국룰임), 그 옆에 숫자 말고 괄호 안에 누적을 봅시다. 제 6 주성분 옆에 0.8564 보여요? 그게 무슨 뜻이냐면 얘가 설명력 85%를 갖기 위해서는 제 6 주성분까지 있어야 한다는 얘기다. 보통 많아야 2~3갠데 쟤는 너무 많아…

 

밑에 있는 건 제 1 주성분을 구성하는 핵심 변수들이라는 얘기. 위에도 썼지만 쟤는 수치형을 좋아해서 범주형 다 빼고 수치형만 돌린거다. 그럼 범주형은요? 더미화하면 돌릴 수는 있는데 결과 개발살나니까 하지 마세요.

 

Scree plot

# 데이터 준비 (제공해주신 로직 유지)
pc_labels = []
for i in range(len(pca.components_)):
    # 최신 Pandas/NumPy 대응을 위해 numeric_scale.columns 사용 권장
    top_feature = numeric_columns[np.argmax(np.abs(pca.components_[i]))]
    pc_labels.append(f"PC{i+1}\n({top_feature})")

exp_var_pca = pca.explained_variance_ratio_

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

# 1. 꺾은선 그래프 (Line + Marker)
sns.lineplot(x=pc_labels, y=exp_var_pca, marker='o', markersize=10, 
                color='purple', linewidth=2.5, label='Explained Variance')

# 2. 각 점 위에 수치 표시 (Annotate)
for i, val in enumerate(exp_var_pca):
    plt.text(i, val + 0.005, f'{val:.2f}', ha='center', va='bottom', fontsize=11)

plt.title('주성분별 설명력 및 대표 변수 (Scree Plot)', fontsize=15)
plt.ylabel('설명 가능한 분산 비율')
plt.xlabel('Principal Components')
plt.xticks(rotation=45)
plt.grid(axis='both', alpha=0.3) # 격자를 가로세로 모두 추가해서 흐름 파악 용이하게 함
plt.tight_layout()
plt.show()

아니 근데 6번까지 필요하다매 왜 제 2주성분에서 완만해지는건데… 저게 뭔데요? 저건 주성분분석 결과를 시각화한 스크리 플롯이라는겁니다. 그 꺾은선이 완만해지는 구간을 엘보라고 하는데, 저 엘보에서 꺾어지는 구간까지 갖다 쓰거든요… 엘보가 PC2야 쟤는…

 

Loading plot

# 로딩 값 추출
loadings = pca.components_.T * np.sqrt(pca.explained_variance_)
loading_df = pd.DataFrame(loadings, index=numeric_columns, columns=[f'PC{i+1}' for i in range(len(pca.components_))])

plt.figure(figsize=(10, 8))
plt.axhline(0, color='black', linewidth=1)
plt.axvline(0, color='black', linewidth=1)

# 주요 변수들만 화살표로 표시 (너무 많으면 복잡하므로 상위 기여 변수 위주)
top_features = loading_df.abs().sum(axis=1).sort_values(ascending=False).head(15).index

for feature in top_features:
    plt.arrow(0, 0, loading_df.loc[feature, 'PC1'], loading_df.loc[feature, 'PC2'], 
                color='r', alpha=0.5, head_width=0.02)
    plt.text(loading_df.loc[feature, 'PC1']*1.15, loading_df.loc[feature, 'PC2']*1.15, 
                feature, color='g', ha='center', va='center', fontsize=10)

plt.xlabel('PC1 (분자의 크기와 지질친화성)')
plt.ylabel('PC2 (극성 및 활성 지표)')
plt.title('변수별 로딩 플롯 (PC1 vs PC2)', fontsize=15)
plt.grid(alpha=0.3)
plt.xlim(-1, 1)
plt.ylim(-1, 1)
plt.tight_layout()
plt.show()

이건 각 칼럼들이 주성분에 미치는 영향이다. 무기 플래그는 별 상관 없는거냐… 아니면 범주형인데 저기 낀거냐… 아무튼, 가로축이 제 1주성분이고 세로축이 제 2주성분이다. 그래서 어떻게 보느냐…

 

1. 화살표의 크기: 해당 요소가 주성분에 얼마나 영향을 미치는가
2. 화살표의 방향: 해당 요소가 주성분과 어떤 상관관계인가(수직이면 노상관)

 

이렇게 보시면 된다. 저기 Heavy Atom 보이시죠? 쟤는 x축이랑 거의 근접한 대신 y축이랑은 거의 수직인데, 제 2 주성분에 1도 영향을 안 끼친다는 얘기다. 반대로 Polar surface area는 제 1 주성분에 1도 영향을 안 끼친다.

 

MCA

mca = prince.MCA(
    n_components=2,
    n_iter=3,
    copy=True,
    check_input=True,
    engine='sklearn',
    random_state=42
)

# 이친구는 데이터프레임만 취급합니다. 
mca = mca.fit(df[categorical_columns])

전처리 위에서 했으니까 걍 돌리십쇼. 근데 쟤는 칼럼 말고 데이터프레임 상태로 줘야 한다.

 

결과 확인하기

mca_coords = mca.row_coordinates(df[categorical_columns])

print(f'Coordinates of the rows (MCA): {len(mca_coords)}')
mca_coords
Coordinates of the rows (MCA): 59
0	1
1	-0.108148	0.57735
2	-0.108148	0.57735
3	-0.108148	0.57735
4	-0.108148	0.57735
5	-0.108148	0.57735
6	-0.108148	0.57735
7	-0.108148	0.57735
8	-0.108148	0.57735
9	-0.108148	0.57735
11	-0.108148	0.57735
12	-0.108148	0.57735
13	-0.108148	0.57735
14	-0.108148	0.57735
15	-0.108148	0.57735
16	-0.108148	0.57735
17	-0.108148	0.57735
18	-0.108148	0.57735
19	-0.108148	0.57735
20	-0.108148	0.57735
21	-0.108148	0.57735
22	-0.108148	0.57735
23	-0.108148	0.57735
24	-0.108148	0.57735
25	-0.108148	0.57735
26	-0.108148	0.57735
27	-0.108148	0.57735
28	-0.108148	0.57735
31	-0.108148	0.57735
32	-0.108148	0.57735
33	-0.108148	0.57735
34	-0.108148	0.57735
35	-0.108148	0.57735
36	-0.108148	0.57735
37	-0.108148	0.57735
38	-0.108148	0.57735
39	-0.108148	0.57735
40	-0.108148	0.57735
41	-0.108148	0.57735
42	-0.108148	0.57735
43	-0.108148	0.57735
44	-0.108148	0.57735
45	-0.108148	0.57735
46	-0.108148	0.57735
47	-0.108148	0.57735
48	-0.108148	0.57735
49	3.082207	0.57735
50	-0.108148	0.57735
51	-0.108148	0.57735
52	-0.108148	0.57735
53	-0.108148	0.57735
54	-0.108148	0.57735
55	-0.108148	0.57735
56	-0.108148	0.57735
57	-0.108148	0.57735
58	-0.108148	0.57735
59	-0.108148	0.57735
60	-0.108148	0.57735
61	3.082207	0.57735
62	-0.108148	0.57735

내가 위에도 썼지만 얘는 MCA를 하기에 그렇게 좋은 데이터가 아님... 범주향 데이터의 대부분이 이름이나 아이디같은거라 다 날렸더니 칼럼 세개 남았거든요? 그래서 좌표고 시각화고 다 개발살났습니다. 시각화는 뭐 걍 스캐터플롯임.

 

FAMD

아마 당신들이 대부분 하게 될지도 모를 분석. 대부분의 데이터는 수치형과 범주형이 섞여있습니다… 퓨어한거 없엉.. 머신러닝용으로 쓰는거면 뭐 0 1 이런 식으로는 표현할듯.

 

# 1. 인덱스 때문에 꼬이는 경우가 많으니 복사본 생성 후 인덱스 초기화
famd_input = df[list(numeric_columns) + list(categorical_columns)].copy().reset_index(drop=True)

# 2. 혹시 모르니 모든 수치형 칼럼을 float로 강제 변환
for col in numeric_columns:
    famd_input[col] = famd_input[col].astype(float)

# 3. FAMD 다시 실행
famd = prince.FAMD(
    n_components=2,
    n_iter=5, # 반복 횟수를 조금 늘려주면 안정적입니다
    engine='sklearn',
    random_state=42
)

# 피팅 시도
try:
    famd = famd.fit(famd_input)
    print("성공! 이제 좌표를 뽑을 수 있습니다.")
except ValueError as e:
    print(f"아직 에러가 나네요: {e}")
    # 만약 또 에러가 나면, 범주형 칼럼이 너무 적어서(3개) 발생할 수 있습니다.

이제 머리 아플 예정이니까 나가서 바람 쐬고 오세요.

 

결과 확인-차원 좌표와 기여도

famd_coords = famd.row_coordinates(famd_input)

famd_coords
component	0	1
0	0.761347	-1.055728
1	6.705050	2.673633
2	-0.433802	-1.758587
3	1.179095	-1.175229
4	-2.277754	0.778783
5	-1.286728	-0.835628
6	0.496398	0.143235
7	2.966788	0.026485
8	-0.907002	-0.248847
9	-1.577409	0.312389
10	0.422361	-0.973843
11	7.586716	2.877328
12	0.208778	-1.126320
13	1.428615	-1.178440
14	-1.562023	0.220391
15	-1.066529	-0.147865
16	0.279641	-1.424065
17	-1.709486	0.209354
18	-0.384484	-0.977643
19	5.725015	0.224686
20	-0.112680	-1.598461
21	2.437106	0.156564
22	-2.229055	0.113714
23	-0.631329	-0.466001
24	1.032493	-1.589061
25	-1.938724	-2.142537
26	0.055179	-3.214237
27	-0.597261	-1.404748
28	-0.289858	-1.695625
29	2.887685	-1.176094
30	3.408455	-0.040365
31	0.422582	-0.839601
32	0.966494	-1.003540
33	-2.570216	1.060088
34	2.334272	2.386978
35	-1.956972	0.799100
36	0.684912	0.122746
37	-1.152568	-0.428253
38	0.871038	-0.004118
39	-1.859993	1.739152
40	2.510349	2.405997
41	-1.043532	0.065698
42	-1.408376	-0.605283
43	-0.993196	-0.338999
44	-2.271696	2.458735
45	-0.996774	-1.663435
46	-4.436712	5.660998
47	-1.282504	0.375607
48	-2.896505	2.160478
49	-1.540748	0.864644
50	-1.309301	0.532292
51	-0.963285	0.327781
52	-0.500237	-0.086713
53	1.107287	1.095981
54	-0.515824	-0.911213
55	2.982609	-0.894089
56	-1.430050	1.204734
57	-2.196237	0.279800
58	-1.131413	-0.272804

악 내눈! 저기 위에 MCA랑 똑같이 쟤도 좌표가 나온다. 그럼 저걸로 차원별 기여도 보나요? 아뇨. 아, 차원이 뭐냐고? 그걸 설명하려면 일단 선형대수학을 배우셔야...

 

# 1. 설명력 (%) - 소수점 첫째 자리까지 + '%' 붙이기
explain_rounded = [f"{v:.1f}%" for v in famd.percentage_of_variance_]

# 2. 원점수 (고윳값) - 소수점 셋째 자리까지 깔끔하게
point_rounded = [round(v, 3) for v in famd.eigenvalues_]

print(f'해당 차원의 설명력: {explain_rounded}')
print(f'해당 차원의 원점수: {point_rounded}')
해당 차원의 설명력: ['32.3%', '13.7%']
해당 차원의 원점수: [np.float64(5.004), np.float64(2.13)]

0차원의 설명력이 32.3%, 1차원의 설명력이 13.7%란다. 밑에는 원점수임.

 

결과 확인-차원별 기여도

# 각 칼럼들이 0차원에 기여하는 순서대로 정렬
df_comp0 = famd.column_contributions_.sort_values(by=0, ascending=False)

# 해서 Top 10만 
df_comp0[[0]].head(10)
component	0
variable	
Heavy Atoms	0.160621
QED Weighted	0.114318
HBA	0.114028
Molecular Weight	0.099723
AlogP	0.088948
Np Likeness Score	0.083103
#RO5 Violations	0.081003
#Rotatable Bonds	0.065774
HBD	0.052470
Aromatic Rings	0.045718

제 0차원에 기여하는 칼럼들의 순위. 이걸 보고 0차원이 뭔지를 유추해야 한다. 얘는 헤비아톰에 분자량에 있는 걸 보니 분자 크기 관련된 차원인 듯 하다.

 

# 각 칼럼들이 1차원에 기여하는 순서대로 정렬
df_comp1 = famd.column_contributions_.sort_values(by=1, ascending=False)

# Top 10
df_comp1[[1]].head(10)
component	1
variable	
Polar Surface Area	0.229303
Bioactivities	0.179545
Targets	0.149604
HBD	0.096879
Max Phase	0.093455
#Rotatable Bonds	0.073360
#RO5 Violations	0.047313
Aromatic Rings	0.044396
HBA	0.029810
AlogP	0.024106

이건 제 1차원에 기여하는 칼럼들의 순위이다. 보니까 반응성이나 생리활성 관련인가배?

 

plt.figure(figsize = (20, 10))

plt.subplot(1, 2, 1)
famd.column_contributions_[0].sort_values().plot(kind='barh', figsize=(8, 6), color='#00498c')
plt.title('Component 0에 대한 변수 기여도')
plt.xlabel('Contribution (%)')

plt.subplot(1, 2, 2)
famd.column_contributions_[1].sort_values().plot(kind='barh', figsize=(8, 6), color='#00498c')
plt.title('Component 1에 대한 변수 기여도')
plt.xlabel('Contribution (%)')
plt.tight_layout()
plt.show()

원래 따로 된 것만 있었는데 한눈에 보려고 합쳤음… 위에 그 0차원 1차원 기여도를 막대그래프로 시각화 한 결과다.

 

Scatter plot

# 나눔스퀘어 폰트는 이미 설정되어 있으니 바로 사용합니다.
plt.figure(figsize=(12, 9))

# 1. 산점도 그리기
scatter = sns.scatterplot(
    x=famd_coords[0], 
    y=famd_coords[1], 
    style=df['Type'], 
    hue=df['Molecular Weight'],
    s=150,           # 점 크기 확대
    alpha=0.7,       # 투명도 조절로 겹친 부분 확인
    edgecolor='w',   # 점 테두리 흰색으로 구분감 상승
    palette='viridis' # 세련된 색감 적용
)

# 2. 축 이름 지어주기 (기여도 분석 결과 반영)
# 설명력(percentage_of_variance_)을 소수점 첫째 자리까지 포함
plt.xlabel(f'Dimension 0: 화합물의 체급 및 복잡도 ({famd.percentage_of_variance_[0]:.1f}%)', fontsize=12)
plt.ylabel(f'Dimension 1: 극성 및 생물학적 활성 프로파일 ({famd.percentage_of_variance_[1]:.1f}%)', fontsize=12)

# 3. 보조 가이드라인 (원점 기준)
plt.axhline(0, color='gray', linestyle='--', alpha=0.3)
plt.axvline(0, color='gray', linestyle='--', alpha=0.3)

# 4. 제목 및 범례 설정
plt.title('FAMD 최종 분석: ChEMBL 화합물의 다차원 지형도', fontsize=16, pad=20)
plt.legend(bbox_to_anchor=(1.05, 1), loc='upper left', title='변수 구분')

plt.grid(True, linestyle=':', alpha=0.4)
plt.tight_layout()
plt.show()

저 위에 좌표 구한 거 있죠? 그걸로 산점도 그린거다. 점 색이랑 모양은 다른 인자로 구별할 수는 있는데, 연속형이시면 범주화 함 하시는게 정신건강에 이롭습니다.

댓글

홈으로 돌아가기

검색 결과

"search" 검색 결과입니다.