EVidenceModeler(EVM)流程做基因组注释

2023-08-23 10:40:45 浏览数 (1)

首先是安装这个流程,试了一下可以使用conda进行安装

代码语言:javascript复制
conda search evidencemodeler

image.png

代码语言:javascript复制
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基因模型 https://www.jianshu.com/p/8c378421af12
  • 3、基因组结构注释 https://www.jianshu.com/p/2cfc7638663d 这个写的很详细

我这里的数据就使用拟南芥的一条染色体,这个数据来源于论文

Chromosome-level assemblies of multiple Arabidopsis genomes reveal hotspots of rearrangements with altered evolutionary dynamics

C24.chr.all.v2.0.fasta的一号染色体

evm 流程第一步用到的命令是

代码语言:javascript复制
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

这里需要的输入数据

  • 染色体序列
  • weights.txt

这里的内容

代码语言:javascript复制
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 这个是转录组的数据

接下来介绍如何得到这些输入文件

首先是对重复序列进行屏蔽

代码语言:javascript复制
RepeatMasker -species Arabidopsis -pa 12 -xsmall -dir repeat chr1.fa

这里拟南芥可以直接指定物种,如果自己研究的物种不是很常见的,需要构建这个重复序列的库

基于同源蛋白的注释

evm_pro.gff3 这个相对简单,这里需要有一个同源蛋白文件和需要注释的染色体基因组,得到gff文件后需要用evm这个流程里的脚本对格式进行转换,gth这个软件的安装,如果braker2这个软件安装好是可以直接用的

代码语言:javascript复制
~/anaconda3/envs/braker2/bin/gth -gff3out -intermediate -duplicatecheck seq -protein proteins.fa -translationtable 1 -genomic chr1.fa.masked -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

基于转录组数据

transcripts.fasta.transdecoder.genome.gff3 这个文件用到转录组数据

代码语言:javascript复制
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

从头预测

代码语言:javascript复制
braker.pl --cores 48 --species=At01 --genome=../repeat/chr1.fa.masked --softmasking --bam=SRR4420296.bam --gff3 --prot_seq=../proteins.fa --prg=gth

这个从头预测的时候也可以指定蛋白证据,但不知道为什么这个生成的结果文件里没有单独的gth结果,不知道是不是需要单独指定某个参数

这一步会生成 augustus.hints.gtf 和 GeneMark-ET/genemark.gtf

也需要做一个格式转换

代码语言:javascript复制
~/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

有了这些结果文件就可以运行开头提到的命令

代码语言:javascript复制
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后怎么处理后面再来更新,因为程序还没有跑完,还没有拿到这个结果

这一步是可以并行的,这个并行怎么用还需要仔细研究一下

推文记录的是自己的学习笔记,内容可能会存在错误,请大家批判着看,欢迎大家指出其中的错误

示例数据和代码可以给推文点赞,然后点击在看,最后留言获取

欢迎大家关注我的公众号

小明的数据分析笔记本

小明的数据分析笔记本 公众号 主要分享:1、R语言和python做数据分析和数据可视化的简单小例子;2、园艺植物相关转录组学、基因组学、群体遗传学文献阅读笔记;3、生物信息学入门学习资料及自己的学习笔记!

0 人点赞