基因家族分析（1）数据准备及软件安装
基因组数据的准备
进行全基因水平的基因家族鉴定之前，需要准备好一套基因组数据，基因组数据可以从公共数据库下载，也可以根据基因组文献提供的地址到指定网站进行下载。一套完整的数据至少包括如下内容：
- 基因组序列文件，fasta 格式
- 基因结构注释文件，gff3 格式
- 所有蛋白质序列文件，fasta 格式
- 所有 cds 序列文件，fasta 格式
由于不同基因组数据库存储的数据格式及命名有各自规范，很多数据下载后不能直接用来做基因家族分析，需要提前处理好，处理原则如下：
1.去除所有可变剪切，一个编码基因保留一个转录本
2.去除 gff3 文件中非编码基因及重复序列等信息下面介绍 JGI、Ensembl 和 NCBI 下载的数据的处理
JGI/phytozome 数据处理
phytozome 是一个收录植物基因组数据的网站，数据整理比较规范，已经提供了去除可变剪切的 cds 和 protein 序列文件。只有 gff3 文件需要过滤处理。
示例数据为拟南芥数据，下载于phytozome13
Athaliana_167_TAIR10.cds_primaryTranscriptOnly.fa #cds序列
Athaliana_167_TAIR10.gene_exons.gff3 #基因结构文件
Athaliana_167_TAIR10.protein_primaryTranscriptOnly.fa #蛋白文件
Athaliana_167_TAIR9.fa #基因组文件
## 提取最长转录本基因ID
awk '$1 ~ /^>/ {print $1}' Athaliana_167_TAIR10.cds_primaryTranscriptOnly.fa | sed 's/^>//' > Ath_mRNA.id
## gff3文件相对于cds和蛋白序列文件，ID部分多了.TAIR10的字符串需要去除掉
sed 's/\.TAIR10//g' Athaliana_167_TAIR10.gene_exons.gff3 > Ath_1.gff3
## 基于mRNA id对gff3文件进行过滤
perl gff_filter_bymRNAID.pl Ath_1.gff3 Ath_mRNA.id geneID_mrnaID.table Ath_final.gff3
## 重命名蛋白序列和cds序列文件及基因组名称，方便后续使用
mv Athaliana_167_TAIR10.cds_primaryTranscriptOnly.fa Ath.cds.fa
mv Athaliana_167_TAIR10.protein_primaryTranscriptOnly.fa Ath.pep.fa
mv Athaliana_167_TAIR9.fa Ath.genome.fa
Ensembl 数据处理
Ensembl 数据库可以下载动物、植物、真菌、细菌等物种基因组数据。数据格式规范，mRNA ID 和 cds ID 基本一致，但没有去除可变剪切的版本，需要自己进行手动处理
示例数据为拟南芥数据，下载自 ensembl.
Arabidopsis_thaliana.TAIR10.47.gff3 # 基因结构文件
Arabidopsis_thaliana.TAIR10.dna.toplevel.fa # 基因组序列文件
Arabidopsis_thaliana.TAIR10.cds.all.fa # cds序列文件
Arabidopsis_thaliana.TAIR10.pep.all.fa # 蛋白序列文件
# 去除gff3文件中ID部分多余字符
cp Arabidopsis_thaliana.TAIR10.47.gff3 Ath.gff3.tmp #复制一份
sed -i 's/=gene:/=/g' Ath.gff3.tmp
sed -i 's/=transcript:/=/g' Ath.gff3.tmp
sed -i 's/=CDS:/=/g' Ath.gff3.tmp
# 基于gff3提取最长cds序列ID，并过滤gff3文件
perl gff_filter_longest.pl Ath.gff3.tmp Ath_gene_mrna_cds.ids Ath_final.gff3
# 提取最长cds ID列表
awk '{print $3}' Ath_gene_mrna_cds.ids > Ath_mRNA.id
##基于最长cds ID信息提取cds和蛋白质序列文件
seqtk subseq Arabidopsis_thaliana.TAIR10.cds.all.fa Ath_mRNA.id > Ath.cds.fasta
seqtk subseq Arabidopsis_thaliana.TAIR10.pep.all.fa Ath_mRNA.id > Ath.pep.fasta
# 基因组文件重命名
mv Arabidopsis_thaliana.TAIR10.dna.toplevel.fa Ath.genome.fasta
NCBI及其它数据库的处理有需要的请私信我
没有cds和蛋白序列的情况
如果没有 cds 和蛋白序列，可以基于 gff 和基因组序列文件使用 gffread进行提取.
gffread Ath_final.gff3 -g Ath.genome.fasta -x Ath.cds.fasta #提取cds序列
gffread Ath_final.gff3 -g Ath.genome.fasta -y Ath.pep.fasta #提取蛋白序列
软件安装
conda安装
用conda安装比对、结构域预测、motif鉴定、进化树构建、多序列比对结果过滤、fasta序列处理工具等等
blast
hmmer
meme
fasttree
trimal
seqkit
gffread
McscanX
JCVI
R包的安装
Peptides #蛋白质等电点和分子量的统计
seqlogo #绘制seqlogo图
pheatmap #绘制热图
msa #多序列比对的R包
windows软件
染色体核型图mapchart
进化树构建 mega
在线软件
进化树美化 https://itol.embl.de/
motif 预测meme MEME - Submission form (meme-suite.org)
基因结构绘制 Gene Structure Display Server 2.0 (gao-lab.org)
顺式作用元件预测 PlantCARE, a database of plant promoters and their cis-acting regulatory elements (ugent.be)
