首先是安裝這個流程往踢,試了一下可以使用conda進行安裝
conda search evidencemodeler
conda install evidencemodeler
安裝好以后很多perl腳本是在 anaconda3/envs/EVM/opt/evidencemodeler-2.1.0/EvmUtils/
這個目錄下
學習這個流程的參考鏈接
1先匪、https://www.zhouxiaozhao.cn/2020/11/26/2020-11-26-annotion(4)/
基因結構注釋(4):整合預測結果2、使用MAKER進行注釋: 訓練SNAP基因模型 http://www.reibang.com/p/8c378421af12
3、基因組結構注釋 http://www.reibang.com/p/2cfc7638663d 這個寫的很詳細
我這里的數(shù)據(jù)就使用擬南芥的一條染色體岸晦,這個數(shù)據(jù)來源于論文
Chromosome-level assemblies of multiple Arabidopsis genomes reveal hotspots of rearrangements with altered evolutionary dynamics
C24.chr.all.v2.0.fasta
的一號染色體
evm 流程第一步用到的命令是
time ~/anaconda3/envs/EVM/opt/evidencemodeler-2.1.0/EvmUtils/evidence_modeler.pl \
--genome ../../repeat/chr1.fa.masked \
--weights /data/myan/raw_data/practice/pan.genome/at.nc/version2.5.2019-10-09/genome.annotation/braker2/evm/weights.txt \
--gene_predictions evm_abinitio.gff3 \
--protein_alignments evm_pro.gff3 \
--transcript_alignments transcripts.fasta.transdecoder.genome.gff3 > evm.out
這里需要的輸入數(shù)據(jù)
染色體序列
weights.txt
這里的內(nèi)容
PROTEIN genomeThreader 4
TRANSCRIPT transdecoder 8
ABINITIO_PREDICTION Augustus 1
ABINITIO_PREDICTION GeneMark.hmm 1
- evm_abinitio.gff3
這個是從頭預測的結果用到的是augustus和genemark 直颅,我這里直接用braker2這個流程先跑一遍,就可以得到這兩個結果文件 將這兩個結果文件合并成了evm_abinitio.gff3
evm_pro.gff3 這個是基于同源蛋白的結果使用的是gth這個程序
transcripts.fasta.transdecoder.genome.gff3 這個是轉錄組的數(shù)據(jù)
接下來介紹如何得到這些輸入文件
首先是對重復序列進行屏蔽
RepeatMasker -species Arabidopsis -pa 12 -xsmall -dir repeat chr1.fa
這里擬南芥可以直接指定物種空猜,如果自己研究的物種不是很常見的,需要構建這個重復序列的庫
基于同源蛋白的注釋
evm_pro.gff3 這個相對簡單,這里需要有一個同源蛋白文件和需要注釋的染色體基因組恨旱,得到gff文件后需要用evm這個流程里的腳本對格式進行轉換辈毯,gth這個軟件的安裝,如果braker2這個軟件安裝好是可以直接用的
~/anaconda3/envs/braker2/bin/gth -gff3out -intermediate -duplicatecheck seq -protein proteins.fa -translationtable 1 -genomic chr1.fa.mask -o homo_protein.gff3 -force
~/anaconda3/envs/EVM/opt/evidencemodeler-2.1.0/EvmUtils/misc/genomeThreader_to_evm_gff3.pl homo_protein.gff3 > evm_pro.gff3
基于轉錄組數(shù)據(jù)
transcripts.fasta.transdecoder.genome.gff3 這個文件用到轉錄組數(shù)據(jù)
hisat2-build ../repeat/chr1.fa.masked chr1 -p 8
hisat2 -p 8 --dta -x chr1 -1 SRR4420296_R1.fastq.gz -2 SRR4420296_R2.fastq.gz -S SRR4420296.sam
samtools sort -@ 8 -O bam -o SRR4420296.bam SRR4420296.sam
stringtie -p 12 -o SRR4420296.gtf SRR4420296.bam
gtf_genome_to_cdna_fasta.pl SRR4420296.gtf ../repeat/chr1.fa.masked > transcripts.fasta
gtf_to_alignment_gff3.pl SRR4420296.gtf > transcripts.gff3
TransDecoder.LongOrfs -t transcripts.fasta
TransDecoder.Predict -t transcripts.fasta
cdna_alignment_orf_to_genome_orf.pl transcripts.fasta.transdecoder.gff3 transcripts.gff3 transcripts.fasta > transcripts.fasta.transdecoder.genome.gff3
從頭預測
braker.pl --cores 48 --species=At01 --genome=../repeat/chr1.fa.masked --softmasking --bam=SRR4420296.bam --gff3 --prot_seq=../proteins.fa --prg=gth
這個從頭預測的時候也可以指定蛋白證據(jù)搜贤,但不知道為什么這個生成的結果文件里沒有單獨的gth結果谆沃,不知道是不是需要單獨指定某個參數(shù)
這一步會生成 augustus.hints.gtf 和 GeneMark-ET/genemark.gtf
也需要做一個格式轉換
~/anaconda3/envs/EVM/opt/evidencemodeler-2.1.0/EvmUtils/misc/augustus_GTF_to_EVM_GFF3.pl augustus.hints.gtf > evm_augustus.gff3
~/anaconda3/envs/EVM/opt/evidencemodeler-2.1.0/EvmUtils/misc/GeneMarkHMM_GTF_to_EVM_GFF3.pl GeneMark-ET/genemark.gtf > evm_genemark.gff3
cat evm_augustus.gff3 evm_genemark.gff3 > evm_abinitio.gff3
有了這些結果文件就可以運行開頭提到的命令
time ~/anaconda3/envs/EVM/opt/evidencemodeler-2.1.0/EvmUtils/evidence_modeler.pl \
--genome ../../repeat/chr1.fa.masked \
--weights /data/myan/raw_data/practice/pan.genome/at.nc/version2.5.2019-10-09/genome.annotation/braker2/evm/weights.txt \
--gene_predictions evm_abinitio.gff3 \
--protein_alignments evm_pro.gff3 \
--transcript_alignments transcripts.fasta.transdecoder.genome.gff3 > evm.out
得到evm.out后怎么處理后面再來更新,因為程序還沒有跑完仪芒,還沒有拿到這個結果
這一步是可以并行的唁影,這個并行怎么用還需要仔細研究一下
推文記錄的是自己的學習筆記,內(nèi)容可能會存在錯誤掂名,請大家批判著看夭咬,歡迎大家指出其中的錯誤
示例數(shù)據(jù)和代碼可以給推文點贊,然后點擊在看铆隘,最后留言獲取
歡迎大家關注我的公眾號
小明的數(shù)據(jù)分析筆記本
小明的數(shù)據(jù)分析筆記本 公眾號 主要分享:1卓舵、R語言和python做數(shù)據(jù)分析和數(shù)據(jù)可視化的簡單小例子;2膀钠、園藝植物相關轉錄組學掏湾、基因組學裹虫、群體遺傳學文獻閱讀筆記;3融击、生物信息學入門學習資料及自己的學習筆記筑公!