0

I have a gridded map shapefile with grid identifier cells (gid).

dput(head(grid))
structure(list(gid = c("183253", "183254", "183255", "183256", 
"183258", "183259"), xcoord = c(6.25, 6.75, 7.25, 7.75, 8.75, 
9.25), ycoord = c(37.25, 37.25, 37.25, 37.25, 37.25, 37.25), 
    col = c("373", "374", "375", "376", "378", "379"), row = c("255", 
    "255", "255", "255", "255", "255"), gwno = c(615L, 615L, 
    615L, 615L, 616L, 616L), country = c("Algeria", "Algeria", 
    "Algeria", "Algeria", "Tunisia", "Tunisia"), km2 = c(3875.43674549, 
    3875.43674549, 3875.43674549, 3875.43674549, 3875.43674549, 
    3875.43674549), geometry = structure(list(structure(list(
        structure(c(6, 6, 6.5, 6.5, 6, 37, 37.5, 37.5, 37, 37
        ), .Dim = c(5L, 2L))), class = c("XY", "POLYGON", "sfg"
    )), structure(list(structure(c(6.5, 6.5, 7, 7, 6.5, 37, 37.5, 
    37.5, 37, 37), .Dim = c(5L, 2L))), class = c("XY", "POLYGON", 
    "sfg")), structure(list(structure(c(7, 7, 7.5, 7.5, 7, 37, 
    37.5, 37.5, 37, 37), .Dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(7.5, 7.5, 
    8, 8, 7.5, 37, 37.5, 37.5, 37, 37), .Dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(8.5, 8.5, 
    9, 9, 8.5, 37, 37.5, 37.5, 37, 37), .Dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg")), structure(list(structure(c(9, 9, 9.5, 
    9.5, 9, 37, 37.5, 37.5, 37, 37), .Dim = c(5L, 2L))), class = c("XY", 
    "POLYGON", "sfg"))), class = c("sfc_POLYGON", "sfc"), precision = 0, bbox = structure(c(xmin = 6, 
    ymin = 37, xmax = 9.5, ymax = 37.5), class = "bbox"), crs = structure(list(
        input = "WGS 84", wkt = "GEOGCRS[\"WGS 84\",\n    DATUM[\"World Geodetic System 1984\",\n        ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n            LENGTHUNIT[\"metre\",1]]],\n    PRIMEM[\"Greenwich\",0,\n        ANGLEUNIT[\"degree\",0.0174532925199433]],\n    CS[ellipsoidal,2],\n        AXIS[\"latitude\",north,\n            ORDER[1],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n        AXIS[\"longitude\",east,\n            ORDER[2],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n    ID[\"EPSG\",4326]]"), class = "crs"), n_empty = 0L)), sf_column = "geometry", agr = structure(c(gid = NA_integer_, 
xcoord = NA_integer_, ycoord = NA_integer_, col = NA_integer_, 
row = NA_integer_, gwno = NA_integer_, country = NA_integer_, 
km2 = NA_integer_), .Label = c("constant", "aggregate", "identity"
), class = "factor"), row.names = c(NA, 6L), class = c("sf", 
"data.frame"))

I have spatial points that represent the sum of fatalities during a period of time.

 dput(head(acled))
structure(list(OBJECTID = c(1, 2, 3, 4, 5, 6), data_id = c(8037001, 
6225870, 8037590, 8037835, 6270986, 8037919), year = c(2010, 
2010, 2009, 2007, 2007, 2007), event_type = c("Battles", "Battles", 
"Battles", "Battles", "Battles", "Battles"), country = c("Angola", 
"Angola", "Angola", "Angola", "Angola", "Angola"), latitude = c(-8.8383, 
-5.5758, -4.7666, -4.7666, -5.5758, -4.7666), longitude = c(13.2344, 
12.1871, 12.55, 12.55, 12.1871, 12.55), fatalities = c(2, 0, 
0, 1, 2, 1), geometry = structure(list(structure(c(13.2344, -8.8383
), class = c("XY", "POINT", "sfg")), structure(c(12.1871, -5.5758
), class = c("XY", "POINT", "sfg")), structure(c(12.55, -4.7666
), class = c("XY", "POINT", "sfg")), structure(c(12.55, -4.7666
), class = c("XY", "POINT", "sfg")), structure(c(12.1871, -5.5758
), class = c("XY", "POINT", "sfg")), structure(c(12.55, -4.7666
), class = c("XY", "POINT", "sfg"))), class = c("sfc_POINT", 
"sfc"), precision = 0, bbox = structure(c(xmin = 12.1871, ymin = -8.8383, 
xmax = 13.2344, ymax = -4.7666), class = "bbox"), crs = structure(list(
    input = "WGS 84", wkt = "GEOGCRS[\"WGS 84\",\n    DATUM[\"World Geodetic System 1984\",\n        ELLIPSOID[\"WGS 84\",6378137,298.257223563,\n            LENGTHUNIT[\"metre\",1]]],\n    PRIMEM[\"Greenwich\",0,\n        ANGLEUNIT[\"degree\",0.0174532925199433]],\n    CS[ellipsoidal,2],\n        AXIS[\"latitude\",north,\n            ORDER[1],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n        AXIS[\"longitude\",east,\n            ORDER[2],\n            ANGLEUNIT[\"degree\",0.0174532925199433]],\n    ID[\"EPSG\",4326]]"), class = "crs"), n_empty = 0L)), sf_column = "geometry", agr = structure(c(OBJECTID = NA_integer_, 
data_id = NA_integer_, year = NA_integer_, event_type = NA_integer_, 
country = NA_integer_, latitude = NA_integer_, longitude = NA_integer_, 
fatalities = NA_integer_), .Label = c("constant", "aggregate", 
"identity"), class = "factor"), row.names = c(NA, 6L), class = c("sf", 
"data.frame"))

Ideally, I would like a final spatialpolygonsdataframe that represents polygons summing the fatalities inside them for each year + the polygon that don't have any fatalities. It would be like the 'union' function between two similar shapefiles.

However, when I try do so with the function point.in.poly(acled, grid) from the spatialEco package I get a spatialpointsdataframe with only the points falling into the corresponding grid polygons. Same with the 'over' function.

Is there a way to get the opposite ? I tried over(grid, ACLED) but it does not give me a polygon per year as there is a total number of 10667 cells in the grid shapefile, in other words, it does not duplicate the polygons to let me distinguish the fatalities for the different years.

1
  • 1
    I've removed the "(or ArcGIS)" from the title because that would get the question closed (multiple questions)
    – Vince
    May 17, 2022 at 22:41

1 Answer 1

1

I found out how to do it.

Two ways:

  1. applying point.in.poly(grid, acled). So x= polygons, y=points.

  2. applying point.in.poly(ACLED, grid) in the right way then merging with the initial grid.

grid.acled <- sp::merge(grid, pts.poly[c("fatalities","event_type", "gid", "year")], by="gid", all.x= TRUE).

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.