微生物组β-多样性——PCoA和可视化

在微生物组测序结果中对于β-多样性的分析一般以PCoA 和 NMDS 为主，并且以 Bray-Curtis 距离最为普遍，因此本文以该距离进行PCoA分析和可视化

（一）数据格式

分别为分析样本的特征表和分组信息

表1. pcoa （特征表）

表2，group（分组信息）

（）9

1.  
   #清除变量
2.  
   rm(list=ls())
3.  
   library(vegan)
4.  
   library(ggplot2)
5.  
   library(ggrepel)
6.  
   #读取数据(特征表)
7.  
   df<- read.delim('pcoa.txt', row.names = 1, sep = '\t', head = TRUE, check.names = FALSE)
8.  
   #读取分组
9.  
   group<-read.delim('group.txt', row.names = 1, sep = '\t', head = TRUE, check.names = FALSE)
10.  
    #数据转置
11.  
    df1<-t(df)
12.  
    #计算距离
13.  
    distance<-vegdist(df1,method='bray')
14.  
    pcoa<- cmdscale(distance,k=(nrow(df1)-1),eig=TRUE)
15.  
    #提取前两个分类解释
16.  
    plot_data<-data.frame({pcoa$point})[1:2]
17.  
    head(plot_data)
18.  
    #前两个分类解释命名
19.  
    names(plot_data)[1:2]<-c('PCoA1','PCoA2')
20.  
    eig=pcoa$eig
21.  
    group1<-group['group']
22.  
    data<-plot_data[match(rownames(group),rownames(plot_data)),]
23.  
    data<-data.frame(group,plot_data)
24.  
    head(data)
25.  
    tail(data)
26.  
    #作图
27.  
    ggplot(data,aes(x=PCoA1,y=PCoA2,shape=group,color=group)) geom_point(alpha=1,size=8) stat_ellipse(level=0.95,size=3) labs(x=paste("PCoA1(",format(100*eig[1]/sum(eig),digits = 4),"%)",sep=""),y=paste("PCoA2(",format(100*eig[2]/sum(eig),digits = 4),"%)",sep="")) geom_vline(aes(xintercept=0),linetype="dotted") geom_hline(aes(yintercept=0),linetype="dotted") theme(panel.background = element_rect(fill='white',colour = 'black'),axis.title.x=element_text(colour = 'black',size=20),axis.title.y = element_text(colour = 'black',size=20),legend.text = element_text(size = 15))

（三） 出图

通过以上代码对微生物组β-多样性以PCoA进行计算并可视化得到如下图

这篇好文章是转载于：学新通技术网