Search code examples
rggplot2tidyversecoordinateslatitude-longitude

How to get latitude and longitude coordinates on same coordinate system as state boundary map


I have some coordinates, say:

library(tidyverse)
library(haven)
library(tidycensus)
library(tigris)

coords <- data.frame(lat = c(38.09720, 36.85298, 31.31517, 21.48344), long = c(-121.38785, -75.97742, -85.85522, -158.03648))

I then get a map of the US:

geo <- get_acs(geography = "state",
               variables = c(x = "B04006_036"),
               year = 2021, 
               geometry = TRUE, 
               keep_geo_vars=TRUE) %>%
    filter(STATEFP!="72")

#to get alaska and hawaii in the picture
geo <- shift_geometry(geo)

I then attempt to plot the state map, with the coordinates overlaid:

ggplot(data = coords) +
    geom_point(aes(x=lat, y=long)) +
    geom_sf(fill = "transparent", color = "gray50", size = 1, data = geo %>% group_by(STATEFP) %>% summarise()) +   
    theme(panel.background = element_rect(fill = 'white')) +
    theme(panel.grid = element_blank(),axis.title = element_blank(),
          axis.text = element_blank(),axis.ticks = element_blank(),
          panel.border = element_blank())

This produces:

enter image description here

However, this doesn't work as it produces a map with all coordinates seemingly in the same place. How can I modify things so that the map and the coordinates are on the same scale?


Solution

  • This is simply a projection error! You need to tell R (more specifically sf) what coordinate reference system your coordinates are in. I will assume they are in the typical WGS 1984 (EPSG:4326). Then we can create an sf object from your coords df and plot with geom_sf.

    library(tidyverse)
    library(haven)
    library(tidycensus)
    library(tigris)
    library(sf)
    
    coords <- data.frame(lat = c(38.09720, 36.85298, 31.31517, 21.48344), long = c(-121.38785, -75.97742, -85.85522, -158.03648)) |>
      st_as_sf(coords = c("long", "lat"), crs = 4326) %>% # specify coordinate columns and coordinate reference system
      shift_geometry() # remember to shift the coords too!
    
    geo <- get_acs(geography = "state",
                   variables = c(x = "B04006_036"),
                   year = 2021, 
                   geometry = TRUE, 
                   keep_geo_vars=TRUE) %>%
        filter(STATEFP!="72")
    
    geo <- shift_geometry(geo)
    
    ggplot() +
        geom_sf(data = coords) + # now you can use geom_sf
        geom_sf(fill = "transparent", color = "gray50", size = 1, data = geo %>% group_by(STATEFP) %>% summarise()) +   
        theme(panel.background = element_rect(fill = 'white')) +
        theme(panel.grid = element_blank(),axis.title = element_blank(),
              axis.text = element_blank(),axis.ticks = element_blank(),
              panel.border = element_blank())
    

    correct plot