[R-sig-Geo] Creating density heatmaps for geographical data

Karl Ove Hufthammer karl at huftis.org
Mon Oct 18 13:43:55 CEST 2010


Michael Sumner wrote:

> Just answering your first question, you could call sm.density for each
> point individually, then trim back to whatever level you want. Sum each
> grid for each point and you might get a better "trim".
> 
> Otherwise, you'd need to hack the density functions in kde2d or
> sm.density, or maybe in KernSmooth to get non-Gaussians, or otherwise to
> ensure the density is only calculated for a smaller local region around
> the point.
> 
> KernSmooth package may well have more control over this stuff.

Thanks for the suggestion. KernSmooth did not have any more options, but I 
finally found a solution I’m mostly happy with. ‘density.ppp’ in ‘spatstat’ 
estimates the intensity function, i.e., the expected number of observations
per area, and accepts weighted observations. It’s possible to specify both
the resolution and the bandwidth, but not the kernel. And 'density.ppp' is 
blazingly fast!

'ppp' objects used by 'spatstat' have an annoying structure, but it’s not 
difficult to convert 'sp' objects to this structure. Be aware of duplicate 
points, though.

Since I use an UTM projection, I get estimates of the expected numbers of 
observations per square meter. Multiplying this with 1000^2 gives the 
expected number per square kilometer, if I have understood things correctly.
I then replace all grid cells with a value < 1 with NA (modifying the ‘v‘ 
variables on the output object directly). Using ‘con2tr’ converts this into
something ‘levelplot’ in ‘lattice’ accepts. Putting this on top of a GSHHS 
map and using 'bpy.colors' basically gives me what I need.

-- 
Karl Ove Hufthammer



More information about the R-sig-Geo mailing list