Augustus操作指南
Augustus指南(Trainning部分)
目录
Augustus指南
官方 Tutorial Index
Input:
Trainning:
Prediction
Output：
可能用到的其他软件
Trainning Augustus
- 对于没有直接基因结构的情况，可以使用Scipio制作基因结构文件
- 将用于Trainning的基因结构信息分为Trainning set和test set
- 为你的物种产生meta参数文件
- 第一次Training
- optimize_augustus.pl的优化
Augustus指南
官方 Tutorial Index
Augustus是一个真核生物基因预测软件，目前有网页服务端和本地版，它基于Hidden-Markov Model（隐马尔科夫链模型HMM）（一个不错的HMM介绍博客）的预测方法，因此需要一个已经研究清楚的物种进行training(学习)之后再对新物种进行prediction（预测），用于trainning的物种应该和需要预测的物种具有较近的亲缘关系
特点：官方介绍
Input:
Trainning:
已知物种序列
核酸序列：fasta格式 ( genome.fa )
蛋白序列：fasta格式 ( proteins.aa )
Prediction
未知基因组序列
核酸序列：fasta格式 ( NewGenome.fa )
Output：
可能用到的其他软件
Scipio
下载地址中下载Scipio Version 1.4，压缩包中包括perl文件和说明文件
Scipio是利用蛋白序列标记出基因组结构的软件，提供核酸序列和蛋白序列后，它可以为你解析出基因的结构。它基于BLAT的序列比对功能。scipio可以很好地跨contigs比对。
用途包括：
检测序列assembly的正确性
识别外显子区域，可变剪切
用已知物种蛋白序列描述未知物种的基因组结构
Scipio依赖：perl中的YAML module
Ubuntu安装yaml-perl命令：
sudo apt-get install libyaml-perl
Scipio依赖：BLAT 下载地址，下载后的文件中有数个可执行文件，请添加进环境变量中，方便Scipio直接调用
Blat，全称 The BLAST-Like Alignment Tool，可以称为"类BLAST 比对工具"，对于DNA序列，BLAT是用来设计寻找95%及以上相似至少40个碱基的序列。对于蛋白序列，BLAT是用来设计寻找80%及以上相似至少20个氨基酸的序列。速度快，共线性输出结果简单易读。Blat 把相关的呈共线性的比对结果连接成为更大的 比对结果，从中也可以很容易的找到 exons 和 introns。
Trainning Augustus
重要提示：如果是第一次training，按照下面所述的步骤可以正常完成，但是如果效果不佳，官方介绍中有一些小经验，本文没有包含这部分内容
Augustus基于Hidden-Markov Model（隐马尔科夫链模型）的预测方法，因此需要一个已经研究清楚的物种进行training(学习)之后再对新物种进行prediction（预测），用于trainning的物种应该和需要预测的物种具有较近的亲缘关系或者就是本物种。
Trainning的实质是用已知的知识调整模型中的种种参数，因此一般已知的知识必须包含各种可能的情况。具体地说，在预测基因的模型中，用于training的数据就必须拥有是基因的序列以及不是基因的序列（基因间隔区，内含子等）。
Trainning的数据会被分为两部分，一部分用于Trainning，另一部分用于检测调整参数后的模型。前者称为training set，后者称为test set。
作为用于training Augustus的已知信息，一般可以有以下几种：
该物种已经存在的基因结构信息（比如：GenBank中的数据）
该物种ESTs信息加上相对于的基因组信息（可使用PASA）
该物种de novo组装的RNA-seq数据加上对应的基因组数据
该物种已知的蛋白序列和对应的基因组数据
近源物种的基因结构信息
已经Trainning过的模型保存参数集后可以反复使用
- 对于没有直接基因结构的情况，可以使用Scipio制作基因结构文件
如果已经有基因结构信息，且文件格式为Genbank格式的文件（*.gb）可以跳过此步。这种文件包含了：
核酸序列
CDS编码区域的序列坐标信息
UTR信息
1.1 Run Scipio
scipio.1.4.pl --blat_output=prot.vs.genome.psl genome.fa proteins.aa > scipio.yaml # takes ~7m
scipio.yaml包含了每个蛋白alignment的细节，可用于下一步提取
输入的基因组5M，耗时约7min
–blat_output=prot.vs.genome.psl的含义就是记录时间。
genome.fa和proteins.aa为fasta格式。
1.2 Extract a GFF file，从scipio.yaml中提取GFF (General Feature Format) 基因结构文件
cat scipio.yaml | yaml2gff.1.4.pl > scipio.scipiogff
scipiogff2gff.pl --in=scipio.scipiogff --out=scipio.gff
cat scipio.yaml | yaml2log.1.4.pl > scipio.log #产生log文件
scipio.gff文件是我们所希望得到的GFF文件，log文件可以帮助检测是否每个蛋白都得到标注。
scipio.gff长这样：
chr2R Scipio CDS 900562 900621 1.000 + 0 transcript_id “392”
chr2R Scipio CDS 904518 904880 1.000 + 0 transcript_id “392”
chr2R Scipio CDS 904940 905131 1.000 + 0 transcript_id “392”
chr2R Scipio CDS 905195 905263 1.000 + 0 transcript_id “392”
chr2R Scipio CDS 3595076 3596041 1.000 + 0 transcript_id “2517”
…
1.3 可选步骤 可视化基因结构（GBrowse）
GBrowse 下载并配置GBrowse，将GFF文件转化成GBrowse可读的文件。
1.4 将GFF文件转化成Genbank格式的文件
Augustus的etrainning学习软件需要输入Genbank格式的文件，这种文件包含了：
核酸序列
CDS编码区域的序列坐标信息
UTR信息
由于大部分的基因结构信息都是GFF或GTF格式，Augustus提供了gff2gbSmallDNA.pl进行格式转化
gff2gbSmallDNA.pl scipio.gff genome.fa 1000 genes.raw.gb
整合gff文件和fasta文件，genes.raw.gb文件中将会有每个基因序列和其1000bp的上下游基因间隔序列。
1.5 去除不适合trainning的基因序列
genes.raw.gb文件已经符合Augustus的要求，但是其中可能包含：
可变剪切的基因
缺失启示密码子的基因
缺失终止密码子的基因
非正常的终止密码子（in-frame stop codons）
这些序列不适合用于软件的trainning，所以最好去除
etraining --species=generic --stopCodonExcludedFromCDS=true genes.raw.gb 2> train.err
得到问题序列文件train.err
–stopCodonExcludedFromCDS=true选项指将gb文件中的终止密码子视为cds中的一部分，这个选项需要视gb文件的来源而定，非scipio来源的文件可能需要设置为false
train.err长这样：
gene 186 transcr. 1 in sequence chr2R_549753-555284: Initial exon has length < 3!
gene 461 transcr. 1 in sequence chr2R_1034318-1036751: Initial Exon doesn’t begin with start codon.
gene 567 transcr. 1 in sequence chr2R_1198437-1201521: Initial Exon doesn’t begin with start codon.
gene 4537 transcr. 1 in sequence chr2R_1354359-1361857: Initial Exon doesn’t begin with start codon.
gene 4783 transcr. 1 in sequence chr2R_1669860-1673241: Terminal exon doesn’t end in stop codon. Variable stopCodonExcludedFromCDS set right?
gene 5161 transcr. 1 in sequence chr2R_2043765-2046183: Single exon gene doesn’t begin with atg codon but with ccc
gene 6319 transcr. 1 in sequence chr2R_3734070-3735386: Initial Exon doesn’t begin with start codon.
gene 3577 transcr. 1 in sequence chr2R_4767472-4770826: Initial Exon doesn’t begin with start codon.
使用下面的命令来筛选掉这些问题序列：
cat train.err | perl -pe ‘s/.in sequence (\S+): ./$1/’ > badgenes.lst
filterGenes.pl badgenes.lst genes.raw.gb > genes.gb
grep -c “LOCUS” genes.raw.gb genes.gb
genes.raw.gb:594 原来序列有594个基因
genes.gb:586 筛选后还剩586个
badgenes.lst为问题基因序列
genes.gb为没有问题的基因，依然是gb格式
- 将用于Trainning的基因结构信息分为Trainning set和test set
随机的将genes.gb文件分为两个部分
randomSplit.pl genes.gb 100
该命令会生成genes.gb.test文件，它包含100条序列。
同时会生成genes.gb.trains文件，包含剩下的序列。
为了满足test的统计学意义，test set必须足够大（100-200 genes），才能起到检测作用。
test set的选择必须满足随机，不能只是选择文件最前面100个或者最后面100个genes
randomSplit.pl是Augustus中包含的软件
3. 为你的物种产生meta参数文件
模型中有两类参数：meta parameters 和 parameters（一般参数）。
前者包括拼接位点模型窗口的大小(the size of the window of the splice site models)、马尔科夫模型顺序(the order of the Markov model)等；
后者则是比如拼接位点模式的分布(the distribution of splice site patterns)、编码区非编码区k-mer概率(the k-mer probabilities of coding and noncoding regions)等。
一般来说meta parameters决定了parameters的算法，而且meta的数量比较少，而parameters则非常多。
Training本质是调整模型中的parameters，但是不会改变meta parameters。
new_species.pl --species=bug
假设我们的物种叫bug
这条命令会在环境变量AUGUSTUS_CONFIG_PATH指定的位置建立一系列文件和文件夹
就像这样：
creating directory /home/mario/augustus/trunk/config/species/bug/ …
creating /home/mario/augustus/trunk/config/species/bug/bug_parameters.cfg …
creating /home/mario/augustus/trunk/config/species/bug/bug_weightmatrix.txt …
creating /home/mario/augustus/trunk/config/species/bug/bug_metapars.cfg …
…
bug_parameters.cfg文件包含了meta parameters和parameters，还包括一些对output文件输出的控制选项。
编辑bug_parameters.cfg文件，将stopCodonExcludedFromCDS选项改为ture
- 第一次Training
4.1 Trainning
这里将进行第一次Training，training所使用的meta parameters是默认值。
etraining --species=bug genes.gb.train
这条命令将会在环境变量AUGUSTUS_CONFIG_PATH指定的位置$AUGUSTUS_CONFIG_PATH/species/bug建立或者更新一些文件，这些文件是exon, intron和intergenic region的parameter文件（training的结果文件）
ls -ort $AUGUSTUS_CONFIG_PATH/species/bug/
展示出与模型相关的各个参数保存文件，
其中：bug_{intron,exon,igenic}.pbl是新生成的文件
4.2 Testing
那么，我们现在就可以使用之前划分出来的test set来检测Trainning的效果
Test set的文件为genes.gb.train
augustus --species=bug genes.gb.test | tee firsttest.out # takes ~1m
这条命令对genes.gb.test中的所有核酸序列信息进行了预测，这里不会使用其中的结构信息的。
之后又将预测信息和gb文件给出的结构信息进行了比较，从而评价trainning的模型。
产生文件为firsttest.out*
grep -A 22 Evaluation firsttest.out
该命令展示了评价部分，输出这样的表格：
******* Evaluation of gene prediction *******
---------------------------------------------
| sensitivity | specificity |
---------------------------------------------|
nucleotide level | 0.975 | 0.89 |
---------------------------------------------/
----------------------------------------------------------------------------------------------------------
| #pred | #anno | | FP = false pos. | FN = false neg. | | |
| total/ | total/ | TP |--------------------|--------------------| sensitivity | specificity |
| unique | unique | | part | ovlp | wrng | part | ovlp | wrng | | |
----------------------------------------------------------------------------------------------------------|
| | | | 115 | 76 | | |
exon level | 511 | 472 | 396 | ------------------ | ------------------ | 0.839 | 0.775 |
| 511 | 472 | | 40 | 5 | 70 | 43 | 3 | 30 | | |
----------------------------------------------------------------------------------------------------------/
----------------------------------------------------------------------------
transcript | #pred | #anno | TP | FP | FN | sensitivity | specificity |
----------------------------------------------------------------------------|
gene level | 118 | 100 | 53 | 65 | 47 | 0.53 | 0.449 |
----------------------------------------------------------------------------/
表格中有三个子表格：
nucleotide level，sensitivity（预测到的百分率），specificity（其中正确的百分率）
exon level， #pred total/unique（预测得到unique外显子总数），#anno total/unique（实际unique外显子总数），TP（正确的预测），FP（假阳性），FN（假阴性）
gene level
100个基因中，预测到53个
83.9%的外显子被预测到
77.5%的外显子预测成功
5. optimize_augustus.pl的优化
脚本optimize_augustus.pl可以通过不断迭代Tranning和Testing的过程，根据模型评价自动修改*_parameters.cfg文件中的meta parameters，使meta parameters的取值最优。
这个脚本的运行会非常费时间，而且最终的效果一般只能提高准确度几个百分点，慎重使用
在本例中Trainning基因组越5M，耗时越1d
使用此脚本augustus和etraining命令必须在环境变量中
5.1 自动调整meta parameters
optimize_augustus.pl --species=bug genes.gb.train # takes ~1d
本句命令最好加上nohup
命令完成后，meta parameters调整完毕。
5.2 retraining Augustus
meta parameters调整好后必须重新Trainning模型，否则没有任何意义。
etraining --species=bug genes.gb.train
同样你可以进行检验
augustus --species=bug genes.test.gb
如果此时你的gene level的sensitivity还是低于20%说明Trainning set不够大，请添加数据。
如果你获得了满意的Trainning结果，请开始prediction吧
Augustus操作指南相关推荐
- 叛乱联机服务器未响应,叛乱沙漠风暴怎么开服 叛乱沙漠风暴开服操作指南详解 安装准备-游侠网...
叛乱沙漠风暴怎么开服?游戏一款多人联机操作游戏,在开服前期要做好相应的准备工作,也就是设置一些选项,这里给大家带来了"xudong162"分享的叛乱沙漠风暴开服操作指南详解,详情一 ...
- oracle读书笔记之Oracle 11g R2 DBA操作指南（前三章）
Oracle 11g R2 DBA操作指南是我学习的第一本讲述oracle数据库知识的书籍,在此记录下自己的学习笔记,对于书中的知识点本人通过实际操作进行验证. 1.Oracle数据库基础 1.1.关 ...
- tcp转串口_PROFIBUS DP与Modbus/TCP网络转换操作指南
>>PROFIBUS DP转Modbus/TCP网关操作指南 摘要: 本文以NT100-RE-DP网关为例,实现Modbus/TCP Client与PROFIBUS DP设备间的通讯,从而 ...
- GitHub使用教程详解——官网操作指南（翻译）
GitHub使用教程详解--官网操作指南(翻译) GitHub使用指南 原文地址:GitHub官网 示例项目:Hello World 十分钟轻松教学 在学习计算机语言编程的过程中创建Hello Wor ...
- 拨号到公司网络操作指南
拨号到公司网络操作指南 转载于:https://blog.51cto.com/itlingm/399066
- 最快捷的Linux命令查询工具来了：「我该怎么做XX」，一句话自动返回操作指南...
贾浩楠 发自 凹非寺 量子位 报道 | 公众号 QbitAI 无论你是Linux萌新用户还是老司机,使用命令行总会遇到不懂的问题. 遇到问题要么去论坛询问.谷歌搜索,要么查自带命令手册(RTFM),过 ...
- 开放应用模型操作指南（一）| 云服务一键接入 OAM 体系
作者 | 邓洪超 阿里云容器平台软件工程师 导读:Open Application Model(OAM)是阿里云联合微软等国际顶级技术团队联合发布的开放应用模型技术.旨在通过全新的应用定义.运维.分 ...
- 计算机文化基础操作考试,(计算机文化基础上机考试操作指南.doc
(计算机文化基础上机考试操作指南 <计算机文化基础>上机考试操作指南 一.考试说明: 本考试为无纸化考试,考生按照监考老师的指引进入考场,在指定的座位参加考试. 二.登录系统 1．用鼠标双 ...
- kubernetest master 节点恢复灾备恢复操作指南
kubernetest master 节点恢复灾备恢复操作指南 本文基本转载别人文章,文末会标明出处. 1. 基本说明 本文档简述了Kubernetes主节点灾备恢复的相关步骤,供在发生k8s m ...
最新文章
- 如何查java当前负载_Linux查看机器负载
- Python 第2周 - Python基础-模块\数据运算
- 为什么越来越多的程序员开始学机器学习？揭秘一个重要原因
- .NET6之MiniAPI(三)：Response
- 教学楼电梯调度需求分析
- getBoundingClientRect使用指南
- 语言余一个负数_Python 中负数取余问题
- 错误信息：org.springframework.beans.factory.NoSuchBeanDefinitionException: No qualifying bean of type
- 结构体别名typedef struct LNode{}LNode,*Linklist疑惑解答
- 北航软件工程考研经验总结
- 微信线下门店二维码扫码支付和退款
- 机电学生写给十年后自己的一封信
- kryo com.esotericsoftware.kryo.KryoException: Buffer underflow.
- jetson nano-mavros-px4
- Egret：一个简单的打砖块游戏
- PAMTRI：用于车辆重新识别的姿势感知多任务学习
- Android SDK 2.1 下载与安装教程
- 【计算机存储单位】字，字节，字符
- Windows11清理C盘空间优化
- i2c协议解读（附读写过程c代码）