使用TransDecoder寻找转录本中的编码区

50 阅读 0 评论 33 点赞

我是靠谱客的博主多情奇异果，最近开发中收集的这篇文章主要介绍使用TransDecoder寻找转录本中的编码区，觉得挺不错的，现在分享给大家，希望可以做个参考。

概述

TransDecoder能够从转录本序列中鉴定候选编码区。这些转录本序列可以来自于Trinity的从头组装，或者来自于Cufflinks或者StringTie的组装结果。

软件安装

从https://github.com/TransDecoder/TransDecoder/releases下载最新版的TransDecoder，以v5.5.0为例

mkdir -p ~/opt/biosoft && cd ~/opt/biosoft
wget https://github.com/TransDecoder/TransDecoder/archive/TransDecoder-v5.5.0.zip
unzip TransDecoder-v5.5.0.zip
mv TransDecoder-TransDecoder-v5.5.0 TransDecoder-v5.5.0

运行TransDecoder

我们从cufflinks或stringtie输出的gtf文件开始分析流程，因此你会有两个输入文件

transcripts.gtf: 记录预测转录本的GTF文件
genome.fasta: 参考基因组序列

第一步: 从GTF文件中提取FASTA序列

~/opt/biosoft/TransDecoder-v5.5.0/util/gtf_genome_to_cdna_fasta.pl transcripts.gtf genome.fasta > transcripts.fasta

第二步: 将GTF文件转成GFF3格式

~/opt/biosoft/TransDecoder-v5.5.0/util/gtf_to_alignment_gff3.pl transcripts.gtf > transcripts.gff3

第三步: 预测转录本中长的开放阅读框, 默认是100个氨基酸，可以用-m修改

~/opt/biosoft/TransDecoder-v5.5.0/TransDecoder.LongOrfs -t transcripts.fasta

第四步: 使用DIAMOND对上一步输出的transcripts.fasta.transdecoder.pep在蛋白数据库中进行搜索，寻找同源证据支持

# 下载数据并解压缩
wget ftp://ftp.uniprot.org/pub/databases/uniprot/current_release/knowledgebase/complete/uniprot_sprot.fasta.gz
gunzip uniprot_sprot.fasta.gz
# 建立索引
diamond makedb --in uniprot_sprot.fasta --db uniprot_sprot.fasta
# BLASTP比对
diamond blastp -d uniprot_sprot.fasta -q transcripts.fasta.transdecoder_dir/longest_orfs.pep --evalue 1e-5 --max-target-seqs 1 > blastp.outfmt6

关于DIAMOND的使用，参考这篇说明DIAMOND: 超快的蛋白序列比对软件

第五步: 预测可能的编码区

~/opt/biosoft/TransDecoder-v5.5.0/TransDecoder.Predict 
    -t transcripts.fasta 
    --retain_blastp_hits blastp.outfmt6

第六步: 生成基于参考基因组的编码区注释文件

~/opt/biosoft/TransDecoder-v5.5.0/util/cdna_alignment_orf_to_genome_orf.pl 
     transcripts.fasta.transdecoder.gff3 
     transcripts.gff3 
     transcripts.fasta > transcripts.fasta.transdecoder.genome.gff3

最终输出文件如下: