PacBio测序表现出了非常强大的优势,比如测序read长、无GC偏好、直接检测变异、直接检测碱基修饰等。尤其在转录中的应用已经相对成熟,下面重点介绍如何去进行序列的比对和去冗余。
关键词:PacBio、minimap2、GMAP、collapse、full length isoforms
因为我们知道Iso-seq测序后经过smrtlink和IsoSeq3软件的处理会得到高质量的全长转录本序列,针对有参考基因组的物种,首先要进行序列的回帖(比对)。
回帖软件有minimap2, GMAP, STAR, BLAT等,重点介绍minimap2和GMAP。
一、比对
1. minimap2
软件链接:https://github.com/lh3/minimap2
该软件支持剪切比对和非剪切比对,所以非常适合转录本的比对,推荐使用2.9及以上版本,且支持不建索引模式的比对。
使用示例如下:
minimap2 -t 30 -ax splice -uf --secondary=no -C5 hg38.fasta hq_isoforms.fasta > hq_isoforms.fasta.sam 2> hq_isoforms.fasta.sam.log
-ax spliced比对模式
--secondary=no 只输出最好的比对结果
2.GMAP
软件链接:http://research-pub.gene.com/gmap/
索引文件必需提前建好(gmap_build命令),推荐使用2018-03-20或更高版本(版本以日期命名)。
使用示例:
gmap -D /gmap_db/ -d hg38 -f samse -n 0 -t 16 --cross-species --max-intronlength-ends 200000 -z sense_force hq_isoforms.fasta > hq_isoforms.fasta.sam 2> hq_isoforms.fasta.sam.log
1和2步骤中,获得了比对文件后(SAM格式)可进一步转换为BAM格式:
samtools view -bS hq_isoforms.fasta.sam > hq_isoforms.fasta.bamsamtools sort hq_isoforms.fasta.bam > hq_isoforms.fasta.sorted.bamsamtools index hq_isoforms.fasta.sorted.bam
二、去冗余
比对文件中记录了低质量比对和同于的基因、转录本异构体,需要进一步进行过滤。
是的,没错,CupCake可以处理。
软件链接:https://github.com/Magdoll/cDNA_Cupcake
直接使用里面的collapse_isoforms_by_sam.py进行处理:
usage: collapse_isoforms_by_sam.py [-h]
[--input INPUT] [--fq] -s SAM -o PREFIX [-c MIN_ALN_COVERAGE] [-i MIN_ALN_IDENTITY] [--max_fuzzy_junction MAX_FUZZY_JUNCTION] [--flnc_coverage FLNC_COVERAGE] [--dun-merge-5-shorter]
比对处理和去冗余:
gmap -D <gmap_db_location> -d hg38 -f samse -n 0 -t 12 -z sense_force hq_isoforms.fastq > hq_isoforms.fastq.sam
sort -k 3,3 -k 4,4n hq_isoforms.fastq.sam > hq_isoforms.fastq.sorted.sam
collapse_isoforms_by_sam.py --input hq_isoforms.fastq --fq -s hq_isoforms.fastq.sorted.sam --dun-merge-5-shorter -o test -c 0.95 -i 0.85
后续可以通过再IGV中观察hq_isoforms.fastq.bam和生成的去冗余后的test.collapsed.gff文件。
另外,
test.group.txt为记录合并冗余后的对应文件,可以知道一共保留了多少个非冗余转录本。test.ignored_ids.txt为没有比对上而被丢弃的转录本编号。
参考:
https://github.com/Magdoll/cDNA_Cupcake/wiki/Best-practice-for-aligning-Iso-Seq-to-reference-genome:-minimap2,-GMAP,-STAR,-BLAT#refgmap
https://github.com/Magdoll/cDNA_Cupcake/wiki/Cupcake-ToFU:-supporting-scripts-for-Iso-Seq-after-clustering-step#collapse