library(keys.lid)
#> SWMM executable not found.
performances <- keys.lid::performances %>%
dplyr::mutate(
scenario_name = kwb.utils::multiSubstitute(strings = .data$scenario_name,
replacements = list("through-trench" = "mulde_rigole",
"through" = "mulde",
"with-berm" = "mit-berme",
"no-berm" = "ohne-berme",
"no-drainmat" = "keine-drainagematte",
"with-drainmat" = "mit-drainagematte",
"no-drainage" = "keine-drainage",
"with-drainage" = "mit-drainage",
"per.hour" = "pro.Stunde",
"extensive" = "extensiv",
"intensive" = "intensiv"))
) %>%
dplyr::group_by(.data$lid_name_tidy,
.data$scenario_name) %>%
dplyr::mutate(scenario_id = dplyr::cur_group_id())
performances_without_lids <- performances %>%
dplyr::filter(.data$lid_name_tidy == "bioretention_cell",
.data$lid_area_fraction == 0,
.data$scenario_id == 1) %>%
dplyr::mutate(lid_name_tidy = "Referenz",
scenario_name = "Regenr\u00FCckhalt ohne LID")
mycolors <- c(rev(RColorBrewer::brewer.pal(name="Reds", n = 6)),
RColorBrewer::brewer.pal(name="Blues", n = 6))
vrr_reference <- performances_without_lids %>%
tidyr::unnest(.data$annual) %>%
ggplot2::ggplot(mapping = ggplot2::aes(x = .data$zone_id,
y = 100*.data$vrr)) +
ggplot2::scale_y_continuous(labels = scales::percent_format(scale = 1)) +
ggplot2::geom_boxplot() +
ggplot2::geom_jitter(mapping = ggplot2::aes(col = factor(.data$year)),
size=2, alpha=0.8) +
ggplot2::scale_color_manual(values = mycolors) +
ggplot2::stat_summary(ggplot2::aes(label=sprintf("%d %%", round(..y..,0))),
fun.y = median,
geom="text",
size = 3,
fontface = "bold",
vjust = -0.5,
) +
ggplot2::labs(title = "Modelleinzugsgebiet (100% Versiegelungsgrad)",
x = "Klimazone",
y = "J\u00E4hrlicher Regenr\u00FCckhalt ohne LIDs",
col = "Jahr") +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="top")
#> Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
#> ℹ Please use `"annual"` instead of `.data$annual`
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
#> Warning: The `fun.y` argument of `stat_summary()` is deprecated as of ggplot2 3.3.0.
#> ℹ Please use the `fun` argument instead.
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
ggplot2::ggsave("vrr_reference.jpg",
plot = vrr_reference,
scale = 1,
width = 7,
height = 5)
#> Warning: The dot-dot notation (`..y..`) was deprecated in ggplot2 3.4.0.
#> ℹ Please use `after_stat(y)` instead.
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
performances_without_lids <- performances %>%
dplyr::ungroup() %>%
dplyr::filter(.data$lid_area_fraction == 0) %>%
tidyr::unnest(.data$annual) %>%
dplyr::rename(vrr_reference = .data$vrr) %>%
dplyr::select(.data$zone_id,
.data$lid_name_tidy,
.data$scenario_id,
.data$year,
.data$vrr_reference
)
#> Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
#> ℹ Please use `"vrr"` instead of `.data$vrr`
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
#> Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
#> ℹ Please use `"zone_id"` instead of `.data$zone_id`
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
#> Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
#> ℹ Please use `"lid_name_tidy"` instead of `.data$lid_name_tidy`
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
#> Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
#> ℹ Please use `"scenario_id"` instead of `.data$scenario_id`
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
#> Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
#> ℹ Please use `"year"` instead of `.data$year`
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
#> Warning: Use of .data in tidyselect expressions was deprecated in tidyselect 1.2.0.
#> ℹ Please use `"vrr_reference"` instead of `.data$vrr_reference`
#> This warning is displayed once every 8 hours.
#> Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
#> generated.
percent_lid_area <- 100
performances_lid <- performances %>%
dplyr::filter(.data$lid_name_tidy == "bioretention_cell" & .data$lid_area_fraction == percent_lid_area/100 |
.data$lid_name_tidy == "permeable_pavement" & .data$lid_area_fraction == percent_lid_area/100 |
.data$lid_name_tidy == "green_roof" & .data$lid_area_fraction == percent_lid_area/100) %>%
tidyr::unnest(.data$annual) %>%
dplyr::left_join(performances_without_lids) %>%
dplyr::mutate(lid_name_tidy = kwb.utils::multiSubstitute(strings = .data$lid_name_tidy,
replacements = list(
"bioretention_cell" = "Versickerungsmulden",
"green_roof" = "Gr\u00FCnd\u00E4cher",
"permeable_pavement" = "Durchl\u00E4ssige Bodenbel\u00E4ge")
))
#> Joining with `by = join_by(zone_id, lid_name_tidy, year, scenario_id)`
plot_lid_performance <- function(lid_name, scenario_ids = NULL) {
if(is.null(scenario_ids)) {
scenario_ids <- performances_lid %>%
dplyr::filter(.data$lid_name_tidy == lid_name) %>%
dplyr::pull(.data$scenario_id) %>%
unique()
}
lid_area_faction <- performances_lid %>%
dplyr::filter(.data$lid_name_tidy == lid_name) %>%
dplyr::pull(lid_area_fraction) %>%
unique()*100
performances_lid %>%
dplyr::filter(.data$lid_name_tidy == lid_name,
.data$scenario_id %in% scenario_ids) %>%
dplyr::mutate(label = sprintf("%s: %s",
.data$lid_name_tidy,
.data$scenario_name)) %>%
ggplot2::ggplot(mapping = ggplot2::aes(x = .data$zone_id,
y = 100*.data$vrr)) +
ggplot2::facet_wrap(~ .data$scenario_name, ncol = 1) +
ggplot2::scale_y_continuous(labels = scales::percent_format(scale = 1)) +
ggplot2::geom_boxplot() +
ggplot2::geom_jitter(mapping = ggplot2::aes(col = factor(.data$year)),
size=2, alpha=0.8) +
ggplot2::scale_color_manual(values = mycolors) +
ggplot2::stat_summary(ggplot2::aes(label=sprintf("%d %%", round(..y..,0))),
fun.y = median,
geom="text",
size = 3,
fontface = "bold",
vjust = -0.5,
) +
ggplot2::labs(title = sprintf("%s (%d %% der Einzugsgebietsfl\u00E4che)",
lid_name,
lid_area_faction),
x = "Klimazone",
y = "J\u00E4hrlicher Regenr\u00FCckhalt mit LIDs",
col = "Jahr") +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="top")
}
lid_names <- unique(performances_lid$lid_name_tidy)
bioretention_cells <- plot_lid_performance(lid_names[1],
scenario_ids = c(1:3, 7))
ggplot2::ggsave(filename = sprintf("vrr_lid_versickerungsmulden_%d-prozent.jpg",
percent_lid_area),
plot = bioretention_cells,
scale = 1,
width = 7,
height = 2.2*5,
limitsize = FALSE)
green_roofs <- plot_lid_performance(lid_names[2],
scenario_ids = 15:18)
ggplot2::ggsave(filename = sprintf("vrr_lid_gruendaecher_%d-prozent.jpg",
percent_lid_area),
plot = green_roofs,
scale = 1,
width = 7,
height = 2.2*5,
limitsize = FALSE)
permeable_pavements <- plot_lid_performance(lid_names[3],
scenario_ids = 19:20)
ggplot2::ggsave(filename = sprintf("vrr_lid_durchlaessige-bodenbelaege_%d-prozent.jpg",
percent_lid_area),
plot = permeable_pavements,
scale = 1,
width = 7,
height = 1.1*5,
limitsize = FALSE)
create_lid_table <- function(percent_lid_area = 50,
group_cols = c("zone_id", "lid_name_tidy", "scenario_id", "scenario_name", "percent_lid_area")) {
performances_lid <- performances %>%
dplyr::filter(.data$lid_name_tidy == "bioretention_cell" & .data$lid_area_fraction == percent_lid_area/100 |
.data$lid_name_tidy == "permeable_pavement" & .data$lid_area_fraction == percent_lid_area/100 |
.data$lid_name_tidy == "green_roof" & .data$lid_area_fraction == percent_lid_area/100) %>%
tidyr::unnest(.data$annual) %>%
dplyr::left_join(performances_without_lids) %>%
dplyr::ungroup() %>%
dplyr::mutate(percent_lid_area = percent_lid_area,
vrr_diff = vrr - vrr_reference,
vrr_diff_perpercent_catchment = vrr_diff/(percent_lid_area/100),
lid_name_tidy = kwb.utils::multiSubstitute(strings = .data$lid_name_tidy,
replacements = list(
"bioretention_cell" = "Versickerungsmulden",
"green_roof" = "Gr\u00FCnd\u00E4cher",
"permeable_pavement" = "Durchl\u00E4ssige Bodenbel\u00E4ge")
)) %>%
dplyr::group_by(dplyr::across(tidyselect::all_of(group_cols))) %>%
dplyr::summarize(vrr_median = median(.data$vrr),
vrr_sd = sd(.data$vrr),
vrr_reference_median = median(.data$vrr_reference),
vrr_reference_sd = sd(.data$vrr_reference),
vrr_diff_median = median(.data$vrr_diff),
vrr_diff_sd = sd(.data$vrr_diff),
vrr_diff_perpercent_catchment_sd = sd(.data$vrr_diff_perpercent_catchment)
) %>%
dplyr::mutate(vrr_diff_perpercent_catchment_median = (.data$vrr_median - .data$vrr_reference_median)/(percent_lid_area/100)) %>%
dplyr::arrange(dplyr::desc(.data$vrr_diff_perpercent_catchment_median))
}
lid_50percent_table <- create_lid_table(percent_lid_area = 50,
group_cols = c("lid_name_tidy", "scenario_id", "scenario_name", "percent_lid_area"))
#> Joining with `by = join_by(zone_id, lid_name_tidy, year, scenario_id)`
#> `summarise()` has grouped output by 'lid_name_tidy', 'scenario_id',
#> 'scenario_name'. You can override using the `.groups` argument.
lid_50percent_table_by_zone <- create_lid_table(percent_lid_area = 50,
group_cols = c("zone_id", "lid_name_tidy", "scenario_id", "scenario_name", "percent_lid_area"))
#> Joining with `by = join_by(zone_id, lid_name_tidy, year, scenario_id)`
#> `summarise()` has grouped output by 'zone_id', 'lid_name_tidy', 'scenario_id',
#> 'scenario_name'. You can override using the `.groups` argument.
lid_100percent_table <- create_lid_table(percent_lid_area = 100,
group_cols = c("lid_name_tidy", "scenario_id", "scenario_name", "percent_lid_area"))
#> Joining with `by = join_by(zone_id, lid_name_tidy, year, scenario_id)`
#> `summarise()` has grouped output by 'lid_name_tidy', 'scenario_id',
#> 'scenario_name'. You can override using the `.groups` argument.
lid_100percent_table_by_zone <- create_lid_table(percent_lid_area = 100,
group_cols = c("zone_id", "lid_name_tidy", "scenario_id", "scenario_name", "percent_lid_area"))
#> Joining with `by = join_by(zone_id, lid_name_tidy, year, scenario_id)`
#> `summarise()` has grouped output by 'zone_id', 'lid_name_tidy', 'scenario_id',
#> 'scenario_name'. You can override using the `.groups` argument.
lid_names <- unique(lid_50percent_table_by_zone$lid_name_tidy)
plot_vrr_per_percent_lidarea <- function(climate_zone_id = 1) {
lid_50percent_table_by_zone %>%
dplyr::filter(.data$zone_id == climate_zone_id) %>%
#dplyr::mutate(label = sprintf("Klimazone %s", .data$zone_id)) %>%
ggplot2::ggplot(mapping = ggplot2::aes(x = forcats::fct_reorder(.data$scenario_name,
.data$vrr_diff_perpercent_catchment_median,
.desc = FALSE),
y = .data$vrr_diff_perpercent_catchment_median,
#col = .data$scenario_name,
fill = .data$lid_name_tidy)) +
ggplot2::scale_y_continuous(labels = scales::percent_format(scale = 1),
breaks = seq(0,1,0.1),
limits = c(0, 1)) +
ggplot2::coord_flip() +
#ggplot2::facet_wrap(~ label, ncol = 1) +
ggplot2::geom_bar(stat="identity",
#color="black",
alpha = 0.5,
position = ggplot2::position_dodge()) +
ggplot2::geom_errorbar(ggplot2::aes(ymin = .data$vrr_diff_perpercent_catchment_median - .data$vrr_diff_perpercent_catchment_sd,
ymax = .data$vrr_diff_perpercent_catchment_median + .data$vrr_diff_perpercent_catchment_sd),
width=.2,
position = ggplot2::position_dodge(.9)
) +
ggplot2::labs(title = sprintf("Klimazone %s", climate_zone_id),
y = "% j\u00E4hrlicher Regenr\u00FCckhalt / % LID-Fl\u00E4chenanteil",
x = "Designszenario",
fill = "LID") +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="top")
}
lid_plots <- lapply(1:5, function(climate_zone_id) { plot_vrr_per_percent_lidarea(climate_zone_id)})
lid_plots <- lapply(1:5, function(climate_zone_id) {
ggplot2::ggsave(filename = sprintf("lid_plot_zone%d.jpg", climate_zone_id),
plot = plot_vrr_per_percent_lidarea(climate_zone_id),
width = 9,
height = 4.5)
})
# pdff <- "lid_plots.pdf"
# mp <- gridExtra::marrangeGrob(lid_plots, nrow=3, ncol=2)
# ggplot2::ggsave(pdff,
# plot = mp,
# width = 14,
# height = 20,
# units = "cm")
#
# ggpubr::ggarrange(lid_plots[[1]] + ggplot2::ggtitle("") + ggpubr::rremove("xy.text"),
# lid_plots[[2]] + ggplot2::ggtitle("") + ggpubr::rremove("xy.text"),
# lid_plots[[3]] + ggplot2::ggtitle("") + ggpubr::rremove("xy.text"),
# lid_plots[[4]] + ggplot2::ggtitle("") + ggpubr::rremove("xy.text"),
# lid_plots[[5]] + ggplot2::ggtitle("") + ggpubr::rremove("xy.text"),
# labels = sprintf("Klimazone %s", 1:5),
# ncol = 2,
# nrow = 3)