目的

一个新的物种（没有参考基因组）用Trinity从头组装得到的转录本fasta文件，需要将其转换成gtf文件

思路

1. 用http://www.bioinformatics.nl/courses/RNAseq/1c%20Assembly.pdf的思路，用StringTie得到gtf
   
2. 用http://seqanswers.com/forums/showthread.php?t=44455的思路

bwa-mem of Trinity.fasta to reference genome
samtools -Sbh | bam2bed | bed12ToBed6 | bedToGenePred stdin stdout | genePredToGtf file stdin assembly.gtf

思路的实操记录

• 1. 用Trinity组装成fasta文件
• 2. cd-hit去掉上述fasta文件的冗余序列
     教程：①教程 | 如何用cd-hit去除冗余序列？ ②无参转录组序列组装 — Trinity
• 3. 利用hisat2-build构建索引 | 教程： RNASEQ分析入门笔记3-HISAT2建立索引并比对序列

  hisat2-build -p 10 test.fasta prefix
  

• 4. 利用以上述组装得到的fasta文件作为参考序列，用hisat2将原始的fastq文件比对到该fasta文件上，得到sam文件

  hisat2 -x prefix -1 raw1.fq.gz -2 raw2.fq.gz -S test.sam &>log.txt
  

• 5. 将sam文件转换成bam文件

  nohup samtools view -bS test.sam > test.bam &
  

• 6. 将bam格式转换成bed6格式

  nohup bamToBed -i ./test.bam > test.bed & # 其实我到了这步已经是bed6文件格式了
  nohup bed12ToBed6 -i test.bed > test.bed &   # 将bed12转换成bed6
  

• 7. 将bed6格式转换成genePred格式

  bedToGenePred ./test.bed ./test.genePred
  

• 8. 将genePred格式转换成gtf格式

  genePredToGtf file test.genePred test.gtf
  

或者直接用管道写成一条命令也行