# R : How to build heatmap with the leaflet package

I read a post about interactive maps with R using the `leaflet` package.

In this article, the author createD a heat map like this:

``````X=cbind(lng,lat)
kde2d <- bkde2D(X, bandwidth=c(bw.ucv(X[,1]),bw.ucv(X[,2])))

x=kde2d\$x1
y=kde2d\$x2
z=kde2d\$fhat
CL=contourLines(x , y , z)

m = leaflet() %>% addTiles()
m %>% addPolygons(CL[[5]]\$x,CL[[5]]\$y,fillColor = "red", stroke = FALSE)
``````

I am not familiar with the `bkde2D` function, so I'm wondering if this code could be generalized to any shapefiles?

What if each node has a specific weight, that we would like to represent on the heat map?

Are there other ways to create a heat map with `leaflet` map in R ?

• bke2d lets you do 2d binning (kernel density estimation) for a set of points (so lng/lat pairs work well). the ks package supports kernel smoothing for data from 1- to 6-dimensions. the akima package can do interpolation (useful when you need a regular grid). it might be worth reading up on the spatial task view for this before attempting to produce something that may not represent the data properly. – hrbrmstr Nov 3 '15 at 18:05
• ok, thanks for the link, I will definitely look this. Actually the bke2d function is not working that well with my data as it works in the example, and I can t figure why. – Felipe Nov 4 '15 at 17:02

Here's my approach for making a more generalized heat map in Leaflet using R. This approach uses `contourLines`, like the previously mentioned blog post, but I use `lapply` to iterate over all the results and convert them to general polygons. In the previous example it's up to the user to individually plot each polygon, so I would call this "more generalized" (at least this is the generalization I wanted when I read the blog post!).

``````## INITIALIZE
library("leaflet")
library("data.table")
library("sp")
library("rgdal")
# library("maptools")
library("KernSmooth")

infile <- "mvthefts.csv"

## Also, clean up variable names, and convert dates
if(!file.exists(infile)){
}
setnames(dat, tolower(colnames(dat)))
setnames(dat, gsub(" ", "_", colnames(dat)))
dat <- dat[!is.na(longitude)]
dat[ , date := as.IDate(date, "%m/%d/%Y")]

## MAKE CONTOUR LINES
## Note, bandwidth choice is based on MASS::bandwidth.nrd()
kde <- bkde2D(dat[ , list(longitude, latitude)],
bandwidth=c(.0045, .0068), gridsize = c(100,100))
CL <- contourLines(kde\$x1 , kde\$x2 , kde\$fhat)

## EXTRACT CONTOUR LINE LEVELS
LEVS <- as.factor(sapply(CL, `[[`, "level"))
NLEV <- length(levels(LEVS))

## CONVERT CONTOUR LINES TO POLYGONS
pgons <- lapply(1:length(CL), function(i)
Polygons(list(Polygon(cbind(CL[[i]]\$x, CL[[i]]\$y))), ID=i))
spgons = SpatialPolygons(pgons)

## Leaflet map with polygons
leaflet(spgons) %>% addTiles() %>%
addPolygons(color = heat.colors(NLEV, NULL)[LEVS])
``````

Here's what you'll have at this point:

``````## Leaflet map with points and polygons
## Note, this shows some problems with the KDE, in my opinion...
## For example there seems to be a hot spot at the intersection of Mayfield and
## Fillmore, but it's not getting picked up.  Maybe a smaller bw is a good idea?

leaflet(spgons) %>% addTiles() %>%
addPolygons(color = heat.colors(NLEV, NULL)[LEVS]) %>%
addCircles(lng = dat\$longitude, lat = dat\$latitude,
radius = .5, opacity = .2, col = "blue")
``````

And this is what the heat map with points would look like:

Here's an area that suggests to me that I need to tune some parameters or perhaps use a different kernel:

``````## Leaflet map with polygons, using Spatial Data Frame
## Initially I thought that the data frame structure was necessary
## This seems to give the same results, but maybe there are some
## advantages to using the data.frame, e.g. for adding more columns
spgonsdf = SpatialPolygonsDataFrame(Sr = spgons,
data = data.frame(level = LEVS),
match.ID = TRUE)
leaflet() %>% addTiles() %>%
color = heat.colors(NLEV, NULL)[spgonsdf@data\$level])
``````
• Scoured the interwebs trying to figure this out and this was by far the best example I found. Plugged it in to my code and it "just worked." Awesome. Thank you! – Jeff Allen Jan 31 '18 at 21:46
• Thanks! I've actually created a repo with several other map examples that might be useful for others github.com/geneorama/wnv_map_demo – geneorama Feb 9 '18 at 17:54
• Thanks for this mini-tutorial. How did you select the `bandwidth` in `bkde2d()` ? – the_darkside Feb 21 '18 at 0:14
• @the_darkside great question. In reality I fiddle with it until I get something I like, I originally developed this map specifically to examine the bandwidth assumptions. In this case I actually used `MASS::bandwidth.nrd(dat\$latitude)` and `MASS::bandwidth.nrd(dat\$longitude)` as the starting points. See `?MASS::kde2d` documentation which links to `bandwith.nrd`. Also see `?KernSmooth::dpik` if you're interested for another approach. – geneorama Feb 21 '18 at 16:02
• if `gridsize = c(100,100)` does that mean there are a total of 10,000 cells? – the_darkside Feb 21 '18 at 17:25

An easy way of creating Leaflet heat maps in R is using the Leaflet.heat plugin. An excellent guide on how to use it can be found here. Hope you find it useful.

Building off of genorama's answer above, you can also convert the output of bkde2D into a raster rather than contour lines, using the fhat values as the raster cell values

``````library("leaflet")
library("data.table")
library("sp")
library("rgdal")
# library("maptools")
library("KernSmooth")
library("raster")

infile <- "mvthefts.csv"

## Also, clean up variable names, and convert dates
if(!file.exists(infile)){
}
setnames(dat, tolower(colnames(dat)))
setnames(dat, gsub(" ", "_", colnames(dat)))
dat <- dat[!is.na(longitude)]
dat[ , date := as.IDate(date, "%m/%d/%Y")]

## Create kernel density output
kde <- bkde2D(dat[ , list(longitude, latitude)],
bandwidth=c(.0045, .0068), gridsize = c(100,100))
# Create Raster from Kernel Density output
KernelDensityRaster <- raster(list(x=kde\$x1 ,y=kde\$x2 ,z = kde\$fhat))

#create pal function for coloring the raster
palRaster <- colorNumeric("Spectral", domain = KernelDensityRaster@data@values)

## Leaflet map with raster
leaflet() %>% addTiles() %>%
colors = palRaster,
opacity = .8) %>%
values = KernelDensityRaster@data@values,
title = "Kernel Density of Points")
``````

This is your output. Note that the low density values still show up as colored in the raster.

We can remove these low density cells with the following :

``````#set low density cells as NA so we can make them transparent with the colorNumeric function
KernelDensityRaster@data@values[which(KernelDensityRaster@data@values < 1)] <- NA

#create pal function for coloring the raster
palRaster <- colorNumeric("Spectral", domain = KernelDensityRaster@data@values, na.color = "transparent")

## Redraw the map
leaflet() %>% addTiles() %>%
colors = palRaster,
opacity = .8) %>%
values = KernelDensityRaster@data@values,
title = "Kernel Density of Points")

``````

Now any raster cell with a value of less than 1 is transparent.

If you want a binned raster, use the colorBin function rather than the colorNumeric function:

``````palRaster <- colorBin("Spectral", bins = 7, domain = KernelDensityRaster@data@values, na.color = "transparent")

## Leaflet map with raster
leaflet() %>% addTiles() %>%
``````gridsize = c(1000,1000)