文献分享：single-cell rna sequencing analysis:astep-by-step overview

今天我们分享篇文献，如上图。

摘要

        由于创新的样品制备和测序技术，单个细胞中的基因表达现在可以在一次实验中测量数千个细胞。自引入以来，单细胞RNA测序(scRNA-seq)方法已经彻底改变了基因组学领域，因为它们通过在单细胞分辨率下探索基因表达谱，为解决细胞异质性创造了前所未有的机会。然而，快速发展的scRNA-seq领域引发了各种分析方法的出现，旨在最大限度地发挥这种新策略的全部潜力。与基于群体的RNA测序方法不同，scRNA测序需要全面的计算工具来解决高数据复杂性，并跟上新兴的单细胞相关挑战。尽管有大量的分析方法，但缺乏普遍的标准化。虽然这反映了该领域的不成熟，但也可能阻碍新来者融入其中。

        在这篇综述中，我们的目标是跨越上述障碍，并提出了四种现成的用于scRNA-seq分析的管道，这些管道很容易被新手获取，可以适应各种生物数据类型。在这里，我们概述了目前可用的单细胞技术，用于细胞分离和文库制备，并提供了一步一步的指南，涵盖了分析scRNA-seq数据的整个规范分析工作流程，包括读取映射、质量控制、基因表达量化、归一化、特征选择、降维和对轨迹推断和差异表达有用的细胞聚类。这样的工作流程指南将为新手和专家用户在分析复杂的scRNA-seq数据集时提供指导，从而进一步扩大单细胞方法在基础科学中的研究潜力，并设想其未来的实施作为该领域的最佳实践。

        关键词:单细胞RNA-seq，实验流程，数据分析教程，计算管道，聚类，Monocle, Seurat, gf-icf, Scanpy

介绍

        在过去的十年中，基于群体的RNA测序方法(又名bulk RNA-seq)在破译包括癌症生物学、发育生物学和细胞稳态在内的广泛领域的全基因组转录组变异方面发挥了重要作用。然而，由于大量RNA-seq数据代表了单个细胞中基因表达的平均值，它可能掩盖了具有最丰富细胞类型或状态的不同亚群的转录趋势(辛普森悖论)。单细胞RNA测序(scRNA-seq)跨越了这一障碍，为在单细胞分辨率下探索基因表达谱提供了前所未有的机会。自2009年首次引入以来，scRNA-seq为揭示复合系统潜在的细胞异质性开辟了新的途径。然而，实际的程序是艰巨的、耗时的、成本高的，并且严重依赖于单一来源的设备。目前，随着高效和低成本技术的出现(表1)，一个典型的实验室工作台足以建立多达数千个细胞的测序文库，因此鼓励使用单细胞技术作为标准程序。

        这些技术进步使得新的细胞类型的发现和细胞动态过程的研究在以前无法实现的空间和时间分辨率，细胞内的变异，如基因相互作用，等位基因表达，和新的RNA加工在分子细胞生物学领域。此外，scRNA-seq成为快速发展的精准医学领域的关键因素。通过scRNA-seq获得的大量新信息有可能重塑我们对发育生物学、基因调控以及健康和疾病中的细胞异质性的理解。        

scRNA-seq的实验室工作流程

        目前，所有scRNA-seq实验室方法都依赖于六个主要步骤(图1):(I)制备有活力的单细胞悬液，(II)评估细胞活力，(III)裂解细胞去除，(IV)个体转录组条形码，(V) cDNA生成，(VI)测序文库生成。在仪器实施方面，最受欢迎的测序平台之一是Illumina®系列，因为它具有成本效益和高质量的输出。该领域的一个相对较新的引入是华大基因测序组合，它甚至可以在单细胞研究中获得等电位测序结果。

        

        图1单细胞RNA测序工作流程。scRNA-seq程序包括六个关键步骤。(1)将样品解离成单细胞悬浮液。(II)由于裂解的细胞可能会对数据产生偏差并造成高噪声干扰，因此必须最大限度地提高输入材料的质量并评估细胞活力。(三)如果活活力低于90%，应采用离心(即密度梯度)或免疫消耗(即FACS或磁分选)过滤死细胞。(IV)根据所选择的技术，以不同的方式捕获和分离单细胞。基于微流体的scRNA-seq技术将单个细胞封装在油包水液滴中，并将独特的引物附着在微粒表面和裂解缓冲液上。然后，每个裂解细胞的mRNA含量被单个引物的poly-A尾部结构域捕获，并用UMI和细胞特异性条形码标记。在这一步骤中可能会发生一些错误，例如在单个液滴中捕获多个细胞或微粒(即多组)，以及亚泊松加载权衡，例如空条形码液滴。(V)然后收集从液滴中捕获的mRNA转录本，进行反转录，(VI)在池中扩增，用于标准测序平台。在文库构建过程中，cDNA分子被标记为样品特异性索引，允许在相同的测序运行中复用不同的捕获。进一步的计算解复用将使用这样的条形码信息来对样本、细胞和转录本进行排序。      

        对单个细胞的转录组进行条形码是所有可用的单细胞方案的关键步骤，并举例说明了与批量rna测序的主要区别。本文提出了两种条形码策略，要么(1)在细胞分离后为每个转录组添加细胞特异性条形码，要么(2)为每个细胞转录组添加不需要物理划分的唯一索引组合(例如split-Seq)。这两种策略都可以进一步划分为具有不同优点和缺点的子类别(表1，补充)。然而，所有scRNA-seq策略都依赖于高质量的输入材料，需要优化任何解离和解冻方案，以最大限度地提高细胞活力。

        在该领域的最新进展中，基于微流体的scRNA-seq技术因其成本效益高、效率高、对保持数据完整性和一致性的适度数据大小要求而受到欢迎。一般来说，微流控技术，如Chromium、inDrop和Drop-seq，依靠细胞、微粒(即微珠)和裂解缓冲液的被动共流，产生油包水滴，从而精确地封装一个细胞和一个微珠。每个细胞的转录内容被捕获并通过附着在单个微粒表面的独特引物扩增。

        这些引物具有相同的基本三元组结构，包括(1)细胞条形码，即单个微粒上所有引物共有的短序列，目的是识别属于同一细胞的所有转录本;(2)唯一分子标识符(UMI)，这是一种分子转录特异性标签，通过识别PCR重复来确保read的完整性[29];(iii) poly-T尾部，用于捕获和扩增每个转录本的30个END。

        理想情况下，每个液滴应该封装一个细胞和一个珠。然而，由于在实践中封装步骤遵循泊松分布，单个液滴内一个头和一个细胞的捕获率遵循双泊松分布。因此，许多基于液滴的方法产生了大量的空液滴和低效的数据组装。

        泊松统计的局限性已经被inDrop和Chromium技术解决了。通过对可变形颗粒的紧密排列，两种方法都测量了亚泊松分布，从而实现了可控制的封装颗粒数量，单粒占用率约为80%。因此，inDrop、Chromium和dopo -seq这三个平台之间的主要区别在于它们各自的捕获效率，这在很大程度上取决于使用的珠粒类型。虽然Drop-seq、inDrops和Chromium捕获了大约5-12%、75%和65%的输入细胞，但它们也分别需要>2 105个、2 103 -104个和>103个输入细胞。

        因此，选择合适的技术是至关重要的，这取决于特定的研究领域和研究要求。当研究高度异质性的样品时，如肿瘤和组织，高通量方法是可取的。然而，在分析低表达基因或对稀有细胞群体进行分类时，高灵敏度策略是最适合的。然而，如果一方面scRNA-seq允许在高分辨率上解剖细胞异质性，它也有两个关键的缺点。首先是低基因检索率，通常每个细胞1-5%的转录本代表高表达基因(每个细胞约5000个基因)，从而导致显著的观察不确定性。这种退出效应导致了高细胞间变异性和低信噪比(SNR)。另一个缺点是scRNA-seq商业技术的成本负担，而非商业平台(inDrop, Drop-seq)需要相当多的操作人员专业知识。因此，该领域的许多实验室仍然无法广泛使用scRNA-seq技术。然而，这些障碍并不会阻碍scRNA-seq扩展到基因组研究前沿之外，因为克服当前的挑战将拓宽医学和生物学研究的未来前景。

  scRNA-seq的计算工作流程

        与之前的全基因组转录组学分析不同，scRNA-seq需要创新的分析工具来解决新出现的单细胞相关挑战，包括大规模数据和由于辍学事件引起的高水平噪声干扰。事实上，目前已有超过600种独立的工具可用于分析和探索单细胞转录组学数据，但由于缺乏通用的标准化，一种新工具加入其中。

        实现标准化管道的挑战源于几个原因，包括该领域的相对不成熟。根据所选择的平台，个别程序步骤的处理可能不同，导致相同输入数据集的下游分析输出不一致。此外，特定分析工具的选择在很大程度上受编程语言偏好(如R或Python)的影响，从而限制了它们的使用范围，使其仅限于使用特定编程语言的狭窄受众。进一步的，可能是最重要的障碍，是需要找到一种通用的分析策略，可以适应各种生物数据类型(细胞系、癌细胞、干细胞等)。然而，由于它们的高度多样性和不同的生物学调查，可能需要特别的计算策略。

        在这篇综述中，我们的目标是解决上述所有挑战，概述了一个标准化的工作流程，将指导读者通过scRNA-seq数据分析的关键步骤，而不考虑特定的工具和不同的生物数据类型。在这里，我们提出了四个现成的计算管道，包括原始计数归一化、特征选择、降维和聚类(图2)。完成这些步骤，用户可以在不丢失任何信息的情况下分析各自的数据。提议的管道涵盖R和Python编程语言，并使用Seurat (R) ， Scanpy (Python) ， Monocle (R)和gf-icf (R)平台，这些平台对于新手来说都很容易访问。

        使用Chromium技术检索的Tabula Muris公共数据集的一个子集演示了采用四种拟议管道的案例研究，用图表和命令行概述了不同的步骤，所有这些都可以在github上获得。由于拟议的管道可能对给定的分析是允许的或过于严格的，我们提供了定制分析设置以满足用户数据要求的指导方针。

        图2单细胞RNA测序计算分析。ScRNA-seq分析包括六个基本步骤，包括原始数据预处理，通过QC协变过滤，归一化，特征选择，线性降维，可视化和聚类:(I)原始读取被处理和量化以生成基因/条形码矩阵。(II)然后对计数基质中的细胞进行过滤，以避免对环境基因表达、凋亡细胞和多胞胎的误解。(III)由于低输入和弱信噪比干扰了分析，因此需要计数读数归一化，然后为下游分析准备数据。(四)选择较少数量的高度可变特征，以实现更快、更准确的过程。(5)基于指定的基因，进行主成分分析，降低数据维数。(VI)聚类和非线性降维步骤利用重要主成分的子集来克服数据噪声。随后，基于它们的PCA评分对细胞进行聚类和可视化。

 4 Raw读取解复用，对齐和表达式量化

        被测序器处理的捕获的转录片段(称为“reads”)被存储为基于文本的FASTQ格式。FASTQ文件包含核苷酸序列读取和它们对应的编码为ASCII(美国信息交换标准代码)字符的质量分数。虽然大多数批量方法适用于全长scRNA-seq数据集的预处理，但30’scRNA-seq协议需要独特的分析工具。30’scRNA-seq原始数据的预处理流程包括三个步骤，(1)将捕获的RNA片段分配到相关样本中，并将其存储在FASTQ文件中(即解复用);(2)将reads与参考基因组比对;(3)量化每个基因的UMI，并将其分配到相关的条形码(即细胞)。最终，每个样本编译成基因/条形码矩阵，可以进一步过滤和分析。

4.1 Demultiplexing

        在这里，我们采用CellRanger平台的常用子流水线，即mkfastq，专门用于预处理从10x-Genomics®平台获得的原始数据。虽然CellRanger为聚类和基因表达分析提供了额外的分析工具，但我们已经将其使用范围缩小到预处理步骤。作为输入，CellRanger mkfastq使用BCL文件形式的原始音序器读取。通过提供示例索引序列，mkfastq将使用示例索引将原始数据解复用到特定于示例的FASTQ文件中。

4.2 Mapping and Expression Quantification

        在量化基因表达之前，首先将原始reads与参考基因组对齐，按基因分组，并分配到其原始细胞条形码。这些步骤既可以通过CellRanger-count通过10 -Genomic®平台检索数据，也可以通过STARsolo工具应用于所有其他协议。这两种工具都需要通过解复用步骤获得的原始FASTQ文件作为输入，并执行:(1)使用预定义的白名单对单元格条形码进行错误校正;(2)利用STAR对准器进行制图;(3)对UMI进行校正和重复删除，最后(4)通过计算每个基因中唯一UMI的数量(即转录本)来定量每个细胞的基因表达。

        通过定位步骤，读取比对将原始序列分配到参考基因组中最合适的位置。虽然比对可以使用转录参考，但最好使用全基因组参考，因为它可以更容易地去除“脱靶”捕获的序列，这些序列不是被迫在转录参考上比对，而是过滤掉(见注释1和2)。接下来，过滤不一致的细胞条形码和umi，以避免数据错误。在此步骤中，在单细胞平台提供的已知细胞条形码序列的预定义列表中验证每个条形码的存在。因此，不兼容的细胞条形码要么被丢弃，要么被单个编辑距离分隔的最丰富的条形码所纠正。同样，CellRanger和STARsolo将评估UMIs的质量，如果它们共享细胞条形码和基因序列，则将单个错配纠正为更高计数的UMI序列。

        CellRanger和STARsolo都输出两个计数矩阵，过滤和未过滤，因此用户可以选择在下游分析中包括哪个。过滤后的计数矩阵由条形码/标识符组成，代表真正的细胞和每个基因的表达水平。与STARsolo不同的是，最后的CellRanger版本(3.1以上)采用了一种称为EmptyDrop的统计方法来区分单元格和空的条形码滴。在这篇评论中，我们将演示如何使用未过滤计数矩阵自主应用EmptyDrop，因为它对STARsolo和Cellranger输出都是常见的。

5质量控制和细胞筛选:如何鉴定活细胞

        目前scRNA-seq的局限性主要与低捕获效率有关，这可能导致技术噪音水平的增加。到目前为止，即使是高度敏感的scRNA-seq方案也会由于裂解或凋亡的细胞而产生一小部分低质量的条形码。因此，在进行下游分析之前，必须过滤掉与活细胞不对应的细胞条形码。这些细胞通常通过检测QC协变量分布中的异常值来识别，并通过阈值过滤(见注释3和4)。该步骤适用于所有scRNA-seq管道，并基于三个QC协变量分布的分析:(1)每个细胞条形码捕获的基因数量;(2)每个条形码中用于鉴定死亡细胞的线粒体读数的比例;(3)每个条形码的唯一umi的数量(即，一个单元的覆盖深度)。

5.1 Identify Empty Barcoded Drops

        当使用基于液滴的技术时，通常会有空液滴，因为细胞被高度稀释以产生单细胞缩放。空滴液可能被游离RNA分子污染，也被称为“环境”RNA，它们来源于细胞裂解，可能被错误地认为是细胞特异性转录物。为避免误导结果，空的条形码滴液不应包括在下游分析中。最近一种识别和过滤空水滴的方法是通过前面提到的由DropletUtils包提供的emptyDrops函数[49]。EmptyDrops是一个函数，用于使用dirichlet -多项式模型测试条形码表达轮廓偏离环境的显着程度。作为输入，它接受一个未经过滤的特征条形码矩阵并返回一个数据帧，其中每个条形码都与通过置换测试获得的p值及其相对FDR校正相关联。对后一个参数设置一个阈值，可以识别与含有细胞的液滴有显著偏差的环境概况，然后将其视为真正的细胞。在这里，我们展示了如何读取从cellRanger count管道生成的数据，并在Tabula Muris数据集的情况下检测空滴(图3a)，其中读取的数据最初是由cellRanger count管道生成的。值得注意的是，由于重要性是通过使用排列来检索的，因此需要设置种子。

5.2 Multiplet Identification

        当两个或多个细胞在一滴中被捕获并被分配到相同的细胞条形码时，就会出现多胞胎。这种错误可能被误解为单个细胞中较高的基因计数。因此，可以通过识别计数深度分布中的异常值来简单地过滤重态。对于由不同样本聚合而成的数据集，具有不同的深度尽管阈值法通常足以识别细胞多胞胎，但最近开发了新的特定工具，提供了更优雅、可能更好的解决方案。

        图3 Tabula Muris数据集上的Cell QC。(a)利用Tabula Muris数据集上的DropUtils R包中的emptyDrop函数检测空液滴。(b)在Tabula Muris数据集的每次独立运行中识别细胞多胞胎。(c) Tabula muris数据集中检测到的基因在细胞中的分布数量。(d) PCA分量作为其在Tabula Muris数据集上解释方差百分比的函数(肘形图)        

5.3 Cells Lysis

        与裂解细胞产生的转录本相关的细胞条形码通常具有计数深度低，检测到的基因很少(图2c)和线粒体读取的高比例。在这种情况下，与细胞质RNA不同，由于线粒体膜未受损，大多数线粒体RNA是保守的。因此，过滤掉超过10%的线粒体相关读取的条形码是可以接受的。然而，在设置阈值时，应该始终考虑数据集的生物学特性，因此线粒体读取可接受百分比的阈值可能会根据研究的生物学模型而变化。对于癌细胞或呼吸或代谢过程增加的特定细胞类型，高水平的线粒体rna是模型本身所固有的(见注6)。

6 Start Working with the Scanpy, Seurat, Monocle, and gf-icf Pipelines

        随着单细胞测序技术的兴起，许多统计方法被开发出来，以解决scRNA-seq分析的不同步骤。Seurat、Monocle 3、Scanpy和gf-icf等不同的工具包将这些独立的方法组合在一起，以提供一个单一的工作流。最流行的基于代码的平台之一是Seurat，它提供了广泛的教程和分析工具。另一个适应良好的平台是Monocle，自首次引入以来，它在很大程度上促进了轨迹推理领域的发展。最新版本Monocle 3提供了伪时间排序和基本的scRNA-seq聚类管道，方便用户使用。Scanpy是该领域的一个相对较新的添加，它允许分析多达一百万个单元格或更多的大型数据集，因为它改进了计算缩放。在这里，我们还测试了最近引入的一种名为gf-icf的方法，该方法基于术语频率逆文档频率(TF-IDF)的数据转换模型，该模型已广泛用于文本挖掘领域，其中稀疏和零膨胀数据很常见。对于下游分析，每个管道使用R或Python编程语言。为了解释输出和理解分析工具的基础，每个步骤都将在所有四个管道中进行检查和比较。

7 Gene Filtering: How to Remove “Noisy” Genes

        scRNA-seq数据集通常包括在数千个细胞中测量的超过25,000个基因，其中许多可能没有信息，因为它们大多含有零计数，应该在开始下游分析之前过滤掉。基因过滤可以通过降低数据的维数和减少多余的零计数来加快数据的处理速度，从而改进数据的归一化步骤和所有的下游分析。通常，定义一个固定的阈值，据此去除在少数细胞中检测到的基因(见注释7和8)。

8数据归一化：如何使基因表达在单个细胞之间具有可比性

        数据归一化解决了由于计数深度变化而产生的不必要的偏差，同时保留了真正的生物差异。从每个细胞中捕获的mRNA数量可能会由于生物变异或在整个scRNA-seq过程中遗传的技术影响而发生变化，包括单细胞制备、文库构建和序列步骤。通过标准化，考虑到每个细胞捕获的mRNA分子的丰度，每个基因的表达被重新调整，以便使基因表达在单个细胞之间具有可比性。估计每个单元格的比例因子的方式主要不同于当前可用的大量归一化方法。

        如上所述，由于生物学和技术原因，scRNA-seq数据通常很稀疏。因此，采用了批量RNA-seq的归一化方法，如TMM 和DESeq，可能受到零通货膨胀的影响。为了解决这一问题，近年来单细胞归一化程序不断发展。然而，目前最常用的scRNA-seq数据归一化方法是百万计数(count per million, CPM)，这是一种继承自bulk RNA-seq的线性全局标度方法。

        大量rna测序。与所研究的生物过程无关的另一个变异源可能来自处理不同批次的样品。当用相同的细胞独立地重复实验时，例如，由不同的操作人员或采样不同的实验时间线时，可能会产生批处理效应。标准的归一化程序不能对批效应进行校正，影响了对真实生物效应的分析。最近已经开发了几种方法来解释scRNA-seq数据中的批处理效应，尽管ComBat，一种最初为微阵列数据开发的方法，在中低复杂度的单细胞实验中也表现良好。

        这里提出的所有四个管道都通过CPM方法解释了规范化步骤。Seurat、Monocle和Scanpy使用CPM的对数变换来减少细胞深度的可变性(见注释9-11)，并使用一些高级选项来重新缩放一些变异源的数据，包括影响细胞周期。使用gf-icf管道，基因通过它们的逆细胞频率被重新缩放，细胞被重新缩放为欧几里得范数等于1 (L2归一化)，以便考虑细胞深度的可变性。

9 Feature Selection: How to Discard “Uninformative” Genes

        一个大规模的scRNA-seq数据集可以很容易地包括在超过10,000个细胞中测量的超过25,000个基因，其中许多基因不具有信息性，因为它们大多含有零计数。特征选择的目的是检测出具有相关生物学信息的基因，排除不具有相关生物学信息的基因。ScRNA-seq数据维度可以保持相当高，即使经过基因过滤，仍然保留了大量的基因(>10,000)。特征选择可以大大加快处理速度，因为它通过过滤“无信息”的基因来降低数据维数。这通常是通过选择有限数量的高可变基因(HVG)来指导进一步分析来完成的。HVG具有很高的信息量，因为它们对数据配置有重大影响，因此可以保持数据的完整性和可再现性。通常，根据测定的大小选择1000-5000 HVG(见注12)。每条管道都有自己独特的方法来检测和选择hgv。使用Scanpy将基因按其平均表达量进行分类，并在每个分类箱中选择方差均值比最高的基因作为hvg。

        另一方面，Seurat首先使用局部多项式回归函数对均值-方差关系进行建模。然后，给定拟合曲线的期望方差和观测均值，对特征值进行标准化，并计算每个基因在所有细胞中的方差(见注13)。与Seurat和Scanpy不同，Monocle不包括这一步，而gf-icf只有在需要识别集群间差异表达基因时才进行特征选择。虽然差异不大，但在gf-icf中，特征选择也是由Chen等人提出的均值/方差关系建模来完成的。gf-icf管道中的特征选择是在归一化步骤中构建的，通过逆细胞频率对基因表达进行重新缩放，并在降维步骤中考虑过滤基因的总数。

10 降维:如何总结和可视化scRNA-seq数据

10.1线性降维：总结scRNA-seq数据

降维的目的是通过最优地保持其关键属性，将数据的复杂性压缩为一个低维空间。降维方法对于scRNA-seq数据的聚类、可视化和总结至关重要。线性降维方法通常作为非线性降维方法的预处理步骤。最流行的线性降维算法是主成分分析（主成分分析）。通常，10-50个重要的主成分被选择，然后用作非线性降维方法的输入。主成分高度表明了数据集中异构性的主要来源。

PCA用于总结数据集的前N个主成分(见注释14)。要使用的主成分的数量通常是通过手动检查弯头图(图3d)来确定的，其中主成分是作为它们所占变异性的函数来绘制的，而要使用的主成分的数量是由观察到“弯头”的点来确定的。还可以使用其他方法，包括jackstraw和每个主成分中主要基因的热图。但是，在选择要使用的重要主成分时，最好选择较高的主成分，以避免信息丢失。

10.2scRNA-seq数据可视化的非线性降维

对scrna数据可视化的降维方法使用了捕获scRNA-seq数据的非线性的方法，避免了表示的过度拥挤（见注15）。两种最常用的方法分别是t-分布随机邻域嵌入（t-SNE）和统一流形近似和投影（UMAP）。t-SNE是一种有效地突出局部数据结构的随机方法低维，表示细胞群为不同的集群。但是，t-SNE不能保持全局结构，因此簇之间的距离毫无意义。UMAP是一种更新的非线性降维技术，它能够保持数据的局部和全局结构，优于实际大规模的t-SNE。对于线性和非线性数据降维存在一些额外的方法，但它超出了本教程的范围，以审查所有现有的方法，而我们更倾向于关注目前被scRNA-seq社区接受的最佳实践和方法。然而，Moon等人对单细胞转录组数据的线性和非线性降维方法进行了详细的回顾。

11 聚类分析:如何识别细胞亚群

        由于转录上不同的细胞群通常对应不同的细胞类型，因此scRNA-seq的一个关键目标在于基于转录相似性鉴定细胞亚群。因此，将细胞组织成组(即集群)可以允许在单个细胞状态下重新检测细胞类型或识别不同的亚群(见注释16)。聚类是一个古老的无监督机器学习问题，其目的是通过知道一组未标记对象的相似性得分(即距离)来确定它们的内在分组。文献中已经提出了大量的距离度量来计算感兴趣对象之间的相似性分数，包括欧几里得距离，余弦距离和基于相关的距离。几种无监督聚类方法已被应用于划分单细胞数据，可进一步分为三类:(1)k-means，(2)分层聚类和(3)社区检测方法。对于单细胞数据分析，所有方法都是在pc降维空间上进行特征选择和数据降维后应用的。然后将识别的细胞簇覆盖到可视化空间上。

        k-means算法使用迭代方法将细胞划分为预定义数量的簇(k)。在每次迭代中，使用欧几里得距离将细胞分配到最近的质心。替代距离，如基于相关或余弦距离，也可用于单细胞数据分析。在每次迭代结束时重新计算质心的位置，由于质心的起始位置是随机选择的，因此k-means算法通常会多次运行。

        虽然速度快，但k-Means需要知道划分细胞的初始簇数(k)，这通常是未知的，必须通过额外的复杂分析来解决。

        分层聚类是一种寻求建立集群层次结构的划分方法，一般分为聚集型和分裂型两种。聚合分层聚类技术遵循“自下而上”的方法，其中最初每个单元代表一个单独的集群，然后逐渐合并相似的集群，直到获得唯一的集群。另一方面，分裂的分层集群遵循“自上而下”的方法，其中所有细胞从单个集群开始，然后逐渐分裂。分层聚类产生一个树形图，其中聚类是通过在预定的距离上切割树获得的，该距离可以使用自举方法启发式地解决。层次聚类在scRNA-seq数据中的应用实例可以在CIDR、SINCERA和producduce中找到。然而，分层聚类方法通常比k-means工作速度慢，并且在大规模scRNA-seq数据集上表现不佳。

        社区检测技术是一种可扩展的聚类方法，它适用于大规模图，可以有效地聚类数十万甚至数百万个单元。根据定义，图g¼(V,E)由节点V(即细胞)和表示细胞对之间相似程度的边的集合组成。这个细胞图可以使用K-最近邻(KNN)算法来构建，该算法应用于pc精简空间，其中每个细胞都连接到它的K个最相似的细胞。然后，任意两个单元之间的边权通过使用它们共享的邻居的比例，通过Jaccard相似性来改进。

        寻找群落意味着将细胞聚集成组，组内的边缘密度高于组间的边缘密度。图的群落结构的一个度量是模块化，即，如果边是随机分布的，则落在给定组内的边的分数减去期望分数。模块化是基于这样一种思想，即随机图不期望具有集群结构。最流行的基于模块化的检测算法是Louvain，它是由PhenoGraph引入的，Seurat、Scanpy和gf-icf也在使用它。

        在运行基于图的聚类时，需要为基于模块化优化的社区检测算法设置分辨率参数。分辨率参数与观测群落的规模相关。特别是，分辨率参数越高，小社区的数量就越多。在我们的管道中，我们将分辨率参数设置为0.5，这通常表示一个很好的权衡。

12差异表达:如何标注细胞群体

        通过差异表达分析识别标记基因(即聚类基因标记)，可以对聚类算法识别的细胞群进行表征和注释。标记基因是通过将每一个细胞群中的细胞与所有其他细胞进行比较来确定的。已经开发了一些差分表达测试方法专门用于处理scRNA-seq数据中dropout元素的存在，包括贝叶斯方法和MAST，但在考虑大规模scRNA-seq数据集时，它们的计算效率不高。因此，更快的测试被用于检测差异表达基因，如Seurat、Scanpy和gf-icf中实现的Wilcoxon秩和测试，而Monocle使用广义加性模型(VGAM)。Seurat、Scanpy和Monocle还提供了其他复杂的测试。一旦确定了每个细胞簇的基因特征，就可以使用基因本体富集分析(GOEA)和基因集富集分析(GSEA)[85]等附加分析来确定每个细胞簇中活跃的生物过程。

13实施管道的效果评价与比较

        为了评估四种管道在识别细胞组方面的性能(图4a-d)，我们通过使用每对簇之间的平均Jaccard系数计算了不同方法产生的簇之间的一致性(图4e)。然后，我们使用从每种方法中检索到的集群对细胞类型进行分层聚类(图4f)，并表明不同的方法产生具有生物学意义的分区。我们还观察到，同一簇中的细胞属于相同的谱系，但具有不同的粒度级别，这可以通过更改用于识别细胞簇的分辨率参数来调整。

14附加分析:如何重建细胞转录动力学

        根据要解决的生物学问题，人们可能会考虑利用其他现有工具进一步研究单细胞数据，这些工具可能提供其他级别的信息。生物机制是高度动态的过程，因此不能总是用离散的方法来很好地描述，比如聚类。细胞可以在环境变化和外部扰动的控制下跨越几种转录状态。因此，为了模拟这样的延续生物系统，包括发育过程，一种新的计算方法，称为轨迹推理，在过去的几年里已经发展起来。这些方法使用在不同时间点收集的生物过程基础细胞群生成的scRNA-seq数据，并尝试沿着进化轨迹计算它们的顺序，这可能具有不同的拓扑结构(即线性，分分叉或甚至更复杂的图结构)。一旦对细胞进行了排序，整个推断轨迹中的基因表达模式就可以用来识别控制细胞命运决定的关键调节基因。我们首先在Monocle中引入了伪时间的概念，将其作为描述发育系统的稳健方法。自从Monocle(现在已经是第三个版本了)以来，可用方法的数量呈指数级增长。最近，一种新提出的推断发育轨迹的方法与其他方法有很大不同，即使用最优运输问题对细胞过程进行建模。有趣的是，迄今为止，已经开发了100多种方法来推断细胞轨迹。一旦轨迹被重建，RNA速度可以叠加到推断的轨迹上，为重建的动力学过程增加方向性。

        

        图4实现管道的细胞聚类、比较和评价。(a-d)由四个实现的管道中的每一个产生的UMAP可视化，其中单元格根据其所属的集群着色。(e)在四个已执行的管道中商定已确定的细胞群。(f)利用每个独立管道的聚类结果对细胞类型进行分层聚类，重构细胞谱系.        

15 讨论和未来的发展方向

        随着scRNA-seq技术的兴起，越来越多的分析方法被引入科学界。尽管分析选项的范围很广，但缺乏标准化导致了很高的进入壁垒。在本综述中，我们提出了四种现成的管道，用于分析scRNAseq数据，可以适应各种生物数据类型。对于新手来说，这些计算管道提供了一个有效而简单的工作流程，包括原始计数的规范化、特征选择、降维和数据聚类。提议的管道涵盖R和Python编程语言，并使用Seurat (R)、Scanpy (Python)、Monocle (R)和gf-icf (R)平台。为了确保数据一致性，具有解释输出的能力是很重要的，因此我们回顾了scRNA-seq分析的关键步骤。我们还强调了指导方针，并提供了标准化的值来过滤和减少数据维度。用户有责任仔细评估其分析的输出，并在必要时调整管道默认设置以适应源数据。此外，随着该领域的迅速发展，本综述可能落后于最新的工具。因此，我们建议将此综述作为基本的工作流程指南，同时与该领域的创新保持一致。

        由于单细胞测序不再局限于转录实验，而是允许捕获其他数据类型，包括DNA, ChIP和ATAC，我们认为未来的管道必须能够应对多组学数据整合。单细胞多组学可以同时获得活细胞各个层面的信息，包括DNA、RNA、蛋白质和表观遗传修饰[90-92]。将这些不同的“组学”信息整合到一个维度中，将使我们对细胞命运调控和表型有更全面的了解。有趣的是，另一项需要高尺度计算工具的新技术是空间转录组学，它可以识别组织的细胞类型和空间组成[93,94]。这种方法可以通过在数据中添加另一个引导维度来帮助提高所调查系统的准确性。通过对细胞进行位置注释，可以在高度异质的系统(如类器官)中精确地聚集不同的亚群，并跟踪它们之间的时空动态。因此，保留空间位置的能力将提供对组织组织、功能和发展的更好感知。

        另一个角度是使用高通量scRNA-seq技术进行个性化医疗。通过scRNAseq筛选不同的细胞类型和组织，以根据患者的个体特征定制适当的药物治疗已经进行了一些努力[95,96]。开发结合机器学习方法的新工具可能会增加精准医学领域的进步，并使其更接近临床应用。我们相信，具有上述特性的创新工具将站在科学的最前沿。

16 Notes

1. 在进行进一步分析之前，建议先进行测序和图谱统计，这些数据通常由用于预处理的生物信息学工具提供。例如，少于70%的条形码相关读取可能表明高水平的环境RNA(由于大量的裂解细胞或组织分离后洗涤不足)。
2. 确定的基因组图谱应该超过总数的80%。
3. 基于qc的离群值检测，即多重滤波和分解细胞滤波，应同时采用这些协变量进行。 4. 过滤离群细胞的阈值应尽可能允许，以避免过度的dropout效应。还可以更进一步一旦进行下游分析，就会进行调整，以更好地解释数据。
4. 由于转录本覆盖率可能在不同的样本之间有所不同，因此为每个样本分别设置阈值是必要的。
5. 在设置阈值时，应考虑数据集的生物学特性，因为增加的呼吸或代谢过程也可能导致高线粒体读取。
6. 7. 所选阈值应尽可能允许，以避免dropout效应或稀有细胞群的移除。
7. 8. 一个可接受的指导方针是将阈值调整为最小的簇大小或在数据集中表达超过1-5%的基因数量。
8. 9. 尽管选择了规范化方法，但数据转换(例如，对数转换)应该始终应用，因为大多数用于下游分析的工具期望得到正态分布的数据。
9. 10. 尽量避免校正生物批次，除非你想推断轨迹，并且这种校正不会掩盖其他感兴趣的生物信息。
10. 11. 对技术协变量和生物协变量进行批量校正时，应同时进行。
11. 12. hvg的选择可能会影响下游分析，尽管已有研究表明，在200到2400 hvg之间的选择并不影响在较低维度(即PCA空间)中的表示[9]。
12. 13. 基于均值和方差的特征选择不能在缩放到零均值和单位方差的数据上进行。
13. 14. 主成分也可用于检查技术协变量对数据的影响[65]，或解决特定基因在数据集中的作用[69]。
14. 15. 非线性降维方法是数据可视化的强大工具，而不是汇总。
15. 16. 下游分析需要汇总数据，例如PCA或扩散图。