1

In R, I am trying to create a choropleth map. I have built a database of businesses, some are part of chains (e.g. McDonalds) and others are independent. I want to calculate how many businesses are within 30km of each point on the map, but treat the different locations of chains as a single business.

For example, if you a point is:

  • 5km from a McDonalds,
  • 10km from Taco Bell
  • 15km from Chick-Fil-A
  • 20km from KFC
  • 25km from McDonalds
  • 35km from Five Guys

The colour will show that there are 4 fast food outlets within 30km.

I am happy to use any R package but I am mostly familiar with tmaps and ggplot2 maps.

At this stage the best approach I can think of is to create polygons for each chain and stack them as transparent layers of the same colour. I don't think this would be very efficient and wouldn't create a very nice looking choropleth.

The other answers I could find were either counting points (e.g https://gis.stackexchange.com/questions/229066/counting-how-many-times-a-point-is-inside-a-set-of-intersecting-polygons-in-r) or for GIS software.

EDIT: I have managed to create a 30km radius from every location of every chain (using rgeos gIntersection). I now have a series of polygons.

To solve my question the additional thing I need to do is create polygons for where:

  • Only one polygon covers the area,
  • Two polygons covers the area,
  • etc.

To try to visual is this I used the answer from https://gis.stackexchange.com/questions/229066/counting-how-many-times-a-point-is-inside-a-set-of-intersecting-polygons-in-r

enter image description here

In the linked question they are trying to count how many polygons cover the numbered points (the image on the right). What I am trying to do is to create the image on the left, where there are polygons of no overlap (1), two overlapping polygons (2) and so on.

Novo88
  • 109
  • 1
  • 7
  • Just a technical question - Are the points you want to draw polygons around different than the businesses' locations? – Justin Landis Dec 07 '20 at 20:39
  • Hi Justin, I basically want a map of New Zealand where it is coloured in based on how many fast food outlets within 30km of each point. I think practically that means building radiuses around each business within each chain (I can do that). Then changing the colour every time there is an overlap of different chains - this I haven't figured out. I thought there might be a neater way. – Novo88 Dec 07 '20 at 22:05
  • 1
    From the comment it sounds like you want to create a grid and tabulate the far food restaurants in each cell. Raster methods would be more appropriate here. – Brian Fisher Dec 08 '20 at 05:17
  • Hi Brian, Can you point me in the direction of where I could find a tutorial or example of how to use Raster methods? – Novo88 Dec 08 '20 at 06:22
  • 1
    Geocomputing in R by Robin Lovelace is a good reference https://geocompr.robinlovelace.net/spatial-class.html#raster-data – Brian Fisher Dec 08 '20 at 19:35
  • I posted an example below. It would be easier to help if you include some sample data so we can see what you are working with. – Brian Fisher Dec 08 '20 at 20:00

1 Answers1

0

I think what you are trying to accomplish would be best approached using a raster approach rather than a chloropleth. To make a chorlopleth, you define a set of (generally irregular) polygons, summarize something within each polygon, then color the polygons based on the attributes. This would be a good approach if you wanted to say how many fast food resteraunts are within each state or county, or how many fast food joints per capita by state.
From your description, however, you are looking for how many fast food joints within a set radius for all points. This is more of a raster question, since you can represent your data on a regular grid.
The raster package is a good start for working with raster data and works well with the sf package.

You need to determine what density you need to accomplish your goal, then use this to determine the resolution of your raster. Once you've got that you can use raster::rasterize() to summarize your (I'm assuming) point data.

I'm assuming you have an object that has the locations of each restaurant, I'll call this object "points".

library(raster)
library(sf)

# create raster template with 30km resolution (I'm assuming your projection is in meters)
raster_template = raster((extent(points), 
                         resolution = 30000,
                         crs = st_crs(points)
                         )
# rasterize your point data
r = rasterize(points, raster_template, fun = "count")

This should create a grid where each cell has the number of points within each 30km cell. You should then be able to plot the raster, but may want to either clip or mask it to just show parts that are within New Zealand

Brian Fisher
  • 1,305
  • 7
  • 17