《生物软件及应用》课程笔记

使用tophat2和cufflinks进行转录组分析

2018-12-06  本文已影响2人  孤久成瘾_1180

处理原始数据和参考基因组数据后,开始比对分析。
将比对所需参考基因组的索引文件和基因组注释文件存放于hg19文件夹,并将sra文件解压至Fastq文件夹。
主要步骤:
1.用Ttimmomatic对fastq数据去接头。

mkdir trim_out #创建一个存放用Trimmomatic去接头的输出文件夹
for i in Fastq/*fastq
do
        echo $i
        a=$(echo $i | cut -d "/" -f5 | cut -d "_" -f1)
        java -jar ~/Biosofts/Trimmomatic-0.38/trimmomatic-0.38.jar SE -phred33 $i ~/trim_out/$a-clean.fastq ILLUMINACLIP:/home/qiujunhui/Biosofts/Trimmomatic-0.38/adapters/TruSeq2-SE.fa:2:30:10
done

2.序列比对需要tophat2软件。
其优点在于tophat2将待测序列与参考基因组比对后,则直接生成bam文件,且生成的bam文件可以直接给cufflinks构建转录本。

mkdir tophat_out #创建一个用tophat比对的输出文件夹
for x in ~/trim_out/*clean.fastq
do
          echo $x
          b=$(echo $x | cut -d "/" -f5 | cut -d "-" -f1)
          mkdir ~/tophat_out/$b
          tophat -p 4 -o ~/tophat_out/$b ~/hg19/genome $x
done            

-p 指定线程
-o 指定输出目录
hg19/genome为bowtie2官网上人类基因组的索引文件

3.利用cufflinks构建转录本。

#mkdir cufflinks_out 创建一个用cufflinks构建转录本的输出文件夹
cd cufflinks_out
mkdir SRR957677
mkdir SRR957678
mkdir SRR957679
mkdir SRR957680
cd ..
cufflinks -g ~/hg19/gencode.v29lift37.annotation.gtf -o ~/cufflinks_out/SRR957677 -p 4 ~/tophat_out/SRR957677/accepted_hits.bam
cufflinks -g ~/hg19/gencode.v29lift37.annotation.gtf -o ~/cufflinks_out/SRR957678 -p 4 ~/tophat_out/SRR957678/accepted_hits.bam
 cufflinks -g ~/hg19/gencode.v29lift37.annotation.gtf -o ~/cufflinks_out/SRR957679 -p 4 ~/tophat_out/SRR957679/accepted_hits.bam
cufflinks -g ~/hg19/gencode.v29lift37.annotation.gtf -o ~/cufflinks_out/SRR957680 -p 4 ~/tophat_out/SRR957680/accepted_hits.bam

4.使用cuffmerge合并转录组。
为了比较不同样品的差异,需要合并实验组和对照组的转录组,cuffmerge不仅能合并两个或多个转录组,还能合并注释过后的基因组信息,从而找到新的基因通过可变剪切提高合并转录组的质量。

mkdir cuffmerge_out #创建一个存放用cuffmerge合并转录本的输出文件夹
#把上一操作得到的transcripts.gtf的绝对路径写到一个文本文件里面
vi accepted.txt
~/home/qiujunhui/cufflinks_out/SRR957677/transcripts.gtf
~/home/qiujunhui/cufflinks_out/SRR957678/transcripts.gtf
~/home/qiujunhui/cufflinks_out/SRR957679/transcripts.gtf
~/home/qiujunhui/cufflinks_out/SRR957680/transcripts.gtf
#保存
cuffmerge -g ~/hg19/gencode.v29lift37.annotation.gtf -o ~/cuffmerge_out -p 4 accpeted.txt

5.使用cuffdiff进行基因表达差异分析

mkdir cuffdiff _out #创建一个存放用cuffdiff进行基因表达差异分析的输出文件夹
cuffdiff -o ~/cuffdiff_out -L lable1,lable2,lable3,lable4 -p 4 -u ~/cuffmerge_out/merged.gtf ~/tophat_out/SRR957677/accepted_hits.bam  ~/tophat_out/SRR957677/accepted_hits.bam 
~/tophat_out/SRR957677/accepted_hits.bam ~/tophat_out/SRR957677/accepted_hits.bam

ps:merged.gtf为上一步用cuffmerge合并的转录组
-L 后面为bam文件的标签,有几个bam文件便取几个标签,本文以4个为例。

cuffdiff输出文件较多,其对每个基因、每个转录片段、每个编码序列、每个基因的不同剪切体进行FPKM、个数及样本间差异进行分析,最后生成机组不同的文件,按照不同的需求,继续后续分析。
cuffdiff计算每个样本中的转录本、初始转录本及基因FPKM:
(1)traisoforms.fpkm_tracking 转录组的FPKM

gens.fpkm_tracking 基因的fpkm

cds fpkm_tracking 编码序列的fpkm

tss_groups.fpkm_tracking 原始转录组的FPKM

(2)Count tracking files

评估每个样本中来自每个 transcript, primary transcript 及gene的fragment数目

(3)Read group tracking

计算在每个repulate中每个transcript, primary transcript和gene的表达量和frage数目

(4)Differential expression test

对于splicing transcript, primary transcripts, genes

and coding sequences.样本之间的表达差异检验

上一篇 下一篇

猜你喜欢

热点阅读