[R-sig-Geo] Overlay, sample, and extract

Edzer Pebesma edzer.pebesma at uni-muenster.de
Mon Nov 12 20:26:11 CET 2012


try:

x = do.call(rbind, lapply(1:length(SpPDF), function(x)
spsample(SpPDF[x,], 10, "random")))
points(x)

On 11/12/2012 08:13 PM, dani jon wrote:
> Dear List,
> 
> I am interested in overlaying a regular polygon grid of 8 by 8 into a
> spatial surface. Essentially, dividing the spatial surface into 64 polygon
> grids and would like to take a sample of 10 points (pixels) with
> corresponding xy coordinate and values for each of the 64 grids.  Here is
> my attempt.
> 
> library(maptools)
> library(raster)
> library(sp)
> library(gstat)
> library(maps)
> 
> ## Creating spatial surface
> 
> xg = seq(1, 100,, 100)
> yg = seq(1, 100,, 100)
> xy = expand.grid(xg, yg)
> names(xy) <- c("x","y")
> g.dummy <- gstat(formula=z~1, locations=~x+y, dummy=T, beta=1,
> model=vgm(psill=0.025,model="
> Exp",range=15), nmax=20)
> yy <- predict(g.dummy, newdata=xy, nsim=1)
> gridded(yy) = ~x+y
> coor.yy <- coordinates(yy[1])
> ytrue <-yy[1]@data
> yyspdf<- SpatialPointsDataFrame(coords=coor.yy, data=ytrue)
> yyspdf[1:10,]
>  coordinates     sim1
> 1       (1, 1) 1.305260
> 2       (2, 1) 1.287697
> 3       (3, 1) 1.219094
> 4       (4, 1) 1.333913
> 5       (5, 1) 1.268508
> 6       (6, 1) 1.263417
> 7       (7, 1) 1.243237
> 8       (8, 1) 1.276870
> 9       (9, 1) 1.231528
> 10     (10, 1) 1.162322
> 
> ## creating regular polygon grid (8x8)
> 
> grd <- GridTopology(c(6.3,6.3), c(12.5,12.5), c(8,8))
> SpP_grd <-as.SpatialPolygons.GridTopology(grd)
> plot(SpP_grd)
> text(coordinates(SpP_grd), sapply(slot(SpP_grd, "polygons"), function(i)
> slot(i, "ID")), cex=0.5)
> ids <- sapply(slot(SpP_grd, 'polygons'), function(i) slot(i, 'ID'))
> trdata <- data.frame(A=rep(c(1,2,3,4,5,6,7,8), 8),
> B=rep(c(1,2,3,4,5,6,7,8), each=8),
> row.names=sapply(slot(SpP_grd, "polygons"), function(i) slot(i, "ID")))
> SpPDF <- SpatialPolygonsDataFrame(SpP_grd, trdata)
> SpPDF at data = cbind(SpPDF at data,ids)
> 
>> SpPDF at data[1:10,]
>     A B ids
> g1  1 1  g1
> g2  2 1  g2
> g3  3 1  g3
> g4  4 1  g4
> g5  5 1  g5
> g6  6 1  g6
> g7  7 1  g7
> g8  8 1  g8
> g9  1 2  g9
> g10 2 2 g10
> 
> ## Overlay
> image(yy[1])
> plot(SpPDF, border="magenta", add=TRUE)
> text(coordinates(SpP_grd), sapply(slot(SpP_grd, "polygons"), function(i)
> slot(i, "ID")), cex=0.5)
> 
> My question is how to overlay yy[1] and SpPDF and take 10 samples within
> each polygon grid (g1-g64) and extract the values and coordinates of
> yy[1]?  So I will have an output with a sample of 640
> 
>  x  y   sim1        ids
> 1  2   1.305260   g1
> 1  3   1.02567     g1
> .
> .
> 95 95 1.3487    g64
> 
> Edzer  suggested something likes this, but it does not exactly  640 samples
> (10*64) as I wanted.
> 
>  x = do.call(rbind, lapply(1:10, function(x) spsample(yy, cellsize = 12.5,
> type = "stratified")))
> 
> Thanks in advance
> Best regards
> Sam
> 
> 	[[alternative HTML version deleted]]
> 
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> 

-- 
Edzer Pebesma
Institute for Geoinformatics (ifgi), University of Münster
Weseler Straße 253, 48151 Münster, Germany. Phone: +49 251
8333081, Fax: +49 251 8339763  http://ifgi.uni-muenster.de
http://www.52north.org/geostatistics      e.pebesma at wwu.de



More information about the R-sig-Geo mailing list