rnaturalearth and sf object

Question

I would like to construct a map with rnaturalearth.

On top of that map I want to have three dots representing three cities, which are connected with lines. Line color corresponds to the value of F (indicates how close those countries are), with the scale color for the map (as it is on the figure)

I have the code for map:

target <- c("Armenia", "Azerbaijan", "Turkey", "Iran", "Georgia", "Russia", "Syria", "Iraq")


countries<- ne_countries(scale = 10,returnclass = 'sf') %>%
  filter(name %in% target)

coordinates <- data.frame(
  region = c("A", "B", "C"),
  longitude = c(39.80949367, 39.65485163, 38.47802923),
  latitude = c(46.7793274, 36.99614676, 43.3863868)
)

map <- ggplot(data=countries)+geom_sf(fill = "lightgoldenrodyellow") +
  xlim(36,49) +
  ylim(36.5,42) + theme_void() +
  theme(panel.background = element_rect(fill = "lightblue2"))

print(map) 

How to plot dots for cities and colored lines?

Thank you!

Answer 1

Your problem shows a bit of a mix of solutions. Thus, I offer 2 Options below. But read up on what you use. For example, how ggplot works, how data frame manipulation are implemented, or what sf does/offers.

For a start:

{ggplot} works with layers that you add on top of each other. This is symbolised by the + operator.

You have produced already the underlying map with your code. For this you worked with the mapdata from {rnaturalearth} and used {sf} for handling the geospatial aspects. For example, this helps a bit with the projection of the map (or you can tap into other functions of the {sf} package.

library(dplyr)    # for data crunching - you use the pipe %>% or now |>
library(ggplot2)  # for plotting
library(sf)       # handling of geospatial data
library(rnaturalearth) # underlying geo map data

# ... your code here ....
map     # calling your base map

This delivers your base map:

Option 1: combining {ggplot} and {sf} ... takes care of geospatial placement

There is one side-aspect of combining {ggplot2} and {sf}. The latter ensures that "ordinary" data frame columns are handled as coordinates (without being a sf-object).
Thus, you could "add" a point layer on too of your base-map.

Note: I also add a text label to show another problem.

map + 
# now let's add a point geom on top of the base map, given by coordinates
  geom_point( data = coordinates
            , mapping = aes(x = longitude, y = latitude)
    # ----- make it big(ger) to see it better
           , size = 5
    ) + 

# -------- we also add a text label to show where we are
geom_text( data = coordinates
    # ----- we displace the label by adding to the longitude
    # ----- note this is for demo purposes, be mindful about 0.5 is decimal (geospatial) degree in this case ... dependent on your zoom/area this jump can differ significantly (but that is another discussion).
         , mapping = aes(x = longitude + 0.5, y = latitude, label = region)
    # ----- arbitrary setting of size to make it "big" and color
    # ----- note: you may use this to have a labelled point by removing the offset
         , size = 8, color = "red"
    )

So we now have (one) point super-imposed on the base-map. Unfortunately, you picked an area that only shows B.

We can now use the same side-effect to plot links between the cities. {ggplot} allows you to plot line-segments defined by a start and end position. Let's create a links data frame:

links <- coordinates |> 
#---------- create end points, i.e. lon2, lat2
   mutate(longitude2 = lead(longitude, default = first(longitude))
         , latitude2 = lead(latitude, default = first(latitude))
#---------- we also add a label for the links and value F
         , name = c("AB", "BC", "CA")
         , F = c(3,6,9)
    ) 
 links
  region longitude latitude longitude2 latitude2 name F
1      A  39.80949 46.77933   39.65485  36.99615   AB 3
2      B  39.65485 36.99615   38.47803  43.38639   BC 6
3      C  38.47803 43.38639   39.80949  46.77933   CA 9

If you plot this .... unfortunately no lines are shown. {ggplot} gives you a warning about missing (points, text, and segments). We have seen this already before ... but it is clear that the area you picked does not allow to show all elements (as they are outside your viewing area from a geo-spatial perspective)

map + 
    #---- add city/ies from above
    geom_point(data = coordinates, aes(x = longitude, y = latitude), size = 5) + 
    geom_text(data = coordinates, aes(x = longitude + 0.5, y = latitude, label = region), size = 8, color = "red") +
    
    #---- add connection links we just have created
    geom_segment(data = links
               , mapping = aes(x = longitude, xend = longitude2  #note end point
                             , y = latitude,  yend = latitude2
    #-------------------- just adding a color to highlight the link value
                             , color = F)
    )

To fix the "missing" points/links issue, check your limits. You have picked latitude values that eliminate A and C.
For this change ylim accordingly:

map2 <- ggplot(data=countries)+geom_sf(fill = "lightgoldenrodyellow") +
         xlim(36,49) +
    # ---------------- set upper latitude to cover all points, i.e. 47
         ylim(36.5,47) + 
         theme_void() + theme(panel.background = element_rect(fill = "lightblue2"))

map2 + 
    #---- add city/ies ---------------------
    geom_point(data = coordinates, aes(x = longitude, y = latitude), size = 5) + 
    geom_text(data = coordinates, aes(x = longitude + 0.5, y = latitude, label = region), size = 8, color = "red") +
    #---- add connection links -------------
    geom_segment(data = links, mapping = aes(x = longitude, xend = longitude2, y = latitude, yend = latitude2, color = F))

To fix the scale, convert your F values into "discrete" variables. {ggplot} will interpret them as continuous. For example, you can turn them into a factor or character when creating the data frame (or inside your plot call, (i) mutate( ... F = as.character(c(3,6,9) ...) or (ii) geom_segment( ... color = as.character(F) ...).

Option 2: using {sf}

As you picked {sf} for the base map, you can also address the problem in {sf}.
The key here is that you have to "cast" your data frames of points and links to so-called sf-objects.

Working with {sf} gives you access to powerful geospatial tools. However, it is not straight forward as you need to read up on other (geospatial) concepts.

For a start:

• standard positions have the geometry POINT
• "lines" are interpreted as geometry LINESTRING

I leave the other geometries out here.

Coercing from a data frame of POINT positions to sf-object is relatively straight forward.

For a LINESTRING we need to reshape the data, as a linestring can conceptually comprise multiple points (and not just the start and end point).
For our case we create a long table of the start and end positions, label then with the name and F value (which are "constant" along the link) by using group_by() and then coerce the linestrings of start- and end-positions.
Work through the code line-by-line if you want to get a better grip on the operations.

library(tidyr) # for pivoting

links_ls <- links |>
    #---- reformat the positions to resemble POINTS ----------
    tidyr::unite(start, longitude, latitude) |>
    tidyr::unite(end, longitude2, latitude2) |>
    #---- turn into a long table
    tidyr::pivot_longer(
        cols = start:end
      , names_to = "start_end"
      , values_to = "coords"
    ) |>
    #------ convert text coordinates back to numeric columns
    tidyr::separate(coords, c("LON", "LAT"), sep = "_")  |>  
    dplyr::mutate_at(vars(LON, LAT), as.numeric) |>
    #------ now cast them to WGS84 positions
    sf::st_as_sf(coords = c("LON", "LAT"), crs = 4326) |>

    #------ we can now cast the grouped points to LINESTRING
    dplyr::group_by(name, F) |>
    dplyr::summarise(.groups = "drop") |>
    sf::st_cast("LINESTRING", crs = 4326)
    

We now have a sf-object for our links (city pair connections).

With this we can go back to {ggplot) and use geom_sf()-layers:

# ----- let's start from our corrected base-map
map2 + 

#------ add a layer of points for our cities
#------ note: increase size to have the text printed inside
  geom_sf(data = cities, size = 8) + 
#------ add the labels to the points, note change color to make
#------ them different from the point color (and thus visible)
  geom_sf_text(data = cities, aes(label = region), color = "white") + 

#------ finally add the connecting links
#------ to force a discrete scale we turn F to character
  geom_sf(data = links_ls, aes(color = as.character(F)))

This yields:

You can now beautify the visualisation by renaming the scale, changing colors, etc.

Good luck.