zoukankan      html  css  js  c++  java
  • OrthoMCL及orthofinder 两种软件进行聚类分析

    OrthoMCL介绍 

    1. OrthoMCL的用途

    基于序列的相似性,OrthoMCL能将一组proteins(比如全基因组的proteins)归类到ortholog groups、in-paralogs groups和co-orthologs。

    2. OrthoMCL-DB

    OrthoMCL-DB包含了很多proteins,这些proteins来自一些已经完全测序的真核或原核生物的基因组。

    OrthoMCL-DB将这些proteins进行了聚类,分成很多的ortholog groups。

    2010.5.31,发布了OrthoMCL-DB第4版,包含 116,536个ortholog groups、1,270,853个proteins、88个真核生物基因组、16个古菌基因组、34个细菌基因组。

    2011.5.31,发布了OrthoMCL-DB第5版,包含 124,740个ortholog groups、1,398,546 个proteins、150个基因组。

    2013年末即将发布OrthoMCL-DB第6版。

    3. OrthoMCL的两种使用方法

    1. OrthoMCL-DB的官网已经将数据中的proteins进行了ortholog的聚类,其网站提供了一个工具,用于接收上传的基因组proteins,再将这些proteins group到相应的ortholog groups中。官网提供的工具Assign your proteins to OrthoMCL Groups用于进行分析。

    2. 如果要对多个基因组的proteomes进行聚类,则可以使用OrthoMCL单机版的软件来进行运算。其用法详见:OrthoMCL的使用。

    4. OrthoMCL算法

    1. 将多个proteomes转换成orthomcl兼容的FASTA文件。

    2. 移除低质量的序列。

    3. All-versus-All BLASTP with 1e-5 cutoff。即使用这些proteomes的protein sequences构建blast数据库,再将所有的这些序列和数据库进行BLASTP比对,取evalue小于1e-5的比对结果。

    4. Filter by percent match length。计算比对结果的percent match length ( 所有hsp中比对上序列的长度之和 / 两条序列中短的那条序列的长度 )。取50%的cutoff值。

    5. 寻找不同物种间potential ortholog pairs(两两物种的protein序列相互是best hits);寻找同一物种内in-paralog pairs(相互之间是better hits,即对于2个序列之中的任意一条序列,和其in-paralog序列之间的evalue值 <= 这条序列和其它物种比对的evalue值).

    6. 根据上一步结果寻找co-ortholog pairs(pairs connected by orhthology and in-paralog,并且pairs之间的evalue值低于1e-5).

    7. 对所有的pairs进行E-values的Normalization,以利于下一步MCL的计算。见下一部分内容,或参考OrthoMCL Algorithm Document。

    8. 将所有的ortholog,in-paralog和co-ortholog pairs,以及它们的标准化后的weight值输入到MCL程序中,来进行聚类分群。MCL documentation

    All-v-all BLAST
    这一步你必须自己做BLAST,即:将你上一步得到的goodProteins.fasta进行多对多的blast,参数建议设置
    -m 8 -F F -b 1000 -v 1000 -a 2
    EXAMPLE:blastall -p blastp -i goodProteins.fasta -d goodProteins.fasta -m 8 -F F -b 1000 -v 1000 -a 2 -o all_VS_all.out.tab
    这一步事实上为MCL提供相似矩阵 

    mcl
    这一步开始对上一步给出的输出文件,进行mcl操作,开始聚类
    EXAMPLE: orthomclSoftware/bin/orthomclDumpPairsFile my_orthomcl_dir/orthomcl.config
    输出文件为mclOutput文件
    EXAMPLE: mcl my_orthomcl_dir/mclInput --abc -I 1.5 -o my_orthomcl_dir/mclOutput
    这里比较重要的参数是-I 具体看mcl文档

    第7步Blast,是整个过程中最关键的一步。有以下2点需要注意:

    1. 数据库中的蛋白质序列数量:在OrthoMCL DB中选取和要分析的物种亲缘关系较近的几个物种的基因组,或下载其它公布的基因组,加上要分析的物种的基因组;使用这些基因组总体的蛋白质序列来构建Blast数据库。如果只是使用要分析的物种的蛋白质序列建数据库,则inparalogs文件中成对的序列实际上是paralogs,数目比真正的inparalogs要多很多。使用所有的OrthoMCL DB中的序列,第5版含150个基因组,信息量太大,不使用几百个核的超算或计算机集群去运行,是很不现实的。

    2. 对数据库中所有的蛋白质序列来使用blast比对到该数据库中得到结果。如果只是对要分析的物种进行Blast,则只能得到inparalogs的信息,而没有orthologs和coorthologs。

    4 分析过程

    4.1 输入文件格式转化

    orthomcl的输入文件为fasta格式的基因或蛋白序列,fasta文件的序列名称要求以样品名开头之后接’|’分隔,之后接每个样品的序列名(如例1),而且样品名和序列名不能有重复。

    命令:orthomclAdjustFasta程序,将fasta文件转换出兼容orthomcl的fasta文件使用命令:(1)orthomclAdjustFasta    A(B,ref)   X1(X2,X3).fa   1,结果输出为A(B,ref).fasta。(单个跑完再合并)。本文生成样品A,B和参考序列ref为例,在compliantFasta文件夹中的序列文件名分别为:A.fa,b.fa,ref.fa。

    例1:

    >A|gene1

    ASSRKSKWQFMGARDAGAKDELRQVYGVSERTESDGAANLIHKLRAINYTLAELGQWCAYKVGQSFLSAL

     >B|contig1

    KDELRQVYGVSERTESD

    4.2 输入文件合并过滤

    使用命令:orthomclFilterFasta  compliantFasta/  10  20。允许的最短的protein长度是10,stop codons最大比例为20%;生成了两个文件(2)goodProteins.fasta和poorProteins.fasta两个文

    4.3 全序列比对

    将上一步的goodProteins.fasta序列进行自身的多序列比对,比对使用软件为blast+,输出结果为all.m8.anno。文件太大可以拆分比对,最后合并
    /share/nas2/genome/bin//blastall -b 1000 -v 1000  -a 2  -p  blastp  -e 1e-5 -F F -d goodProteins.fasta  -i goodProteins.fasta.div1/goodProteins.fasta.f2.106.seq   -o /goodProteins.fasta.div1/goodProteins.fasta.f2.106.seq.blast   -m 8

    cat goodProteins.fasta.f2.*.seq.blast  >(3)all_VS_all.out.tab  还可以去除重复(一列,二列)最后获得 (4)all_VS_all.result

    4.4  导入比对结果
    将比对结果导入mysql数据库,包含以下几个步骤:
    A. 将比对结果转化为规定格式,命名为similarSequences.txt,命令为:(5)orthomclBlastParser all_VS_all.result seq > similarSequences.txt
    B. 将similarSequences.txt导入到数据库中,命令为:orthomclLoadBlast orthomcl.config.template similarSequences.txt

    4.5 寻找paired蛋白
    输入为数据库中的表SimilarSequences,和数据库的空表InParalog, Ortholog, CoOrtholog tables;输出为对这些空表的操作,命令为:orthomclPairs orthomcl.config.template orthomcl_pairs.log cleanup=no。

    4.6 将数据从mysql导出

    生成(6)mcllnput文件和pairs目录。这个目录包含三个文件:
    ortholog.txt, coortholog.txt, inparalog.txt。
    每一个文件有三列: proteinA, protein B, their normalized score (See the Orthomcl Algorithm Document)。

    命令为:orthomclDumpPairsFiles orthomcl.config.template。

    4.7 使用mcl对paired蛋白聚类
    命令为:mcl mclInput --abc -I 1.5 -o (7) mclOutput。

    4.8 对结果编号
    命令为:orthomclMclToGroups gf 1 < mclOutput > (8)groups.txt。家族名为gf_1,gf_2,gf_3...,格式如图2 。

    论到直系同源基因分析的时候,大部分教程都是介绍OrthoMCL,这是2003年发表的一个工具,目前的引用次数已经达到了3000多,但这个软件似乎在2013年之后就不在更新,而且安装时还需要用到MySQL(GitHub上有人尝试从MySQL转到sqlite)。

    而OrthoFinder则是2015年出现的软件,目前已有400多引用。该软件持续更新,安装更加友好,因此我决定使用它来做直系同源基因的相关分析。

    OrthoFinder: solving fundamental biases in whole genome comparisons dramatically improves orthogroup inference accuracy提到,它的优点就是比其他的直系同源基因组的推断软件准确,并且速度还快。

    此外他还能分析所提供物种的系统发育树,将基因树中的基因重复事件映射到物种树的分支上,还提供了一些比较基因组学中的统计结果。

    OrthoFinder的分析过程

    OrthoFinder的分析过程分为如下几步:

    1. BLAST all-vs-all搜索。使用BLASTP以evalue=10e-3进行搜索,寻找潜在的同源基因。(除了BLAST, 还可以选择DIAMOND和MMSeq2)
    2. 基于基因长度和系统发育距离对BLAST bit得分进行标准化。
    3. 使用RBNHs确定同源组序列性相似度的阈值
    4. 构建直系同源组图(orthogroup graph),用作MCL的输入
    5. 使用MCL对基因进行聚类,划分直系同源组



    来源:

    http://www.chenlianfu.com/?p=1839

    https://orthomcl.org/orthomcl/?rm=proteomeUploadForm 

    https://www.jianshu.com/p/16e0bbb2ba19

    https://blog.csdn.net/ych99ych99/article/details/17436367?utm_source=blogxgwz0

    https://www.plob.org/article/6174.html

    https://www.cnblogs.com/renping/p/10126229.html 

  • 相关阅读:
    港湾交换机-交换机配置
    华为交换机-SNMP配置
    H3C交换机-SNMP配置
    CISCO交换机-SNMP配置
    docker容器与宿主交互数据
    MYSQL总览
    mysql数据类型优化
    compact 创建一个包含变量名为数组的键和它们的值为数组的值的数组
    php单例模式
    完形填空格式
  • 原文地址:https://www.cnblogs.com/bio-mary/p/12818863.html
Copyright © 2011-2022 走看看