R provides a rich environment for cartographic plotting. In the following, we want to have a look at the following plotting library that supports the visualization of geographical data as maps:

Maps with ggplot2

The powerful ggplot2 library is capable of producing quite nice-looking maps. In order to showcase these capacities, we will visualize the location of weather stations maintained by Deutscher Wetterdienst (DWD) (German Weather Service). The data was downloaded from the DWD data portal on April 21, 2017.

First, before plotting any data, we must prepare our data set. Hence, we start by importing the data set and assigning it a proper name. In this example, we use the read.csv2() function to load the data.

dwd <- read.csv2("https://userpage.fu-berlin.de/soga/data/raw-data/DWD.csv",
  stringsAsFactors = FALSE, sep = ","
)

The dwd data set consists of 599 rows, each of them representing a particular weather station in Germany, and 22 columns, each corresponding to a variable/feature related to that particular weather station. For the purpose of this tutorial, we are only interested in the variables LAT and LON, corresponding to the geographic coordinates of each particular weather station in the data set, as well as the record length in years (RECORD.LENGTH). Hence, we subset our data set based on these variables. Further, we make sure that we exclude all missing values.

## subset data
dwd <- dwd[, c("LAT", "LON", "RECORD_LENGTH")]
## drop missing values
dwd <- dwd[complete.cases(dwd), ]

dwd$LAT <- as.numeric(dwd$LAT)
dwd$LON <- as.numeric(dwd$LON)
nrow(dwd)
## [1] 599
head(dwd)
##       LAT     LON RECORD_LENGTH
## 1 47.8413  8.8493            55
## 2 50.7827  6.0941           160
## 3 52.9335  8.2370            45
## 4 48.2156  8.9784            30
## 5 48.6159 13.0506            64
## 6 52.4853  7.9126            55

For visualization we categorize the variable RECORD.LENGTH into five groups:

To achieve this task, we apply the mutate() function provided by the dplyr package. If you get an error message, type install.packages('dplyr') to install the package.

library(dplyr)
dwd <- dwd %>%
  mutate(RECORD_LENGTH_CATEGORY = cut(RECORD_LENGTH,
    breaks = c(-Inf, 10, 30, 60, 90, Inf),
    labels = c(
      "very short (<10)",
      "short (10-30)",
      "middle (30-60)",
      "long (60-90)",
      "very long (>90)"
    )
  ))

Now we are ready to plot a map to visualize the spatial distribution of the DWD weather stations. For this, we rely on the terra package and on the ggplot2 package as well as the geodata package.

library(terra)
library(ggplot2)
library(geodata)

If you get an error message that one or more packages are not found, use the command install.packages(c('terra', 'ggplot2', 'geodata')) to install the required packages.

To provide a geographical context for the location of our DWD weather stations within our map, we want to add the administrative boundaries of Germany on a federal-state scale. Therefore, we use the gadm() function provided by the geodata package to obtain the associated polygon files:

germany <- gadm(country = "Germany", level = 1, path = ".")

In the next step, we build a ggplot2() graph, layer by layer.

germany_sf <- sf::st_as_sf(germany)
DWD_sf <- sf::st_as_sf(dwd, coords=c('LON', 'LAT'), crs = crs(germany_sf))
ggplot() + 
  geom_sf(data = germany_sf, 
          colour = "grey10",
          fill = "#fff7bc") + 
  geom_sf(data = DWD_sf, 
          aes(col = RECORD_LENGTH_CATEGORY), 
          alpha = 0.5, size = 1.5) + 
  coord_sf() +
  theme_bw() +
  xlab("Longitude") +
  ylab("Latitude") +
  labs(color = "Record length") +
  ggtitle("DWD Weather Stations")

mapview package

The mapview package provides functions to very quickly and conveniently create interactive visualizations of spatial data. It allows quick interactive plotting to examine and visually investigate both aspects of spatial data, the geometries and their attributes. Visit the project’s website for more information and a quick tour of the functionality of this awesome package.

Our goal is to localize GeoCampus Lankwitz, Berlin, Germany and plot the location onto an interactive map.

First, we have to determine the geographical location of GeoCampus Lankwitz. Therefore, we make use of the geocode_OSM() function in the tmaptools package. This function allows us to geocode a location (find its latitude and longitude) using OpenStreetMap Nominatim. If you have not already, install the tmaptools package using install.packages("tmaptools").

Let us retrieve the coordinates for the GeoCampus Lankwitz using the geocode_OSM() function as a data frame. We then extract only the lat and lon coordinates from the result.

library(tmaptools)
## Warning: Paket 'tmaptools' wurde unter R Version 4.3.3 erstellt
geocampus <- geocode_OSM(q = "Malteserstrasse 74-100, Berlin", as.data.frame = TRUE)
geocampus_latlon <- geocampus[, c("lat", "lon")]
geocampus_latlon
##        lat      lon
## 1 52.42625 13.35803

Nice!

In the next step, we transform the geocampus date frame object into a simple feature object. Therefore, we apply the st_as_sf() function from the sf package. In addition, we provide that object with a spatial coordinate reference system. In our example we use the EPSG identifier \(4326\), which corresponds to the World Geodetic System (WGS84).

library(sf)
geocampus_sf <- st_as_sf(
  x = geocampus_latlon, coords = c("lon", "lat"),
  crs = 4326
)
geocampus_sf
## Simple feature collection with 1 feature and 0 fields
## Geometry type: POINT
## Dimension:     XY
## Bounding box:  xmin: 13.35803 ymin: 52.42625 xmax: 13.35803 ymax: 52.42625
## Geodetic CRS:  WGS 84
##                    geometry
## 1 POINT (13.35803 52.42625)

Once we have a spatial feature, plotting it using the mapview() function from the same-titled mapview package is very easy. Note that the mapview() function has many additional arguments (type help(mapview) into the R console), allowing for customization. Additionally, we can use the leafpop package to embed tables, images or graphs in pop-ups of our interactive map (Feel free to scroll over the blue marker in the map and click on it!).

library(mapview)
img <- "http://www.osa.fu-berlin.de/geographie_bsc/_medien/bilder_geocampus/bild_010/010_546.jpeg"

mapview(geocampus_sf,
  map.types = "Esri.WorldImagery",
  popup = leafpop:::popupIframe(img),
  label = "GeoCampus Lankwitz"
)

leaflet package

Another great package for producing interactive maps in R is the leaflet package. For more information and a thorough introduction, we recommend to review its website Leaflet for R.

We will reproduce the example above using the mapview package.

In order to get started, we install the leaflet package by typing install.packages("leaflet"). Again, we retrieve the coordinates for the GeoCampus Lankwitz using the geocode_OSM() function from the tmaptools package. In this case, the extracted lat and lon coordinates are already enough to plot the map, and we do not need to convert them to a spatial feature.

library(tmaptools)

geocampus <- geocode_OSM(q = "Malteserstrasse 74-100, Berlin", as.data.frame = TRUE)
geocampus_latlon <- geocampus[, c("lat", "lon")]
geocampus_latlon
##        lat      lon
## 1 52.42625 13.35803

We then create a map widget using the leaflet() function. Note that leaflet maps are commonly produced using the pipe operator %>% of the magrittr package. It makes our code much more readable and allows us to gradually add elements to our map.

library(leaflet)

img <- "http://www.osa.fu-berlin.de/geographie_bsc/_medien/bilder_geocampus/bild_010/010_546.jpeg"

geocampus_latlon %>%
  leaflet() %>%
  addProviderTiles("Esri.WorldImagery") %>%
  addMarkers(~lon, ~lat, label = "GeoCampus Lankwitz", popup = paste0("<img src = ", img, ">"))

Citation

The E-Learning project SOGA-R was developed at the Department of Earth Sciences by Kai Hartmann, Joachim Krois and Annette Rudolph. You can reach us via mail by soga[at]zedat.fu-berlin.de.

Creative Commons License
You may use this project freely under the Creative Commons Attribution-ShareAlike 4.0 International License.

Please cite as follow: Hartmann, K., Krois, J., Rudolph, A. (2023): Statistics and Geodata Analysis using R (SOGA-R). Department of Earth Sciences, Freie Universitaet Berlin.