[R] 2D density tophat
Aaron Robotham
asgr at st-andrews.ac.uk
Fri Dec 5 16:51:13 CET 2008
In case anyone other than me was interested, a pretty efficient
circular tophat can be made using the fields function fields.rdist.near:
CircHat=function (x, y, h=1, gridres = c((max(x)-min(x))/25,(max(y)-
min(y))/25), lims = c(range(x), range(y)),density=FALSE)
{
require(fields)
nx <- length(x)
ny <- length(y)
n=c(1+(lims[2]-lims[1])/gridres[1],1+(lims[4]-lims[3])/gridres[2])
if (length(y) != nx)
stop("data vectors must be the same length")
if (any(!is.finite(x)) || any(!is.finite(y)))
stop("missing or infinite values in the data are not allowed")
if (any(!is.finite(lims)))
stop("only finite values are allowed in 'lims'")
gx <- seq(lims[1], lims[2], by = gridres[1])
gy <- seq(lims[3], lims[4], by = gridres[2])
fullgrid=expand.grid(gx,gy)
if (missing(h))
h <- c(bandwidth.nrd(x), bandwidth.nrd(y))
temp
=
table
(fields
.rdist
.near
(as
.matrix
(fullgrid
),as.matrix(cbind(x,y)),mean.neighbor=ceiling(length(x)*pi*h^2/
((lims[2]-lims[1])*(lims[4]-lims[3]))),delta=h)$ind[,1])
pad=rep(0,length(gx)*length(gy))
pad[as.numeric(names(temp))]=as.numeric(temp)
z <- matrix(pad, length(gx), length(gy))
if(density){z=z/(nx*pi*h^2)}
return=list(x = gx, y = gy, z = z)
}
It works in more or less the same way as kde2d, but by default it
returns counts not densities
Aaron
On 1 Dec 2008, at 11:46, Aaron Robotham wrote:
> Hello R users,
>
> I have successfully created a square (or more generally, rectangular)
> tophat smoothing routine based on altering the already available
> KDE2D. I would be keen to implement a circular tophat routine also,
> however this appears to be much more difficult to write efficiently (I
> have a routine, but it's very slow). I tried to create one based on
> using crossdist (in spatstat) to create a distance matrix between my
> data and the sampling grid, but it doesn't take a particularly large
> amount of data (or hi-res grid) for memory to be a big problem. The 2D
> density routines I have been able to find either don't support a
> simple tophat, or don't use the absolute distances between the
> sampling grid and the data. Should anyone know of more general 2D
> density routines that might support circular tophats, or know of a
> simple and efficient method for creating them, I would be very
> grateful.
>
> Thanks for your time,
>
> Aaron
>
> PS: I tried sending this on Friday originally, but as far as I know
> that didn't work, so should another post appear from me asking the
> same thing I apologise in advance.
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help at 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.
More information about the R-help
mailing list