ChIP-seq数据比对注意事项

测试开头

今天是生信星球陪你的第616天

大神一句话，菜鸟跑半年。我不是大神，但我可以缩短你走弯路的半年~

就像歌儿唱的那样，如果你不知道该往哪儿走，就留在这学点生信好不好~

这里有豆豆和花花的学习历程，从新手到进阶，生信路上有你有我！

豆豆写于2020.5.6
首先是ChIP-seq分析的前言介绍部分：
1：了解ChIP-seq的实验流程
2：继续了解ChIP-seq
3：关于ChIP-seq的实验对照与偏差来源
4：ChIP-seq的实验设计补充
5：ChIP-seq数据库及实战数据介绍
然后开始实战部分：
6：ChIP-seq计算资源准备与实战数据下载
7：ChIP-seq数据质控和过滤
这一次将介绍不同的比对方法和参数设置注意事项

比对的目的很简单，就是“追根溯源”，看看我们过滤后的reads能匹配到基因组的什么位置。

面临的挑战有：

基因组很大（比如人类的有30亿碱基对，也就是60亿碱基）
测序的reads数量很多（30M就是三千万条reads）
测序read和基因组之间并非严丝合缝，中间肯定有不匹配的碱基（可能是测序错误，也有可能是就是和参考基因组不同，存在变异位点），那么到底存在几个mismatch才算匹配成功呢？

1 介绍比对方法和工具

ChIP-seq属于DNA测序，得到的reads应使用：连续短序列比对工具（contiguous short read mappers），一般采用的方法是：种子序列定位及延伸（seed-and-extend）。

第一：通过扫描参考基因组序列，对参考基因组序列建立哈希表，将序列分成一定长度的小片段（k-mer），这种小片段也被称之为种子。根据算法，seed可以完全和基因组匹配或者允许存在一些错配（mismatch）
第二：在目标序列中查找和种子序列相同的片段并标记，以这些标记点为锚点向左右最大限度延伸，并且中间不能有gap
第三：将不合条件的舍弃，符合条件的结果将输出保存。

一篇中文文章：四种常用的生物序列比对软件比较 （http://swxxx.alljournals.cn/html/2016/1/20160110.html）

不同软件的算法在比对层面（sensitivity and precision）和运算层面（消耗的时间和内存）各异；并且选择的seed长度不同，比对和运算结果也不同，一般来讲，短的seed可以提高比对的敏感度，而长的seed可以提高比对速度；大部分算法都会给比对结果提供一个质量值（在SAM/BAM的结果有体现）

关于不同的工具，在这一篇：Ye, H., Meehan, J. et al. (2015). Alignment of Short Reads: A Crucial Step for Application of Next-Generation Sequencing Data in Precision Medicine（https://www.ncbi.nlm.nih.gov/pmc/articles/PMC4695832/）中有介绍，ChIP-seq使用比较多的是Bowtie2，它兼顾了准确性和速度，消耗内存也较小，还支持三代测序的比对；除此以外，BWA作为ENCODE计划使用的比对工具表现也很好。

2 参数设置与注意事项

2.1 错配 mismatch

由于测序存在一定的错误率，还存在一些SNV（single nucleotide variants）的存在，一些reads不会完美比对到参考基因组。但为了不丢失这些测序数据，比对时就可以运行一定数量的错配存在。

关于测序错误率：

一代测序：测序读长长，能达到800-1K bp，且测序用时短，测序准确度高准确性高达99.999%
Illumina二代测序：低于0.1%
三代测序：错误率12-15%

关于SNV：

A single-nucleotide variant (SNV) is a variation in a single nucleotide without any limitations of frequency and may arise in somatic cells. A somatic single-nucleotide variation (e.g., caused by cancer) may also be called a single-nucleotide alteration.

错配设定的阈值取决于样本类型和进行的实验，大部分比对软件都支持以下方式：

使用固定的数字：直接指定每次比对允许的错配数
使用相对的比例：根据read 长度设定错配出现的比例

如果是存在变异的细胞（如癌细胞）进行的ChIP-seq或者有些非模式物种参考基因组质量比较低，那么可以允许更多的错配数；另外还要根据不同测序平台的错误率进行调整这个阈值

2.2 多重比对 Multi-mapping

多重比对就是指reads同时比对到了基因组的多个位置，并且比对质量都不错，尤其在短序列比对中比较明显。一般这种模糊比对的来源是基因组上的重复区域，例如人类基因组上有超过10%的是Alu元件；另外在植物中，串联重复事件的发生频率比较高，例如拟南芥的全基因组重复也会形成重复区域

Alu元件：Alu element是人类基因组中一组散在分布的相关序列，每个长约300bp，由两段130bp的重复序列与一段31bp的间隔序列构成，两侧各有一段6~20bp的正向重复序列，不含反转录酶的编码序列。单个成员的每个末端上有Alu（藤黄节杆菌的缩写）限制酶的切割位点，并由此命名。Alu元件是人类基因组中丰度最高的转座元件

处理这种多重比对，比较保守的办法就是指保留唯一比对的reads；如果要考虑这部分多重比对的reads，有时会将所有的多重比对reads（这样会导致比对结果的条目超过reads数）或者随机选择一条reads作为结果。

ChIP-seq的DNA片段一般在200bp左右，如果有许多reads唯一比对到了较短的重复区域，那么这个结合位点依然能够被捕获到。如果认为的确存在这种情况，那么可以再加上与重复序列数据库Repbase(https://www.girinst.org/server/RepBase/index.php)的比对，一般能得到更高的覆盖度并且对重复区域可以更准确地统计【但这个数据库需要所在机构的订阅】

注释基因组重复序列的RepeatMasker 这款软件就是基于了这个数据库：RepeatMasker makes use of curated libraries of repeats and currently supports Dfam ( profile HMM library derived from Repbase sequences ) and Repbase, a service of the Genetic Information Research Institute.

2.3 其他参数

基因组版本

大部分参考基因组都有许多版本，一般的建议是使用最新版。参考基因组中一般只包含常见的染色体就好（例如人类的chr1-22,X,Y,M），一些未完全组装的scaffold区域可以不要，避免对比对产生干扰。

如果研究物种没有参考基因组，可以试着比对到近缘物种基因组。不过这种情况下，可以指定允许更高的容错能力（调高mismatch参数），毕竟同源基因组还是存在很多的不同。
另外，还可以根据ChIP-seq的reads自己拼接（就像无参转录组的操作），感兴趣的可以看看：He, X., Cicek, A.E. et al. (2015). De novo ChIP-seq analysis