[R-sig-Geo] kriging inside borders

[Roger Bivand]

On Mon, 23 May 2005, Chloe ARCHINARD wrote:

> Here is in attached file the map plotted by 
> 
> > plot(x, col="grey")
> 
> > points(xp[is.na(res)], col="red")
> 
> > points(xp[!is.na(res)], col="blue")
> 
> It seems to be OK! When I use krige.conv to make kriging, in fact it
> makes prediction only inside each polygon cell. I explain me, it make
> one prediction per cell and the result of kriging is a mosaic of colors
> and not a map with gradients. I can't use contour() and I use :

> Plot(x,y,kc$predict) #x,y are coordinates for xx

> (If it is not clear I can send you the maps.) but thank you for your
> help, I will try something else.

The attachment is unfortunately in *.emf format, which can only be opened 
on Windows. A PNG graphic would have been enough, I think.

As far as your data are concerned, I think another possibility is to 
respect the fact that the rings were actually an incomplete grid. I guess 
that the prediction locations are the centroids of the grid cells. These 
suit the "sp" package class of SpatialPixels very well:

> shapefile <- read.shape("xdf.shp")
Shapefile Type: Polygon   # of Shapes: 4 
> x <- as.SpatialRings.Shapes(shapefile$Shapes, IDs=1:length(shapefile$Shapes))
> pts <- getSRSringsLabptSlots(x)
> pts
     [,1] [,2]
[1,]  1.5  0.5
[2,]  0.5  1.5
[3,]  2.5  1.5
[4,]  1.5  2.5

> check.locations(pts)
     [,1] [,2]
[1,]  1.5  0.5
[2,]  0.5  1.5
[3,]  2.5  1.5
[4,]  1.5  2.5

Make your prediction:

> kr.out <- krige.conv(..., locations=pts, ...)

Then use the "sp" magic to attach the points (which are regularly spaced) 
to the predictions:

res <- SpatialPixelsDataFrame(SpatialPoints(pts), list(z=kr.out$predict, 
z.se=krige.var))

should work. Try print(res) to confirm that it is a 
SpatialPixelsDataFrame. Then spplot(res, "z") will plot the predictions - 
the help page for spplot shows how to set the breaks and colour palette. 
This avoids the overlay step, and will work provided the points are spaced 
as a (rectangular) grid, even though many grid cells are missing.

If you want to share graphics, PNG - device png() - are much smaller. If 
you can put them on a web server, that's even better.

Hope this helps,

Roger



> Chloé 
> 
> 

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no