導(dǎo)讀
-
菌群分類數(shù)據(jù)的可視化,包括豐度在個體之間以以及組間的一個比較怎么生動形象肋杖,一直是我想解決的問題永罚,今天看到一篇文章上面用的是ggplot2包做的一個面積圖,主要利用到了里面的geom_area()函數(shù)宋舷,這樣的圖有一個好處就是可以看菌群豐度隨著采樣時間點的變化而變化。那么我在想瓢姻,如果橫坐標是個體號呢祝蝠,怎么想辦法也能像下圖(圖1)這樣展示。所以我想通過自己平時下載的數(shù)據(jù)來重現(xiàn)這個圖。好了绎狭,話不多說细溅,開始作圖!
作圖流程
文件準備
- 這里用到的是R包
Rhea
里構(gòu)建的數(shù)據(jù)坟岔,用的是科水平,數(shù)據(jù)如圖(圖2)如下:
這個表的行名為科水平分類名摔桦,列名為樣品名社付,但是表格的具體形式不一定要和我的一致,行列隨意邻耕,反正可以隨意轉(zhuǎn)置
R代碼
rm(list = ls())
# 加載系列包
library(reshape2)
library(ggplot2)
library(phyloseq)
library(magrittr)
# 加載數(shù)據(jù)
setwd("D:\\元基因組\\3. 16S 擴增子下游測序分析&其他分析方法集錦\\Rhea_amplicon_analysis\\4.Taxonomic-Binning")
otu_file<-"Families.all.tab"
otu_table <- read.table (otu_file,
check.names = FALSE,
header = TRUE,
dec = ".",
sep = "\t",
comment.char = "")
otu_table <- otu_table[1:10,] # 為了展示方便,這里只取前面10個科
otu_table[,2:ncol(otu_table)] <- sweep(otu_table[,2:ncol(otu_table)],2, colSums(otu_table[,2:ncol(otu_table)]),'/')*100
col <- c("f__Bacteroidaceae" = "#E5CD34", "f__Deferribacteraceae" = "#6C6011",
"f__Desulfovibrionaceae" = "#2B5B51", "f__Enterobacteriaceae" = "#1D9376", "f__Erysipelotrichaceae" = "#A3E4D5",
"f__Lachnospiraceae" = "#4D2599", "f__Lactobacillaceae" = "#1E0748", "f__Peptostreptococcaceae" = "#9C74E3",
"f__Porphyromonadaceae" = "#6B2B0B", "f__Rikenellaceae" = "#E46628")
# 將原始數(shù)據(jù)表進行轉(zhuǎn)換兄世,將寬矩陣轉(zhuǎn)為長矩陣
otu_table_1 <- melt(otu_table)
names(otu_table_1)[1] <- "family"
otu_table_1$variable <- as.factor(otu_table_1$variable)
p <- ggplot(otu_table_1,aes(x=variable, y = value)) +
geom_area(aes(group=family, fill = family),position = "stack") + ##這里要添加group和fill,不然個體之間連接不成面積圖啼辣,而是和條形圖差不多。并且要選擇stack是堆砌圖御滩。
#facet_wrap(~vessel) +
theme_minimal() +
scale_fill_manual(values = col) +
#scale_x_datetime(date_labels="Day %d") +
ylab("Proportions") +
#guides(fill = element_text(size=16))+
theme(panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
legend.position="bottom",
axis.title.x = element_blank(),
axis.text.x = element_text(angle=45,color="black",vjust = 0.95,hjust = 0.95,size=12),
axis.title.y = element_text(size=16),
axis.text.y = element_text(size=14),
strip.text.x = element_text(size=18),
legend.text = element_text(size=14),
legend.title = element_text(size=16))
ggsave("family.tiff", height=8, width=20, units="in")
結(jié)果展示
-
結(jié)果展示如下圖所示(圖3)