李榮遠 龍法寧

1(廣西師范大學(xué)計算機科學(xué)與工程學(xué)院 廣西 桂林 541004)

2(玉林師范學(xué)院計算機科學(xué)與工程學(xué)院 廣西 玉林 537000)

0 引 言

鑒定細胞類型及亞型成為scRNA-seq重要的應(yīng)用之一,大量的半監(jiān)督、無監(jiān)督聚類方法被開發(fā)出來。基于聚類的分類方法假定聚類中所有細胞均屬于同一類型,因此可以進行集群標(biāo)記。但這種假設(shè)通常是錯誤的,集群中除了主要細胞類型外,通常還包含少量占比的多種細胞類型[1]。常用的聚類算法如k-means、hierarchical clustering需要設(shè)置相應(yīng)的類別數(shù),其中類別數(shù)的設(shè)置對聚類結(jié)果影響較大。因此,一種無須先進行聚類就可對每個細胞進行分類的方法可解決該問題。ScRNA-seq測序技術(shù)由于本身技術(shù)噪聲、批次效應(yīng)等問題,導(dǎo)致下游分析困難,且單細胞表達譜維度較高,至少上萬維,普通方法難以分辨細胞類型。超圖神經(jīng)網(wǎng)絡(luò)能較好地處理大規(guī)模多模態(tài)數(shù)據(jù)集[2]。基于此,本文嘗試將特征工程結(jié)合超圖神經(jīng)網(wǎng)絡(luò)應(yīng)用于單細胞測試數(shù)據(jù)集上來證明方法的有效性。

1 相關(guān)工作

繼《基因組計劃》后,Regev等[3]提出《人類細胞圖譜計劃》,該計劃旨在繪制人類眾多細胞類型及狀態(tài)。計劃描述出人類每個細胞,探索先前未知細胞,如研究細胞類型、細胞之間關(guān)系和細胞組成成分等,因此多種單細胞RNA測序技術(shù)(ScRNA-seq)蜂擁而至。目前主流的微流控技術(shù)有10X Genomics和微孔板技術(shù)[4]等,其測序通量較高,能對單細胞全長測序,基因檢測率較高。ScRNA-seq能測序出單細胞整個轉(zhuǎn)錄組的基因表達,方便各種方法進行細胞分類,識別細胞亞型,探索細胞的異質(zhì)性[5]。

目前單細胞主流研究方向主要有:轉(zhuǎn)錄組學(xué)、基因組學(xué)和軌跡推斷[6]等分支。本文主要研究單細胞轉(zhuǎn)錄組方向,目前面臨一些挑戰(zhàn),單個細胞轉(zhuǎn)錄狀態(tài)的所有表征能全面了解細胞內(nèi)RNA的相互作用,但scRNA-seq觀測值零值較多,給定的基因中沒有唯一的分子標(biāo)識符。例如,缺失率(Dropout)通常被描述為scRNA-seq數(shù)據(jù)中的零值,但是該項觀測通常將兩種不同零值類型混為一起。一是歸因于噪聲,其基因已經(jīng)表達但未被測序技術(shù)檢測到的零值;二是歸因于生物學(xué)上真正的零值。不建議將Dropout作為觀察值為零的總稱[7-8],這些零值歸因于技術(shù)限制:可能是生物變異、細胞裂解和人工操作等因素。較多零值使高維數(shù)據(jù)變得稀疏,稀疏度取決于所用的scRNA-seq平臺、測序深度和基因自身表達水平。ScRNA-seq數(shù)據(jù)的稀疏性可能會阻礙下游分析,難以正確建模或處理,因此需要進一步開發(fā)新方法。M3Drop[9]實現(xiàn)兩種針對零值的處理方法。第一種方法適合于全轉(zhuǎn)錄組測序協(xié)議產(chǎn)生的數(shù)據(jù);第二種適合數(shù)據(jù)集符合負(fù)二項分布ZINB(零膨脹負(fù)二項式)模型NBDrop,它能提取出具有較高零值的數(shù)據(jù)特征。兩種基于Dropout的特征選擇方法NBDrop和M3Drop均比基于方差的特征選擇方法好,但不適合小樣本量或高噪聲數(shù)據(jù)。

特征工程在總體預(yù)測性能中起關(guān)鍵作用,包括質(zhì)量控制(刪除細胞、基因)、規(guī)范化處理、基因選擇(選高變異基因)或降維,對聚類或分類模型起關(guān)鍵作用。文獻[10]提出3種基因選擇方法及14種聚類算法分析。其中SC3[11]和Seurat[12]表現(xiàn)出最好結(jié)果,但SC3運行大數(shù)據(jù)集耗時較長,效率較低。Seurat在處理scRNA-sq數(shù)據(jù)集上效果較好,其中特征工程中選出高度變異的基因是其關(guān)鍵。文獻[1]提出scPred方法,這是一種新的可推廣方法,該方法結(jié)合基于降維的特征選擇和機器學(xué)習(xí)的預(yù)測方法,對單細胞進行高度準(zhǔn)確的分類。但該方法局限于一些特定數(shù)據(jù)集,對多種單細胞測序平臺兼容效果差,分類準(zhǔn)確率低。

基于超圖的思想在單細胞數(shù)據(jù)集方面應(yīng)用較少,SAME[13]使用超圖將多種方法進行聚類集成,將聚類標(biāo)簽作為輸入,從而構(gòu)建一種共識機制。SAFE[14]集成四種最先進的聚類方法:SC3、CIDR、Seurat和t-SNE+k-means,該方法運行時間開銷較大。

基于標(biāo)記的單細胞分類必不可少,超圖學(xué)習(xí)[15]在處理大規(guī)模數(shù)據(jù)集有較好效果。超圖處理復(fù)雜數(shù)據(jù)更靈活,HGNN根據(jù)超邊卷積學(xué)習(xí)數(shù)據(jù)之間的相關(guān)性,有效地進行傳統(tǒng)的超圖學(xué)習(xí)。超高維多模態(tài)數(shù)據(jù)符合這種模型處理[16],scRNA-seq數(shù)據(jù)集正是一種超高維數(shù)據(jù)集。

綜上所述,單細胞多核超圖分類整體流程如圖1所示。

圖1 單細胞多核超圖分類整體流程

(1) 數(shù)據(jù)預(yù)處理。收集多種平臺數(shù)據(jù)集,將單細胞數(shù)據(jù)經(jīng)過預(yù)處理后,選出高度變異的基因(High Variable Gene,HVG)。

(2) 多核KNN圖構(gòu)建。計算細胞之間各相似距離,根據(jù)距離構(gòu)建各細胞K近鄰圖,將多核KNN圖合并為超圖。

(3) 構(gòu)建超圖神經(jīng)網(wǎng)絡(luò)學(xué)習(xí)HGNN(Hypergraph Neural Network)。通過點-邊-點特征表示,構(gòu)建超圖學(xué)習(xí)分類器。

(4) 實驗對比驗證。采用多種數(shù)據(jù)類型進行驗證三種方法。主流分類方法有scPred、HGNN和HVG-MHGNN(Multi-kernel Hypergraph Neural Network based on High Variable Gene)。

2 scRNA-seq的HVG-MHGNN模型構(gòu)建

2.1 scRNA-seq數(shù)據(jù)集整理

ScRNA-seq測序平臺多樣,為驗證實驗方法,收集了多種類型數(shù)據(jù)。數(shù)據(jù)集主要來源conquer數(shù)據(jù)庫、10X、GEO、ArrayExpress。Trapnell、TrapnellTCC和Petropoulos來源于conquer。整理多種類型數(shù)據(jù)集工作量較大,一些作者只提供原始ATGC序列數(shù)據(jù)集,一般都是上百GB,普通平臺難以分析。本次收集作者處理過的基因表達值,最終形成count統(tǒng)計量類型數(shù)據(jù)。由于一些數(shù)據(jù)集作者標(biāo)簽隱藏在數(shù)據(jù)庫中,需要查看作者原始論文人工提取標(biāo)簽,以方便驗證。如Petropoulos數(shù)據(jù)集從88個人類植入前胚胎中獲得了1 529個單細胞RNA-seq,觀察胚胎從3天到第7天發(fā)育情況。其他數(shù)據(jù)集如表1所示。

表1 多種數(shù)據(jù)類型scRNA-seq數(shù)據(jù)集

2.2 數(shù)據(jù)預(yù)處理

單細胞數(shù)據(jù)預(yù)處理已開發(fā)多種優(yōu)秀程序包,如scater質(zhì)量控制,Seurat選出高變異基因(HVG)。Scater實現(xiàn)刪除低質(zhì)量細胞及基因(零值、ERCC和pink-in),Seurat能選出高度變異的基因。這些預(yù)處理能使單細胞數(shù)據(jù)降低一定維度,但一般也接近上千維,接著用經(jīng)典的PCA降維,scVI通過深度自編碼器學(xué)習(xí)達到降維作用,對不同的數(shù)據(jù)集有一定效果。本實驗為體現(xiàn)超圖在單細胞數(shù)據(jù)上的表現(xiàn)效果,只使用Seurat提取高度變異基因,進行超圖學(xué)習(xí)。

2.3 相似距離

細胞之間的距離是聚類的核心關(guān)鍵,歐氏距離被廣泛應(yīng)用于各分類和聚類方法中,但歐氏距離有一定的局限性,細胞聚類效果并不好。Pearson和Spearmam相關(guān)性對細胞分類或聚類效果較好,其范圍在[-1,1]之間,值越大,相關(guān)性越高,與距離成反比。本次使用式(1)和式(2)計算細胞之間的Pearson和Spearman相似距離。單一距離具有一定的偶然性,本次將三者進行結(jié)合,對構(gòu)造超圖具有明顯的有效性。

(1)

(2)

2.4 超圖神經(jīng)網(wǎng)絡(luò)

圖2 細胞超圖表示

(1) 構(gòu)建頂點和邊之間關(guān)系。

超圖頂點和邊之間的關(guān)系如式(3)所示,如果節(jié)點與節(jié)點之間相連接,用1表示,否則為0。

(3)

(2) 構(gòu)建超圖目標(biāo)函數(shù)。

考慮超圖每個頂點的下游分類問題,每個頂點的標(biāo)簽應(yīng)該能夠應(yīng)用到超圖結(jié)構(gòu)中,整個學(xué)習(xí)的目標(biāo)函數(shù)用式(4)描述。

(4)

式中:Remp(f)表示監(jiān)督的損失;f(·)表示分類函數(shù);Ω(f)表示超圖的規(guī)范化,定義如式(5)所示。

(5)

(3) 通過超邊卷積獲取特征Y。

超圖卷積由兩個子模塊組成:頂點卷積子模塊和超邊卷積子模塊。頂點卷積將頂點特征集合到上邊緣,然后上邊緣卷積將相鄰的上邊緣特征集合到形心頂點。采用文獻[15]的方法提取卷積后的特征如式(6)所示。

(6)

式中:H表示節(jié)點到邊關(guān)聯(lián)矩陣;X表示節(jié)點的特征;W、Θ代表學(xué)習(xí)參數(shù)。

(4) 構(gòu)建超圖神經(jīng)網(wǎng)絡(luò)分類器。

將多模態(tài)數(shù)據(jù)劃分訓(xùn)練集和測試集,根據(jù)多模態(tài)數(shù)據(jù)之間復(fù)雜的相關(guān)性構(gòu)建多個超邊結(jié)構(gòu)群,接著對超邊群進行相連得到超邊關(guān)聯(lián)矩陣H,將關(guān)聯(lián)矩陣H和節(jié)點特征輸入到HGNN,得到節(jié)點輸出標(biāo)簽。HGNN實現(xiàn)點到邊再到點的特征轉(zhuǎn)換,有效提取超圖上的高階相關(guān)性,通過多層不斷學(xué)習(xí)經(jīng)過Softmax得到預(yù)測標(biāo)簽。

2.5 多核超圖的構(gòu)建

在本文的基因?qū)ο蠓诸惾蝿?wù)中,N個可視對象數(shù)據(jù)的特征可以表示為X=[x1,x2,…,xn]T。本文建立超圖根據(jù)兩個特征之間歐幾里得、Pearson和Spearman距離來計算d(xi,xj)。每個頂點代表一個細胞對象,每個超邊是由一個頂點和它的K個最近鄰構(gòu)成,總共有N個超邊,每個超邊包含K+1頂點。因此,得到單個矩陣H∈RN×N,H中有N×(K+1)項等于1,其他的等于0;n個多核矩陣合并為H∈RN×nN。基因數(shù)據(jù)以圖形結(jié)構(gòu)組織,每個超邊是通過連接一個頂點和它們的鄰居節(jié)點來構(gòu)建鄰接關(guān)系。n個核得到n×N個超邊和H∈RN×nN。

3 實驗與結(jié)果分析

3.1 數(shù)據(jù)預(yù)處理

基因數(shù)據(jù)預(yù)處理,是分類準(zhǔn)確的關(guān)鍵,大部分處理scRNA-seq數(shù)據(jù)的步驟包括刪除少量異常細胞、刪除表達值為0的基因、0值占比大于一定量的基因,這些操作對一些方法有一定效果。本實驗直接對原始數(shù)據(jù)使用公認(rèn)較好的Seurat方法選高度變異的基因。

具體步驟:創(chuàng)建Seurat對象(Create SeuratObject),不刪除細胞和基因。

規(guī)范化處理如式(7)所示。

X=log(X+1)

(7)

式中:X是表示基因表達矩陣,X值可以是原始counts或FPKM、TPM值;X+1是為了防止表達值為0時取log出錯。

根據(jù)離散度值(Dispersion)選出前10%基因(HVG),如圖3所示。

圖3 高度變異基因

選用Pollen數(shù)據(jù)集,其中基因23 730個,選擇高度變異的2 373個基因,如圖3灰點所示。并標(biāo)記出排名前十的基因,如Spike1、HBG2等基因。

3.2 超圖神經(jīng)網(wǎng)絡(luò)學(xué)習(xí)

將基因數(shù)據(jù)集分兩類,訓(xùn)練集70%,測試集30%。根據(jù)細胞兩兩之間的歐氏距離、Pearson和Spearman建立超圖:其中每個頂點代表細胞,每個超邊由細胞與細胞之間的K(K=10,15,20)個最近鄰連接。將預(yù)處理數(shù)據(jù)直接輸入超圖神經(jīng)網(wǎng)絡(luò)進行學(xué)習(xí),其中HGNN為兩層,使用Softmax生成預(yù)測標(biāo)簽。

3.3 實驗結(jié)果對比

實驗使用不同測序技術(shù)平臺產(chǎn)生的8個數(shù)據(jù)集進行驗證,數(shù)據(jù)來源及描述見表1。其中scPred方法中,作者使用閾值為0.9能有較高準(zhǔn)確率,但在本文數(shù)據(jù)集中,閾值設(shè)置0.9準(zhǔn)確率較低。為提高準(zhǔn)確率,本次實驗降低閾值,設(shè)置為0.7。所有HGNN訓(xùn)練次數(shù)設(shè)置為600,學(xué)習(xí)率0.001。實驗結(jié)果如表2所示。

表2 多種算法在不同數(shù)據(jù)集上準(zhǔn)確率的比較(%)

其中HGNN表示原始數(shù)據(jù)集直接經(jīng)過HGNN處理;MHGNN為3種距離進行多核合并得出結(jié)果(其中Retina數(shù)據(jù)集較大,使用3種距離合并后進行超圖學(xué)習(xí),需占用大量內(nèi)存,導(dǎo)致內(nèi)存溢出,這里只選擇歐氏距離和pearson距離);HVG-MHGNN是本文的方法,先選出高度變異的基因,再進行多核超圖學(xué)習(xí)(其中Sala數(shù)據(jù)集準(zhǔn)確率較低,可能不適合分類研究,或者其先驗知識標(biāo)注標(biāo)簽不對)。根據(jù)四個實驗結(jié)果,HVG-MHGNN準(zhǔn)確率在8個數(shù)據(jù)集中的5個最高,平均準(zhǔn)確率最高。

實驗平臺采用WinServer 2019,處理器為E5-2620v4,2.10 GHz,內(nèi)存為96 GB。多種算法在不同數(shù)據(jù)集上運行時間性能比較如表3所示,單位為秒。其中scPred使用R語言,HGNN采用Python的PyTorch。實驗中由于使用語言及多距離合并方式不同。scPred只能與HGNN能進行實驗對比,R語言在處理大型數(shù)據(jù)集時效率較低,Python讀取數(shù)據(jù)及運行速度較快,較大數(shù)據(jù)集應(yīng)選擇Python處理更合理。MHGNN和HVG-MHGNN實驗結(jié)果對比顯示,MHGNN經(jīng)過特征選擇后(HVG-MHGNN)運行效率較高。

表3 多種算法在不同數(shù)據(jù)集上運行時間性能比較 單位:s

4 結(jié) 語

多種單細胞RNA測序技術(shù)為識別細胞類型和細胞異質(zhì)性提供便利。本文提出一種基于特征選擇(HVG)多核超圖神經(jīng)網(wǎng)絡(luò)分類方法(HVG-MHGNN),對于超高維數(shù)據(jù)集,經(jīng)過特征工程,然后根據(jù)多距離視角合并再進行HGNN特征表示學(xué)習(xí)對細胞分類準(zhǔn)確率有較大提升。方法的合理性可從以下幾點思考:(1) 主流的分類方法scPred已證明在生物學(xué)和臨床場景中有很高準(zhǔn)確性,但HVG-MHGNN與scPred對比,在相同數(shù)據(jù)集上分類準(zhǔn)確率更高。(2) 多視角數(shù)據(jù)在單細胞測序成本上和采樣技術(shù)上有一定限制。從單視角數(shù)據(jù)多距離的角度考慮單細胞測序數(shù)據(jù)聚類有一定的合理性,本實驗結(jié)果已經(jīng)證明MHGNN比HGNN準(zhǔn)確率較高。(3) MHGNN在情感分析和鏈路預(yù)測已經(jīng)有一定的應(yīng)用場景,在單細胞測序數(shù)據(jù)上應(yīng)用較少。

基于此,通過多種實驗平臺不同數(shù)據(jù)集證明了HVG-MHGNN方法的有效性。該方法雖在準(zhǔn)確率有所提升,但在處理較大數(shù)據(jù)集時,超圖學(xué)習(xí)需要大量內(nèi)存,這也是圖學(xué)習(xí)的一個缺點。從多視角收集單細胞數(shù)據(jù)是今后研究的一種主流趨勢,這更能給臨床帶來診斷的參考意義。后續(xù)將繼續(xù)研究超圖如何在單細胞轉(zhuǎn)錄組數(shù)據(jù)集進行聚類,以及如何處理大型數(shù)據(jù)集。