美文网首页RNA-seq注释和富集
如何获取GO和KEGG注释、描述信息

如何获取GO和KEGG注释、描述信息

作者: 生信小白杀手 | 来源:发表于2022-07-25 15:35 被阅读0次

    在非模式生物中,我们想要获取相应物种的GO或KEGG注释信息,或者描述信息的渠道很少。数据分布不集中,需要我们花费很长时间去收集,整合,才能得到相应的信息去做功能富集分析。下面以小麦为例,我们如何去收集功能富集前的准备工作。

    一、功能注释

    对一些未知功能的序列做初步的功能注释,通常的方法就是先下载一些功能数据库(NR、NT、uniprot、GO、KEGG),然后将序列进行比对,从而预测功能。但是这些数据库通常比较庞大,下载慢,处理起来也比较复杂。因此,我们一般采用在线做功能注释。有些网站例如eggNOG-MAPPERPANNZERg:Profiler等。我认为在线功能注释网站中比较快的是eggNOG-MAPPER,PANNZER专注做功能注释约有10年了,也是比较专业的,大家也可以多在几个网站中做,然后比较一下结果。取自己想要的,合适的结果。

    下面就以eggNOG-MAPPER网站为例,做以下序列的功能注释:

    示例数据取自小麦中国春v2.1版本的高可信蛋白序列。一般我们以蛋白序列作为功能注释的输入文件,示例文件名:Example_pep.fa


    image.png
    image.png

    上面为填写后的网页,我一般采用默认参数,然后提交。在邮箱中会看到下面界面。


    image.png
    点击后会转到下面界面
    image.png
    开始后只需等待就可以啦
    image.png

    当运行完成之后会给你邮箱发送一个任务完成邮件。
    任务完成啦,然后下载结果。


    image.png
    一般下载 out.emapper.annotations和out.emapper.annotations.xlsx两个结果就够了,其中out.emapper.annotations为后面处理数据的输入文件。
    image.png
    我们来看一下结果文件里面有什么内容
    image.png
    文件内容这么乱,想得到我们需要的输入文件,还需要进行进一步整理,下面就使用TBtools软件进行整理,确实很方便。使用GO&KEGG模块里面的eggNOG-mappper Helper功能进行处理。
    image.png
    输入文件为out.emapper文件
    image.png
    在目录下就会产生每个功能注释库所对应的蛋白质ID
    image.png

    out.emapper.annotations.GO文件里,有GO与蛋白质ID对应关系,这样我们就初步得到我们所需要的文件啦。


    image.png
    image.png
    image.png
    对于out.emapper.annotations.KEGG_Pathway文件中,我们可以用以下命令简单处理:
    grep map out.emapper.annotations.KEGG_Pathway.txt > out.emapper.annotations.KEGG_map.tx
    grep ko out.emapper.annotations.KEGG_Pathway.txt > out.emapper.annotations.KEGG_ko.txt
    

    这样我们就能提取出只有map号和KO号的文件啦。


    image.png
    image.png

    到此我们前景基因的功能注释就完成啦。

    二、获取该物种所对应的GO、KEGG注释

    获取物种与GO号所对应的关系的方法有以下几种方式:

    1、通过GO官网进行下载http://current.geneontology.org/annotations/index.html或者http://current.geneontology.org/products/pages/downloads.html下载对应物种的注释信息。很可惜,官网里面并没有小麦的注释信息。

    2、 从COA项目中下载ftp://ftp.ebi.ac.uk/pub/databases/GO/goa/,但是也没有小麦的注释信息。

    image.png
    3、从NCBI中下载数据ftp.ncbi.nih.gov
    image.png
    下载gene2go.gz文件,打开后根据Tax id进行筛选,小麦的Tax id是4565,很可惜也没有。
    4、从Bioconductor 获取,虽然AnnotationHub包提供了常见物种的注释信息,却没有小麦的注释信息。
    5、将小麦参考基因组蛋白质序列与uniprot蛋白质库进行blast,然后得到GO号。由于uniprot库下载速度慢,并且比较大,耗时耗力。我采用的是下面的方法得到背景基因。
    6、将参考蛋白质序列在eggNOG-mapper在线网站进行比对,然后对比对结果进行处理,这样就能得到GO号以及KEGG的K num、map号、ko号等等。
    image.png

    三、GO、KEGG描述信息

    GO描述信息可以在http://current.geneontology.org/ontology/go-basic.obo得到,是所有物种的描述信息,不影响做GO富集分析。进入下面页面后右击点击另存为,保存文件。然后提取文件中信息即可得到简化的描述信息。

    image.png
    使用下面python进行处理一下:
    with open("go-basic.obo","r") as file:
        lib={}
        for line in file:
            line=line.strip()
            col_name=line.split(":")[0]
            if col_name == "id":
                id=line.split(" ",maxsplit=1)[1]
                lib[id]=""
            if col_name == "name":
                name=line.split(" ",maxsplit=1)[1]
                lib[id]=lib[id]+"@"+name
            if col_name == "namespace":
                namespace=line.split(" ",maxsplit=1)[1]
                lib[id]=lib[id]+"@"+namespace
    out=open("GO_basic_Description.txt","a+")
    out.write("Class"+"\t"+"GO_IDs"+"\t"+"Description"+"\n")
    for key in lib.keys():
        go_id=key
        go_name=lib[key].split("@")[1]
        go_namespace=lib[key].split("@")[2]
        if go_namespace == "molecular_function":
            go_namespace="MF"
            out.write(go_namespace+"\t"+go_id+"\t"+go_name+"\n")
        if go_namespace == "biological_process":
            go_namespace="BP"
            out.write(go_namespace+"\t"+go_id+"\t"+go_name+"\n")
        if go_namespace == "cellular_component":
            go_namespace="CC"
            out.write(go_namespace+"\t"+go_id+"\t"+go_name+"\n")
    

    得到下面的文件 GO_basic_Description.txt:

    image.png
    KEGG描述信息同样也可以在KEGG官网中得到。参考网址https://www.kegg.jp/kegg/catalog/org_list.html
    image.png
    然后找到小麦(Triticumaestivum)缩写为taes。点击raes进入小麦专栏
    image.png
    点击Brite hierarchy
    image.png
    然后在点击KEGG Orthology
    image.png
    然后在download栏中选择一种格式进行下载
    image.png
    json格式
    image.png
    KEG格式
    image.png
    发现上面两个格式不是我们想要的,还需要进行处理,显得有点麻烦。
    下面推荐一种较为简便的方式,使用R包“KEGGREST”,这个包就比较简便,人性化。
    下面代码就是提取相应的信息。只需要修改输出文件路径即可。
    #获取KEGG数据库信息
    #加载包
    library(KEGGREST)
    ##查看KEGG数据库包含的数据
    listDatabases()#
    ##获取pathway(所有物种)数据集中的数据
    pathway<- keggList("pathway")
    head(pathway)
    #转换数据集,导出数据集
    pathway_data<-as.data.frame(pathway)
    write.table(pathway_data,"<path>/KEGG_pathway_allspacies_database.txt",row.names = T,col.names = F,sep = "\t")
    ##对单个(小麦)数据库进行物种的选择
    taes_pathway <-keggList("pathway","taes")#taes是小麦的缩写
    taes_pathway_data<-as.data.frame(taes_pathway)
    write.table(taes_pathway_data,"<path>/KEGG_pathway_wheat_database.txt",row.names = T,col.names = F,sep = "\t")
    ##获取KO(所有基因)数据集中的数据
    Ko<-keggList("ko")
    Ko_data<-as.data.frame(Ko)
    write.table(Ko_data,"<path>/KEGG_KO_allspacies_database.txt",row.names = T,col.names = F,sep = "\t")
    ##获取KO(小麦)数据集中的数据
    taes_Ko<-keggList("ko","taes")
    taes_Ko_data<-as.data.frame(taes_Ko)
    write.table(taes_Ko_data,"<path>/KEGG_KO_wheat_database.txt",row.names = T,col.names = F,sep = "\t")
    

    得到文件如下:
    KEGG_pathway_allspacies_database.txt
    KEGG_pathway_wheat_database.txt
    KEGG_KO_allspacies_database.txt
    KEGG_KO_wheat_database.txt
    这个包不仅可以提取map号、ko号还可以提取其他的东西。下图就是所有可提的部分,只要修改以下keggList("pathway","taes")参数即可。


    image.png

    得到KEGG的描述信息后需要用一个python脚本简单处理一下,输入文件为KEGG_pathway_allspacies_database.txt、KEGG_KO_allspacies_database.txt

    #将得到的KO、pathway文件进行处理,得到可以富集的描述信息文件
    with open("KEGG_KO_allspacies_database.txt","r") as file:
        out=open("KEGG_KO_allspacies_description.txt","a+")
        for line in file:
            line=line.strip()
            if "[EC:" in line.split("\t",maxsplit=1)[1] :
                desc=line.split("\t",maxsplit=1)[1].rsplit("[EC:",maxsplit=1)[0].split(";")[1].strip('"')
                id=line.split("\t",maxsplit=1)[0].split(":")[1][:-1]
                out.write(id+"\t"+desc+"\n")
            else:
                if ";" in line.split("\t",maxsplit=1)[1]:
                    out.write(line.split("\t",maxsplit=1)[0][:-1].split(":")[1]+"\t"+line.split("\t",maxsplit=1)[1].split(";")[1].strip('"')+"\n")
                else:
                    out.write(line.split("\t",maxsplit=1)[0][:-1].split(":")[1]+"\t"+line.split("\t",maxsplit=1)[1].strip('"')+"\n")
    with open("KEGG_pathway_allspacies_database.txt") as file1:
        out1=open("KEGG_pathway__allspacies_description.txt","a+")
        for line1 in file1:
            line1=line1.strip()
            map_id=line1.split("\t") [0].strip('"').split(":")[1]
            descript=line1.split("\t") [1].strip('"')
            out1.write(map_id+"\t"+descript+"\n")
    

    输出文件为KEGG_KO_allspacies_description.txt、KEGG_pathway__allspacies_description.txt


    image.png
    image.png

    这样我们GO、KEGG富集前准备文件就完成啦!

    相关文章

      网友评论

        本文标题:如何获取GO和KEGG注释、描述信息

        本文链接:https://www.haomeiwen.com/subject/payuirtx.html