多元表型與基因型的全基因組關(guān)聯(lián)研究中的統(tǒng)計方法

勾建偉,劉應安,夏業(yè)茂

?

多元表型與基因型的全基因組關(guān)聯(lián)研究中的統(tǒng)計方法

勾建偉,劉應安,夏業(yè)茂

南京林業(yè)大學理學院, 江蘇 南京 210037

復雜疾病往往需要多元表型變量共同刻畫，然而，常規(guī)的全基因組關(guān)聯(lián)研究僅是檢驗單個表型，該策略往往由于忽略某些表型的信息而降低檢驗復雜疾病與基因關(guān)聯(lián)性的效能。多元表型與基因型的全基因關(guān)聯(lián)研究可以更有效揭示復雜疾病與基因的關(guān)聯(lián)，同時對現(xiàn)有的統(tǒng)計方法提出挑戰(zhàn)。本文綜述了現(xiàn)有的處理多元表型與基因型的全基因組關(guān)聯(lián)的統(tǒng)計方法，討論了多元因變量多重回歸模型框架下的回歸系數(shù)與方差逆矩陣同時選擇的組合稀疏方法的前景。該組合稀疏方法既能充分利用多元表型變量的相關(guān)性信息，增強變量選擇的準確性，又能提供基因變異和表型關(guān)聯(lián)程度的度量指標。同時基于多元因變量的稀疏方法可以推廣到基因組學數(shù)據(jù)的整合分析中。

多元表型; 基因型; 全基因組關(guān)聯(lián)性;多元統(tǒng)計分析

傳統(tǒng)的全基因組關(guān)聯(lián)研究(Genome-wide association studies，GWAS)是通過考察單核苷酸多態(tài)性(Single nucleotide polymorphism, SNP)位點與單個表型變量之間的關(guān)系。然而，刻畫復雜疾病的量往往是錯綜復雜的，絕非僅僅單個表型變量就能夠充分描述。因此，簡單的用單一測量表型指標定義復雜疾病往往缺乏代表性。例如，身體質(zhì)量指數(shù)(Body Mass Index，BMI)并不能全面地代表心腦血管疾病相關(guān)的肥胖特征，因為它只是粗略地測量了在給定體表面積的平均體重，并沒有展示脂肪的分布。研究表明BMI、腰圍、臀圍這三個表型決定的體型更能反映心腦血管病風險[1]。因此，基于BMI、腰圍、臀圍的三個表型變量特征的GWAS與傳統(tǒng)基于單個表型變量BMI的GWAS相比更合理。另一方面，從數(shù)量遺傳學的觀點看，同一種復雜疾病表型變量不僅受多個不同基因的影響，同一基因的變異也不只影響到一個表型變量，往往會同時影響到多個表型變量，這些表型變量之間往往具有相關(guān)性?？梢?，復雜疾病往往是由多個表型變量共同表現(xiàn)出來，同時檢驗多個表型變量比檢驗一個表型變量來分析疾病的相關(guān)性的效能要高。然而在傳統(tǒng)GWAS背景下，通常將多元表型變量信息綜合為一個得分(Score)單變量來進行基因關(guān)聯(lián)分析。該策略往往會忽略某些對疾病發(fā)生、發(fā)展、轉(zhuǎn)歸有影響的表型信息，同時表型變量數(shù)據(jù)可能是二分類數(shù)據(jù)、有序數(shù)據(jù)、或者連續(xù)型數(shù)據(jù)以及三者間的任意混合，此時如何構(gòu)造得分變量也面臨挑戰(zhàn)。傳統(tǒng)的GWAS雖然已經(jīng)取得豐碩成果，但距離完全揭開疾病的遺傳密碼相差甚遠，已有研究[2]提出表型變量的復雜性、測量偏差及表型變量選取的不恰當都會大大削弱基因的關(guān)聯(lián)性研究，造成信息損失?；趩蝹€表型與基因的關(guān)聯(lián)性分析難以全面揭示遺傳與復雜疾病的聯(lián)系。因此，我們需要通過多元表型與基因型的全基因關(guān)聯(lián)研究，充分利用眾多相關(guān)表型變量，探測出與疾病相關(guān)聯(lián)的基因變異位點，從而為疾病的發(fā)病機制研究提供依據(jù)。

1 處理多元表型與基因型的全基因組關(guān)聯(lián)研究的現(xiàn)有統(tǒng)計方法

傳統(tǒng)的檢驗多元表型與基因相關(guān)性的方法是逐個對表型單變量與基因關(guān)聯(lián)性進行檢驗，然后采用多重檢驗矯正。然而同一疾病的多元表型變量往往具有相關(guān)性，多重檢驗校正會降低效能。已有的研究利用多元表型變量間的相關(guān)性，從多元表型變量中提取類似潛變量因子的策略進行降維，然后再進行關(guān)聯(lián)性分析，現(xiàn)有的統(tǒng)計方法概括如下。

1.1 多元因變量的主成分分析方法

通過對因變量進行主成分析，消除因變量之間的多重相關(guān)性[3]。對多表型的GWAS中，大多數(shù)研究[4,5]是采用前幾個因變量的主成分與基因位點進行線性回歸分析。然而，對于選取主成分的個數(shù)沒有嚴格標準。已有研究[6]通過模擬不同相關(guān)性的表型變量，比較了不同情形下的多元因變量的主成分分析方法，并得出僅采用少數(shù)幾個因變量主成分進行相關(guān)分析往往會降低檢驗的效能這一結(jié)論。同時，這類多元因變量的主成分方法僅僅綜合了因變量的信息，對于超高維的基因自變量并沒有給出有效的降維分析策略。

1.2 典型相關(guān)分析方法

分別提取自變量集與因變量集的最大主成分，通過兩個主成分的相關(guān)關(guān)系推測自變量集與因變量集之間的相關(guān)關(guān)系[7]。由于典型相關(guān)分析也采用了主成分提取的思路，因此典型相關(guān)成分也擁有主成分降維的性質(zhì)。已有研究[8-10]采用典型相關(guān)分析，對多表型變量與單個基因位點或者多個基因位點的相關(guān)性進行檢驗。然而，典型相關(guān)分析只是通過優(yōu)化樣本之間的協(xié)方差準則來求解典型變量，并沒有考慮原始樣本在高維空間中分布的結(jié)構(gòu)，當GWAS這種相對于變量數(shù)目來說樣本無法充分大的數(shù)據(jù)結(jié)構(gòu)下，往往受樣本影響，不太穩(wěn)定。

1.3 偏最小二乘回歸方法

偏最小二乘方法是多元因變量對多自變量的回歸分析方法[11]。該方法考慮到自變量與因變量的相關(guān)程度，確保提取的成分對因變量解釋能力強。已有研究[12]利用偏最小二乘回歸方法對多元表型變量與多個位點集的相關(guān)性進行檢驗。然而，該方法也存在一些缺點，無法給出模型參數(shù)檢驗統(tǒng)計量的分布，只能采用置換檢驗(Permutation Test)得到檢驗的值。置換檢驗運算對于高維數(shù)據(jù)運算緩慢，無法應用在超高維的全基因組水平上，只能用于后期的驗證性研究。同時，偏最小二乘提取出潛變量后，僅執(zhí)行假設檢驗而不能參數(shù)估計，因此無法給出基因位點與表型變量相關(guān)性大小的度量，從而無法評價基因變異與疾病相關(guān)性大小。

1.4 非參數(shù)檢驗方法

若干用于多元表型GWAS的非參數(shù)檢驗方法也被提出，Zhang等[13]提出了廣義Kendall’s秩檢驗用于檢驗單個基因位點與多表型變量的關(guān)聯(lián)性。朱文圣等[14]提出改進的廣義Kendall’s秩檢驗，研究了包含協(xié)變量調(diào)整下的基因與多表型變量的相關(guān)性檢驗。Van der Sluis等[15]提出組合值法，對每個表型變量關(guān)聯(lián)分析的值提出一個組合值，可以有效的處理多元表型與基因關(guān)聯(lián)性分析。

1.5 其它多元因變量的方法

基本是上述所有方法的不同擴展或改進，上面提到的方法主要是將基因位點作為一個集合(SNPs集)，檢驗SNPs集與元多表型變量的相關(guān)性，適用于變量維數(shù)和樣本數(shù)相差不大的情形，然后通過多重檢驗應用到全基因組。對于高維自變量數(shù)據(jù)結(jié)構(gòu)的全基因關(guān)聯(lián)研究，另外一條思路是降維。通常的做法就是在上面的方法上增加懲罰函數(shù)，提出對應的稀疏方法。例如，已有研究[16]提出偏最小二乘方法與1懲罰結(jié)合的稀疏偏最小二乘方法，對多元表型變量與基因位點集的相關(guān)性給出了檢驗。該方法被應用到多元表型的若干GWAS中[17,18]。

上述所有方法不預先設定模型結(jié)構(gòu)形式，所選方法都由所觀測到的樣本數(shù)據(jù)所決定，因此具有較大的適應性，在降低建模偏差方面有較大優(yōu)勢，但也存在明顯的局限和不足：首先，變量的維數(shù)較高時，上述方法所涉及的方差及協(xié)方差往往是奇異的，給數(shù)值計算與理論性質(zhì)帶來了難度。其次，上述方法對基因變異的相關(guān)性都是采用的檢驗，給出檢驗的值，但很難給出基因關(guān)聯(lián)效應的估計值，無法合理解釋基因變異影響性狀表型的大小。從理論和實際應用角度講，好的多元表型與基因關(guān)聯(lián)分析方法應該具備以下幾個特點：能夠給出更快速更有效的估計值，提供度量基因變異和表型變量關(guān)聯(lián)大小的客觀指標；能夠有效給出表型變量與基因位點回歸參數(shù)的漸進分布；能夠?qū)z傳數(shù)據(jù)的潛在模型提出更好的解釋，揭示疾病和基因變異之間的復雜關(guān)聯(lián)性。

2 稀疏多元因變量多重回歸方法的前景和意義

稀疏回歸或者懲罰回歸模型是在損失函數(shù)的基礎(chǔ)上增加懲罰函數(shù)，通過調(diào)整參數(shù)的選擇，將不顯著參數(shù)壓縮為零，對顯著參數(shù)進行很小壓縮或者不壓縮，并且給出該參數(shù)的估計。稀疏回歸具有以下優(yōu)勢：允許變量的維數(shù)超過樣本總數(shù)，通過變量選擇選出顯著的變量；能夠在進行變量選擇的同時，給出模型參數(shù)的優(yōu)良估計。從而，稀疏回歸方法在高維數(shù)據(jù)分析中有著良好的表現(xiàn)，是近年來高維數(shù)據(jù)領(lǐng)域中十分熱門的課題，比如：LASSO[19]、LARS[20]、Elastic net[21]、SCAD[22,23]、自適應LASSO[24]、成組LASSO[25]、自適應成組LASSO[26]、兩水平懲罰方法[27]以及圖LASSO方法[28]。

另一方面，多元因變量多重線性回歸模型建立了多元因變量與多自變量之間如下的線性關(guān)系：=+其中表示×的因變量矩陣，表示×的設計矩陣，是×的未知系數(shù)矩陣，表示×的誤差量矩陣，并假定的個行向量獨立同分布于元正態(tài)分布(0;?)。

其中Ω=Σ-1，稱為精度矩陣，刻畫了元因變量的兩兩條件相關(guān)性。

模型中待估的回歸參數(shù)刻畫預測變量與因變量的線性關(guān)系，精度矩陣表示多元因變量的兩兩條件相關(guān)性。當變量維數(shù)特別大時，假設僅有一部分自變量與因變量相關(guān)，其余變量都是噪聲變量，此時系數(shù)矩陣應該是稀疏的；當因變量維數(shù)特別大，若干因變量之間的相關(guān)性應該比較小，此時精度矩陣Ω應該是稀疏的。當回歸系數(shù)矩陣和精度矩陣都滿足稀疏性時，可以通過參數(shù)的懲罰函數(shù)，實現(xiàn)參數(shù)的稀疏性約束。稀疏的多元因變量多重線性回歸的模型估計轉(zhuǎn)變?yōu)椋?/p>

同時，稀疏多元因變量多重回歸模型不僅可以實現(xiàn)變量選擇，還可以給出回歸系數(shù)矩陣和精度矩陣的估計?；貧w系數(shù)的估計值可以反映自變量和因變量的相關(guān)性；利用精度矩陣與高斯圖解模型的關(guān)系，我們可以根據(jù)精度矩陣刻畫出因變量的相關(guān)關(guān)系圖結(jié)構(gòu)，進一步深入理解復雜疾病的遺傳機制。

綜上，在多元因變量多重線性回歸框架下發(fā)展和應用懲罰回歸方法，將同時具備上面提到的三個特點。稀疏的多元因變量多重線性回歸將為基于多元表型與基因型的復雜疾病關(guān)聯(lián)研究提供更加靈敏和強有力的統(tǒng)計分析工具。

3 討論

針對多元表型與基因型數(shù)據(jù)的復雜疾病關(guān)聯(lián)分析，借助組合稀疏回歸模型，提出在多元因變量多重回歸的框架下進行組合懲罰回歸的設想，根據(jù)不同的數(shù)據(jù)結(jié)構(gòu)選擇和構(gòu)造不同的懲罰函數(shù)1與2。稀疏的多元因變量多重回歸模型中回歸系數(shù)矩陣和精度矩陣Ω同時約束的最優(yōu)解求解往往不易求。現(xiàn)有的方法主要是分為兩類：令2(Ω)=0的降秩方法與令1()=0的高維協(xié)方差選擇的圖LASSO(graph LASSO)方法。這兩類方法無法同時考慮自變量和因變量進行變量選擇的高維數(shù)據(jù)結(jié)構(gòu)。我們可以借助懲罰似然與懲罰條件似然的關(guān)系，探索同時約束回歸系數(shù)矩陣和精度矩陣Ω的最優(yōu)解算法。

針對現(xiàn)有多元表型關(guān)聯(lián)性分析檢驗方法統(tǒng)計量分布未知和實際應用受限等問題，需要給出稀疏回歸模型參數(shù)估計理論框架，包括研究估計值的無偏性、相合性以及推導漸進分布；結(jié)合坐標下降優(yōu)化算法，給出參數(shù)計算快速有效的算法。方法比較、估計值的理論研究（無偏性，相合性以及漸進分布）、數(shù)值模擬和實例分析等方面系統(tǒng)探索新提出的方法，為多元表型的GWAS數(shù)據(jù)和基因組學數(shù)據(jù)的整合分析提供實用、高效和靈敏的關(guān)聯(lián)性分析工具其順利實施不但能夠豐富和發(fā)展關(guān)聯(lián)性分析理論，也對進一步深入理解疾病遺傳基礎(chǔ)和解釋遺傳缺失具有重要意義。

[1] Wells JC, Treleaven P, Cole TJ. BMI compared with 3-dimensional body shape: the UK National Sizing Survey[J]. Am J. Clin. Nutr., 2007(85):419-425

[2] Sluis SVD, Verhage M, Posthuma D,. Phenotypic complexity, measurement bias, and poor phenotypic resolution contribute to the missing heritability problem in genetic association studies[J]. Plos One, 2010,5(11):e13929

[3] Cadima J, Jolliffe IT. Loading and correlations in the interpretation of principle components[J]. Journal of Applied Statistics, 1995(22):203-214

[4] Liu F, Van DLF, Schurmann C,. A genome-wide association study identifies five loci influencing facial morphology in Europeans[J]. Plos genetics, 2012,8(9):e1002932

[5] Zhang F, Guo X, Wu S,. Genome-wide pathway association studies of multiple correlated quantitative phenotypes using principle component analyses[J]. Plos One, 2012,7(12):e53320

[6] Aschard H, Vilhjálmsson BJ, Greliche N,. Maximizing the power of principal-component analysis of correlated phenotypes in genome-wide association studies[J]. The American Journal of Human Genetics, 2014,94(5):662-676

[7] Hotelling H. Relations between two sets of variates[J]. Biometrika, 1936,28(3/4):321-377

[8] Galesloot TE, Van SK, Kiemeney LA,. A comparison of multivariate genome-wide association methods[J]. Plos One, 2014,9(4):e95923

[9] Seoane JA, Campbell C, Day INM,. Canonical Correlation Analysis for Gene-Based Pleiotropy Discovery[J]. Plos Computational Biology, 2014,10(10):e1003876

[10] Tang CS, Ferreira MAR. A gene-based test of association using canonical correlation analysis[J]. Bioinformatics, 2012,28(6):845-850

[11] Geladi P, Kowalski BR. Partial least-squares regression: a tutorial[J]. Analytica chimica acta, 1985,185(86):1-17

[12] Xue F, Li S, Luan J,. A latent variable partial least squares path modeling approach to regional association and polygenic effect with applications to a human obesity study[J]. Plos One, 2012,7(2):e31927

[13] Zhang H, Liu CT, Wang X. An association test for multiple traits based on the generalized Kendall’s tau[J]. Journal of the American Statistical Association, 2010,105(490):473-481

[14] Zhu W, Jiang Y, Zhang H. Nonparametric covariate-adjusted association tests based on the generalized Kendall's Tau.[J]. Journal of the American Statistical Association, 2012,107(497):1-11

[15] Van dSS, Posthuma D, Dolan CV. TATES: efficient multivariate genotype-phenotype analysis for genome-wide association studies[J]. Plos genetics, 2013,9(1):e1003235

[16] Chun H, Ballard DH, Cho J,Identification of association between disease and multiple markers via sparse partial least‐squares regression[J]. Genetic epidemiology, 2011,35(6):479-486

[17] Allen GI, Peterson C, Vannucci M,. Regularized partial least squares with an application to NMR spectroscopy[J]. Statistical Analysis and Data Mining the ASA Data Science Journal, 2013,6(4):302-314

[18] Rinnan A, Andersson M, Ridder C,. Recursive weighted partial least squares (rPLS): an efficient variable selection method using PLS[J]. Journal of Chemometrics, 2014,28(5):439-447

[19] Tibshirani R. Regression shrinkage and selection via the lasso[J]. Journal of the Royal Statistical Society: Series B (Methodological), 1996,58(1):267-288

[20] Efron B, Hastie T, Johnstone I, et al. Least angle regression[J]. Institute of Mathematical Statistics, 2004,32(2):407-451

[21] Zou H, Hastie T. Regularization and variable selection via the elastic net[J]. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 2005,67(2):301-320

[22] Fan J, Li R. Variable Selection via Nonconcave Penalized Likelihood and its Oracle Properties[J]. Journal of the American Statistical Association, 2001,96(456):1348-1360

[23] Kim Y, Choi H, Oh HS. Smoothly clipped absolute deviation on high dimensions[J]. Journal of the American Statistical Association, 2008,103(484):1665-1673

[24] Zou H. The Adaptive Lasso and Its Oracle Properties[J]. Journal of the American Statistical Association, 2006,101(476):1418-1429

[25] Yuan M, Lin Y. Model selection and estimation in regression with grouped variables[J]. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 2006,68(1):49-67

[26] Wang H, Leng C. A note on adaptive group lasso[J]. Computational Statistics & Data Analysis, 2008,52(12):5277-5286

[27] Seetharaman I. Consistent bi-level variable selection via composite group bridge penalized regression[D].Kansas, USA: Kansas State Univesity, 2013

[28] Fang Y, Wang R, Dai B,. Graph-based learning via auto-grouped sparse regularization and kernelized extension[J]. Transactions on Knowledge & Data Engineering, 2015,27(1):142-154

[29] Yuan M. Dimension reduction and coefficient estimation in multivariate linear regression[J]. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 2007,69(3):329-346

Statistical Method of Multiple Phenotype and Genotype in Genome-wide Association Study

GOU Jian-wei, LIU Ying-an, XIA Ye-mao

210037,

The complex disease is often identified with multiple traits. However, the common genome-wide association studies only test the association between a single phenotype trait and a large number of single nucleotide polymorphisms (SNPs), Which will lead to loss in statistical power when multivariate phenotype are involved. The genome-wide association study of multiple phenotype and genotype can more effectively reveal the association between complex diseases and genes, and challenge existing statistical methods. This paper reviewed the existing statistical methods for dealing with genome-wide associations between multivariate phenotype and genotype, and discussed the prospects of combined sparse methods for simultaneous selection of regression coefficients and variance inverse matrices under the framework of multivariate dependent multivariate regression models. The combined sparse method not only improved the accuracy of variable selection by utilizing joint information among multivariate phenotype, but also offered a measure of the association between genotype and phenotype. At the same time, the method could extend the proposed sparse multivariate regression methods to the integrative analysis of genomic data.

Multiple phenotype; genotype; genome-wide association; multivariate statistical analysis

O212

A

1000-2324(2018)05-0906-05

10.3969/j.issn.1000-2324.2018.05.038

2017-09-02

2017-10-13

國家自然科學基金(11471161);江蘇省高?；?15KJB110010)

勾建偉(1982-),女,博士,講師,研究方向為統(tǒng)計基因組學. E-mail:gjw1983@139.com