首页
学习
活动
专区
圈层
工具
发布
首页
学习
活动
专区
圈层
工具
MCP广场
社区首页 >专栏 >整合单细胞和空转数据多种方法之R包semla

整合单细胞和空转数据多种方法之R包semla

作者头像
生信菜鸟团
发布于 2024-04-25 10:26:40
发布于 2024-04-25 10:26:40
73602
代码可运行
举报
文章被收录于专栏:生信菜鸟团生信菜鸟团
运行总次数:2
代码可运行

我注意到这个R包也提供了一个整合单细胞和空转数据的算法NNLS(https://ludvigla.github.io/semla/articles/cell_type_mapping_with_NNLS.html)。因此,本文对这个算法进行测试。

一 . NNLS代码实战

Step0.R包安装

github在https://github.com/ludvigla/semla

官方教程在:https://ludvigla.github.io/semla/articles/cell_type_mapping_with_NNLS.html

安装:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
remotes::install_github("ludvigla/semla")
Step1.加载示例数据
代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
library(readr)
library(semla)
library(tibble)
library(dplyr)
library(Seurat)
library(purrr)
library(patchwork)
library(ggplot2)
#install.packages('RcppML')  
library(RcppML)
代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
# 空转
brain_st_cortex <- read_rds("./Rawdata/brain_st_cortex.rds")

# 单细胞
brain_sc <- read_rds("./Rawdata/brain_sc.rds")

## Rename the cells/spots with syntactically valid names
brain_st_cortex <- RenameCells(brain_st_cortex, new.names=make.names(Cells(brain_st_cortex)))
brain_sc <- RenameCells(brain_sc, new.names=make.names(Cells(brain_sc)))

可视化空转数据:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
## Visualize the ST data
SpatialDimPlot(brain_st_cortex)

image-20231105144510038

可视化单细胞数据:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
## Visualize the scRNA-seq data
DimPlot(brain_sc, label = T, label.size = 4.5, group.by = "cell_type")

image-20231105144526583

Step2. 标准流程

这里,我们将对10x Visium数据和单细胞数据应用相同的对数归一化处理流程。这里将高变基因的数量设置得非常高,因为稍后我们将对单细胞数据中的高变基因与10x Visium数据中的高变基因取交集。

由于NNLS算法运行非常快,因此不需要像常规处理那样只取2000个高变基因。相反,我们可以使用在单细胞和10x Visium数据中存在交集的所有基因。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
# Normalize data and find variable features for Visium data
brain_st_cortex <- brain_st_cortex |>
  NormalizeData() |>
  FindVariableFeatures(nfeatures = 10000)

# Normalize data and run vanilla analysis to create UMAP embedding
brain_sc <- brain_sc |>
  NormalizeData() |>
  FindVariableFeatures() |> 
  ScaleData() |> 
  RunPCA() |> 
  RunUMAP(reduction = "pca", dims = 1:30)

# Rerun FindVariableFeatures to increase the number before cell type deconvolution
brain_sc <- brain_sc |> 
  FindVariableFeatures(nfeatures = 10000)
Step3. Run NNLS

使用RunNNLS()函数,输入单细胞数据和空转数据,指定单细胞细胞类型的标签名称,即可运行NNLS反卷积:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
DefaultAssay(brain_st_cortex) <- "Spatial"

ti <- Sys.time()
brain_st_cortex <- RunNNLS(object = brain_st_cortex, 
                            singlecell_object = brain_sc, 
                            groups = "cell_type",
                           minCells_per_celltype = 10)
# ── Predicting cell type proportions ──
# 
# ℹ Fetching data from Seurat objects
# →   Filtering out features that are only present in one data set
# →   Kept 2961 features for deconvolution
# ℹ Preparing data for NNLS
# →   Downsampling scRNA-seq data to include a maximum of 50 cells per cell type
# →   Kept 15 cell types after filtering
# →   Calculating cell type expression profiles
# ℹ Predicting cell type proportions with NNLS for 15 cell types
# ℹ Returning results in a new 'Assay' named 'celltypeprops'
# ℹ Setting default assay to 'celltypeprops'
# ✔ Finished
sprintf("RunNNLS completed in %s seconds", round(Sys.time() - ti, digits = 2))
#"RunNNLS completed in 0.91 seconds"

1秒钟不到就运行结束了。

注意:这个函数默认会过滤掉小于10个细胞的细胞类型,我们可以指定minCells_per_celltype = 参数进行设置。

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
# Check available cell types
rownames(brain_st_cortex)
# [1] "Astro"   "Endo"    "L2/3 IT" "L4"      "L5 IT"   "L5 PT"   "L6 CT"   "L6 IT"   "L6b"     "Lamp5"   "NP"      "Pvalb"  
# [13] "Sncg"    "Sst"     "Vip" 
Step4. 可视化
代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
# Plot selected cell types
DefaultAssay(se_brain_spatial) <- "celltypeprops"

selected_celltypes <- c('L6 CT', 'L5 IT', 'L6b', 'L2/3 IT', 'L4', 'Astro')
plots <- lapply(seq_along(selected_celltypes), function(i) {
  SpatialFeaturePlot(brain_st_cortex,features = selected_celltypes[1]) +
    scale_fill_gradientn(colours = RColorBrewer::brewer.pal(n = 9, name = "Spectral") %>% rev())+
  plot_layout(guides = "collect") & 
  theme(legend.position = "right", legend.margin = margin(b = 50),
        legend.text = element_text(angle = 0),
        plot.title = element_blank())
  }) %>%  setNames(nm = selected_celltypes)

wrap_plots(plots, ncol=4)

image-20240419182522039

这里我们对比一下cell2location的结果【整合单细胞和空转数据多种方法之Cell2location】:

image-20231228230033246

就肉眼来看,结果还是挺一致的。

这里使用NNLS和cell2location的结果做一个基于spot的相关性分析:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
library(ggpubr)
cell2locat_res = read.csv("./Rawdata/st_cell2location_res.csv", row.names = 1, check.names = F)
head(cell2locat_res)
table(row.names(cell2locat_res) == colnames(brain_st_cortex))
# TRUE 
# 1075 
  • 例如L6 CT细胞:
代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
input.plot = data.frame(L6CT_cell2locat = cell2locat_res$`q05cell_abundance_w_sf_L6 CT`,
                        L6CT_NNSL = as.numeric(brain_st_cortex@assays$celltypeprops["L6 CT",]*100))

ggscatter(input.plot, cor.method = "p",
          x = "L6CT_cell2locat", y = "L6CT_NNSL",
          title = "L6 CT cells",
          add = "reg.line", conf.int = TRUE, 
          add.params = list(fill = "lightgray"))+
  stat_cor(method = "p",
           label.sep = "\n",
           label.y.npc = "top",
           label.x = 10, 
           label.y = 1)

image-20240419203539900

  • L2/3 IT细胞:
代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
input.plot = data.frame(L23ICT_cell2locat = cell2locat_res$`q05cell_abundance_w_sf_L2/3 IT`,
                        L23ICT_NNSL = as.numeric(brain_st_cortex@assays$celltypeprops["L2/3 IT",]* 100))

ggscatter(input.plot, cor.method = "p",
          x = "L23ICT_cell2locat", y = "L23ICT_NNSL",
          title = "L2/3 IT cells",
          add = "reg.line", conf.int = TRUE, 
          add.params = list(fill = "lightgray"))+
  stat_cor(method = "p",
           label.sep = "\n",
           label.y.npc = "top",
           label.x = 20, 
           label.y = 2)

image-20240419203624188

  • L4细胞:
代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
input.plot = data.frame(L4_cell2locat = cell2locat_res$q05cell_abundance_w_sf_L4,
                        L4_NNSL = as.numeric(brain_st_cortex@assays$celltypeprops["L4",]*100))

ggscatter(input.plot, cor.method = "p",
          x = "L4_cell2locat", y = "L4_NNSL",
          title = "L4 cells",
          add = "reg.line", conf.int = TRUE, 
          add.params = list(fill = "lightgray"))+
  stat_cor(method = "p",
           label.sep = "\n",
           label.y.npc = "top",
           label.x = 16, 
           label.y = 1)

image-20240419203656828

  • .丰度较低的星形胶质细胞:
代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
library(ggpubr)
input.plot = data.frame(Astro_cell2locat = cell2locat_res$q05cell_abundance_w_sf_Astro,
                        Astro_NNSL = as.numeric(brain_st_cortex@assays$celltypeprops["Astro",]*100))

ggscatter(input.plot, cor.method = "p",
          x = "Astro_cell2locat", y = "Astro_NNSL",
          title = "Astro cells",
          add = "reg.line", conf.int = TRUE, 
          add.params = list(fill = "lightgray"))+
  stat_cor(method = "p",
           label.sep = "\n",
           label.y.npc = "top",
           label.x = 7, 
           label.y = 1)

image-20240419203731108

尽管两种算法的相关性非常好,但是可以看到NNSL计算得到的每种细胞类型的结果区间(区间比较大)均大于cell2location算法的结果区间。

我们还可以使用semla的内置函数MapMultipleFeatures()在同一张slide里可视化多种细胞类型的分布情况(cell2location也有类似功能),但是需要先将seurat对象转换为semla需要的对象,这点我在【空转可视化R包semla(一)入门介绍】介绍过:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
semla.data <- UpdateSeuratForSemla(brain_st_cortex)
# Plot multiple features
MapMultipleFeatures(semla.data, 
                    pt_size = 2, max_cutoff = 0.99,
                    override_plot_dims = TRUE, 
                    colors = c("#332288", "#88CCEE", "#44AA99", "#117733", "#DDCC77", "#CC6677"),
                    features = selected_celltypes) +
  plot_layout(guides = "collect")

image-20240419183205407

semla还提供了细胞类型的共定位分析可视化函数。通过计算不同细胞类型之间的成对相关性,我们可以了解哪些细胞类型经常出现在相同的spot上:

代码语言:javascript
代码运行次数:0
运行
AI代码解释
复制
cor_matrix <- FetchData(semla.data, selected_celltypes) |> 
  mutate_all(~ if_else(.x<0.1, 0, .x)) |>  # Filter lowest values (-> set as 0)
  cor()

diag(cor_matrix) <- NA
max_val <- max(cor_matrix, na.rm = T)
cols <- RColorBrewer::brewer.pal(7, "RdYlBu") |> rev(); cols[4] <- "white"
pheatmap::pheatmap(cor_matrix, 
                   breaks = seq(-max_val, max_val, length.out = 100),
                   color=colorRampPalette(cols)(100),
                   cellwidth = 14, cellheight = 14, 
                   treeheight_col = 10, treeheight_row = 10, 
                   main = "Cell type correlation\nwithin spots")

image-20240419183405325

三. 总结

总的来说,NNLS解卷积算法运行速度非常快。R包semla还提供了一些额外的可视化方案,但NNLS解卷积算法的精确度尚待进一步探讨(特别是需要和其他解卷积算法进行benchmark分析)。

本文参与 腾讯云自媒体同步曝光计划,分享自微信公众号。
原始发表:2024-04-19,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 生信菜鸟团 微信公众号,前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体同步曝光计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
暂无评论
推荐阅读
编辑精选文章
换一批
整合单细胞和空转数据多种方法之Cell2location
除了我们上期介绍过的CellTrek算法【整合单细胞和空转数据多种方法之CellTrek】,还有非常多的算法用于整合单细胞和空转数据,如何快速系统的了解更多的主流整合算法呢?这么多种算法我们又应该选择哪种呢?答案是查阅综述+大量实践!
生信菜鸟团
2024/01/04
8.3K0
整合单细胞和空转数据多种方法之Cell2location
整合单细胞和空转数据多种方法之CellTrek
CellTrek发表于2022年的Nature Biotechnology,题为《Spatial charting of single-cell transcriptomes in tissues》。CellTrek可以结合单细胞和空间转录组数据准确地定位组织内单个细胞的位置,并构建空间细胞图谱。gitHub在https://github.com/navinlabcode/CellTrek
生信菜鸟团
2023/11/07
3K1
整合单细胞和空转数据多种方法之CellTrek
R语言之可视化①③散点图+拟合曲线目录
gene2) Pearson's product-moment correlation data: data gene1 and data$gene2 t = 2.4858, df = 395, p-value = 0.01334 95 percent confidence interval: 0.02600102 0.21984192 cor 0.1241053
用户1359560
2018/12/14
2.4K0
R语言之可视化①③散点图+拟合曲线目录
两次单细胞差异分析后的结果进行相关性散点图绘制
也就是说,它并不是拿两次差异分析各自统计学显著的基因的交集去绘图,而是把在两次差异分析至少有一次是统计学显著的基因拿过去。(说起来一样的绕口,让我们看看后面的代码)
生信技能树jimmy
2022/01/17
3.3K0
两次单细胞差异分析后的结果进行相关性散点图绘制
SPOTlight || 用NMF解卷积空间表达数据
空间解析基因表达谱是理解组织组织和功能的关键。然而,目前空间转录组分析技术(Spatial Transcriptomics,ST)尚未达到单细胞分辨率,往往需要结合单细胞RNA测序(scRNA-seq)信息来反褶积(或解卷积,Deconvolute )空间数据集。SPOTlight利用这两种数据类型的优势,能够将ST与scRNA-seq数据集成,从而推断出复杂组织中细胞类型和状态的位置。SPOTlight基于一个种子的非负矩阵因子分解回归(Seeded NMF regression ),使用细胞类型标记基因和非负最小二乘(NNLS)初始化,随后去卷积ST捕获位置(spot)。在作者的文章中,在示例数据人类胰腺癌中,成功地将患者切片划分为健康和癌区,并进一步精细绘制正常和肿瘤细胞状态。SPOTlight 流程如下:
生信技能树
2021/10/22
1.6K0
SPOTlight || 用NMF解卷积空间表达数据
Giotto|| 空间表达数据分析工具箱 Seurat 新版教程:分析空间转录组数据(上) Seurat 新版教程:分析空间转录组数据(下) scanpy教程:空间转录组数据分析 10X Visium:空间转录组样本制备到数据分析 空间信息在空间转录组中的运用 定量免疫浸润在单细胞研究中的应用
生信技能树jimmy
2021/01/12
3K0
SPOTlight || 用NMF解卷积空间表达数据
空转 | CellChat-V2,揭秘空间转录组数据的细胞通讯分析
仍然使用空转 | 结合scRNA完成空转spot注释(Seurat Mapping) & 彩蛋(封面的空转主图代码)推文中的空转数据进行示例展示。
生信补给站
2023/11/09
5.6K0
空转 | CellChat-V2,揭秘空间转录组数据的细胞通讯分析
空转|CARD-结合scRNA解决空间转录组spot注释,还能增强空间精度?!
前面已经介绍过了Seurat 空转 | 结合scRNA完成空转spot注释(Seurat Mapping) & 彩蛋(封面的空转主图代码)和 SPOTlight空转 | 我,SPOTlight,用解卷积,解决空间转录组spot注释! 联合单细胞进行空间转录组spot注释的方法,本文介绍下20202年发表于NBT的文献Spatially informed cell type deconvolution for spatial transcriptomics 的CARD方法。
生信补给站
2023/08/25
1.9K0
空转|CARD-结合scRNA解决空间转录组spot注释,还能增强空间精度?!
整合单细胞数据和Bulk数据的多种方法(二):R包Scissor
上个帖子介绍了整合单细胞数据和Bulk数据的最新方法scAB,整合单细胞数据和Bulk数据的多种方法(一):R包scAB。除了scAB以外,衔接单细胞数据和Bulk测序数据的方法还有Scissor, scPrognosis 和DEGAS 等工具。本系列文章计划依次介绍这些算法,提供给大家更多的选择。本文将继续介绍Scissor算法。
生信技能树
2023/02/27
4.5K0
整合单细胞数据和Bulk数据的多种方法(二):R包Scissor
空转 | 结合scRNA完成空转spot注释(Seurat Mapping) & 彩蛋(封面的空转主图代码)
前面空间转录组|没有单细胞数据如何做空转spot “注释”?文献和代码都给你!介绍了不结合单细胞进行空转spot的注释方式,文献中更多还是结合单细胞转录组数据进行spot注释 。整合scRNA-seq和空间转录组数据主要有(1)映射(Mapping)和(2)去卷积(Deconvolution)两种方法。本文介绍下使用Seurat的实现Mapping结合的方式 ,主要是将基于scRNA的细胞亚型定位到HPRI图谱上,后续会介绍SPOTlight ,CARD 等Deconvolution 的方式。
生信补给站
2023/08/25
4.2K2
空转 | 结合scRNA完成空转spot注释(Seurat Mapping) &  彩蛋(封面的空转主图代码)
空转反卷积方法哪家好:来看看主流的10种算法大比拼
空间转录组学技术主要分为两类:基于成像的方法和基于测序的方法。基于成像的方法(如smFISH、MERFISH和osmFISH)通过直接成像单个RNA分子,提供RNA表达水平的定量信息及其在单细胞中的空间定位。基于测序的方法(如10× Genomics Visium平台和Spatial Transcriptomics平台)则通过在带有位置条形码的逆转录引物上进行测序和计算重建,捕获组织样本中的基因表达。在相关文献中,“空间转录组学”(ST)是一个广义概念,而“Spatial Transcriptomics”则指特定的技术平台。
生信技能树
2025/07/08
1270
空转反卷积方法哪家好:来看看主流的10种算法大比拼
Seurat4.0系列教程22:空间转录组的分析、可视化与整合
Seurat4.0系列教程告一段落,但这决不是终点。这个系列教程只是给大家打开一扇窗,知道Seurat4.0有这些功能可用,少走弯路,起到一个抛砖引玉的作用,后续还要自己深入研究。不要像我当初入门单细胞之时,为了找到整合方法耗费大量时间及不必要的金钱。君子生非异也,善假于物也。学,然后知不足。加油吧!少年!用好单细胞测序及分析这个技术,为人类健康研究做出自己的贡献!不足之处,恳请批评指正!
生信技能树jimmy
2021/08/20
11.1K0
空间单细胞|基于图像的空间数据分析(2)
在这篇指南[1]中,我们介绍了Seurat的一个新扩展功能,用以分析新型的空间解析数据,将重点介绍由不同成像技术生成的三个公开数据集。
数据科学工厂
2024/07/05
6830
空间单细胞|基于图像的空间数据分析(2)
空转 | 我,SPOTlight,用解卷积,解决空间转录组spot注释!
前面介绍过整合scRNA-seq和空间转录组数据主要有(1)映射(Mapping)和(2)去卷积(Deconvolution)两种方法。前面空转 | 结合scRNA完成空转spot注释(Seurat Mapping) & 彩蛋(封面的空转主图代码)介绍了使用Seurat Mapping的方式进行spot注释,本文介绍一种经典的解卷积方法-SPOTlight。
生信补给站
2023/08/25
2.5K0
空转 | 我,SPOTlight,用解卷积,解决空间转录组spot注释!
空间转录组|没有单细胞数据如何做空转spot “注释”?文献和代码都给你!
在上文空间转录组|Load10X_Spatial函数修改适配多形式数据 + 空转标准流程得到空转数据T0的降维聚类结果后,面临和单细胞一样的注释问题。空间细胞类型注释主要分2大类:不联合单细胞数据 或者 联合单细胞数据。
生信补给站
2023/08/25
2.7K0
空间转录组|没有单细胞数据如何做空转spot “注释”?文献和代码都给你!
整合单细胞数据和Bulk数据的多种方法(一):R包scAB
最近刷到NAR最新的一篇算法文章《scAB detects multiresolution cell states with clinical significance by integrating single-cell genomics and bulk sequencing data》,DOI10.1093/nar/gkac1109。通过整合单细胞基因组学(包括scRNA-seq和sc-ATAC-seq数据)和Bulk RNA 测序数据及表型信息,在这里作者提出了一种称为 scAB 的新算法。
生信技能树
2023/02/27
4K1
整合单细胞数据和Bulk数据的多种方法(一):R包scAB
Seurat新版教程:分析空间转录组数据(下)
Seurat提供了两个工作流程来识别与组织空间位置相关的分子特征。第一种是根据组织内预先标注的解剖区域进行差异表达,这种差异表达可以通过非监督聚类或先验知识来确定。这种策略在这种情况下有效,因为上面的集群显示出明显的空间差异。
生信技能树jimmy
2020/03/27
3.1K0
知识积累----空转转录因子TF活性的计算框架
分析基于spot的ST数据集,以估计spot特异性TF活性,并识别与特定细胞类型、空间结构域、病理区域和配体/受体相关的TF。大多数基于spot的ST分析方法的一个关键限制是缺乏单细胞分辨率,因为空间RNA-seq测量结合了多个细胞。因此,推断TFs的点特异性活性来自多个细胞。为了克服这一限制,使用反卷积方法,通过计算整合空间RNA-seq与从scRNA-seq谱中获得的细胞类型的参考转录组特征来估计每个spot的细胞类型比例,通过线性回归进行TF和细胞类型之间的关联分析,通过反卷积方法从细胞类型比例单独预测每个TF的活性,以确定其活性在病理区域和空间域中不同的TF,并确定其基因表达与邻近位置的TF相关的配体/受体。
追风少年i
2024/09/05
1810
知识积累----空转转录因子TF活性的计算框架
方法分享--Spotiphy获取单细胞精度的空间数据(visium)
追风少年i
2025/03/14
970
方法分享--Spotiphy获取单细胞精度的空间数据(visium)
使用ggpubr包的stat_cor函数一步到位绘制相关性散点图并且添加统计学指标
再比如前面笔记两次单细胞差异分析后的结果进行相关性散点图绘制提到的两次差异分析结果的对比,就使用了ggpubr包的ggscatter函数绘制了相关性散点图:
生信技能树
2023/02/27
2.1K0
使用ggpubr包的stat_cor函数一步到位绘制相关性散点图并且添加统计学指标
推荐阅读
相关推荐
整合单细胞和空转数据多种方法之Cell2location
更多 >
领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档