拟时分析

生信喵 发表于 2022-8-8 02:34:22

构建细胞分化轨迹
一、Trajectory analysis 轨迹推断
   轨迹推断（Trajectory Inference，TI），也称为细胞分析轨迹（differentiation trajectories）基于单细胞转录组数据，利用模型预测细胞分化过程。
   由于分析结果是基于某个组织一瞬间的群体单细胞数据预测得到的，并非真实时间维度上的分化结果，因此也被称为伪时序分析 (pseudotime analysis)，也可以叫做拟时分析。
   拟时分析通过分析关键基因的表达模式，将所有细胞按照发育时间的先后排布在拟时间轴上，模拟发育过程中的细胞分化过程。通过对细胞轨迹的分析，我们可以挖掘出细胞分化过程中经历的细胞类型变化、伴随发育过程变化的动态变化基因、祖细胞不同的分化命运等与生命发育息息相关的信息。

二、monocle3 进行细胞轨迹推断
文档：https://cole-trapnell-lab.github.io/monocle3/docs/trajectories/
#1 读取数据
getwd()
setwd('/home/xhs/jyxy/16.rscrna/')
expression_matrix <- readRDS("./monocle3/packer/packer_embryo_expression.rds")
cell_metadata <- readRDS("./monocle3/packer/packer_embryo_colData.rds")
gene_annotation <- readRDS("./monocle3/packer/packer_embryo_rowData.rds")

cds <- new_cell_data_set(expression_matrix,cell_metadata = cell_metadata,gene_metadata = gene_annotation)

#2 数据处理以及去除批次效应
cds <- preprocess_cds(cds, num_dim = 50)
cds <- align_cds(cds, alignment_group = "batch", residual_model_formula_str = "~ bg.300.loading + bg.400.loading + bg.500.1.loading + bg.500.2.loading + bg.r17.loading + bg.b01.loading + bg.b02.loading")

#3 数据降维以及结果可视化
cds <- reduce_dimension(cds)
#设置group_label_size = 3,cell_size = 0.7调整文本和点的大小
plot_cells(cds, label_groups_by_cluster=FALSE,color_cells_by = "cell.type",group_label_size = 3,cell_size = 0.7)

#基因沿着轨迹变化
ciliated_genes <- c("che-1","hlh-17","nhr-6","dmd-6","ceh-36","ham-1")
plot_cells(cds,genes=ciliated_genes,label_cell_groups=FALSE,show_trajectory_graph=FALSE)

#4 聚类
cds <- cluster_cells(cds)
plot_cells(cds, color_cells_by = "partition")

#5 构建Learn the trajectory graph
cds <- learn_graph(cds)
plot_cells(cds,color_cells_by = "cell.type",label_groups_by_cluster=FALSE,label_leaves=FALSE,label_branch_points=FALSE)

#6排序 Order the cells in pseudotime
plot_cells(cds,color_cells_by = "embryo.time.bin",label_cell_groups=FALSE,label_leaves=TRUE,label_branch_points=TRUE, graph_label_size=1.5)

#7 选取root跟节点
cds <- order_cells(cds)

plot_cells(cds,
      color_cells_by = "pseudotime",
      label_cell_groups=FALSE,
      label_leaves=FALSE,
      label_branch_points=FALSE,
      graph_label_size=1.5)

#8 自定义函数，确定根节点
get_earliest_principal_node <- function(cds, time_bin="130-170"){
cell_ids <- which(colData(cds)[, "embryo.time.bin"] == time_bin)
closest_vertex <-
cds@principal_graph_aux[["UMAP"]]$pr_graph_cell_proj_closest_vertex
closest_vertex <- as.matrix(closest_vertex)
root_pr_nodes <-
igraph::V(principal_graph(cds)[["UMAP"]])$name[as.numeric(names
                                                         (which.max(table(closest_vertex))))]
root_pr_nodes
}
cds <- order_cells(cds, root_pr_nodes=get_earliest_principal_node(cds))
plot_cells(cds,color_cells_by = "pseudotime",label_cell_groups=FALSE,label_leaves=FALSE,label_branch_points=FALSE,graph_label_size=1.5)

#9 选取树枝Subset cells by branch
cds_sub <- choose_graph_segments(cds)

#10 绘制3D 轨迹图
cds_3d <- reduce_dimension(cds, max_components = 3)
cds_3d <- cluster_cells(cds_3d)
cds_3d <- learn_graph(cds_3d)
cds_3d <- order_cells(cds_3d, root_pr_nodes=get_earliest_principal_node(cds))

cds_3d_plot_obj <- plot_cells_3d(cds_3d, color_cells_by="partition")

页: [1]

生信人's Archiver

拟时分析