R学习 - 箱线图 小提琴图

箱线图


箱线图是能同时反映数据统计量和整体分布，又很漂亮的展示图。在2014年的Nature Method上有2篇Correspondence论述了使用箱线图的好处和一个在线绘制箱线图的工具。就这样都可以发两篇Nature method，没天理，但也说明了箱线图的重要意义。
下面这张图展示了Bar plot、Box
plot、Volin plot和Bean plot对数据分布的反应。从Bar plot上只能看到数据标准差或标准误不同；Box plot可以看到数据分布的集中性不同；Violin plot和Bean plot展示的是数据真正的分布，尤其是对Biomodal数据的展示。
Boxplot从下到上展示的是最小值，第一四分位数(箱子的下边线)、中位数(箱子中间的线)、第三四分位数(箱子上边线)、最大值，具体解读参见http://mp.weixin.qq.com/s/t3UTI_qAIi0cy1g6ZmHtwg。


    stat_smooth
    stat_smooth
    stat_smooth
    stat_smooth
    stat_smooth
  • Nature Method文章

    一步步解析箱线图绘制

    假设有这么一个基因表达矩阵，第一列为基因名字，后面几列为样品名字，想绘制下样品中基因表达的整体分布。
    profile="Name;2cell_1;2cell_2;2cell_3;4cell_1;4cell_2;4cell_3;zygote_1;zygote_2;zygote_3
A;4;6;7;3.2;5.2;5.6;2;4;3
B;6;8;9;5.2;7.2;7.6;4;6;5
C;8;10;11;7.2;9.2;9.6;6;8;7
D;10;12;13;9.2;11.2;11.6;8;10;9
E;12;14;15;11.2;13.2;13.6;10;12;11
F;14;16;17;13.2;15.2;15.6;12;14;13
G;15;17;18;14.2;16.2;16.6;13;15;14
H;16;18;19;15.2;17.2;17.6;14;16;15
I;17;19;20;16.2;18.2;18.6;15;17;16
J;18;20;21;17.2;19.2;19.6;16;18;17
L;19;21;22;18.2;20.2;20.6;17;19;18
M;20;22;23;19.2;21.2;21.6;18;20;19
N;21;23;24;20.2;22.2;22.6;19;21;20
O;22;24;25;21.2;23.2;23.6;20;22;21"
    读入数据并转换为ggplot2需要的长数据表格式(经过前面几篇的练习，这应该都很熟了)
    profile_text <- read.table(text=profile,header=T, row.names=1, quote="",sep=";", check.names=F)
#在melt时保留位置信息
# melt格式是ggplot2画图最喜欢的格式
#好好体会下这个格式，虽然多占用了不少空间，但是确实很方便
library(ggplot2)
library(reshape2)
data_m <- melt(profile_text)
head(data_m)
    variable value
    12cell_14
    22cell_16
    32cell_18
    42cell_110
    52cell_112
    62cell_114
    像往常一样，就可以直接画图了。
    # variable和value为矩阵melt后的两列的名字，内部变量, variable代表了点线的属性，value代表对应的值。
p <- ggplot(data_m, aes(x=variable,y=value),color=variable) +
geom_boxplot() +
theme(axis.text.x=element_text(angle=50,hjust=0.5,vjust=0.5)) +
theme(legend.position="none")
p
#图会存储在当前目录的Rplots.pdf文件中，如果用Rstudio，可以不运行dev.off()
dev.off()
    箱线图出来了，看上去还可以，再加点色彩。
    # variable和value为矩阵melt后的两列的名字，内部变量, variable代表了点线的属性，value代表对应的值。
p <- ggplot(data_m, aes(x=variable,y=value),color=variable) +
geom_boxplot(aes(fill=factor(variable))) +
theme(axis.text.x=element_text(angle=50,hjust=0.5,vjust=0.5)) +
theme(legend.position="none")
p
#图会存储在当前目录的Rplots.pdf文件中，如果用Rstudio，可以不运行dev.off()
dev.off()

    再看看Violin plot
    # variable和value为矩阵melt后的两列的名字，内部变量, variable代表了点线的属性，value代表对应的值。
p <- ggplot(data_m, aes(x=variable,y=value),color=variable) +
geom_violin(aes(fill=factor(variable))) +
theme(axis.text.x=element_text(angle=50,hjust=0.5,vjust=0.5)) +
theme(legend.position="none")
p
#图会存储在当前目录的Rplots.pdf文件中，如果用Rstudio，可以不运行dev.off()
dev.off()

    还有Jitter plot (这里使用的是ggbeeswarm包)
    library(ggbeeswarm)
#为了更好的效果，只保留其中一个样品的数据
# grepl类似于Linux的grep命令，获取特定模式的字符串
data_m2 <- data_m[grepl("_3",data_m$variable),]
# variable和value为矩阵melt后的两列的名字，内部变量, variable代表了点线的属性，value代表对应的值。
p <- ggplot(data_m2, aes(x=variable,y=value),color=variable) +
geom_quasirandom(aes(colour=factor(variable)))+
theme_bw() + theme(panel.grid.major =element_blank(),
panel.grid.minor = element_blank(),legend.key=element_blank()) +
theme(legend.position="none")
#也可以用geom_jitter(aes(colour=factor(variable)))代替geom_quasirandom(aes(colour=factor(variable)))
#但个人认为geom_quasirandom给出的结果更有特色
ggsave(p,filename="jitterplot.pdf", width=14, height=8,units=c("cm"))

    绘制单个基因(A)的箱线图

    为了更好的展示效果，下面的矩阵增加了样品数量和样品的分组信息。
    profile="Name;2cell_1;2cell_2;2cell_3;2cell_4;2cell_5;2cell_6;4cell_1;4cell_2;4cell_3;4cell_4;4cell_5;4cell_6;zygote_1;zygote_2;zygote_3;zygote_4;zygote_5;zygote_6
A;4;6;7;5;8;6;3.2;5.2;5.6;3.6;7.6;4.8;2;4;3;2;4;2.5
B;6;8;9;7;10;8;5.2;7.2;7.6;5.6;9.6;6.8;4;6;5;4;6;4.5"
profile_text <- read.table(text=profile,header=T, row.names=1, quote="",sep=";", check.names=F)
data_m = data.frame(t(profile_text['A',]))
data_m$sample = rownames(data_m)
#只挑选显示部分
# grepl前面已经讲过用于匹配
data_m[grepl('_[123]', data_m$sample),]
    Asample
    2cell_14.02cell_1
    2cell_26.02cell_2
    2cell_37.02cell_3
    4cell_13.24cell_1
    4cell_25.24cell_2
    4cell_35.64cell_3
    zygote_12.0 zygote_1
    zygote_24.0 zygote_2
    zygote_33.0 zygote_3
    获得样品分组信息(这个例子比较特殊，样品的分组信息就是样品名字下划线前面的部分)
    #可以利用strsplit分割，取出其前面的字符串
# R中复杂的输出结果多数以列表的形式体现，在之前的矩阵操作教程中
#提到过用str函数来查看复杂结果的结构，并从中获取信息
group =unlist(lapply(strsplit(data_m$sample,"_"), function(x) x[1]))
data_m$group = group
data_m[grepl('_[123]', data_m$sample),]
    Asamplegroup
    2cell_14.02cell_12cell
    2cell_26.02cell_22cell
    2cell_37.02cell_32cell
    4cell_13.24cell_14cell
    4cell_25.24cell_24cell
    4cell_35.64cell_34cell
    zygote_12.0 zygote_1 zygote
    zygote_24.0 zygote_2 zygote
    zygote_33.0 zygote_3 zygote
    如果没有这个规律，也可以提到类似于下面的文件，指定样品所属的组的信息。
    sampleGroup_text="Sample;Group
zygote_1;zygote
zygote_2;zygote
zygote_3;zygote
zygote_4;zygote
zygote_5;zygote
zygote_6;zygote
2cell_1;2cell
2cell_2;2cell
2cell_3;2cell
2cell_4;2cell
2cell_5;2cell
2cell_6;2cell
4cell_1;4cell
4cell_2;4cell
4cell_3;4cell
4cell_4;4cell
4cell_5;4cell
4cell_6;4cell"
#sampleGroup =read.table(text=sampleGroup_text,sep="\t",header=1,check.names=F,row.names=1)
#data_m <- merge(data_m, sampleGroup,by="row.names")
#会获得相同的结果，脚本注释掉了以免重复执行引起问题。
    矩阵准备好了，开始画图了(小提琴图做例子，其它类似)
    #调整下样品出现的顺序
data_m$group <- factor(data_m$group,levels=c("zygote","2cell","4cell"))
# group和A为矩阵中两列的名字，group代表了值的属性，A代表基因A对应的表达值。
#注意看修改了的地方
p <- ggplot(data_m, aes(x=group,y=A),color=group) +
geom_violin(aes(fill=factor(group))) +
theme(axis.text.x=element_text(angle=50,hjust=0.5,vjust=0.5)) +
theme(legend.position="none")
p
#图会存储在当前目录的Rplots.pdf文件中，如果用Rstudio，可以不运行dev.off()
dev.off()

    长矩阵绘制箱线图

    常规矩阵绘制箱线图要求必须是个方正的矩阵输入，而有时想比较的几个组里面检测的值数目不同。比如有三个组，GrpA组检测了6个病人，GrpB组检测了10个病人，GrpC组是12个正常人的检测数据。这时就很难形成一个行位检测值，列为样品的矩阵，长表格模式就适合与这种情况。
    long_table <- "Grp;Value
GrpA;10
GrpA;11
GrpA;12
GrpB;5
GrpB;4
GrpB;3
GrpB;2
GrpC;2
GrpC;3"
long_table <-read.table(text=long_table,sep="\t",header=1,check.names=F)
p <- ggplot(data_m, aes(x=Grp,y=Value),color=Grp) +
geom_violin(aes(fill=factor(Grp))) +
theme(axis.text.x=element_text(angle=50,hjust=0.5,vjust=0.5)) +
theme(legend.position="none")
p
    长表格形式自身就是常规矩阵melt后的格式，这种用来绘制箱线图就很简单了，就不举例子了。

    箱线图-一步绘制

    绘图时通常会碰到两个头疼的问题：
    1.有时需要绘制很多的图，唯一的不同就是输入文件，其它都不需要修改。如果用R脚本，需要反复替换文件名，繁琐又容易出错。(R也有命令行参数，不熟，有经验的可以尝试下)
    2.每次绘图都需要不断的调整参数，时间久了不用，就忘记参数怎么设置了；或者调整次数过多，有了很多版本，最后不知道用哪个了。
    为了简化绘图、维持脚本的一致，我用bash对绘图命令做了一个封装，通过配置修改命令行参数，生成相应的绘图脚本，然后再绘制。
    首先把测试数据存储到文件中方便调用。数据矩阵存储在boxplot.normal.data、sampleGroup和boxplot.melt.data文件中(TAB键分割，内容在文档最后。如果你手上有自己的数据，也可以拿来用)。
    使用正常矩阵默认参数绘制箱线图
    # -f:指定输入的矩阵文件，第一列为行名字，第一行为header
列数不限，列名字不限；行数不限，行名字默认为文本
sp_boxplot.sh -f boxplot.normal.data
    箱线图出来了，但有点小乱。
    # -f:指定输入的矩阵文件，第一列为行名字，第一行为header
列数不限，列名字不限；行数不限，行名字默认为文本
# -P: none,去掉legend(uppercase P)
# -b: X-axis旋转45度
# -V: TRUE绘制小提琴图
sp_boxplot.sh -f boxplot.normal.data -Pnone -b 45 -V TRUE

    绘制单个基因的小提琴图加抖动图
    # -q:指定某一行的名字，此处为基因名，绘制基因A的表达图谱
# -Q:指定样本分组，绘制基因A在不同样品组的表达趋势
# -F Group: sampleGroup中第二列的名字，指代分组信息，根据需要修改
# -J TRUE:绘制抖动图jitter plot
# -L:设置X轴样品组顺序
# -c TRUE -C "'red', 'pink',
'blue'":指定每个箱线图的颜色
sp_boxplot.sh -f boxplot.normal.data -q A-Q sampleGroup -F Group -V TRUE -J TRUE -L "'zygote','2cell','4cell'"-c TRUE -C "'red', 'pink', 'blue'" -P none

    使用melted矩阵默认参数绘箱线图
    # -f:指定输入文件
# -m TRUE:指定输入的矩阵为meltedformat
# -d Expr：指定表达值所在的列
# -F Rep:指定子类所在列，也就是legend
# -a Group：指定X轴分组信息
# -j TRUE: jitter plot
sp_boxplot.sh -f boxplot.melt.data -m TRUE-d Expr -F Rep -a Group-j TRUE

    #如果没有子类，则-a和-F指定为同一值
# -R TRUE:旋转boxplot
sp_boxplot.sh -f boxplot.melt.data -m TRUE-d Expr -a Group -F Group -J TRUE -R TRUE

    参数中最需要注意的是引号的使用：
    外层引号与内层引号不能相同 凡参数值中包括了空格，括号，逗号等都用引号括起来作为一个整体。
    为了推广，也为了激起大家的热情，如果想要sp_boxplot.sh脚本的，还需要劳烦大家动动手，转发此文章到朋友圈，并留言索取。
    也希望大家能一起开发，完善功能。
    #boxplot.normal.data
Name2cell_12cell_22cell_32cell_42cell_52cell_64cell_14cell_24cell_34cell_44cell_54cell_6zygote_1zygote_2zygote_3zygote_4zygote_5zygote_6
A4675863.25.25.63.67.64.8243242.5
B68971085.27.27.65.69.66.8465464.5
C81011912107.29.29.67.611.68.8687686.5
D1012131114129.211.211.69.613.610.881098108.5
E12141513161411.213.213.611.615.612.8101211101210.5
F14161715181613.215.215.613.617.614.8121413121412.5
G15171816191714.216.216.614.618.615.8131514131513.5
H16181917201815.217.217.615.619.616.8141615141614.5
I17192018211916.218.218.616.620.617.8151716151715.5
J18202119222017.219.219.617.621.618.8161817161816.5
L19212220232118.220.220.618.622.619.8171918171917.5
M20222321242219.221.221.619.623.620.8182019182018.5
N21232422252320.222.222.620.624.621.8192120192119.5
O22242523262421.223.223.621.625.622.8202221202220.5

    #boxplot.melt.data
GeneSampleGroupExprRep
Azygote_1zygote21
Azygote_2zygote42
Azygote_3zygote33
Azygote_4zygote24
Azygote_5zygote45
Azygote_6zygote2.56
A2cell_12cell41
A2cell_22cell62
A2cell_32cell73
A2cell_42cell54
A2cell_52cell85
A2cell_62cell66
A4cell_14cell3.21
A4cell_24cell5.22
A4cell_34cell5.63
A4cell_44cell3.64
A4cell_54cell7.65
A4cell_64cell4.86

    #sampleGroup
SampleGroup
zygote_1zygote
zygote_2zygote
zygote_3zygote
zygote_4zygote
zygote_5zygote
zygote_6zygote
2cell_12cell
2cell_22cell
2cell_32cell
2cell_42cell
2cell_52cell
2cell_62cell
4cell_14cell
4cell_24cell
4cell_34cell
4cell_44cell
4cell_54cell
4cell_64cell

    Reference

    *http://blog.genesino.com//2017/06/R-Rstudio

    关注生信宝典，几千人一起学生信


​

最后编辑于 ：
©著作权归作者所有,转载或内容合作请联系作者
  • 人面猴
    序言：七十年代末，一起剥皮案震惊了整个滨河市，随后出现的几起案子，更是在滨河造成了极大的恐慌，老刑警刘岩，带你破解...
    沈念sama阅读 212,222评论 6赞 493
  • 死咒
    序言：滨河连续发生了三起死亡事件，死亡现场离奇诡异，居然都是意外死亡，警方通过查阅死者的电脑和手机，发现死者居然都...
    沈念sama阅读 90,455评论 3赞 385
  • 救了他两次的神仙让他今天三更去死
    文/潘晓璐 我一进店门，熙熙楼的掌柜王于贵愁眉苦脸地迎上来，“玉大人，你说我怎么就摊上这事。” “怎么了？”我有些...
    开封第一讲书人阅读 157,720评论 0赞 348
  • 道士缉凶录：失踪的卖姜人
    文/不坏的土叔 我叫张陵，是天一观的道长。 经常有香客问我，道长，这世上最难降的妖魔是什么？ 我笑而不...
    开封第一讲书人阅读 56,568评论 1赞 284
  • 港岛之恋（遗憾婚礼）
    正文 为了忘掉前任，我火速办了婚礼，结果婚礼上，老公的妹妹穿的比我还像新娘。我一直安慰自己，他们只是感情好，可当我...
    茶点故事阅读 65,696评论 6赞 386
  • 恶毒庶女顶嫁案：这布局不是一般人想出来的
    文/花漫 我一把揭开白布。 她就那样静静地躺着，像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上，一...
    开封第一讲书人阅读 49,879评论 1赞 290
  • 城市分裂传说
    那天，我揣着相机与录音，去河边找鬼。 笑死，一个胖子当着我的面吹牛，可吹牛的内容都是我干的。 我是一名探鬼主播，决...
    沈念sama阅读 39,028评论 3赞 409
  • 双鸳鸯连环套：你想象不到人心有多黑
    文/苍兰香墨 我猛地睁开眼，长吁一口气：“原来是场噩梦啊……” “哼！你这毒妇竟也来了？” 一声冷哼从身侧响起，我...
    开封第一讲书人阅读 37,773评论 0赞 268
  • 万荣杀人案实录
    序言：老挝万荣一对情侣失踪，失踪者是张志新（化名）和其女友刘颖，没想到半个月后，有当地人在树林里发现了一具尸体，经...
    沈念sama阅读 44,220评论 1赞 303
  • 护林员之死
    正文 独居荒郊野岭守林人离奇死亡，尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 36,550评论 2赞 327
  • 白月光启示录
    正文 我和宋清朗相恋三年，在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 38,697评论 1赞 341
  • 活死人
    序言：一个原本活蹦乱跳的男人离奇死亡，死状恐怖，灵堂内的尸体忽然破棺而出，到底是诈尸还是另有隐情，我是刑警宁泽，带...
    沈念sama阅读 34,360评论 4赞 332
  • 日本核电站爆炸内幕
    正文 年R本政府宣布，位于F岛的核电站，受9级特大地震影响，放射性物质发生泄漏。R本人自食恶果不足惜，却给世界环境...
    茶点故事阅读 40,002评论 3赞 315
  • 男人毒药：我在死后第九天来索命
    文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹，春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 30,782评论 0赞 21
  • 一桩弑父案，背后竟有这般阴谋
    文/苍兰香墨 我抬头看了看天上的太阳。三九已至，却和暖如春，着一层夹袄步出监牢的瞬间，已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 32,010评论 1赞 266
  • 情欲美人皮
    我被黑心中介骗来泰国打工， 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留，地道东北人。 一个月前我还...
    沈念sama阅读 46,433评论 2赞 360
  • 代替公主和亲
    正文 我出身青楼，却偏偏与公主长得像，于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子，可洞房花烛夜当晚...
    茶点故事阅读 43,587评论 2赞 350

推荐阅读更多精彩内容

  • R学习-线图
    线图 线图是反映趋势变化的一种方式，其输入数据一般也是一个矩阵。 单线图 假设有这么一个矩阵，第一列为转录起始位点...
    生信宝典阅读 917评论 0赞 0
  • 男人戒烟后……
    很多男性都知道吸烟有害健康，嚷嚷着戒烟的人很多，实际上成功的人却很少，很多人都处于“吸了戒、戒了吸”的状态，这样一...
    商佳阅读 649评论 0赞 0
  • 【iOS】Objective-C 对象模型
    元类(Meta Class) 元类存储着一个类的所有类方法，当我们向一个对象发送消息时，runtime会在这个对象...
    cod_mm阅读 734评论 0赞 1
  • 子玉（2）
    二 子玉的生活很平静，是从什么时候开始改变的呢？后来的子玉仔细回想，感觉像是生活突然就来了一道雷，把自己的生活劈的...
    雨羽yuyu阅读 401评论 0赞 0
  • Android一个APK多个入口（多个桌面图标）的实现
    前言   Android应用一般都是一个APK一个桌面图标，但有时候我们需要实现一个APK在桌面上有多个图标（比如...
    张明云阅读 8,752评论 1赞 36