转录组分析流程|TransDecoder预测转录本的开放阅读框(二)

您所在的位置:网站首页 转录组unigene和transcript 转录组分析流程|TransDecoder预测转录本的开放阅读框(二)

转录组分析流程|TransDecoder预测转录本的开放阅读框(二)

2023-07-19 08:50| 来源: 网络整理| 查看: 265

使用TransDecoder预测CDS

TransDecoder按照其官网的说明,主要用于识别转录本序列中的潜在的编码区域,也就是预测CDS。转录本可以由RNA-Seq数据通过Trinity组装来的,也可以由RNA-Seq比对到参考基因组上构建的转录本。 最新版本的TransDecoder可在此处找到。

TransDecoder识别可能的编码区域是基于以下几个标准:

在转录本序列中发现一个最小长度的开放阅读框(ORF)

类似于GeneID软件计算的对数似然分数>0

当ORF在第一个阅读框中得分时,与其他5个阅读框中的得分相比,上述编码得分最高

如果发现一个候选ORF被另一个候选ORF的坐标完全封装,则报告较长的ORF。然而,一个转录本可以报告多个ORF(允许操纵子、嵌合体等)

建立/训练/使用PSSM来完善起始密码子预测。

可选:假定肽与噪声截止分数以上的Pfam结构域匹配

Step 1: 提取长开放阅读框

TransDecoder.LongOrfs -t target_transcripts.fasta

默认情况下,TransDecoder.LongOrfs将识别至少100个氨基酸长的ORF。您可以通过’-m’参数降低此值,但可以知道,使用更短的最小长度标准,误报ORF预测的比率会急剧增加。

Step 2: (optional)

通过blast或pfam搜索鉴定与已知蛋白具有同源性的ORF。 参考下面的将同源性搜索包括为ORF保留标准部分。

Step 3: 预测可能的编码区域

TransDecoder.Predict -t target_transcripts.fasta [ homology options ] 最终的文件可以在当前目录找到,也就是后缀为.pep, .cds, .gff3和.bed的文件

一般来说,我们会使用TransDecoder对无参转录组的拼接结果序列预测其CDS,所以我们可以先将拼接序列用BLAST比对nr以及swissprot蛋白数据库,然后提取其比对上的同源序列的位置来识别CDS,最后再通过TransDecoder的第一步和第三步来预测那些未比对上的序列的CDS。

以有参考基因组的转录结果GTF文件预测编码区域

此处的过程与上述过程相同,不同之处在于,我们必须首先生成一个与转录本序列相对应的fasta文件,最后,我们重新计算GFF3格式的基因组注释文件

使用基因组和transcripts.gtf文件构建transcript fasta文件,如下所示:util/gtf_genome_to_cdna_fasta.pl transcripts.gtf test.genome.fasta > transcripts.fasta

接下来,将成绩单结构GTF文件转换为alignment-GFF3格式的文件(之所以这样做,是因为我们的进程对gff3而不是对gtf文件进行操作-没什么大的影响)。像这样将gtf转换为alignment-gff3:util/gtf_to_alignment_gff3.pl transcripts.gtf > transcripts.gff3

运行上述过程以生成最佳候选ORF预测:TransDecoder.LongOrfs -t transcripts.fasta (optionally, identify peptides with homology to known proteins) TransDecoder.Predict -t transcripts.fasta [ homology options ]

最后,生成基于基因组的编码区域注释文件:util/cdna_alignment_orf_to_genome_orf.pl \ transcripts.fasta.transdecoder.gff3 \ transcripts.gff3 \ transcripts.fasta > transcripts.fasta.transdecoder.genome.gff3

输出文件说明

创建了一个工作目录(例如transcripts.transdecoder_dir /)来运行和存储管道的中间部分,其中包含:

1. longest_orfs.pep: 满足最小长度标准的所有ORF,无论编码潜力如何。

2. longest_orfs.gff3: 目标转录本中所有ORF的位置

3. longest_orfs.cds: 所有检测到的ORF的核苷酸编码序列

4. longest_orfs.cds.top_500_longest: top 500最长的ORF,用于训练编码序列的马尔可夫模型

5. hexamer.scores: 每个k-mer的对数似然分数(编码/随机)

6. longest_orfs.cds.scores: 6个阅读框中每个ORF的对数似然和分数

7. longest_orfs.cds.scores.selected: 根据评分标准选择的ORF的加入(如顶部所述)

01. transcripts.fasta.transdecoder.pep: 最终候选ORF的肽序列;删除较长ORF中的所有较短候选项

02. transcripts.fasta.transdecoder.cds: 最终候选ORF编码区的核苷酸序列

03. transcripts.fasta.transdecoder.gff3: 最终选定ORF的目标转录本中的位置

04. transcripts.fasta.transdecoder.bed: 描述ORF位置的bed格式文件,最适合使用GenomeView或IGV查看

最终输出文件如下:

*.pep (是最终的候选ORF编码的蛋白序列)*.cds (最终预测的CDS序列)*.gff3 (是表示ORF的目标转录本的位置)*.bed (用于后期的IGV可视化,以BED格式存放ORF位置信息) 一键化脚本

也有作者发布了脚本,只需要准备参考基因组的数据库位置信息以及相关软件即可自动进行分析,得出结果,非常方便快捷。 具体有需要脚本的后续可以留言~~

2021.12.26更新,之前没有讲相关脚本,最近看到有朋友留言,我在简单写一下后续脚本操作内容

首先脚本下载地址:

Bitbucket​​  

git clone https://bitbucket.org/yangya/phylogenomic_dataset_construction.git

Bitbucket​​下载到本地之后

下载需要作为参考的 cds序列 和在一起做参考,使用python2 运行 刚下载下来的文件夹中 的脚本,具体位置:phylogenomic_dataset_construction/scripts/transdecoder_wrapper.py

python2 /home/lixingze/software/phylogenomic_dataset_construction/scripts/transdecoder_wrapper.py Usage: python transdecoder_wrapper.py transcripts num_cores strand(stranded or non-stranded) output_dir # transcripts替换为前面组装得到的Trinity.fasta # num_cores 为运行线程数 # strand 选项二选一 stranded or non-stranded # output_dir 自定义输出文件夹

我实际运行的命令:

python2 /home/lixingze/software/phylogenomic_dataset_construction/scripts/transdecoder_wrapper.py Trinity.fasta 10 stranded transcoder

还有更多内容可以看上面或者留言中的链接中作者所介绍的,下一期内容点击这里​​​​​​​



【本文地址】


今日新闻


推荐新闻


CopyRight 2018-2019 办公设备维修网 版权所有 豫ICP备15022753号-3