Skip to content

Extract values with over()

6 messages · ASANTOS, Edzer Pebesma, Robert J. Hijmans +1 more

Original
ASANTOS
# 
Dear colleges,

      I try to extract values with over () function , but not working 
with my script, when I make:

require(sp)
#Create random points
x <- runif(n=500,min=0, max=20)
y <- runif(n=500,min=0, max=20)
xy<-cbind(x,y)
plot(x=x,y=y,axes=F, xlab="", ylab="",
      xlim=c(-1,20),ylim=c(0,21))
pt0 <- c(0.5,0.5)##First point of grid
##Create subarea
limx<-c(2.5,2.5,18,18,2.5)
limy<-c(2.5,18,18,2.5,2.5)
lim=cbind(limx,limy)
lim<-as.matrix(lim)
##Convert in spatial objects of sp package
#Random points
pontos <- SpatialPoints(cbind(x,y))
pontos<- SpatialPoints(list(x,y))
pontos<- SpatialPoints(data.frame(x,y))
#Create the grid
grid0= GridTopology(c(0.5,0.5), c(1,1), c(20,20))
grid<-SpatialGrid(grid = grid0)
##Subarea
limite=SpatialPolygons(list(Polygons(list(Polygon(lim)),"lim")))
##For view
plot(limite)#Subarea
points(pontos,col='red')#Random points
fullgrid(grid)=FALSE
clip<-overlay(grid,limite)##Select the grid inside the subarea
sel.grid <-grid[!is.na(clip),]##Removings Na's that represents the grid 
outside the subarea
points(sel.grid,col='blue')# New Grid
over(pontos,sel.grid)

was expected that the over() give the frequency of points inside each 
parcel, but I don't know what is the problem,
any suggestion,
Thanks
Edzer Pebesma
#
On 11/21/2011 08:47 PM, ASANTOS wrote:
it only extracts values when there are values to extract, please read
the documentation.
I wonder where you got this expectation from. If the documentation gave
you this expectation, please let me know.

table() can compute frequencies from indices.
Paulo Justiniano Ribeiro Jr
# 
Dear Edzer and Alexandre

wouldn't the following work as requested?

fullgrid(grid) <- FALSE
class(grid)  ## SpatialPixels

ap3 <- as.SpatialPolygons.SpatialPixels(grid)
table(over(pontos, ap3))



Paulo Justiniano Ribeiro Jr
LEG (Laboratorio de Estatistica e Geoinformacao)
Universidade Federal do Parana
Caixa Postal 19.081
CEP 81.531-990
Curitiba, PR  -  Brasil
Tel: (+55) 41 3361 3573
VOIP: (+55) (41) (3361 3600) 1053 1066
Fax: (+55) 41 3361 3141
e-mail: paulojus AT  ufpr  br
http://www.leg.ufpr.br/~paulojus
On Mon, 21 Nov 2011, Edzer Pebesma wrote:

            

      
      
      
    

        

      
      
    

  
  


  


      

    

    

      
      

        


  

        

      Edzer Pebesma
    

    

      
      #
    

  


  

      
On 11/21/2011 09:24 PM, Paulo Justiniano Ribeiro Jr wrote:

            

      
      
      
    

        
Paulo, that is a step in the right direction, as

            

      
      
      
    

        
1 2 3 4 5
2 1 3 1 2

gives frequencies for each cell, but this would also be obtained by the
grid itself:

            

      
      
      
    

        
1 2 3 4 5
2 1 3 1 2

However, as

            

      
      
      
    

        
[1] 298

            

      
      
      
    

        
[1] 400

shows we have empty grid cells. but this gives gaps (implicit 0 values)
for empty cells in grid. The simpler alternative is to do the over with
returnList=TRUE, which returns ALL points in each (grid cell, now
polygon) in ap3:

            

      
      
      
    

        
[[1]]
[1] 416 459

[[2]]
[1] 385

[[3]]
[1]  17 132 351

[[4]]
[1] 387

[[5]]
[1] 311 322

taking the length of each vector gives the frequencies (INCLUDING the
zero frequencies that table does NOT give you):

n = sapply(over(ap3, pontos, returnList=TRUE), length)
x = SpatialPixelsDataFrame(grid, data.frame(n=n))
spplot(x, sp.layout=list("sp.points", pontos), at=0:6-.5)

Without doubt, there are other ways of doing this, for instance with
package raster.

Thanks for helping. It seems to be hard to make complex things simple,
and I'd be happy to hear suggestions for making this even simpler.

            

      
      
      
    

        

      
      
    

        

  
    

      
      
    

  


  


      

    

    
1 day later

    

      
      

        
      

    

    

      
      

        


  

        

      Paulo Justiniano Ribeiro Jr
    

    

      
      #
    

  


  

      
Edzer
thanks for pointing this

My complete solution goes a bit different without using returnList
and taking advantege of the named vector returned by over()

##
require(sp)
set.seed(12)
xy <- cbind(runif(n=500,min=0, max=20),
             runif(n=500,min=0, max=20))
## points
pontos <- SpatialPoints(xy)
## grid of "quadrats"
## i. as SpatialGrid
grid <- SpatialGrid(grid = GridTopology(c(0.5,0.5), c(1,1), c(20,20)))
## ii. as SpatialPixels
gridPix <- grid
fullgrid(gridPix) <- FALSE
## iii. as SpatialPolygons
gridPol <- as.SpatialPolygons.SpatialPixels(gridPix)

PtGr <- table(over(pontos, grid));length(PtGr)
POS <- as.numeric(names(PtGr))

N <- numeric(length(grid))
N[POS] <- PtGr
N

plot(gridPol)
plot(pontos, add=T, pch=19, cex=0.5)
plot(grid, pch=as.character(N), add=T, cex=0.6, col=2)



Paulo Justiniano Ribeiro Jr
LEG (Laboratorio de Estatistica e Geoinformacao)
Universidade Federal do Parana
Caixa Postal 19.081
CEP 81.531-990
Curitiba, PR  -  Brasil
Tel: (+55) 41 3361 3573
VOIP: (+55) (41) (3361 3600) 1053 1066
Fax: (+55) 41 3361 3141
e-mail: paulojus AT  ufpr  br
http://www.leg.ufpr.br/~paulojus

        
On Tue, 22 Nov 2011, Edzer Pebesma wrote: