Orthomcl使用说明

结果文件处理：/home/liqiang/my-uxe/10/orthomcl2species.pl

                            /home/liqiang/my-uxe/10/orthogenes_statistics.pl

1.Orthomcl是用来干嘛的？
  我们都知道，我们在注释完一些基因组后，会拿到大量的蛋白，而这些蛋白，孤零零的放在那，你不觉得别扭吗？除了常规的pfam、interpro做完了，再做点什么呢？聚个类，找个家族，看看家族水平上的进化什么的貌似也蛮不错，今天就介绍一下如何利用orthomcl来做蛋白家族吧。
OrthoMCL groups proteins into “ortholog groups.” That name is a little misleading
because the groups contain proteins related by:
orthology (recent descent)
in-paralogy (recent duplication)
co-orthology (recent descent and duplication).
Orthomcl会将你所给的各个种的蛋白使用mcl的方法聚到一起，具有很强的鲁棒性
 
2.去哪下载？
用最新版本的orthomclSoftware-v2.0.3.tar.gz吧
http://orthomcl.org/common/downloads/software/v2.0/
 
3.安装
这个就不说了，解压即可
 
4.所需辅助模块
（1）UNIX：OrthMCL Pairs program 仅仅在UNIX上进行过测试
           MCL program也仅仅和UNIX兼容
 
（2）BLAST：我们推荐NCBI BLAST，有两个原因，从理论上来讲，鼻祖；从实际运行商来讲，NCBI BLAST支持 
     tab分割的output格式，这个我们将在下面的第7步用到.
 
（3）Database：从V1.4版，Orthomcl就开始用数据库的方式管理数据，这样大大的加强和提升了运行效率，在
     这你可以选择mysql或者oracle，下文我选择mysql来进行阐述.
 
（4）Hardware：推荐：memory至少4G，disk至少100G
 
（5）Perl：需要standare perl和DBI
 
（6）MCL program：参照 http://www.micans.org/mcl/sec_description1.html
 
（7）Time：时间复杂度：all-v-all Balst在500个cpu的cluster上大概运行3天
                       使用orthmclPaire找pairs大概需要16h
                       MCL找groups大概需要2h
 
5.详细过程
（1）安装和配置相关的数据库文件：这个靠你自己了。
 
（2）安装mcl
从这里下载 http://www.micans.org/mcl/src/mcl-latest.tar.gz
 
安装手册: http://www.micans.org/mcl/sec_description1.html
（3）安装Orthmcl：
tar -zxvf orthomclSoftware-v2.0.3.tar.gz
解压完之后，你将会看见这样的文件结构
  orthomclSoftware/
    bin/
     ...
    doc/
     UserGuide.txt
     orthomcl.config.template
    lib/
bin下面就是下面你要做各步骤的一些脚本
将doc下的orthomcl.config.template配置文件模版拷出来，修改里面的参数配置
# this config assumes a mysql database named 'orthomcl'.  adjust according
# to your situation.
dbVendor=mysql #你所选择的数据库，在这我选择使用mysql
dbConnectString=dbi:mysql:orthomcl:node69:3310 #设置你使用的数据库和hostname及其使用端口，默认是
                                        3307，在这由于我的服务器上这些端口都被占了，所以我选择3310
dbLogin=xiaenhua #你mysql的用户名
dbPassword=xeh #密码
similarSequencesTable=SimilarSequences #下面都是中间产生的各种表
orthologTable=Ortholog
inParalogTable=InParalog
coOrthologTable=CoOrtholog
interTaxonMatchView=InterTaxonMatch
percentMatchCutoff=50  #Coverage cutoff值  这里选择50%的Coverage，视你自己而定
evalueExponentCutoff=-5  #blast 筛选的e-value  用过blast的都不默认
oracleIndexTblSpc=NONE
（4）orthomclInstallSchema
这一部分就是将你刚才配置的config文件，对mysql进行配置，建立在你所create的database下，建立一些表，Note：在做这步前，请先在你的mysql中新建一个数据库，如create database orthomcl，下面我就使用这个数据库来操作数据。
 EXAMPLE: orthomclSoftware/bin/orthomclInstallSchema my_orthomcl_dir/orthomcl.config my_orthomcl_dir/install_schema.log
 
（5）orthomclAdjustFasta
该步将会将你的pep文件转换为orthmcl所要求的文件，其实也就是一个改写的过程，格式：
EXAMPLE: orthomclSoftware/bin/orthomclAdjustFasta hsa Homo_sapiens.NCBI36.53.pep.all.fa 1
注：hsa表示你的种名
    hsa Homo_sapiens.NCBI36.53.pep.all.fa你的蛋白序列
    1 表示将在你的每个蛋白名前加上种名，并且用|隔开
 
（6）orthomclFilterFasta
这一步将会对你刚才改写的蛋白文件进行过滤，去除长度小于XX（自己设定），stop coden所占百分比的序列
格式：
EXAMPLE: orthomclSoftware/bin/orthomclFilterFasta my_orthomcl_dir/compliantFasta 10 20
注：10 pep序列长度不低于10
    20 stop coden所占的百分比不高于20%
 
（7）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提供相似矩阵 
 
（8）orthomclBlastParser
将上一步得到的blast比对结果进行解析，使用我们开始设好的阈值进行筛选，e-value：1e-5 ；Coverage：50%
EXAMPLE: orthomclSoftware/bin/orthomclBlastParser my_blast_results my_orthomcl_dir/compliantFasta >> my_orthomcl_dir/similarSequences.txt
 
（9）orthomclLoadBlast 
这一步，将我们刚才解析好的blast结果导入到mysql中，便于下面的数据操作
EXAMPLE: orthomclSoftware/bin/orthomclLoadBlast my_orthomcl_dir/orthomcl.config my_orthomcl_dir/similarSequences.txt
在这需要提供我们先前配置好的config文件
 
（10）orthomclPairs
这一步，将在database中SimilarSequences表中的数据，进行pairs的运算，产生三个表格存在mysql
  - PotentialOrthologs table
  - PotentialInParalogs table
  - PotentialCoOrthologs table
EXAMPLE: orthomclSoftware/bin/orthomclPairs my_orthomcl_dir/orthomcl.config my_orthomcl_dir/orthomcl_pairs.log cleanup=no
 
（11）orthomclDumpPairsFiles
这一步，将数据库中pairs表进行处理，生成mclInput文件和另外一个文件夹pairs，在这个pairs中，包含着这些蛋白之间的关系，格式如下：
   - protein A
   - protein B
   - their normalized score (See the Orthomcl Algorithm Document).
 
（12）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文档
 
（13）orthomclMclToGroups
将mcl的输出结果转换为groups.txt
在这个文件中，每一行表示一个家族
EXAMPLE：orthomclMclToGroups my_prefix 1 < mclOutput > groups.txt
注：my_prefix 指定在groups.txt中每个家族的前缀，如：GF_ 则在groups.txt中，每个家族以GF_开始
    1 表示家族从1开始编码
 
到这一步位置，你的蛋白家族就算是完工了，groups.txt中的格式如下：
GF_1: r5|r5_37273 r5|r5_7773 r5|r5_39887 r5|r5_42220 r5|r5_40330 r5|r5_37989 r5|r5_40295
GF_2: r5|r5_4187 r5|r5_37986 r5|r5_8138 r5|r5_38203 r5|r5_3914 r5|r5_9613 r5|r5_9656 r5
GF_3: .....
 
每行代表一个家族，找到了这些家族，下面就可以开始你的其他分析了