朱 淵， 何瑞瑞， 劉 源， 朱華慶*， 李 棟,2),*

(1)安徽醫(yī)科大學(xué)基礎(chǔ)醫(yī)學(xué)院， 合肥 230032；2)河北大學(xué)生命科學(xué)學(xué)院， 河北 保定 071002；3)軍事醫(yī)學(xué)研究院生命組學(xué)研究所，北京蛋白質(zhì)組研究中心，蛋白質(zhì)組學(xué)國(guó)家重點(diǎn)實(shí)驗(yàn)室， 北京 102206)

免疫系統(tǒng)中，不同的免疫功能由高度特化的不同類型的細(xì)胞執(zhí)行，各種細(xì)胞在免疫系統(tǒng)中的交互作用類似于一個(gè)社會(huì)網(wǎng)絡(luò)，而細(xì)胞因子則是連接這一網(wǎng)絡(luò)的重要分子。細(xì)胞因子是由免疫細(xì)胞和某些非免疫細(xì)胞經(jīng)刺激而合成、分泌的一類具有廣泛生物學(xué)活性的小分子多肽蛋白質(zhì)或糖蛋白質(zhì)。在免疫系統(tǒng)中，免疫細(xì)胞通過(guò)交換細(xì)胞因子實(shí)現(xiàn)細(xì)胞間的動(dòng)態(tài)通訊，從而產(chǎn)生各種適應(yīng)環(huán)境的反應(yīng)。細(xì)胞因子包括白細(xì)胞介素(interleukin, IL)[1, 2]、干擾素(interferon, IFN)、腫瘤壞死因子(tumor necrosis factor, TNF)[3]、集落刺激因子(colony stimulating factor, CSF)、趨化因子(chemokines)[4]和生長(zhǎng)因子(growth factor)[5]等，以自分泌、旁分泌、內(nèi)分泌3種方式作用于自身細(xì)胞或其他細(xì)胞，通過(guò)結(jié)合相應(yīng)的受體誘導(dǎo)細(xì)胞凋亡、調(diào)控細(xì)胞發(fā)育和分化、調(diào)控機(jī)體免疫應(yīng)答、介導(dǎo)炎癥反應(yīng)和促進(jìn)組織修復(fù)[6]。眾多細(xì)胞通過(guò)細(xì)胞因子在機(jī)體內(nèi)相互促進(jìn)或相互制約，形成極其復(fù)雜的細(xì)胞因子調(diào)節(jié)網(wǎng)絡(luò)。

目前，已發(fā)展了多種實(shí)驗(yàn)方法用于發(fā)現(xiàn)細(xì)胞-細(xì)胞因子的信息傳遞網(wǎng)絡(luò)[7, 8]，但這些方法通常存在周期長(zhǎng)、設(shè)備要求高和成本高等缺陷，有必要發(fā)展生物信息學(xué)方法助力細(xì)胞因子-細(xì)胞相互作用的發(fā)現(xiàn)。為此，Kveler等[9]發(fā)展了一種文本挖掘的方法，從廣泛疾病條件下的所有可用的PubMed摘要中提取340種細(xì)胞類型與140種細(xì)胞因子的關(guān)系，建立了目前規(guī)模最大的細(xì)胞-細(xì)胞因子相互作用數(shù)據(jù)集immuneXpresso (iX)。進(jìn)而通過(guò)聚類的方式來(lái)系統(tǒng)地預(yù)測(cè)細(xì)胞和細(xì)胞因子之間的相互作用關(guān)系。然而，iX數(shù)據(jù)集以及[9]文中預(yù)測(cè)的839對(duì)細(xì)胞-細(xì)胞因子的結(jié)果，同整個(gè)的免疫細(xì)胞-細(xì)胞因子相互作用空間(iX數(shù)據(jù)集提供，包括345種免疫細(xì)胞和143種細(xì)胞因子)相比仍然是極為有限，有必要發(fā)展更為高效的預(yù)測(cè)模型，進(jìn)一步推進(jìn)細(xì)胞-細(xì)胞因子相互作用的研究。細(xì)胞與細(xì)胞因子相互作用的預(yù)測(cè)本質(zhì)是一種關(guān)聯(lián)關(guān)系的研究，變分圖自編碼器(VGAE)作為一種無(wú)監(jiān)督學(xué)習(xí)框架，被廣泛應(yīng)用于鏈接預(yù)測(cè)任務(wù)[10-12]，其有望在細(xì)胞-細(xì)胞因子相互作用的預(yù)測(cè)中發(fā)揮重要作用。

VGAE將變分自編碼器(variational auto-encoder, VAE)中的變分思想引入到圖，充分利用圖在表征和計(jì)算領(lǐng)域的優(yōu)點(diǎn)[13]。圖是現(xiàn)實(shí)生活中廣泛存在的一類數(shù)據(jù)[14]，生物醫(yī)學(xué)領(lǐng)域中的分子網(wǎng)絡(luò)、化合物分子等都可以用圖來(lái)表示。利用VGAE中的圖卷積神經(jīng)網(wǎng)絡(luò)(graph convolutional neural network, GCN)能夠很好地融合圖的網(wǎng)絡(luò)拓?fù)涮卣骱凸?jié)點(diǎn)屬性特征[15, 16]。除此之外，自編碼器因其較強(qiáng)的特征學(xué)習(xí)能力、訓(xùn)練速度快、更少的中間參數(shù)等優(yōu)點(diǎn)而被廣泛使用[17]。

因此，本文提出一種基于VGAE預(yù)測(cè)細(xì)胞-細(xì)胞因子相互作用的深度學(xué)習(xí)模型——DeepCKI。此模型利用細(xì)胞因子組成的蛋白質(zhì)-蛋白質(zhì)相互作用(protein-protein interaction, PPI)網(wǎng)絡(luò)和蛋白質(zhì)不同類型的特征構(gòu)成圖，通過(guò)GCN學(xué)習(xí)圖的節(jié)點(diǎn)以及節(jié)點(diǎn)屬性之間的內(nèi)在規(guī)律和更加深層次的語(yǔ)義特征。利用已知的細(xì)胞-細(xì)胞因子相互作用和學(xué)習(xí)到節(jié)點(diǎn)嵌入向量訓(xùn)練深度神經(jīng)網(wǎng)絡(luò)分類器，進(jìn)而預(yù)測(cè)新的細(xì)胞-細(xì)胞因子之間的關(guān)聯(lián)。

1 材料與方法

1.1 細(xì)胞-細(xì)胞因子相互作用的數(shù)據(jù)集

本文從iX數(shù)據(jù)庫(kù)(http://www.immunexpresso.org)下載了文本挖掘獲得的人類細(xì)胞-細(xì)胞因子相互作用數(shù)據(jù)集。該數(shù)據(jù)集包括細(xì)胞作用于細(xì)胞因子以及細(xì)胞因子作用于細(xì)胞兩大類(Table 1)。其中，細(xì)胞因子作用于細(xì)胞可進(jìn)一步細(xì)分為：細(xì)胞因子正調(diào)控細(xì)胞(cytokine-cell+)、細(xì)胞因子負(fù)調(diào)控細(xì)胞(cytokine-cell-)、細(xì)胞因子調(diào)控細(xì)胞(cytokine-cell)。匯總所有數(shù)據(jù)集并刪除重復(fù)項(xiàng)，得到3 345對(duì)細(xì)胞-細(xì)胞因子相互作用。

Table 1 Details of cell-cytokine interaction datasets

1.2 圖的概念與定義

本文用圖來(lái)描述蛋白質(zhì)相互作用網(wǎng)絡(luò)。圖(graph)是由節(jié)點(diǎn)(node)和連接這些節(jié)點(diǎn)的邊(edge)組成的數(shù)據(jù)結(jié)構(gòu)，圖定義為G=(V, E)。其中，V表示節(jié)點(diǎn)的集合，E表示節(jié)點(diǎn)之間相連邊的集合。通用的圖表示是一個(gè)五元組：G(V, E, A, X, D)。其中，AN×N代表圖的鄰接矩陣，XN×F代表節(jié)點(diǎn)的特征矩陣，DN×N代表度矩陣，N和F分別代表節(jié)點(diǎn)的數(shù)量和節(jié)點(diǎn)特征向量的維度。

1.3 用于預(yù)測(cè)細(xì)胞-細(xì)胞因子相互作用的生物學(xué)證據(jù)

1.3.1 蛋白質(zhì)-蛋白質(zhì)相互作用網(wǎng)絡(luò) 本文從STRING(v11.0)[18]數(shù)據(jù)庫(kù)下載PPI(9606.protein.links.v11.0.txt.gz)數(shù)據(jù)。細(xì)胞-細(xì)胞因子相互作用數(shù)據(jù)集劃分為不同的類型。根據(jù)每種類型數(shù)據(jù)集中的細(xì)胞因子，篩選2個(gè)蛋白質(zhì)都是細(xì)胞因子且“combined score”大于300的高可靠的PPI，然后構(gòu)建鄰接矩陣AN×N表示蛋白質(zhì)之間的連接關(guān)系，而如果使用STRING數(shù)據(jù)庫(kù)中經(jīng)過(guò)實(shí)驗(yàn)驗(yàn)證的PPI，最多只能整合55對(duì)PPI，這么小的規(guī)模無(wú)法構(gòu)建有效的預(yù)測(cè)體系，所以，采用了全部證據(jù)來(lái)源的相互作用。蛋白質(zhì)vi和vj之間有相互作用，那么對(duì)應(yīng)的鄰接矩陣的元素Aij=1，否則Aij=0。鄰接矩陣對(duì)角元素通常設(shè)為0。

1.3.2 蛋白質(zhì)序列 從Swiss-Prot[19]數(shù)據(jù)庫(kù)中下載所有人類蛋白質(zhì)及其序列數(shù)據(jù)。采用聯(lián)合三聯(lián)體(conjoint triad,CT)[20]編碼將蛋白質(zhì)序列轉(zhuǎn)換為固定維度的嵌入向量，而包含模糊氨基酸(B、O、J、U、X、Z)的蛋白質(zhì)被刪去。常見(jiàn)的20種氨基酸根據(jù)偶極子和側(cè)鏈體積的差異被分成7類(Table 2)，所有同一類的氨基酸都被認(rèn)為是相同的。CT編碼嵌入向量的維度為7 × 7 × 7 = 343。任意3個(gè)連續(xù)的氨基酸作為1個(gè)三聯(lián)體單位，嵌入向量第i維值為三聯(lián)體在蛋白質(zhì)序列中出現(xiàn)的頻率。

Table 2 Classification of amino acids according to their dipoles and volumes of the side chains

CT定義為：

V=[n0,n1,…,nq-1]

其中,ni是每種三聯(lián)體在蛋白質(zhì)序列中出現(xiàn)的頻率，q等于343。

1.3.3 蛋白質(zhì)的亞細(xì)胞定位、結(jié)構(gòu)域和功能注釋 從Swiss-Prot數(shù)據(jù)庫(kù)中下載蛋白質(zhì)的亞細(xì)胞定位(subcellular localization)、結(jié)構(gòu)域(domain)和GO(gene ontology)注釋信息，采用詞袋模型(bag-of-words)[21]編碼上述數(shù)據(jù)，以蛋白質(zhì)的亞細(xì)胞定位為例，假設(shè)本文獲得n個(gè)不同的亞細(xì)胞定位，則亞細(xì)胞定位被編碼為一個(gè)長(zhǎng)度為n的二進(jìn)制向量，向量中的每個(gè)元素表示該蛋白質(zhì)是否被標(biāo)注了亞細(xì)胞定位。對(duì)于無(wú)任何亞細(xì)胞定位注釋的蛋白質(zhì)，它被表示為一個(gè)全為0的長(zhǎng)度為n的向量。采用詞袋模型編碼蛋白質(zhì)結(jié)構(gòu)域和GO條目時(shí)，為了避免維度災(zāi)難和降低復(fù)雜性，只使用在數(shù)據(jù)集中出現(xiàn)5次以上的結(jié)構(gòu)域條目(term)，以及10次以上的GO注釋。

1.4 DeepCKI模型架構(gòu)

DeepCKI整體架構(gòu)如Fig.1所示。模型主要分為3個(gè)模塊：輸入特征轉(zhuǎn)換模塊、VGAE模塊和DNN分類器模塊。3個(gè)模塊具體介紹如下：

Fig.1 The framework of DeepCKI model (A) Input feature conversion module. This module takes input as the human PPI network downloaded from the STRING database and protein amino acid sequences from the Swiss-Prot database. The PPI network is converted to the format of the adjacency matrix A, and a node attributes matrix X encoded by different protein features. (B) VGAE module. VGAE is an encoder-decoder approach. The encoder is two-layer graph convolutional network and the decoder is a dot product decoder. We obtained the mean embedding vector μ of cytokines for downstream classification tasks. (C) DNN module. The mean vector μ and cell-cytokines label matrix as input for the DNN classification, and the DNN classifier outputs the probabilities of the cell-cytokines interactions

1.4.1 輸入特征轉(zhuǎn)換模塊 如Fig.1A所示，此模塊以 STRING 數(shù)據(jù)庫(kù)中下載的人類PPI網(wǎng)絡(luò)和來(lái)自Swiss-Prot數(shù)據(jù)庫(kù)的蛋白質(zhì)氨基酸序列為輸入，PPI網(wǎng)絡(luò)被轉(zhuǎn)換為鄰接矩陣A的格式。CT方法編碼蛋白質(zhì)氨基酸序列，詞袋模型編碼亞細(xì)胞定位、結(jié)構(gòu)域和GO 注釋，轉(zhuǎn)換后的不同維度的二進(jìn)制向量作為特征矩陣X。

1.4.2 VGAE模塊 如Fig.1B所示，VGAE是一種無(wú)監(jiān)督特征提取方法，它基于網(wǎng)絡(luò)結(jié)構(gòu)和節(jié)點(diǎn)特征生成細(xì)胞因子的潛在表示。其一般包含編碼器和解碼器2個(gè)部分，基本思路是：編碼器采用2層GCN結(jié)構(gòu)，學(xué)習(xí)節(jié)點(diǎn)隱含向量表示的均值(μ)和方差(σ)的分布，然后從正態(tài)分布中采樣得到節(jié)點(diǎn)的向量表示zi，解碼器利用學(xué)習(xí)到的潛變量Z通過(guò)內(nèi)積的方式重構(gòu)鄰接矩陣A。編碼器與解碼器的詳細(xì)構(gòu)造如下：

1.4.2.1 編碼器

VGAE的編碼器由GCN組成。它以鄰接矩陣A和特征矩陣X作為輸入，生成潛在變量Z作為輸出。首先，2層GCN定義為：

(1)

Z=μ+σ⊙

(2)

q(zi|X,A) =N(zi|μi,diag(σi2))

(3)

1.4.2.2 解碼器

VGAE解碼器將使用學(xué)習(xí)到隱變量zi的內(nèi)積來(lái)重構(gòu)鄰接矩陣：

(4)

1.4.2.3 損失函數(shù)

L=Eq(Z|X,A)[logp(A|Z)-

KL(q(Z|X,A)||p(Z)]

(5)

1.4.3 DNN分類器模塊 正如Fig.1(C)所示，DNN分類模塊由3部分組成：輸入層(input layers)、隱藏層(hidden layers)與輸出層(output layer)。DNN分類器的輸入是每個(gè)細(xì)胞因子經(jīng)過(guò)GCN編碼后提取的均值向量μ和已知的細(xì)胞-細(xì)胞因子相互作用關(guān)系。隱藏層與輸出層包含功能神經(jīng)元，神經(jīng)元之間的“連接權(quán)”，以及每個(gè)功能神經(jīng)元的閾值會(huì)在學(xué)習(xí)過(guò)程中根據(jù)訓(xùn)練數(shù)據(jù)進(jìn)行調(diào)整。DNN分類器是一個(gè)多標(biāo)簽分類模型，輸出為m×n形狀的概率矩陣，每一行表示某個(gè)細(xì)胞因子與n個(gè)細(xì)胞發(fā)生相互作用的概率，每一列表示m個(gè)細(xì)胞因子與某類細(xì)胞發(fā)生相互作用的概率?；诟怕示仃嚺c細(xì)胞-細(xì)胞因子構(gòu)成的標(biāo)簽矩陣，計(jì)算出各個(gè)閾值下預(yù)測(cè)結(jié)果的真陽(yáng)性率和假陽(yáng)性率，進(jìn)而繪制n條ROC(Receiver Operating Characteristic)曲線，每條曲線對(duì)應(yīng)一類細(xì)胞與m個(gè)細(xì)胞因子的二分類結(jié)果。DeepCKI代碼的編寫(xiě)參考了Graph2GO[22]文章。

1.5 預(yù)測(cè)模型的性能評(píng)估

1.5.1 k-倍交叉驗(yàn)證 為了驗(yàn)證模型的效果，本文將所有已知細(xì)胞-細(xì)胞因子相互作用的數(shù)據(jù)集數(shù)據(jù)分為k份，k-1份作為訓(xùn)練集，1份作為驗(yàn)證集，依次輪換訓(xùn)練集和驗(yàn)證集k次。交叉驗(yàn)證可降低由于訓(xùn)練集和驗(yàn)證集單次劃分而導(dǎo)致的偶然性，提高泛化能力。

1.5.2 ROC曲線 ROC的全稱是“受試者工作特征”曲線[23]，該曲線能夠反映在不同閾值下敏感度與特異度之間的關(guān)系。一個(gè)理想的預(yù)測(cè)模型應(yīng)該同時(shí)具有較高的真陽(yáng)性率(true positive rate, TPR)和較低的假陽(yáng)性率(False Positive Rate, FPR)。ROC曲線圖中，每個(gè)點(diǎn)以對(duì)應(yīng)的FPR值為橫坐標(biāo)，以TPR值為縱坐標(biāo)。ROC曲線下面積(area under the curve of ROC, AUC)衡量模型優(yōu)劣的一種評(píng)價(jià)指標(biāo)。本文使用機(jī)器學(xué)習(xí)的Scikit-learn[24]包中的roc_curve和auc函數(shù)繪制ROC曲線，并計(jì)算AUC值。

1.5.3 精確率、召回率和F1-score 真陽(yáng)性(true positive, TP)指的是將陽(yáng)性樣本正確預(yù)測(cè)為陽(yáng)性樣本；假陽(yáng)性 (false positive, FP)指的是將陰性樣本錯(cuò)誤預(yù)測(cè)為陽(yáng)性樣本；真陰性(true negative, TN)指的是將陰性樣本正確預(yù)測(cè)為陰性樣本；假陰性(false negative, FN)指的是將陽(yáng)性樣本錯(cuò)誤預(yù)測(cè)為陰性樣本。本文同時(shí)使用精確率(precision)、召回率(recall)和F1-socre評(píng)估模型的預(yù)測(cè)性能。具體計(jì)算公式如下所示：

(6)

(7)

(8)

2 結(jié)果與討論

2.1 細(xì)胞-細(xì)胞因子預(yù)測(cè)模型DeepCKI的建立

本文通過(guò)整合細(xì)胞因子的PPI網(wǎng)絡(luò)和蛋白質(zhì)序列特征，建立了細(xì)胞-細(xì)胞因子預(yù)測(cè)模型——DeepCKI，該模型包含3個(gè)模塊，其中輸入特征轉(zhuǎn)換模塊實(shí)現(xiàn)蛋白質(zhì)相互作用網(wǎng)絡(luò)和蛋白質(zhì)屬性特征到編碼向量的轉(zhuǎn)換；VGAE模塊使用Adam優(yōu)化器進(jìn)行迭代地更新神經(jīng)網(wǎng)絡(luò)權(quán)重，學(xué)習(xí)率lr等于0.001，迭代次數(shù)epochs等于60，2層GCN隱藏層神經(jīng)元個(gè)數(shù)分別是800和400；DNN模塊的隱藏層為3層，每層神經(jīng)元的個(gè)數(shù)分別為1 024、512、256，在訓(xùn)練過(guò)程中對(duì)每層的輸入數(shù)據(jù)加1個(gè)批標(biāo)準(zhǔn)化處理(batch normalization)，隱藏層的激活函數(shù)為L(zhǎng)eakyRelu，添加Dropout防止模型過(guò)擬合，Adam優(yōu)化器實(shí)現(xiàn)對(duì)模型的優(yōu)化，二進(jìn)制交叉熵?fù)p失(binary crossentropy)作為損失函數(shù)，輸出層使用sigmod激活函數(shù)完成最終的多標(biāo)簽分類任務(wù)。本文利用已知細(xì)胞-細(xì)胞因子相互作用數(shù)據(jù)集訓(xùn)練DeepCKI模型，訓(xùn)練的目標(biāo)是不斷調(diào)節(jié)模型的超參數(shù)，提高模型預(yù)測(cè)的準(zhǔn)確率，降低損失函數(shù)的值。DeepCKI基于Tensorflow深度學(xué)習(xí)框架，代碼和數(shù)據(jù)集都已上傳到https://github.com/zhuyuan804/DeepCKI/tree/master。

2.2 特征選擇對(duì)DeepCKI模型預(yù)測(cè)效果的影響

為了驗(yàn)證不同類型的特征對(duì)模型預(yù)測(cè)性能的影響，進(jìn)而尋找最合適的特征或特征組合。本文使用蛋白質(zhì)序列、亞細(xì)胞定位、結(jié)構(gòu)域、GO 功能4種蛋白質(zhì)屬性特征來(lái)訓(xùn)練模型，從中選擇與細(xì)胞因子信息傳遞功能最相關(guān)的特征。其中，不同特征的嵌入向量以首尾橫向拼接的方式實(shí)現(xiàn)特征的組合。5倍交叉驗(yàn)證(5-fold cross-validation)評(píng)估模型預(yù)測(cè)效果，ROC曲線下面積評(píng)估蛋白質(zhì)不同屬性特征對(duì)DeepCKI模型預(yù)測(cè)性能的影響。

由Fig.2A可知，在只使用單類型特征的情況下，蛋白質(zhì)序列作為特征時(shí)模型的ROC曲線下面積等于0.8701，優(yōu)于亞細(xì)胞定位、結(jié)構(gòu)域、GO功能。在Fig.2(B)中，將序列特征分別與亞細(xì)胞定位和結(jié)構(gòu)域進(jìn)行組合，模型預(yù)測(cè)性能略低于只使用序列信息。雖然GO和序列的組合對(duì)模型的預(yù)測(cè)效果略有提升，但可以忽略不計(jì)?？傊A(yù)測(cè)結(jié)果說(shuō)明，蛋白質(zhì)序列特征比其他類型特征包含更多的信息量，多種特征的組合并不能顯著提升模型的預(yù)測(cè)性能，而僅使用序列特征就能達(dá)到,較好的預(yù)測(cè)效果。因此，本文僅使用蛋白質(zhì)序列特征對(duì)細(xì)胞-細(xì)胞因子相互作用進(jìn)行預(yù)測(cè)。

Fig.2 Prediction model performance comparison against using different protein attributes (A) Individual protein attributes sequence, subcellular location, protein domains and GO. (B) Combined protein attributes. sequence and location, sequence and domain, location and domain, sequence, location and domain, GO and sequence

2.3 DeepCKI模型與兩類經(jīng)典深度學(xué)習(xí)模型預(yù)測(cè)效果的比較

VAE與DNN是兩類經(jīng)典的深度學(xué)習(xí)算法，可以用于鏈接預(yù)測(cè)或分類任務(wù)。VAE無(wú)法讀入圖，DNN僅是對(duì)細(xì)胞-細(xì)胞因子相互作用數(shù)據(jù)進(jìn)行簡(jiǎn)單的分類訓(xùn)練，并不能捕捉細(xì)胞因子之間的內(nèi)在相似性，而DeepCKI可以用于圖，因此，本文建立的DeepCKI模型在預(yù)測(cè)細(xì)胞-細(xì)胞因子相互作用時(shí)，將會(huì)比上述2種經(jīng)典深度學(xué)習(xí)模型更有優(yōu)勢(shì)。為了驗(yàn)證這一假設(shè)，本文利用相同的去除重復(fù)的3 345對(duì)細(xì)胞-細(xì)胞因子相互作用數(shù)據(jù)集對(duì)3個(gè)模型的預(yù)測(cè)性能進(jìn)行比較，并利用5倍交叉驗(yàn)證評(píng)估模型預(yù)測(cè)性能。

從Fig.3的結(jié)果可知，基于變分圖自編碼器的DeepCKI模型取得了最優(yōu)的結(jié)果(AUC=0.8701)，遠(yuǎn)高于VAE模型(AUC=0.7637)和DNN模型(AUC=0.6960)。DeepCKI模型預(yù)測(cè)性能得到顯著提高的原因是因?yàn)閂AGE相對(duì)于VAE而言增加了對(duì)圖特征的提取，編碼器的GCN結(jié)構(gòu)從輸入的網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)和蛋白質(zhì)序列中捕獲細(xì)胞因子包含的深層次信息。結(jié)果可知，將VGAE學(xué)習(xí)到細(xì)胞因子的隱含表示均值向量μ輸入到DNN分類器，比簡(jiǎn)單使用DNN進(jìn)行分類更為有效。

Fig.3 Comparison of DeepCKI, VAE and DNN Each model was evaluated using 5-fold cross-validation. DeepCKI has the highest area under ROC curve (AUC=0.8701) among the three models

2.4 DeepCKI模型在不同類型細(xì)胞-細(xì)胞因子數(shù)據(jù)集上的表現(xiàn)

細(xì)胞-細(xì)胞因子之間具有多種類型的關(guān)聯(lián)關(guān)系(詳情可查看Table 1)，其中包括：細(xì)胞正調(diào)控細(xì)胞因子(例如，巨噬細(xì)胞和T細(xì)胞產(chǎn)生IL-10促進(jìn)淋巴細(xì)胞性脈絡(luò)叢腦膜炎病毒(LCMV)克隆感染的慢性[25])、細(xì)胞因子正調(diào)控細(xì)胞(例如IL-2在體外能有效誘導(dǎo)T細(xì)胞擴(kuò)增[26])、細(xì)胞因子負(fù)調(diào)控細(xì)胞(例如IL-21和IL-12抑制宮頸癌患者調(diào)節(jié)性T細(xì)胞的分化[27])、細(xì)胞因子調(diào)控細(xì)胞(例如IL-7, IL-4, IL-6，IL-10作用于B細(xì)胞[28])。為了考察DeepCKI模型預(yù)測(cè)不同類型細(xì)胞-細(xì)胞因子預(yù)測(cè)的能力，本文使用Table 1中不同類型的細(xì)胞-細(xì)胞因子相互作用數(shù)據(jù)對(duì)DeepCKI模型進(jìn)行訓(xùn)練，采用AUC、精確率、召回率和F1-score共4種評(píng)估指標(biāo)評(píng)估模型預(yù)測(cè)性能。由Fig.4可知，雖然在擁有最多PPI數(shù)據(jù)的cytokin-cell得不到最高的召回率和F1-score，但他的精確率和F1-score值最高。模型在4種類型數(shù)據(jù)集上的AUC值均在0.8 以上，證明了模型具有較好的預(yù)測(cè)能力。此外，本文發(fā)現(xiàn)在,4類相互作用中，已知PPI越多，該類AUC值越高，表明數(shù)據(jù)集越完整，算法可以學(xué)習(xí)到更多的信息，預(yù)測(cè)更加準(zhǔn)確。本文使用cell-cytokine+數(shù)據(jù)集將該模型與上述2種深度學(xué)習(xí)模型(VAE和DNN)，以及17種來(lái)源于Pykeen[29]包中基于知識(shí)圖譜用來(lái)預(yù)測(cè)節(jié)點(diǎn)之間關(guān)系的嵌入方法。結(jié)果正如Fig.5所示，DeepDKI具有最優(yōu)的預(yù)測(cè)性能，進(jìn)一步證明了模型的高效預(yù)測(cè)性能。

Fig.4 AUC, precision, recall, and F1-Score values of DeepCKI model against four different datasets during 5-fold cross validation The AUC values of the model on the four types of data sets are all above 0.8, which proves that the model has good predictive ability

Fig.5 The performance of the DeepCKI and 19 other methods on the cell-cytokine interaction datasets These methods include two deep learning models (VAE and DNN). Seventeen embedding methods from the pykeen package are used to predict relationships between nodes. DeepDKI has the best predictive performance

2.5 DeepCKI模型預(yù)測(cè)細(xì)胞-細(xì)胞因子網(wǎng)絡(luò)

為了進(jìn)一步驗(yàn)證DeepCKI模型的預(yù)測(cè)能力，本文使用3 345對(duì)細(xì)胞-細(xì)胞因子相互作用數(shù)據(jù)集對(duì)DeepCKI模型進(jìn)行訓(xùn)練，在全蛋白質(zhì)組范圍內(nèi)進(jìn)行細(xì)胞-細(xì)胞因子的關(guān)聯(lián)預(yù)測(cè)。刪除與已知數(shù)據(jù)集重疊的部分，本文預(yù)測(cè)獲得12 410對(duì)細(xì)胞-細(xì)胞因子相互作用(結(jié)果未顯示)。利用Cytoscape[30]開(kāi)源軟件對(duì)預(yù)測(cè)分值排名前100的細(xì)胞-細(xì)胞因子相互作用進(jìn)行可視化展示，結(jié)果正如Fig.5展示，涉及21種細(xì)胞與41種細(xì)胞因子，其中有36對(duì)在最新發(fā)表的文獻(xiàn)中被證實(shí)，表明DeepCKI具有發(fā)現(xiàn)新的細(xì)胞-細(xì)胞因子相互作用的關(guān)系的能力。

結(jié)果正如Fig.6所示，巨噬細(xì)胞(macrophage)和T細(xì)胞具有最多的細(xì)胞因子數(shù)量。文獻(xiàn)證明，巨噬細(xì)胞存在于大部分組織中，是一類具有吞噬作用的天然免疫細(xì)胞，參與細(xì)胞碎片和病原體的識(shí)別、吞噬和降解，在炎癥和宿主防御中以及疾病進(jìn)展中發(fā)揮至關(guān)重要的作用[31]。當(dāng)巨噬細(xì)胞受到刺激時(shí)，通過(guò)釋放TNF、IL-1、IL-6、IL-20、趨化因子(C-X-C motif)配體9 (chemokine (C-X-C motif) ligand 9, CXCL9)等細(xì)胞因子發(fā)揮重要調(diào)節(jié)作用。其中，TNF是引起感染性休克的主要細(xì)胞因子之一，在下丘腦中TNF刺激促腎上腺皮質(zhì)釋放激素的釋放，抑制食欲，誘導(dǎo)發(fā)燒[32]。IL-1是一個(gè)多效應(yīng)的細(xì)胞因子，在炎癥期間IL-1可刺激肝產(chǎn)生急性期蛋白(acute phase proteins)，并作用于中樞神經(jīng)系統(tǒng)，誘導(dǎo)發(fā)熱和前列腺素分泌[33]。IL-6具有促炎和抗炎的雙向調(diào)節(jié)功能，影響從免疫到組織修復(fù)和新陳代謝的調(diào)節(jié)過(guò)程[32]。IL-20可以通過(guò)調(diào)節(jié)脂肪生成和巨噬細(xì)胞失調(diào)參與肥胖[32, 34]，促進(jìn)白細(xì)胞和上皮細(xì)胞之間的細(xì)胞通訊[35]。巨噬細(xì)胞來(lái)源的CXCL9為免疫檢查點(diǎn)阻斷后的抗腫瘤免疫應(yīng)答所必需[36]。T細(xì)胞是免疫系統(tǒng)中核心組分及效應(yīng)細(xì)胞，具有多種細(xì)胞亞型，包括CD8+T細(xì)胞，CD4+ T細(xì)胞以及T調(diào)節(jié)細(xì)胞等。CD8+T細(xì)胞可將細(xì)胞毒素釋放到受感染的細(xì)胞，導(dǎo)致細(xì)胞死亡。腫瘤微環(huán)境中的CD8+T細(xì)胞可生成IL-2、IL-12和IFNγ，靶向殺死腫瘤細(xì)胞。CD4+ T細(xì)胞分泌的CCL2可召集免疫調(diào)節(jié)細(xì)胞在感染部位集聚。調(diào)節(jié)性T細(xì)胞在活化時(shí)可分泌免疫調(diào)節(jié)因子，例如IL-10、轉(zhuǎn)化生長(zhǎng)因子β和IL-35。不同細(xì)胞因子也參與調(diào)節(jié)T細(xì)胞的功能，CXCL1通過(guò)調(diào)節(jié)T細(xì)胞功能來(lái)幫助宿主防御微生物敗血癥[37]，CCL8誘導(dǎo)活化的T細(xì)胞趨化[37, 38]，MIF(macrophage migration inhibitory factor)抑制T細(xì)胞活化[39]。細(xì)胞因子IL23A可以同最多類型的細(xì)胞關(guān)聯(lián)。一方面，細(xì)胞因子IL-23A可以由多種免疫細(xì)胞分泌，例如巨噬細(xì)胞[40]，CD4陽(yáng)性T細(xì)胞(CD4-positive T cell)[41]和B細(xì)胞[42]。另一方面，IL-23A可以調(diào)控多種細(xì)胞的功能，例如可使極化T細(xì)胞向不同的效應(yīng)功能表型分化，參與活化T細(xì)胞增殖的正向調(diào)節(jié)，抑制調(diào)節(jié)性T細(xì)胞活性，導(dǎo)致T細(xì)胞依賴性結(jié)腸炎[43]。Fig.5網(wǎng)絡(luò)的分析結(jié)果與文獻(xiàn)報(bào)道一致，其一定程度上證明了DeepCKI模型預(yù)測(cè)結(jié)果的可靠性。

Fig.6 Insights into the top 100 predicted cell-cytokine interaction network Orange and purple nodes, respectively, stand for cell and cytokine. Edges represent cell and cytokine interactions. Solid lines indicate that the predicted cell-cytokine interactions have been validated in the literature, and dotted lines indicate that they have not been validated. " (+)” denotes cytokines that promote cells, e.g. induces"; "(-)" indicates that cytokines inhibit cells, e.g. "decreased"; "(·)" means cytokines act on cells, a neutral regulation, e.g. "correlate”; (*) indicates that cells secrete cytokines

了解免疫細(xì)胞間網(wǎng)絡(luò)的通訊對(duì)理解疾病中的免疫反應(yīng)至關(guān)重要。本文構(gòu)建的預(yù)測(cè)細(xì)胞-細(xì)胞因子相互作用的深度學(xué)習(xí)模型加速推進(jìn)對(duì)細(xì)胞間通訊的系統(tǒng)研究。該模型采用圖結(jié)構(gòu)設(shè)計(jì)，以細(xì)胞因子組成的鄰接矩陣A和蛋白質(zhì)序列組成的特征矩陣X作為輸入，通過(guò)整合網(wǎng)絡(luò)節(jié)點(diǎn)和節(jié)點(diǎn)特征生成潛在表示，以DNN分類器輸出預(yù)測(cè)的細(xì)胞-細(xì)胞因子相互作用。通過(guò)對(duì)不同蛋白質(zhì)特征進(jìn)行篩選，發(fā)現(xiàn)僅使用序列作為初始特征時(shí)模型的預(yù)測(cè)性能最優(yōu)，說(shuō)明蛋白質(zhì)序列特征比其他類型特征包含更多的信息量。此外，具有GCN結(jié)構(gòu)的DeepCKI模型性能，優(yōu)于無(wú)GCN結(jié)構(gòu)的變分自編碼器和不經(jīng)過(guò)編碼直接使用特征信息進(jìn)行分類的深度神經(jīng)網(wǎng)絡(luò)，說(shuō)明GCN能從圖中學(xué)習(xí)到節(jié)點(diǎn)和邊的內(nèi)在規(guī)律，以及節(jié)點(diǎn)屬性蘊(yùn)含更加深層次的特征信息。此外，該模型在不同類型數(shù)據(jù)集的訓(xùn)練中均表現(xiàn)出了魯棒性和有效性，并且新預(yù)測(cè)的細(xì)胞和細(xì)胞因子，通過(guò)查閱文獻(xiàn)證明其與多種疾病的發(fā)生發(fā)展密切相關(guān)。

綜上結(jié)果可以推測(cè)，DeepCKI模型具有發(fā)現(xiàn)新的細(xì)胞-細(xì)胞因子相互作用的能力，有助于為大規(guī)模的細(xì)胞-細(xì)胞因子的實(shí)驗(yàn)研究提供一定的理論指導(dǎo)。但模型尚不能實(shí)現(xiàn)對(duì)細(xì)胞-細(xì)胞因子激活/抑制關(guān)系及方向性的預(yù)測(cè)。未來(lái)，可在此方面進(jìn)一步拓展模型的性能，同時(shí)對(duì)模型的可解釋性進(jìn)行深入研究。