addcontours {MQMF}

addcontours simplifies adding contours to an xy plot of points

Description

addcontours is used to add contours to a dense plot of xy points such as might be generated when conducting an analysis of the the uncertainty associated with a stock assessment, or other analysis using a bootstrap, a Bayesian MCMC, or even using asymptotic errors and sampling from a multi-variate normal distribution. addcontours first uses the kde2d function from the MASS package to translate the density of points into 2-D kernal densities, and then searches through the resulting densities for those points that would identify approximate contours. Finally it calls the contour function to add the identified approximate contours to the xy plot.

Usage

addcontours(
  xval,
  yval,
  xrange,
  yrange,
  ngrid = 100,
  contval = c(0.5, 0.9),
  lwd = 1,
  col = 1
)

Arguments

xval

the vector of x-axis values, one half of the data pairs

yval

the vector of y-axis values, the other half of the data

xrange

the range of x-axis data included in the graph

yrange

the range of y-axis data included in the graph

ngrid

the number of subdivisions by which to split the data along each axis; defaults to 100

contval

the contour values, defaults to those containing 50 and 90 percent i.e. c(0.5, 0.9)

lwd

the width of the contour lines, defaults=1

col

the col of the contour lines, default=1

Value

nothing but it does add contours to a plot of points

Examples

 library(mvtnorm)
 library(MASS)
 data(abdat)
 param <- log(c(r= 0.42,K=9400,Binit=3400,sigma=0.05)) 
 bestmod <- nlm(f=negLL,p=param,funk=simpspm,indat=abdat, 
                hessian=TRUE,logobs=log(abdat$cpue),
                typsize=magnitude(param),iterlim=1000)
 optpar <- bestmod$estimate
 vcov <- solve(bestmod$hessian)      # solve inverts matrices
 columns <- c("r","K","Binit","sigma")
 N <- 1000  # the contours improve as N increases; try 5000
 mvnpar <- matrix(exp(rmvnorm(N,mean=optpar,sigma=vcov)),
                  nrow=N,ncol=4,dimnames=list(1:N,columns))
 xv <- mvnpar[,"K"]
 yv <- mvnpar[,"r"]
 # plotprep(width=6,height=5,newdev = FALSE)
 plot(xv,yv,type="p") # use default 0.5 and 0.9 contours.
 addcontours(xv,yv,range(xv),range(yv),lwd=2,col=2)
 points(mean(xv),mean(yv),pch=16,cex=1.5,col=2)