張偉偉， 吳志強(qiáng)

(東華理工大學(xué) 理學(xué)院，江西 南昌 330013)

癌癥是一類由于細(xì)胞分裂和凋亡機(jī)制失常而導(dǎo)致的疾病，通常表現(xiàn)為惡性腫瘤(Hanahan et al., 2011)。差異表達(dá)的基因可能直接或間接關(guān)系到癌癥的發(fā)生和易感性。差異表達(dá)基因的檢測能夠從基因?qū)哟谓沂景┌Y的發(fā)生機(jī)理并且能夠更好地探索癌癥治療的新方法(Wu et al., 2013)。

腫瘤樣本的差異表達(dá)分析通常是將腫瘤樣本和正常樣本直接進(jìn)行比對。然而，臨床上獲得的腫瘤組織往往包含不同比例的正常細(xì)胞，如果不考慮正常細(xì)胞的混入，直接進(jìn)行差異甲基化/基因表達(dá)或樣本聚類等下游分析，勢必會導(dǎo)致結(jié)果出現(xiàn)偏差甚至錯誤。在表觀遺傳層面上，針對腫瘤純度導(dǎo)致差異甲基化分析出現(xiàn)偏差的問題，Zheng等(2017)提出了腫瘤純度校正的差異甲基化分析模型 InfiniumDMC。Zhang等(2017)提出了一種考慮腫瘤細(xì)胞純度的腫瘤樣本亞型聚類方法 InfiniumClust。同樣，腫瘤樣本的“不純”也會使腫瘤樣本在基因上的均值估計產(chǎn)生偏差，進(jìn)而導(dǎo)致差異表達(dá)分析結(jié)果出現(xiàn)偏差。為此Aran等(2015)將腫瘤純度作為線性“協(xié)變量”加入回歸模型，盡管這樣做得到了更多差異表達(dá)的基因。但是，通過分析發(fā)現(xiàn)，腫瘤純度與基因差異之間的關(guān)系是乘性的而非原來認(rèn)為的線性關(guān)系。忽視腫瘤純度，或者將腫瘤純度作為“協(xié)變量”加入回歸模型都會使差異分析結(jié)果出現(xiàn)偏差。

為解決上述問題，本文在正態(tài)分布假設(shè)下構(gòu)建了腫瘤樣本基因表達(dá)量的混合概率模型，通過廣義最小二乘法和Wald檢驗來確定每個基因在癌癥和正常細(xì)胞之間的差異性，并通過TCGA數(shù)據(jù)分析來說明所提模型的優(yōu)越性。

1 考慮腫瘤細(xì)胞純度的差異表達(dá)基因模型

基因的原始數(shù)據(jù)不服從正態(tài)分布，經(jīng)過log2變換后數(shù)據(jù)的正態(tài)性會變強(qiáng)。圖1以肺腺癌為例，圖1a顯示肺腺癌正常樣本變換前后基因正態(tài)分布擬合檢驗的統(tǒng)計量，每個點代表一個基因，顯然變換后的統(tǒng)計量變得更加小，特別是對于那么違反正態(tài)假設(shè)的基因，因此log2變換有助于數(shù)據(jù)正態(tài)化。其次，對變換后的正常樣本進(jìn)行正態(tài)分布擬合檢驗，圖1b顯示大部分p值均分分布，這表明正態(tài)分布假設(shè)對于絕大多數(shù)基因是成立的。雖然有一小部分基因具有較小的p值，但這些基因可能是正常樣本中不同細(xì)胞的差異表達(dá)基因。因此，在以下模型假設(shè)中可以將每個基因在正常樣本中的表達(dá)量假設(shè)成一個單一的正態(tài)分布。

圖1 正常樣本上所有基因正態(tài)分布檢驗的p值和統(tǒng)計量Fig.1 P-values and statistics for the normal distribution of all genes in normal samples

Yis′ =(1-λs)Xis+λsYis

=(1-λs)Xis+λs(Xis+δis)

=Xis+λsδis

令Zi=[Xi1,Xi2,…,Xin0,Yi1′,Yi2′…,Yin1′]。其中前n0項為正常樣本的基因表達(dá)量，后n1項為腫瘤樣本的基因表達(dá)量。上述問題變?yōu)榫€性模型Zis=mi+asμi+s,s=1,2,…,n0+n1。模型的參數(shù)可以通過廣義最小二乘法來求解，令

把腫瘤純度作為線性模型的一個實驗設(shè)計因子引入，會降低腫瘤樣本的內(nèi)部方差，因此，從理論上來說對于差異表達(dá)基因的估計會有所提升。選取 InfiniumPurify估計的腫瘤細(xì)胞純度數(shù)據(jù)，并將所提方法命名為LinReg，通過與limma (Ritchie et al., 2015)在TCGA的14種癌癥類型(正常樣本數(shù)量超過10個的癌癥類型)上的差異分析比較來說明所提方法的優(yōu)越性。

2 TCGA腫瘤樣本差異表達(dá)基因分析結(jié)果

圖2 兩種方法在不同癌癥類型中各指標(biāo)比對圖Fig.2 Comparison diagram of the two methods in different cancer types

首先，把每一種癌癥類型的腫瘤樣本隨機(jī)平均分成兩份，將這兩份腫瘤樣本分別和正常樣本利用limma和LinReg進(jìn)行差異表達(dá)基因分析(圖2)。在同樣的顯著性指標(biāo)下(FDR < 0.001)，對每一種方法將這兩組探測到的差異表達(dá)基因取交集，并且來比較一下交集中所包含的基因的數(shù)量。圖2a顯示對于絕大多數(shù)TCGA的癌癥類型，線性模型比limma得到了更多的交集基因，表明線性方法相比較limma具有更好的探測出診斷生物標(biāo)志物的性能。接下來比較這兩種方法在同樣的顯著性指標(biāo)下(FDR < 0.001)得到的差異表達(dá)的基因數(shù)量。如圖2b所示，對于絕大多數(shù)TCGA的癌癥類型，線性模型比limma得到了更多的差異表達(dá)的基因，表明線性方法敏感性更強(qiáng)。雖然不同癌癥類型具有不同的病理學(xué)基礎(chǔ)和差異表達(dá)的基因，但從整體上看，不同癌癥類型應(yīng)具有一些相同點，這些共有的性質(zhì)會帶來癌癥類型之間較高的統(tǒng)計量相關(guān)性。圖2c中列舉了不同癌癥類型與其他癌癥類型的平均統(tǒng)計量相關(guān)性，結(jié)果同樣表明線性方法比limma的相關(guān)性更強(qiáng)，說明其更具一致性。最后，來分析一下差異表達(dá)基因的生物學(xué)意義。選取FDR最小的前2 000個基因，將這2 000個基因與KEGG的“PATHWAY_IN_CANCER”上的328個基因取交集，來看一下這14種癌癥上兩種方法的交集基因個數(shù)。圖2d顯示線性方法在BLCA，BRCA，ESCA，KICH，KIRC，LIHC，LUSC，PRAD，STAD和UCEC這10種癌癥上具有更多的交集中的基因，其中LinReg在UCEC上的交集基因的個數(shù)幾乎是limma的5倍。這些結(jié)果說明線性方法得到的基因更具有生物學(xué)意義，這些基因?qū)τ诎┌Y診斷和治療具有重要意義。

圖3 只被線性方法探測到的差異表達(dá)的基因與腫瘤純度的關(guān)系Fig.3 The relationship between differentially expressed genes only detected by linear method and tumor purities

為了更好地理解腫瘤純度在差異表達(dá)基因分析中所起的作用，本文任意選取了一個只能被LinReg探測到差異的基因來研究該問題。圖3左面面板顯示該基因log2變換后在腫瘤樣本和正常樣本的表達(dá)量分布，該基因的均值在兩組之間是沒有差異的，limma檢驗的p值為0.872；中間的面板顯示該基因在腫瘤樣本上log2變換后的表達(dá)量和腫瘤樣本純度的相關(guān)系數(shù)，相關(guān)系數(shù)達(dá)-0.67；右面面板顯示校正腫瘤純度后，兩組之間均值有明顯差異，利用LinReg進(jìn)行檢驗發(fā)現(xiàn)p=9.41×10-33。這些結(jié)果說明腫瘤樣本間較大的組內(nèi)方差是由于腫瘤樣本的純度引起的，因此利用癌旁的正常樣本進(jìn)行線性提純后，發(fā)現(xiàn)該基因在兩組之間具有明顯差異，這個結(jié)果進(jìn)一步說明了差異分析時校正腫瘤純度的必要性。

3 結(jié)論

本文針對腫瘤樣本“不純”導(dǎo)致差異分析結(jié)果出現(xiàn)偏差的問題，構(gòu)建了腫瘤樣本基因表達(dá)量的混合概率模型。該模型是把腫瘤純度作為影響差異表達(dá)基因分析的因素引入差異分析，相比直接使用腫瘤樣本作為癌癥細(xì)胞樣本與正常樣本直接比對的方法更加合理，有效避免了因為腫瘤樣本“不純”而導(dǎo)致差異表達(dá)分析出現(xiàn)的偏差。同樣腫瘤純度也影響著基因共表達(dá)網(wǎng)絡(luò)和共甲基化網(wǎng)絡(luò)，計劃將

在進(jìn)一步工作中構(gòu)建腫瘤純度校正的基因共表達(dá)網(wǎng)絡(luò)和共甲基化網(wǎng)絡(luò)，并分析其與傳統(tǒng)的共表達(dá)/共甲基化網(wǎng)絡(luò)在功能模塊分析中的差異。