R语言绘制密度提琴图_vlambda技术博客

vlambda
2020-06-09

R语言绘制密度提琴图

就叫它“密度提琴图”吧

某日，老大发问，有这样的一张图，有什么思路解决它吗。

（来源文献“A genome-based model for adjusting radiotherapy dose (GARD): a retrospective, cohort-based study”）

好了，第一眼看去，提琴图嘛，对不对？

其它的都还好说，貌似唯一难点就是怎样根据“提琴宽度”分配渐变色了。所谓提琴宽度，也就是对应数值的频数分布（密度）。

然后就找方法，怎么在提琴图中，加入渐变色。当然最好是R，最好是ggplot2。

然后试了半天，关键词用了个遍，比方说“gradient violin”、“density violin”、“histogram violin”、“heatmap violin”之类的，毛线都没看到......

唯一找到一个比较像的示例，本已喜出望外，然后发现它只是按分位数为提琴图上色，我的心情嘞......

（ https://stackoverflow.com/questions/22278951/combining-violin-plot-with-box-plot ）

然后一上午过去了……

中午吃完饭回来，突然想到，密度？等等，我记得好像有个ggplot2的拓展包，可以画“半边的密度图”（看了下文就知道我为什么一开始这样称呼它了）。好吧，人家学名叫“核密度图”。

然后找到那个包那个函数。

  
    
    
  
   
     
     
   library(ggplot2)
   
     
     
   
library(ggridges)
   
     
     
   
 
   
     
     
   
#详情 ?geom_density_ridges_gradient，这是文档中下方示例
   
     
     
   
ggplot(lincoln_weather, aes(x = `Mean Temperature [F]`, y = `Month`, fill = ..x..)) +
   
     
     
   
geom_density_ridges_gradient(scale = 3, rel_min_height = 0.01) +
   
     
     
   
scale_x_continuous(expand = c(0.01, 0)) +
   
     
     
   
scale_y_discrete(expand = c(0.01, 0)) +
   
     
     
   
scale_fill_gradientn(colours = terrain.colors(10)) +
   
     
     
   
labs(title = 'Temperatures in Lincoln NE in 2016') +
   
     
     
   
theme_ridges(font_size = 13, grid = TRUE) + theme(axis.title.y = element_blank())

有没有有联想到什么呢？对吧，这个核密度图，好像就是那个“密度提琴图”的一半形状，对不对？

那么，如果我先画张正方向的核密度图，再来张反方向的，一组合，是不是就刚好拼成一张“密度提琴图”了？

趁着感觉来了，开整。

  
    
    
  
   
     
     
   #还是文档中的示例数据，简化一下方便展示
   
     
     
   
#“正方向”
   
     
     
   
ggplot(lincoln_weather, aes(x = `Mean Temperature [F]`, y = `Month`, fill = ..x..)) +
   
     
     
   
geom_density_ridges_gradient(scale = 0.5, rel_min_height = 0.01) +
   
     
     
   
scale_fill_gradientn(colours = terrain.colors(10))
   
     
     
   
 
   
     
     
   
#“反方向”
   
     
     
   
ggplot(lincoln_weather, aes(x = `Mean Temperature [F]`, y = `Month`, fill = ..x..)) +
   
     
     
   
geom_density_ridges_gradient(scale = -0.5, rel_min_height = 0.01) +
   
     
     
   
scale_fill_gradientn(colours = terrain.colors(10))

没错，就是这样。然后改改颜色、背景、图例啥的，组合在一起，就可以了。

  
    
    
  
   
     
     
   #“密度提琴图”
   
     
     
   
library(ggplot2)
   
     
     
   
library(ggridges)
   
     
     
   
library(grid)
   
     
     
   
 
   
     
     
   
#计算中位数值
   
     
     
   
med <- aggregate(lincoln_weather$`Mean Temperature [F]`, by = list(lincoln_weather$`Month`), FUN = median)
   
     
     
   
 
   
     
     
   
#“正方向”
   
     
     
   
p1 <- ggplot(lincoln_weather, aes(x = `Mean Temperature [F]`, y = `Month`, fill = ..density..)) +
   
     
     
   
geom_density_ridges_gradient(scale = 0.3, gradient_lwd = NA, rel_min_height = NA,size = NA) +
   
     
     
   
scale_fill_gradientn(colours = colorRampPalette(c('#663490', '#5054A3', '#37A3D2', '#3EB897', '#6EBB46', '#C0D857'))(30)) +
   
     
     
   
theme(panel.grid = element_blank(), panel.background = element_rect(fill = 'transparent', color = 'black'), plot.background = element_blank()) +
   
     
     
   
labs(x = 'Mean Temperature', y = 'Month', fill = '')
   
     
     
   
 
   
     
     
   
#“反方向”
   
     
     
   
p2 <- ggplot(lincoln_weather, aes(x = `Mean Temperature [F]`, y = `Month`, fill = ..density..)) +
   
     
     
   
geom_density_ridges_gradient(scale = -0.3, gradient_lwd = NA, rel_min_height = NA,size = NA) +
   
     
     
   
scale_fill_gradientn(colours = colorRampPalette(c('#663490', '#5054A3', '#37A3D2', '#3EB897', '#6EBB46', '#C0D857'))(30)) +
   
     
     
   
theme(panel.grid = element_blank(), panel.background = element_blank(), plot.background = element_blank()) +
   
     
     
   
labs(x = 'Mean Temperature', y = 'Month', fill = '') +
   
     
     
   
#把中位数以红点的形式添在图中
   
     
     
   
annotate('text', label = '.', x = med$x, y = med$Group.1, size = 10, colour = 'red', vjust = -0.02)
   
     
     
   
 
   
     
     
   
grid.newpage()
   
     
     
   
print(p1, vp = viewport(x = 0.5, y = 0.5))
   
     
     
   
print(p2, vp = viewport(x = 0.5, y = 0.5))