擴增子統計繪圖1箱線圖：Alpha多樣性
阿新 • • 發佈：2017-08-23
ggplot2 highlight hole ima pre 運行 提取 posit row 繪制Alpha多樣性線箱圖
繪圖和統計全部為R語言，建議復制代碼，在Rstuido中運行，並設置工作目錄為存儲之前分析結果文件的result目錄
# 運行前，請在Rstudio中菜單欄選擇“Session - Set work directory -- Choose directory”，彈窗選擇之前分析目錄中的result文件夾
# 安裝相關軟件包，如果末安裝改為TRUE運行即可安裝
if (FALSE){
source("https://bioconductor.org/biocLite.R")
biocLite(c("ggplot2"))
}
# 加載相關軟件包
library("ggplot2") # load related packages
# 讀入實驗設計和Alpha多樣性值
design = read.table("design.txt", header=T, row.names= 1, sep="\t")
alpha = read.table("alpha.txt", header=T, row.names= 1, sep="\t")
# 以Observed OTU為例進行可視化和統計分析，其它指數將observed_otus替換為shannon, chao1, PD_whole_tree即可計算
# 合並Alpha指數與實驗設計
index = cbind(alpha, design[match(rownames(alpha), rownames(design)), ])
# 繪圖代碼、預覽、保存PDF
p = ggplot(index, aes(x=genotype, y=observed_otus, color=genotype))+
geom_boxplot(alpha=1, outlier.size=0, size=0.7, width=0.5, fill="transparent") +
geom_jitter( position=position_jitter(0.17), size=1, alpha=0.7)+
labs(x="Groups", y="observed_otus index")
p
ggsave(paste("alpha_observed_otus.pdf", sep=""), p, width = 5, height = 3)
# 統計組間是否顯著差異
# anova對指數與分組統計
observed_otus_stats <- aov(observed_otus ~ genotype, data = index)
# 使用TukeyHSD對組間進行檢驗，效正pvalue
Tukey_HSD_observed_otus <- TukeyHSD(observed_otus_stats, ordered = FALSE, conf.level = 0.95)
# 結果中提取需要的結果
Tukey_HSD_observed_otus_table <- as.data.frame(Tukey_HSD_observed_otus$genotype)
# 預覽結果
Tukey_HSD_observed_otus_table
# 保存結果到文件，按Pvaule值由小到大排序
write.table(Tukey_HSD_observed_otus_table[order(Tukey_HSD_observed_otus_table$p, decreasing=FALSE), ], file="alpha_observed_otus_stats.txt",append = FALSE, quote = FALSE, sep="\t",eol = "\n", na = "NA", dec = ".", row.names = TRUE,col.names = TRUE)
各組間的統計結果如下：主要看最後一列p adj(Adjust P-value)是否顯著，本文數據不顯著
diff lwr upr p adj OE-KO -7.52380952380952 -24.480725165752 9.43310611813294 0.515429907536906 WT-KO -6.11111111111111 -21.9728532782553 9.75063105603303 0.604309699204896 WT-OE 1.4126984126984 -15.5442172292441 18.3696140546409 0.976169656924344
擴增子統計繪圖1箱線圖：Alpha多樣性