R：microtable包alpha多样性计算及箱线图绘图

rm(list = ls())
setwd("C:\\Users\\Administrator\\Desktop\\microtable") #设置工作目录
library(microeco)
library(magrittr)
feature_table <- read.table('feature_table_g.txt', header = TRUE, row.names = 1, sep = "\t") #特征表
sample_table <- read.table('sample_table_g.txt', header = TRUE, row.names = 1, sep = "\t") #样本表
tax_table <- read.table('tax_table_g.txt', header = TRUE, row.names = 1, sep = "\t", fill = TRUE) #分类表

dataset <- microtable$new(sample_table = sample_table,
                          otu_table = feature_table, 
                          tax_table = tax_table)
dataset$tidy_dataset() 
dataset$rarefy_samples(sample.size = 1000000)

dataset$cal_alphadiv(PD = FALSE) #计算 alpha 多样性

t1 <- trans_alpha$new(dataset = dataset, group = "Subgroup", by_group = "Type") #设置数据与样本分组
head(t1$data_stat)

#t1$cal_diff(method = "KW")
t1$cal_diff(method = "wilcox")

#head(t1$res_diff)
write.table(t1$res_diff, file = "alpha_diversity_differences.txt", sep = "\t", quote = FALSE, row.names = FALSE)

###############################################################
library(ggplot2)
library(ggpubr)
library(dplyr)
library(multcompView)
index <- read.table('Simpson_index.txt', header = TRUE, row.names = 1)

# 绘制箱线图
p <- ggboxplot(index, x = "DAS", y = "Simpson", color = "Gene", width = 0.75,
               size = 1.5,  # 直接在ggboxplot中设置箱线图的线条粗细
               add = "jitter",
               add.params = list(size = 5, alpha = 0.7)) +
  theme_minimal(base_size = 16) +
  labs(x = "", y = "Simpson") +  # 设置y轴标题
  theme(legend.position = "top",
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        axis.line = element_line(colour = "black", linewidth = 1.5),
        axis.text.x = element_text(size = 22, face = "plain", color = "black", angle = 0, hjust = 0.5, family = "TNR"),
        axis.text.y = element_text(size = 22, face = "plain", color = "black", family = "TNR"),
        axis.title.y = element_text(size = 36, face = "plain", color = "black", family = "TNR"),  # 自定义Y轴标题文字大小
        legend.title = element_blank(),  # 移除图例标题
        axis.ticks = element_line(color = "black"),
        axis.ticks.length = unit(0.25, "cm"),
        legend.key.size = unit(1.2, "cm"),
        legend.text = element_text(size = 18, family = "TNR")) +
  scale_color_manual(values = c("B73" = "#8FC9E2", "Mo17" = "#ECC97F"))

# 手动添加显著性标记
dodge_width <- 0.375
p <- p + annotate("segment", x = 1 - dodge_width, xend = 1 + dodge_width, y = 0.99, yend = 0.99, size = 1.2) +
  annotate("text", x = 1, y = 0.993, label = "*", size = 10, vjust = 0) +
  annotate("segment", x = 2 - dodge_width, xend = 2 + dodge_width, y = 0.99, yend = 0.99, size = 1.2) +
  annotate("text", x = 2, y = 0.993, label = "**", size = 10, vjust = 0) +
  annotate("segment", x = 3 - dodge_width, xend = 3 + dodge_width, y = 0.99, yend = 0.99, size = 1.2) +
  annotate("text", x = 3, y = 0.993, label = "*", size = 10, vjust = 0) +
  annotate("segment", x = 4 - dodge_width, xend = 4 + dodge_width, y = 0.99, yend = 0.99, size = 1.2) +
  annotate("text", x = 4, y = 0.993, label = "*", size = 10, vjust = 0)
# 显示图表
print(p)
# 保存图形
ggsave("boxplot.png", plot = p, width = 10, height = 8, dpi = 1200, bg = "white")