复现NC图表---配对散点小提琴的绘制

看到一篇NC文章绘制的小提请图，文章提供了数据和代码，文章的绘制方式是用的函数，一般直接使用也会出错，这里我们解析下绘图的内部机制，一步一步可能有助于理解为什么ggplot等绘图函数的一些参数。****注释代码及数据已上传群文件！

图片

整体而言，小提琴的绘制是依赖vioplot这个包的。

BiocManager::install("vioplot")
library(vioplot)
library(dplyr)
setwd("D:/KS科研分享与服务_公众号文章/vioplot配对连线小提琴图")
df <- read.csv("df.csv", header = T)

之后我们为排队数据的分组添加上颜色，并设置一列坐标，主要是为了抖动点的绘制，这样有一个动态的感觉，否则只有一个坐标，点就是竖着排列的。

cols <- data.frame(Disease.state=c('Patient','Relative'),
                   color = c('#E69F00', "#55B4EA"))
df_merge <- full_join(df, cols, by = "Disease.state")

# colnames(df)# [1] "sample.id"     "DonorID"       "Family.ID"     "Disease.state" "Richness"  

df_merge$dummy <- 1
df_merge$dummy[df_merge$Disease.state != "Patient" ] <- 2
df_merge$dummy  <- df_merge$dummy + sample(-100 : 100, nrow(df_merge), replace = T)/2000

然后接可以绘制基本的小提琴图了，计算两组平均值并添加。

var.pat <- df$Richness[df$Disease.state == "Patient"]
var.rel <- df$Richness[df$Disease.state == "Relative"]


vioplot(var.pat,  var.rel, drawRect = F,
        names = c("Patients", "Relatives"),
        col=unique(df_merge$color),
        ylim = c(0,240))

segments(.7, mean(var.pat),1.3, mean(var.pat), lwd = 2)
segments(1.7, mean(var.rel), 2.3, mean(var.rel), lwd = 2)

添加散点和配对连线。

添加散点和配对连线。

map.pat <- df_merge[df_merge$Disease.state == "Patient", ]
map.rel <- df_merge[df_merge$Disease.state == "Relative", ]
families <- sort(unique(df_merge$Family.ID))

for(i in 1: length(families)){
  #分别获取两个分组数据
  map.pat.i <- map.pat[map.pat$Family.ID ==  families[i], ]
  map.rel.i <- map.rel[map.rel$Family.ID ==  families[i], ]
  #获取分组每个ID在不同组中的值var.pat.i <- var.pat[map.pat$Family.ID ==  families[i] ]
  var.rel.i <- var.rel[map.rel$Family.ID ==  families[i] ]

  segments(y0 = var.pat.i,  x0 =  map.pat.i$dummy,
                 y1 =  var.rel.i, x1 =  map.rel.i$dummy,
                 col = 8)

  points(map.pat.i$dummy, var.pat.i, pch = 21, bg = 8)
  points(map.rel.i$dummy, var.rel.i, pch = 21, bg = 8)

}

最后添加上标题等就完成了。

title(ylab ="Richness", cex=1.5)
#显著性检验
res <- t.test(df_merge$Richness ~ df_merge$Disease.state, paired = T)
res
segments(1, 210,2,210, lwd = 2)
text(x=1.5, y=220, labels="Pvalue = 8.324e-09")

**其实用到的函数大多是R的基本函数，我想是可以对一些绘图有理解的。此外，这里作图很不方便，但是步骤是全的。感兴趣的小伙伴可以将其修改，封装为一个通用的函数也是不错的选择。参考我们之前的帖子：https://mp.weixin.qq.com/s?__biz=Mzg5OTYzMzY5Ng==&mid=2247486346&idx=1&sn=5aa734906d6045046346564d73283627&chksm=c0510cc5f72685d3f05f8fef2fcdbd5905c8fa49287263aff909316926a1fd21eb280d21f145&scene=21#wechat_redirect

最后，觉得分享对你有帮助，点个赞，分享一下呗！