制作RDP数据库的USEARCH版本用于扩增子物种注释

制作RDP数据库的USEARCH版本用于扩增子物种注释

RDP是在扩增子分析中物种注释的常用的数据库之一，USEARCH、VSEARCH是常用的扩增子分析工具。RDP数据库训练集 trainset 18(https://sourceforge.net/projects/rdp-classifier/)仅提供了QIIME使用的数据库格式，

而USEARCH作者Robert Edger提供的RDP数据库(http://www.drive5.com/sintax)仅为几年前的transet 16 (rdp_16s_v16_sp.fa.gz)。

我联系让作者更新一下，居然回复是商业软件申请RDP数据库的授权使用比较麻烦，你们自己制作一下也非常简单。好吧，只能自己动手，丰衣足食了。

数据下载

在 https://sourceforge.net/projects/rdp-classifier/ 下载QIIME格式的序列和物种注释 RDPClassifier_16S_trainsetNo18_QiimeFormat.zip

目前最新的为RDP 18，2020年更新

数据的解压和基本统计

# 解压
unzip RDPClassifier_16S_trainsetNo18_QiimeFormat.zip
cd RDPClassifier_16S_trainsetNo18_QiimeFormat
# 统计序列，21195 seqs, 30M
head -n2 RefOTUs.fa
seqkit stat RefOTUs.fa
# 统计物种注释，3界，67门，120纲，210目，545科，3240属和8682种
head Ref_taxonomy.txt
for i in `seq 1 7`; do
cut -f2 Ref_taxonomy.txt| sed 's/; /\t/g' | cut -f $i |sort|uniq -c|sort -k1,1nr |wc -l
done

RDP训练集18版，包括3界，67门，120纲，210目，545科，3240属和8682种。

观察QIIME(输入)和USEARCH(输出)格式

输入文件：整合前序列和物种注释

# 序列
>AJ000684
gaacgctggcggcgtgcttaacacatgcaagtcgaacggaaaggtctcttcggagatactcgagt

# 物种注释
AJ000684 k__Bacteria; p__Actinobacteria; c__Actinobacteria; o__Mycobacteriales; f__Mycobacteriaceae; g__Mycobacterium; s__heidelbergense

输出文件：整合后样式——注释+序列(参考USEARCH的RDP16: rdp_16s_v16_sp.fa)

>X71857_S000021696;tax=d:Bacteria,p:Firmicutes,c:Clostridia,o:Clostridiales,f:Clostridiaceae_1,g:Clostridium_sensu_stricto,s:Clostridium_puniceum;
TGATCCTGGCTCAGGACGAACGCTGGCGGCGTGCTTAACACATGCAAGTCGAGCGATGAAGTTCCTTCGGGAATGGATTA

制作USEARCH物种注释数据库单文件格式

本质上只需要两行命令即可实现格式的输入。首先用sed调整物种注释格式与目标一致，其次使用awk追加至序列文件的ID后面。结果一定要用head检查，确保符合预期。代码如下：

# 修改注释为uearch格式
sed 's/^/>/;s/k__/tax=d:/;s/__/;/;s/; /,/g;s/$/;/' Ref_taxonomy.txt > Ref_taxonomy_usearch.txt
# 预览
head -n3 Ref_taxonomy_usearch.txt

# 添加注释至序列ID
awk 'NR==FNR{a[$1]=$2}NR>FNR{print $0"\t"a[$1]}' Ref_taxonomy_usearch.txt RefOTUs.fa|sed 's/\t$//;s/\t/;/'|less -S > rdp_16s_v18_sp.fa
# 预览结果，应与上面示例的输出文件完全一致，序列字母大小写对结果无影响
head -n4 rdp_16s_v18_sp.fa|cat -A

最终结果下载链接：https://github.com/YongxinLiu/EasyMicrobiome/blob/main/usearch/rdp_16s_v18_sp.fa.gz ，点击页面中的Download下载。

赞 (0)

相关推荐