快捷導(dǎo)航

R語(yǔ)言實(shí)現(xiàn)PCA主成分分析圖的示例代碼

更新時(shí)間：2022年04月29日 08:32:34 作者：watermark's

主成分分析（Principal?Component?Analysis，PCA）是一種無(wú)監(jiān)督的數(shù)據(jù)降維方法，通過主成分分析可以盡可能保留下具備區(qū)分性的低維數(shù)據(jù)特征。本文將用R語(yǔ)言實(shí)現(xiàn)PCA主成分分析圖，需要的可以參考一下

簡(jiǎn)介

主成分分析（Principal Component Analysis，PCA）是一種無(wú)監(jiān)督的數(shù)據(jù)降維方法，通過主成分分析可以盡可能保留下具備區(qū)分性的低維數(shù)據(jù)特征。主成分分析圖能幫助我們直觀地感受樣本在降維后空間中的分簇和聚合情況，這在一定程度上亦能體現(xiàn)樣本在原始空間中的分布情況，這對(duì)于只能感知三維空間的人類來(lái)說(shuō)，不失為一種不錯(cuò)的選擇。

再舉個(gè)形象的栗子，假如你是一本養(yǎng)花工具宣傳冊(cè)的攝影師，你正在拍攝一個(gè)水壺。水壺是三維的，但是照片是二維的，為了更全面的把水壺展示給客戶，你需要從不同角度拍幾張圖片。下圖是你從四個(gè)方向拍的照片：

PCA形象解釋說(shuō)明

第一張圖里水壺的背面可以看到，但是看不到前面。

第二張圖是拍前面，可以看到壺嘴，這張圖可以提供了第一張圖缺失的信息，但是壺把看不到了。

第三張俯視圖既可以看到壺嘴，也可以看到壺把，但是無(wú)法看出壺的高度。

第四張圖是你打算放進(jìn)目錄的，水壺的高度，頂部，壺嘴和壺把都清晰可見。

PCA的設(shè)計(jì)理念與此類似，它可以將高維數(shù)據(jù)集映射到低維空間的同時(shí)，盡可能的保留更多變量。

開始作圖

使用 R 語(yǔ)言能做出像 SIMCA-P 一樣的 PCA 圖嗎？

答案是肯定的，使用 R 語(yǔ)言不僅能做出像 SIMCA-P 一樣的 PCA 圖，還能做出比 SIMCA-P 更好看的圖，而且好看的上限僅取決于個(gè)人審美風(fēng)格。

1. PCA 分析圖本質(zhì)上是散點(diǎn)圖

主成分分析圖 = 散點(diǎn)圖 + 置信橢圓，散點(diǎn)的橫縱坐標(biāo)對(duì)應(yīng) PCA 的第一主成分、第二主成分。

library(ggplot2)

# 數(shù)據(jù)準(zhǔn)備
data = subset(iris, select = -Species)
class = iris[["Species"]]

# PCA
pca = prcomp(data, center = T, scale. = T)
pca.data = data.frame(pca$x)
pca.variance = pca$sdev^2 / sum(pca$sdev^2)

ggplot(pca.data, aes(x = PC1, y = PC2)) +
    geom_point(size = 3) +
    geom_hline(yintercept = 0) + 
    geom_vline(xintercept = 0) +
    stat_ellipse(aes(x = PC1, y = PC2), linetype = 2, size = 0.5, level = 0.95) + 
    theme_bw()

散點(diǎn)圖+橢圓

2. 為不同類別著色

接下來(lái)想給散點(diǎn)加上分類顏色：

library(ggplot2)

# 數(shù)據(jù)準(zhǔn)備
data = subset(iris, select = -Species)
class = iris[["Species"]]

# PCA
pca = prcomp(data, center = T, scale. = T)
pca.data = data.frame(pca$x)
pca.variance = pca$sdev^2 / sum(pca$sdev^2)

ggplot(pca.data, aes(x = PC1, y = PC2, color = class)) +
    geom_point(size = 3) +
    geom_hline(yintercept = 0) + 
    geom_vline(xintercept = 0) +
    stat_ellipse(aes(x = PC1, y = PC2), linetype = 2, size = 0.5, level = 0.95) + 
    theme_bw()

給不同類別樣本著色

顏色是加上了，但是橢圓咋變成了 3 個(gè)？

原來(lái)是 stat_ellipse 函數(shù)默認(rèn)對(duì)每個(gè)類別的數(shù)據(jù)計(jì)算自己的置信區(qū)間。如何對(duì)多類樣本只計(jì)算一個(gè)置信區(qū)間呢？查看 stat_ellipse 的幫助文檔：

inherit.aesdefault TRUE, If FALSE, overrides the default aesthetics, rather than combining with them. This is most useful for helper functions that define both data and aesthetics and shouldn't inherit behaviour from the default plot specification,

原來(lái)是 stat_ellipse 函數(shù)默認(rèn)會(huì)繼承 ggplot 中的 aes 設(shè)置，如果希望 stat_ellipse 使用自己的 aes 設(shè)置，需要將參數(shù) inherit.aes 設(shè)置為 FALSE。

library(ggplot2)

# 數(shù)據(jù)準(zhǔn)備
data = subset(iris, select = -Species)
class = iris[["Species"]]

# PCA
pca = prcomp(data, center = T, scale. = T)
pca.data = data.frame(pca$x)
pca.variance = pca$sdev^2 / sum(pca$sdev^2)

ggplot(pca.data, aes(x = PC1, y = PC2, color = class)) +
    geom_point(size = 3) +
    geom_hline(yintercept = 0) + 
    geom_vline(xintercept = 0) +
    stat_ellipse(aes(x = PC1, y = PC2), linetype = 2, size = 0.5, level = 0.95, inherit.aes = FALSE) + 
    theme_bw()

給不同類別樣本著色

3. 樣式微調(diào)

接下來(lái)對(duì)樣式進(jìn)行微調(diào)：為不同類別樣本自定義著色，添加 x 軸、y 軸標(biāo)題，添加 title：

library(ggplot2)

# 數(shù)據(jù)準(zhǔn)備
data = subset(iris, select = -Species)
class = iris[["Species"]]

# PCA
pca = prcomp(data, center = T, scale. = T)
pca.data = data.frame(pca$x)
pca.variance = pca$sdev^2 / sum(pca$sdev^2)

# 自定義顏色
palette = c("mediumseagreen", "darkorange", "royalblue")

ggplot(pca.data, aes(x = PC1, y = PC2, color = class)) +
    geom_point(size = 3) +
    geom_hline(yintercept = 0) + 
    geom_vline(xintercept = 0) +
    stat_ellipse(aes(x = PC1, y = PC2), linetype = 2, size = 0.5, level = 0.95, inherit.aes = FALSE) + 
    theme_bw() +
    scale_color_manual(values = palette) +
    theme(panel.grid.major.x = element_blank(),
          panel.grid.minor.x = element_blank(),
          panel.grid.major.y = element_blank(),
          panel.grid.minor.y = element_blank()) +
    labs(x = paste0("PC1: ", signif(pca.variance[1] * 100, 3), "%"), 
         y = paste0("PC2: ", signif(pca.variance[2] * 100, 3), "%"), 
         title = paste0("PCA of iris")) +
    theme(plot.title = element_text(hjust = 0.5))