CODEWITHSUNDEEP
rpolygon
ianaj
ianaj

Reputation: 99

Aggregating Polygon Data in a Buffer in R

I have an excel file with gas stations and another shape file that has polygons (counties) and the corresponding population per county.

I've converted the excel file into a geospatial one using the st_as_sf command. The county shape file has polygons of all counties in the city and per county has a column which indicates the population in the county.

I was able to create a buffer around the gas stations and using map view able to plot the county, buffer and gas station together.

My question is how many person live within x meters (my buffer) from this gas station. In the picture the blue areas are the buffers around the gas station and the layer below shows the different counties coloured by "Population". Searching around has focused on counting points in a polygon but my query is a bit different.

I essentially want to pull from the population column in my county shape file and say 400 ppl live in the buffer of gas station 1 etc etc..for all gas stations.

picture of gas stations and buffers in blue

#create a buffer around fuel stations 1000 metres

petrol_buffer <- st_buffer(stations_geo,1000) %>% 
    rename(fuel_station_buffer=geometry) %>% 
     st_as_sf() %>% 
      mutate(buffer_id =1:29)

#create centroids of communities so that i can see which centroids fall within the fuel station buffer

   population_geo <- population_geo %>% 
    mutate(centroids= st_centroid(geometry))

#determine which centroids fall into fuel station buffers

in_buffer<-    st_intersection(population_geo$centroids,  petrol_buffer$fuel_station_buffer) %>% 
 as_tibble() %>% 
 rename(in_buffer_geo=geometry) %>% 
  st_as_sf(  crs = 4326) 

  in_buffer <- in_buffer %>% 
    st_join(population_geo, by='centroids')

#this is where im stuck i need to find a way to link this back to my fuel station file so i can know which fuel stations are being referenced by the 196 communities

Upvotes: 0

Views: 602

Answers (1)

Grzegorz Sapijaszko
Grzegorz Sapijaszko

Reputation: 3604

In example below I will use tigris package to get county geometries (instead of shapefile) and create few buffers (as gas station). Please note, that population density is an example, not real one.

library(tigris)
library(sf)
LA <- counties(state = "LA", cb = TRUE)
LA <- LA |>
  mutate(area = as.double(st_area(geometry)/10^6)) |>
  mutate(density = 200000/area) |>
  st_as_sf()
plot(LA["density"])

Let's create sample buffers -- in our case just two, for simplicity:

buffers <- st_sample(LA, 2) |>
  st_buffer(dist = 50000) |>
  st_as_sf() |>
  mutate(ID = 1:2)
  
plot(buffers$x, lwd = 3, add = TRUE)

Now, we have to find those counties (and their parts) which intersects with our buffers. The easiest way is just call st_intersection() function which in this case works like filter and returns the list of geometries, then we have to calculate the area/number of people and summarize it by our buffer ID.

LA |>
 st_intersection(buffers) |>
  mutate(aa = as.double(st_area(geometry)/10^6)) |>
  mutate(how_many_peple = aa * density) |>
  group_by(ID) |>
  summarise(sum(how_many_peple))

Simple feature collection with 2 features and 2 fields
[...]
# A tibble: 2 × 3
     ID `sum(how_many_peple)`                                                                                     geometry
  <int>                 <dbl>                                                                                <POLYGON [°]>
1     1               823361. ((-91.95035 30.28296, -91.94902 30.28298, -91.94413 30.28305, -91.94413 30.28882, -91.939...
2     2               292827. ((-93.29358 30.05181, -93.29358 30.05433, -93.28855 30.05446, -93.28855 30.05734, -93.286...

Please note area is in km^2 and density in number of people/km^2.

Upvotes: 0

Related Questions