数据基础架构

2020-03-30 14:59:36 浏览数 (2)

好的书籍是人类进步的阶梯,但有些人却找不到优秀的阶梯,为此我们开设了书籍翻译这个栏目,作为你学习之路的指路明灯;分享国内外优秀书籍,弘扬分享精神,做一个知识的传播者。

希望大家能有所收获!

序言

通过Bioconductor进行单细胞分析(一)

本教程是翻译自 Orchestrating Single-Cell Analysis with Bioconductor(https://osca.bioconductor.org/), 由于前三章是序言和软件介绍等无关紧要的内容,所以我们从第四章数据基础框架开始

第四章 数据基础架构

使用Bioconductor软件包的优势之一是它们利用了通用的数据基础架构,从而使分析可以在各种软件包之间实现互操作。此外,要使此基础架构强大且可扩展,需要付出大量工程努力。在这里,我们将详细描述

使用Bioconductor软件包的优势之一是它们利用了通用的数据基础架构,从而使分析可以在各种软件包之间实现互操作。此外,要使此基础架构强大且可扩展,需要付出大量工程努力。在这里,我们将详细描述SingleCellExperiment对象(或简写成sce`形式),以描述如何构造,在下游分析中使用该对象以及如何存储各种类型的主数据和元数据。

4.1

先决条件

按以下方式安装SingleCellExperiment包:

代码语言:javascript复制
BiocManager::install('SingleCellExperiment')

此外,还会用到scaterscran和CRAN包uwot。安装命令如下:

代码语言:javascript复制
BiocManager::install(c('scater', 'scran', 'uwot'))

对于本节,需要加载的是SingleCellExperiment包:

代码语言:javascript复制
library(SingleCellExperiment)

4.2

SingleCellExperiment Class

图1 Overview of the SingleCellExperiment class object

4.2.1 Primary Data:The `assays` Slot

SingleCellExperimentsce)对象是基于在Bioconductor的单细胞分析应用的基础。该sce对象是一个S4对象(https://adv-r.hadley.nz/s4.html),与R中其他可用的方法相比,它本质上为数据的构造和访问提供了一种更为形式化的方法。S4对象最大最重要的用处就是可以安全有效检查数据并且能够扩展数据内容。

假设的sce是一艘船,sceslot可被认为是货箱——可以作为sce中独立的实体。并且每个slot都以自己的格式保存数据。可以理解为水果和砖头需要不同的货箱。在sce中,某些slot需要数字矩阵,而其他slot可能需要数据框。

要构建基本sce对象,我们需要的是一个slot:

  • assays slot:包含基本数据,例如列表中的counts,列表中的每个条目均采用矩阵格式,其中行对应于特征(基因),列对应于样本(cell)(图1A,蓝色框

让我们从生成3个细胞,10个基因的简单数据开始:

代码语言:javascript复制
counts_matrix <- data.frame(cell_1 = rpois(10, 10), 
                    cell_2 = rpois(10, 10), 
                    cell_3 = rpois(10, 30))
rownames(counts_matrix) <- paste0("gene_", 1:10)
counts_matrix <- as.matrix(counts_matrix) # must be a matrix object!

现在可以构建我们的第一个SingleCellExperiment对象,使用已定义的构造函数SingleCellExperiment()。请注意,我们以命名列表的形式提供数据,并且列表的每个条目都是一个矩阵。在这里,我们将counts_matrix命名为counts

代码语言:javascript复制
sce <- SingleCellExperiment(assays = list(counts = counts_matrix))

要查看对象,我们可以简单地在控制台输入sce以查看一些相关信息,这些信息将显示可供我们使用的各种slot的概述(可能有也可能没有数据)。

代码语言:javascript复制
sce
## class: SingleCellExperiment 
## dim: 10 3 
## metadata(0):
## assays(1): counts
## rownames(10): gene_1 gene_2 ... gene_9 gene_10
## rowData names(0):
## colnames(3): cell_1 cell_2 cell_3
## colData names(0):
## reducedDimNames(0):
## spikeNames(0):
## altExpNames(0):

要访问我们刚刚提供的计数数据,我们可以执行以下任一操作:

  • assay(sce, "counts") 这是最通用的方法,我们可以提供测定的名称作为第二个参数。
  • counts(sce)与上述相同,但适用于"counts"
代码语言:javascript复制
counts(sce)
##         cell_1 cell_2 cell_3
## gene_1      10      5     37
## gene_2       7     11     30
## gene_3       5     12     34
## gene_4      12     10     27
## gene_5      12     11     32
## gene_6      12      8     28
## gene_7       7     13     22
## gene_8      16     11     32
## gene_9      17     10     29
## gene_10     10      8     29
代码语言:javascript复制
## assay(sce, "counts") ## same as above in this special case

4.2.2 Extending the `assays ` Slot

使assay slot特别强大的原因是它可以容纳主要数据的多种表示形式。这对于存储数据的原始版本和规范化版本特别有用。我们可以使用scranscater包来计算初始主数据的规范化和用对数转换表示初始数据。

请注意,您需要将结果附加到sce对象上。

代码语言:javascript复制
sce <- scran::computeSumFactors(sce)
sce <- scater::normalize(sce)

再次查看该对象,我们看到这些函数添加了一些新条目:

代码语言:javascript复制
sce
## class: SingleCellExperiment 
## dim: 10 3 
## metadata(1): log.exprs.offset
## assays(2): counts logcounts
## rownames(10): gene_1 gene_2 ... gene_9 gene_10
## rowData names(0):
## colnames(3): cell_1 cell_2 cell_3
## colData names(0):
## reducedDimNames(0):
## spikeNames(0):
## altExpNames(0):

具体来说,我们看到该assays slot 已由两个条目组成:counts(初始数据)和logcounts(规范化数据)。与counts相似,该logcounts名称也是一个可以直接访问的特殊名称logcounts(sce)

代码语言:javascript复制
logcounts(sce)
##         cell_1 cell_2 cell_3
## gene_1    4.06   3.25   4.45
## gene_2    3.58   4.31   4.16
## gene_3    3.14   4.43   4.33
## gene_4    4.31   4.18   4.02
## gene_5    4.31   4.31   4.25
## gene_6    4.31   3.87   4.07
## gene_7    3.58   4.54   3.74
## gene_8    4.70   4.31   4.25
## gene_9    4.79   4.18   4.12
## gene_10   4.06   3.87   4.12
代码语言:javascript复制
## assay(sce, "logcounts") ## same as above

要查看中所有可用的检测方法sce,我们可以输入:

代码语言:javascript复制
assays(sce)
## List of length 2
## names(2): counts logcounts

虽然上面的功能演示了sceassays可以自动添加,但是在某些情况下,我们可能希望执行自己的计算并将结果保存到assaysslot中。让我们附加偏移的数据,在原来的基础上 100

代码语言:javascript复制
counts_100 <- assay(sce, "counts")   100
assay(sce, "counts_100") <- counts_100 # assign a new entry to assays slot

然后,我们可以使用访问器assays()(请注意,这是复数!)来查看到目前为止assay在slot中的所有条目。请注意,要查看所有测定,我们使用复数assays()访问器,并使用单数访问器assay()检索单个测定条目(作为矩阵),并提供想要检索的测定名称。

代码语言:javascript复制
assays(sce)
## List of length 3
## names(3): counts logcounts counts_100

这些条目也可以在默认视图中看到sce

代码语言:javascript复制
sce
## class: SingleCellExperiment 
## dim: 10 3 
## metadata(1): log.exprs.offset
## assays(3): counts logcounts counts_100
## rownames(10): gene_1 gene_2 ... gene_9 gene_10
## rowData names(0):
## colnames(3): cell_1 cell_2 cell_3
## colData names(0):
## reducedDimNames(0):
## spikeNames(0):
## altExpNames(0):

图1B(深蓝色框)assays以图形方式表示了slot的这种扩展,显示了将logcounts矩阵添加到assays slot中。

以类似的方式,sce许多slot可通过分配而扩展,从而允许与面向单细胞的生物导体封装之外的功能互操作所需的多种定制功能。

4.2.3 Column (Meta)Data: `colData` Slot

为了进一步注释我们的sce对象,第一个也是最有用的信息是添加描述主要数据列的元数据,例如实验的样本或cell 。此数据输入到colData slot:

  • colData slot:提供描述样本(cell)的元数据,以data.frame(或DataFrame)形式保存,其中行对应于cell,列对应于样本(cell)元数据功能(例如id,批处理,作者等)(图1A ,橙色框)。

因此,让我们为cell提供一些元数据,从批处理变量开始,其中cell1和2在批处理1中,而cell3在批处理2中。

代码语言:javascript复制
cell_metadata <- data.frame(batch = c(1, 1, 2))
rownames(cell_metadata) <- paste0("cell_", 1:3)

现在,我们可以采用两种方法将cell_metadata追加到我们现有的sce,或者通过SingleCellExperiment()构造函数从头开始并从一开始就提供。现在,我们将从头开始,但还将展示如何添加数据:

代码语言:javascript复制
## From scratch:
sce <- SingleCellExperiment(assays = list(counts = counts_matrix),
                           colData = cell_metadata)

## Appending to existing object (requires DataFrame() coercion)
## colData(sce) <- DataFrame(cell_metadata)

与相似assays,我们可以看到colData默认填充到sce

代码语言:javascript复制
sce
## class: SingleCellExperiment 
## dim: 10 3 
## metadata(0):
## assays(1): counts
## rownames(10): gene_1 gene_2 ... gene_9 gene_10
## rowData names(0):
## colnames(3): cell_1 cell_2 cell_3
## colData names(1): batch
## reducedDimNames(0):
## spikeNames(0):
## altExpNames(0):

通过colData()访问数据:

代码语言:javascript复制
colData(sce)
## DataFrame with 3 rows and 1 column
##            batch
##        <numeric>
## cell_1         1
## cell_2         1
## cell_3         2

最后,一些软件包会自动添加到colDataslot中,例如,该scater中的calculateQCMetrics()函数,该函数会附加许多质量控制数据。在这里,我们显示附加了质量控制指标的前五列。

代码语言:javascript复制
sce <- scater::calculateQCMetrics(sce)
colData(sce)[, 1:5]
## DataFrame with 3 rows and 5 columns
##            batch is_cell_control total_features_by_counts log10_total_features_by_counts
##        <numeric>       <logical>                <integer>                      <numeric>
## cell_1         1           FALSE                       10               1.04139268515822
## cell_2         1           FALSE                       10               1.04139268515822
## cell_3         2           FALSE                       10               1.04139268515822
##        total_counts
##           <integer>
## cell_1          108
## cell_2           99
## cell_3          300
4.2.3.1 Using `colData` for Subsetting

colData最常用的功能是提取子集。一种简单的访问colData的方式是通过使用$运算符,这是访问colDataslot内变量的捷径:

代码语言:javascript复制
sce$batch
## [1] 1 1 2
代码语言:javascript复制
## colData(sce)$batch # same as above

如果只需要批次1中的cell,则可以按如下方式对sce对象进行取子集(我们可以在上取子集,是因为在这里是通过cell/样本进行过滤的)。

代码语言:javascript复制
sce[, sce$batch == 1]
## class: SingleCellExperiment 
## dim: 10 2 
## metadata(0):
## assays(1): counts
## rownames(10): gene_1 gene_2 ... gene_9 gene_10
## rowData names(7): is_feature_control mean_counts ... total_counts log10_total_counts
## colnames(2): cell_1 cell_2
## colData names(10): batch is_cell_control ... pct_counts_in_top_200_features
##   pct_counts_in_top_500_features
## reducedDimNames(0):
## spikeNames(0):
## altExpNames(0):

4.2.4 Feature Metadata: `rowData`/`rowRanges`

最后,这些行还具有自己的元数据slot,用于存储与sce对象功能有关的信息:

  • rowData slot:包含data.frameDataFrame)格式的数据,该数据描述与主要数据行相对应的方面(图1A,绿色框)。

此外,还有一个特殊的slot,适用于具有基因组坐标的特征:

  • rowRanges slot:以GRangesList(其中每个条目均为GenomicRanges格式)的数据形式描述染色体,开始和结束坐标,特征(基因,基因组区域)等。

这两种方法都能通过各自的存取访问,rowRanges()rowData()。在我们的例子中,rowRanges(sce)产生一个空列表:

代码语言:javascript复制
rowRanges(sce) # empty
## GRangesList object of length 10:
## $gene_1
## GRanges object with 0 ranges and 0 metadata columns:
##    seqnames    ranges strand
##       <Rle> <IRanges>  <Rle>
##   -------
##   seqinfo: no sequences
## 
## $gene_2
## GRanges object with 0 ranges and 0 metadata columns:
##    seqnames    ranges strand
##       <Rle> <IRanges>  <Rle>
##   -------
##   seqinfo: no sequences
## 
## $gene_3
## GRanges object with 0 ranges and 0 metadata columns:
##    seqnames    ranges strand
##       <Rle> <IRanges>  <Rle>
##   -------
##   seqinfo: no sequences
## 
## ...
## <7 more elements>

上一节中调用calculateQCMetrics(sce)的结果保存在sce对象的rowDataslot中,如下所示(为简便起见,仅显示了前三列):

代码语言:javascript复制
rowData(sce)[, 1:3]
## DataFrame with 10 rows and 3 columns
##         is_feature_control      mean_counts log10_mean_counts
##                  <logical>        <numeric>         <numeric>
## gene_1               FALSE 17.3333333333333  1.26324143477458
## gene_2               FALSE               16  1.23044892137827
## gene_3               FALSE               17  1.25527250510331
## gene_4               FALSE 16.3333333333333  1.23888208891514
## gene_5               FALSE 18.3333333333333  1.28630673884327
## gene_6               FALSE               16  1.23044892137827
## gene_7               FALSE               14  1.17609125905568
## gene_8               FALSE 19.6666666666667  1.31527043477859
## gene_9               FALSE 18.6666666666667  1.29373075692248
## gene_10              FALSE 15.6666666666667  1.22184874961636

以类似于colDataslot的方式,可以在创建SingleCellExperiment对象时在开始时提供此类功能元数据,我们将其留给读者作为练习。

4.2.4.1 Subsetting by Rows

sce在特征/基因级别上将对象细分,我们可以通过提供数字索引或名称向量来进行类似于其他R对象的行设置操作:

代码语言:javascript复制
sce[c("gene_1", "gene_4"), ]
## class: SingleCellExperiment 
## dim: 2 3 
## metadata(0):
## assays(1): counts
## rownames(2): gene_1 gene_4
## rowData names(7): is_feature_control mean_counts ... total_counts log10_total_counts
## colnames(3): cell_1 cell_2 cell_3
## colData names(10): batch is_cell_control ... pct_counts_in_top_200_features
##   pct_counts_in_top_500_features
## reducedDimNames(0):
## spikeNames(0):
## altExpNames(0):
代码语言:javascript复制
## sce[c(1, 4), ] # same as above in this case

4.2.5 Size Factors Slot: `sizeFactors`

简而言之,我们已经通过调用scran::computeSumFactors(sce)(添加了一个sizeFactors slot)

  • sizeFactorsslot:在数字向量中包含有关样本/细胞归一化因子的信息,该因子用于产生归一化的数据(图1B,棕色框
代码语言:javascript复制
sce <- scran::computeSumFactors(sce)
sce <- scater::normalize(sce)
sizeFactors(sce)
## [1] 0.639 0.586 1.775

4.3

A Brief Recap: From `se` to `sce`

到目前为止,我们已经涵盖了assays(原始数据),colData(样本元数据),rowData/ rowRanges(特征元数据)和SingleCellExperiment的sizeFactorsslot。

需要注意的是,SingleCellExperiment类是从SummarizedExperimentse)类导出,并继承assayscolDatarowData/ rowRangesslot。大多数SummarizedExperiment功能保留在SingleCellExperiment中。这使与之配合使用的现有方法在SingleCellExperiment对象上类似地工作。

那么SingleCellExperiment又有什么新消息呢?对于我们的讨论,最重要的更改是添加了一个名为reducedDims的新slot。

4.4

The `reducedDims` Slot

reducedDimsslot是一个新增功能,专门用于存储通过PCA,tSNE,UMAP等方法获得的原始数据的降维表示。

  • reducedDimsslot:包含数字matrix条目的列表,这些条目描述了降维的原始数据表示,因此行代表原始数据的列(样本/cell),而列则代表维度

最重要的是,就像assaysslot一样,该reducedDimsslot可以容纳许多条目的列表。因此,它可以在reducedDimsslot内保存给定数据集的PCA,TSNE和UMAP表示形式。

在我们的示例中,我们可以使用scater中的函数runPCA()计算数据的PCA。我们看到,sce现在显示了一个新的reducedDim,并且通过访问器reducedDim()对来自的规范化数据logcounts(sce)产生了运行PCA的结果。

代码语言:javascript复制
sce <- scater::runPCA(sce)
reducedDim(sce, "PCA")
##           PC1    PC2
## cell_1 -1.126 -0.309
## cell_2  0.991 -0.455
## cell_3  0.135  0.764
## attr(,"percentVar")
## [1] 71.9 28.1

由此,我们还可以使用scater函数runTSNE()来计算tSNE表示形式,并且可以在默认视图中sce 看到它:

代码语言:javascript复制
sce <- scater::runTSNE(sce, perplexity = 0.1)
## Perplexity should be lower than K!
代码语言:javascript复制
reducedDim(sce, "TSNE")
##         [,1]  [,2]
## cell_1  1231 -5557
## cell_2  4197  3859
## cell_3 -5428  1699

我们可以通过访问器reducedDims()查看reducedDimsslot中所有条目的名称(请注意,该名称是复数形式,与reducedDim()不同:

代码语言:javascript复制
reducedDims(sce)
## List of length 2
## names(2): PCA TSNE

现在,假设我们要尝试使用一种不同的维算法,该算法尚未实现与SingleCellExperiment的直接兼容性,例如,使用scater。为了适应这种情况(或者,当我们想直接对数据本身而不是通过包装器运行降维方法时),我们可以直接添加到reducedDimsslot中。这类似于我们之前assays使用自定义条目counts_100扩展slot的方式。

下面,我们展示如何直接通过uwot包中的umap(),而不是通过scater关联包装函数runUMAP()来运行,并保存中间结果,然后将它们添加到sce我们先前拥有的对象中。

代码语言:javascript复制
u <- uwot::umap(t(logcounts(sce)), n_neighbors = 2)
reducedDim(sce, "UMAP_uwot") <- u

reducedDim(sce, "UMAP_uwot")
##          [,1]   [,2]
## cell_1  0.526 -0.402
## cell_2 -0.266  0.632
## cell_3 -0.260 -0.230
## attr(,"scaled:center")
## [1]  9.12 15.20

现在,当我们查看访问reducedDims()器输出时,我们还可以看到其条目:

代码语言:javascript复制
reducedDims(sce)
## List of length 3
## names(3): PCA TSNE UMAP_uwot

4.5

The `metadata` Slot

一些分析产生的结果不适合上述slot。值得庆幸的是,有一个slot仅用于这种类型的混乱数据,实际上它可以容纳任何类型的数据,只要它在命名列表中即可:

  • metadata slot:一个命名的条目列表,列表中的每个条目都可以是您想要的任何内容

例如,假设我们有一些喜欢的基因,例如高度可变的基因,我们希望将其保存在内部sce以便以后使用。我们可以简单地通过如下所示将其附加到元数据slot来做到这一点:

代码语言:javascript复制
my_genes <- c("gene_1", "gene_5")
metadata(sce) <- list(favorite_genes = my_genes)
metadata(sce)
## $favorite_genes
## [1] "gene_1" "gene_5"

同样,我们可以通过$运算符附加更多信息:

代码语言:javascript复制
your_genes <- c("gene_4", "gene_8")
metadata(sce)$your_genes <- your_genes
metadata(sce)
## $favorite_genes
## [1] "gene_1" "gene_5"
## 
## $your_genes
## [1] "gene_4" "gene_8"

4.6

About Spike-Ins

您可能已经注意到sce默认视图中带有spikeNamesSingleCellExperiment对象对通过掺入(ERCC)对照 的实验数据预留了通道。我们将此留给有兴趣的读者,以在SingleCellExperiment入门小插图中(https://bioconductor.org/packages/release/bioc/vignettes/SingleCellExperiment/inst/doc/intro.html)了解更多信息。

4.7

Recording Analyses in SingleCellExperiment

在随后的部分中,我们将显示一个示例工作流,该工作流使用SingleCellExperiment对象作为其基础,并且与上面的SingleCellExperiment类的演练类似,将连续追加新条目以保存分析结果。因此,SingleCellExperiment可以以此方式用作分析记录。这使得它对于协作特别有用,因为可以通过图形用户界面(例如iSEE)来传输,分析甚至可视化对象。

4.8

结论

Bioc-verse单细胞相关的软件包,与SingleCellExperiment class可以连接。正是这种与SingleCellExperiment的连接,使得许多这些程序包在scRNA-seq分析的整个过程中都易于互操作和模块化。此外,它允许任何人以SingleCellExperiment为基础,建立自己的scRNA-SEQ分析工具,促进发展。

0 人点赞