0

I have 2 shapefiles: data (spatial point dataframe) and polys (polygons dataframe). I want to do an overlap but is seems that it does not work.

Here are data and polys:

> data
class       : SpatialPointsDataFrame 
features    : 12527 
extent      : 10.20075, 20.6108, 54.08669, 57.75905  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs 
variables   : 3
names       :    timestamp_pretty,     timestamp,     imo 
min values  : 01/04/2006 00:00:55, 1143849655232, 9048392 
max values  : 30/04/2006 23:59:36, 1146441576823, 9191541 
> polys
class       : SpatialPolygonsDataFrame 
features    : 436375 
extent      : 4210000, 5441000, 3395000, 4813000  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=laea +lat_0=52 +lon_0=10 +x_0=4321000 +y_0=3210000 +ellps=GRS80 +units=m +no_defs 
variables   : 2
names       :    Id, Count 
min values  :     0,     0 
max values  : 99999,     9

to manage the overlapping, I use

proj4string(data) <- proj4string(polys)                            # to confirm the same reference system
inside <- !is.na(over(data, as(polys, "SpatialPolygons")))    # overlapping shape file and data

and then mean(inside) to check the average of points in polys.

But there is nothing happening, mean is always 0. I used this many times before and it always work, I guess it does not work not because the extent of the 2 sph files are different.

Is there a way to edit this?

  • It is not enough to set the CRS via proj4string(data) <- proj4string(polys). You need to transform one of the data to ensure identical CRS. data <- spTransform(data, CRS(proj4string(polys))) should do. – TimSalabim Jan 16 '18 at 9:09
  • @TimSalabim But following the description of data and polys, the CRS is already the same one, no? – Floni Jan 16 '18 at 9:27
  • 1
    Looking at the coordinates we see that even though the CRS is the same, the coordinates of data seem to be in longitude & latitude. Hence, none of the points are located inside any of the polygons. is the data excerpt you show above before or after you set the proj4string? If this is the raw data, then you can try to set proj4string(data) <- CRS("+init=epsg:4326") before the transformation. If the raw data are ok, then you need to do the transformation INSTEAD of setting the CRS via proj4string. – TimSalabim Jan 16 '18 at 9:38
  • @ TimSalabim, the data excerpt is after ´proj4string(data) <- proj4string(polys)´. I will continue to try... – Floni Jan 16 '18 at 11:16
  • Ok, then just replace proj4string(data) <- proj4string(polys) with data <- spTransform(data, CRS(proj4string(polys))) and things should work. – TimSalabim Jan 16 '18 at 12:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.