Ellipses of calculation errors in R

Question

I am trying to build an U-Pb isochron, which is basically an xy plot, with error ellipses around each datapoint to indicate the accuracy of this data point. Here is an example of what I'm trying to achieve, with the only difference being that I have much more data:

I have symmetrical + and - errors for each data point in the spaces x and y that I would like to use to dictate the shape of the error ellipse.

Below is the code that I have written so far. This computes my data, but instead of error ellipses, I currently have errors.

# Select data file. fname<-file.choose() # Rename imported data file. upbiso <- read.table(file=fname, header= TRUE) # Attaches data file as default data source. attach(upbiso) # Reads and display column headers in console. names(upbiso) # Sets margins around plot. par(mar=c(5,5,4,2)) # Plots 238U/206Pb vs 207Pb/206Pb with superscript notation for isotopic masses, xlim and ylim sets min and max limit for axes, las rotates y axis labels, cex.lab increases the size of the axis labels. plot(UPb, PbPb, xlab = expression({}^238*"U/"*{}^206*"Pb"), ylab = expression({}^207*"Pb/"*{}^206*"Pb"), xlim = c(0,2500),ylim = c(0, 1), las=1, cex.lab = 1.5) # Opens Oceanographic package to run errorbars function and runs 1 sigma percent error bars for xy co-ordinates. library(oce) errorbars(UPb,PbPb,UPbErrP,PbPbErrP, percent=T) 

How can I swap errors for error ellipses?

Here is a link from Google docs to my data, which is in .txt format.

https://docs.google.com/file/d/0B75P9iT4wwlTNUpQand2WVRWV2s/edit?usp=sharing

The columns are UPb , UPbErrP (error in UPb with 1 sigma%), UPbErrAbs (absolute error in UPb), and then again for PbPb data.

If you need me to clarify anything, just let me know and I will do my best.

Answer 1

A few months ago, I wrote a small function to draw ellipses to answer another question . On the contrary, this is not much, we can achieve something useful for your problem, I think.

 ellipse <- function(a,b,xc,yc,...){ # a is the length of the axis parallel to the x-axis # b is the length of the axis parallel to the y-axis # xc and yc are the coordinates of the center of the ellipse # ... are any arguments that can be passed to function lines t <- seq(0, 2*pi, by=pi/100) xt <- xc + a*cos(t) yt <- yc + b*sin(t) lines(xt,yt,...) } plot(UPb, PbPb, pch=19, xlab = expression({}^238*"U/"*{}^206*"Pb"), ylab = expression({}^207*"Pb/"*{}^206*"Pb"), xlim = c(0,2500),ylim = c(0, 1), las=1, cex.lab = 1.5) apply(upbiso, 1, function(x)ellipse(a=x[2]*x[1]/100, b=x[5]*x[4]/100, xc=x[1], yc=x[4], col="red"))