I created a relative simple map for a systematic review to represent where the included studies have been performed. The ‘rworldmap’ package was perfect for me. However, in a more recent systematic review, we would combined countries or cities in the same map. Indeed, authors could analyse the data from an (inter)national survey or collect data in one or several cities. It was challenging for (a DIY in R like) me. Consequently, I share my R code with fake data.

Useful packages

library (rworldmap)
library (maps)
library (tidyverse)

Data preparation for countries related data
To identify the “ISO” for each countries included in the “rworldmap” packages.
There is 2 examples.

df_contries_freq <- data.frame(Country = c(
  "Argentina", "Australia", "Austria","Belarus", "Belgium", "Bulgaria", "Canada", "Chile", "Colombia", "Costa Rica", "Croatia","Cyprus", "Czech Rep.", "Denmark", "Dominican Rep.","Estonia", "Fiji", "Finland", "France","Georgia", "Germany", "Greece", "Hungary", "India", "Indonesia", "Iran", "Ireland", "Israel","Italy", "Japan", "Jordan", "Latvia", "Lithuania", "Luxembourg", "Malaysia", "Mexico", "Morocco","Netherlands", "New Zealand", "Norway", "Philippines", "Poland", "Portugal", "Romania", "Russia", "Saudi Arabia",
  "Slovakia", "Slovenia", "Solomon Islands", "South Africa", "South Korea", "Spain", "Sri Lanka", "Sweden", "Switzerland","Thailand", "Turkey", "Ukraine", "United Kingdom", "United States", "China"),
  Freq_countries = c( 1, 1, 1, 9, 1, 1, 7, 1, 1, 7, 1, 18, 1, 1, 1, 2, 1, 8, 1, 5, 3, 1, 1, 6, 1, 1, 4, 8, 1, 1, 11, 1, 11, 1, 13, 1, 1, 17, 1, 1, 3, 1, 7, 1, 1, 6, 3, 1, 4, 1, 1, 1, 1, 1, 1, 16, 1, 3, 1, 4, 2 ))

df_contries_freq_iso <- left_join(df_contries_freq, iso, by =c("Country"))

df_contries_freq02 <- data.frame(Country = c("East Timor","Guam","Malta","Serbia","Singapore","United Arab Emirates","Vietnam"), Freq_countries = c(7,1,4,1,6,1,8), ISO3V10 = c("TLS","GUM","MLT","SRB","SGP","ARE","VNM"))

df_contries_freq_iso2 <- full_join(df_contries_freq_iso, df_contries_freq02, by =c("Country", "Freq_countries","ISO3V10"))

Data preparation for cities with “world.cities” data from “maps” package

plot_cities <- subset(world.cities, 
                     name %in% c("Beijing","Porto Alegre","Boston","Rotterdam","Tokyo","Bangkok", "Boston","Taipei","Seattle","Niamey","Maebashi","Tsukuba","Rochester","New York","Arnhem", "Groningen","Essen","Bochum","Mülhe", "Cambridge","Sydney","Shanghai","Phoenix","Tucson","Sacramento","Framingham","Minneapolis","New York","Pittsburgh","Paris","Montpellier","Marseille","Lille","Caen")
                     & country.etc %in% c("China","Brazil","USA","Netherlands","Japan","Thailand","USA","Taiwan","Niger","Japan","Japan","USA","USA","Netherlands","Netherlands","Germany","Germany","Germany","USA","Australia","China","USA","USA","USA","USA","USA","USA","USA","France","France", "France", "France", "France"))

To prepare the worldwide map without Antartica and with the Viridis palette (easier to read by those with colorblindness)

Map_for_review <- joinCountryData2Map(df_contries_freq_iso2, joinCode = "ISO3", nameJoinColumn = "ISO3V10")
Map_for_review <- subset(Map_for_review, continent != "Antarctica")

colourPalette <- viridis::viridis(7, direction = -1) #viridis palette

A <- mapCountryData(Map_for_review, nameColumnToPlot="Freq_countries", missingCountryCol = gray(.8), catMethod = "logFixedWidth", colourPalette = colourPalette, addLegend=FALSE)
points(plot_cities$long, plot_cities$lat, pch =  21, col = "black", bg="white", cex=1)
do.call( addMapLegend, c(A, legendWidth=0.5, legendMar = 8))

UPDATE : I also found another useful tutorial about harvest plots for systematic review