A simple custom R contour plotter

Base R and other contributed R packages have a wide range of functions for graphing contour plots and surfaces.

In working on a project I needed to create a type of contour plot, but in which I could be sure to represent each cell in the plot with its own color without interpolation.

I wasn’t sure exactly what function would meet that need, so I decided to make my own. (It was pretty easy.)

surfplot<-function(Z,x,y,...){
  if(hasArg(xlab)) xlab<-list(...)$xlab
  else xlab<-"x"
  if(hasArg(ylab)) ylab<-list(...)$ylab
  else ylab<-"y"
  if(hasArg(cex.axis)) cex.axis<-list(...)$cex.axis
  else cex.axis<-1
  if(hasArg(digits)) digits<-list(...)$digits
  else digits<-2
  if(hasArg(title)) title<-list(...)$title
  else title<-"log(prob)"
  if(hasArg(cols)) col<-list()$cols
  else cols<-hcl.colors(10)
  if(hasArg(min_z)) min_z<-list(...)$min_z
  else min_z<-min(Z)
  if(hasArg(max_z)) max_z<-list(...)$max_z
  else max_z<-max(Z)
  if(hasArg(asp)) asp<-list(...)$asp
  else asp<-NULL
  bp<-diff(range(x))*0.05
  plot(NA,xlim=c(min(x),max(x)+2*bp),
    ylim=c(min(y),max(y)+0.05*diff(range(y))),
    bty="n",xlab=xlab,ylab=ylab,axes=FALSE,asp=asp)
  axis(1,at=seq(min(x),max(x),length.out=5),
    cex.axis=cex.axis)
  axis(2,at=seq(min(y),max(y),length.out=5),
    cex.axis=cex.axis)
  cols<-colorRampPalette(cols)(length(y))
  levs<-seq(min_z,max_z,length.out=length(y))
  dx<-mean(x[2:length(x)]-x[1:(length(x)-1)])
  dy<-mean(y[2:length(y)]-y[1:(length(y)-1)])
  for(i in 1:length(x)) for(j in 1:length(y)){
    COL<-cols[which.min(abs(Z[i,j]-levs))]
    polygon(x[i]+c(-dx/2,dx/2,dx/2,-dx/2),
      y[j]+c(-dy/2,-dy/2,dy/2,dy/2),col=COL,
      border=FALSE)
  }
  for(i in 1:length(cols)){
    polygon(max(x)+bp+c(-bp/2,bp/2,bp/2,-bp/2),
      y[i]+c(-dy/2,-dy/2,dy/2,dy/2),col=cols[i],
      border=FALSE)
  }
  polygon(max(x)+bp+c(-bp/2,bp/2,bp/2,-bp/2),
    c(min(y)-dy/2,min(y)-dy/2,max(y)+dy/2,max(y)+dy/2))
  axis(4,at=seq(min(y),max(y),length.out=5),
    labels=round(seq(min_z,max_z,length.out=5),
      digits),las=1,cex.axis=0.6,
    pos=max(x)+1.5*bp)
  text(x=max(x)+bp,y=max(y),title,pos=3,
    cex=0.7)
}

For fun, let’s do a bivariate normal probability density.

x1<-x2<-seq(-2,2,length.out=200)
d<-matrix(NA,length(x1),length(x2))
v<-matrix(c(1,0.7,0.7,1),2,2)
for(i in 1:length(x1)) for (j in 1:length(x2))
  d[i,j]<-mnormt::dmnorm(c(x1[i],x2[j]),varcov=v)
par(mar=c(4.1,4.1,1.1,1.1))
surfplot(d,x1,x2,xlab=expression(x[1]),
  ylab=expression(x[2]),title="density")

The behavior of this function is probably closest to graphics::filled.contour in base R.

data(volcano)
par(mar=c(2.6,2.6,1.1,1.1))
filled.contour(volcano)

For fun, let’s try our custom function with the same data:

data(volcano)
par(mar=c(2.6,2.6,1.1,1.1))
cols<-colorRampPalette(c("#FFFFED","orange","darkred"))(100)
surfplot(volcano,
  seq(0,1,length.out=nrow(volcano)),
  seq(0,1,length.out=ncol(volcano)),
  cols=cols,title="elevation")

Ha! For guessing at the color palette, that’s closer than I expected.