查看: 10187|回复: 3

[R语言] sc-RAN-seq 数据分析||Seurat新版教程: Integration and Label Transfer

[复制链接]

迅猛龙

Rank: 8Rank: 8

主题
161
注册时间
2015.11.14
在线时间
78 小时

发表于 2019.9.26 16:37:22 | 显示全部楼层 |阅读模式
如Stuart, ButlerComprehensive Integration of Single-Cell Data所述。Seurat v3引入了集成多个单细胞数据集的新方法。这些方法的目的是识别存在于不同数据集中的共享细胞状态(shared cell states),即使它们是从不同的个体、实验条件、技术甚至物种中收集来的。
我们的方法旨在首先识别数据集对之间的“锚(anchors)”。这些代表了个体细胞之间成对的对应关系(每个数据集中有一个),我们假设它们源自相同的生物状态。然后,这些“锚”用于协调数据集,或将信息从一个数据集传输到另一个数据集。下面,我们将演示集成分析的多种应用,并介绍2019年文章的中描述的之外的新功能(cell,Comprehensive Integration of Single-Cell Data)。为了帮助指引用户,我们简要介绍以下功能:
  • 标准工作流程
描述标准Seurat v3集成工作流,并将其应用于集成(跨不同技术)收集的多个人类胰岛数据集。我们还演示了如何使用Seurat v3作为分类器,将集群标签传输到新收集的数据集中。我们向新用户推荐这个。
  • SCTransform
描述v3集成工作流的一个修改,以便应用于使用我们的新规范化方法SCTransform进行规范化的数据集。我们将此方法应用于与前面描述的相同的胰岛数据集,并集成来自8种不同技术的人类PBMC数据集 eight different technologies,作为人类细胞图谱的系统技术基准。
我们向熟悉SCTransform归一化方法的高级用户推荐这个方法。您可以在我们最近的预印本(preprint)中阅读有关SCTransform的更多信息,并了解如何将其应用于单独的方法中的单个数据集( vignette)。
  • Reference-based (基于参照)
描述v3集成工作流的修改,其中将数据集的子集(或单个数据集)作为“参考”(reference)列出。这种方法可以显著提高速度,特别是在需要集成大量数据集的情况下。我们将此方法应用于上面描述的8个PBMC数据集,并观察到相同的结果,尽管处理时间大大减少。
对于正在集成许多数据集并希望提高速度的用户,我们推荐使用这个方法。
  • Reciprocal PCA
    描述v3集成工作流的一个修改,其中使用互反PCA(reciprocal PCA)代替正则相关分析来减少锚点查找中使用的维度。这种方法可以提高处理大型数据集的速度和效率。
    我们建议用户在处理大量数据集或细胞时,如使用许多实验条件、复制或患者的实验设计,寻找速度/内存方面的改进。然而,这个工作流可能很难对齐高度不同的样本(例如跨物种,或跨模式,集成)。对于“turbo”模式,可以考虑结合“基于参照”的方式整合,如下所示。


标准工作流程
在这个例子工作流中,我们演示了我们最近在论文中引入的两种新方法,Comprehensive Integration of Single Cell Data:
  • 将多个不同的scna -seq数据集组装到一个参考数据集中
  • 将细胞类型标签从参考数据集转移到新的查询数据集
在本例中,我们选择了通过四种技术(CelSeq (GSE81076) CelSeq2 (GSE85241)、Fluidigm C1 (GSE86469)和SMART-Seq2 (E-MTAB-5061)生成的人类胰岛细胞数据集。为了方便起见,我们通过SeuratData包分发这个数据集。
新方法的代码在Seurat v3中实现。您可以使用install.packages从CRAN下载和安装。
除了新的方法之外,Seurat v3还包含了许多旨在改进Seurat对象和用户交互的改进。为了帮助用户熟悉这些更改,我们为常见任务编写了一个命令备忘单( command cheat sheet)。
载入包和数据:
[AppleScript] 纯文本查看 复制代码
library(Seurat)
packageVersion("Seurat")

>[1] ‘3.1.0’

library(SeuratData)
#InstallData("pbmc3k")
#InstallData("panc8")
#install.packages("C:/Users/Administrator/Desktop/panc8.SeuratData_3.0.2.tar.gz", repos = NULL, type = "source")
#install.packages("C:/Users/Administrator/Desktop/pbmcsca.SeuratData_3.0.0.tar.gz", repos = NULL, type = "source")
[AppleScript] 纯文本查看 复制代码
> library(panc8.SeuratData)
> data("panc8")
> panc8
An object of class Seurat 
34363 features across 14890 samples within 1 assay 
Active assay: RNA (34363 features)
> summary(panc8)
Length  Class   Mode 
     1 Seurat     S4 
> citation('panc8.SeuratData')

在出版物中使用程序包时引用‘panc8.SeuratData’:

  Satija Lab (2019). panc8.SeuratData: Eight Pancreas
  Datasets Across Five Technologies. R package version
  3.0.2.

A BibTeX entry for LaTeX users is

  @Manual{,
    title = {panc8.SeuratData: Eight Pancreas Datasets Across Five Technologies},
    author = {Satija Lab},
    year = {2019},
    note = {R package version 3.0.2},
  }

Warning message:
In citation("panc8.SeuratData") :
  程序包‘panc8.SeuratData’里的DESCRIPTION文件中没有日期这一域
要构造参考数据集,我们将在各个数据集之间标识“锚”。首先,我们将组合的对象拆分为一个列表,每个数据集作为一个元素。
[AppleScript] 纯文本查看 复制代码
> pancreas.list <- SplitObject(panc8, split.by = "tech")
> summary(pancreas.list)
           Length Class  Mode
celseq     1      Seurat S4  
celseq2    1      Seurat S4  
smartseq2  1      Seurat S4  
fluidigmc1 1      Seurat S4  
indrop     1      Seurat S4  
在找到锚之前,我们执行标准的预处理(log-normalization),并为每个锚分别识别变量特性。注意,Seurat v3实现了一种改进的基于方差稳定转换(“vst”)的变量特征选择方法。
[AppleScript] 纯文本查看 复制代码
pancreas.list <- pancreas.list[c("celseq", "celseq2", "fluidigmc1", "smartseq2")]
for (i in 1:length(pancreas.list)) {
  pancreas.list[[i]] <- NormalizeData(pancreas.list[[i]], verbose = FALSE)
  pancreas.list[[i]] <- FindVariableFeatures(pancreas.list[[i]], selection.method = "vst",  nfeatures = 2000, verbose = FALSE)
}
整合3个胰岛细胞数据集
接下来,我们使用findintegrationanchor函数来标识锚,该函数接受Seurat对象的列表(list)作为输入。在这里,我们将三个对象集成到一个参考数据集中(稍后我们将使用第四个对象)。
[AppleScript] 纯文本查看 复制代码
> reference.list <- pancreas.list[c("celseq", "celseq2", "smartseq2")]
> pancreas.anchors <- FindIntegrationAnchors(object.list = reference.list, dims = 1:30)
Computing 2000 integration features
Scaling features for provided objects
  |++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed = 02s
Finding all pairwise anchors
  |                                                  | 0 % ~calculating  Running CCA
Merging objects
Finding neighborhoods
Finding anchors
    Found 3514 anchors
Filtering anchors
    Retained 2761 anchors
Extracting within-dataset neighbors
  |+++++++++++++++++                                 | 33% ~25s          Running CCA
Merging objects
Finding neighborhoods
Finding anchors
    Found 3500 anchors
Filtering anchors
    Retained 2728 anchors
Extracting within-dataset neighbors
  |++++++++++++++++++++++++++++++++++                | 67% ~12s          Running CCA
Merging objects
Finding neighborhoods
Finding anchors
    Found 6174 anchors
Filtering anchors
    Retained 4561 anchors
Extracting within-dataset neighbors
  |++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed = 49s
我们在这里使用所有默认参数来标识锚,包括数据集的“维数”(30;您可以尝试在较大范围内更改此参数,例如在10到50之间)。
然后我们将这些锚传递给IntegrateData函数,该函数返回一个Seurat对象。
返回的对象将包含一个新的Assay,其中包含一个针对所有细胞的完整(或“批量校正”)表达矩阵,使它们能够被联合分析。
[AppleScript] 纯文本查看 复制代码
> pancreas.integrated <- IntegrateData(anchorset = pancreas.anchors, dims = 1:30)
Merging dataset 1 into 2
Extracting anchors for merged samples
Finding integration vectors
Finding integration vector weights
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Integrating data
Merging dataset 3 into 2 1
Extracting anchors for merged samples
Finding integration vectors
Finding integration vector weights
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Integrating data
[AppleScript] 纯文本查看 复制代码
> head([url=mailto:pancreas.integrated@meta.data]pancreas.integrated@meta.data[/url])
        orig.ident nCount_RNA nFeature_RNA   tech replicate
D101_5        D101   4615.810         1986 celseq    celseq
D101_7        D101  29001.563         4209 celseq    celseq
D101_10       D101   6707.857         2408 celseq    celseq
D101_13       D101   8797.224         2964 celseq    celseq
D101_14       D101   5032.558         2264 celseq    celseq
D101_17       D101  13474.866         3982 celseq    celseq
        assigned_cluster celltype dataset
D101_5              <NA>    gamma  celseq
D101_7              <NA>   acinar  celseq
D101_10             <NA>    alpha  celseq
D101_13             <NA>    delta  celseq
D101_14             <NA>     beta  celseq
D101_17             <NA>   ductal  celseq
运行IntegrateData后,Seurat对象将包含一个使用集成表达矩阵的Assay 。注意,原始值(未更正的值)仍然存储在“RNA”测试中的对象中,因此您可以来回切换。
然后我们可以使用这个新的矩阵进行下游分析和可视化。在这里,我们缩放集成的数据,运行PCA,并使用UMAP可视化结果。集成的数据集按细胞类型聚类,而不是按技术分群。
[AppleScript] 纯文本查看 复制代码
library(ggplot2)
library(cowplot)
# switch to integrated assay. The variable features of this assay are automatically
# set during IntegrateData
DefaultAssay(pancreas.integrated) <- "integrated"

# Run the standard workflow for visualization and clustering
?ScaleData
pancreas.integrated <- ScaleData(pancreas.integrated, verbose = FALSE)

pancreas.integrated <- RunPCA(pancreas.integrated, npcs = 30, verbose = FALSE)
pancreas.integrated <- RunUMAP(pancreas.integrated, reduction = "pca", dims = 1:30)
[AppleScript] 纯文本查看 复制代码
Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
This message will be shown once per session
20:55:25 Read 5683 rows and found 30 numeric columns
20:55:25 Using Annoy for neighbor search, n_neighbors = 30
20:55:26 Building Annoy index with metric = cosine, n_trees = 50
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
20:55:28 Writing NN index file to temp file C:\Users\ADMINI~1\AppData\Local\Temp\Rtmpo5X81p\file2870430b3367
20:55:28 Searching Annoy index using 1 thread, search_k = 3000
20:55:30 Annoy recall = 100%
20:55:31 Commencing smooth kNN distance calibration using 1 thread
20:55:31 Found 2 connected components, falling back to 'spca' initialization with init_sdev = 1
20:55:31 Initializing from PCA
20:55:31 PCA: 2 components explained 44.22% variance
20:55:32 Commencing optimization for 500 epochs, with 252460 positive edges
0%   10   20   30   40   50   60   70   80   90   100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
20:55:59 Optimization finished
[AppleScript] 纯文本查看 复制代码
p1 <- DimPlot(pancreas.integrated, reduction = "umap", group.by = "tech")
p2 <- DimPlot(pancreas.integrated, reduction = "umap", group.by = "celltype", label = TRUE, 
              repel = TRUE) + NoLegend()
plot_grid(p1, p2)
使用参考数据集对细胞类型进行分类
Seurat v3还支持将参考数据(或元数据)投影到查询对象上。虽然许多方法是一致的(这两个过程都是从识别锚开始的),但数据映射(data transfer)和整合之间有两个重要的区别:
  • 在数据映射中,Seurat不纠正或修改查询表达式数据。
  • 在数据映射中,Seurat有一个选项(默认设置)来将参考的PCA结构投射到查询上,而不是学习与CCA的联合结构。我们通常建议在scna -seq数据集之间投影数据时使用此选项。
找到锚之后,我们使用TransferData函数根据参考数据(参考数据的细胞类型标签向量)对查询数据集的细胞进行分类。TransferData返回一个带有预测id和预测分数的矩阵,我们可以将其添加到查询元数据中(query metadata)
[AppleScript] 纯文本查看 复制代码
> pancreas.query <- pancreas.list[["fluidigmc1"]]
> pancreas.anchors <- FindTransferAnchors(reference = pancreas.integrated, query = pancreas.query, 
+ dims = 1:30)
Performing PCA on the provided reference using 2000 features as input.
Projecting PCA
Finding neighborhoods
Finding anchors
Found 953 anchors
Filtering anchors
Retained 885 anchors
Extracting within-dataset neighbors
> predictions <- TransferData(anchorset = pancreas.anchors, refdata = pancreas.integrated$celltype, 
+ dims = 1:30)
Finding integration vectors
Finding integration vector weights
0% 10 20 30 40 50 60 70 80 90 100%
[----|----|----|----|----|----|----|----|----|----|
**************************************************|
Predicting cell labels
> pancreas.query <- AddMetaData(pancreas.query, metadata = predictions)
因为我们有来自完整集成分析的原始标签注释,所以我们可以评估预测的细胞类型注释与完整引用的匹配程度。在这个例子中,我们发现在细胞类型分类上有很高的一致性,超过97%的细胞被正确标记。
[AppleScript] 纯文本查看 复制代码
pancreas.query$prediction.match <- pancreas.query$predicted.id == pancreas.query$celltype
table(pancreas.query$prediction.match)

FALSE TRUE 
18 620 
为了进一步验证这一点,我们可以检查特定胰岛细胞群的一些典型细胞类型标记(cell type markers)。请注意,即使这些细胞类型中的一些只由一两个细胞(例如epsilon细胞),我们仍然能够正确地对它们进行分类。
[AppleScript] 纯文本查看 复制代码
> table(pancreas.query$predicted.id)

acinar activated_stellate alpha beta delta 
21 17 248 258 22 
ductal endothelial epsilon gamma macrophage 
33 13 1 17 1 
mast schwann 
2 5 
> VlnPlot(pancreas.query, c("REG1A", "PPY", "SST", "GHRL", "VWF", "SOX10"), group.by = "predicted.id")


本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
我的意中人是一个盖世英雄,有一天他会踩着七彩祥云来娶我……
回复

使用道具 举报

迅猛龙

Rank: 8Rank: 8

主题
161
注册时间
2015.11.14
在线时间
78 小时

 楼主| 发表于 2019.9.26 16:37:23 | 显示全部楼层
本帖最后由 犬夜× 于 2019.9.26 16:46 编辑

SCTransform
在前面的选项卡中,我们演示了如何在使用标准规范化(log-normalization)对每个数据集进行预处理之后集成数据集。在这里,我们修改工作流,以利用改进的预处理和规范化工作流:SCTransform。您可以在我们最近的预印本中阅读有关SCTransform的更多信息,并了解如何将其应用于单个数据集。
  • 创建要集成的Seurat对象列表
  • 对每个数据集分别执行SCTransform规范化
  • 在对象列表上运行PrepSCTIntegration函数
  • 集成数据集,并进行联合分析
首先,设置Seurat对象列表,分别对每个对象运行SCTransform:
[AppleScript] 纯文本查看 复制代码
m(list=ls())
library(Seurat)
library(SeuratData)
library(ggplot2)
options(future.globals.maxSize = 4000 * 1024^2)
data("panc8")
pancreas.list <- SplitObject(panc8, split.by = "tech")
pancreas.list <- pancreas.list[c("celseq", "celseq2", "fluidigmc1", "smartseq2")]
for (i in 1:length(pancreas.list)) {
pancreas.list[[i]] <- SCTransform(pancreas.list[[i]], verbose = FALSE)
}
接下来,选择下游集成的特性,并运行PrepSCTIntegration,这将确保计算了所有必需的Pearson残差。
[AppleScript] 纯文本查看 复制代码
?SelectIntegrationFeatures
pancreas.features <- SelectIntegrationFeatures(object.list = pancreas.list, nfeatures = 3000)
pancreas.list <- PrepSCTIntegration(object.list = pancreas.list, anchor.features = pancreas.features, verbose = FALSE)
接下来,识别锚并集成数据集。命令与标准工作流相同,但请确保设置了normalization.method = 'SCT':
[AppleScript] 纯文本查看 复制代码
pancreas.anchors <- FindIntegrationAnchors(object.list = pancreas.list, normalization.method = "SCT", 
anchor.features = pancreas.features, verbose = FALSE)
pancreas.integrated <- IntegrateData(anchorset = pancreas.anchors, normalization.method = "SCT", 
verbose = FALSE)
现在继续对集成数据集进行下游分析(即可视化、聚类)。命令与标准工作流相同,但在集成之后不运行ScaleData函数。您可以看到,在集成之后,细胞按照它们的生物细胞类型(已经预先注释过)分组,而不是按照它们的技术分组。
[AppleScript] 纯文本查看 复制代码
pancreas.integrated <- RunPCA(pancreas.integrated, verbose = FALSE)
pancreas.integrated <- RunUMAP(pancreas.integrated, dims = 1:30)


plots <- DimPlot(pancreas.integrated, group.by = c("tech", "celltype"), combine = FALSE)
plots <- lapply(X = plots, FUN = function(x) x + theme(legend.position = "top") + guides(color = guide_legend(nrow = 3, 
byrow = TRUE, override.aes = list(size = 3))))
CombinePlots(plots)


为了进一步演示这个工作流,我们接下来将它应用于来自8种不同技术的一系列人类PBMC数据集,这些数据集由human Cell Atlas作为系统的技术基准生成。数据是单细胞门户下载的。为了方便起见,我们通过SeuratData包分发了这个数据集。
[AppleScript] 纯文本查看 复制代码
data("pbmcsca")
pbmc.list <- SplitObject(pbmcsca, split.by = "Method")
for (i in names(pbmc.list)) {
pbmc.list[[i]] <- SCTransform(pbmc.list[[i]], verbose = FALSE)
}
pbmc.features <- SelectIntegrationFeatures(object.list = pbmc.list, nfeatures = 3000)
pbmc.list <- PrepSCTIntegration(object.list = pbmc.list, anchor.features = pbmc.features)
pbmc.anchors <- FindIntegrationAnchors(object.list = pbmc.list, normalization.method = "SCT", 
anchor.features = pbmc.features)
pbmc.integrated <- IntegrateData(anchorset = pbmc.anchors, normalization.method = "SCT")

pbmc.integrated <- RunPCA(object = pbmc.integrated, verbose = FALSE)
pbmc.integrated <- RunUMAP(object = pbmc.integrated, dims = 1:30)
[AppleScript] 纯文本查看 复制代码
plots <- DimPlot(pbmc.integrated, group.by = c("Method", "CellType"), combine = FALSE)
plots <- lapply(X = plots, FUN = function(x) x + theme(legend.position = "top") + guides(color = guide_legend(nrow = 4, 
byrow = TRUE, override.aes = list(size = 2.5))))
CombinePlots(plots)
请再次注意,细胞是按照共享的生物亚型分组的,而不是按照技术分组的。在这种情况下,作者分别在每个数据集中标注“CellType”,你可以看到这与我们的综合分析是一致的:
[AppleScript] 纯文本查看 复制代码
DimPlot(pbmc.integrated, group.by = "CellType", split.by = "Method", ncol = 3)
然而,我们的综合分析揭示了这些广泛的单元类型之间更精细的细分,因为组合的数据集具有更大的统计能力。我们没有在集成的数据集上进行完整的聚类,但是请注意,我们可以观察CD4+ (Memory/Naive)、CD8(多个细胞毒性群体)和B细胞(多个发育阶段)的亚群。下面,我们在集成可视化上可视化原始度量的表达式(注意,您可以添加split.by = 'Method‘ 为每种技术单独绘制这些图)"
[AppleScript] 纯文本查看 复制代码
DefaultAssay(pbmc.integrated) <- "RNA"
# Normalize RNA data for visualization purposes
pbmc.integrated <- NormalizeData(pbmc.integrated, verbose = FALSE)
FeaturePlot(pbmc.integrated, c("CCR7", "S100A4", "GZMB", "GZMK", "GZMH", "TCL1A"))
Reference-based
在这里,我们对Seurat集成工作流进行了额外的修改,我们将其称为“基于参照”集成。在前面的工作流中,我们在所有对数据集之间标识锚。虽然这使得数据集在下游集成中具有同等的权重,但它也可能成为计算密集型的。例如,当集成10个不同的数据集时,我们执行45种不同的成对比较。
作为一种替代方法,我们在这里介绍了指定一个或多个数据集作为集成分析的“参照”的可能性,其余的指定为“查询”数据集。在此工作流中,我们不确定查询数据集对之间的锚,从而减少了比较的数量。例如,当将10个数据集与一个指定为参照的数据集集成时,我们只执行9个比较。基于参照的集成可以应用于log-normalized 或sctransform规范化数据集。
在实践中,我们经常在两种方法之间观察到非常相似的结果,但是计算时间大大减少。我们使用与前一个选项卡相同的PBMC分析来演示这一点,输出高度一致。因此,我们向集成大量数据集并希望提高计算效率的用户推荐这种工作流。
这种方法要求用户预先选择数据集作为参照。在本例中,我们使用了10X v3数据集,因为数据集同时包含高细胞数和高灵敏度。用户可以在选择引用时应用类似的标准。此外,您可以在下面的工作流中更改此设置(例如,尝试选择10X v2或Drop-seq数据集作为参考),在最终结果中只有非常小的差异。
工作流程如下所示。命令与前一个例子相同,只是我们在findintegrationanchor函数中指定了一个参照参数。数据最初是从 the Broad Single Cell Portal下载的,但是为了方便起见,我们通过SeuratData包分发该数据集。
[AppleScript] 纯文本查看 复制代码
library(Seurat)
library(SeuratData)
library(ggplot2)
InstallData("pbmcsca")
data("pbmcsca")
pbmc.list <- SplitObject(pbmcsca, split.by = "Method")
for (i in names(pbmc.list)) {
pbmc.list[[i]] <- SCTransform(pbmc.list[[i]], verbose = FALSE)
}
pbmc.features <- SelectIntegrationFeatures(object.list = pbmc.list, nfeatures = 3000)
pbmc.list <- PrepSCTIntegration(object.list = pbmc.list, anchor.features = pbmc.features)

# This command returns dataset 5. We can also specify multiple refs. (i.e. c(5,6))
reference_dataset <- which(names(pbmc.list) == "10x Chromium (v3)")
[AppleScript] 纯文本查看 复制代码
pbmc.anchors <- FindIntegrationAnchors(object.list = pbmc.list, normalization.method = "SCT", 
anchor.features = pbmc.features, reference = reference_dataset)
pbmc.integrated <- IntegrateData(anchorset = pbmc.anchors, normalization.method = "SCT")

pbmc.integrated <- RunPCA(object = pbmc.integrated, verbose = FALSE)
pbmc.integrated <- RunUMAP(object = pbmc.integrated, dims = 1:30)
[AppleScript] 纯文本查看 复制代码
plots <- DimPlot(pbmc.integrated, group.by = c("Method", "CellType"), combine = FALSE)
plots <- lapply(X = plots, FUN = function(x) x + theme(legend.position = "top") + guides(color = guide_legend(nrow = 4, 
byrow = TRUE, override.aes = list(size = 2.5))))
CombinePlots(plots)
[AppleScript] 纯文本查看 复制代码
DefaultAssay(pbmc.integrated) <- "RNA"
# Normalize RNA data for visualization purposes
pbmc.integrated <- NormalizeData(pbmc.integrated, verbose = FALSE)
FeaturePlot(pbmc.integrated, c("CCR7", "S100A4", "GZMB", "GZMK", "GZMH", "TCL1A"))
Reciprocal PCA
对于非常大的数据集,CCA(Canonical Correlation Analysis )有时可能成为计算上的卡点步骤。在这个工作流中,我们引入了另一种降维形式——reciprocal PCA,以确定寻找锚的有效空间。在使用reciprocal PCA确定任意两个数据集之间的锚点时,我们将每个数据集投影到另一个数据集的主元分析空间中,并根据相同的互邻域要求约束锚点。所有下游集成步骤保持不变,我们能够“纠正”(或协调)数据集。
如下所示,工作流程包括以下步骤:
  • 创建要集成的Seurat对象列表
  • 对每个数据集分别执行规范化、特征选择和缩放
  • 对列表中的每个对象运行PCA
  • 集成数据集,并进行联合分析
总的来说,我们在这两个工作流中观察到惊人的相似结果,但是当使用reciprocal PCA时,计算时间和内存显著减少。然而,如果数据集高度分散(例如,跨模态映射或跨物种映射),则只能使用一小部分特性来促进集成,并且您可以使用CCA观察到更好的结果。
对于包含许多数据集的大型研究,我们还建议将reciprocal PCA与基于引用的集成或SCTransform标准化相结合(请参阅前面的详细信息)。
对于这些例子,我们将使用来自人类细胞图谱的“免疫细胞图谱”数据,可以在这里找到。首先,我们将集成这个数据集的随机下采样版本,然后是整合的数据集。随机下采样的版本包含40K个细胞(每个供体5K个),可以通过我们新的SeuratData包获得。
[AppleScript] 纯文本查看 复制代码
library(Seurat)
library(SeuratData)
# Get data from seurat-data
InstallData(pkgs = "bm40k.SeuratData")
在获取数据后,首先进行标准归一化和变量特征选择。
[AppleScript] 纯文本查看 复制代码
library(ggplot2)
# set up future for parallelization
library(future)
library(future.apply)
plan("multiprocess", workers = 4)
options(future.globals.maxSize = 20000 * 1024^2)
# load in the data
data(hcabm40k)
bm40k.list <- SplitObject(hcabm40k, split.by = "orig.ident")
bm40k.list <- lapply(X = bm40k.list, FUN = function(x) {
x <- NormalizeData(x, verbose = FALSE)
x <- FindVariableFeatures(x, verbose = FALSE)
})
接下来,选择用于下游集成的特性,并在列表中的每个对象上运行PCA,这是运行对等的PCA工作流所必需的。
[AppleScript] 纯文本查看 复制代码
features <- SelectIntegrationFeatures(object.list = bm40k.list)
bm40k.list <- lapply(X = bm40k.list, FUN = function(x) {
x <- ScaleData(x, features = features, verbose = FALSE)
x <- RunPCA(x, features = features, verbose = FALSE)
})
由于这个数据集同时包含男性和女性,我们将选择一个男性和一个女性(BM1和BM2)在基于引用的工作流中使用。我们通过检测XIST基因的表达来确定供体性别。接下来,识别锚并集成数据集。命令与标准工作流相同,但请确保设置reduce = 'rpca':
[AppleScript] 纯文本查看 复制代码
anchors <- FindIntegrationAnchors(object.list = bm40k.list, reference = c(1, 2), reduction = "rpca", 
dims = 1:30)
bm40k.integrated <- IntegrateData(anchorset = anchors, dims = 1:30)
现在继续对集成数据集进行下游分析(即可视化、聚类)。命令与标准工作流相同。
[AppleScript] 纯文本查看 复制代码
bm40k.integrated <- ScaleData(bm40k.integrated, verbose = FALSE)
bm40k.integrated <- RunPCA(bm40k.integrated, verbose = FALSE)
bm40k.integrated <- RunUMAP(bm40k.integrated, dims = 1:30)
DimPlot(bm40k.integrated, group.by = "orig.ident")
为了进一步演示这个工作流,我们接下来将它应用于完整的数据集。虽然我们不通过SeuratData提供此功能,但您可以在HCA数据门户网站上找到counts矩阵。下面的命令与40k集成相同,只是我们将维度增加到了50,以反映增加的细胞数量和多样性。
[AppleScript] 纯文本查看 复制代码
bm280k.data <- Read10X_h5("../data/ica_bone_marrow_h5.h5")
bm280k <- CreateSeuratObject(counts = bm280k.data, min.cells = 100, min.features = 500)
[AppleScript] 纯文本查看 复制代码
bm280k.list <- SplitObject(bm280k, split.by = "orig.ident")
bm280k.list <- future_lapply(X = bm280k.list, FUN = function(x) {
x <- NormalizeData(x, verbose = FALSE)
x <- FindVariableFeatures(x, verbose = FALSE)
})
features <- SelectIntegrationFeatures(object.list = bm280k.list)
bm280k.list <- future_lapply(X = bm280k.list, FUN = function(x) {
x <- ScaleData(x, features = features, verbose = FALSE)
x <- RunPCA(x, features = features, verbose = FALSE)
})
[AppleScript] 纯文本查看 复制代码
anchors <- FindIntegrationAnchors(object.list = bm280k.list, reference = c(1, 2), reduction = "rpca", 
dims = 1:50)
bm280k.integrated <- IntegrateData(anchorset = anchors, dims = 1:50)
bm280k.integrated <- ScaleData(bm280k.integrated, verbose = FALSE)
bm280k.integrated <- RunPCA(bm280k.integrated, verbose = FALSE)
bm280k.integrated <- RunUMAP(bm280k.integrated, dims = 1:50)
DimPlot(bm280k.integrated, group.by = "orig.ident")


本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有帐号?立即注册

x
我的意中人是一个盖世英雄,有一天他会踩着七彩祥云来娶我……
回复 支持 反对

使用道具 举报

钵水母

Rank: 3Rank: 3

主题
4
注册时间
2017.4.23
在线时间
29 小时

发表于 2019.9.26 18:12:57 | 显示全部楼层
真的很受用
回复 支持 反对

使用道具 举报

帝王蝶

Rank: 4

主题
0
注册时间
2018.1.22
在线时间
127 小时

发表于 2019.9.27 10:15:49 | 显示全部楼层
lalalalal!
回复

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

快速回复 返回顶部 返回列表