This report is automatically generated with the R
package knitr
(version 1.39
)
by .
--- title: 'NPDES Allocation Comparison' author: "Mercury Program and Basin Planning Unit" date: "2/15/2023" output: html_document: code_folding: show toc: TRUE toc_float: TRUE toc_depth: 3 runtime: shiny assets: css: - "http://fonts.googleapis.com/css?family=Raleway:300" - "http://fonts.googleapis.com/css?family=Oxygen" --- <style> body{ font-family: 'Oxygen', sans-serif; font-size: 16px; line-height: 24px; } h1,h2,h3,h4 { font-family: 'Raleway', sans-serif; } .container { width: 1250px; } h3 { background-color: #D4DAEC; text-indent: 50px; } h4 { text-indent: 75px; margin-top: 35px; margin-bottom: 5px; } </style> ```{r setup, echo=TRUE, warning=FALSE, message=FALSE, fig.width=9.5} knitr::opts_chunk$set(echo=TRUE, warning=FALSE, message=FALSE, fig.width=9.5) ``` ```{r Libraries, echo=FALSE} library(kableExtra) # better formatting of tables library(ggbeeswarm) # for geom_beeswarm in ggplot - stops over plotting library(janitor) # for clean_names() library(shiny) # Had issue trying to set WD with Shiny in R project, reset working directory of rproj wd <- rstudioapi::getActiveProject() setwd(wd) source("R Functions/functions_estimate NDDNQ values.R") source("R Functions/functions_QA data.R") source("R Functions/ggplot_theme_border.r") source("R Functions/ggplot_force_panelsizes.r") file <- function(file_name){ paste0(wd, "/Reeval_Allocations/04_Percent Allocations/", file_name) } ``` ## Load Data ```{r, eval=T, echo=T} # Change chunck option to "eval=TRUE" if changes were made to "Tables 8.4a-h MeHg Inputs Per Subarea.xlsx" # e.g. Mean concentrations were added to WWTP or MS4 facilities # Idea is to read "Tables 8.4a-h MeHg Inputs Per Subarea.xlsx" and concatenate sheets into one dataframe # how to: https://readxl.tidyverse.org/articles/readxl-workflows.html # Will likely need to clean up spreadsheet and condense multiple header rows into one row & Add "NA" for missing values # Code will need to specify range to read for each sheet sheets <- file("Tables 8.4a-h MeHg Inputs Per Subarea.xlsx") %>% excel_sheets() %>% set_names() sheets <- grep("summary|north|south", sheets, ignore.case=T, value=T, invert=T) # do not select summary sheets or Yolo ByPass North & South # specify Cell Ranges for sheets # Central Delta, Marsh Creek, MokelumneCosumnes Rivers, Sacramento River, San Joaquin River, West Delta, & Yolo Bypass - Combined ranges <- list("A1:O20", "A1:O10", "A1:O10", "A1:O20", "A1:O29", "A1:O12", "A1:O24") # Concatenate sheets allocations_npdes <- map2_df( sheets, ranges, ~ read_xlsx(file("Tables 8.4a-h MeHg Inputs Per Subarea.xlsx"), sheet = .x, range = .y, na=c("NA", "depends on tributary", "depends on discharge location", "depends on location", "Not in 2010 TMDL Staff Report")), .id = "Subarea" ) %>% clean_names() %>% filter(grepl("wwtp|ms4", me_hg_source, ignore.case=T)) %>% write_csv(file("output/WWTP & MS4 Allocation Complilation.csv")) ``` ### Add Future Growth Allocations to assigned Facilities ```{r, echo=T} # Load Formatted Data allocations_npdes <- read_csv(file("output/WWTP & MS4 Allocation Complilation.csv")) %>% mutate(subarea_graph = case_when(grepl("Central", subarea) ~ "Central\nDelta", grepl("Marsh", subarea) ~ "Marsh\nCreek", grepl("Moke", subarea) ~ "Moke/Cos\nRivers", grepl("Sac", subarea) ~ "Sacramento\nRiver", grepl("San", subarea) ~ "San Joaquin\nRiver", grepl("West", subarea) ~ "West\nDelta", grepl("Yolo", subarea) ~ "YB\nCombined", T ~ subarea), subarea_table = case_when(grepl("Central", subarea) ~ "Central Delta", grepl("Marsh", subarea) ~ "Marsh Creek", grepl("Moke", subarea) ~ "Moke/Cos Rivers", grepl("Sac", subarea) ~ "Sacramento River", grepl("San", subarea) ~ "San Joaquin River", grepl("West", subarea) ~ "West Delta", grepl("Yolo", subarea) ~ "YB Combined", T ~ subarea)) %>% arrange(me_hg_source, subarea, specific_source) allocations_npdes %>% kbl(align='lcccccl', caption='WWTP & MS4 Individual WLAs') %>% kable_paper('hover', full_width=T) %>% kable_styling(fixed_thead = T) allocations_facility_sum <- allocations_npdes %>% mutate(specific_source = sub(" \\(allowable future growth)", "" , specific_source)) %>% # remove " (allowable future growth)" to add facility loads together group_by(me_hg_source, subarea_graph, subarea_table, specific_source) %>% summarise(across(ends_with("_g_yr"), ~sum(.x, na.rm=T)), .groups="drop") %>% replace (.==0, NA_real_) %>% # Replace 0's with NA so 0 values aren't graphed mutate(dmcp_review_compliance = ifelse(dmcp_review_median_annual_me_hg_load_g_yr <= dmcp_review_me_hg_allocation_g_yr, 1, 0), tmdl_2010_compliance = ifelse(tmdl_2010_staff_report_average_annual_me_hg_load_g_yr <= tmdl_2010_staff_report_me_hg_allocation_g_yr, 1, 0), dmcp_review_average_2010reduction = ifelse(!is.na(tmdl_2010_staff_report_average_annual_me_hg_load_g_yr), ifelse(dmcp_review_average_annual_me_hg_load_g_yr <= tmdl_2010_staff_report_average_annual_me_hg_load_g_yr, 1, 0), NA_real_), dmcp_review_average_2010compliance = ifelse(!is.na(tmdl_2010_staff_report_me_hg_allocation_g_yr), ifelse(dmcp_review_average_annual_me_hg_load_g_yr <= tmdl_2010_staff_report_me_hg_allocation_g_yr, 1, 0), NA_real_), dmcp_review_average_2010either = ifelse(dmcp_review_average_2010reduction == 1 | dmcp_review_average_2010compliance == 1, 1, 0)) %>% mutate(#WWTP Name Changes # specific_source=sub(" WWTP", "", specific_source), specific_source=sub("City of | Vocational Institution| CSD", "", specific_source), specific_source=sub("Regional County", "Reg Co", specific_source), specific_source=sub("Sacramento", "Sac", specific_source), #MS4 Name Changes specific_source=sub("San Joaquin", "SJ", specific_source), specific_source=sub("Community Services District", "CSD", specific_source), specific_source=sub("Sac Stormwater Quality Partnership", "SSQP", specific_source)) writexl::write_xlsx(allocations_facility_sum, path=file('output/allocations_facility_summary.xlsx')) allocations_facility_sum %>% kbl(align='lcccccl', caption='WWTP & MS4 Total WLAs (Summed with Future Growth)') %>% kable_paper('hover', full_width=T) %>% kable_styling(fixed_thead = T) ``` ### DMCP Review & 2010 Allocation comparision ```{r} allocation_comparison <- allocations_facility_sum %>% select(specific_source, me_hg_source, subarea_table, tmdl_2010_staff_report_average_annual_me_hg_load_g_yr, tmdl_2010_staff_report_me_hg_allocation_g_yr, dmcp_review_median_annual_me_hg_load_g_yr, dmcp_review_me_hg_allocation_g_yr ) %>% rename(Facility = specific_source, `Facility Type` = me_hg_source, Subarea = subarea_table, `2010 TMDL Avg Annual Load (g/yr)` = tmdl_2010_staff_report_average_annual_me_hg_load_g_yr, `2010 TMDL Allocation (g/yr)` = tmdl_2010_staff_report_me_hg_allocation_g_yr, `DMCP Review Median Annual Load (g/yr)` = dmcp_review_median_annual_me_hg_load_g_yr, `DMCP Review Allocation (g/yr)` = dmcp_review_me_hg_allocation_g_yr ) %>% # filter(!is.na(`2010 TMDL Allocation (g/yr)`)) %>% mutate(`2010 Load Reduction (g/yr)`= ifelse( `2010 TMDL Avg Annual Load (g/yr)` - `2010 TMDL Allocation (g/yr)` > 0, `2010 TMDL Avg Annual Load (g/yr)` - `2010 TMDL Allocation (g/yr)`, 0)) %>% mutate(`DMCP Review Load Reduction (g/yr)` = ifelse(`DMCP Review Median Annual Load (g/yr)` - `DMCP Review Allocation (g/yr)` > 0, `DMCP Review Median Annual Load (g/yr)` - `DMCP Review Allocation (g/yr)`, 0)) %>% # mutate(`DMCP Review Load Reduction (g/yr)` = ifelse(!is.na(`2010 TMDL Allocation (g/yr)`), ifelse(`DMCP Review Median Annual Load (g/yr)` - `DMCP Review Allocation (g/yr)` > 0, `DMCP Review Median Annual Load (g/yr)` - `DMCP Review Allocation (g/yr)`, 0), NA_real_)) %>% arrange(desc(`Facility Type`), `DMCP Review Load Reduction (g/yr)`, Facility) %>% add_row(Facility = "Cumulative Total", summarize(., across(contains("(g/yr)"), ~sum(.x, na.rm=T)))) %>% mutate(`2010 TMDL % Reduction Needed of Load to Meet Allocation` = round((`2010 TMDL Avg Annual Load (g/yr)` - `2010 TMDL Allocation (g/yr)`)/`2010 TMDL Avg Annual Load (g/yr)`, 2), `DMCP Review % Reduction Needed of Load to Meet Allocation` = round((`DMCP Review Median Annual Load (g/yr)` - `DMCP Review Allocation (g/yr)`)/`DMCP Review Median Annual Load (g/yr)`, 2)) writexl::write_xlsx(allocation_comparison, path=file('output/allocation_comparison.xlsx')) allocation_comparison %>% kbl(align='lcccccl', caption='WWTP & MS4 Total WLAs (Summed with Future Growth)') %>% kable_paper('hover', full_width=T) %>% kable_styling(fixed_thead = T) ``` ## Setup Graph Theme ```{r} colors_graphs <- RColorBrewer::brewer.pal(n = 11, name = 'RdYlBu') # "#A50026" "#D73027" "#F46D43" "#FDAE61" "#FEE090" "#FFFFBF" "#E0F3F8" "#ABD9E9" "#74ADD1" "#4575B4" "#313695" used_colors <- c("DMCP Review\nProposed WLA" = colors_graphs[9], "DMCP Review\nMedian Load"=colors_graphs[3], "2010 TMDL\nWLA"=colors_graphs[9], "2010 TMDL\nAvg Load"=colors_graphs[11], "DMCP Review\nAvg Load"=colors_graphs[3]) base_theme <- theme_set(theme_light() + theme(text = element_text(size=14, family="sans"), axis.title=element_text(size=12, face="bold"), axis.title.x = element_blank(), axis.text.x = element_text(size=10.5, angle = 45, hjust = 1), panel.grid.major.y=element_line(size=1), panel.grid.minor=element_line(size=0.5), legend.key.height=unit(1.5, "cm"), legend.text = element_text(size = 14), legend.title = element_blank() ) ) ``` ## Idea to Replace Long labels with Alphabetic ```{r} # https://stackoverflow.com/questions/30614105/create-abbreviated-legends-manually-for-long-x-labels-in-ggplot2 ``` ## DMCP Review Load & Allocation Comparison ```{r} dmcp_allocations <- allocations_facility_sum %>% select(me_hg_source:dmcp_review_me_hg_allocation_g_yr) %>% rename(`DMCP Review\nMedian Load` = dmcp_review_median_annual_me_hg_load_g_yr, `DMCP Review\nProposed WLA` = dmcp_review_me_hg_allocation_g_yr) %>% pivot_longer(cols=c("DMCP Review\nMedian Load", "DMCP Review\nProposed WLA"), names_to="DMCP Review", values_to="Methylmercury Load (g/yr)") ``` ### WWTP Compliance ```{r} allocation_wwtp <- dmcp_allocations %>% filter(grepl("wwtp", me_hg_source, ignore.case=T)) allocation_wwtp %>% kbl(align='lcccccl', caption='WWTP Total WLAs (Summed with Future Growth)') %>% kable_paper('hover', full_width=T) %>% kable_styling(fixed_thead = T) # wwtp_wla_plot <- allocation_wwtp %>% # ggplot() + # scale_fill_manual(values = used_colors[unique(dmcp_allocations$`DMCP Review`)]) zoom_value <- 0.65 # Graph All wwtp_wla_all <- allocation_wwtp %>% ggplot() + scale_fill_manual(values = used_colors[unique(dmcp_allocations$`DMCP Review`)]) + # Gray Rectangle Annotation to highlight Zoom Region annotate(geom="rect", xmin=-Inf, xmax=Inf, ymin=0, ymax=zoom_value, fill=alpha("grey", 0.3)) + # Full Scale of all data points geom_point(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=`DMCP Review`), size=3, colour="black", pch=21) + # Axis padding scale_y_continuous(expand = expansion(mult = c(0.01, 0.03))) + facet_grid(~ subarea_graph, scales="free", space="free") + # add minor grid lines and adjust space between plots theme(axis.text.x = element_text(angle = 59), panel.border = theme_border(type=c("right", "left"), colour="grey52", size = 0.6), panel.spacing.x = unit(0, "lines")) # Graph Zoom of WWTPs less than < 1 g/yr wwtp_wla_zoom <- allocation_wwtp %>% ggplot() + scale_fill_manual(values = used_colors[unique(dmcp_allocations$`DMCP Review`)]) + # Gray Rectangle Annotation to highlight Zoom Region annotate(geom="rect", xmin=-Inf, xmax=Inf, ymin=0, ymax=zoom_value, fill=alpha("grey88", 0.3)) + # Zoom in on data points less than 1 g/yr # geom_point(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=`DMCP Review`), # size=3, colour="black", pch=21) + geom_beeswarm(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=`DMCP Review`), size=3, colour="black", pch=21, cex=3, method="swarm") + # Axis padding scale_y_continuous(limits=c(0,zoom_value), expand = expansion(mult = c(0.01, 0))) + facet_grid(~ subarea_graph, scales="free", space="free") + # add minor grid lines and add space between plots theme(axis.text.x = element_text(angle = 59), panel.border = theme_border(type=c("right", "left"), colour="grey52", size = 0.6), panel.spacing.x = unit(0, "lines")) # View Plots wwtp_wla_all wwtp_wla_zoom # Save Graphs ggsave(filename = file("output/Fig 8.1a - DMCP Review Load & Allocation Comparison_WWTP_full.png"), plot = wwtp_wla_all, width = 9, height = 5.75, units = "in", dpi = 300) ggsave(filename = file("output/Fig 8.1b - DMCP Review Load & Allocation Comparison_WWTP.zoom.png"), plot = wwtp_wla_zoom, width = 9, height = 5.75, units = "in", dpi = 300) ``` ### MS4 Compliance Graph ```{r} allocation_ms4 <- dmcp_allocations %>% filter(grepl("ms4", me_hg_source, ignore.case=T)) allocation_ms4 %>% kbl(align='lcccccl', caption='WWTP Total WLAs (Summed with Future Growth)') %>% kable_paper('hover', full_width=T) %>% kable_styling(fixed_thead = T) allocation_ms4 %>% ggplot() + # geom_point(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=`DMCP Review`), # size=3, colour="black", pch=21) + geom_beeswarm(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=`DMCP Review`), size=3, colour="black", pch=21, cex=3.5, method="swarm") + scale_fill_manual(values = used_colors[unique(dmcp_allocations$`DMCP Review`)]) + scale_y_continuous(expand = expansion(mult = c(0.01, 0.025))) + facet_grid(. ~ subarea_graph, scales = "free", space = "free") + force_panelsizes(cols = c(.45, .25, .45, .5, 0.9, .35, .5)) + theme(axis.text.x = element_text(angle = 52), panel.border = theme_border(type=c("right", "left"), colour="grey52", size = 0.6), panel.spacing.x = unit(0, "lines")) ggsave(filename = file("output/Fig 8.3 - DMCP Review Load & Allocation Comparison_MS4.png"), # plot = predictionGraph, width = 9, height = 5.75, units = "in", dpi = 300) ``` ## Median Load & Allocation Compared to Mean Load & Allocation in 2010 TMDL Staff Report ```{r} dmcp_2010_comparison <- allocations_facility_sum %>% filter(!is.na(tmdl_2010_staff_report_me_hg_allocation_g_yr)) %>% # filter for facilities that were in 2010 TMDL select(me_hg_source:specific_source, dmcp_review_median_annual_me_hg_load_g_yr, dmcp_review_me_hg_allocation_g_yr, tmdl_2010_staff_report_average_annual_me_hg_load_g_yr, tmdl_2010_staff_report_me_hg_allocation_g_yr) %>% rename(`DMCP Review\nMedian Load` = dmcp_review_median_annual_me_hg_load_g_yr, `DMCP Review\nWLA` = dmcp_review_me_hg_allocation_g_yr, `2010 TMDL\nAvg Load` = tmdl_2010_staff_report_average_annual_me_hg_load_g_yr, `2010 TMDL\nWLA` = tmdl_2010_staff_report_me_hg_allocation_g_yr) %>% pivot_longer(cols=c("DMCP Review\nMedian Load", "DMCP Review\nWLA", "2010 TMDL\nAvg Load", "2010 TMDL\nWLA"), names_to="Comparison", values_to="Methylmercury Load (g/yr)") ``` ### WWTP Comparison ```{r} compare_wwtp <- dmcp_2010_comparison %>% filter(grepl("wwtp", me_hg_source, ignore.case=T)) compare_wwtp_plot <- compare_wwtp %>% ggplot() + scale_fill_manual(values = used_colors[unique(compare_wwtp$`Comparison`)]) zoom_value <- 2.7 # Graph All compare_wwtp_all <- compare_wwtp_plot + # Gray Rectangle Annotation to highlight Zoom Region annotate(geom="rect", xmin=-Inf, xmax=Inf, ymin=0, ymax=zoom_value, fill=alpha("grey", 0.3)) + # Full Scale of all data points geom_point(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=`Comparison`), size=3, colour="black", pch=21) + # Axis padding scale_y_continuous(expand = expansion(mult = c(0.01, 0.03))) + facet_grid(~ subarea_graph, scales="free", space="free") + # add minor grid lines and adjust space between plots theme(axis.text.x = element_text(angle = 56), panel.border = theme_border(type=c("right", "left"), colour="grey52", size = 0.6), panel.spacing.x = unit(0, "lines")) # Graph Zoom of WWTPs less than < 1 g/yr compare_wwtp_zoom <- compare_wwtp_plot + # Gray Rectangle Annotation to highlight Zoom Region annotate(geom="rect", xmin=-Inf, xmax=Inf, ymin=0, ymax=zoom_value, fill=alpha("grey88", 0.3)) + # Zoom in on data points less than 1 g/yr # geom_point(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=`Comparison`), # size=3, colour="black", pch=21) + geom_beeswarm(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=Comparison), size=3, colour="black", pch=21, cex=3.4, method="swarm") + # Axis padding scale_y_continuous(limits=c(0,zoom_value), expand = expansion(mult = c(0.01, 0))) + facet_grid(~ subarea_graph, scales="free", space="free") + # add minor grid lines and add space between plots theme(axis.text.x = element_text(angle = 56), axis.title.y = element_text(margin = margin(t = 0, r = 15, b = 0, l = 0)), # adds space after y-axis title so plot gridlines line up with compare_wwtp_all panel.border = theme_border(type=c("right", "left"), colour="grey52", size = 0.6), panel.spacing.x = unit(0, "lines")) # View Plots compare_wwtp_all compare_wwtp_zoom # Save Graphs ggsave(filename = file("output/Fig 8.Xa - Compare DMCP Review Median to 2010 TMDL_WWTP_full.png"), plot = compare_wwtp_all, width = 9, height = 5.75, units = "in", dpi = 300) ggsave(filename = file("output/Fig 8.Xb - Compare DMCP Review Median to 2010 TMDL_WWTP.zoom.png"), plot = compare_wwtp_zoom, width = 9, height = 5.75, units = "in", dpi = 300) ``` ### MS4 Comparison ```{r} compare_ms4 <- dmcp_2010_comparison %>% filter(grepl("ms4", me_hg_source, ignore.case=T)) # Look into geom_dotplot if overplotting is an issue compare_ms4 %>% ggplot() + # geom_point(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=Comparison), # size=3, colour="black", pch=21) + geom_beeswarm(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=Comparison), size=3, colour="black", pch=21, cex=3.4, method="swarm") + scale_fill_manual(values = used_colors[unique(dmcp_2010_comparison$Comparison)]) + scale_y_continuous(expand = expansion(mult = c(0.01, 0.025))) + facet_grid(. ~ subarea_graph, scales = "free", space = "free") + force_panelsizes(cols = c(.5, .2, .35, .6, .25, .25)) + theme(axis.text.x = element_text(angle = 61), panel.border = theme_border(type=c("right", "left"), colour="grey52", size=0.6), panel.spacing.x = unit(0, "lines")) ggsave(filename = file("output/Fig 8.X - Compare DMCP Review Median to 2010 TMDL_MS4.png"), # plot = predictionGraph, width = 9, height = 5.75, units = "in", dpi = 300) ``` ## Mean Load Compared to Mean Load & Allocation in 2010 TMDL Staff Report ```{r} dmcp_2010_comparison <- allocations_facility_sum %>% filter(!is.na(tmdl_2010_staff_report_me_hg_allocation_g_yr)) %>% # filter out facilities that were in 2010 TMDL select(me_hg_source:specific_source, dmcp_review_average_annual_me_hg_load_g_yr:tmdl_2010_staff_report_me_hg_allocation_g_yr) %>% rename(`DMCP Review\nAvg Load` = dmcp_review_average_annual_me_hg_load_g_yr, `2010 TMDL\nAvg Load` = tmdl_2010_staff_report_average_annual_me_hg_load_g_yr, `2010 TMDL\nWLA` = tmdl_2010_staff_report_me_hg_allocation_g_yr) %>% pivot_longer(cols=c("DMCP Review\nAvg Load", "2010 TMDL\nAvg Load", "2010 TMDL\nWLA"), names_to="Comparison", values_to="Methylmercury Load (g/yr)") ``` ### WWTP Comparison ```{r} compare_wwtp <- dmcp_2010_comparison %>% filter(grepl("wwtp", me_hg_source, ignore.case=T)) compare_wwtp_plot <- compare_wwtp %>% ggplot() + scale_fill_manual(values = used_colors[unique(compare_wwtp$`Comparison`)]) zoom_value <- 2.7 # Graph All compare_wwtp_all <- compare_wwtp_plot + # Gray Rectangle Annotation to highlight Zoom Region annotate(geom="rect", xmin=-Inf, xmax=Inf, ymin=0, ymax=zoom_value, fill=alpha("grey", 0.3)) + # Full Scale of all data points geom_point(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=`Comparison`), size=3, colour="black", pch=21) + # Axis padding scale_y_continuous(expand = expansion(mult = c(0.01, 0.03))) + facet_grid(~ subarea_graph, scales="free", space="free") + # add minor grid lines and adjust space between plots theme(axis.text.x = element_text(angle = 56), panel.border = theme_border(type=c("right", "left"), colour="grey52", size = 0.6), panel.spacing.x = unit(0, "lines")) # Graph Zoom of WWTPs less than < 1 g/yr compare_wwtp_zoom <- compare_wwtp_plot + # Gray Rectangle Annotation to highlight Zoom Region annotate(geom="rect", xmin=-Inf, xmax=Inf, ymin=0, ymax=zoom_value, fill=alpha("grey88", 0.3)) + # Zoom in on data points less than 1 g/yr # geom_point(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=`Comparison`), # size=3, colour="black", pch=21) + geom_beeswarm(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=Comparison), size=3, colour="black", pch=21, cex=3.4, method="swarm") + # Axis padding scale_y_continuous(limits=c(0,zoom_value), expand = expansion(mult = c(0.01, 0))) + facet_grid(~ subarea_graph, scales="free", space="free") + # add minor grid lines and add space between plots theme(axis.text.x = element_text(angle = 56), axis.title.y = element_text(margin = margin(t = 0, r = 15, b = 0, l = 0)), # adds space after y-axis title so plot gridlines line up with compare_wwtp_all panel.border = theme_border(type=c("right", "left"), colour="grey52", size = 0.6), panel.spacing.x = unit(0, "lines")) # View Plots compare_wwtp_all compare_wwtp_zoom # Save Graphs ggsave(filename = file("output/Fig 8.2a - Compare DMCP Review Avg to 2010 TMDL_WWTP_full.png"), plot = compare_wwtp_all, width = 9, height = 5.75, units = "in", dpi = 300) ggsave(filename = file("output/Fig 8.2b - Compare DMCP Review Avg to 2010 TMDL_WWTP.zoom.png"), plot = compare_wwtp_zoom, width = 9, height = 5.75, units = "in", dpi = 300) ``` ### MS4 Comparison ```{r} compare_ms4 <- dmcp_2010_comparison %>% filter(grepl("ms4", me_hg_source, ignore.case=T)) # Look into geom_dotplot if overplotting is an issue compare_ms4 %>% ggplot() + # geom_point(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=Comparison), # size=3, colour="black", pch=21) + geom_beeswarm(aes(x=specific_source, y=`Methylmercury Load (g/yr)`, fill=Comparison), size=3, colour="black", pch=21, cex=3.4, method="swarm") + scale_fill_manual(values = used_colors[unique(dmcp_2010_comparison$Comparison)]) + scale_y_continuous(expand = expansion(mult = c(0.01, 0.025))) + facet_grid(. ~ subarea_graph, scales = "free", space = "free") + force_panelsizes(cols = c(.5, .2, .35, .6, .25, .25)) + theme(axis.text.x = element_text(angle = 61), panel.border = theme_border(type=c("right", "left"), colour="grey52", size=0.6), panel.spacing.x = unit(0, "lines")) ggsave(filename = file("output/Fig 8.4 - Compare DMCP Review Avg to 2010 TMDL_MS4.png"), # plot = predictionGraph, width = 9, height = 5.75, units = "in", dpi = 300) ```
## Error: <text>:17:1: unexpected '<' ## 16: --- ## 17: < ## ^
The R session information (including the OS info, R version and all packages used):
sessionInfo()
## R version 4.2.2 (2022-10-31 ucrt) ## Platform: x86_64-w64-mingw32/x64 (64-bit) ## Running under: Windows 10 x64 (build 22621) ## ## Matrix products: default ## ## locale: ## [1] LC_COLLATE=English_United States.utf8 LC_CTYPE=English_United States.utf8 ## [3] LC_MONETARY=English_United States.utf8 LC_NUMERIC=C ## [5] LC_TIME=English_United States.utf8 ## ## attached base packages: ## [1] stats graphics grDevices utils datasets methods base ## ## loaded via a namespace (and not attached): ## [1] compiler_4.2.2 magrittr_2.0.3 tools_4.2.2 stringi_1.7.6 highr_0.9 ## [6] knitr_1.39 stringr_1.4.0 xfun_0.31 evaluate_0.15
Sys.time()
## [1] "2023-12-26 14:14:11 PST"