目的
一個(gè)新的物種（沒(méi)有參考基因組）用Trinity從頭組裝得到的轉(zhuǎn)錄本fasta文件牲尺，需要將其轉(zhuǎn)換成gtf文件
思路
- 用http://www.bioinformatics.nl/courses/RNAseq/1c%20Assembly.pdf的思路陪腌，用StringTie得到gtf
- 用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
思路的實(shí)操記錄
- 用
Trinity組裝成fasta文件
- 用
- cd-hit去掉上述fasta文件的冗余序列
教程：①教程 | 如何用cd-hit去除冗余序列缆瓣？ ②無(wú)參轉(zhuǎn)錄組序列組裝 — Trinity
- cd-hit去掉上述fasta文件的冗余序列
-
- 利用
hisat2-build構(gòu)建索引 | 教程： RNASEQ分析入門(mén)筆記3-HISAT2建立索引并比對(duì)序列
hisat2-build -p 10 test.fasta prefix - 利用
-
- 利用以上述組裝得到的fasta文件作為參考序列，用
hisat2將原始的fastq文件比對(duì)到該fasta文件上材泄，得到sam文件
hisat2 -x prefix -1 raw1.fq.gz -2 raw2.fq.gz -S test.sam &>log.txt - 利用以上述組裝得到的fasta文件作為參考序列，用
-
- 將sam文件轉(zhuǎn)換成bam文件
nohup samtools view -bS test.sam > test.bam & -
- 將
bam格式轉(zhuǎn)換成bed6格式
nohup bamToBed -i ./test.bam > test.bed & # 其實(shí)我到了這步已經(jīng)是bed6文件格式了 nohup bed12ToBed6 -i test.bed > test.bed & # 將bed12轉(zhuǎn)換成bed6 - 將
-
- 將
bed6格式轉(zhuǎn)換成genePred格式
bedToGenePred ./test.bed ./test.genePred - 將
-
- 將
genePred格式轉(zhuǎn)換成gtf格式
genePredToGtf file test.genePred test.gtf - 將
或者直接用管道寫(xiě)成一條命令也行
