正在更新。英文是真的看不习惯。
参考:使用MAKER进行基因注释(高级篇之AUGUSTUS模型训练) - 简书 (jianshu.com)
Training AUGUSTUS (wisc.edu)
Augustus 进行基因注释 - 斩毛毛 - 博客园 (cnblogs.com)
首先是我曾经理解上的错误:
1.“若存在已经被训练的物种(augustus --species=help查看,或者augustus的config/species路径下有文件夹列表),则直接使用一下代码进行预测基因,”
说实话,我一开始理解错了。我以为是使用其他近缘物种(比如我的近缘物种是rice)也可以对我自己新测的物种参考。但是他们想说的是,如果你测序的是水稻,那你可以直接用水稻的augustus训练好的文件。
2.“b. 这些基因的基因结构一定要足够的准确。”
我一开始就更不理解了。我都找你来训练了,我怎么能确认我的基因结构是一定准确的?大概意思需要转录组RNA数据就可以?
总结训练集要求:基因之间不重复(只要一个基因的一个转录本),至少随机100-200个基因,要准确。
理论的事我就不多说了。我直接运行算了。
我是已经用genewise进行过了同源注释。用了四个单子叶的做近缘物种,genewise.gff的BUSCO还都挺高(水稻的低一些),如果没有RNA数据用同源注释结果应该也是可以的。
水稻genewise的BUSCO结果 同源注释的gffBUSCO结果
同源注释的gffBUSCO结果 同源注释的gffBUSCO结果
准备好gff之后,运行一个prepare脚本后,可以开始训练了。
#!bin/bash
genome=genome.fa #基因组文件
species=Hb21 #我自己设置的物种名
unset PERL5LIB; export PATH=/useful/perl-5.30.2/bin:$PATH
export PATH=/share/app/blat-319/blat:$PATH
source activate augustus
export AUGUSTUS_CONFIG_PATH=/01.bin/Augustus/augustus_config #公共软件无法写入新文件夹,可以设置自己的路径
#perl /bin/perfect_gene/perfect_gene.pl --sco 99 --start 0 --stop 1 $genome ../Oryza_sativa/Oryza_sativa.IRGSP-1.0.pep.all.fa.genewise.gff ../Brachypodium_distachyon/Brachypodium_distachyon.Brachypodium_distachyon_v3.0.pep.all.fa.genewise.gff ../Musa_acuminata/Musa_acuminata_v2.pep.all.fa.genewise.gff ../Sorghum_bicolor/Sorghum_bicolor.Sorghum_bicolor_NCBIv3.pep.all.fa.genewise.gff
#上面这个命令应该是用同源注释做训练集。
perl /autoAugTrain.pl --genome=$genome --species=$species --trainingset=genome.fa.gff.nr.gff --cpus 30
上面这个perfect_gene.pl脚本是别人写的,我不好分享。内容大概是:(待补充)
/autoAugTrain.pl这个好像在augustus的script文件夹下
最后会在你设置的AUGUSTUS_CONFIG_PATH生成这些文件:
#在AUGUSTUS_CONFIG_PATH/species/Hb21下
Hb21_weightmatrix.txt
Hb21_metapars.cfg
Hb21_metapars.utr.cfg
Hb21_metapars.cgp.cfg
Hb21_parameters.cfg.orig1
Hb21_parameters.cfg
Hb21_intron_probs.pbl
Hb21_exon_probs.pbl
Hb21_igenic_probs.pbl
Hb21_exon_probs.pbl.withoutCRF
Hb21_igenic_probs.pbl.withoutCRF
Hb21_intron_probs.pbl.withoutCRF
#运行的脚本路径autoAugTrain/training下
utr
training.gff
training.gb
training.gb.train
training.gb.test
training.gb.onlytrain
training.gb.train.test
train.err
train.out
optimize.out
tmp_opt_Hb21
train.withoutCRF.err
train.withoutCRF.out
test
#
就是你的augustus比标准数据库里多了你训练的物种。
然后我发现我只有训练集,没有测试集对这个进行测试。而且根据参考(avrilomics: Training the Augustus gene-finding software),他说默认会进行5次参数的优化,我这只有一次Hb21_parameters.cfg.orig1
。所以我想再运行几次,再生成一个测试集(使用训练集的测试结果会虚高)。
如果要用RNA结果训练augustus...RNA注释是有两种策略,,一种是使用HISAT2 + StringTie先比对再组装, 一种是从头组装,然后使用PASA将转录本比对到基因组上。(基因结构注释(3):转录组预测 - 简书 (jianshu.com))
RNA下机数据经过对基因组index,
samtools faidx genome.fasta
bwa index genome.fasta
java -jar picard/2.23.8/picard.jar CreateSequenceDictionary R=genome.fasta O=genome.fasta.dict
去接头
conda activate rna
java -jar /01.Software/Trimmomatic-0.39/trimmomatic-0.39.jar PE -threads 16 -phred33 RNA/*85_1.fq.gz RNA/*85_2.fq.gz /cleandata/*85_1.paired.fq.gz /cleandata/*85_1.unpaired.fq.gz /cleandata/*85_2.paired.fq.gz /cleandata/*85_2.unpaired.fq.gz ILLUMINACLIP:/01.Software/Miniconda/envs/rna/adapter/adapter.fa:2:35:4:12:true LEADING:3 TRAILING:3 SLIDINGWINDOW:5:15 MINLEN:50 2>trimming.log
#我的conda下环境rna安了很多处理RNA的软件。
#路径是简化的。trimmomatic运行命令可去查找别的教程。
#ILLUMINACLIP的fa文件和测序方联系获得。就是你接头的序列。trimmomatic软件的目录下也自带了一些接头。trimmomatic参数参考其他教程吧,不赘述。
之后运行
#!usr/bin/bash
hisat2-build -p 4 genome.fasta genome.fasta
hisat2 -p 4 --max-intronlen 500000 --sensitive --dta --dta-cufflinks --phred33 --no-discordant --no-mixed -x genome.fasta -1 cleandata/*85_1.paired.fq.gz -2 cleandata/*85_2.paired.fq.gz -S genome_rna.sam
samtools view -bF 4 -S genome_rna.sam -b -o genome_rna.bam
samtools sort genome_rna.bam -o genome_rna.bam.sort
stringtie genome_rna.bam.sort -p 1 -o genome_rna.transcript.gtf
xzg是用的PASA,虽然我用的stringtie,但是估计也要回到PASA,因为最后要用EVM整合。我问了chatGPT有关这两种方法的区别,回答如下(警告:AI自动生成结果,不一定对,注意甄别)
“PASA和StringTie是两个常用的生物信息学分析软件,用于转录组数据的分析和注释。它们有一些区别,包括其用途、功能和特点。
PASA(Program to Assemble Spliced Alignments)是一个用于转录组注释和组装的软件。它的主要功能是利用已知的基因组和转录组测序数据,通过比对和组装转录本来提供基因结构注释。PASA能够将已知的基因组注释与转录组测序数据进行比对,通过多轮迭代的方法,生成准确的转录本注释,并提供详细的转录本结构信息。
相比之下,StringTie是一个转录组装和表达定量分析工具。它可以将转录组测序数据进行拼接并组装成转录本,同时对每个转录本进行定量表达分析。StringTie具有高效的组装算法,能够有效地处理复杂的基因结构和低表达水平的转录本。它还能够生成基因和转录本的表达水平估计,提供丰富的定量信息。
以下是PASA和StringTie之间的一些区别:
用途:PASA主要用于转录组注释和组装,旨在提供准确的基因结构注释。StringTie主要用于转录组装和表达定量分析,旨在提供转录本的组装和表达量估计。
功能:PASA具有注释转录本结构的功能,能够比对和组装转录本,并生成准确的注释结果。StringTie具有高效的组装算法,能够处理复杂的基因结构,并提供表达量估计和可视化分析。
灵活性:PASA相对较为复杂,需要多轮迭代和大量的参数调整来生成准确的注释结果。StringTie具有简单的工作流程和默认参数,适用于一般的转录组装和表达量分析。
可扩展性:PASA可以与其他转录组装软件和注释工具进行集成,如BLAT、Exonerate和BLAST等。StringTie可以与其他转录组装和表达量分析软件进行集成,如Cufflinks和Ballgown等。
综上所述,PASA和StringTie在转录组分析中有不同的用途和功能。选择哪个软件取决于具体的研究目标和需求。如果需要准确的转录本注释和结构信息,PASA可能是一个更好的选择。如果主要关注转录组装和表达量分析,StringTie可能更适合。”
网友评论