[R] A function for plotting a boxplot with added dot and bars (formean and SE) - please help improve my code

Wed Aug 12 14:23:58 CEST 2009

Two solutions with errorbars

library(ggplot2)
#option 1
lowCI <- function(x){mean(x) + qt(0.025, length(x) - 1) * sd(x) /
sqrt(length(x))}
highCI <- function(x){mean(x) + qt(0.975, length(x) - 1) * sd(x) /
sqrt(length(x))}
ggplot(mtcars, aes(x = factor(round(wt)), y = mpg, colour = factor(am)))
+ geom_boxplot() + geom_point(stat = "summary", fun.y = mean, position =
position_dodge(width = 0.75)) + geom_errorbar(stat = "summary", fun.ymin
= lowCI, fun.ymax = highCI, position = position_dodge(width = 0.75),
linetype = 2)

#option2
mtcars$fwt <- factor(round(mtcars$wt))
MeanCI <- ddply(mtcars, c("fwt", "am"), function(z){
	with(z, mean(mpg) + qt(c(low = 0.025, mean = 0.5, high = 0.975),
length(mpg) - 1) * sd(mpg) / sqrt(length(mpg)))
})
ggplot(mtcars, aes(x = fwt, y = mpg, colour = factor(am))) +
geom_boxplot() + geom_point(data = MeanCI, aes(y = mean), position =
position_dodge(width = 0.75)) + geom_errorbar(data = MeanCI, aes(ymin =
low, ymax = high), position = position_dodge(width = 0.75), linetype =
2)

------------------------------------------------------------------------
----
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature
and Forest
Cel biometrie, methodologie en kwaliteitszorg / Section biometrics,
methodology and quality assurance
Gaverstraat 4
9500 Geraardsbergen
Belgium
tel. + 32 54/436 185
Thierry.Onkelinx at inbo.be
www.inbo.be

To call in the statistician after the experiment is done may be no more
than asking him to perform a post-mortem examination: he may be able to
say what the experiment died of.
~ Sir Ronald Aylmer Fisher

The plural of anecdote is not data.
~ Roger Brinner

The combination of some data and an aching desire for an answer does not
ensure that a reasonable answer can be extracted from a given body of
data.
~ John Tukey

________________________________

Van: Tal Galili [mailto:tal.galili at gmail.com] 
Verzonden: woensdag 12 augustus 2009 12:48
Aan: ONKELINX, Thierry
CC: r-help at r-project.org
Onderwerp: Re: [R] A function for plotting a boxplot with added dot and
bars (formean and SE) - please help improve my code

Thanks Thierry, 
I didn't know about that, it's always nice to see different ways of
doing things.
Question: How would you add the bars ?

Hope to see more interesting replies,
Tal

On Wed, Aug 12, 2009 at 1:29 PM, ONKELINX, Thierry
<Thierry.ONKELINX at inbo.be> wrote:

	This can be done much easier and transparent with ggplot2

	library(ggplot2)
	ggplot(mtcars, aes(x = factor(round(wt)), y = mpg, colour =
factor(am)))
	+ geom_boxplot() + geom_point(stat = "summary", fun.y = "mean",
position
	= position_dodge(width = 0.75))

	HTH,

	Thierry

------------------------------------------------------------------------
	----
	ir. Thierry Onkelinx
	Instituut voor natuur- en bosonderzoek / Research Institute for
Nature
	and Forest
	Cel biometrie, methodologie en kwaliteitszorg / Section
biometrics,
	methodology and quality assurance
	Gaverstraat 4
	9500 Geraardsbergen
	Belgium
	tel. + 32 54/436 185
	Thierry.Onkelinx at inbo.be
	www.inbo.be

	To call in the statistician after the experiment is done may be
no more
	than asking him to perform a post-mortem examination: he may be
able to
	say what the experiment died of.
	~ Sir Ronald Aylmer Fisher

	The plural of anecdote is not data.
	~ Roger Brinner

	The combination of some data and an aching desire for an answer
does not
	ensure that a reasonable answer can be extracted from a given
body of
	data.
	~ John Tukey

	-----Oorspronkelijk bericht-----
	Van: r-help-bounces at r-project.org
[mailto:r-help-bounces at r-project.org]
	Namens Tal Galili
	Verzonden: woensdag 12 augustus 2009 11:48
	Aan: r-help at r-project.org
	Onderwerp: [R] A function for plotting a boxplot with added dot
and bars
	(formean and SE) - please help improve my code

	Hello people,

	A while back I wanted to plot boxplots with interactions that
will have
	a dot for the mean of the sample + bars for the SE.
	After searching for some code, I found something that did it for
one
	level, but couldn't find something that will allow for
interactions the
	way the original boxplot does.

	After playing with the original code, I found a way for allowing
the
	boxplot code to introduce interactions to it. The price for that
was
	that I was forced to use a little different syntax for the
function then
	that of the original boxplot. In order to use the function, one
must use
	lists for the function input arguments, and until now I didn't
find a
	way for doing the same task from a formula input.
	So for example, the original boxplot will be written like this:
	boxplot(y ~ A*B)
	Where as my function will look like this:
	boxplot.2(y , list(A,B) )

	In this e-mail I am giving away the code for:
	1) helping out others searching for this solution. and
	2) in hope to have more experienced R programmers come by and
improve on
	this code (by, for example, removing the for loop in it, or
allowing to
	use a formula instead of a list input)

	Here is the code:

	boxplot.2 <- function(fo.head, list.fo.tail = list(1),
print.mean = T,
	plot.CI = T, add.mean.sd.to.boxplot.names = T, plot.round.factor
= 2
	,...) {
	 require(gplots)
	 tmp   <- split(fo.head ,  list.fo.tail)
	 means <- sapply(tmp, mean)
	 stdev <- sqrt(sapply(tmp, var))  #IS right - because of the
sqrt !!!
	<- sqrt(sapply(tmp, var))  # was var, I changed it to sd
	 n     <- sapply(tmp,length)
	 ciw   <- qt(0.975, n-1) * stdev / sqrt(n)
	 old.names = attributes(tmp)$names
	 length.of.names = length(old.names)
	 new.names = as.list(old.names)
	 for(i in c(1:length.of.names))
	 {
	   new.names[[i]] <- paste(old.names[i], " (",
	round(means[i],plot.round.factor) ,",",
round(ciw[i],plot.round.factor)
	,")", sep = "")
	 }

	 if(add.mean.sd.to.boxplot.names)
	      {
	          if(length(list.fo.tail) == 1)
	          {
	            sub.text = paste("mean:", round(means,3),
	"(SE:",round(ciw,3), " ; N:" ,n, ")")
	            boxplot(tmp, xlab = sub.text, cex.axis = min(1,
	mean(c(1,12/
	max(nchar(new.names)) ))  ), ...)
	          } else {
	            boxplot(tmp, names = new.names, cex.axis = min(1,
	mean(c(1,12/
	max(nchar(new.names)) ))  ), ...)
	          }
	 } else {  boxplot(tmp, ...) } # else - make boxplot without
them

	 # adding the points and the bars
	 points(means, col = 'red', pch = 19)
	 if(require(gplots) & plot.CI) {
	 plotCI(x=means, uiw=ciw, col="white", barcol="blue",
	        xaxt="n" , add = T)
	        }

	 if(print.mean) {
	     small.mean.and.se.table <-
rbind(round(means,3),round(ciw,3), n)
	     rownames(small.mean.and.se.table) <- c("means:", "SE's:",
"N")
	     print(t(small.mean.and.se.table))
	   }
	} # <-- end boxplot3

	# Here is a small example:
	data(mtcars)
	head(mtcars)
	attach(mtcars)
	boxplot.2(mpg , list(round(wt), am) , data = mtcars, las = 2)

	Hope this will benefit others, and also hope that others will
improve on
	this code and give it back to the community, Tal galili

	--
	----------------------------------------------

	My contact information:
	Tal Galili
	Phone number: 972-50-3373767
	FaceBook: Tal Galili
	My Blogs:
	http://www.r-statistics.com/
	http://www.talgalili.com
	http://www.biostatistics.co.il

	       [[alternative HTML version deleted]]

	______________________________________________
	R-help at r-project.org mailing list
	https://stat.ethz.ch/mailman/listinfo/r-help
	PLEASE do read the posting guide
	http://www.R-project.org/posting-guide.html
	and provide commented, minimal, self-contained, reproducible
code.

	Dit bericht en eventuele bijlagen geven enkel de visie van de
schrijver weer
	en binden het INBO onder geen enkel beding, zolang dit bericht
niet bevestigd is
	door een geldig ondertekend document. The views expressed in
this message
	and any annex are purely those of the writer and may not be
regarded as stating
	an official position of INBO, as long as the message is not
confirmed by a duly
	signed document.

-- 
----------------------------------------------

My contact information:
Tal Galili
Phone number: 972-50-3373767
FaceBook: Tal Galili
My Blogs:
http://www.r-statistics.com/
http://www.talgalili.com
http://www.biostatistics.co.il

Dit bericht en eventuele bijlagen geven enkel de visie van de schrijver weer 
en binden het INBO onder geen enkel beding, zolang dit bericht niet bevestigd is
door een geldig ondertekend document. The views expressed in  this message 
and any annex are purely those of the writer and may not be regarded as stating 
an official position of INBO, as long as the message is not confirmed by a duly 
signed document.