映射至基因组mapping

1、映射至基因组(Mapping)第一步的工作是比对(alignment)o对于RNA-seq的比对，从来都不是一件容易 的事情。其难点如下：1. 没有很好的比对模板。现在的比对模板都是基因组模板，而不是真正的转 录组模板，也就是说，这对本来就不是很长的短序来说，它很有可能是界 于两个exon之间。我们在比对junction的时候，一般还是假设它如果没 能在基因组模板中找到合适的位置的时候，才考虑它是否是界于junction 上。这种人为的假设可能并不准确。2. SNPs,碱基插入，删除，错配，或者质量不高的测序结果，从模板至比 对序列本身，都存在着比基因比对更为复杂的问题。3. 短序可能会有多个

2、100%的匹配位点。4. 有些基因组可能需要庞大的内存空间。为了解决最后一个问题，人们使用了很多办法，但基本上都会基于事先建立的引 索库。即所谓"启发式"比对(heuristic match)o首先使用一定长度的(通常是11 个碱基)的序列做为索引用的关键字，在匹配这一索引字之后，就很大程度地缩 小了其需要匹配的模板范围。但是这一办法的问题在于不容易解决问题2中的 空格，错配问题。所以在很多软件使用时，会要求人工确认高保真区，以及最高 允许2?3个错配。现在比较快的“启发式”比对主要有两种算法，一种是哈希表(hashtable), 一种是 BW Ik缩转换(Burrows

3、Wlieeler transform, BWT)。前者速度快，但是对内存要 求比后者要高。对于问题3, 般而言，大部分软件使用的办法是只保留一个匹配位点，其中, 有些是只保留第一个匹配位点，有些是按照概率分布选取保留的位点。当然，前 面已经提到过，可以使用paired-end read来尽量避免问题3的出现。对于问题1,可以使用外显子库来确定junction readso有两种办法，一种是依 靠已知的外显子库来构建，另一种办法就是依据已经匹配好的短序来构建外显子 库(de novo assembly of transcriptome) 后者的不足是运算量大，对测序覆盖范 围要求高，最好是使用p

4、aired-end readso还有人发现，对于ploy(A)的处理会减少不能映身的短序数。比如，Pickrell et al. 就发现，对于46bp的Illumina reads, 87%的短序可以映射至模板，7%可以映 射至junction libraryc如果对那些不能映射的短序，将在头或者尾含有的超过连 续4个的A或者T去除，就可以得到约0.005%的映射。11综合评价 (Summarizing mapped reads)这一步，主要是基本于不同水平(外显子水平，转录水平，或者基因水平)进行 统讣。最简单的办法就是统计落在每个外显上的短序数。但是有研究表明，很多 (可能超过15%)的短

5、序会落在外显子两侧，这会影响统计的结果。另一种办 法就是统会落在内显子区域的短序数。无论如何，即使是基因水平的综合评价，也还是有其它的一些问题。比如 overlapping的基因的统计。比如junction的统计。标准化(Normalization)标准化对于样品内及样品间的比较而言是非常重要的。标准化被分为两类，样品 内及样品间(between- and within-library)。样品内标准化使得在同一样品内不得基因之间的表达差异变得有意义。最常用到 的一个办法就是使用落在同一基因内的短序数除以单位基因长度。比较常用的 单位是 RPKM (reads per kilobase of e

6、xon model per million mapped reads)。但是 这一方法也受到样品制备和测序方法的干扰。而对于样品间标准化，最简单而直接的办法使用短序总数来平衡表达量。然而短 序总数受测序深度的干扰，而且单个基因的短序数与实际的表达量并不一定会 呈线性比较关系。人们乂使用四分位(quantile normlization)标准化的办法。但是 有研究说这一办法并没有实际的价值。还有提出使用对数分布法则(power law distributions)来进行样品间标准化。但没有研究对这一处理方式进行验证。差异表达(Differential expression)差异表达分析的最终LI

7、的是将那些差异表达的基因(外显子等等)从海量数据中 提取出来。最终的结果显示一般来说是表格化的，这一表格按照一定的规则排序, 让人们能够尽可能简单地拿到想要的结果。山于RNA-seq结果的离散性，人们一般都会使用统计模型来拟合实验得到的结 果。一般而言，RNA-seq的结果是比较附合伯松分布(poisson distribution) W 这一结果得到了单通道Illumina GA测序结果的实验验证。但是，伯松分布分析 结果常常在多组重复的样品间带来较高的假阳性，因为它低佔了生物取样的样品 间误差。所以RNA-seq如何设 置重复是一个很重要的问题。为了平衡重复样品 所带来的误差，人们使用 了

8、 serial analysis of gene expression (SAGE) data。现有的软件一般都是针对较为简单的实验设讣的。而对于复朵的实验设讣，比如 说成对样品，时间依赖样品等等，还没有专门的，较好的解决方案。大多数都使 用edgeR的线性模型来进行分析。22后期系统生物学分析简单地讲，前景是广阔的，但U前为止手段还是比较有限的，基本上就是GO分 析。在转录组高通量测序数据分析过程中，主要就是转录本的构建（拼接）及表达水 平的衡量，其次还包括可变剪接之类的转录后修饰等研究分析。因此，做过 RNA-seq转录组测序数据分析的童鞋们一定使用过或至少听说过Cufflinks这个 软

9、件。Cufflinks是加利福尼亚大学伯克利分校数学和计算机生物实验室，llj LiorPachter 领导的StevenSalzberg团队，和马里兰大学生物信息和计算机生物中心的 Steven Salzberg小组，以及加州理工学院的Barbara Wold实验室联合作用的结 果。Cufflinks利用Top hat比对的结果（alignments）来组装转录本，估汁这些 转录本的丰度，并且检测样本间的差异表达及可变剪接。这个软件其实是个套装, 包括四个部分分别命名为：cufflinkscuffcompare、cuffmerge 及 cuffdiff。tophat->Cufflink

10、s->cuffdiff转录组分析大致流程如下：第一步，利用tophat/bowtie比对结果（bam格式）及参考基因组构建转录本, 最终的转录本是以gtf格式保存的。第二步.Cuffcompare主要是对两个或多个转录本集合中转录本相似悄况的比较, 例如将笫一步构建出的转录本与ENSEMBL数据库中的转录本进行比较，评佔转 录本构建情况，此外，根据构建的转录本与已知ENSEMBL数据库中的转录本的 相对位置定义了一系列分类，例如内含子区域、反义、基因间区域转录本等等近 10种分类。第三步，cuffmerge是将多个转录本集合合并成一套转录本集合，例如将在多个 组织样本中构建的多套转录本合

11、并成一套转录本，cuffinerge能够很好地完成去 除冗余。第四步，cuffdiff衡量两个或多个样本间差异表达的基因，例如癌症与正常组织 33间差异表达的转录本，此外还能衡量差异可变剪接体。至此,转录本测序常规数据分析基本结束，接下来进行实验验证或深入数据分析。0K,如果对cufflinks的流程理解没问题的话，现在问题来了！ ！在tophat->cufflinks -> cuffdiff流程中，第一步对每一个样本会产生一个的gtf 格式的转录本，因此若有n个样本就会有n个gtf.但是在对n个样本之间做差 异表达分析时，由于各样本之间的转录本集合不相同。在做cuffdiff时，

12、却要求 输入一个gtf，以便对该gtf文件中转录本在不同样本间进行差异分析。因此需 要将n个样本的gtf文件合并成一个。问题是，cufflinks有cuffcompare和 cuffhierge两个合并转录本的方法，这是应该选用那个？这两个有什么区别？个人理解及经验分享：这是所有初做NGS分析的很因惑的问题，事实上，对于多个样本构建的多套转 录本，如何得到统一的一套转录本，有下面三种方法：1在转录本拼接之前，把各样本的比对bam文件合并，然后用合并的bam跑 cufflinks2每个样本的比对bam文件分别单独跑cufflinks，各样本的转录本构建后，再 用cuffcompare合并为一套转

13、录本。3每个样本的比对bam文件分别单独跑cufflinks，各样本的转录本构建后，再 用cuffmerge合并为一套转录本。这三种做法区别在于：第一种方法流程相对简单，所有的工作都抛给cufflinks 一人完成，你都不需要 知道cuffmergecuffcompare的用法。貌似是种完美解决方案。但很大的问题 是：cufflinks能处理得了最终合并的bam吗？对于小物种的样本还可以，但对 于人，若是7, 8个样本合成的bam, cufflinks吃不销！ 44笫二和三种方法是类似的，都是在保留可变剪切结构的询提下，将转录本合并。 不同的是，cuffcompare只有A、B两条转录本结构相同的时候，才将A、B合并。 而cuffmerge是A、B某些部分互相overlap,就将它俩合并。事实上，cuffmerge 再做合并的时候，是把overlap的transfrag重新调用了 cufflinks,合成一个