3.0 转录本比对组装HISAT2+StringTie+cuff

HISAT2

下载参考基因组：Sbicolor_454_v3.0.1.fa sbi301（Phytozome）

#建立索引：
hisat2-build -p 8 Sbicolor_454_v3.0.1.fa sbi301
#比对：
hisat2 -p 16 -x ./sbi301/sbi301 -1 CK18c_11.fq.gz -2 CK18c_22.fq.gz -S CK18.sam
hisat2 -p 16 -x ./sbi301/sbi301 -1 T18_rmrna_1.fq.gz -2 T18_rmrna_2.fq.gz -S T18.sam
转为bam格式：samtools sort -@ 5 -o T18.bam T18.sam

StringTie

参考文章：

https://www.bioinfo-scrounger.com/archives/284/
https://blog.csdn.net/hill_night/article/details/44829965

2020-09-25 (3).png

下载参考注释：Phytozome Sbicolor_454_v3.1.1.gene_exons.gff3.gz

#转格式：
gffread Sbicolor_454_v3.1.1.gene_exons.gff3.gz -T -o sbi311gx.gtf

#组装转录本：
stringtie -p 16 -G ~/sbi311gx.gtf -o T18.gtf T18.bam
stringtie -p 16 -G ~/sbi311gx.gtf -o CK18.gtf CK18.bam
stringtie -p 16 -G ~/sbi311gx.gtf -o CK69.gtf CK69.bam
stringtie -p 16 -G ~/sbi311gx.gtf -o T69.gtf T69.bam

#四样本合并，得到非冗余转录本集
 stringtie --merge -p 8 \
-o stringtie_merge.gtf \
CK18.gtf T18.gtf CK69.gtf T69.gtf

#以stringtie_merge.gtf 为参考注释重新组装
stringtie -e -B -p 16 -G ~/4_stringtie/stringtie_merge.gtf -o T18.gtf ~/T18.bam
#-B代表输出Ballgown可读文件 

Cuffcompare

将各样品的注释文件与参考基因组+参考注释比较，得到新转录本（lncRNA)

• 筛选条件：
  1）FPKM>=0.5 2）coverage >1 3）Length > 200 4）分类为i j o u x
  输入文件：cufcomp.CK18.gtf.tmap 等
  结果文件：class.CK18等 （内容为注释，格式gtf）
• 提取序列
  awk命令，比对转录本id
  以class.CK18 从 CK18.gtf 中提取，文件名为CK18class.gtf
  gffread： 注释文件 、 参考基因组.fa 提取序列
  结果文件：CK18class.fa （G：备份数据）

cuffcompare -r ~/sbi311gx.gtf -s ~/Sbicolor_454_v3.0.1.fa -o cufcomp ../T18.gtf 
awk '{if($7>=0.5 && $10 > 1 && $11 >200) print $0}' cufcomp.T18.gtf.tmap > filter.T18
awk '{if ($3=="u" || $3=="x" || $3=="i" || $3=="j" || $3=="o"){print $0}}' filter.T18 > class.T18
awk 'NR==FNR{a[$5]=1} NR!=FNR{if(a[substr($12,2,length($12)-3)]==1) print $0}' class.T18 T18.gtf > T18class.gtf
#从参考基因组中提取序列
gffread T18class.gtf -g ~/Sbicolor_454_v3.0.1.fa -w T18class.fa