HISAT2+STRINGTIE分析转录组数据

一、构建基因组索引

#! /bin/sh
hisat2-build ref.fa ref

二、序列比对，将clean data 比对到参考基因组。
创建files目录，files下创建samples.txt文件，将样本名称按行写在位于/files目录下的samples.txt文件中。

#! /bin/sh
#（将样本名称按行写在位于/files目录下的samples.txt文件中）
# run some hisat2 alignments

while read line
do
    reads1=${line}_1.fq.gz
    reads2=${line}_2.fq.gz
    hisat2 -p 16 --dta -x ../genome/ref -1 ~/data/cleandata/qingcai_zise/rna-seq/$reads1 -2 ~/data/cleandata/qingcai_zise/rna-seq/$reads2 -S ${line}.sam
    ~/usr/bin/samtools_1.9 sort -@ 16 -o ${line}.bam ${line}.sam
done   <../files/samples.txt

三、 组装转录本并定量表达基因

#! /bin/bash
while read line
do
    stringtie -p 16 -G ../files/gene.gff3 -o ${line}.gtf -l ${line} ../aligment/${line}.bam
done < ../files/samples.txt

四、专录本合并并比较
在files目录下新建mergelist.txt文件，文件写入上一步生成的每个gtf完整路径。

/sevzone/home/rna_seq/assemble/Greenvegetables_Z113-01-T01_good_p.gtf
/sevzone/home/rna_seq/assemble/Greenvegetables_Z113-01-T02_good_p.gtf
/sevzone/home/rna_seq/assemble/Greenvegetables_Z113-01-T03_good_p.gtf
/sevzone/home/rna_seq/assemble/Greenvegetables_Z113-01-T04_good_p.gtf
/sevzone/home/rna_seq/assemble/Greenvegetables_Z113-01-T05_good_p.gtf
/sevzone/home/rna_seq/assemble/Greenvegetables_Z113-01-T06_good_p.gtf

合并GTF

stringtie --merge -p 8 -G ../files/gene.gff3 -o stringtie_merged.gtf ../files/mergelist.txt

gffcompare -r ../files/gene.gff3 -G -o merged_stringtie_merged.gtf stringtie_merged.gtf

五、重新组装转录本并估算基因表达丰度，并为ballgown创建读入文件
新建目录expression，基于新合并的GTF计算表达量。

#! /bin/bash
while read line
do
    stringtie -e -B -p 8 -G ../assemble/merged_stringtie_merged.gtf.annotated.gtf -o ./${line}/${line}.gtf ../aligment/${line}.bam

done < ../files/samples.txt