monocle 2拟时序分析

monocle做拟时序分析首先要构建CDS需要3个矩阵:expr.matrix、pd、fd,其次将Seurat中的对象转换为monocle识别的对象。然后选择想要做拟时序依据的基因就可以了,如果已知开始和结束的细胞,将过程开始时收集的细胞与结束时收集的细胞简单地进行比较,并找到差异表达的基因,做拟时序依据的基因,根据时间点的差异分析选择基因通常非常有效,但是如果我们没有时间序列数据,可以选择离散度和表达量高的基因。

library(monocle)  packageVersion("monocle")  #monocle构建CDS需要3个矩阵:expr.matrix、pd、fd  # 将Seurat中的对象转换为monocle识别的对象  #cds <- importCDS(GetAssayData(seurat.object))  #选择做拟时序的亚群  Mono_tj<-subset(seurat.object, idents = c(1,2,4,6,7))    Mono_matrix<-as(as.matrix(GetAssayData(Mono_tj,slot = "counts")), 'sparseMatrix')  #构建featuredata,一般featuredata需要两个col,一个是gene_id,一个是gene_short_name,row对应counts的rownames  feature_ann<-data.frame(gene_id=rownames(Mono_matrix),gene_short_name=rownames(Mono_matrix))  rownames(feature_ann)<-rownames(Mono_matrix)  #  Mono_fd<-new("AnnotatedDataFrame", data = feature_ann)  #  #Seurat object中的@meta.data一般会存放表型相关的信息如cluster、sample的来源、group等,所以选择将metadata转换为phenodata  sample_ann<[email protected]  #rownames(sample_ann)<-colnames(Mono_matrix)    Mono_pd<-new("AnnotatedDataFrame", data =sample_ann)  #build new cell data set  Mono.cds<-newCellDataSet(Mono_matrix,phenoData =Mono_pd,featureData =Mono_fd,expressionFamily=negbinomial.size())    #查看phenodata、featuredata  head(pData(Mono.cds))  head(fData(Mono.cds))  #预处理  Mono.cds <- estimateSizeFactors(Mono.cds)  Mono.cds <- estimateDispersions(Mono.cds)  #筛选基因,这里可以根据自己的需要筛选特定的基因  disp_table <- dispersionTable(Mono.cds)  unsup_clustering_genes <- subset(disp_table, mean_expression >= 0.1)  Mono.cds <- setOrderingFilter(Mono.cds, unsup_clustering_genes$gene_id)  #用DDRtree 进行降维分析  Mono.cds <- reduceDimension(    Mono.cds,    max_components = 2,    method = 'DDRTree')  #计算psudotime值  Mono.cds <- orderCells(Mono.cds)  head(pData(Mono.cds))    plot_cell_trajectory(Mono.cds,cell_size = 1)

image.png

plot_cell_trajectory(Mono.cds, color_by = "Pseudotime")

image.png

plot_cell_trajectory(Mono.cds, color_by = "seurat_clusters",cell_size = 1)

image.png

欢迎关注~

参考:http://cole-trapnell-lab.github.io/monocle-release/docs/