查看原文
其他

同一数据多变量分组的boxplot?

2017-10-31 Y叔 biobabble

小密圈的问题,不是三两句话可以说明白的事情,必须要写文来解答,上一次写文是《听说你还不会画热图》,里面正好吐槽了某知乎大V的「除了ggplot2之外其它都是鸡肋」,这次正好也可以再次呼应一下。

这个图明显是R的base graphics做的,图是可圈可点的,能做出这图来,也已经是告别了只会用plot的低级趣味。这其实是4个图拼起来的,第一个图只是多了个y轴而已,它们画起来是一样的,那就是只画boxplot,不画x和y轴,(你可能会说不是有x轴?),x轴是后面再加上去的,而且加x轴的时候,不写labels,只有线条没有文本(你可能又会说明明有文本!),因为axis这个函数只支持labels要么是水平的,要么是垂直的,旋转某个角度是不支持的,所以labels是额外再打上去的。这里一张小小的图,门道还是挺多的。

set.seed(2017-10-30) d <- data.frame(riskScore = abs(rnorm(100)),                BMI = sample(1:2, 100, replace=T),                stage = sample(1:2, 100, replace=T),                age = sample(1:2, 100, replace=T),                gender = sample(1:2, 100, replace=T)) head(d) ##     riskScore BMI stage age gender ## 1 0.008282743   1     2   1      1 ## 2 0.499375414   1     1   2      1 ## 3 0.188257548   1     1   1      1 ## 4 0.330772189   2     1   2      1 ## 5 0.790797457   1     2   1      1 ## 6 1.943465449   1     1   2      1

先搞一个数据集,都只是随机数,纯粹是为了演示用而已。下面我将定义一个myboxplot,它画boxplot,不带x和y轴,然后加x轴不带labels,再额外打labels,因为提问者的图中还有pvalue,我顺道把pvalue也整合进这个myboxplot里去,可以用pvalue=NULL来关掉这个功能。

关于箱式图,可以参考我之前写的《boxplot》,而这里用到T检验,可以参考《什么是T检验》。

myboxplot <- function(x, data, col = NULL, xlab, pvalue="auto") {    boxplot(x, data, axes = FALSE, col = col)    axis(1, at = 1:2, labels =FALSE)    text(1:2, y=par()$usr[3]-0.08*(par()$usr[4]-par()$usr[3]),         srt=60, xpd=T, adj=1, labels = xlab)    if (pvalue == "auto") {        pvalue <- round(t.test(x, data=data)$p.value, 3)    }    if (!is.null(pvalue)) {        plab <- paste("p =", pvalue)        text(1.5, y = par()$usr[4]*1.05, xpd=T, label=plab, col=col)    } }

万事具备,有函数,有数据,我们先初始化画4个column,然后你只要调用myboxplot,分4次画4个图,就大功告成了,第一个图的时候,把y轴给加上。

layout(t(1:4)) par(oma=c(2, 4, 4, 0), mar=c(5,2,1,1), cex=1) myboxplot(riskScore~age, data=d, col='red', xlab=c("age < 60", "age > 60")) axis(2, las=1) myboxplot(riskScore~gender, data=d, col='green', xlab=c("Male", "Female")) myboxplot(riskScore~stage, data=d, col='blue', xlab=c("pStage 1-2", "pStage 1-2")) myboxplot(riskScore~BMI, data=d, col='cyan', xlab=c("BMI < 24", "BMI > 24"))

假如我们想要用ggplot2来画,该怎么搞?首先毫无意外,要把数据整理成ggplot2喜欢的样子,我定义一个convert函数专门来搞这个数据:

convert <- function(d) {    lapply(2:ncol(d), function(i) {        d2 <- d[, c(1,i)]        d2$type = colnames(d2)[2]        colnames(d2) = c("riskScore", "category", "type")        return(d2)    }) %>% do.call('rbind', .) } dd <- convert(d) head(dd) ##     riskScore category type ## 1 0.008282743        1  BMI ## 2 0.499375414        1  BMI ## 3 0.188257548        1  BMI ## 4 0.330772189        2  BMI ## 5 0.790797457        1  BMI ## 6 1.943465449        1  BMI

然后就可以直接ggplot来画了,这里蛋疼的是颜色不是我们想要的那种用type来上色,如果你指定用type,那么不好意思,不会用category分组,如果你指定group = factor(category)呢?又不好意思了,不会让type来分组了,也就是说你画出来的是按category分的两个box,而不会有不同type是不同组数据的切分了,这就是ggplot2蛋疼之外,语法太高级,以至于有些情况没办法以它的语法表达的时候,是非常困难的。

library(ggplot2) ggplot(dd, aes(type, riskScore, fill=factor(category))) + geom_boxplot()

当然可以通过分面来补救:

ggplot(dd, aes(type, riskScore, group=factor(category), fill=type)) +    geom_boxplot() + facet_grid(.~type, scales = "free_x")

分面也不好啊,分面的strip text和x axis text重了,当然这个x axis text不是我们想要的。

如果我不想用分面呢?你就得用另外的变量去欺骗它,让它分好组,比如这里我用color=factor(category),这样会给boxplot的外框加颜色,但这个颜色不是我们想要用的,category分1和2,在不同的type里意义是不一样的。第二点,它的legend也不是我们想要的,所以这里又需要额外的设置了,我们要指定颜色统一,要去掉legend。

ggplot(dd, aes(type, riskScore, color=factor(category), fill=type)) + geom_boxplot() +    scale_color_manual(values=rep('black',2), guide=FALSE)

这个图就像模像样了,和上面分面的其实差不多,有一个共同点,x axis text不是我们想要的,怎么改?还是illustrator吧,带着ggplot2的枷锁改起来可费劲了。

当然画图嘛,功夫一半在画图上,另一半在于对数据的操作,既要用ggplot2,又要用得爽,你得有70岁的觉悟,「七十而从心所欲,不逾矩」,在矩(俗称枷锁)之下,从心所欲,关键还是对数据和作图系统的理解。这数据还得再变一下,我们不分组了,不就画8个boxplot么,把不同category的1和2全部换掉,换成不同的变量,然后画8个box就完事了。

bmi = c("BMI < 24", "BMI > 24") stage = c("pStage 1-2", "pStage 3-4")     age = c("age < 60", "age > 60") gender = c("Male", "Female") d$BMI = bmi[d$BMI] d$stage = stage[d$stage] d$age = age[d$age] d$gender = gender[d$gender] dd = convert(d) dd$category = factor(dd$category, levels=c(age, gender, stage, bmi))     p1 = ggplot(dd, aes(category, riskScore, fill=type)) + geom_boxplot() +        theme(axis.text.x = element_text(angle=60, vjust=1, hjust=1)) p2 = p1 + facet_grid(.~type, scales="free_x") cowplot::plot_grid(p1, p2, ncol=2)    

数据使然,我们很容易就想到这有4组,每一组有两类,我们要分组,枷锁就来了,后面要改细节,标x axis text,简直是恶梦,除非你放弃治疗用illustrator。如果我们能够放开分组的概念,出王八拳,倒是豁然开朗。这里额外强调一点的是base graphics的作图,容易理解,符合直觉这一块,还是很厉害的,hadley wikham也说ggplot2是试图结合base + lattice的优点。多学一点「外语」对理解和应用「母语」是有帮助的。第二点,放开套路,多试试王八拳,能把对方打趴下的拳,就是好拳!

您可能也对以下帖子感兴趣

文章有问题?点此查看未经处理的缓存