You can follow the next steps to achieve your goal:
- Transform your data as a
data.frame
object
- Convert the
data.frame
to SpatialPointsDataFrame
(SPDF) object
- Add a Coordinate Reference System (CRS) to SPDF objects
- Build a circular buffer around the points of the second dataset with the given radius
- Query which points of the first dataset are inside the buffer polygon
- Subset and sum of data by year
Try the commented R
code below.
- Goal: get the sum of the "newusdco" for the points in the corresponding area of the second dataset in the same years.
Step 1:
### Load libraries
library('sp') # for handling spatial objects
library('rgeos') # for buffers geometry
library('mapview') # for interactive map viewing in R (so nice!)
### Example data
# Data from question
dataset1 <- data.frame("eventid" = c("112130.01", "112087.01", "115838.01", "115839.01", "115839.02"),
"year" = c(1989, 1989, 1989, 1989, 1989),
"lat" = c(-8.838333, -15.333333, -8.833333, -15.332072, -16.510087),
"long" = c(13.234444, 15, 14.5, 12.660025, 15.435258),
"newusdco" = c(201289, 1143688, 1276358, 230453.667, 230453.667))
# Data from question
dataset2 <- data.frame("id" = c("70", "70", "70"),
"year" = c(1989, 1990, 1991),
"latitude" = c(12.19, 12.10, 11.93),
"longitude" = c(38.78, 39.03, 38.04),
"radius" = c(400, 300, 400),
"conflictarea" = c(455918, 280999, 444130),
"conflictsize" = c(530, 530, 530),
"conflictterritory" = c("Ethiopia", "Ethiopia", "Ethiopia"))
# Data not from question: created data with the aim of give an example to match points inside buffer polygons of dataset2
dataset1b <- data.frame("eventid" = c("112130.01", "112087.01", "115838.01", "115839.01", "115839.02", "112130.01", "115839.01"),
"year" = c(1989, 1989, 1989, 1989, 1989, 1990, 1990),
"lat" = c(12.191, 12.101, 11.931, 12.189, 12.990, 12.1911, 12.1891),
"long" = c(38.781, 39.031, 38.041, 38.779, 39.029, 38.781, 38.7791),
"newusdco" = c(201289, 1143688, 1276358, 230453.667, 230453.667, 201289, 230453.667))
Step 2:
### Convert to SpatialPointsDataFrame (SPDF) object
# Establish coordinates
coordinates(dataset1) <- c("long", "lat")
coordinates(dataset2) <- c("longitude", "latitude")
coordinates(dataset1b) <- c("long", "lat")
Step 3:
# Add Coordinate Reference System to SPDF objects
wgs84 <- "+proj=longlat +datum=WGS84 +no_defs +ellps=WGS84 +towgs84=0,0,0"
proj4string(dataset1) <- CRS(wgs84)
proj4string(dataset2) <- CRS(wgs84)
proj4string(dataset1b) <- CRS(wgs84)
### Plot data
spplot(dataset1, zcol = "newusdco", main = "dataset1 plot")
spplot(dataset2, zcol = "radius", main = "dataset2 plot")
# Interactive map
mapview(dataset1, color = "#CC0000") + mapview(dataset2, color = "#73D216") + mapview(dataset1b, color = "#75507B")

- Answer to the first problem "aggregate the data to the lat/long+radius data"
Step 4:
# Find UTM projection for SpatialObject
# Reference: https://stackoverflow.com/questions/9186496/determining-utm-zone-to-convert-from-longitude-latitude
# Also check: https://epsg.io/
long2UTM <- function(long) {
(floor((long + 180)/6) %% 60) + 1
}
long2UTM(coordinates(dataset2)[,"longitude"]) # UTM zone 37
# UTM
projUTM2 <- "+proj=utm +zone=37 +north +datum=WGS84 +units=m +no_defs"
dataset2UTM <- spTransform(dataset2, CRS(projUTM2))
# Make polygons from dataset2 circle buffer of radius = radius (in meters)
dataset2Buffer <- gBuffer(spgeom = dataset2UTM, byid = TRUE, width = dataset2UTM$radius)
# plot
mapview(dataset1, color = "#CC0000") + mapview(dataset2, color = "#73D216") + mapview(dataset1b, color = "#75507B") + mapview(dataset2Buffer, color = "#8AE234")

Step 5:
# Spatial Query: points in polygon
# Check which points of dataset1 are over dataset2
dataset2Buffer <- spTransform(dataset2Buffer, CRS(wgs84)) # Transform to WGS84
pointsInPolygons1 <- sp::over(x = dataset1, y = dataset2Buffer, returnList = TRUE) # Empty: no points over buffer
pointsInPolygons1b <- sp::over(x = dataset2Buffer, y = dataset1b, returnList = TRUE) # 6 points over buffer
$`1`
eventid year newusdco
1 112130.01 1989 201289.0
4 115839.01 1989 230453.7
6 112130.01 1990 201289.0
7 115839.01 1990 230453.7
$`2`
eventid year newusdco
2 112087.01 1989 1143688
$`3`
eventid year newusdco
3 115838.01 1989 1276358
- Answer to the second problem: "get the sum of 'newusdco' in the conflict site in the same years"
Step 6:
results <- lapply(X = pointsInPolygons1b, FUN = function(x) aggregate(newusdco ~ year, FUN = sum, data = x)
)
$`1`
year newusdco
1 1989 431742.7
2 1990 431742.7
$`2`
year newusdco
1 1989 1143688
$`3`
year newusdco
1 1989 1276358