Hello Roger, Thanks for the suggestions.
I finally managed to do it using the output of kde2d - The code is pasted below. Actually this made me realize that the outcome of kde2d can be quite influenced by outliers if a boundary box is not given (try running the code without the boundary box, e.g.lims, and you will see). library(MASS) x1 = rnorm(10000) x2 = rnorm(10000) nBins=200 z1 = kde2d(x1,x2,n=nBins, lims=c(-4,4,-4,4)) x1.cut = cut(x1, seq(z1$x[1], z1$x[nBins],len=(nBins+1) )) x2.cut = cut(x2, seq(z1$y[1], z1$y[nBins],len=(nBins+1) )) xy.cuts = data.frame(x1.cut,x2.cut, ord=1:(length(x1.cut)) ) density = data.frame( X=rep(factor(levels(x1.cut)),rep(nBins,nBins) ), Y=rep(factor(levels(x2.cut)),nBins) , Z= as.vector(z1$z)) xy.density = merge( xy.cuts, density, by=c(1,2), sort=FALSE, all.x=TRUE) xy.density = xy.density[order(x=xy.density$ord),] xy.density$Z[is.na(xy.density$Z)]=0 Z.lim = quantile(xy.density$Z, prob=c(0.1)) to.plot = xy.density$Z > Z.lim[1] par(mfrow=c(1,2)) plot(x1,x2, xlim=c(-3,3) ,ylim=c(-3,3)) plot(x1[to.plot], x2[to.plot], xlim=c(-3,3),ylim=c(-3,3)) On 19 November 2010 21:52, Roger Koenker <rkoen...@illinois.edu> wrote: > >> Also I wouldn't be surprised if there is a trick with quantile that >> escapes my mind. > > I would be surprised... this is all very dependent on how you do the > density estimation, but you might look at contourLines and then try > to use in.convex.hull() from tripack, but beware that things get more > complicated if the density is multi-modal. > > You might also look at one of the packages that do "data depth" > sorts of things. > > Roger Koenker > rkoen...@illinois.edu > > > > > > ______________________________________________ R-help@r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.
