基因組數(shù)據(jù)的準(zhǔn)備
進(jìn)行全基因水平的基因家族鑒定之前菇篡,需要準(zhǔn)備好一套基因組數(shù)據(jù)异希,基因組數(shù)據(jù)可以從公共數(shù)據(jù)庫下載袭景,也可以根據(jù)基因組文獻(xiàn)提供的地址到指定網(wǎng)站進(jìn)行下載狠角。一套完整的數(shù)據(jù)至少包括如下內(nèi)容:
- 基因組序列文件号杠,fasta 格式
- 基因結(jié)構(gòu)注釋文件,gff3 格式
- 所有蛋白質(zhì)序列文件丰歌,fasta 格式
- 所有 cds 序列文件姨蟋,fasta 格式
由于不同基因組數(shù)據(jù)庫存儲的數(shù)據(jù)格式及命名有各自規(guī)范,很多數(shù)據(jù)下載后不能直接用來做基因家族分析动遭,需要提前處理好芬探,處理原則如下:
1.去除所有可變剪切,一個(gè)編碼基因保留一個(gè)轉(zhuǎn)錄本
2.去除 gff3 文件中非編碼基因及重復(fù)序列等信息下面介紹 JGI厘惦、Ensembl 和 NCBI 下載的數(shù)據(jù)的處理
JGI/phytozome 數(shù)據(jù)處理
phytozome 是一個(gè)收錄植物基因組數(shù)據(jù)的網(wǎng)站偷仿,數(shù)據(jù)整理比較規(guī)范,已經(jīng)提供了去除可變剪切的 cds 和 protein 序列文件宵蕉。只有 gff3 文件需要過濾處理酝静。
示例數(shù)據(jù)為擬南芥數(shù)據(jù),下載于phytozome13
Athaliana_167_TAIR10.cds_primaryTranscriptOnly.fa #cds序列
Athaliana_167_TAIR10.gene_exons.gff3 #基因結(jié)構(gòu)文件
Athaliana_167_TAIR10.protein_primaryTranscriptOnly.fa #蛋白文件
Athaliana_167_TAIR9.fa #基因組文件
## 提取最長轉(zhuǎn)錄本基因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文件進(jìn)行過濾
perl gff_filter_bymRNAID.pl Ath_1.gff3 Ath_mRNA.id geneID_mrnaID.table Ath_final.gff3
## 重命名蛋白序列和cds序列文件及基因組名稱别智,方便后續(xù)使用
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 數(shù)據(jù)處理
Ensembl 數(shù)據(jù)庫可以下載動物、植物稼稿、真菌薄榛、細(xì)菌等物種基因組數(shù)據(jù)。數(shù)據(jù)格式規(guī)范让歼,mRNA ID 和 cds ID 基本一致敞恋,但沒有去除可變剪切的版本,需要自己進(jìn)行手動處理
示例數(shù)據(jù)為擬南芥數(shù)據(jù)谋右,下載自 ensembl.
Arabidopsis_thaliana.TAIR10.47.gff3 # 基因結(jié)構(gòu)文件
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 #復(fù)制一份
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和蛋白質(zhì)序列文件
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及其它數(shù)據(jù)庫的處理有需要的請私信我
沒有cds和蛋白序列的情況
如果沒有 cds 和蛋白序列,可以基于 gff 和基因組序列文件使用 gffread進(jìn)行提取.
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安裝比對、結(jié)構(gòu)域預(yù)測啸蜜、motif鑒定坑雅、進(jìn)化樹構(gòu)建、多序列比對結(jié)果過濾衬横、fasta序列處理工具等等
blast
hmmer
meme
fasttree
trimal
seqkit
gffread
McscanX
JCVI
R包的安裝
Peptides #蛋白質(zhì)等電點(diǎn)和分子量的統(tǒng)計(jì)
seqlogo #繪制seqlogo圖
pheatmap #繪制熱圖
msa #多序列比對的R包
windows軟件
染色體核型圖mapchart
進(jìn)化樹構(gòu)建 mega
在線軟件
進(jìn)化樹美化 https://itol.embl.de/
motif 預(yù)測meme MEME - Submission form (meme-suite.org)
基因結(jié)構(gòu)繪制 Gene Structure Display Server 2.0 (gao-lab.org)
順式作用元件預(yù)測 PlantCARE, a database of plant promoters and their cis-acting regulatory elements (ugent.be)