利用monocle3分析单细胞数据

2022-10-25 19:42:39 浏览数 (1)

背景

开始monocle3案例数据学习。

文档：https://cole-trapnell-lab.github.io/monocle3/docs/clustering/

一、读入数据

代码语言：javascript复制
#1 读入数据
expression_matrix <- readRDS("monocle3/celegans/cao_l2_expression.rds")
cell_metadata <- readRDS("monocle3/celegans/cao_l2_colData.rds")
gene_annotation <- readRDS("monocle3/celegans/cao_l2_rowData.rds")

# Make the CDS object
cds <- new_cell_data_set(expression_matrix,
                         cell_metadata = cell_metadata,
                         gene_metadata = gene_annotation)

二、数据处理

数据处理包括表达数据标准化和批次效应的去除，对数据进行标准化使用 preprocess_cds函数，相当于 seurat 中 NormalizeData ScaleData RunPCA。

代码语言：javascript复制
cds <- preprocess_cds(cds, num_dim = 100)

num_dim 是降维的维数，默认为 100，这里是根据各组成部分解释的差异散点趋势图，也可以自己调整，建议安装默认维数进行降维。

代码语言：javascript复制
plot_pc_variance_explained(cds)

monocle3 数据处理

三、降维以及可视化

3.1 UMAP 方法

代码语言：javascript复制
cds <- reduce_dimension(cds) 
plot_cells(cds)
plot_cells(cds, color_cells_by="cao_cell_type")

#绘制特定基因在细胞中表达
plot_cells(cds, genes=c("cpna-2", "egl-21", "ram-2", "inos-1"))
代码语言：javascript复制
#可以使用多线程，前面读取好cds数据重新preprocess
cds <- preprocess_cds(cds, num_dim = 100)
cds <- reduce_dimension(cds, cores = 12) 
plot_cells(cds, color_cells_by="cao_cell_type")
#和单线程的结果看起来差别不大，但是提醒是有差别的。

3.2 tSNE 方法

代码语言：javascript复制
cds <- reduce_dimension(cds, reduction_method="tSNE")
plot_cells(cds, reduction_method="tSNE", color_cells_by="cao_cell_type", group_label_size=3.5,cell_size = 0.8)

3.3 去除批次效应

代码语言：javascript复制
plot_cells(cds, color_cells_by="plate", label_cell_groups=FALSE)
cds <- align_cds(cds, num_dim = 100, alignment_group = "plate")
cds <- reduce_dimension(cds)
plot_cells(cds, color_cells_by="plate", label_cell_groups=FALSE)

3.4 细胞聚类

代码语言：javascript复制
cds <- cluster_cells(cds, resolution=1e-5)
plot_cells(cds)
代码语言：javascript复制
#使用PAGA算法进行聚类
plot_cells(cds, color_cells_by="partition", group_cells_by="partition")

#根据细胞类型进行涂色
plot_cells(cds, color_cells_by="cao_cell_type")

#更改细胞标注策略，不对聚类的进行标注。
plot_cells(cds, color_cells_by="cao_cell_type", label_groups_by_cluster=FALSE)

四、 寻找 marker 基因

代码语言：javascript复制
#寻找每个cluster的marker基因
marker_test_res <- top_markers(cds, group_cells_by="partition", 
                               reference_cells=1000, cores=8)


#根据pseudo_R2方法对标记基因进行排序，取出每个cluster排名最高的基因
library(dplyr)
top_specific_markers <- marker_test_res %>% filter(fraction_expressing >= 0.10) %>% group_by(cell_group) %>% top_n(1, pseudo_R2)
top_specific_marker_ids <- unique(top_specific_markers %>% pull(gene_id))

#绘制每个cluster中不同marker基因的表达量
plot_genes_by_group(cds, top_specific_marker_ids, group_cells_by="partition",
                    ordering_type="maximal_on_diag",max.size=3)

#选取更多marker基因进行绘图
top_specific_markers <- marker_test_res %>%  filter(fraction_expressing >= 0.10) %>% group_by(cell_group) %>%  top_n(3, pseudo_R2)
top_specific_marker_ids <- unique(top_specific_markers %>% pull(gene_id))
plot_genes_by_group(cds,top_specific_marker_ids,group_cells_by="partition", ordering_type="cluster_row_col", max.size=3)

五、 对细胞进行注释

代码语言：javascript复制
#根据marker基因对细胞进行鉴定。
colData(cds)$assigned_cell_type <- as.character(partitions(cds))
#为每个cluster重新命名。
colData(cds)$assigned_cell_type <- dplyr::recode(colData(cds)$assigned_cell_type,
                                                 "1"="Body wall muscle",
                                                 "2"="Germline",
                                                 "3"="Motor neurons",
                                                 "4"="Seam cells",
                                                 "5"="Sex myoblasts",
                                                 "6"="Socket cells",
                                                 "7"="Marginal_cell",
                                                 "8"="Coelomocyte",
                                                 "9"="Am/PH sheath cells",
                                                 "10"="Ciliated neurons",
                                                 "11"="Intestinal/rectal muscle",
                                                 "12"="Excretory gland",
                                                 "13"="Chemosensory neurons",
                                                 "14"="Interneurons",
                                                 "15"="Unclassified eurons",
                                                 "16"="Ciliated neurons",
                                                 "17"="Pharyngeal gland cells",
                                                 "18"="Unclassified neurons",
                                                 "19"="Chemosensory neurons",
                                                 "20"="Ciliated neurons",
                                                 "21"="Ciliated neurons",
                                                 "22"="Inner labial neuron",
                                                 "23"="Ciliated neurons",
                                                 "24"="Ciliated neurons",
                                                 "25"="Ciliated neurons",
                                                 "26"="Hypodermal cells",
                                                 "27"="Mesodermal cells",
                                                 "28"="Motor neurons",
                                                 "29"="Pharyngeal gland cells",
                                                 "30"="Ciliated neurons",
                                                 "31"="Excretory cells",
                                                 "32"="Amphid neuron",
                                                 "33"="Pharyngeal muscle")
#重新绘图
plot_cells(cds, group_cells_by="partition", color_cells_by="assigned_cell_type")

每个 cluster 重新命名之后进行绘图

写在最后：有时间我们会努力更新的。大家互动交流可以前去论坛，地址在下面，复制去浏览器即可访问，弥补下公众号没有留言功能的缺憾。原地址暂未启用（bioinfoer.com）。

代码语言：javascript复制
sx.voiceclouds.cn

有些板块也可以预设为大家日常趣事的分享等，欢迎大家来提建议。

数据处理 com marker

0 人点赞