目錄
主成分方差解釋率計算
通常，求得了PCA降維后的特征值，我們就可以繪圖，但各個維度的方差解釋率沒有得到，就無法獲得PC坐標的百分比。
有些工具的結果是提供了維度標准差的，如ggbiplot繪圖時，直接會給你算出各個坐標的方差解釋率。但我覺得這類工具繪圖遠不如ggplot本身，此時，就需要自己計算。
當理解了PCA的原理和含義后，就比較容易得到。網上一大堆，這里不介紹。
以ggbiplot數據為例，並將算出結果與之比較。
if(!require(devtools))
install.packages("devtools")
library(devtools)
if(!require(ggbiplot))
install_github("vqv/ggbiplot")
library(ggbiplot)
data(wine)
pca <- prcomp(wine, scale. = TRUE)
ggbiplot(pca,
# groups = wine.class,
ellipse = TRUE, circle = TRUE,
obs.scale = 1, var.scale = 1) +
scale_color_discrete(name = '') +
theme(legend.direction = 'horizontal', legend.position = 'top')
R自帶函數prcomp的結果中得到PCA的5個對象結果，其中包含了標准差（sdev）和特征向量（x）。
> str(pca)
List of 5
$ sdev : num [1:13] 2.169 1.58 1.203 0.959 0.924 ...
$ rotation: num [1:13, 1:13] -0.14433 0.24519 0.00205 0.23932 -0.14199 ...
..- attr(*, "dimnames")=List of 2
.. ..$ : chr [1:13] "Alcohol" "MalicAcid" "Ash" "AlcAsh" ...
.. ..$ : chr [1:13] "PC1" "PC2" "PC3" "PC4" ...
$ center : Named num [1:13] 13 2.34 2.37 19.49 99.74 ...
..- attr(*, "names")= chr [1:13] "Alcohol" "MalicAcid" "Ash" "AlcAsh" ...
$ scale : Named num [1:13] 0.812 1.117 0.274 3.34 14.282 ...
..- attr(*, "names")= chr [1:13] "Alcohol" "MalicAcid" "Ash" "AlcAsh" ...
$ x : num [1:178, 1:13] -3.31 -2.2 -2.51 -3.75 -1.01 ...
..- attr(*, "dimnames")=List of 2
.. ..$ : NULL
.. ..$ : chr [1:13] "PC1" "PC2" "PC3" "PC4" ...
- attr(*, "class")= chr "prcomp"
手動計算方差解釋率：
> pca$sdev^2/sum(pca$sdev^2)*100
#注意平方
[1] 36.1988481 19.2074903 11.1236305 7.0690302 6.5632937 4.9358233
[7] 4.2386793 2.6807489 2.2221534 1.9300191 1.7368357 1.2982326
[13] 0.7952149
可看出，前兩個主成分與圖中一致。當然如果沒有標准差結果，我們也可以根據特征向量計算出來：
> sdev<- apply(pca$x,2,sd)
> sdev
PC1 PC2 PC3 PC4 PC5 PC6 PC7
2.1692972 1.5801816 1.2025273 0.9586313 0.9237035 0.8010350 0.7423128
PC8 PC9 PC10 PC11 PC12 PC13
0.5903367 0.5374755 0.5009017 0.4751722 0.4108165 0.3215244
繪圖示例
一個示例，可在此基礎上進一步優化。如樣本要再分組，可加shape。
ggplot(data=data.frame(pca$x), aes(PC1,PC2)) +
stat_ellipse(aes(fill=wine.class),type="norm",geom="polygon",alpha=0.2,color=NA)+
geom_point(size=2)+
# scale_size(guide=FALSE)+
scale_color_manual(values = col)+
geom_vline(xintercept = 0,linetype="dotted")+
geom_hline(yintercept = 0,linetype="dotted")+
labs(x=paste0("PC1", sprintf("(%0.2f%%)",100*pca$sdev[1]^2/sum(pca$sdev^2))),
y=paste0("PC2", sprintf("(%0.2f%%)",100*pca$sdev[2]^2/sum(pca$sdev^2))))
