Please have a look at the reprex at the end of the plot. My goal (which I achieve step by step) is to construct a panel plot combining different color maps of European countries while removing in the map the inner border in Cyprus (no political statement on my side, I have just to do it). I somehow achieve this, but in a rather cumbersome way, so I wonder if someone can suggest some improvements.
Thanks!
library(tidyverse)
library(rnaturalearth)
#> Support for Spatial objects (`sp`) will be deprecated in {rnaturalearth} and will be removed in a future release of the package. Please use `sf` objects with {rnaturalearth}. For example: `ne_download(returnclass = 'sf')`
library(sf)
#> Linking to GEOS 3.11.1, GDAL 3.6.2, PROJ 9.1.1; sf_use_s2() is TRUE
## data I want to plot in a facet plot.
df_plot <- structure(list(expenditure_year = c(2019, 2019, 2019, 2019, 2019,
2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019,
2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019, 2019,
2020, 2020, 2020, 2020, 2020, 2020, 2020, 2020, 2020, 2020, 2020,
2020, 2020, 2020, 2020, 2020, 2020, 2020, 2020, 2020, 2020, 2020,
2020, 2020, 2020, 2020, 2020, 2021, 2021, 2021, 2021, 2021, 2021,
2021, 2021, 2021, 2021, 2021, 2021, 2021, 2021, 2021, 2021, 2021,
2021, 2021, 2021, 2021, 2021, 2021, 2021, 2021, 2021, 2021, 2022,
2022, 2022, 2022, 2022, 2022, 2022, 2022, 2022, 2022, 2022, 2022,
2022, 2022, 2022, 2022, 2022, 2022, 2022, 2022, 2022, 2022, 2022,
2022, 2022, 2022, 2022), iso2 = c("IT", "FI", "FR", "NL", "CY",
"BE", "SE", "GR", "DE", "LT", "BG", "PT", "CZ", "IE", "EE", "HU",
"SK", "ES", "PL", "SI", "DK", "AT", "RO", "LU", "MT", "LV", "HR",
"IT", "FI", "FR", "NL", "CY", "BE", "SE", "GR", "DE", "LT", "BG",
"PT", "CZ", "IE", "EE", "HU", "SK", "ES", "PL", "SI", "DK", "AT",
"RO", "LU", "MT", "LV", "HR", "IT", "FI", "FR", "NL", "CY", "BE",
"SE", "GR", "DE", "LT", "BG", "PT", "CZ", "IE", "EE", "HU", "SK",
"ES", "PL", "SI", "DK", "AT", "RO", "LU", "MT", "LV", "HR", "IT",
"FI", "FR", "NL", "CY", "BE", "SE", "GR", "DE", "LT", "BG", "PT",
"CZ", "IE", "EE", "HU", "SK", "ES", "PL", "SI", "DK", "AT", "RO",
"LU", "MT", "LV", "HR"), percentage_share = c(0.450577656119157,
0.92248580410076, 0.975584667105617, 0.409803641820049, 0.539818102589104,
0.96123699068389, 0.826519143748971, 0.570729391308522, 1.51184843024544,
1.56004652853805, 0.682394691276564, 0.529918189934815, 1.51537052534534,
0.340282452391896, 1.23411076526779, 1.78839893466253, 0.648651854236538,
0.474837572951868, 1.11702742041526, 0.858524606698324, 1.43464767790405,
0.505813688224417, 0.655039437261184, 0.310317708208196, 1.6984935879487,
1.07006270258955, 1.33269223306193, 2.02652222755559, 1.52387828834052,
2.27427055972995, 0.875130252470089, 1.50963072632275, 1.43394831922113,
1.02805626070114, 4.31903491607162, 3.71236026065522, 2.00271869460953,
1.43822454206211, 1.63660567657213, 2.74250293931477, 0.407182979009171,
1.69199343784178, 3.12677589492696, 1.76946859138244, 0.88979342454491,
5.16611685495998, 3.54047898922094, 2.48037431884307, 1.99660745861322,
1.23774421769849, 0.879255567282404, 4.24216627476823, 2.70305053222405,
2.66814646113911, 2.11099250992334, 1.59620461257039, 2.69079496650438,
1.50012485828527, 2.13084813620244, 1.29750144508215, 1.04199044469184,
3.03512300799337, 3.17970791856142, 1.58793142828813, 1.68604027084116,
1.42812640590373, 3.62657019057752, 0.721667603152213, 1.41210754275081,
3.61800356104344, 2.22916452630624, 1.11319571460128, 1.81329721467997,
3.00851645408846, 1.56776732632244, 2.46506682189924, 1.58129951074308,
0.813434327102068, 3.99050063429371, 2.81167234901301, 2.20079792507412,
1.41830583744773, 0.944648699957193, 1.71235447267469, 1.03518425244823,
0.538183029124816, 1.02736324514123, 0.893949584792747, 1.37693432013489,
2.22128961439947, 1.13396128434949, 1.01149050905179, 0.950283072433532,
1.43451796238923, 0.487096865306793, 0.92012157362591, 2.08815527711486,
0.884735492296972, 1.27188526690518, 0.927269440893475, 0.84679781968768,
1.39310961513098, 1.45742166676803, 1.46526610855142, 0.665295824788144,
1.50079360929557, 1.00760277849241, 0.839755449738563), iso3 = c("ITA",
"FIN", "FRA", "NLD", "CYP", "BEL", "SWE", "GRC", "DEU", "LTU",
"BGR", "PRT", "CZE", "IRL", "EST", "HUN", "SVK", "ESP", "POL",
"SVN", "DNK", "AUT", "ROM", "LUX", "MLT", "LVA", "HRV", "ITA",
"FIN", "FRA", "NLD", "CYP", "BEL", "SWE", "GRC", "DEU", "LTU",
"BGR", "PRT", "CZE", "IRL", "EST", "HUN", "SVK", "ESP", "POL",
"SVN", "DNK", "AUT", "ROM", "LUX", "MLT", "LVA", "HRV", "ITA",
"FIN", "FRA", "NLD", "CYP", "BEL", "SWE", "GRC", "DEU", "LTU",
"BGR", "PRT", "CZE", "IRL", "EST", "HUN", "SVK", "ESP", "POL",
"SVN", "DNK", "AUT", "ROM", "LUX", "MLT", "LVA", "HRV", "ITA",
"FIN", "FRA", "NLD", "CYP", "BEL", "SWE", "GRC", "DEU", "LTU",
"BGR", "PRT", "CZE", "IRL", "EST", "HUN", "SVK", "ESP", "POL",
"SVN", "DNK", "AUT", "ROM", "LUX", "MLT", "LVA", "HRV"), country = c("Italy",
"Finland", "France", "Netherlands", "Cyprus", "Belgium", "Sweden",
"Greece", "Germany", "Lithuania", "Bulgaria", "Portugal", "Czechia",
"Ireland", "Estonia", "Hungary", "Slovakia", "Spain", "Poland",
"Slovenia", "Denmark", "Austria", "Romania", "Luxembourg", "Malta",
"Latvia", "Croatia", "Italy", "Finland", "France", "Netherlands",
"Cyprus", "Belgium", "Sweden", "Greece", "Germany", "Lithuania",
"Bulgaria", "Portugal", "Czechia", "Ireland", "Estonia", "Hungary",
"Slovakia", "Spain", "Poland", "Slovenia", "Denmark", "Austria",
"Romania", "Luxembourg", "Malta", "Latvia", "Croatia", "Italy",
"Finland", "France", "Netherlands", "Cyprus", "Belgium", "Sweden",
"Greece", "Germany", "Lithuania", "Bulgaria", "Portugal", "Czechia",
"Ireland", "Estonia", "Hungary", "Slovakia", "Spain", "Poland",
"Slovenia", "Denmark", "Austria", "Romania", "Luxembourg", "Malta",
"Latvia", "Croatia", "Italy", "Finland", "France", "Netherlands",
"Cyprus", "Belgium", "Sweden", "Greece", "Germany", "Lithuania",
"Bulgaria", "Portugal", "Czechia", "Ireland", "Estonia", "Hungary",
"Slovakia", "Spain", "Poland", "Slovenia", "Denmark", "Austria",
"Romania", "Luxembourg", "Malta", "Latvia", "Croatia")), row.names = c(NA,
-108L), class = c("tbl_df", "tbl", "data.frame"))
ww_ini <- ne_countries(scale = "medium",
type = 'countries',
returnclass = "sf")
bb <- ne_download(type = "wgs84_bounding_box", category = "physical",
returnclass = "sf")
### select a subset of the data for a single year to plot
df_plot_sel <- df_plot |>
filter(expenditure_year==2022)
ww_plot <- ww_ini |>
select(adm0_a3) |>
left_join(y=df_plot_sel, by=c("adm0_a3"="iso3"))
## start with a map plot without colors
gpl <- ggplot(data = ww_plot) +
geom_sf( aes(fill=percentage_share), col = "black", lwd = 0.3 )+
coord_sf(xlim=c(-20,45), ylim=c(30, 73) ) +
theme_minimal()
gpl
#### The above is OK, but I need to remove the boundary inside Cyprus (no political statement whatsoever from my side).
## See https://stackoverflow.com/questions/68672575/r-rnaturalearth-and-sf-remove-a-single-border-from-map
cyprus <- ww_ini |>
select(adm0_a3) |>
filter(adm0_a3 %in% c("CYP", "CYN")) |>
mutate(adm0_a3 = "CYP") |> ## I renamed as Cyprus my version
group_by(adm0_a3) |>
summarise()
ww_ini_cyp <- ww_ini |>
select(adm0_a3) |>
filter(!adm0_a3 %in% c("CYP", "CYN")) |>
bind_rows(cyprus)
gpl_cyp <- ggplot(data = ww_ini_cyp) +
geom_sf( col = "black", lwd = 0.3 )+
coord_sf(xlim=c(-20,45), ylim=c(30, 73) ) +
theme_minimal()
### and this is OK. I now have a map without the inner border in Cyprus
gpl_cyp
### Now let us replot a colormap with Cyprus without border
ww_plot_cyp <- ww_ini_cyp |>
## select(adm0_a3) |>
left_join(y=df_plot_sel, by=c("adm0_a3"="iso3"))
gpl_cyp <- ggplot(data = ww_plot_cyp) +
geom_sf( aes(fill=percentage_share), col = "black", lwd = 0.3 )+
coord_sf(xlim=c(-20,45), ylim=c(30, 73) ) +
theme_minimal()
gpl_cyp ## and this works
## Now let us try a facet plot
ww_plot_cyp_facet <- ww_ini_cyp |>
left_join(y=df_plot, by=c("adm0_a3"="iso3"))
gpl_cyp_facet <- ggplot(data = ww_plot_cyp_facet) +
geom_sf( aes(fill=percentage_share), col = "black", lwd = 0.3 )+
facet_wrap( ~expenditure_year, nrow = 2, scales = "fixed")+
coord_sf(xlim=c(-20,45), ylim=c(30, 73) ) +
theme_minimal()
gpl_cyp_facet ## and this is a disaster. How do I fix it?
## The problem is that the expenditure_year is missing for adm0_a3 codes
## I just need to create it for the all the adm0_a3 codes
years <- df_plot|>
pull(expenditure_year) |>
unique()
iso_list <- ww_ini$adm0_a3 |>
unique()
df <- crossing(years, iso_list)
ww_plot_cyp_facet_fix <- ww_ini_cyp |>
left_join(df, by=c("adm0_a3"="iso_list")) |>
left_join(df_plot, by=c("adm0_a3"="iso3",
"years"="expenditure_year"))
gpl_cyp_facet_fix <- ggplot(data = ww_plot_cyp_facet_fix) +
geom_sf( aes(fill=percentage_share), col = "black", lwd = 0.3 )+
facet_wrap( ~years, nrow = 2, scales = "fixed")+
coord_sf(xlim=c(-20,45), ylim=c(30, 73) ) +
theme_minimal()
gpl_cyp_facet_fix ## and this works, but if feels very cumbersome. Any
## suggestion is appreciated.
sessionInfo()
#> R version 4.3.2 (2023-10-31)
#> Platform: x86_64-pc-linux-gnu (64-bit)
#> Running under: Debian GNU/Linux 12 (bookworm)
#>
#> Matrix products: default
#> BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.11.0
#> LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.11.0
#>
#> locale:
#> [1] LC_CTYPE=en_GB.UTF-8 LC_NUMERIC=C
#> [3] LC_TIME=en_GB.UTF-8 LC_COLLATE=en_GB.UTF-8
#> [5] LC_MONETARY=en_GB.UTF-8 LC_MESSAGES=en_GB.UTF-8
#> [7] LC_PAPER=en_GB.UTF-8 LC_NAME=C
#> [9] LC_ADDRESS=C LC_TELEPHONE=C
#> [11] LC_MEASUREMENT=en_GB.UTF-8 LC_IDENTIFICATION=C
#>
#> time zone: Europe/Brussels
#> tzcode source: system (glibc)
#>
#> attached base packages:
#> [1] stats graphics grDevices utils datasets methods base
#>
#> other attached packages:
#> [1] sf_1.0-14 rnaturalearth_0.3.4 lubridate_1.9.3
#> [4] forcats_1.0.0 stringr_1.5.0 dplyr_1.1.3
#> [7] purrr_1.0.2 readr_2.1.4 tidyr_1.3.0
#> [10] tibble_3.2.1 ggplot2_3.4.4 tidyverse_2.0.0
#>
#> loaded via a namespace (and not attached):
#> [1] s2_1.1.4 utf8_1.2.4 generics_0.1.3
#> [4] class_7.3-22 KernSmooth_2.23-22 stringi_1.7.12
#> [7] lattice_0.22-5 hms_1.1.3 digest_0.6.33
#> [10] magrittr_2.0.3 evaluate_0.23 grid_4.3.2
#> [13] timechange_0.2.0 fastmap_1.1.1 jsonlite_1.8.7
#> [16] e1071_1.7-13 DBI_1.1.3 httr_1.4.7
#> [19] fansi_1.0.5 scales_1.2.1 cli_3.6.1
#> [22] rlang_1.1.2 units_0.8-4 munsell_0.5.0
#> [25] reprex_2.0.2 withr_2.5.2 yaml_2.3.7
#> [28] tools_4.3.2 tzdb_0.4.0 colorspace_2.1-0
#> [31] vctrs_0.6.4 R6_2.5.1 proxy_0.4-27
#> [34] lifecycle_1.0.4 classInt_0.4-10 fs_1.6.3
#> [37] pkgconfig_2.0.3 pillar_1.9.0 gtable_0.3.4
#> [40] glue_1.6.2 Rcpp_1.0.11 rnaturalearthdata_0.1.0
#> [43] xfun_0.41 tidyselect_1.2.0 knitr_1.45
#> [46] farver_2.1.1 htmltools_0.5.7 labeling_0.4.3
#> [49] rmarkdown_2.25 wk_0.9.0 compiler_4.3.2
#> [52] sp_2.1-1
Created on 2023-11-24 with reprex v2.0.2
If you can switch from Natural Erath data to Eurostat GISCO shapes (accessible through {giscoR}
package), you'd get just one shape for Cyprus.
When joining sf object and a data.frame, consider using right_join()
to get sf object as a result while keeping all records of df_plot
.
library(sf)
#> Linking to GEOS 3.11.2, GDAL 3.6.2, PROJ 9.2.0; sf_use_s2() is TRUE
library(dplyr)
library(ggplot2)
# using Eurostat GISCO shapes and ETRS89 / ETRS-LAEA projection
countries_sf <- giscoR::gisco_get_countries(year = "2020", epsg = "3035")
bb_3035 <- st_bbox(c(xmin = -20, ymin = 30, xmax = 45, ymax = 73), crs = "WGS84") |>
st_as_sfc() |>
st_transform(st_crs(countries_sf))
countries_sf <- st_crop(countries_sf, bb_3035)
# we want to have sf object as a first argument while joining
# to df_plot so all df_plot records would remain, hence the right_join
countries_sf |>
right_join(df_plot, by = join_by(ISO3_CODE == iso3)) |>
ggplot() +
geom_sf(data = countries_sf) +
geom_sf(aes(fill = percentage_share), col = "black", linewidth = 0.3) +
facet_wrap(~expenditure_year, nrow = 2) +
coord_sf(expand = FALSE) +
theme_minimal() +
theme(legend.position = "bottom")
Created on 2023-11-24 with reprex v2.0.2