2

I am looking for a way to represent density points using heatmaps in R for printed media (e.g. using Rmarkdown to generate a PDF).

The map should have a basemap (OpenStreetMap or Stammen) and instead of presenting all points I want them to be presented as a density heatmap overlayed.

I know there are several ways to create heatmaps in R, like leaflet.extras's addHeatmap() for HTML maps. Unfortunately, they can't be exported to PDF (or other formats different to html).

ggmap, on the contrary, is aimed for printed maps, but unfortunately, it has two main problems: 1) it requires google API key (I do not want to do do it, neither I want to use google services) 2) it is quite cumbersome, as it requires knowing bounding box beforehand and, since I have to do plenty of maps, I need an automated way to get that bounding box (ideally without requiring other packages, such as osmdata's getbb()).

tmap would be a great solution, as it allows re-using same code for printed and interactive maps, but I haven't found any documentation about generating heatmaps.

Is there any way to achieve what I am looking for?

ccamara
  • 892
  • 7
  • 17
  • 2
    What do you understand as a "heatmap"? To me its a kernel-smoothed density of points. You can compute that as a raster using several techniques, then all you need to do is plot a raster using tmap. Look at how to do point density estimation. – Spacedman Dec 05 '20 at 12:21

1 Answers1

4

You can use ggplot to create a kernel density map. Have a look at an old question on kernel density maps for some discussion on what 'heatmap' means in several contexts.

# Setup ----
library(ggplot2)

Reproducible data ----


obs <- data.frame(x = 35 + (rnorm(100)/5), y = -2.5 + (rnorm(100)/5),
                  count = rbinom(100, 10, 0.4))
obs <- rbind(obs, 
             data.frame(x = 34 + (rnorm(100)/5), y = -3 + (rnorm(100)/5),
                        count = rbinom(100, 10, 0.4))
             )


Plot ----


ggplot(obs, aes(x = x, y = y)) +
  stat_density_2d(aes(fill = ..density..), geom = "raster", contour = FALSE) +
  scale_fill_distiller(palette = 'OrRd', direction = 1) +
  scale_x_continuous(expand = c(0, 0)) +
  scale_y_continuous(expand = c(0, 0)) +
  theme(legend.position = 'none') +
  geom_point(cex = 0.5)

enter image description here

For a tiled basemap you can use the stamen maps:

library(ggmap)
osm <- get_stamenmap(bbox = c(min(obs$x), min(obs$y), max(obs$x), max(obs$y)))
g <- ggmap(osm)
g +   stat_density_2d(data = obs, aes(fill = ..density..), geom = "raster", contour = FALSE, , alpha = 0.5) + coord_cartesian() +
  scale_fill_distiller(palette = 'BuPu', direction = 1) +
  theme(legend.position = 'none')

enter image description here

Simbamangu
  • 14,773
  • 6
  • 59
  • 93
  • thank you for the solution and disambiguation! I love the idea of only using ggplot, with no further package. However, I will explore how to add a basemap to your code before accepting it. – ccamara Dec 08 '20 at 09:09
  • Could you edit your question and describe what else you need on the map? – Simbamangu Dec 08 '20 at 14:27
  • The additional code I added should solve it - though it's annoying to have to use the coord_cartesian to reproject – Simbamangu Dec 09 '20 at 14:36