Comparing Hector integrated vs non integrated DOECLIM

1 Set Up

library(data.table)
library(dplyr)
library(tidyr)
library(ggplot2)

BASE_DIR <- here::here()

theme_set(theme_bw())


data <- rbind(read.csv("~/Desktop/hector_doeclim.csv", stringsAsFactors = FALSE), 
              read.csv("~/Desktop/hector_land-tracking.csv", stringsAsFactors = FALSE)) %>%  
    dplyr::filter(grepl(pattern = "ssp", x = scenario))

1.1 Objective

How much does Hector output change by when DOECLIM is integrated vs. fully integrated.

2 Compare Hector Outputs

2.1 Temperature

data %>% 
    dplyr::filter(variable %in% c("Tgav", "Tgav_land", "Tgav_ocean_air", "Tgav_ocean_ST")) %>% 
    ggplot(aes(year, value, color = version, linetype = version, groupby = interaction(variable, version, scenario))) + 
    geom_line() + 
    facet_wrap("variable", scales = "free") + 
    theme_bw() + 
    theme(legend.position = "bottom") + 
    labs(y = NULL, x = NULL)

data %>% 
    dplyr::filter(variable %in% c("Temp_HL", "Temp_LL")) %>% 
    ggplot(aes(year, value, color = version, linetype = version, groupby = interaction(variable, version, scenario))) + 
    geom_line() + 
    facet_wrap("variable", scales = "free") + 
    theme_bw() + 
    theme(legend.position = "bottom") + 
    labs(y = NULL, x = NULL)

Not surprising that we see a large change here since we are swtiching from tas to tos data so the fact that the HL and LL temp is now lower is not surprising.

2.2 Carbon Fluxes

data %>% 
    dplyr::filter(variable %in% c("atm_ocean_flux", "atm_ocean_flux_LL", "atm_ocean_flux_HL")) %>% 
    ggplot(aes(year, value, color = version, linetype = version, groupby = interaction(variable, version, scenario))) + 
    geom_line() + 
    facet_wrap("variable", scales = "free") + 
    theme_bw() + 
    theme(legend.position = "bottom") + 
    labs(y = NULL, x = NULL)

2.3 Ocean Vars

data %>% 
    dplyr::filter(variable %in% c("pH_HL", "pH_LL", "CO3_HL", "CO3_LL", "PCO2_HL", "PCO2_LL")) %>% 
    ggplot(aes(year, value, color = version, linetype = version, groupby = interaction(variable, version, scenario))) + 
    geom_line() + 
    facet_wrap("variable", scales = "free") + 
    theme_bw() + 
    theme(legend.position = "bottom") + 
    labs(y = NULL, x = NULL)

data %>% 
    dplyr::filter(variable %in% c("carbon_HL", "carbon_LL", "carbon_IO", "carbon_DO")) %>% 
    ggplot(aes(year, value, color = version, linetype = version, groupby = interaction(variable, version, scenario))) + 
    geom_line() + 
    facet_wrap("variable", scales = "free") + 
    theme_bw() + 
    theme(legend.position = "bottom") + 
    labs(y = NULL, x = NULL)

2.4 Terrestiral Carbon

data %>% 
    dplyr::filter(variable %in% c("npp")) %>% 
    ggplot(aes(year, value, color = version, linetype = version, groupby = interaction(variable, version, scenario))) + 
    geom_line() + 
    facet_wrap("variable", scales = "free") + 
    theme_bw() + 
    theme(legend.position = "bottom") + 
    labs(y = NULL, x = NULL)

3 Percent Difference

data %>% 
    tidyr::spread(version, value) %>% 
    dplyr::mutate(percent = 100 * (`full doeclim integration` - `land-tracking`)/ `land-tracking`) -> 
    percent_change_data 

percent_change_data %>% 
    dplyr::group_by(variable, scenario) %>%  
    dplyr::summarise(value = mean(percent)) %>% 
    tidyr::spread(scenario, value) %>% 
    knitr::kable(digits = 3, caption = "(doeclim - land tracking) / land tracking")

`summarise()` has grouped output by 'variable'. You can override using the `.groups` argument.

(doeclim - land tracking) / land tracking

variable	ssp119	ssp126	ssp245	ssp370	ssp434	ssp460	ssp534-over	ssp585
atm_ocean_flux	5.901	5.427	4.598	4.839	5.593	5.063	5.426	4.806
atm_ocean_flux_HL	0.734	0.828	0.870	0.979	0.826	0.965	0.930	1.090
atm_ocean_flux_LL	6.040	6.892	3.754	3.785	-12.513	4.261	7.906	3.677
carbon_DO	0.005	0.005	0.004	0.006	0.005	0.006	0.006	0.007
carbon_HL	0.042	0.047	0.049	0.058	0.050	0.057	0.055	0.065
carbon_IO	0.003	0.004	0.003	0.004	0.004	0.004	0.004	0.004
carbon_LL	0.053	0.058	0.058	0.066	0.062	0.067	0.066	0.071
CO3_HL	-2.358	-2.470	-2.587	-2.807	-2.499	-2.707	-2.626	-2.999
CO3_LL	-4.242	-4.344	-4.462	-4.658	-4.376	-4.555	-4.482	-4.834
DIC_HL	0.042	0.047	0.049	0.058	0.050	0.057	0.055	0.065
DIC_LL	0.053	0.058	0.058	0.066	0.062	0.067	0.066	0.071
npp	0.057	0.060	0.062	0.063	0.062	0.063	0.066	0.064
PCO2_HL	0.020	0.028	0.046	0.047	0.004	0.044	0.041	0.053
PCO2_LL	0.076	0.088	0.112	0.124	0.084	0.111	0.109	0.138
pH_HL	-0.016	-0.017	-0.019	-0.021	-0.016	-0.020	-0.019	-0.023
pH_LL	-0.014	-0.015	-0.017	-0.019	-0.015	-0.017	-0.016	-0.021
Temp_HL	-21.335	-21.482	-21.563	-21.813	-21.549	-21.733	-21.617	-21.948
Temp_LL	-5.555	-5.638	-5.702	-5.864	-5.679	-5.802	-5.748	-5.988
Tgav	2.774	2.750	1.721	2.704	2.740	2.713	2.737	2.687
Tgav_land	-0.007	-0.033	-0.821	-0.082	-0.043	-0.072	-0.044	-0.099
Tgav_ocean_air	1.475	1.451	-0.211	1.407	1.442	1.416	1.438	1.390
Tgav_ocean_ST	1.475	1.451	-0.211	1.407	1.442	1.416	1.438	1.390

3.1 Temperature

percent_change_data %>%
    dplyr::filter(variable %in% c("Tgav", "Tgav_land", "Tgav_ocean_air", "Tgav_ocean_ST")) %>% 
    ggplot(aes(year, percent, color = scenario)) + 
    geom_hline(yintercept = 10, color = "grey", size = 1) + 
    geom_hline(yintercept = -10, color = "grey", size = 1) + 
    geom_line() + 
    coord_cartesian(ylim = c(-50, 50)) + 
    facet_wrap("variable") + 
    labs(caption = "grey lines at -10 & 10%") + 
    labs(x = NULL, y = "%")

percent_change_data %>%
    dplyr::filter(variable %in% c("Temp_HL", "Temp_LL")) %>% 
    ggplot(aes(year, percent, color = scenario)) + 
    geom_hline(yintercept = 10, color = "grey", size = 1) + 
    geom_hline(yintercept = -10, color = "grey", size = 1) + 
    geom_line() + 
    facet_wrap("variable") + 
    labs(caption = "grey lines at -10 & 10%") + 
    labs(x = NULL, y = "%")

3.2 Carbon Fluxes

percent_change_data %>%
    dplyr::filter(variable %in% c("atm_ocean_flux", "atm_ocean_flux_LL", "atm_ocean_flux_HL")) %>% 
    ggplot(aes(year, percent, color = scenario)) + 
    geom_hline(yintercept = 10, color = "grey", size = 1) + 
    geom_hline(yintercept = -10, color = "grey", size = 1) + 
    geom_line() + 
    facet_wrap("variable") + 
    coord_cartesian(ylim = c(-75, 75)) +
    labs(caption = "grey lines at -10 & 10%") + 
    labs(x = NULL, y = "%")

3.3 Ocean Vars

percent_change_data %>%
    dplyr::filter(variable %in% c("pH_HL", "pH_LL", "CO3_HL",
                                  "CO3_LL", "PCO2_HL", "PCO2_LL")) %>% 
    ggplot(aes(year, percent, color = scenario)) + 
    geom_hline(yintercept = 10, color = "grey", size = 1) + 
    geom_hline(yintercept = -10, color = "grey", size = 1) + 
    geom_line() + 
    facet_wrap("variable") + 
    labs(caption = "grey lines at -10 & 10%") + 
    labs(x = NULL, y = "%")

percent_change_data %>%
    dplyr::filter(variable %in% c("carbon_HL", "carbon_LL", "carbon_IO", "carbon_DO")) %>% 
    ggplot(aes(year, percent, color = scenario)) + 
    geom_hline(yintercept = 10, color = "grey", size = 1) + 
    geom_hline(yintercept = -10, color = "grey", size = 1) + 
    geom_line() + 
    facet_wrap("variable") + 
    labs(caption = "grey lines at -10 & 10%") + 
    labs(x = NULL, y = "%")

3.4 Terrestiral Carbon

percent_change_data %>%
    dplyr::filter(variable %in% c("npp")) %>% 
    ggplot(aes(year, percent, color = scenario)) + 
    geom_hline(yintercept = 10, color = "grey", size = 1) + 
    geom_hline(yintercept = -10, color = "grey", size = 1) + 
    geom_line() + 
    facet_wrap("variable") + 
    labs(caption = "grey lines at -10 & 10%") + 
    labs(x = NULL, y = "%")