上海锐翌生物科技有限公司

服务电话:021-51001612

邮箱:support@realbio.cn

技术课堂

使用OrthoMCL软件进行基因家族分析
发布时间:2018-03-13 11:21   点击率:

在基因组学分析当中,我们在对reads进行组装和预测之后会得到样品的基因序列,这些基因序列除了可进行后续的基因数据库注释分析外,还可以进行基因家族分析。将样品与参考基因组的基因序列进行基因家族分析便可分析样品间或样品与参考基因组间的相似性与差异性。今天我们要介绍的OrthoMCL软件就可以简单便捷地进行基因家族分析。


1. 数据库的配置
OrthoMCL的分析需要先行建立mysql账户并建立相应的数据库。关于mysql用户的创建我们不在此进行介绍,我们以已经建立好的账户(账户名user,密码123456)为例进行操作。
A. 在linux环境下输入mysql -uuser -p123456进入mysql界面;
B. 输入create database orthomcl;建立一个名为orthomcl的空数据库用以存放分析时的中间文档;
C. 输入\q退出mysql界面。


2. 分析过程

2.1 软件下载
OrthoMCL的分析需要OrthoMCL软件本体和mcl软件。
OrthoMCL软件下载地址为:http://orthomcl.org/common/downloads/software/,解压缩后,其中包含文件夹:bin、config、doc、lib四个文件夹。
其中mcl软件下载地址为:http://www.micans.org/mcl/src/mcl-latest.tar.gz;下载后使用:’./configure && make && make install’安装即可。
OrthoMCL软件路径(示例)为:~/orthomclSoftware-v2.0.9,将命令export PATH=~/orthomclSoftwarev2.0.9/bin/:$PATH写入~/.bashrc,之后输入命令source ~/.bashrc。


2.2 配置OrthoMCL软件
进入~/orthomclSoftware-v2.0.9路径下,输入:cp~/orthomclSoftwarev2.0.9/doc/OrthoMCLEngine/Main/orthomcl.config.template ~/example,将~/orthomcl.config.template拷贝到工作目录(以~/example为例)中,该文件为OrthoMCL的配置文件,以使用mysql数据库为例,其中的内容如下:
dbVendor=mysql   #使用的数据库为mysql
dbConnectString=dbi:mysql:orthomcl   #使用之前建立的名为orthomcl的数据库
dbLogin=user    #创建的用户名
dbPassword=123  #密码
similarSequencesTable=SimilarSequences #
orthologTable=Ortholog
inParalogTable=InParalog
coOrthologTable=CoOrtholog
interTaxonMatchView=InterTaxonMatch
percentMatchCutoff=50
evalueExponentCutoff=-5
oracleIndexTblSpc=NONE


2.3 输入文件
orthomcl的输入文件为fasta格式的基因或蛋白序列,fasta文件的序列名称要求以样品名开头之后接’|’分隔,之后接每个样品的序列名(如图1),而且样品名和序列名不能有重复。本文以样品A,B和参考序列ref为例,序列文件名分别为:A.fa,b.fa,ref.fa。

图1 fasta格式


2.4 处理输入文件

创建输入文件夹命名为seq,进入到seq文件夹输入:cd sed,之后将所有样品文件存放到一起命名为all.fa,命令为cat a.fa b.fa ref.fa > all.fa。
为保证序列格式的正确使用orthomclAdjustFasta程序,将fasta文件转换出兼容orthomcl的fasta文件使用命令:orthomclAdjustFasta tem all.fa 1,结果输出为tem.fasta。
之后对生成的tem.fasta进行过滤,使用命令:orthomclFilterFasta seq/ 10 20。允许的最短的protein长度是10,stop codons最大比例为20%;生成了两个文件goodProteins.fasta和poorProteins.fasta两个文件。


2.5 全序列比对
本文中已经包含参考序列所以可以直接将上一部的goodProteins.fasta序列进行自身的多序列比对。若为输入参考序列可将样品序列与OrthoMCL官网提供的参考序列存放到同一文件再进行多序列比对,多序列比对使用软件为blast+,输出结果为all.m8.anno。


2.6 导入比对结果
将比对结果导入mysql数据库,包含以下几个步骤:
A. 将orthomcl.blastout中的注释行去掉,生成新的文件blastresult,否则下一个命令会出错,命令为grep -P "^[^#]" all.m8.anno > all.m8
B. 将比对结果转化为规定格式,命名为similarSequences.txt,命令为:orthomclBlastParser all.m8 seq > similarSequences.txt
C. 由于evalue值为0时不识别故将0改为1e-181,命令为:perl -p -i -e 's/0\t0/1\t-181/' similarSequences.txt
D. 将similarSequences.txt导入到数据库中,命令为:orthomclLoadBlast orthomcl.config.template similarSequences.txt


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


2.8 将数据从mysql导出
自动生成了文件mclInput和文件夹pairs;文件夹中包含3个文件coorthologs.txt,inparalogs.txt,orthologs.txt。命令为:orthomclDumpPairsFiles orthomcl.config.template。


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


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

 
锐翌原创文章,未经授权严禁转载。
 图2 输出格式