我认为需要自己写,那是一段异常艰难的时光。
学,然后知不足。
与scRNA-seq参考整合
为了促进Slide-seq数据集的细胞类型注释,我们利用了由Saunders *,Macosko *等人制作的现有小鼠单细胞RNA-seq海马数据集。数据可在此处作为已处理的Seurat对象下载,而原始计数矩阵可在DropViz网站上获得。
ref <- readRDS("../data/mouse_hippocampus_reference.rds")
本文的原始注释在Seurat对象的metadata 中提供。这些注释以几种“resolutions”提供,从大类(ref$class
)到亚群ref$subcluster
。出于本教程的目的,我们将对细胞类型注释(ref$celltype
)进行修改,使之达到良好的平衡。
我们将从运行Seurat标签转移方法开始,以预测每个bead的主要细胞类型。
anchors <- FindTransferAnchors(reference = ref, query = slide.seq, normalization.method = "SCT",
npcs = 50)
predictions.assay <- TransferData(anchorset = anchors, refdata = ref$celltype, prediction.assay = TRUE,
weight.reduction = slide.seq[["pca"]], dims = 1:50)
slide.seq[["predictions"]] <- predictions.assay
然后,我们可以可视化一些主要类型的预测分数。
DefaultAssay(slide.seq) <- "predictions"
SpatialFeaturePlot(slide.seq, features = c("Dentate Principal cells", "CA3 Principal cells", "Entorhinal cortex",
"Endothelial tip", "Ependymal", "Oligodendrocyte"), alpha = c(0.1, 1))
slide.seq$predicted.id <- GetTransferPredictions(slide.seq)
Idents(slide.seq) <- "predicted.id"
SpatialDimPlot(slide.seq, cells.highlight = CellsByIdentities(object = slide.seq, idents = c("CA3 Principal cells",
"Dentate Principal cells", "Endothelial tip")), facet.highlight = TRUE)
识别空间变异基因
正如Visium教程中提到的那样,我们可以通过两种通用方法来识别空间变异基因:预先标注的解剖区域之间的差异表达检测或基因特征对空间位置的依赖性的统计信息。
在这里,我们演示后者,通过在FindSpatiallyVariableFeatures()中设置method = 'moransi'。 Moran's I计算总体空间自相关性,并给出一个统计量(类似于相关系数),该统计量用于衡量要素对空间位置的独立性。这使我们能够根据基因表达的变化对基因进行排名。为了便于快速估算此统计信息,我们实施了一种基本的分级策略,该策略将在Slide-seqpuck 上绘制一个矩形网格,并对每个分级中的基因和位置取平均值。x和y方向上的bins 分别由
x.cuts和
y.cuts参数控制。此外,虽然不是必需的,但安装可选
Rfast2软件包(
install.packages('Rfast2')),将通过更有效的方式来减少运行时间。
DefaultAssay(slide.seq) <- "SCT"
slide.seq <- FindSpatiallyVariableFeatures(slide.seq, assay = "SCT", slot = "scale.data", features = VariableFeatures(slide.seq)[1:1000],
selection.method = "moransi", x.cuts = 100, y.cuts = 100)
现在,我们可视化由Moran's I识别的前6个基因的表达。
SpatialFeaturePlot(slide.seq, features = head(SpatiallyVariableFeatures(slide.seq, selection.method = "moransi"),
6), ncol = 3, alpha = c(0.1, 1), max.cutoff = "q95")