【单细胞技术贴】空间转录组与单细胞转录组的整合分析（下篇）

专题介绍：单细胞RNA-seq被评为2018年重大科研进展，但实际上这是老技术。2015年，商品化单细胞RNA测序流程已经建立，成果发表在Cell上。今年井喷式发文章，关注点那么高，是因为最近这项技术全面商品化了。

SPOTlight 的目标是提供一种工具，能够对包含细胞混合物的每个捕获位置中存在的细胞类型和细胞类型比例进行解卷积，最初是为 10X 的 Visium – 空间转录组学技术开发的，它可用于返回混合物的所有技术细胞。SPOTlight 基于通过NMFreg模型为每种细胞类型查找主要的topic profile signatures，并找到最适合想要解卷积的点的组合。其分析流程如下图所示：

该软件为GitHub上面开源，可以从GitHub上下载SPOTlight（https://github.com/MarcElosua/SPOTlight）然后本地安装，或者直接通过网址安装，具体安装命令如下：

install.packages("devtools")
devtools::install_github("https://github.com/MarcElosua/SPOTlight")#如果网速不快的话，可以先下载，然后本地安装

使用

library(RCTD)
library(Matrix)
library(Hmisc)
library(annotools)
library(data.table)
library(presto)
library(Seurat)
library(SeuratWrappers)
require(dplyr)
library(data.table)
library(SeuratData)
library(gt)
library(SPOTlight)
library(igraph)
library(RColorBrewer)

准备数据

这里的单细胞转录组数据为已经鉴定好细胞类型的单细胞数据的seurat对象，一般为Seurat4.0以上版本，细胞类型的列名为cell_annotation，一般细胞类型可以通过命令：sc@meta.data$cell_annotation进行查看，如果细胞类型中有未知细胞类型比如none，则需要删掉这部分细胞，不然会影响解卷积分析结果。然后对单细胞的进行marker基因分析，这里用的是RunPrestoAll函数，这个函数效率较高，与FindAllMarkers函数的功能一样，但是效率差不多是快10倍以上。

set.seed(123)###单细胞转录组数据准备
sc<-readRDS('seurat.celltype.rds')#一般为鉴定好细胞类型的seurat对象，并且sc@meta.data矩阵中含有列名为cell_annotation列
Idents(sc)<-gsub(' |/|-','_',sc@meta.data[as.vector(names(Idents(sc))),'cell_annotation'])
sc<-subset(sc,idents=as.vector(setdiff(as.vector(unique(Idents(sc))),c('none')))) ####删掉未知细胞类型
cluster_markers_all<-RunPrestoAll(sc,only.pos = TRUE,slot = "data",verbose = TRUE)###空间转录组数据的准备
rds_file='spatial.rds'
spatial <- readRDS(rds_file)
DefaultAssay(spatial) <- "RNA"

空间转录组的解卷机

解卷积函数参数较多，这里必选参数为se_sc 、counts_spatial、clust_vr、cluster_markers，主要为输入的Seurat对象、空间转录组原始表达矩阵counts、细胞类型注释信息列名、marker基因。其他参数可以根据需要进行调整，一般不需要调整。

###解卷积分析
spotlight_ls <- spotlight_deconvolution(se_sc = sc,  counts_spatial = spatial@assays$RNA@counts,  clust_vr = "cell_annotation", # Variable in sc_seu containing the cell-type annotation  cluster_markers = cluster_markers_all, # Dataframe with the marker genes  cl_n = 100, # number of cells per cell type to use  hvg = 3000, # Number of HVG to use  ntop = NULL, # How many of the marker genes to use (by default all)  transf = "uv", # Perform unit-variance scaling per cell and spot prior to factorzation and NLS  method = "nsNMF", # Factorization method  min_cont = 0 # Remove those cells contributing to a spot below a certain threshold   )###一般2500多个spots差不多不到10min左右

结果展示和说明

prefix='TEST'nmf_mod <- spotlight_ls[[1]]decon_mtrx <- spotlight_ls[[2]]#topic结果展示 Next we can take a look at the how the individual topic profiles of each cell within each cell-type behave.Here we expect that all the cells from the same cell type show a similar topic profile distribution, if not there might be a bit more substructure in that cluster and we may only be capturing one or the other.h <- NMF::coef(nmf_mod[[1]])rownames(h) <- paste("Topic", 1:nrow(h), sep = "_")topic_profile_plts <- SPOTlight::dot_plot_profiles_fun(  h = h,  train_cell_clust = nmf_mod[[2]])pdf(paste(prefix,'topic.celltype.pdf',sep='_'),w=10,h=7)p<-topic_profile_plts[[2]] + ggplot2::theme(  axis.text.x = ggplot2::element_text(angle = 90),   axis.text = ggplot2::element_text(size = 12))print(p)dev.off()pdf(paste(prefix,'topic.celltype.every.pdf',sep='_'),w=45,h=12)p<-topic_profile_plts[[1]] + theme(axis.text.x = element_text(angle = 90),axis.text = element_text(size = 12))print(p)dev.off()

这里主要是退topic结果进行展示：

其实这里的topic相当于通过NMF进行聚类，将相似的细胞聚类在一起，一个topic相当于一个聚类，上图展示了每个Topic对象的细胞类型，如果细胞类型与topic能够一一对象，这样的结果就最好不过，如果不对应，可能细胞分的不够细。

上图展示了每个topic对应细胞类型，一般理想情况下，一个topic对应一种细胞类型，当然也不是完全对应，可能有部分细胞处于其他topic。上图展示了B细胞类型在不同topic中的分布情况，从结果可以看出，大部分细胞属于B细胞。

下面对每个topic的基因进行展示webshot包，而且这里为便于基因筛选，这里输出了动态网页表格，这里需要具体结果过程如下：

library("webshot")#webshot::install_phantomjs() devtools::install_github("wch/webshot")basis_spotlight <- data.frame(NMF::basis(nmf_mod[[1]]))colnames(basis_spotlight) <- unique(stringr::str_wrap(nmf_mod[[2]], width = 30))dtable<-basis_spotlight %>%  dplyr::arrange(desc(Macrophage)) %>%  round(., 5) %>%   DT::datatable(., filter = "top",caption = paste('Table 1: The most important genes for each topic .',prefix,sep=' ')) 
html <- paste(prefix,'dtable.html',sep='_')saveWidget(dtable, html)

上述结果输出了html文件，可以打个html文件，然后对表观进行筛选，具体的结果如下图所示：

上述动态表格是可以筛选，或者排序，这里值的大小没有实际的意义，值越大代表对细胞类型的贡献度越大。

# This is the equivalent to setting min_cont to 0.04decon_mtrx_sub <- decon_mtrx[, colnames(decon_mtrx) != "res_ss"]decon_mtrx_sub[decon_mtrx_sub < 0.08] <- 0decon_mtrx <- cbind(decon_mtrx_sub, "res_ss" = decon_mtrx[, "res_ss"])rownames(decon_mtrx) <- colnames(anterior)
decon_df <- decon_mtrx %>%  data.frame() %>%  tibble::rownames_to_column("barcodes")
anterior@meta.data <- anterior@meta.data %>%  tibble::rownames_to_column("barcodes") %>%  dplyr::left_join(decon_df, by = "barcodes") %>%  tibble::column_to_rownames("barcodes")cell_types_all <- colnames(decon_mtrx)[which(colnames(decon_mtrx) != "res_ss")]
SPOTlight::spatial_scatterpie(se_obj = anterior,                              cell_types_all = cell_types_all,                              img_path = "sample_data/spatial/tissue_lowres_image.png",                              pie_scale = 0.4)

上图展示了解卷积的每个spots中细胞类型比例的饼图（每个饼代表一个spots，饼图中不同颜色代表不同的细胞类型），可以通过上图看出，大部分的spots并不是一种细胞类型，是由多种细胞类型混合而成，因此解卷积的方法进行空间转录组的方法是是否合适。当然也可以展示单独某个细胞类型的结果，如下图所示：

SPOTlight::spatial_scatterpie(se_obj = anterior,                              cell_types_all = cell_types_all,                              img_path = "sample_data/spatial/tissue_lowres_image.png",                              cell_types_interest = "L6b",                              pie_scale = 0.8)

上图展示了含有L6b细胞类型的结果，也就是说只要spots中含有L6b细胞类型，其结果都将在上图进行展示。

总结

单细胞转录组与空间转录组学的整合可以产生组织中细胞亚群的高分辨率图谱。空间转录组的技术正在迅速发展，但没有一种单一的空间转录组学技术适合所有应用。根据所提出的生物学问题，可以设计实验方法以将任何空间转录组学方法与单细胞转录组相结合。除了开发增强方法之外，仔细选择用于整合此类数据的算法也至关重要，因为尚不存在能够以单细胞分辨率、单细胞转录组深度和全转录组覆盖空间解析组织的空间转录组学方法。这种综合方法开始在空间上绘制发育和疾病中的特定细胞亚群，并阐明这些细胞群协同塑造组织表型的机制。

这里介绍了两种解卷积的方法进行单细胞转录组和空间转录组整合分析，当然，能够用于解卷积的方法还有多种方法比如SpatialDWLS、stereoscope等等，以后可能也有更多的方法，不同的方法有不同的优缺点，比如这里的RCTD方法有点在于可以很好的跨平台和细胞亚群定位等等，老师可以根据自己的研究需求，挑选合适的方法。就个人而言，比较喜欢SPOTlight的展示方法，该结果展示，能够直接的展示每个spots中不同类型的细胞比例情况，结果较为清晰明确。比如下图所示，图A展示了单细胞的细胞类型，图B展示了空间转录组中每个spots中不同细胞类型的比例，该图就较为明显的看出主要细胞类型分布，当然RCTD结果也可以用此方法进行展示，不过需要对脚本或者可视化函数进行较大的修改。

针对不同的方法，在SpatialDWLS中进行了比较，下次再对此文献进行解读，老师可以自己去看一下SpatialDWLS的文献，通过不同的方法比较，发现RCTD和SpatialDWLS方法都不错。具体结果如下：

图A为SpatialDWLS整体流程，图B/C/D为不同方法的比较，从结果来说，SpatialDWLS、RCTD、SPOTlight结果准确性不错，SPOTlight可能运行效率低一些，其他的相差不大。

展望

随着空间转录组和单细胞的转录组的技术发展，今后二者将会越来越普及，将成为科研界常用的工具。比如单细胞转录组技术，在前几年，我们还想着10x的单细胞平台有多方便和便宜，随着这几年的技术发展，特别是今年的海量单细胞技术的推出，一次性使用可以同时检测70多万的细胞，这将大大改变科研的方式，使得今后的单细胞图谱更为精确，单细胞研究越来越准确。

空间转录组也是同样的，现在价格偏贵，这会大大的限制其应用，而且空间转录组具体能够给我们带来什么样的结果？这些结果对解决生物学问题具体有什么样的帮助？这都需要我们考虑的问题，就目前而言，空间转录组提供了一个空间信息，因此应用较大的可能是空间的细胞通讯和发育轨迹，这是纯分析角度来说，具体其他研究问题，可能会带来更精确的结果。

目前空间转录组还存在较多的问题，将来可能随着技术的发展，目前存在的问题都会一一的解决。比如目前空间转录组的spots是混合细胞，远远达不到单细胞水平，这是将在较长一个时间内没法达到的水平，因此这会一定角度的限制了空间转录组的结果的准确性，这也是为什么目前需要将空间转录组与单细胞转录组结合分析的缘故，虽然可能有部分技术能到亚细胞级别，但技术不成熟，属于实验室中技术，没法进行推广；另外，目前10xVisum空间转录组的玻片较小，6.5×6.5mm，如果切片较大的样品，将远远不能覆盖完整的组织样品，不过最近华大基因推出的空间转录组的玻片到cm级别，这也是其优势之一；此外空间转录的算法和工具并不成熟，目前虽然有很多的空间转录组的工具，但是目前的工具都还处于初始阶段，没有一个完全成熟的工具；再者空间转录组目前的成本偏贵，这是一个较为重要的限制因素。

总而言之，空间转录组技术目前存在一定的局限，但是随着技术的发展，其分辨率、算法工具、成本将都会有较大的提升，空间转录组技术将会给科研界带来新的视角，为人类的临床、药物开发、伴随诊断、生物标志物开发提供强有力的支持（https://pages.10xgenomics.com/visium-clinical-translation-network.html）。

参考文献：
[1] Longo S K ,  Guo M G ,  Ji A L , et al. Integrating single-cell and spatial transcriptomics to elucidate intercellular tissue dynamics[J]. Nature Reviews Genetics.
[2] V Marx. Method of the Year: spatially resolved transcriptomics[J]. Nature Methods, 2021, 18(1):9-14.
[3] Asp M , Bergenstrhle J , Lundeberg J . Spatially Resolved Transcriptomes—Next Generation Tools for Tissue Exploration[J]. BioEssays, 2020, 42(10).
[4] Cable D M , Murray E , Zou L S , et al. Robust decomposition of cell type mixtures in spatial transcriptomics[J]. Nature Biotechnology, 2021:1-10.
[5] Dong R ,  Yuan G C . SpatialDWLS: accurate deconvolution of spatial transcriptomic data[J]. Genome Biology, 2021.