基于因果推斷方法檢驗mRNA在年齡和低級別膠質(zhì)瘤中的中介效應(yīng)*

王碧玨 郭志飛 楊海濤 李 治 王菊平 曹紅艷，6△ 周立業(yè)△

【提 要】 目的 構(gòu)建年齡、mRNA表達(dá)水平和低級別膠質(zhì)瘤(LGG)患者結(jié)局之間的因果中介推斷模型，識別介導(dǎo)年齡和LGG結(jié)局的中介因子，為LGG的治療和預(yù)后提供方法指導(dǎo)。方法 以TCGA數(shù)據(jù)庫中LGG數(shù)據(jù)為例，采用CIT和VanderWeele方法建立因果中介模型，識別年齡和LGG結(jié)局之間的介導(dǎo)基因。結(jié)果 識別出11個在年齡和LGG結(jié)局之間起中介作用的基因，中介間接效應(yīng)的RR值均大于1，表明這11個基因的表達(dá)異常可能是LGG不良結(jié)局的危險因素。結(jié)論 基于因果推斷的方法能夠識別介導(dǎo)年齡和低級別膠質(zhì)瘤結(jié)局的重要基因，為LGG的預(yù)后治療提供了潛在的預(yù)后標(biāo)志物和治療靶點。

人腦膠質(zhì)瘤多起源于神經(jīng)膠質(zhì)細(xì)胞，屬于神經(jīng)上皮源性腫瘤，是最常見的原發(fā)性顱內(nèi)惡性腫瘤，發(fā)病率(3～6.4)/10萬，分為低級別膠質(zhì)瘤(low-grade gliomas，LGG)(WHOⅠ、Ⅱ級)和高級別膠質(zhì)瘤(high-grade gliomas，HGG)(WHOⅢ、Ⅳ級)[1-2]。LGG約占膠質(zhì)瘤的15%～25%，具有侵襲性生長、放化療抵抗等惡性生物學(xué)特征，并且有向高級別膠質(zhì)瘤進(jìn)展的傾向[3-4]。LGG的發(fā)生發(fā)展是基因、環(huán)境共同調(diào)控的結(jié)果，年齡≥40歲是LGG患者不良預(yù)后的重要因素，不可逆及持續(xù)增長的年齡如何影響疾病一直是研究的熱點，但并未得到一致的結(jié)論[5]?；虮磉_(dá)在特定蛋白的合成中調(diào)控細(xì)胞增殖與分化，參與LGG癌變的發(fā)生發(fā)展，是LGG預(yù)后和靶向治療的重要候選生物標(biāo)志物[6]。年齡和基因表達(dá)均與LGG預(yù)后有關(guān)聯(lián)，故假設(shè)年齡、基因表達(dá)和LGG預(yù)后之間存在因果中介關(guān)系，即年齡可能調(diào)控基因表達(dá)影響LGG的預(yù)后。

因果推斷方法根據(jù)推斷理論的不同分為檢驗法和模型法，其中常用的檢驗法和模型法分別為因果推斷檢驗(causal inference test，CIT)[7]方法和VanderWeel[8]中介模型。本文以TCGA(The Cancer Genome Atlas)數(shù)據(jù)為例，聯(lián)合使用CIT因果中介方法和VanderWeele模型構(gòu)建以年齡(以40歲為二分類的分類標(biāo)準(zhǔn))為暴露、基因表達(dá)為中介、LGG患者預(yù)后為結(jié)局的因果推斷模型，探討年齡、mRNA表達(dá)水平和LGG患者結(jié)局之間的關(guān)系，識別中介因子。所識別的中介基因?qū)長GG患者的靶向治療提供理論依據(jù)，有利于改善患者預(yù)后，提高生命質(zhì)量。

材料與方法

1.數(shù)據(jù)資料來源與整理

使用TCGA Assembler軟件在TCGA網(wǎng)站下載LGG數(shù)據(jù)[9-10]，從中選取包含臨床和mRNA表達(dá)數(shù)據(jù)的患者共515例，刪除零表達(dá)值所占比例≥30%的基因，最終mRNA有17236個變量，結(jié)局為截止至隨訪時間是否死亡。

2.方法

(1)CIT

Millstein(2009)提出了判斷因果關(guān)系的CIT方法，將因果推斷問題視為一條數(shù)學(xué)條件的鏈，所謂“因果關(guān)系”，即中介的差異在一定比例上決定了結(jié)局的某些差異，即使這個比例很小。CIT方法的適用條件為：①連續(xù)型或離散型暴露變量，且不限制為單變量，也可是一組變量；②連續(xù)型潛在中間變量；③連續(xù)型或二分類結(jié)局變量[11-12]。CIT方法的基本思想是，引入暴露L，潛在中間變量G和結(jié)局變量T，三者需要滿足以下四個條件：①L與T相關(guān)，②G與T|L相關(guān)，③L與G|T相關(guān)，④L與T|G相互獨(dú)立。以年齡為L，基因表達(dá)為G，生存結(jié)局為T，四個條件對應(yīng)的模型如下：

(1)

(2)

(3)

(4)

H0：γ1=0，H1：γ1≠0，

H0：βi=0，H1：βi≠0，

H0：γ3=0，H1：γ3≠0，

H0：γ4≠0，H1：γ4=0.

CIT方法為基于似然的假設(shè)檢驗方法，若以上四個條件均滿足，則推斷潛在中間變量是結(jié)局的原因，即對每一個條件計算P值，將所有P值中的最大值作為整體檢驗的結(jié)果。

(2)VanderWeele模型

(5)

(6)

將總效應(yīng)分解為限制性直接效應(yīng)(controlled direct effects，CDE)、自然直接效應(yīng)(natural direct effects，NDE)和自然間接效應(yīng)(natural indirect effects，NIE)。通過計算直接效應(yīng)和間接效應(yīng)的RR值，進(jìn)一步估計處理因素對結(jié)局的直接效應(yīng)和間接效應(yīng)。間接效應(yīng)的計算見公式(7)，利用bootstrap法重復(fù)1000次計算RR值95%置信區(qū)間。

RRNIE=exp(θm×βx+θinteraction×βx)

(7)

(3)統(tǒng)計方法實現(xiàn)

應(yīng)用R軟件進(jìn)行統(tǒng)計分析，CIT第四步在cit包實現(xiàn)，bootstrap在boot包實現(xiàn)，顯著水平α=0.05。

結(jié) 果

1.CIT方法分析結(jié)果

(1)第一步分析結(jié)果

第一步校正性別協(xié)變量后，研究年齡與LGG的相關(guān)性，P值為0.000316，RR值為0.124，95%CI為(0.060，0.192)，表明年齡與LGG有關(guān)聯(lián)，即年齡對LGG患者的結(jié)局有直接效應(yīng)。

(2)第二步分析結(jié)果

第二步在校正性別協(xié)變量及控制年齡基礎(chǔ)上，研究LGG與每個基因表達(dá)水平的相關(guān)性。以每個基因的表達(dá)為自變量、LGG為應(yīng)變量擬合logistic回歸模型(2)，對檢驗水準(zhǔn)α進(jìn)行Bonferroni校正后，有557個(3.23%)基因與LGG結(jié)局有關(guān)聯(lián)。圖1為第二步驟的Mahattan圖，橫坐標(biāo)表示基因在23條染色體上的位置，縱坐標(biāo)為-log(P)尺度。所得有意義的基因?qū)⑦M(jìn)入到第三步的分析。

(3)第三步分析結(jié)果

第三步在校正性別協(xié)變量及控制LGG結(jié)局的基礎(chǔ)上，研究每個基因的表達(dá)水平與年齡之間的相關(guān)性。具體過程是以每個基因的mRNA表達(dá)水平為應(yīng)變量、年齡為自變量來擬合線性回歸模型(3)，用所得P值表示第j個基因的表達(dá)水平與年齡的關(guān)聯(lián)性，對檢驗水準(zhǔn)α進(jìn)行Bonferroni校正后，有25個(4.49%)基因與年齡有關(guān)聯(lián)。圖2為第三步驟的Mahattan圖。所得有意義的基因?qū)⑦M(jìn)入到第四步的分析中。

(4)第四步分析結(jié)果

第四步在校正性別協(xié)變量及控制基因表達(dá)水平基礎(chǔ)上，研究LGG與年齡之間是否相互獨(dú)立。具體過程是以LGG為應(yīng)變量、年齡為自變量來擬合線性回歸模型(4)。該步驟利用R軟件的“cit”包進(jìn)行分析，可得到P值和基于置換的FDR值(q值)，q值小于0.05表示控制基因的表達(dá)水平后LGG與年齡之間相互獨(dú)立。有11個(44%)基因的表達(dá)水平有統(tǒng)計學(xué)意義，分別為ARL3、DPY19L1、E2F7、FAM171A1、IGFBP2、LAMA2、MRPS16、PER3、PLAT、RANBP17和RASL10A。因此年齡可能通過調(diào)節(jié)這11個基因的表達(dá)水平來影響LGG結(jié)局。11個基因的位置信息和CIT方法四個條件的P值見表1。

圖1 17236個基因的Mahattan圖

圖2 557個基因的Mahattan圖

表1 CIT方法分析結(jié)果及基因的位置信息

2.VanderWeele分析結(jié)果

CIT方法識別出11個在年齡和LGG間起中介作用的基因，采用VanderWeele模型對中介的間接效應(yīng)(NIE)進(jìn)行估計，并利用bootstrap方法計算RR值的95%置信區(qū)間。圖3是11個基因的表達(dá)水平分別在年齡和LGG結(jié)局中的間接效應(yīng)的RR值。

結(jié)果表明11個基因在年齡和LGG之間的效應(yīng)均有意義，且均大于1，說明這11個基因的表達(dá)異?？赡苁荓GG不良結(jié)局的危險因素。

圖3 11個基因的mRNA表達(dá)水平分別在年齡和LGG之間的間接效應(yīng)大小

討 論

本研究基于因果推斷方法檢驗mRNA在年齡和LGG之間的中介效應(yīng)，年齡通過調(diào)控11個基因的表達(dá)，進(jìn)而影響LGG的預(yù)后，為LGG患者的治療提供了潛在的預(yù)后標(biāo)志物和治療靶點。

ARL3在膠質(zhì)瘤中表達(dá)下調(diào)[13]，其通過參與免疫微環(huán)境和血管生成，進(jìn)而影響膠質(zhì)瘤的生長、侵襲和轉(zhuǎn)移。DPY19L1在膠質(zhì)瘤中表達(dá)下調(diào)，有可能通過調(diào)節(jié)大腦皮層中谷氨酸能神經(jīng)元的徑向遷移而作用于膠質(zhì)瘤[14-15]。LAMA2基因編碼層粘連蛋白α2鏈，是一種重要的細(xì)胞外基質(zhì)蛋白，有研究發(fā)現(xiàn)層粘連蛋白可能通過促進(jìn)腦室下區(qū)少突膠質(zhì)細(xì)胞祖細(xì)胞的存活來調(diào)節(jié)產(chǎn)后少突膠質(zhì)細(xì)胞的產(chǎn)生[16]。因此，這些基因可能是膠質(zhì)瘤的預(yù)后標(biāo)志物和治療靶點。E2F7、PER3和RANBP17基因均為與細(xì)胞增殖、細(xì)胞周期有關(guān)的因子，其中，E2F7在膠質(zhì)瘤中高表達(dá)，可能通過與CDK1的協(xié)同作用參與膠質(zhì)瘤的發(fā)生發(fā)展[17]。E2F7的過表達(dá)可能通過上調(diào)Beclin-1誘導(dǎo)自噬而導(dǎo)致LGG的不良預(yù)后和結(jié)局[18]，因此提示E2F7可能是神經(jīng)膠質(zhì)瘤的不良預(yù)后指標(biāo)。PER3是人類Period基因家族中最重要的一員，可通過調(diào)控生物節(jié)律等重要生物通路而參與腫瘤的發(fā)生[19]。目前已發(fā)現(xiàn)PER3在乳腺癌[20]、頭頸部鱗狀細(xì)胞癌[21]等癌癥中表達(dá)異常，且與腫瘤浸潤程度、復(fù)發(fā)率和預(yù)后有較大關(guān)系。但其在LGG中的作用仍需要進(jìn)一步的驗證。RANBP17與細(xì)胞外基質(zhì)和細(xì)胞黏附功能相關(guān)，是膠質(zhì)瘤獨(dú)立的預(yù)后因素[22]，且RANBP17的表達(dá)水平與年齡呈負(fù)相關(guān)[23]，因此可以推斷RANBP17可能是在年齡和膠質(zhì)瘤之間發(fā)揮作用的中介因子。IGFBP2和RASL10A為抑癌基因。IGFBP2在膠質(zhì)瘤中通常為過表達(dá)[24]，且IGFBP2的瞬時表達(dá)對膠質(zhì)瘤進(jìn)展和生存有顯著影響，抑制IGFBP2也能抑制腫瘤細(xì)胞的擴(kuò)散[25]。LGG中的IDH突變促進(jìn)腫瘤抑制因子的信號傳導(dǎo)，同時抑制IGFBP2來改善患者的生存[26]，提示IGFBP2是膠質(zhì)瘤進(jìn)展的驅(qū)動因素和預(yù)后因素，是治療膠質(zhì)瘤的可行治療靶點。RASL10A又稱RRP22，其具有腫瘤抑制作用[27]，可以抑制腫瘤細(xì)胞生長，促進(jìn)非依賴性細(xì)胞死亡，降低侵襲性，抑制軟瓊脂中膠質(zhì)瘤細(xì)胞生長，從而起到抑癌功能[28]。

在11個基因中，F(xiàn)AM171A1、PLAT和MRPS16與膠質(zhì)瘤的關(guān)系尚不清楚，但有研究發(fā)現(xiàn)他們與其他癌癥有關(guān)。FAM171A1是一種調(diào)節(jié)人體細(xì)胞形狀和侵襲性生長的物質(zhì)，其表達(dá)上調(diào)可能會促進(jìn)癌細(xì)胞的侵襲性增長[29]。PLAT基因在乳腺癌和卵巢癌中發(fā)揮作用[30]。MRPS16是已知的與胚胎或未產(chǎn)期致死的基因[31]，截斷型MRPS16變體可能是胎兒死亡的潛在原因。本研究通過嚴(yán)格的因果推斷方法檢驗了其在年齡和LGG結(jié)局間的影響，可為LGG今后的研究提供新的生物學(xué)假設(shè)。

本研究僅控制了性別協(xié)變量，CIT方法在處理混雜因素上具有穩(wěn)健性，可彌補(bǔ)潛在混雜因素干擾的不足。VanderWeele中介模型在暴露和中介變量存在交互作用時，其效應(yīng)估計也較為穩(wěn)健。本文未研究年齡和基因表達(dá)間的交互作用，今后將進(jìn)一步探討。LGG預(yù)后隨訪數(shù)據(jù)的結(jié)局變量包括生存時間和是否死亡，但目前未見針對生存分析的CIT方法，今后將進(jìn)一步研究生存分析的CIT方法，提高模型的識別能力。

綜上所述，本研究基于因果推斷的方法識別出介導(dǎo)年齡和LGG患者結(jié)局的重要基因，年齡通過調(diào)控其基因表達(dá)水平進(jìn)而對LGG結(jié)局產(chǎn)生影響，為膠質(zhì)瘤患者的治療提供了重要的方法參考。