01

概況

單細(xì)胞轉(zhuǎn)錄組測(cè)序（scRNA-seq）技術(shù)能夠?qū)?xì)胞群中的每一個(gè)細(xì)胞進(jìn)行大規(guī)模的全轉(zhuǎn)錄組分析。它的核心分析是將單細(xì)胞聚類，以揭示細(xì)胞亞型，并根據(jù)細(xì)胞之間的關(guān)系推斷細(xì)胞譜系。本文綜述了在過(guò)去幾年間發(fā)展起來(lái)的，用于單細(xì)胞轉(zhuǎn)錄組分析中聚類的機(jī)器學(xué)習(xí)和統(tǒng)計(jì)方法，重點(diǎn)介紹了如何將一些常見(jiàn)的聚類方法，如層次聚類、基于圖的聚類、混合模型、k-means、集成學(xué)習(xí)、神經(jīng)網(wǎng)絡(luò)和基于密度的聚類等加以調(diào)整及應(yīng)用，從而解決單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)分析中的獨(dú)特挑戰(zhàn)，例如低表達(dá)基因的缺失，轉(zhuǎn)錄本的不均勻覆蓋，以及由技術(shù)偏差和不相關(guān)的混雜生物變異所帶來(lái)的細(xì)胞標(biāo)記的失真。我們?cè)u(píng)價(jià)了標(biāo)準(zhǔn)化、dropouts推測(cè)以及降維等預(yù)處理步驟如何提高聚類效果。此外，我們還將介紹一些能夠?qū)r(shí)間序列樣本和多個(gè)細(xì)胞群進(jìn)行聚類并且檢測(cè)罕見(jiàn)細(xì)胞類型的新方法。最后，本文對(duì)部分開(kāi)發(fā)用于單細(xì)胞轉(zhuǎn)錄組聚類分析的軟件進(jìn)行了實(shí)驗(yàn)和比較，以評(píng)估其性能和效率，為未來(lái)的數(shù)據(jù)分析提供一定的指導(dǎo)和方向。

02

介紹

細(xì)胞的轉(zhuǎn)錄組分析可以捕捉基因的表達(dá)活性，從而揭示細(xì)胞的身份和功能。在傳統(tǒng)的bulk-RNA測(cè)序中，轉(zhuǎn)錄組是通過(guò)從生物樣本中收集的大量細(xì)胞轉(zhuǎn)錄水平的平均值來(lái)測(cè)量的，這些平均后的表達(dá)值被用于基因共表達(dá)模塊的識(shí)別和樣本聚類。由于忽略了單個(gè)細(xì)胞的特性，這些傳統(tǒng)的方法無(wú)法在單細(xì)胞分辨率上研究重要的生物學(xué)問(wèn)題，如細(xì)胞在早期發(fā)育過(guò)程中的不同功能角色、復(fù)雜組織中的不同細(xì)胞類型和細(xì)胞譜系關(guān)系。目前，scRNA-seq技術(shù)已廣泛用于量化單個(gè)細(xì)胞中的mRNA水平。在單細(xì)胞轉(zhuǎn)錄組的實(shí)驗(yàn)操作中，使用不同的捕獲方法（如FACS，F(xiàn)luidigm C1，microdroplet microfluidics）分離單細(xì)胞，然后對(duì)RNA進(jìn)行逆轉(zhuǎn)錄并擴(kuò)增測(cè)序。單細(xì)胞轉(zhuǎn)錄組的應(yīng)用已經(jīng)帶來(lái)了重要的生物學(xué)見(jiàn)解和發(fā)現(xiàn)，例如，對(duì)癌癥中腫瘤異質(zhì)性的理解。

細(xì)胞聚類是單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)分析中識(shí)別細(xì)胞亞群結(jié)構(gòu)的必要步驟，然而目前仍然存在一些挑戰(zhàn)。首先，由細(xì)胞的自身特征（如細(xì)胞所處周期階段、細(xì)胞大?。┖图夹g(shù)（捕獲方法、捕獲效率、PCR擴(kuò)增、測(cè)序深度等）引入的技術(shù)噪音和偏差。這些噪音和偏差將導(dǎo)致轉(zhuǎn)錄組的基因覆蓋極度不均勻，從而造成零覆蓋區(qū)域和dropouts的產(chǎn)生。另外，當(dāng)一個(gè)隊(duì)列的多個(gè)樣本同時(shí)進(jìn)行分析時(shí)，樣本間的技術(shù)偏差和變異將會(huì)主導(dǎo)細(xì)胞的聚類，導(dǎo)致細(xì)胞群體的形成更偏向于不同樣本來(lái)源而非細(xì)胞類型，即批次效應(yīng)。

在本文中，我們回顧了最近發(fā)展的用于提升單細(xì)胞轉(zhuǎn)錄組聚類效果或其相關(guān)的統(tǒng)計(jì)和機(jī)器學(xué)習(xí)方法。這些方法涉及：（1）用于基因表達(dá)值的標(biāo)準(zhǔn)化、dropouts推測(cè)、數(shù)據(jù)降維以及細(xì)胞特異Marker鑒定的數(shù)據(jù)預(yù)處理方法；（2）傳統(tǒng)的聚類算法，包括基于劃分的聚類、層次聚類、混合模型、基于圖的聚類、基于密度的聚類、神經(jīng)網(wǎng)絡(luò)、集成聚類和近鄰傳播聚類等；（3）在時(shí)間序列樣本和多個(gè)批次的細(xì)胞群中進(jìn)行聚類并檢測(cè)罕見(jiàn)細(xì)胞類型的新方法。我們還討論了單細(xì)胞轉(zhuǎn)錄組聚類分析中的幾個(gè)重要方面，包括細(xì)胞間相似性度量，特征值提取和單細(xì)胞聚類結(jié)果的評(píng)估。此外，我們對(duì)十多個(gè)軟件包進(jìn)行了比較，以評(píng)估它們?cè)诖笠?guī)模單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)集上的聚類性能和效率。最后，我們對(duì)聚類分析中存在的一些挑戰(zhàn)進(jìn)行了討論。

03

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

在單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)的聚類分析中，數(shù)據(jù)預(yù)處理對(duì)于減少技術(shù)變異和噪聲（如捕獲效率低、擴(kuò)增偏差、GC含量、總RNA含量和測(cè)序深度的差異等）以及建庫(kù)和測(cè)序過(guò)程中產(chǎn)生的dropouts至關(guān)重要。高維的基因表達(dá)矩陣通常需要經(jīng)過(guò)標(biāo)準(zhǔn)化及降維映射到低維空間中，一些計(jì)算方法還利用到統(tǒng)計(jì)學(xué)和數(shù)學(xué)方法來(lái)解決dropouts事件。

標(biāo)準(zhǔn)化

原始的單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)通常從兩個(gè)層面進(jìn)行標(biāo)準(zhǔn)化：細(xì)胞的標(biāo)準(zhǔn)化和基因的標(biāo)準(zhǔn)化。細(xì)胞的標(biāo)準(zhǔn)化是為了消除擴(kuò)增偏差和其他細(xì)胞特異性的效應(yīng)，可以通過(guò)常用的reads計(jì)數(shù)標(biāo)準(zhǔn)化方法實(shí)現(xiàn)，如FPKM、RPKM、TPM等?；赨MI建庫(kù)的實(shí)驗(yàn)方案，理論上已經(jīng)避免了與擴(kuò)增或測(cè)序深度相關(guān)的誤差，因?yàn)楸幌嗤琔MI標(biāo)記的reads只會(huì)統(tǒng)計(jì)一次。然而，由于測(cè)序文庫(kù)通常是不飽和的，標(biāo)準(zhǔn)化對(duì)于該類型的數(shù)據(jù)也是有效的。細(xì)胞標(biāo)準(zhǔn)化的另一個(gè)方法是使用“spike-in”，它的基本思想是，由技術(shù)原因帶來(lái)的誤差對(duì)于內(nèi)外源基因的影響是相同的。另外，使用對(duì)數(shù)轉(zhuǎn)換進(jìn)行原始計(jì)數(shù)值的處理也非常常見(jiàn)。

基因標(biāo)準(zhǔn)化的目的是為了防止一些高表達(dá)基因主導(dǎo)了分析。常用的基因標(biāo)準(zhǔn)化方法如，在PCA中包含的z-score標(biāo)準(zhǔn)化。從過(guò)往的經(jīng)驗(yàn)中可以看到，基因的標(biāo)準(zhǔn)化可以提高算法的收斂和聚類效果。值得注意的是，數(shù)據(jù)的標(biāo)準(zhǔn)化處理將會(huì)使其失去原本基因表達(dá)的相對(duì)尺度，并且由于表達(dá)值的平移，造成表達(dá)矩陣變得不那么稀疏，這可能會(huì)影響到大規(guī)模數(shù)據(jù)集的聚類結(jié)果。

在SINCERA包中，對(duì)基因的標(biāo)準(zhǔn)化方法即是z-score，對(duì)細(xì)胞的標(biāo)準(zhǔn)化則是使用截尾均值（Trimmed mean）。一些工具會(huì)執(zhí)行更為特殊的標(biāo)準(zhǔn)化。例如，BISCUIT通過(guò)學(xué)習(xí)代表技術(shù)誤差的參數(shù)，在聚類過(guò)程中進(jìn)行迭代標(biāo)準(zhǔn)化；RaceID將每個(gè)細(xì)胞內(nèi)的總表達(dá)計(jì)數(shù)標(biāo)準(zhǔn)化到所有細(xì)胞表達(dá)計(jì)數(shù)的中位值。

此外，如果基因或者細(xì)胞顯現(xiàn)出極低的表達(dá)信號(hào)（基因表達(dá)值過(guò)低或者細(xì)胞表達(dá)基因過(guò)少），通常會(huì)將其移除，因?yàn)樗鼈兺碇摷傩盘?hào)。在不同的研究中，為去除低表達(dá)基因和細(xì)胞建立了不同的閾值，這主要根據(jù)分析中囊括的細(xì)胞和基因的數(shù)量而有所不同。例如，scVDMC對(duì)PBMC樣本的處理中，表達(dá)值低于3的基因和總表達(dá)計(jì)數(shù)值小于200的細(xì)胞都將被去除。

雖然基因和細(xì)胞的標(biāo)準(zhǔn)化在目前大多數(shù)的單細(xì)胞數(shù)據(jù)分析流程中是常見(jiàn)的，但關(guān)于其對(duì)聚類結(jié)果的影響仍存在一些爭(zhēng)論。一項(xiàng)研究的分析表明，基于bulk的標(biāo)準(zhǔn)化方法在單細(xì)胞上的應(yīng)用可能會(huì)對(duì)其分析產(chǎn)生嚴(yán)重的不良后果，例如在聚類前進(jìn)行的高變基因的檢測(cè)。相同的，也有研究表明，通過(guò)中位數(shù)或者“spike-in”進(jìn)行標(biāo)準(zhǔn)化無(wú)法解決dropouts存在的問(wèn)題，反而可能消除每種細(xì)胞類型特有的生物隨機(jī)性，這兩者都會(huì)導(dǎo)致潛在的細(xì)胞類型的不恰當(dāng)聚類或表征。

通過(guò)下面的例子，我們可以認(rèn)識(shí)到標(biāo)準(zhǔn)化的重要性。



Figure 1. 巨噬細(xì)胞群t-SNE圖 來(lái)自Zilionis等人數(shù)據(jù)集的巨噬細(xì)胞群t-SNE圖。(A)依據(jù)總計(jì)數(shù)值上色；(B)依據(jù)基因S100A9原始計(jì)數(shù)值上色；(C)依據(jù)標(biāo)準(zhǔn)化后的S100A9的表達(dá)值上色。

從圖1A,B很容易看出，S100A9的原始表達(dá)值與總計(jì)數(shù)高度相關(guān)，兩個(gè)圖的中心區(qū)域計(jì)數(shù)和表達(dá)量都較低，而外圍區(qū)域計(jì)數(shù)和表達(dá)量較高。我們能得出的唯一結(jié)論是，當(dāng)細(xì)胞中捕獲的轉(zhuǎn)錄本總量增加時(shí)，S100A9轉(zhuǎn)錄本的數(shù)量也會(huì)增加。這顯然沒(méi)什么意義。而在圖1C中，經(jīng)過(guò)標(biāo)準(zhǔn)化后的S100A9表達(dá)值與總計(jì)數(shù)之間似乎沒(méi)有相關(guān)性。我們可以說(shuō)，S100A9表達(dá)的差異不依賴于測(cè)序深度等技術(shù)噪音，而應(yīng)該來(lái)自(主要)生物因素。

Dropout

單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)中一個(gè)重要的技術(shù)誤差被稱為“dropouts”。Dropout事件是指在反轉(zhuǎn)錄過(guò)程中由于缺失或轉(zhuǎn)錄本表達(dá)過(guò)低而導(dǎo)致基因未表達(dá)的錯(cuò)誤定量。先前的研究也表明，簡(jiǎn)單的數(shù)據(jù)標(biāo)準(zhǔn)化并不能解決該問(wèn)題。因此，一些聚類算法中包含了特定的機(jī)制以矯正dropouts。例如，Seurat通過(guò)跨細(xì)胞的基因共表達(dá)模式，在聚類前進(jìn)行標(biāo)記基因的挑選。

另外，也可以通過(guò)計(jì)算配對(duì)相似性來(lái)估算dropouts。CIDR便是在聚類前進(jìn)行缺失值的填補(bǔ)。首先分析單細(xì)胞中可能出現(xiàn)的dropouts，識(shí)別每個(gè)細(xì)胞中的候選dropout基因，計(jì)算每個(gè)基因的dropout率；然后使用候選基因的dropout率來(lái)估算表達(dá)水平，即當(dāng)dropout事件以高概率被識(shí)別時(shí)，檢測(cè)算法會(huì)從其它細(xì)胞的表達(dá)譜中對(duì)該基因的表達(dá)值進(jìn)行填充；最后，利用矯正后的值計(jì)算細(xì)胞間的不相似度，進(jìn)行層次聚類。Seurat和SNN-Cliq是基于共享最近鄰SNN來(lái)度量細(xì)胞相似性。已經(jīng)證明，在稀疏的高維數(shù)據(jù)中，SNN考慮到周圍的近鄰數(shù)據(jù)點(diǎn)，更適合應(yīng)用于存在dropouts的聚類分析。

在一個(gè)更復(fù)雜的概率圖模型中，BISCUIT明確估計(jì)了每個(gè)細(xì)胞中的基因表達(dá)，以及通過(guò)數(shù)據(jù)分布和先驗(yàn)分布估算的代表技術(shù)和生物學(xué)變異的參數(shù)。其中，代表著未觀測(cè)到的基因真實(shí)表達(dá)水平的隨機(jī)變量被引入圖模型中并通過(guò)吉布斯抽樣來(lái)估算表達(dá)值。

降維

降維通常用于將高維基因表達(dá)數(shù)據(jù)投射到低維空間，使分析聚焦于低維空間中的相關(guān)信號(hào)，從而更好地實(shí)現(xiàn)數(shù)據(jù)的可視化、聚類分析等，幫助進(jìn)行生物學(xué)解釋。當(dāng)維數(shù)大于樣本數(shù)時(shí)，降維還有助于解決樣本不足的統(tǒng)計(jì)學(xué)問(wèn)題。許多降維方法已經(jīng)應(yīng)用于單細(xì)胞轉(zhuǎn)錄組聚類算法，包括PCA、多維尺度變換(MDS)、t分布、隨機(jī)近鄰嵌入(t-SNE)、典型相關(guān)分析(CCA)、潛在狄利克雷分布(LDA)以及嵌入其他模型的降維等等。

PCA：將原本數(shù)據(jù)點(diǎn)映射到與協(xié)方差矩陣的最大特征值相關(guān)聯(lián)的特征向量（即主成分），以保留原始數(shù)據(jù)中的大部分方差。例如，pcaReduce在聚類前將表達(dá)矩陣映射到一個(gè)含有K-1個(gè)主成分的空間中；SC3使用PCA和拉普拉斯變換應(yīng)用于距離矩陣以獲得一致性矩陣并進(jìn)行層次聚類。此外，在聚類之后，PCA也被廣泛應(yīng)用于二維或三維的數(shù)據(jù)可視化。PCA是一種基于假設(shè)數(shù)據(jù)為高斯分布的線性投影方法，為了捕捉數(shù)據(jù)中的非線性結(jié)構(gòu)，可以使用核主成分分析與非線性核映射相結(jié)合。

MDS：也稱為主坐標(biāo)分析(PCoA)。MDS將數(shù)據(jù)點(diǎn)映射到低維空間，通過(guò)最小化所有配對(duì)數(shù)據(jù)點(diǎn)的原始空間中的距離與投影空間中的距離之間的差值，從而在低維嵌入保持原始高維空間中的數(shù)據(jù)點(diǎn)之間的距離。CIDR便是使用MDS來(lái)計(jì)算細(xì)胞的不相似矩陣。MDS的優(yōu)點(diǎn)是在低維空間中保持原始的成對(duì)距離，易于實(shí)現(xiàn)非線性特征嵌入。然而，MDS不能擴(kuò)展到大規(guī)模數(shù)據(jù)，因?yàn)楸仨氂?jì)算成對(duì)距離來(lái)最小化目標(biāo)函數(shù)。

t-SNE：是一種將距離轉(zhuǎn)換為概率的方法。t-SNE構(gòu)造一個(gè)與原始空間及映射后的低維空間中數(shù)據(jù)點(diǎn)之間的相似性相關(guān)的概率分布，然后最小化兩個(gè)分布之間的Kullback-Leibler散度。t-SNE被廣泛應(yīng)用于單細(xì)胞數(shù)據(jù)分析中的數(shù)據(jù)可視化。

CCA：是一種基于互協(xié)方差矩陣的降維方法。給定兩個(gè)或多個(gè)數(shù)據(jù)集，該方法查找每個(gè)數(shù)據(jù)集的映射，以最大化數(shù)據(jù)集之間的相關(guān)性。在單細(xì)胞轉(zhuǎn)錄組的數(shù)據(jù)分析中，CCA通常用于不同來(lái)源樣本的整合，如Seurat（圖2）。

Figure 2. Seurat CCA數(shù)據(jù)整合示意圖

LDA：該方法最初是在自然語(yǔ)言處理中提出的。LDA假設(shè)一個(gè)文檔（document）是通過(guò)如下方法生成的：首先從具有狄利克雷先驗(yàn)的話題（topic）的多項(xiàng)分布中對(duì)話題進(jìn)行抽樣，然后對(duì)文檔中的單詞（word）進(jìn)行抽樣，這些單詞的多項(xiàng)式分布是基于每個(gè)話題的狄利克雷先驗(yàn)條件。然后，每個(gè)文檔都可以在包含k個(gè)話題的低維空間中表示。cellTree使用LDA學(xué)習(xí)“topics”作為潛在特征來(lái)表示細(xì)胞，其中“words”是受所選的潛在特征制約的基因表達(dá)水平。LDA的生成過(guò)程產(chǎn)生了一組可解釋的潛在特征。

相似度及核函數(shù)

在許多聚類方法的計(jì)算過(guò)程中，不是使用降維的方法，而是通過(guò)核函數(shù)或相似度函數(shù)來(lái)計(jì)算單個(gè)細(xì)胞之間的配對(duì)相似性進(jìn)行聚類。核函數(shù)策略將從N × M表達(dá)矩陣中計(jì)算獲得N × N相似矩陣，以期望通過(guò)核映射或相似函數(shù)在隱式特征映射空間中減少原始特征空間中的差異（如果使用有效的核函數(shù)）。SNN-cliq和Seurat使用SNN作為相似圖。cellTree在用LDA找到的話題直方圖上通過(guò)卡方找到細(xì)胞間的距離。DTWscore利用時(shí)間序列樣本為每個(gè)基因找到細(xì)胞對(duì)之間的動(dòng)態(tài)時(shí)間規(guī)整（DTW）距離，以選擇高度可變的基因，其中DTW距離是基于兩個(gè)時(shí)間序列在最佳規(guī)整路徑上的比對(duì)計(jì)算的?；赥CC的聚類使用細(xì)胞間的Jensen-Shannon距離作為譜聚類或近鄰傳播聚類的輸入。SIMLR結(jié)合多個(gè)核來(lái)學(xué)習(xí)得到細(xì)胞相似矩陣，并使用秩約束和圖擴(kuò)散來(lái)解決dropouts問(wèn)題。

大多數(shù)其他方法使用更標(biāo)準(zhǔn)的相似性函數(shù)或距離函數(shù)。BackSPIN，DendroSplit，ICGS和SINCERA在層次聚類策略中使用Pearson相關(guān)來(lái)尋找最佳分割點(diǎn)。GiniClust和RaceID也分別使用相關(guān)性矩陣進(jìn)行DBSCAN和k-means聚類。參考成分分析（RCA）計(jì)算單個(gè)細(xì)胞和參考細(xì)胞之間的表達(dá)譜之間的相關(guān)性，作為聚類的新特征，以最小化技術(shù)差異和批次效應(yīng)。SC3使用斯皮爾曼、皮爾森和歐氏距離來(lái)計(jì)算細(xì)胞間的配對(duì)相似性或距離以獲得一致性矩陣。

審核編輯：劉清