Concentration vs Emission: What’s the difference anways?

Objective

What is the difference between the emission driven and concentration driven runs? In theory we expect the emission driven runs with the carbon cycle turned on to have feedback mechanimss impacting the temperature and co2 otuput. Corinne brought up the idea, does it even really matter? Is the difference 5 ppm or 100 ppm? Are the ESMs and Hector senstive enough to the carbon cycle for the emission driven runs to be different from the concentration driven runs? Could the climate paramterizations do an adequate job of representing the ESMs?

1. Is there a difference in the emission vs concentration driven runs for the ESMs and Hector? / How well does the concenration climate paramterizations do with emulating the carbon cycle?
2. Solving for beta at fixed q10 values.

1. Is there a difference in the emission vs concentration driven runs for the ESMs and Hector?

# What models have esmission driven results? 
cmip_individual %>%  
    filter(grepl('esm', experiment)) %>% 
    pull(model) %>% 
    unique() -> 
    emission_driven_models
# Import the co2 concentration for the concentration driven runs (from the hector inputs). 
# Format the data frame so that it has the model name and ensemble information.
read.csv(system.file('input/constraints/rcp85_co2ppm.csv', package = 'hector'), skip = 1) %>% 
    rename(year = Date, value = Ca_constrain) %>% 
    mutate(experiment = if_else(year <= 2005, 'historical', 'rcp85'), 
           variable = 'co2') %>% 
    mutate(join = 1) %>% 
    left_join(tibble(join = 1, 
                     ensemble = c("r1i1p1", "r2i1p1", "r3i1p1")) %>% 
                         left_join(tibble(join = 1, 
                                   model = emission_driven_models),  by = 'join'), 
              by = 'join') %>% 
    select(-join) -> 
    concentration_co2_values
# Import the CMIP ESM data and make sure that the concentration results does not include co2. 
cmip_individual %>% 
    filter(model %in% emission_driven_models) %>% 
    filter(year <= 2100 & experiment %in% c('historical', 'rcp45')) %>% 
    filter(variable != 'co2') %>% 
    bind_rows(concentration_co2_values) -> 
    concenctraion_rlst
# Import the CMIP ESM emisssion driven results. 
cmip_individual %>% 
    filter(model %in% emission_driven_models) %>% 
    filter(year <= 2100 & experiment %in% c('esmHistorical', 'esmrcp45')) %>% 
    bind_rows(concenctraion_rlst) %>% 
     mutate(driven = if_else(grepl('esm', experiment), 'emissions', 'concentration')) %>%  
    select(year, model, ensemble, variable, value, driven) %>%  
    distinct() -> 
    cmip_conc_emiss
# Calculate the difference between the emission driven anf the concentration driven runs. 
cmip_conc_emiss %>% 
    tidyr::spread(driven, value) %>% 
    na.omit %>% 
    mutate(dif = emissions - concentration) -> 
    diff_emis_conc

How do the ESM emission driven and concentration runs differ from one another?

diff_emis_conc %>%  
    ggplot(aes(concentration, emissions, color = model)) + 
    geom_point() + 
    geom_abline(intercept = 0, slope = 1, col = 'black') +
    facet_wrap('variable', scales = 'free') + 
    labs(y = 'emission driven ESM', 
         x = 'concentration driven ESM')

The black line show the 1:1 relationship between the x and y axis, so if the emission driven and concentration driven runs are fairly similar to one another then we would expect the pattern to be centered around the black line like we see in heatflux and tas. But it does look like the emission driven runs have a higher co2 concentration fo at least 2 of the models. CanESM2 does not seem to be particuarly senstive to the carbon cycle feedbacks.

What is the difference bewteen the runs? What does the noise look like of the conc and emission driven runs look like?

diff_emis_conc %>% 
    group_by(variable, model) %>% 
    summarise(min_dif = min(dif), 
              mean_dif = mean(dif), 
              max_dif = max(dif), 
              sd_dif = sd(dif)) %>% 
    ungroup()

How well does the Emission Hector calibrated with climate best fits models the emission ESM output?

How does the temp and co2 data for Hector calibrated with the results from the bestFit_conc_temp_heat2100 compare with the emission driven results?

paths     <- list.files(DIR, full.names = TRUE) 
to_import <- paths[which(grepl(paste(emission_driven_models, collapse = '|'), paths))]
lapply(to_import, function(input){
    object <- readRDS(input)
    
    data.frame(matrix(object$optim_rslt$par, nrow = 1, dimnames = list(NULL, names(object$optim_rslt$par)))) %>%  
        mutate(model = gsub(pattern = 'conc_|.rds', replacement = '', x = basename(input)))
}) %>% 
    bind_rows -> 
    params_df

Compare Hector emissions and ESM emission driven run results. How well does the concenration climate paramterizations do with emulating the carbon cycle?

path <- system.file('input/hector_rcp85.ini', package = 'hector')
core <- newcore(path)
split(x = params_df, f = params_df$model) %>% 
    lapply(function(input){params <- input[which(names(input) != 'model')]
parameterize_core(params = params, core = core)
reset(core)
run(core, runtodate = 2100)
fetchvars(core, 1850:2100, vars = c(HEAT_FLUX(), GLOBAL_TEMP(), ATMOSPHERIC_CO2())) %>% 
    mutate(experiment = if_else(year < 2005, 'esmHistorical', 'esmrcp85'), 
           model = input[['model']]) %>%  
    mutate(scenario = 'hector') %>% 
    mutate(variable = if_else(variable == 'Tgav', 'tas', variable)) %>% 
    mutate(variable = if_else(variable == 'Ca', 'co2', variable))
}) %>%  
    bind_rows() -> 
    emission_hector
emission_hector %>% 
    bind_rows(cmip_individual %>%  
                  mutate(scenario = 'cmip')) %>% 
    filter(grepl('esm', experiment)) -> 
    emission_hector_output 
emission_hector_output %>% 
    ggplot(aes(year, value, color = scenario, linetype = model)) + 
    geom_line() + 
    facet_grid(variable ~ model, scales = 'free') 

How much does the carbon cycle impact Hector output?

bind_rows(emission_hector_output %>%  mutate(run = 'emission'), 
          concentration_hector_output %>% mutate(run = 'concentration')) %>% 
    filter(grepl('hector', scenario)) -> 
    to_plot 
to_plot %>% 
    ggplot(aes(year, value, color = run, group = experiment)) + 
    geom_line() + 
    facet_grid(variable~model, scale = 'free')

to_plot %>% 
    select(year, variable, value, model, run) %>% 
    spread(run, value) %>%  
    ggplot(aes(concentration, emission, color = model)) +
    geom_point() + 
    geom_abline(intercept = 0, slope = 1) +
    facet_wrap('variable', scales = 'free') + 
      labs(y = 'emission driven ESM', 
         x = 'concentration driven ESM')

I guess I do not know if this sort of difference really seems like it is that notable to me. What do you guys think?

2. Solving for beta at fixed q10 values.

list.files(DIR2, '.rds', full.names = TRUE) %>% 
    lapply(function(input){readRDS(input)[['hector_output']]}) %>% 
    bind_rows() %>% 
    mutate(beta = signif(beta, digits = 3), 
           q10 = signif(q10, digits = 3), 
           beta_q10 = paste0(beta, '_', q10)) %>% 
    filter(scenario == 'esmrcp85') %>%  
    mutate(variable = if_else(variable == 'Tgav', 'tas', variable)) %>%  
    mutate(variable = if_else(variable == 'Ca', 'co2', variable))-> 
    to_plot

How much do the beta and q10 values range?

summary(to_plot$beta)

    Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
0.000000 0.003237 0.057250 0.061660 0.108500 0.165000 

The beta values that optim solved for do not span the range we see in the DOE PI meeting mcmc PDFs… which I think it curious, very curios.

summary(to_plot$q10)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
  1.000   1.998   3.000   3.000   4.003   5.000 

ggplot() + 
    geom_line(data = cmip_individual %>% filter(grepl('esm', experiment) & model == 'CESM1-BGC') %>%  
                  mutate(beta_q10 = 'cmip'),  
              aes(year, value, color = beta_q10), 
              color = 'grey', size = 2) +
    geom_line(data = to_plot, aes(year, value, color = beta_q10), alpha = 0.5) +
    facet_wrap('variable', scale = 'free') + 
    theme_bw(base_size = 18) + 
        theme(legend.position = 'none') 

We see a large range in the potential q10 and beta values, and they results are fairly similar to one another, they fall within the range of the cmip data especially for the tas data. For the co2 data it is some what hard to see but the solved hector runs with the various beta and q10 values do fall on top of the esm cmip data. Which confirms our suspicion that beta and q10 are capable of trading off with one another.