Skip to content

Instantly share code, notes, and snippets.

@andrewheiss
Created April 13, 2021 02:43
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save andrewheiss/1d6e922fdca4798904781f6fc7b2497a to your computer and use it in GitHub Desktop.
Save andrewheiss/1d6e922fdca4798904781f6fc7b2497a to your computer and use it in GitHub Desktop.
library(tidyverse)
library(sf)
library(albersusa)
library(tidygeocoder)

usa_map <- usa_sf("longlat") %>% 
  filter(!(name %in% c("Alaska", "Hawaii")))

some_cities <- tribble(
  ~name,
  "Fort Worth, Texas", 
  "Houson, Texas",    
  "Atlanta, Georgia"                           
)

# Geocode the cities
lat_longs <- some_cities %>%
  geocode(name, method = "osm", lat = latitude , long = longitude)

# It worked! But not quite—it guessed way wrong for Houston
lat_longs
#> # A tibble: 3 x 3
#>   name              latitude longitude
#>   <chr>                <dbl>     <dbl>
#> 1 Fort Worth, Texas     32.8    -97.3 
#> 2 Houson, Texas         53.5     -2.11
#> 3 Atlanta, Georgia      33.7    -84.4

# This thing works a lot better with actual addresses instead of general cities


# Turn those numeric coordinates into actual GIS data
cities <- lat_longs %>% 
  st_as_sf(coords = c("longitude", "latitude"), crs = 4326)
cities
#> Simple feature collection with 3 features and 1 field
#> Geometry type: POINT
#> Dimension:     XY
#> Bounding box:  xmin: -97.33275 ymin: 32.75318 xmax: -2.113108 ymax: 53.53647
#> Geodetic CRS:  WGS 84
#> # A tibble: 3 x 2
#>   name                          geometry
#> * <chr>                      <POINT [°]>
#> 1 Fort Worth, Texas (-97.33275 32.75318)
#> 2 Houson, Texas     (-2.113108 53.53647)
#> 3 Atlanta, Georgia  (-84.39026 33.74899)

# Plot!
ggplot() +
  # Add state layer
  geom_sf(data = usa_map, fill = "#FF4136", color = "white", size = 0.5) +
  # Add cities layer with points
  geom_sf(data = cities, size = 3) +
  # Add cities layer with labels
  geom_sf_label(data = cities, aes(label = name),
                nudge_y = 80000, fontface = "bold") +  # Move the labels up a bunch
  # Use the Albers projection
  coord_sf(crs = st_crs("ESRI:102003")) +  
  theme_void()  # Blank theme

Created on 2021-04-12 by the reprex package (v1.0.0)

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment