单细胞数据挖掘实战：文献复现（二）批量创建Seurat对象及质控
单细胞数据挖掘实战：文献复现（七）MG 和 Mo/MΦ 评分
这篇文献的第四个结果就是观察了marker基因在Hom-MG、 Act-MG 和 Mo/MΦ 细胞中的表达情况，结果主要体现在Fig. 4中。
一、将cell_type由三类细分为四类
- Hom-MG： h.Microglia <- ctrl_Microglia
- Act-MG: a.Microglia <- tumor_Microglia
- Mo/MΦ: Macrophages
- BAM
table(Idents(seu_object))
# MG Mo/MΦ BAM
#31959 5624 1680
#选出"MG","Mo/MΦ"
seu_object <- subset(seu_object, idents = c("MG","Mo/MΦ"))
table(Idents(seu_object))
# MG Mo/MΦ
#31959 5624
#将Microglia细分为"h.Microglia"和"a.Microglia"
Idents(seu_object) <- seu_object$seurat_clusters
seu_object$cell_type_4_groups <- plyr::mapvalues(Idents(seu_object),
from=c(0:4,6:11,13,14,16,18,19),
to = c("h.Microglia", "h.Microglia", "a.Microglia", "h.Microglia", "Macrophages", "h.Microglia",
"h.Microglia", "Macrophages","a.Microglia", "a.Microglia", "Macrophages","a.Microglia",
"h.Microglia", "Macrophages" ,"Macrophages", "Macrophages"))
# Figure 4a
# left panel
DimPlot(seu_object, group.by = "cell_type_4_groups")
1.png
二、Act-MG 和 Mo/MΦ 中差异上调基因的表达水平
数据处理
Idents(seu_object) <- seu_object$cell_type_4_groups
#markers_ActMG_vs_HomMG
markers_ActMG_vs_HomMG <- FindMarkers(object = seu_object, ident.1 = "a.Microglia",
ident.2 = "h.Microglia", only.pos = TRUE, min.pct = 0.25,
logfc.threshold = 0.25)
markers_ActMG_vs_HomMG$gene <- rownames(markers_ActMG_vs_HomMG)
#markers_MoM_vs_ActMG
markers_MoM_vs_ActMG <- FindMarkers(object = seu_object, ident.1 = "Macrophages",
ident.2 = "a.Microglia", only.pos = TRUE, min.pct = 0.25,
logfc.threshold = 0.25)
markers_MoM_vs_ActMG$gene <- rownames(markers_MoM_vs_ActMG)
画图
Fig. 4b
#画图
genes <- unique(c(markers_ActMG_vs_HomMG$gene, markers_MoM_vs_ActMG$gene))
gene_expression_data <- GetAssayData(object = seu_object, slot = "data")
gene_expression_data <- as.data.frame(t(gene_expression_data[genes, ]))
common_genes <- intersect(markers_ActMG_vs_HomMG$gene, markers_MoM_vs_ActMG$gene)
markers_ActMG_only <- markers_ActMG_vs_HomMG$gene[!(markers_ActMG_vs_HomMG %in% common_genes)]
markers_MoM_only <- markers_MoM_vs_ActMG$gene[!(markers_MoM_vs_ActMG$gene %in% common_genes)]
genes <- unique(c(markers_ActMG_vs_HomMG$gene, markers_MoM_vs_ActMG$gene))
cell_types_selected <- c("h.Microglia", "a.Microglia", "Macrophages")
names(cell_types_selected) <- c("h.Microglia", "a.Microglia", "Macrophages")
genes_mean_expr <- sapply(cell_types_selected, function(cell_type) {
Matrix::rowMeans(seu_object@assays$RNA@data[genes,
colnames(seu_object)[seu_object$cell_type_4_groups == cell_type]],
na.rm = T)
})
colnames(genes_mean_expr)[1:3]<-c("Hom_MG", "Act_MG", "MoMphi")
genes_mean_expr <- as.data.frame(genes_mean_expr)
genes_mean_expr$gene <- rownames(genes_mean_expr)
labels<-c("Ly6c2", "Ccl5", "Ly6i", "Lyz2", "Lgals3","Ifitm2",
"Tgfbi","Tmsb10", "Il1rn","Ass1","Ifitm3","Il1b",
"Irf7","Ccr2", "H2-Aa","H2-Ab1", "H2-Eb1","Cd74","Ifit3",
"Il18bp", "Mif", "Apoe", "Stat1", "Ccl12", "H2-D1", "Ccl4", "Ccl3", "Ly86")
genes_mean_expr$color <- 0
genes_mean_expr[genes_mean_expr$gene %in% markers_ActMG_only, "color"]<-"Act-MG"
genes_mean_expr[genes_mean_expr$gene %in% markers_MoM_only, "color"]<-"MoM"
genes_mean_expr[genes_mean_expr$gene %in% common_genes, "color"]<-"common"
genes_mean_expr$color<-factor(genes_mean_expr$color, levels=c("Act-MG", "MoM", "common"))
genes_mean_expr_labeled <- genes_mean_expr[labels, ]
col_Macro<-"#FABF00"
col_ActMG<-"#2F8EA1"
pdf(file = "fig4_scat.pdf",width = 10, height = 10)
ggplot(genes_mean_expr, aes(x=Act_MG, y=MoMphi))+
geom_jitter(aes(fill=color),shape=21, color="white", alpha=0.7, size=4)+
geom_text_repel(data=genes_mean_expr_labeled, aes(label=gene), nudge_y=0.2, size=5,
direction="both")+
geom_abline(intercept=0, slope=1)+
scale_fill_manual(values=c(col_ActMG, col_Macro, "black"))+
xlim(0,5)+
ylim(0,5)+
xlab("Act-MG")+
ylab("MoMphi")+
coord_fixed()+
theme_bw(base_size = 18)+
theme(panel.grid = element_blank())
dev.off()
2.png
Fig. 4c
Hom-MG 与 Act-MG 和 Act-MG 与 Mo/MΦ 中前 25 个上调基因的表达情况
markers_MoM_vs_ActMG <- markers_MoM_vs_ActMG[order(-markers_MoM_vs_ActMG$avg_log2FC),]
markers_ActMG_vs_HomMG <- markers_ActMG_vs_HomMG[order(-markers_ActMG_vs_HomMG$avg_log2FC),]
genes_for_heatmap <- unique(c(markers_MoM_vs_ActMG$gene[1:25], markers_ActMG_vs_HomMG$gene[1:25]))
genes_mean_expr_heatmap <- genes_mean_expr[genes_for_heatmap, ]
mat_breaks <- seq(0, abs(max(genes_mean_expr_heatmap[,1:3])), length=51)
library(pheatmap)
pheatmap(
mat = genes_mean_expr_heatmap[, 1:3],
color = colorRampPalette(rev(c("#810f7c", "#8856a7", "#8c96c6", "#b3cde3", "#edf8fb")))(length(mat_breaks) - 1),
breaks = mat_breaks,
border_color = NA,
cluster_cols = F,
cluster_rows = F,
show_colnames = TRUE,
show_rownames = TRUE,
treeheight_col = 0,
treeheight_row = 0,
gaps_row = 25,
drop_levels = TRUE,
fontsize = 10,
angle_col = 0,
main = "Top upregulated genes in Act-MG and MoM"
)
3.png
往期单细胞数据挖掘实战：
