#' Plot Quantities of Interest in a Zelig-fashion
#'
#' Various graph generation for different common types of simulated results from
#' Zelig
#' @usage simulations.plot(y, y1=NULL, xlab="", ylab="", main="", col=NULL, line.col=NULL, axisnames=TRUE)
#' @param y A matrix or vector of simulated results generated by Zelig, to be
#' graphed.
#' @param y1 For comparison of two sets of simulated results at different
#' choices of covariates, this should be an object of the same type and
#' dimension as y. If no comparison is to be made, this should be NULL.
#' @param xlab Label for the x-axis.
#' @param ylab Label for the y-axis.
#' @param main Main plot title.
#' @param col A vector of colors. Colors will be used in turn as the graph is
#' built for main plot objects. For nominal/categorical data, this colors
#' renders as the bar color, while for numeric data it renders as the background
#' color.
#' @param line.col A vector of colors. Colors will be used in turn as the graph is
#' built for line color shading of plot objects.
#' @param axisnames a character-vector, specifying the names of the axes
#' @return nothing
#' @author James Honaker
simulations.plot <-function(
y, y1=NULL,
xlab="", ylab="",
main="",
col=NULL,
line.col=NULL,
axisnames=TRUE
) {
## Univariate Plots ##
if(is.null(y1)){
if (is.null(col))
col <- rgb(100,149,237,maxColorValue=255)
if (is.null(line.col))
line.col <- "black"
# Character
if (is.character(y)) {
# Try to cast y as integers, note that this is not always possible for the
# general case of characters
newy <- tryCatch(
as.numeric(y),
warning = function (w) NULL,
error = function (e) NULL
)
# If:
# newy is not NULL (can be cast as a numeric) AND
# newy is actually a collection of integers (not just numeric)
# Then:
# we can tabulate (so sick)
if (!is.null(newy) && all(as.integer(y) == y)) {
# Create a sequence of names
nameseq <- paste("Y=", min(newy):max(newy), sep="")
# Set the heights of the barplots.
# Note that tablar requires that all out values are greater than zero.
# So, we subtract the min value (ensuring everything is at least zero)
# then add 1
bar.heights <- tabulate(newy - min(newy) + 1) / length(y)
# Barplot with (potentially) some zero columns
output <- barplot(
bar.heights,
xlab=xlab, ylab=ylab, main=main, col=col[1],
axisnames=axisnames, names.arg=nameseq
)
}
# Otherwise, we stick with old-style tables
else {
y <- if (is.null(levels(y)))
factor(y)
else
factor(y, levels = levels(y))
bar.heights <- table(y)/length(y)
bar.names <- paste("Y=", names(bar.heights), sep="")
output <- barplot(
bar.heights,
xlab=xlab, ylab=ylab, main=main, col=col[1],
axisnames=axisnames, names.arg=bar.names
)
}
}
## Numeric
else if(is.numeric(y)){
den.y <- density(y)
output <- plot(den.y, xlab=xlab, ylab=ylab, main=main, col=line.col[1])
if(!identical(col[1],"n")){
polygon(den.y$x, den.y$y, col=col[1])
}
}
## Comparison Plots ##
}
else{
## Character - Plot and shade a matrix
if(is.character(y) & is.character(y1) & length(y)==length(y1) ){
newy<-trunc(as.numeric(y))
newy1<-trunc(as.numeric(y1))
yseq<-min(c(newy,newy1)):max(c(newy,newy1))
nameseq<- paste("Y=",yseq,sep="")
n.y<-length(yseq)
colors<-rev(heat.colors(n.y^2))
lab.colors<-c("black","white")
comp<-matrix(NA,nrow=n.y,ncol=n.y)
for(i in 1:n.y){
for(j in 1:n.y){
flag<- newy==yseq[i] & newy1==yseq[j]
comp[i,j]<-mean(flag)
}
}
old.pty<-par()$pty
old.mai<-par()$mai
par(pty="s")
par(mai=c(0.3,0.3,0.3,0.1))
image(z=comp, axes=FALSE, col=colors, zlim=c(min(comp),max(comp)),main=main )
locations.x<-seq(from=0,to=1,length=nrow(comp))
locations.y<-locations.x
for(m in 1:n.y){
for(n in 1:n.y){
text(x=locations.x[m],y=locations.y[n],labels=paste(round(100*comp[m,n])), col=lab.colors[(comp[m,n]> ((max(comp)-min(comp))/2) )+1])
}
}
axis(side=1,labels=nameseq, at=seq(0,1,length=n.y), cex.axis=1, las=1)
axis(side=2,labels=nameseq, at=seq(0,1,length=n.y), cex.axis=1, las=3)
box()
par(pty=old.pty,mai=old.mai)
## Numeric - Plot two densities on top of each other
}else if(is.numeric(y) & is.numeric(y1)){
if(is.null(col)){
col<-c("blue","red")
}else if(length(col)<2){
col<-c(col,col)
}
if(is.null(col)){
semi.col.x <-rgb(142,229,238,150,maxColorValue=255)
semi.col.x1<-rgb(255,114,86,150,maxColorValue=255)
col<-c(semi.col.x,semi.col.x1)
}else if(length(col)<2){
col<-c(col,col)
}
den.y<-density(y)
den.y1<-density(y1,bw=den.y$bw)
all.xlim<-c(min(c(den.y$x,den.y1$x)),max(c(den.y$x,den.y1$x)))
all.ylim<-c(min(c(den.y$y,den.y1$y)),max(c(den.y$y,den.y1$y)))
output<-plot(den.y,xlab=xlab,ylab=ylab,main=main,col=col[1],xlim=all.xlim,ylim=all.ylim)
par(new=TRUE)
output<-plot(den.y1,xlab=xlab,ylab=ylab,main="",col=col[2],xlim=all.xlim,ylim=all.ylim)
if(!identical(col[1],"n")){
polygon(den.y$x,den.y$y,col=col[1])
}
if(!identical(col[1],"n")){
polygon(den.y1$x,den.y1$y,col=col[2])
}
}
}
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.