rebuild cdc obesity maps

you find a short script here

or you can download the images

Visualizing fMRI

Analyze

• first read data (download from R for medicine and biology companion website: http://www.jblearning.com/catalog/9780763758080/)
• extract one point in time
• here you see two snapshots from the resulting 3d object (opengl has to be installed)
• plot using rgl

library("oro.dicom")
library(rgl)
library("misc3d")
library("AnalyzeFMRI")
x <- f.read.analyze.volume("finger.img")
volume <- x[ , , , 1]
contour3d(volume, 1:64, 1:64,1.5*(1:25), lev=c(500,700,800),alpha=c(0.2,0.5,1),color=c("white","red","green"))

another map style

ggmap

• just another (more fancy) map style
• for the preliminaries look here

map <- get_cloudmademap(bbox=c(left=12.25,bottom=51.285,
                        right=12.5,top=51.47),
                        maptype = 58461,
                        api_key = api_key,
                        zoom=13)
ggmap(map) + geom_point(aes(x=lon,y=lat),colour="red",data=gc) 

plot addresses as points on map using ggmap

ggmap

• given addresses
• find long and lat
• get a map of the city
• plot point at location of the addresses

library(ggmap)
library(stringr)
mystreets <- read.table("streets.txt",header=F,sep="|")
head(mystreets) ## addresses of physicians in Leipzig

                          V1
1  Rosa-Luxemburg-Str. 20-30
2            Marktstraße 2-6
3       Engelsdorfer Str. 21
4               Jahnallee 59
5            Johannisplatz 1
6         Hugo-Krone-Platz 9

## adds city and country to get complete addresses
ads <- paste(mystreets$V1, ', Leipzig, Germany',sep='')

## get rid of space at begin and end of the strings
ads <- str_trim(ads)


## get long and lat from google so you 
## need to have access to the internet
gc <- geocode(ads)

## get data of leipzig (limits of the map)
leipzig <- geocode("leipzig",output="more")


## get map from cloudmademap, also possible googlemaps etc
## for cloudemademap you need additionaly a api key
map <- get_cloudmademap(bbox=c(left=leipzig$west,
                          bottom=leipzig$south,
                          right=leipzig$east,
                          top=leipzig$north),
                        maptype = 997,api_key = api_key,zoom=13)

## plot all together
ggmap(map) + geom_point(aes(x=lon,y=lat),colour="red",data=gc)  

axis tick labels with comma as decimal point (german style)

german commas - commas as decimal point

• how to format the axis tick label, as example we use formatting numbers in german style (comma as dec)
• first we have to write a little function using the format function from the base package (which makes it very easy)

## the formatter function
germ_commas <- function(x,...){
  format(x, dec="," , trim = TRUE, scientific = FALSE, ...)
}

## create example data
x <- y <- rnorm(10)
z <- data.frame(x=x,y=y)

## create the plot using the function above inside the scales
ggplot(z,aes(x=x,y=y)) +
  geom_point() +
  scale_x_continuous(labels=germ_commas,breaks=seq(-2.5,2.5,by=0.25)) +
  scale_y_continuous(labels=germ_commas,breaks=seq(-2.5,2.5,by=0.25)) +
  opts(axis.text.x = theme_text(size=13,angle=90)) +
  opts(axis.text.y = theme_text(size=13))

Transformations axis vs values

log transformation: axis vs values

• sometimes it is useful to transform the values
• in ggplot you have the choice between transforming the values or the axis
• here is an example to compare both
• the data set movie is part of the ggplot package, additional used are the gridExtra package (grid.arrange), the scales package (labels formatter "dollar")
• in the beginning a cut the year var in 15 even spaced intervals
• on the lefthandside is the plot with the transformation of the axis, on the righthandside with transformation of values

movies$years <- cut(movies$year,breaks=15)

p <- ggplot(movies, aes(x=years,y=budget,colour=years)) +
## bar filling
  geom_boxplot(fill="black") + 
## add mean points
  stat_summary(fun.y="mean",geom="point") + 
 ## customize tick labels
  opts(axis.text.x = theme_text(angle=270,hjust=0,vjust=1)) +
## get rid of the legend
  opts(legend.position="none") 

## log transformation of the y axis
p1 <- p + scale_y_continuous(labels=dollar,trans="log")

## compare with
## here we do the transformation while mapping:
p2 <- ggplot(movies, aes(x=years,y=log(budget),colour=years)) +
  geom_boxplot(fill="black") +
  stat_summary(fun.y="mean",geom="point") +
  opts(axis.text.x = theme_text(angle=270,hjust=0,vjust=1)) +
  opts(legend.position="none")

## plot both graphs on one page
grid.arrange(p1,p2,ncol=2)

Dotplot with facets

• x contains variables riskgroup, aod (for age of death) and sex

ggplot(x,aes(x=riskgroup,y=aod,colour=riskgroup)) + 
  geom_dotplot() + 
  facet(~sex) +
  xlab("Risk Group") +
  ylab("Age of Death")