Landing

Column

About

This digital tool provides interactive figures and data tables for the following paper

Humpenöder, F., Popp, A., Schleussner, C.-F., Orlov, A., Windisch, M.G., Menke, I., Pongratz, J., Havermann, F., Thiery, W., Luo, F., et al. (2022). Overcoming global inequality is critical for land-based mitigation in line with the Paris Agreement. Nat Commun 13, 7453. 10.1038/s41467-022-35114-7.

https://www.nature.com/articles/s41467-022-35114-7

Abstract

Transformation pathways for the land sector in line with the Paris Agreement depend on the assumption of globally implemented greenhouse gas (GHG) emission pricing, and in some cases also on inclusive socio-economic development and sustainable land-use practices. In such pathways, the majority of GHG emission reductions in the land system is expected to come from low- and middle-income countries, which currently account for a large share of emissions from agriculture, forestry and other land use (AFOLU). However, in low- and middle-income countries the economic, financial and institutional barriers for such transformative changes are high. Here, we show that if sustainable development in the land sector remained highly unequal and limited to high-income countries only, global AFOLU emissions would remain substantial throughout the 21st century. Our model-based projections highlight that overcoming global inequality is critical for land-based mitigation in line with the Paris Agreement. While also a scenario purely based on either global GHG emission pricing or on inclusive socio-economic development would achieve the stringent emissions reductions required, only the latter ensures major co-benefits for other Sustainable Development Goals, especially in low- and middle-income regions.

Further information

Contact for interactive figures: Florian Humpenöder ()

Interactive figures are realized with the R packages ggplot2, patchwork and ggiraph.

The website is generated from an Rmd file using the R flexdashboard package.

Column

Information on the LAMACLIMA project

LAMACLIMA – LAnd MAnagement for CLImate Mitigation and Adaptation

The project aims to investigate how changes in land cover and land management can help to meet the mitigation and adaptation objectives of the Paris Agreement, as well as the Sustainable Development Goals.

https://climateanalytics.org/projects/lamaclima/

Overview

Column

Figure 6

Column

Information

Figure 6: AFOLU GHG emissions over time for five scenarios. a) shows results at global level over time. b) shows results at regional level for 2020, 2050 and 2100. Black solid lines in a) and dots in b) show the net effect across AFOLU GHG emissions. Grey solid lines in a) show median values across 1.5°C compatible illustrative mitigation pathways (LD, Ren, SP) from the IPCC AR6 Scenarios Database24. CO2 includes emissions from land-use change such as deforestation and conversion of non-forest ecosystems, as well as removals from re/afforestation and natural succession on abandoned agricultural land. CH4 includes emissions from enteric fermentation, animal waste management and rice cultivation. N2O includes emissions from agricultural soils (fertilizer application) and animal waste management. N2O and CH4 emissions have been converted into CO2 equivalents using IPCC AR6 GWP100 factors of 273 and 27, respectively.

Data Fig6a

AFOLU GHG emissions (Gt CO2 eq per year)

Data Fig6b

AFOLU GHG emissions (Gt CO2 eq per year)

Breakdown

Column

Figure 7

Column

Information

Figure 7: AFOLU GHG emission breakdown by source. Data is shown at global level for five scenarios. a) shows CO2 emissions and removals from land-use change and management. Carbon losses consist of emissions from deforestation, conversion of non-forest ecosystems, drained peatlands and wood harvest. Carbon gains consist of carbon storage in wood products, re/afforestation (Aff CO2-price and Aff NDC), timber plantations and regrowth of natural vegetation (secondary forests and other natural land). The black line shows the net effect of carbon losses and carbon gains at global level. b) shows CH4 emissions from agriculture. c) shows N2O emissions from agriculture. Further details on AFOLU GHG emission sources and sinks are provided in the methods section and in Table S1.

Data Fig7a

CO2 emissions from land-use change (Gt CO2 per year)

Data Fig7b

CH4 emissions from agriculture (Mt CH4 per year)

Data Fig7c

N2O emissions from agriculture (Mt N2O per year)

Synergies & Trade-offs

Column

Figure 3

Column

Information

Figure 3: Summary of results for five scenarios and six indicators (mapped to SDGs) for 2050 at regional and global level. a) CO2 emissions and removals from land-use change and management, b) N2O emissions from agriculture, c) CH4 emissions from agriculture, d) change in forest area without plantations, e) nitrogen fixation and f) agricultural water use. Black dots show the global net effect for each indicator and scenario. Blue vertical lines show the respective value in year 2020. Percentage labels in panels b, c, e and f show the relative change of the global indicator level for each scenario compared to the Global-Inequality scenario. For panels a and d it is not meaningful to show percentage changes because these indicators can, by definition, turn net-zero (in case of divergent regional dynamics) or net-negative at global level.

Data Fig3a

CO2 emissions from land-use change (Gt CO2 eq per year)

Data Fig3b

Annual N2O emissions from agriculture (Gt CO2 eq per year)

Data Fig3c

Annual CH4 emissions from agriculture (Gt CO2 eq per year)

Data Fig3d

Forest area without plantations (Change in Mha wrt to 2020)

Data Fig3e

Nitrogen fixation (Mt N per year)

Data Fig3f

Agricultural water use (km3 per year)

Overview

Column

Figure 4

Column

Information

Figure 4: Land-use change of major land types over time compared to 2020 for five scenarios. a) show results at global level over time. b) shows results at regional level for 2050 and 2100. Cropland includes food, non-food, and feed crops. Bioenergy includes 2nd generation bioenergy (fast growing grasses and trees such as miscanthus and poplar). Pasture includes rangeland and managed pasture areas. Timber refers to timber plantations for wood production. Aff NDC includes prescribed re/afforestation according to National Determined Contributions (NDCs) towards the objectives of the Paris Agreement. Aff CO2-Price refers to endogenous re/afforestation based on the underlying scenario-specific CO2 price trajectory. Secondary forest includes modified and regrown forest, e.g., following wood harvest or cropland abandonment. Primary forest is intact forest without signs of human intervention. Urban land includes built-up area. Other natural land is a residual category, which includes among others non-forest ecosystems, deserts, and shrublands.

Data Fig4a

Land-use change (Change in Mha wrt to 2020)

Data Fig4b

Land-use change (Change in Mha wrt to 2020)

Productivity

Column

Figure 5

Column

Information

Figure 5: Land-use intensity factor τ. Data is shown at regional level for five scenarios. The τ factor reflects the degree of crop yield amplification caused by human activities. A duplication of τ implies a doubling of crop yields under fixed environmental conditions. See Figure S10 for validation data.

Data Fig5

Land-use intensity factor τ

Water & Nitrogen

Column

Figure 8

Column

Information

Figure 8: Change in agricultural water use, nitrogen fixation and forest area compared to 2020 for five scenarios at regional level for the years 2050 and 2100. The black dot indicates the net effect at global level.

Data Fig8

Scenario Narratives

Column

Figure 1

Column

Information

Figure 1: Key socio-economic drivers and developments. Data is shown at regional and global level for the two main scenarios Global-Inequality (same for Global-GHGPrice and Global-EnvirProt) and Global-Sustainability (same for Global-SustDemand). a) Population, b) Per-capita income, c) Per-capita calorie intake, d) Total demand for crops (including food and feed) and livestock products in million-ton dry matter per year, e) Prevalence of underweight (body mass index < 18.5) and f) Prevalence of obesity (BMI > 30). Details on prevalence of underweight and obesity can be found in Table S1.

Data Fig1a

Population (billion people)

Data Fig1b

Income (USD05 PPP cap-1 yr-1)

Data Fig1c

Per-capita calorie intake (kcal cap-1 day-1)

Data Fig1d

Total demand for crops and livestock products (Mt DM yr-1)

Data Fig1e

Prevalence of underweight BMI < 18.5 (million people)

Data Fig1f

Prevalence of obesity BMI > 30 (million people)

Overview

Column

Figure 2

Column

Information

Figure 2: Per-capita calorie supply. The sum of crops, livestock products, fish and secondary products reflects food intake, which together with food waste sums up to food supply. Percentage values indicate the share of food waste in total food supply. Data is shown at regional level for the two main scenarios Global-Inequality (same for Global-GHGPrice and Global-EnvirProt) and Global-Sustainability (same for Global-SustDemand).

Data Fig2

Per-capita calorie supply (kcal cap-1 day-1)

Details

Column

Figure S3

Column

Information

Figure S3: Per-capita calorie intake. a) shows regional data for 2020. The height of each rectangle shows per-capita kcal intake, the width shows the population of the region, so that the area of the rectangles refers to the total calorie intake for each region. b) shows EAT-Lancet recommendations (planetary health diet), aggregated to global level. In the Global-Sustainability scenario, all regions converge towards the EAT-Lancet recommendation by 2050.

---
title: "LAMACLIMA"
output: 
  flexdashboard::flex_dashboard:
    storyboard: false
    orientation: columns
    self_contained: false
    social: menu
    source_code: embed

---

```{r include=FALSE}
library(DT)
library(data.table)
library(knitr)
library(htmlwidgets)
opts_chunk$set(message=FALSE)

#source plotting function
#source("LAMACLIMA_functions.R")

library(luplot)
library(ggplot2)
library(luscale)
library(data.table)
library(ggrepel)
library(patchwork)
library(quitte)
library(RColorBrewer)
library(mip)
library(ggiraph)
library(ggiraphExtra)
library(scales)
library(writexl)
#source("theme_my.R")
theme_my <- function(base_size = 11, base_family = "",rotate_x=FALSE, panel.spacing=3) {
  txt <- element_text(size = base_size, colour = "black", face = "plain")
  bold_txt <- element_text(size = base_size, colour = "black", face = "bold")
  
  theme_bw(base_size = base_size, base_family = base_family) +
    theme(
      legend.key = element_blank(), 
      strip.background = element_rect(color="black",fill="grey95"),
      axis.text.x = if(rotate_x) element_text(angle = 90, hjust = 1, vjust = 0.5) else element_text(angle = 0),
      axis.title.x = element_text(vjust=0),
      axis.title.y = element_text(margin = margin(t = 0, r = 5, b = 0, l = 0)),
      axis.title.y.right = element_text(margin = margin(t = 0, r = 0, b = 0, l = 5)),
      panel.grid.minor.x = element_blank(),
      panel.grid.major.x = element_blank(),
      
      panel.spacing.x = unit(panel.spacing, "mm"),
      panel.spacing.y = unit(panel.spacing, "mm"),
      
      text = txt, 
      plot.title = bold_txt, 
      
      axis.title = bold_txt, 
      axis.text = txt, 
      
#      legend.title = bold_txt, 
#      legend.text = txt,
      legend.text = element_text(margin = margin(r = 10), hjust = 0,size = base_size, colour = "black", face = "plain"),
      legend.title = element_text(margin = margin(r = 10), hjust = 0,size = base_size, colour = "black", face = "bold")
      )
  # +
  #    theme(legend.position="bottom",legend.box = "vertical",legend.title.align = 0)
}


set.seed(42)

### Contact
# Florian Humpenoeder
# Potsdam Institute for Climate Impact Research
# humpenoeder@pik-potsdam.de
# September 2022

if(sys.nframe() == 0L){
  #set revision prefix
  rev <- "HR_LAMA91"

  #load data as data.table
  rep <- as.data.table(readRDS(paste0("report_WP4_",rev,".rds")))

  #load validation data as data.table
  hist <- as.data.table(readRDS("validation2.rds"))

  #Rename GLO to World
  rep[region=="GLO",region:="World"]

  #split up the scenario name by "_"
  rep[, c("scenario") := gsub(paste0(rev,"_"),"",scenario)]
  rep$scenario <- factor(rep$scenario)

  #mapping for regional aggregation
  map <- data.frame(matrix(nrow=12,ncol=2))
  names(map) <- c("region_class","region")
  map[1,] <- c("OECD90+EU","EUR")
  map[2,] <- c("OECD90+EU","USA")
  map[3,] <- c("OECD90+EU","CAZ")
  map[4,] <- c("OECD90+EU","JPN")
  map[5,] <- c("OECD90+EU","NEU")
  map[6,] <- c("ASIA","CHA")
  map[7,] <- c("ASIA","IND")
  map[8,] <- c("LAM","LAM")
  map[9,] <- c("SSA","SSA")
  map[10,] <- c("ASIA","OAS")
  map[11,] <- c("ROW","MEA")
  map[12,] <- c("ROW","REF")
  map[13,] <- c("World","World")
  rep <- merge(rep,map)

  #order of regions
  reg_sub_order <- c("OECD90+EU","ASIA","LAM","SSA","ROW")
  reg_sel <- c("World","EUR","LAM","SSA")
}


##### Main figures

### Fig 1 Main Drivers

Fig1 <- function(rep,reg_sub_order) {
  scen_sel <- rev(c("Global-Sustainability","Global-Inequality"))
  names(scen_sel) <- rev(c("Global-Sustainability","Global-Inequality"))
  
  var <- c("Population")
  title <- c("Population")
  unit <- c("billion people")
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- droplevels(b)
  pop <- b
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel)
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  b[,value:=value/1000]
  b$unit <- unit
  # b_glo <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region == "World",]
  # b_glo <- droplevels(b_glo)
  # b_glo[,value:=value/1000]
  # b_glo$region <- NULL
  # b_glo <- setcolorder(b_glo,names(b))
  
  # b$period <- as.factor(b$period)
  # b_glo$period <- as.factor(b_glo$period)
  p <- ggplot(b,aes(y=value,x=period))+facet_wrap(vars(scenario),ncol = 2)+theme_my(rotate_x = TRUE)+xlab(NULL)
  p <- p + geom_bar_interactive(aes(tooltip = paste0(region_class," ",period,"\n",formatC(value,big.mark = ",",format="fg")," ",unit),data_id=interaction(period),fill=region_class),position="stack",stat = "identity",width=5)
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World ",x,"\n",formatC(y,big.mark = ",",format="fg")," billion people")),data_id=period))
  p <- p + scale_fill_brewer_interactive("Region",palette = "Set2",direction = "-1")
  p1 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+ylab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(0,2,4,6,8,10))# + labs(caption = paste(Sys.Date()))
  
  #Data for table
  tableColOrder <- c("Scenario","Region","Year","Value")
  b <- b[,c("scenario","region_class","period","value")]
  b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(b,c("scenario","region_class","period","value"),tableColOrder)
  Fig1a <- copy(b)
  
  var <- c("Income")
  title <- c("Income")
  unit <- expression(bold('USD05 PPP cap'^{-1} ~ 'yr'^{-1}))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- cbind(b,pop$value)
  b <- b[,.(value=weighted.mean(value,V2)),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  p <- ggplot(b,aes(y=value,x=period))+facet_wrap(vars(scenario),ncol = 2)+theme_my(rotate_x = TRUE)+xlab(NULL)
  p <- p+geom_line_interactive(aes(color=region_class,group=region_class,data_id=period),size=1)
  p <- p+geom_point_interactive(aes(tooltip = paste0(region_class," ",period,"\n",formatC(value,big.mark = ",",format="fg")," USD/cap/yr"),data_id=interaction(period),color=region_class,group=region_class),size=2)
  p <- p + scale_color_brewer("Region",palette = "Set2",direction = "-1")
  p2 <- p+theme(legend.position = "none")+guides(color="none")+ylab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(0,2,4,6,8,10))# + labs(caption = paste(Sys.Date()))
  
  #Data for table
  tableColOrder <- c("Scenario","Region","Year","Value")
  b <- b[,c("scenario","region_class","period","value")]
  b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(b,c("scenario","region_class","period","value"),tableColOrder)
  Fig1b <- copy(b)
  
  var <- c("Nutrition|Calorie Intake")
  title <- c("Per-capita calorie intake")
  unit <- expression(bold('kcal cap'^{-1} ~ 'day'^{-1}))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- cbind(b,pop$value)
  b <- b[,.(value=weighted.mean(value,V2)),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  b$variable <- title
  #b <- setcolorder(b,names(Fig1a))
  p <- ggplot(b,aes(y=value,x=period))+facet_wrap(vars(scenario),ncol = 2)+theme_my(rotate_x = TRUE)+xlab(NULL)
  p <- p+geom_line_interactive(aes(color=region_class,group=region_class,data_id=period),size=1)
  p <- p+geom_point_interactive(aes(tooltip = paste0(region_class," ",period,"\n",formatC(value,big.mark = ",",format="fg")," kcal/cap/day"),data_id=interaction(period),color=region_class,group=region_class),size=2)
  p <- p + scale_color_brewer("Region",palette = "Set2",direction = "-1")
  p3 <- p+theme(legend.position = "none")+guides(color="none")+ylab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(0,2,4,6,8,10))# + labs(caption = paste(Sys.Date()))
  
  #Data for table
  tableColOrder <- c("Scenario","Region","Year","Value")
  b <- b[,c("scenario","region_class","period","value")]
  b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(b,c("scenario","region_class","period","value"),tableColOrder)
  Fig1c <- copy(b)
  
  var <- c("Demand|++|Crops","Demand|++|Livestock products")
  title <- c("Total demand for crops and livestock products")
  unit <- expression(bold('Mt DM yr'^{-1}))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  b$variable <- title
  #b <- setcolorder(b,names(Fig1a))
  p <- ggplot(b,aes(y=value,x=period))+facet_wrap(vars(scenario),ncol = 2)+theme_my(rotate_x = TRUE)+xlab(NULL)
  p <- p + geom_bar_interactive(aes(tooltip = paste0(region_class," ",period,"\n",formatC(value,big.mark = ",",format="fg")," ",unit),data_id=interaction(period),fill=region_class),position="stack",stat = "identity",width=5)
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World ",x,"\n",formatC(y,big.mark = ",",format="fg")," Mt DM/yr")),data_id=period))
  p <- p + scale_fill_brewer("Region",palette = "Set2",direction = "-1")
  p4 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+ylab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(0,2,4,6,8,10))# + labs(caption = paste(Sys.Date()))
  
  #Data for table
  tableColOrder <- c("Scenario","Region","Year","Value")
  b <- b[,c("scenario","region_class","period","value")]
  b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(b,c("scenario","region_class","period","value"),tableColOrder)
  Fig1d <- copy(b)
  
  var <- c("SDG|SDG02|Prevalence of underweight")
  title <- c("Prevalence of underweight (BMI < 18.5) - SDG 2")
  unit <- c("million people")
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  #b <- cbind(b,pop$value)
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  b$variable <- title
  b$unit <- unit
  #b <- setcolorder(b,names(Fig1a))
  p <- ggplot(b,aes(y=value,x=period))+facet_wrap(vars(scenario),ncol = 2)+theme_my(rotate_x = TRUE)+xlab(NULL)
  p <- p + geom_bar_interactive(aes(tooltip = paste0(region_class," ",period,"\n",formatC(value,big.mark = ",",format="fg")," ",unit),data_id=interaction(period),fill=region_class),position="stack",stat = "identity",width=5)
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World ",x,"\n",formatC(y,big.mark = ",",format="fg")," ",unit)),data_id=period))
  p <- p + scale_fill_brewer("Region",palette = "Set2",direction = "-1")
  p5 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+ylab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(0,2,4,6,8,10))# + labs(caption = paste(Sys.Date()))
  
  #Data for table
  tableColOrder <- c("Scenario","Region","Year","Value")
  b <- b[,c("scenario","region_class","period","value")]
  b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(b,c("scenario","region_class","period","value"),tableColOrder)
  Fig1e <- copy(b)
  
  var <- c("SDG|SDG03|Prevalence of obesity")
  title <- c("Prevalence of obesity (BMI > 30) - SDG 3")
  unit <- c("million people")
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- cbind(b,pop$value)
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  b$variable <- title
  b$unit <- unit
  #b <- setcolorder(b,names(Fig1a))
  p <- ggplot(b,aes(y=value,x=period))+facet_wrap(vars(scenario),ncol = 2)+theme_my(rotate_x = TRUE)+xlab(NULL)
  p <- p + geom_bar_interactive(aes(tooltip = paste0(region_class," ",period,"\n",formatC(value,big.mark = ",",format="fg")," ",unit),data_id=interaction(period),fill=region_class),position="stack",stat = "identity",width=5)
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World ",x,"\n",formatC(y,big.mark = ",",format="fg")," ",unit)),data_id=period))
  p <- p + scale_fill_brewer("Region",palette = "Set2",direction = "-1")
  p6 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+ylab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(0,2,4,6,8,10))# + labs(caption = paste(Sys.Date()))
  b[period %in% c(2020,2100),.(value=round(sum(value),0)),by = .(model,scenario,unit,period)]
  
  #Data for table
  tableColOrder <- c("Scenario","Region","Year","Value")
  b <- b[,c("scenario","region_class","period","value")]
  b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(b,c("scenario","region_class","period","value"),tableColOrder)
  Fig1f <- copy(b)
  
  #writexl::write_xlsx(list(Fig1a = Fig1a,Fig1b = Fig1b,Fig1c = Fig1c,Fig1d = Fig1d,Fig1e = Fig1e,Fig1f = Fig1f),"Fig1_drivers.xlsx")
  
  # combined <- (p1 + p2) / (p3 + p4) / (p5 + p6) + plot_annotation(tag_levels = 'a')
  # combined <- combined + plot_layout(guides = "collect",nrow = 3,heights = c(1,1,1)) & theme(legend.position = "bottom")
  # ggsave(filename = "Fig1_drivers.pdf",combined,width = 8,height = 6,scale=1.4)
  # ggsave(filename = "Fig1_drivers.png",combined,width = 8,height = 6,scale=1.4)
  #return(list(plot=p1,data=NULL,dataH12=NULL))
  
  combined <- girafe(code = print((p1 + p2) / (p3 + p4) / (p5 + p6) + plot_annotation(tag_levels = 'a') + plot_layout(guides = "collect",nrow = 3,heights = c(1,1,1)) & theme(legend.position = "bottom")),
    options = list(
      opts_sizing(rescale = TRUE, width = .95),
      opts_hover(girafe_css(css = "fill:NULL;cursor:pointer;", area = "fill:NULL;")),
      opts_hover_inv(css = "opacity:0.1;"),
      opts_selection(girafe_css(css = "fill:NULL;stroke:black"), only_shiny = FALSE, type="multiple",selected = NULL),
      opts_tooltip(css = "background-color:white;padding:5px;border-radius:2px;border: black 1px solid;color:black;")
    ),
    width_svg = 12,
    height_svg = 8
  )
  
  return(list(plot=combined,Fig1a=Fig1a,Fig1b=Fig1b,Fig1c=Fig1c,Fig1d=Fig1d,Fig1e=Fig1e,Fig1f=Fig1f))
}

Fig2 <- function(rep,reg_sub_order) {
  ### Fig2 Food waste
  scen_sel <- c("Global-Sustainability","Global-Inequality")
  names(scen_sel) <- c("Global-Sustainability","Global-Inequality")

  var <- c("Population")
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- droplevels(b)
  pop <- b
  
  title <- "Food waste"
  var <- rev(c("Nutrition|Calorie Intake|+|Crops","Nutrition|Calorie Intake|+|Livestock products","Nutrition|Calorie Intake|+|Fish","Nutrition|Calorie Intake|+|Secondary products","SDG|SDG12|Food waste"))
  names(var) <- rev(c("Crops","Livestock products","Fish","Secondary products","Food waste"))
  unit <- expression(bold('kcal' ~ 'capita'^{-1} ~ 'day'^{-1}))
  b <- rep[variable %in% var & region!="World" & scenario %in% scen_sel & period >= 2020 & period <= 2100,]
  b <- cbind(b,pop$value)
  b <- b[,.(value=weighted.mean(value,V2)),by = .(region_class,model,scenario,variable,unit,period)]
  b$region_class <- factor(b$region_class,levels = c(reg_sub_order))
  b$scenario <- factor(b$scenario,levels = scen_sel, labels=names(scen_sel))
  b$variable <- factor(b$variable,levels = var,labels = names(var))
  b <- b[period %in% c(2020,2050,2100),]
  b[period==2020,value:=value[scenario=="Global-Sustainability"],by = .(region_class,model,variable,unit,period)]
  b[period==2020,scenario:="All scenarios"]
  b <- unique(b)
  
  b[,label2:=round((value[variable=="Food waste"]/sum(value))*100,1),by=.(region_class,scenario,model,period)]
  b[variable!="Food waste",label2:=NA]
  b[,labelpos:=sum(value)-value[variable=="Food waste"]/2,by=.(region_class,scenario,model,period)]
  b[variable!="Food waste",labelpos:=NA]
  b[variable=="Food waste",]

  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)

  p <- ggplot(b,aes(x=value,y=scenario))+facet_grid(vars(period),vars(region_class),space = "free_y",scales = "free_y")+theme_my(panel.spacing = 3,rotate_x = 0)+xlab(unit)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=variable,data_id=interaction(variable),tooltip=paste0(region_class," ",period,"\n",formatC(value,big.mark = ",",format="fg")," kcal/cap/day")),position="stack",stat = "identity")+theme(legend.position = "bottom")
  p <- p + scale_fill_manual("Type",values = c("#c51b8a","#74c476","#2b8cbe","#756bb1","#d95f0e"))+guides(fill=guide_legend("Type",ncol=5,title.position = "left", byrow = TRUE,reverse=TRUE))
  p <- p + geom_text_interactive(aes(label=paste0(label2," %"),x=labelpos,data_id=interaction(variable),tooltip=paste0(region_class," ",period,"\n",formatC(value,big.mark = ",",format="fg"))),angle=90,color="white")

  #Data for table
  tableColOrder <- c("Scenario","Region","Year","Value")
  b <- b[,c("scenario","region_class","period","value")]
  b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(b,c("scenario","region_class","period","value"),tableColOrder)
  Fig2 <- copy(b)
  
  combined <- girafe(ggobj = p,
                     options = list(
                       opts_sizing(rescale = TRUE, width = .95),
                       opts_hover(girafe_css(css = "fill:NULL;cursor:pointer;", area = "fill:NULL;")),
                       opts_hover_inv(css = "opacity:0.1;"),
                       opts_selection(girafe_css(css = "fill:NULL;stroke:black"), only_shiny = FALSE, type="multiple",selected = NULL),
                       opts_tooltip(css = "background-color:white;padding:5px;border-radius:2px;border: black 1px solid;color:black;")
                     ),
                     width_svg = 12,
                     height_svg = 6
  )

  return(list(plot=combined,Fig2=Fig2))
  
  # writexl::write_xlsx(list(Fig2 = Fig2),"Fig2_food_waste.xlsx")
  # 
  # ggsave(filename = "Fig2_food_waste.pdf",p,width = 8,height = 4,scale=1.4)
  # ggsave(filename = "Fig2_food_waste.png",p,width = 8,height = 4,scale=1.4)

}

Fig3 <- function(rep,reg_sub_order) {
  ### Fig 3 Summary

  years <- c(2050)
  scen_sel <- rev(c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability"))
  names(scen_sel) <- rev(c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability"))

  ##Fig3a
  var <- c("Emissions|CO2|Land|+|Land-use Change")
  title <- expression(bold('Annual CO'[2] ~ 'emissions land-use change & land management (SDG 13)'))
  unit <- expression(bold('Gt CO'[2] ~ 'eq' ~ 'yr'^{-1}))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- b[,.(value=sum(value)/1000),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  Fig3a <- copy(b)
  Fig3a$variable <- "Annual CO2 emissions land-use change"
  Fig3a$unit <- "Gt CO2eq per year"
  tableColOrder <- c("Scenario","Region","Year","Value")
  Fig3a <- Fig3a[,c("scenario","region_class","period","value")]
  Fig3a[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig3a,c("scenario","region_class","period","value"),tableColOrder)
  
  b[,value2020:=value[period==2020],by=.(region_class,model,scenario,variable,unit)]
  b <- b[period %in% years,]
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)
  p <- ggplot(b,aes(y=as.factor(scenario)))+facet_wrap(vars(period),ncol = 3)+theme_my(rotate_x = FALSE)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(positive,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p+geom_bar_interactive(aes(x=negative,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(negative,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p + scale_fill_brewer("Region",palette = "Set2",direction = "-1")
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World 2050\n",formatC(x,big.mark = ",",format="fg")," Gt CO2eq per year")),data_id=period,group = scenario,x=value))
  #p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = "point", mapping = aes(group = scenario,x=value))
  p <- p + geom_vline(xintercept = 0)#+scale_x_continuous(breaks = c(-6000,-3000,0,3000,6000))
  p <- p + geom_vline(data = b[scenario==names(scen_sel)[1],],aes(xintercept = sum(value2020)),linetype="solid",colour="blue")
  p <- p + geom_text(data = b[scenario==names(scen_sel)[1],],aes(x=sum(value2020), label="\n2020", y=1), colour="blue", angle=90, hjust=0.3)
  p1 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+xlab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(-500,-250,0,250,500))# + labs(caption = paste(Sys.Date()))

  ##Fig3b
  var <- c("Emissions|N2O|Land|+|Agriculture")
  title <- expression(bold('Annual' ~ N[2]*O ~ 'emissions from agriculture (SDG 13)'))
  unit <- expression(bold('Gt CO'[2] ~ 'eq' ~ 'yr'^{-1}))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- b[,.(value=sum(value)*273/1000),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  Fig3b <- copy(b)
  Fig3b$variable <- "Annual N2O emissions from agriculture (SDG 13)"
  Fig3b$unit <- "Gt CO2eq per year"
  tableColOrder <- c("Scenario","Region","Year","Value")
  Fig3b <- Fig3b[,c("scenario","region_class","period","value")]
  Fig3b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig3b,c("scenario","region_class","period","value"),tableColOrder)
  
  b[,value2020:=value[period==2020],by=.(region_class,model,scenario,variable,unit)]
  b <- b[period %in% years,]
  bb<-copy(b)
  bb <- bb[,.(value=sum(value)),by = .(model,scenario,variable,unit,period)]
  bb[,label:=round((value/value[scenario=="Global-Inequality"]-1)*100),by=.(model,variable,unit,period)]
  bb[,label:= ifelse(label > 0,paste0("+",label,"%"),paste0(label,"%"))]
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)

  p <- ggplot(b,aes(y=as.factor(scenario)))+facet_wrap(vars(period),ncol = 3)+theme_my(rotate_x = FALSE)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(positive,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p+geom_bar_interactive(aes(x=negative,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(negative,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p + scale_fill_brewer("Region",palette = "Set2",direction = "-1")
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World 2050\n",formatC(x,big.mark = ",",format="fg")," Gt CO2eq per year")),data_id=period,group = scenario,x=value))
  p <- p + geom_vline(data = b[scenario==names(scen_sel)[1],],aes(xintercept = sum(value2020)),linetype="solid", colour="blue")
  p <- p + geom_text(data = b[scenario==names(scen_sel)[1],],aes(x=sum(value2020), label="\n2020", y=1), colour="blue", angle=90, hjust=0.3)
  p <- p + geom_label_repel(data=bb,aes(x=value,label = label),direction = "x",nudge_x = 0.05,size=3.5)
  p <- p+scale_x_continuous(guide = guide_axis(check.overlap = TRUE),expand = expansion(mult = c(0.05,0.1)))
  p2 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+xlab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(-500,-250,0,250,500))# + labs(caption = paste(Sys.Date()))

  ##Fig3c
  var <- c("Emissions|CH4|Land|+|Agriculture")
  title <- expression(bold('Annual' ~ CH[4] ~ 'emissions from agriculture (SDG 13)'))
  unit <- expression(bold('Gt CO'[2] ~ 'eq' ~ 'yr'^{-1}))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- b[,.(value=sum(value)*27/1000),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  Fig3c <- copy(b)
  Fig3c$variable <- "Annual CH4 emissions from agriculture (SDG 13)"
  Fig3c$unit <- "Gt CO2eq per year"
  tableColOrder <- c("Scenario","Region","Year","Value")
  Fig3c <- Fig3c[,c("scenario","region_class","period","value")]
  Fig3c[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig3c,c("scenario","region_class","period","value"),tableColOrder)
  
  b[,value2020:=value[period==2020],by=.(region_class,model,scenario,variable,unit)]
  b <- b[period %in% years,]
  bb<-copy(b)
  bb <- bb[,.(value=sum(value)),by = .(model,scenario,variable,unit,period)]
  bb[,label:=round((value/value[scenario=="Global-Inequality"]-1)*100),by=.(model,variable,unit,period)]
  bb[,label:= ifelse(label > 0,paste0("+",label,"%"),paste0(label,"%"))]
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)

  p <- ggplot(b,aes(y=as.factor(scenario)))+facet_wrap(vars(period),ncol = 3)+theme_my(rotate_x = FALSE)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(positive,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p+geom_bar_interactive(aes(x=negative,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(negative,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p + scale_fill_brewer("Region",palette = "Set2",direction = "-1")
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World 2050\n",formatC(x,big.mark = ",",format="fg")," Gt CO2eq per year")),data_id=period,group = scenario,x=value))
  p <- p + geom_vline(data = b[scenario==names(scen_sel)[1],],aes(xintercept = sum(value2020)),linetype="solid", colour="blue")
  p <- p + geom_text(data = b[scenario==names(scen_sel)[1],],aes(x=sum(value2020), label="\n2020", y=1), colour="blue", angle=90, hjust=0.3)
  p <- p + geom_label_repel(data=bb,aes(x=value,label = label),direction = "x",nudge_x = 0.05,size=3.5)
  p <- p+scale_x_continuous(guide = guide_axis(check.overlap = TRUE),expand = expansion(mult = c(0.05,0.1)))
  p3 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+xlab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(-500,-250,0,250,500))# + labs(caption = paste(Sys.Date()))

  ##Fig3d
  ### add forest area without timber to rep object
  var <- c("Resources|Land Cover|Forest|Managed Forest|+|NPI/NDC","Resources|Land Cover|Forest|Managed Forest|+|Afforestation","Resources|Land Cover|Forest|Natural Forest|+|Secondary Forest","Resources|Land Cover|Forest|Natural Forest|+|Primary Forest")
  b <- rep[variable %in% var,]
  save <- names(b)
  b <- b[,.(value=sum(value)),by = .(region_class,region,model,scenario,unit,period)]
  b[,variable:="Forest area without plantations"]
  b[,unit:="Change in Mha compared to 2020"]
  setcolorder(b,save)
  b <- droplevels(b)
  b[, value:=value-value[period==2020], by = .(region_class,region,model,scenario,variable,unit)]
  rep <- rbind(rep,b)

  var <- c("Forest area without plantations")
  title <- expression(bold('Forest area without plantations (SDG 15)'))
  unit <- expression(bold('Change in Mha compared to 2020'))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  Fig3d <- copy(b)
  Fig3d$variable <- "Forest area without plantations (SDG 15)"
  Fig3d$unit <- "Change in Mha compared to 2020"
  tableColOrder <- c("Scenario","Region","Year","Value")
  Fig3d <- Fig3d[,c("scenario","region_class","period","value")]
  Fig3d[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig3d,c("scenario","region_class","period","value"),tableColOrder)
  
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)
  b <- b[period %in% years,]
  p <- ggplot(b,aes(y=as.factor(scenario)))+facet_wrap(vars(period),ncol = 3)+theme_my(rotate_x = FALSE)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(positive,big.mark = ",",format="fg")," Change in Mha compared to 2020"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p+geom_bar_interactive(aes(x=negative,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(negative,big.mark = ",",format="fg")," Change in Mha compared to 2020"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p + scale_fill_brewer("Region",palette = "Set2",direction = "-1")
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World 2050\n",formatC(x,big.mark = ",",format="fg")," Change in Mha compared to 2020")),data_id=period,group = scenario,x=value))
  p <- p + geom_vline(xintercept = 0,linetype="dotted")#+scale_x_continuous(breaks = c(-6000,-3000,0,3000,6000))
  p4 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+xlab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(-500,-250,0,250,500))# + labs(caption = paste(Sys.Date()))

  ##Fig3e
  var <- c("SDG|SDG15|Industrial and intentional biological fixation of N")
  title <- expression(bold('Nitrogen fixation (SDG 15)'))
  unit <- expression(bold("Mt N yr"^{-1}))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  Fig3e <- copy(b)
  Fig3e$variable <- "Nitrogen fixation (SDG 15)"
  Fig3e$unit <- "Mt N per year"
  tableColOrder <- c("Scenario","Region","Year","Value")
  Fig3e <- Fig3e[,c("scenario","region_class","period","value")]
  Fig3e[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig3e,c("scenario","region_class","period","value"),tableColOrder)
  
  b[,value2020:=value[period==2020],by=.(region_class,model,scenario,variable,unit)]
  b <- b[period %in% years,]
  bb<-copy(b)
  bb <- bb[,.(value=sum(value)),by = .(model,scenario,variable,unit,period)]
  bb[,label:=round((value/value[scenario=="Global-Inequality"]-1)*100),by=.(model,variable,unit,period)]
  bb[,label:= ifelse(label > 0,paste0("+",label,"%"),paste0(label,"%"))]
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)

  p <- ggplot(b,aes(y=as.factor(scenario)))+facet_wrap(vars(period),ncol = 3)+theme_my(rotate_x = FALSE)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(positive,big.mark = ",",format="fg")," Mt N per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p+geom_bar_interactive(aes(x=negative,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(negative,big.mark = ",",format="fg")," Mt N per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p + scale_fill_brewer("Region",palette = "Set2",direction = "-1")
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World 2050\n",formatC(x,big.mark = ",",format="fg")," Mt N per year")),data_id=period,group = scenario,x=value))
  p <- p + geom_vline(data = b[scenario==names(scen_sel)[1],],aes(xintercept = sum(value2020)),linetype="solid",colour="blue")
  p <- p + geom_text(data = b[scenario==names(scen_sel)[1],],aes(x=sum(value2020), label="\n2020", y=1), colour="blue", angle=90, hjust=0.3)
  p <- p + geom_label_repel(data=bb,aes(x=value,label = label),direction = "x",nudge_x = 0.05,size=3.5)
  p <- p+scale_x_continuous(guide = guide_axis(check.overlap = TRUE),expand = expansion(mult = c(0.05,0.15)))
  p5 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+xlab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(-500,-250,0,250,500))# + labs(caption = paste(Sys.Date()))
  b[,.(value=round(sum(value),0)),by = .(model,scenario,variable,unit,period)]
  b[,.(value=round(sum(value),0)),by = .(region_class,model,scenario,variable,unit,period)]

  ##Fig3f
  var <- c("SDG|SDG06|Agricultural water use")
  title <- expression(bold('Agricultural water use (SDG 6)'))
  unit <- expression(bold(km^3 ~ "yr"^{-1}))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >=2020 & period <= 2100 & region != "World",]
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(reg_sub_order))
  Fig3f <- copy(b)
  Fig3f$variable <- "Agricultural water use (SDG 6)"
  Fig3f$unit <- "km3 per year"
  tableColOrder <- c("Scenario","Region","Year","Value")
  Fig3f <- Fig3f[,c("scenario","region_class","period","value")]
  Fig3f[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig3f,c("scenario","region_class","period","value"),tableColOrder)
  

  b[,value2020:=value[period==2020],by=.(region_class,model,scenario,variable,unit)]
  b <- b[period %in% years,]
  bb<-copy(b)
  bb <- bb[,.(value=sum(value)),by = .(model,scenario,variable,unit,period)]
  bb[,label:=round((value/value[scenario=="Global-Inequality"]-1)*100),by=.(model,variable,unit,period)]
  bb[,label:= ifelse(label > 0,paste0("+",label,"%"),paste0(label,"%"))]
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)

  p <- ggplot(b,aes(y=as.factor(scenario)))+facet_wrap(vars(period),ncol = 3)+theme_my(rotate_x = FALSE)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(positive,big.mark = ",",format="fg")," km3 per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p+geom_bar_interactive(aes(x=negative,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(negative,big.mark = ",",format="fg")," km3 per year"),data_id=interaction(scenario)),position="stack",stat = "identity")
  p <- p + scale_fill_brewer("Region",palette = "Set2",direction = "-1")
  p <- p + geom_vline(data = b[scenario==names(scen_sel)[1],],aes(xintercept = sum(value2020)),linetype="solid",colour="blue")
  p <- p + geom_text(data = b[scenario==names(scen_sel)[1],],aes(x=sum(value2020), label="\n2020", y=1), colour="blue", angle=90, hjust=0.3)
  p <- p + geom_label_repel(data=bb,aes(x=value,label = label),direction = "x",nudge_x = 1,box.padding = 0.05,force=5,size=3.5)
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World 2050\n",formatC(x,big.mark = ",",format="fg")," km3 per year")),data_id=period,group = scenario,x=value))
  p <- p+scale_x_continuous(guide = guide_axis(check.overlap = TRUE),expand = expansion(mult = c(0.05,0.15)))
  p6 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Region",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+xlab(unit)+ggtitle(title)#+scale_x_continuous(breaks = c(-500,-250,0,250,500))# + labs(caption = paste(Sys.Date()))
  
  combined <- girafe(code = print(p1 / p2 / p3 / p4 / p5 / p6 + plot_annotation(tag_levels = 'a')  + plot_layout(guides = "collect",ncol=2,byrow = FALSE) & theme(legend.position = "bottom")),
                     options = list(
                       opts_sizing(rescale = TRUE, width = .95),
                       opts_hover(girafe_css(css = "fill:NULL;cursor:pointer;", area = "fill:NULL;")),
                       opts_hover_inv(css = "opacity:0.1;"),
                       opts_selection(girafe_css(css = "fill:NULL;stroke:black"), only_shiny = FALSE, type="multiple",selected = NULL),
                       opts_tooltip(css = "background-color:white;padding:5px;border-radius:2px;border: black 1px solid;color:black;")
                     ),
                     width_svg = 12,
                     height_svg = 8
  )
  #combined
  return(list(plot=combined,Fig3a=Fig3a,Fig3b=Fig3b,Fig3c=Fig3c,Fig3d=Fig3d,Fig3e=Fig3e,Fig3f=Fig3f))
  
# 
#   writexl::write_xlsx(list(Fig3a=Fig3a,Fig3b=Fig3b,Fig3c=Fig3c,Fig3d=Fig3d,Fig3e=Fig3e,Fig3f=Fig3f),"Fig3_summary.xlsx")
# 
#   combined <- p1 / p2 / p3 / p4 / p5 / p6 + plot_annotation(tag_levels = 'a')
#   combined <- combined + plot_layout(guides = "collect",ncol=2,byrow = FALSE) & theme(legend.position = "bottom")
#   ggsave(filename = "Fig3_summary.pdf",combined,width = 12,height = 7,scale=1.2)
#   ggsave(filename = "Fig3_summary.png",combined,width = 12,height = 7,scale=1.2)
  
} 

Fig4 <- function(rep,reg_sub_order) {
  
  ### Fig 4 Land

  unit <- "Change in Mha compared to 2020"
  
  ##Fig4a
  b <- rep[region == "World" & period <= 2100,]
  var <- c("Resources|Land Cover|+|Cropland","Resources|Land Cover|Cropland|+|Bioenergy crops","Resources|Land Cover|+|Pastures and Rangelands","Resources|Land Cover|Forest|Managed Forest|+|Plantations","Resources|Land Cover|Forest|Managed Forest|+|NPI/NDC","Resources|Land Cover|Forest|Managed Forest|+|Afforestation","Resources|Land Cover|Forest|Natural Forest|+|Secondary Forest","Resources|Land Cover|Forest|Natural Forest|+|Primary Forest","Resources|Land Cover|+|Other Land","Resources|Land Cover|+|Urban Area")
  names(var) <- c("Cropland","Bioenergy","Pasture","Timber","Aff NDC","Aff CO2-Price","Secondary Forest","Primary Forest","Other Natural","Urban")
  scen_sel <- c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability")
  names(scen_sel) <- c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability")
  b <- b[scenario %in% scen_sel & variable %in% var,]
  b <- droplevels(b)
  b$variable <- factor(b$variable,levels = rev(var),labels = names(rev(var)))
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b[, value2:=value[variable=="Cropland"]-value[variable=="Bioenergy"], by = .(model,scenario,region,period)]
  b[variable=="Cropland", value:=value2, by = .(model,scenario,region,period)]
  b$value2 <- NULL
  b <- b[period >= 2020,]
  b[, value:=value-value[period==2020], by = .(model,scenario,region,variable)]
  b$unit <- unit
  Fig4a <- copy(b)
  tableColOrder <- c("Scenario","Region","Year","Variable","Value")
  Fig4a <- Fig4a[,c("scenario","region_class","period","variable","value")]
  Fig4a[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig4a,c("scenario","region_class","period","variable","value"),tableColOrder)
  
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)
  
  p <- ggplot(b,aes(x=period))+facet_grid(vars(region),vars(scenario))+theme_my(panel.spacing = 3,rotate_x = 90)+ylab(unit)
  p <- p + geom_hline(yintercept = 0,linetype="dotted")
  p <- p+geom_bar_interactive(aes(y=positive,fill=variable,tooltip = paste0(scenario,"\n",variable,"\n",region_class," ",period,"\n",formatC(positive,big.mark = ",",format="fg")," Mha"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p+geom_bar_interactive(aes(y=negative,fill=variable,tooltip = paste0(scenario,"\n",variable,"\n",region_class," ",period,"\n",formatC(negative,big.mark = ",",format="fg")," Mha"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p + scale_fill_manual("Land type",values=rev(c("chocolate4","brown3", "#E6AB02","#6a51a3","#9e9ac8","#ae017e","#66A61E","darkgreen" , "#41b6c4","#08589e")))
  p1 <- p+theme(legend.position = "bottom")+guides(fill=guide_legend(ncol=5,title.position = "left", byrow = TRUE,reverse=TRUE))+xlab(NULL)
  
  ##Fig4b
  var <- c("Resources|Land Cover|+|Cropland","Resources|Land Cover|Cropland|+|Bioenergy crops","Resources|Land Cover|+|Pastures and Rangelands","Resources|Land Cover|Forest|Managed Forest|+|Plantations","Resources|Land Cover|Forest|Managed Forest|+|NPI/NDC","Resources|Land Cover|Forest|Managed Forest|+|Afforestation","Resources|Land Cover|Forest|Natural Forest|+|Secondary Forest","Resources|Land Cover|Forest|Natural Forest|+|Primary Forest","Resources|Land Cover|+|Other Land","Resources|Land Cover|+|Urban Area")
  names(var) <- c("Cropland","Bioenergy","Pasture","Timber","Aff NDC","Aff CO2-Price","Secondary Forest","Primary Forest","Other Natural","Urban")
  scen_sel <- rev(c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability"))
  names(scen_sel) <- rev(c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability"))
  b <- rep[variable %in% var & scenario %in% scen_sel & region != "World" & period %in% c(2020,2050,2100),]
  b <- droplevels(b)
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b$variable <- factor(b$variable,levels = var,labels = names(var))
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,variable,unit,period)]
  b <- droplevels(b)
  b[variable %in% c("Afforestation","NDC"),value:=sum(value),by = .(region_class,model,scenario,unit,period)]
  b <- b[!variable == "NDC",]
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = c(reg_sub_order))
  b[, value2:=value[variable=="Cropland"]-value[variable=="Bioenergy"], by = .(model,scenario,region_class,period)]
  b[variable=="Cropland", value:=value2, by = .(model,scenario,region_class,period)]
  b$value2 <- NULL
  b[, value:=value-value[period==2020], by = .(region_class,model,scenario,variable,unit)]
  b <- b[period != 2020,]
  b$unit <- unit
  Fig4b <- copy(b)
  tableColOrder <- c("Scenario","Region","Year","Variable","Value")
  Fig4b <- Fig4b[,c("scenario","region_class","period","variable","value")]
  Fig4b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig4b,c("scenario","region_class","period","variable","value"),tableColOrder)
  
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)

  p <- ggplot(b,aes(y=scenario))+facet_grid(vars(period),vars(region_class),scales = "free",space="free")+theme_my(rotate_x = FALSE)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=variable,tooltip = paste0(scenario,"\n",variable,"\n",region_class," ",period,"\n",formatC(positive,big.mark = ",",format="fg")," Mha"),data_id=interaction(scenario)),position="stack",stat = "identity",width=0.95)
  p <- p+geom_bar_interactive(aes(x=negative,fill=variable,tooltip = paste0(scenario,"\n",variable,"\n",region_class," ",period,"\n",formatC(negative,big.mark = ",",format="fg")," Mha"),data_id=interaction(scenario)),position="stack",stat = "identity",width=0.95)
  p <- p + scale_fill_manual("Land type",values=c("chocolate4","brown3", "#E6AB02","#6a51a3","#9e9ac8","#ae017e","#66A61E","darkgreen","#41b6c4","#08589e"))
  p <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("Land type",ncol=5,title.position = "left", byrow = TRUE,reverse=F))+xlab(unit)+scale_x_continuous(guide = guide_axis(check.overlap = TRUE),breaks= c(-400,-200,0,200,400)) + theme(panel.spacing.x = unit(6, "mm"))# + labs(caption = paste(Sys.Date()))+scale_x_continuous(guide = guide_axis(check.overlap = TRUE),breaks= pretty_breaks())
  p2 <- p

  # writexl::write_xlsx(list(Fig4a = Fig4a,Fig4b = Fig4b),"Fig4_land.xlsx")
  # 
  # combined <- p1 / p2 + plot_annotation(tag_levels = 'a')
  # combined <- combined + plot_layout(guides = "collect",nrow = 2,heights = c(1.1,0.9)) & theme(legend.position = "bottom")
  # ggsave(filename = "Fig4_land.pdf",combined,width = 8,height = 6,scale=1.3)
  # ggsave(filename = "Fig4_land.png",combined,width = 8,height = 6,scale=1.3)
  
  combined <- girafe(code = print(p1 / p2 + plot_annotation(tag_levels = 'a') + plot_layout(guides = "collect",nrow = 2,heights = c(1.1,0.9)) & theme(legend.position = "bottom")),
                     options = list(
                       opts_sizing(rescale = TRUE, width = .95),
                       opts_hover(girafe_css(css = "fill:NULL;cursor:pointer;", area = "fill:NULL;")),
                       opts_hover_inv(css = "opacity:0.1;"),
                       opts_selection(girafe_css(css = "fill:NULL;stroke:black"), only_shiny = FALSE, type="multiple",selected = NULL),
                       opts_tooltip(css = "background-color:white;padding:5px;border-radius:2px;border: black 1px solid;color:black;")
                     ),
                     width_svg = 13,
                     height_svg = 8.5
  )
  
  return(list(plot=combined,Fig4a=Fig4a,Fig4b=Fig4b))
}

Fig5 <- function(rep,reg_sub_order) {
  
  ### Fig5 TAU
  
  scen_sel <- c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability")
  names(scen_sel) <- c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability")
  
  var <- c("Resources|Land Cover|+|Cropland")
  crop <- rep[scenario %in% scen_sel & variable %in% var & period >=1995 & period <= 2100 & region != "World",]
  crop_hist <- hist[variable==var & region!="World" & period >= 1975  & period <= 2007 & model == "FAO_crop_past",]
  
  var <- c("Productivity|Landuse Intensity Indicator Tau")
  title <- c("Productivity|Landuse Intensity Indicator Tau")
  unit <- expression(bold("Land-use intensity factor "*tau))
  b <- rep[scenario %in% scen_sel & variable %in% var & period >= c(1995) & region != "World",]
  b <- cbind(b,crop$value[crop$period %in% b$period])
  b <- b[,.(value=weighted.mean(value,V2)),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = reg_sub_order)
  Fig5 <- copy(b)
  tableColOrder <- c("Scenario","Region","Year","Value")
  Fig5 <- Fig5[,c("scenario","region_class","period","value")]
  Fig5[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig5,c("scenario","region_class","period","value"),tableColOrder)
  
  p <- ggplot(b,aes(y=value,x=period))+facet_wrap(vars(scenario),ncol = 5)+theme_my(rotate_x = TRUE)+xlab(NULL)
  p <- p+geom_line_interactive(aes(color=region_class,data_id=period),size=1.3)
  p <- p+geom_point_interactive(aes(tooltip = paste0(region_class," ",period,"\n",formatC(value,big.mark = ",",format="fg")),data_id=interaction(period),color=region_class,group=region_class),size=2)
  p <- p + scale_color_brewer("Region",palette = "Set2",direction = "1")
  p <- p+theme(legend.position = "bottom")+ylab(unit)#+ggtitle(title)#+scale_x_continuous(breaks = c(0,2,4,6,8,10))# + labs(caption = paste(Sys.Date()))
  
  
  #Data for table
  tableColOrder <- c("Scenario","Region","Year","Value")
  b <- b[,c("scenario","region_class","period","value")]
  b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(b,c("scenario","region_class","period","value"),tableColOrder)
  Fig2 <- copy(b)
  
  combined <- girafe(ggobj = p,
                     options = list(
                       opts_sizing(rescale = TRUE, width = .95),
                       opts_hover(girafe_css(css = "fill:NULL;cursor:pointer;", area = "fill:NULL;")),
                       opts_hover_inv(css = "opacity:0.1;"),
                       opts_selection(girafe_css(css = "fill:NULL;stroke:black"), only_shiny = FALSE, type="multiple",selected = NULL),
                       opts_tooltip(css = "background-color:white;padding:5px;border-radius:2px;border: black 1px solid;color:black;")
                     ),
                     width_svg = 12,
                     height_svg = 5
  )
  
  return(list(plot=combined,Fig5=Fig5))
  
  # writexl::write_xlsx(list(Fig5 = Fig5),"Fig5_tau.xlsx")
  # 
  # ggsave(filename = "Fig5_tau.pdf",p1,width = 10,height = 3.5,scale=1.2)
  # ggsave(filename = "Fig5_tau.png",p1,width = 10,height = 3.5,scale=1.2)

}
Fig6 <- function(rep,reg_sub_order) {
  ### Fig 6 GHG emissions

  ##Fig6a
  b <- rep[region == "World" & period <= 2100,]
  unit <- expression(bold('Gt CO'[2] ~ 'eq yr'^{-1}))#"Gt CO2eq per Year"
  var <- c("Emissions|CO2|Land|+|Land-use Change","Emissions|CH4|Land|+|Agriculture","Emissions|N2O|Land|+|Agriculture")
  names(var) <- c('CO[2]','CH[4]','N[2]*O')
  scen_sel <- c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability")
  names(scen_sel) <- c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability")
  b <- b[scenario %in% scen_sel & variable %in% var,]
  b <- droplevels(b)
  b$variable <- factor(b$variable,levels = rev(var),labels = names(rev(var)))
  b$variableName <- factor(b$variable,levels = names(rev(var)),labels = c("CO2","CH4","N2O"))
  b$scenario <- factor(b$scenario,levels = scen_sel,labels=names(scen_sel))
  b[ variable=="CH[4]", value:=value*27, by = .(model,scenario,region,period)]
  b[ variable=="N[2]*O", value:=value*273, by = .(model,scenario,region,period)]
  b[, value:=value/1000, by = .(model,scenario,region,period,variable)]
  b <- b[period >= 2020,]
  b$unit <- "Gt CO2eq per year"
  Fig6a <- copy(b)
  levels(Fig6a$variable) <- c("N2O","CH4","CO2")
  tableColOrder <- c("Scenario","Region","Year","Variable","Value")
  Fig6a <- Fig6a[,c("scenario","region_class","period","variable","value")]
  Fig6a[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig6a,c("scenario","region_class","period","variable","value"),tableColOrder)
  
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)
  p <- ggplot(b,aes(x=period))+facet_grid(vars(region),vars(scenario))+theme_my(panel.spacing = 3,rotate_x = 90)+ylab(unit)
  p <- p + geom_hline(yintercept = 0,linetype="dotted")
  p <- p+geom_bar_interactive(aes(y=positive,fill=variable,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variableName,"\n",formatC(positive,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p+geom_bar_interactive(aes(y=negative,fill=variable,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variableName,"\n",formatC(negative,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p + scale_fill_manual("AFOLU emission type",values=rev(c("#1b9e77","#d95f02","#7570b3")),labels=function(x) parse(text=x))
  p <- p + stat_summary(fun = "sum", colour = "black", size = 1, geom = "line", mapping = aes(group = scenario,y=value))
  p <- p+theme(legend.position = "bottom")+guides(fill=guide_legend(ncol=4,title.position = "left", byrow = TRUE,reverse=TRUE))+xlab(NULL)

  AR6 <- as.data.table(read.csv("IPCC_AR6_IMP_C1_AFOLU_Emissions.csv"))
  Fig6a_IPCC_AR6_raw <- copy(AR6)
  AR6[Variable=="Emissions|CH4|AFOLU","Value"] <- AR6[Variable=="Emissions|CH4|AFOLU","Value"]*27
  AR6[Variable=="Emissions|N2O|AFOLU","Value"] <- AR6[Variable=="Emissions|N2O|AFOLU","Value"]*273/1000
  AR6[,"Value"] <- AR6[,"Value"]/1000
  AR6 <- AR6[,.(Value=sum(Value)),by = .(Pathway,Year)]
  AR6 <- AR6[,.(Value=median(Value)),by = .(Year)]
  AR6 <- AR6[order(Year),]
  names(AR6) <- c("period","value")
  Fig6a_IPCC_AR6_median <- copy(AR6)
  Fig6a_IPCC_AR6_median$unit <- "Gt CO2eq per year"

  dat_text <- data.frame(
  #  label = c("grey line shows \nmedian values accross \n1.5°C compatible \nillustrative mitigation \npathways from \nthe IPCC AR6 WG3"),
    label = c("grey line indicates \n1.5°C compatible \nnet AFOLU \nGHG emissions","black line shows \nnet AFOLU \nGHG emissions"),
    scenario = c("Global-Inequality","Global-GHGPrice"),
    y     = c(-5,9),
    x     = c(2060,2065)
  )
  dat_text$scenario <- factor(dat_text$scenario)
  p <- p+geom_line(data=AR6,aes(y=value),color="grey70",linetype="solid",size=1)+geom_label(dat_text,mapping = aes(x = x, y = y, label = label),size=4)
  p1 <- p

  ##Fig6b
  b <- rep[region != "World" & period %in% c(2020,2050,2100),]
  var <- c("Emissions|CO2|Land|+|Land-use Change","Emissions|CH4|Land|+|Agriculture","Emissions|N2O|Land|+|Agriculture")
  names(var) <- c('CO[2]','CH[4]','N[2]*O')
  scen_sel <- rev(c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability"))
  names(scen_sel) <- rev(c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability"))
  b <- b[scenario %in% scen_sel & variable %in% var,]
  b <- droplevels(b)
  b$variable <- factor(b$variable,levels = var,labels = names(var))
  b$scenario <- factor(b$scenario,levels = scen_sel,labels=names(scen_sel))
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,variable,unit,period)]
  b$variableName <- factor(b$variable,levels = names(var),labels = c("CO2","CH4","N2O"))
  b$region_class <- factor(b$region_class,levels = c(reg_sub_order))
  b <- droplevels(b)
  b[ variable=="CH[4]", value:=value*27, by = .(model,scenario,region_class,period)]
  b[ variable=="N[2]*O", value:=value*273, by = .(model,scenario,region_class,period)]
  b[, value:=value/1000, by = .(model,scenario,region_class,period,variable)]
  b[period==2020,scenario:="All scenarios"]
  b <- unique(b)
  b$unit <- "Gt CO2eq per year"
  Fig6b <- copy(b)
  levels(Fig6b$variable) <- c("CO2","CH4","N2O")
  tableColOrder <- c("Scenario","Region","Year","Variable","Value")
  Fig6b <- Fig6b[,c("scenario","region_class","period","variable","value")]
  Fig6b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig6b,c("scenario","region_class","period","variable","value"),tableColOrder)
  
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)
  unit <- expression(bold('Gt CO'[2] ~ 'eq yr'^{-1}))#"Gt CO2eq per Year"
  p <- ggplot(b,aes(y=scenario))+facet_grid(vars(period),vars(region_class),scales = "free",space = "free")+theme_my(rotate_x = FALSE)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=variable,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variableName,"\n",formatC(positive,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(scenario)),position="stack",stat = "identity",width=0.95)
  p <- p+geom_bar_interactive(aes(x=negative,fill=variable,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variableName,"\n",formatC(negative,big.mark = ",",format="fg")," Gt CO2eq per year"),data_id=interaction(scenario)),position="stack",stat = "identity",width=0.95)
  p <- p + scale_fill_manual("AFOLU emission type",values=c("#1b9e77","#d95f02","#7570b3"),labels=function(x) parse(text=x))+geom_vline(xintercept = 0,linetype="dotted")
  p <- p + stat_summary(fun = "sum", colour = "black", size = 1, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World 2050\n",formatC(x,big.mark = ",",format="fg")," Gt CO2eq per year")),data_id=scenario,group = scenario,x=value))
  #p <- p + stat_summary(fun = "sum", colour = "black", size = 1, geom = "point", mapping = aes(group = scenario,x=value))
  p <- p+theme(legend.position = "bottom")+guides(fill=guide_legend("AFOLU emission type",ncol=4,title.position = "left", byrow = TRUE,reverse=F))+xlab(unit)+scale_x_continuous(guide = guide_axis(check.overlap = TRUE),expand = expansion(mult = c(0.05,0.1)),breaks=pretty_breaks())#breaks= function(x) seq(round(min(x)/0.5)*0.5, round(max(x)/0.5)*0.5, by = 1)#+scale_x_continuous(breaks = c(-2,0,2,4))# + labs(caption = paste(Sys.Date()))
  p2 <- p

  combined <- girafe(code = print(p1 / p2 + plot_annotation(tag_levels = 'a') + plot_layout(guides = "collect",nrow = 2,heights = c(0.9,1.1)) & theme(legend.position = "bottom")),
                     options = list(
                       opts_sizing(rescale = TRUE, width = .95),
                       opts_hover(girafe_css(css = "fill:NULL;cursor:pointer;", area = "fill:NULL;")),
                       opts_hover_inv(css = "opacity:0.1;"),
                       opts_selection(girafe_css(css = "fill:NULL;stroke:black"), only_shiny = FALSE, type="multiple",selected = NULL),
                       opts_tooltip(css = "background-color:white;padding:5px;border-radius:2px;border: black 1px solid;color:black;")
                     ),
                     width_svg = 13,
                     height_svg = 8.5
  )
  
  return(list(plot=combined,Fig6a=Fig6a,Fig6b=Fig6b))
  
  # writexl::write_xlsx(list(Fig6a = Fig6a,Fig6a_IPCC_AR6_median = Fig6a_IPCC_AR6_median,Fig6a_IPCC_AR6_raw=Fig6a_IPCC_AR6_raw,Fig6b=Fig6b),"Fig6_emis_GHG.xlsx")
  # 
  # combined <- p1 / p2 + plot_annotation(tag_levels = 'a')
  # combined <- combined + plot_layout(guides = "collect",nrow = 2,heights = c(0.9,1.1)) & theme(legend.position = "bottom")
  # ggsave(filename = "Fig6_emis_GHG.pdf",combined,width = 8,height = 7,scale=1.3)
  # ggsave(filename = "Fig6_emis_GHG.png",combined,width = 8,height = 7,scale=1.3)
  
}

Fig7 <- function(rep,reg_sub_order) {
  ### Fig7 Emission Breakdown
  
  ##Fig7a
  scen_sel <- c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability")
  names(scen_sel) <- c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability")
  
  var <- c("Emissions|CO2|Land|Land-use Change|+|Gross LUC",
           "Emissions|CO2|Land|Land-use Change|+|Peatland",
           "Emissions|CO2|Land|Land-use Change|Wood Harvest|Short Lived Products",
           "Emissions|CO2|Land|Land-use Change|Wood Harvest|Storage",
           "Emissions|CO2|Land|Land-use Change|Regrowth|CO2-price AR",
           "Emissions|CO2|Land|Land-use Change|Regrowth|NPI_NDC AR",
           "Emissions|CO2|Land|Land-use Change|Regrowth|Timber Plantations",
           "Emissions|CO2|Land|Land-use Change|Regrowth|Secondary Forest",
           "Emissions|CO2|Land|Land-use Change|Regrowth|Other Land")
  
  names(var) <- c("Gross Land-use Change",
                  "Managed Peatland",
                  "Wood Harvest",
                  "Wood Storage Pool",
                  "Managed Aff CO2-price",
                  "Managed Aff NDC",
                  "Managed Timber Plantations",
                  "Regrowth Secondary Forest",
                  "Regrowth Other Natural Land")
  
  var_colors <- c("darkgreen",
                  "#E6AB02",
                  "chocolate4",
                  "#08589e",
                  "#ae017e",
                  "#9e9ac8",
                  "#6a51a3",
                  "#66A61E",
                  "#41b6c4")
  
  title <- expression(bold('CO'[2] ~ 'emissions and removals from land-use change and management'))
  unit <- expression(bold('Gt CO'[2] ~ 'yr'^{-1}))#"Gt CO2eq per Year"
  b <- rep[variable %in% var & region=="World" & scenario %in% scen_sel & period >= 2005,]
  b[,value:=value/1000]
  b$scenario <- factor(b$scenario,levels = scen_sel, labels=names(scen_sel))
  b$variable <- factor(b$variable,levels = var, labels=names(var))
  b <- droplevels(b)
  Fig7a <- copy(b)
  tableColOrder <- c("Scenario","Region","Year","Variable","Value")
  Fig7a <- Fig7a[,c("scenario","region_class","period","variable","value")]
  Fig7a[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig7a,c("scenario","region_class","period","variable","value"),tableColOrder)
  
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)
  p <- ggplot(b,aes(x=period))+facet_grid(rows = NULL,vars(scenario))+theme_my(panel.spacing = 3,rotate_x = 90)+ylab(unit)+ggtitle(title)
  p <- p+geom_bar_interactive(aes(y=positive,fill=variable,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variable,"\n",formatC(positive,big.mark = ",",format="fg")," Gt CO2 per year"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p+geom_bar_interactive(aes(y=negative,fill=forcats::fct_rev(variable),tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variable,"\n",formatC(negative,big.mark = ",",format="fg")," Gt CO2 per year"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p + scale_fill_manual("Carbon source / sink",values=var_colors)
  p <- p+theme(legend.position = "right")+guides(fill=guide_legend(ncol=1,title.position = "top", byrow = TRUE,reverse=FALSE))+xlab(NULL)
  p <- p + geom_hline(yintercept = 0,linetype="dotted")
  p <- p + stat_summary(fun = "sum", colour = "black", size = 1, geom = "line", mapping = aes(group = scenario,y=value))+scale_x_continuous(breaks = c(2005,2025,2050,2075,2100))
  p1 <- p
  
  ##Fig7b
  var <- c(
    "Emissions|CH4|Land|Agriculture|+|Animal waste management",
    "Emissions|CH4|Land|Agriculture|+|Enteric fermentation",
    "Emissions|CH4|Land|Agriculture|+|Rice",
    "Emissions|CH4|Land|Peatland|+|Managed")
  
  names(var) <- c("Animal Waste Management",
                  "Enteric Fermentation",
                  "Rice",
                  "Managed Peatland")
  
  var_colors <- c("#377eb8",
                  "#f781bf",
                  "#999999",
                  "#E6AB02")
  
  title <- expression(bold('CH'[4] ~ 'emissions from agriculture'))
  unit <- expression(bold('Mt CH'[4] ~ 'yr'^{-1}))#"Gt CO2eq per Year"
  b <- rep[variable %in% var & region=="World" & scenario %in% scen_sel & period >= 2005,]
  b$scenario <- factor(b$scenario,levels = scen_sel, labels=names(scen_sel))
  b$variable <- factor(b$variable,levels = var, labels=names(var))
  b <- droplevels(b)
  Fig7b <- copy(b)
  tableColOrder <- c("Scenario","Region","Year","Variable","Value")
  Fig7b <- Fig7b[,c("scenario","region_class","period","variable","value")]
  Fig7b[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig7b,c("scenario","region_class","period","variable","value"),tableColOrder)
  
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)
  p <- ggplot(b,aes(x=period))+facet_grid(rows = NULL,vars(scenario))+theme_my(panel.spacing = 3,rotate_x = 90)+ylab(unit)+ggtitle(title)
  p <- p+geom_bar_interactive(aes(y=positive,fill=variable,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variable,"\n",formatC(positive,big.mark = ",",format="fg")," Mt CH4 per year"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p+geom_bar_interactive(aes(y=negative,fill=forcats::fct_rev(variable),tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variable,"\n",formatC(negative,big.mark = ",",format="fg")," Mt CH4 per year"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p + scale_fill_manual("Methane source",values=var_colors)
  p <- p+theme(legend.position = "right")+guides(fill=guide_legend(ncol=1,title.position = "top", byrow = TRUE,reverse=FALSE))+xlab(NULL)
  p <- p + geom_hline(yintercept = 0,linetype="dotted")
  p <- p + scale_x_continuous(breaks = c(2005,2025,2050,2075,2100))
  p2 <- p
  
  ##Fig7c
  var <- c(  "Emissions|N2O|Land|Agriculture|+|Animal Waste Management",
             "Emissions|N2O|Land|Agriculture|+|Agricultural Soils",
             "Emissions|N2O|Land|Peatland|+|Managed")
  
  names(var) <- c("Animal Waste Management",
                  "Agricultural Soils",
                  "Managed Peatland")
  
  var_colors <- c("#377eb8",
                  "#bf5b17",
                  "#E6AB02")
  
  title <- expression(bold(N[2]*O ~ 'emissions from agriculture'))
  unit <- expression(bold('Mt' ~ N[2]*O ~ 'yr'^{-1}))
  b <- rep[variable %in% var & region=="World" & scenario %in% scen_sel & period >= 2005,]
  b$scenario <- factor(b$scenario,levels = scen_sel, labels=names(scen_sel))
  b$variable <- factor(b$variable,levels = var, labels=names(var))
  b <- droplevels(b)
  Fig7c <- copy(b)
  tableColOrder <- c("Scenario","Region","Year","Variable","Value")
  Fig7c <- Fig7c[,c("scenario","region_class","period","variable","value")]
  Fig7c[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig7c,c("scenario","region_class","period","variable","value"),tableColOrder)
  
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)
  p <- ggplot(b,aes(x=period))+facet_grid(rows = NULL,vars(scenario))+theme_my(panel.spacing = 3,rotate_x = 90)+ylab(unit)+ggtitle(title)
  p <- p+geom_bar_interactive(aes(y=positive,fill=variable,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variable,"\n",formatC(positive,big.mark = ",",format="fg")," Mt N2O per year"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p+geom_bar_interactive(aes(y=negative,fill=forcats::fct_rev(variable),tooltip = paste0(scenario,"\n",region_class," ",period,"\n",variable,"\n",formatC(negative,big.mark = ",",format="fg")," Mt N2O per year"),data_id=interaction(period)),position="stack",stat = "identity",width=5)
  p <- p + scale_fill_manual("Nitrous oxide source",values=var_colors)
  p <- p+theme(legend.position = "right")+guides(fill=guide_legend(ncol=1,title.position = "top", byrow = TRUE,reverse=FALSE))+xlab(NULL)
  p <- p + geom_hline(yintercept = 0,linetype="dotted")
  p <- p + scale_x_continuous(breaks = c(2005,2025,2050,2075,2100))
  p3 <- p
  
  combined <- girafe(code = print(p1 + p2 + p3 + plot_annotation(tag_levels = 'a') + plot_layout(guides = "keep",ncol = 1,heights = c(2,1,1))),
                     options = list(
                       opts_sizing(rescale = TRUE, width = .95),
                       opts_hover(girafe_css(css = "fill:NULL;cursor:pointer;", area = "fill:NULL;")),
                       opts_hover_inv(css = "opacity:0.1;"),
                       opts_selection(girafe_css(css = "fill:NULL;stroke:black"), only_shiny = FALSE, type="multiple",selected = NULL),
                       opts_tooltip(css = "background-color:white;padding:5px;border-radius:2px;border: black 1px solid;color:black;")
                     ),
                     width_svg = 13,
                     height_svg = 8.5
  )
  
  return(list(plot=combined,Fig7a=Fig7a,Fig7b=Fig7b,Fig7c=Fig7c))
  
  # combined <- p1 + p2 + p3 + plot_annotation(tag_levels = 'a')# + labs(caption = paste(Sys.Date()))
  # combined <- combined + plot_layout(guides = "keep",ncol = 1,heights = c(2,1,1))
  # 
  # writexl::write_xlsx(list(Fig7a = Fig7a,Fig7b = Fig7b,Fig7c = Fig7c),"Fig7_emis_brakedown.xlsx")
  # 
  # ggsave(filename = "Fig7_emis_brakedown.pdf",combined,width = 10,height = 7,scale=1.3)
  # ggsave(filename = "Fig7_emis_brakedown.png",combined,width = 10,height = 7,scale=1.3)
  
}

Fig8 <- function(rep,reg_sub_order) {
  
  ### Fig 8 Indicators Details
  
  scen_sel <- rev(c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability"))
  names(scen_sel) <- rev(c("Global-Inequality","Global-EnvirProt","Global-GHGPrice","Global-SustDemand","Global-Sustainability"))
  
  years <- c(2020,2050,2100)
  
  ### add forest area without timber to rep object
  var <- c("Resources|Land Cover|Forest|Managed Forest|+|NPI/NDC","Resources|Land Cover|Forest|Managed Forest|+|Afforestation","Resources|Land Cover|Forest|Natural Forest|+|Secondary Forest","Resources|Land Cover|Forest|Natural Forest|+|Primary Forest")
  b <- rep[variable %in% var,]
  save <- names(b)
  b <- b[,.(value=sum(value)),by = .(region_class,region,model,scenario,unit,period)]
  b[,variable:="Forest area without plantations"]
  b[,unit:="Change in Mha compared to 2020"]
  setcolorder(b,save)
  b <- droplevels(b)
  b[, value:=value-value[period==2020], by = .(region_class,region,model,scenario,variable,unit)]
  rep <- rbind(rep,b)
  
  var <- c("Resources|Water|Withdrawal|Agriculture","SDG|SDG15|Industrial and intentional biological fixation of N","Forest area without plantations")
  names(var) <- c("atop('Agricultural water use (SDG 6)','km'^3*' yr'^-1*' (change compared to 2020)')","atop('Nitrogen fixation (SDG 15)','Mt N yr'^-1*' (change compared to 2020)')","atop('Forest area without plantations (SDG 15)',' Mha (change compared to 2020)')")
  
  b <- rep[scenario %in% scen_sel & variable %in% var & period %in% years & region != "World",]
  b <- b[,.(value=sum(value)),by = .(region_class,model,scenario,variable,unit,period)]
  b$scenario <- factor(b$scenario,levels = scen_sel,labels = names(scen_sel))
  b <- droplevels(b)
  b$region_class <- factor(b$region_class,levels = rev(c(reg_sub_order,"World")))
  b$variable <- factor(b$variable,levels = var,labels = names(var))
  b[,value:=value-value[period==2020],by=.(region_class,model,scenario,variable,unit)]
  b <- b[period != 2020,]
  Fig8 <- copy(b)
  tableColOrder <- c("Scenario","Region","Year","Variable","Unit","Value")
  Fig8 <- Fig8[,c("scenario","region_class","period","variable","unit","value")]
  Fig8[,value:=formatC(value,big.mark = ",",format="fg")]
  setnames(Fig8,c("scenario","region_class","period","variable","unit","value"),tableColOrder)
  
  levels(Fig8$Variable) <- c("Ag. water use","Nitrogen fixation","Forest area")
  b$positive <- ifelse(b$value >= 0, b$value, 0)
  b$negative <- ifelse(b$value < 0, b$value, -1e-36)
  
  p <- ggplot(b,aes(y=as.factor(scenario)))+facet_grid(vars(period),vars(variable),scales = "free",space="free_y",labeller = label_parsed)+theme_my(rotate_x = FALSE)+ylab(NULL)
  p <- p+geom_bar_interactive(aes(x=positive,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(positive,big.mark = ",",format="fg")),data_id=interaction(scenario)),position="stack",stat = "identity",width = 0.7)
  p <- p+geom_bar_interactive(aes(x=negative,fill=region_class,tooltip = paste0(scenario,"\n",region_class," ",period,"\n",formatC(negative,big.mark = ",",format="fg")),data_id=interaction(scenario)),position="stack",stat = "identity",width = 0.7)
  # p <- p+geom_bar(aes(x=positive,fill=region_class),position="stack",stat = "identity",width = 0.7)
  # p <- p+geom_bar(aes(x=negative,fill=region_class),position="stack",stat = "identity",width = 0.7)
  p <- p + geom_vline(xintercept = 0,linetype="dotted")
  p <- p + scale_fill_brewer("Scenario",palette = "Set2",direction = "-1")
  p <- p + stat_summary(fun = "sum", colour = "black", size = 2, geom = ggiraph:::GeomInteractivePoint, mapping = aes(tooltip = after_stat(paste0("World 2050\n",formatC(x,big.mark = ",",format="fg")," Gt CO2eq per year")),data_id=period,group = scenario,x=value))
  #p <- p + stat_summary(fun = "sum", colour = "black", size = 1, geom = "point", mapping = aes(group = scenario,x=value))
  p <- p+theme(legend.position = "bottom")+theme(strip.text.x = element_text(size = 10))+guides(fill=guide_legend("Scenario",nrow=1,title.position = "left", byrow = TRUE,reverse=TRUE))+xlab(NULL)#+scale_y_discrete(expand = expansion(mult = c(0.1,0.1)))+scale_x_continuous(expand = expansion(mult = c(0.1,0.1)))#+scale_x_continuous(breaks = c(-500,-250,0,250,500))# + labs(caption = paste(Sys.Date()))
  
  combined <- girafe(ggobj = p,
                     options = list(
                       opts_sizing(rescale = TRUE, width = .95),
                       opts_hover(girafe_css(css = "fill:NULL;cursor:pointer;", area = "fill:NULL;")),
                       opts_hover_inv(css = "opacity:0.1;"),
                       opts_selection(girafe_css(css = "fill:NULL;stroke:black"), only_shiny = FALSE, type="multiple",selected = NULL),
                       opts_tooltip(css = "background-color:white;padding:5px;border-radius:2px;border: black 1px solid;color:black;")
                     ),
                     width_svg = 12,
                     height_svg = 6
  )
  
  return(list(plot=combined,Fig8=Fig8))
  
  
  # writexl::write_xlsx(list(Fig8 = Fig8),"Fig8_indicator_detail.xlsx")
  # 
  # ggsave(filename = "Fig8_indicators_detail.pdf",p,width = 10,height = 5,scale=1.2)
  # ggsave(filename = "Fig8_indicators_detail.png",p,width = 10,height = 5,scale=1.2)
}

FigS3 <- function(rep,reg_sub_order) {
  
  #### Supplementary Information Figures
  
  ### FigS3 EAT-Lancet Diet
  
  scen_sel <- c("Global-Sustainability")
  
  var <- c("Nutrition|Calorie Intake")
  EAT_glo <- rep[scenario %in% scen_sel & variable %in% var & period == 2050 & region == "World",]
  
  var <- c("Population")
  pop <- rep[scenario %in% scen_sel & variable %in% var & period == 2020 & region != "World",]
  
  var <- c(
    "Nutrition|Calorie Intake|Crops|+|Cereals",
    "Nutrition|Calorie Intake|Crops|+|Oil crops",
    "Nutrition|Calorie Intake|Crops|Other crops|+|Fruits Vegetables Nuts",
    "Nutrition|Calorie Intake|Crops|Other crops|+|Potatoes",
    "Nutrition|Calorie Intake|Crops|Other crops|+|Pulses",
    "Nutrition|Calorie Intake|Crops|Other crops|+|Tropical roots",
    "Nutrition|Calorie Intake|Livestock products|+|Dairy",
    "Nutrition|Calorie Intake|Livestock products|+|Eggs",
    "Nutrition|Calorie Intake|Livestock products|+|Monogastric meat",
    "Nutrition|Calorie Intake|Livestock products|+|Poultry meat",
    "Nutrition|Calorie Intake|Livestock products|+|Ruminant meat",
    "Nutrition|Calorie Intake|Secondary products|+|Alcoholic beverages",
    "Nutrition|Calorie Intake|Secondary products|+|Brans",
    "Nutrition|Calorie Intake|Secondary products|+|Molasses",
    "Nutrition|Calorie Intake|Secondary products|+|Oils",
    "Nutrition|Calorie Intake|Secondary products|+|Sugar",
    "Nutrition|Calorie Intake|+|Fish")
  
  names(var) <- gsub("Nutrition\\|Calorie Intake\\|","",var)
  names(var) <- gsub("\\|Other crops","",names(var))
  names(var) <- gsub("\\+\\|","",names(var))
  
  fill <- c(brewer.pal(6,"Oranges"),brewer.pal(5,"Purples"),brewer.pal(5,"Greens"),"#3182bd")
  
  b <- rep[scenario %in% scen_sel & variable %in% var & period == 2020 & region != "World",]
  b <- merge(b,pop[,c("model","scenario","region","value")],by = c("model","scenario","region"),)
  setnames(b, "value.x", "value")
  setnames(b, "value.y", "pop")
  
  b$variable <- factor(b$variable,levels = var,labels = names(var))
  b <- b[order(-ave(b$value, b$region, FUN=sum),model,scenario,variable,),]
  
  b$ymax <- ave(b$value, b$region, FUN=cumsum)
  b$ymin <- b$ymax - b$value
  
  # Calculate the future positions on the x axis of each bar (left border, central position, right border)
  b$right <- ave(b$pop, b$variable, FUN=cumsum)
  b$left <- b$right - b$pop
  b$mean <- b$left + (b$right - b$left) / 2
  mag <- b
  
  p <- ggplot(mag, aes(ymin = 0)) +
    geom_rect_interactive(aes(xmin=left, xmax = right, ymax = ymax, ymin = ymin, fill = variable,data_id=variable,tooltip=paste0(region,"\n",variable,"\n",formatC(ymax-ymin,big.mark = ",",format="fg")," kcal/cap/day","\n",formatC(right-left,big.mark = ",",format="fg")," million people")),color="white") +
    xlab("Population (million)")+scale_fill_manual("Calorie Source",values=fill, guide = guide_legend(reverse = TRUE)) +ylab("kcal per capita per day")
  
  p1 <- p + geom_hline(data=EAT_glo,aes(yintercept=value),color="red",size=2,alpha=0.5) + geom_text(data=EAT_glo,aes(7756,value, vjust = -1,hjust="inward"),label="EAT-Lancet recommendation",fontface = "bold") + geom_label_repel(data=mag[variable==names(var[length(var)])],aes(x=mean,y=ymax,label=region),size=3,direction = "y",max.overlaps = 10,nudge_y=10)+theme_classic()+theme(legend.position="right")+scale_x_continuous(expand = c(0,0),limits = c(0,7820),breaks=c(1000,2000,3000,4000,5000,6000,7000,7756))+scale_y_continuous(limits = c(0,3000),expand = c(0,0),breaks = c(0,500,1000,1500,2000,2500,3000))
  p1 <- p1 + geom_text(x=4000,y=3000,vjust="inward",hjust=0.5,label="Calorie intake in 2020",fontface = "bold")
  
  b <- rep[scenario %in% scen_sel & variable %in% var & period == 2050 & region == "World",]
  b$variable <- factor(b$variable,levels = var,labels = names(var))
  b <- b[order(-ave(b$value, b$region, FUN=sum),model,scenario,variable,),]
  
  b$ymax <- ave(b$value, b$region, FUN=cumsum)
  b$ymin <- b$ymax - b$value
  
  b$region <- "EAT-Lancet"
  
  EAT <- b
  
  p <- ggplot(EAT, aes(xmin = 0,xmax=100)) +
    geom_rect_interactive(aes(ymax = ymax, ymin = ymin, fill = variable,data_id=variable,tooltip = paste0("World","\n",variable,"\n",formatC(ymax-ymin,big.mark = ",",format="fg")," kcal/cap/day")),color="white") +
    xlab("")+scale_fill_manual("Calorie Source",values=fill, guide = guide_legend(reverse = TRUE)) +# + scale_fill_brewer("Calorie Source",palette = "Paired",direction = "-1") +
    ylab(NULL) + geom_label_repel(data=EAT[variable==names(var[length(var)])],aes(x=50,y=ymax,label=region),size=2.5,direction = "y",max.overlaps = 10,nudge_y = 10)+theme_classic()+theme(legend.position="right")+scale_x_continuous(expand = c(0,0))+scale_y_continuous(limits = c(0,3000),expand = c(0,0),breaks = c(0,500,1000,1500,2000,2500,3000))
  p2 <- p + theme(axis.text.x=element_blank(),axis.ticks.x=element_blank(),axis.text.y=element_blank(),axis.ticks.y=element_blank(),axis.line = element_blank())
  
  combined <- girafe(code = print(p1 + p2 + plot_annotation(tag_levels = 'a') + plot_layout(guides = "collect",nrow = 1,heights = c(1,1),widths = c(0.9,0.1)) & theme(legend.position = "right")),
                     options = list(
                       opts_sizing(rescale = TRUE, width = .95),
                       opts_hover(girafe_css(css = "fill:NULL;cursor:pointer;", area = "fill:NULL;")),
                       opts_hover_inv(css = "opacity:0.1;"),
                       opts_selection(girafe_css(css = "fill:NULL;stroke:black"), only_shiny = FALSE, type="multiple",selected = NULL),
                       opts_tooltip(css = "background-color:white;padding:5px;border-radius:2px;border: black 1px solid;color:black;")
                     ),
                     width_svg = 11,
                     height_svg = 5
  )
  
  return(list(plot=combined))
  
  # combined <- p1 + p2 + plot_annotation(tag_levels = 'a')
  # combined <- combined + plot_layout(guides = "collect",nrow = 1,heights = c(1,1),widths = c(0.9,0.1)) & theme(legend.position = "right")
  # ggsave(filename = "FigS3_kcal_pop.pdf",combined,width = 11,height = 6,scale=1)
  # ggsave(filename = "FigS3_kcal_pop.png",combined,width = 11,height = 6,scale=1)
  
}


#load data
convert=FALSE
if(convert) {
  rev <- "HR_LAMA91"
  #load data as data.table
  rep <- as.data.table(readRDS(paste0("report_WP4_",rev,".rds")))
  #split up the scenario name by "_". You have to adjust this for your scenario name pattern
  rep[, c("scenario") := gsub(paste0(rev,"_"),"Global-",scenario)]
  rep[, c("scenario") := gsub(paste0("-Inequality-"),"-",scenario)]
  rep$scenario <- factor(rep$scenario)
  saveRDS(rep,"report.rds",compress = "xz")

  hist <- read.quitte("validation.mif")
  hist <- as.data.table(hist)
  hist <- na.omit(hist)
  hist <- droplevels(hist)
  saveRDS(hist,"validation2.rds")
}
rep <- readRDS("report.rds")
rep[,varunit:=paste(variable,unit,sep=":")]
write.csv(unique(rep$varunit),"vars_magpie.csv",row.names = F,quote = F)
rep[,varunit:=NULL]

hist <- readRDS(paste0("validation2.rds"))
hist[,varunit:=paste(model,scenario,variable,unit,sep=":")]
write.csv(unique(hist$varunit),"vars_hist.csv",row.names = F,quote = F)
hist[,varunit:=NULL]

#mapping for regional aggregation
map <- read.csv("mapping_regions.csv")
rep <- merge(rep,map)
hist <- merge(hist,map)

#region order
reg_sub_order <- unique(map$region_class)[-6]

#scen selection and order
scen_sel <- c("Global-Inequality","Global-Sustainability")

```

Landing
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### About {data-height=750}

**This digital tool provides interactive figures and data tables for the following paper**

Humpenöder, F., Popp, A., Schleussner, C.-F., Orlov, A., Windisch, M.G., Menke, I., Pongratz, J., Havermann, F., Thiery, W., Luo, F., et al. (2022). Overcoming global inequality is critical for land-based mitigation in line with the Paris Agreement. Nat Commun 13, 7453. 10.1038/s41467-022-35114-7.

https://www.nature.com/articles/s41467-022-35114-7

**Abstract**

Transformation pathways for the land sector in line with the Paris Agreement depend on the assumption of globally implemented greenhouse gas (GHG) emission pricing, and in some cases also on inclusive socio-economic development and sustainable land-use practices. In such pathways, the majority of GHG emission reductions in the land system is expected to come from low- and middle-income countries, which currently account for a large share of emissions from agriculture, forestry and other land use (AFOLU). However, in low- and middle-income countries the economic, financial and institutional barriers for such transformative changes are high. Here, we show that if sustainable development in the land sector remained highly unequal and limited to high-income countries only, global AFOLU emissions would remain substantial throughout the 21st century. Our model-based projections highlight that overcoming global inequality is critical for land-based mitigation in line with the Paris Agreement. While also a scenario purely based on either global GHG emission pricing or on inclusive socio-economic development would achieve the stringent emissions reductions required, only the latter ensures major co-benefits for other Sustainable Development Goals, especially in low- and middle-income regions. 

### Further information {data-height=250}
Contact for interactive figures: Florian Humpenöder (humpenoeder@pik-potsdam.de)

Interactive figures are realized with the R packages [ggplot2](https://ggplot2.tidyverse.org), [patchwork](https://patchwork.data-imaginist.com/index.html) and [ggiraph](https://davidgohel.github.io/ggiraph/index.html).

The website is generated from an Rmd file using the R [flexdashboard](https://pkgs.rstudio.com/flexdashboard/) package.


Column {data-width=300}
-----------------------------------------------------------------------

### Information on the LAMACLIMA project

LAMACLIMA – LAnd MAnagement for CLImate Mitigation and Adaptation

The project aims to investigate how changes in land cover and land management can help to meet the mitigation and adaptation objectives of the Paris Agreement, as well as the Sustainable Development Goals. 

https://climateanalytics.org/projects/lamaclima/

Overview {data-navmenu="GHG emissions"}
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### Figure 6

```{r echo=FALSE}
x <- Fig6(rep,reg_sub_order)
x[["plot"]]
```


Column {data-width=300 .tabset .tabset-fade}
-----------------------------------------------------------------------

### Information

Figure 6: AFOLU GHG emissions over time for five scenarios. a) shows results at global level over time. b) shows results at regional level for 2020, 2050 and 2100. Black solid lines in a) and dots in b) show the net effect across AFOLU GHG emissions. Grey solid lines in a) show median values across 1.5°C compatible illustrative mitigation pathways (LD, Ren, SP) from the IPCC AR6 Scenarios Database24. CO2 includes emissions from land-use change such as deforestation and conversion of non-forest ecosystems, as well as removals from re/afforestation and natural succession on abandoned agricultural land. CH4 includes emissions from enteric fermentation, animal waste management and rice cultivation. N2O includes emissions from agricultural soils (fertilizer application) and animal waste management. N2O and CH4 emissions have been converted into CO2 equivalents using IPCC AR6 GWP100 factors of 273 and 27, respectively.

### Data Fig6a

**AFOLU GHG emissions (Gt CO~2~ eq per year)**

```{r echo=FALSE}
datatable(x[["Fig6a"]], options = list(pageLength = 100))
```

### Data Fig6b

**AFOLU GHG emissions (Gt CO~2~ eq per year)**

```{r echo=FALSE}
datatable(x[["Fig6b"]], options = list(pageLength = 100))
```

Breakdown {data-navmenu="GHG emissions"}
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### Figure 7

```{r echo=FALSE}
x <- Fig7(rep,reg_sub_order)
x[["plot"]]
```


Column {data-width=300 .tabset .tabset-fade}
-----------------------------------------------------------------------

### Information

Figure 7: AFOLU GHG emission breakdown by source. Data is shown at global level for five scenarios. a) shows CO2 emissions and removals from land-use change and management. Carbon losses consist of emissions from deforestation, conversion of non-forest ecosystems, drained peatlands and wood harvest. Carbon gains consist of carbon storage in wood products, re/afforestation (Aff CO2-price and Aff NDC), timber plantations and regrowth of natural vegetation (secondary forests and other natural land). The black line shows the net effect of carbon losses and carbon gains at global level. b) shows CH4 emissions from agriculture. c) shows N2O emissions from agriculture. Further details on AFOLU GHG emission sources and sinks are provided in the methods section and in Table S1.

### Data Fig7a

**CO~2~ emissions from land-use change (Gt CO~2~ per year)**

```{r echo=FALSE}
datatable(x[["Fig7a"]], options = list(pageLength = 100))
```

### Data Fig7b

**CH~4~ emissions from agriculture (Mt CH~4~ per year)**

```{r echo=FALSE}
datatable(x[["Fig7b"]], options = list(pageLength = 100))
```

### Data Fig7c

**N~2~O emissions from agriculture (Mt N~2~O per year)**

```{r echo=FALSE}
datatable(x[["Fig7c"]], options = list(pageLength = 100))
```


Synergies & Trade-offs
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### Figure 3

```{r echo=FALSE}
x <- Fig3(rep,reg_sub_order)
x[["plot"]]
#htmlwidgets::saveWidget(x[["plot"]],"Fig3.html")
```


Column {data-width=300 .tabset .tabset-fade}
-----------------------------------------------------------------------

### Information

Figure 3: Summary of results for five scenarios and six indicators (mapped to SDGs) for 2050 at regional and global level. a) CO2 emissions and removals from land-use change and management, b) N2O emissions from agriculture, c) CH4 emissions from agriculture, d) change in forest area without plantations, e) nitrogen fixation and f) agricultural water use. Black dots show the global net effect for each indicator and scenario. Blue vertical lines show the respective value in year 2020. Percentage labels in panels b, c, e and f show the relative change of the global indicator level for each scenario compared to the Global-Inequality scenario. For panels a and d it is not meaningful to show percentage changes because these indicators can, by definition, turn net-zero (in case of divergent regional dynamics) or net-negative at global level.

### Data Fig3a

**CO~2~ emissions from land-use change (Gt CO~2~ eq per year)**

```{r echo=FALSE}
datatable(x[["Fig3a"]], options = list(pageLength = 100))
```

### Data Fig3b

**Annual N~2~O emissions from agriculture (Gt CO~2~ eq per year)**

```{r echo=FALSE}
datatable(x[["Fig3b"]], options = list(pageLength = 100))
```

### Data Fig3c

**Annual CH~4~ emissions from agriculture (Gt CO~2~ eq per year)**

```{r echo=FALSE}
datatable(x[["Fig3c"]], options = list(pageLength = 100))
```

### Data Fig3d

**Forest area without plantations (Change in Mha wrt to 2020)**

```{r echo=FALSE}
datatable(x[["Fig3d"]], options = list(pageLength = 100))
```

### Data Fig3e

**Nitrogen fixation (Mt N per year)**

```{r echo=FALSE}
datatable(x[["Fig3e"]], options = list(pageLength = 100))
```

### Data Fig3f

**Agricultural water use (km^3^ per year)**

```{r echo=FALSE}
datatable(x[["Fig3f"]], options = list(pageLength = 100))
```

Overview {data-navmenu="Land"}
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### Figure 4

```{r echo=FALSE}
x <- Fig4(rep,reg_sub_order)
x[["plot"]]
```


Column {data-width=300 .tabset .tabset-fade}
-----------------------------------------------------------------------

### Information

Figure 4: Land-use change of major land types over time compared to 2020 for five scenarios. a) show results at global level over time. b) shows results at regional level for 2050 and 2100. Cropland includes food, non-food, and feed crops. Bioenergy includes 2nd generation bioenergy (fast growing grasses and trees such as miscanthus and poplar). Pasture includes rangeland and managed pasture areas. Timber refers to timber plantations for wood production. Aff NDC includes prescribed re/afforestation according to National Determined Contributions (NDCs) towards the objectives of the Paris Agreement. Aff CO2-Price refers to endogenous re/afforestation based on the underlying scenario-specific CO2 price trajectory. Secondary forest includes modified and regrown forest, e.g., following wood harvest or cropland abandonment. Primary forest is intact forest without signs of human intervention. Urban land includes built-up area. Other natural land is a residual category, which includes among others non-forest ecosystems, deserts, and shrublands.

### Data Fig4a

**Land-use change (Change in Mha wrt to 2020)**

```{r echo=FALSE}
datatable(x[["Fig4a"]], options = list(pageLength = 100))
```

### Data Fig4b

**Land-use change (Change in Mha wrt to 2020)**

```{r echo=FALSE}
datatable(x[["Fig4b"]], options = list(pageLength = 100))
```

Productivity {data-navmenu="Land"}
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### Figure 5

```{r echo=FALSE}
x <- Fig5(rep,reg_sub_order)
x[["plot"]]
```


Column {data-width=300 .tabset .tabset-fade}
-----------------------------------------------------------------------

### Information

Figure 5: Land-use intensity factor τ. Data is shown at regional level for five scenarios. The τ factor reflects the degree of crop yield amplification caused by human activities. A duplication of τ implies a doubling of crop yields under fixed environmental conditions. See Figure S10 for validation data.

### Data Fig5

**Land-use intensity factor τ**

```{r echo=FALSE}
datatable(x[["Fig5"]], options = list(pageLength = 100))
```

Water & Nitrogen
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### Figure 8

```{r echo=FALSE}
x <- Fig8(rep,reg_sub_order)
x[["plot"]]
```


Column {data-width=300 .tabset .tabset-fade}
-----------------------------------------------------------------------

### Information

Figure 8: Change in agricultural water use, nitrogen fixation and forest area compared to 2020 for five scenarios at regional level for the years 2050 and 2100. The black dot indicates the net effect at global level.

### Data Fig8

```{r echo=FALSE}
datatable(x[["Fig8"]], options = list(pageLength = 100))
```

Scenario Narratives
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### Figure 1

```{r echo=FALSE}
#, fig.width = 12, fig.height= 8
x <- Fig1(rep,reg_sub_order)
x[["plot"]]
x[["plot"]] <- girafe_options(x[["plot"]],opts_sizing(rescale = FALSE))
#htmlwidgets::saveWidget(x[["plot"]],"Fig1.html",)
#https://stackoverflow.com/questions/71695258/how-to-resize-html-widget-using-savewidget-in-htmlwidgets-r-reopened-question
```


Column {data-width=300 .tabset .tabset-fade}
-----------------------------------------------------------------------

### Information

Figure 1: Key socio-economic drivers and developments. Data is shown at regional and global level for the two main scenarios Global-Inequality (same for Global-GHGPrice and Global-EnvirProt) and Global-Sustainability (same for Global-SustDemand). a) Population, b) Per-capita income, c) Per-capita calorie intake, d) Total demand for crops (including food and feed) and livestock products in million-ton dry matter per year, e) Prevalence of underweight (body mass index < 18.5) and f) Prevalence of obesity (BMI > 30). Details on prevalence of underweight and obesity can be found in Table S1.

### Data Fig1a

**Population (billion people)**

```{r echo=FALSE}
datatable(x[["Fig1a"]], options = list(pageLength = 100))
```

### Data Fig1b

**Income (USD05 PPP cap^-1^ yr^-1^)**

```{r echo=FALSE}
datatable(x[["Fig1b"]], options = list(pageLength = 100))
```

### Data Fig1c

**Per-capita calorie intake (kcal cap^-1^ day^-1^)**

```{r echo=FALSE}
datatable(x[["Fig1c"]], options = list(pageLength = 100))
```

### Data Fig1d

**Total demand for crops and livestock products (Mt DM yr^-1^)**

```{r echo=FALSE}
datatable(x[["Fig1d"]], options = list(pageLength = 100))
```

### Data Fig1e

**Prevalence of underweight BMI < 18.5 (million people)**

```{r echo=FALSE}
datatable(x[["Fig1e"]], options = list(pageLength = 100))
```

### Data Fig1f

**Prevalence of obesity BMI > 30 (million people)**

```{r echo=FALSE}
datatable(x[["Fig1f"]], options = list(pageLength = 100))
```

Overview {data-navmenu="Food System"}
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### Figure 2

```{r echo=FALSE}
x <- Fig2(rep,reg_sub_order)
x[["plot"]]
```


Column {data-width=300 .tabset .tabset-fade}
-----------------------------------------------------------------------

### Information

Figure 2: Per-capita calorie supply. The sum of crops, livestock products, fish and secondary products reflects food intake, which together with food waste sums up to food supply. Percentage values indicate the share of food waste in total food supply. Data is shown at regional level for the two main scenarios Global-Inequality (same for Global-GHGPrice and Global-EnvirProt) and Global-Sustainability (same for Global-SustDemand). 

### Data Fig2

**Per-capita calorie supply (kcal cap^-1^ day^-1^)**

```{r echo=FALSE}
datatable(x[["Fig2"]], options = list(pageLength = 100))
```

Details {data-navmenu="Food System"}
=======================================================================

Column {data-width=700}
-----------------------------------------------------------------------

### Figure S3

```{r echo=FALSE}
x <- FigS3(rep,reg_sub_order)
x[["plot"]]
```


Column {data-width=300 .tabset .tabset-fade}
-----------------------------------------------------------------------

### Information

Figure S3: Per-capita calorie intake. a) shows regional data for 2020. The height of each rectangle shows per-capita kcal intake, the width shows the population of the region, so that the area of the rectangles refers to the total calorie intake for each region. b) shows EAT-Lancet recommendations (planetary health diet), aggregated to global level. In the Global-Sustainability scenario, all regions converge towards the EAT-Lancet recommendation by 2050.