1. Augustus training
首选，需要有至少 200 个完整基因模型的数据。 例如： 使用 genome-guided 方法进行 trinity 有参考基因组的 de novo 组装；再使用 PASA 将组装出来的 inchworm 序列比对到基因组； 再提取完整基因模型数据，得到文件 trainingSet_CompleteBest.gff3 。
然后，将 GFF3 文件转换为 GeneBank 文件：
$ /opt/biosoft/augustus-3.0.3/scripts/gff2gbSmallDNA.pl trainingSet_CompleteBest.gff3 ../../genome.fasta 50 trainingSetComplete.gb
使用 genebank 格式的文件预先进行一次 Augustus training，一般会得到一些错误提示
$ etraining --species=generic --stopCodonExcludedFromCDS=false trainingSetComplete.gb 2> train.err
根据错误提示，提取出有错误的基因模型
$ cat train.err | perl -ne 'print "$1\n" if /in sequence (\S+):/' > badlist
去除错误的基因模型，得到能用于 training 的基因模型
$ /opt/biosoft/augustus-3.0.3/scripts/filterGenes.pl badlist trainingSetComplete.gb > training.gb
在进行 Augustus training 之前，最好保证这些基因模型两两之间在氨基酸水平的 identity < 70% ：
先获取这些基因模型的 Proteins 序列
$ /opt/biosoft/augustus-3.0.3/scripts/gbSmallDNA2gff.pl training.gb > training.gff2
$ perl -ne ‘print “$1\n” if /gene_id \”(\S+?)\”/’ training.gff2 | uniq > trainSet.lst
$ perl extract_genes.pl trainingSet_CompleteBest.gff3 trainSet.lst > training.gff3
$ /opt/biosoft/EVM_r2012-06-25/EvmUtils/gff3_file_to_proteins.pl training.gff3 genome.fasta prot > training.protein.fasta
$ perl -p -i -e ‘s/(>\S+).*/$1/’ training.protein.fasta
$ perl -p -i -e ‘s/\*//’ training.protein.fasta
对这些 proteins 序列构建 blast 数据库，并将 proteins 序列比对到此数据库
$ makeblastdb -in training.protein.fasta -dbtype prot -title training.protein.fasta -parse_seqids -out training.protein.fasta
$ blastp -db training.protein.fasta -query training.protein.fasta -out training.protein.fasta.out -evalue 1e-5 -outfmt 6 -num_threads 8
提取 blast 结果中 identity >= 70% 的比对信息，identity 高的 proteins 序列仅保留一条
$ grep -v -P “\t100.00\t” training.protein.fasta.out | perl -ne ‘split; print if $_[2] >= 70’ > blast.out
$ perl delete_high_identity_proteins_in_training_gff3.pl training.protein.fasta blast.out training.gff3 > training.new.gff3
获得去除了冗余的基因模型
$ /opt/biosoft/augustus-3.0.3/scripts/gff2gbSmallDNA.pl training.new.gff3 genome.fasta 50 training.gb
进行 Augustus Training without CRF
初始化本物种的 HMM 文件
$ /opt/biosoft/augustus-3.0.3/scripts/new_species.pl --species=my_species
如果有 RNA-Seq 数据，获得能用于 training 的基因模型一般会有好几千个。我们将基因模型随机分成两份： 第一份 300 个基因，用于检测 training 的精确性； 另外一份含有更多的基因，用于进行 Augustus training。
$ /opt/biosoft/augustus-3.0.3/scripts/randomSplit.pl training.gb 300
使用大份的 traing.gb.train 进行 Augustus Training
$ etraining --species=my_species training.gb.train >train.out
根据输出结果 train.out 来修正参数文件 species_parameters.cfg 中终止密码子的频率
$ tag=$(perl -ne 'print "$1" if /tag:\s+\d+\s+\((\S+)\)/' train.out)
$ perl -p -i -e "s#/Constant/amberprob.*#/Constant/amberprob $tag#" /opt/biosoft/augustus-3.0.3/config/species/lentinula_edodes/lentinula_edodes_parameters.cfg
$ taa=$(perl -ne 'print "$1" if /taa:\s+\d+\s+\((\S+)\)/' train.out)
$ perl -p -i -e "s#/Constant/ochreprob.*#/Constant/ochreprob $taa#" /opt/biosoft/augustus-3.0.3/config/species/lentinula_edodes/lentinula_edodes_parameters.cfg
$ tga=$(perl -ne 'print "$1" if /tga:\s+\d+\s+\((\S+)\)/' train.out)
$ perl -p -i -e "s#/Constant/opalprob.*#/Constant/opalprob $tga#" /opt/biosoft/augustus-3.0.3/config/species/lentinula_edodes/lentinula_edodes_parameters.cfg
根据 training 的结果，进行第一遍精确性评估
$ augustus --species=my_species training.gb.test > test.1.out
再将 training.gb.train 分成两份： 第一份 800 个基因，剩下基因为另外一份。这两份基因都用于 Optimizing meta parameters of AUGUSTUS
$ /opt/biosoft/augustus-3.0.3/scripts/randomSplit.pl training.gb.train 800
使用上面的 2 份基因模型，进行 Augustus training 的优化
$ /opt/biosoft/augustus-3.0.3/scripts/optimize_augustus.pl --species=my_species --rounds=5 --cpus=16 --kfold=16 training.gb.train.test --onlytrain=training.gb.onlytrain --metapars=/opt/biosoft/augustus-3.0.3/config/species/my_species/my_species_metapars.cfg > optimize.out
按如上参数，则程序会对 my_species_metapars.cfg 文件中的 28 个参数进行优化，总共优化 5 轮或有一轮找不到可优化的参数为止。每进行一个参数的优化时： 将 training.gb.train.test 文件中 800 个基因随机分成 16 等份，取其中 15 等份和 training.gb.onlytrain 中的基因模型一起进行 training，剩下的 1 等份用于精确行评估； 要对每个等份都进行一次精确性评估；使用 16 个 CPU 对 16 个等份并行进行 training 和 精确性评估，得到平均的精确值；优化的每个参数会分别 3~4 个值，每个值得到一个 training 的精确值，对参数的多个设定值进行比较，找到最佳的值。
此外， training 的精确值的算法： accuracy value = (3*nucleotide_sensitivity + 2*nucleotide_specificity + 4*exon_sensitivity + 3*exon_specificity + 2*gene_sensitivity + 1*gene_specificity) / 15 。
优化参数完毕后，需要再此使用 training.gb.train 进行一遍 training
$ etraining --species=my_species training.gb.train
再进行第二遍的精确性评估，一般进行优化后，精确性会有较大幅度的提高
$ augustus --species=my_species training.gb.test > test.2.withoutCRF.aout
进行 Augustus Training with CRF(Conditional Random Field)
在进行 Training with CRF 之前，先备份非 CRF 的参数文件
$ cd /opt/biosoft/augustus-3.0.3/config/species/species/
$ cp species_exon_probs.pbl lentinula_edodes_exon_probs.pbl.withoutCRF
$ cp species_igenic_probs.pbl lentinula_edodes_igenic_probs.pbl.withoutCRF
$ cp species_intron_probs.pbl lentinula_edodes_intron_probs.pbl.withoutCRF
$ cd -
在 training 的时候，加入 --CRF 参数，这样进行 training 比较耗时
$ etraining --species=my_species training.gb.train --CRF=1 1>train.CRF.out 2>train.CRF.err
再次进行精确行评估
$ augustus --species=my_species training.gb.test > test.2.CRF.out
比较 CRF 和 非 CRF 两种情况下的精确性。一般情况下，CRF training 的精确性要高些。若 CRF training 的精确行低些，则将备份的参数文件还原回去即可。
2. Training UTR parameters for Augustus
Training UTR 有利于利用 exonpart hint。 当使用 exonpart hint 时，对基因结构预测有显著提高。
首先，需要获得用于 training 的基因模型，这些基因模型需要同时带有 5’UTR 和 3’UTR 。
$ mkdir utr
$ cd utr
从用于 Augustus Training 的数据中提取同时带有 5'UTR 和 3'UTR 的基因模型。
$ perl -ne 'print "$2\t$1\n" if /.*\t(\S+UTR)\t.*ID=(\S+).utr/' ../training.new.gff3 | sort -u | perl -ne 'split; print "$_[0]\n" if ($g eq $_[0]); $g = $_[0];' > bothutr.lst
$ perl -e 'open IN, "bothutr.lst"; while () {chomp; $keep{$_}=1} $/="//\n"; while (<>) { if (/gene=\"(\S+?)\"/ && exists $keep{$1}) {print} }' ../training.gb.test > training.gb.test
$ perl -e 'open IN, "bothutr.lst"; while () {chomp; $keep{$_}=1} $/="//\n"; while (<>) { if (/gene=\"(\S+?)\"/ && exists $keep{$1}) {print} }' ../training.gb.train > training.gb.train
使用 traing.gb.train 中的基因模型进行 Training UTR parameters
$ /opt/biosoft/augustus-3.0.3/scripts/randomSplit.pl training.gb.train 400
$ mv training.gb.train.train training.gb.onlytrain
$ optimize_augustus.pl --species=my_species --cpus=16 --rounds=5 --kfold=16 --UTR=on --trainOnlyUtr=1 --onlytrain=training.gb.onlytrain --metapars=/opt/biosoft/augustus-3.0.3/config/species/lentinula_edodes/lentinula_edodes_metapars.utr.cfg training.gb.train.test > optimize.out
进行 without CRF 和 with CRF 的精确性比较，选取其中精确性较高的参数文件
$ etraining --species=species --UTR=on training.gb.train
$ augustus --species=species --UTR=on training.gb.test > test.withoutCRF.out
$ etraining --species=species --UTR=on training.gb.train --CRF=1
$ augustus --species=species --UTR=on training.gb.test > test.CRF.out
3. Creating hints from RNA-Seq data with Tophat
先使用 tophat 将 RNA-Seq 数据比对到屏蔽了重复序列的基因组，得到 bam 文件。对 bam 文件进行转换，得到 intron 和 exonpart hints。
得到 intron 的 hints 信息
$ /opt/biosoft/augustus-3.0.3/bin/bam2hints --in=tophat.bam --out=hints.intron.gff --maxgenelen=30000 --intronsonly
得到 exonpart 的 hints 信息
$ /opt/biosoft/augustus-3.0.3/bin/bam2wig tophat.bam tophat.wig
$ cat tophat.wig | /opt/biosoft/augustus-3.0.3/scripts/wig2hints.pl --width=10 --margin=10 --minthresh=2 --minscore=4 --prune=0.1 --src=W --type=ep --UCSC=unstranded.track --radius=4.5 --pri=4 --strand="." > hints.exonpart.gff
cat hints.intron.gff hints.exonpart.gff > hints.rnaseq.gff
4. Creating hints from Protiens with Exonerate
要先得到临近物种的 pritein 序列，然后使用 tblastn 将这些 protein 序列比对到基因组，再将相似性较高的序列使用 exonerate 比对到基因组，对 exonerate 结果进行分析，得到 hint 信息。
将 protein 序列比对到屏蔽了重复序列和转录组序列匹配区域的基因组。这样得到的 hints 主要位于转录组中未表达的基因处，以免和转录组的 hint 有冲突，或得到连接 2 个相邻基因的错误 hint。
$ makeblastdb -in genome.estMasked.fa -dbtype nucl -title genome -parse_seqids -out genome
$ blast.pl tblastn genome proteins.fasta 1e-8 24 blast 5
对 blast 的结果进行分析，挑选相似性高的 pritein 序列用于 exonerate 分析
$ perl tblastn_xml_2_exonerate_commands.pl --exonerate_percent 50 --exonerate_maxintron 20000 --cpu 24 blast.xml proteins.fasta genome.estMasked.fa
将 exonerate 结果转换为 hint 文件
$ /opt/biosoft/augustus-3.0.3/scripts/exonerate2hints.pl --in=exonerate.out --maxintronlen=20000 --source=P --minintron=31 --out=hints.protein.gff
5. Creating hints from RepeatMasker output
在转座子重复区，不倾向于存在编码区，该区域可以作为 nonexonpart hint。将 RepeatMakser 结果转换为此类 hint，有力于 exon 的预测。这样比使用屏蔽了重复序列的基因组进行 Augustus gene prediction 要好。
不对 Simple_repeat, Low_complexity 和 Unknown 这 3 类进行屏蔽，因为这些区域存在 CDS 的可能性相对较大。
$ grep -v -P "position|begin|Unknown|Simple_repeat|Low_complexity" genome.repeat.out | perl -pe 's/^\s*$//' | perl -ne 'chomp; s/^\s+//; @t = split(/\s+/); print $t[4]."\t"."repmask\tnonexonpart\t".$t[5]."\t".$t[6]."\t0\
