首都醫科大學數學教研室(100069) 華 琳 鄭衛英 張金旺 劉 紅 閆 巖
同胞對連鎖分析的三種統計方法比較*
首都醫科大學數學教研室(100069) 華 琳 鄭衛英 張金旺 劉 紅△閆 巖
*:北京市教委科技發展計劃面上項目(KM200910025006)
△通訊作者:劉紅
在很多物種中,無論是低等的酵母,還是高等的人類,許多基因的表達水平都展示出豐富的自然變異。而目前,很多科研工作都集中于研究疾病發展不同階段的基因表達水平,較少研究健康人群基因表達的自然變化。事實上,健康人群個體間的遺傳變異常常和人類復雜疾病的易感性及復雜疾病性狀有密切聯系。而人類基因組中廣泛存在并被廣泛應用的第三代基因遺傳標記,即SNP標記,是一種能夠穩定遺傳的早期突變,研究者可以通過對SNP的相關分析和高密度的SNP圖譜來定位一系列復雜疾病的相關基因。
我們通過分析健康人群淋巴母細胞基因表達的自然變化,發現可能存在影響基因表達的標記位點,從而提示這些標記位點可能對遺傳、進化及復雜疾病的易感性有重要的生物學意義。我們利用三種統計分析方法,逐步判別分析、逐步logistic回歸和同胞對Haseman-Elston回歸,將Centre d'Etude du Polymorphisme Humain(CEPH)家族淋巴母細胞的基因表達譜數據和SNP基因型數據結合在一起。我們將1號染色體上的HSPA6的基因表達水平作為數量性狀,利用同胞對的H-E(Haseman-Elston)線性回歸方法對50個SNP標記進行連鎖定位,同時對基因表達水平進行重新設置,使之成為分類變量,分別采用逐步判別分析、逐步logistic回歸分析方法作連鎖分析,從而發現影響HSPA6表達水平的SNP遺傳標記。我們采用matlab 6.5編寫程序實現數據轉換,采用著名數量遺傳學家Dr.Robert C.Elston和他的研究小組開發的遺傳流行病學的統計分析系統 SAGE 5.2的 SIBPAL〔1〕模塊實現Haseman-Elston回歸。采用SAS軟件實現逐步判別分析、廣義可加模型和逐步logistic回歸。
同胞對的Haseman-Elston線性回歸法是根據同胞對間的表型差異與同胞對共享來自同一祖先的標記等位基因IBD(identity-by-descent)個數平均值估計的關系模型進行連鎖定位的。
設一個家系中s個同胞的第j對同胞的HSPA6表達值為x1j,x2j,同胞對的基因表達差異:

其中 β0為截距,βi和 di(i=1,2,…,m)分別為 m 個遺傳標記中第i個標記的加性和顯性遺傳方差,(i=1,2,…,m)分別為第 i個標記位點上同胞對共享IBD個數平均值估計值和同胞對共享兩個等位基因的概率估計值,ck是反應協變量yk的函數f(yk)對表達值的影響系數,e為隨機誤差。
當只考慮一個標記時,模型則退化為

我們采用單標記回歸法(即模型2),分別對CEPH 家族的4 個家系(1333,1340,1341,1345)54 個個體1號染色體的50個SNP標記作連鎖分析,估計回歸參數并進行統計學檢驗。檢驗假設為:H0:θ=(無連鎖);H1:θ<(存在連鎖),若 P <0.05,則可判斷這一SNP標記與HSPA6的基因表達水平存在連鎖。
逐步判別分析是在每一步選一個判別能力最強的指標變量進入判別函數,直至判別函數中所有指標的判別能力都有意義,而判別函數外的所有指標都不能納入。
設有m個指標變量,分別用x1,x2,…,xm來表示,有G個類別。具體步驟為:
(1)計算總離差陣 T,類內離差陣 W,T-1,W-1,行列式|T|,|W|及 Wilks統計量

其中U越小,判別能力越強。
(2)逐步剔除無統計學意義的指標變量
假設有p個指標變量已被選入判別函數。對于未入選的m-p個指標變量中的任何一個指標變量xr,我們將W和T矩陣分解為

若 F >Fα(G-1,N-p-G),則認為第 r個指標變量的判別能力有統計學意義。
我們將1號染色體上的HSPA6基因的表達水平計算均值。重新設置表達變量,如果個體的表達值大于均值,設為1,小于均值設為0。對98對同胞對作分析,如果同胞對的兩個個體表達變量均為0,我們作為G1類;如果兩個體表達變量一個為0,一個為1,我們作為G2類;如果兩個個體表達變量均為1,則作為G3類。用每個SNP標記的IBD值作為指標變量。
將98對同胞對中每個同胞對的HSPA6基因表達水平作差并取絕對值,即

將d'作為分類變量,將50個SNP標記的IBD作為協變量,作逐步logistic回歸分析。
首先采用Matlab 6.5軟件編寫程序將表達譜數據和SNP基因型數據轉化為SAGE5.2可識別的輸入文件。由于SIBPAL模塊需要有Geneibd文件,所以應先通過Geneibd模塊產生Geneibd文件。產生Geneibd文件需要有四個輸入文件,即家系資料文件、參數文件、等位基因位置文件和由Freq模塊生成的等位基因頻率文件。將生成的Geneibd數據文件通過我們編寫的程序輸入到SAS程序中,分別作逐步判別分析和逐步logistic回歸分析。我們將四種統計分析方法的結果進行比較,結果見表1。
從表1中可以看出逐步判別分析識別出6個SNP標記,逐步logistic回歸識別出3個SNP標記,Haseman-Elston線性回歸識別出2個SNP標記。從連鎖分析的統計結果來看,顯然逐步判別分析要優于其他方法,逐步logistic回歸和Haseman-Elston線性回歸的效果相似。三種方法均發現SNP標記位點rs5556161和rs991191有意義,提示這兩個標記位點附近可能存在影響HSPA6基因表達水平的性狀位點。基因HSPA6位于1號染色體(chr1q23),將此基因映射到 GO(Gene Ontology)數據庫,在生物學過程中分別為GO:6457((protein folding)和 GO:6986(response to unfolded protein),說明該基因與蛋白質的展開折疊相關。分子功能為 GO:166(nucleotide binding)和 GO:5524(ATP binding),說明該基因的分子功能是核酸結合與ATP結合。SNP rs5556161和SNP rs991191可能影響該基因的表達水平,從而影響基因功能。

表1 三種統計方法的連鎖分析結果比較(98對同胞對)
另一方面,我們可以通過看逐步判別分析過程中F統計量的變化來分析基因間的交互作用〔4〕(表2)。

表2 逐步判別分析中6個SNP標記的F統計量變化
從表 2中,我們可以觀察到,在第二步中,rs5556161的F值從6.10下降到0.50,而 rs1511687的F值從2.86上升到8.92,說明 rs991191可能與rs5556161和rs1511687存在交互作用。同樣,在第三步中,rs5556161的 F值從 0.50上升到 7.40,rs1511687的F值從8.92下降到2.74,rs2209698的F值從0.02上升到3.68,說明這三個SNP位點均可能與rs1924761存在交互作用。
本文采用了三種統計分析方法作同胞對的連鎖分析。它們同時均發現了2個SNP位點rs5556161和rs991191有統計學意義,說明這兩個標記位點附近可能存在影響HSPA6基因表達水平的性狀位點。但比較而言,逐步判別分析要優于其他兩種方法。本文采用的數據是健康人群淋巴母細胞的基因表達數據。而在實際應用中,這三種方法均可用于復雜疾病的連鎖分析。由于同胞對Haseman-Elston回歸要求性狀為數量性狀,而很多復雜疾病的數量性狀常常很難表達,逐步判別分析和逐步logistic回歸分析恰好可以將復雜性狀用分類變量表示,進而再作連鎖分析。因此,對于數量性狀位點QTL或表達數量性狀位點eQTL,適用于采用H-E回歸進行連鎖分析。而對于狀態變量(疾病或正常),可選擇采用逐步logistic回歸和逐步判別分析方法。逐步判別分析和逐步logistic回歸有很多類似的地方。但是,如果研究的位點較多(如全基因組的連鎖分析),會使得logistic回歸中變量大量增加,從而導致結果的不準確性,應進行多重校正。相比較而言,逐步判別分析在分析較多位點時,統計效能要高些。另一方面,如果我們想分析位點或基因的交互作用,通過逐步判別分析中SNP標記的F統計量變化,就可以初步發現這種交互作用,當然這種交互作用還需要大樣本的數據加以驗證。而且,逐步判別分析法還有一個優點,就是不需要預先了解基因表達值和各個SNP標記之間的關系,因此它不會象一般的線性回歸那樣對模型要求很高〔5〕。總之,SNP數據的統計挖掘技術作為一種新的連鎖分析方法必將會在復雜疾病的連鎖分析中起著重要作用。
1.SAGE Statistical Analysis for Genetic Epidemiology.Release 5.2.0
2.楊興云,張瑞杰,宮濱生,等.心血管疾病遺傳連鎖分析及發病風險因子的研究.數理醫藥學雜志,2004,17:485-488.
3.Li X,Rao SQ,Katly LM,et al.Genetic mapping of complex discrete human diseases by discriminant analysis.Progress in Natural Science,2002,12:431-437.
4.Cuo Z,Li X,Rao SQ,et al.Multivariate sibpair linkage analysis of longitudinal phenotypes by three stepwise analysis approaches.BMC Genetics,2003,4:17.
5.Rao SQ,Li L,Li X,et al.Genetic linkage analysis of longitudinal hypertension phenotypes using three summarymeasures.BMC Genetics,2003,4:24-31.