少肌癥相關表型單變量與雙變量GWAS 分析比較

海 榮 張 壘 高學文 (內(nèi)蒙古自治區(qū)人民醫(yī)院干部保健所，呼和浩特 010020)

瘦體重是少肌癥重要的復雜性狀。其性狀值受到多個微效基因以及環(huán)境因素的影響，進行關聯(lián)研究是為探索少肌癥分子病理提供新線索的有效方法。Deng 等［1］研究顯示，在美國白種人12q(Log of Odds Score，LOD=1.79)、17p(LOD=1.64)、5q(LOD=1.59)和7p(LOD=1.52)處可能存在顯著影響瘦體重變化的數(shù)量性狀位點。最近，Zhao［2］等報道5q(LOD=3.50)、15q(LOD=2.72)和7q(LOD=2.79)處可能存在顯著影響瘦體重變化的數(shù)量性狀位點。Deng 等［3］報道，6p(LOD=3.46)可能存在顯著影響瘦體重變化的數(shù)量性狀位點。Liu等［4］采用全基因組關聯(lián)研究方法鑒定出與瘦體重顯著相關的TRHR 基因。

AAM 是指女性的初次月經(jīng)的年齡，月經(jīng)初潮是女性性成熟的重要標志，也是女性具有生殖能力的生物學信號［5］。AAM 的提早和延遲對女性身體健康有很大的影響。AAM 是由多種遺傳因素和環(huán)境因素決定［6，7］，最近的研究表明遺傳因素對AAM 有重要決定作用。從臨床上來看，研究影響AAM 的遺傳因素具有重要意義。

多變量分析是一種能同時考慮所有表型的好方法，目前認為在增加統(tǒng)計效力和參數(shù)估計精確性上比單變量分析更有優(yōu)勢。過去的多變量分析主要是采用雙變量連鎖分析為主，在當前疾病基因定位研究的前沿，連鎖分析顯然已經(jīng)相對“過時”，目前，雙變量連鎖分析已經(jīng)成功運用到疾病基因定位領域，但是關于雙變量關聯(lián)分析，研究得很少。到目前為止，運用基于GEE 理論進行多變量的關聯(lián)分析方法只有Lange 等［8］發(fā)展的基于家系資料的多變量關聯(lián)分析、最近Liu 等［9］發(fā)展的基于群體的關于兩個混雜性狀(一個質(zhì)量性狀、一個數(shù)量性狀)的關聯(lián)分析和Pei 等［10］發(fā)展的基于群體的兩個數(shù)量性狀的單體型關聯(lián)分析。然而，關于這三種方法的運用還很少，說明雙變量關聯(lián)研究還沒有引起足夠的重視。本研究為進一步深入研究影響瘦體重的潛在因素采用了雙變量全基因組關聯(lián)研究方法，首次對白種人女性人群的LBM 與AAM 進行關聯(lián)分析。

1 材料與方法

1.1 樣本信息 樣本來自美國奧馬哈及堪薩斯城及周邊城市，驗證樣本來自中國西安及周邊地區(qū)的漢族人群。一般信息見表1。

1.2 主要的試劑及大型儀器 由華大基因與上海理工大學系統(tǒng)生物醫(yī)學中心提供。

1.3 表型測量 測量體重、身高、腰圍、臀圍以及前臂長。一般指標測量結束后采10 ml 血4℃冰箱保存。然后采用Hologic DEXA4500 掃描儀測量四肢，全身瘦體重以及全身脂肪量。

1.4 DNA 提取 采用Gentra systems Minneapolis，MN，試劑盒提取，操作按說明進行。

1.5 SNP 分型 采用了Affymetrix 人類基因組SNP6.0 分型芯片。

1.6 質(zhì)量控制 基因分型后，必須對分型結果進行質(zhì)量檢測以剔除部分不合格的數(shù)據(jù)。

1.7 統(tǒng)計學處理

1.7.1 人群分層 由于全基因組關聯(lián)分析大多采用隨機樣本。隨機樣本的一個最大問題就是樣本在構成上可能存在分層現(xiàn)象，不同祖先的人可能具有不同的遺傳背景，而這些遺傳背景很可能會干擾關聯(lián)分析結果導致假陽性或假陰性。Structure2.2 與EIGENSTRAT 軟件分析的人群分層分析結果見圖1。

1.7.2 SNP 與LBM 與AAM 的雙變量關聯(lián)分析 首先，我們采用SPSS18.0 的逐步回歸法檢驗了年齡、年齡的平方、身高是否對LBM 具有顯著影響，然后以顯著的因素(P ＜0.05)為自變量對LBM 進行線性校正。校正后的殘差使用Minitab(Minitab Inc.，State College，PA，USA)中 的 Kolmogorov-Smirnov 法對其正態(tài)性進行檢驗，對AAM 也采用Minitab 軟件檢驗是否正態(tài)分布。在全基因組關聯(lián)研究中，由于SNP 的數(shù)量很大，多重比較是一個很難解決的問題。為防止發(fā)生假陽性或假陰性結果，我們采用Bonferroni 校正法對原始P 值進行校正。一般認為校正后P 值應達到整體全基因組水平上顯著的。其次，我們采用R 軟件(http://www.rproject.org)進行單、雙變量關聯(lián)分析。

2 結果

2.1 四肢LBM 分析結果 全基因研究樣本和驗證樣本的基本信息列于表1 和表2 中，單變量分析發(fā)現(xiàn)SNP rs1860547 和rs11030746 全基因組水平上不顯著(初潮年齡P=0.08)，具體信息見表3，在白種人研究樣本中SNP rs1860547 和rs11303746 在雙變量GWA 分析中全基因組水平非常顯著(P=2.67E-14.14)，在中國人驗證樣本中也達到了顯著水平(P值＜0.05)，具體信息見表4。

2.2 研究樣本和驗證樣本群體分層分析結果 群體分層是引起假陽性或假陰性關聯(lián)的一個重要干擾因素。因此，需要對樣本進行遺傳背景評估。我們使用EIGENSTRAT 軟件對研究群體和驗證樣本分析其遺傳結構，兩個樣本人群分析結果見圖1。

2.3 LBM 與AAM 雙變量全基因組關聯(lián)分析結果 研究發(fā)現(xiàn)雙變量GWA分析SNPrs1860547 和rs11030746 與LBM 與AAM 顯著關聯(lián)。單變量GWA 分析發(fā)現(xiàn)以上兩個SNP 與LBM 密切相關。在genecard 網(wǎng)站中查詢發(fā)現(xiàn) SNP rs1860547 和rs11030746 左側與右側發(fā)現(xiàn)3 個顯著關聯(lián)基因KCNA1、KCNA4、KCNA5，見表3、4。

表1 研究樣本和驗證樣本基本信息Tab.1 Basic characteristics of the study sample

表2 四肢LBM 分析結果Tab.2 Appendicular LBM analysis results

圖1 群體分層Fig.1 Population stratification

表3 單變量GWA 分析信息Tab.3 Univariate GWA analysis information

表4 雙變量GWA 分析信息Tab.4 Bivariate GWA analysis information

3 討論

本研究采用雙變量全基因組關聯(lián)分析方法，首次在白種人群中開展了LBM 和AAM 關聯(lián)分析。共得2 個SNP 得 到 了 驗 證(P ＜0.05)即 SNP rs1860547 和SNP rs11030746，前者左側47 513 bp處有KCNA1 基因，在右側78 149 bp 處有KCNA5 基因，在后者的右側395 214 bp 處有KCNA4 基因，據(jù)報道KCNA(鉀電壓閥門通道，混合器相關亞家族)基因與骨骼肌神經(jīng)系統(tǒng)疾病密切相關。鉀離子通道(KCN)家族主要的功能是可以通過調(diào)控細胞內(nèi)pH值和離子濃度來維持正常細胞體積及細胞內(nèi)生物分子活性所需的離子濃度范圍。在細胞信號傳導中具有重要功能，包括調(diào)節(jié)胰島素分泌、平滑肌收縮、上皮細胞電解質(zhì)轉(zhuǎn)運和心率等［11］。它主要分布于心臟、骨骼肌、血管、神經(jīng)、氣管、胃腸道、腎臟、血液和內(nèi)分泌腺體等。根據(jù)其結構和功能特點分為4 種，即電壓敏感型、受體偶聯(lián)型、鈣敏感型及其余類型等，共50 多個亞型。其中電壓敏感型的鉀離子通道又稱電壓依賴性鉀離子通道(Voltage-gated Potassium channels，Kv)，是鉀離子通道超家族中的重要成員［12］。2003 年，國際上通過Kv 命名法將鉀離子通道家族分為3 大類，其中Kv1-Kv6 和Kv8-Kv9 為一個家族;Kv7 為一個家族;Kv10-Kv12 為一個家族;對應的基因以KCN 命名(圖2)。例如電壓門控通道鉀離子蛋白α 第1、4、5 亞基蛋白分別命名為Kv1.1、Kv1.4、Kv1.5，編碼它們的基因分別命名為KCNA1、KCNA4、KCNA5 基因［13］。其中，KCNA 基因編碼的Kvl(Shaker)包含最多亞型(見圖2)。2007 年Hoegg 等［14］采用最大似然法詳細繪制了KCNA 基因家族樹狀圖(圖2)。

目前，所有離子通道病中骨骼肌離子通道病是了解最為清楚的，已發(fā)現(xiàn)十幾種與各種離子通道病變有關的骨骼肌疾病。Browne［15］第一次報道KCNA1 基因的表達會影響發(fā)作性共濟失調(diào)的發(fā)病進程。最近研究發(fā)現(xiàn)此基因的缺失型突變可以引起神經(jīng)性肌肉強直和良性家族性癲癇［16］。Wipff 等［17］研究發(fā)現(xiàn)特發(fā)性肺動脈高壓患者的KCNA5 基因變異會降低KCNA5 表達，使肺動脈平滑肌細胞過度增殖和凋亡減少，從而與特發(fā)性肺動脈高壓相關聯(lián)。Olson 等［18］首次報道了由于KCNA5 基因突變所導致的Kv1.5 功能缺失，是復極受損和心房纖顫新的危險因素。Bouchard 等［19］對高加索人群采用基因組掃描技術進行了運動能力相關基因研究時發(fā)現(xiàn)KCNA4 基因通過對心肌電信號傳導來調(diào)節(jié)影響心臟收縮能力從而進一步影響運動耐力。有研究發(fā)現(xiàn)Kv1.7 通道蛋白(由KCNA7 基因編碼)是骨骼肌細胞膜上低氧誘導的一個分子調(diào)節(jié)元件。當骨骼肌缺氧時，可以通過Kv1.7 通道蛋白的結構的變化進行調(diào)節(jié)。從以上幾位學者的研究結果中發(fā)現(xiàn)KCNA1、KCNA4、KCNA5 基因與心肌和平滑肌功能密切相關。KCNA7 基因與骨骼肌耐缺氧能力有關。目前沒有 發(fā) 現(xiàn) 有 關KCNA1、KCNA4、KCNA5 基 因 與骨骼肌相關性報道。我們的研究第一次發(fā)現(xiàn)KCNA1、KCNA4、KCNA5 基因與LBM 和AAM 有關。

圖2 KCN 家族最大似然法繪制圖［14］Fig.2 Maximum likelihood tree of KCN gene family［14］

骨骼肌具有一套復雜的結構，骨骼肌纖維一般為長圓柱形。一般長約1～40 mm，直徑10～100 μm。骨骼肌纖維表面為肌膜，每條肌纖維周圍均有一薄層結締組織稱為肌內(nèi)膜，由數(shù)條至數(shù)十條肌纖維集合成肌束，肌束外有較厚的結締組織稱為肌束膜，由許多肌束組成一塊肌肉，其表面的結締組織稱肌外膜，即深筯膜。以上各結締組織中均有豐富的血管和神經(jīng)。

根據(jù)骨骼肌結構我們推測，骨骼肌各結締組織中豐富的脊髓運動神經(jīng)元伸出的軸突至外周支配肌纖維，信號從神經(jīng)到肌膜再到肌肉收縮裝置的整個過程，離子通道均起了重要作用。因此我們進一步推測脊髓運動神經(jīng)元以及骨骼肌細胞膜上存在很多由KCNA1、KCNA4、KCNA5 基因編碼的相應的通道蛋白Kv1.1、Kv1.4、Kv1.5。KCNA1、KCNA4、KCNA5 基因突變等遺傳變異引起相應的通道蛋白Kv1.1、Kv1.4、Kv1.5 結構變化而無法正常調(diào)節(jié)電信號從神經(jīng)到肌膜再到肌肉收縮裝置。從而使脊髓運動神經(jīng)元功能性失神經(jīng)，進而影響肌纖維代謝和收縮。Brown 等［20］研究發(fā)現(xiàn)功能性失神經(jīng)引起的肌肉代謝和運動神經(jīng)元營養(yǎng)缺失可以引起肌肉量的減少而導致少肌癥。在圖2 中我們也可以觀察到KCNA1、KCNA4、KCNA5 基因處于同一個亞家族，而且KCNA1 和KCNA5 基因位置都是12p13 區(qū)，而且它們之間的距離比較近。因此，我們推測KCNA1、KCNA4、KCNA5 功能具有一定的相似性。

到目前為止未見有關KCNA1、KCNA4、KCNA5與AAM 相關性報道。從青春期開始，卵巢內(nèi)的卵子陸續(xù)發(fā)育成熟并排出。與此同時卵巢分泌雌激素和孕激素，促使子宮內(nèi)膜增厚和血管增生，為受精卵在子宮內(nèi)發(fā)育成熟和種植創(chuàng)造條件。排出的卵子如果沒有受精，卵巢的雌激素和孕激素的分泌會很快減少，引起子宮內(nèi)膜組織壞死脫落，血管破裂出血。脫落的子宮內(nèi)膜碎片連同血液一起由陰道排出稱為月經(jīng)初潮。那么KCNA1、KCNA4、KCNA5 以什么途徑來影響初潮年齡呢?因此，根據(jù)我們研究的結果和子宮卵巢結構來進行推測和假設KCNA1、KCNA4、KCNA5 基因是如何影響初潮年齡。在前面的內(nèi)容中提到鉀離子通道是迄今為止分布最廣、亞型最多、作用最復雜的一類離子通道，也存在于內(nèi)分泌腺體等部位的細胞。卵巢是女性重要的性腺，其功能與AAM 密切相關。由KCNA1、KCNA4、KCNA5基因編碼的相應通道蛋白Kv1.1、Kv1.4、Kv1.5 可能廣泛分布于卵巢腺體細胞膜上，當這些基因發(fā)生突變時相應的通道蛋白結構也發(fā)生變化，進而影響雌激素和孕激素的分泌信號的傳遞而影響子宮內(nèi)膜增厚和血管增生。因此，KCNA1、KCNA4、KCNA5 基因突變的個體和此三個基因正常的個體初潮年齡時間可能不同。

［1］Deng HW，Deng H，Liu YJ，et al.Genome wide linkage scan for quantitative-trait loci for obesity Phenotypes［J］.Am J Hum Genet，2002，70(5):1138-1151.

［2］Zhao LJ，Xiao P，Liu YJ，et al.A genome-wide linkage scan for quantitative trait loci underlying obesity related phenotypes in 434 Caucasian families［J］.Hum Genet，2007，121(1):145-148.

［3］Deng FY，Xiao P，Lei SF，et al.Bivariate whole genome linkage analysis for femoral neck geometric parameters and total body lean mass［J］.J Bone Miner Res，2007，22(6):808-816.

［4］Liu XG，Tan LJ，Lei SF，et al.Genome-wide association and replication studies identified TRHR as an important gene for lean body mass［J］.Am J Hum Genet，2009，84(3):418-423.

［5］Chen CT，F(xiàn)ernández-Rhodes L，Brzyski RG，et al.Replication of loci influencing ages at menarche and menopause in Hispanic women:the Women's Health Initiative SHARe Study［J］.Hum Mol Genet，2012，21(6):1419-1432.

［6］Rees M.The age of menarche［J］.ORGYN，1995，(4):2-4.

［7］Treloar SA，Martin NG.Age at menarche as a fitness trait:nonadditive genetic variance detected in a large twin sample［J］.Am J Hum Genet，1990，47(1):137-148.

［8］Lange C，Silverman EK，Xu X，et al.A multivariate family-based association test using generalized estimating equations:FBAT-GEE［J］.Biostatistics，2003，4 (2):195-206.

［9］Liu J，Pei Y，Papasian CJ，et al.Bivariate association analyses for the mixture of continuous and binary traits with the use of extended generalized estimating equations［J］.Genet Epidemiol，2009，33(3):217-227.

［10］Pei YF，Zhang L，Liu J，et al.Multivariate association test using haplotype trend regretssion［J］.Ann Hum Genet，2009，73(4):456-464.

［11］Coetzee WA，Amarillo Y，Chiu J，et al.Molecular diversity of K+channels［J］.Ann NY Acad Sci，1999，868:233-285.

［12］Wickenden A.K channels as therapeutic drug targets［J］.Pharmacol Ther，2002，94(1-2):157-182.

［13］Catterall WA，Chandy KG，Clapham DE.International Union of pharmacology:Approaches to the Nomenclature of Voltage-Gated Ion Channels［J］.Pharmacol Rev，2003，55(4):573-574.

［14］Hoegg S，Meyer A.Phylogenomic analyses of KCNA gene clusters in vertebrates:why do gene clusters stay intact?［J］.BMC Evol Biol，2007，7(139):3-12.

［15］Browne DL，Gancher ST，Nutt JG，et al.Episodic ataxia/myokymia syndrome is associated with point mutations in the human potassium channel gene，KCNA1［J］.Nat Genet，1994，8(2):136-140.

［16］Ishida S，Sakamoto Y，Nishio T，et al.Kcna1-mutant rats dominantly display myokymia，neuromyotonia and spontaneous epileptic seizures［J］.Brain Res，2012，1435:154-166.

［17］Wipff J，Dieudé P，Guedj M，et al.Association of a KCNA5 gene polymorphism with systemic sclerosis-associated pulmonary arterial hypertension in the European Caucasian population［J］.Arthritis Rheum，2010，62(10):3093-3100.

［18］Olson TM，Alekseev AE，Liu XK，et al.Kv1.5 channelo pathy due to KCNA5 loss of function mutation causes human atrial fibrillation［J］.Hum Mol Genet，2006，15(14):2185-2191.

［19］Bouchard C，Rankinen T，Chagnon YC，et al.Genomic scan for maximal oxygen uptake and its response to training in the HERITAGE Family Study［J］.J Appl Physiol，2000，88(2):551-559.

［20］Brown WF.A method for estimating the number of motor units in thenar muscles and the changes in motor unit count with ageing［J］.J Neurol Neurosurg Psychiatry，1972，35(6):845-852.