## Lecture 4: Graphics ## Author: Adrian Waddell ## Date: November 6, 2012 ## Scatterplots x <- c(5,2,5,6,2,6,8,3,5,6) y <- c(3,7,5,6,3,3,9,6,3,4) dev.new(width=5,height=6,pointsize=18) plot(x,y) ## pch plot(x,y,pch=16) dev.new() par(mai=c(1,1,1,1)/2) xg <- rep(1:5,5) yg <- rep(1:5,each=5) plot(xg, yg, pch=1:25, axes=FALSE, ylab='', xlab='', main="pch", xlim=c(0,6), ylim=c(0,6), bg='#E2307E') box() text(xg,yg,labels=1:25,pos=1, offset=0.8) ?text plot(x,y,pch=1:16) plot(x,y,pch='$') ## Colors plot(x,y,pch=16,col='forestgreen') grep('green',colors(),fixed=TRUE,value=TRUE) cols <- colors()[sample(1:length(colors()),length(x))] plot(x,y, pch=16, col=cols, bg='magenta') orange <- rgb(245,104,52, maxColorValue=255) orange <- '#F56834' # Hexadecimal representation plot(1,1,pch=16,col=orange) ??hex as.hexmode(c(245,104,52)) ?hsl ?hsv ## RColorBrewer install.packages("RColorBrewer") library(RColorBrewer) example(RColorBrewer) help(package="RColorBrewer") dev.new(pointsize=18) par(mai=c(0.2,1,0.2,0)) ?par display.brewer.all() display.brewer.all(type="div") display.brewer.all(type="seq") display.brewer.all(type="qual") ## You do not need to understand this at this point install.packages('ggplot2') library(ggplot2) cols <- brewer.pal(11, name="PiYG") df <- data.frame(y= rep(0,11), x=1:11) df ggplot(df,aes(x,y,fill=factor(x)))+ geom_raster()+ scale_fill_manual(values=cols) + theme(axis.ticks = element_blank(), axis.text = element_blank())+ xlab('')+ylab('') cols2 <- colorRampPalette(brewer.pal(11, name="PiYG"))(40) df2 <- data.frame(y= rep(0,40), x=1:40) ggplot(df2,aes(x,y,fill=factor(x)))+ geom_raster()+ scale_fill_manual(values=cols2) + theme(axis.ticks = element_blank(), axis.text = element_blank())+ xlab('')+ylab('') plot(1:40,1:40,pch=19, col=cols2) ## Magnifying plot symbols ?par plot(x,y,pch=16,col="#E2307E",cex=3) plot(x,y,pch=16,col="#E2307E",cex=3,xlim=c(0,10), ylim=c(0,10)) plot(x,y,pch=16,col="#E2307E",cex=3,xlim=c(0,10), ylim=c(0,10), xlab="age in months", ylab="# doctor visits", main="Babies") plot(x,y,pch=16,cex=12,col=rgb(226,48,126,alpha=60,maxColorValue=255)) as.hexmode(60) plot(x,y,pch=16,col="#E2307E3C",cex=12) ## plot piece by piece ?points dev.new(width=5,height=6,pointsize=18) par(mai=c(.6,.1,.2,.1)) plot(0,0,type='n',xlim=c(0,1), ylim=c(0,1),axes=FALSE,xlab='',ylab='') points(.5,.7,pch=1,cex=16) lines(c(.5,.5),c(.8,.7),lwd=3) points(c(.4,.6),c(.8,.8), pch=21,cex=5,col='black',bg='#8AD3CD') points(c(.4,.6)+.02,c(.7,.7) + .1, pch=16,cex=2,col='black') polygon(c(.4,.6,.6,.4),c(.4,.4,.2,.2)+.15,col='#938A86') segments(c(.45,.55),c(.35,.35),c(.38,.62),c(.1,.1),lwd=25,col='#C4D5A5') segments(c(.41,.59),c(.5,.5),c(.2,.8),c(.6,.6),lwd=25,col='#F8C850') text(.8,.95,"This is Freddie") box(col='#A885B9',lwd=3) -axis(side=1) ## Motorcycles library(AER) data(Journals) data(package = "AER") data(DoctorVisits) head(DoctorVisits) with(DoctorVisits,ftable(visits,private,gender)) data(MotorCycles) help("MotorCycles") is.ts(MotorCycles) MotorCycles ?rect bright <- '#F8C850' dark <- '#F8EBCA' dev.new(width=12,height=5,pointsize=12) par(mar=c(5, 6, 4, 2) + 0.1) plot(MotorCycles, main=" Motor Cycles in The Netherlands", type='b',pch=16,xlab="Year",ylab="# of Motor Cycles\n in thousands") abline(v = c(1950,1960,1970,1980,1990)) rect(par("usr")[1],par("usr")[3],1950,par("usr")[4],col = bright, lwd=0) rect(1950,par("usr")[3],1960,par("usr")[4],col = dark, lwd=0) rect(1960,par("usr")[3],1970,par("usr")[4],col = bright, lwd=0) rect(1970,par("usr")[3],1980,par("usr")[4],col = dark, lwd=0) rect(1980,par("usr")[3],1990,par("usr")[4],col = bright, lwd=0) rect(1990,par("usr")[3],par("usr")[2],par("usr")[4],col = dark, lwd=0) box() lines(MotorCycles, lwd=2) points(MotorCycles, pch=19) ## Boxplot data(SmokeBan) head(SmokeBan) dim(SmokeBan) help(SmokeBan) dev.new(pointsize=12) with(SmokeBan,boxplot(age~smoker)) with(SmokeBan,boxplot(age~smoker+ban)) with(SmokeBan,boxplot(age~smoker+ban+education)) par(mar=c(5, 10, 4, 2) + 0.1) with(SmokeBan,boxplot(age~smoker+ban+education,horizontal=TRUE, las=1)) with(SmokeBan,boxplot(age~education+ban+smoker,horizontal=TRUE, las=1)) with(SmokeBan,boxplot(age~education+ban+smoker,horizontal=TRUE, las=1, at = c(1:5,7:11,14:18,20:24))) rect(par("usr")[1],19,par("usr")[2],par("usr")[4],col="#F7D7E5",lwd=0) rect(par("usr")[1],12.5,par("usr")[2],19,col="#F0ACCA",lwd=0) rect(par("usr")[1],6,par("usr")[2],12.5,col="#EDE4F2",lwd=0) rect(par("usr")[1],par("usr")[3],par("usr")[2],6,col="#A885B9",lwd=0) with(SmokeBan,boxplot(age~education+ban+smoker,horizontal=TRUE, las=1, at = c(1:5,7:11,14:18,20:24),add=TRUE, pch=19, col=rep(c('#F56834','#F8C850','#8E9A7A','#8AD3CD','#28457C'),4))) ## Histogram hist(Journals$price) hist(Journals$price, nclass=50) hist(Journals$price, nclass=50, freq=FALSE) ??lognormal with(Journals,curve(dlnorm(x,meanlog=mean(log(price)),sdlog=sd(log(price)), log=FALSE),from=0,to=2000, add=TRUE,lwd=4,col="#E2307E")) ## Kernel Density Estimation ?density lines(density(Journals$price, bw=30, kernel="gaussian"), col="#28457C", lwd=4) with(Journals,curve(dlnorm(x,meanlog=mean(log(price)),sdlog=sd(log(price)), log=FALSE),from=0,to=2000, add=TRUE,lwd=4,col="#E2307E")) ?par par(mfrow=c(2,3)) for(bw in c(50,80,100,200,300,500)) { hist(Journals$price, nclass=50, freq=FALSE, main=paste("bandwidth=",bw,sep='')) with(Journals,curve(dlnorm(x,meanlog=mean(log(price)),sdlog=sd(log(price)), log=FALSE),from=0,to=2000, add=TRUE,lwd=4,col="#E2307E")) lines(density(Journals$price, bw=bw, kernel="gaussian"), col="#28457C", lwd=4) } install.packages("MASS") library("MASS") ?density ?truehist truehist(Journals$price,h=50,x0=0,xlim=c(0,2000)) ## Curve curve(sin(x),from=0,to=2*pi,xaxt='n', lwd=4, col="#E2307E") axis(1, at=c(0, pi/2, pi, 3/2*pi, 2*pi), labels=expression(0,frac(pi,2),pi,frac(3*pi,2),2*pi),mgp=c(3, 2, 0)) abline(h = 0) ## Mathematical Annotation in R ?expression ?plotmath ## Distributions ?dnorm ??'Binomial Distribution' curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E") text(2,0.3, expression(paste('f(x)=',frac(1,sqrt(2*pi)),e^{-frac(1,2)*x^2})), cex=2) p <- c(pnorm(q=3,mean=0,sd=1),.75,.5,.25) cols <- c('#BE8ED5','#DBBDEA','#E3D5EA','#E9E7EA') for(i in 1:4){ q <- qnorm(p=p[i],mean=0,sd=1) area <- curve(dnorm(x, mean=0, sd=1), from=-3,to=q, add=TRUE, type='n') polygon(x=c(area$x,q,-3),y=c(area$y,0,0), col=cols[i]) } curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="black", add=TRUE) text(c(-1.1,-.35,.35,1.1), rep(0.1,4),rep('25%',5)) ## Explain d curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E", main="dnorm") text(2,0.3, expression(paste('f(x)=',frac(1,sqrt(2*pi)),e^{-frac(1,2)*x^2})), cex=2) ## Explain q curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E", main="qnorm", axes=FALSE) box() q <- qnorm(p=.38,mean=0,sd=1) area <- curve(dnorm(x, mean=0, sd=1), from=-3,to=q, add=TRUE, type='n') polygon(x=c(area$x,q,-3),y=c(area$y,0,0), col='#F8C850') text(-1,0.1,'p') curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E", main="dnorm", add=TRUE) axis(1,at=q,"qnorm(p)", cex=3) ## explain p curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E", main="pnorm", axes=FALSE) box() q <- qnorm(p=.22,mean=0,sd=1) area <- curve(dnorm(x, mean=0, sd=1), from=-3,to=q, add=TRUE, type='n') polygon(x=c(area$x,q,-3),y=c(area$y,0,0), col='#F8C850') text(-1.25,0.1,'pnorm(q)') axis(1,at=q,"q", cex=3) curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E", add=TRUE) ## Explain r hist(rnorm(1000, mean=0, sd=1), freq=FALSE, ylim = c(0,0.6),nclass=30) curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E", add=TRUE) ## Save Images setwd('~/Desktop/') ?pdf() pdf(file="fig01.pdf",width=10,height=5) hist(rnorm(1000, mean=0, sd=1), freq=FALSE, ylim = c(0,0.6),nclass=30) curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E", add=TRUE) dev.off() pdf(file="fig02.pdf",width=10,height=5,pointsize=20) hist(rnorm(1000, mean=0, sd=1), freq=FALSE, ylim = c(0,0.6),nclass=30) curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E", add=TRUE) dev.off() pdf(file="fig03.pdf",width=10,height=5,pointsize=5) hist(rnorm(1000, mean=0, sd=1), freq=FALSE, ylim = c(0,0.6),nclass=30) curve(dnorm(x, mean=0, sd=1), from=-3,to=3,lwd=4, col="#E2307E", add=TRUE) dev.off() ?dev.new ## Maps install.packages('maps') help(package='maps') library(maps) ?map map("state", interior = FALSE) map("state", boundary = FALSE, lty = 2, add = TRUE) # ggplot2 install.packages("ggplot2") ## lattice install.packages("lattice") ## rggobi ## see installation instructions http://www.ggobi.org/ install.packages("rggobi") ## RnavGraph ## See installation instructions in vignette ## http://cran.r-project.org/web/packages/RnavGraph/vignettes/RnavGraph.pdf install.packages("RnavGraph","RnavGraphImageData") ## Touching up graphics with inkscape or Adobe illustrator #setwd("~/Desktop/") #svg(file="sincurve.svg") setwd("~/Desktop/") pdf(file="sincurve.pdf") ## SVG or PDF curve(sin(x),from=0,to=2*pi,xaxt='n', lwd=4, col="#E2307E") axis(1, at=c(0, pi/2, pi, 3/2*pi, 2*pi), labels=expression(0,frac(pi,2),pi,frac(3*pi,2),2*pi),mgp=c(3, 2, 0)) abline(h = 0) dev.off()