生信技能树单细胞数据挖掘笔记(5)

生信技能树单细胞数据挖掘笔记(1)：数据读入

生信技能树单细胞数据挖掘笔记(2)：创建Seurat对象并进行质控、筛选高变基因并可视化

生信技能树单细胞数据挖掘笔记(3)：降维与聚类

生信技能树单细胞数据挖掘笔记(4)：其他分析（周期判断、double诊断、细胞类型注释）

生信技能树单细胞数据挖掘笔记(5)：轨迹分析

轨迹分析

通过构建细胞间的变化轨迹来重塑细胞随时间的变化过程。

##4、拟时轨迹分析----
library(monocle)
data <- as(as.matrix(scRNA@assays$RNA@counts), 'sparseMatrix')
# count矩阵
pd <- new('AnnotatedDataFrame', data = scRNA@meta.data)
# meta表转成特定格式
fData <- data.frame(gene_short_name = row.names(data), row.names = row.names(data))
fd <- new('AnnotatedDataFrame', data = fData)
# 基因名表转成特定格式
mycds <- newCellDataSet(data,
                        phenoData = pd,
                        featureData = fd,
                        expressionFamily = negbinomial.size())
#expressionFamily参数用于指定表达矩阵的数据类型，有几个选项可以选择：
#稀疏矩阵用negbinomial.size()，FPKM值用tobit()，logFPKM值用gaussianff()
##save(mycds, file = "../../tmp/mycds_raw.Rdata")
##rm(list = ls())
##load("../../tmp/mycds_raw.Rdata")
#library("monocle")
mycds <- estimateSizeFactors(mycds)
mycds <- estimateDispersions(mycds, cores=4, relative_expr = TRUE) #!
#完成数据导入和预处理后，就可以考虑选择特定基因代表细胞的发育特征

读入maker基因

mycds <- estimateSizeFactors(mycds)
mycds <- estimateDispersions(mycds, cores=4, relative_expr = TRUE) #!
#完成数据导入和预处理后，就可以考虑选择特定基因代表细胞的发育特征
#这里可以选取我们之前挑选的marker gene
load("../markergene.Rdata")
markers.gene <- all.markers$gene
mycds <- setOrderingFilter(mycds, markers.gene)

#降维
mycds <- reduceDimension(mycds, max_components = 2, method = 'DDRTree') #!
#耗时，耗内存

细胞排序

#降维
mycds <- reduceDimension(mycds, max_components = 2, method = 'DDRTree') #!
#耗时，耗内存
#排序
mycds <- orderCells(mycds)
#save(mycds,file = "../../tmp/mycds_reduced.Rdata")
load("../../tmp/mycds_reduced.Rdata")

p1 <- plot_cell_trajectory(mycds, color_by = "seurat_clusters")
ggsave("../../out/3.4trajectory_1.pdf", plot = p1) 
p2 <- plot_cell_trajectory(mycds, color_by = "State")
ggsave("../../out/3.4trajectory_2.pdf", plot = p2) 

by_cluster
by_state

---

转载请注明：周小钊的博客>>>生信技能树单细胞数据挖掘笔记(5)