Seurat新版教程：分析空间转录组数据（下）
分享是一种态度
作者 | 周运来
男，
一个长大了才会遇到的帅哥，
稳健，潇洒，大方，靠谱。
一段生信缘，一棵技能树，
一枚大型测序工厂的螺丝钉，
一个随机森林中提灯觅食的津门旅客。
空间变量特征的识别
Seurat提供了两个工作流程来识别与组织空间位置相关的分子特征。第一种是根据组织内预先标注的解剖区域进行差异表达，这种差异表达可以通过非监督聚类或先验知识来确定。这种策略在这种情况下有效，因为上面的集群显示出明显的空间差异。
1de_markers <- FindMarkers(brain, ident.1 = 4, ident.2 = 6)
2
3|++++++++++++++++++++++++++++++++++++++++++++++++++| 100% elapsed=37s
4
5SpatialFeaturePlot(object = brain, features = rownames(de_markers)[1:3], alpha = c(0.1, 1), ncol = 3)
6
在FindSpatiallyVariables中实现的另一种方法是在没有预注释的情况下搜索显示空间模式的特性。默认的方法(method = 'markvariogram ')受到 Trendsceek,的启发，后者将空间转录组学数据建模为标记点过程，并计算一个' variogram '，它识别其表达水平取决于其空间位置的基因。更具体地说，这个过程计算伽玛(r)值，测量两个点之间一定的“r”距离的相关性。默认情况下，我们在这些分析中使用的r值为‘5’，并且只计算可变基因的这些值(其中的变异是独立于空间位置计算的)，以节省时间。
现在，我们可视化的表达前6个特征确定了这一措施。
1top.features <- head(SpatiallyVariableFeatures(brain, selection.method = "markvariogram"), 6)
2SpatialFeaturePlot(brain, features = top.features, ncol = 3, alpha = c(0.1, 1))
3
可视化解剖区域的子集
与单细胞对象一样，您可以对该对象进行子集设置，以将重点放在数据的子集上。在这里，我们大致划分了额叶皮质。这个过程也促进了这些数据与下一节的皮层scRNA-seq数据集的整合。首先，我们取集群的一个子集，然后根据精确的位置进一步细分。设置好亚组后，我们可以在完整图像或裁剪图像上看到皮质细胞。
1cortex <- subset(brain, idents = c(1, 2, 3, 5, 6, 7))
2# now remove additional cells, use SpatialDimPlots to visualize what to remove
3# SpatialDimPlot(cortex,cells.highlight = WhichCells(cortex, expression = image_imagerow > 400 |
4# image_imagecol < 150))
5cortex <- subset(cortex, image_imagerow > 400 | image_imagecol < 150, invert = TRUE)
6cortex <- subset(cortex, image_imagerow > 275 & image_imagecol > 370, invert = TRUE)
7cortex <- subset(cortex, image_imagerow > 250 & image_imagecol > 440, invert = TRUE)
81p1 <- SpatialDimPlot(cortex, crop = TRUE, label = TRUE)
2p2 <- SpatialDimPlot(cortex, crop = FALSE, label = TRUE, pt.size.factor = 1, label.size = 3)
3plot_grid(p1, p2)
4
与单细胞数据关联分析（空间细胞类型定义）
在~50um时，visium检测的斑点将包含多个细胞的表达谱。对于越来越多可用scRNA-seq数据的系统，用户可能有兴趣“反卷积”每个空间体素，以预测单元类型的底层组成。在准备这篇文章的过程中，我们测试了各种各样的脱卵方法和整合方法（decovonlution and integration methods），使用的是来自Allen研究所的参考scRNA-seq数据集（reference scRNA-seq dataset），其中包含了大约14000只成年小鼠的皮层细胞分类，并使用SMART-Seq2 protocol 生成。
我们一致认为，使用集成方法(与反褶积方法相反)可以获得更好的性能，这可能是因为空间和单细胞数据集的噪声模型本质上是不同的，而集成方法的特殊设计是为了对这些差异具有鲁棒性。因此我们应用“锚”的集成工作流一如最近在 Seurat v3介绍的，使注释的概率从一个reference 数据query 数据集转移。因此,我们遵循这里介绍的标签转移流程,利用sctransform正常化,但预测新方法被开发来完成这项任务。
我们首先加载数据( here提供下载，我只能假装我下载成功了_)，预处理scRNA-seq 参考数据集，然后执行标签传输。该过程为每个spot输出每个scRNA-seq派生类的概率分类。我们将这些预测添加到Seurat对象中作为新的试验。
1allen_reference <- readRDS("~/Downloads/allen_cortex.rds")
2# note that setting ncells=3000 normalizes the full dataset but learns noise models on 3k cells
3# this speeds up SCTransform dramatically with no loss in performance
4library(dplyr)
5allen_reference <- SCTransform(allen_reference, ncells = 3000, verbose = FALSE) %>% RunPCA(verbose = FALSE) %>%
6 RunUMAP(dims = 1:30)
7
8# After subsetting, we renormalize cortex
9cortex <- SCTransform(cortex, assay = "Spatial", verbose = FALSE) %>% RunPCA(verbose = FALSE)
10# the annotation is stored in the 'subclass' column of object metadata
11DimPlot(allen_reference, group.by = "subclass", label = TRUE)
1anchors <- FindTransferAnchors(reference = allen_reference, query = cortex, normalization.method = "SCT")
2predictions.assay <- TransferData(anchorset = anchors, refdata = allen_reference$subclass, prediction.assay = TRUE,
3 weight.reduction = cortex[["pca"]])
4cortex[["predictions"]] <- predictions.assay现在，我们得到了每个班每个点的预测分数。在额叶皮层区域特别有趣的是层状兴奋性神经元。在这里，我们可以区分这些神经元亚型的不同顺序层，例如:
image
根据这些预测分数，我们还可以预测位置受到空间限制的细胞类型。我们使用基于标记点过程的相同方法来定义空间变量特征，但使用细胞类型预测评分作为“标记”，而不是使用基因表达。
1cortex <- FindSpatiallyVariableFeatures(cortex, assay = "predictions", features = rownames(cortex),
2 r.metric = 5, slot = "data")
3top.clusters <- head(SpatiallyVariableFeatures(cortex), 4)
4SpatialPlot(object = cortex, features = top.clusters, ncol = 2)
5
最后，我们证明我们的整合程序能够恢复已知的神经元和非神经元亚群的空间定位模式，包括层兴奋性亚群、第1层星形胶质细胞和皮层灰质。
1SpatialFeaturePlot(cortex, features = c("Astro", "L2/3 IT", "L4", "L5 PT", "L5 IT", "L6 CT", "L6 IT",
2 "L6b", "Oligo"), pt.size.factor = 1, ncol = 2, crop = FALSE, alpha = c(0.1, 1))
处理Seurat中的多个slices
这个老鼠大脑的数据集包含了另一个与大脑另一半相对应的切片。这里我们读取它并执行相同的初始规范化。
1brain2 <- LoadData("stxBrain", type = "posterior1")
2brain2 <- SCTransform(brain2, assay = "Spatial", verbose = FALSE)为了处理同一个Seurat对象中的多个片，我们提供了merge函数。其实之前就有的呀。
1brain.merge <- merge(brain, brain2)这样就可以对底层的RNA表达数据进行联合降维和聚类。
1DefaultAssay(brain.merge) <- "SCT"
2VariableFeatures(brain.merge) <- c(VariableFeatures(brain), VariableFeatures(brain2))
3brain.merge <- RunPCA(brain.merge, verbose = FALSE)
4brain.merge <- FindNeighbors(brain.merge, dims = 1:30)
5brain.merge <- FindClusters(brain.merge, verbose = FALSE)
6brain.merge <- RunUMAP(brain.merge, dims = 1:30)
7最后，可以在单个UMAP图中共同可视化数据。SpatialDimPlot和SpatialFeaturePlot将默认将所有片绘制为列，将groupings/features 绘制为行。
1DimPlot(brain.merge, reduction = "umap", group.by = c("ident", "orig.ident"))
1SpatialDimPlot(brain.merge)
1SpatialFeaturePlot(brain.merge, features = c("Hpca", "Plp1"))
image
我们要感谢Nigel Delaney和Stephen Williams对Seurat 分析空间数据代码的有益反馈和贡献。
教程官网(https://satijalab.org/seurat/v3.1/spatial_vignette.html)