Skip to content
Prev 139919 / 398502 Next

Smoothing z-values according to their x, y positions

Dear Bert,

Thanks for your reply - I indeed saw a lot of functions using:
help.search("smooth")

The problem is that most seem to not be very appropriate to what I'd like, or
they seem extremely complicated (e.g. gma). I am probably missing
something as I don't
see how to use Loess. From my poor understanding, it seems to be for 2D data.
Here I want to smooth the third "z" component.

In the meantime I adapted a function from Greg Warnes, which sort of
does the job (although the smoothing is not very nice).
So I'd be very happy to learn how do do something similar with loess!

x = runif(1000)
y = runif(1000)
z = x+y
hist2d.el(x,y,z, nbins=10)

## -- This function was adapted from Greg Warnes
hist2d.el <- function( x,y=NULL, z=NULL, nbins=200, my.cuts = NULL,
same.scale=FALSE, na.rm=TRUE, show=TRUE, col=c("black",
heat.colors(12)), ... )
  {
    if(is.null(y))
      {
        if(ncol(x) != 2) stop("If y is ommitted, x must be a 2 column matirx")
        y <- x[,2]
        x <- x[,1]
      }

    if(length(nbins)==1)
      nbins <- rep(nbins,2)

    nas <- is.na(x) | is.na(y)

    if(na.rm)
      {
        x <- x[!nas]
        y <- y[!nas]
      }
    else
      stop("missinig values not permitted if na.rm=FALSE")

    if(!is.null(my.cuts)){

      nbins = rep(length(my.cuts)-1,2)
      x.cuts=my.cuts
      y.cuts=my.cuts

    } else {

      if(same.scale)
        {
          x.cuts <- seq( from=min(x,y), to=max(x,y), length=nbins[1]+1)
          y.cuts <- seq( from=min(x,y), to=max(x,y), length=nbins[2]+1)
        }
      else
        {
          x.cuts <- seq( from=min(x), to=max(x), length=nbins[1]+1)
          y.cuts <- seq( from=min(y), to=max(y), length=nbins[2]+1)
        }
    }

    print(nbins)

    index.x <- cut( x, x.cuts, include.lowest=TRUE)
    index.y <- cut( y, y.cuts, include.lowest=TRUE)

    m <- matrix( 0, nrow=nbins[1], ncol=nbins[2],
                dimnames=list( levels(index.x),
                               levels(index.y) ) )

    m.sum <- matrix( 0, nrow=nbins[1], ncol=nbins[2],
                dimnames=list( levels(index.x),
                               levels(index.y) ) )

    for( i in 1:length(index.x) )
      m[ index.x[i], index.y[i] ] <-  m[ index.x[i], index.y[i] ] + 1

    for( i in 1:length(index.x) )
      m.sum[ index.x[i], index.y[i] ] <-  m.sum[ index.x[i], index.y[i] ] + z[i]

    m[which(m==0)]=1
    m = m.sum/(m)

    xvals <- x.cuts[1:nbins[1]]
    yvals <- y.cuts[1:nbins[2]]

    if(show)
      image( xvals,yvals, m, col=col,...)

                                        #
invisible(list(counts=m,x=xvals,y=yvals))
    invisible(m)
  }
On 19/03/2008, Bert Gunter <gunter.berton at gene.com> wrote: