準(zhǔn)備數(shù)據(jù)禾酱，即OTU table，原始的或者rare的均可资厉。
下面是繪圖代碼：
```
library(tidyverse)
library(ggthemes)
library(ggsci)
#加載數(shù)據(jù)
otutab <- read_tsv("otutab.txt")
# 定義rank函數(shù)（需要排除為0的OTU）
myrank <- function(x){
? res <- length(x)-rank(x, ties.method = "first")+1
? res[x==0] <- NA
? return(res)
}
# 生成rank矩陣
otutab_rank <- map_dfr(otutab, myrank)
rownames(otutab_rank) <- rownames(otutab)
# 重塑數(shù)據(jù)
otutab2 <- otutab %>% rownames_to_column("FeatureID") %>% gather("Group", "Count", -FeatureID)
otutab_rank2 <- otutab_rank %>% rownames_to_column("FeatureID") %>%
? gather("Group", "Rank", -FeatureID)
# 合并otutab數(shù)據(jù)和ranks數(shù)據(jù)
data <- left_join(otutab2, otutab_rank2) %>% group_by(Group) %>%
? mutate(Abundance = Count/sum(Count)) %>% na.omit()
# 繪制 Rank-Abundance泄朴，注意y軸數(shù)據(jù)進(jìn)行對(duì)數(shù)變化
ggplot(data, aes(x = Rank, y = Abundance, color = Group)) + geom_line() +
? scale_y_log10() +
?scale_color_jco(name = NULL) +
? labs(title="Rank-Abundance Curves") + ylab("% Abundance")+
? theme_bw() +
? theme(panel.grid=element_blank(), plot.title = element_text(lineheight=2.5, face="bold",hjust=0.5))
```
數(shù)據(jù)概覽和結(jié)果
原始數(shù)據(jù)otutab:
otutab_rank:
繪圖數(shù)據(jù)data:
Rank-Abundance曲線
