Control for Population Structure and Relatedness for Binary Traits in Genetic Association Studies via Logistic Mixed Models

Han Chen; Chaolong Wang; Matthew P. Conomos; Adrienne M. Stilp; Zilin Li; Tamar Sofer; Adam A. Szpiro; Wei Chen; John M. Brehm; Juan C. Celedón; Susan Redline; George J. Papanicolaou; Timothy A. Thornton; Cathy C. Laurie; Kenneth Rice; Xihong Lin

doi:10.1016/j.ajhg.2016.02.012

美国人类遗传学杂志。2016年4月7日；98（4）：653–666。

2016年3月24日在线发布。数字对象标识：10.1016/j.ajhg.2016.02.012

预防性维修识别码：项目经理4833218

PMID：27018471

基于Logistic混合模型的遗传关联研究中种群结构和二元性状相关性的控制

韩晨（Han Chen）,^1,⁸ 王朝龙,^1,^2,⁸ 马修·P。科诺莫斯,^三阿德里安·M·。斯蒂尔普,^三李子林,^1,⁴ 塔马尔·索弗,^三亚当A。斯匹诺,^三魏晨,⁵ 约翰·M·。布雷姆,⁵ 胡安·C。塞莱多恩,⁵ 苏珊·雷德兰,⁶ 乔治·J·。帕帕尼古劳,⁷ 蒂莫西A。桑顿,^三凯西·C。劳里,^三肯尼思·赖斯,^三和林喜红^1,^∗

韩晨

¹美国马萨诸塞州波士顿哈佛T.H.Chan公共卫生学院生物统计学系，邮编02115

查找文章依据韩晨（Han Chen）

王朝龙

¹美国马萨诸塞州波士顿哈佛T.H.Chan公共卫生学院生物统计学系，邮编02115

²新加坡基因组研究所计算与系统生物学，新加坡138672

查找文章依据王朝龙

马修·P。科诺莫斯

^三美国华盛顿大学生物统计系，西雅图，WA 98195

查找文章依据马修·P。科诺莫斯

阿德里安·M·。Stilp（静止）

^三美国华盛顿大学生物统计系，西雅图，WA 98195

查找文章依据阿德里安·M·。斯蒂尔普

李子林

¹美国马萨诸塞州波士顿哈佛T.H.Chan公共卫生学院生物统计学系，邮编02115

⁴清华大学数学系，北京100084

查找文章依据李子林

塔马尔·索弗

^三华盛顿大学生物统计学系，西雅图，华盛顿，98195，美国

查找文章依据塔马尔·索弗

亚当A。斯匹诺

^三美国华盛顿大学生物统计系，西雅图，WA 98195

查找文章依据亚当A。斯匹诺

魏晨

⁵美国宾夕法尼亚州匹兹堡市匹兹堡大学UPMC匹兹堡儿童医院儿科肺内科、过敏和免疫科

查找文章依据魏晨

约翰·M·。布雷姆

⁵美国宾夕法尼亚州匹兹堡市匹兹堡大学UPMC匹兹堡儿童医院儿科肺内科、过敏和免疫科

查找文章依据约翰·M·。布雷姆

胡安·C。塞莱多恩

⁵美国宾夕法尼亚州匹兹堡市匹兹堡大学UPMC匹兹堡儿童医院儿科肺内科、过敏和免疫科

查找文章依据胡安·C。塞莱多恩

苏珊·雷德兰

⁶美国马萨诸塞州波士顿布莱根妇女医院医学和神经病学系睡眠和昼夜节律紊乱科02115

查找文章依据苏珊·雷德兰

乔治·J·。帕帕尼古劳

⁷美国马里兰州贝塞斯达NIH国家心脏、肺和血液研究所心血管科学部预防和人口科学项目，邮编：20892

查找文章依据乔治·J·。帕帕尼古劳

蒂莫西A。桑顿

^三美国华盛顿大学生物统计系，西雅图，WA 98195

查找文章依据蒂莫西A。桑顿

凯西·C。劳里

^三美国华盛顿大学生物统计系，西雅图，WA 98195

查找文章依据凯西·C。劳里

肯尼思·赖斯

^三华盛顿大学生物统计学系，西雅图，华盛顿，98195，美国

查找文章依据肯尼思·赖斯

林喜红

¹美国马萨诸塞州波士顿哈佛T.H.Chan公共卫生学院生物统计学系，邮编02115

查找文章依据林喜红

作者信息文章注释版权和许可信息 PMC免责声明

关联数据

补充资料: 文件S1。图S1–S12和表S1和S2
毫米c1.pdf（400万）
GUID:8D16ABEA-0F91-49E8-A3FF-A27D2F2DE975
文件S2。条款和补充数据
mmc2.pdf格式（460万）
制导：A7105D4F-10B6-4164-9D3B-D5F53DEFC4CC

摘要

线性混合模型（LMM）广泛应用于全基因组关联研究（GWAS），以解释连续性状和二元性状的群体结构和相关性。由于LMM无法控制哮喘GWAS中的I型错误，这是一个二元特征，我们表明，当人口分层导致违反LMM的恒残留方差假设时，LMM通常不适合分析二元特征。为了克服这个问题，我们开发了一种计算效率高的logistic混合模型方法，用于二进制性状的全基因组分析，即广义线性混合模型关联检验（GMMAT）。该方法每GWAS拟合一次logistic混合模型，并在二元性状和个体遗传变异之间无关联的无效假设下进行分数测试。我们在模拟研究和实际数据分析中表明，在分析各种研究设计中的二元性状时，GMMAT有效地控制了种群结构和相关性。

介绍

群体分层是遗传关联研究中混淆的主要来源。¹随着计算效率算法的发展，线性混合模型（LMM）在全基因组关联研究（GWAS）中越来越流行，用于控制种群分层以及家族或隐秘关联性。²^,^三^,⁴^,⁵^,⁶^,⁷^,⁸^,⁹^,¹⁰然而，在人类遗传学中，GWAS通常使用二进制特征进行；例如，病例对照研究对疾病受试者（病例受试者）和无病受试者进行抽样，并回顾性收集其基因型和暴露数据。队列研究前瞻性地跟踪一组健康受试者，收集他们随后对二元疾病结局的状态评估。尽管LMM在二进制性状GWAS中的使用很受欢迎，但它们对普遍无效假设的依赖似乎被忽视了。具体来说，正如通常使用的那样，LMM假设性状具有恒定的剩余方差，而在存在协变量的情况下，二元性状通常会违反这一假设。因此，我们将表明，在存在群体分层的情况下，将LMM拟合到二元性状可能无法控制I型错误率，并产生不正确的p值估计。LMM的典型使用也忽略了病例对照研究中的有偏抽样，这导致了有偏效应估计。

我们对二元特征LMM的研究是由西班牙裔社区健康研究/拉丁裔研究（HCHS/SOL）中的二元特征哮喘分析引发的，最初使用具有三种随机效应的LMM来解释其多阶段抽样设计中的遗传相关性以及家庭和街区组成员。对于哮喘，血统是西班牙裔/拉丁裔人群中已知的遗传关联混杂因素，在HCHS/SOL中，哮喘病例的比例因血统群体而异（例如，波多黎各人为25.6%，南美洲为4.1%，见表1).¹¹尽管通过包括主成分（PC）对祖先进行了回归调整¹²^,¹³而祖先群体的指标作为协变量，LMM对该性状的结果似乎无效(图1)与所有其他祖先群体相比，具有次要等位基因的变异体在波多黎各遗传中分别较不常见/更常见，显示出明显的保守性/反保守性。

在单独的窗口中打开

图1

HCHS/SOL中哮喘GWAS分析的关联检验p值的分位数-分位数图

（A）所有SNP。

（B）第1类：波多黎各人与非波多黎各人为预期方差比小于0.8的SNP。

（C）第2类：波多黎各人与非波多黎各人为预期方差比率在0.8至1.25之间的SNP。

（D）第3类：波多黎各人与非波多黎各人为预期方差比大于1.25的SNP。

缩写如下：LMM，使用LMM对组合样品进行联合分析；LMM meta是一种反向变量加权固定效应meta分析方法，用于结合分别分析波多黎各人和非波多黎各人民的LMM结果。

表1

样本质量控制后HCHS/SOL中六个西班牙裔/拉丁裔背景组的哮喘病例受试者和对照受试者人数

组	案例受试者	对照受试者	样本大小	受试者比例	性状方差
中美洲	55	1,173	1,228	4.5%	0.043
古巴	182	1,722	1904年	9.6%	0.086
多米尼加	99	933	1,032	9.6%	0.087
墨西哥人	172	4,189	4,361	3.9%	0.038
波多黎各人	478	1,390	1,868	25.6%	0.190
南美	34	795	829	4.1%	0.039

在单独的窗口中打开

性状方差计算为病例和对照受试者比例的乘积。

根据公布的病例对照GWAS，¹⁴^,¹⁵^,¹⁶^,¹⁷这种令人担忧的行为似乎并不为人所知。二进制特征LMM的用户已经向Pirinen等人发出了呼吁。，⁹世卫组织表示，当遗传变异的影响大小很小并且没有人口分层时，线性模型与logistic模型在p值计算方面非常接近。我们将表明，在存在人口分层的情况下，在队列和病例对照研究中，将LMM与二元性状进行拟合，会因异方差性而导致GWAS中的I型错误率不正确，即模型中的剩余方差相对于协变量的不同值不是恒定的。¹⁸与线性回归一样，LMM假设同方差，即无论协变量值如何，剩余方差都是常数(图2).¹⁹

在单独的窗口中打开

图2

二元性状的真均值-方差关系和线性模型假设的常数均值-方差的关系，以HCHS/SOL哮喘数据为例

对于具有平均值的二进制特征π，其方差为π(1 −π)，随平均值变化。这种异方差性可以通过逻辑回归得到适当解释。线性模型不恰当地假设二元性状的方差不随平均值而变化，并且是一个常数（同方差）。例如，波多黎各人的二元特征（哮喘状态）的方差远远大于其他五个人口的方差，因为波多黎各人口的哮喘疾病比例远远高于其他人口。由人口分层引起的这种异方差性导致LMM计算的p值可能不正确，但使用GMMAT的logistic混合模型适当考虑了这一点。

最近，人们提出了用于病例对照研究的责任阈值线性混合模型。²⁰^,²¹拟合这些模型需要了解基础研究人群中的疾病流行率以及准确的遗传力和责任估计，这在实践中可能无法获得或可能难以确定，尤其是当疾病流行率在不同亚人群之间存在差异时。例如，在HCHS/SOL中，不同祖先背景组的哮喘患病率不同；¹¹在存在人群分层的病例对照研究中，当潜在亚群未知且疾病流行率在潜在亚群之间存在差异时，很难指定可用于分析的疾病流行率。我们将证明，在中等或强人口分层的情况下，这些方法可能无法控制I型错误率。

为了解决这些问题，我们建议使用logistic混合模型²²在二元性状的遗传关联研究中考虑群体分层和相关性，同时自然地考虑到这些性状的非恒定方差。由于逻辑混合模型的计算成本高于LMM和常规逻辑模型，因此它们尚未用于大规模GWAS。尽管SAS PROC GLIMMIX允许使用遗传关系矩阵拟合逻辑混合模型，但对于大规模GWAS来说，这是不可行的，因为其计算负担与拟合基因组中的大量逻辑混合模型有关，每个变体一个。我们开发了广义线性混合模型关联测试（GMMAT），该测试要求每个GWAS在零假设下只拟合一次logistic混合模型，因此在计算上对大规模GWAS是可行的。

具体地说，GMMAT首先适用于零logistic混合模型，该模型仅包括作为固定效应的协变量，包括祖先PC，但不包括个体遗传变异，并使用随机效应来解释固定效应PC未捕获的剩余人口分层以及相关度。在拟合该零模型时，GMMAT使用惩罚拟似然（PQL）²²以及计算效率高的平均信息受限最大似然（AI-REML）算法。⁶^,²³这种拟合零模型对于GWAS中的所有遗传变异都是相同的。然后GMMAT进行分数测试²⁴针对每个个体的遗传变异来测试其与二元性状的关联。与拟合零logistic混合模型的成本相比，每个遗传变异的分数测试的计算成本可以忽略不计，因此该过程的计算速度很快，甚至对于大规模GWAS也是如此。

作为一种完整的建模方法，GMMAT解释了性状的二元性，特别是在存在协变量时的非恒定方差，因此在存在人口分层和相关性的情况下正确控制了I型错误率。GMMAT还可以进行复杂的抽样设计，如分层设计，并通过合并多种随机效应，在抽样个体之间实现未观察到的共享环境效应。

材料和方法

Logistic混合模型和分数测试

对于单变量测试，我们考虑以下logistic混合模型：

罗吉特(π_我) = X（X）_我α+G公司_我β+ b条_我,

哪里π_我= P（P）(年_我= 1|X（X）_我, G公司_我, b条_我)是受试者的二元表型（例如疾病状态）的概率我，取决于其协变量、基因型和随机效应b条_我,X（X）_我是1×第页主题协变量的行向量我,α是一个第页×1固定协变量效应的列向量，包括截距，G公司_我是受试者遗传变异的基因型我、和β是基因型效应。我们假设 $b条 \sim N个 (0, \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})$ 是一个n个×1随机效应的列向量，其中τ_k个是方差分量参数和V（V）_k个已知的n个×n个相关矩阵。当方差分量的数量K（K）= 1,V（V）₁通常是从大量遗传变异中估计的遗传关系矩阵。我们允许多重随机效应来解释复杂的抽样设计，例如分层设计和环境协方差结构。二元表型年_我假设独立于随机效应b条.

执行零假设的分数测试H（H）₀:β=0，我们需要拟合零logistic混合模型，这对于所有遗传变异都是一样的

罗吉特(π_i0（i0）) = X（X）_我α+ b条_我,

（方程式1）

哪里π_我0= P（P）(年_我= 1|X（X）_我, b条_我). 我们拟合方程1使用惩罚准似然（PQL）方法。²²具体来说，让W公司= d日我一克{v（v）_我0}，其中v（v）_我0= π_我0(1 − π_我0), $X（X） = {(\begin{matrix} {X（X）}_{1}^{T型} & {X（X）}_{2}^{T型} & \dots & {X（X）}_{n个}^{T型} \end{matrix})}^{T型}$ 成为n个×第页包含截距的协变量矩阵，以及 $\tilde{年}$ 是包含组件的“工作向量” ${\tilde{年}}_{我} = {X（X）}_{我} α + {b条}_{我} + {v（v）}_{我 0}^{- 1} (年_{我} - π_{我 0})$ .在零假设下H（H）₀:β=0，我们迭代拟合工作线性混合模型 $\tilde{年} = X（X） α + b条 + ε$ 哪里ε ∼ N个(0, W公司⁻¹). 我们使用计算效率高的AI-REML算法⁶^,²³估计 ${\hat{τ}}_{k个}$ 获得方差分量估计后，固定效应的估计 $\hat{α}$ 和随机效应 $\hat{b条}$ 可以更新，后面跟着工作矢量 $\tilde{年}$ 。这一过程一直持续到趋同。

在每次迭代中，我们基于Cholesky分解执行矩阵反演，其复杂性为O（运行）(n个^三)和复杂度为的矩阵乘法O（运行）(pn码²+K（K）²n个²+第页²n个)，其中n个是样本量，第页是协变量的数量，以及K（K）是方差分量的数量。因此，将逻辑混合模型与K（K）已知关联矩阵是O（运行）(在里面^三+智能功率网络²+iK公司²n个²+知识产权²n个)，其中我是达到收敛所需的迭代次数。

的分数H（H）₀:β=0是 $T型 = {G公司}^{T型} (年 - {\hat{π}}_{0})$ ，其中 $G公司 = {(\begin{matrix} {G公司}_{1} & {G公司}_{2} & \dots & {G公司}_{n个} \end{matrix})}^{T型}$ 是n个×1基因型柱载体， $年 = {(\begin{matrix} 年_{1} & 年_{2} & \dots & 年_{n个} \end{matrix})}^{T型}$ 是n个×1列结果向量，以及 ${\hat{π}}_{0}$ 是下拟合值的向量H（H）₀，这对所有SNP都是一样的。分数的估计方差为 $V（V）一第页 (T型 | {H（H）}_{0}) = {G公司}^{T型} \hat{P（P）} G公司$ 在零假设下，其中 $\hat{P（P）} = {\hat{Σ}}^{- 1} - {\hat{Σ}}^{- 1} X（X） {({X（X）}^{T型} {\hat{Σ}}^{- 1} X（X）)}^{- 1} {X（X）}^{T型} {\hat{Σ}}^{- 1}$ 、和 $\hat{Σ} = {\hat{W公司}}^{- 1} + \sum_{k个 = 1}^{K（K）} {\overset{图6}{τ}}_{k个} {V（V）}_{k个}$ 每个遗传变异的测试涉及向量-矩阵-向量乘法和两个向量的内积计算，因此分数测试步骤的复杂性为O（运行）(量子数²+量子数)，其中q个是测试的遗传变异的总数。此外，如有必要，可以轻松并行分数测试步骤。我们使用C++库Armadillo执行矩阵和向量计算。关于拟合logistic混合模型和执行分数测试的更多详细信息，请参见附录A.

HCHS/SOL基因分型与统计分析

同意遗传研究的HCHS/SOL受试者从血液中提取DNA。这些DNA样本在Illumina Microarray Services的SOL HCHS Custom 15041502 B3阵列（由Papanicolaou、Rotter和Taylor设计和开发的定制内容）上进行基因分型。该阵列包括Illumina Omni 2.5M阵列（HumanOmni2.5-8v1-1）和为HCHS/SOL选择的其他自定义内容，包括祖先信息标记、美洲印第安人特有的变体以及已知的GWAS点击。²⁵Illumina Microarray Services、LA Biomed和SOL遗传分析中心（GAC）进行了质量控制，以生成推荐的SNP和样品级质量过滤器。检查样本的注释性别与遗传性别、染色体总体异常、亲缘关系和群体结构、漏接率、批量效应、重复样本不一致和孟德尔错误。在SNP水平上，检查Hardy-Weinberg平衡、次要等位基因频率（MAF）、重复探针不一致性和失访率。

研究参与者是通过多阶段调查抽样程序招募的，其中个人是在家庭内进行抽样的，而家庭是在街区组内进行抽样。²⁶该研究包括12803名拉美裔/拉丁裔个体的遗传数据，这些个体属于六个不同的拉美裔或拉丁裔背景群体。²⁵HCHS/SOL研究由参与机构的机构审查委员会批准，并获得所有参与者的书面知情同意。标准质量控制，与Laurie等人描述的类似。，²⁷用于筛选SNP和质量较差的样本。此外，删除了缺少哮喘诊断信息的样本和主成分分析中确定为祖先异常值的56个样本，剩下11222个样本用于分析。我们还过滤了整个样本中MAF小于5%的SNPs，得到了一组1299221个待分析的常染色体SNPs。在中美洲、古巴、多米尼加、墨西哥、波多黎各和南美血统的每个遗传背景组的分析参与者中，哮喘病例的比例分别为4.5%、9.6%、9.6%，3.9%、25.6%和4.1%(表1). 波多黎各群体明显是其他群体中的局外人，因此我们将其与所有其他背景群体的集合进行了比较，后者的合并比例为5.8%。为了说明抽样设计产生的潜在共享环境影响，除了标准化遗传关系矩阵外，还包括代表家庭和街区组成员的相关矩阵，线性和逻辑混合模型共有三种随机影响。前五名个人电脑以及拉美裔/拉丁裔背景群体指标被用于调整两种模型中的血统。其他固定效应协变量包括现场中心、性别、年龄、吸烟量、卷烟包装年数和抽样重量（即，根据抽样调查设计对样本进行不成比例选择的权重²⁶). 以这种方式将抽样权重作为固定效应处理，可以有效地调整设计变量的边际效应。²⁸

模拟研究

我们进行合并模拟²⁹从20×20网格中为总共8000名创始人生成基因型，其中有1000000个独立SNP(图S1)用相邻细胞之间的迁移率模拟空间上连续的种群（每个细胞20个个体）M（M）=10代表欧洲人口结构。³⁰^,³¹为了模拟后代队列的基因型，我们首先从原始队列中为20×20网格中的每个细胞取样（不替换）10对父母，然后使用基因删除算法为每个家庭模拟两个孩子，³²从而在后代队列中获得8000个个体。我们将这两个队列结合起来，得到16000个同时具有人口结构和关联性的个体。对于个人j个在家里我，成为病例的概率π_我j个计算依据为

罗吉特(π_ij公司) = α₀+ α₁Z轴_我+ 一_我j个,

哪里Z轴_我=1，如果是家庭我来自左上角的10×10网格，并且Z轴_我否则=0。参数α₀和α₁在左上角的高危人群组中，疾病流行率为0.28，无随机效应，在低危人群组中为0.05，无随机影响。随机效应模拟为

一_{我} = (\begin{matrix} 一_{我 1} \\ 一_{我 2} \\ 一_{我 三} \\ 一_{我 4} \end{matrix}) 至 N个 ((\begin{matrix} 0 \\ 0 \\ 0 \\ 0 \end{matrix}), σ^{2} (\begin{matrix} 1 & 0 & 0.5 & 0.5 \\ 0 & 1 & 0.5 & 0.5 \\ 0.5 & 0.5 & 1 & 0.5 \\ 0.5 & 0.5 & 0.5 & 1 \end{matrix})),

其中方差分量参数σ²设置为2。我们从原始和后代组合队列中随机抽取10000名个体，形成一项模拟队列研究，并使用625583个MAF大于5%的创始群体遗传变异计算其标准化遗传关系矩阵。我们删除了家族指标，并使用遗传关系矩阵比较了线性和logistic混合模型。我们将前十个祖先PC作为协变量纳入了这两个模型。我们分析了MAF大于5%的常见遗传变异。在罕见疾病的模拟病例对照研究中，我们选择α₀和α₁因此，左上角的高危人群组的患病率为0.045，无随机效应，而低危人群组的发病率为0.005。我们从原始和后代队列中随机抽取1667名病例受试者和8333名对照受试者，形成一项病例对照研究，总样本量为10000，病例对照比例为1:5。我们进行了与上述相同的分析。

结果

HCHS/SOL哮喘GWAS数据分析

在基于人群的HCHS/SOL队列研究中，我们使用GMMAT比较了LMM和logistic混合模型的结果，以分析医生诊断的哮喘。

图1A显示了哮喘分析中每种方法的总分位数（QQ）图的行为。LMM的膨胀结果是由于违反了其常数剩余方差（同方差）假设。对于二元性状，群体分层影响了群体特异性均值（队列研究中的疾病流行率）和性状方差；疾病流行率接近0.5的人群差异较大(表1和图2). 线性模型和LMM中假设的均值-方差关系对于二元性状来说是错误的。尽管可以通过调整种群或祖先PC来解释种群结构的混淆，但不同种群之间的不平等二进制性状方差并没有在线性模型和LMM中得到适当的建模。当具有较大二元性状方差（即较高疾病风险）的群体也具有较高的MAF，从而具有较大的基因型方差时，基于LMM的无遗传关联标准检验统计量的方差被低估。这会导致测试统计中出现通货膨胀，反之亦然。

为了证明这一点，我们根据波多黎各人和非波多黎各人民的基因型差异，将基因组中所有常见SNP分为三类：（1）具有预期差异比率的SNP，即Var（SNP）=2M（M）A类F类(1 − M（M）A类F类)波多黎各人（高危人群）比非波多黎各（低危人群）低0.8；（2）波多黎各人与非波多黎各人为变量（SNP）之比在0.8至1.25之间的SNP；（3）波多黎各人与非波多黎各人为Var（SNP）比值大于1.25的SNP。在这一分类中，1类SNPs（n=144815[11%]）在波多黎各人（高风险）中的MAF明显低于非波多黎各人（低风险）。2类SNP在两组中具有相似的MAF（n=982805[76%]）。3类SNP（n=171601[13%]）在波多黎各人（高风险）中的MAF明显高于非波多黎各（低风险）。LMM和GMMAT在2类SNP中表现良好(图1C） ●●●●。LMM的问题对于1类和3类SNP更为明显：图1B显示了1类SNP的保守p值，以及图1D显示了预期的3类SNP的反保守p值。这表明，对于HCHS/SOL哮喘GWAS中至少24%的SNP，LMM计算的p值可能不正确。在删除两组相关个体后，通过分别分析波多黎各人和非波多黎各人民进行的荟萃分析提高了LMM的表现，但结果仍然没有很好地针对1类和3类SNP进行校准(图1B和1D），可能是由于组内人口分层引起的异方差问题。相比之下，GMMAT在所考虑的所有SNP类别中表现良好。

模拟研究

我们还在零假设下进行了广泛的模拟研究，以在各种GWAS病例对照和队列研究设计中说明这一问题，从而排除了低分子肝素在HCHS/SOL中的膨胀结果是由多基因效应引起的可能性。我们将分析局限于MAF大于5%的常见SNP。首先，在一项对10000名具有隐性关联的个体进行的模拟队列研究中，我们从空间连续人口图中模拟了疾病流行率分别为28%（人口1）和5%（人口2）的两个人口组。选择疾病流行率模拟HCHS/SOL中波多黎各和非波多黎各人的哮喘疾病比例(表1). 在LMM和GMMAT中，我们调整了前十名祖先PC作为固定效应协变量。图3A包括3200个零模拟数据集的结果，每个数据集有625583个SNP和10000个相关受试者。LMM的I类错误率为1.26×10⁻⁷标称5×10⁻⁸水平，与5.0×10相比⁻⁸用于GMMAT。请注意，由于包括了20多亿个p值，因此该水平的I类错误率估计值是准确的，标准误差为5×10⁻⁹.

在单独的窗口中打开

图3

10000名相关个体的模拟队列研究

在无遗传关联的零假设下，将来自3200个模拟重复的关联测试p值的分位数-分位数图与625583个常见SNP进行组合，得到超过20亿个零p值。

（A）所有SNP。

（B）第1类：人群1（高风险）与人群2（低风险）的预期方差比小于0.8的SNP。

（C）第2类：人群1（高风险）与人群2（低风险）的预期方差比率在0.8和1.25之间的SNP。

（D）第3类：人群1（高风险）与人群2（低风险）的预期方差比大于1.25的SNP。

以HCHS/SOL为例，我们将所有常见SNP分为三类：（1）群体1（高风险）与群体2（低风险）的预期方差之比小于0.8的SNP；（2）群体1的期望方差与群体2的期望方差之比在0.8和1.25之间的SNP；（3）群体1的期望方差与群体2的期望方差之比大于1.25的SNP。第1类包括在种群1中MAF显著低于种群2（24%）的SNP。第二类包括两个群体中具有可比MAF的SNP（58%）。第3类包括群体1中MAF明显高于群体2（18%）的SNPs。对于2类SNP，LMM和GMMAT都能正确控制I类错误率(图3C） ●●●●。然而，尽管对祖先个人计算机进行了调整，LMM还是降低了第一类SNP的I型错误率(图3B）以及夸大了第3类SNP的I类错误率(图3D） ●●●●。正如我们由HCHS/SOL研究激发的模拟研究所证明的那样，当不同人群的患病率差异导致不同人群之间的二元特征方差存在较大差异时，LMM往往无法控制中等或强人口分层情况下的I型错误率。然而，在使用GWAS中所有SNP的整个QQ图中，这个问题并不总是很明显，因为典型GWAS中的SNP数量通常在10个范围内⁵和10⁶，置信区间相对较宽，通货膨胀和通缩很可能在整个QQ图中相互掩盖。

为了说明这一点，图S2显示了仅一个模拟数据集的结果，该数据集具有625583个常见SNP和10000名受试者，与图3。LMM中的p值有轻微膨胀迹象图S2A、但对于1类SNP，结果较为保守（n=151206[24%]）(图S2B）和3类SNP的防伪(图S2D）（n=111455[18%]）。这些结果表明，尽管通货膨胀和通缩在整个QQ图中可能并不明显，但LMM的p值在该分析中可能对至少42%的SNP不正确。

接下来，我们模拟了一项病例对照研究，总样本量为10000，病例对照比为1:5，来自与我们的模拟队列研究相同的两个人群，存在潜在相关性。与队列研究相比，病例对照研究的样本过多。人群1（高危）和人群2（低危）的患病率分别为4.5%和0.5%。图S3A表明，根据p值的QQ图，GMMAT工作良好。基因组控制通货膨胀因子表明，LMM和GMMAT都具有适当的中间p值。然而，标称水平为5×10时的I类错误率⁻⁸为1.75×10⁻⁷对于LMM和4.6×10⁻⁸分别用于GMMAT。此外，当我们如前所述将SNP分为三类时，我们观察到第1类SNP出现了强烈的通货紧缩(图S3B）以及第3类SNP的强劲通货膨胀(图S3D）如模拟队列研究所示，从LMM分析中得出。额外的模拟研究表明，在存在人口分层的情况下，在没有隐性关联的队列研究中，LMM也存在同样的问题(图S4). 尽管通过LMM和荟萃分析分别分析每个群体，可以适当解释离散群体的混淆(图S5)，这种策略在存在连续种群结构的混杂情况下表现不佳(图S6). 此外，LMM在对具有中度到重度人口分层的无关个体进行的病例对照研究中效果不佳(图S7). 这些表明，尽管LMM在基于人群和基于家庭的遗传关联研究中都得到了广泛应用，但由于模型中的表型变异指定错误，LMM通常不适用于二元性状，可能产生错误的p值（参见附录B详细信息）。

我们还进行了额外的模拟来比较GMMAT和ROADTRIPS，后者对病例对照数据中的关联性进行了回顾性测试，³³以及最近开发的作为表型（LEAP）方法的责任估计器。LEAP适用于解释病例对照确定的责任混合模型。²⁰在与所述相同的病例对照设置中图S7，由于未能充分考虑人口分层，ROADTRIPS没有得到很好的校准。研究发现，LEAP的QQ图总体表现良好，但第1类SNP的I型错误率过高，无论前十大祖先PC是否作为协变量进行了调整(图S8).

接下来，我们通过改变两个人群组的病例对照抽样比例，模拟了疾病流行率为1%的病例对照研究。我们首先考虑了两个人群中的平衡病例和控制。具体来说，当两组病例对照比均为1:1时，没有人群分层，ROADTRIPS、LEAP、LMM和GMMAT都能适当控制I型错误率(图S9). 当病例对照比不同但在两个人群中翻转时（人群1为4:1，人群2为1:4），存在人群分层。然而，在这种情况下，由于二元性状的方差在两个种群中是相同的，LMM和GMMAT都表现良好，然而，ROADTRIPS在尾部夸大了I型错误率，LEAP在总体QQ图上有适当的结果，但在两个群体中MAF不同的SNP中夸大或缩小了I型错率(图S10).

我们还考虑了两个群体的病例对照比率不同的情况（群体1中为25:2，群体2中为25:48），这导致群体1中二元性状的方差小于群体2。ROADTRIPS和LEAP对前十名PC进行了调整，结果显示QQ图总体上出现了通货膨胀，没有协变量的LEAP和LMM使第1类SNP的I型错误率膨胀，而第3类SNP则使I型错误比率膨胀(图S11). 就所需的计算资源而言，与GMMAT相比，LEAP需要20倍以上的内存来分析样本量为10000的病例对照研究。

我们还进行了模拟研究，以评估GMMAT在估计遗传变异优势比方面的性能。我们发现，与大样本（n=10000）中的真实参数值相比，在零假设和替代假设下拟合logistic混合模型时，GMMAT在估计优势比方面存在微小偏差。此外，就绝对偏差而言，性能似乎与使用MACAU相似，³⁴最近开发的基于马尔可夫链蒙特卡罗的二项式计数数据混合模型方法，而GMMAT比值比估计的变量较小(图S12).

在没有人口分层的情况下，LMM和GMMAT具有可比较的功效，但它们的功效都不如logistic回归(表S1). 在存在人群分层的情况下，LMM对高风险人群中MAF低于低风险人群的因果遗传变异的作用小于GMMAT(表S1)，由于此类变体的保守型I错误率控制，如所示图S4B。

计算速度和内存使用

我们将GMMAT包与SAS PROC GLIMMIX在计算速度和内存使用方面进行了基准测试。为了拟合样本量为2000的零模型，当用一个方差分量拟合一个逻辑混合模型时，GMMAT所需时间少于SAS PROC GLIMMIX所需时间的1.5%，而用三个方差分量所需时间则少于0.6%(表S2)，得出的数值结果至少保留到小数点后四位。平均而言，对于随机效应的一个方差分量，通常用于解释GWAS中的遗传相关性，SAS PROC GLIMMIX需要大约28分钟才能在Intel Xeon E5-2690 CPU（2.90 GHz）的单核上拟合零模型，而使用GMMAT大约需要22秒。除了遗传相关性外，SAS PROC GLIMMIX还需要约1.2小时的时间，而GMMAT需要约22秒的时间，因为随机效应的三个方差分量说明了复杂的抽样设计，GMMAT对同一核心上的1000000个遗传变异进行得分测试需要约14分钟，而无需并行化。实际上，不同变体的分数测试可以很容易地在一个计算集群中并行运行，测试每个遗传变体只需要大约0.8毫秒。

GMMAT在分析2000个人时需要不到1GB的内存。样本量为10000，使用GMMAT对一个方差分量的零模型进行拟合大约需要18分钟，对三个方差分量进行拟合需要34分钟，对单个核心上的1000000个遗传变异进行得分测试大约需要3.6小时（测试每个遗传变异大约需要13毫秒）(表S2). SAS PROC GLIMMIX报告，当指定100 GB内存来适应10000个个体的空模型时，内存不足，而GMMAT要求一个方差分量少于14 GB，三个方差分量小于21 GB。

讨论

在分析GWAS中的二元性状时，我们建议使用logistic混合模型来校正人口分层和相关性。提出的GMMAT对具有二元性状的队列和病例对照GWAS中的遗传关联进行了计算效率高的分数测试。我们证明了GMMAT在控制I型错误率方面是有效的。相反，即使将祖先项作为协变量包含在内，在存在人口分层的情况下，将LMM应用于二元性状可能会导致不正确的I型错误率，特别是当群体具有异质性疾病风险或病例对照比，导致不同的二元性状方差时。在这种情况下，LMM只有在被测试关联的遗传变异体的MAF在所有人群中大致相同时，即在人口结构没有或弱混杂的情况下才近似有效。标准LMM下的同方差假设很重要，但在以前使用LMM进行二元性状的遗传关联研究中，这一假设在很大程度上被忽视了。尽管人们普遍认为LMM可以用来解释连续性状和二元性状的群体分层，但我们的结果表明，LMM在存在群体分层和关联性的情况下会导致不正确的I型错误率和p值。最近，Conomos等人。²⁵结果表明，连续性状也可能违反同方差假设，并且在这种情况下，考虑群体特定残差的模型优于标准LMM。对于离散群体中的二元性状，当每个群体中的群体结构没有混杂或混杂较弱时，我们已经表明，元分析是一种有效的方法，可以将使用LMM分别分析每个群体的结果结合起来(图S5). 此外，我们可以使用逻辑回归分析同质无关样本，这比LMM和GMMAT更有效(表S1). 然而，实际上，在大规模的遗传关联研究中，尤其是混合群体中，没有组内混杂或组内混杂较弱的群体可能并不总是被明确定义。当组内人群分层仍然存在时，我们还表明，使用LMM的元分析方法与标准LMM方法类似，仍然存在校准错误(图S6). 此外，通过比较特定群体的二元特征方差和MAF来检验LMM的有效性通常是困难的。相反，GMMAT在拟合logistic混合模型时不需要进行这种检查。

实际上，QQ图已广泛用于GWAS中的模型诊断。我们的结果表明，良好的QQ图不足以识别由于模型错误指定而导致的无效结果。具体来说，QQ图可以显示LMM正确控制二进制特征的I型错误率(图S2A），但这可能只是反映了不同SNP子集（1类和3类）的保守p值和非保守p值之间的平衡(图S2B和S2D）。需要注意的是，这些错误（不正确的p值）不会“抵消”；虚假的重要结果，即i类错误，不会通过省略真实信号（即II类错误）而“纠正”。更普遍地说，仅仅查看所有SNP的QQ图可能会导致错误的信念，即分析（例如，对二进制性状使用标准LMM）是合适的，而这些结果在很大程度上可能存在严重缺陷。

在对病例进行过抽样的病例对照研究中，LMM的I型错误率不正确，这是由于不同人群抽样方案中的病例对照比率不相等导致的人口分层导致的二元性状方差不相等，即使所有亚群中的疾病流行率相同，且总病例对照比为1:1(图S11). 由于ROADTRIPS目前不允许进行协变量调整，因此它在存在中度到重度人口分层的情况下不起作用。³⁵最近开发的负债混合模型²⁰^,²¹需要准确估计潜在总体人群中的疾病流行率、遗传率和易感性，这在实际存在人口分层和未知亚群的情况下可能很难获得。这些模型通常适用于没有或弱人口分层的情况，例如，在不同人口之间的病例对照比率相同的情况下，但在中等或强人口分层情况下，例如，当不同人口之间病例对照比率不同时，可能无法控制I型错误率。除了那些解释基因相关性的因素外，它们目前也无法处理多种随机效应，例如六氯环己烷/甲醛中的家庭和群体。责任阈值混合线性模型（LTMLM）方法仅适用于基于人群的病例对照研究设计，没有混淆因素，相关性低，因为它不能调整协变量或处理家庭数据。²¹LEAP方法可以容纳协变量并适用于家庭数据，但包含协变量带来了技术和统计挑战。²⁰GMMAT提供了一种灵活的方法，不需要知道疾病流行率或遗传力和责任估计，并提供有效的p值，同时适当控制I型错误率。

最近，Song等人。³⁶提出了一种基因型条件关联检验，用于解释关联检验中的群体结构。然而，正如作者所指出的，他们的方法没有考虑到家庭或隐秘的关系。相比之下，我们的logistic混合模型方法更灵活，可以在基于人群和基于家庭的队列和病例对照研究中，以及复杂的抽样设计中，同时考虑人口结构和相关性（如我们的HCHS/SOL哮喘示例所示）。因此，除了人口分层之外，我们的方法还可以应用于更广泛的遗传关联研究，包括家庭数据、隐性关联性、未观察到的共享环境影响和非随机抽样研究设计，而无需以不同的方式对其建模。

我们提供了一个开源的R包GMMAT，用于拟合逻辑混合模型并在GWAS中执行基于分数的测试。在广义线性模型的一般框架中，该软件包还可以应用于其他类型的连续和离散性状³⁷允许不同的链接功能和不同的均值-方差关系。此外，从同一疾病的不同研究中获得的得分统计数据可以很容易地合并到荟萃分析中。³⁸

最近提出的MACAU算法³⁴使用基于马尔可夫链蒙特卡罗的方法，在二项式计数数据的混合模型中实现Wald测试。然而，通过logistic混合模型对大量单个变量进行Wald和似然比测试，目前对于中大样本GWAS来说在计算上是不切实际的，需要进一步研究开发高效算法。与线性混合模型不同，logistic混合模型可以通过在替代假设下拟合模型来直接用于估计优势比。通过为感兴趣的候选遗传变异体子集拟合替代logistic混合模型来估计优势比在计算上是可行的。如果未来可以解决计算问题，逻辑混合模型也可以用于GWAS中的风险预测。

致谢

这项工作得到了赠款P01 CA134294、R35 CA197449和R01 HL113338（给H.C.、Z.L.和X.L.）、K99 HL130593（给H.C）和R37 CA076404（给C.W.、Z.L和X.L）的支持。Z.L.还获得了中国奖学金委员会的部分奖学金。我们感谢西班牙裔社区健康研究/拉丁裔研究（HCHS/SOL）的参与者和工作人员对本研究的贡献。HCHS/SOL的基线检查是一项合作研究，由国家心脏、肺和血液研究所（NHLBI）与北卡罗来纳大学（N01-HC65233）、迈阿密大学（N01-HC65234）、阿尔伯特·爱因斯坦医学院（N01-HC65235）、西北大学（N01-HC65236）签订的合同支持，和圣地亚哥州立大学（N01-HC65237）。以下研究所、中心和办公室通过向NHLBI转移资金，为HCHS/SOL的第一阶段做出了贡献：国家少数民族健康和健康差异研究所、国家耳聋和其他沟通障碍研究所、全国牙科和颅面研究所（NIDCR）、，国家糖尿病、消化和肾脏疾病研究所、国家神经疾病和中风研究所以及NIH膳食补充剂办公室。华盛顿大学遗传分析中心得到了NHLBI和NIDCR合同（HHSN268201300005C AM03和MOD03）的支持。基因分型工作得到了NHLBI（HSN26220/20054C）、国家转化科学促进中心临床转化科学研究所（UL1TR000124）和NIDDK糖尿病研究中心（DK063491）的支持。HCHS/SOL出版物委员会已对本手稿进行了审查，以确定其科学内容以及数据解释与之前HCHS/SSOL出版物的一致性。

笔记

发布日期：2016年3月24日

脚注

补充数据包括12张图和2张表，可在以下网站上找到http://dx.doi.org/10.1016/j.ajhg.2016.02.012.

附录A：GMMAT的推导

广义线性混合模型

以下推导基于广义线性混合模型（GLMMs），当链接函数为logit且分散参数固定为1时，logistic混合模型是GLMMs的特例。在单变量测试的背景下，我们考虑以下GLMM

η_我= 克(μ_我) = X（X）_我α+G公司_我β+ b条_我,

哪里X（X）_我是1×第页主题协变量的行向量我,α是一个第页×1固定协变量效应的列向量，包括截距，G公司_我是受试者感兴趣的遗传变异的基因型我、和β是固定的基因型效应。我们假设 $b条 \sim N个 (0, \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})$ 是一个n个×1随机效应列向量，τ_k个是方差分量参数，τ是一个K（K）×1列向量τ_k个、和V（V）_k个已知的n个×n个矩阵。我们还假设，给定随机效应b条，结果年_我与平均值有条件独立E类(年_我|b条) = μ_我和方差 $V（V）一第页 (年_{我} | b条) = ϕ 一_{我}^{- 1} v（v） (μ_{我})$ ，其中ϕ是分散参数（对于二进制和泊松数据ϕ= 1),一_我是已知重量，以及v（v）（‧）是方差函数。线性预测器η_我是条件平均值的单调函数μ_我通过链接功能η_我= 克(μ_我). 对于二进制特征年_我,μ_我= π_我= P（P）(年_我= 1|X（X）_我, G公司_我, b条_我)是受试者二元结果（如疾病状态）的概率我.

对于主题我，给定随机效应的拟似然b条是

q个 我_{我} (α, β; b条) = \int_{年_{我}}^{μ_{我}} \frac{一_{我} (年_{我} - μ)}{ϕ v（v） (μ)} d日 μ .

对数积分拟似然函数(α, β, ϕ, τ)是

q个 我 (α, β, ϕ, τ) = 日志 \int 经验 {\sum_{我 = 1}^{n个} q个 我_{我} (α, β; b条)} {(2 π)}^{- \frac{n个}{2}} {| \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个} |}^{- \frac{1}{2}} \times 经验 {- \frac{1}{2} {b条}^{T型} {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} b条} d日 b条 .

（方程式A1）

让

（f） (b条) = \sum_{我 = 1}^{n个} q个 我_{我} (α, β; b条) - \frac{1}{2} {b条}^{T型} {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} b条,

我们可以使用拉普拉斯方法来近似n个-量纲积分

\int 经验 {（f） (b条)} d日 b条 \approx {(2 π)}^{\frac{n个}{2}} {| - {（f）}^{″} (\tilde{b条}) |}^{- \frac{1}{2}} 经验 {（f） (\tilde{b条})},

因此方程式A1变为

q个 我 (α, β, ϕ, τ) = - \frac{1}{2} 日志 | \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个} | - \frac{1}{2} 日志 | - {（f）}^{″} (\tilde{b条}) | + （f） (\tilde{b条}),

（方程式A2）

哪里

\tilde{b条} = \underset{b条}{参数 最大值} （f） (b条)

是的解决方案（f）^′(b条) = 0.

的一阶偏导数q个我_我(α, β; b条)关于b条是

\frac{\partial q个 我_{我}}{\partial b条} = \frac{\partial q个 我_{我}}{\partial μ_{我}} \frac{\partial μ_{我}}{\partial η_{我}} \frac{\partial η_{我}}{\partial b条} = \frac{一_{我} (年_{我} - μ_{我})}{ϕ v（v） (μ_{我})} \frac{1}{克^{'} (μ_{我})} {Z轴}_{我}^{T型},

哪里Z轴_我是1×n个指标向量，以便b条_我= Z轴_我b条, $我_{n个} = (\begin{matrix} {Z轴}_{1}^{T型} & {Z轴}_{2}^{T型} & \dots & {Z轴}_{n个}^{T型} \end{matrix})$ ，二阶导数为

\frac{\partial^{2} q个 我_{我}}{\partial b条 \partial {b条}^{T型}} = - \frac{一_{我} {Z轴}_{我}^{T型} {Z轴}_{我}}{ϕ v（v） (μ_{我}) {[克^{'} (μ_{我})]}^{2}} - \frac{一_{我} (年_{我} - μ_{我}) {v（v）}^{'} (μ_{我}) {Z轴}_{我}^{T型} {Z轴}_{我}}{ϕ {v（v）}^{2} (μ_{我}) {[克^{'} (μ_{我})]}^{2}} - \frac{一_{我} (年_{我} - μ_{我}) 克^{″} (μ_{我}) {Z轴}_{我}^{T型} {Z轴}_{我}}{ϕ v（v） (μ_{我}) {[克^{'} (μ_{我})]}^{三}} .

对于规范链接函数，v（v）(μ_我)克^′(μ_我)=1，后两项变为0。让

Δ= d日我一克{克^′(μ_我)},

W公司 = d日 我 一 克 {\frac{一_{我}}{ϕ v（v） (μ_{我}) {[克^{'} (μ_{我})]}^{2}}},

然后是方程式A2类成为

\begin{matrix} q个 我 (α, β, ϕ, τ) = - \frac{1}{2} 日志 | \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个} | - \frac{1}{2} 日志 | \sum_{我 = 1}^{n个} \frac{一_{我} {Z轴}_{我}^{T型} {Z轴}_{我}}{ϕ v（v） (μ_{我}) {[克^{'} (μ_{我})]}^{2}} + {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} | + \sum_{我 = 1}^{n个} q个 我_{我} (α, β; \tilde{b条}) - \frac{1}{2} {\tilde{b条}}^{T型} {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} \tilde{b条} \\ = - \frac{1}{2} 日志 | \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个} W公司 + 我 | + \sum_{我 = 1}^{n个} q个 我_{我} (α, β; \tilde{b条}) - \frac{1}{2} {\tilde{b条}}^{T型} {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} \tilde{b条} . \end{matrix}

（方程式A3）

我们假设权重矩阵W公司相对于条件平均值变化缓慢（继布雷斯洛和克莱顿之后²²)，这是

\frac{\partial W公司}{\partial μ_{我}} \approx 0,

然后我们求方程的导数A3号:

\frac{\partial q个 我 (α, β, ϕ, τ)}{\partial α} = \sum_{我 = 1}^{n个} \frac{一_{我} (年_{我} - μ_{我})}{ϕ v（v） (μ_{我})} \frac{1}{克^{'} (μ_{我})} {X（X）}_{我}^{T型} = {X（X）}^{T型} W公司 Δ (年 - μ),

\frac{\partial q个 我 (α, β, ϕ, τ)}{\partial β} = \sum_{我 = 1}^{n个} \frac{一_{我} (年_{我} - μ_{我})}{ϕ v（v） (μ_{我})} \frac{1}{克^{'} (μ_{我})} {G公司}_{我} = {G公司}^{T型} W公司 Δ (年 - μ),

（方程式A4）

\begin{matrix} \frac{\partial q个 我 (α, β, ϕ, τ)}{\partial b条} = \sum_{我 = 1}^{n个} \frac{一_{我} (年_{我} - μ_{我})}{ϕ v（v） (μ_{我})} \frac{1}{克^{'} (μ_{我})} {Z轴}_{我}^{T型} - {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} b条 \\ = W公司 Δ (年 - μ) - {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} b条 . \end{matrix}

在零假设下H（H）₀:β=0，如果ϕ和τ众所周知，我们共同选择 $\hat{α} (ϕ, τ)$ 和 $\hat{b条} (ϕ, τ)$ 最大化方程式A3号，然后 $\overset{图6}{b条} (ϕ, τ) = \tilde{b条} (\overset{图6}{α} (ϕ, τ), β = 0)$ 因为 $\tilde{b条}$ 最大化（f）(b条)对于给定的(α, β). 定义工作向量 $\tilde{年}$ 包含元素 ${\tilde{年}}_{我} = η_{我} + 克^{'} (μ_{我}) (年_{我} - μ_{我})$ ，解决方案

{\begin{matrix} {X（X）}^{T型} W公司 Δ (年 - μ) = 0 \\ W公司 Δ (年 - μ) = {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} b条 \end{matrix}

可以写为系统的解决方案

[\begin{matrix} {X（X）}^{T型} W公司 X（X） & {X（X）}^{T型} W公司 \\ W公司 X（X） & {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} + W公司 \end{matrix}] [\begin{matrix} α \\ b条 \end{matrix}] = [\begin{matrix} {X（X）}^{T型} W公司 \tilde{年} \\ W公司 \tilde{年} \end{matrix}] .

让 $Σ = {W公司}^{- 1} + \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个}$ ,P（P）= Σ⁻¹ − ∑⁻¹X（X）(X（X）^T型Σ⁻¹X（X）)⁻¹X（X）^T型Σ⁻¹，然后

{\begin{matrix} \hat{α} = {({X（X）}^{T型} Σ^{- 1} X（X）)}^{- 1} {X（X）}^{T型} Σ^{- 1} \tilde{年} \\ \hat{b条} = (\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个}) Σ^{- 1} (\tilde{年} - X（X） \overset{图6}{α}) \end{matrix}

是使方程最大化的解

A3号。我们注意到

\tilde{年} - \hat{η} = \tilde{年} - X（X） \hat{α} - \hat{b条} = {我 - (\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个}) Σ^{- 1}} (\tilde{年} - X（X） \hat{α}) = {W公司}^{- 1} Σ^{- 1} (\tilde{年} - X（X） \hat{α}) = {W公司}^{- 1} P（P） \tilde{年} .

方差分量参数的估计

继布雷斯洛和克莱顿之后，²²我们忽视了W公司在τ并使用Pearson chi-square统计来近似偏差

- 2 ϕ \sum_{我 = 1}^{n个} q个 我_{我} (α, β; b条) = \sum_{我 = 1}^{n个} - 2 \int_{年_{我}}^{μ_{我}} \frac{一_{我} (年_{我} - μ)}{v（v） (μ)} d日 μ \approx \sum_{我 = 1}^{n个} \frac{一_{我} {(年_{我} - μ_{我})}^{2}}{v（v） (μ_{我})} .

然后是方程式A3号最大值变为

\begin{matrix} q个 我 (\hat{α} (ϕ, τ), β = 0, ϕ, τ) \approx - \frac{1}{2} 日志 | \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个} W公司 + 我 | - \frac{1}{2} \sum_{我 = 1}^{n个} \frac{一_{我} {(年_{我} - {\overset{图6}{μ}}_{我})}^{2}}{ϕ v（v） ({\hat{μ}}_{我})} - \frac{1}{2} {\hat{b条}}^{T型} {(\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个})}^{- 1} \hat{b条} \\ = - \frac{1}{2} 日志 | \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个} W公司 + 我 | - \frac{1}{2} {(年 - \hat{μ})}^{T型} Δ W公司 Δ (年 - \hat{μ}) - \frac{1}{2} {(\tilde{年} - X（X） \hat{α})}^{T型} Σ^{- 1} (\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个}) Σ^{- 1} (\tilde{年} - X（X） \hat{α}) \\ = - \frac{1}{2} 日志 | Σ W公司 | - \frac{1}{2} {(\tilde{年} - \hat{η})}^{T型} W公司 (\tilde{年} - \hat{η}) - \frac{1}{2} {(\tilde{年} - X（X） \overset{图6}{α})}^{T型} Σ^{- 1} (\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个}) Σ^{- 1} (\tilde{年} - X（X） \hat{α}) \\ = - \frac{1}{2} 日志 | W公司 | - \frac{1}{2} 日志 | Σ | - \frac{1}{2} {\tilde{年}}^{T型} P（P） {W公司}^{- 1} P（P） \tilde{年} - \frac{1}{2} {\tilde{年}}^{T型} P（P） (\sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个}) P（P） \tilde{年} \\ = c - \frac{1}{2} 日志 | Σ | - \frac{1}{2} {\tilde{年}}^{T型} P（P） \tilde{年} . \end{matrix}

类似地，限制最大似然（REML）版本为

q个 我_{R（右）} (\hat{α} (ϕ, τ), β = 0, ϕ, τ) = c_{R（右）} - \frac{1}{2} 日志 | \sum | - \frac{1}{2} 日志 | {X（X）}^{T型} Σ^{- 1} X（X） | - \frac{1}{2} {\tilde{年}}^{T型} P（P） \tilde{年} .

让 ${V（V）}_{0} = d日我一克 {一_{我}^{- 1} v（v） (μ_{我}) {[克^{'} (μ_{我})]}^{2}} = ϕ^{- 1} {W公司}^{- 1}$ ，然后 $Σ = ϕ {V（V）}_{0} + \sum_{k个 = 1}^{K（K）} τ_{k个} {V（V）}_{k个}$ 和的一阶导数 $q个我_{R（右）} (\hat{α} (ϕ, τ), β = 0, ϕ, τ)$ 关于ϕ和τ是

\frac{\partial q个 我_{R（右）} (\hat{α} (ϕ, τ), β = 0, ϕ, τ)}{\partial ϕ} = \frac{1}{2} {\tilde{年} P（P） {V（V）}_{0} P（P） \tilde{年} - t吨 第页 (P（P） {V（V）}_{0})},

\frac{\partial q个 我_{R（右）} (\hat{α} (ϕ, τ), β = 0, ϕ, τ)}{\partial τ_{k个}} = \frac{1}{2} {\tilde{年} P（P） {V（V）}_{k个} P（P） \tilde{年} - t吨 第页 (P（P） {V（V）}_{k个})} .

我们定义平均信息⁶^,²³矩阵人工智能包含以下条目

A类 我_{ϕ ϕ} = \frac{1}{2} \tilde{年} P（P） {V（V）}_{0} P（P） {V（V）}_{0} P（P） \tilde{年},

A类 我_{ϕ τ_{k个}} = \frac{1}{2} \tilde{年} P（P） {V（V）}_{0} P（P） {V（V）}_{k个} P（P） \tilde{年},

A类 我_{τ_{k个} τ_{我}} = \frac{1}{2} \tilde{年} P（P） {V（V）}_{k个} P（P） {V（V）}_{我} P（P） \tilde{年} .

让θ是要估计的方差分量参数，当ϕ ≠ 1,θ= (ϕ, τ)、和人工智能是一个(K（K）+ 1)×(K（K）+1）矩阵。对于二进制和泊松数据，ϕ= 1,θ=========================================================τ、和人工智能是一个K（K）×K（K）矩阵仅包含A类我_{τ_k个τ_我}.

我们使用以下算法拟合空GLMM：

1
用拟合广义线性模型τ=0并获得 ${\hat{α}}^{(0)}$ 和工作矢量 ${\tilde{年}}^{(0)}$ ;
2
使用 $θ^{(0)} = V（V）一第页 ({\tilde{年}}^{(0)}) / K（K）$ （如果ϕ=1）或 $θ^{(0)} = V（V）一第页 ({\tilde{年}}^{(0)}) / (K（K） + 1)$ （如果ϕ≠1）作为的初始值θ;
三。
对于每个 $k个 = 0,1, \dots, K（K）$ ，更新θ使用 $θ_{k个}^{(1)} = θ_{k个}^{(0)} + 2 {n个}^{- 1} {θ_{k个}^{(0)}}^{2} (\partial q个我_{R（右）} (θ^{(0)}) / \partial θ_{k个})$ ;
4
使用 ${\tilde{年}}^{(1)} = {\tilde{年}}^{(0)}$ 作为 $\tilde{年}$ 和更新θ⁽²⁾= θ⁽¹⁾+ {A类我⁽¹⁾}⁻¹(∂q个我_R（右）(θ⁽¹⁾)/∂θ)；
5
计算 ${\hat{α}}^{(2)}$ 和 ${\overset{图6}{b条}}^{(2)}$ 使用 ${\tilde{年}}^{(1)}$ 和θ⁽²⁾;
6
更新 ${\tilde{年}}^{(2)}$ 使用 ${\hat{α}}^{(2)}$ 和 ${\hat{b条}}^{(2)}$ ;
7
重复步骤4-6，直到 $2 最大值 {| {\hat{α}}^{(我)} - {\hat{α}}^{(我 - 1)} | / (| {\hat{α}}^{(我)} | + | {\hat{α}}^{(我 - 1)} |), | {\hat{θ}}^{(我)} - {\hat{θ}}^{(我 - 1)} | / (| {\hat{θ}}^{(我)} | + | {\hat{θ}}^{(我 - 1)} |)} \leq$ 公差。

分数测试

一次 $(\hat{α}, \hat{ϕ}, \hat{τ})$ 是在零假设下估计的H（H）₀:β=0，分数测试可以通过评估方程来构建A4（A4）在 $(\hat{α}, β = 0, \hat{ϕ}, \hat{τ})$ ，即

T型 = \frac{\partial q个 我 (\hat{α}, β = 0, \hat{ϕ}, \hat{τ})}{\partial β} = {G公司}^{T型} \hat{W公司} \hat{Δ} (年 - \hat{μ}) = {G公司}^{T型} \hat{W公司} (\tilde{年} - \hat{η}) = {G公司}^{T型} \overset{图6}{P（P）} \tilde{年} .

其在零假设下的方差为

V（V） 一 第页 (T型 | {H（H）}_{0}) = E类 {\frac{\partial q个 我 (\hat{α}, β = 0, \hat{ϕ}, \hat{τ})}{\partial β} \frac{\partial q个 我 (\hat{α}, β = 0, \hat{ϕ}, \hat{τ})}{\partial β^{T型}}} = E类 ({G公司}^{T型} \overset{图6}{P（P）} \tilde{年} {\tilde{年}}^{T型} \hat{P（P）} G公司) = {G公司}^{T型} \hat{P（P）} G公司,

最后一个等式成立是因为 $\hat{P（P）} \hat{Σ} \hat{P（P）} = \hat{P（P）}$ .

附录B其他模拟研究

人口分层的无关个体

我们进行了额外的模拟研究，以比较存在人口分层的无关个体的LMM和GMMAT。我们使用了聚合模型²⁹从20×20网格中模拟16000个无亲缘关系个体和1000000个遗传变异的基因型(图S1)具有相邻细胞间迁移率的空间连续种群（每个细胞40个个体）M（M）=10以模拟欧洲的人口结构。³⁰^,³¹对于个人我，成为病例的概率π_我计算依据为

罗吉特(π_我) = α₀+ α₁Z轴_我,

哪里Z轴_我=1（如果是个人）我来自左上角的10×10网格（人口1）Z轴_我否则=0（总体2）。参数α₀和α₁选择1组（高危）和2组（低危）的患病率分别为0.28和0.05。请注意，当使用线性链接函数时，平均值模型没有被误定，因此我们的模拟设置不支持逻辑模型。我们随机抽取了10000名受试者，对常见疾病进行了模拟队列研究(图S4)并使用MAF大于5%的625504个遗传变异计算其标准化遗传关系矩阵。我们调整了两种模型中排名前十的PC，并将3200个零模拟副本的结果合并。我们将我们的方法与使用LMM分别分析人群1和人群2的策略进行了比较，然后将结果合并到元分析中（我们称之为LMM元分析），使用100个零模拟重复(图S5). 此外，我们考虑了连续的人口分层，其中Z轴_我对于每个单元格，是其行和列坐标的最小值图S1，范围从0到19，以及α₀被选中的人群中的疾病流行率Z轴_我=0为0.02α₁= 0.2 (图S6).

我们还模拟了一项罕见疾病的病例对照研究(图S7). 我们选择了α₀和α₁因此，1号人群（高危人群）的患病率为0.045，2号人群（低危人群）的发病率为0.005。我们随机抽取1667名病例受试者和8333名对照受试者进行病例对照研究，总样本量为10000，病例对照比为1:5。我们进行了与上述相同的分析，并结合了3200个空模拟副本的结果。

与现有方法的比较

我们还比较了GMMAT和ROADTRIPS³³以及最近开发的负债（probit）混合模型。²⁰^,²¹ROADTRIPS不考虑协变量。因为LTMLM²¹不允许协变量调整，我们使用LEAP²⁰并将疾病流行率设置为0.015，这是高风险和低风险人群的合并流行率。我们比较了无协变量的ROADTRIPS和LEAP，以及作为协变量的前十名PC的LEAP与LMM和GMMAT。模拟设置与中的相同图S7，但只显示了一个模拟复制的p值(图S8).

此外，我们模拟了三个病例对照环境，在两个人群组中使用相同的流行率但不同的病例对照抽样方案图S1两组的患病率均设为0.01。我们将ROADTRIPS和LEAP与LMM和GMMAT进行了比较，并对前十名PC进行了调整。

在第一组中，我们从人群1中随机抽取1250名病例受试者和1250名对照受试者，从人群2中随机抽取3750名病例受试者和3750名对照受试者。这是一项平衡的病例对照研究，在两组人群中进行平衡设计。由于病例对照比相同，因此没有人群分层(图S9).

在第二个设置中，我们从人口1中随机抽取2000名病例受试者和500名对照受试者，从人口2中随机抽取1500名病例受试者和6000名对照受试者。这是一项不平衡的病例对照研究，在两个群体中，设计不平衡，但二元性状方差相等。人口分层是由两个人口组中不同的病例对照比率（人口1为4:1，人口2为1:4）而不是不同的疾病流行率创建的(图S10).

在第三组中，我们从1组中随机抽取2500名病例受试者和200名对照受试者，从2组中随机选取2500名病例和4800名对照受试者。这是一项平衡的病例对照研究，在两个群体中具有不平衡的设计和不相等的二元性状方差。由于两组的患病率均为0.01，因此没有高风险或低风险组，但群体1为低二进制性状方差组，方差为0.0686，群体2为高二进制性状方差，方差为0.2252(图S11).

具有遗传效应的模拟

我们首先进行了模拟研究，以评估GMMAT在人口分层情况下估计遗传变异优势比的性能。我们使用的基因型数据与我们对具有人口分层的无关个体进行的零模拟中使用的数据相同。对于个人我，成为案例的概率π_我计算依据为

罗吉特(π_我) = α₀+ α₁Z轴_我+ βG公司_我,

哪里Z轴_我=1（如果是个人）我来自人口1，Z轴_我否则=0（总体2）。G公司_我是因果遗传变异的附加编码基因型。参数α₀和α₁选择的人群1的疾病流行率为0.28，人群2的患病率为0.05G公司_我=0。β选择的优势比为1.0、1.1、1.25、1.5至2.0。总样本量为10000。我们比较了GMMAT和MACAU计算的优势比估计值和p值(图S12). 然后，我们比较了LMM和GMMAT在识别人群1中MAF低于人群2的因果变量方面的功效(表S1).

我们还对1667名病例受试者和8333名无人群分层的对照受试者进行了病例对照研究。我们假设患有G公司_我=0，并在人群1和人群2中使用相同的抽样方案。我们比较了logistic回归、LMM和GMMAT在识别因果变量方面的功效(表S1).

Web资源

此处显示的数据URL如下：

犰狳，http://arma.sourceforge.net/
GMMAT公司，http://www.hsph.harvard.edu/xlin/software.html#gmmat
澳门，http://www.xzlab.org/software.html
理学硕士，http://home.uchicago.edu/rhudson1/source/mksamples.html
道路，http://faculty.washington.edu/athornt/software/ROADTRIPS2（网址：http://faculty.washington.edu/athornt/software/ROADTRIPS2）/

补充数据

文件S1。图S1–S12以及表S1和S2：

单击此处查看。^{（4.0M，pdf格式）}

文件S2。文章和补充数据：

单击此处查看。^{（460万，pdf）}

工具书类

1Lander E.S.、Schork N.J.复杂性状的遗传解剖。科学。1994;265:2037–2048.[公共医学][谷歌学者]

2Aulchenko Y.S.，de Koning D.J.，Haley C.使用混合模型和回归的全基因组快速关联：基于全基因组系谱的数量性状位点关联分析的快速简单方法。遗传学。2007;177:577–585. [PMC免费文章][公共医学][谷歌学者]

三。Kang H.M.、Zaitlen N.A.、Wade C.M.、Kirby A.、Heckerman D.、Daly M.J.、Eskin E.模型生物关联图中种群结构的有效控制。遗传学。2008;178:1709–1723. [PMC免费文章][公共医学][谷歌学者]

4Kang H.M.、Sul J.H.、Service S.K.、Zaitlen N.A.、Kong S.Y.、Freimer N.B.、Sabatti C.、Eskin E.用于解释全基因组关联研究中样本结构的方差分量模型。自然遗传学。2010;42:348–354. [PMC免费文章][公共医学][谷歌学者]

5Zhang Z.、Ersoz E.、Lai C.Q.、Todhunter R.J.、Tiwari H.K.、Gore M.A.、Bradbury P.J.、Yu J.、Arnett D.K.、Ordovas J.M.、Buckler E.S.适用于全基因组关联研究的混合线性模型方法。自然遗传学。2010;42:355–360. [PMC免费文章][公共医学][谷歌学者]

6Yang J.、Lee S.H.、Goddard M.E.、Visscher P.M.GCTA：全基因组复杂性状分析工具。Am.J.Hum.遗传学。2011;88:76–82. [PMC免费文章][公共医学][谷歌学者]

7Lippert C.、Listgarten J.、Liu Y.、Kadie C.M.、Davidson R.I.、Heckerman D.FaST全基因组关联研究线性混合模型。自然方法。2011;8:833–835.[公共医学][谷歌学者]

8Zhou X.，Stephens M.关联研究的全基因组高效混合模型分析。自然遗传学。2012;44:821–824. [PMC免费文章][公共医学][谷歌学者]

9Pirinen M.，Donnelly P.，Spencer C.C.A.在大规模数据集上使用线性混合模型进行高效计算，并应用于遗传研究。附录申请。斯达。2013;7:369–390. [谷歌学者]

10Loh P.R.、Tucker G.、Bulik-Sullivan B.K.、Vilhjálmsson B.J.、Finucane H.K.、Salem R.M.、Chasman D.I.、Ridker P.M.、Neale B.M.、Berger B.高效的贝叶斯混合模型分析增加了大型队列中的关联力。自然遗传学。2015;47：284–290。 [PMC免费文章][公共医学][谷歌学者]

11Barr R.G.、Avilés-Santa L.、Davis s.M.、Aldrich T.、Gonzalez F.、Henderson A.G.、Kaplan R.C.、LaVange L.、Liu K.、Loredo J.s.西班牙裔移民中的肺部疾病和年龄：西班牙人社区健康研究/拉丁裔研究（HCHS/SOL）的结果美国J.Respir。批评。护理医学。2016年；193:386–395. [PMC免费文章][公共医学][谷歌学者]

12Price A.L.、Patterson N.J.、Plenge R.M.、Weinblatt M.E.、Shadick N.A.、Reich D.主成分分析纠正了全基因组关联研究中的分层。自然遗传学。2006;38:904–909.[公共医学][谷歌学者]

13Tucker G.、Price A.L.和Berger B.提高GWAS的能力，避免人口分层与PC-Select的混淆。遗传学。2014;197:1045–1049. [PMC免费文章][公共医学][谷歌学者]

14Fingerlin T.E.、Murphy E.、Zhang W.、Peljto A.L.、Brown K.K.、Steele M.P.、Loyd J.E.、Cosgrove G.P.、Lynch D.、Groshong S.全基因组关联研究确定了肺纤维化的多个易感基因座。自然遗传学。2013;45:613–620. [PMC免费文章][公共医学][谷歌学者]

15Cortes A.、Hadler J.、Pointon J.P.、Robinson P.C.、Karaderi T.、Leo P.、Cremin K.、Pryce K.、Harris J.、Lee S.、国际强直性脊柱炎遗传学联合会（IGAS）澳大利亚-安哥拉-美国脊椎关节炎联合会（TASC）Française d’Etude Génétique des Spondylarthritites（GFEGS）Nord-Trondelag Health Study（HUNT）加拿大脊椎关节炎研究联合会（SPARCC）Wellcome Trust Case Control Consortium 2（WTCCC2）通过免疫相关基因座的高密度基因分型确定强直性脊柱炎的多种风险变体。自然遗传学。2013;45:730–738. [PMC免费文章][公共医学][谷歌学者]

16Fakiola M.、Strange A.、Cordell H.J.、Miller E.N.、Pirinen M.、Su Z.、Mishra A.、Mehrotra S.、Monteiro G.R.、Band G.、LeishGEN Consortium。Wellcome Trust Case Control Consortium 2 HLA-DRB1-HLA-DQA1 HLA-II类区域的常见变异与内脏利什曼病易感性相关。自然遗传学。2013;45:208–213. [PMC免费文章][公共医学][谷歌学者]

17Liu J.Z.、Hov J.R.、Folseraas T.、Ellinghaus E.、Rushbrook S.M.、Doncheva N.T.、Andreassen O.A.、Weersma R.K.、Weismüller T.J.、Eksteen B.、UK-PSCC财团。国际PSC研究小组。国际IBD遗传学联合会免疫相关疾病区域的密集型基因分型确定了9个新的原发性硬化性胆管炎风险位点。自然遗传学。2013;45:670–675. [PMC免费文章][公共医学][谷歌学者]

18Huber，P.J.（1967年）。非标准条件下极大似然估计的行为。第五届伯克利数理统计与概率研讨会论文集，第1卷：统计学（加州伯克利：加利福尼亚大学出版社），221–233。

19Jarque C.M.，Bera A.K.回归残差的正态性、同方差性和序列独立性的有效检验。经济。莱特。1980;6:255–259. [谷歌学者]

20Weissbrod O.，Lippert C.，Geiger D.，Heckerman D.在已确定的病例对照研究中，准确的责任估计可以提高能力。自然方法。2015;12:332–334.[公共医学][谷歌学者]

21Hayeck T.J.、Zaitlen N.A.、Loh P.R.、Vilhjalmsson B.、Pollack S.、Gusev A.、Yang J.、Chen G.B.、Goddard M.E.、Visscher P.M.带有病例对照确定校正的混合模型增加了关联力。Am.J.Hum.遗传学。2015;96:720–730. [PMC免费文章][公共医学][谷歌学者]

22Breslow N.E.，Clayton D.G.广义线性混合模型中的近似推断。美国期刊。统计协会。1993;88:9–25. [谷歌学者]

23Gilmour A.R.，Thompson R.，Cullis B.R.平均信息REML：线性混合模型方差参数估计的有效算法。生物计量学。1995;51:1440–1450. [谷歌学者]

24Rao C.R.关于几个参数的统计假设的大样本检验及其在估计问题中的应用。数学。程序。外倾角。菲洛斯。Soc公司。1948;44：50–57。 [谷歌学者]

25Conomos M.P.、Laurie C.A.、Stilp A.M.、Gogarten S.M.、McHugh C.P.、Nelson S.C.、Sofer T.、Fernández-Rowdes L.、Justice A.E.、Graff M.美国西班牙裔/拉丁裔人群的遗传多样性和关联研究：在西班牙人社区健康研究/拉丁族研究中的应用。Am.J.Hum.遗传学。2016年；98:165–184. [PMC免费文章][公共医学][谷歌学者]

26Lavange L.M.、Kalsbeek W.D.、Sorlie P.D.、Avilés-Santa L.M.，Kaplan R.C.、Barnhart J.、Liu K.、Giachello A.、Lee D.J.、Ryan J.西班牙裔社区健康研究/拉丁裔研究中的样本设计和队列选择。流行病学年鉴。2010;20:642–649. [PMC免费文章][公共医学][谷歌学者]

27Laurie C.C.、Doheny K.F.、Mirel D.B.、Pugh E.W.、Bierut L.J.、Bhangale T.、Boehm F.、Caporaso N.E.、Cornelis M.C.、Edenberg H.J.、GENEVA调查员全基因组关联研究中基因型数据的质量控制和质量保证。基因。流行病。2010;34:591–602. [PMC免费文章][公共医学][谷歌学者]

28Pfeffermann D.复杂调查数据建模：为什么建模？为什么这是一个问题？我们如何处理它。调查。Methodol公司。2011;37:115–136. [谷歌学者]

29Hudson R.R.在Wright-Fisher中性遗传变异模型下生成样本。生物信息学。2002;18:337–338.[公共医学][谷歌学者]

30Mathieson I.，McVean G.空间结构种群中罕见和常见变异的差异混淆。自然遗传学。2012;44:243–246. [PMC免费文章][公共医学][谷歌学者]

31Wang C.，Zhan X.，Bragg-Gresham J.，Kang H.M.，Stambolian D.，Chew E.Y.，Branham K.E.，Heckenlavy J.，Fulton R.，Wilson R.K.，FUSION Study基于序列的关联研究的祖先估计和种群分层控制。自然遗传学。2014;46:409–415. [PMC免费文章][公共医学][谷歌学者]

32MacCluer J.W.、VandeBerg J.L.、Read B.、Ryder O.A.通过计算机模拟进行谱系分析。动物园生物。1986;5：147–160。 [谷歌学者]

33Thornton T.、McPeek M.S.ROADTRIPS：部分或完全未知人群和系谱结构的病例对照关联测试。Am.J.Hum.遗传学。2010;86:172–184. [PMC免费文章][公共医学][谷歌学者]

34Lea A.J.，Tung J.，Zhou X.亚硫酸盐测序数据中鉴定差异DNA甲基化的灵活、有效的二项式混合模型。公共科学图书馆-遗传学。2015;11：e1005650。 [PMC免费文章][公共医学][谷歌学者]

35Price A.L.、Zaitlen N.A.、Reich D.、Patterson N.全基因组关联研究中人口分层的新方法。Nat.Rev.基因。2010;11:459–463. [PMC免费文章][公共医学][谷歌学者]

36Song M.，Hao W.，Storey J.D.任意结构群体的遗传关联测试。自然遗传学。2015;47:550–554. [PMC免费文章][公共医学][谷歌学者]

37McCullagh P.，Nelder J.A.第二版。查普曼和霍尔/CRC；1989.广义线性模型。[谷歌学者]

38Lee S.、Abecasis G.R.、Boehnke M.、Lin X.Rare-variant关联分析：研究设计和统计检验。Am.J.Hum.遗传学。2014;95：5–23。 [PMC免费文章][公共医学][谷歌学者]

文章来自美国人类遗传学杂志由以下人员提供美国人类遗传学学会