拟时序分析就是差异分析的细节剖析-linux运维-PHP中文网

拟时序分析就是差异分析的细节剖析

蓮花仙者

发布： 2025-07-09 11:42:01

原创

398人浏览过

关于大样本量差异分析是否可以转为拟时序分析，以及两个分组的差异分析仅仅是上下调的问题，以下是经过伪原创处理后的文章：

许多朋友在后台表示对单细胞数据分析中的拟时序分析感到困惑。恰逢最近看到了一篇清晰展现拟时序分析重要性的文献，这里与大家分享。这篇文献完美地展示了为什么单纯的差异分析不足以揭示全部信息，而拟时序分析则是对差异分析的深入剖析。

这篇发表在NATURE COMMUNICATIONS | (2021) 的文章：《CD177 modulates the function and homeostasis of tumor-infiltrating regulatory T cells》，链接是：https://www.php.cn/link/1ee9bee2c7227c35ac1ca90f2e4fb172：

研究开始时使用了13,433个外周血单核细胞（PB）和12,239个肿瘤浸润细胞（TI），通过降维聚类分群后，根据FOXP3和CD25 (IL2RA)定位到了Treg亚群，分别是160个PB和574个TI Treg细胞。

拟时序分析就是差异分析的细节剖析降维聚类分群

如果直接对这两个分组进行差异分析，可以识别出273个差异表达基因（DEGs）（Log fold-change > 1, adjusted p-value < 0.05）。

此外，作者在自己的ccRCC单细胞矩阵以及一个公共数据集HCC中也进行了类似的差异分析，并筛选出共有基因：

拟时序分析就是差异分析的细节剖析差异基因及其交集

这样的差异分析虽然进行了交集筛选，但仍遗漏了许多细节，仅得到基因的上下调信息，而忽略了每个基因在两个单细胞亚群中具体的渐变趋势。

由于细胞数量并不多，运行Monocle 2算法，可以清楚地看到上下调的变化趋势，甚至发现隐藏的变化模式：

首先是拟时序分析的轨迹流形图，展示了ccRCC中Treg细胞的轨迹，使用Monocle 2算法。实线和虚线代表由表达谱定义的不同细胞轨迹/命运。

拟时序分析就是差异分析的细节剖析拟时序分析结果

虽然拟时序分析只展示了上述的一个图，但具体的代码步骤还是有一定难度的。以下是可直接使用的代码示例，我们以SeuratData包中的pbmc3k数据集为例，主要是将Seurat包的对象转换为monocle中的单细胞对象：

代码语言：javascript

library(SeuratData) # 加载seurat数据集
getOption('timeout')
options(timeout=10000)
#InstallData("pbmc3k")
data("pbmc3k")
sce <- as.SingleCellExperiment(pbmc3k[, 1:10, ])
sc_cds <- importCDS(sce, import_all = TRUE)
cds <- estimateSizeFactors(cds)
cds <- estimateDispersions(cds)
cds <- detectGenes(cds, min_expr = 0.1)
cds <- setOrderingFilter(cds, ordering_genes = row.names(subset(fData(cds), num_cells_expressed >= 10)))
cds <- reduceDimension(cds, max_components = 2, method = 'DDRTree')
cds <- orderCells(cds)
save(cds, file = 'input_cds.Rdata')

登录后复制

确保前面的步骤无误后，接下来就可以运行拟时序分析的主程序：

代码语言：javascript

rm(list=ls())
options(stringsAsFactors = F)
library(monocle)
library(Seurat)
load(file = 'input_cds.Rdata')
# 接下来很重要，到底是看哪个性状的轨迹
colnames(pData(cds))
table(pData(cds)$Cluster)
table(pData(cds)$Cluster, pData(cds)$celltype)
plot_cell_clusters(cds, 1, 2)
## 我们这里并不能使用 monocle的分群
# 还是依据前面的 seurat分群, 其实取决于自己真实的生物学意图
pData(cds)$Cluster = pData(cds)$celltype
table(pData(cds)$Cluster)
Sys.time()
diff_test_res <- differentialGeneTest(cds, fullModelFormulaStr = "~Cluster")
save(diff_test_res, cds, file = 'output_of_phe2_monocle.Rdata')

登录后复制

有了上面的output_of_phe2_monocle.Rdata文件后，就是拟时序分析的结果，接下来就可以进行各种可视化。首先看看每个细胞的时序信息：

LimeSurvey在线问卷管理系统

LimeSurvey是一款在线问卷管理系统,具有问卷的设计、修改、发布、回收和统计等多项功能。同时它也是一个开源软件,其最新版本的软件包可以完全免费获取和使用。它集成了调查程序开发、调查问卷的发布以及数据收集等功能,使用它,用户不必了解这些功能的编程细节。网上收集的调查数据可以导出多种文件格式以便分析，例如 spss数据格式 *.dat文件。

198

查看详情

代码语言：javascript

## 后面是对前面的结果进行精雕细琢
rm(list=ls())
options(stringsAsFactors = F)
library(Seurat)
library(gplots)
library(ggplot2)
library(monocle)
library(dplyr)
load(file = 'output_of_phe2_monocle.Rdata')
cds = my_cds_subset
colnames(pData(cds))
table(pData(cds)$State, pData(cds)$Cluster)
library(ggsci)
p1 = plot_cell_trajectory(cds, color_by = "Cluster") + scale_color_nejm()
p1
ggsave('trajectory_by_cluster.pdf')
plot_cell_trajectory(cds, color_by = "celltype")
p2 = plot_cell_trajectory(cds, color_by = "Pseudotime")
p2
ggsave('trajectory_by_Pseudotime.pdf')
p3 = plot_cell_trajectory(cds, color_by = "State") + scale_color_npg()
p3
ggsave('trajectory_by_State.pdf')
library(patchwork)
p1 + p2 / p3

登录后复制

可视化结果如下：

拟时序分析就是差异分析的细节剖析拟时序分析的多元化结果