大家對主成分分析(principal components analysis, PCA) 都很熟悉，但是今天我們來介紹下主坐标分析(principal coordinate analysis, PCoA)。那麼這兩個差了個o字母具體有什麼差別？首先PCA是常用的降維算法；利用線性變換，将資料變換到一個新的坐标系統中;然後再利用降維的思想，使得任何資料投影的第一大方差在第一個坐标(稱為第一主成分)上，第二大方差在第二個坐标(第二主成分)上。這種降維的思想首先減少資料集的維數，同時還保持資料集的對方差貢獻最大的特征，最終使資料直覺呈現在二維坐标系。PCoA主要是探索資料相似度或者相異度可視化方法。可呈現研究資料相似性或差異性的可視化坐标，是一種非限制性的資料降維分析方法，可用來研究樣本群落組成的相似性或相異性。其實通俗的講，PCA主要是基于原始資料矩陣的降維；PCoA主要是基于樣本的原始資料計算出來的距離矩陣的降維。如果樣本數目比較多，而物種數目比較少，那肯定首選PCA;如果樣本數目比較少，而物種數目比較多，那肯定首選PCoA。
接下來我們看下在R中如何去實作，首先安裝ape包和vegan包，聯合使用才能達到最終的目的。包的安裝我們就不贅述了，其在CRAN平台，直接install.packages()。
首先是資料的導入，我們利用vegan自帶的資料dune。具體的資料集的構成大家可以直接在包的資訊中去看。接下來我們首先基于dune資料構造距離矩陣，需要用到的函數vegdist。
vegdist(x, method="bray",binary=FALSE, diag=FALSE, upper=FALSE,
na.rm = FALSE, ...)
其中主要的參數：
Method 支援了目前大部分的距離計算函數。"manhattan", "euclidean", "canberra","clark", "bray", "kulczynski","jaccard", "gower", "altGower","morisita", "horn", "mountford", "raup","binomial", "chao", "cao" or"mahalanobis"。
Diag 是否顯示對角線的值。
Upper 是否顯示對角線以上的值
library(vegan)
data(dune)
data(dune.env)
dune.dist #構造距離矩陣。
接下來就是利用ape中的pcoa函數擷取PCOA分析結果。當然也可以應用我們R自帶的函數cmdscale。
pcoa(D, correction="none",rn=NULL)
其中主要參數：
D 不用多說就是距離矩陣
Correction 主要指需不需進行校正。
Rn 主要是指的每個樣本的名稱标簽，可以進行自定義
接下來我們看下執行個體：
res=pcoa(dune.dist,correction="none")
head(res$value)
其中主要的值是特征值的一些相關的轉換的值。
head(res$vectors)
其中主要是和PCA中主成分類似的柱坐标的值，進行了排序展示，一般選擇前兩個繪制二維可視散點圖。
biplot(res)#可視化PCOA 的結果
為了進一步完善我們的可靠性，我們還可以利用vegan中的ANOSIM相似性分析是一種非參數檢驗，用來檢驗組間(兩組或多組)差異是否顯著大于組内差異，進而判斷分組是否有意義。
anosim(x, grouping, permutations = 999, distance = "bray", strata = NULL, parallel = getOption("mc.cores"))
其中的參數我們不再贅述，直接看下執行個體：
dune.ano <- with(dune.env, anosim(dune.dist, Management))
summary(dune.ano)
plot(dune.ano)
至此，我們的PCOA的分析過程可以實作，那麼如何優化我們輸出的可視化圖像，我們需要用到ggplot2這個包可以對我們的值進行更加友好的可視化。直接看執行個體：
library(ggplot2)
P= dune.ano$signif
R=round( dune.ano$statistic,3)
data=res$vectors[,1:2]
data=data.frame(data)
colnames(data)=c('x','y')
eig=as.numeric(res$value[,1])
type=dune.env$Management
data=cbind(data,type)
p =ggplot(data, aes(x=x, y=y, color=type))+
geom_point(alpha=.7, size=2) +
labs(x=paste("PCoA 1 (",format(100* eig[1] / sum(eig), digits=4), "%)", sep=""),
y=paste("PCoA 2 (", format(100*eig[2] / sum(eig), digits=4), "%)", sep=""),
title="PCoA")
p +theme_classic()+stat_ellipse(type ="t", linetype = 2)+
annotate("text",x=0.004,y=-0.15,parse=TRUE,size=4,label=paste('R:',R),family="serif",font,colour="darkred")+
annotate("text",x=0,y=-0.17,parse=TRUE,size=4,label=paste('p:',P),family="serif",font,colour="darkred")
歡迎大家互相學習！