Tax4Fun: 16S 预测群落代谢功能实例

这一期我们用下载后的数据进行代谢功能预测,为分析自己的数据做好准备。

将OTUs丰度表转成标准格式

对于Tax4Fun软件,其输入为OTU丰度表,且各OTU的注释信息必须来自SILVA数据库。此外Tax4Fun目前具有3个版本,分别对应了3个SILVA物种注释数据集版本,安装注释的时候最好二者相互对应。所以,为了预测结果的准确性,首先需要保证自己的OTUs注释信息来自于SILVA数据库的注释。如果自己手里得到的OTU注释信息并非来自SILVA数据库,比如来自Greengenes、RDP等,那就需要使用SILVA数据库进行重新注释,且最好选择相应的Tax4Fun版本。

以下示例选择Tax4Fun的SILVA123参考数据集进行细菌群落功能预测,因此我们也最好使用SILVA数据库中的SILVA123数据集对OTU进行物种注释。

使用QIIME为OTU添加注释

首先,对于所有非SILVA数据库注释的OTU物种信息,都可以可以使用这些OTU的代表序列重新进行注释。这里使用到:1、无注释信息的OTU丰度表(otu_table.txt),2、各OTU代表序列文件(otu.fasta)。使用QIIME分析平台完成注释功能。

我们使用各OTU的代表序列与SILVA123数据库中的已知物种序列进行比对,找到每个OTU序列的最相似序列,这样得到得到每个OTU所对应的物种分类单元(OTU注释信息)。同时将OTU丰度表转化为便于QIIME识别的biom格式,并将各OTU注释信息添加至OTU表的最后一列。最终得到具有SILVA123物种注释信息的OTU丰度表格。

备注:以下QIIME每个步骤所得的中间结果均放置在网盘附件“qiime_result”中,可在文末查看。

##注意:以下步骤均在 Linux环境下中完成,运行程序或脚本均来自QIIME自带命令#在 otu_table.tsv 开头添加一行“# Constructed from biom file”,以便将 otu_table.tsv 转为 qiime 可识别样式cp otu_table.txt otu_table.tsvsed -i '1i# Constructed from biom file' otu_table.tsv #将 otu_table.tsv 转换为 otu_table.biom,需要提前安装QIIME软件biom convert -i otu_table.tsv -o otu_table.biom --table-type="OTU table" --to-json
# 下载SILVA123数据库,网址:http://tax4fun.gobics.dewget  http://tax4fun.gobics.de/SilvaSSURef_123_NR.zipunzip SilvaSSURef_123_NR.zip
#OTU 注释,输出结果 otu_table_tax_assignments.txt 即为注释文件assign_taxonomy.py -i otu.fasta -r SILVA_123_SSURef_Nr99_tax_silva.fasta -t SILVA_123_SSURef_Nr99_tax_silva.tax -o ./ #添加 otu 注释信息至 biom 格式的 otu 表(otu_table.biom )的最后一列,并将列名称命名为 taxonomybiom add-metadata -i otu_table.biom --observation-metadata-fp otu_tax_assignments.txt -o otu_table.silva.biom --sc-separated taxonomy --observation-header OTUID,taxonomy  #otu_table.silva.biom 转换为 otu_table.silva.tsvbiom convert -i otu_table.silva.biom -o otu_table.silva.tsv --to-tsv --header-key taxonomy

根据以上步骤,最终得到带有注释信息的OTU丰度表,文件名为:otu_table.silva.tsv,第一行为“# Constructed from biom file”;最后一列为OTU注释信息,其余行列即为各OTU名称、各样本名称、以及各OTU在各样本中的丰度分布。

otu_table.silva.tsv 文件截图
Tax4Fun预测落功能

接下来我们正式运行Tax4Fun进行功能预测,将得到KEGG功能基因及代谢通路的丰度。

Tax4Fun会根据OTUs丰度及注释信息,在其数据库中找到对应的物种类别,根据16S拷贝数将物种数据标准化后,统计所对应物种类别其本身所具有的功能基因即所含代谢通路的数量,预测细菌群落中所具有的功能基因丰度以及代谢通路丰度。因此我们可知,OTU注释结果中,若注释信息越详细(分类水平越细,越接近种(species)水平),则理论上功能注释结果的准确度越高。

首先进入R环境,在Linux下输入R。

##加载 Tax4Fun 包library(Tax4Fun) ##Tax4Fun 默认流程#此处使用“otu_table.silva.tsv”作为输入文件QIIMESingleData <- importQIIMEData('otu_table.silva.tsv')## SILVA123 需要提前下载,下载命令 wget http://tax4fun.gobics.de/Tax4Fun/ReferenceData/SILVA123.zipfolderReferenceData <- 'SILVA123' #功能基因丰度预测Tax4FunOutput <- Tax4Fun(QIIMESingleData, folderReferenceData, fctProfiling = TRUE, refProfile = 'UProC', shortReadMode = TRUE, normCopyNo = TRUE)tax4fun_gene <- as.data.frame(t(Tax4FunOutput$Tax4FunProfile))gene <- rownames(tax4fun_gene)write.table(cbind(gene, tax4fun_gene), 'tax4fun.gene.txt', row.names = FALSE, sep = 't', quote = FALSE) #代谢通路(KEGG 代谢第三层级)丰度预测Tax4FunOutput <- Tax4Fun(QIIMESingleData, folderReferenceData, fctProfiling = FALSE, refProfile = 'UProC', shortReadMode = TRUE, normCopyNo = TRUE)tax4fun_pathway <- as.data.frame(t(Tax4FunOutput$Tax4FunProfile))pathway <- rownames(tax4fun_pathway)write.table(cbind(pathway, tax4fun_pathway), 'tax4fun.pathway.txt', row.names = FALSE, sep = 't', quote = FALSE)
上面是我自己运行的命令,在R环境下

上面两个流程里面,一个重要参数“fctProfiling”可分别设定为“TRUE”和“FALSE”:

  • 当“fctProfiling = TRUE”时,将输出功能基因丰度预测结果;
  • 当“fctProfiling = FALSE”时,将输出代谢通路丰度预测结果。

所以我们运行两次,分别得到功能基因和代谢通路的结果,其他参数一般情况下无需更改,默认即可。

对于功能基因丰度预测结果文件名字为:tax4fun.gene.txt,内容截图如下图所示,第1列为各功能基因的ID(可根据此ID在KEGG数据库中检索以查看其详细功能信息等)、名称及其对应的酶类型(酶命名规则可参考https://enzyme.expasy.org/enzyme_details.html),第2列及之后为各样本中各功能基因的相对丰度。

tax4fun.gene.txt文件部分图

同样,代谢通路丰度预测结果文件名字为:tax4fun.pathway.txt,内容如下图所示,第1列为各KEGG代谢通路KO编号(可根据此ID在KEGG数据库中检索以查看其详细功能信息等)及其名称,第2列及之后为各样本中各代谢通路的相对丰度。注意,此处的代谢通路为KEGG第三层级代谢通路结果。

查询每个KO的具体功能


上面两个文件就是Tax4Fun的预测内容,怎么知道每个KO或者pathway的具体详细的信息呢?需要到KEGG官网查看,网址为:https://www.genome.jp/kegg/pathway.html

例如,查询代谢通路“Glycolysis / Gluconeogenesis”了解详情,即“ko00010”,过程和结果如下所示。

 

 

 

点击后下图就是K00010的详细信息。

 

之后就是差异分析了,今天先介绍到这里,下期为大家呈现功能基因或者代谢通路的差异分析。

参考文献:

Aßhauer K P, Wemheuer B, Daniel R, et al. Tax4Fun: predicting functional profiles from metagenomic 16S rRNA data[J]. Bioinformatics, 2015, 31(17): 2882-2884.

http://wap.sciencenet.cn/blog-3406804-1162514.html

示例数据下载地址

百度网盘:链接:https://pan.baidu.com/s/1gLH0gyfkQ36a3yoYQjBOeQ 提取码:q0t6

蓝奏云:https://www.lanzous.com/b00z9evpa 密码:d487

生物信息学

用SnpSift过滤VCF文件

2020-3-1 23:58:04

生物信息学

R语言实战(18) 处理缺失数据的高级方法

2020-3-2 19:14:31

声明 本网站部分文章源于互联网,出于传递更多信息和学习之目的转载,并不保证内容正确或赞同其观点。
如转载稿涉及失效、版权等问题,请立即联系管理员;我们会予以修改、删除相关文章,请留言反馈
Notice: When your legal rights are being violated, please send an email to: [email protected]
搜索