科研视角----空间转录组NMF分析的生物学解读

原创

追风少年i

发布于 2024-12-18 14:03:30

4390

作者，Evil Genius

上一篇脚本更新---NMF在空间转录组的运用我们分享了空间转录组如何进行NMF分析，但是光分析还不够，分析结果如何进行生物学解读才是重点。

这一篇我们就是来解读分析结果。比较简单的内容。

关键点1、使用非负矩阵分解（NMF）将空间数据解卷积为30个“因子”（描述相关基因、途径和细胞类型）。每个因子都由共变基因的表达驱动。因子显示了不同细胞类型和结构的特征，包括混合细支气管上皮细胞类型（因子1 (F1)），平滑肌细胞（F10）和血浆和B细胞（F6）。

这个地方因子的分布如何关联细胞类型的？需要大家考虑其中的运用。而且这个地方也说明了一点，NMF针对空间的细胞注释有帮助。

因子活性揭示通路和细胞相互作用

不同样本间的因子分布差异，与不同样本的生物学功能差异有关。与HC相比，IPF中有11个因子更为普遍。这些因子与重要的IPF细胞形态和过程相关，包括ECM相关通路，以及纤维化（F4和F14）或IPF“蜂窝”形成（F5和F21）的重叠区域。F5显示树突状细胞和巨噬细胞标记物，而F21显示表达MUC5B的气道上皮特征，F21的特征可能反映了先前鉴定的MUC5B+、BPIFB1+和SCGB3A1+ ipf相关细胞群。F9是IPF肺特有的，其特征是与氧化应激、炎症、ECM重塑和血管改变相关的基因。此外，它特别存在于密度较小的纤维化IPF样本中（“B1”），提示参与早期纤维化过程。在ECM/纤维化相关因子中，F4和F14在IPF组织中都有明显的活性，可能反映了病理性重塑。F14活性与推断的KRT5−/KRT17+ AbBa细胞、肌成纤维细胞和最近描述的HAS1-hi成纤维细胞亚型5的细胞类型密度相关，特别是在IPF样本中。

因子如何与细胞类型相关联？最后我们附上分析的代码

准备好东西，1、基础分析的rds文件；2、单细胞空间联合分析的细胞矩阵。

library(STutility)
library(harmony)
library(RColorBrewer)
library(patchwork)
library(dplyr)
library(tidyr)
library(purrr)
library(broom)
library(ggpmisc)

adata = readRDS('test.rds')
####读取单细胞空间联合分析的矩阵，推荐cell2location
cell_anno <- read.csv('联合矩阵',row.names = 1)
head(cell_anno)

##### Add metadata to se obj ####
adata <- AddMetaData(adata, cell_anno)

细胞与因子的相关性

###### Cell~Factor correlation ######

###### Cell~Factor correlation ######
groups_use <- c("All", "Healthy", "IPF_1", "IPF_2", "IPF_3", "IPF_4")   ###选取样本

factor_cell_cor_data_group <- setNames(lapply(groups_use, function(g){
  if (g=="All") {
    factor_cell_cor_subset <- factor_cell_cor_data
  } else if (g=="Healthy") {
    factor_cell_cor_subset <- subset(factor_cell_cor_data, condition == "control")
  } else {
    factor_cell_cor_subset <- subset(factor_cell_cor_data, subject_alias == g)
  }
  factor_cell_density_hm <- cor(factor_cell_cor_subset[,c(new_c2l_names_filtered, factor_names)])
  diag(factor_cell_density_hm) <- 0
  factor_cell_density_hm <- as.data.frame(factor_cell_density_hm)
  factor_cell_density_hm <- factor_cell_density_hm[factor_names, new_c2l_names_filtered]
  }),
  nm = groups_use)

saveRDS(factor_cell_cor_data_group, file.path(DIR_OBJ_C2L, "hs_visium_A_nmf_1-30_factor_cell_cor_data_group.rds"))


# Plot pheatmap
pal_length <- 11
ph_colors <- RColorBrewer::brewer.pal(pal_length, "RdBu") %>% rev()
abs_max_val <- lapply(factor_cell_cor_data_group, max) %>% unlist() %>% max()

pdf(file = file.path(DIR_FIG_C2L, "hs_visium_c2l_res_cell_type_A-NMF30_factor_cor_heatmap_per_subject.pdf"), width = 10, height = 10, useDingbats = F)
for(g in groups_use){
  factor_cell_density_hm <- factor_cell_cor_data_group[[g]]

  if(g=="all"){abs_max_val <- max(cell_density_hm)}
  ph_breaks <- c(seq(-abs_max_val, 0, length.out=ceiling(pal_length/2) + 1),
                 seq(abs_max_val/pal_length, abs_max_val, length.out=floor(pal_length/2)))

  print(
    pheatmap::pheatmap(t(factor_cell_density_hm),
                       cellwidth = 10,
                       cellheight = 10,
                       color = ph_colors,
                       breaks = ph_breaks,
                       main = g)
  )
}
dev.off()