A tutorial on conducting genome‐wide association studies : Quality control and statistical analysis
A tutorial on conducting genome‐wide association studies : Quality control and statistical analysis
最近在深入学习全基因组关联分析，本篇文献是一份很基础的GWAS GuideLine，每位想入门GWAS的朋友，可以首先看看这篇手把手教GWAS的文献。（强烈建议大家阅读文献的英文原文，并跑跑作者的代码！）
Introduction
本篇文章主要是提供一份guideline来帮助大家进行GWAS相关的数据分析。
GWAS研究的目的是识别等位基因频率随表型的变化而系统变化中的单核苷酸多态性（SNP） 【例如精神分裂症人群和健康人群之间的病例对照等等】。 识别与表型相关联的SNP可以为揭示造成不同表型的生物学机制提供新的思路。得益于技术的进步，现在我们可以全面的调查分布在整个基因组中的snp对生物体的影响。
目前，GWAS方法在探究snp对精神病症的影响方面已经有很多成功的例子，例如：精神分裂症，精神紊乱等等。这些所提到的GWAS研究，在很大的程度上依赖于我们对人类基因组遗传结构的深度了解，这些对人类基因组遗传知识来源于HapMap计划((http://hapmap.ncbi.nlm.nih.gov/; Gibbs et al., 2003)和千人基因组计划（1000 genomes， http://www.1000genomes.org/）。HapMap计划描述了常见的snp人类DNA序列中的模式。千人基因组计划则是对常见的snp和稀少的snp在基因组水平上绘制map图
由于单个个体snp对整个GWAS结果的影响还是比较小，所以精神病领域的研究者们开发了整合snp效果的方法模型来提高GWAS结果的准确性。本篇文章中我们会主要聚焦多基因风险评分(polygenic risk score, PRS)分析。PRS方法将多个SNP的效应整合在一起，来对疾病进行预测和评分
。。。 暂时省略
本篇文章是一份Guidline，为哪些没有遗传背景，但是想将遗传学方法应用在自己研究领域的研究者提供帮助。
文章的主体纲要如下：
- 首先，我们会展示如何对基因型数据进行严格的质量控制，以保证后续GWAS的成功进行。这一部分内容我们还会讨论如何选用合适的方法来处理不同种族之间的异质性所引起的问题。
- 第二，我们会介绍常见的用于关联SNP和表型性状的统计模型，并阐述如何控制潜在的干扰因子。
- 第三，我们会展示如何进行一个多基因风险评分分析。
需要注意的是，我们这篇文章没有写如何进行基因型填充，所以对基因型填充有需要的同学，可以看看文章 van Leeuwen and colleagues (2015). 本文所有的代码都存放在Github上 https://github.com/MareesAT/GWA_tutorial/ 。
准备软件
- PLINK（Version 1.09） ，下载地址：https://www.cog‐genomics.org/plink/1.9/.
- R语言
- Linux环境
- 其它处理SNP数据的软件推荐：Genabel，SNPTEST
1.数据格式准备
因为我们主要是用PLINK来完成GWAS分析，所以我们需要认真了解一下输入PLINK软件的数据格式。
对于PLINK来说，它既可以处理文本格式的文件，也可以处理二进制格式的文件。但是大文本的文件处理起来十分消耗计算资源，所以我们一般 推荐使用二进制格式的输入文件 。
- 文本格式的PLINK数据包括两份文件
.ped文件和.map文件.ped文件包含个体信息（例如个体标识符ID，性别等等）以及他们的基因型信息.map文件包含遗传标记的信息（染色体号，snp号等等）
- 二进制格式的PLINK数据则包括三份文件
.bed文件，.fam文件和.bim文件.bed文件含有 每个个体的识别符（ID）和每个个体相对应的基因型.fam文件含有 个体信息（例如性别之类的）.bim文件含有 遗传标记的信息（染色体号，snp号等等）
例如目前要做一个精神狂躁症的研究，那么 .bed文件 主要包括的就是病人们和健康对照组的基因型 。 .fam文件 则需要写入和这项研究相关的数据。（例如家庭关系，性别，疾病症状等等） 。.bim文件 就要提供snp的物理位置等等信息。如果我们需要进行协变量分析，那么往往要加入第四份文件，这份协变量文件包含了每个协变量的值。
2. PLINK的基础使用命令
‐‐file {your_file}输入文本格式文件‐‐bfile {your_file}输入二进制格式文件--assoc关联分析,这一步会对每个SNP和研究者感兴趣的性状进行卡方检验--out {outfile}输出文件
更多命令可见官网：http://zzz.bwh.harvard.edu/plink/
3. 遗传数据的质控
任何GWAS研究中，很重要的一步就是数据的质控，如果没有有效的质量控制，那么GWAS研究得出的结论也没有意义。
数据中的错误可以由多个因素引起，例如：DNA样本的低质量，DNA在芯片上杂交的时候不是很成功，检测基因型的芯片探针有误，样本混杂或者出现了污染。
由于质控环节出现了问题，Sebastiani et al. (2010)等人在Science上发表的文章也被撤回了，所以我们一定要注意数据的 质量控制 ！
3.1 使用HapMap数据进行数据模拟
为了使用真实的数据解释GWAS中的每个步骤，我们通过HapMap计划的公共数据集，模拟出了一个样本数为207的数据集（按照二进制文件的格式来模拟的）。对于本篇指南而言，为了构建同源数据集，所以我们只选了来自北欧和西欧的犹他州的居民。
因为我们选取的模拟数据集，相对于HapMap计划来说较小，所以模拟数据集中的遗传效应会比我们做常规的复杂性状的遗传研究中的遗传效应更大。很重要的一点是，大样本集（至少1000个样本，但大多数是上万个样本）需要检测复杂性状的遗传风险因素。
我们对HapMap数据集进行模拟得到的性状存放在：https://github.com/MareesAT/GWA_tutorial/ (1_QC_GWAS.zip) 中
3.2 数据质量控制步骤概览
表格中给出了进行质量控制的7个步骤，其中质控的阈值需要根据研究的问题来动态调整。
质控（过滤SNP和去掉不好的样本）的这七个步骤基于以下内容：
SNP missingness个体信息和SNP的缺失- 命令行： --geno, --mind
- 功能描述：对于 – geno， 去除在大部分受试者中都缺失的snp； 对于 --mind，去除基因型缺失较大的个体。
- 阈值设置：我们推荐在第一次过滤SNP和个体的时候选用一个较为宽松的阈值（0.2; >20%），这将过滤掉具有非常高水平缺失的SNP和个体。然后，再使用一个更为严格的阈值来进行过滤(0.02) 。
需要注意的是必须先进行SNP过滤，才能进行个体过滤
sex discrepancy受试者信息填写的性别和遗传性别不一致(见性别差异)- 命令行：–check-sex
- 功能描述：根据X染色体杂合度/纯合率检查数据集中记录的个体性别与遗传性别之间的差异【就是看记录的性别信息和实际的性别是否有出入】
- 阈值设置：这一步骤可以帮助我们发现样本是否发生了混淆。如果多个受试者都存在这种差异，则应该仔细检测数据。男性的X染色体纯合度估计值应该 > 0.8，女性的这个值应该 < 0.2
minor allele frequency (MAF)最小等位基因频率- 命令行：–maf
- 功能描述：仅选取高于设定MAF阈值的SNP
- 阈值设置：具有低MAF的SNP是罕见的，因此缺乏检测SNP-表型关联的能力。 这些SNP也更容易出现基因分型错误。 MAF阈值应取决于您的样品
大小，较大的样本可以使用较低的MAF阈值。 分别对于大（N = 100.000）对中等样品（N = 10000），通常使用0.01和0.05作为MAF阈值。
Hardy–Weinberg equilibrium (HWE)哈迪-温伯格平衡(HWE)偏离- 命令行：–hwe
- 功能描述： 排除偏离Hardy-Weinberg平衡的标记
- 阈值设置：hwe是检验基因型分型是否发生错误的常见指标。（他也可能指示发生了进化选择）。对于二元性状来说（例如性别）我们建议排除病例组中HWE - P值小于1e-10的标记，和控制组中 HWE - P值小于1e -6的标记。对于数量性状来说，我们推荐HWE - P值小于1e-6
heterozygosity rate杂合率- 命令行：见脚本
- 功能描述：去除具有较高或者较低杂合率的个体
- 阈值设置：偏差可能表明样品污染，近亲繁殖。我们建议去除偏离样本的杂合率平均值±3 SD的个体。
relatedness亲缘性- 命令行： – genome, --min
- 功能描述：-- genome 用于计算两两样本直接的血缘同一性 ； --min 设置阈值并创建具有高于所选阈值的相关性的个体列表。意味着可以检测到与例如pi-hat> 0.2（即，第二度亲属）相关的受试者
- 阈值设置： 使用独立的SNP进行此项分析，并且将它限制在常染色体中。隐秘相关性可能会干扰关联分析。 如果您有基于家庭的样本（例如，父母后代），则不需要删除相关对，但统计分析应考虑家庭相关性。 然而，对于基于人群的样本，我们建议使用0.2的pi-hat阈值，这与文献一致
population stratification群体分层- 命令行：–genome, --cluster --mds-plot k
- 功能描述： - -genome 用于计算两两样本直接的血缘同一性 。‐‐cluster ‐‐mds‐plot k 基于IBS 产生数据中的K维子结构
- 阈值设置： 需要使用独立的SNP进行此项分析，并且将它限制在常染色体中。K是需要我们定义的维数（通常为10）。 这是质量控制的一个重要步骤，由多个程序组成，但出于完整性的原因，我们简要地在表中叙述这一步骤。 该步骤将在“控制群体分层”部分中更详细地描述。
我们的在线教程提供了手把手的按照上述1-7步进行质量控制的资源 https://github.com/MareesAT/GWA_tutorial/ (1_QC_GWAS.zip + 2_Population_stratification.zip)
如果一个样本包括了多个种群（例如 非洲，亚洲，欧洲），我们建议将这些种群分别分开后，在各自的种群中进行关联分析，并且使用合适的方法来整合这些结果（例如可以使用meta分析）
4. 群体结构分层
GWAS系统偏差的一个重要来源是人口结构分层，如专栏1所述。已经表明，即使是单一种族人口中的微妙程度的人口分层，也会引起偏差。
因此，测试和控制种群分层的存在，是一个必不可少的质量控制（QC）的步骤。目前有很多方法可以用来校正群体结构分层(Price, Zaitlen, Reich, & Patterson, 2010)。在本篇文章中，我们以PLINK软件中包含的**多维缩放法（MDS）**来举例说明 如何进行群体结构分层。
多维缩放法（MDS） 计算样本中任何一对个体之间共享的等位基因的全基因组平均比例，以生成每个个体的遗传变异的定量指标（组分）【这句英文确实没看懂】 。 我们可以通过绘制单个组分的分数，来探究是否存在遗传上彼此比预期更相似的个体组。举例：在一个包含了亚洲和欧洲人种的遗传研究中，MDS计算就会显示亚洲人之间的遗传相似性会高于欧洲人。
通过这一步可以有效的找到那些偏离整体目标群体的游离样本。这些游离样本我们在之后都会去除。
5. 关联分析
在进行了质量控制，并计算了MDS组分后，数据就可以用于接下来的关联检验了。我们需要通过性状的特质，以及研究的情况来选择合适的统计检验模型。在接下来的指南中，我们会提供适用于二元特征的各种类型关联的脚本（例如，酒精依赖患者 VS 健康人群） 或者数量性状 （例如，每周消费的酒精饮料数量）。
PLINK提供一个自由度的等位基因检验，在这个检验中，性状的值或者是二元性状取对数后的值，会随着风险等位基因数量的变化增加或者减少。另外，非累加检验是必要的，例如基因型关联检测，显性基因功能测验，隐形基因功能测验。但是，非累加检验并不是应用的很广泛，因为利用统计模型检验非累加效应实际效果是比较差的。加入基于R的插件后，PLINK中还可以使用更多复杂的分析。
5.1 二元结果测量
在PLINK软件中，snp和二元性状（值1 =未受影响，值2 =受影响; 0和-9表示缺失;这些PLINK中的默认选项是可以更改的）的关联检验，可以使用参数 --assoc 或者 --logistic 来完成。
--assoc参数，可以执行卡方检验，需要注意的是卡方检验不允许包含协变量--logistic参数，可以执行回归分析，允许包含协变量
--logistic 参数比 --assoc 参数更加灵活，但是它消耗的计算时间也更多。
5.2 数量性状
在PLINK软件中，snp和数量性状之间的关联检验可以使用参数 --assoc 或者 --linear 来完成。
当PLINK检测到数量性状时（例如值超过了1，2，0 或者是缺失），--assoc参数会自动的使用t检验的渐进版本来比较两个均值。这个选项不允许使用协变量。
而 --linear 参数把每个独立的SNP作为预测因子，执行线性回归分析。和--logistic选项类似，--linear选项可以使用协变量，但是执行速度比-assoc选项要慢一些
5.3 多重检验校正
现代的基因型分型芯片，可以同时对400万个标记进行基因分型，从而大大增加了检验的次数，产生多重检验的负担。SNP基因型填充也会增加关联检验的次数。各种模拟都表明，无论研究的实际SNP密度如何，当基因组范围的显著性设置在5x10-8时，在欧洲人群队列上，我们可以完全控制全基因组中独立SNP的数量。
由于个体的遗传多样性，当检测非洲人群队列时，阈值会设置的更严格（基本接近1x10-8）。【Hoggart, Clark, De Lorio, Whittaker, & Balding, 2008】
目前有三种方法来检测全基因组范围的显著性：Bonferroni校正，Benjamini-Hochberg错误发现率（FDR）和置换检验（permutation testing）
6. PRS分析
待完善
1万+
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
