跟着Cell学单细胞转录组分析(十):Monocle2拟时分析演

上节完成了拟时分析，剩下的内容就比较简单了，只需要可视化。我们可视化一些一般文章中出现的图。用到的文件是上节分析得到的cds文件。图的主题可以结合ggplot2进行修饰。

先做一个细胞群的谱系分化图。从这个图可以看出我们关注的细胞分化轨迹。

library(ggsci)
plot_cell_trajectory(cds, color_by = "Cluster")  + scale_color_nejm()

图片

plot_cell_trajectory(cds, color_by = "State")  + scale_color_npg()

图片

plot_cell_trajectory(cds, color_by = "Pseudotime")

image.png

如果细胞轨迹全部在一起，很难看出不同细胞状态在分支上的位置，这时，我们可以将每个状态单独画出来，看起来比较清晰。

plot_cell_trajectory(cds, color_by = "State") +
  facet_wrap(~State, nrow = 1)

图片

处理细胞谱系拟时可视化，我们还关注分化轨迹过程中基因的情况。选定关注的基因，看看其在拟时中的表达。

pData(cds)$TGFBR2 = log2( exprs(cds)['TGFBR2',]+1)

图片

通过拟时基因表达模式聚类。

cds$id <- rownames(cds)
library(dplyr)
cds %>% arrange(qval) %>% head(10) %>% select(id) -> gene_to_cluster
gene_to_cluster <- gene_to_cluster$id
my_pseudotime_cluster <- plot_pseudotime_heatmap(cds[gene_to_cluster,],
                                                 num_clusters = 3,
                                                 cores = 8,
                                                 show_rownames = TRUE)

image.png

BEAM进行统计分析。

BEAM_res <- BEAM(my_cds_subset, branch_point = 1, cores = 8)
BEAM_res <- BEAM_res[order(BEAM_res$qval),]
BEAM_res <- BEAM_res[,c("gene_short_name", "pval", "qval")]
head(BEAM_res)
table(BEAM_res$qval < 1e-4)
plot_genes_branched_heatmap(my_cds_subset[row.names(subset(BEAM_res, qval < 1e-4)),],
                            branch_point = 1,
                            num_clusters = 4,
                            cores = 8,
                            use_gene_short_name = TRUE,
                            show_rownames = TRUE)

图片

拟时分析的内容很丰富，也很多，在不同的研究中有不同的意义，这里只是简单展示了几种常见的可视化结果，对于结果的解读，应用，还需要结合具体的生物学意义，通过推断与生物学背景结合，才能让分析彰显意义。