❝1、《KS科研分享与服务》公众号有QQ交流群,进入门槛是20元(完全是为了防止白嫖党,请理解),请考虑清楚。群里有免费推文的注释代码和示例数据(终身拥有),没有付费内容,群成员福利是购买单个付费内容半价!
2、《KS科研分享与服务》微信VIP群只针对购买打包代码的小伙伴(公众号所有付费内容合集)!微信群不是单独的,是对于打包的人答疑解惑和交流的平台、群成员专享视频教程,帖子提前发布,以及其他更多福利!
点击:→ 加入微信vip群:2024-2025《KS科研分享与服务》付费内容打包集合
3、需进QQ群或者打包代码入微信VIP的小伙伴请添加作者微信了解,请备注目的,除此之外请勿添加,谢谢!
详情请联系作者:
❞
Macrophage <- subset(human_data, celltype=='Macrophage')
diff <- FindMarkers(Macrophage, ident.1 = "GM", ident.2 = "BM",
group.by = "group", logfc.threshold = 0.25,min.pct = 0.25)
diff$gene <- rownames(diff)
diff$group <- ""
diff$group <- ifelse(diff$avg_log2FC>0,"up",'down')
#富集分析,我们这里就以KEGG为例子
group <- data.frame(gene=diff$gene,group=diff$group)#分组情况
#gene转化为ID
Gene_ID <- bitr(diff$gene, fromType="SYMBOL", toType="ENTREZID", OrgDb="org.Hs.eg.db")
#构建文件并分析
data <- merge(Gene_ID,group,by.x='SYMBOL',by.y='gene')
diff_KEGG <- compareCluster(ENTREZID~group,
data=data,
fun = "enrichKEGG",#函数选择什么定义什么分析
pAdjustMethod = "BH",
pvalueCutoff = 0.01,
qvalueCutoff = 0.01,
organism= "hsa")#物种
#将gene ID转化为gene symbol
diff_KEGG = setReadable(diff_KEGG,OrgDb = "org.Hs.eg.db", keyType = "ENTREZID")
#获取富集分析表格文件
diff_KEGG <- diff_KEGG@compareClusterResult
diff_KEGG <- diff_KEGG %>%
group_by(group) %>%
top_n(n = 5, wt = -qvalue)
#排序
diff_KEGG$group <- factor(diff_KEGG$group, levels = c("up","down"))
# 使用排序索引重新排列数据框
diff_KEGG <- diff_KEGG[order(diff_KEGG$group), ]
#terms因子顺序
diff_KEGG$Description <- factor(diff_KEGG$Description, levels = diff_KEGG$Description)
#展示的基因,我们选择每个terms展示5个基因,实际情况可以展示自己关注的基因
diff_KEGG$geneID <- sapply(strsplit(diff_KEGG$geneID , "/"), function(x) paste(x[1:5], collapse = "/"))
ggplot(diff_KEGG, aes(x = -log10(qvalue), y = rev(Description), fill = group))+
geom_bar(stat = "identity", width = 0.5)+
geom_text(aes(x=0.1,y=rev(Description),label = Description),size=3.5, hjust =0)+
theme_classic()+
theme(axis.text.y = element_blank(),
axis.ticks.y = element_blank(),
axis.title.y = element_text(colour = 'black', size = 12),
axis.line = element_line(colour = 'black', linewidth =0.5),
axis.text.x = element_text(colour = 'black', size = 10),
axis.ticks.x = element_line(colour = 'black'),
axis.title.x = element_text(colour = 'black', size = 12),
legend.position = "none")+
scale_x_continuous(expand = c(0,0))+
scale_fill_manual(values = c("#CB5640","#65B0C6"))+
geom_text(data = diff_KEGG,
aes(x = 0.1, y = rev(Description), label = geneID, color = group),
size = 4,
fontface = 'italic',
hjust = 0,
vjust = 2.3)+
scale_color_manual(values = c("#CB5640","#65B0C6"))+
scale_y_discrete(expand = c(0.1,0))+
labs(title = "Enrichment of genes",
y=c("Down Up"))
联系客服