看到一篇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
最后,觉得分享对你有帮助,点个赞,分享一下呗!