孟德爾隨機化研究潰瘍性結腸炎與甲狀腺功能減退發生風險的因果關系

花胤 王曉燕 王振 辛永寧 劉守勝

1南京醫科大學青島臨床醫學院（青島市市立醫院）（山東青島 266011）；青島市市立醫院2國際門診，3感染性疾病科，4臨床研究中心（山東青島 266011）

潰瘍性結腸炎（ulcerative colitis，UC）是一種病程反復的炎癥性腸?。╥nflammatory bowel diseases，IBD），常由復雜的環境、免疫以及遺傳因素共同引起，其從結腸遠端開始發病并向近端擴展直至全結腸［1］。作為一種免疫介導的炎癥性腸病，雖然UC主要影響消化道，但患者常出現腸道以外癥狀，包括肌肉骨骼、皮膚、眼部、肝膽、血液系統及內分泌疾病等［2］。甲狀腺功能減退（hypothyroidism，HT）作為一種內分泌疾病，常由一種自身免疫疾病橋本病引起，表現為免疫細胞攻擊甲狀腺，最終引起甲狀腺的功能喪失導致HT 的發生［3-5］。有研究［6-8］發現UC 患者發生甲狀腺疾病包括HT 的風險較正常個體更高。然而另外的研究［9］則發現IBD 患者與正常個體之間的HT 發病風險無明顯差異。近期的一項回顧性研究［10］發現，IBD 會降低患者發生HT 的風險，并不推薦IBD 患者在缺少臨床癥狀的情況下實施甲狀腺超聲檢查。因此，目前雖已有較多研究報道了UC 與HT 之間的關系［11］，但結論仍然具有爭議，并且既往研究大多是靠觀察性研究得到的關聯結論，可能受到某些混雜因素及反向因果關系的干擾。孟德爾隨機化（mendelian randomization，MR）作為一種疾病之間因果關系的推斷方法在遺傳病學領域廣泛應用，MR 通過遺傳變異作為工具變量來推斷暴露變量能否影響結局，可以避免觀察性研究中可能出現的缺點。因此，本研究擬通過兩樣本孟德爾隨機化（two sample mendelian randomization，TSMR）分析來確定UC 與HT 之間發病風險的相關性。

1 資料與方法

1.1 研究原理 孟德爾隨機化的基本原理是通過工具變量來分析暴露與結局關聯的因果關系。其分析的準確性需要滿足下述3 個核心假設：（1）關聯性假設：工具變量需要與暴露因素強相關；（2）排他性假設：工具變量只能通過暴露因素對結局產生影響；（3）獨立性假設：工具變量需要和結局及混雜因素沒有關聯［12］（圖1）。

圖1 孟德爾隨機化的3 個核心假設Fig.1 Key assumptions of Two-Sample MR Analysis

1.2 數據來源 通過搜索全基因組關聯研究（genome wide association study，GWAS）。分別獲得UC 和HT 的GWAS 匯總數據，其中UC 樣本數據來自芬蘭生物銀行（https：//www.finngen.fi/fi），共計211 498個樣本（觀察組：1 198；對照組：210 300），共計16 380 453 個SNP 位點。HT 樣本數據來自IEU open GWAS（https：//gwas.mrcieu.ac.uk/datasets/ebi-a-GCST90018862/），此數據庫樣本量為410 141（觀察組：30 155；對照組：379 986），共計24 138 872 個SNP 位點（表1）。

1.3 工具變量的篩選 根據以下條件篩選暴露的工具變量：（1）選取與暴露的關聯強度達到全基因組顯著關聯水平（P＜ 5 × 10-8）的SNPs；（2）排除SNPs 之間的連鎖不平衡影響（linkage disequilibrium，LD），設置連鎖不平衡參數（r2）閾值為0.01，遺傳距離為1 000 kb，從而確保SNP的獨立性；（3）排除位于回文序列中的工具變量；（4）使用PhenoScanner數據庫剔除與混雜因素相關的遺傳工具變量，重點排除表型與結局有相關意義的SNP。

1.4 統計分析 分別采用逆方差加權法（inverse variance weighted，IVW）、Egger 回歸法（MR-Egger）、中位數加權法（weighted median）開展TSMR，并采用森林圖（forest plot），散點圖（scatter plot）可視化上述方法的OR值（odds ratio）及95%CI（confidence interval）。IVW 法通過Wald估計值獲得加權平均值，是本研究的主要分析工具，其主要原理是將每個遺傳變異看作獨立的工具變量，并通過薈萃分析合并工具變量的因果效應估值，當所有的遺傳變異都是有效的工具變量時，IVW 法是因果效應的最有效評價工具。MR-Egger 回歸法與IVW 法類似，但可以通過回歸截距評估工具變量的水平多效性。中位數加權法對因果效應的評估最少只需要50%的有效工具變量，因此能在存在無效工具變量的情況下，提供穩定的因果效應值。關于敏感性分析，本研究利用MR-PRESSO（MR Pleiotropy Residual Sum and Outlier，MR-PRESSO）法評價是否存在離群的SNP 位點［13］，采用MREgger 的截距值評估水平多效性，采用Cochran′s Q檢測及漏斗圖（funnel test）評估SNPs 中的異質性，采用留一法分析（leave one out）分析是否存在單個SNP 能夠影響因果評估的結果，以上統計分析均通過開源軟件R（版本4.3.1）及R包“TwoSampleMR”（版本0.5.7）進行。以P＜ 0.05 為差異有統計學意義。

1.5 工具變量評價 工具變量的強度通過統計量F計算，其公式為F=R2× (N-k- 1)/(1 -R2) ×k。其中R2代表工具變量解釋的變異比例，N為暴露樣本的樣本量，k為工具變量的數量，而R2=(2 ×EAF×(1 -EAF) ×beta2)/(N×SE2)，EAF為次要等位基因頻率，beta為等位基因效應值，SE為標準誤差［14］。

2 結果

2.1 UC 對HT 的MR 分析

2.1.1 遺傳工具變量確定 本研究選擇與UC 相關具有全基因組顯著性（P＜ 5 × 10-8）的遺傳變量SNP 位點進行匯總，篩選出同時滿足假設1、2、3 的SNP，得到4 個符合假設的SNP。利用Pheno Scanner 數據庫去除與結局HT 具有相關意義的SNP（n= 0），通過GWAS 提取結局HT 的數據，獲取以上4 個SNP 與結局的顯著關系。合并暴露與結局數據集，并剔除了存在連鎖不平衡和回文結構的SNP，最終共納入了4 個SNP 作為工具變量進行MR 分析。

2.1.2 TSMR 研究結果 本研究的TSMR 結果提示UC 可降低HT 的發生風險。IVW 法提供了主要的因果證據（OR= 0.957，95%CI：0.924 ～ 0.990，SE= 0.017，P＜ 0.05），三種MR 統計方法的結果見表2，MR-egger 法的結果（OR= 0.936，95%CI：0.924 ～ 0.990，SE= 0.032，P＞ 0.05）雖然不顯著但OR方向同樣支持IVW 法的結果。Weighted median法的結果（OR= 0.944，SE= 0.021，P＜ 0.001）其OR方向也支持了IVW 法的因果判斷。因此可知UC對HT 具有負向因果關系，TSMR 分析的散點圖和森林圖見圖2、3。

表2 TSMR 分析結果Tab.2 Results of the TSMR analysis

圖2 兩樣本孟德爾隨機化散點圖Fig.2 Scatter Plot of Two-Sample MR Analysis

圖3 兩樣本孟德爾隨機化森林圖Fig.3 Forest Plots of Two-Sample MR Analysis

2.1.3 敏感性分析 Cochran′s Q 檢驗評估異質性Q= 2.566，P= 0.463，表明納入的SNPs 無明顯異質性，漏斗圖也進一步驗證了無異質性的結果（圖4）。對于水平多效性，MR-Egger 法截距值為0.011 7（P= 0.523），表明納入的SNPs無明顯水平多效性，說明工具變量并不通過暴露以外的途徑影響結局。MR-PRESSO 法未發現任何異常的SNP（P= 0.548）。留一法檢驗顯示去除任意SNP 后結果仍然穩定（圖5）。

圖4 兩樣本孟德爾隨機化漏斗圖Fig.4 Funnel plot of Two-Sample MR Analysis

圖5 “留一法”敏感性分析結果Fig.5 The results of sensitivity analysis of “leave-one-out”

2.1.4 工具變量評價 利用R2及F值的計算公式得出，其中4 個SNP 的統計量F值均＞ 10（31.07 ～65.43），表明本研究的工具變量均為強工具變量。

2.2 HT 對UC 的反向因果分析

2.2.1 遺傳工具變量確定 同樣選擇與HT相關具有全基因組顯著性（P＜ 5 × 10-8）的遺傳變量SNP位點進行匯總，篩選出同時滿足假設1、2、3 的SNP，得到70 個符合假設的SNP。利用PhenoScanner 數據庫去除與結局UC 具有相關意義的SNP（n= 0），通過GWAS 提取結局UC 的數據，獲取以上70 個SNP 與結局的顯著關系。合并暴露與結局數據集，并剔除了2 個存在回文結構的SNP（rs2412976、rs2921053），并根據MR-PRESSO 結果剔除1 個離群的SNP（rs9271365），最終共納入了67 個SNP 作為工具變量進行MR 分析。

2.2.2 TSMR 研究結果 本研究的TSMR 結果提示HT 對UC 不存在反向的因果關系。IVW 法（P=0.521），MR-Egger 法（P= 0.231），中位數加權法（P= 0.544）的結果均不提示存在顯著因果關系（P＞ 0.05）。

2.2.3 敏感性分析 同樣采用IVW 法評估異質性Q= 78.831，P= 0.133，表明納入的SNPs 無明顯異質性，漏斗圖也進一步加強了無異質性的結果。MR-Egger 法截距值為-0.014（P= 0.306），表明納入的SNPs 無明顯水平多效性，說明工具變量并不通過暴露以外的途徑影響結局。此外留一法檢驗示去除任意SNP 后結果仍然穩定。

2.2.4 工具變量評價 反向因果關系統計量F均＞ 10（30.03 ～ 0523.14），都為強工具變量。

3 討論

在本研究中，我們通過TSMR 評價了UC 和HT之間潛在的雙向因果關系。結果證實了UC 對HT的負向因果關系，表明并不推薦對無癥狀的UC 患者開展額外的HT 檢查。此外我們識別了4 個聯系UC 與HT 的SNP 位點。這些發現有助于加深對UC 與HT 之間潛在的病理機制理解。

在UC 對HT 的MR 分析結果中我們發現，MREgger 法得到了與另兩種方法不一致的結果。有研究認為MR-Egger 法只能用作來評估工具變量的敏感性分析，而不能代替IVW 法的結果，并且本研究的質量評估中，工具變量并不存在異質性和水平多效性［15］。所以盡管MR-egger 法存在不一致的結果，我們仍然判斷在UC 和HT 之間存在負向的因果關系。

既往的觀察性研究［7-8，10，16-20］已經發現了一些UC 和HT 之間的關系。DORE 等［10］在一項回顧性病例對照研究指出，UC 患者與健康個體相比患有HT 的風險更低（OR=0.33；95%CI：0.17 ～ 0.66）。RICART 等［19］也發現IBD 患者患有自身免疫性甲狀腺疾病的風險較低。此外POORAN 等［20］在一項橫斷面研究中發現，同屬于IBD 的克羅恩病患者發生甲亢的風險與正常個體相似，但HT 的風險更低。這些研究的結果與我們的孟德爾分析結果一致，即UC 患者與HT 之間存在負向因果關系。雖然也有一些文獻指出UC 會增加HT 的風險，例如KAPPELMAN 等［7］的研究中發現兒童UC 患者發生HT 的風險更高，MAHMOUD 等［8］發現IBD 患者中患有HT 的頻率更高。但是BAR YEHUDA 等［21］2019 年發表的一項納入12 625 例IBD 患者的IBD與自身免疫性疾病的研究表明，IBD 患者與大部分自身免疫性疾病相關如：風濕性關節炎、過敏性哮喘、銀屑病等，但這其中并不包括甲狀腺疾病。針對本研究UC 對HT 負向的因果關系解釋，有研究指出針對IBD 的免疫抑制治療可能抑制了自身免疫性HT 的產生，PASCHOU 等［22］發現IBD患者經過抗腫瘤壞死因子（anti-TNF）如阿達利姆單抗（adalimumab），英夫利昔單抗（infliximab）治療后，與一般治療法相比甲狀腺自身抗體呈下降趨勢。

本研究的優勢主要在于首次運用了孟德爾隨機化法研究了UC 與HT 之間的因果關系，并且我們還進行了反向TSMR 分析，以驗證可能的反向因果關系。但是在我們的研究中仍然存在一些局限性。第一，本研究的選用的數據集是匯總類型的數據，無法通過年齡分層或者性別來進行亞組分析。第二，數據樣本集包括了歐洲人群，這限制了對非歐洲人群的應用。需要進一步的研究驗證對非歐洲人群的適用性。第三，本研究只納入了4 個工具變量，雖然都為強工具變量，但存在工具偏倚的可能，因此需要謹慎地解釋研究結果。在未來的研究中，可納入更多的遺傳標記作為工具變量，并進一步探索環境和生活方式因素對UC與HT 之間的因果關系，以及多種混雜因素對UC與HT 之間因果關系的影響。

【Author contributions】HUA Yin and WANG Zhen participated in the planning process and wrote the manuscript draft.WANG Xiaoyan performed the data analyses.XIN Yongning and LIU Shousheng interpreted the results of the data analyses and approved the final manuscript.All authors read and approved the final manuscript as submitted.

【Conflict of interest】The authors declare no conflict of interest.