點擊上方「藍字」關注我

寫在前面

• 前面我們知道了降維分析

• 學習了PCoA分析

• 今天學習PCA分析...

PCA(Principal Component Analysis)，即主成分分析方法，是一種使用最廣泛的數據降維算法。PCA的主要思想是將n維特征映射到k維上，這k維是全新的正交特征也被稱為主成分，是在原有n維特征的基礎上重新構造出來的k維特征。PCA的工作就是從原始的空間中順序地找一組相互正交的坐標軸，新的坐標軸的選擇與數據本身是密切相關的。其中，第一個新坐標軸選擇是原始數據中方差最大的方向，第二個新坐標軸選取是與第一個坐標軸正交的平面中使得方差最大的，第三個軸是與第1,2個軸正交的平面中方差最大的。依次類推，可以得到n個這樣的坐標軸。通過這種方式獲得的新的坐標軸，我們發現，大部分方差都包含在前面k個坐標軸中，后面的坐標軸所含的方差幾乎為0。于是，我們可以忽略余下的坐標軸，只保留前面k個含有絕大部分方差的坐標軸。事實上，這相當于只保留包含絕大部分方差的維度特征，而忽略包含方差幾乎為0的特征維度，實現對數據特征的降維處理。

思考：我們如何得到這些包含最大差異性的主成分方向呢？

答案：事實上，通過計算數據矩陣的協方差矩陣，然后得到協方差矩陣的特征值特征向量，選擇特征值最大(即方差最大)的k個特征所對應的特征向量組成的矩陣。這樣就可以將數據矩陣轉換到新的空間當中，實現數據特征的降維。

?PCoA(Principal Co-ordinates Analysis)分析即主坐標分析。它與PCA類似，通過一系列的特征值和特征向量進行排序后，選擇主要排在前幾位的特征值，找到距離矩陣中最主要的坐標，結果是數據矩陣的一個旋轉，它沒有改變樣本點之間的相互位置關系，只是改變了坐標系統。如何取舍？

在微生物分析中我們會基于beta多樣性分析得到的距離矩陣，進行PCA和PCoA分析，具體距離矩陣可見：beta多樣性矩陣。PCA是基于樣本的相似矩陣(如歐式距離)來尋找主成分，而PCoA是基于相異距離矩陣(歐式距離以外的其他距離，包括binary_jaccard ,bray_curtis ,unweighted_unifrac和weighted_unifrac距離)來尋找主坐標。

在分析的過程中PCA和PCoA分析都會用到降維的思想，但是在降維的過程中必然會造成數據損失，多數情況下，我們在做降維處理的時候，期望維數越低越好，這樣我們就可以最大程度地保真原始數據。PCA基于物種豐度矩陣就意味著PCA分析的矩陣維度等于物種數目。同樣的道理，PCoA基于樣本間的距離矩陣就意味著PCoA分析的矩陣維度與樣本數目相關。因此，如果樣本數目比較多，而物種數目比較少，那肯定首選PCA；如果樣本數目比較少，而物種數目比較多，那肯定首選PCoA。

上期我們已經介紹了PCoA，今天我們主講PCA。

發表在2019年代的Microbome (IF>10.0)上的一篇文章Comparative metagenomic and metatranscriptomic analyses reveal the breed effect on the rumen microbiome and its associations with feed efficiency in beef cattle，利用宏基因組和宏轉錄組的關聯分析揭示不同牛品種對瘤胃微生物菌群結構及飼料利用率的影響。作者在分析品種對瘤胃細菌(屬水平)和古菌(種)的影響時候，做了一個PCoA分析，結果如下圖：

但是如果我們不是分析微生物多樣性數據，而是分析其他環境變量，能不能同樣做出此圖呢？當然可以，但是我們用到的吧不是PCoA，而是PCA分析了。話不多說，我們開始實戰。

逐準備數據

在Excel里面輸入我們的數據，行為組，檢測變量為列，切記一定要遵循此規則，不然后面會出現bug。

加載包計算PC1

、PC2以及PC3

library(openxlsx)mydatalibrary(vegan)#PCA主要功能包

把mydata轉化為矩陣，注意第一列為字符，我們要排除第一列

A

Acolnames(A) rownames(A) rownames(A)library(vegan)env##PCA 排序#環境變量需要標準化，詳情 ?rdapca_env #I 型標尺summary(pca_env, scaling = 1)#II 型標尺#summary(pca_env, scaling = 2)#各主成分(PCA軸)特征根，特征根是固定的，和標尺選擇無關pca_eig #除以特征根總和，即可得各主成分(PCA軸)解釋量pca_exp 以上就是我們后面作圖要用到的PC1、PC2以及PC3

---

計算坐標尺度

#提取對象排序坐標，以 I 型標尺為例，以前兩軸為例#scores() 提取樣方坐標site.scaling1 #或者在 summary 中提取site.scaling1 site.scaling1#write.table(site.scaling1, 'site.scaling1.txt', col.names = NA, sep = '\t', quote = FALSE)

提取變量排序坐標，以 II 型標尺為例，以前兩軸為例

#scores() 提取環境變量坐標

env.scaling2

#或者 summary 中提取

env.scaling2

env.scaling2

#write.table(env.scaling2, 'env.scaling2.txt', col.names = NA, sep = '\t', quote = FALSE)

以上就是計算標尺坐標

---

ggplot2作圖

我們用Excel打開輸出的pca_site.txt文件。加入我們的組別信息：我們共有三組(group1所示)，個數如group2所示：

我們讀取pca-sitepca_site

pca_site

library(scatterplot3d)#加載包

scatterplot3d(pca_site$PC1, pca_site$PC2, pca_site$PC3, pch=16,color = rep(c('red2', 'purple2','blue2'), c(3,6,5)),xlab = paste('PCA1: 23.14%'), ylab = paste('PCA2: 13.21%'), zlab = paste('PCA3: 10.02%'))

c(3,6,5)為樣本數量，其中樣本數量必須與顏色一致

是不是很簡單易操作，還在猶豫什么趕緊動動手操作吧！

? ? ? ? ? ? 快開學啦

點贊關注我

總結

以上是生活随笔為你收集整理的pca降维的基本思想_R语言进行PCA分析的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。