文献名：Spatial charting of single-cell transcriptomes in tissues
期刊：Nature Biotechnology
IF：33.1
发表日期：21 March 2022
单细胞 RNA 测序 (scRNA-seq) 方法可以分析单细胞的转录组，但不能保留空间信息。相反，空间转录组学 (ST) 分析可以描绘组织切片中的空间区域，但没有单细胞基因组分辨率。为了解决这个问题，已经设计了诸多计算方法如（cell2location、RCTD、SPOTlight）等来将 ST spot解卷积为不同细胞类型的比例。然而，空间去卷积方法仅限于推断每个spot的细胞类型比例，无法实现单细胞分辨率。此外，去卷积方法将细胞类型进一步解析为反映不同生物学功能的更细粒度的“细胞状态”（表达程序）的能力有限。最后，大多数反卷积方法只能预测分类标签，而不能以空间分辨率推断连续的细胞信息（例如，谱系轨迹、基因特征、连续表型）等。
而CellTrek可以根据 scRNA-seq 和 ST 数据将单个细胞直接映射回组织切片中的空间坐标。这种方法提供了一种不同于 ST 反卷积的新模式，能够更灵活、更直接地研究具有空间地形的单细胞数据。CellTrek 工具包还提供了两个下游分析模块，包括用于空间共定位分析的 SColoc 和用于空间共表达分析的 SCoexp。
CellTrek的计算原理
首先CellTrek将ST和scRNA-seq数据集成到共享的特征空间中。CellTrek使用ST数据训练多元随机森林（RF）模型，使用共享降维特征来预测空间坐标。为了提高空间分辨率，CellTrek在ST数据上引入了空间非线性插值方法。然后将训练的模型应用于共嵌入的数据，推导出一个RF距离矩阵，该矩阵测量每个ST点和空间坐标监督的单个细胞之间的表达相似度。基于上述的RF距离矩阵，CellTrek使用阈值化后的相互最近邻（MNN）生成稀疏点小区图。最后，CellTrek从相邻点传输细胞的空间坐标。为了提高兼容性，CellTrek可以接受从其他方法（例如，novoSpaRc20）计算出的任何细胞位置概率/距离矩阵，作为细胞空间图的输入。
CellTrek工作原理
操作方法
•1.首先在通过Github链接下载示例数据，这里以小鼠大脑数据为示例：
￮scRNA: https://www.dropbox.com/s/ruseq3necn176c7/brain_sc.rds?dl=0
￮ST:https://www.dropbox.com/s/azjysbt7lbpmbew/brain_st_cortex.rds?dl=0
•2.加载数据
导入数据后，使用
Golibrary(CellTrek)library(Seurat)library(ggplot2)library(dplyr)library(plotly)library(psycModel)library(pagedown)library(visNetwork)library(viridis)scdata<-readRDS('brain_sc.rds')stdata<-readRDS('brain_st_cortex.rds')#细胞ID中横线等特殊字符转为点stdata <- RenameCells(stdata, new.names=make.names(Cells(stdata)))scdata <- RenameCells(scdata, new.names=make.names(Cells(scdata)))
•3.共嵌入scRNA与ST数据，将单细胞RNA数据与空间转录组数据通过Seurat包中TransferData函数进行映射，寻找锚点基因后，将共有的锚点数据合并为新的数据后再进行降维聚类得到共嵌入的降维结果。
Gobrain_traint <- CellTrek::traint(st_data=stdata, sc_data=scdata, sc_assay='RNA', cell_names='cell_type')p3<-DimPlot(brain_traint, group.by = "type")ggsave(p3,file='type_DimPlot.pdf')
scRNA与ST共嵌入降维结果
•4.使用 ST 数据训练multivariate random forests (RF) model，以使用共嵌入的降维结果特征预测空间坐标。引入非线性插值法增强 ST spot数据以增加空间分辨率。然后将训练后的模型应用于共嵌入数据以导出 RF 距离矩阵，该矩阵测量 ST spot和由空间坐标监督的单个细胞之间的表达相似性。基于 RF 距离矩阵，在阈值过滤后使用相互最近邻 (MNN) 生成稀疏spot细胞图，最后从相邻spot传输细胞的空间坐标。为了提高兼容性，CellTrek 可以接受从其他方法（例如 novoSpaRc）计算的任何细胞位置概率/距离矩阵作为细胞空间图表的输入。
Gobrain_celltrek <- CellTrek::celltrek(st_sc_int=brain_traint, int_assay='traint', sc_data=scdata, sc_assay = 'RNA',reduction='pca', intp=T, intp_pnt=5000, intp_lin=F, nPCs=30, ntree=1000,dist_thresh=0.55, top_spot=5, spot_n=5, repel_r=20, repel_iter=20, keep_model=T)$celltreksource('celltrek_visplot.R')#celltrek_df必须包含coord_x,coord_y,id_new以及其他做着色的分类型/连续型列，coord_x,coord_y来源于brain_celltrek@images$anterior1@coordinatescelltrek_df<-brain_celltrek@meta.data %>% dplyr::select(coord_x, coord_y, cell_type:id_new)#rdata中的图片img<-brain_celltrek@images$anterior1@image#scale_fac是低分辨率的缩放因子scale_fac<-brain_celltrek@images$anterior1@scale.factors$lowresplot<-celltrek_visplot(celltrek_df, img, scale_fac, color_var='cell_type',pnt_size=2.8)plotly::save_image(plot, file = "plotly_plot.pdf")
CellTrek结果
•5.细胞共定位分析，SColoc基于CellTrek结果图，可以选择三种方法即 Kullback-Leibler 散度 (KL)、Delaunay 三角测量 (DT) 和 K-最近邻距离 (KD)等方法计算细胞型空间相异矩阵，并使用最小生成树 (MST)生成树结构。在自举样本上重复执行这些步骤，计算出相异矩阵的一致性矩阵 或 MST，生成最终的单元类型空间图表示
细胞共定位原理
Go#以谷氨酸能神经元细胞类型为例，建议去掉细胞数少于20的细胞类型，细胞分类不能包含'X1','X2'以及'1','2'等字段，会与代码冲突glut_cell <- c('L2/3 IT', 'L4', 'L5 IT', 'L5 PT', 'NP', 'L6 IT', 'L6 CT', 'L6b')names(glut_cell) <- make.names(glut_cell)brain_celltrek_glut <- subset(brain_celltrek, subset=cell_type %in% glut_cell)brain_celltrek_glut$cell_type <- factor(brain_celltrek_glut$cell_type, levels=glut_cell)brain_sgraph_KL <- CellTrek::scoloc(brain_celltrek_glut, col_cell='cell_type', use_method='KL', eps=1e-50)## 从图中提取最小生成树（MST）结果brain_sgraph_KL_mst_cons <- brain_sgraph_KL$mst_consrownames(brain_sgraph_KL_mst_cons) <- colnames(brain_sgraph_KL_mst_cons) <- glut_cell[colnames(brain_sgraph_KL_mst_cons)]## We then extract the metadata (including cell types and their frequencies)brain_cell_class <- brain_celltrek@meta.data %>% dplyr::select(id=cell_type) %>% table %>% as.data.framesource('scoloc_vis.R')plot2<-scoloc_vis(brain_sgraph_KL_mst_cons, meta_data=brain_cell_class,edge_val_cutoff=0.3,node_col_column='id')visSave(plot2, file = "network.html")system('wkhtmltopdf network.html network.pdf')
细胞共定位结果
•6.为了研究不同的表达程序是否分布在不同的区域，作者开发了 SCoexp，它利用 CellTrek 坐标来检测目标细胞内的共表达基因模块。对于感兴趣的细胞，基于 CellTrek 图，SCoexp 首先根据其空间距离使用 RBF 计算空间核矩阵。接下来，基于空间核矩阵和细胞基因表达矩阵，SCoexp计算空间加权基因共表达。然后使用 CC（一致性聚类） 或 WGCNA 鉴定基因模块。对于识别出的共表达模块，可以计算模块活性分数并将其映射回 CellTrek 坐标。
细胞共表达分析原理
Go## 选择感兴趣的细胞类型brain_celltrek_l5 <- subset(brain_celltrek, subset=cell_type=='L5 IT')brain_celltrek_l5@assays$RNA@scale.data <- matrix(NA, 1, 1)brain_celltrek_l5$cluster <- gsub('L5 IT VISp ', '', brain_celltrek_l5$cluster)p4<-DimPlot(brain_celltrek_l5, group.by = 'cluster')ggsave(p4,file='DimPlot_interst.pdf')## 选择2000个高变基因并且排除线粒体、核糖体、高0表达的基因brain_celltrek_l5 <- FindVariableFeatures(brain_celltrek_l5)vst_df <- brain_celltrek_l5@assays$RNA@meta.features %>% data.frame %>% mutate(id=rownames(.))nz_test <- apply(as.matrix(brain_celltrek_l5[['RNA']]@data), 1, function(x) mean(x!=0)*100)hz_gene <- names(nz_test)[nz_test<20]mt_gene <- grep('^Mt-', rownames(brain_celltrek_l5), value=T)rp_gene <- grep('^Rpl|^Rps', rownames(brain_celltrek_l5), value=T)vst_df <- vst_df %>% dplyr::filter(!(id %in% c(mt_gene, rp_gene, hz_gene))) %>% arrange(., -vst.variance.standardized)feature_temp <- vst_df$id[1:2000]brain_celltrek_l5_scoexp_res_cc <- CellTrek::scoexp(celltrek_inp=brain_celltrek_l5, assay='RNA', approach='cc', gene_select = feature_temp, sigm=140, avg_cor_min=.4, zero_cutoff=3, min_gen=40, max_gen=400)brain_celltrek_l5_k <- rbind(data.frame(gene=c(brain_celltrek_l5_scoexp_res_cc$gs[[1]]), G='K1'),data.frame(gene=c(brain_celltrek_l5_scoexp_res_cc$gs[[2]]), G='K2')) %>%magrittr::set_rownames(.$gene) %>% dplyr::select(-1)pdf('scoexpheatmap.pdf')pheatmap::pheatmap(brain_celltrek_l5_scoexp_res_cc$wcor[rownames(brain_celltrek_l5_k), rownames(brain_celltrek_l5_k)],clustering_method='ward.D2', annotation_row=brain_celltrek_l5_k, show_rownames=F, show_colnames=F,treeheight_row=10, treeheight_col=10, annotation_legend = T, fontsize=8,color=viridis(10), main='L5 IT spatial co-expression')dev.off()## 模块打分结果brain_celltrek_l5 <- AddModuleScore(brain_celltrek_l5, features=brain_celltrek_l5_scoexp_res_cc$gs, name='CC_', nbin=10, ctrl=50, seed=42)p5<-FeaturePlot(brain_celltrek_l5, grep('CC_', colnames(brain_celltrek_l5@meta.data), value=T))ggsave(p5,file='FeaturePlot_interst.pdf')p6<-SpatialFeaturePlot(brain_celltrek_l5, grep('CC_', colnames(brain_celltrek_l5@meta.data), value=T))ggsave(p6,file='SpatialFeaturePlot_interst.pdf')
感兴趣细胞中基因共表达分析结果
模块打分降维图结果
模块打分空间图结果
参考资料
https://www.nature.com/articles/s41587-022-01233-1
https://github.com/navinlabcode/CellTrek
CytoTRACE：细胞分化潜能分析|生信开发实战
Visium HD数据分析之Bin2Cell |生信开发实战
ROC分析介绍|生信开发实战
如果使用find_circ来鉴定circRNA|生信开发实战
本文系联川生物公众号原创文章，未经授权禁止转载，侵权必究！
扫描下方二维码
点分享
点点赞
点在看
