复现一篇Nature medicine的火山图，如下：

（reference：Sex-dependent APOE4 neutrophil–microglia interactions drive cognitive impairment in Alzheimer’s disease）

话不多说，看代码（直播视频：哔哩哔哩https://space.bilibili.com/471040659?spm_id_from=333.1007.0.0）：

setwd('D:\\KS项目\\公众号文章\\复现火山图')

#Loading packages
library(Seurat)
library(ggplot2)
library(ggrepel)
library(cowplot)
library(grid)
#load data & prepare data for plot
uterus <- readRDS("D:/KS项目/公众号文章/RNA速率分析/uterus.rds")
SMC <- subset(uterus, celltype=='Smooth muscle cells')
EH <- FindMarkers(uterus, ident.1 = "EEC", ident.2 = "HC",logfc.threshold = 0, min.pct = 0.2, group.by = 'orig.ident')#差异基因分析用来构建演示作图数据
AH <- FindMarkers(uterus, ident.1 = "AEH", ident.2 = "HC",logfc.threshold = 0, min.pct = 0.2, group.by = 'orig.ident')#差异基因分析用来构建演示作图数据


EH$gene <- rownames(EH)#添加gene列
AH$gene <- rownames(AH)#添加gene列

EH$log_pval <- -log10(EH$p_val_adj)#p值取对数
AH$log_pval <- -log10(AH$p_val_adj)#p值取对数



AH$log_pval <- -AH$log_pval #这一组的log_pval变成负值，主要是为了图的倒转结合

#sig gene for label，or you can chose any genes if you want
sig_EH <- EH[EH$p_val_adj <= 0.01 & abs(EH$avg_log2FC)>=1.2,]
sig_AH <- AH[AH$p_val_adj <= 0.01 & abs(AH$avg_log2FC)>=1.2,]


#merge data
df <- rbind(EH,AH)


#设置自己需要的渐变背景
# rect1 <- rasterGrob(alpha(c('#FFF7F3', '#FDE0DD', '#FCC5C0', '#FA9FB5', '#F768A1','#DD3497', '#7A0177'),0.5),
#                     width=unit(1,"npc"), height = unit(1,"npc"), interpolate = TRUE)
# rect2<- rasterGrob(alpha(c('#fff7fb','#ece7f2','#d0d1e6','#a6bddb','#74a9cf','#3690c0','#0570b0','#045a8d'),0.5), 
#                    width=unit(1,"npc"), height = unit(1,"npc"), interpolate = TRUE) 
# rect3<- rasterGrob(c("#FFFFD9", "#EDF8B1", "#C7E9B4", "#7FCDBB", "#41B6C4"), width=unit(1,"npc"), height = unit(1,"npc"), interpolate = TRUE) 
# 
# rect4<- rasterGrob(c("#DC6F58","#F7B698","#FAE7DC","#E1EDF3","#A7CFE4"),width=unit(1,"npc"), height = unit(1,"npc"), interpolate = TRUE) 


rect1 <- rasterGrob(matrix(rev(alpha(c('#FFF7F3', '#FDE0DD', '#FCC5C0', '#FA9FB5', '#F768A1','#DD3497', '#7A0177'),0.5)),nrow = 1),
                    width=unit(1,"npc"), height = unit(1,"npc"), interpolate = TRUE)
rect2<- rasterGrob(matrix(alpha(c('#fff7fb','#ece7f2','#d0d1e6','#a6bddb','#74a9cf','#3690c0','#0570b0','#045a8d'),0.5),nrow = 1), 
                   width=unit(1,"npc"), height = unit(1,"npc"), interpolate = TRUE) 
rect3<- rasterGrob(matrix(rev(c("#FFFFD9", "#EDF8B1", "#C7E9B4", "#7FCDBB", "#41B6C4")),nrow = 1), width=unit(1,"npc"), height = unit(1,"npc"), interpolate = TRUE) 

rect4<- rasterGrob(matrix(rev(c("#DC6F58","#F7B698","#FAE7DC","#E1EDF3")),nrow = 1),width=unit(1,"npc"), height = unit(1,"npc"), interpolate = TRUE) 



#plot
p = ggplot(df, aes(avg_log2FC,log_pval))+
  annotation_custom(rect1,xmin=-Inf,xmax=-0.3,ymin=20,ymax=Inf)+
  annotation_custom(rect2, xmin=0.3,xmax=Inf,ymin=20,ymax=Inf)+
  annotation_custom(rect3,xmin=-Inf,xmax=-0.3,ymin=-20,ymax=-Inf)+
  annotation_custom(rect4, xmin=0.3,xmax=Inf,ymin=-20,ymax=-Inf)+
  geom_hline(aes(yintercept=20), color = "#999999", linetype="dashed", size=1) +#添加横线
  geom_hline(aes(yintercept=-20), color = "#999999", linetype="dashed", size=1) +#添加横线
  geom_vline(aes(xintercept=-0.3), color = "#999999", linetype="dashed", size=1) + #添加纵线
  geom_vline(aes(xintercept=0.3), color = "#999999", linetype="dashed", size=1) + #添加纵线
  geom_point_rast(stroke = 0.5, size=2, shape=21, fill="grey")+  
  geom_point_rast(data = df[abs(df$log_pval) >= 20  &  abs(df$avg_log2FC) >=0.3,], 
             aes(fill=avg_log2FC),stroke = 0.5, size=3, shape=21)+##显著基因特别标注，size比非显著的大一点，颜色用logFC表示
  scale_fill_gradient2(limits=c(-max(df$avg_log2FC), max(df$avg_log2FC)), low="blue", mid="whitesmoke", high = "red", na.value = 'blue')+ #密度颜色设置，并限定范围
  geom_text_repel(data = sig_EH[sig_EH$avg_log2FC >0,], aes(label=gene),
                  nudge_y = 80,
                  nudge_x = 2,
                  color = 'black', size = 4,fontface = 'italic',
                  min.segment.length = 0,
                  max.overlaps = Inf,
                  box.padding = unit(0.5, 'mm'),
                  point.padding = unit(0.5, 'mm')) + #label需要标注的基因
  geom_text_repel(data = sig_EH[sig_EH$avg_log2FC <0,], aes(label=gene),
                  nudge_y = 80,
                  nudge_x = -2,
                  color = 'black', size = 4,fontface = 'italic',
                  min.segment.length = 0,
                  max.overlaps = Inf,
                  box.padding = unit(0.5, 'mm'),
                  point.padding = unit(0.5, 'mm')) + #label需要标注的基因
  geom_text_repel(data = sig_AH[sig_AH$avg_log2FC >0,],aes(label=gene),
                  color = 'black',size = 4, fontface = 'italic',
                  nudge_y = -80,
                  nudge_x = 2,
                  min.segment.length =0,
                  max.overlaps = Inf,
                  point.padding = unit(0.5, 'mm'))+#label需要标注的基因
  geom_text_repel(data = sig_AH[sig_AH$avg_log2FC <0,],aes(label=gene),
                  color = 'black',size = 4, fontface = 'italic',
                  nudge_y = -80,
                  nudge_x = -2,
                  min.segment.length =0,
                  max.overlaps = Inf,
                  point.padding = unit(0.5, 'mm'))+#label需要标注的基因
  theme_classic()+ #设置主题
  theme(aspect.ratio = 1,
        axis.text.x=element_text(colour="black",size =12),
        axis.text.y=element_text(colour="black",size =12),
        axis.ticks=element_line(colour="black"),
        axis.title = element_blank(),
        plot.margin=unit(c(0,0,1,0),"line"),
        legend.direction = 'horizontal',
        legend.position = 'top',
        legend.justification=c(1,2),
        legend.key.width=unit(0.5,"cm"),
        legend.key.height = unit(0.3, "cm"),
        legend.title.position = 'top')+ #修改主题标签文字等等
  geom_rect(aes(xmin =-Inf, xmax = Inf, ymin = 20, ymax = 400),
            fill = "transparent", color = "red", size = 1.5,linetype="dashed")+
  geom_rect(aes(xmin =-Inf, xmax = Inf, ymin = -20, ymax = -400),
            fill = "transparent", color = "blue", size = 1.5,linetype="dashed")+
  scale_y_continuous(expand = c(0,0),breaks = c(-200,0,200), labels = c(200,0,200))#修改y轴标签


p


#add other label
ggdraw(xlim = c(0, 1), ylim = c(0,1.1))+  # 设置绘图区域的界限
  draw_plot(p,x = 0, y =0) +  # 添加主图（热图）
  draw_line(x = c(0.3,0.8), y = c(0.01,0.01),lineend = "round",
            size =1, col = "black",  # 添加箭头
            arrow=arrow(angle = 15, length = unit(0.1,"inches"),type = "closed"))+
  draw_text(text = "Log2FC", size = 12,
            x = 0.2, y = 0.02,color="black",fontface = "italic")+
  draw_line(x = c(0.1,0.1), y = c(0.2,0.85),
            lineend = "round",size =1, col = "black",  # 添加箭头
            arrow=arrow(angle = 15, length = unit(0.1,"inches"),type = "closed"))+
  draw_text(text = "-Log(P)", size = 12,angle=90,
            x = 0.1, y = 0.15,color="black",fontface = "italic")+
  draw_text(text = "EEC vs HC", size = 12,angle=90,x = 0.9, y = 0.7,color="black",fontface = "bold")+
  draw_text(text = "AEC vs HC", size = 12,angle=90,x = 0.9, y = 0.3,color="black",fontface = "bold")