Cufflinks基于reads与参考基因组的比对来进行全长转录本的组装，不需要额外提供gff/gtf这种描述gene/transcript的feature文件，软件会生成gff文件。
1. 比对
第一步仍然是建索引比对，这里用到tophat2（以bowtie2为搜索引擎）。
bowtie2-build -f Arabidopsis_thaliana.TAIR10.dna.toplevel.fa TAIR10.bybowtie
mv Arabidopsis_thaliana.TAIR10.dna.toplevel.fa TAIR10.bybowtie.fa #和一些比对软件（比如（二）中介绍的bowtie2）不同的是，tophat比对需要参考基因组跟索引有相同的前缀；而像bowtie2只需要利用索引即可，所以没有这个要求
$ ls
TAIR10.bybowtie.1.bt2 TAIR10.bybowtie.2.bt2 TAIR10.bybowtie.3.bt2 TAIR10.bybowtie.4.bt2 TAIR10.bybowtie.fa TAIR10.bybowtie.rev.1.bt2 TAIR10.bybowtie.rev.2.bt2
for i in `seq 2 7`
do
tophat2 -r 230 -p 8 -o ./SRR328680${i} ~/learn_rnaseq/mRNA-seq/bowtie2/ref/TAIR10.bybowtie \ #r表示insert size
~/learn_rnaseq/mRNA-seq/data/SRR328680${i}_1.fastq.gz \
~/learn_rnaseq/mRNA-seq/data/SRR328680${i}_2.fastq.gz
done
#输出文件，其中accepted_hits.bam后面会用到
$ ls -t
unmapped.bam logs/ accepted_hits.bam junctions.bed deletions.bed insertions.bed align_summary.txt prep_reads.info
2. 组装
2.1 官方流程
Cufflinks的流程图，可以看出当有多个处理（或重复）时，流程建议分开跑tophat和cufflinks，然后再把分别组装得到的转录本合并。
Cufflinks在组装的时候可以提供额外的gff文件用于指导组装过程。
for i in {2..7}
do
cd ~/learn_rnaseq/mRNA-seq/map/map_tophat2/SRR328680${i}
cufflinks -p 2 -o ./creat/ ./accepted_hits.bam &
done
之后每个creat文件夹下面会生成4个文件
genes.fpkm_tracking isoforms.fpkm_tracking skipped.gtf transcripts.gtf
接下来将这6个样本的transcripts.gtf合并
$ cat list
SRR3286802_transcripts.gtf
SRR3286803_transcripts.gtf
SRR3286804_transcripts.gtf
SRR3286805_transcripts.gtf
SRR3286806_transcripts.gtf
SRR3286807_transcripts.gtf
nohup cuffmerge -p 8 list -s ~/learn_rnaseq/mRNA-seq/bowtie2/ref/TAIR10.bybowtie.fa -g ~/learn_rnaseq/mRNA-seq/ref/new.GCF_000001735.4_TAIR10.1.gtf &
#之后生成一个merged_asm文件夹，merged.gtf为合并后的转录本gtf文件
$ ls merged_asm/
logs/ merged.gtf
之后做差异表达，结果相当丰富
nohup cuffdiff ./merged_asm/merged.gtf \
SRR3286802_accepted_hits.bam,SRR3286803_accepted_hits.bam,SRR3286804_accepted_hits.bam \
SRR3286805_accepted_hits.bam,SRR3286806_accepted_hits.bam,SRR3286807_accepted_hits.bam &
网友评论