MR‐Egger 回歸在孟德?tīng)栯S機(jī)化分析中的應(yīng)用

徐藝耘 劉振球 樊 虹 張 欣 施婷婷 吳 聲 張鐵軍，2△

（1復(fù)旦大學(xué)公共衛(wèi)生學(xué)院流行病學(xué)教研室 上海 200032；2復(fù)旦大學(xué)義烏研究院 義烏 322000）

孟德?tīng)栯S機(jī)化（Mendelian randomization，MR）是借助遺傳變異（genetic variation）作為工具變量（instrumental variable，IV）來(lái)推斷暴露因素與結(jié)局之間因果關(guān)聯(lián)的方法，有效避免了反向因果關(guān)聯(lián)和潛在混雜因素導(dǎo)致的偏倚。隨著全基因組關(guān)聯(lián)研究（genome-wide association study，GWAS）的增多，孟德?tīng)栯S機(jī)化分析廣泛應(yīng)用于觀察性流行病學(xué)研究中。傳統(tǒng)孟德?tīng)栯S機(jī)化分析方法，如逆方差加權(quán)法（inverse-variance weighted，IVW），假定所有的遺傳變異均滿(mǎn)足工具變量的3 個(gè)核心假設(shè)[1］：（1）工具變量與暴露因素之間有強(qiáng)相關(guān)（關(guān)聯(lián)性假設(shè)）；（2）工具變量與混雜因素?zé)o關(guān)聯(lián)（獨(dú)立性假設(shè)）；（3）工具變量只能通過(guò)暴露因素對(duì)結(jié)局產(chǎn)生作用，不能通過(guò)其他途徑影響結(jié)局發(fā)生（排他性假設(shè)）。當(dāng)工具變量存在多效性時(shí)，因果效應(yīng)的估計(jì)則會(huì)有偏[2］。MR-Egger 回歸是在IVW 的基礎(chǔ)上修正而來(lái)的一種基于匯總數(shù)據(jù)的多工具變量孟德?tīng)栯S機(jī)化方法。與IVW 不同的是，該方法僅需滿(mǎn)足工具變量多效性效應(yīng)獨(dú)立于工具變量與暴露因素之間的關(guān)聯(lián)（instrument strength independent of direct effect，InSIDE）假設(shè)和無(wú)測(cè)量誤差（no measurement error，NOME）假設(shè)，不如工具變量的3 個(gè)核心假設(shè)要求嚴(yán)格。同時(shí)，該方法既能檢測(cè)多效性又能校正多效性偏倚，所以在以多個(gè)遺傳變異作為工具變量的研究中，采用MR-Egger 回歸能夠保持孟德?tīng)栯S機(jī)化方法的有效性。 本研究將對(duì)MR-Egger 回歸進(jìn)行闡述，并用此方法探索高密度脂蛋白膽固醇（highdensity lipoprotein cholesterol，HDL-C）與冠狀動(dòng)脈疾?。╟oronary artery disease，CAD）之間的因果關(guān)聯(lián)。

資料和方法

基本原理MR-Egger 回歸是由Bowden 等[3］于2015年提出，廣泛運(yùn)用于采用多個(gè)單核苷酸多態(tài)性（single nucleotide polymorphism，SNP）作為工具變量進(jìn)行因果推斷的孟德?tīng)栯S機(jī)化分析中，尤其適用于遺傳變異存在定向多效性的情況。該方法的核心在于加權(quán)線性回歸時(shí)考慮截距項(xiàng)的存在，利用截距項(xiàng)來(lái)衡量工具變量間平均多效性的大小，斜率則是因果效應(yīng)的無(wú)偏估計(jì)。MR-Egger 回歸的效用取決于兩個(gè)前提：第一，InSIDE 假設(shè)要求SNP 對(duì)暴露因素的影響?yīng)毩⒂谒鼈儗?duì)結(jié)局的多效性影響，弱化了排他性假設(shè)；第二，NOME假設(shè)要求工具變量與暴露因素的關(guān)聯(lián)中不存在測(cè)量誤差，I2統(tǒng)計(jì)量可以評(píng)估是否滿(mǎn)足NOME 假設(shè)，如果不滿(mǎn)足，則會(huì)產(chǎn)生弱工具變量偏倚[4］。

記J 個(gè)工具變量為Gj（j=1，2，…，J），暴露因素為X，結(jié)局為Y，工具變量Gj 對(duì)暴露因素X 和結(jié)局Y的效應(yīng)估計(jì)值分別為β?Xj和β?Yj及其對(duì)應(yīng)的標(biāo)準(zhǔn)誤為和，遺傳變異與結(jié)局之間的關(guān)聯(lián)可表示為：

其中，αj 表示工具變量Gj 對(duì)結(jié)局Y 的直接影響，當(dāng)αj≠0 時(shí)，工具變量具有多效性；θ 表示暴露因素X 對(duì)結(jié)局Y 的影響；θβXj則為工具變量Gj 通過(guò)暴露因素X 對(duì)結(jié)局Y 的間接影響[5］。

對(duì)于滿(mǎn)足工具變量所有假設(shè)的單一遺傳變異Gj（αj = 0），暴露因素X 對(duì)結(jié)局Y 的因果效應(yīng)可被一致性地估計(jì)為關(guān)聯(lián)估計(jì)的比值，比值及相應(yīng)的標(biāo)準(zhǔn)誤分別為：。

由于MR-Egger 回歸并沒(méi)有限制截距項(xiàng)為0，回歸方程可表示為：

當(dāng)InSIDE 假設(shè)和NOME 假設(shè)均滿(mǎn)足時(shí)，θ是因果效應(yīng)的無(wú)偏估計(jì)[6］：

若滿(mǎn)足InSIDE 假設(shè)但不滿(mǎn)足NOME 假設(shè)，因果效應(yīng)估計(jì)值則會(huì)有偏：

其中，s2為測(cè)量誤差，s2= 0 時(shí)才滿(mǎn)足NOME 假設(shè)；I2統(tǒng)計(jì)量是在Meta 分析中評(píng)估研究間異質(zhì)性的工具，的大小可用來(lái)估計(jì)，代表SNP-暴露關(guān)聯(lián)的真實(shí)方差，取值范圍為0～1，越接近1，MR-Egger回歸的估計(jì)越接近真實(shí)因果效應(yīng)。其中，Cochran’s Q 統(tǒng)計(jì)量是評(píng)估工具變量間異質(zhì)性的工具，是SNP-暴露關(guān)聯(lián)的平均值。

IVW 法是由Burgess 等[2］于2013年提出的一種加權(quán)線性回歸模型。與MR-Egger 回歸不同，該方法限制回歸直線通過(guò)原點(diǎn)，將兩個(gè)或多個(gè)SNP 的效應(yīng)值進(jìn)行整合，因果效應(yīng)估計(jì)值為：

實(shí)例分析HDL-C 資料來(lái)自全球脂質(zhì)遺傳學(xué)（Global Lipids Genetics Consortium，GLGC）數(shù)據(jù)庫(kù)（http：//csg. sph. umich. edu/willer/public/lipids2013/），該數(shù)據(jù)庫(kù)中有關(guān)基因位點(diǎn)與HDL-C的GWAS 研究結(jié)果發(fā)表于2013年，該研究樣本量為94595 人，分析了2418527 個(gè)與HDL-C 相關(guān)的SNP 位 點(diǎn)[7］。 與HDL-C相關(guān)的SNP篩選標(biāo)準(zhǔn)：（1）達(dá)到全基因組統(tǒng)計(jì)顯著性水平，即P＜5×10-8；（2）去除連鎖不平衡（linkage disequilibrium，LD）的SNPs，排除標(biāo)準(zhǔn)為r2＞0.01，kb= 5000，篩選出125個(gè)與HDL-C 相關(guān)且相互獨(dú)立的SNP，提取每個(gè)SNP 的rs 編號(hào)、效應(yīng)等位基因（effect allele，EA）和次要等位基因（non-effect allele，non-EA）、β 系數(shù)、P值和標(biāo)準(zhǔn)誤。 CAD 數(shù)據(jù)來(lái)源于2018年Van Der Harst 等[8］發(fā)表的GWAS 統(tǒng)計(jì)結(jié)果，該研究將CARDIoGRAMplusC4D[Coronary Artery Disease Genome wide Replication and Meta-analysis（CARDIoGRAM） plus The Coronary Artery Disease（C4D）Genetics］的數(shù)據(jù)與英國(guó)生物銀行（UK Biobank，UKBB）的CAD 數(shù)據(jù)進(jìn)行GWAS分析，共包括547216 人，其中病例組122733 人，對(duì)照組424528 人，分析了69033 個(gè)與CAD 相關(guān)的SNP位點(diǎn)。 設(shè)置最小等位基因頻率MAF=0.01，從GLGC 數(shù)據(jù)庫(kù)篩選的SNP 在CAD 來(lái)源的GWAS研究中均存在，經(jīng)等位基因?qū)R后，rs11637365、rs3790106、rs4986970、rs6589581 和rs964184 由于存在回文結(jié)構(gòu)被剔除，最終有120 個(gè)SNP 作為工具變量分析HDL-C 與CAD 的因果關(guān)聯(lián)。上述研究人群均為歐洲血統(tǒng)，具有相似的年齡和性別，符合兩樣本孟德?tīng)栯S機(jī)化分析的要求。 利用R4.0.3 中的TwoSampleMR 包進(jìn)行兩樣本孟德?tīng)栯S機(jī)化分析，檢驗(yàn)水準(zhǔn)α=0.05。 因果效應(yīng)分析使用IVW 法和MR-Egger 回歸，對(duì)比兩者差異；敏感性分析包括利用MR-Egger 回歸進(jìn)行多效性檢驗(yàn)（pleiotropy test）和逐個(gè)剔除檢驗(yàn)（leave-one-out sensitivity test），“l(fā)eave-one-out”是指逐一剔除SNP，分別計(jì)算剩下SNP 的合并效應(yīng)，如果剔除某個(gè)SNP 后其他SNP估計(jì)出來(lái)的MR 結(jié)果和總結(jié)果相差很大，說(shuō)明MR結(jié)果對(duì)該SNP 是敏感的。

結(jié)果

本研究針對(duì)HDL-C 和CAD 選取有效SNP 作為IV，經(jīng)篩選確定了120 個(gè)SNP 作為IV，表1 列出了前10 個(gè)SNP 的相關(guān)信息。采用MR-Egger 回歸和IVW 法估計(jì)因果效應(yīng)（表2）。 IVW 結(jié)果支持HDL-C 和CAD之間存在因果關(guān)系（OR=0.82，95%CI：0.75～0.89），且HDL-C 每增加1 個(gè)標(biāo)準(zhǔn)差，CAD 風(fēng)險(xiǎn)降低18%。而MR-Egger 法不支持（OR=0.96，95%CI：0.83～1.11），即HDL-C 對(duì)CAD 的發(fā)生風(fēng)險(xiǎn)沒(méi)有影響。多效性檢驗(yàn)結(jié)果顯示，MR-Egger回歸的截距為-0.01，P=0.008，說(shuō)明因果分析結(jié)果會(huì)受到水平多效性的影響。在水平多效性存在的情況下，傳統(tǒng)孟德?tīng)栯S機(jī)化方法IVW 的估計(jì)出現(xiàn)了偏差，MR-Egger 回歸則能識(shí)別并校正多效性，給出接近真實(shí)值的因果估計(jì)。敏感性分析結(jié)果見(jiàn)圖1，所有線條均在0 的左側(cè)，說(shuō)明無(wú)論去除哪個(gè)SNP 都不會(huì)對(duì)結(jié)果產(chǎn)生根本影響，即孟德?tīng)栯S機(jī)化的結(jié)果穩(wěn)健。

表1 工具變量SNP 信息Tab 1 Information of the instrumental variable SNP

表2 因果效應(yīng)估計(jì)Tab 2 Estimation of causal effects

圖1 逐個(gè)剔除檢驗(yàn)Fig 1 Leave‐one‐out sensitivity analysis

討論

早期孟德?tīng)栯S機(jī)化研究?jī)A向于使用單一遺傳變異，并集中在單一研究人群中的特定危險(xiǎn)因素與疾病的關(guān)聯(lián)研究。近年GWAS 中大量增長(zhǎng)的基因型-表型關(guān)聯(lián)導(dǎo)致眾多暴露因素的大量遺傳變異被識(shí)別出來(lái)[9］。許多遺傳變異被認(rèn)為具有多效性效應(yīng)，要證明排他性假設(shè)的成立幾乎不可能。 MREgger 回歸區(qū)別于傳統(tǒng)孟德?tīng)栯S機(jī)化方法之處在于放寬了SNP 之間不存在水平多效性的要求，它假設(shè)在基因-暴露關(guān)聯(lián)和基因變異對(duì)結(jié)果的直接影響之間沒(méi)有相關(guān)性，是比更嚴(yán)格的排他性標(biāo)準(zhǔn)更弱的要求，即使所有工具變量是無(wú)效的，MR-Egger 回歸也能得到準(zhǔn)確的因果效應(yīng)估計(jì)[10］。因此，在工具變量存在多效性的情況下，MR-Egger 回歸是一個(gè)更好的選擇。

本研究借助GLGC 數(shù)據(jù)庫(kù)和Van Der Harst 的GWAS 統(tǒng)計(jì)結(jié)果，提取與HDL-C 和CAD 均密切相關(guān)且相互獨(dú)立的SNP 作為工具變量，利用MREgger 回歸和IVW 法進(jìn)行兩樣本孟德?tīng)栯S機(jī)化分析，估計(jì)HDL-C 和CAD 的因果關(guān)聯(lián)并進(jìn)行比較。IVW 結(jié)果顯示兩者之間具有因果關(guān)系，HDL-C 對(duì)于CAD 有保護(hù)作用，但I(xiàn)VW 法假設(shè)截距為0，未考慮所納入的工具變量的基因多效性；而MR-Egger回歸并不局限于截距為0，是在IVW 基礎(chǔ)上的修正，同時(shí)考慮基因多效性對(duì)結(jié)果的影響，因此給出了并不一致的估計(jì)結(jié)果，即HDL-C 與CAD 之間不存在因果關(guān)聯(lián)，與以往孟德?tīng)栯S機(jī)化研究結(jié)果相吻合[11］。隨后的逐一剔除檢驗(yàn)顯示MR-Egger 結(jié)果更具有穩(wěn)健性。因此，MR-Egger 回歸在孟德?tīng)栯S機(jī)化分析中有較好的應(yīng)用價(jià)值。

與傳統(tǒng)孟德?tīng)栯S機(jī)化方法相比，MR-Egger 回歸的偏倚更小，而且因果無(wú)效假設(shè)的排斥率更接近名義上的5%[3］。MR-Egger 回歸的局限性是它在因果估計(jì)方面的統(tǒng)計(jì)效力較低。如果沒(méi)有證據(jù)表明IVW 和MR-Egger 之間存在系統(tǒng)性差異，那么通常采用IVW，雖然IVW 可能不那么穩(wěn)健，但在此情況下，IVW 估計(jì)值的精度更高。有研究者指出MREgger 回歸應(yīng)被視為一種敏感性分析，用來(lái)檢驗(yàn)是否違反工具變量假設(shè)，而不是作為上述更常見(jiàn)方法的替代[12］。實(shí)際應(yīng)用時(shí)應(yīng)結(jié)合多種方法的結(jié)果進(jìn)行綜合考量，如果多種方法的結(jié)果一致，則結(jié)論更為可靠，從而更好地應(yīng)用于傳統(tǒng)觀察性流行病學(xué)研究，為病因推斷提供更強(qiáng)有力的證據(jù)，并有利于為未來(lái)疾病預(yù)防的干預(yù)目標(biāo)設(shè)定優(yōu)先級(jí)。

作者貢獻(xiàn)聲明徐藝耘 查閱文獻(xiàn)，論文構(gòu)思、撰寫(xiě)和修訂，數(shù)據(jù)統(tǒng)計(jì)分析。劉振球，樊虹，張欣，施婷婷，吳聲 結(jié)果解釋?zhuān)撐男薷摹?張鐵軍論文構(gòu)思、指導(dǎo)和修訂。

利益沖突聲明所有作者均聲明不存在利益沖突。