比较基因组分析，首先需要选定要比较的物种，然后下载这些物种的基因组数据。
基因组数据准备
- 数据准备
1.基因结构注释文件， gff/gtf格式
2.蛋白序列文件，fasta格式
3.cds序列文件，fasta格式 - 数据来源
1.自己组装的基因组
2.数据库下载 - 处理原则
1.一个基因保留一个转录本
2.ID要一致
数据库
- GI-Phytozome下载植物基因组
https://phytozome-next.jgi.doe.gov/
提供去除可变剪切的cds和pep
Gff/gtf文件需要过滤
- Ensembl 基因组数据库
植物： http://plants.ensembl.org/index.html， 79 个记录
后生动物：http://metazoa.ensembl.org/index.htm ，112个记录
原生生物：http://protists.ensembl.org/index.html ，237个记录
真菌：http://fungi.ensembl.org/index.html ， 1,014 个记录
细菌： http://bacteria.ensembl.org/index.html, 44,048 个记录
蛋白序列和cds序列均包含所有转录本，需要过滤
gff文件包含所有转录本，需要过滤
NCBI下载基因组
https://www.ncbi.nlm.nih.gov/genome/
有些物种只有基因组序列，没有上传基因注释结果
基因组的染色体名称为NCBI的accession编号
蛋白序列和cds序列均包含所有转录本，需要过滤
蛋白序列ID和基因ID无法直接对应，需要借助gff文件物种特异数据库
拟南芥： https://www.arabidopsis.org/
苦荞： http://www.mbkbase.org/Pinku1/
水稻： http://rice.plantbiology.msu.edu/
...
由于可变剪切的存在，一条基因仅保留一条最长的转录本。蛋白、CDS序列ID需要和注释文件一致。
提取最长转录本
Phytozome 提供最长转录本数据，可以直接下载。当使用Ensembl 下载的数据时，可以参考以下脚本提取最长cds。
# 基因注释文件 Arabidopsis_thaliana.TAIR10.47.gff3
# cds序列文件 Arabidopsis_thaliana.TAIR10.cds.all.fa
# 蛋白序列文件 Arabidopsis_thaliana.TAIR10.pep.all.fa
# 基因组序列文件 Arabidopsis_thaliana.TAIR10.dna.toplevel.fa
## 去除gff3文件中ID部分多余字符
cp Arabidopsis_thaliana.TAIR10.47.gff3 Ath.gff3
sed -i 's/=gene:/=/g' Ath.gff3
sed -i 's/=transcript:/=/g' Ath.gff3
sed -i 's/=CDS:/=/g' Ath.gff3
## 基于gff3提取最长cds序列ID，并过滤gff3文件
perl ./gff_longest.pl Ath.gff3 Ath_id Ath_longest.gff3
## 提取最长的cds ID
awk '{print $2}' Ath_id > Ath_longest_id
## 基于最长的cds提取cds和蛋白质序列文件
seqtk subseq Arabidopsis_thaliana.TAIR10.cds.all.fa Ath_longest_id > Ath_longest.cds.fasta
seqtk subseq Arabidopsis_thaliana.TAIR10.pep.all.fa Ath_longest_id > Ath_longest.pep.fasta
## 如果没有cds和蛋白文件，也可以基于过滤后gff3文件从genome里提取cds并翻译成蛋白
gffread Ath_longest.gff3 \
-g Arabidopsis_thaliana.TAIR10.dna.toplevel.fa \ #基因组序列
-x Ath_longest.cds.fasta \ #输出cds序列
-y Ath_longest.pep.fasta #输出蛋白序列
