주성분 분석 (Principal Component Analysis)
상관관계가 있는 변수들을 선형결합하여 변수를 축약하는 기법. 요인 분석의 한 종류.
첫번째 주성분이 전체 변동을 가장 많이 설명하고, 두번째 주성분은 첫번째 주성분과 상관성이 낮아 첫번째 주성분이 설명하지 못하는 나머지 변동을 정보의 손실없이 가장 많이 설명할 수 있도록 변수들을 조합.
1. 목적
- 여러 변수들 간에 내재하는 상관관계, 연관성을 이용해 소수의 주성분으로 차원을 축소
- 다중공선성이 존재하는 경우, 상관성이 적은 주성분으로 변수들을 축소하여 모형 개발에 활용
- 주성분분석을 통해 차원을 축소한 후 군집분석을 수행하면 결과와 연산속도를 개선할 수 있음
- 다량의 센서 데이터를 주성분분석으로 차원 축소한 후 시계열로 분포나 추세의 변화를 분석하여 고장 징후를 사전에 파악
2. 공분산 행렬 / 상관계수 행렬
주성분 분석의 문제는 척도에 영향을 받는다는 점이다. 변수들의 선형결합을 유도할 때 분산을 이용하기 때문에 결과적으로 공분산 행렬로부터 유도된 주성분은 측정단위의 크기에 좌우된다.
설문조사처럼 모든 변수들이 같은 수준으로 점수화가 된 경우에는 공분산 행렬을 사용해도 되지만 변수들의 scale이 서로 많이 다른 경우에는 값이 큰 특정 변수가 전체적인 경향을 좌우하기 때문에 상관계수 행렬을 사용하여 추출해야 한다.
3. 주성분 선택법
주성분분석 결과에서 누적기여율(cumulative proportion)이 85%이상인 주성분까지 선택
screen plot에서 고유값이 수평을 유지하기 전단계로 주성분의 수를 선택
library(MVA)
df <- read.csv("data/open_closed.csv")
head(df)
## mechanics_C vectors_C algebra_O analysis_O statistics_O
## 1 77 82 67 67 81
## 2 63 78 80 70 81
## 3 75 73 71 66 81
## 4 55 72 63 70 68
## 5 63 63 65 70 63
## 6 53 61 72 64 73
library(psych)
pairs.panels(df)
# 설문조사처럼 같은 scale 점수화가 된 경우에는 공분산행렬 사용해도 되지만
# 변수들의 scale이 많이 다른 경우 특정 변수가 전체적인 경향을 좌우하기 때문에
# 상관계수 행렬을 사용하여 분석하는 것이 좋다.
pca_c <- prcomp(df) # default : 공분산 행렬을 이용한 선형 결합
summary(pca_c)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5
## Standard deviation 26.2105 14.2166 10.1856 9.19948 5.67039
## Proportion of Variance 0.6191 0.1821 0.0935 0.07627 0.02898
## Cumulative Proportion 0.6191 0.8013 0.8948 0.97102 1.00000
pca_c
## Standard deviations:
## [1] 26.210490 14.216577 10.185642 9.199481 5.670387
##
## Rotation:
## PC1 PC2 PC3 PC4 PC5
## mechanics_C -0.5054457 -0.74874751 0.2997888 -0.296184264 -0.07939388
## vectors_C -0.3683486 -0.20740314 -0.4155900 0.782888173 -0.18887639
## algebra_O -0.3456612 0.07590813 -0.1453182 0.003236339 0.92392015
## analysis_O -0.4511226 0.30088849 -0.5966265 -0.518139724 -0.28552169
## statistics_O -0.5346501 0.54778205 0.6002758 0.175732020 -0.15123239
pca_r <- prcomp(df, scale = T) # scale : 상관계수행렬 이용
summary(pca_r)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5
## Standard deviation 1.7835 0.8600 0.66706 0.62281 0.49658
## Proportion of Variance 0.6362 0.1479 0.08899 0.07758 0.04932
## Cumulative Proportion 0.6362 0.7841 0.87310 0.95068 1.00000
# 주성분 갯수의 선택
# (1) 주성분의 누적 중요도(Cumulative Proportion)가 70 ~ 90 % 사이에서 선택.
# (2) 평균 고유값보다 작은 고유값을 갖는 주성분을 버림.
# 상관계수 행렬 이용시 평균 분산(표준편차) = 1 이므로 Standard deviation이 1 또는 0.7 보다 작은 것은 버림.
# (3) Scree diagram에서 기울기가 완만해지는 시점, 즉 팔꿈치에 대응하는 지점까지 선택.
pca_r
## Standard deviations:
## [1] 1.7835302 0.8599836 0.6670571 0.6228101 0.4965788
##
## Rotation:
## PC1 PC2 PC3 PC4 PC5
## mechanics_C -0.3996045 -0.6454583 0.62078249 -0.1457865 -0.1306722
## vectors_C -0.4314191 -0.4415053 -0.70500628 0.2981351 -0.1817479
## algebra_O -0.5032816 0.1290675 -0.03704901 -0.1085987 0.8466894
## analysis_O -0.4569938 0.3879057 -0.13618182 -0.6662561 -0.4221885
## statistics_O -0.4382444 0.4704545 0.31253342 0.6589164 -0.2340223
# pc1 : 전체적으로 잘하고 못하는 학생들을 분리
# pc2 : mechanics 는 잘하지만 statistics는 못하는 특성의 학생들 분리
# scaling 했기 때문에 각 주성분의 분산 = 1, 전체 총분산 = 5 (변수 5개)
pca_r$sdev # 고유값 = 각 주성분의 분산(Standard deviation의 제곱)
## [1] 1.7835302 0.8599836 0.6670571 0.6228101 0.4965788
pca_r$rotation[, 1] # 고유벡터 = 각 주성분의 Rotation 값
## mechanics_C vectors_C algebra_O analysis_O statistics_O
## -0.3996045 -0.4314191 -0.5032816 -0.4569938 -0.4382444
pc1 <- pca_r$x[, 1]
pc2 <- pca_r$x[, 2]
cor(pc1, pc2) # pc1, pc2의 공분산은 0 에 수렴
## [1] -8.547928e-17
df_s <- scale(df)
summary(df_s)
## mechanics_C vectors_C algebra_O analysis_O
## Min. :-2.2277 Min. :-3.16354 Min. :-3.35087 Min. :-2.5383
## 1st Qu.:-0.5121 1st Qu.:-0.65345 1st Qu.:-0.52728 1st Qu.:-0.7364
## Median : 0.1456 Median : 0.03112 Median :-0.05669 Median : 0.1562
## Mean : 0.0000 Mean : 0.00000 Mean : 0.00000 Mean : 0.0000
## 3rd Qu.: 0.5888 3rd Qu.: 0.71569 3rd Qu.: 0.62568 3rd Qu.: 0.6951
## Max. : 2.1757 Max. : 2.38908 Max. : 2.76690 Max. : 1.5708
## statistics_O
## Min. :-1.9302
## 1st Qu.:-0.6553
## Median :-0.1337
## Mean : 0.0000
## 3rd Qu.: 0.5328
## Max. : 2.2424
rot1 <- pca_r$rotation[, 1]
rot1
## mechanics_C vectors_C algebra_O analysis_O statistics_O
## -0.3996045 -0.4314191 -0.5032816 -0.4569938 -0.4382444
plot(df_s%*%rot1, pca_r$x[, 1]) # 완벽한 상관관계를 보인다.
screeplot(pca_r, type = "l") # 떨어지는 각도가 완만해지는 2까지 주성분으로 선택.
biplot(pca_r)
# 행렬도(biplot)
# - 원변수와 PC간의 관계를 그래프로 표현
# - 화살표는 원변수와 PC의 상관계수. PC와 평행할수록 해당 PC에 큰 영향.
# - 화살표 벡터의 길이가 원변수의 분산을 표현
# Example
df <- heptathlon # 여자 철인 7종 경기
df
## hurdles highjump shot run200m longjump javelin
## Joyner-Kersee (USA) 12.69 1.86 15.80 22.56 7.27 45.66
## John (GDR) 12.85 1.80 16.23 23.65 6.71 42.56
## Behmer (GDR) 13.20 1.83 14.20 23.10 6.68 44.54
## Sablovskaite (URS) 13.61 1.80 15.23 23.92 6.25 42.78
## Choubenkova (URS) 13.51 1.74 14.76 23.93 6.32 47.46
## Schulz (GDR) 13.75 1.83 13.50 24.65 6.33 42.82
## Fleming (AUS) 13.38 1.80 12.88 23.59 6.37 40.28
## Greiner (USA) 13.55 1.80 14.13 24.48 6.47 38.00
## Lajbnerova (CZE) 13.63 1.83 14.28 24.86 6.11 42.20
## Bouraga (URS) 13.25 1.77 12.62 23.59 6.28 39.06
## Wijnsma (HOL) 13.75 1.86 13.01 25.03 6.34 37.86
## Dimitrova (BUL) 13.24 1.80 12.88 23.59 6.37 40.28
## Scheider (SWI) 13.85 1.86 11.58 24.87 6.05 47.50
## Braun (FRG) 13.71 1.83 13.16 24.78 6.12 44.58
## Ruotsalainen (FIN) 13.79 1.80 12.32 24.61 6.08 45.44
## Yuping (CHN) 13.93 1.86 14.21 25.00 6.40 38.60
## Hagger (GB) 13.47 1.80 12.75 25.47 6.34 35.76
## Brown (USA) 14.07 1.83 12.69 24.83 6.13 44.34
## Mulliner (GB) 14.39 1.71 12.68 24.92 6.10 37.76
## Hautenauve (BEL) 14.04 1.77 11.81 25.61 5.99 35.68
## Kytola (FIN) 14.31 1.77 11.66 25.69 5.75 39.48
## Geremias (BRA) 14.23 1.71 12.95 25.50 5.50 39.64
## Hui-Ing (TAI) 14.85 1.68 10.00 25.23 5.47 39.14
## Jeong-Mi (KOR) 14.53 1.71 10.83 26.61 5.50 39.26
## Launa (PNG) 16.42 1.50 11.78 26.16 4.88 46.38
## run800m score
## Joyner-Kersee (USA) 128.51 7291
## John (GDR) 126.12 6897
## Behmer (GDR) 124.20 6858
## Sablovskaite (URS) 132.24 6540
## Choubenkova (URS) 127.90 6540
## Schulz (GDR) 125.79 6411
## Fleming (AUS) 132.54 6351
## Greiner (USA) 133.65 6297
## Lajbnerova (CZE) 136.05 6252
## Bouraga (URS) 134.74 6252
## Wijnsma (HOL) 131.49 6205
## Dimitrova (BUL) 132.54 6171
## Scheider (SWI) 134.93 6137
## Braun (FRG) 142.82 6109
## Ruotsalainen (FIN) 137.06 6101
## Yuping (CHN) 146.67 6087
## Hagger (GB) 138.48 5975
## Brown (USA) 146.43 5972
## Mulliner (GB) 138.02 5746
## Hautenauve (BEL) 133.90 5734
## Kytola (FIN) 133.35 5686
## Geremias (BRA) 144.02 5508
## Hui-Ing (TAI) 137.30 5290
## Jeong-Mi (KOR) 139.17 5289
## Launa (PNG) 163.43 4566
pairs.panels(df)
# 달리기는 숫자가 작은 것이 좋은 것이기 때문에 반대로 변환.
df$hurdles <- with(df, max(hurdles) - hurdles)
df$run200m <- with(df, max(run200m) - run200m)
df$run800m <- with(df, max(run800m) - run800m)
df
## hurdles highjump shot run200m longjump javelin
## Joyner-Kersee (USA) 3.73 1.86 15.80 4.05 7.27 45.66
## John (GDR) 3.57 1.80 16.23 2.96 6.71 42.56
## Behmer (GDR) 3.22 1.83 14.20 3.51 6.68 44.54
## Sablovskaite (URS) 2.81 1.80 15.23 2.69 6.25 42.78
## Choubenkova (URS) 2.91 1.74 14.76 2.68 6.32 47.46
## Schulz (GDR) 2.67 1.83 13.50 1.96 6.33 42.82
## Fleming (AUS) 3.04 1.80 12.88 3.02 6.37 40.28
## Greiner (USA) 2.87 1.80 14.13 2.13 6.47 38.00
## Lajbnerova (CZE) 2.79 1.83 14.28 1.75 6.11 42.20
## Bouraga (URS) 3.17 1.77 12.62 3.02 6.28 39.06
## Wijnsma (HOL) 2.67 1.86 13.01 1.58 6.34 37.86
## Dimitrova (BUL) 3.18 1.80 12.88 3.02 6.37 40.28
## Scheider (SWI) 2.57 1.86 11.58 1.74 6.05 47.50
## Braun (FRG) 2.71 1.83 13.16 1.83 6.12 44.58
## Ruotsalainen (FIN) 2.63 1.80 12.32 2.00 6.08 45.44
## Yuping (CHN) 2.49 1.86 14.21 1.61 6.40 38.60
## Hagger (GB) 2.95 1.80 12.75 1.14 6.34 35.76
## Brown (USA) 2.35 1.83 12.69 1.78 6.13 44.34
## Mulliner (GB) 2.03 1.71 12.68 1.69 6.10 37.76
## Hautenauve (BEL) 2.38 1.77 11.81 1.00 5.99 35.68
## Kytola (FIN) 2.11 1.77 11.66 0.92 5.75 39.48
## Geremias (BRA) 2.19 1.71 12.95 1.11 5.50 39.64
## Hui-Ing (TAI) 1.57 1.68 10.00 1.38 5.47 39.14
## Jeong-Mi (KOR) 1.89 1.71 10.83 0.00 5.50 39.26
## Launa (PNG) 0.00 1.50 11.78 0.45 4.88 46.38
# delete outlier : 모든 종목에서 이상치 보이는 Launa (PNG) 삭제. hurdles = 0 (불참).
df2 <- df[df$hurdles > 0, ]
pairs.panels(df2)
# 주성분 분석
pca <- prcomp(df2[, -8], scale = T) # score 변수 제외
summary(pca)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6
## Standard deviation 2.0793 0.9482 0.9109 0.68320 0.54619 0.33745
## Proportion of Variance 0.6177 0.1284 0.1185 0.06668 0.04262 0.01627
## Cumulative Proportion 0.6177 0.7461 0.8646 0.93131 0.97392 0.99019
## PC7
## Standard deviation 0.26204
## Proportion of Variance 0.00981
## Cumulative Proportion 1.00000
pca
## Standard deviations:
## [1] 2.0793370 0.9481532 0.9109016 0.6831967 0.5461888 0.3374549 0.2620420
##
## Rotation:
## PC1 PC2 PC3 PC4 PC5
## hurdles -0.4503876 0.05772161 -0.1739345 0.04840598 -0.19889364
## highjump -0.3145115 -0.65133162 -0.2088272 -0.55694554 0.07076358
## shot -0.4024884 -0.02202088 -0.1534709 0.54826705 0.67166466
## run200m -0.4270860 0.18502783 0.1301287 0.23095946 -0.61781764
## longjump -0.4509639 -0.02492486 -0.2697589 -0.01468275 -0.12151793
## javelin -0.2423079 -0.32572229 0.8806995 0.06024757 0.07874396
## run800m -0.3029068 0.65650503 0.1930020 -0.57418128 0.31880178
## PC6 PC7
## hurdles 0.84665086 -0.06961672
## highjump -0.09007544 0.33155910
## shot -0.09886359 0.22904298
## run200m -0.33279359 0.46971934
## longjump -0.38294411 -0.74940781
## javelin 0.07193437 -0.21108138
## run800m -0.05217664 0.07718616
plot(pca, type = "l") # PC2까지 주성분으로 선택
plot(df2$score, pca$x[,1])
cor(df2$score, pca$x[,1])
## [1] -0.9931168
# PC1과 score는 -1 의 상관관계. 전 종목에서 잘하고 못하는 선수 구분
biplot(pca)
# pc1 축과 평행일수록 상관관계가 높고 roatation 값이 높다.
# pc2 축과 수직인 변수들은 별로 영향력이 없음을 의미한다.
# pc2 : 높이뛰기는 못하지만 장거리를 잘하는 선수들이라고 할 수 있다.
biplot(pca, cex = 0.8, choices = c(1,3)) # pc1 & pc3
qqnorm(pca$x[,1])
qqline(pca$x[,1])
# Practice 2
# Bulls : 경매시장에서 거래된 76마리의 2살 이하 황소의 특성과 거래가격(SalePr)에 관한 자료
# SalePr와 Breed 변수를 제외한 7개의 변수를 사용해 주성분분석을 시행하여 아래의 질문에 답하시오.
# 공분산 행렬과 상관계수 행렬을 사용하여 각각 분석하고 비교하시오.
bulls <- read.csv("data/bulls.csv", header = T)
head(bulls)
## Breed SalePr YrHgt FtFrBody PrctFFB Frame BkFat SaleHt SaleWt
## 1 1 2200 51.0 1128 70.9 7 0.25 54.8 1720
## 2 1 2250 51.9 1108 72.1 7 0.25 55.3 1575
## 3 1 1625 49.9 1011 71.6 6 0.15 53.1 1410
## 4 1 4600 53.1 993 68.9 8 0.35 56.4 1595
## 5 1 2150 51.2 996 68.6 7 0.25 55.0 1488
## 6 1 1225 49.2 985 71.4 6 0.15 51.4 1500
bullsV7 <- bulls[, -c(1,2)] # Breed, SalePr 변수 제거
head(bullsV7)
## YrHgt FtFrBody PrctFFB Frame BkFat SaleHt SaleWt
## 1 51.0 1128 70.9 7 0.25 54.8 1720
## 2 51.9 1108 72.1 7 0.25 55.3 1575
## 3 49.9 1011 71.6 6 0.15 53.1 1410
## 4 53.1 993 68.9 8 0.35 56.4 1595
## 5 51.2 996 68.6 7 0.25 55.0 1488
## 6 49.2 985 71.4 6 0.15 51.4 1500
library(psych)
pairs.panels(bullsV7) # BkFat - 한쪽에 치우져 있어서 변환이 필요한 것처럼 보인다.
cor(bullsV7$YrHgt, bullsV7$BkFat) # -0.34
## [1] -0.344277
cor(bullsV7$YrHgt, log(bullsV7$BkFat)) # -0.41
## [1] -0.4063381
cor(bullsV7$SaleWt, bullsV7$BkFat) # 0.21
## [1] 0.2075349
cor(bullsV7$SaleWt, log(bullsV7$BkFat)) # 0.15
## [1] 0.1543811
# 로그 변환 했으나 correlation 증가가 미미하므로 변환하지 않고 분석을 진행하는 것으로 한다.
# 1. 각 주성분의 표준편차와 그 주성분을 계산하는데 사용된 rotation값을 찾으시오.
# (1) 공분산 행렬 이용
pca_cov <- prcomp(bullsV7)
pca_cov
## Standard deviations:
## [1] 143.45596038 69.81887125 2.33005787 1.82107340 0.68471173
## [6] 0.27212808 0.06722679
##
## Rotation:
## PC1 PC2 PC3 PC4
## YrHgt -5.887328e-03 -0.0096800709 -0.286337289 -0.608787152
## FtFrBody -4.870470e-01 -0.8726966457 0.034277115 0.003226954
## PrctFFB -8.526499e-03 -0.0291964492 -0.904388519 0.425174911
## Frame -3.111988e-03 -0.0048861100 -0.133266834 -0.311194400
## BkFat -6.919922e-05 0.0004925452 0.018864084 0.005278296
## SaleHt -9.329509e-03 -0.0085770135 -0.284214793 -0.593037047
## SaleWt -8.732589e-01 0.4871927200 -0.004846824 0.005597435
## PC5 PC6 PC7
## YrHgt -0.5355689528 -0.5097273178 0.0245917521
## FtFrBody -0.0004437402 -0.0004566049 -0.0002530995
## PrctFFB -0.0083876301 0.0103890723 0.0142927590
## Frame -0.3905733600 0.8552041268 -0.0379840767
## BkFat -0.0119061237 0.0437862261 0.9987777777
## SaleHt 0.7485979836 0.0823314748 0.0138200628
## SaleWt -0.0026647979 -0.0003410092 -0.0002556156
# (2) 상관계수 행렬 이용
pca <- prcomp(bullsV7, scale = T)
pca
## Standard deviations:
## [1] 2.0299502 1.1563431 0.8610357 0.6491727 0.4310521 0.3827563 0.2169256
##
## Rotation:
## PC1 PC2 PC3 PC4 PC5
## YrHgt -0.4499313 -0.042790217 -0.41570891 0.1133565 -0.06587066
## FtFrBody -0.4123256 0.129836547 0.45029241 0.2474787 0.71934339
## PrctFFB -0.3555618 -0.315507785 0.56827313 0.3147874 -0.57936738
## Frame -0.4339569 0.007728211 -0.45234503 0.2428179 -0.14299538
## BkFat 0.1867048 0.714719363 -0.03873196 0.6181171 -0.16023789
## SaleHt -0.4528538 0.101315086 -0.17665043 -0.2157694 0.10953536
## SaleWt -0.2699470 0.600514834 0.25331192 -0.5824327 -0.29054729
## PC6 PC7
## YrHgt 0.07223418 -0.77492612
## FtFrBody 0.17706072 -0.01776760
## PrctFFB -0.12780009 0.00239740
## Frame 0.43414400 0.58233705
## BkFat -0.20801720 -0.04244214
## SaleHt -0.79928778 0.23672329
## SaleWt 0.27656055 -0.04703601
# 주성분의 표준편차(Standard deviation) = 고유값
# 주성분을 계산하는데 사용된 rotation = 고유벡터
# 2. 적절한 주성분의 개수를 선택하고 근거를 설명하시오.
# (1) 공분산 행렬 이용
plot(pca_cov, type = "l")
# plot에서 제3 주성분부터 기울기가 완만하게 변하고 그 이후로는 거의 의미가 없기 때문에
# 제3 주성분까지 3개의 주성분을 선택한다.
summary(pca_cov)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6
## Standard deviation 143.4560 69.8189 2.33006 1.82107 0.68471 0.2721
## Proportion of Variance 0.8082 0.1914 0.00021 0.00013 0.00002 0.0000
## Cumulative Proportion 0.8082 0.9996 0.99985 0.99998 1.00000 1.0000
## PC7
## Standard deviation 0.06723
## Proportion of Variance 0.00000
## Cumulative Proportion 1.00000
# 각 주성분의 중요도(Proportion of Variance) 누적합인 Cumulative Proportion을 보면
# 제1, 제2 주성분 만으로도 전체 데이터의 99.96%를 설명할 수 있다.
# 그러므로 2개의 주성분만 선택한다.
# (2) 상관계수 행렬 이용
plot(pca, type = "l")
# plot에서는 기울기가 완만하게 변하는 팔꿈치 부분을 명확하게 정하기 쉽지 않다.
# 제2 ~ 제5 주성분의 기울기가 비슷해 보인다.
summary(pca)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 2.0300 1.1563 0.8610 0.6492 0.43105 0.38276 0.21693
## Proportion of Variance 0.5887 0.1910 0.1059 0.0602 0.02654 0.02093 0.00672
## Cumulative Proportion 0.5887 0.7797 0.8856 0.9458 0.97235 0.99328 1.00000
# 각 주성분의 중요도 누적합인 Cumulative Proportion을 보면
# 제3 주성분까지 선택할 경우 전체 데이터 88.56% 에 대한 설명력을 갖는다.
# 또한 제4 주성분의 고유값(Standard deviation)이 0.7 보다 작기 때문에 제3 주성분까지 3개의 주성분을 선택한다.
# 3. 각 주성분의 rotation값을 표와 그래프를 사용해 비교하고 주성분의 의미를 해석하시오.
a1 <- pca$rotation[,1]
a1
## YrHgt FtFrBody PrctFFB Frame BkFat SaleHt
## -0.4499313 -0.4123256 -0.3555618 -0.4339569 0.1867048 -0.4528538
## SaleWt
## -0.2699470
bulls_x <- scale(bullsV7) %*% a1
plot(bulls_x, pca$x[, 1], ylab = "PC1", xlab = "scale(bulls) %*% rotation")
cor(bulls_x, pca$x[, 1]) # cor = 1
## [,1]
## [1,] 1
# 7개의 각 변수값에 rotation 값을 행렬곱한 결과가 주성분의 값이다.
# 즉, 주성분의 rotation값은 각 변수에 대해 다음과 같은 설명력을 갖는다.
# PC1 = (-0.45 * YrHgt) + (-0.41 * FtFrBody) + (-0.36 * PrctFFB) +
# (-0.43 * Frame) + (0.19 * BkFat) + (-0.45 * SaleHt) + (-0.27 * SaleWt)
par(mfrow=c(1,2))
barplot(pca$rotation[,1], col = rainbow(8), ylim = c(-0.6,0.4), las = 2, main = "PC1")
abline(h = -0.3, col="blue")
barplot(pca$rotation[,2], col = rainbow(8), ylim = c(-0.4,0.8), las = 2, main = "PC2")
abline(h = 0.3, col="blue")
par(mfrow=c(1,1))
# 4. 행렬도를 사용해 원변수와 주성분의 관계, 원변수 간의 상관관계, 특이한 관측치의 존재 유무 등을 파악하고 설명하시오.
# (1) 공분산 행렬 이용
biplot(pca_cov)
## Warning in arrows(0, 0, y[, 1L] * 0.8, y[, 2L] * 0.8, col = col[2L], length
## = arrow.len): zero-length arrow is of indeterminate angle and so skipped
head(bullsV7)
## YrHgt FtFrBody PrctFFB Frame BkFat SaleHt SaleWt
## 1 51.0 1128 70.9 7 0.25 54.8 1720
## 2 51.9 1108 72.1 7 0.25 55.3 1575
## 3 49.9 1011 71.6 6 0.15 53.1 1410
## 4 53.1 993 68.9 8 0.35 56.4 1595
## 5 51.2 996 68.6 7 0.25 55.0 1488
## 6 49.2 985 71.4 6 0.15 51.4 1500
# 제1, 제2 주성분 만으로도 전체 데이터의 99.96%를 설명할 수 있지만
# 변수들의 단위가 다르기 때문에 결국 큰 값을 가진 FtFrBody와 SaleWt 두 변수가 주성분을 좌우하게 된다.
# 그러므로 변수간 단위가 크게 차이나는 경우 scaling이 적용되는 상관계수 행렬을 사용해야 한다.
# (2) 상관계수 행렬 이용
pca
## Standard deviations:
## [1] 2.0299502 1.1563431 0.8610357 0.6491727 0.4310521 0.3827563 0.2169256
##
## Rotation:
## PC1 PC2 PC3 PC4 PC5
## YrHgt -0.4499313 -0.042790217 -0.41570891 0.1133565 -0.06587066
## FtFrBody -0.4123256 0.129836547 0.45029241 0.2474787 0.71934339
## PrctFFB -0.3555618 -0.315507785 0.56827313 0.3147874 -0.57936738
## Frame -0.4339569 0.007728211 -0.45234503 0.2428179 -0.14299538
## BkFat 0.1867048 0.714719363 -0.03873196 0.6181171 -0.16023789
## SaleHt -0.4528538 0.101315086 -0.17665043 -0.2157694 0.10953536
## SaleWt -0.2699470 0.600514834 0.25331192 -0.5824327 -0.29054729
## PC6 PC7
## YrHgt 0.07223418 -0.77492612
## FtFrBody 0.17706072 -0.01776760
## PrctFFB -0.12780009 0.00239740
## Frame 0.43414400 0.58233705
## BkFat -0.20801720 -0.04244214
## SaleHt -0.79928778 0.23672329
## SaleWt 0.27656055 -0.04703601
biplot(pca)
# BkFat, SaleWt를 제외한 모든 변수가 제1 주성분에 대해 비슷한 정도로 기여를 하고 있다.
# BkFat, SaleWt는 제2 주성분을 특징짓는 변수이다.
# BkFat 변수는 제1 주성분의 특징과 반대되는 성향을 보인다.
# 즉 BkFat 값이 큰 경우 SaleHt, YrHgt, Frame 등의 값은 작은 경향을 보인다.
bulls[c(26, 44), ] # BkFat 값이 큰 소
## Breed SalePr YrHgt FtFrBody PrctFFB Frame BkFat SaleHt SaleWt
## 26 1 1800 47.7 944 66.5 5 0.40 53.3 1556
## 44 5 975 48.6 936 65.3 5 0.35 51.4 1550
bulls[c(63, 57), ] # PrctFFB 값이 큰 소
## Breed SalePr YrHgt FtFrBody PrctFFB Frame BkFat SaleHt SaleWt
## 63 8 1825 53.0 1055 76.8 8 0.10 56.7 1526
## 57 8 1925 52.7 1141 78.5 7 0.15 55.6 1572
bulls[c(26, 44, 63, 57), ]
## Breed SalePr YrHgt FtFrBody PrctFFB Frame BkFat SaleHt SaleWt
## 26 1 1800 47.7 944 66.5 5 0.40 53.3 1556
## 44 5 975 48.6 936 65.3 5 0.35 51.4 1550
## 63 8 1825 53.0 1055 76.8 8 0.10 56.7 1526
## 57 8 1925 52.7 1141 78.5 7 0.15 55.6 1572
# 5. 첫 두개의 주성분을 사용해 산점도를 그리고 Breed를 서로 다른 색깔과 기호로 표시하시오.
# 주성분에 의해 다른 종의 황소를 구분할 수 있는가? 이상점이 있는가? 있다면 어떤 특성을 가진 소인가?
bulls$Breed <- factor(bulls$Breed)
summary(bulls)
## Breed SalePr YrHgt FtFrBody PrctFFB
## 1:32 Min. : 975 Min. :47.20 Min. : 841.0 Min. :64.90
## 5:17 1st Qu.:1375 1st Qu.:49.17 1st Qu.: 932.5 1st Qu.:68.60
## 8:27 Median :1550 Median :50.35 Median : 990.5 Median :70.85
## Mean :1742 Mean :50.52 Mean : 995.9 Mean :70.88
## 3rd Qu.:1900 3rd Qu.:51.73 3rd Qu.:1039.2 3rd Qu.:72.25
## Max. :4600 Max. :54.80 Max. :1383.0 Max. :81.40
## Frame BkFat SaleHt SaleWt
## Min. :5.000 Min. :0.1000 Min. :49.40 Min. :1285
## 1st Qu.:6.000 1st Qu.:0.1500 1st Qu.:52.83 1st Qu.:1474
## Median :6.000 Median :0.1500 Median :54.30 Median :1538
## Mean :6.316 Mean :0.1967 Mean :54.13 Mean :1555
## 3rd Qu.:7.000 3rd Qu.:0.2500 3rd Qu.:55.50 3rd Qu.:1648
## Max. :8.000 Max. :0.5000 Max. :59.60 Max. :1904
par(mfcol = c(1,2))
biplot(pca)
plot(pca$x[,1], pca$x[,2], xlab = "PC1", ylab = "PC2",
pch = as.numeric(bulls$Breed), col = as.numeric(bulls$Breed))
text(pca$x[,1], pca$x[,2], labels = as.character(bulls$Breed),
cex = 0.7, pos = 3, col = as.numeric(bulls$Breed))
par(mfcol = c(1,1))
# 그 plot을 biplot과 같이 보면, 제1 주성분으로 8번종의 소는 구별해낼 수 있지만,
# 1번과 5번 종의 소는 비슷한 특성으로 보이며 섞여 있다.
# 하지만 1번종은 BkFat, SaleWt이 5번 종보다 큰 경향을 보인다.
plabels <- paste(as.character(bulls$Breed), rownames(bulls))
plot(pca$x[,1], pca$x[,2], xlab = "PC1", ylab = "PC2", xlim = c(-8,4), ylim = c(-4,7),
pch = as.numeric(bulls$Breed), col = as.numeric(bulls$Breed))
text(pca$x[,1], pca$x[,2], labels = rownames(bulls),
cex = 0.7, pos = 3, col = as.numeric(bulls$Breed))
legend(-2, 7, c(1,5,8), col = c("black", "red", "green"), fill = c("black", "red", "green"))
lambda <- pca$sdev * sqrt(nrow(pca$x))
Rot <- t(t(pca$rotation)*lambda)
arrows(rep(0,nrow(pca$rotation)), rep(0,nrow(pca$rotation)), Rot[,1], Rot[,2], col = "grey")
text(Rot[,1:2], rownames(Rot), col = "blue")
# 이상점(outlier)
# 16번 소는 BkFat에서 특이하게 큰 값을 보인다.
# 51번 소는 FtFrBody에서 특이하게 큰 값을 보인다.
bulls[16, ]
## Breed SalePr YrHgt FtFrBody PrctFFB Frame BkFat SaleHt SaleWt
## 16 1 2300 49.6 975 68.2 6 0.5 52.9 1842
bulls[51, ]
## Breed SalePr YrHgt FtFrBody PrctFFB Frame BkFat SaleHt SaleWt
## 51 8 1450 53.3 1383 81.4 8 0.2 59.6 1904
par(mfcol=c(1,2))
boxplot(bulls$BkFat, xlab = "BkFat")
boxplot(bulls$FtFrBody, xlab = "FtFrBody")
par(mfcol=c(1,1))
# 6. 첫 주성분을 사용해 Q-Q plot을 그리고 해석하시오.
qqnorm(pca$x[,1])
qqline(pca$x[,1])
shapiro.test(pca$x[,1]) # p-value > 0.05 이므로 정규분포를 따른다고 할 수 있다.
##
## Shapiro-Wilk normality test
##
## data: pca$x[, 1]
## W = 0.96961, p-value = 0.06496