1. 程式人生 > >PCA、PCoA、NMDS、Anosim R繪圖

PCA、PCoA、NMDS、Anosim R繪圖

基於R包 

PCA R繪圖 vegan

### 輸入檔案的格式處理
>head(dataT)
            sp1    sp2    sp3    sp4    ...
sample1    相對丰度
sample2
sample3
  ...

>dataTD=decostand(data,"hell")  ## 資料的標準化“採用total標準化以後再取平方根”
>pca1=rda(dataTD)          ## RDA冗餘分析
>pc1=c(pca1$CA$eig/sum(pca1$CA$eig))[1]*100  ## 計算第一主成分
>pc2=c(pca1$CA$eig/sum(pca1$CA$eig))[2]*100  ## 計算第二主成分
##繪製散點圖
>plot(pca1,display="si",scaling=1,type="n", main="")
>points(pca1, dis="si", scaling=1,col=c("#C1E168", "#C1E168", "#C1E168", "#FD9347" ),pch=c(15, 15, 15, 17 ),cex=1)
##legend("bottomright", legend=levels(groups$g1), col=mycol, =shape,bty="n",cex=0.8) ##新增圖例
> ordispider(pca1,groups = groups$group,display = "si", scaling=1,col = mycol)
#新增輔助線,ordispider把專案組合至它們的(加權)類質心,ordiellpse畫出類標準差、標準誤差或者置信區域的橢圓。

PS:
head(groups)
samples group
sample1    A    
sample2    A
sample3    B
 ...

PCoA

dataTB=vegdist(data, method="bray") ### 計算矩陣距離
dataTBpcoA=cmdscale(dataTB)

NMDS

nmds1=metaMDS(dataTD,distance = "bray", k = 2,trymax=20)

CCA

cca1=cca(dataTD)

 Anosim

  ar = anosim(dataT,groups$group,permutations = 999, distance = "bray")
  summary(ar)
  dc = data.frame(dis=ar$dis.rank, class=ar$class.vec)