基因家族生信分析

基因家族生信分析
一、什么是基因家族
概念：是来源于同一个祖先，有一个基因通过基因重复而产生两个或更多的拷
贝而构成的一组基因，他们在结构和功能上具有明显的相似性，编码相似的蛋白质产物。

划分：
按功能划分：把一些功能类似的基因聚类，形成一个家族。

按照序列相似程度划分:一般将同源的基因放在一起认为是一个家族。

1.常见基因家族：
WRKY基因家族：是植物前十大蛋白质基因家族之一，大量研究表明，WRKY 基因家族的许多成员参与调控植物的生长发育，形态建成与抗病虫。

NBS-LRR抗病基因家族：是植物中最大类抗病基因家族之一。

MADS-BOX基因家族：是植物体内的重要转录因子，它们广泛地调控着植物的生长、发育和生殖等过程。

在植物中参与花器官的发育，开花时间的调节，在果实，根，茎，叶的发育中都起着重要的作用。

热激蛋白70家族（HSP70）是一类在植物中高度保守的分子伴侣蛋白，在细胞中协助蛋白质正确折叠。

二、基因家族分析流程：
●利用蛋白保守域结构提取号在Pfam数据库提取其隐马尔科夫模型矩
阵文件（*.hmm）
●在数据库（Ensemble 、JGI、NVBI）下载你所需要的物种的基因组数
据（*.fa,*.gff）
●在虚拟机中Bio-Linux中的hummsearch程序，用隐马尔科夫模型矩
阵文件在蛋白序列文件中搜索含有该保守结构域的蛋白
●将蛋白序列导入MEGA软件构建进化树（可以阐明成员之间系统进化
关系，从进化关系上揭示其多样性）
●利用MEME搜索蛋白质的保守结构域
利用MEME搜索基因家族成员的motif可以揭示基因家族在物种内的多样化及其功能，如果他们都含有相同的motif表明其功能具有
相似性，如果部分家族成员含有其他不同的motif,很可能这些成员有
其他特异功能，或者可以归分为一个亚族
●绘制基因染色体位置图
从*.gff文件中抽取我们搜索到的基因位置信息，_v2.0/在线绘
制基因染色体位置图
通过染色体位置分布，可以了解基因主要分布字哪条染色体上，及是
否能形成基因簇（被认为是通过重组与错配促进基因交流）
●基因结构分析
从gff文件中抽取基因的结构信息，绘制转录本结构图。

●计算串联重复基因的Ka,Ks
1.首先将筛选到的基因的cds序列进行多序列对比，筛选
identity > 75%，tength大于对比的两条序列中较长的那条的长
度的75%，将筛选到的基因分别用clustalw进行比对，比对结果
导入KsKs_Calculster计算Ka，Ks、
Ka/ks比，计算核苷酸的非同义替代（ka）与核苷酸的同义替代
（ks）的平均速率。

2.Ka/ks比值<1表明：通过纯化选择降低了氨基酸变化的速率；比
值=1表示中性选择；比值>1,表明这些基因可能已经收到积极选
择，有利于适应性遗传，这些受正向选择的基因将作为以后的研
究重点。

软件的安装
从图片中获得进入NCBI-blast官网复制blast-linux版本的链接
在Linux终端
1.blast的安装
#wget blast链接
#tar xvfz 文件名解压缩文件
# mv 解压缩文件/root/local/app
# mv 解压缩文件blast
# vi .bashrc
#在最后一行添加export $PATH=/root/local/app/blast/bin:$PATH 并保存退出#source .bashrc 运行
#blastp -version 查看是否安装成功。

2．hummer的安装
#yum install -y wget //安装wget
#wget hmmer源码链接
#tar -zxvf hmmer-3.2.1
# vi .bashrc
#(在最末端添加的语句) PATH=$PATH:~/biosoft/hmmer-
#yum install -y gcc
#./configure
#make
#make check
#make install
#which hmmsearch 查看是否安装成功。

3.perl的安装
#wget 源代码链接
# tar xvfz perl- 解压缩
#cd perl-5.28.1
#./configure
#make
#make install 安装完成。

3.bioperl 的安装
#wget -O - | bash
#perlbrew install-cpanm
#/root/perl5/perlbrew/bin/cpanm Bio::Perl
三、具体操作：
1.保守域结构分析
下载蛋白保守结构域文件、cds、cDNA、gff注释文件和隐马尔科夫矩阵模型。

以拟南芥为例：
下载完成后，需要将文件传到Linux系统上进行分析：
打开虚拟机输入ip a将虚拟机IP连接到Xshell上，在Xshell上进行操作，将文件通过xftp（同样需要连接IP）传到Linux系统上，然后进行解压。

(一个命令解压多个文件:gunzip*.gz)
接下来用hummsearch寻找含有该蛋白保守结构域的蛋白及核酸序列
安装hummsearch
yum install -y wget //安装wget
#wget hmmer源码链接
#tar -zxvf hmmer-3.2.1
# vi .bashrc
#(在最末端添加的语句) PATH=$PATH:~/biosoft/hmmer-
#yum install -y gcc
#./configure
#make
#make check
#make install
#which hmmsearch 查看是否安装成功。

解压文件
移动到APP目录下面
在app目录下面新建文件夹 mkdir hmmer
将hmmer- mmove -v c:/hmmer-3.2.1 c:/hmmer
删除安装包
打开文字编辑器
vi ~/.bashrc
在文字编辑器里最后一行添加以上内容
安装好wget
如果make check出现错误百度用以下方法解决
出现complete安装完成
#source ~/.bashrc
#which hmmsearch
至此hmmer安装完成。

虚拟机操作：
1.导入下载好的文件；
2.hmmsearch --cut_tc –domtblout NB-ARC.txt NB-ARC.hmm Arabidopsis_thaliana.TAIR10
可以用editplus打开.txt文件
3. perl domain_xulie.pl 结果文件.txt 蛋白序列文件domain.fa 1e-20
4.clustalw进行多序列比对，得到aln文件和dnd文件。

5.hmmbuild 拟南芥特异的hmm模型文件domain.aln
6.hmmsearch –cut_tc –domtblout result.txt newhmm文件蛋白质序列文件
7.在Excel中，根据特定的evalue进行筛选，并对第一列进行去重复，得到第一列去重复的id，保存为id.txt
8. 用perl脚本根据id提取序列
Per get_fa_by_id.pl id.txt 蛋白序列wenjain >结果输出文件
可以根据的得到的序列文件进行后续的构建进化树、motif分析等。

2.搜索基因家族成员的MOTIF
2.1 需要准备的文件
1. 拟南芥NBS基因蛋白质序列
2. 蛋白保守结构域的隐马尔科夫模型矩阵文件
2.2 MOTIF的搜索
使用meme软件
命令：meme nbs_pep.fa -protein -oc nbs_motif -nostatus -maxsize 600000 -mor anr -nmotifs 10 -minw 6 -maxw 50
搜索结果存放在nbs_motif文件夹中。

文件夹中的eps文件可以用AI打开编辑，可以另存为png或jpg格式，也可打开网页版，也可用tbtools软件打开，下载motif在基因上的位置信息。

3. 绘制基因在染色体上的位置图
3.1 需要准备的文件
1. 拟南芥NBS基因id
2. 拟南芥的注释文件（gff3文件）
3. 拟南芥基因组长度
4.1在线绘图工具：MapGene2Chrom
4.2 samtools faidx 拟南芥. 可得到拟南芥. 该文件包括各个染色体，染色体长度。

4.3 对基因的id文件在Excel中进行分列，去重复处理。

4.4 使用处理过的id文件，对拟南芥的注释文件进行筛选
使用perl脚本得到基因在染色体上的位置。

命令：perl get_gene_gff.pl -in1 基因的id文件-in2 拟南芥gff3文件-out 新文件名称
4.5 新文件存放的是基因在染色体上的位置
4.6在在线文件MapGene2chrom 中，将基因在染色体上的位置信息文件复制到，input1框中，在input2中粘入samtools得到的fai文件。

4.绘制转录本的结构图
4.1 需要准备的文件
1. 拟南芥NBS基因转录本id（通过家族成员鉴定得到的蛋白id文件）
2.拟南芥基因的注释文件（gtf文件）
3.在线绘图工具：Gene Structure Display Server 2.0
http：//
4.2 具体方法
1. 准备gtf文件：输入命令：gffread gff3注释文件-T -o 输出文件（gtf文件）
2.editplus 打开gtf文件，去除”transcript:”
3. 使用perl脚本提取拟南芥转录本结构信息：命令：perl get_gtf.pl -in1 拟南芥转录本id文件-in2 gtf文件-out 输出文件（nbs_gtf.txt）
4.通过在线绘图工具，进行绘图。

5.筛选出串联重复基因
5.1准备文件
1.拟南芥NBS基因CDS序列
串联重复基因筛选标准【（a）length of alignable sequence covers>75% of longer
gene,and(b) similarity of aligned regions >75%】
参考文献：Extent of gene duplication in the genomes of Drosophila, nematode, and yeast.
2.由于筛选时产生的文件较多，因此创建新的目录：mkdir 新目录
3.用editplus 打开家族成员的id文件，对转录本id进行处理，使一个基因只拿一个转录本。

4.把id复制到Excel，首先排序处理，然后进行分列，然后以第一列删除重复值。

最后将第一列和第二列进行合并。

将处理好的id导入Linux。

5.使用perl脚本提取cds序列：命令：perl get_fa_by_id.pl id文件拟南芥cds 序列文件> cds.fa
6.使用blast软件筛选串联重复基因
6.1建立目标序列的数据库：makeblastdb -in cds.fa -dbtype nucl -title cds.fa 6.2 进行多序列比对：blastn -query cds.fa -db cds.fa -evalue 1e-20 -outfmt 9 -out result.txt
6.3 用editplus打开
6.4 得到cds序列的长度，使用samtools工具建立索引：命令：samtools faidx cds.fa
6.5 用perl脚本对result.txt进行筛选，perl KAKS_SHAIXUAN.pl -in1 -in2 result.txt -out cleanresult.txt
6.6 用editplus打开，将内容复制到Excel，在id后插入一列用公式：if （A1>B1,A1&B1,B1&A1）。

然后全选，以第C列删除重复值。

并保存到新的文件中，并导入到Linux中
7. 计算串联重复基因的KaKs。

7.1准备文件
1.串联重复基因的CDS序列文件
7.2 将成对的串联重复序列保存在一个文件中。

、
方法一：复制需要找到的序列的id，在editplus中按ctrl+F搜索，找到后复制粘贴到一个文件中。

方法二：首先将成对的id保存在同一个文件，导入到Linux中，在Linux中，利
用perl脚本提取序列：perl get_fa_by_id.pl 新的id文件cds.fa文件>id1.fa 7.3 计算KaKS
1.计算之前需要使用CLUSTAW对序列进行比较。

可获得id1.aln。

2.使用KaKs_calculator工具将id1.aln文件转换成id1.axt文件
命令：axtvenvertor id1.aln id1.axt
3.计算KaKs，输入命令：KaKs_calculstor -i id.axt -o id1_kaks.txt
4.如果报错，则把两条序列长度保持一致。

依此将所有的串联重复基因对，进行计算。

四．基因家族成员的鉴定（未知隐马尔科夫模型）
1.鉴定测略
•在NCBI数据库中尽量多下载几个物种的需要鉴定的蛋白保守结构域序列，以及所要研究物种的所有蛋白序列
•在虚拟机中本地建库，并进行blast
建库命令：makeblastdb -in 研究物种的蛋白序列文件-dbtype prot -title 库名称•进行序列比对
命令：blastp -query 下载的多个物种序列文件-db 库名称-evalue 1e-10 -outfmt 6 -out 结构域.blast
•使用sed 命令去除表头和结尾得到新的new结构域. blast
•Awk ‘{print$1}’ new结构域. Blast |less 可查看打印的结果
• Awk ‘{print$1}’ new结构域. Blast > id.txt
•cat id.txt|sort|uniq > idd.txt去重复
•Perl get_fa_by_id.pl idd.txt去重复蛋白序列wenjain >结果输出文件
•在Pfam或者NCBI的cdd中搜索检查是否有相关蛋白结构域。

之后再进行motif分析
一些命令及软件应用说明
（参考一些视频资料）
hmmsearch使用说明
用途：利用蛋白保守结构域的隐马尔科夫模型搜索蛋白序列中具有该保守结构域的蛋白
用法：hmmsearch --cut_tc --domtblout result.txt *.hmm
说明：result.txt是输出的结果文件，*.hmm在pfam数据库下载的模型，
hmmbuild使用说明
用途：利用clustalw比对生成的aln文件构建蛋白保守结构域的隐马尔科夫模型
用法：hmmbuild new.hmm domain.aln
说明：new.hmm 是结果文件也就是构建的蛋白保守结构域的隐马尔科夫模型，domain.aln是clustalw比对生成的aln文件
domain_xulie.pl脚本使用说明
用途：提取hmmsearch搜索结果中蛋白序列中保守结构域的序列，用于构建新的物种特异的蛋白保守结构域的隐马尔科夫模型
用法：perl domain_xulie.pl（脚本不在使用目录下要写全路径）hmmoutfile domain.fasta E-value
说明：hmmoutfile 是hmmsearch搜索结果文件
domain.fasta是结果存放文件也就是蛋白序列中保守结构域的序列，E-value 是提取序列时设定的E值
get_fa_by_id.pl使用说明
用途：通过ID号获取其相应的基因或蛋白序列
用法：perl perl get_fa_by_id.pl id.txt cds.fastat >id_cds.fasta
说明：id.txt是包含你的ID的文件，cds.fasta是你丛数据库中下载的包含所有cds序列的文件，id_cds.fasta是输出文件内容是ID对应的序列
samtools faidx
用途：提取fasta文件信息
用法：samtools faidx *.fa
说明：输入文件是fasta文件，自动生成输出目录*.fa.fai,结果的fai文件第一列是你输入的fasta文件的ID第二列是其序列长度
Gffread 使用说明
用途：将基因组注释文件gff3转化成基因的注释文件gtf
用法：gffread my.gff3 -T -o my.gtf
说明：my.gff3是输入文件基因组注释文件，my.gtf是输出文件是基因的注释文件
Get_gene_gff.pl
用途：想要绘制基因的染色体位置图必须要拿到基因在染色体上的具体信息，该脚本就是从总的gff文件抽去你需要的基因的信息如：所在染色体，起始终止位置等信息
用法：Perl Get_gene_gff.pl -in1 gene_id.txt -in2 my.gff3 -out gene_location.txt 说明：gene_id.txt 是第一个输入文件基因的ID文件，my.gff3是第二个输入文件是物种基因组所有蛋白序列，gene_location.txt 是结果输出文件
Get_gtf.pl脚本使用说明
用途：从基因注释文件gtf文件中提取转录本的结构信息
用法：perl get_gtf.pl -in1 id.txt -in2 gene.gtf -out structure.txt
说明：id.txt 是第一个输入文件是id文件，gene.gtf是第二个输入文件是基因注释文件gtf文件structure.txt是结果输出文件存放着转录本的结构信息
KaKs_shaixuan.pl 使用说明
用途：在多序列比对结果文件中筛选identity大于75%，比对上的序列长度大于对比的两条序列中最长序列的长度的75%
用法：perl KaKs_shaixuan.pl -in1 cds.fai -in2 result.txt -out shaixuan.txt
说明：cds.fai 是samtools faidx 对cds的序列文件fasta作用，生成的文件，result.txt 是拿cds进行多序列比对得到的结果，shaixuan.txt 是筛选后的结果存放的文本
aln文件转化axt文件
命令：
/home/manager/share/KaKs_Calculator2.0/KaKs_Calculator2.0/src/AXTConv ertor tong.aln tong.axt
说明：红色部分是KaKs_Calculator软件存放的位置（可变部分，放在那个路径就写全路径），蓝色部分是AXTConvertor工具在KaKs_Calculator的路径是不变的，tong.aln 文件是输入文件是序列比对生成的文件，tong.axt是结果文件，用来计算KAKS。

（注：KaKs_Calculator软件存放在参考资料部分，是压缩文件，解压即可使用，不需要安装）
计算KAKS
命令：
/home/manager/share/KaKs_Calculator2.0/KaKs_Calculator2.0/bin/Linux/Ka Ks_Calculator -i tong.axt -o kaks1.txt
说明：tong.axt 文件是输入文件，kaks1.txt是结果文件，存放同源基因的kaks circos.pl 脚本使用说明
用途：生成绘制圈图的配置文件
用法：perl circos.pl --chr Arabidopsis_
--circle link.txt --type link --circle text.txt --type text --od circose1
说明：Arabidopsis_，link.txt 存放同源基因对的位置信息，text.txt 存放同源基因的位置信息，od 是结果输出的文件夹
绘制circos图
命令：/biosoft/circos/circos-0.69-6/bin/circos -conf circos.conf
说明：circos.conf是circos.pl 脚本结果目录下的conf文件，
评语：。