1 文章信息題目：Current best practices in single-cell RNA-seq analysis: a tutorial 發(fā)表日期：2019年6月19日 雜志：Mol Syst Biol 文章在：https://www./doi/10.15252/msb.20188746 DOI：https:///10.15252/msb.20188746  圖1 2 摘要單細(xì)胞領(lǐng)域日新月異，大量的工具被開(kāi)發(fā)出來(lái)，但很難去判斷是否好用，而且如何組建一個(gè)分析流程是一個(gè)難點(diǎn)。本文將詳細(xì)介紹單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)分析的步驟，包括預(yù)處理（質(zhì)控、歸一化標(biāo)準(zhǔn)化、數(shù)據(jù)矯正、挑選基因、降維）以及細(xì)胞和基因?qū)用娴南掠畏治觥２⑶易髡邔⒄麄€(gè)流程應(yīng)用在了一個(gè)公共數(shù)據(jù)集作為展示（詳細(xì)說(shuō)明在：https://www.github.com/theislab/single-cell-tutorial），目的是幫助新入坑用戶建立一個(gè)知識(shí)體系，已入坑用戶更新知識(shí)體系。 3 前言
現(xiàn)在已經(jīng)可以利用scRNA研究斑馬魚(yú)、青蛙、渦蟲(chóng)的細(xì)胞異質(zhì)性(Briggs et al, 2018; Plass et al, 2018; Wagner et al, 2018) ，重新理解以前的細(xì)胞群體，但這個(gè)領(lǐng)域面臨的一個(gè)問(wèn)題就是沒(méi)有成熟的標(biāo)準(zhǔn)化流程。標(biāo)準(zhǔn)化之路的困難有：大量分析方法和工具的誕生（截止2019.3.7 已經(jīng)有385種工具）、爆炸式增長(zhǎng)的數(shù)據(jù)量(Angerer et al, 2017; Zappia et al, 2018)。另外根據(jù)不同研究目的，各種分支也突顯，例如在細(xì)胞分化過(guò)程中預(yù)測(cè)細(xì)胞命運(yùn)(La Manno et al, 2018)。在我們眼界大開(kāi)的同時(shí)，分析流程標(biāo)準(zhǔn)化就變得更加困難。 在未來(lái)分析流程標(biāo)準(zhǔn)化之路上，困難還會(huì)存在于技術(shù)整合層面。比如現(xiàn)在大量的scRNA工具都是用R和Python寫(xiě)的，跨平臺(tái)分析需求在增長(zhǎng)，而對(duì)編程語(yǔ)言的喜好也決定了工具的選擇。很多好用的分析工具將自己限制在用各自的編程語(yǔ)言開(kāi)發(fā)的環(huán)境中，例如Seurat、Scater、Scanpy。 接下來(lái)，就一起看看作者列出了哪些他認(rèn)為比較好的軟件和流程吧 先上一個(gè)scRNA分析總體流程圖：  圖2 4 預(yù)處理和可視化4.1 首先看一下實(shí)驗(yàn)過(guò)程比較詳細(xì)的介紹可以看：Ziegenhain et al (2017); Macosko et al (2015); Svensson et al (2017). 原文描述的關(guān)鍵點(diǎn)是：
感覺(jué)原文描述的還沒(méi)有illumina給出的詳細(xì)，那么就看看illumina的圖文并茂版：
 圖3  圖4 原始測(cè)序數(shù)據(jù)要經(jīng)過(guò)處理得到表達(dá)矩陣，注意這里有兩種表述方式：molecular counts (count matrices) 【也即是使用UMI的】和 read counts (read matrices)，取決于是否使用UMI。而作者介紹的流程中，默認(rèn)使用 count matrices，除非readmatrices和 count matrices得到的結(jié)果存在差異，才會(huì)特別介紹read matrices
原始數(shù)據(jù)處理工具主要有：CellRanger、indrops、SEQC、zUMIs 它們主要做了這么幾件事：
得到的矩陣行是轉(zhuǎn)錄本，列是barcodes【這里用barcodes而不是直接叫細(xì)胞，是因?yàn)椴煌?xì)胞的reads也可能屬于同一個(gè)barcode =》如果出現(xiàn)一孔/液滴多細(xì)胞（doublet情況），那么barcode在多個(gè)細(xì)胞都是一樣的】當(dāng)然也會(huì)出現(xiàn)有barcode但實(shí)際沒(méi)有細(xì)胞的情況（一個(gè)孔/液滴沒(méi)有細(xì)胞即droplet，但這個(gè)孔/液滴也會(huì)賦予barcode）
關(guān)于10X實(shí)驗(yàn)環(huán)節(jié)，可以看我之前寫(xiě)的：https://mp.weixin.qq.com/s/0DEybX7GnuDFhfY1uj9t9A  圖5 4.2 質(zhì)控在正式分析之前，先要確定barcode是不是對(duì)應(yīng)真正的細(xì)胞（上面已經(jīng)了解了barcode和細(xì)胞的關(guān)系），也就是進(jìn)行Cell QC，主要考慮三個(gè)因素（這幾個(gè)因素也就是現(xiàn)在流程中常用的過(guò)濾指標(biāo)）：
先看圖A：其中這個(gè)小的直方圖就是把count depth小于4000的放大，這里設(shè)定了一個(gè)閾值1500，也就是一個(gè)barcode中至少有1500的表達(dá)量 圖B：每個(gè)細(xì)胞中包含的基因數(shù)直方圖?？梢钥吹綑M坐標(biāo)有一個(gè)小的峰在400附近，這里設(shè)定的閾值是700 圖C：依舊是看count depth。從高到低排列count depth值，可以過(guò)濾一些空的液滴（empty droplets），看到從”肘部“也就是縱坐標(biāo)1500左右開(kāi)始迅速下降 圖D：看線粒體比例。如果占比很高并且細(xì)胞類(lèi)型不是線粒體特別豐富的那種（如心肌細(xì)胞），可能說(shuō)明這個(gè)細(xì)胞本身的基因數(shù)不多并且總體表達(dá)量也不高  圖6 以上三個(gè)指標(biāo)固然重要，但如果只關(guān)注其中某一個(gè)，也會(huì)產(chǎn)生誤導(dǎo)作用，所以作者建議看問(wèn)題一定要全面，并且要把數(shù)據(jù)和生物學(xué)知識(shí)結(jié)合起來(lái)。作者舉了個(gè)例子：比如線粒體表達(dá)量相對(duì)較高的細(xì)胞也可能參與了呼吸過(guò)程。細(xì)胞總體表達(dá)量低或者基因數(shù)量少，也可能是因?yàn)楫?dāng)時(shí)取的細(xì)胞處于靜止；細(xì)胞表達(dá)量很高，也可能因?yàn)楸旧砑?xì)胞體積就比較大。的確，細(xì)胞與細(xì)胞之間的總表達(dá)量還是存在較大差異的。未來(lái)也許QC會(huì)提供更多的選擇。 除了檢查細(xì)胞完整度，QC還要進(jìn)行轉(zhuǎn)錄本層面上的檢查。原始的count矩陣一般包含超過(guò)20000個(gè)基因。這里一般要根據(jù)在細(xì)胞中有表達(dá)的數(shù)量進(jìn)行過(guò)濾，但這個(gè)閾值要根據(jù)總體細(xì)胞數(shù)和預(yù)計(jì)的分群情況來(lái)靈活調(diào)整。比如有的細(xì)胞類(lèi)型本身就數(shù)量比較少（也許就50個(gè)），那么如果我們要設(shè)定”在少于50個(gè)細(xì)胞中有表達(dá)的基因“這種條件，那么可能會(huì)丟失那些總共就50個(gè)細(xì)胞中的marker基因，最終導(dǎo)致鑒定的細(xì)胞亞群會(huì)缺失。 質(zhì)控的目的就是給下游提供更高質(zhì)量的數(shù)據(jù)，但一開(kāi)始誰(shuí)也不知道這個(gè)質(zhì)量高不高，只能先進(jìn)行下游分析，看看結(jié)果（比如細(xì)胞分群結(jié)果）再判斷。尤其是針對(duì)異質(zhì)性高的細(xì)胞群體
小結(jié)
4.3 歸一化/標(biāo)準(zhǔn)化
表達(dá)矩陣中的每個(gè)count值都表示成功的細(xì)胞捕獲、成功的反轉(zhuǎn)錄、成功的測(cè)序。但即使是相同類(lèi)型的細(xì)胞，它們的count depth（也就是每個(gè)細(xì)胞的全部表達(dá)量）也會(huì)有變化，變化的來(lái)源就在于上面說(shuō)的那三步。因此在比較兩個(gè)細(xì)胞時(shí)，任何差異都可能由于實(shí)驗(yàn)測(cè)序誤差產(chǎn)生，而不是真的生物學(xué)差異。歸一化就是解決這個(gè)問(wèn)題，它把要比較的兩個(gè)count值根據(jù)各自身處的環(huán)境求出一個(gè)相對(duì)豐度，也就是放在了一個(gè)水平上考慮，減少實(shí)驗(yàn)測(cè)序誤差，突出更多的生物學(xué)差異。 最常用的歸一化方法就是：count depth scaling，也稱(chēng)為counts per million（CPM），這個(gè)方法常用于bulk轉(zhuǎn)錄組，它會(huì)根據(jù)每個(gè)細(xì)胞的總表達(dá)量計(jì)算一個(gè) size factor ，然后對(duì)其中各個(gè)基因表達(dá)量進(jìn)行normalize。
單細(xì)胞測(cè)序中使用的歸一化方法由于細(xì)胞種類(lèi)和基因錯(cuò)綜復(fù)雜，有人就在bulk的基礎(chǔ)上進(jìn)行了改動(dòng)。例如：Weinreb et al (2018) 先排除了表達(dá)量超過(guò)總體5%的基因，然后再計(jì)算size factor，主要是預(yù)防少量極高表達(dá)量基因的存在；Scran包有個(gè)pooling‐based size factor estimation方法，允許更高的細(xì)胞異質(zhì)性存在；另外Scran包在批次矯正和差異分析環(huán)節(jié)也比其他歸一化方法表現(xiàn)更好(Buttner et al, 2019)。 在單細(xì)胞RNA測(cè)序領(lǐng)域，目前有三種常用方法：其一是以10x Genomics為代表的微滴(droplet-based)測(cè)序；其二是以Namocell為代表的PCR板(plate-based)測(cè)序；其三是以BD Rhapsody為代表的微孔(micro-well-based)測(cè)序。就測(cè)序長(zhǎng)度來(lái)說(shuō)，Smart-seq/C1和Smart-seq2基于full length的測(cè)序方案，CEL-seq2, Drop-seq, MARS-seq, SCRBseq是基于UMI的測(cè)序方案。 不能指望某一種方法適用于所有類(lèi)型的scRNA數(shù)據(jù)，(Cole et al, 2019)就發(fā)現(xiàn)不同的歸一化方法對(duì)于不同類(lèi)型數(shù)據(jù)集表現(xiàn)不同，使用scone工具可以幫助選擇合適的方法。 一般在歸一化后，數(shù)據(jù)都會(huì)變成 使用log轉(zhuǎn)換的一個(gè)好處就是：讓數(shù)據(jù)更加集中，減少數(shù)據(jù)的偏斜度，從而近似于許多下游分析工具對(duì)數(shù)據(jù)為正態(tài)分布的假設(shè)（盡管scRNA數(shù)據(jù)并不是真正的符合正態(tài)分布），比如在差異表達(dá)分析和批次矯正環(huán)節(jié) 小結(jié)
4.4 數(shù)據(jù)矯正與整合數(shù)據(jù)矯正的對(duì)象種技術(shù)和生物因素都有，例如：不同批次、捕獲失?。╠ropout）、不同細(xì)胞周期。這些在之前的歸一化中沒(méi)有被矯正，但這些差異因素都可能會(huì)后面的分析產(chǎn)生影響，它們現(xiàn)在都是導(dǎo)致差異的”嫌疑人“之一。這里要做的就是把這些差異來(lái)源去掉（Regressing out 《=》【專(zhuān)門(mén)查的詞典】 同義詞partialling out ：剔除） 4.4.1 首先是生物因素最常見(jiàn)的生物矯正因素就是：轉(zhuǎn)錄組中的細(xì)胞周期信息。簡(jiǎn)單一點(diǎn)的方式就像Scanpy和Seurat對(duì)細(xì)胞周期評(píng)分進(jìn)行簡(jiǎn)單線性回歸；復(fù)雜點(diǎn)的方式就像scLVM和f‐scLVM。用來(lái)計(jì)算細(xì)胞周期分?jǐn)?shù)的marker基因可以從文獻(xiàn)中獲得 (Macosko et al, 2015)。另外，這些方法還能用來(lái)去除其他已知的生物因素，例如線粒體基因表達(dá)量（可以作為細(xì)胞應(yīng)激的標(biāo)記）。 需要注意的是：
4.4.2 然后是技術(shù)因素最常見(jiàn)的技術(shù)矯正因素就是：樣本測(cè)序深度、批次、噪音。 去除測(cè)序深度的影響，可以促進(jìn)軌跡推斷算法的表現(xiàn)，因?yàn)樗枰诩?xì)胞之間找變化的路徑，只要放在同一水平才能看到更準(zhǔn)確的總體表達(dá)高低。 批次的來(lái)源可能是：細(xì)胞捕獲的時(shí)期不同、文庫(kù)制備使用的芯片不同、測(cè)序使用的lane不同。由此產(chǎn)生的效應(yīng)存在于多個(gè)層面:一次實(shí)驗(yàn)中各個(gè)細(xì)胞群之間、同一實(shí)驗(yàn)室中進(jìn)行的不同實(shí)驗(yàn)之間、或來(lái)自不同實(shí)驗(yàn)室的數(shù)據(jù)集之間。這里主要介紹第一種和最后一種情況：
看一下Combat矯正前后的差別：其中顏色表示不同樣本  圖7 去噪也是矯正的一種類(lèi)型。單細(xì)胞數(shù)據(jù)的一個(gè)特點(diǎn)就是含有許多噪音來(lái)源，其中一個(gè)就是dropout。一些工具就用來(lái)推斷dropout，用適當(dāng)?shù)谋磉_(dá)量來(lái)替代0，例如：MAGIC、DCA、scVI、SAVER、scImpute。去噪可以提高基因間相關(guān)性的估計(jì)。這一步可以和歸一化、批次矯正及其他下游分析整合起來(lái)，例如基于Python的scVI工具。但任何方法都可能導(dǎo)致矯正過(guò)度或不足。 4.4.3 小結(jié)
4.5 挑選基因、降維、可視化人類(lèi)的scRNA數(shù)據(jù)中可能會(huì)包含25000個(gè)基因，但其中許多基因并非能提供有用信息，還有很多基因表達(dá)量直接為0。即使在QC階段去掉這些表達(dá)量為0的基因，一個(gè)單細(xì)胞數(shù)據(jù)的基因空間依然會(huì)有超過(guò)15000個(gè)維度（一個(gè)基因表示一個(gè)維度），因此需要降低維度 4.5.1 首先挑選基因就是挑那些真正”具有情報(bào)價(jià)值“的基因，也就是會(huì)數(shù)據(jù)變化起作用的基因。因此我們這里會(huì)挑選名為HVG的基因，也就是highly variable genes。根據(jù)數(shù)據(jù)集的復(fù)雜程度不同，HVGs一般會(huì)有1000-5000個(gè)（如下圖就對(duì)不同數(shù)據(jù)集的HVGs做了個(gè)統(tǒng)計(jì)）  圖8 之前有研究表明，HVGs數(shù)量從200到2400，它們降維后的表現(xiàn)差不多（Klein et al (2015)，作者建議先盡量多選一些HVGs。 比較流行的挑選HVGs的方法有Scanpy和Seurat，而且最好是在去除技術(shù)因素后挑選，避免因?yàn)榕?、測(cè)序等因素導(dǎo)致錯(cuò)誤挑選HVG。當(dāng)然還有其他挑選的方法，看Yip et al (2018). 4.5.2 接著降維挑出來(lái)HVGs后，就是降維了，力求在最少的維度中捕捉到最多的數(shù)據(jù)特征。 常用的降維方法：A-F分別是：PCA、t-SNE、diffusion maps、UMAP、ForceAtlas2（force‐directed graph）、Variance explained by the first 31 principal components (PCs)。關(guān)于單細(xì)胞數(shù)據(jù)的降維方法，詳細(xì)可以看：Moon et al (2018)  圖9 其中兩個(gè)應(yīng)用比較廣的方法是：PCA（Pearson, 1901）和diffusion maps (Coifman et al, 2005) 【diffusion maps 于2015年在單細(xì)胞領(lǐng)域走紅 Haghverdi et al (2015) 】
4.5.3 最后可視化可視化一般使用非線性降維的方法。最常用的就是2008年提出的t-SNE（ t‐distributed stochastic neighbour embedding）。t-SNE的一個(gè)特性就是關(guān)注局部而忽視整體，因此帶來(lái)的一個(gè)影響就是：可視化結(jié)果可能夸大了細(xì)胞群之間的差異，忽略了這些細(xì)胞群之間的潛在聯(lián)系 另外，使用t-SNE的一大難點(diǎn)就是 除了t-SNE，還有2018年推出的UMAP和SPRING可以用，在缺乏明確的生物學(xué)問(wèn)題時(shí)，可以用UMAP作為不錯(cuò)的數(shù)據(jù)探索。 小結(jié)
4.6 「總結(jié)」 預(yù)處理的各個(gè)階段作者貼心將預(yù)處理比作5種類(lèi)型數(shù)據(jù)的處理： 原始數(shù)據(jù)（raw data）、歸一化數(shù)據(jù)（normalized data）、矯正后的數(shù)據(jù)（corrected data）、挑選后的數(shù)據(jù)（feature‐selected data）、降維后的數(shù)據(jù)（dimensionality‐reduced data） 這5個(gè)階段又分成3個(gè)層次：
其中每個(gè)步驟適時(shí)調(diào)整，例如單一批次的數(shù)據(jù)集，就可以跳過(guò)矯正批次這一步  圖10 5 下游分析之細(xì)胞層面下游分析的目的是解釋生物問(wèn)題，例如根據(jù)表達(dá)量將細(xì)胞劃分成不同的類(lèi)型；相似細(xì)胞間表達(dá)量的微小變化也會(huì)體現(xiàn)連續(xù)的分化路徑；基因表達(dá)量之間的相關(guān)性可能與基因共表達(dá)有關(guān)... 下游分析也是有細(xì)胞層面和基因?qū)用妫?/p>
 圖11
5.1 細(xì)胞分群5.1.1 先是：分群方法
將細(xì)胞分群基本就是任何單細(xì)胞分析的必經(jīng)之路。群的劃分就是根據(jù)細(xì)胞中基因表達(dá)譜的相似性，表達(dá)譜的相似性是由于歐幾里得距離量度決定的，而距離量度又是利用的降維的數(shù)據(jù)。一般有兩種方法計(jì)算：clustering algorithms、community detection methods
5.1.2 然后是：分群后的注釋這個(gè)過(guò)程主要是基因?qū)用娴牟僮?，為每個(gè)cluster找marker gene（也就是能代表這個(gè)cluster的基因，而這個(gè)基因又和已知的細(xì)胞類(lèi)型有關(guān)）。任何的分群算法和參數(shù)設(shè)置都會(huì)將一整團(tuán)細(xì)胞分成多個(gè)群，但這些群是否真的有意義，就要靠這一步來(lái)和生物背景結(jié)合起來(lái)。 我們希望看到的是存在很多類(lèi)型的細(xì)胞，來(lái)說(shuō)明細(xì)胞異質(zhì)性的問(wèn)題，但這里關(guān)于細(xì)胞類(lèi)型這個(gè)定義還是存在爭(zhēng)議。首先，細(xì)胞類(lèi)型的劃分怎樣算是清楚，對(duì)于一些人來(lái)說(shuō)，”T cells“這個(gè)名稱(chēng)可以叫一個(gè)細(xì)胞類(lèi)型，但還有人認(rèn)為，必須繼續(xù)深入，像”CD4+ T cells“、”CD8+ T cells“才叫細(xì)胞類(lèi)型；另外，即使是同一種細(xì)胞類(lèi)型的細(xì)胞也會(huì)有不同的發(fā)育狀態(tài)，因此它們也會(huì)顯示不同的分群結(jié)果。但不管如何，它們都是當(dāng)時(shí)細(xì)胞的一種身份（identity）
因此，我們將分群的結(jié)果稱(chēng)為不同身份的細(xì)胞（cell identities）會(huì)比不同類(lèi)型的細(xì)胞（cell types）要好一些【即每個(gè)亞群可能并不是真的不同類(lèi)型細(xì)胞，只是顯示了此時(shí)此刻的細(xì)胞身份】 對(duì)于不同細(xì)胞身份的注釋?zhuān)陙?lái)也隨之細(xì)胞圖譜的研究而加速，例如小鼠腦細(xì)胞圖譜 (Zeisel et al, 2018) 、人類(lèi)細(xì)胞圖譜 (Regev et al, 2017)的發(fā)現(xiàn)，產(chǎn)生了許多參考數(shù)據(jù)庫(kù)。在缺乏相關(guān)背景的情況下，我們可以借用數(shù)據(jù)庫(kù)中已發(fā)現(xiàn)的細(xì)胞marker 基因套入我們的細(xì)胞，幫助判斷細(xì)胞身份。需要注意：通常使用的細(xì)胞表面marker基因在細(xì)胞身份鑒定方面存在局限性(Tabula Muris Consortium et al, 2018)
圖12
利用差異分析，分成兩組：某個(gè)cluster中的細(xì)胞、數(shù)據(jù)集中其余全部的細(xì)胞。然后重點(diǎn)關(guān)注這個(gè)cluster中上調(diào)的基因，因?yàn)閙arker基因一般具有更強(qiáng)的表達(dá)作用。差異分析也會(huì)使用簡(jiǎn)單的統(tǒng)計(jì)檢驗(yàn)，例如Wilcoxon rank‐sum test、t-test，將基因的差異大小排個(gè)序，選出排名靠前的基因來(lái)作為marker基因
將數(shù)據(jù)集中選出的marker基因和參考數(shù)據(jù)集進(jìn)行比對(duì)，統(tǒng)計(jì)方法可以是：enrichment tests、the Jaccard index、other overlap statistics 參考數(shù)據(jù)集可以是網(wǎng)頁(yè)工具： www.mousebrain.org、 http:///，可以將選出的marker基因在參考數(shù)據(jù)集中進(jìn)行可視化，幫助判斷這個(gè)marker基因是什么細(xì)胞身份
細(xì)胞分群、分群注釋、重分群、重注釋...這個(gè)循環(huán)很耗費(fèi)時(shí)間。自動(dòng)化注釋方法加快了這個(gè)過(guò)程，例如scmap (Kiselev et al, 2018b) 、Garnett (preprint: Pliner et al, 2019) ，但這樣的方法有利有弊。自動(dòng)化提高了速度，但相比手動(dòng)注釋也降低了靈活性。畢竟自動(dòng)化工具使用的參考數(shù)據(jù)集中可能并不包含我們數(shù)據(jù)中的這樣細(xì)胞。因此，有自動(dòng)化工具也不能完全拋棄手動(dòng)挑選，尤其針對(duì)大型數(shù)據(jù)集中多種多樣的細(xì)胞。自動(dòng)化的過(guò)程可以先幫我們粗略地給細(xì)胞加個(gè)標(biāo)記，如果有需要，我們可以繼續(xù)手動(dòng)對(duì)這種細(xì)胞繼續(xù)劃分子細(xì)胞。對(duì)于小型數(shù)據(jù)集或者缺乏參考基因集的，手動(dòng)注釋就足夠了。 5.1.3 注意
5.1.4 細(xì)胞分群衍生——細(xì)胞組成分析（Compositional analysis）就像上面的圖12中的C圖，顯示的是近端(上圖)和遠(yuǎn)端(下圖)腸上皮區(qū)域的細(xì)胞身份組成圖（顏色越深細(xì)胞密度越大）。研究細(xì)胞組成的變化也是一個(gè)新方向，例如沙門(mén)氏菌感染已被證明會(huì)增加小鼠腸上皮細(xì)胞的比例 (Haber et al, 2017)。 這個(gè)分析既需要足夠多的細(xì)胞數(shù)量來(lái)推斷各個(gè)cluser的占比，又需要足夠的樣本數(shù)量來(lái)證明是單純一個(gè)樣本得cluster數(shù)量這樣變還是總體都會(huì)這樣變。相關(guān)的分析工具還沒(méi)有太多，未來(lái)的開(kāi)發(fā)可能會(huì)借鑒單細(xì)胞質(zhì)譜流式（mass cytometry）或者是宏基因組分析【單細(xì)胞與宏基因組的結(jié)合...】 5.2 軌跡分析5.2.1 軌跡推斷Trajectory inference軌跡推斷就是為了找到不同細(xì)胞身份、分化或者生物過(guò)程中漸進(jìn)式非同步的變化，構(gòu)建出的一個(gè)動(dòng)態(tài)模型。它認(rèn)為單細(xì)胞數(shù)據(jù)實(shí)際上就是一個(gè)連續(xù)過(guò)程中的快照（snapshot），這個(gè)過(guò)程可以通過(guò)在細(xì)胞空間中尋找最小化相鄰細(xì)胞間轉(zhuǎn)錄變化的路徑來(lái)重建
圖13 2014年Monocle和Wanderlust先推出了軌跡推斷，之后誕生的分析方法更加豐富，它們?cè)诮Ｂ窂降膹?fù)雜性上有所不同，從簡(jiǎn)單的linear or bifurcating(分叉) trajectories，到復(fù)雜的graphs, trees, or multifurcating(多叉) trajectories。Saelens et al, 2018)進(jìn)行過(guò)軌跡推斷方法的比較，結(jié)論是沒(méi)有一種方法對(duì)所有類(lèi)型的軌跡推斷有效，應(yīng)該根據(jù)預(yù)期軌跡的復(fù)雜度來(lái)選擇。不過(guò)，Slingshot在簡(jiǎn)單軌跡推斷中優(yōu)于其他方法(Street et al, 2018) 。如果期望得到更復(fù)雜的軌跡，PAGA值得推薦。軌跡推斷是一個(gè)不確定的過(guò)程，可以用多種方法來(lái)進(jìn)行佐證。
5.2.2 基因表達(dá)量的動(dòng)態(tài)變化在擬時(shí)序（pseudotime）中變化的基因描述了軌跡，這組與軌跡相關(guān)的基因有望包含調(diào)控建模過(guò)程的基因，可以用來(lái)識(shí)別潛在的生物過(guò)程。 目前很少有專(zhuān)門(mén)分析基因表達(dá)動(dòng)態(tài)變化的工具。BEAM將Monocle的軌跡推斷整合進(jìn)來(lái)，允許檢測(cè)在軌跡分支過(guò)程中相關(guān)基因的動(dòng)態(tài)變化。另外還有LineagePulse (https://github.com/YosefLab/LineagePulse)考慮了dropout技術(shù)噪音但還在開(kāi)發(fā)中。 下面這樣的圖在Slingshot的幫助文檔就有提及：https:///packages/release/bioc/vignettes/slingshot/inst/doc/vignette.html 【4.1：Identifying temporally expressed genes】 Slingshot基因表達(dá)量的動(dòng)態(tài)變化 5.2.3 細(xì)胞亞穩(wěn)態(tài)分析 Metastable states
擬時(shí)序分析會(huì)展示出不同階段細(xì)胞數(shù)量的多少。假設(shè)細(xì)胞以無(wú)偏的方式采樣，其中軌跡中的稠密區(qū)域就表示轉(zhuǎn)錄時(shí)首選的方案。當(dāng)把軌跡理解為一條時(shí)間線時(shí)（例如在發(fā)育這個(gè)時(shí)間線），這些密集的區(qū)域可能代表細(xì)胞的亞穩(wěn)態(tài)，可以結(jié)合擬時(shí)間坐標(biāo)來(lái)繪制直方圖，找到這些亞穩(wěn)態(tài)【因此看到B圖中很多種狀態(tài)，但C中直方圖認(rèn)為這幾個(gè)密集的區(qū)域才屬于亞穩(wěn)態(tài)】 Metastable states 5.2.4 整合分群與軌跡分析
將分群的結(jié)果當(dāng)成節(jié)點(diǎn)（node），將軌跡當(dāng)成節(jié)點(diǎn)之間的橋梁（edge），所以將動(dòng)靜數(shù)據(jù)結(jié)合在了一起。利用partition‐based graph abstraction（PAGA）這個(gè)工具就能得到類(lèi)似下面這個(gè)圖。
整合分群與軌跡分析 6 下游分析之基因?qū)用?/h3>之前都是對(duì)細(xì)胞進(jìn)行分析，但細(xì)胞中的基因分析會(huì)提供更多的信息。例如差異表達(dá)分析、基因集分析和基因調(diào)控網(wǎng)絡(luò)推斷，不是表面上研究細(xì)胞異質(zhì)性，而是基于異質(zhì)性探索基因表達(dá)相關(guān)的原因 6.1 差異表達(dá)分析
這個(gè)方法也是常規(guī)bulk轉(zhuǎn)錄組中經(jīng)常做的。不過(guò)單細(xì)胞相比于bulk轉(zhuǎn)錄組的一個(gè)優(yōu)勢(shì)就是：可以深入一個(gè)層次，原來(lái)bulk只是看一塊組織的平均表達(dá)量，但現(xiàn)在經(jīng)過(guò)分群后，能得到一塊組織中各種各樣的亞群，再結(jié)合差異分析，對(duì)理解異質(zhì)性問(wèn)題更有幫助。 雖然都是朝著一個(gè)方向前進(jìn)，但單細(xì)胞和bulk轉(zhuǎn)錄組的差異分析方法還是不同的。
但最近(Soneson & Robinson, 2018)研究表明，基于大批量的差異分析，bulk分析方法的性能與最好的單細(xì)胞分析方法相當(dāng)。當(dāng)bulk方法進(jìn)行改進(jìn)，加入基因權(quán)重分析后，表現(xiàn)要好于單細(xì)胞原有工具。例如：bulk差異分析工具DESeq2/EdgeR + ZINB‐wave工具估算的權(quán)重。 不過(guò)，bulk差異分析工具的性能雖然好，但是計(jì)算的效率很難提升。畢竟單細(xì)胞數(shù)據(jù)樣本數(shù)量越來(lái)越多，程序跑的時(shí)間長(zhǎng)短也成了衡量工具優(yōu)劣的重要因素。單細(xì)胞工具MAST脫穎而出。在單個(gè)數(shù)據(jù)集的小范圍比較中，完勝bulk和其他單細(xì)胞方法(Vieth et al, 2017)。而且MAST比bulk方法快了10到100倍 (Van den Berge et al, 2018) 。 小結(jié)
6.2 基因集分析
例如差異分析我們往往能得到上千基因，為了比較方便解讀，一般會(huì)把有共同特性的基因歸為一組，然后檢查我們歸類(lèi)的可靠性 【grouping the genes into sets based on shared characteristics and testing whether these characteristics are overrepresented in the candidate gene list.】 我們一般關(guān)注基因在生物過(guò)程（biological processes, BP）中的富集，可以使用MSigDB、GO、KEGG pathway、Reactome數(shù)據(jù)庫(kù) 另外，單細(xì)胞中的一個(gè)新進(jìn)展就是利用成對(duì)基因標(biāo)簽進(jìn)行配體受體分析（ ligand–receptor analysis）
配體-受體成對(duì)標(biāo)簽可以從：CellPhoneDB數(shù)據(jù)庫(kù)獲得，然后用來(lái)解釋cluster之間高表達(dá)基因的聯(lián)系 例如，利用celltalker 就可以做 Celltalker分析 6.3 基因調(diào)控網(wǎng)絡(luò) gene regulatory network (GRN)
方法例如：SCONE、PIDC、SCENIC (Single-Cell rEgulatory Network Inference and Clustering)，但發(fā)展還不是很完善，推斷的調(diào)控關(guān)系不是很穩(wěn)定【謹(jǐn)慎使用】 7 分析平臺(tái)現(xiàn)在開(kāi)發(fā)了很多平臺(tái)，整合了一套分析流程，有基于R的(McCarthy et al, 2017; Butler et al, 2018) ，python的 (Wolf et al, 2018)，本地的(Patel, 2018; preprint: Scholz et al, 2018) ，網(wǎng)頁(yè)版帶可視化的(Gardeux et al, 2017; Zhu et al, 2017) Seurat是使用最廣泛的，Scater在QC和預(yù)處理中表現(xiàn)優(yōu)異；除此以外，基于Python的scanpy也逐漸發(fā)展起來(lái)，它對(duì)于大量細(xì)胞的標(biāo)準(zhǔn)化方面表現(xiàn)不錯(cuò) 如果不使用命令行，可視化界面也有，只不過(guò)用戶只能跑人家已經(jīng)寫(xiě)好的腳本，操作靈活性不足。這樣的平臺(tái)更多的用處是在可視化探索上，例如Granatum、ASAP。未來(lái) Human Cell Atlas（HCA）會(huì)在數(shù)據(jù)可視化探索上迅速發(fā)展： https://www./data-sharing 8 結(jié)語(yǔ)8.1 作者的結(jié)語(yǔ)作者把流程測(cè)試和說(shuō)明都放在了：https://github.com/theislab/single-cell-tutorial 感興趣的可以跟著走一遍，比較一下不同的工具。作者希望這一篇能代表單細(xì)胞領(lǐng)域目前發(fā)展的一個(gè)最新動(dòng)向。他也提到，新方法層出不窮，本文介紹的大量的方法是經(jīng)過(guò)實(shí)踐比較、驗(yàn)證過(guò)的。目前可用的方法不管是運(yùn)行效率還是易用性可能都不如最新開(kāi)發(fā)的方法，但要注意：新方法在未被大量驗(yàn)證之前都需小心使用。而且新方法一般都是針對(duì)單個(gè)層面（比如降維、分群、軌跡推斷等），大體的分析流程基本固定了。 未來(lái)整合深度學(xué)習(xí)和單細(xì)胞多組學(xué)是兩個(gè)重要的發(fā)展方向，流程化運(yùn)行更是趨勢(shì)。 隨著文庫(kù)制備和測(cè)序技術(shù)的進(jìn)步，未來(lái)的單細(xì)胞平臺(tái)必將可以處理多種類(lèi)型數(shù)據(jù)：DNA甲基化、蛋白豐度等等。 8.2 劉小澤的結(jié)語(yǔ)
三天的時(shí)間，基本每天都會(huì)花半天時(shí)間在閱讀這篇綜述上。從第一眼看到它的文章邏輯，就感覺(jué)：嗯是它，沒(méi)錯(cuò)了！連午覺(jué)都不想睡了。 一開(kāi)始想強(qiáng)迫自己看下去，沒(méi)想到，越看越精彩。尤其是將整個(gè)流程和自己的知識(shí)結(jié)合起來(lái)，就看得比較順暢。為了更加易讀，我在其中加了很多注釋?zhuān)ㄖ白约簩?xiě)的一些推文和網(wǎng)上一些好的資源，可以幫助梳理知識(shí)點(diǎn)。 最后，希望看完本文對(duì)你有幫助??！ 歡迎關(guān)注我們的公眾號(hào)～_～ Welcome to our bioinfoplanet! |
|
|
