My original post on this topic was far more detailed but I got very little response so my second attempt was much more focussed. Here it is to see if you are able to offer more help:
Sent: 08 May 2009 12:46 To: r-sig-geo@stat.math.ethz.ch Subject: [R-sig-Geo] Calculating Risk at point locations (withoutrasterising) I would appreciate advice on how to do a spatial statistics task: detecting points at a high risk of Cases (and cluster analysis). I have semi-randomly located vector point data with 2 attributes (Cases and Controls) at each point, i.e. marked point data. I want to apply spatial smoothing as I hypothesize that the underlying process is such that clusters of cases are of interest (although I also know that the process operates on a small enough scale that isolated points with high numbers of cases are also of interest). Having read the book by Bivand et al "Applied Spatial Data Analysis with R" I was following the approach of an inhomogeneous Poisson process analysis (p.175) to calculate a risk image of log(cases/controls), which incidentally gives a few challenges with divide by zeros. This risk image (I used data class ppp in R library SpatStat as I had marked point data) is a raster interpolated image, even though I originated with vector point data. However, the aim of my analysis is to score each of the (vector) points with its risk value, i.e. give each of my original vector points a new attribute called Risk. For my application I am not interested in what happens in between my point locations (as the point locations have a meaning). So my problem is how to get the Risk value at the original point locations. My current approach involves converting the points to raster, so I lose the precise point location and I have a large data set so in practice cannot set a very fine resolution. So my key question is: please, can someone offer me advice on how to get the Risk value at the original point locations? I wonder if I can avoid raster and do all processing in a vector map (perhaps by interpolation or spline approximation). This would mean I can avoid resampling a raster map back to the original vector point locations as occasionally the locations are so close (less than resolution) this will create inaccuracy. (Background: I have not done much spatial statistics before so am trying to keep it simple; I have used GRASS; not much experience with R). When I originally import the map into R it has class SpatialPointsDataFrame. Below are a couple of code examples: > Case <- smooth.ppp(Cases,75, eps=50) # smooth my Case data with a kernel width of 75m; 50m pixel dimension in the resulting smoothed ppp. > plot(eval.im( log((Case+9e-6)/(Control+1e-5)) + 0.3521755)) # to get the Risk map; the 9e6- and 1e-5 are how I avoid divide by zero type problems. -----Original Message----- From: Paul Hiemstra [mailto:p.hiems...@geo.uu.nl] Sent: 12 May 2009 10:37 To: Elliott,MR,Mike,BMK6 R Cc: r-sig-geo@stat.math.ethz.ch Subject: Re: [R-sig-Geo] Interpolation / smoothing of points (without raster grids) mike.elli...@openreach.co.uk wrote: > Hello - can anyone suggest functions/packages that allow interpolation > / smoothing of point data (without using raster grids)? i.e. which > give the results at the original point locations. > > Many thanks, Mike Elliott. > > _______________________________________________ > R-sig-Geo mailing list > R-sig-Geo@stat.math.ethz.ch > https://stat.ethz.ch/mailman/listinfo/r-sig-geo > If you could tell us more in detail what you want to do, we could offer you more guidance. To add to the promotion of packages that we wrote ourselves (edzer :)), you can also use automap, a package for automatic interpolation based on gstat and sp (two other packages). For cross-validation (if this is what you need) the following script will do the trick: library(automap) # Load the data data(meuse) coordinates(meuse) = ~x+y data(meuse.grid) gridded(meuse.grid) = ~x+y # Perform cross-validation kr.cv = autoKrige.cv(log(zinc)~1, meuse, model = c("Exp")) kr_dist.cv = autoKrige.cv(log(zinc)~sqrt(dist), meuse, model = c("Exp")) kr_dist_ffreq.cv = autoKrige.cv(log(zinc)~sqrt(dist)+ffreq, meuse, model = c("Exp")) # Compare the results compare.cv(kr.cv, kr_dist.cv, kr_dist_ffreq.cv) compare.cv(kr.cv, kr_dist.cv, kr_dist_ffreq.cv, bubbleplots = TRUE) compare.cv(kr.cv, kr_dist.cv, kr_dist_ffreq.cv, bubbleplots = TRUE, col.names = c("OK","UK1","UK2")) compare.cv(kr.cv, kr_dist.cv, kr_dist_ffreq.cv, bubbleplots = TRUE, col.names = c("OK","UK1","UK2"), plot.diff = TRUE) automap is available from CRAN. cheers and good luck, Paul -- Drs. Paul Hiemstra Department of Physical Geography Faculty of Geosciences University of Utrecht Heidelberglaan 2 P.O. Box 80.115 3508 TC Utrecht Phone: +3130 274 3113 Mon-Tue Phone: +3130 253 5773 Wed-Fri http://intamap.geo.uu.nl/~paul _______________________________________________ R-sig-Geo mailing list R-sig-Geo@stat.math.ethz.ch https://stat.ethz.ch/mailman/listinfo/r-sig-geo