查看原文
其他

基于microeco包进行PCoA分析!!!

王志山 科研后花园 2023-09-08

       microeco 包是基于 R 语言的一个生态学分析包,专门用于分析微生物组数据。该包提供了各种功能用于处理、分析和可视化各种微生物组数据,如 Alpha 多样性、Beta 多样性、NMDS 分析、LEfSe 分析和 co-occurrence 分析等。本期内容是展示如何基于microeco 包进行PCoA分析!


加载所需R包


rm(list=ls())#clear Global Environmentsetwd('D:/桌面/PCoA')#设置工作路径#加载包library(microeco)library(ggplot2)

        microeco包是一个微生物组分析的集成包,依赖包很多,所以安装也比较繁琐,具体安装方式大家可以去该包的官网进行查看,这里小编不作过多赘述!





加载数据


与此前介绍PCoA推文的数据是一致的,需要加载所需的OTU特征表和分组信息表!

#OTU特征表otu <- read.table(file="otu.txt",sep="\t",header=T,check.names=FALSE ,row.names=1)#读入分组文件group <- read.table("group.txt", sep='\t', header=T)





构造microtable


microeco包对于数据处理的方式与phyloseq包是类似的,这里我们需要将加载的OTU特征表与分组信息构造成microtable,具体如下:

data <- microtable$new(sample_table = group, otu_table = otu)#是否使所有文件信息统一data





PCoA分析


##计算beta多样性data$cal_betadiv()#unifrac = F表示不计算unifrac指标##距离计算beta <- trans_beta$new(dataset = data, #数据 group = "group", #分组 measure = "bray")#距离计算方法## PCoA分析beta$cal_ordination(ordination = "PCoA")head(beta$res_ordination)#查看数据# $model# $correction# [1] "none" "1" # # $note# [1] "No correction was applied to the negative eigenvalues"# # $values# Eigenvalues Relative_eig Rel_corr_eig Broken_stick Cum_corr_eig Cumul_br_stick# 1 0.748331060 0.529789796 0.478455549 0.29289683 0.4784555 0.2928968# 2 0.341284734 0.241616551 0.223988859 0.19289683 0.7024444 0.4857937# 3 0.129024261 0.091344247 0.091293348 0.14289683 0.7937378 0.6286905# 4 0.080089708 0.056700453 0.060701710 0.10956349 0.8544395 0.7382540# 5 0.056516568 0.040011571 0.045964865 0.08456349 0.9004043 0.8228175# 6 0.037601664 0.026620541 0.034140134 0.06456349 0.9345445 0.8873810# 7 0.017746969 0.012564176 0.021727890 0.04789683 0.9562724 0.9352778# 8 0.014433485 0.010218356 0.019656452 0.03361111 0.9759288 0.9688889# 9 0.006879324 0.004870298 0.014933937 0.02111111 0.9908627 0.9900000# 10 0.000000000 0.000000000 0.009137255 0.01000000 1.0000000 1.0000000# 11 -0.002393079 -0.001694209 0.000000000 0.00000000 1.0000000 1.0000000# 12 -0.017009084 -0.012041781 0.000000000 0.00000000 1.0000000 1.0000000# # $vectors# Axis.1 Axis.2 Axis.3 Axis.4 Axis.5 Axis.6# A_1 -0.25702536 0.204025208 0.105983433 -0.10079152 0.080021609 -0.029949218# A_4 -0.23313882 0.228981753 -0.007621535 0.12499984 -0.071976959 0.054397737# A_2 -0.22609227 -0.007240737 -0.099493180 0.04868246 0.007904816 0.040525817# A_3 -0.27403139 0.179610380 -0.135946135 -0.04596674 0.018890164 -0.042951015# C_7 0.23577220 0.131910626 0.203035051 -0.01190617 -0.017126534 -0.054687225# B_3 -0.08382270 -0.176991999 0.095687362 -0.10395148 0.010210922 0.119568446# C_5 0.45884845 0.114568653 -0.123729611 -0.11772510 -0.114145063 0.007494586# B_1 -0.09512579 -0.218164092 0.039757723 -0.01256771 -0.046923698 0.023167756# B_4 -0.07889908 -0.222735764 0.006804823 0.01945288 0.007497136 -0.081002220# B_2 -0.08675324 -0.228844682 -0.043921776 0.02759774 -0.051815585 -0.071971759# C_6 0.30043563 0.039100543 0.077922268 0.15278039 0.016631182 0.015454261# C_8 0.33983237 -0.044219890 -0.118478422 0.01939540 0.160832009 0.019952833# Axis.7 Axis.8 Axis.9# A_1 -0.021052300 -0.012930257 -0.004179052# A_4 -0.005481859 0.030047950 0.006552689# A_2 -0.033474280 0.026965546 0.019183055# A_3 0.043128229 -0.039391721 -0.013978077# C_7 -0.018361015 0.041384984 -0.000367667# B_3 0.059849357 0.015382290 -0.001888682# C_5 -0.001954602 -0.005636821 0.011633491# B_1 -0.082924339 -0.058122599 -0.010497898# B_4 0.026633906 -0.003898877 0.058770209# B_2 0.018200342 0.040231685 -0.049193013# C_6 0.041880012 -0.058922771 -0.009404508# C_8 -0.026443450 0.024890592 -0.006630548# # $trace# [1] 1.412506# # attr(,"class")# [1] "pcoa"# # $scores# PCo1 PCo2 PCo3 samples group# A_1 -0.25702536 0.204025208 0.105983433 A_1 A# A_4 -0.23313882 0.228981753 -0.007621535 A_4 A# A_2 -0.22609227 -0.007240737 -0.099493180 A_2 A# A_3 -0.27403139 0.179610380 -0.135946135 A_3 A# C_7 0.23577220 0.131910626 0.203035051 C_7 C# B_3 -0.08382270 -0.176991999 0.095687362 B_3 B# C_5 0.45884845 0.114568653 -0.123729611 C_5 C# B_1 -0.09512579 -0.218164092 0.039757723 B_1 B# B_4 -0.07889908 -0.222735764 0.006804823 B_4 B# B_2 -0.08675324 -0.228844682 -0.043921776 B_2 B# C_6 0.30043563 0.039100543 0.077922268 C_6 C# C_8 0.33983237 -0.044219890 -0.118478422 C_8 C# # $eig# PCo1 PCo2 PCo3 # 53.0 24.2 9.1



可视化


p <- beta$plot_ordination(plot_color = "group", plot_type = c("point", "centroid","ellipse"), point_size = 4, point_alpha = 0.5, ellipse_chull_fill = T, centroid_segment_alpha = 0.9, centroid_segment_size = 1, centroid_segment_linetype = 1)+ theme_bw()+ theme(panel.grid = element_blank(), axis.title.x=element_text(size=12), axis.title.y=element_text(size=12,angle=90), axis.text.y=element_text(size=10), axis.text.x=element_text(size=10), axis.text = element_text(size=8,colour = "black"))+ geom_vline(xintercept = 0, linetype = 2) + geom_hline(yintercept = 0, linetype = 2)p

     可以看到,microeco包的可视化方式是基于ggplot2包的,所以对于图形的修改也就很方便,只需要通过添加或者修改特定参数即可实现既定目的。



组间差异性的检验及添加


这里我们使用Adonis方法进行检验:

#提取作图数据df <- beta$res_ordination$scores

#计算组间差异性Adonis <- adonis2(beta$use_matrix~group,data=df, distance = "bray", # 距离算法 permutations = 999) # 排列次数Adonis

#手动添加显著性信息p+geom_text(aes(x=0.5,y=-0.23),label="p<0.001",size=5,color="red")






参考:https://chiliubio.github.io/microeco_tutorial/

好看你就

点点

您可能也对以下帖子感兴趣

文章有问题?点此查看未经处理的缓存