其他
R语言绘制物种Rank-abundance曲线
在微生物16S/ITS/18S测序分析报告中，我们经常可以看到这样的图，称为Rank-abundance曲线。
Rank Abundance曲线是将样品中的OTUs（可以粗略理解为“物种”）按相对丰度（或者包含的序列数目）由大到小排序得到对应的排序编号，再以OTUs的排序编号为横坐标，OTUs中的相对丰度（也可用该等级OTU中序列数的相对百分含量）为纵坐标，将这些点用折线连接，即绘制得到Rank Abundance曲线，它可直观的反映样品中物种的丰富度和均匀度。在水平方向上，物种的丰富度由曲线的宽度来反映，物种的丰富度越高，曲线在横轴上的跨度越大；在垂直方向上，曲线的平滑程度，反映了样品中物种的均匀程度，曲线越平缓，物种分布越均匀。
下图中，横坐标代表各OTU的排序编号，纵坐标代表各OTU的相对丰度，样本曲线的延伸终点的横坐标位置为对应样本的测序数量。如果曲线越平滑下降表明样本的物种多样性越高，而曲线快速陡然下降表明样本中的优势菌群所占比例很高，多样性较低。
本篇白鱼小编继续带大家绘制这种Rank-abundance曲线。示例数据、R脚本等，已上传至百度盘（提取码488s）。https://pan.baidu.com/s/1HInewdVJ-T15hPHxI8gmYg 网盘文件“otu_table.txt”为某16S扩增子测序所得OTU丰度表格，可以理解为物种丰度表。表中每一列为一个样本，每一行为一种OTU，交叉区域为每种OTU在各样本中的丰度。
Rank-abundance曲线绘制
otu <- read.delim('otu_table.txt', row.names = 1, sep = '\t', stringsAsFactors = FALSE, check.names = FALSE)
otu <- t(otu)
#统计，BiodiversityR 包 rankabundance() 实现 OTU 排序
#虽然自己排序也很简单，但还是导个包统计省事点……
library(BiodiversityR)
otu_relative <- otu / rowSums(otu) #转化为相对丰度
rank_dat <- data.frame()
for (i in rownames(otu_relative)) {
rank_dat_i <- data.frame(rankabundance(subset(otu_relative, rownames(otu_relative) == i), digits = 6))[1:2]
rank_dat_i$sample <- i
rank_dat <- rbind(rank_dat, rank_dat_i)
}
rank_dat <- subset(rank_dat, abundance != 0)
#ggplot2 作图，更好的可视化效果，还请自行修改 ggplot2 作图细节
library(ggplot2)
ggplot(rank_dat, aes(rank, log(abundance, 10), color = sample)) +
geom_line() +
scale_colour_manual(limits = c('a1', 'a2', 'a3', 'a4', 'a5', 'a6'), values = c('orange', 'purple', 'green', 'blue', 'red', 'gray40')) +
labs(x = 'OTUs rank', y = 'Relative adundance (%)', color = NULL) +
theme(panel.grid = element_blank(), panel.background = element_rect(fill = 'transparent', color = 'black'), legend.key = element_rect(fill = 'transparent')) +
scale_y_continuous(breaks = 0:-5, labels = c('100', '10', '1', '0.1', '0.01', '0.001'), limits = c(-5, 0))
事实上，BiodiversityR包也自带了一些函数，如rankabunplot()、rankabuncomp()等，可以自动绘制Rank-abundance曲线。不过貌似效果比较一般。
#BiodiversityR 包自带的作图函数#例如
rankabuncomp(otu, y = data.frame(name = rownames(otu)), type = 'l', factor = 'name', scale = 'logabun', legend = FALSE)