Need code for plotting occupancy estimates in R

[hwaddle]

I am looking for a way to plot occupancy estimates with CI's in Program R. I can't figure out how to get the error bars since my data are not X-Y (the X is categorical, like species or year).

I would like to make plots similar to the ones MARK makes for survival estimates. Does anyone have a way to make plots like this in R?

Any code or hints would be greatly appreciated.

Hardin Waddle

[bacollier]

There is a wealth of literature and websites devoted to building graphics in R; including several email lists and forums you could search on for graphic advice in R.

I assume that your statement 'error bars since my data are not x-y' means that you have a bunch of occupancy estimates and SE (which you should be able to use to derive CI's) for each species or year and you just need to figure out how to plot them with a categorical x-axis.

Here is a simple example that should get you started.

data.one<-c(0.27, 0.21, 0.16, 0.14, 0.13, 0.09)
id<-c(1, 2, 3, 4, 5, 6)
#note--barplot2 is part of package-gregmisc so you might need to load it
# Barplot with Confidence intervals
lower<-c(data.one*0.85)
upper<-c(data.one*1.15)
barplot <- barplot2(data.one, beside=T, names.arg=id, plot.ci=T, ci.l=lower, ci.u=upper, col=c("white"), ylim=c(0, 0.35))

[hwaddle]

Well I admit I still have a lot to learn about making plots in R, but I have spent all day with books and online sources trying to solve this. The above code is fine for bar graphs, but I really wanted just a point at the mean of the estimate and error bars on it like you see in the MARK plots.

plotCI() seems like it would work, except it only wants X-Y data and I can't seem to trick it into accepting my categorical X data.

[cooch]

Well I admit I still have a lot to learn about making plots in R, but I have spent all day with books and online sources trying to solve this. The above code is fine for bar graphs, but I really wanted just a point at the mean of the estimate and error bars on it like you see in the MARK plots.

plotCI() seems like it would work, except it only wants X-Y data and I can't seem to trick it into accepting my categorical X data.

This really should be submitted to the various R discussion lists, since it has nothing in particular to do with PRESENCE, or occupancy analysis, but is a technical question on how to use R to handle plotting results. There are hundreds of sources of information out there on using R. You really should RTFM (where T = 'their') first.

[dhewitt]

Apologies in advance, but this comes up enough with RMark that I think it's worth giving people a simple start.

If you have estimates for a series of years and your data frame looks like this:

Year Estimate Lower Upper
2002 0.4743 0.2821 0.8462
2003 0.7767 0.641 0.9487
2004 0.7293 0.4872 0.9487
2005 0.5947 0.2821 0.9231
2006 0.5093 0.2051 0.8974

The plotCI() function in the plotrix package works well for this. Something like:

# Read in the data
ests <- read.table("estimates.txt", header=T)
# Eliminate white space (depends on output type)
oldpar <- par(mar=c(4.5, 5.5, 1, 1))
# Plot the estimates and CIs, but suppress axes
plotCI(ests$Year, ests$Estimate, yaxt="n", xaxt="n",
ui=ests$Upper, li=ests$Lower, gap=0.02, sfrac=0.01,
lwd=2, cex=1.6, pch=21, pt.bg="blue",
ylab="", xlab="", ylim=c(0, 1), cex.axis=1.5)
# Add the x-axis
axis(1, at=ests$Year, las=2, lwd=1.5, cex.axis=1.5)
# Add the y-axis
axis(2, las=1, lwd=1.5, cex.axis=1.5)
# Make the outline of the box a little thicker
box(lty=1, lwd=2)
# Add the y-axis label
mtext(side=2, line=4, "Probability of occupancy", cex=1.8)
# Reset the graphical parameters
par(oldpar)

If your x-axis categories are not years (integer variables), you'll have to do some more tweaking. One easy option is to coerce your factor (e.g., habitat types A, B, C, and D) into integers with as.integer() wrapped around the variable in the plotCI() call. Then the x-axis call above could be modified coarsely to make it what you want... e.g.,

axis(1, at=c(1, 2, 3, 4), lwd=1.5, cex.axis=1.5,
labels=c("A", "B", "C", "D"))