AR6 Climate assessment on reporting

[gabriel-abrahao]

Purpose of this PR

This PR adds the option to run AR6's climate-assessment as part of the normal output.R procedure, under the option ar6Climate. If one wants it to be run automatically, it can be added to cfg$output. This includes the probabilistic assessment of temperature change, compatible with the AR6 standards.

By default, the whole 600 member MAGICC v7.5.3 ensemble is run, which can take a while. Only certain percentiles are output, with each percentile in a different variable. Although the original output is yearly, which can be found inside a subfolder climate-temp, only data for the REMIND timesteps are added to the standard MIF (REMIND_generic***). All new variables have the AR6 climate diagnostics| prefix. The original temperature variable, Temperature|Global Mean, is still calculated with the old (v6) MAGICC. The comparable variable in the new assessment is AR6 climate diagnostics|Surface Temperature (GSAT)|MAGICCv7.5.3|XXth percentile, with separate variables for 10 different percentiles.

The subfolder climate-temp includes the full output of climate-assessment, including files describing temperature exceedence probabilities and the consequent climate classification of the scenario (C1, C2, etc).

Type of change

• New feature

Checklist:

• My code follows the coding etiquette
• I performed a self-review of my own code
• I explained my changes within the PR, particularly in hard-to-understand areas
• I checked that the in-code documentation is up-to-date
• I adjusted the reporting in remind2 where it was needed
• I adjusted forbiddenColumnNames in readCheckScenarioConfig.R in case the PR leads to deprecated switches
• All automated model tests pass (FAIL 0 in the output of make test)
• The changelog CHANGELOG.md has been updated correctly

[gabriel-abrahao]

There are still two big pain points that should probably be addressed before this is merged. Could someone familiar with Python and our venv setup take a look at this please? Maybe @pfuehrlich-pik ?

1. The script makes use of scripts in climate-assessment. However, the scripts folder does not seem to be included when installing the package with pip. I've opened an issue at the package regarding this, but maybe there's a simpler solution. Either way, we need to assure that the version of the scripts and the version of the package used is the one installed in the venv. climate-assessment is already part of requirements.txt, but in my tests it seemed to be using the one in my own personal environment. Currenlty the script works, but I've hardcoded a scripts folder that definitely shouldn't be there.

2. The script is pretty slow as it is, but it can be substantially sped up if more processors are available simply by changing this line. However, how to determine and/or set up this within the current output.R setup? Maybe @LaviniaBaumstark or @fbenke-pik know?

[pfuehrlich-pik]

Hi Gabriel, I don't really know much about python venvs, Mika was our expert here. However, it seems like you are using run_harm_inf.py and run_clim.py which are both quite trivial scripts. They import climate_assessment.cli and then call a single function from there. This should also be possible from the pip-installed package (or is it called module in python?), because that includes cli.py (it seems like this subfolder is included for pip-installations). So maybe you don't actually need the scripts, but can import cli directly and call the functions yourself?

[orichters]

Nice work, Gabriel, thanks…

My first set of comments ;)

I have a preference for not calling it "ar6climate", because in AR6 different climate models were used, but rather MAGICC7.5.3 or something like that.

[fbenke-pik]

Concerning this line
Sys.setenv(MAGICC_WORKER_NUMBER=1)

You should be able get this number by acessing "SLURM_NTASKS" on the cluster. A failsafe solution could look like this:
workerNumber <- as.numeric(Sys.getenv("SLURM_NTASKS", unset = 1))

There might be other useful environment variables available. Please check out this tutorial: https://www.pik-potsdam.de/en/institute/about/it-services/hpc/user-guides/slurm

And currently, via output.R, only tasks = 1 can be selected.

But of course, if this is the case, then SLURM_NTASKS will always be 1.

[orichters]

@gabriel-abrahao: Merge this to improve slurm options in output.R: gabriel-abrahao#1

[orichters]

I noted some minor things that have changed since Gabriel started working on it and that should be kept.

[orichters]

Because I created a new merge conflict, I have fixed the stuff that I noted.

[orichters]

Some results for /p/tmp/oliverr/remind-smallfix/output/SSP2EU-NPi_2024-02-20_15.40.25, median temperatures (so 50.0th percentile)

   period MAGICC6 MAGICC7.5.3
    <int>   <dbl>       <dbl>
 1   2005   0.881       0.868
 2   2010   0.895       0.976
 3   2015   0.999       1.12
 4   2020   1.13        1.23
 5   2025   1.27        1.37
 6   2030   1.40        1.48
 7   2035   1.55        1.62
 8   2040   1.68        1.75
 9   2045   1.81        1.86
10   2050   1.93        2.00
11   2055   2.03        2.10
12   2060   2.13        2.22
13   2070   2.30        2.46
14   2080   2.45        2.66
15   2090   2.58        2.83
16   2100   2.68        2.97
17   2110   2.77       NA
18   2130   2.90       NA
19   2150   3.02       NA

[orichters]

Comparison for NGFS3, so 2022, where I have both MAGICC6 from REMIND and MAGICC7.5.3 as run by IIASA at hand:

   scenario         period MAGICC6 MAGICC7.5.3
   <fct>             <int>   <dbl>       <dbl>
 1 Current Policies   2010   0.953       0.976
 2 Current Policies   2020   1.16        1.24
 3 Current Policies   2030   1.40        1.49
 4 Current Policies   2040   1.65        1.76
 5 Current Policies   2050   1.92        2.04
 6 Current Policies   2060   2.18        2.29
 7 Current Policies   2070   2.41        2.54
 8 Current Policies   2080   2.63        2.77
 9 Current Policies   2090   2.84        2.98
10 Current Policies   2100   3.04        3.20
11 NDCs               2010   0.953       0.976
12 NDCs               2020   1.16        1.24
13 NDCs               2030   1.41        1.50
14 NDCs               2040   1.64        1.73
15 NDCs               2050   1.83        1.90
16 NDCs               2060   1.99        2.04
17 NDCs               2070   2.12        2.16
18 NDCs               2080   2.23        2.26
19 NDCs               2090   2.32        2.33
20 NDCs               2100   2.40        2.39
21 Net Zero 2050      2010   0.953       0.976
22 Net Zero 2050      2020   1.16        1.24
23 Net Zero 2050      2030   1.42        1.48
24 Net Zero 2050      2040   1.54        1.56
25 Net Zero 2050      2050   1.58        1.55
26 Net Zero 2050      2060   1.58        1.50
27 Net Zero 2050      2070   1.57        1.47
28 Net Zero 2050      2080   1.57        1.43
29 Net Zero 2050      2090   1.56        1.40
30 Net Zero 2050      2100   1.55        1.38

click to expand code!

library(quitte)
library(tidyr)
library(dplyr)

data <- as.quitte("/p/projects/remind/users/oliverr/data/NGFS3.4_IAM_data_R5.xlsx")
vars <- c("Temperature|Global Mean", "AR6 climate diagnostics|Surface Temperature (GSAT)|MAGICCv7.5.3|50.0th Percentile")
scens <- c("Net Zero 2050", "Current Policies", "NDCs")
t <- seq(2010, 2100, 10)

d <- data %>%
  filter(.data$variable %in% vars, .data$scenario %in% scens, .data$model %in% "REMIND NGFS3", .data$period %in% t) %>%
  mutate(variable = gsub(vars[[1]], "MAGICC6", .data$variable, fixed = TRUE)) %>%
  mutate(variable = gsub(vars[[2]], "MAGICC7.5.3", .data$variable, fixed = TRUE)) %>%
  mutate(unit = "K") %>%
  quitteSort() %>%
  df_variation() %>%
  pivot_wider(names_from = "variable") %>%
  droplevels() %>%
  print(n = 200)

[orichters]

I did another test using NGFS4, so 2023, data. I took a REMIND NGFS run from 2023, ran the ar6climate reporting script, and compare it to the AR6 assessment I have received from IIASA:

   period REMIND_MAGICC6 REMIND_MAGICC7.5.3 IIASA_MAGICC7.5.3
 1   2010          0.956              0.976             0.976
 2   2020          1.16               1.23              1.23
 3   2030          1.40               1.46              1.48
 4   2040          1.64               1.71              1.73
 5   2050          1.89               1.93              1.96
 6   2060          2.10               2.13              2.16
 7   2070          2.29               2.34              2.35
 8   2080          2.48               2.52              2.53
 9   2090          2.65               2.71              2.69
10   2100          2.81               2.87              2.84

There is some slight variation that I don't understand. But overall looks good.

click to expand code!

library(quitte)
library(tidyr)
library(dplyr)

files <- c("/p/tmp/oliverr/remind-smallfix/output/C_h_cpol-rem-6/REMIND_generic_C_h_cpol-rem-6.mif",
           "/p/projects/remind/users/oliverr/data/NGFS4_IAM_data_S5.xlsx")
data <- as.quitte(files)

models <- c("REMIND NGFS4 S5", "REMIND")
vars <- c("Temperature|Global Mean", "AR6 climate diagnostics|Surface Temperature (GSAT)|MAGICCv7.5.3|50.0th Percentile")
scens <- c("Current Policies", "C_h_cpol-rem-6")
t <- seq(2010, 2100, 10)

d <- data %>%
  filter(model %in% models) %>%
  filter(scenario %in% scens) %>%
  filter(variable %in% vars) %>%
  filter(period %in% t) %>%
  mutate(unit = factor("K")) %>%
  mutate(variable = gsub(vars[[1]], "MAGICC6", variable, fixed = TRUE)) %>%
  mutate(variable = gsub(vars[[2]], "MAGICC7.5.3", .data$variable, fixed = TRUE)) %>%
  mutate(model = factor(ifelse(model == models[[1]], "NGFS", "REMIND"))) %>%
  mutate(variable = paste(model, variable, sep = "_")) %>%
  mutate(model = factor(models[[2]])) %>%
  mutate(scenario = factor(scens[[1]])) %>%
  df_variation() %>%
  pivot_wider(names_from = variable)
print(d, n = 30)

[orichters]

Code is fine from my side. Still might be worth investigating why the results slightly differ, but this is a big step forward, thank you @tonnrueter and @gabriel-abrahao.

@LaviniaBaumstark: Can you have a look as well, please, mostly whether you consider the last table good enough in terms of results.

[LaviniaBaumstark]

I would see no reason for not merging. Differences in REMIND_MAGICC6 vs. REMIND_MAGICC7.5.3 is no suprice. Why REMIND_MAGICC7.5.3 vs. IIASA_MAGICC7.5.3 are not the same is maybe related to some inconcistend parameterisation but nothing we need to clarify before merging. We are clode enough.

[LaviniaBaumstark]

but: please run make test

[orichters]

but: please run make test

Please run make test but temporarily and locally out-comment this line before. Is pointless without it. And merge develop into this branch, if not, it will fail on memory issues.

[gabriel-abrahao]

I agree with Lavinia here, the differences in Oli's tests don't seem large enough to justify not merging. I've put a lot more work into making sure everything is correct in this version here than in the previous one, so if there's an actual mistake, it's more likely that the differences came from fixing it. That said, a lot of segfaults and overflows usually happen in some individual MAGICC parametrizations. The IIASA people told me they are aware of it but didn't really give me an explanation yet, so maybe it's related to that. At the end of the day, that's the official AR6 pipeline, so if there's a problem, everyone else also has it, and they are aware.

[tonnrueter]

Thanks y'all for chiming in! I'll run the tests and merge asap