其他
基于microeco包进行PCoA分析!!!
microeco 包是基于 R 语言的一个生态学分析包,专门用于分析微生物组数据。该包提供了各种功能用于处理、分析和可视化各种微生物组数据,如 Alpha 多样性、Beta 多样性、NMDS 分析、LEfSe 分析和 co-occurrence 分析等。本期内容是展示如何基于microeco 包进行PCoA分析!
rm(list=ls())#clear Global Environment
setwd('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)
microeco包对于数据处理的方式与phyloseq包是类似的,这里我们需要将加载的OTU特征表与分组信息构造成microtable,具体如下:
data <- microtable$new(sample_table = group,
otu_table = otu)#是否使所有文件信息统一
data
##计算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
这里我们使用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")
好看你就
点点
我