Gráficos de densidad y puntos basados en mapas en R, usando ggmap

Ya he proporcionado ejemplos sobre cómo crear mapas de calor (es decir, diagramas de densidad) usando deckgl y Leaflet en R. En esta publicación, quiero proporcionar un ejemplo de cómo visualizar los atributos espaciales de un conjunto de datos usando el paquete ggmap en R. Empiezo por cargando los paquetes que querré usar para mi análisis.

#devtools::install_github("dkahle/ggmap", ref = "tidyup") 
# ya que actualmente ggmap no está en CRAN
library(ggmap)
library(ggplot2)
library(dplyr)
library(gridExtra)

Querré usar la base de datos predeterminada disponible en R para fines de demostración en este ejemplo de codificación. Por lo tanto, proporciono un vistazo de ese conjunto de datos mostrando sus entradas principales.

head(crime)

##                      time     date hour premise            offense  beat
## 82729 2010-01-01 07:00:00 1/1/2010    0     18A             murder 15E30
## 82730 2010-01-01 07:00:00 1/1/2010    0     13R            robbery 13D10
## 82731 2010-01-01 07:00:00 1/1/2010    0     20R aggravated assault 16E20
## 82732 2010-01-01 07:00:00 1/1/2010    0     20R aggravated assault  2A30
## 82733 2010-01-01 07:00:00 1/1/2010    0     20A aggravated assault 14D20
## 82734 2010-01-01 07:00:00 1/1/2010    0     20R           burglary 18F60
##           block    street type suffix number   month    day
## 82729 9600-9699   marlive   ln      -      1 january friday
## 82730 4700-4799 telephone   rd      -      1 january friday
## 82731 5000-5099  wickview   ln      -      1 january friday
## 82732 1000-1099   ashland   st      -      1 january friday
## 82733 8300-8399    canyon           -      1 january friday
## 82734 9300-9399     rowan   ln      -      1 january friday
##                       location           address       lon      lat
## 82729    apartment parking lot   9650 marlive ln -95.43739 29.67790
## 82730 road / street / sidewalk 4750 telephone rd -95.29888 29.69171
## 82731        residence / house  5050 wickview ln -95.45586 29.59922
## 82732        residence / house   1050 ashland st -95.40334 29.79024
## 82733                apartment       8350 canyon -95.37791 29.67063
## 82734        residence / house     9350 rowan ln -95.54830 29.70223

Notará: El conjunto de datos ya contiene coordenadas de longitud y latitud para todas las entradas de datos. Esta es la propiedad espacial de nuestro conjunto de datos.

A continuación, proporciono un ejemplo de cómo se pueden «extraer» mosaicos de mapas base del paquete ggmap. En el siguiente código recortado, construyo los mosaicos del mapa base para EE. UU.

# trazar un mapa base de ggmap
us <- c(left = -125, bottom = 25.75, right = -67, top = 49)
map <- get_stamenmap(us, zoom = 5, maptype = "toner-lite",legend="none")
plot(map)

Entonces la función «get_stamenmap» es del paquete ggmap.

Ahora estamos listos para crear un primer gráfico, basado en las propiedades espaciales de nuestro conjunto de datos. A continuación, muestro la distribución de las escenas del crimen de asesinato, según las coordenadas proporcionadas por el conjunto de datos del «crimen». La función «qmplot» es del paquete ggmap.

scatterplot_murder <- qmplot(x=lon,y=lat,data=filter(crime,offense=="murder"),legend="none",color=I("darkred"))
plot(scatterplot_murder)

A continuación, dibujaré un mapa de calor (es decir, una gráfica de densidad). Para esto, necesitaré especificar el parámetro «geom» en la función «qmplot» a «polígono». Además, necesito utilizar las funciones «stat_density_2d» y «scale_fill_gradient2». La estimación de densidad se basa en la estimación de densidad de kernel 2D. Se calcula mediante la función «stat_density_2d».

# crea otros tipos de gráficos con el paquete ggmap
densityplot_murder <- qmplot(x=lon, y=lat, 
                             data = filter(crime,offense=="murder"), 
                             geom = "blank",
                             maptype = "toner-background", 
                             darken = .7, 
                             legend = "topright") + stat_density_2d(aes(fill = ..level..), 
                  geom = "polygon", 
                  alpha = .5,
                  color = NA) + scale_fill_gradient2(low = "blue", 
                       mid = "green", 
                       high = "red")

plot(densityplot_murder)

En este ejemplo, la visualización aún no es perfecta y podría mejorarse aún más. Las formas de hacerlo serían, por ejemplo, ajustando el cálculo de estimación de densidad.

Asegúrese de revisar mis otras publicaciones sobre visualización de datos espaciales en R.