scCB2 1.8.0
基于液滴的单细胞RNA测序(scRNA-seq)是一种强大的和广泛使用的方法，用于分析单个细胞的全基因组基因表达。目前基于液滴的商业技术(如10X Genomics)利用凝胶珠，每个凝胶珠都含有由珠特异性条形码组成的寡核苷酸指数，结合独特的分子标识符(UMIs)和寡核苷酸dt标签来启动聚腺苷酸RNA。感兴趣的单细胞在微流控芯片的一个通道中与试剂结合，在另一个通道中与凝胶珠结合，形成乳液中的凝胶珠，或GEMs。在凝胶珠被溶解之前，寡核苷酸指数在每个GEM反应囊泡内结合聚腺苷酸化RNA，将结合的寡核苷酸释放到溶液中进行反转录。根据设计，每个结果的cDNA分子包含一个UMI和一个gem特异性条形码，理想情况下，标记来自单个细胞的mRNA，但实际情况往往不是这样。为了在基于液滴的单细胞RNA-seq实验中区分真实细胞与背景条形码，我们介绍了CB2而且scCB2，其对应的R包。
CB2通过引入聚类步骤来扩展EmptyDrops方法，该步骤将相似的条形码分组，然后进行统计测试，以识别具有不同背景的表情分布的组。虽然在许多设置中都有优势，但用户将从注意到这一点中受益CB2不检测双胞或多胞，因此，一些高计数的鉴定可能由两个或多个细胞组成。识别多胞胎的方法可能在应用后被证明是有用的CB2．同样重要的是要注意，任何区分细胞与背景条形码的方法在识别低质量细胞时在技术上是正确的，因为受损细胞表现出与背景不同的表达谱。具体来说，线粒体基因在受损细胞中的表达通常较高。这些细胞通常对下游分析不感兴趣，因此应该被移除。中的GetCellMat函数scCB2可用于此目的。
从Bioconductor安装:
如果(!requireNamespace("BiocManager"， quiet = TRUE)) install.packages("BiocManager")::install("scCB2")QuickCB2是一个一体机功能，应用CB2对10x Cell Ranger原始数据，并获得一个由CB2在默认设置下识别的真实细胞矩阵。通过指定AsSeurat = TRUE，返回一个Seurat对象，以便用户可以直接应用Seurat管道进行下游分析。
用法:
library(scCB2) #如果原始数据在一个目录中有三个独立的文件#并且您想要在默认的1%控制FDR: RealCell <- QuickCB2(dir = "/path/to/raw/data/directory") #如果原始数据是HDF5格式并且#您想要一个Seurat对象在默认的FDR阈值下:RealCell_S <- QuickCB2(h5file = "/path/to/raw/data/HDF5"， AsSeurat = TRUE)下面的示例说明了它是如何工作的以及最终的输出是什么样的详细步骤．
计算速度与输入数据集的大小和结构有关。作为参考，scCB2在使用6个核心的中型数据集(约30,000个基因和8000个细胞调用后)上运行用时不到10分钟。
对于希望将CB2输出单元格矩阵保存为10x格式的用户(例如“matrix。mtx”、“条形码。Tsv”和“genes.tsv”)，有现成的软件包可以提供帮助。例如在包装中DropletUtils：
DropletUtils::write10xCounts(path = "/path/to/save/data"， x = RealCell)目前应用最广泛的基于液滴的协议是10x Chromium。我们的包提供了直接读取10x Cell Ranger输出文件的功能，并生成一个可以读入r的特征-条形码计数矩阵在这里．
我们的包包含一个10倍数据的小子集，mbrainSub，对应的前5万个条形码1k脑细胞来自E18小鼠．
我们首先生成10x的输出文件mbrainSub，然后使用内置函数读取它。
库(scCB2)库(summarizeexperiment)数据(mbrainSub)数据。dir <- file.path(tempdir()，"CB2_example") DropletUtils::write10xCounts(data. path)mbrainSub, version = "3") list.files(data.dir)##[1] "条形码.tsv.gz" "features.tsv.gz" "matrix.mtx.gz"对于<3的Cell Ranger版本，来自10x Cell Ranger输出的原始数据包含“条形码”。tsv”、“基因。Tsv”和“matrix.mtx”。对于Cell Ranger版本>=3，输出文件为" barcode .tsv.gz "， " features.tsv.gz "和" matrix.mtx.gz "。现在我们将这些文件读入R，并与原始数据矩阵进行比较。
mbrainSub_2 <- Read10xRaw(data.dir) same (mbrainSub_2, mbrainSub_2)##[1]真如果原始数据不是来自10x Chromium管道，用户可以手动创建逐条形码特征计数矩阵，其中行表示基因，列表示条形码。基因id和条形码id应该是唯一的。计数矩阵的格式可以是稀疏矩阵或标准矩阵。
CB2以及其他类似方法的关键参数是背景截断，它将条形码分为两组:(1)最可能是背景空液滴的小条形码;(2)其余的条形码可以是背景的，也可以是单元格的，有待测试。这些小条形码将用于估计背景分布，从而指导从背景中识别细胞。对背景分布进行无偏估计是至关重要的。
默认情况下，背景截止值设置为100，这意味着所有UMI总计数小于或等于100的条形码都用于估计背景分布。根据经验，这种设置在大多数真实世界的数据集中工作得很好。然而，对于具有特殊结构的数据集，或检测到的单元格数量意外地低或高，建议重新评估背景截断的选择。
适当的背景截断线应该足够大，以包含足够的背景信息，但也不能太大，以免将真实的单元格错误地放入背景组中。基于经验知识，我们建议(1)将超过90%的条形码放入背景，或(2)将超过10%的UMI计数放入背景。这保证了我们有足够的信息对背景截止进行无偏估计。从100开始，满足任一条件的最小临界值就是推荐临界值。
check_cutoff <- CheckBackgroundCutoff(mbrainSub) check_cutoff$summary_table## cutoff n_bg_bcs n_bg_counts prop_bg_bcs prop_bg_counts ## 1 100 49714 97013 0.99428 0.09659609 ## 2 150 49791 106163 0.99582 0.10570677 0.99644 0.11109750 ## 4 250 49846 116966 0.99692 0.11646334 ## 5 300 49861 121121 0.99722 0.12060049 ## 6 350 49861 124047 0.99740 0.12351391 ## 8 450 49873 126456 0.99746 0.12591256 ## 9 500 49879 127901 0.99758 0.12735135check_cutoff recommended_cutoff美元## [1] 100在本例中，默认背景截止100是合适的，因为它将超过90%的条形码放入背景中，以及超过10%的UMI计数放入背景中。一般来说，我们建议总是检查背景截止。
主要功能CB2FindCell接受原始计数矩阵作为输入，并返回实际单元格、测试统计信息和p值。现在我们应用CB2FindCell在mbrainSub，将FDR控制在0.01水平(默认)，假设总计数小于或等于100的条形码均为背景空液滴(默认)，采用2核并行计算(默认:2)。详细信息请参见CB2FindCell ?．
CBOut <- CB2FindCell(mbrainSub, FDR_threshold = 0.01, lower = 100, Ncores = 2)时差1.27602分钟str(assay(CBOut)) #细胞基质##正式类“dgCMatrix”[包“矩阵”]有6个插槽## ..@ I: int[1:274054] 62 469 538 663 758 903 1251 1266 1396 1532…# # . .@ p: int[1:175] 0 219 2503 4608 6288 6512 6796 7194 7335 7789…# # . .@ Dim: int [1:2] 27998 174 ## ..@ Dimnames: 2 ## .. ..的列表$: chr [1:27998] "Xkr4" "Gm1992" "Gm37381" "Rp1"…## .. ..$: chr [1:174] " aaactggtgctcttc -1" "AAACGGGAGCCACGTC-1" "AAACGGGAGCGAGAAA-1" " aaacgggcagtttac -1"…# # . .@ x : num [1:274054] 1 1 1 1 1 1 1 1 1 1 ... ## ..@ factors : list()str(metadata(CBOut)) #测试统计数据，p值等## $ ClusterStat:'data.frame': 5 obs四个变量:## ..$ corr: num[1:5] 0.6715 0.6773 0.0509 0.38 0.5547 ##…$ count: num [1:5] 260 315 271 113 4679 ## ..$ pval: num [1:5] 0.297702 0.102897 0.000999 0.003996 0.000999 ## ..$ padj: num [1:5] 0.2977 0.12862 0.0025 0.00666 0.0025 ## $ Cluster:5 ## ..$: chr [1:85] "AACACGTAGGCAAAGA-1" " acagccgcaaacgga -1" "AACACGTGTATAGGTA-1" "AACTCAGTCTCGGACG-1"…# # . .$: chr [1:41] "AACTTTCAGGCCCGTT-1" " aactttcgccgtga -1" " aactttcgtaccgta -1" " aaccatggttgccctt -1"…# # . .$: chr [1:3] " aaggcaggtcttgcg -1" "ACATCAGCACATTAGC-1" "AAGACCTAGACTGGGT-1" ## ..$ : chr [1:6] "AATCCAGAGATCTGAA-1" "AACCATGTCTCGCATC-1" "AACACGTTCTCTAGGA-1" "AACTCCCGTGCTTCTC-1" ... ## ..$ : chr [1:23] "ACCAGTAAGTGCGATG-1" "AAGTCTGAGGACAGCT-1" "ACATACGCAAGCCGTC-1" "AAGGCAGGTCTAGAGG-1" ... ## $ BarcodeStat:'data.frame': 203 obs. of 4 variables: ## ..$ logLH: num [1:203] -762 -770 -886 -454 -9992 ... ## ..$ count: num [1:203] 337 304 312 207 5717 ... ## ..$ pval : num [1:203] 0.9416 0.4256 0.0006 1 0.0001 ... ## ..$ padj : num [1:203] 1 0.663139 0.001137 1 0.000205 ... ## $ background : Named num [1:11121] 8 3 2 0 2 2 0 1 3 1 ... ## ..- attr(*, "names")= chr [1:11121] "Mrpl15" "Lypla1" "Tcea1" "Rgs20" ...如果读者对输出的测试信息不感兴趣，GetCellMat能直接提取出真实的细胞矩阵吗CB2FindCell输出。它还提供了一种过滤选项，根据线粒体基因表达的比例去除破碎的细胞。现在我们应用GetCellMat在CBOut，过滤出线粒体比例大于0.25的细胞(默认值:1，不过滤)。
GetCellMat(CBOut, MTfilter = 0.25) str(RealCell)##正式类“dgCMatrix”[包“矩阵”]有6个插槽## ..@ I: int[1:273036] 62 469 538 663 758 903 1251 1266 1396 1532…# # . .@ p: int[1:172] 0 219 2503 4608 6288 6512 6796 7194 7335 7789…# # . .@ Dim: int [1:2] 27998 171 ## ..@ Dimnames: 2 ## .. ..的列表$: chr [1:27998] "Xkr4" "Gm1992" "Gm37381" "Rp1"…## .. ..$: chr [1:171] " aaactggtgctcttc -1" "AAACGGGAGCCACGTC-1" "AAACGGGAGCGAGAAA-1" " aaacgggcagtttac -1"…# # . .@ x : num [1:273036] 1 1 1 1 1 1 1 1 1 1 ... ## ..@ factors : list()后CB2预处理后，真实的细胞矩阵仍然是矩阵格式，因此可以直接用于下游统计分析。例如，如果我们想使用修拉管道，我们可以轻松地创建一个Seurat对象使用
SeuratObj <- Seurat::CreateSeuratObject(counts = RealCell, project = "mbrain_example"一个Seurat 27998类的对象，在1个测定法中横跨171个样品##活性测定:RNA(27998个特征，0个可变特征)在默认参数下，可以直接使用all-in-one功能QuickCB2从10倍原始数据中得到真实的单元格矩阵。
RealCell_Quick <- QuickCB2(dir = data. properties)目录:Ncores = 2时差1.278609分钟str (RealCell_Quick)##正式类“dgCMatrix”[包“矩阵”]有6个插槽## ..@ I: int[1:273436] 62 469 538 663 758 903 1251 1266 1396 1532…# # . .@ p: int[1:17 17] 0 219 2503 4608 6288 6512 6796 7194 7335 7789…# # . .@ Dim: int [1:2] 27998 170 ## ..@ Dimnames: 2 ## .. ..的列表$: chr [1:27998] "Xkr4" "Gm1992" "Gm37381" "Rp1"…## .. ..$: chr [1:170] " aaactggtgctcttc -1" "AAACGGGAGCCACGTC-1" "AAACGGGAGCGAGAAA-1" " aaacgggcagtttac -1"…# # . .@ x : num [1:273436] 1 1 1 1 1 1 1 1 1 1 ... ## ..@ factors : list()现在可以进行规范化和聚类等下游分析了。Seurat教程示例:https://satijalab.org/seurat/v3.1/pbmc3k_tutorial.html
sessionInfo ()## R版本4.2.1(2022-06-23)##平台:x86_64-pc-linux-gnu(64位)##运行在Ubuntu 20.04.5 LTS ## ##矩阵产品:默认## BLAS: /home/biocbuild/bbs-3.16-bioc/R/lib/libRblas。/home/biocbuild/bbs-3.16-bioc/R/lib/libRlapack。所以## ## locale: ## [1] LC_CTYPE=en_US。UTF-8 LC_NUMERIC= c# # [3] LC_TIME=en_GB LC_COLLATE= c# # [5] LC_MONETARY=en_US。utf - 8 LC_MESSAGES = en_US。UTF-8 ## [7] LC_PAPER=en_US。UTF-8 LC_NAME= c# # [9] LC_ADDRESS=C lc_phone = c# # [11] LC_MEASUREMENT=en_US。utf - 8 LC_IDENTIFICATION = C附加基本包:# # # # # # [1]stats4统计图形grDevices跑龙套数据集方法# # # # # #[8]基地其他附加包:# # [1]SummarizedExperiment_1.28.0 Biobase_2.58.0 # # [3] GenomicRanges_1.50.0 GenomeInfoDb_1.34.0 # # [5] IRanges_2.32.0 S4Vectors_0.36.0 # # [7] BiocGenerics_0.44.0 MatrixGenerics_1.10.0 # # [9] matrixStats_0.62.0 scCB2_1.8.0 # # [11] knitr_1.40 BiocStyle_2.26.0 # # # #通过加载一个名称空间(而不是附加):[11] foreach_1.5.2 htmltools_0.5.3 ## [13] fansi_1.0.3 magrittr_2.0.3 ## [15] tensor_1. 1.5 cluster_2.1.4 ## [17] doParallel_1.0.17 ROCR_1.0-11 # [19] limma_3.54.0 globals_0.16.1 ## [21] r.l uts_2.12.1 spatstat.sparse_3.0-0 ## [23] colorspace_2.0-3 ggrepel_0.9.1 ## [25] xfun_0.34 dplyr_1.0.10 ## [27] RCurl_1.98-1.9jsonlite_1.8.3 # # [29] spatstat.data_3.0-0 progressr_0.11.0 # # [31] survival_3.4-0 zoo_1.8-11 # # [33] iterators_1.0.14 glue_1.6.2 # # [35] polyclip_1.10-4 gtable_0.3.1 # # [37] zlibbioc_1.44.0 XVector_0.38.0 # # [39] leiden_0.4.3 DelayedArray_0.24.0 # # [41] DropletUtils_1.18.0 Rhdf5lib_1.20.0 # # [43] future.apply_1.9.1 SingleCellExperiment_1.20.0 # # [45] HDF5Array_1.26.0 abind_1.4-5 # # [47] scales_1.2.1 DBI_1.1.3 # # [49] edgeR_3.40.0 spatstat.random_3.0-0 # # [51] miniUI_0.1.1.1 Rcpp_1.0.9 # # [53]viridisLite_0.4.1 xtable_1.8-4 # # [55] spatstat.core_2.4-4 reticulate_1.26 # # [57] dqrng_0.3.0 htmlwidgets_1.5.4 # # [59] httr_1.4.4 RColorBrewer_1.1-3 # # [61] ellipsis_0.3.2 Seurat_4.2.0 # # [63] ica_1.0-3 pkgconfig_2.0.3 # # [65] R.methodsS3_1.8.2 scuttle_1.8.0 # # [67] uwot_0.1.14 deldir_1.0-6 # # [69] sass_0.4.2 locfit_1.5 - 9.6 # # [71] utf8_1.2.2 reshape2_1.4.4 # # [73] tidyselect_1.2.0 rlang_1.0.6 # # [75] later_1.3.0 munsell_0.5.0 # # [77] tools_4.2.1 cachem_1.0.6 # # [79] cli_3.4.1 generics_0.1.3## [81] ggridges_0.5.4 evaluate_0.17 ## [83] matrix_1.5 .4.1 fastmap_1.1.0 ## [85] goftest_1.2-3 yaml_2.3.6 ## [87] fitdistrplus_1.1-8 purrr_0.3.5 ## [91] pbapply_1.5-0 future_1.28.0 ## [93] sparseMatrixStats_1.10.0 mime_0.12 ## [95] R.oo_1.25.0 compiler_4.2.1 ## [97] plotly_4.10.0 png_0.1-7 ## [99] spatstat.utils_3.0-1 tibble_3.1.8 ## [101] bslib_0.4.0 stringi_1.7.8 ## [103] rgeos_0.5-9 lattice_0.20-45 ## [105] Matrix_1.5-1 vctrs_0.5.0 ## # [107] pillar_1.8.1[109] rhdf5filters_1.10.0 BiocManager_1.30.19 ## [111] spatstat.geom_3.0-3 lmtest_0.9-40 ## [113] jquerylib_0.1.4 RcppAnnoy_0.0.20 ## [115] data.table_1.14.4 cowplot_1.1.1 ## [117] bitops_1.0-7 irlba_2.3.5.1 ## [119] httpuv_1.6.6 patchwork_1.1.2 ## [121] R6_2.5.1 bookdown_0.29 ## [123] promises_1.2.0.1 gridExtra_2.3 ## [125] KernSmooth_2.23-20 parallelly_1.32.1 ## [129] codetools_0.2-18 MASS_7.3-58.1 ## [129] assertthat_0.2.1 rhdf5_2.42.0 ## [113] SeuratObject_4.1.2 sctransform_0.3.5 ## [133] GenomeInfoDbData_1.2.9 mgcv_1.8-41 ## [135] parallel_4.2.1 rpart_4.1.19 ## [137] grid_4.2.1 beachmat_2.14.0 ## [139] tidyr_1.2.1 rmarkdown_2.17 ## [141] DelayedMatrixStats_1.20.0 Rtsne_0.16 ## [143] shiny_1.7.3倪，Z.，陈，S.， Brown, J.， & Kendziorski, C.(2020)。CB2提高了基于液滴的单细胞RNA测序数据的细胞检测能力。生物工程学报，21(1)，137。https://doi.org/10.1186/s13059-020-02054-8