原文鏈接：
https://www.embopress.org/doi/10.15252/msb.20188746
主編評語
這篇文章最好的地方不只在于推薦了工具雁乡，提供了一套分析流程，更在于詳細(xì)介紹了有哪些方法可用嘴拢、什么時候用哪些方法和有什么注意事項胜榔。這也是我們易生信培訓(xùn)過程中廣泛討論的問題胳喷。
Abstract
單細(xì)胞RNA-seq使研究者能夠以前所未有的分辨率研究基因表達(dá)圖譜。這一潛力吸引著更多科研工作者應(yīng)用單細(xì)胞分析技術(shù)解決研究問題苗分。隨著可用的分析工具越來越多，如何組合成一個最新最好的數(shù)據(jù)分析流程也越來越難牵辣。我們詳細(xì)闡述了一個典型的單細(xì)胞轉(zhuǎn)錄組分析各個步驟的細(xì)節(jié)和注意事項摔癣，包括預(yù)處理（質(zhì)控、標(biāo)準(zhǔn)化、數(shù)據(jù)校正择浊、特征選擇戴卜、降維）和細(xì)胞/基因水平的下游分析等∽裂遥基于獨立的比較研究投剥，我們?yōu)槊恳徊蕉纪扑]了當(dāng)前最好的方法和操作建議。隨后把這些最好的工具整合成一個分析流程并應(yīng)用于一套公共數(shù)據(jù)集的分析以演示其效果担孔。案例研究具體可見https://www.github.com/theislab/single-cell-tutorial江锨。這篇綜述為這個領(lǐng)域的新人提供了一份學(xué)習(xí)單細(xì)胞分析的指南，并且也能幫助老用戶更新他們的分析流程糕篇。
背景
近年來啄育，單細(xì)胞RNA測序（scRNA-seq）大大提高了我們對生物系統(tǒng)的了解。我們已經(jīng)能夠在研究斑馬魚拌消、青蛙和渦蟲（zebrafish, frogs and planaria）等生物細(xì)胞異質(zhì)性的同時發(fā)現(xiàn)先前未知的細(xì)胞群體挑豌。這項技術(shù)的巨大潛力激勵了計算生物學(xué)家開發(fā)了一系列分析工具。盡管開發(fā)者為了確保單個工具的可用性付出了巨大的努力墩崩，但是由于該領(lǐng)域的相對不成熟氓英，對于單細(xì)胞數(shù)據(jù)分析的新手來說，入門的障礙是缺少一份標(biāo)準(zhǔn)指南鹦筹。在本文中铝阐，我們提供了一份scRNA-seq分析的參考教程，并概述了當(dāng)前的最佳實踐方案盛龄，為將來的分析標(biāo)準(zhǔn)化奠定了基礎(chǔ)饰迹。
分析標(biāo)準(zhǔn)化面臨的挑戰(zhàn)來源于越來越多的可用分析方法（截至2019年3月7日有385種工具）和數(shù)據(jù)集規(guī)模爆炸性的提高。使得我們一直在尋找新的方法來分析處理我們的數(shù)據(jù)余舶。例如啊鸭，最近已經(jīng)有方法可以預(yù)測細(xì)胞分化過程中的命運選擇。盡管分析工具的不斷改進(jìn)有助于產(chǎn)生新的科學(xué)推論匿值，但它也使分析流程的標(biāo)準(zhǔn)化變得更為復(fù)雜赠制。
標(biāo)準(zhǔn)化的另一個挑戰(zhàn)在于軟件技術(shù)方面。用于scRNA-seq數(shù)據(jù)的分析工具是用不同的編程語言編寫的 - 最主要的是R和Python挟憔。盡管跨編程語言的支持越來越多钟些，但使用的編程語言確實影響了對分析工具的選擇。諸如Seurat绊谭，Scater或Scanpy等常用工具提供了集成環(huán)境來開發(fā)流程并包含大量分析工具政恍。然而，出于維護(hù)的需要达传， 這些平臺限制了它們只能使用各自的編程語言開發(fā)的工具篙耗。引申開來迫筑，當(dāng)前可用的scRNA-seq分析教程也存在語言限制，而且許多是圍繞這3大平臺進(jìn)行講述的 （易生信2020年最新的單細(xì)胞課程會同時提供R和Python版本的最新分析流程）：
-
R and bioconductor tools:
https://github.com/drisso/bioc2016singlecell宗弯；
https://hemberg-lab.github.io/scRNA.seq.course/; https://www.ncbi.nlm.nih.gov/pubmed/27909575?dopt=Abstract脯燃； Scanpy: https://scanpy.readthedocs.io/en/stable/tutorials.html蒙保；
考慮到上述挑戰(zhàn)辕棚，我們選擇概述當(dāng)前分析過程中每一步的最佳實踐和獨立于編程語言的常用工具，而不是評估一個完整的分析流程邓厕。我們將帶領(lǐng)讀者完成scRNA-seq分析流程的各個步驟逝嚎，介紹每個步驟的最優(yōu)方法，并討論分析過程中會遇到的陷阱和當(dāng)前未解決的問題邑狸。當(dāng)然由于一些工具比較新懈糯，并且缺乏系統(tǒng)比較，最好的工具很難決定单雾，我們列出了流行的可用工具赚哗。列出的工具從基因count矩陣開始到潛在的分析終點。在我們的Github上有結(jié)合了已建立的當(dāng)前最佳實踐流程的案例研究：https：//github.com/theislab/single-cell-tutorial/硅堆。我們在實際的示例工作流程中應(yīng)用了當(dāng)前的最佳實踐來分析公共數(shù)據(jù)集屿储。該分析流程使用Jupyter notebook和rpy2整合了R和Python的工具。結(jié)合現(xiàn)有的文檔渐逃，這已經(jīng)可以被視作一個可接受的分析流程模板了够掠。
圖1. 典型單細(xì)胞轉(zhuǎn)錄組分析流程。
原始測序數(shù)據(jù)預(yù)處理和比對后獲得Gene count矩陣茄菊，作為本分析流程的開始疯潭。這些結(jié)果圖是Count矩陣采用最佳實戰(zhàn)分析流程進(jìn)行預(yù)處理和下游分析獲得的。數(shù)據(jù)來源于Haber et al 2017腸上皮細(xì)胞文章面殖。
預(yù)處理和可視化
原始測序數(shù)據(jù)經(jīng)過處理得到分子計數(shù)矩陣（count matrix）竖哩，或者reads count （讀數(shù)矩陣）。這取決于單細(xì)胞文庫構(gòu)建方案中是否包含唯一分子標(biāo)識符（UMI脊僚， unique molecular identifiers）（參考下面的Box1. Hemberg-lab單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)分析（六）- 構(gòu)建表達(dá)矩陣相叁，UMI介紹）。原始數(shù)據(jù)處理流程如Cell Ranger (10X單細(xì)胞測序分析軟件:Cell ranger辽幌、從拆庫到定量), indrops, SEQC, zUMIs負(fù)責(zé)測序數(shù)據(jù)質(zhì)控增淹，確定reads來源的細(xì)胞 (barcodes) （這一步也稱為demultiplexing）和mRNA分子 （生信寶典注：單細(xì)胞測序不只獲得mRNA，更準(zhǔn)確說是帶poly-A尾巴的RNA乌企，包括mRNA和lncRNA等）虑润、基因組比對和定量。獲得的reads或count矩陣的行數(shù)等于barcodes的數(shù)目加酵，列數(shù)等于基因數(shù)目（生信寶典注：也可能反過來拳喻，行為基因梁剔，列為barcodes）。這里使用術(shù)語barcodes而不是cell舞蔽，因為分配給相同barcode的所有reads可能并不只是來源于同一細(xì)胞。barcode可能會錯誤地標(biāo)記多個細(xì)胞（doublet）码撰，也可能不會標(biāo)記任何細(xì)胞（空液滴/孔） (Hemberg-lab單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)分析（四）- 文庫拆分和細(xì)胞鑒定)渗柿。
盡管reads data和 count data 的測量噪聲級別不同 （生信寶典注：基于UMI的數(shù)據(jù)，獲得的是分子計數(shù)脖岛，count data朵栖，噪聲更低），但在分析流程中的處理步驟是相同的柴梆。為簡單起見陨溅，在本教程中，我們將數(shù)據(jù)稱為計數(shù)矩陣（count data）绍在。在reads和計數(shù)矩陣的結(jié)果不同的地方门扇，將會特別提到reads矩陣。
Box1：scRNA-seq實驗流程的關(guān)鍵要素
從生物樣品到單細(xì)胞數(shù)據(jù)需要多步實驗操作偿渡。典型的工作流程包含單細(xì)胞解離臼寄、單細(xì)胞分離、文庫構(gòu)建和測序溜宽。我們在這里簡要介紹這些步驟吉拳。不同protocols的更詳細(xì)的解釋和比較可參考下面三篇文章：
Ziegenhain et al (2017)：https://www.embopress.org/doi/10.15252/msb.20188746#msb188746-bib-0159
Macosko et al (2015)：https://www.embopress.org/doi/10.15252/msb.20188746#msb188746-bib-0083
Svensson et al (2017)：https://www.embopress.org/doi/10.15252/msb.20188746#msb188746-bib-0126
- 單細(xì)胞實驗的輸入材料通常是生物組織樣品。第一步适揉，單細(xì)胞解離：消化組織產(chǎn)生單細(xì)胞懸液留攒。為了分別分析每個細(xì)胞中的
mRNA，必須分離單細(xì)胞嫉嘀。根據(jù)實驗方案不同炼邀，單細(xì)胞分離的方式也有所不同。
- 基于平板的技術(shù)將細(xì)胞分到到板上的孔中吃沪。
- 基于液滴的方法則依賴于微流體液滴捕獲單個細(xì)胞汤善。
- 在這兩種情況下，都可能出現(xiàn)一些問題票彪，如多個細(xì)胞一起被捕獲（
doublets or multiplets）红淡、非活細(xì)胞被捕獲或根本沒有細(xì)胞被捕獲（空液滴/孔）〗抵基于液滴的方法需要通過低的輸入細(xì)胞濃度來保持低的doublets率在旱，因此空液滴是特別常見的 （生信寶典注：一般beads和細(xì)胞的輸入比例是20:1）。
每個孔或液滴均包含必要的試劑以裂解細(xì)胞膜并進(jìn)行文庫構(gòu)建 (生信寶典注：植物單細(xì)胞就要注意了推掸，需要提前去除細(xì)胞壁)桶蝎。文庫構(gòu)建包括捕獲細(xì)胞內(nèi)mRNA驻仅、反轉(zhuǎn)錄為cDNA分子并進(jìn)行擴(kuò)增等過程栏尚。因為文庫構(gòu)建時每個細(xì)胞是獨立的丰泊，所以每個細(xì)胞的mRNA也就特異的標(biāo)記了孔特異性或液滴特異性細(xì)胞barcode。此外恨诱，許多實驗方案還使用唯一分子標(biāo)識符（UMI）標(biāo)記捕獲的RNA分子胜茧。一般在測序之前需要先擴(kuò)增細(xì)胞cDNA以增加其被檢測的可能性粘优。但微量擴(kuò)增更容易引入PCR偏好性。UMI使我們能夠區(qū)分測到的reads是來源于mRNA分子的不同擴(kuò)增拷貝還是來源于獨立的mRNA分子呻顽，從而可以進(jìn)行更準(zhǔn)確的定量雹顺。
每個細(xì)胞單獨構(gòu)建的cDNA文庫都帶有cell barcode和/或UMI（取決于protocol），后續(xù)將這些文庫混合在一起測序廊遍。測序產(chǎn)生的reads數(shù)據(jù)進(jìn)行質(zhì)量控制 (Hemberg-lab單細(xì)胞轉(zhuǎn)錄組數(shù)據(jù)分析（三）- 原始數(shù)據(jù)質(zhì)控)嬉愧，根據(jù)其barcodes序列分組（demultiplexing），并且進(jìn)行后續(xù)比對定量喉前。對基于UMI的protocols没酣，reads的數(shù)據(jù)可以進(jìn)一步demultiplexed以得到捕獲的mRNA分子的計數(shù)（count data）。也就是本套流程的起始輸入數(shù)據(jù)卵迂。
質(zhì)控
在分析單細(xì)胞基因表達(dá)數(shù)據(jù)之前四康，我們必須確保所有barcode都對應(yīng)于有效細(xì)胞 (viable cells，有活力的細(xì)胞)狭握。質(zhì)控有3個指標(biāo)：
測到的轉(zhuǎn)錄本分子總數(shù)
測到的基因總數(shù)
來源于線粒體基因的轉(zhuǎn)錄本所占比例
質(zhì)控就是檢查這3個指標(biāo)的分布中是否存在異常峰并設(shè)置閾值去除闪金。這些異常的barcodes可能對應(yīng)于死細(xì)胞、細(xì)胞膜破損的細(xì)胞或doublets论颅。比如哎垦，如果某個barcode對應(yīng)的樣品測到的分子總數(shù)低、檢測到的基因數(shù)少恃疯、線粒體基因所占比例高漏设，則表明該樣品可能存在細(xì)胞膜破損導(dǎo)致細(xì)胞質(zhì)RNA漏出，只有線粒體中的RNA保留了下來今妄。相反郑口，如果某個barcode對應(yīng)的樣品有異常高的總分子數(shù)和檢測到的基因數(shù)，則有可能這個樣品包含2個或以上細(xì)胞（doublets）盾鳞。因此犬性，高的總分子數(shù)通常用來過濾潛在的doublets。3個最近發(fā)表的doublets檢測工具提供了更好的解決方案：DoubletDecon, Scrublet, Doublet Finder腾仅。
分開考慮這三個QC變量中的任何一個都可能導(dǎo)致對細(xì)胞狀態(tài)的錯誤解讀乒裆。例如，線粒體計數(shù)相對較高的細(xì)胞可能是呼吸過程比較活躍 （生信寶典注：如心臟細(xì)胞總mRNA的30%是線粒體基因推励，具體見對一篇單細(xì)胞RNA綜述的評述：細(xì)胞和基因質(zhì)控參數(shù)的選擇）鹤耍。同樣肉迫，其他QC變量也具有相應(yīng)的生物學(xué)意義「寤疲總分子數(shù)和/或基因數(shù)低的細(xì)胞可能是處于靜息狀態(tài)的細(xì)胞群體喊衫；總分子數(shù)和/或基因數(shù)高的細(xì)胞也可能是細(xì)胞自身體積較大。實際上杆怕，細(xì)胞之間的總分子數(shù)自身也可能有很大差異（具體見Github上的案例研究）格侯。因此，在做出是否過濾的單閾值決策時财著，應(yīng)聯(lián)合考慮3個QC變量，并且應(yīng)將這些閾值設(shè)置的盡可能寬松撑碴，以避免無意間濾除有效的細(xì)胞群撑教。將來，考慮多元QC依賴的過濾模型可能會提供更敏感的QC選擇方式醉拓。
包含不同類型細(xì)胞混合體的數(shù)據(jù)集的每個QC指標(biāo)可能會呈現(xiàn)多個分布聚集峰伟姐。例如，圖2D顯示了具有不同QC分布的兩個細(xì)胞群亿卤。如果未執(zhí)行先前的過濾步驟（請注意愤兵，Cell Ranger也會執(zhí)行細(xì)胞QC），則應(yīng)僅將最低總分子數(shù)和基因數(shù)的barcode視為無效細(xì)胞排吴。另一個閾值設(shè)定準(zhǔn)則是考慮所選閾值過濾掉的細(xì)胞比例秆乳。比如過濾高分子總數(shù)細(xì)胞時，過濾掉的細(xì)胞比例不應(yīng)超過預(yù)期的doublet率钻哩。
除了檢查細(xì)胞的完整性外屹堰，還必須在轉(zhuǎn)錄本水平上執(zhí)行QC步驟。原始計數(shù)矩陣通常包含超過20,000個基因街氢〕都可以通過濾除在多數(shù)細(xì)胞中不表達(dá)的基因 (通常認(rèn)為這些基因?qū)?xì)胞異質(zhì)性分析貢獻(xiàn)不大)，大大減少該數(shù)目珊肃。設(shè)置此閾值的準(zhǔn)則是使用感興趣的最小細(xì)胞亞群大小荣刑，并為dropout事件留出一些余地 (生信寶典注：比最小細(xì)胞亞群大小數(shù)字再小一點)。例如伦乔，濾除掉在少于20個細(xì)胞中表達(dá)的基因可能會使鑒定少于20個細(xì)胞的細(xì)胞亞群變得困難厉亏。對于dropout比率高的數(shù)據(jù)集，此閾值也可能會使不太大的細(xì)胞亞群難以被檢測到烈和。閾值的選擇應(yīng)隨待分析數(shù)據(jù)集中的細(xì)胞數(shù)量和既定的下游分析而定叶堆。
另外可以直接對計數(shù)矩陣 (count matrix)執(zhí)行進(jìn)一步的質(zhì)量控制。輸入型基因表達(dá)（Ambient gene expression）是指某個barcode檢測到的轉(zhuǎn)錄本是源自其他細(xì)胞在建庫前發(fā)生破裂而釋放到細(xì)胞懸液中的mRNA斥杜。這些增加的外來計數(shù)會影響下游分析虱颗，如Marker基因鑒定或其他差異表達(dá)分析過程沥匈。由于基于液滴的scRNA-seq數(shù)據(jù)集中存在大量空液滴，因此可以通過空液滴建模分析細(xì)胞懸液中的RNA構(gòu)成和豐度來校正這一影響忘渔。最近開發(fā)的SoupX使用這種方法直接校正count matrix高帖。另外，在下游分析中直接忽略這些有強(qiáng)影響的輸入型基因也是處理這個問題的一個實用方法畦粮。
質(zhì)控是用于保證下游分析時數(shù)據(jù)質(zhì)量足夠好散址。由于不能先驗地確定什么是足夠好的數(shù)據(jù)質(zhì)量，所以只能基于下游分析結(jié)果（例如宣赔，細(xì)胞簇注釋）來對其進(jìn)行判斷预麸。因此，在分析數(shù)據(jù)時可能需要反復(fù)調(diào)整參數(shù)進(jìn)行質(zhì)量控制 (生信寶典注：反復(fù)分析儒将，多次分析吏祸，是做生信的基本要求。這也是為啥需要上易生信培訓(xùn)班而不是單純依賴公司的分析钩蚊。)贡翘。從寬松的QC閾值開始并研究這些閾值的影響，然后再執(zhí)行更嚴(yán)格的QC砰逻，通常是有益的鸣驱。這種方法對于細(xì)胞種類混雜的數(shù)據(jù)集特別有用，在這些數(shù)據(jù)集中蝠咆，特定細(xì)胞類型或特定細(xì)胞狀態(tài)可能會被誤解為低質(zhì)量的異常細(xì)胞踊东。在低質(zhì)量數(shù)據(jù)集中，可能需要嚴(yán)格的質(zhì)量控制閾值刚操。數(shù)據(jù)集的質(zhì)量可以通過實驗質(zhì)量控制指標(biāo)來確定 (see Appendix Supplementary Text S2递胧，簡單來說就是原始數(shù)據(jù)測序質(zhì)量(FastQC的結(jié)果、序列比對率和比對到已經(jīng)注釋的基因區(qū)的reads比率赡茸、實際檢測到的細(xì)胞數(shù)和預(yù)估細(xì)胞數(shù)是否吻合)缎脾。在QC閾值迭代優(yōu)化過程中，要避免數(shù)據(jù)挑選 (data peeking)占卧。QC閾值不應(yīng)用于改善統(tǒng)計檢驗的結(jié)果遗菠。相反，可以根據(jù)數(shù)據(jù)集可視化和聚類中QC變量的分布來評估QC選取的閾值是否合理华蜒。
圖2. 質(zhì)控變量的分布圖 (小鼠腸道上皮細(xì)胞數(shù)據(jù)集)
(A) 每個細(xì)胞檢測到的總分子數(shù)的分布用直方圖展示辙纬。右上角的子圖展示了總分子數(shù)小于4,000的部分的分布。因為在count數(shù)為1200處有個峰叭喜，所以閾值設(shè)置為了1,500贺拣。(B) 每個細(xì)胞檢測到的基因數(shù)的直方圖展示。在400個基因處有個小峰（噪音峰），閾值設(shè)置為700譬涡。（生信寶典注：這閾值設(shè)置的隨意性也沒誰了~~~）(C) 每個細(xì)胞檢測到的總分子數(shù)從高到底繪制rank plot闪幽，類似于Cell Ranger檢測并過濾空液滴的log-log plot。在總分子數(shù)為1500時涡匀，存在一個快速降低的拐點盯腌，則1500為篩選閾值。(D) 同時展示檢測到的基因數(shù)（縱軸）陨瘩、總分子數(shù)（橫軸）和線粒體基因的比例 (顏色)腕够。線粒體基因只在低基因數(shù)和低總分子數(shù)的細(xì)胞中比例高。這些細(xì)胞已經(jīng)被上面設(shè)置的總分子數(shù)和基因數(shù)閾值過濾掉了舌劳。質(zhì)控參數(shù)聯(lián)合可視化顯示更低的基因數(shù)篩選閾值可能也是足夠的帚湘。
陷阱和建議：
- 通過可視化檢測到的基因數(shù)量、總分子數(shù)和線粒體基因的表達(dá)比例的分布中的異常峰來執(zhí)行QC甚淡。聯(lián)合考慮這3個變量大诸，而不是單獨考慮一個變量。
- 盡可能設(shè)置寬松的QC閾值材诽；如果下游聚類無法解釋時再重新設(shè)定嚴(yán)格的QC閾值。
- 如果樣品之間的QC變量分布不同（存在多個強(qiáng)峰）恒傻，則需要考慮樣品質(zhì)量差異脸侥，應(yīng)按照
Plasschaert et al. (2018)的方法為每個樣品分別確定QC閾值。
標(biāo)準(zhǔn)化
計數(shù)矩陣中的每個數(shù)值代表細(xì)胞中一個mRNA分子被成功捕獲盈厘、逆轉(zhuǎn)錄和測序（Box 1）睁枕。由于每個操作步驟固有的可變性，即便同一個細(xì)胞測序兩次獲得的計數(shù)深度也可能會有所不同沸手。因此外遇，當(dāng)基于原始計數(shù)數(shù)據(jù)比較細(xì)胞之間的基因表達(dá)時，得到的差異可能來自于技術(shù)原因契吉。Normalization可以通過調(diào)整計數(shù)數(shù)據(jù) (scaling count data)等解決這一問題跳仿，以獲得細(xì)胞之間可比的相對基因表達(dá)豐度。
普通轉(zhuǎn)錄組分析有很多可用的標(biāo)準(zhǔn)化方法捐晶。盡管其中一些方法已應(yīng)用于scRNA-seq分析菲语，但單細(xì)胞數(shù)據(jù)特有的變異特征（例如technical dropouts: zero counts due to sampling ）催生了scRNA-seq特異性標(biāo)準(zhǔn)化方法的發(fā)展。(什么惑灵？你做的差異基因方法不合適山上？)
最常用的標(biāo)準(zhǔn)化方法是測序深度標(biāo)準(zhǔn)化，也稱為“每百萬計數(shù)”或CPM normalization英支。該方法來自普通轉(zhuǎn)錄組表達(dá)分析佩憾，使用每個細(xì)胞的測序深度作為size factor對計數(shù)數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化。CPM標(biāo)準(zhǔn)化假設(shè)數(shù)據(jù)集中的所有細(xì)胞最初都包含相等數(shù)量的mRNA分子，并且計數(shù)深度差異來源于技術(shù)問題妄帘。這一方法的變體有：把size factor放大10^6倍如CPM楞黄，或size factor乘以數(shù)據(jù)集中所有細(xì)胞的測序深度的中位數(shù)。downsampling方法也基于同樣的假設(shè)寄摆，從原始數(shù)據(jù)中隨機(jī)抽取預(yù)設(shè)的等量reads以保證所有細(xì)胞測序深度相同谅辣。downsampling方法在扔掉一些數(shù)據(jù)的同時增加了technical dropout rates，而CPM或其它全局的標(biāo)準(zhǔn)化方式則不會影響婶恼。因此downsampling方法可能更真實地反應(yīng)了相同測序深度下不同細(xì)胞的表達(dá)譜特征桑阶。
由于單細(xì)胞數(shù)據(jù)集通常由大小和分子數(shù)不同的異質(zhì)細(xì)胞群體組成，因此通常需要更復(fù)雜的標(biāo)準(zhǔn)化方法勾邦。例如蚣录，Weinreb et al對CPM算法進(jìn)行了擴(kuò)展，在計算size factors時排除在任何細(xì)胞中總計數(shù)大于5％的基因眷篇。這一方法屏蔽掉少數(shù)高表達(dá)基因?qū)傮w表達(dá)變化的影響萎河。軟件包Scran的pooling‐based size factor方法對細(xì)胞異質(zhì)性的影響處理更好。首先把細(xì)胞合并到一起避免technical dropout效應(yīng)蕉饼，然后基于基因表達(dá)的線性回歸模型估算size factor虐杯。這一方法允許細(xì)胞有少于50%的差異表達(dá)基因，并且在不同的測試評估研究中這一標(biāo)準(zhǔn)化方法都表現(xiàn)最好昧港。評估發(fā)現(xiàn)擎椰，scran比其他標(biāo)準(zhǔn)化方法對后續(xù)批次校正和差異表達(dá)分析效果更好。并且在scran作者小范圍的比較中也展現(xiàn)出這個方法穩(wěn)定性比較好创肥。
CPM, high‐count filtering CPM和scran 使用線性全局縮放對計數(shù)數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化达舒。還有非線性標(biāo)準(zhǔn)化方式可以處理更多的混雜因素的影響。許多此類方法涉及對計數(shù)數(shù)據(jù)進(jìn)行參數(shù)化建模叹侄。例如巩搏， Mayer et al 使用技術(shù)相關(guān)協(xié)變量（例如測序深度和每個基因的計數(shù)）作為模型參數(shù)對count data擬合負(fù)二項模型。模型的殘差作為標(biāo)準(zhǔn)化后的基因表達(dá)值趾代。這種方法可以在標(biāo)準(zhǔn)化的同時校正技術(shù)差異和生物來源的差異（例如批次校正或?qū)?xì)胞周期影響的校正）贯底。已表明，非線性標(biāo)準(zhǔn)化方法優(yōu)于全局標(biāo)準(zhǔn)化方法撒强，尤其是在有較強(qiáng)批次效應(yīng)影響的情況下丈甸。因此，非線性標(biāo)準(zhǔn)化方法特別適用于基于板的scRNA-seq數(shù)據(jù)尿褪，因為不同板之間往往會產(chǎn)生批次效應(yīng)睦擂。此外，基于板的數(shù)據(jù)與基于液滴的數(shù)據(jù)相比杖玲，每個細(xì)胞的測序深度變化更大顿仇。理論上非線性標(biāo)準(zhǔn)化方法或諸如downsampling的放法更適合于板的數(shù)據(jù) (plate)，但仍需要進(jìn)行比較研究以確認(rèn)這一觀點。在本教程中臼闻，我們傾向于將標(biāo)準(zhǔn)化和數(shù)據(jù)校正（批處理校正鸿吆、噪聲校正等）步驟分開，以強(qiáng)調(diào)數(shù)據(jù)的不同處理階段述呐。因此惩淳，我們專注于全局縮放標(biāo)準(zhǔn)化方法。
我們不能期望一個標(biāo)準(zhǔn)化方法適用于所有類型的scRNA-seq數(shù)據(jù)乓搬。例如思犁，Vieth et al 表明，reads data和count data適用不同的模型进肯。確實激蹲，Cole et al 發(fā)現(xiàn)不同的標(biāo)準(zhǔn)化方法只對適合的數(shù)據(jù)集表現(xiàn)最佳，并認(rèn)為應(yīng)使用他們的scone工具對數(shù)據(jù)集進(jìn)行評估以選擇適用的標(biāo)準(zhǔn)化方法江掩。此外学辱，scRNA-seq技術(shù)可分為全長和3'富集方法。全長方案的數(shù)據(jù)可能會受益于考慮了基因長度的標(biāo)準(zhǔn)化方法环形，而3'富集測序數(shù)據(jù)則不然策泣。TPM歸一化是全長scRNA-seq數(shù)據(jù)常用的歸一化方法，它來自bulk RNA-seq分析抬吟。
標(biāo)準(zhǔn)化是對細(xì)胞計數(shù)數(shù)據(jù)進(jìn)行縮放處理以使其在細(xì)胞之間可比萨咕，也可以在基因?qū)用鎸蛴嫈?shù)進(jìn)行歸一化 (scale)以便于基因內(nèi)部進(jìn)行直接比較∞志基因歸一化是指一個基因減去其在所有樣品表達(dá)的均值然后除以其在所有樣品表達(dá)值的標(biāo)準(zhǔn)差任洞。歸一化后蓄喇，這個基因在所有樣品表達(dá)值均值為0发侵，用單位方差形式表示其表達(dá)值。歸一化后妆偏，所有基因在下游分析時權(quán)重是一樣的刃鳄。（生信寶典注：歸一化后，原始基因的表達(dá)豐度信息就沒了钱骂，換成了無標(biāo)度的標(biāo)準(zhǔn)差的表示叔锐。）是否對基因進(jìn)行歸一化目前尚無達(dá)成共識。盡管流行Seurat教程通常應(yīng)用gene scaling见秽，但Slingshot方法的作者在其教程中選擇了不對基因進(jìn)行scaling愉烙。兩種選擇的爭議點在：所有基因不論表達(dá)高低在進(jìn)行下游分析時權(quán)重一致 (生信寶典注：scale后所有基因權(quán)重一致)，還是基因表達(dá)量的絕對值對下游分析也有貢獻(xiàn) (生信寶典注：突出高表達(dá)基因?qū)ο掠蔚母笥绊?/em>) 解取。為了從數(shù)據(jù)中保留盡可能多的生物相關(guān)信息步责，在本教程中，我們選擇不進(jìn)行歸一化 （scaling genes）這一步操作。
標(biāo)準(zhǔn)化后蔓肯，數(shù)據(jù)矩陣通常進(jìn)行log(x+1)轉(zhuǎn)換遂鹊。此轉(zhuǎn)換具有三個重要作用。
首先蔗包，對數(shù)轉(zhuǎn)換后的表達(dá)式值之間的差值可對應(yīng)于對數(shù)轉(zhuǎn)換后的倍數(shù)變化秉扑，這是衡量基因表達(dá)變化的常用方法。
其次调限，對數(shù)轉(zhuǎn)換可減輕（但不能消除）單細(xì)胞數(shù)據(jù)的均值-方差關(guān)系 (
mean-variance relationship) （生信寶典注：均值-方差關(guān)系展示數(shù)據(jù)在特征空間的分布關(guān)系舟陆。方差越大數(shù)據(jù)分布越廣，后續(xù)采用線性回歸算法時效果越差旧噪。）吨娜。最后，對數(shù)轉(zhuǎn)換可以減少數(shù)據(jù)的偏態(tài)分布淘钟，從而使數(shù)據(jù)近似于正態(tài)分布宦赠，更符合許多下游分析工具對數(shù)據(jù)分布的假設(shè)要求。
盡管scRNA-seq數(shù)據(jù)實際上不是對數(shù)正態(tài)分布的米母，但這三個效果使對數(shù)轉(zhuǎn)換成為一種粗略但有用的工具勾扭。這一有用性在下游差異表達(dá)分析或批次校正時有更好的體現(xiàn)。但是應(yīng)該注意的是铁瞒，數(shù)據(jù)的對數(shù)轉(zhuǎn)換會在數(shù)據(jù)中引入虛假的差異表達(dá)結(jié)果妙色。而且如果size factor在組間差異更大時影響尤其明顯。
陷阱和建議：
- 我們建議使用
scran來標(biāo)準(zhǔn)化非全長數(shù)據(jù)集慧耍。另一種選擇是通過scone評估標(biāo)準(zhǔn)化方法未辆，尤其是對基于板的數(shù)據(jù)集。全長scRNA-seq數(shù)據(jù)可以使用bulk轉(zhuǎn)錄組的方法校正基因長度的影響蟹腾。
- 基因表達(dá)歸一化（
scale）到均值為0和單位方差尚無共識私杜。我們傾向于不進(jìn)行scale操作。
- 標(biāo)準(zhǔn)化后的數(shù)據(jù)應(yīng)進(jìn)行
log(x+1)轉(zhuǎn)換泌豆，以使數(shù)據(jù)更符合正態(tài)分布定庵，滿足下游分析方法的初始假設(shè)。
數(shù)據(jù)校正和整合
如前所述踪危，數(shù)據(jù)標(biāo)準(zhǔn)化是去除實驗過程中隨機(jī)性的影響 (count sampling)蔬浙。但是，標(biāo)準(zhǔn)化后的數(shù)據(jù)可能仍然包含有與研究不相干的因素帶來的影響贞远。數(shù)據(jù)校正的目的就是進(jìn)一步去除技術(shù)因素和非關(guān)注的生物學(xué)混雜因素畴博，例如批次、dropout或細(xì)胞周期不同帶來的影響 (如何火眼金睛鑒定那些單細(xì)胞轉(zhuǎn)錄組中的混雜因素)蓝仲。需要注意的是這些混淆因素并不總是需要校正俱病。相反蜜唾，考慮在分析中移除哪些因素取決于下游分析目的。我們建議分別考慮對生物學(xué)和技術(shù)混雜因素 (covariates)的校正庶艾，因為它們用于不同的目的并代表不同的分析挑戰(zhàn)袁余。
去除與研究不相干的生物因素的影響
盡管校正實驗中技術(shù)因素的影響可能對揭示潛在的生物學(xué)信號至關(guān)重要，但對不關(guān)注的生物學(xué)因素的校正可用來挑選出特定的感興趣的目標(biāo)生物學(xué)信息咱揍。最常見的生物因素校正是消除細(xì)胞周期對轉(zhuǎn)錄組中基因表達(dá)的影響颖榜。可以使用Scanpy和Seurat對每個細(xì)胞的細(xì)胞周期評分進(jìn)行簡單的線性回歸校正或通過應(yīng)用了更復(fù)雜的混合模型的專用程序包如scLVM或f-scLVM進(jìn)行校正煤裙。用于計算細(xì)胞周期評分的標(biāo)記基因列表可在文獻(xiàn)中獲取 (Seurat亮點之細(xì)胞周期評分和回歸)掩完。這些方法還可用于校正其他已知的生物因素的影響，例如線粒體基因表達(dá)硼砰，通常是細(xì)胞應(yīng)激的標(biāo)識且蓬。
在校正特定的生物因素影響之前，應(yīng)考慮幾個方面题翰。首先恶阴，校正生物學(xué)因素影響并不總是有助于解釋scRNA-seq數(shù)據(jù)。雖然消除細(xì)胞周期影響可以改善對發(fā)育軌跡的推斷豹障，但細(xì)胞周期信號也可以提供有意義的生物學(xué)信息冯事。例如，可以根據(jù)細(xì)胞周期評分確定增殖細(xì)胞群 (在Github的案例研究中有講)血公。同樣昵仅，必須根據(jù)具體研究的對象理解生物學(xué)信號。假設(shè)生物過程發(fā)生在同一有機(jī)體內(nèi)累魔，則這些過程之間可能存在依賴性摔笤。因此，校正一個過程的影響可能會無意間掩蓋另一個過程的信號垦写。另外吕世，也有研究者提出細(xì)胞大小對轉(zhuǎn)錄組的影響通常也歸因于細(xì)胞周期不同。因此梯澜，通過標(biāo)準(zhǔn)化或?qū)Ｓ霉ぞ撸ɡ?code>cgCorrect）在校正細(xì)胞大小的同時也部分校正了scRNA-seq數(shù)據(jù)中細(xì)胞周期的影響寞冯。
去除技術(shù)影響
用于移除不相干生物因素影響的回歸模型也可應(yīng)用于校正實驗技術(shù)因素的影響渴析。單細(xì)胞數(shù)據(jù)中最突出的技術(shù)影響因素是測序深度和批次晚伙。盡管標(biāo)準(zhǔn)化可以使得細(xì)胞之間的基因定量數(shù)據(jù)可比，但測序深度的影響仍保留在數(shù)據(jù)中俭茧。這種測序深度效應(yīng)既可以是生物學(xué)自身的影響咆疗，也可能是技術(shù)帶來的影響。例如母债，細(xì)胞大小可能不同午磁，因此mRNA分子數(shù)量也可能不同尝抖。而且，標(biāo)準(zhǔn)化后技術(shù)帶來的計數(shù)影響可能依然存在迅皇，因為沒有辦法推斷由于實驗過程中的隨機(jī)性帶來的未檢測到的基因的表達(dá)昧辽。對測序深度進(jìn)行回歸校正可以提高軌跡推斷算法的性能，該算法依賴于查找細(xì)胞之間的過渡（具體見Github上的案例研究）登颓。在校正多個協(xié)變量（干擾因素）（例如細(xì)胞周期和測序深度）時搅荞，應(yīng)在單個步驟中對所有協(xié)變量（干擾因素）執(zhí)行回歸，以解決協(xié)變量（干擾因素）之間的依賴性框咙。
消除測序深度影響的另一個策略是使用更嚴(yán)格的標(biāo)準(zhǔn)化程序咕痛，例如downsampling或非線性歸一化方法 (見前面標(biāo)準(zhǔn)化部分)。這些方法可能更適用于基于板的scRNA-seq數(shù)據(jù)集喇嘱，因為細(xì)胞之間較大的測序深度差異可能掩蓋細(xì)胞之間的異質(zhì)性茉贡。
批次效應(yīng)和數(shù)據(jù)整合
當(dāng)將細(xì)胞分組操作時可能會帶來批次效應(yīng)（DESeq2差異基因分析和批次效應(yīng)移除），比如不同芯片上的細(xì)胞者铜、不同測序通道中的細(xì)胞或在不同時間點收集的細(xì)胞都?xì)w類于不同的組腔丧。實驗操作過程中細(xì)胞所經(jīng)歷的不同環(huán)境可能會影響轉(zhuǎn)錄組的測量結(jié)果或甚至影響細(xì)胞自身的轉(zhuǎn)錄變化。所產(chǎn)生的影響存在多個層面：同一實驗不同的細(xì)胞組作烟、同一實驗室的不同實驗或不同實驗室的數(shù)據(jù)集之間悔据。在這里，我們把第一種情況與后面兩種情況區(qū)分開俗壹。校正同一實驗中樣品或細(xì)胞之間的批次效應(yīng)是bulk RNA測序批次效應(yīng)的一種經(jīng)典方案科汗。我們將其與整合來自多個實驗的數(shù)據(jù)（稱為數(shù)據(jù)整合）區(qū)分開。通常批次效應(yīng)校正使用線性方法绷雏，而非線性方法則用于數(shù)據(jù)整合头滔。
最近對經(jīng)典批次校正方法的比較顯示，ComBat在中低到中復(fù)雜度的單細(xì)胞實驗中也表現(xiàn)良好 (收藏 Combat開發(fā)者涎显、北大生信平臺” 單細(xì)胞分析坤检、染色質(zhì)分析” 視頻和PPT分享)。ComBat構(gòu)建了基因表達(dá)的線性模型期吓，其中批次貢獻(xiàn)在數(shù)據(jù)的均值和方差中均得到校正（圖3）早歇。與采用什么計算方法無關(guān)，批次校正的最佳方法是通過巧妙的實驗設(shè)計來避免存在不同批次讨勤。通過合并不同實驗條件和樣品的細(xì)胞一起進(jìn)行后續(xù)操作可以避免批次效應(yīng)箭跳。使用細(xì)胞標(biāo)記等策略或根據(jù)基因組變異信息，可以對實驗中合并的細(xì)胞進(jìn)行拆分潭千。（生信寶典注：最不建議的實驗設(shè)計方式是對照組樣品是一個批次谱姓，處理組樣品是一個批次；如果這么做刨晴，將不能確定差異是來源于批次屉来，還是真的存在生物差異路翻。）
圖3. 批次校正前后UMAP可視化結(jié)果展示
細(xì)胞按其來源的樣本上色。在批次校正前可以看到顏色區(qū)分很明顯（細(xì)胞所來源的樣品對細(xì)胞分群影響大）茄靠，表明批次影響比較大茂契。使用ComBat校正批次后，同一個顏色分布比較散了慨绳，細(xì)胞所來源的樣品對細(xì)胞分群影響變小了账嚎。
與批次校正相比，數(shù)據(jù)整合面臨的另一個挑戰(zhàn)是數(shù)據(jù)集之間的組成差異 (compositional differences)儡蔓。估計批次效應(yīng)時郭蕉，ComBat使用同一批次的所有細(xì)胞來擬合批次校正參數(shù)。這種方法將批次效應(yīng)與數(shù)據(jù)集之間非共有的細(xì)胞類型或狀態(tài)之間的生物學(xué)差異混淆喂江。數(shù)據(jù)整合方法（Cell 深度| 一套普遍適用于各類單細(xì)胞測序數(shù)據(jù)集的錨定整合方案）召锈，例如典型相關(guān)分析（Canonical Correlation Analysis, CCA），相互最近鄰（Mutual Nearest Neighbours, MNN), Scanorama, RISC, scGen, LIGER, BBKNN和Harmony已經(jīng)開發(fā)用于解決這個問題获询。盡管數(shù)據(jù)整合方法也可應(yīng)用于簡單的批次校正處理涨岁，但是鑒于非線性數(shù)據(jù)整合方法會增大自由度，使用時需要注意過度校正問題吉嚣。例如梢薪，在較簡單的批次校正設(shè)置中，ComBat的表現(xiàn)優(yōu)于MNN（Buttner et al.尝哆，2019）秉撇。需要進(jìn)一步進(jìn)行數(shù)據(jù)整合和批次校正方法之間的比較研究，才可以評估這些方法的適用范圍秋泄。
缺失值填充
技術(shù)數(shù)據(jù)校正的另一種類型是缺失值填充（也稱為降噪或插補, denoising or imputation）琐馆。單細(xì)胞轉(zhuǎn)錄組的數(shù)據(jù)包含各種噪聲。這種噪音的一個特別突出的來源是dropout恒序。推斷dropouts事件瘦麸，用推斷出的合適的表達(dá)值替換這些零以減少數(shù)據(jù)集中的噪聲成為幾種最新工具的目標(biāo) (MAGIC, DCA, scVI, SAVE, scImpute)。已證明進(jìn)行缺失值填充可改善基因與基因相關(guān)性的估計歧胁。此外滋饲，這一步也可以與標(biāo)準(zhǔn)化、批次校正和其他下游分析整合喊巍，就像在scVI工具中實現(xiàn)的那樣屠缭。盡管大多數(shù)數(shù)據(jù)校正方法都將標(biāo)準(zhǔn)化后的數(shù)據(jù)作為輸入，但是某些缺失值填充方法是基于預(yù)期的負(fù)二項噪聲分布玄糟，需要基于原始計數(shù)矩陣進(jìn)行操作勿她。在應(yīng)用缺失值填充時袄秩，應(yīng)考慮到?jīng)]有一種方法是完美的阵翎。因此逢并，任何方法都可能會對數(shù)據(jù)中的噪聲進(jìn)行過高校正或校正不足。確實郭卫，已有報道表明缺失值填充可能引入錯誤的相關(guān)信號砍聊。鑒于在實際應(yīng)用中難以評估缺失值填充是否得當(dāng)，用戶選擇是否應(yīng)用這一方法也是很大的挑戰(zhàn)贰军。當(dāng)前缺失值填充方法是否能拓展應(yīng)用到大數(shù)據(jù)集還是一個問題玻蝌。鑒于這些考慮，目前尚無關(guān)于應(yīng)如何使用缺失值填充的共識词疼。謹(jǐn)慎的方法是僅在視覺展示數(shù)據(jù)時使用缺失值填充俯树，而非在探索性數(shù)據(jù)分析過程中基于填充的數(shù)據(jù)作出推論或假設(shè)。全面的實驗驗證在這里尤為重要贰盗。
陷阱和建議：
- 僅在進(jìn)行軌跡推斷和校正的生物學(xué)混雜因素不影響感興趣的生物學(xué)過程時才校正這些因素的影響许饿。
- 如果校正的話，所有因素同時校正而不是分別校正技術(shù)和非關(guān)注的生物因素變量舵盈。
- 基于板的數(shù)據(jù)集預(yù)處理時需要校正
count數(shù)的影響陋率，建議采用非線性標(biāo)準(zhǔn)化方法或downsampling方法進(jìn)行標(biāo)準(zhǔn)化。
- 當(dāng)批次之間的細(xì)胞類型和狀態(tài)組成一致時秽晚，建議通過
ComBat執(zhí)行批次校正瓦糟。
- 數(shù)據(jù)整合和批次校正應(yīng)通過不同的方法進(jìn)行。數(shù)據(jù)整合工具可能會過度校正簡單的批次效應(yīng)赴蝇。
- 用戶需要對只在缺失值填充后才能發(fā)現(xiàn)的信號格外注意菩浙。探索性分析時最好不進(jìn)行缺失值填充操作。
特征選擇句伶、降維和可視化
人單細(xì)胞RNA-seq數(shù)據(jù)集可包含多達(dá)25,000個基因的表達(dá)值芍耘。對于一個給定的scRNA-seq數(shù)據(jù)集，其中有許多基因都不能提供有用信息熄阻，并且大多只包含零計數(shù)斋竞。即使在QC步驟中濾除了這些零計數(shù)基因后，單細(xì)胞數(shù)據(jù)集的特征空間也可能超過15,000個維度（即還會剩余15,000多基因）秃殉。為了減輕下游分析工具的計算負(fù)擔(dān)坝初、減少數(shù)據(jù)中的噪聲并方便數(shù)據(jù)可視化，可以使用多種方法來對數(shù)據(jù)集進(jìn)行降維钾军。
特征選擇
scRNA-seq數(shù)據(jù)集降維的第一步通常是特征選擇鳄袍。在此步驟中，對數(shù)據(jù)集基因進(jìn)行過濾僅保留對數(shù)據(jù)的變異性具有信息貢獻(xiàn)的基因（在數(shù)據(jù)中變異大的基因）吏恭。這些基因通常被定義為高變化基因（HVG拗小，highly variable genes）。根據(jù)任務(wù)和數(shù)據(jù)集的復(fù)雜性樱哼，通常選擇1,000到5,000個HVG用于下游分析哀九。Klein et al.的初步結(jié)果表明剿配，下游分析對HVG的數(shù)量不太敏感。在HVG數(shù)量從200到2,400之間選擇不同的數(shù)目時阅束，評估顯示PCA結(jié)果相差不大呼胚。基于此結(jié)果息裸，我們寧愿選擇更多的HVG用于下游分析 （err on the side of：可以借鑒的英語短語）蝇更。
圖EV1. 不同大小數(shù)據(jù)集選擇HVG基因數(shù)量的分布展現(xiàn)。（數(shù)據(jù)來源于手動整理的最近發(fā)表的多篇scRNA-seq文章呼盆，顏色代表不同的技術(shù)方式）
在Scanpy和Seurat中都實現(xiàn)了一種簡單而流行的選擇HVG的方法年扩。在這里，基因按其均值表達(dá)進(jìn)行分組访圃，將每個組內(nèi)方差/均值比最高的基因選為每個分組的HVG常遂。該算法在不同軟件中輸入不同，Seurat需要原始count data挽荠；Cell Ranger需要對數(shù)轉(zhuǎn)換的數(shù)據(jù)克胳。最佳地，應(yīng)在技術(shù)等干擾因素校正之后選擇HVG圈匆，以避免選擇僅由于例如批次效應(yīng)等引入的高可變基因漠另。Yip et al.綜述了其他HVG選擇方法。
特征選擇后跃赚，可以通過專用的降維算法進(jìn)一步對單細(xì)胞表達(dá)矩陣進(jìn)行降維笆搓。這些算法將表達(dá)式矩陣映射到低維空間中，同時以盡可能少的維數(shù)捕獲數(shù)據(jù)中所有的信息纬傲。鑒于單細(xì)胞RNA測序數(shù)據(jù)固有的低維性特征满败，這一方法是合適的。也就是說叹括，細(xì)胞表達(dá)圖譜構(gòu)成的生物形態(tài) (biological manifold)可以使用遠(yuǎn)少于基因數(shù)目的維度信息來展示算墨。降維旨在找到這些維度。
降維有兩個主要目標(biāo)：可視化和信息匯總（summarization）汁雷【秽郑可視化是嘗試在二維或三維空間最優(yōu)地展示數(shù)據(jù)集。降維后的維度值就是數(shù)據(jù)在新的空間進(jìn)行可視化如繪制散點圖時的坐標(biāo)值侠讯。信息匯總沒有規(guī)定輸出的維數(shù)挖藏；但更高的維數(shù)對表示原有數(shù)據(jù)的差異越來越不重要（生信寶典注：可以理解為PCA中各個主成分對于原始數(shù)據(jù)差異的解釋依次降低）。匯總技術(shù)可通過計算數(shù)據(jù)的固有維數(shù)來將數(shù)據(jù)降維到基本組成（主）成分厢漩，從而有助于下游分析膜眠。雖然不應(yīng)使用二維可視化數(shù)據(jù)來匯總數(shù)據(jù)集，但匯總方法獲得的降維數(shù)據(jù)可用來可視化數(shù)據(jù)。另外宵膨，專用的可視化技術(shù)通常會更好地展示數(shù)據(jù)的原始結(jié)構(gòu)和變異性 （還在用PCA降維架谎？快學(xué)學(xué)大牛最愛的t-SNE算法吧, 附Python/R代碼）。
降維后的維度是通過對基因表達(dá)向量進(jìn)行線性或非線性組合生成的 (PCA主成分分析實戰(zhàn)和可視化 附R代碼和測試數(shù)據(jù))柄驻。特別是在非線性情況下狐树，降維后的數(shù)據(jù)難以解釋其生物含義焙压。圖4中顯示了一些常用降維方法的示例應(yīng)用鸿脓。隨著可供選擇的方法的增加，詳細(xì)討論這些方法已超出了本教程的范圍涯曲。相反野哭，我們簡要概述了可能有助于用戶在常見的降維方法之間進(jìn)行選擇時需要考慮的因素。Moon et al. 提供了有關(guān)單細(xì)胞分析降維的更詳細(xì)的綜述 (Manifold learning‐based methods for analyzing single‐cell RNA‐sequencing data. Curr Opin Syst Biol 7: 36–46)幻件。
圖4. scRNA‐seq數(shù)據(jù)常用的可視化方法拨黔。
(A) PCA, (B) t‐SNE, (C) diffusion maps, (D) UMAP and (E) A force‐directed graph layout via ForceAtlas2. 細(xì)胞根據(jù)測序深度上色 (F) 前31個主成分解釋的原始數(shù)據(jù)的差異。圖中的拐點( “elbow”)用于選擇下游分析相關(guān)的主成分绰沥，位于PCs 5 和PCs 7之間篱蝇。
主成分分析（PCA）和diffusion maps是兩種常用的降維方法，在單細(xì)胞分析中也很流行徽曲。主成分分析是一種線性方法零截，通過最大化每個可能維度中捕獲的殘差 (residual variance)來進(jìn)行降維。盡管PCA不能像非線性方法那樣在更少的維度捕獲原始數(shù)據(jù)更多的信息秃臣，但是它是許多當(dāng)前可用的聚類或軌跡推斷工具的基礎(chǔ)涧衙。實際上，PCA通常用作非線性降維方法的預(yù)處理步驟奥此。通常弧哎，PCA通過其前N個主成分來代表原始數(shù)據(jù)集，其中N可以通過elbow算法（參見圖4F）或基于置換檢驗的jackstraw方法確定稚虎。PCA簡單線性化的優(yōu)勢是：降維空間中的距離在該空間的所有區(qū)域具有一致的解釋撤嫩。因此，我們可以將感興趣的統(tǒng)計量與主成分進(jìn)行關(guān)聯(lián)分析蠢终，以評估其重要性非洲。例如，可以將主成分投影到技術(shù)干擾協(xié)變量上蜕径，以研究QC两踏、數(shù)據(jù)校正和標(biāo)準(zhǔn)化過程的效果（Buttner et al.，2019）兜喻，或顯示基因在數(shù)據(jù)集中的重要性 (PCA bi-plot)梦染。diffusion map是一種非線性數(shù)據(jù)降維技術(shù)。由于diffusion component強(qiáng)調(diào)數(shù)據(jù)的轉(zhuǎn)換，因此主要用于諸如細(xì)胞分化之類的連續(xù)過程帕识。通常泛粹，每個diffusion component（即diffusion map 維度）突出顯示不同細(xì)胞群體的異質(zhì)性。
可視化
可視化時一般使用非線性降維方法（圖4）肮疗。scRNA-seq數(shù)據(jù)可視化的最常見的降維方法是t‐SNE （ t‐distributed stochastic neighbour embedding） （還在用PCA降維晶姊？快學(xué)學(xué)大牛最愛的t-SNE算法吧, 附Python/R代碼）。t‐SNE的維度著重于以犧牲全局結(jié)構(gòu)為代價來保留局部相似性 (生信寶典注：PCA則是盡可能多的保留全局差異)伪货。因此们衙，這些可視化可能會夸大細(xì)胞群體之間的差異，而忽略群體之間的潛在聯(lián)系 (t‐SNE dimensions focus on capturing local similarity at the expense of global structure. Thus, these visualizations may exaggerate differences between cell populations and overlook potential connections between these populations)碱呼。另一個困難是對參數(shù)perplexity parameter的選擇蒙挑，因為t-SNE圖會因為這個參數(shù)值不同而顯示出明顯不同的分簇數(shù)。t‐SNE的常見替代方法是UMAP (Uniform Approximation and Projection method)或基于圖的工具愚臀，例如SPRING忆蚀。UMAP和SPRING的力導(dǎo)向布局算法ForceAtlas2(*force‐directed layout algorithm*)可以說是數(shù)據(jù)潛在拓?fù)浣Y(jié)構(gòu)展示最好的方法（Wolf et al，2019）姑裂。在此比較中馋袜，使UMAP與眾不同的是它的速度快和能應(yīng)用于更大規(guī)模數(shù)據(jù)的能力（Becht et al，2018）舶斧。因此欣鳖，在沒有特定生物學(xué)問題限制的情況下，我們將UMAP視為探索性數(shù)據(jù)可視化分析的最佳實踐捧毛。此外观堂，UMAP還可以把數(shù)據(jù)降維到二維以上的新數(shù)據(jù)嗜侮。盡管我們不知道UMAP`在數(shù)據(jù)匯總中的任何應(yīng)用头谜，但它可能是PCA的一個合適的替代方案眷蜈。
在細(xì)胞水平上經(jīng)典可視化方法的替代方法是PAGA (partition‐based graph abstraction）萨脑。事實證明竿秆，此工具可以充分展示數(shù)據(jù)的拓?fù)浣Y(jié)構(gòu)绍绘，同時使用聚簇生成粗粒度的可視化圖像没卸。結(jié)合以上任何一種可視化方法戳晌，PAGA可以生成粗粒度的可視化圖像 （coarse‐grained visualizations）泞辐，從而可以簡化單細(xì)胞數(shù)據(jù)的解釋笔横，尤其是在測序細(xì)胞量大或整合了大量細(xì)胞的情況下。
陷阱和建議：
- 我們建議根據(jù)數(shù)據(jù)集的復(fù)雜程度選擇
1,000至5,000個高可變基因 (HVG)咐吼。
- 當(dāng)將基因表達(dá)值歸一化為
均值為0和單位方差時吹缔，或者將模型擬合的殘差用作標(biāo)準(zhǔn)化表達(dá)值時，不能使用基于基因表達(dá)均-方差的特征選擇方法锯茄。因此厢塘，在選擇HVG之前茶没，必須考慮要執(zhí)行哪些預(yù)處理。
- 信息匯總 (
summarization晚碾，類比于挑選重要的主成分)和可視化應(yīng)使用不同的降維方法抓半。
- 我們建議使用UMAP進(jìn)行探索性分析可視化；使用PCA做為通用數(shù)據(jù)降維方法格嘁；
diffusion maps可以在軌跡推斷時替代PCA笛求。
- UMAP+PAGA是可視化特別復(fù)雜的數(shù)據(jù)集的合適替代方法。
數(shù)據(jù)預(yù)處理步驟
雖然我們在上面按順序概述了scRNA‐seq中常見的預(yù)處理步驟糕簿，但下游分析通常需要不同級別的預(yù)處理數(shù)據(jù)探入，建議根據(jù)下游應(yīng)用調(diào)整預(yù)處理的各個步驟。為了向新用戶說明這種情況冶伞，我們將預(yù)處理分為五個數(shù)據(jù)處理階段：（i）原始數(shù)據(jù)新症，（ii）標(biāo)準(zhǔn)化后的數(shù)據(jù)步氏，（iii）校正后的數(shù)據(jù)响禽，（iv）特征選擇后的數(shù)據(jù)，以及（v）降維后的數(shù)據(jù)荚醒。數(shù)據(jù)處理的這些階段可以歸類為三個預(yù)處理層：測量的數(shù)據(jù)芋类，校正的數(shù)據(jù)和降維的數(shù)據(jù)。細(xì)胞和基因的質(zhì)量控制篩選是必須執(zhí)行的步驟界阁，因此在此處略過侯繁。盡管層的順序代表了scRNA-seq分析的典型工作流程，但也可以跳過某一步或在處理階段的順序中稍作更改泡躯。例如贮竟，單批次數(shù)據(jù)集可能不需要批次校正。在表1中较剃，我們總結(jié)了預(yù)處理數(shù)據(jù)的每一層所適用的下游處理程序咕别。
表1
表1中的預(yù)處理階段分為三組：測量的數(shù)據(jù)，校正后的數(shù)據(jù)和降維后的數(shù)據(jù)写穴。我們將測量的數(shù)據(jù)定義為原始數(shù)據(jù)和處理后但保留了表達(dá)值為零的基因的數(shù)據(jù)惰拱。應(yīng)用細(xì)胞特異的size factor的全局標(biāo)準(zhǔn)化方法即使在log(x+1)轉(zhuǎn)換后仍保留零表達(dá)值。相比之下啊送，校正后的數(shù)據(jù)會去除不需要的變化信息進(jìn)而替換零表達(dá)值偿短。校正后的數(shù)據(jù)層表示數(shù)據(jù)的“最干凈”版本，它是數(shù)據(jù)代表的生物信號的最近似值馋没。我們稱最終預(yù)處理層為降維后的數(shù)據(jù)（reduced data）昔逗。此數(shù)據(jù)層強(qiáng)調(diào)數(shù)據(jù)的主要差異，可以使用降維后的特征集來描述篷朵。
前述特性決定了預(yù)處理數(shù)據(jù)對特定下游應(yīng)用的適用性勾怒。作為最后的預(yù)處理階段，降維后的數(shù)據(jù)是廣泛應(yīng)用的數(shù)據(jù)層。但是控硼，差異表達(dá)分析只能在基因空間內(nèi)進(jìn)行生物學(xué)解釋泽论，而無法（或無法完全）體現(xiàn)在降維后的數(shù)據(jù)中。降維后數(shù)據(jù)的優(yōu)勢在于生物數(shù)據(jù)信息的匯總(summarization)和影響生物數(shù)據(jù)的噪聲的降低卡乾。因此翼悴，降維后的數(shù)據(jù)應(yīng)用于后續(xù)的探索性方法（如可視化、鄰域圖推斷幔妨、聚類）以及計算復(fù)雜的下游分析工具（如軌跡推斷）鹦赎。實際上，許多軌跡推斷方法在工具本身中都包含了降維功能误堡。
單個基因的表達(dá)譜只能在基因空間中進(jìn)行比較古话，這一信息存在于測量的數(shù)據(jù)和校正后的數(shù)據(jù)中。表達(dá)譜可以進(jìn)行視覺和統(tǒng)計學(xué)比較锁施。我們認(rèn)為視覺比較和統(tǒng)計比較應(yīng)該在不同的數(shù)據(jù)層上進(jìn)行陪踩。對于目測基因表達(dá)，校正的數(shù)據(jù)是最合適的悉抵，有助于結(jié)果解釋肩狂。如果對原始數(shù)據(jù)進(jìn)行視覺比較，則要求用戶牢記數(shù)據(jù)中的偏差以準(zhǔn)確解釋結(jié)果姥饰。但是傻谁，此處應(yīng)分別考慮技術(shù)和生物協(xié)變量的校正數(shù)據(jù)。雖然對生物協(xié)變量的校正可能會增加特定生物信號的強(qiáng)度列粪，但其對潛在的生物學(xué)意義表示可能不甚準(zhǔn)確审磁，并且會掩蓋可能相關(guān)的其他生物信號。因此岂座，經(jīng)過生物學(xué)不相關(guān)因素校正的數(shù)據(jù)主要適用于專注于特定生物學(xué)過程（例如軌跡推斷方法）的分析工具态蒂。
基因表達(dá)的統(tǒng)計比較需要基于測量數(shù)據(jù)層。去噪掺逼、校正批次或其他噪聲源的方法總不是太完美吃媒。因此，數(shù)據(jù)校正方法不可避免地會對數(shù)據(jù)進(jìn)行過度或不足校正吕喘，因此會以意想不到的方式改變至少某些基因表達(dá)譜的方差赘那。基因表達(dá)的統(tǒng)計檢驗依賴于將背景方差評估為數(shù)據(jù)噪聲的空模型氯质。隨著數(shù)據(jù)校正趨于減少背景變化（圖EV2）募舟，背景變化被數(shù)據(jù)校正方法過度校正的基因?qū)⒏锌赡鼙辉u估為顯著差異表達(dá)基因。此外闻察，某些數(shù)據(jù)校正方法（例如ComBat）將不吻合于實驗設(shè)計的表達(dá)信號定義為噪聲拱礁，隨后將其從數(shù)據(jù)中刪除琢锋。這一基于實驗設(shè)計的優(yōu)化方法除了可能造成數(shù)據(jù)中的噪聲被低估還會導(dǎo)致對效應(yīng)大小 （effect size）的高估。鑒于這些考慮呢灶，與校正后的數(shù)據(jù)相比吴超，基于測量數(shù)據(jù)進(jìn)行差異檢驗是一種更為保守的方法。使用測量的數(shù)據(jù)時鸯乃，可以并且應(yīng)該在差異統(tǒng)計分析模型中考慮技術(shù)協(xié)變量 （參考DESeq2差異基因分析和批次效應(yīng)移除）鲸阻。
圖EV2. 批次校正和降噪后變異系數(shù)的變化
負(fù)值代表數(shù)據(jù)校正后變異系數(shù)降低。第一行展示的是ComBat校正(A) mouse intestinal epithelium (mIE) 和 (B) mouse embryonic stem cell (mESC) 數(shù)據(jù)變異系數(shù)的變化缨睡。第二行展示的是DCA降噪后(C) mIE 和 (D) mESC 數(shù)據(jù)變異系數(shù)的變化鸟悴。
最近一篇單細(xì)胞差異表達(dá)分析方法的文章也支持以上觀點, 它只使用原始或標(biāo)準(zhǔn)化后的矩陣作為輸入 (Soneson & Robinson, 2018)。這篇研究中標(biāo)準(zhǔn)化后的數(shù)據(jù)僅限于全局標(biāo)準(zhǔn)化方法奖年。然而當(dāng)前可用的非線性標(biāo)準(zhǔn)化方法模糊了標(biāo)準(zhǔn)化和數(shù)據(jù)校正的界限（具體見標(biāo)準(zhǔn)化部分）细诸。這樣標(biāo)準(zhǔn)化之后的數(shù)據(jù)已不適合用于差異基因分析。
陷阱和建議：
- 原始測量的數(shù)據(jù)用于差異基因統(tǒng)計檢驗分析陋守；
- 校正后的數(shù)據(jù)用于數(shù)據(jù)的可視化比較震贵；
- 降維后的數(shù)據(jù)用于其他基于數(shù)據(jù)形態(tài)（
biological data manifold）的下游分析。
下游分析
經(jīng)過預(yù)處理后嗅义，我們稱之為下游分析的方法指應(yīng)用于生物學(xué)發(fā)現(xiàn)并描述潛在的生物學(xué)系統(tǒng)的方法屏歹。通過將可以解釋的模型擬合到數(shù)據(jù)中獲得相應(yīng)的結(jié)論隐砸，比如有相似基因表達(dá)譜的細(xì)胞群代表一個細(xì)胞簇之碗、相似細(xì)胞之間基因表達(dá)的微小變化指示連續(xù)（分化）軌跡；或具有相關(guān)表達(dá)趨勢的基因指示共調(diào)控等季希。
下游分析可分為細(xì)胞水平和基因水平的方法褪那，如圖5所示。細(xì)胞水平分析通常著重于兩種結(jié)構(gòu)的描述：簇和軌跡式塌。這些結(jié)構(gòu)又可以在細(xì)胞和基因水平上進(jìn)行分析博敬，即簇分析和軌跡分析方法。廣義地講峰尝，簇分析方法試圖將細(xì)胞分類為離散的多個組來解釋數(shù)據(jù)的異質(zhì)性偏窝。相反，在軌跡分析中武学，數(shù)據(jù)被視為細(xì)胞連續(xù)變化動態(tài)過程的一個個快照祭往，軌跡分析方法研究了這一連續(xù)變化過程。
圖5. 下游分析方法總覽（藍(lán)色背景的方法是基因水平的分析方法）火窒。
在此硼补，我們在詳細(xì)描述獨立于這些細(xì)胞結(jié)構(gòu)進(jìn)行的基因水平分析之前，先描述細(xì)胞水平和基因水平的聚類和軌跡分析工具熏矿。
聚類分析
聚類已骇。將細(xì)胞聚類成簇通常是任何單細(xì)胞分析的第一個中間結(jié)果离钝。聚類成簇使我們可以推斷成員細(xì)胞的身份。簇是通過基于細(xì)胞基因表達(dá)譜的相似性將細(xì)胞分組得到的褪储。表達(dá)譜相似性是通過對將降維的數(shù)據(jù)進(jìn)行距離度量確定的卵渴。相似性評分的一個常見示例是在主成分降維后的表達(dá)空間上計算的歐氏距離。存在兩種根據(jù)這些相似性分?jǐn)?shù)生成細(xì)胞簇的方法：聚類算法和基于圖的社群檢測方法 （community detection methods）鲤竹。
聚類是直接基于距離矩陣的經(jīng)典無監(jiān)督機(jī)器學(xué)習(xí)問題奖恰。通過最小化簇內(nèi)距離或在降維后的表達(dá)空間中鑒定密集區(qū)域，將細(xì)胞分配給簇宛裕。流行的k均值聚類算法通過確定聚類質(zhì)心并將細(xì)胞分配給最近的聚類質(zhì)心瑟啃，將細(xì)胞分為k個類。質(zhì)心位置經(jīng)過迭代優(yōu)化揩尸。這種方法需要輸入預(yù)期的簇數(shù)蛹屿，但這通常是未知的，必須進(jìn)行啟發(fā)式校準(zhǔn)岩榆。k均值在單細(xì)胞數(shù)據(jù)中的應(yīng)用因所使用的距離度量而異错负。標(biāo)準(zhǔn)距離度量是使用歐氏距離，替代距離包括余弦相似度勇边、基于相關(guān)的距離度量或SIMLR方法犹撒，該方法使用Gaussian kernels為每個數(shù)據(jù)集計算距離度量。最近的比較表明粒褒，與k均值一起使用或作為Gaussian kernels分析的基礎(chǔ)時识颊，基于相關(guān)的距離可能會勝過其他距離度量方法。（基因共表達(dá)聚類分析及可視化）
社群檢測方法是一種圖分類算法奕坟，依賴于單細(xì)胞數(shù)據(jù)的網(wǎng)絡(luò)圖表示祥款。圖是使用K最近鄰方法（KNN圖）獲得的。細(xì)胞在圖中表示為節(jié)點月杉。每個細(xì)胞與其K個最相似的細(xì)胞相連刃跛，通常對主成分降維后的數(shù)據(jù)計算歐氏距離作為相似性度量。根據(jù)數(shù)據(jù)集的大小苛萎，K通常設(shè)置為5到100個最近的鄰居桨昙。獲得的KNN結(jié)果圖捕獲了表達(dá)圖譜的基礎(chǔ)拓?fù)浣Y(jié)構(gòu)（Wolf et al., 2019）。表達(dá)譜相似的細(xì)胞集合對應(yīng)于圖的緊密連接區(qū)域腌歉，然后使用社群檢測方法檢測圖中這些緊密連接區(qū)域蛙酪。社群檢測方法通常比聚類方法要快，因為只有相鄰的細(xì)胞對會被視為屬于同一群集究履。因此滤否，這種方法大大減少了可能的細(xì)胞簇的搜索空間。
在開創(chuàng)性的PhenoGraph方法（Levine et al.最仑，2015）出現(xiàn)之后藐俺，聚類單細(xì)胞數(shù)據(jù)集的標(biāo)準(zhǔn)方法已演變成多分辨率模塊優(yōu)化算法 (multi‐resolution modularity optimization）炊甲，即在單細(xì)胞KNN圖中應(yīng)用Louvain算法。此方法是在Scanpy和Seurat單細(xì)胞分析平臺中應(yīng)用的默認(rèn)聚類方法欲芹。評估發(fā)現(xiàn)卿啡，這一方法應(yīng)用于單細(xì)胞RNA測序數(shù)據(jù)以及流式細(xì)胞和質(zhì)譜數(shù)據(jù)分析時優(yōu)于其他聚類方法。從概念上講菱父，Louvin算法將組內(nèi)細(xì)胞連接數(shù)高于預(yù)期的一組細(xì)胞視作一個簇颈娜。(生信寶典注：Louvain algorithm算法采用迭代方式，先計算每個點加入到其鄰居節(jié)點所在社區(qū)帶來的模塊性評估收益浙宜，然后將其加入收益最大的鄰居節(jié)點所在社區(qū)官辽，對所有未歸類點重復(fù)進(jìn)行這一過程，直至結(jié)果穩(wěn)定粟瞬。然后再把每個社區(qū)作為一個節(jié)點同仆，重復(fù)上一步。) 優(yōu)化后的模塊鑒定函數(shù)包括一個分辨率參數(shù)裙品，該參數(shù)允許用戶確定簇的近似大小俗批。通過對KNN圖取子集，還可以只聚類一部分特殊的細(xì)胞簇市怎。這種子集聚簇模式可以允許用戶識別某個細(xì)胞類型簇內(nèi)的細(xì)胞狀態(tài)岁忘，但也可能會發(fā)現(xiàn)僅由數(shù)據(jù)中的噪聲帶來的聚類模式。
陷阱和建議：
- 我們建議在單細(xì)胞KNN圖上執(zhí)行Louvain社群檢測算法獲得聚類簇区匠。
- 細(xì)胞聚簇不一定使用單一分辨率進(jìn)行干像。對特定細(xì)胞簇進(jìn)行細(xì)化聚類是鑒定數(shù)據(jù)集中亞結(jié)構(gòu)的有效方法，可以發(fā)現(xiàn)細(xì)胞亞群辱志。
細(xì)胞簇注釋
在基因水平上蝠筑，對每個簇的分析是鑒定其標(biāo)記基因 (Marker genes)。這些所謂的標(biāo)記基因代表了細(xì)胞簇的特征揩懒，用于給細(xì)胞簇一個有生物學(xué)意義的標(biāo)簽（Celaref | 單細(xì)胞測序細(xì)胞類型注釋工具）。該標(biāo)簽表示集群中細(xì)胞的身份挽封。由于任何聚類算法都會聚類出細(xì)胞簇已球，因此聚類獲得的生物簇的準(zhǔn)確性只能通過其生物學(xué)注釋進(jìn)行衡量 (生信寶典注：這也是前面和易生信課程中反復(fù)強(qiáng)調(diào)的，細(xì)胞過濾時標(biāo)準(zhǔn)盡量松一些辅愿，根據(jù)聚類結(jié)果回看之前的參數(shù)設(shè)置是否合理智亮。)。
盡管把單細(xì)胞數(shù)據(jù)中鑒定到的簇歸類為一個有生物意義的細(xì)胞類型是一個誘人的結(jié)果点待，但是細(xì)胞身份的確定不是那么簡單（Wagner et al., 2016阔蛉；Clevers et al., 2017）。首先癞埠，并不總是清楚細(xì)胞類型定義的精確程度状原。例如聋呢，“T細(xì)胞”可能是某些人滿意的細(xì)胞類型標(biāo)記，但其他人可能會在數(shù)據(jù)集中尋找T細(xì)胞亞型并區(qū)分CD4+T細(xì)胞和CD8+T細(xì)胞 颠区。另外削锰，同一細(xì)胞類型的不同狀態(tài)可能會被分到不同的細(xì)胞簇”侠常基于上述原因器贩，最好使用術(shù)語cell identities而不是cell types。在對細(xì)胞進(jìn)行聚類和注釋之前朋截，用戶必須確定要注釋到的細(xì)胞信息精確水平蛹稍，從而確定合適的聚簇分辨率。
識別和注釋細(xì)胞簇依賴于外部信息部服，這些信息描述了每個細(xì)胞類群預(yù)期的表達(dá)譜稳摄。隨著小鼠腦圖集（Zeisel et al., 2018）、人類細(xì)胞圖集（Regev et al., 2017）和其他工作的繼續(xù)饲宿，可用的參考數(shù)據(jù)庫越來越多厦酬。這些數(shù)據(jù)庫極大地方便了細(xì)胞簇注釋。在沒有相關(guān)參考數(shù)據(jù)庫的情況下瘫想，可以通過將細(xì)胞聚類簇的Marker基因與文獻(xiàn)中的Marker基因進(jìn)行比較來注釋細(xì)胞身份（請參閱Github上的案例研究）或通過直接對文獻(xiàn)報道的標(biāo)記基因在自己研究的數(shù)據(jù)中進(jìn)行表達(dá)值可視化來確定細(xì)胞身份（圖6B）仗阅。應(yīng)當(dāng)注意的是，后一種方法將用戶限制于普通轉(zhuǎn)錄組數(shù)據(jù)鑒定出的經(jīng)典細(xì)胞類型国夜，而不是對細(xì)胞身份 (cell identities)的理解减噪。而且，研究表明常用的細(xì)胞表面標(biāo)志基因在定義細(xì)胞身份上是效果有限的（Tabula Muris Consortium et al., 2018）车吹。
圖6. 聚類分析結(jié)果
(A) UMAP可視化Louvain算法鑒定的細(xì)胞簇筹裕，并標(biāo)記注釋的細(xì)胞類型名稱。(B) 通過細(xì)胞類型Marker基因的表達(dá)標(biāo)記細(xì)胞簇所代表的細(xì)胞類型窄驹，如干細(xì)胞 (Slc12a2), 腸上皮細(xì)胞 (Arg2), 杯狀細(xì)胞 (Tff3) 和潘氏細(xì)胞 (Defa24)朝卒。灰色代表低表達(dá)乐埠，紅色代表高表達(dá)抗斤。Marker基因也可能在其它細(xì)胞簇有表達(dá)，比如Tff3在杯狀細(xì)胞和潘氏細(xì)胞都有表達(dá)丈咐，只是在杯狀細(xì)胞更高瑞眼。(C) 腸上皮細(xì)胞區(qū)域細(xì)胞類型組成熱圖（上面是近端區(qū)域，下面是遠(yuǎn)端區(qū)域）棵逊。相對細(xì)胞密度高的區(qū)域用深紅色表示伤疙。
有兩種方法可以使用參考數(shù)據(jù)庫信息注釋細(xì)胞簇：使用計算出的Marker基因或使用完整的基因表達(dá)譜×居埃可以通過對目標(biāo)簇的細(xì)胞和數(shù)據(jù)集中的所有其他細(xì)胞進(jìn)行差異表達(dá)（DE）分析來鑒定標(biāo)記基因集徒像。通常黍特，我們關(guān)注在目標(biāo)簇中上調(diào)的基因。由于預(yù)期標(biāo)記基因表達(dá)變化幅度較大厨姚，因此通常使用簡單的統(tǒng)計檢驗（例如Wilcoxon秩和檢驗或t檢驗）就可以對兩組細(xì)胞的基因進(jìn)行差異檢驗分析衅澈。根據(jù)統(tǒng)計檢驗結(jié)果，選擇排名最高的基因視為標(biāo)記基因谬墙〗癫迹可以通過比較自己的數(shù)據(jù)中鑒定出的標(biāo)記基因和來自參考數(shù)據(jù)集的標(biāo)記基因?qū)?xì)胞簇進(jìn)行注釋。一些在線工具如www.mousebrain.org (Zeisel et al, 2018) 或http://dropviz.org/ (Saunders et al, 2018)允許用戶在參考數(shù)據(jù)集中可視化用戶鑒定出的Marker基因從而方便細(xì)胞簇的注釋拭抬。
鑒定標(biāo)記基因時應(yīng)注意兩個方面部默。首先，用于鑒定標(biāo)記基因的P值基于以下假設(shè)：即鑒定出的細(xì)胞簇是有生物學(xué)意義的造虎。如果細(xì)胞簇鑒定過程中存在不確定性傅蹂，則必須在統(tǒng)計檢驗中考慮簇分配與標(biāo)記基因鑒定之間的關(guān)系。這個關(guān)系起因于鑒定細(xì)胞簇和標(biāo)記基因都是基于同一套基因表達(dá)數(shù)據(jù)算凿。差異基因檢測的零假設(shè)（null hypothesis）是兩組細(xì)胞整體基因的表達(dá)值具有相同的分布份蝴。然而，由于這兩個聚類組是基于基因表達(dá)變化的聚類結(jié)果得到的氓轰，其基因表達(dá)譜從本質(zhì)上肯定存在差異婚夫。因此即使在應(yīng)用splatter隨機(jī)生成的數(shù)據(jù)的聚類結(jié)果中也可以鑒定出顯著的Marker基因 (Zappia et al, 2017)。為了在聚類數(shù)據(jù)中獲得合理的顯著性度量署鸡，可以使用置換檢驗減少聚類步驟帶來的影響案糙。補充介紹S3中對這一檢驗有詳細(xì)描述。最近的一個差異表達(dá)工具也專門解決了這個問題（Preprint：Zhang et al靴庆，2018）时捌。
在當(dāng)前的研究中，P值通常會被夸大炉抒，這可能導(dǎo)致對標(biāo)記基因數(shù)量的高估奢讨。但是，基于P值對基因進(jìn)行排序則不受影響端礼。假設(shè)聚類結(jié)果是有生物學(xué)意義的禽笑，那么排名最高的標(biāo)記基因仍將是最佳候選標(biāo)記基因。首先蛤奥，我們可以通過可視化展示粗略地查驗獲得的標(biāo)記基因。我們強(qiáng)調(diào)僚稿，特別是通過無監(jiān)督聚類方法定義細(xì)胞聚類簇時凡桥，會導(dǎo)致夸大的P值。當(dāng)改為通過單個基因的表達(dá)確定細(xì)胞簇的身份時蚀同，P值可以解釋為相對于其他基因的期望值缅刽。盡管很常見啊掏，但是這種單變量的聚類簇注釋方法除了在特定情況下不建議使用（例如，β細(xì)胞中的胰島素或紅細(xì)胞中的血紅蛋白）衰猛。其次迟蜜，標(biāo)記基因在數(shù)據(jù)集中將一個細(xì)胞簇與其他細(xì)胞簇區(qū)分開，不僅依賴于細(xì)胞簇鑒定結(jié)果啡省，還依賴于數(shù)據(jù)集的組成。如果數(shù)據(jù)集組成不能準(zhǔn)確表示背景基因表達(dá)，則檢測到的標(biāo)記基因?qū)⑵蛟谄渌?xì)胞中未檢測到的基因茅诱。特別是在細(xì)胞多樣性低的數(shù)據(jù)集中計算標(biāo)記基因時籍嘹，必須特別考慮這一方面。
最近结序，自動細(xì)胞簇注釋已可用障斋。通過直接將注釋好的參考簇的基因表達(dá)譜與單個細(xì)胞進(jìn)行比較，使用諸如scmap（Kiselev et al., 2018b）或Garnett（Preprint：Pliner et al., 2019）之類的工具可以在參考集和自己的數(shù)據(jù)集之間比較注釋徐鹤。因此垃环，這些方法可以同時執(zhí)行細(xì)胞注釋和細(xì)胞簇鑒定，而無需數(shù)據(jù)驅(qū)動的細(xì)胞聚類步驟返敬。但由于不同實驗條件下細(xì)胞類型和狀態(tài)組成不同（Segerstolpe et al.遂庄，2016;Tanay＆Regev，2017）救赐，基于參考數(shù)據(jù)的聚類不應(yīng)取代數(shù)據(jù)驅(qū)動的聚類過程涧团。
聚類、聚類注釋经磅、重新聚類或子聚類以及重新注釋的迭代是很耗時的過程泌绣。自動化的細(xì)胞簇注釋方法大大加快了此過程。但是预厌，自動和手動方法各自存在優(yōu)點和局限性阿迈，速度的提高與靈活性的降低并存，因此很難有一種全局適應(yīng)的方法轧叽。如前所述苗沧，參考數(shù)據(jù)集很難與正在研究的數(shù)據(jù)集包含完全一樣的細(xì)胞類型。因此炭晒，不應(yīng)放棄用于手動注釋標(biāo)記基因的計算待逞。特別是對于包含許多細(xì)胞簇的大型數(shù)據(jù)集，當(dāng)前的最佳實踐是兩種方法（自動+手動）的組合网严。為了提高速度识樱，可以使用自動化細(xì)胞類型注釋工具粗略地標(biāo)記細(xì)胞并確定可能需要鑒定細(xì)胞子群的簇。隨后，應(yīng)針對細(xì)胞簇計算標(biāo)記基因怜庸，并將其與參考數(shù)據(jù)集或文獻(xiàn)中的已知標(biāo)記基因集進(jìn)行比較当犯。對于較小的數(shù)據(jù)集和缺少參考數(shù)據(jù)庫的數(shù)據(jù)集，手動注釋就足夠了割疾。
陷阱和建議：
- 不要使用標(biāo)記基因的P-value值來驗證細(xì)胞簇身份嚎卫，尤其是當(dāng)檢測到的標(biāo)記基因無益于注釋細(xì)胞簇時。P值可能會被夸大宏榕。
- 請注意拓诸，由于數(shù)據(jù)集的細(xì)胞類型和狀態(tài)組成，同一細(xì)胞簇的標(biāo)記基因可能在數(shù)據(jù)集之間有所不同担扑。
- 如果存在相關(guān)的參考數(shù)據(jù)集恰响，我們建議綜合使用自動聚類注釋和基于數(shù)據(jù)的標(biāo)記基因的手動注釋來注釋細(xì)胞類型。
細(xì)胞組成分析 （Compositional analysis）涌献。在細(xì)胞層面胚宦，我們可以根據(jù)其組成結(jié)構(gòu)分析聚類數(shù)據(jù)。細(xì)胞組成數(shù)據(jù)分析圍繞不同樣品落入每個細(xì)胞簇的細(xì)胞比例進(jìn)行分析燕垃。這一比例可能在疾病狀態(tài)下發(fā)生變化枢劝。例如，研究表明沙門氏菌感染會增加小鼠腸上皮區(qū)域腸上皮細(xì)胞（enterocytes）的比例（Haber et al.卜壕，2017）您旁。
要研究單細(xì)胞數(shù)據(jù)的細(xì)胞組成變化，需要足夠的細(xì)胞數(shù)量來穩(wěn)健地評估細(xì)胞簇的比例轴捎，并需要足夠的樣本數(shù)量來評估細(xì)胞簇組成中的背景變化鹤盒。由于合適的數(shù)據(jù)集直到最近才可用，因此專用工具尚待開發(fā)侦副。在上述小鼠研究中侦锯，細(xì)胞身份計數(shù)使用泊松分布建模，使用實驗條件作為協(xié)變量和檢測到的細(xì)胞總數(shù)作為偏移量秦驯。在此尺碰，可以對回歸系數(shù)進(jìn)行統(tǒng)計檢驗，以評估特定細(xì)胞群體的出現(xiàn)頻率是否發(fā)生了顯著變化译隘。但是亲桥，對同一數(shù)據(jù)集中其他細(xì)胞簇的統(tǒng)計檢驗分析并非彼此獨立。如果一個細(xì)胞聚類簇的比例發(fā)生變化固耘，那么所有其他細(xì)胞聚類簇的比例也會發(fā)生變化题篷。因此，使用該模型無法評估整體構(gòu)成是否發(fā)生了顯著變化厅目。在沒有專用工具的情況下悼凑，細(xì)胞構(gòu)成數(shù)據(jù)的可視化展示會給不同樣品之間細(xì)胞組成變化提供有用的信息（圖6C）偿枕。該領(lǐng)域的未來發(fā)展可能會借鑒質(zhì)譜流式細(xì)胞技術(shù)（Mass Cytometry）（例如Tibshirani et al., 2002璧瞬；Arvaniti＆Claassen., 2017户辫；Lun et al., 2017；Weber et al., 2018）或微生物組的方法（Gloor et al., 2017）嗤锉，組成數(shù)據(jù)分析在這些領(lǐng)域已經(jīng)受到了更多關(guān)注渔欢。
陷阱和建議：
- 統(tǒng)計檢驗分析時需要考慮到樣本之間的細(xì)胞簇比例的變化是互相依賴的（非獨立的）
Trajectory analysis
軌跡推斷。細(xì)胞多樣性不能通過離散的分類系統(tǒng)（例如細(xì)胞聚類）充分描述瘟忱。驅(qū)動觀察到的細(xì)胞異質(zhì)性發(fā)展的生物進(jìn)程是一個連續(xù)過程（Tanay＆Regev奥额，2017）。因此访诱，為了捕獲細(xì)胞身份之間的過渡狀態(tài)垫挨、不同的分化分支或生物學(xué)功能的漸進(jìn)式非同步變化，我們需要動態(tài)的基因表達(dá)模型触菜。這類方法稱為軌跡推斷（trajectory inference九榔，TI）。
軌跡推斷方法將單細(xì)胞數(shù)據(jù)視為連續(xù)過程的一個個快照（教你如何定義新亞群 | 在單細(xì)胞水平上解析人肝硬化的纖維化微環(huán)境）涡相。這一過程通過最小化相鄰細(xì)胞之間的轉(zhuǎn)錄改變構(gòu)建細(xì)胞空間的轉(zhuǎn)換路徑（圖7A和B）哲泊。這些路徑上的細(xì)胞排序由偽時間變量 （pseudotime variable）描述。雖然此變量是基于距離根細(xì)胞的轉(zhuǎn)錄距離計算的催蝗，但通常被解釋為發(fā)育時間的代名詞（Moignard et al.切威，2015; Haghverdi et al.，2016; Fischer et al.丙号，2018; Griffiths et al.先朦，2018）。
圖7. 小鼠腸上皮細(xì)胞軌跡分析結(jié)果展示
(A) Slingshot鑒定出的遠(yuǎn)端和近端腸上皮細(xì)胞分化軌跡犬缨。遠(yuǎn)端譜系細(xì)胞按pseudotime從紅色到藍(lán)色上色喳魏。數(shù)據(jù)集中其它細(xì)胞用灰色點展示。(B) PCA空間中Slingshot 軌跡和細(xì)胞簇比較展示遍尺。簇名字簡寫如下：cuEP—enterocyte progenitors; Imm. Ent.—immature enterocytes; Mat. Ent.—mature enterocytes; Prox.—proximal; Dist.—distal. (C) 遠(yuǎn)端腸上皮細(xì)胞軌跡中每個pseudotime對應(yīng)的細(xì)胞密度分布截酷。柱子按每個pseudotime 區(qū)間主要的細(xì)胞簇上色。(D) 數(shù)據(jù)集投射到UMAP空間的抽象圖展示乾戏。每個簇用不同顏色點展示迂苛。出現(xiàn)在其它軌跡中的簇特別標(biāo)記處理用于比較分析。TA代表短暫增殖下?lián)?(transit amplifying cells)鼓择。(E)R包GAM展示腸上皮細(xì)胞中基因表達(dá)隨pseudotime的動態(tài)變化三幻。
自Monocle（Trapnell et al. 2014）和Wanderlust（Bendall et al. 2014）建立了TI （trajectory inference）領(lǐng)域以來，可用的TI方法數(shù)量激增呐能。當(dāng)前可用的TI方法的差別在于構(gòu)建的發(fā)育軌跡模型拓?fù)浣Y(jié)構(gòu)復(fù)雜性不同念搬，從簡單的線性軌跡或二分支軌跡到復(fù)雜樹形軌跡抑堡、多分支軌跡或組合多種拓?fù)浣Y(jié)構(gòu)軌跡。在最近對TI方法進(jìn)行的全面比較中（Saelens et al.朗徊，2018）發(fā)現(xiàn)沒有一種單獨的方法可以在所有類型的軌跡分析中都表現(xiàn)最優(yōu) （NBT|45種單細(xì)胞軌跡推斷方法比較首妖，110個實際數(shù)據(jù)集和229個合成數(shù)據(jù)集）。相反爷恳，應(yīng)根據(jù)預(yù)期軌跡的復(fù)雜性選擇TI（軌跡分析）方法有缆，研究比較表明Slingshot（Street et al.，2018）在簡單軌跡分析如從線性軌跡到二分支和多分軌跡表現(xiàn)最佳温亲。如果預(yù)期數(shù)據(jù)對應(yīng)更復(fù)雜的軌跡棚壁，作者建議使用PAGA（Wolf et al.，2019）栈虚。如果知道精確的軌跡模型袖外，則可以選擇使用更特定的方法來提高性能（Saelenset al.，2018）魂务。通常曼验，應(yīng)使用多種方法來確定評估推斷出的軌跡，以避免方法偏差头镊。
在典型的分析流程中蚣驼，軌跡推斷（TI）方法應(yīng)用于降維后的數(shù)據(jù)。如果使用的TI工具自帶了降維功能相艇，則基于校正后的數(shù)據(jù)進(jìn)行分析颖杏。由于通常在細(xì)胞內(nèi)同時發(fā)生多種生物學(xué)過程，因此消除其他生物過程的影響對鑒定預(yù)期軌跡可能很有用坛芽。例如留储，T細(xì)胞在成熟過程中可能會經(jīng)歷細(xì)胞周期轉(zhuǎn)換（Buettner et al.，2015）咙轩。此外获讳，由于幾種性能最好的TI方法依賴于聚類后的數(shù)據(jù)，因此TI通常在聚類之后執(zhí)行活喊。軌跡中的細(xì)胞簇可能表示穩(wěn)態(tài)或亞穩(wěn)態(tài)細(xì)胞（請參見“亞穩(wěn)態(tài)”丐膝；圖7B和C）。隨后钾菊，可以將RNA velocities （RNA速度帅矗，或RNA表達(dá)動力學(xué)）疊加到軌跡上確認(rèn)發(fā)育方向（La Manno et al.，2018）煞烫。（生信寶典注：新生轉(zhuǎn)錄本成熟過程中需要進(jìn)行剪接操作浑此。對于一個穩(wěn)定表達(dá)的基因，總會在細(xì)胞中找到存在一定比例的未剪接的非成熟RNA形式滞详，用于補充老的轉(zhuǎn)錄本的降解凛俱。如果一個基因剛被激活紊馏，短時間內(nèi)將會有高比例的未成熟轉(zhuǎn)錄本。相反蒲犬，當(dāng)一個基因被抑制時朱监，轉(zhuǎn)錄過程會早于轉(zhuǎn)錄本降解過程而被抑制，未成熟轉(zhuǎn)錄本的比例會降低暖哨。因此對于細(xì)胞中每個基因赌朋，未剪接的mRNA相對于剪接的mRNA的比例（RNA velocity）可以推斷瞬時表達(dá)動力學(xué)，進(jìn)一步推演組織內(nèi)發(fā)生的細(xì)胞轉(zhuǎn)變篇裁。https://www.nature.com/articles/d41586-018-05882-8）
推斷的軌跡不一定要完全對應(yīng)生物發(fā)育過程。首先赡若，推斷的軌跡僅表示轉(zhuǎn)錄相似性达布。很少有TI方法在其模型中包括不確定性評估（Griffiths et al., 2018）。因此逾冬，需要更多的信息來驗證是否確實捕獲了生物過程黍聂。這些信息可以來源于干擾實驗、推斷的調(diào)控基因動力學(xué)以及RNA velocity數(shù)據(jù)的支持等身腻。
陷阱和建議：
- 我們建議使用
Saelens et al.（2018）的綜述（NBT|45種單細(xì)胞軌跡推斷方法比較产还，110個實際數(shù)據(jù)集和229個合成數(shù)據(jù)集）作為指南 。
- 推斷的軌跡不需要完全對應(yīng)生物過程嘀趟， 應(yīng)該收集更多的證據(jù)來解釋軌跡脐区。
基因表達(dá)動力學(xué)。推斷的軌跡也可能是轉(zhuǎn)錄噪聲擬合的結(jié)果她按，一個排除方式是在基因水平上對鑒定出的軌跡進(jìn)行進(jìn)一步分析牛隅。在整個偽時間范圍內(nèi)連續(xù)平穩(wěn)變化的基因是決定軌跡的關(guān)鍵基因，可用于識別潛在的生物學(xué)過程酌泰。此外媒佣，這一組與軌跡相關(guān)的基因應(yīng)該包含調(diào)控相應(yīng)生物進(jìn)程的基因。調(diào)節(jié)基因可以幫助我們理解該生物過程是如何和為什么被觸發(fā)的陵刹，并且可能是潛在的藥物靶標(biāo)（Gashaw et al., 2012）默伍。
早期尋找軌跡相關(guān)基因的方法是沿軌跡在細(xì)胞簇之間進(jìn)行差異基因分析（DE）。現(xiàn)在通過將基因表達(dá)與偽時間 (pusedotime)信息進(jìn)行回歸來鑒定隨軌跡顯著變化的基因衰琐。為了保證表達(dá)值隨偽時間變量平滑變化也糊，首先通過擬合樣條曲線 (fit a spline)或其他局部回歸方法（例如loess）對偽時間進(jìn)行平滑處理〉舛回歸框架的區(qū)別有兩點：一是其假設(shè)的噪聲模型不同显设，二是構(gòu)建的表達(dá)-偽時間函數(shù)的類別不同。通過基因?qū)螘r間的依賴性進(jìn)行模型篩選獲得潛在的調(diào)控基因辛辨。這種基于偽時間回歸進(jìn)行的DE基因測試也會受到軌跡推斷方法的影響捕捂，就像不同細(xì)胞簇之間的差異基因分析會受到聚類方法的影響一樣（參見“細(xì)胞簇注釋”部分）瑟枫。因此，在此步驟中獲得的P值不應(yīng)視為顯著性的評估指攒。
當(dāng)前幾乎沒有專用于分析基因表達(dá)時間動力學(xué)的分析工具慷妙。BEAM是整合于Monocle 軌跡推斷方法中的一個工具，可以分析分支特異的基因表達(dá)動力學(xué)允悦。除了這個膝擂，用戶可以選擇使用LineagePulse (https://github.com/YosefLab/LineagePulse)，這個工具考慮了dropout噪音隙弛，但仍然在開發(fā)中架馋。或者基于標(biāo)準(zhǔn)R包或Limma包寫作自己的分析測試工具全闷。在Slingshot`的使用幫助中可以找到一個這樣應(yīng)用的例子供參考叉寂。
由于可用的工具很少，因此尚不能確定研究基因時間表達(dá)動力學(xué)的最佳方法总珠。使用上述所有方法肯定可以對基因表達(dá)動力學(xué)進(jìn)行探索性研究屏鳍。將來，高斯過程 (Gaussian processes)可能會提供一個natural model來研究基因的時間動力學(xué)局服。此外钓瞭，對調(diào)節(jié)模塊整體而不是單個基因進(jìn)行統(tǒng)計檢驗分析可能會提高信噪比并更方便生物學(xué)解釋。
亞穩(wěn)態(tài)淫奔。軌跡的細(xì)胞水平分析是指研究偽時間軸上的細(xì)胞密度山涡。假設(shè)以無偏方式對細(xì)胞進(jìn)行隨機(jī)采樣，軌跡中細(xì)胞密集的區(qū)域表示優(yōu)選的轉(zhuǎn)錄狀態(tài)搏讶。當(dāng)將軌跡解釋為時間過程時佳鳖，這些密集區(qū)域可能代表例如發(fā)育中的亞穩(wěn)態(tài)（Haghverdi et al., 2016）。我們可以通過繪制偽時間坐標(biāo)的直方圖來找到這些亞穩(wěn)態(tài)（圖7C）媒惕。
細(xì)胞水平聯(lián)合分析系吩。聚類和軌跡推斷代表單細(xì)胞數(shù)據(jù)的兩個不同視角，并可以在粗粒度(coarse‐grained)的圖形展示中進(jìn)行統(tǒng)一妒蔚。用點代表單細(xì)胞簇穿挨，軌跡作為簇之間相連的邊，可以同時展示數(shù)據(jù)的靜態(tài)和動態(tài)屬性肴盏。這一聯(lián)合分析在PAGA工具中首先提出科盛。PAGA提出一個統(tǒng)計模型進(jìn)行細(xì)胞簇互作分析，將細(xì)胞相似度高于期望值的兩個簇之間用線連接菜皂。在最近的綜述中 (Saelens et al, 2018)贞绵，PAGA相比于其它軌跡分析方法整體表現(xiàn)更優(yōu)。PAGA是唯一一個討論到的可以處理不相連的拓?fù)浣Y(jié)構(gòu)和包含環(huán)形結(jié)構(gòu)的復(fù)雜軌跡圖的一個工具恍飘。這一特征使得PAGA在可視化整個數(shù)據(jù)集并進(jìn)行探索分析時很有用榨崩。
基因水平分析谴垫。到目前為止，雖然我們專注于表征細(xì)胞結(jié)構(gòu)的基因水平分析方法母蛛，但是單細(xì)胞數(shù)據(jù)的基因水平分析有更廣泛的內(nèi)容翩剪。差異表達(dá)分析、基因集分析和基因調(diào)控網(wǎng)絡(luò)分析都可以直接研究數(shù)據(jù)中的分子信號彩郊。這些方法不是描述細(xì)胞異質(zhì)性前弯，而是把這種異質(zhì)性作為理解基因表達(dá)差異的背景。
差異表達(dá)分析秫逝。表達(dá)數(shù)據(jù)分析的一個常見問題是兩個實驗條件之間是否有基因發(fā)生了差異表達(dá)恕出。差異基因分析源于大量細(xì)胞基因表達(dá)分析（Scholtens＆von Heydebreck，2005）筷登，已經(jīng)有過很好的描述剃根。相對于bulk差異分析，單細(xì)胞的優(yōu)勢是可以基于細(xì)胞簇進(jìn)行細(xì)胞異質(zhì)性評估前方。分析某個類型細(xì)胞簇在不同實驗條件下的轉(zhuǎn)錄響應(yīng)（Kang et al.，2018）廉油。雖然用于解決同樣的問題惠险，普通和單細(xì)胞差異基因分析工具在方法上差別很大。普通轉(zhuǎn)錄組的差異基因分析方法解決的是少量樣品中精確評估基因表達(dá)變化的問題抒线，而單細(xì)胞數(shù)據(jù)中不存在這一問題班巩。另一方面，單細(xì)胞數(shù)據(jù)包含特有的技術(shù)噪音如dropout和高細(xì)胞變異性 (Hicks et al, 2017; Vallejos et al, 2017)嘶炭。這些影響因素在設(shè)計單細(xì)胞專用的差異分析方法時都要有考慮抱慌。然而，在最近一個大規(guī)模的差異基因分析方法評估中發(fā)現(xiàn)普通轉(zhuǎn)錄組差異基因檢測工具與表現(xiàn)最好的單細(xì)胞差異分析工具性能相當(dāng) (Soneson & Robinson, 2018)眨猎。進(jìn)一步地抑进，當(dāng)在普通轉(zhuǎn)錄組的差異分析工具中引入基因權(quán)重對單細(xì)胞數(shù)據(jù)建模時，這些工具的表現(xiàn)優(yōu)于同類的單細(xì)胞差異分析工具 (Van den Berge et al, 2018)睡陪∷律基于這個比較分析，DESeq2 (Love et al, 2014) 和EdgeR (Robinson et al, 2010) 加ZINB‐wave (Risso et al, 2018)估計出的權(quán)重信息表現(xiàn)最好兰迫。不過仍需加權(quán)的普通轉(zhuǎn)錄組差異分析的獨立比較研究進(jìn)一步確認(rèn)這個結(jié)果信殊。
轉(zhuǎn)錄組差異分析工具的加權(quán)建模后效果的改善是以犧牲計算效率為代價的。隨著單細(xì)胞數(shù)據(jù)中細(xì)胞量越來越大汁果，算法的運行時間在方法選擇中越來越重要涡拘。因此單細(xì)胞工具MAST被視作加權(quán)普通轉(zhuǎn)錄組差異分析工具的替代方法。MAST應(yīng)用hurdle model消除dropout的影響据德，并構(gòu)建基因表達(dá)相對于實驗條件和技術(shù)協(xié)變量的變化模型鳄乏。這是在前述評估研究中表現(xiàn)最好的單細(xì)胞差異基因分析方法 (Soneson & Robinson, 2018)跷车，并且在一個單一數(shù)據(jù)集的小范圍研究中表現(xiàn)優(yōu)于其它普通或單細(xì)胞差異基因分析方法 (Vieth et al, 2017)。MAST相比于加權(quán)的普通轉(zhuǎn)錄組差異分析工具有10-100倍的速度提升汞窗，還可以使用limma–voom模型進(jìn)一步獲得10倍的速度提升姓赤。雖然limma是普通轉(zhuǎn)錄組基因差異分析方法，但limma–voom與MAST性能相當(dāng)仲吏。
因為用于差異基因分析的是標(biāo)準(zhǔn)化后但未進(jìn)行技術(shù)變量和非相關(guān)生物變量校正的數(shù)據(jù)不铆，因此在分析過程中考慮這些混雜變量對獲得穩(wěn)定的結(jié)果是關(guān)鍵的。通常差異基因分析工具都允許用戶考慮混雜變量（confounders）裹唆，但是把哪些變量納入模型需要謹(jǐn)慎評估誓斥。比如，在多數(shù)單細(xì)胞數(shù)據(jù)集中许帐，樣本和實驗處理條件是互斥的劳坑，因為幾乎不可能對一個樣品進(jìn)行多個實驗處理。如果我們在模型中同時考慮樣品和實驗條件兩個協(xié)變量成畦，跟這些協(xié)變量相關(guān)的變化則不會被明確確定距芬。因此在考慮不同實驗條件下的差異分析時，通常不在模型中引入樣品信息作為協(xié)變量循帐。當(dāng)校正多個確定的批次分類變量時框仔，可視化展示協(xié)變量的混雜影響將會變得更困難。在這個情況下拄养，需要先確認(rèn)模型設(shè)計矩陣是否滿秩 （the model design matrix is full rank离斩，full rank是指模型設(shè)計矩陣中的每一列都不能通過對其它列的線性組合得來，即它們之間不存在線性相關(guān)）瘪匿。即便一些設(shè)計矩陣不是滿秩跛梗，差異基因分析工具通常也不會給出警告而是繼續(xù)運行。這時獲得的結(jié)果將可能不是預(yù)期的分析方向棋弥。
我們這兒描述的場景中核偿，實驗條件協(xié)變量是在實驗設(shè)計中決定的。因此在同一簇內(nèi)基于這一協(xié)變量的差異基因分析是獨立于聚類過程的嘁锯。這將基于實驗條件的差異基因分析和基于細(xì)胞簇的差異基因分析區(qū)別開來宪祥。基于實驗條件的差異基因檢測獲得的P-values是顯著性檢驗統(tǒng)計指標(biāo)家乘，需要進(jìn)一步做多重檢驗校正蝗羊。為了降低多重檢驗校正的壓力，不感興趣的基因通常需要預(yù)先排除掉仁锯。雖然假基因和非編碼基因也可以提供有用信息耀找，但在分析中也通常會被忽略掉。
陷阱和建議：
- DE測試不應(yīng)基于校正后的數(shù)據(jù)（去噪，批次校正等）野芒，而應(yīng)基于原始測量數(shù)據(jù)并在模型中引入干擾協(xié)變量蓄愁。
- 用戶不應(yīng)依賴差異基因檢測工具校正帶有混雜協(xié)變量的模型。設(shè)計矩陣模型應(yīng)當(dāng)保證滿秩狞悲。
- 我們建議使用
MAST或limma進(jìn)行差異基因分析撮抓。
基因集分析。基因水平的分析通常會獲得一長串難以解釋其生物含義的基因列表荸恕。比如乖酬，處理組和對照樣品通常會有數(shù)以千計差異表達(dá)基因。我們可以通過基因共有的特征對基因進(jìn)行分組并檢驗這些特征在候補基因列表中的出現(xiàn)是否顯著富集 (overrepresented)以便更好的解釋結(jié)果融求。
不同用途的基因集信息可以在校正過的基因注釋數(shù)據(jù)庫中獲取咬像。為了解釋差異基因結(jié)果，通常按共有的生物進(jìn)程對基因進(jìn)行分組生宛。生物進(jìn)程信息存儲于MSigDB县昂、Gene Ontology，通路數(shù)據(jù)存儲于KEGG和Reactome數(shù)據(jù)庫陷舅。基因列表的功能注釋富集分析有多種工具可用七芭，具體見Huang et al (2009)和 Tarca et al (2013)的綜述。
分析領(lǐng)域富集分析的一個新操作是基于配對基因標(biāo)記的配體-受體分析。通過受體和共軛配體之間的表達(dá)關(guān)系推斷細(xì)胞簇之間的互作。CellPhoneDB提供了受體-配體對信息，并用統(tǒng)計模型解釋簇之間的高表達(dá)基因是否可以配對 (Zepp et al, 2017; Zhou et al, 2017; Cohen et al, 2018; Vento‐Tormo et al, 2018)戴尸。
基因調(diào)控網(wǎng)絡(luò)拂盯。基因不是獨立發(fā)揮作用的。相反撒璧，基因的表達(dá)水平是由與其他基因和小分子之間的復(fù)雜調(diào)控決定的。揭示這些調(diào)控作用是基因調(diào)控網(wǎng)絡(luò)（GRN）推斷方法的目標(biāo)（SCENIC | 從單細(xì)胞數(shù)據(jù)推斷基因調(diào)控網(wǎng)絡(luò)和細(xì)胞類型）。
基因調(diào)控網(wǎng)絡(luò)推斷是基于對基因共表達(dá)如相關(guān)性砚尽、互信息（mutual information）或回歸模型分析進(jìn)行測量（Chen＆Mar，2018）辉词。如果在考慮了其他所有基因都作為混雜因素后必孤，兩個目標(biāo)基因之間依然存在共表達(dá)關(guān)系則可以推斷這兩個基因之間存在因果調(diào)控關(guān)系∪鹛桑基因調(diào)控關(guān)系鑒定與軌跡相關(guān)的調(diào)控基因檢測相關(guān)敷搪。確實兴想，幾種單細(xì)胞GRN推論方法使用的機(jī)械微分方程模型中都引入了軌跡信息（Inferring gene regulatory relationships is related to the detection of trajectory‐associated regulatory genes. Indeed, several single‐cell GRN inference methods use trajectories with mechanistic differential equation models）（Ocone et al., 2015; Matsumoto et al., 2017）。
雖然存在專門針對scRNA-seq數(shù)據(jù)開發(fā)的GRN推斷方法（SCONE：Matsumoto et al., 2017; PIDC：Chan et al., 2017; SCENIC：Aibar et al., 2017）赡勘，但最近的比較分析表明普通和單細(xì)胞轉(zhuǎn)錄組的方法在這些數(shù)據(jù)表現(xiàn)都不好（Chen＆Mar嫂便，2018）。GRN推斷方法仍可提供識別生物過程的因果調(diào)節(jié)子的有用信息闸与，但我們建議謹(jǐn)慎使用這些方法毙替。
陷阱和建議：
- 用戶應(yīng)警惕推斷的調(diào)控關(guān)系中的不確定性〖钙基因調(diào)控關(guān)系豐富的基因模塊比調(diào)控信息稀疏的模塊更可靠蔚龙。
分析平臺。單細(xì)胞分析流程是多個獨立開發(fā)的工具的集合映胁。為了促進(jìn)這些工具之間的數(shù)據(jù)傳遞木羹，單細(xì)胞工具開發(fā)時都考慮使用一致的數(shù)據(jù)格式。這些工具為構(gòu)建分析流程提供了基礎(chǔ)解孙。目前可用的分析平臺有基于R（McCarthy et al., 2017; Butler et al., 2018）或Python（Wolf et al., 2018）的命令行坑填，或具有圖形用戶界面（GUI）的本地應(yīng)用程序（Patel，2018; Scholz et al., 2018）或Web服務(wù)（Gardeux et al., 2017; Zhu et al., 2017）弛姜。Zhu et al.（2017）和Zappia et al.（2018）對這些分析平臺進(jìn)行了綜述脐瑰。
在命令行平臺中，Scater（McCarthy et al.廷臼，2017）和Seurat（Butler et al.苍在，2018）可以輕松地與R Bioconductor項目提供的各種分析工具進(jìn)行交互（Huber et al.，2015）（讓你的單細(xì)胞數(shù)據(jù)動起來荠商！|iCellR(二)）寂恬。Scater在質(zhì)量控制和數(shù)據(jù)預(yù)處理方面具有特殊優(yōu)勢，而Seurat可以說是最受歡迎和最全面的分析平臺莱没，提供了大量分析模塊和教程初肉。scanpy是一個基于Python的不斷發(fā)展的新的單細(xì)胞分析平臺，該平臺可以擴(kuò)展應(yīng)用于分析更大規(guī)模的單細(xì)胞數(shù)據(jù)集饰躲。它可以整合Python寫作的很多工具牙咏，尤其是在機(jī)器學(xué)習(xí)中流行的工具。
圖形用戶界面工具可以幫助非專業(yè)用戶構(gòu)建單細(xì)胞分析工作流嘹裂。通常會基于固定的分析流程引導(dǎo)用戶妄壶，這些分析流程有助于簡化分析，但同時也限制了用戶的靈活性焦蘑，在探索性分析時特別有用盯拱。諸如Granatum（Zhu et al.，2017）和ASAP（Gardeux et al.，2017）之類的平臺在集成的工具方面有所不同狡逢，其中Granatum包含了更多的方法宁舰。作為Web服務(wù)器，這兩個平臺都直接可用奢浑，但是其不同的計算框架將限制它們擴(kuò)展到大規(guī)模數(shù)據(jù)集的能力蛮艰。ASAP的測試數(shù)據(jù)集只包含92個細(xì)胞∪副耍基于Web的GUI平臺的替代產(chǎn)品是FASTGenomics（Scholz et al., 2018）壤蚜、iSEE（Rue-Albrecht et al., 2018）、IS-CellR（Patel徊哑，2018）和Granatum等軟件包袜刷，它們在本地服務(wù)器上運行。這些平臺和GUI封裝工具可以隨著本地電腦計算力的增強(qiáng)進(jìn)行性能擴(kuò)展莺丑。將來著蟹，Human Cell Atlas（https://www.humancellatlas.org/data-sharing）的不斷發(fā)展將促進(jìn)功能更強(qiáng)大的可應(yīng)用于大規(guī)模細(xì)胞數(shù)據(jù)的可視化探索工具的開發(fā)。
結(jié)論與展望
我們綜述了典型的scRNA‐seq分析流程的各個步驟梢莽，并在案例研究教程（https://www.github.com/theislab/single-cell-tutorial）中實現(xiàn)了這些步驟萧豆。本教程旨在通過對當(dāng)前可用分析方法的比較研究確定當(dāng)前最佳分析實踐流程。雖然匯總單個最佳實踐工具不能保證流程最佳昏名，但我們希望我們的分析流程能夠代表單細(xì)胞分析領(lǐng)域最新技術(shù)現(xiàn)狀的一個掠影涮雷。它為新來者提供了進(jìn)入該領(lǐng)域的合適切入點，并希望能有助于人類細(xì)胞圖譜的scRNA-seq分析流程構(gòu)建（Regev et al.轻局，2018）洪鸭。應(yīng)該注意的是，方法比較必然落后于最新的方法開發(fā)仑扑。因此卿嘲，我們提到了可能尚未被獨立評估的新工具。隨著新工具和更好工具的進(jìn)一步發(fā)展和比較研究夫壁，文中推薦的每個單獨的工具也需要一直更新，但是有關(guān)數(shù)據(jù)處理階段的一般注意事項應(yīng)保持不變 （生信寶典注：工具一直在變沃疮，但核心原理不變盒让。來易生信，掌握核心知識司蔬。）邑茄。
特別受關(guān)注并且可能改變現(xiàn)有分析流程的兩個開發(fā)方向是深度學(xué)習(xí)分析和單細(xì)胞多組學(xué)整合分析。由于深度學(xué)習(xí)可以很容易拓展到大規(guī)模數(shù)據(jù)分析俊啼，已經(jīng)廣泛用于計算機(jī)視覺到自然語言處理等多個領(lǐng)域肺缕，而且在基因組領(lǐng)域應(yīng)用也越來越多。機(jī)器學(xué)習(xí)首先應(yīng)用于scRNA-seq的降維和降噪分析領(lǐng)域，例如scVis：Ding et al.同木，2018; scGen：Preprint：Lotfollahi et al.浮梢，2018; DCA：Eraslan et al.，2019）彤路。最近秕硝，深度學(xué)習(xí)已形成完善的分析流程用于擬合數(shù)據(jù)、對數(shù)據(jù)進(jìn)行去噪并在模型框架內(nèi)執(zhí)行下游分析如聚類和差異表達(dá)等（scVI：Lopez et al., 2018）洲尊。在這一分析方法中远豺，可以在保留數(shù)據(jù)變化的準(zhǔn)確估計的同時，把噪聲和批次效應(yīng)的影響納入下游統(tǒng)計分析模型中坞嘀。諸如此類的整合建模方法有可能取代當(dāng)前雜糅了多個工具的分析流程躯护。
隨著單細(xì)胞組學(xué)技術(shù)的進(jìn)步，對組學(xué)分析流程開發(fā)的需求將會增長（Tanay＆Regev丽涩，2017）棺滞。未來的單細(xì)胞分析平臺將必須能夠處理不同的數(shù)據(jù)源，例如DNA甲基化（Smallwood et al., 2014）内狸，染色質(zhì)可及性（Buenrostro et al., 2015）或蛋白質(zhì)豐度數(shù)據(jù)（Stoeckius et al., 2017）检眯，以及能對這些數(shù)據(jù)進(jìn)行整合分析。如果是這個目的昆淡，將不能只使用我們流程起點的單個reads或計數(shù)矩陣作為輸入锰瘸。而且，整合RNA velocity分析的平臺已經(jīng)在適應(yīng)多組學(xué)（多模態(tài)）數(shù)據(jù)結(jié)構(gòu)昂灵，輸入的是未剪接和剪接的reads數(shù)據(jù)（La Manno et al避凝，2018）。單細(xì)胞多組學(xué)整合可以通過一致性聚類方法眨补、多組學(xué)因子分析（Argelaguet et al., 2018）或多組學(xué)基因調(diào)控網(wǎng)絡(luò)推斷（Colomé-Tatché＆Theis管削，2018）來實現(xiàn)。具有這些功能的分析流程將是未來開發(fā)的方向撑螺。我們預(yù)計這種多組學(xué)分析流程將建立在我們?yōu)閟cRNA-seq分析奠定的基礎(chǔ)上含思。
文章解讀：苑曉梅
文章校對：生信寶典
推薦閱讀
參考文獻(xiàn)：
Luecken MD, Theis FJ. Current best practices in single-cell RNA-seq analysis:a tutorial. Mol Syst Biol. 2019 Jun 19;15(6):e8746\. doi: 10.15252
