文章目录
STAR for ENCODE Transcriptome
STAR 全称是Spliced Transcripts Alignment to a Reference (STAR),直白的说,就是用后缀数组算法做转录组的比对(Spliced Transcripts == 剪切过的转录本)。号称 aligning to the human genome 550 million 2*76 bp paired-end reads per hour on a modest 12-core server, 突出一个快字。做 Spliced Transcripts 数据分析有两个难点,一是 reads 包含mismatch,insertions and deletions,做到精准 alignment 比较困难。二是转录本是来自不连续的基因组区域的,具体点解释,rna-seq 拿到的RNA是基因组上不同位置的序列剪切,拼接来的,现在需要把 reads 往回拼,需要处理一下。STAR的诞生主要是为了处理ENCODE 转录组transcriptome项目中产生的 80 billion Illumina reads。
STAR 的算法
STAR 主张将非连续的序列直接 align 到参考基因组上。主要分为 searching step和 Clustering, stitching and scoring 两步。
1. Seed search
算法核心为 Maximal Mappable Prefix (MMP),类似 Maximal Exact (Unique) Match 的最大匹配算法。假设有一条read R,坐标 location i,参考基因组 G。首先把 R 通过后缀树算法在 G 上定位,找到把可能性排序,由于R 可能由于剪切拼接而来,R 的片段可能会匹配到 G 的不同位置。
比如图(a)里面的这个Read,MMP1,MMP2 分别匹配到了 genome 上的两个位置。首先,从read 的第一个碱基开始比对。很明显,图上的这个 read 是 match 不到连续 genome 上去的,但是,我们可以找到剪切位点。然后,MMP算法再从剪切位点开始 match。因为MMP 算法是从没有match的位置开始算的,所有速度非常快,这是 MMP 算法区别于其他算法的特色。这个算法并不考虑剪切拼接点是否合理,也没有引入参考的junction database。找文字匹配的方法是基于uncompressed suffix arrays (SAs) 后缀法。除了splice junction(前文说的挑战中的第二个),关于基因的 mismatches & indels(前文说的挑战中的第一个),MMP 会插入锚点,并进行 matching。如果是一整段 read 的序列都不行(可能是 poly A,library adapter 或者单纯就是poor sequencing)就trim掉。此外,MMP search 是双向的。
2. Clustering, stitching and scoring
在这个阶段,STAR 把 read 比对上的基因组片段 seeds 组装(stitching together)起来,建立alignments。STAR 提出一个最优程序,就是通过限制锚点 anchors 优化 seeds 的集中流程。把这些seeds集中在一起,我们就大概知道这条 read 上的序列都来自基因组那些位置了。关于 paired-end read STAR 有额外的处理流程。如果 genome 的一个 window 装不下某条 read,STAR会找更多的window 直到 cover 整条 reads,这个功能适合融合基因等状况。对于错配、缺失啥的,标准评分如下。以这个标准,选分最高的 read 的 alignments。
STAR 的用法:
构建 index
STAR --runThreadN 6 --runMode genomeGenerate --genomeDir ~/project/pengG4/rtstop/index_hg38_2/ --genomeFastaFiles ~/project/pengG4/rtstop/index_hg38_2/hg38.fa --sjdbGTFfile ~/project/pengG4/rtstop/gtf/hg38_ensemble.gtf --sjdbOverhang 100
May 11 09:32:55 ..... started STAR run
May 11 09:32:55 ... starting to generate Genome files
May 11 09:34:57 ..... processing annotations GTF
May 11 09:35:44 ... starting to sort Suffix Array. This may take a long time...
May 11 09:36:28 ... sorting Suffix Array chunks and saving them to disk...
May 11 11:45:46 ... loading chunks from disk, packing SA...
May 11 11:52:40 ... finished generating suffix array
May 11 11:52:40 ... generating Suffix Array index
May 11 11:59:40 ... completed Suffix Array index
May 11 11:59:41 ..... inserting junctions into the genome indices
May 11 12:07:21 ... wr