14.4 Define metropolitan areas

We define metropolitan areas here as pixels of 20 km2 inhabited by more than 500,000 people

pop_agg <- terra::aggregate(reclass$pop, fact = 20, fun = sum, na.rm = TRUE)

pop_agg <- pop_agg[pop_agg > 500000, drop = FALSE] 

summary(pop_agg)

##       pop         
##  Min.   : 510335  
##  1st Qu.: 586917  
##  Median : 782416  
##  Mean   : 771034  
##  3rd Qu.: 823811  
##  Max.   :1204870  
##  NA's   :1431

metros <- pop_agg |> 
  terra::patches(directions = 8) |>
  terra::as.polygons() |>
  sf::st_as_sf()

metros$names = c("Hamburg", "Berlin", "Düsseldorf", "Leipzig",
                 "Frankfurt am Main", "Nürnberg", "Stuttgart", "München")

metros_points <- st_centroid(metros)

## Warning: st_centroid assumes attributes are constant over geometries

tm_shape(pop_agg/1000) +
  tm_raster(palette = "GnBu",
            title = "Number of people\n(in 1,000)") +
  tm_shape(ger) +
  tm_borders() +
  tm_shape(metros) +
  tm_borders(col = "gold", lwd = 2) +
  tm_shape(metros_points) +
  tm_text(text = "names", ymod = 0.6, shadow = TRUE, size = 0.75,
          fontface = "italic") +
  tm_layout(legend.outside = TRUE, legend.outside.size = 0.3)

We could also reverse geocode for the city names

metro_names <- sf::st_centroid(metros, of_largest_polygon = TRUE) |>
  tmaptools::rev_geocode_OSM(as.data.frame = TRUE) |>
  select(city, town, state)

## Warning: st_centroid assumes attributes are constant over geometries

metro_names <- dplyr::mutate(metro_names, city = ifelse(is.na(city), town, city))

metro_names

##                city    town               state
## 1           Hamburg    <NA>                <NA>
## 2            Berlin    <NA>                <NA>
## 3           Velbert Velbert Nordrhein-Westfalen
## 4           Leipzig    <NA>             Sachsen
## 5 Frankfurt am Main    <NA>              Hessen
## 6          Nürnberg    <NA>              Bayern
## 7         Stuttgart    <NA>   Baden-Württemberg
## 8           München    <NA>              Bayern