前言

单细胞的拟时序分析是利用算法模拟不同细胞之间的分化趋势，目前已经是一种相对常规化的分析手段。通过拟时序分析，可以模拟不同细胞之间的潜在分化联系，因此可以起到动态预测的作用。最近的一些文章开始利用模拟的分化路线关联基因表达和功能通路的变化，这样也能动态反应不同细胞之间功能状态的变化趋势。

今天Immugent就来用代码实操的方式，来介绍如何使用SCPA将功能通路和拟时序分化路线联合分析。

代码实操

library(SCPA)
library(Seurat)
library(tidyverse)
library(magrittr)
library(dyno)
library(ComplexHeatmap)
library(circlize)

naive_cd4 <- load("naive_cd4.rds")
DimPlot(naive_cd4, split.by = "Hour")

图片

naive_cd4 <- subset(naive_cd4, idents = "Treg", invert = T)

df <- as.matrix(naive_cd4[["RNA"]]@data)
var_genes <- names(sort(apply(df, 1, var), decreasing = TRUE))[1:1000]

counts <- Matrix::t(as(as.matrix(naive_cd4@assays$RNA@counts[var_genes,]), 'sparseMatrix'))
expression <- Matrix::t(as(as.matrix(naive_cd4@assays$RNA@data[var_genes,]), 'sparseMatrix'))

dataset_n4 <- wrap_expression(expression = expression,
                              counts = counts)

图片

model_n4 <- infer_trajectory(dataset_n4, method = ti_slingshot(), verbose = T)

plot_dimred(model_n4, 
            "pseudotime", 
            pseudotime = calculate_pseudotime(model_n4), 
            hex_cells = F,
            plot_trajectory = T, 
            size_cells = 1, alpha_cells = 0.8) + 
  theme(aspect.ratio = 1)

plot_dimred(model_n4, 
            grouping = group_onto_nearest_milestones(model_n4), 
            hex_cells = F,
            plot_trajectory = T, 
            size_cells = 1, alpha_cells = 0.8) + 
  theme(aspect.ratio = 1)

图片

mile_group <- data.frame(group_onto_nearest_milestones(model_n4)) %>%
  set_colnames("milestone") %>%
  rownames_to_column("cell")
  
naive_cd4$milestone <- mile_group$milestone
cd4_pseudo <- list()
for (i in 1:max(mile_group$milestone)) {
  cd4_pseudo[[i]] <- seurat_extract(naive_cd4, meta1 = "milestone", value_meta1 = i)
}

pathways <- "combined_metabolic_pathways.csv"

cd4_metabolism <- compare_pathways(samples = cd4_pseudo, 
                                   pathways = pathways)
                                   
cd4_metabolism <- cd4_metabolism %>%
  data.frame() %>%
  select(Pathway, qval) %>%
  column_to_rownames("Pathway")

col_hm <- colorRamp2(colors = c("white", "red"), breaks = c(0, max(mstone_out)))

Heatmap(t(cd4_metabolism),
        name = "Qvalue",
        col = col_hm,
        border = T,
        rect_gp = gpar(col = "white", lwd = 0.1),
        heatmap_height = unit(2, "cm"),
        show_column_dend = F,
        show_row_names = F,
        show_column_names = F)                               

图片

plot_rank(scpa_out = cd4_metabolism, 
          pathway = "hallmark_gly",
          base_point_size = 2.5,
          highlight_point_size = 3)

图片

小结

虽然拟时序分析得到的是一个连续的值，但是我们在实际分析时可以根据需要将其分成不同阶段，如：早期，中期和晚期。但是单纯的以时间分期是没有太大的意义，需要联合不同功能定义出的细胞亚群。就像T细胞可以分成naive, effector以及memory阶段一样，这样就赋予了拟时间实际的生物学意义。在此基础上，算出不同阶段之间功能通路的差异，我们就能得知每个阶段主导的功能通路。

好啦，本期分享到此就结束了，我们下期再会~~