写了个函数rrnDBcorrectOTU：根据16S拷贝数校正OTU---生物信息知识分享

为了给以上两个小目标的达成一些奖励，这两天不务正业了一下，写了一个函数打包上传到了Github上，没错，就是rrnDBcorrectOTU。

rrnDB数据库2014年发表在Nucleic Acids Research，收集了NCBI上细菌和古菌的16S rRNA拷贝数。

由于很多物种有超过一个的16S rRNA拷贝数，PCR扩增的时候16S rRNA拷贝数多的物种扩增出的序列更多，测序数据会偏离真实的物种数量。

目前一些文章已经利用rrnDB数据库校正得到的OTU。

关于rrnDB数据库的介绍及用法，可参见刘永鑫老师的博客文章：

https://blog.csdn.net/woodcorpse/article/details/84437097

我做的就是输入OTU表，RDP classifer的结果及rrnDB数据库，可以得到OTU、物种信息及对应拷贝数信息的整合结果；以及利用rrnDB矫正之后的OTU。

目前的代码写的比较低端且臃肿，但是已经可用了~

虽然我现在也用不上，只图快乐~

rrnDB数据库可以从这里下载~

https://rrndb.umms.med.umich.edu/static/download/

rrnDB格式如下，参考了一下以往文献，基本都利用mean，即平均拷贝数进行计算。

冲冲冲！

#安装 devtools::install_github("Listen-Lii/rrnDBcorrectOTU") library(rrnDBcorrectOTU) #三个输入文件 #OTU otu = read.table(file="otu.tabular",header=T,row.names=1,sep="t") #RDP classifer 文件 classifer = read.table(file="Classifier.txt",header=T,row.names=1,sep="t") #下载的rrnDB数据库 rrnDB = read.csv(file="rrnDB-5.6_pantaxa_stats_RDP.tsv",header=T,sep="t") #简单粗暴出结果 result = rco(otu,classifer,rrnDB) #这个表是OTU、物种信息及对应拷贝数信息的整合结果 result$whole.res #这个表是矫正之后的OTU result$correct.table

result$whole.res

第一列是物种，第二列是物种的分类水平，第三列为该物种在该分类水平的平均拷贝数，注意值为-1表示rrnDB数据库中没有该物种的拷贝数信息。第四列开始为样本。

result$correct.table

校正后的OTU表

END

---来自腾讯云社区的---生物信息知识分享