如何理解GWAS中Manhattan plot和QQ plot所传递的信息
如何理解GWAS中Manhattan plot和QQ plot所传递的信息
全文2,123字，阅读6分钟。配图来源：GWAS Catalog
----/ START /----
在GWAS研究中，Manhattan plot和QQ plot是最常画的两类图，它们可以把跟研究的性状（比如，基因型和身高）显著相关的基因位点清晰地展现出来，不少读者朋友应该都懂得如何画这样的图，但我想应该不是每个人都能够真正知道其中所蕴含的道理。
Manhattan plot（曼哈顿图）比较简单，它是把GWAS分析之后所有SNP位点的p-value在整个基因组上从左到右依次画出来。并且，为了可以更加直观地表达结果，通常都会将p-value转换为-log10(p-value)。这样的话，基因位点-log10(p-value)在Y轴的高度就对应了与表型性状或者疾病的关联程度，关联度越强（即，p-value越低）就越高。而且，一般而言，由于连锁不平衡（LD）关系的原因，那些在强关联位点周围的SNP也会跟着显示出类似的信号强度，并依次往两边递减。由于这个原因，我们在曼哈顿图上就会看到一个个整齐的信号峰（如下图红色部分）。而这些峰所处的位置一般也是整个研究中真正关心的地方。
GWAS研究中的曼哈顿图
GWAS研究中，p-value阈值一般要在10^-6次方甚至10^-8次方以下，也就说曼哈顿图中Y轴大于6甚至大于8的那些SNP位点才是比较值得研究的，不过事无绝对，有些时候也要看你的实际数据表现。 另外，曼哈顿图这个名字的由来实际上源自纽约曼哈顿区高楼大厦夜间的灯景在河面上的倒影（如下图）。
曼哈顿夜景
Q-Q plot（QQ图）虽然所用的数据和上面曼哈顿图的一样，但是它要表达的信息比起曼哈顿图来要丰富得多，而且在这两个图中更加能够体现GWAS结果好坏的是QQ plot——它是GWAS研究中更加重要的质控图。这也是我在这篇文章主要讨论的内容。
其实，一直以来QQ plot是统计学分析中的常用图，在1968年Wilk.M.B的这篇文章（doi:10.1093/biomet/55.1.1）就提出了如何绘制这样的图已经它的用途。
QQ plot全称是quantile-quantile plot，也就是分位图，是一种通过比较两个概率分布的分位数从而实现对两个概率分布进行比较的概率图方法(在统计学上较常用)。之所以可以这样做的原因是，如果两个概率分布相同，那么它们的分位数也应该相同或者重叠在同一条直线上。
在GWAS分析中，当我们通过曼哈顿图看到某些SNP和表型性状（或者疾病）有着很强的相关信号（比如，p-value < 10^-6甚至10^-8）时，依然不能直接认为这些位点就与表型显著相关的。这是因为基因组上基因位点的突变通常有两个来源：
- 第一是自然选择（Selection），我这里所说的自然选择不仅指达尔文在《进化论》中所描述的物竞天择，还指所有对物种适应性有影响作用的“力量”，比如高辐射环境、疾病、病毒等，这也是我们在GWAS研究中真正关心的突变；
- 第二是遗传漂变（genetics drift），它是一种比较随机的基因组突变而且数量也不少，虽然也是物种演化的一种重要力量，但是由于它的突变都比较随机，目前认为它与环境的变迁没有必然联系，但也会在某些时候，有些随机的突变带来了生存优势，便会在种群中显示出它的作用。但绝大多数情况下，对于已经在群体中稳定存在的性状而言，并不认为它们有明显的作用，所以GWAS研究是不关心这一类突变的，我们要把它们全部排除掉。如果你发现自己得到的结果全部是这样的变异的话，那么，应该重新考虑一下如何重新设计这个分析，包括是否应该增加样本量以及想办法排除技术错误以及干扰因素等方面，或者也可能它们之间就是没有关系。
强关联信号遗传漂变的存在，在GWAS上会给我们带来一个问题，那就是无法直观地把它们识别并排除出去（很多时候是很难直接在曼哈顿图上发现的），甚至你都判断不了自己的研究是否只是充满了这类无效信息。那么现在的问题就是，我们应该如何做才能够有效地判断出这个研究中所获得关联结果确实是与表型性状或者疾病相关的呢？
这就是我们要用上QQ-plot的地方了。在GWAS分析里面，QQ-plot的纵轴是SNP位点的p-value值（这是实际得到的结果，observed），与曼哈顿图一样也是表示为 -log10(p-value)；横轴是则是均匀分布的概率值（这是Expecte的结果），同样也是换算为-log10。横轴的这个概率值是如何计算的呢？实际上，它就是均匀分布的分位数——至于为何要用均匀分布而不是其它的分布，我会在下一段文字中进行详述。分位数的个数与GWAS研究的SNP位点数是一一对应的。比如我们研究中使用了5百万个基因位点，那么分位数的个数也是5百万个，从1/5000000,2/5000000,3/5000000,...一直往下排直到5000000/5000000，当然都是转换为-log10，然后与GWAS p-value一起作图而成（如下图）。
GWAS研究的QQ图
得到QQ plot之后，如果通过它来协作判断我们的GWAS结果到底是好还是坏呢？
严格来讲，这里其实是不应该用好坏来形容的，应该用是否与表型性状相关来形容。
判断的秘密就在横轴为什么要用均匀分布而不是选择其它分布上。这是因为均匀分布恰好可以用来近似描述基因组上的随机漂变现象。如果表型性状并非真的受自然选择所左右，那么你应该会看到GWAS p-value的分布和均匀分布的结果将集中在一条直线上，如果不是那么就应该能够看到相互分离的情况，特别是p-value越低的时候分离程度就越高，QQ-plot会翘起来（这是因为GWAS的零假设就是与随机突变相比没有区别）。
而且，我们知道基因组上的随机漂变是一定存在的，所以一定会有位点与随机漂变相关，特别是是在p-value比较大的位点看起来就应该和随机漂变重叠，这就表现在QQ-plot的前半部分里。这位点的分布会和均匀分布重叠！而且，比较好的结果是，当p-value < 10^-3时，GWAS结果开始与均匀分布出现快速分离——也就是说，自然选择的力量明显地显示出来了，使得结果在群体中快速摆脱随机性，最后看到一个高高翘起的QQ-plot。这时基本就可以断定，我们所研究的表型和基因型之间是存在着显著相关的自然选择作用的。
这也是我们在评估一个GWAS研究时最基本的一个判断。
----/ END /----
※ ※ ※
你还可以读
。