Plot Polygon / Rectangle with Leaflet R

Get Ready

I want to use polygons or rectangles to create a gridded heat map with leaflet in R. In this gridded heatmap, I want to show water stress as a function of precipitation and potential evaporation.

To do this, I created a dataset with this precipitation/potential evaporation index (from 2014-05-15 to 2014-09-15) in New York City (40.7127, -74.0059) using 10*15 and 1 arc minute grid cell. The visualization will show different levels of water stress on give area.

To get this data, use this gist with reproducible data frame.



In [1]:

    
source("https://gist.githubusercontent.com/yizhexu/a0f5b4806061005c2d93/raw/010788b1f0f55e87481d44f70d1a9429328f1c3d/data.R")

This gist gives two data sets: baseline_index and coordinates. Now print a few rows of data.



In [2]:

    
baseline_index[1:3]









    Out[2]:





	0.604648997169154
	0.604661799596444
	0.604674620472379



In [3]:

    
coors[1:3,]









    Out[3]:





latitude longitude

	124 40.62937 -74.12257
	248 40.64603 -74.12257
	372 40.6627 -74.12257



In [4]:

    
length(baseline_index)
dim(coors)









    Out[4]:




150






    Out[4]:





	150
	2

Plot Polygon

The first step to plot polygon, is to use known coordincates (coors) to generate coordinates that can be used to create a polygon. I wrote a function to do this:



In [5]:

    
library(sp)
create_poly <- function(latitude, longitude, x, size = 1) {
  lat_m <- latitude[x] - size/2/60
  lat_p <- latitude[x] + size/2/60
  lng_m <- longitude[x] - size/2/60
  lng_p <- longitude[x] + size/2/60

  p <- Polygon(matrix(c(lng_m, lat_m, lng_p, lat_m, lng_p, lat_p, lng_m, lat_p), ncol = 2, byrow = TRUE))
  Polygons(list(p),x)
}

This create_poly function will create a 1*1 polygon named x using a given set of coordinates as the center of the polygon.

Now use this function to create polygons for coors coordinates.



In [6]:

    
sps <- SpatialPolygons(lapply(1:150, function(i) {
  create_poly(coors$latitude, coors$longitude, i)}
))
proj4string(sps) = CRS("+init=epsg:3857")



In [7]:

    
summary(sps)









    Out[7]:





Object of class SpatialPolygons
Coordinates:
        min       max
x -74.13090 -73.88090
y  40.62103  40.80437
Is projected: TRUE 
proj4string : [+init=epsg:3857]

Now add the baseline_index data to each polygon



In [8]:

    
data <- data.frame(index = baseline_index, row.names = row.names(sps))
spdf <- SpatialPolygonsDataFrame(sps, data)
summary(spdf)









    Out[8]:





Object of class SpatialPolygonsDataFrame
Coordinates:
        min       max
x -74.13090 -73.88090
y  40.62103  40.80437
Is projected: TRUE 
proj4string : [+init=epsg:3857]
Data attributes:
     index       
 Min.   :0.4550  
 1st Qu.:0.5165  
 Median :0.5783  
 Mean   :0.5632  
 3rd Qu.:0.6047  
 Max.   :0.6769

Plot this on a map using leaflet package



In [9]:

    
library(leaflet)
pal <- colorQuantile(
  palette = "YlOrRd",
  domain = spdf@data$index, n = 5
)

leaf_poly <- leaflet() %>%
  addTiles() %>%
  addPolygons(data = spdf, fillOpacity = 0.5, smoothFactor = 0.5,
              color = ~pal(spdf@data$index),
              stroke = FALSE) %>%
  addLegend(pal = pal, values = spdf@data$index, opacity = 1)



In [10]:

    
library(htmlwidgets)
saveWidget(leaf_poly, file="leaf_poly.html")

Plot Rectangle

Another way I tried is to plot this collection of plots as rectangles. This method seems a lot cleaner because instead of assemble polygons before hand, just plot out rectangles - my polygons are rectangles anyway.



In [11]:

    
create_rectangle <- function(coors, x, size = 1) {
  lat_m <- coors$latitude[x] - size/2/60
  lat_p <- coors$latitude[x] + size/2/60
  lng_m <- coors$longitude[x] - size/2/60
  lng_p <- coors$longitude[x] + size/2/60
  c(lng_m, lat_p, lng_p, lat_m)
}

list <- do.call(rbind, lapply(1:150, function(i) {
  create_rectangle(coors, i)}
))

pal <- colorQuantile(
  palette = "YlOrRd",
  domain = spdf@data$index, n = 5
)

leaf_rect <- leaflet() %>%
  addTiles() %>%
  addRectangles(lng1 = list[, 1], lat1 = list[, 2],
                lng2 = list[, 3], lat2 = list[, 4],
                fillColor = pal(baseline_index), fillOpacity = 0.5, stroke = FALSE) %>%
  addLegend(pal = pal, values = spdf@data$index, opacity = 1)



In [12]:

    
saveWidget(leaf_rect, file="leaf_rect.html")



In [ ]: