Question

I am using the beanplot function of the beanplot package, and I cannot figure out a way to use the wd= parameter with a good result.

What I want

  • display beanplots whose width depends on the size of the sample.
  • (and now) understand what is the use of 'wd=' and how to use it.

So far, when I try to use the wd parameter, as a list it gives me an error, as a vector, it gives me something weird (wd seems to be multiplied to the density estimates values)

Example

library(beanplot)
set.seed(2000)
Test <- data.frame(
  x=rnorm(30), 
  f1 = factor(c(rep('A', 10), rep('B',20))), 
  f2=factor(c('M','F'))
  )

beanplot(x~f1,Test,
         col=list('orange','yellow'),
         wd=c(1:2/2),
         boxwex = 1
)

beanplot1

Was it helpful?

Solution

Okay, after a bit of personal trial/error, I got my two answers:

First, after looking into the code, wd is not supposed to be used that way, and does not support to have several values (unlike 'col='). It seems that 'wd=' is not meant to be used directly for specifying the beanplot width in a normal use. Indeed, when no 'wd' (or wd=NA) is provided, wd is computed from 'maxwidth=', but normalized to the highest value of all the density functions. Hence it seems more appropriate to specify 'maxwidth=' if one wants to control the width of the beans.

Second, I wrote some code to actually achieve what I want. It may be messy, please fell free to enhance it. This example can be used for any variation of parameter that is not supported by the function beanplot natively.

You will notice that the code shows an example of why 'wd=' is still important, since we want all the width values of the bean to be normalized one time, for all the density plots.

col_list = list('orange','yellow')
wd_fac = aggregate(x~f1,Test,length)$x # get the size relative to the number of points in each bean
wd_fac <- wd_fac/max(wd_fac)

par(mfrow=c(1,2))
## Normal plot
beanplot(x~f1,Test, col=col_list)
## Our plot
bp <- beanplot(x~f1,Test,what=c(T,F,F,F)) # plots the line and frames + get the general parmeters
sapply(1:length(levels(Test$f1)),
       function(X){
         beanplot(subset(Test, f1 == levels(f1)[X])$x,
                  col=col_list[X],
                  bw = bp$bw, # we may want to keep the bandwidth
                  wd= bp$wd*wd_fac[X], 
                  at=X,
                  what=!c(T,F,F,F),
                  add = T)}
)

beanplot

Licensed under: CC-BY-SA with attribution
Not affiliated with StackOverflow
scroll top