[R-sig-Geo] kriging inside borders

[Roger Bivand]

On Fri, 20 May 2005, Chloe ARCHINARD wrote:

> My exact commands are :
> library(maptools)
> shapefile=read.shape('my.shp')
> x=Map2poly(shapefile)
> 
> and the length is :
> length(mappolys)
> [1] 615

OK. So what I think you have is a grid of cells, possibly with cells 
missing around the edges and maybe in the middle too. The splancs polygon 
cannot have holes, by the way. I think that what you could do is to choose 
the locations= values you need within the polygon grid first, and not use 
borders=. I've tried to replicate your situation, and have made a 
shapefile with rectangular grid cells with uneven edges and a hole. Then I 
did:

library(sp)
library(maptools)
shapefile <- read.shape("<your file>.shp")
x <- as.SpatialRings.Shapes(shapefile$Shapes, IDs=1:length(shapefile$Shapes))
xp <- SpatialPoints(cbind(runif(50, 0, 3), runif(50, 0, 3)))

# here you do xp <- SpatialPoints(<your location points matrix>)

res <- overlay(x, xp)

# overlay() sets the value of res to the ID of the polygon the point falls 
# in, or to NA if the point does not fall in any polygon.

plot(x, col="grey")
points(xp[is.na(res)], col="red")
points(xp[!is.na(res)], col="blue")

and you can choose the points for the prediction locations as:

<my prediction locations> <- xp[!is.na(res)]

Please check on the plot to see if overlay divides up the points correctly 
first, just to be sure. If this works, it is will be a nice example of how 
the coordination of spatial data functions in the "sp" package can make 
things easier (if not, any reports on difficulties with "sp" are always 
welcome).

Roger

>  But i'm not sure that it's the good polylist object i need because it's
> a polygon grid in fact.

> Thank you for your help
> Chloé
> -----Message d'origine-----
> De : Roger Bivand [mailto:Roger.Bivand at nhh.no] 
> Envoyé : jeudi 19 mai 2005 19:07
> À : Chloe ARCHINARD
> Cc : r-sig-geo at stat.math.ethz.ch
> Objet : Re: [R-sig-Geo] kriging inside borders
> 
> On Thu, 19 May 2005, Chloe ARCHINARD wrote:
> 
> 
> > Hello, I'm using geoR package to compute kriging with the function
> 
> > krige.conv() and I used maptools (read.shape () + Map2poly () )  and
> 
> > splancs (borders) packages to include a object of class polylist for the
> 
> > attribute borders.
> 
> 
> >  But I've got this message of error when 'borders' is in krige.conv () :  
> 
> > krige.conv: there are no prediction locations inside the borders.  My
> 
> > polylist is not good? The coordinates seam to be OK when I draw a map
> 
> > with this polylist object. And if I use 'borders' in the contour
> 
> > function : message error is :
> 
> 
> > Error in as.data.frame.default(borders) : can't coerce class "polylist"
> 
> > into a data.frame
> 
> 
> Firstly, it does help if the exact commands are given, rather than your 
> 
> choice of them. I think the main issue here is that a polylist object is a 
> 
> list of matrices with attributes, not a matrix. If your polylist object 
> 
> is:
> 
> 
> x <- Map2poly(read.shape("my.shp"))
> 
> 
> what is length(x)? If it is 1, you may be able to say x[[1]] to extract 
> 
> the first matrix in the list, but even then it may not work. If you can 
> 
> show us the commands you are giving, and the length of your polylist 
> 
> object, it should be possible to find something that works. 
> 
> 
> 
> > If someone could help me, thanks a lot!
> 
> 
> > Chloé Archinard
> 
> > 
> 
> > 
> 
> 
> -- 
> 
> 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
> 
> 
> 
> -- 
> passerelle antivirus du campus CNRS de Montpellier
> --
> 
> 
> 
> 

-- 
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