使用R语言进行perMANVOA分析NMDS结果多个分组之间显著性

最新推荐文章于 2023-02-23 19:57:53 发布

Nan_White_

最新推荐文章于 2023-02-23 19:57:53 发布

阅读量1.5k

点赞数

分类专栏： R 文章标签： python 开发语言

本文链接：https://blog.csdn.net/xin996ya/article/details/127309265

版权

R 专栏收录该内容

2 篇文章 1 订阅

订阅专栏

library(vegan)
otu_table.csv <- read.csv("C:/Users/baixin/Desktop/data_nmds.csv",header = T,row.names = 1,check.names=FALSE)
#check.names=FALSE，防止把数值列名转换成添加字母的形式,row.name=1,规定首列为行名，header = T首行为列名
#View(otu_table.csv)#Seeing the data by the type of table.以表格形式查看
#otu_table.csv#输出台查看
otu_transition <- data.frame(t(otu_table.csv))#将OTU 丰度表转换为需要的格式，otu正常为首行表示样本名称，首列表示多个样本检测指标名称
#View(otu_transition)
#otu_transition
#样本分组文件
group <- read.csv("C:/Users/baixin/Desktop/groups.csv")#读入分组文件
#View(group)

##PERMANOVA 分析（所有分组间比较，即整体差异）
#（2）若是使用 OTU 丰度表，则需要在计算时指定所依据的距离类型，这里依然使用 Bray-Curtis 距离
adonis_result_otu_transition <- adonis2(otu_transition ~ site, group, permutations = 999, distance = 'bray',by="margin")   
#这条命令的详情和上述命令所表示的信息一致，表中必须全为数值，不能有一行或列全为0，by="margin"表示变量排列顺序对判断无关
#adonis_result_otu_transition



#可选输出
otuput <- data.frame(adonis_result_otu_transition, check.names = FALSE, stringsAsFactors = F)
if (otuput[1,'Pr(>F)'] <= 0.001) otuput[1, 'SIG'] <- '***' else if (otuput[1,'Pr(>F)'] <= 0.01) otuput[1, 'SIG'] <- '**' else if (otuput[1,'Pr(>F)'] <= 0.05) otuput[1, 'SIG'] <- '*'
otuput

write.csv(otuput,"C:/Users/baixin/Desktop/otuput.csv")



##PERMANOVA 分析（使用循环处理，进行小分组间比较，如两组间）
#推荐使用 OTU 丰度表作为输入数据，每次筛选分组后重新计算样本距离，避免由于样本数减少可能导致的距离变动而造成误差
group_name <- unique(group$site)#判定一共分为几组，得到组名附值给group_name
#group_name#输出“b”，“a”，“c”
adonis_result_two <- NULL
for (i in 1:(length(group_name) - 1)) {#从范围中选择i和j
  for (j in (i + 1):length(group_name)) {
    group_ij <- subset(group, site %in% c(group_name[i], group_name[j]))#选择两个分组
    otu_ij <- otu_transition[group_ij$yangdian, ]#group_ij$yangdian将行名符合group_ij的数据提取出来赋值给otu_ij
    adonis_result_otu_ij <- adonis2(otu_ij~site, group_ij, permutations = 999, distance = 'bray',by="margin")#随机置换检验 999 次
    adonis_result_two <- rbind(adonis_result_two,
                               c(paste(group_name[i], group_name[j], sep = '/'), 
                                 '以Bray-Curtis比较', unlist(data.frame(adonis_result_otu_ij, check.names = FALSE)[1, ])))
  }
}
#View(otu_ij)
#可选标记 “*” 显著性，nrow(adonis_result_two)表示adonis_result_two数量
#'Pr(>F)'需要在输出控制台复制，不同版本可能表达不同,
adonis_result_two <- data.frame(adonis_result_two,check.names = FALSE, stringsAsFactors = F)
for (i in 1:nrow(adonis_result_two)) {
  if (adonis_result_two[i, 'Pr(>F)'] <= 0.001) adonis_result_two[i, 'SIG'] <- '***'
  else if (adonis_result_two[i, 'Pr(>F)'] <= 0.01) adonis_result_two[i, 'SIG'] <- '**'
  else if (adonis_result_two[i, 'Pr(>F)'] <= 0.05) adonis_result_two[i, 'SIG'] <- '*'
}#此段代码只能放在上一行转置代码后
adonis_result_two
write.csv(adonis_result_two ,"C:/Users/baixin/Desktop/output_two.csv")