基于稀疏自編碼器的細(xì)胞類型反卷積模型SMCTD 設(shè)計與實現(xiàn)

朱賢振 李雪玲

摘要：單細(xì)胞 RNA 測序已成為研究生物學(xué)重要特征的強大高分辨率工具。然而，其測序條件苛刻，價格成本高昂。目前細(xì)胞類型反卷積能夠很好地解決這些限制問題，SMCTD（Sparse Model Cell Type Deconvolution） 使用稀疏自編碼器優(yōu)化TAPE（Tissue-AdaPtive autoEncoder） ，使其在直腸癌和PBMC模擬數(shù)據(jù)上預(yù)測細(xì)胞類型比列具有更高的靈敏度、準(zhǔn)確性和整體性能，同時在預(yù)測細(xì)胞類型特異性基因表達上表現(xiàn)更優(yōu)。

關(guān)鍵詞：單細(xì)胞測序；細(xì)胞類型反卷積；深度學(xué)習(xí)；稀疏自編碼器；一致性相關(guān)系數(shù)

中圖分類號：TP311 文獻標(biāo)識碼：A

文章編號：1009-3044（2024）11-0009-04

隨著二代深度測序（NGS） 、單細(xì)胞轉(zhuǎn)錄組測序（scRNA-seq） 、空間轉(zhuǎn)錄組（Spatial Transcriptomics） 技術(shù)、細(xì)胞類型反卷積算法的發(fā)展，為整合單細(xì)胞測序數(shù)據(jù)和大規(guī)模 bulk 基因表達譜，研究癌癥微環(huán)境中的細(xì)胞組成和基因表達，提供了重要手段。然而，現(xiàn)有的細(xì)胞類型反卷積方法[1]的準(zhǔn)確度和解析顆粒度有很大提升空間，發(fā)展基于深度機器學(xué)習(xí)的細(xì)胞類型反卷積算法，能夠?qū)崿F(xiàn)并加速大范圍內(nèi)高通量臨床數(shù)據(jù)的精確分析[2]。

TAPE[3]是一種連接批量RNA-seq和單細(xì)胞RNAseq[4]的深度學(xué)習(xí)方法，可在短時間內(nèi)實現(xiàn)精確的反卷積。通過構(gòu)建可解釋的解碼器并在獨特的方案下進行訓(xùn)練，TAPE 可以自適應(yīng)地預(yù)測細(xì)胞類型分?jǐn)?shù)和細(xì)胞類型特異性基因表達。與多個數(shù)據(jù)集上的流行方法相比，TAPE 在細(xì)胞類型水平上具有更好的整體性能和相當(dāng)?shù)臏?zhǔn)確性。此外，它在不同細(xì)胞類型中更穩(wěn)健、更快、更靈敏，可以提供具有生物學(xué)意義的預(yù)測。然而TAPE是基于傳統(tǒng)自編碼器的深度學(xué)習(xí)模型，需要對所有數(shù)據(jù)進行處理降維，容易過擬合，盡管傳統(tǒng)自編碼器可以學(xué)習(xí)數(shù)據(jù)的特征表示，但并不保證這些特征是有意義、可解釋的。模型可能會學(xué)習(xí)到捕捉數(shù)據(jù)中的一些噪聲或冗余信息，而不是真正有用的特征。

稀疏自編碼器相比自編碼器特征的學(xué)習(xí)與選擇能力更好，稀疏自編碼器通過強制編碼層的神經(jīng)元保持較低的激活率，鼓勵網(wǎng)絡(luò)只激活最重要的特征神經(jīng)元[5]。這可以使更有意義、更具區(qū)分性的特征被學(xué)習(xí)和保留，從而提高模型對數(shù)據(jù)的表征能力。也能降低過擬合的風(fēng)險： 通過在編碼層施加稀疏性約束，稀疏自編碼器可以降低過擬合的風(fēng)險。這是因為模型被迫僅僅選擇最重要的特征，而不會過度適應(yīng)訓(xùn)練數(shù)據(jù)中的噪聲或不重要的變化。更好的泛化能力： 由于稀疏自編碼器傾向于學(xué)習(xí)更有意義的特征，它們通常能夠更好地泛化到新的、未見過的數(shù)據(jù)，從而提高模型的泛化能力。還擁有更高效的特征表示：稀疏自編碼器可以學(xué)習(xí)數(shù)據(jù)的更緊湊、更高級的表示。這種表示可以更好地捕捉數(shù)據(jù)中的關(guān)鍵模式和結(jié)構(gòu)，從而在后續(xù)的任務(wù)（如分類、聚類等）中表現(xiàn)更好。

1 細(xì)胞類型反卷積模型構(gòu)建

細(xì)胞類型反卷積是能夠?qū)⒋笠?guī)模bulk基因表達譜反卷積得到bulk中細(xì)胞類型比例和細(xì)胞類型特異性基因的一種方法，傳統(tǒng)細(xì)胞類型反卷積需要制作簽名矩陣，簽名矩陣是基于一組特定基因的表達模式構(gòu)建的。這些特定基因通常與某個生物學(xué)特性、狀態(tài)或功能相關(guān)聯(lián)。通過分析大量單細(xì)胞數(shù)據(jù)，可以從中提取出這些特定基因的表達模式，并將這些模式整合到一個矩陣中，即簽名矩陣[1]。非常煩瑣，最近幾年興起的深度學(xué)習(xí)細(xì)胞類型反卷積能夠自主學(xué)習(xí)bulk基因表達譜中特征，無須構(gòu)建簽名矩陣即可進行反卷積。這些方法中，TAPE是反卷積性能最出色的之一，比同樣使用深度學(xué)習(xí)的Scaden [2]更穩(wěn)健，與傳統(tǒng)細(xì)胞類型反卷積方法CIBERSORTx [1]性能旗鼓相當(dāng)。而TAPE 使用自編碼器的缺點也讓我們思考用更優(yōu)秀的網(wǎng)絡(luò)構(gòu)建模型來提升性能。

1.1 模型設(shè)置

如圖1所示，我們引入稀疏自編碼器，與傳統(tǒng)自編碼器相比，引入稀疏性目標(biāo)值和稀疏性權(quán)重，稀疏性目標(biāo)值是一個預(yù)設(shè)的小數(shù)，表示隱藏層神經(jīng)元的期望平均激活度，通常是0.05或者更小。稀疏性權(quán)重是一個常數(shù)，表示稀疏性懲罰項在總代價函數(shù)中的權(quán)重，用來控制稀疏性的強度。圖1中B 表示輸入的bulk基因表達數(shù)據(jù)，C 表示經(jīng)過稀疏自編碼器編碼再解碼重現(xiàn)的bulk基因表達數(shù)據(jù)，X 表示通過編碼器得到的細(xì)胞類型比例。圖1左半部分表示為編碼器，是一個回歸模型，負(fù)責(zé)將高維bulk基因表達數(shù)據(jù)映射到低維的細(xì)胞類型比例數(shù)據(jù)。相反，圖1右半部分可以根據(jù)細(xì)胞類型比例數(shù)據(jù)重建bulk基因表達數(shù)據(jù)。

1.2 模型訓(xùn)練

我們先預(yù)設(shè)稀疏性目標(biāo)值和稀疏性權(quán)重，稀疏性目標(biāo)值代表期望的神經(jīng)元平均激活度，稀疏性權(quán)重用來控制稀疏性強度。然后在代碼中定義了KL散度函數(shù)（Kullback-Leibler Divergence） ，使用KL散度函數(shù)和稀疏性目標(biāo)值和權(quán)重計算稀疏性懲罰損失[5]。

我們使用大約5 000個bulk樣本進行訓(xùn)練。使用預(yù)測細(xì)胞類型比例和真實細(xì)胞類型比例之間的MAE （平均絕對誤差）與稀疏性懲罰損失的和來優(yōu)化編碼器的參數(shù)，并使用重構(gòu)bulk數(shù)據(jù)和原始bulk數(shù)據(jù)之間的MAE 與稀疏性懲罰損失的和來優(yōu)化解碼器和編碼器。

1.3 模型預(yù)測

使用反卷積模型進行預(yù)測需要預(yù)先準(zhǔn)備單細(xì)胞參考數(shù)據(jù)，行為細(xì)胞類型，列為基因名稱，文件為TXT 格式。為了使反卷積結(jié)果更精準(zhǔn)更具有生物學(xué)意義，單細(xì)胞參考數(shù)據(jù)要與需要預(yù)測的bulk 數(shù)據(jù)為同一組織的，并且擁有相同的細(xì)胞類型。

需要預(yù)測的bulk 數(shù)據(jù)需要指定分隔符，行為樣本名稱，列為基因名稱，數(shù)據(jù)類型最好為“counts”，若使用“TPM”或者“FPKM”格式須自備基因長度文件使數(shù)據(jù)最終以counts格式運行在程序上。

我們把準(zhǔn)備好的單細(xì)胞參考數(shù)據(jù)和bulk數(shù)據(jù)輸入程序，將模式選擇為“overall”，然后選擇合適的數(shù)據(jù)類型及基因長度文件，自適應(yīng)參數(shù)選擇為“True”或者“False”。如果是“ True”，那么它將會預(yù)測輸出簽名矩陣，反之，則返回空值。等待程序運行完成，會得到預(yù)測的細(xì)胞類型比例數(shù)據(jù)（行為樣本，列為細(xì)胞類型）及選擇可得到的簽名矩陣。

2 細(xì)胞類型反卷積模型性能比較

由于公共數(shù)據(jù)庫中同一樣本中既測bulk數(shù)據(jù)又測單細(xì)胞數(shù)據(jù)的少之又少，所以為了精準(zhǔn)測出反卷積模型的性能，因此有必要進行偽bulk數(shù)據(jù)測試進行估計。偽bulk數(shù)據(jù)是通過具有基本事實（預(yù)定義的細(xì)胞類型比例）的單細(xì) 胞基因表達數(shù)據(jù)在計算機中生成的。也就是說，偽bulk數(shù)據(jù)是許多單細(xì)胞基因表達數(shù)據(jù)的總和。我們將使用TAPE中的偽bulk模擬程序模擬bulk。

2.1 GSE176078單細(xì)胞數(shù)據(jù)模擬預(yù)測比較

首先，我們預(yù)設(shè)真實細(xì)胞類型比例，再將從GEO 中下載的乳腺癌數(shù)據(jù)集GSE176078 [6]的樣本作為參考的單細(xì)胞數(shù)據(jù)生成模擬bulk對五種細(xì)胞類型進行反卷積性能預(yù)測。設(shè)定兩個參數(shù)指標(biāo)：MAE（平均絕對誤差）和CCC（一致性相關(guān)系數(shù)）[7]，MAE是對每個細(xì)胞類型的預(yù)測值與其對應(yīng)的實際值之間的絕對差值進行求和，然后對所有數(shù)據(jù)點的絕對差值求平均值，數(shù)值越小性能越好。CCC是評價細(xì)胞類型比例預(yù)測值與真實值之間的一致性的指標(biāo)，CCC值越接近1代表性能越好。

最終結(jié)果如圖2所示，我們的模型SMCTD（稀疏自編碼器）在誤差方面要比TAPE和Scaden都低，在CCC方面也要比TAPE更出色，略遜于Scaden，綜合兩方面來看，SMCTD是三者中反卷積性能最出色的。

接下來測試模型在細(xì)胞類型增加的情況下的性能，我們繼續(xù)用GSE176078模擬bulk，這次bulk包含14 種細(xì)胞亞型，例如：Monocyte、Fibroblasts、NK cells 等。評價指標(biāo)同上。

如圖3所示，結(jié)果表明，在這種情況下，所有方法都出現(xiàn)性能下降的情況，但這些方法的MAE與預(yù)測五種細(xì)胞類型情形中的 MAE相當(dāng)，這表明這些方法可以預(yù)測接近真實值的值。同時SMCTD也是這種情況下誤差最低的算法，其CCC值也為三種方法第二高，說明其性能是三者間最出色的。

2.2 PBMC 單細(xì)胞數(shù)據(jù)模擬預(yù)測比較

接下來我們使用10X Genomics官網(wǎng)的PBMC（外周血單個核細(xì)胞）單細(xì)胞數(shù)據(jù)[8]模擬bulk（其中包含七種細(xì)胞類型）進行反卷積，評價指標(biāo)同上。

最終結(jié)果如圖4 所示，可以發(fā)現(xiàn)在PBMC 數(shù)據(jù)上Scaden 表現(xiàn)是最出色的，我們的模型SMCTD 雖然比Scaden 略差，但要比TAPE在MAE、CCC值兩個方面都要更優(yōu)秀。

綜上所述，我們發(fā)現(xiàn)模型SMCTD在預(yù)測細(xì)胞類型方面有著不遜于TAPE、Scaden的性能，雖然這三種方法在細(xì)胞類型增多的情況下都會出現(xiàn)性能下降，這種情況是可以預(yù)見的，所以解決這一問題也是未來研究的方向之一。

2.3 組織適應(yīng)性細(xì)胞類型特異性基因表達預(yù)測

SMCTD不僅可以預(yù)測細(xì)胞類型分?jǐn)?shù)，同樣可以自適應(yīng)地預(yù)測細(xì)胞類型特異性基因表達。也就是說，SMCTD 只需要模擬數(shù)據(jù)來訓(xùn)練，如果給出相應(yīng)的bulk RNAseq數(shù)據(jù)，它可以預(yù)測細(xì)胞類型特異性的基因表達。此功能使SMCTD 能夠剖析不同細(xì)胞類型中的bulk基因表達，并發(fā)現(xiàn)不同細(xì)胞類型中一些潛在差異表達的基因。

我們測試了預(yù)測的細(xì)胞類型特異性 GEP 的正確性。為了測試這一點，我們測量了每種細(xì)胞類型的預(yù)測基因表達值與從單細(xì)胞RNA-seq獲得的原始基因表達值之間的一致性（圖5） 。這里，bulk 數(shù)據(jù)用GSE176078單細(xì)胞數(shù)據(jù)模擬生成，而單細(xì)胞數(shù)據(jù)是乳腺癌癌單數(shù)據(jù)集GSE176078。由于在訓(xùn)練階段使用Log2 和 MinMaxScaler（） 將輸入的 RNA-seq 數(shù)據(jù)轉(zhuǎn)換為 0-1 值，因此按細(xì)胞類型分組的基因表達值的總和也以這種方式轉(zhuǎn)換以進行比較與預(yù)測的相對基因表達值。

由圖5、圖6和表1可得，SMCTD在四種細(xì)胞類型上的基因表達預(yù)測一致性相關(guān)系數(shù)都要高于TAPE，僅在免疫細(xì)胞中表現(xiàn)不佳，考慮到其他四種細(xì)胞類型的良好一致性，這種失真可能是由個體差異引起的。圖表中顯示的一致性證明SMCTD正確預(yù)測了細(xì)胞類型特異性基因表達，為進一步的基因表達分析奠定了基礎(chǔ)。

3 結(jié)束語

本文優(yōu)化了一個細(xì)胞類型反卷積模型TAPE，使用了稀疏自編碼器作為模型基礎(chǔ)，增強了模型的稀疏性，提高了模型性能，使細(xì)胞類型反卷積模型的誤差降低，且提高了相關(guān)性，在我們測試的兩類數(shù)據(jù)和多種細(xì)胞類型上均能體現(xiàn)出。此外，在預(yù)測細(xì)胞類型特異性基因表達上SMCTD也比TAPE在大多數(shù)細(xì)胞類型上的結(jié)果更準(zhǔn)確。

但模型在細(xì)胞類型過多的情況下表現(xiàn)下降，這是未來要攻克的方向之一。在預(yù)測細(xì)胞類型特異性基因表達方面，SMCTD同樣有一定的優(yōu)化上升空間。

參考文獻：

[1] NEWMAN A M， STEEN C B， LIU C L， et al. Determining cell type abundance and expression from bulk tissues with digital cytometry [J]. Nat Biotechnol， 2019， 37（7）： 773.

[2] MENDEN K， MAROUF M， OLLER S， et al. Deep learningbased cell composition analysis from tissue expression profiles[J]. Science Advances，2020，6（30）：eaba2619.

[3] CHEN Y S，WANG Y X，CHEN Y L，et al. Deep autoencoder for interpretable tissue-adaptive deconvolution and cell-typespecific gene analysis[J]. Nature Communications，2022，13（1）：6735.

[4] WANG Z，GERSTEIN M，SNYDER M. RNA-Seq：a revolutionary tool for transcriptomics[J]. Nature Reviews Genetics，2009，10（1）：57-63.

[6] WU S Z，AL-ERYANI G，RODEN D L，et al. A single-cell and spatially resolved atlas of human breast cancers[J]. Nature Genetics，2021，53（9）：1334-1347.

[5] NG A. Sparse autoencoder[R]. CS294A Lecture Notes， 2011：72.

[6] WU S Z，AL-ERYANI G，RODEN D L，et al. A single-cell and spatially resolved atlas of human breast cancers[J]. Nature Genetics，2021，53（9）：1334-1347.

[7] LIN L I. A concordance correlation coefficient to evaluate reproducibility[J]. Biometrics，1989，45（1）：255-268.

[8] CHEMISTRY） K P F A H D V. single cell gene expression dataset by cell ranger 2. 1. 0 [DS]. 10X Genomics， 2017，

【通聯(lián)編輯：李雅琪】