---
title: "The energy and carbon inequality corridor for a 1.5 degree compatible and just Europe"
author:
  - Ingram S. Jaccard:
      email: jaccard@pik-potsdam.de
      institute: [PIK]
      correspondence: false
  - Peter-Paul Pichler:
      email: pichler@pik-potsdam.de
      institute: [PIK]
      correspondence: false
  - Johannes Többen:
      email: toebben@pik-potsdam.de
      institute: [PIK, GWS]
      correspondence: false
  - Helga Weisz:
      email: weisz@pik-potsdam.de
      institute: [PIK, HU]
      correspondence: false
institute:
  - PIK: Social Metabolism and Impacts, Potsdam Institute for Climate Impact Research, Member of the Leibniz Association, PO Box 60 12 03, Potsdam, 14412, Germany
  - HU: Department of Cultural History & Theory and Department of Social Sciences, Humboldt University Berlin, Unter den Linden 6, Berlin, 10117, Germany
  - GWS: Gesellschaft für Wirtschaftliche Strukturforschung (GWS) mbH, Heinrichstraße 30, 49080 Osnabrück, Germany
output: 
    bookdown::word_document2:
      fig_caption: yes
      reference_docx: "../templates/template.docx" # Insert path for the DOCX file
      pandoc_args:
      - --lua-filter=../templates/scholarly-metadata.lua
      - --lua-filter=../templates/author-info-blocks.lua
      - --lua-filter=../templates/pagebreak.lua
bibliography: references.bib
csl: "../templates/vancouver.csl" # Insert path for the bib-style
keywords: |
  keyword 1; keyword 2; keyword 3
highlights: |
  These are the highlights. 
abstract: |
---

<!-- This is the format for text comments that will be ignored during renderings. Do not put R code in these comments because it will not be ignored. -->

<!-- With the following code you can access and display values from the yaml header above. -->

Keywords: `r rmarkdown::metadata$keywords`

Highlights: `r rmarkdown::metadata$highlights`

Abstract: The call for a decent life for all within planetary limits poses a dual challenge: Provide all people with the essential resources needed to live well and, collectively, to not exceed the source and sink capacity of the biosphere to sustain human societies. In this paper, we examined the space of possible distributions of household energy and carbon footprints that satisfy both minimal energy requirements for a decent living and maximum supply of decarbonized energy to achieve the 1.5 degree target in 2050 for the populations of 28 European countries. Estimates for a range of minimum energy requirements for a decent life, as well as estimates for the maximum available energy supply, were taken from the 1.5 degree scenario literature. The maximum available energy supply is determined by how quickly the emission intensity of the global energy system can be reduced (expansion of renewables, efficiency improvements, etc.) so that the global energy mix remains within the specified emission reduction pathways. [not finished yet]

<!-- The following code chunk defines some general settings how code chunks should behave. -->

```{r setup, echo = FALSE, include = FALSE, message = FALSE}
knitr::opts_chunk$set(
  collapse = TRUE,
  warning = FALSE,
  message = FALSE,
  echo = FALSE,
  comment = "#>",
  fig.path = "../figures/",
  dpi = 300
)

if (!require("pacman")) install.packages("pacman")
pacman::p_load(tidyverse,
               janitor,
               here,
               wbstats,
               ISOcodes,
               viridis,
               hrbrthemes,
               wesanderson,
               glue,
               ggridges,
               patchwork,
               kableExtra,
               readxl,
               flextable)

pal <- wes_palette("Cavalcanti1", 5, type = "discrete")
extrafont::loadfonts()

options(scipen=999)

```

```{r libraries}

library(here)
```

```{r functions}
# load data wrangling functions
source(here("analysis", "R", "wrangler_functions.R"))

```

```{r load-data, include=FALSE}

## load result data for EU deciles
eu_q_count = 10

# summary countries aggregated by country quintiles and eu ntile
dat_country_summary_by_cquint_and_euntile = get_country_summary_by_cquint_and_euntile(eu_q_count)
# pivot to long format for plotting and attach readable indicator names
cols_ex = c("year", "iso2", "quint", "eu_q_rank")
pdat_country_summary_by_cquint_and_euntile =
  pivot_results_longer_adorn(dat_country_summary_by_cquint_and_euntile, cols_ex)

# summary of countries by EU quantile without sectoral resolution
dat_country_summary_by_eu_ntile = get_country_summary_by_eu_ntile(eu_q_count)
# pivot to long format for plotting and attach readable indicator names
cols_ex = c("year", "iso2", "eu_q_rank")
pdat_country_summary_by_eu_ntile = 
  pivot_results_longer_adorn(dat_country_summary_by_eu_ntile, cols_ex)

# summary of countries by country quintile with aggregate sectoral resolution
dat_sector_summary_by_country_quintile = get_sector_summary_by_country_quintile(eu_q_count)
# pivot to long format for plotting and attach readable indicator names
cols_ex = c("year", "iso2", "quint", "eu_q_rank", "sector_agg_id")
pdat_sector_summary_by_country_quintile =
  pivot_results_longer_adorn(dat_sector_summary_by_country_quintile, cols_ex)

# summary of eu ntile with aggregate sectoral resolution
dat_sector_summary_by_eu_ntile = get_sector_summary_by_eu_ntile(eu_q_count)
# pivot to long format for plotting and attach readable indicator names
cols_ex = c("year", "eu_q_rank", "sector_agg_id")
pdat_sector_summary_by_eu_ntile =
  pivot_results_longer_adorn(dat_sector_summary_by_eu_ntile, cols_ex)

```

<!-- The actual document text starts here: -->

# Introduction

Decarbonizing the energy system in accordance with the Paris Accord requires a deep transformation of both the supply and the demand side [@grubler_low_2018]. On both sides, however, necessary transformation is restricted by different factors. On the supply side, there exist economic and physical upper limits of how much energy can be provided from renewable sources on the one hand, and how much CO2 removal infrastructure is used to compensate for remaining emissions from fossil fuels on the other. On the demand side [@creutzig_towards_2018], by contrast, there are lower limits to how much energy is minimally required for a decent life [@grubler_low_2018 @millward-hopkins_providing_2020], depending on different assumptions about production-consumption infrastructures and service provision [@creutzig_towards_2019], as well as the prevalent social ideas about what constitutes decent living [@rao_energy_2019 @millward-hopkins_providing_2020]. Maximum possible energy supply and minimum necessary energy demand describe the corridor in which the simultaneous achievement of climate targets and a decent living for all is possible and, at the same time, restricts the distribution of available energy services among the population. If this dual objective is taken seriously in European climate policy, then there are practical limits to how unequal the society of the future can be, which go beyond the purely political. In fact, a limited energy supply creates an obvious, if rarely acknowledged, zero-sum game where energetic over-consumption by some has to be compensated by less consumption by others. 

The average household energy footprint of European citizens was around 170 GJ per capita in 2015 [@eurostat_eurostat_nodate-3 @stadler_exiobase_2018] and the household carbon footprint around 7 tonnes CO2eq per capita in 2015 [@eurostat_eurostat_nodate-4]. However, the differences in average household energy and carbon footprints are large within and between different regions in Europe. Energy footprints ranged from less than 100 GJ per capita to over 300 GJ per capita [@oswald_large_2020], and carbon footprints from below 2.5 tonnes CO2eq per capita to 55 tonnes CO2eq per capita [@ivanova_unequal_2020]. Depending on the assumptions of different global mitigation scenarios, the average footprints likely need to be reduced to between 15.7 and 100 GJ per capita [@grubler_low_2018 @millward-hopkins_providing_2020], or 0.5 and 2.1 tCO2eq per capita [@akenji_1.5-degree_2019] by 2050, respectively.

In this paper, we assess under what conditions European energy inequality is compatible with the achievement of global climate goals and a decent standard of living, taking both inequality within and between European countries into account. To this end, we first construct household energy and carbon footprints for harmonized European expenditure deciles in 2015, combining data from EUROSTAT's Household Budget Survey (HBS) with the Environmentally-Extended Multi-Regional Input-Output (EE-MRIO) model EXIOBASE. We analyze the distribution of energy and carbon intensities across European expenditure deciles and consumption categories, and compare this current structure to a hypothetical situation where all European deciles use the best technology available in Europe. Finally we examine how the energy inequality across European household expenditure deciles would need to change to achieve the dual goal of climate protection and a decent standard of living for all. 

While the European Green Deal already recognizes that inequalities in income, energy infrastructure, energy consumption and greenhouse gas emissions lead to different responsibilities and capacities in achieving the emission savings targets, a quantification of the attainable corridor for 1.5 compatible and just transition in Europe is missing in the literature.

# Materials and methods

## Income-stratified national household energy and carbon footprints

We used the EE-MRIO model EXIOBASE for 2015 (version3, industry-by-industry) [@stadler_exiobase_2018] and the European household budget survey (HBS) macro-data from EUROSTAT for 2015 [@eurostat_database_nodate] to calculate income-stratified national household energy and carbon footprints (together denoted as environmental footprints in this paper). The EUROSTAT HBS publishes mean household expenditure by income quintile, in purchasing power standard (PPS), by COICOP consumption category, country and year. We choose EXIOBASE as the EE-MRIO for this study because of its European focus, with nearly all countries in the EUROSTAT HBS also found as stand-alone countries in EXIOBASE (see SI - table 5), its detailed environmental extension data, and its year coverage. 

To integrate HBS data into EXIOBASE we created correspondence tables between the EXIOBASE sectors and the matching COICOP consumption categories used in HBS (see SI, p xx for details). We then used the relative shares of the COICOP consumption categories of each income quintile in the HBS to decompose the matching EXIOBASE national household final demand expenditure per sector and per income quintile. Using standard input-output techniques (see SI) we calculated ‘total’ (i.e. indirect supply chain) energy use and carbon intensities per EXIOBASE sector and multiplied them with the income-stratified EXIOBASE national household expenditure, to estimate the supply chain part of national household energy and carbon footprints by national income quintile. 

We used the energy use extensions ‘gross total energy use’ from EXIOBASE, which converts final energy consumption in the IEA energy balance data from the territorial to residence principle following SEEA energy accounting [@stadler_exiobase_2018] and the EXIOBASE GHG emission extensions CO2, CH4, N2O, SF6, HFCs and PFCs, from combustion, non-combustion, agriculture and waste, but not land-use change. Direct household energy use and carbon emissions are included in the environmental footprints. 

## European household expenditure deciles

To calculate European household expenditure deciles we first ranked the national income quintiles from the HBS of 28 European countries (in total 140 national quintiles) according to their mean expenditure in PPS and aggregated the result to 10 European expenditure groups. For brevity we call them expenditure deciles in the rest of the paper. Our coverage of European countries is limited to those included in both the EUROSTAT HBS data and EXIOBASE in 2015. This resulted in a country sample that includes the non-EU members Norway and Turkey, but excludes the EU members Italy and Luxembourg. 

## Units of analysis

The unit of analysis for our energy and carbon footprint calculations is the household. We normalized our results to average adult equivalent per household and per national decile because this is how the EUROSTAT HBS publishes its data. The first adult in the household is given a weight of 1.0, each adult thereafter 0.5, and each child 0.3 [@eurostat_description_2016]. 

For our calculations of attainable corridors for achieving the dual goal of climate protection and a decent standard of living for all, we adjusted the total per capita results from published 1.5 scenarios to adult equivalents in order to better compare with our environmental footprint estimates (see SI pp xx for details). Data on decarbonization scenarios, Minimum final energy use for a decent living are from Grubler et al. (2018) [@grubler_low_2018] and Millward-Hopkins et al. (2020) [@millward-hopkins_providing_2020], maximum final energy use compatible with the 1.5 degree target is from the IIASA scenario database [@riahi_shared_2017 @gea_gea_nodate].

As inequality measure we use the 10:10 ratio, i.e. the expenditure or the environmental footprint of the top European expenditure decile divided by that of the bottom European expenditure decile. Thus, an expenditure 10:10 ratio of 5 means that one adult equivalent in the top decile spent 5 times more on average than one adult equivalent in the bottom decile. 

## Computing maximum permissible inequality

Based on an hypothetical current best technology distribution across European household expenditure deciles, for each value combination of maximum energy supply from [xx] scenarios [@riahi_shared_2017 @gea_gea_nodate] and minimum energy use requirements from [@grubler_low_2018 @millward-hopkins_providing_2020], the maximum permissible inequality was calculated as a 10:10 ratio using the formula [insert formula]. The remaining global emissions budget to achieve the 1.5 degree target from the scenarios was allocated in proportion to population (equal per capita allocation). All data and procedures are described in detail in the supplementary information (SI).

# Results and discussion

## Environmental footprints are less unequal than expenditure levels

```{r ntiles-total}

p1 = pdat_country_summary_by_eu_ntile %>%
  filter(year == 2015, 
         indicator == "total_fd_me") %>%
  group_by(eu_q_rank) %>%
  summarise(value = sum(value)*0.000001,
            eu_ntile_name = first(eu_ntile_name)) %>%
  ggplot(aes(x=eu_ntile_name, y=value)) +
    geom_col(position = position_dodge(), fill=pal[1]) +
    #theme_ipsum() +
    theme_minimal() +
    theme(text=element_text(family="Liberation Sans Narrow")) +
    labs(x="", y="Expenditure (trn€)") +
    theme(axis.text.x = element_text(angle = 90)) +
    scale_x_discrete(labels = c("D01","D02","D03","D04","D05","D06","D07","D08","D09","D10")) 

p2 = pdat_country_summary_by_eu_ntile %>%
  filter(year == 2015, 
         indicator == "total_energy_use_tj") %>%
  group_by(eu_q_rank) %>%
  summarise(value = sum(value)*0.000001,
            eu_ntile_name = first(eu_ntile_name)) %>%
  ggplot(aes(x=eu_ntile_name, y=value)) +
    geom_col(position = position_dodge(), fill=pal[1]) +
    #theme_ipsum() +
    theme_minimal() +
    theme(text=element_text(family="Liberation Sans Narrow")) +
    labs(x="", y="Energy footprint (EJ)") +
    theme(axis.text.x = element_text(angle = 90)) +
    scale_x_discrete(labels = c("D01","D02","D03","D04","D05","D06","D07","D08","D09","D10")) 

p3 = pdat_country_summary_by_eu_ntile %>%
  filter(year == 2015, 
         indicator == "total_co2eq_kg") %>%
  group_by(eu_q_rank) %>%
  summarise(value = sum(value)*0.000000001,
            eu_ntile_name = first(eu_ntile_name)) %>%
  ggplot(aes(x=eu_ntile_name, y=value)) +
    geom_col(position = position_dodge(), fill=pal[1]) +
    #theme_ipsum() +
    theme_minimal() +
    theme(text=element_text(family="Liberation Sans Narrow")) +
    labs(x="", y="Carbon footprint (MtCO2eq)") +
    theme(axis.text.x = element_text(angle = 90)) +
    scale_x_discrete(labels = c("D01","D02","D03","D04","D05","D06","D07","D08","D09","D10")) 


p_top = p1 + p2 + p3

```

```{r ntiles-intensity-violin}


p1 = dat_country_summary_by_eu_ntile %>%
  filter(year == 2015) %>%
  ggplot(aes(x=factor(eu_q_rank), y=pe_co2eq_kg)) +
  geom_violin(aes(weight=total_fd_me), fill=pal[1], color=pal[1], alpha=0.5) +
  geom_point( alpha=0.3) +
  geom_segment(data=dat_country_summary_by_eu_ntile %>%
                 filter(year == 2015) %>%
                 group_by(eu_q_rank) %>%
                 summarise(pe_co2eq_kg = weighted.mean(pe_co2eq_kg,total_fd_me)),
               aes(y=pe_co2eq_kg, yend=pe_co2eq_kg, x=eu_q_rank-0.3, xend=eu_q_rank+0.3), size=1.5) +
    #theme_ipsum() +
    theme_minimal() +
    theme(text=element_text(family="Liberation Sans Narrow")) +
    labs(x="", y="Carbon intensity per expenditure (kgCO2eq/€)") +
    theme(axis.text.x = element_text(angle = 90)) +
    scale_x_discrete(labels = c("D01","D02","D03","D04","D05","D06","D07","D08","D09","D10")) 

p2 = dat_country_summary_by_eu_ntile %>%
  filter(year == 2015) %>%
  ggplot(aes(x=factor(eu_q_rank), y=pe_energy_use_mj)) +
  geom_violin(aes(weight=total_fd_me), fill=pal[1], color=pal[1], alpha=0.5) +
  geom_point( alpha=0.3) +
  geom_segment(data=dat_country_summary_by_eu_ntile %>%
                 filter(year == 2015) %>%
                 group_by(eu_q_rank) %>%
                 summarise(pe_energy_use_mj = weighted.mean(pe_energy_use_mj,total_fd_me)),
               aes(y=pe_energy_use_mj, yend=pe_energy_use_mj, x=eu_q_rank-0.3, xend=eu_q_rank+0.3), size=1.5) +
    #theme_ipsum() +
    theme_minimal() +
    theme(text=element_text(family="Liberation Sans Narrow")) +
    labs(x="", y="Energy intensity per expenditure (MJ/€)") +
    theme(axis.text.x = element_text(angle = 90)) +
    scale_x_discrete(labels = c("D01","D02","D03","D04","D05","D06","D07","D08","D09","D10")) 

dat3 = dat_country_summary_by_eu_ntile %>%
  filter(year == 2015) %>%
  mutate(intensity_e_c = total_co2eq_kg*0.001/total_energy_use_tj)

p3 = dat3 %>%
  ggplot(aes(x=factor(eu_q_rank), y=intensity_e_c)) +
  geom_violin(aes(weight=total_energy_use_tj), fill=pal[1], color=pal[1], alpha=0.5) +
  geom_point( alpha=0.3) +
  geom_segment(data=dat3 %>%
                 filter(year == 2015) %>%
                 group_by(eu_q_rank) %>%
                 summarise(intensity_e_c = weighted.mean(intensity_e_c,total_energy_use_tj)),
               aes(y=intensity_e_c, yend=intensity_e_c, x=eu_q_rank-0.3, xend=eu_q_rank+0.3), size=1.5) +
    #theme_ipsum() +
    theme_minimal() +
    theme(text=element_text(family="Liberation Sans Narrow")) +
    labs(x="", y="Carbon intensity per energy (gCO2eq/TJ)") +
    theme(axis.text.x = element_text(angle = 90)) +
    scale_x_discrete(labels = c("D01","D02","D03","D04","D05","D06","D07","D08","D09","D10")) 


p_bottom = p2 + p1 + p3

#p3 & theme(plot.margin =  margin(1, 1, 1, 1, "mm"))
```

```{r , fig.width=12, fig.height=8}

a = p_top / p_bottom + plot_annotation(tag_levels = 'a') +
  plot_layout(guides = 'collect')  & 
  theme(plot.margin = unit(c(0.25,0.25,0.25,0.25), "cm"),
        legend.position = 'bottom',
        axis.title.y = element_text(size=13, hjust = 0.5),
        axis.text.x = element_text(size = 12),
        axis.text.y = element_text(size = 12),
        legend.text = element_text(size=12),
        legend.title = element_text(size=13))

ggsave(here("analysis", "figures", "figure1.pdf"), device=cairo_pdf)

```

```{r values-in-text1}

# values in text

## inequality
exp = pdat_country_summary_by_eu_ntile %>%
  filter(year == 2015, 
         indicator == "total_fd_me") %>%
  group_by(eu_q_rank) %>%
  summarise(value = sum(value)*0.000001,
            eu_ntile_name = first(eu_ntile_name))

exp_10_10 = round((exp %>% filter(eu_q_rank == 10))$value/(exp %>% filter(eu_q_rank == 1))$value,digits = 1)

energy = pdat_country_summary_by_eu_ntile %>%
  filter(year == 2015, 
         indicator == "total_energy_use_tj") %>%
  group_by(eu_q_rank) %>%
  summarise(value = sum(value)*0.000001,
            eu_ntile_name = first(eu_ntile_name))

energy_10_10 = round((energy %>% filter(eu_q_rank == 10))$value/(energy %>% filter(eu_q_rank == 1))$value,digits = 1)

co2eq = pdat_country_summary_by_eu_ntile %>%
  filter(year == 2015, 
         indicator == "total_co2eq_kg") %>%
  group_by(eu_q_rank) %>%
  summarise(value = sum(value)*0.000000001,
            eu_ntile_name = first(eu_ntile_name))

co2eq_10_10 = round((co2eq %>% filter(eu_q_rank == 10))$value/(co2eq %>% filter(eu_q_rank == 1))$value,digits = 1)

## total per decile

exp_bottom_decile = round((exp %>% filter(eu_q_rank == 1))$value, digits = 1)

exp_top_decile = round((exp %>% filter(eu_q_rank == 10))$value, digits = 1)

energy_bottom_decile = round((energy %>% filter(eu_q_rank == 1))$value, digits = 1)

energy_top_decile = round((energy %>% filter(eu_q_rank == 10))$value, digits = 1)

co2eq_bottom_decile = round((co2eq %>% filter(eu_q_rank == 1))$value, digits = 1)

co2eq_top_decile = round((co2eq %>% filter(eu_q_rank == 10))$value, digits = 1)


## per adult equivalent per decile

aeu = pdat_country_summary_by_eu_ntile %>%
  filter(year == 2015,
         indicator == "total_adult_eq") %>%
  group_by(eu_q_rank) %>%
  summarise(value = sum(value),
            eu_ntile_name = first(eu_ntile_name)) 

exp_pae = exp %>%
  rename(total_fd_me = value) %>%
  left_join(aeu, by = c("eu_q_rank", "eu_ntile_name")) %>%
  mutate(pae_fd_e = (total_fd_me/value)*1000000000000)

fd_pae_bottom_decile = round((exp_pae %>% filter(eu_q_rank == 1))$pae_fd_e, digits = 0)

fd_pae_top_decile = round((exp_pae %>% filter(eu_q_rank == 10))$pae_fd_e, digits = 0)

energy_pae = energy %>%
  rename(total_energy_use_tj = value) %>%
  left_join(aeu, by = c("eu_q_rank", "eu_ntile_name")) %>%
  mutate(pae_energy_use_gj = (total_energy_use_tj/value)*1000000000)

energy_pae_bottom_decile = round((energy_pae %>% filter(eu_q_rank == 1))$pae_energy_use_gj, digits = 1)

energy_pae_top_decile = round((energy_pae %>% filter(eu_q_rank == 10))$pae_energy_use_gj, digits = 1)

co2eq_pae = co2eq %>%
  rename(total_co2eq_kg = value) %>%
  left_join(aeu, by = c("eu_q_rank", "eu_ntile_name")) %>%
  mutate(pae_co2eq_t = (total_co2eq_kg/value)*1000000)

co2eq_pae_bottom_decile = round((co2eq_pae %>% filter(eu_q_rank == 1))$pae_co2eq_t, digits = 1)

co2eq_pae_top_decile = round((co2eq_pae %>% filter(eu_q_rank == 10))$pae_co2eq_t, digits = 1)


## intensities

mean_energy_intens = dat_country_summary_by_eu_ntile %>%
                 filter(year == 2015) %>%
                 group_by(eu_q_rank) %>%
                 summarise(pe_energy_use_mj = weighted.mean(pe_energy_use_mj,total_fd_me))

mean_energy_intens_bottom_decile = round((mean_energy_intens %>% filter(eu_q_rank == 1))$pe_energy_use_mj, digits = 1)

mean_energy_intens_top_decile = round((mean_energy_intens %>% filter(eu_q_rank == 10))$pe_energy_use_mj, digits = 1)

mean_co2eq_of_energy_intens = dat_country_summary_by_eu_ntile %>%
  filter(year == 2015) %>%
  mutate(intensity_e_c = total_co2eq_kg*0.001/total_energy_use_tj) %>%
  group_by(eu_q_rank) %>%
  summarise(intensity_e_c = weighted.mean(intensity_e_c,total_energy_use_tj))

mean_co2eq_of_energy_intens_bottom_decile = round((mean_co2eq_of_energy_intens %>% filter(eu_q_rank == 1))$intensity_e_c, digits = 1)

mean_co2eq_of_energy_intens_top_decile = round((mean_co2eq_of_energy_intens %>% filter(eu_q_rank == 10))$intensity_e_c, digits = 1)

```

Increasing expenditure generally translated into larger environmental footprints across European expenditure deciles (Figure 1a). However, the energy and carbon inequality was much lower than the expenditure inequality (Figure 1b and 1c). The top decile divided by the bottom decile (the 10:10 ratio) was `r exp_10_10` for expenditure, `r energy_10_10` for energy and `r co2eq_10_10` for carbon (Figure 1a-c). Total expenditure ranged from `r exp_bottom_decile` trn€ to `r exp_top_decile` trn€ between bottom and top decile, or `r fd_pae_bottom_decile`€ to `r fd_pae_top_decile`€ per adult equivalent (ae), the energy footprint from `r energy_bottom_decile` EJ to `r energy_top_decile` EJ (or `r energy_pae_bottom_decile` GJ/ae to `r energy_pae_top_decile` GJ/ae), and the carbon footprint from `r co2eq_bottom_decile` MtCO2eq to `r co2eq_top_decile` MtCO2eq (or `r co2eq_pae_bottom_decile` tCO2eq/ae to `r co2eq_pae_top_decile` tCO2eq/ae).

The reason for this is evident from Figures 1d-f. Both the energy intensity of consumption, measured as energy use per € expenditure (d), and the carbon intensity of energy, measured as carbon footprint per energy footprint (f), decreased from bottom to top expenditure decile. The population weighted average energy intensity of consumption decreased from `r mean_energy_intens_bottom_decile` MJ/€ in the bottom decile to less than half (`r mean_energy_intens_top_decile` MJ/€) in the top decile. Likewise, the carbon intensity of energy was higher in the bottom decile (`r mean_co2eq_of_energy_intens_bottom_decile` gCO2eq/TJ) compared to the top decile (`r mean_co2eq_of_energy_intens_top_decile` gCO2eq/TJ). The carbon intensity of consumption in Figure 1f combines the effects of the intensities displayed in Figure 1d and 1e. Across all population weighted intensities per deciles, the variance in the lower four deciles is much higher (Figure 1d-f). 

```{r figure1, out.width="98%", fig.cap="Household expenditure and environmental footprints and intensities across European expenditure deciles. Total expenditures (a), energy footprint (b), and carbon footprint (c) per decile. Energy intensity of consumption as energy footprint per expenditure (d), carbon intensity of consumption as carbon footprint per expenditure (e), and carbon intensity of energy as carbon footprint per energy footprint (f)."}
knitr::include_graphics(here::here("analysis", "figures", "figure1.pdf"))
```

The different intensities of household consumption across European expenditure deciles can be attributed to a combination of two plausible causes: first, the composition of consumption baskets could systematically differ according to the level of household expenditure. Second, the energy and carbon intensity within individual consumption sectors could systematically differ according to the level of household expenditure.

```{r , fig.width=8, fig.height=2.5}

pdat_basket = get_sector_summary_by_eu_ntile_direct(eu_q_count) %>%
  ungroup() %>%
  filter(year==2015) %>%
  left_join(read_csv(here("analysis/data/derived/sectors_agg_method1_ixi.csv")), by="sector_agg_id") %>% 
  group_by(eu_q_rank) %>%
  mutate(value = total_fd_me/sum(total_fd_me)*100) %>%
  select(five_sectors, eu_q_rank, 
         value) %>%
  mutate(indicator = "Expenditure share (%)")

pdat_int_co2eq = get_sector_summary_by_eu_ntile_direct(eu_q_count) %>%
  ungroup() %>%
  filter(year==2015) %>%
  left_join(read_csv(here("analysis/data/derived/sectors_agg_method1_ixi.csv")), by="sector_agg_id") %>% 
  mutate(value = (total_co2eq_kg)/(total_fd_me*1000000)) %>%
  select(five_sectors, eu_q_rank, 
         value) %>%
  mutate(indicator = "Carbon intensity (kgCO2eq/€)")

pdat_int_energy = get_sector_summary_by_eu_ntile_direct(eu_q_count) %>%
  ungroup() %>%
  filter(year==2015) %>%
  left_join(read_csv(here("analysis/data/derived/sectors_agg_method1_ixi.csv")), by="sector_agg_id") %>% 
  mutate(value = (total_energy_use_tj)/(total_fd_me)) %>%
  select(five_sectors, eu_q_rank,
         value) %>%
  mutate(indicator = "Energy intensity (MJ/€)")

library(viridis)
pdat = pdat_int_co2eq %>%
  bind_rows(pdat_int_energy) %>%
  bind_rows(pdat_basket) 

pdat$facet = factor(pdat$indicator, levels = c("Expenditure share (%)", "Energy intensity (MJ/€)", "Carbon intensity (kgCO2eq/€)"))

#a = ggplot(pdat, aes(x=value, y=five_sectors, color=eu_q_rank)) + 
#  #geom_point(shape=18, alpha=0.75) +
#  geom_text(label="I", alpha=0.75, fontface="bold") +
#  scale_color_viridis(option = "A", end = 0.8,
#                      direction = -1, 
#                      breaks = c(1, 5, 10), 
#                      name="European\ndecile") +
#  facet_wrap(~facet, scales="free_x") +
#  theme_minimal() +
#  labs(x="",y="") +
#  theme(text=element_text(size=13, family="Liberation Sans Narrow"),
#        panel.spacing = unit(2.5, "lines")) 

a = ggplot(pdat, aes(x=value, y=five_sectors, color=eu_q_rank)) +
  geom_text(label="I", alpha=0.75, fontface="bold") +
  scale_colour_viridis_b(option = "A", end = 0.8,
                         direction = -1,
                         breaks = c(2,3,4,5,6,7,8,9),
                         #guide_legend(),
                         guide = guide_colorsteps(show.limits=TRUE),
                         name="European\ndecile") +
  facet_wrap(~facet, scales="free_x",
             strip.position="bottom") +
  theme_minimal() +
  labs(x="",y="") +
  theme(text=element_text(size=13, family="Liberation Sans Narrow"),
        panel.spacing = unit(2.5, "lines"),
        strip.placement="outside") 

ggsave(here("analysis", "figures", "figure2.pdf"), device=cairo_pdf)
```

```{r figure2, out.width="100%", fig.cap="Sectoral expenditure shares and carbon intensities of European expenditure deciles. Share of expenditure per final demand sector of total spending per decile in percent (a) and carbon intensity per final demand sector in kgCO2eq/€."}
knitr::include_graphics(here::here("analysis", "figures", "figure2.pdf"))
```

```{r values-in-text2}

exp_share_housing_bottom_decile = round((pdat_basket %>% filter(eu_q_rank == 1, five_sectors == "housing"))$value, digits = 1)

exp_share_housing_top_decile = round((pdat_basket %>% filter(eu_q_rank == 10, five_sectors == "housing"))$value, digits = 1)

int_co2eq_housing_bottom_decile = round((pdat_int_co2eq %>% filter(eu_q_rank == 1, five_sectors == "housing"))$value, digits = 1)

int_co2eq_housing_top_decile = round((pdat_int_co2eq %>% filter(eu_q_rank == 10, five_sectors == "housing"))$value, digits = 1)

exp_share_services_bottom_decile = round((pdat_basket %>% filter(eu_q_rank == 1, five_sectors == "services"))$value, digits = 1)

exp_share_services_top_decile = round((pdat_basket %>% filter(eu_q_rank == 10, five_sectors == "services"))$value, digits = 1)

```

Our results show that both of these factors play a role (Figure 2). Lower-income households, on average, spend larger shares of their expenditure in the housing sector. The bottom decile spent an average of `r exp_share_housing_bottom_decile`% of their household expenditure on housing, while the top decile spent `r exp_share_housing_top_decile`%. The housing sector stands out with a carbon intensity of consumption more than 3 times higher in the bottom decile (`r int_co2eq_housing_bottom_decile` kgCO2eq/€) than in the top decile (`r int_co2eq_housing_top_decile` kgCO2eq/€). At the same time, housing is by far the most carbon intensive sector and has the highest variance in carbon intensity among expenditure deciles. Overall, with increasing decile, the shares of mobility and services expenditures increase and the shares of housing, food and goods decrease. Households in the top decile spend about `r exp_share_services_top_decile`% in the service sector, which has the lowest carbon intensity of all consumption sectors, compared to `r exp_share_services_bottom_decile`% in the bottom decile. 

The tendency for energy and carbon intensity to decrease with increasing affluence has been reported for the global level [@hubacek_global_2017] between countries and also within Europe [@sommer_carbon_2017]. Our results show that the four lowest European expenditure deciles make up 80% to 100% of the population in Poland, Romania, Bulgaria and the Czech Republic, while less than 20% of the population in the higher-income European countries (Scandinavia, Germany, France, Austria, the Netherlands, Belgium, the UK, and Ireland) are in the lowest European expenditure deciles. Note here that this does not imply that there are no high-income households in Eastern Europe. Our analysis is based on average expenditure data from national income quintiles. This aggregation cuts off the lower and higher tails of the respective national expenditure distributions (see SI - Supplementary Note and Map). 

The high intensities in the bottom four European expenditure deciles can be attributed in large part to more inefficient and dirtier domestic energy supply and demand technologies for heating and electricity generation in Poland, Bulgaria, the Czech Republic, and Romania. Poland alone was responsible for about 40% of total coal combustion for heat production in Europe in 2015 [@eurostat_eurostat_nodate-2], and had a higher average intensity of carbon per MJ of heat delivered than both Europe and the world [@werner_international_2017]. We did not account here for subsidies which could also have attributed to high energy and carbon intensities (see SI limitations, pp xx). 

## Inequality across final consumption sectors

```{r , fig.width=8, fig.height=5}

pdat = get_sector_summary_by_eu_ntile_direct(eu_q_count) %>%
  ungroup() %>%
  filter(year==2015) %>%
  left_join(read_csv(here("analysis/data/derived/sectors_agg_method1_ixi.csv")), by="sector_agg_id") %>% 
  select(five_sectors, eu_q_rank, contains("pae_co2eq"), contains("pae_energy")) %>%
  pivot_longer(cols=-c(five_sectors, eu_q_rank), names_to="indicator", values_to="value") %>%
  mutate(indicator_type = if_else(str_detect(indicator, "co2eq"), 
                                  "tCO2eq per adult eq", 
                                  "Energy GJ per adult eq")) %>%
  filter(indicator != "pae_co2eq_t", 
         indicator != "pae_energy_use_gj")
  
  
pal <- wes_palette("Zissou1", 4, type = "discrete")
pal = pal[c(2,1,3,4)]
  
p1 = ggplot(pdat %>% filter(indicator_type == "tCO2eq per adult eq"), 
       aes(x=(eu_q_rank), y=value, fill=indicator)) +
  geom_col(position = position_stack()) +
  facet_wrap(~five_sectors, ncol = 5) +
  scale_fill_manual(values = pal, name="Location", labels=c("Direct", "Domestic", "Europe", "non-Europe")) +
  scale_x_continuous(breaks = c(1,2,3,4,5,6,7,8,9,10), 
                     labels = c("D01","D02","D03","D04","D05","D06","D07","D08","D09","D10")) +
  #theme_ipsum() +
  theme_minimal() +
  theme(text=element_text(family="Liberation Sans Narrow")) +
  labs(x="", y="Carbon footprint tCO2eq per adult eq") +
  theme(legend.position = "bottom", 
          axis.text.x = element_text(angle = 90, size = 12),
          axis.text.y = element_text(size = 12),
          strip.text = element_text(size = 12),
          axis.title.y = element_text(size = 10))


check = pdat %>% filter(indicator_type == "Energy GJ per adult eq") %>%
  group_by(indicator) %>%
  summarise(value = sum(value))

p2 = ggplot(pdat %>% filter(indicator_type == "Energy GJ per adult eq"), 
       aes(x=(eu_q_rank), y=value, fill=indicator)) +
  geom_col(position = position_stack()) +
  facet_wrap(~five_sectors, ncol = 5) +
  scale_fill_manual(values = pal, name="Location", labels=c("Direct", "Domestic", "Europe", "non-Europe")) +
  scale_x_continuous(breaks = c(1,2,3,4,5,6,7,8,9,10), 
                     labels = c("D01","D02","D03","D04","D05","D06","D07","D08","D09","D10")) +
  #theme_ipsum() +
  theme_minimal() +
  theme(text=element_text(family="Liberation Sans Narrow")) +
  labs(x="", y="Energy footprint GJ per adult eq") +
  theme(legend.position = "bottom", 
          axis.text.x = element_text(angle = 90, size = 12),
          axis.text.y = element_text(size = 12),
          strip.text = element_text(size = 12),
          axis.title.y = element_text(size = 10))

a = p1 / p2 + plot_layout(guides = "collect") & theme(legend.position = 'bottom', 
                                                  plot.margin =  margin(1, 1, 1, 1, "mm"),
                                                  panel.spacing = unit(c(1, 1, 1, 1), "mm"))


hm = pdat %>% filter(indicator_type == "Energy GJ per adult eq") %>%
  group_by(eu_q_rank) %>%
  summarise(value = sum(value))

ggsave(here("analysis", "figures", "figure3.pdf"))
```

```{r values-in-text3}

co2eq_per_sector = pdat %>%
  filter(indicator_type == "tCO2eq per adult eq") %>%
  group_by(five_sectors, eu_q_rank) %>%
  summarise(value = sum(value))

mobility_co2eq_10_10 = round((co2eq_per_sector %>% filter(eu_q_rank == 10, five_sectors == "mobility"))$value/
                                (co2eq_per_sector %>% filter(eu_q_rank == 1, five_sectors == "mobility"))$value,digits = 0)

goods_co2eq_10_10 = round((co2eq_per_sector %>% filter(eu_q_rank == 10, five_sectors == "goods"))$value/
                                (co2eq_per_sector %>% filter(eu_q_rank == 1, five_sectors == "goods"))$value,digits = 0)

services_co2eq_10_10 = round((co2eq_per_sector %>% filter(eu_q_rank == 10, five_sectors == "services"))$value/
                                (co2eq_per_sector %>% filter(eu_q_rank == 1, five_sectors == "services"))$value,digits = 1)

food_co2eq_10_10 = round((co2eq_per_sector %>% filter(eu_q_rank == 10, five_sectors == "food"))$value/
                                (co2eq_per_sector %>% filter(eu_q_rank == 1, five_sectors == "food"))$value,digits = 1)

energy_per_sector = pdat %>%
  filter(indicator_type == "Energy GJ per adult eq") %>%
  group_by(five_sectors, eu_q_rank) %>%
  summarise(value = sum(value))

mobility_energy_10_10 = round((energy_per_sector %>% filter(eu_q_rank == 10, five_sectors == "mobility"))$value/
                                (energy_per_sector %>% filter(eu_q_rank == 1, five_sectors == "mobility"))$value,digits = 0)

goods_energy_10_10 = round((energy_per_sector %>% filter(eu_q_rank == 10, five_sectors == "goods"))$value/
                                (energy_per_sector %>% filter(eu_q_rank == 1, five_sectors == "goods"))$value,digits = 0)

services_energy_10_10 = round((energy_per_sector %>% filter(eu_q_rank == 10, five_sectors == "services"))$value/
                                (energy_per_sector %>% filter(eu_q_rank == 1, five_sectors == "services"))$value,digits = 1)

food_energy_10_10 = round((energy_per_sector %>% filter(eu_q_rank == 10, five_sectors == "food"))$value/
                                (energy_per_sector %>% filter(eu_q_rank == 1, five_sectors == "food"))$value,digits = 1)

```

The final consumption sectors (housing, mobility, food, goods, and services) contributed very differently to the total environmental footprint of European households in 2015 (Figure 3). On average, housing and mobility are the two largest sectors, accounting for nearly two thirds of both the energy and carbon footprints. However, there are big differences between the sectors when looking at the respective contributions of each expenditure decile. For housing there is very little difference between deciles in both the energy and the carbon footprint. The bottom four deciles even have higher carbon footprints from housing than most top deciles, which can be explained by the extreme differences in intensity shown in Figure 2. Mobility was the most unequal sector, with footprints in the top decile 10 times higher than the bottom decile, corroborating findings in [@ivanova_quantifying_2020] and [@oswald_large_2020]. Goods was the second most unequal final consumption sector (10:10 ratios around `r goods_energy_10_10` for both footprints), followed by services (10:10 ratios of `r services_energy_10_10` for energy and `r services_co2eq_10_10` for carbon) and then food (10:10 ratios of `r food_energy_10_10` for both footprints).

```{r figure3, out.width="100%", fig.cap="Energy and carbon footprints by final consumption sector and European expenditure decile in 2015, further broken down by emission source location."}
knitr::include_graphics(here::here("analysis", "figures", "figure3.pdf"))
```

```{r values-in-text4}

co2eq_per_source = pdat %>%
  filter(indicator_type == "tCO2eq per adult eq") %>%
  group_by(five_sectors, indicator) %>%
  summarise(value = sum(value)) %>%
  mutate(share = value/sum(value))

housing_co2eq_direct = round(((co2eq_per_source %>% filter(five_sectors == "housing", 
                                                          indicator == "pae_co2eq_direct_t"))$share)*100, digits = 0)

energy_per_source = pdat %>%
  filter(indicator_type == "Energy GJ per adult eq") %>%
  group_by(five_sectors, indicator) %>%
  summarise(value = sum(value)) %>%
  mutate(share = value/sum(value))

housing_energy_direct = round(((energy_per_source %>% filter(five_sectors == "housing", 
                                                          indicator == "pae_energy_use_direct_gj"))$share)*100, digits = 0)

```

The geographical source of the household energy and carbon footprints also varies with consumption categories (Figure 3). The housing footprint was almost entirely domestic, with `r housing_co2eq_direct`% and `r housing_energy_direct`% respectively coming from direct household energy use and carbon emissions from heating and cooling, and the rest embedded primarily along the domestic supply chain. The mobility footprint, on the other hand, was around one fourth non-European. The majority of the mobility footprint, above 60%, came from vehicle fuel, either directly from household, or indirectly, i.e. embedded along household's supply chains. The goods footprint was mostly non-European, while services and food were both around one third non-European. These results suggest that proposed future carbon border-adjustment mechanisms [@european_commission_communication_2019] will especially impact the goods and mobility footprints of the higher deciles, and to a lesser extent the food and services footprints.

# Counterfactual: a 1.5°C compatible Europe

Global 1.5°C compatible decarbonisation scenarios achieve a similar climate outcome with different assumptions about the transformation of energy supply and demand, from renewable capacity, deployment of carbon-capture-and-storage (CCS), to socio-technological demand transformation. Table 1 shows some final energy use results for the year 2050 from six different decarbonisation scenarios, already adjusted from total GJ/capita to household GJ/adult equivalent. The original total GJ/capita scenario results are from different world regions (OECD, West EU, Global North, and Global), depending on the regional disaggregation of the publicly available scenario results, and so should not be interpreted as perfectly comparable with each other. For the purposes of our study, however, we are simply interested in the range of scenario results within which to situate our household footprint results, presented below in the ‘Inequality in a 1.5°C compatible Europe’ section and Figure 5.  

```{r}

df_scenario_info = read_excel(here("analysis/data/raw/scenarios.xlsx"), sheet="overview") %>%
  select(scenario, type,fe_gj_aeu = final_energy_gj_per_aeu_2050,
         ccs_required = primary_energy_fossil_w_ccs2050_ej,
         description) %>%
  arrange(desc(fe_gj_aeu)) %>%
  mutate(fe_gj_aeu = round(fe_gj_aeu),
         ccs_required = round(ccs_required)) %>%
  select(scenario,type,fe_gj_aeu)

flextable(df_scenario_info) %>%
  autofit() %>%
  set_header_labels(scenario = "scenario", 
                    type = "type",
                    fe_gj_aeu = "final energy in 2050: household GJ/adult equivalent") %>%
  set_caption("Decarbonisation scenarios.") %>%
  align(j=3, align = "center") %>%
  align(align = "left", part = "header") %>%
  width(width = 2.1) 
```

The various global supply side scenarios (SSP1-1.9, SSP2-1.9, GEA efficiency, IEA ETP B2DS)[@riahi_shared_2017 @gea_gea_nodate @grubler_low_2018] envisage household European energy use falling from the 2015 level of 92 EJ to around 21-31 EJ by 2050, equivalent to a per household reduction from a current average of 250 GJ to 64-94 GJ per adult equivalent. The differences in energy use in 2050 in the scenarios reflect different model assumptions about the rate of expansion of renewable energy and CCS capacity. These scenarios rely on substantial amounts of CCS, which is still a fairly speculative technology, and we therefore interpret them as ranges for the upper limits of 1.5°C-compatible energy supply [@riahi_shared_2017 @gea_gea_nodate].

It is even more difficult to determine a lower limit for the minimum amount of energy needed for a decent life. Such a lower limit depends strongly on the prevalent socio-cultural idea of what constitutes a decent life, and, perhaps even more strongly, on the physical infrastructure available to deliver this life. The two global demand side scenarios LED [@grubler_low_2018] and DLE [@millward-hopkins_providing_2020] that attempt to define such a limit conclude that, in principle, a very low energy footprint, between 16-53 GJ per household adult equivalent, could be sufficient. However, these scenarios rely on socio-technological transformations on a scale that, especially at the lower end, far exceed the current political discourse on the subject. These scenarios are 1.5°C compatible without resorting to any CCS but they all implicitly (LED) [@grubler_low_2018] or explicitly (DLE) [@millward-hopkins_providing_2020] assume near complete equality of consumption across the population. To put these low energy demand numbers in perspective, the average energy footprint in our sample is 250 per adult equivalent in 2015, about a factor 5 above the high estimate. Households in the bottom European expenditure decile, which almost entirely fell within the Eurostat definition of severe material deprivation [@eurostat_living_nodate], still had an energy footprint of 130 GJ per adult equivalent in 2015 (roughly 80 GJ/capita), factor 2.5 above the high estimate. 

Based on these two constraints, the upper limit on the supply side and the lower limit on the demand side, it is possible to make a generalized estimate of how much inequality in the distribution of energy consumption is numerically possible, if at the same time global warming is to be kept below 1.5°C above pre-industrial levels and a decent life for all is to be made possible. Before we can make this evaluation, we must take into account the existing large differences in the technological efficiency of energy provision (Figure 2). These differences will be adjusted in the next step. 

## Current empirical best technology per sector

```{r}
pal <- wes_palette("Zissou1", 4, type = "discrete")
pal = pal[c(2,1,3,4)]

pdat_int_energy = get_sector_summary_by_eu_ntile_direct(eu_q_count) %>%
  ungroup() %>%
  filter(year==2015) %>%
  left_join(read_csv(here("analysis/data/derived/sectors_agg_method1_ixi.csv")), 
            by="sector_agg_id") %>% 
  mutate(intensity_energy = (total_energy_use_tj)/(total_fd_me)) %>%
  select(five_sectors, eu_q_rank, 
         intensity_energy) %>%
  filter(eu_q_rank == 10) %>%
  select(-eu_q_rank)

pdat_final_demand = get_sector_summary_by_eu_ntile_direct(eu_q_count) %>%
  ungroup() %>%
  filter(year==2015) %>%
  left_join(read_csv(here("analysis/data/derived/sectors_agg_method1_ixi.csv")), 
            by="sector_agg_id") %>%
  left_join(pdat_int_energy, by="five_sectors") %>%
  mutate(total_energy_use_tj_new = (total_fd_me)*intensity_energy) %>%
  mutate(total_energy_use_tj_diff = total_energy_use_tj-total_energy_use_tj_new) %>%
  select(five_sectors, eu_q_rank, 
         total_energy_use_tj_diff, total_energy_use_tj_new) %>%
  group_by(eu_q_rank) %>%
  summarise(total_energy_use_tj_new = sum(total_energy_use_tj_new)*0.000001,
            total_energy_use_tj_diff = sum(total_energy_use_tj_diff)*0.000001) %>%
  pivot_longer(-c(eu_q_rank), names_to = "indicator", values_to = "value")


p1 = ggplot(pdat_final_demand, aes(x=eu_q_rank, y=value, fill=indicator, alpha=indicator)) +
  geom_col(position = position_stack()) +
  scale_fill_manual(values=c(pal[1], pal[2]), 
                    labels=c("2015", "Best technology"), name="Energy\nfootprint") +
  scale_alpha_manual(values=c(0.3,1),labels=c("2015", "Best technology"), name="Energy\nfootprint")+
  scale_x_continuous(breaks = c(1,2,3,4,5,6,7,8,9,10), 
                     labels = c("D01","D02","D03","D04","D05","D06","D07","D08","D09","D10")) +
  labs(y="Energy footprint (EJ)", x="") +
  theme_minimal() +
  coord_flip() +
  theme(legend.position="bottom")

```

```{r}

pdat_energy_country = get_sector_summary_by_country_quintile_direct() %>%
  ungroup() %>%
  filter(year==2015) %>%
  left_join(read_csv(here("analysis/data/derived/sectors_agg_method1_ixi.csv")), 
            by="sector_agg_id") %>%
  left_join(pdat_int_energy, by="five_sectors") %>%
  mutate(total_energy_use_tj_new = (total_fd_me)*intensity_energy) %>%
  #mutate(total_energy_use_tj_save = total_energy_use_tj_new/total_energy_use_tj*100) %>%
  select(iso2, 
         total_energy_use_tj_new, total_energy_use_tj) %>%
  group_by(iso2) %>%
  summarise(total_energy_use_tj_save = 100 - sum(total_energy_use_tj_new)/sum(total_energy_use_tj)*100)

```

```{r include=FALSE}

library(rworldmap)
library(ggthemes)



df_country_mean = pdat_energy_country %>%
  mutate(iso2 = if_else(iso2 == "EL", "GR", iso2)) %>%
  mutate(iso2 = if_else(iso2 == "UK", "GB", iso2)) %>%
  left_join(ISO_3166_1 %>% select(iso2 = Alpha_2, iso3 = Alpha_3), by="iso2")


quantile_rank_map = joinCountryData2Map(df_country_mean, 
                                        joinCode = "ISO_A3", 
                                        nameJoinColumn = "iso3")
quantile_rank_map_poly = fortify(quantile_rank_map) #extract polygons 

quantile_rank_map_poly = merge(quantile_rank_map_poly, quantile_rank_map@data, by.x="id", by.y="ADMIN", all.x=T)
quantile_rank_map_poly = quantile_rank_map_poly %>% arrange(id, order)

```

```{r dle-country, fig.width=8}
map = ggplot() + 
  geom_polygon(data = quantile_rank_map_poly, aes(long, lat, group = group), color="black") +
  geom_polygon(data = quantile_rank_map_poly, aes(long, lat, group = group, 
                                        fill=total_energy_use_tj_save)) + 
  #scale_fill_gradient(name="Mean decile \nrank", low="skyblue", high = "skyblue4", na.value = "#EEEEEE") +
  scale_fill_viridis(direction = -1, discrete = F, na.value = "white", name="Energy\nsavings (%)") +
  theme_map() +
  labs(x="",y="")+
  coord_map("bonne", lat0 = 50,xlim = c(-9, 40), ylim = c(38, 68), clip="on") +
  theme(legend.position = c(0.85, 0.6))
  
```

```{r fig.height=7, fig.width=11}

a = p1 + map + plot_annotation(tag_levels = 'a') +
  theme(plot.margin = unit(c(0.25,0.25,0.25,0.25), "cm"),
        legend.position = 'bottom',
        axis.title.y = element_text(size=13, hjust = 0.5),
        legend.text = element_text(size=12),
        legend.title = element_text(size=13))

ggsave(here("analysis", "figures", "figure4.pdf"))
```

```{r figure4, out.width="100%", fig.cap="Energy savings through sectoral best current technology by expenditure deciles a) and country b)."}
knitr::include_graphics(here::here("analysis", "figures", "figure4.pdf"))
```

Our results show that in 2015, higher-income people in higher-income countries had access to the most energy-efficient energy services across all final consumption sectors (Figure 2). Since we are interested in the numerically possible inequality in the distribution of actual consumption of goods and services in the next section, these efficiency differences must first be adjusted. In practice, this corresponds, for example, to the need for large-scale investments in the technical efficiency of heat, electricity and hot water supply, especially in Eastern Europe [@bianco_understanding_2019]. Improving technical efficiency is already a major part of the European Union (EU) platform, and new transition funds for lower-income countries, whether public or private under a Green Deal framework, need to be appropriately targeted, and at an appropriately large scale, to reduce the high intensities of consumption in the lower deciles [@european_commission_communication_2019 @european_commission_european_2020]. Figure 4 shows the energy footprint savings per decile (Fig. 4a) that would have occurred in 2015 if all deciles had the same efficiency per final consumption sector as the top decile. Around 17 EJ would have been saved in total, and the energy footprint of the bottom decile would have been nearly half its 2015 value. Fig. 4b shows saved energy per country, with Eastern European countries especially saving large proportions of their 2015 footprint, over 60% for Bulgaria and Estonia for example.

## Inequality in a 1.5°C compatible Europe

```{r }

pdat_int_energy = get_sector_summary_by_eu_ntile_direct(eu_q_count) %>%
  ungroup() %>%
  filter(year==2015) %>%
  left_join(read_csv(here("analysis/data/derived/sectors_agg_method1_ixi.csv")), 
            by="sector_agg_id") %>% 
  mutate(intensity_energy = (total_energy_use_tj)/(total_fd_me)) %>%
  select(five_sectors, eu_q_rank, 
         intensity_energy) %>%
  filter(eu_q_rank == 10) %>%
  select(-eu_q_rank)

pdat_final_demand = get_sector_summary_by_eu_ntile_direct(eu_q_count) %>%
  ungroup() %>%
  filter(year==2015) %>%
  left_join(read_csv(here("analysis/data/derived/sectors_agg_method1_ixi.csv")), 
            by="sector_agg_id") %>%
  left_join(pdat_int_energy, by="five_sectors") %>%
  mutate(total_energy_use_tj_new = (total_fd_me)*intensity_energy) %>%
  mutate(total_energy_use_tj_diff = total_energy_use_tj-total_energy_use_tj_new) %>%
  select(eu_q_rank,total_energy_use_tj_new) %>%
  group_by(eu_q_rank) %>%
  summarise(total_energy_use_tj_new = sum(total_energy_use_tj_new)) %>%
  mutate(pae_energy_use_tj = total_energy_use_tj_new/33417583,
         pae_energy_use_gj = pae_energy_use_tj*1000)

df_energy_deciles = pdat_final_demand %>%
  select(eu_q_rank, pae_energy_use_gj)

ineq_curr = df_energy_deciles$pae_energy_use_gj[10]/df_energy_deciles$pae_energy_use_gj[1]

```

```{r }

library(readxl)

df_scenario_info = read_excel(here("analysis/data/raw/scenarios.xlsx"), sheet="overview") %>%
  select(scenario, fe_gj_aeu = final_energy_gj_per_aeu_2050,
         ccs_required = primary_energy_fossil_w_ccs2050_ej,
         description) %>%
  arrange(fe_gj_aeu) %>%
  mutate(fe_gj_aeu = round(fe_gj_aeu),
         ccs_required = round(ccs_required)) 

mea = c(min(df_scenario_info$fe_gj_aeu),max(df_scenario_info$fe_gj_aeu))
mer = c(min(df_scenario_info$fe_gj_aeu),53)

c_mean_mea = round((27+mea[2])/2)
c_mean_mer = round((mer[1]+mer[2])/2)

```

```{r , eval = FALSE}

# run once to save file. If changing scenarios included or input data, need to re-run and save 'scenarios_fine.rds'

# vectorized function that returns scaled quantiles given 
#quantile index column and column with quantile averages
#qidx: quantile index
#qavg_to_scale: column to scale
#first_target: target value of first quantile
#mean_target: target mean of scaled quantiles
scaled_quantiles <- function(.data, 
                          qidx, 
                          qavg_to_scale, 
                          first_target, 
                          mean_target) {
  
  # cumbersomely extract current quantile mean
  mean_current = .data %>%
    ungroup() %>%
    summarise(mean_cur = first(mean({{qavg_to_scale}}))) %>%
    pull(mean_cur)
    
  # cumbersomely extract current first wuantile value
  first_current = .data %>%
    ungroup() %>%
    arrange({{qidx}}) %>%
    summarise(first_cur = first({{qavg_to_scale}})) %>%
    pull(first_cur)
  
  df_tmp = .data %>%
    mutate(tmp = {{qavg_to_scale}}*mean_target/mean_current)
  
  first_tmp = df_tmp$tmp[1]
  
  df_tmp = df_tmp %>%
    mutate(scaled = mean_target-(mean_target-tmp) * (mean_target-first_target)/(mean_target-first_tmp)) %>%
    select({{qidx}}, scaled) %>%
    mutate(v_mean = mean_target, v_first = first_target)
}

## run once to save file
df_all = NULL
for (min_energy in seq(from=mer[1], to=mer[2], by=0.25)) {
  for (mean_energy in seq(from=mea[1], to=mea[2], by=0.25)) {
    if (min_energy <= mean_energy) {
      df_all = df_all %>%
      bind_rows(df_energy_deciles %>%
                  scaled_quantiles(eu_q_rank, pae_energy_use_gj, min_energy, mean_energy))
    }
  }
}

saveRDS(df_all, here("analysis/data/derived/scenarios_fine.rds"))
```

```{r , fig.width=7, fig.height=5.5}

round_by = 10

df_all = readRDS(here("analysis/data/derived/scenarios_fine.rds")) %>%
  filter(eu_q_rank %in% c(1,10)) %>%
  group_by(v_mean, v_first) %>%
  summarise(ratio = last(scaled)/first(scaled)) %>%
  mutate(bin_ratio = if_else((ratio*100)%%round_by > round_by*0.5, 
                             ratio*100+(round_by-(ratio*100)%%round_by),
                             ratio*100-(ratio*100)%%round_by)) %>%
  group_by(bin_ratio, v_first) %>%
  summarise(v_mean = mean(v_mean)) %>%
  mutate(bin_ratio = bin_ratio*0.01)

df_scenario = df_all %>%
  filter(v_mean %in% df_scenario_info$fe_gj_aeu) %>%
  left_join(df_scenario_info, by=c("v_mean"="fe_gj_aeu"))

library(wesanderson)

a = df_all %>%
  ggplot(aes(x=v_first, y=bin_ratio, fill=v_mean)) +
  geom_tile() +
  geom_hline(yintercept = ineq_curr, alpha=0.8, color="grey", linetype=2) +
  geom_line(data=df_scenario, aes(color=scenario, group=scenario)) +
  annotate(geom="text", x=max(df_all$v_first)-5.7,y=ineq_curr+0.6,label = "Current (2015) 10:10 ratio") +
  #scale_fill_viridis("Mean energy\navailable") +
  scale_fill_gradient("Mean energy\navailable (GJ/ae)",
                      low=wes_palette("Chevalier1")[3], 
                      high = wes_palette("Rushmore1")[4]) +
  scale_color_manual(values=wes_palette("Darjeeling1")) +
  theme_minimal() +
  labs(x="Minimum energy requirement (GJ/ae)", y="Maximum energy inequality (10:10 ratio)", color = "Scenario")+
  theme(text=element_text(family="Liberation Sans Narrow"),
        axis.text.x = element_text(size = 13),
        axis.text.y = element_text(size = 13)) +
  scale_y_continuous(breaks = c(2.5,5,7.5,10,12.5)) # +
  #theme_ipsum()

ggsave(here("analysis", "figures", "figure5.pdf"))
```

Based on this counterfactual distribution of the energy footprint using homogeneous supply technologies, we can now scale down energy use across European expenditure deciles to meet supply constraints and, where necessary, "squeeze" the distribution to not undershoot minimum energy use requirements in any decile (Figure 5). 

Both the DLE and LED scenarios satisfy energy demand for decent living and are compatible with the 1.5 degree target without resorting to CCS technologies [@millward-hopkins_providing_2020 @grubler_low_2018]. The DLE scenario explicitly envisions absolute global inequality (10:10 ratio of 1) in consumption, except for small differences in required energy consumption based on climatic and demographic factors, as well as differences in population density [@millward-hopkins_providing_2020]. The LED scenario does not explicitly discuss distributional aspects beyond giving different final energy values for the Global North (53GJ/ae) and the Global South (27GJ/ae) [@grubler_low_2018]. However, due to the bottom-up construction of this demand scenario, these values can be interpreted as estimates for the minimum required energy use. 

The descriptions of the energy supply scenarios do not include specific details about how the energy footprints are distributed within the population. The energy savings here are achieved primarily through efficiency improvements, and perhaps also generally assumed demand reductions. 

```{r figure5, out.width="70%", fig.align="center", fig.cap="Mean energy available for Europe in decarbonisation scenarios, positioned in option space between a range of minimum energy requirements and range of associated maximum inequality. All expenditure deciles have 'best technology' already."}
knitr::include_graphics(here::here("analysis", "figures", "figure5.pdf"))
```

The colored curves in Figure 5 represent constant average household energy footprints according to the different scenarios. The slopes of the curves connect different assumptions about minimal energy for a decent living (on the x-axis) to the corresponding energy inequality that is consistent with the average energy availability. It is clear from Figure 5 that at current inequality levels, only the two scenarios with heavy CCS deployment [add scenario acronym] and GEA efficiency are possible and only if we assume in addition extremely low minimum energy use requirements (below 27 GJ/ae). This 27 GJ/ae is roughly the value the low-energy demand (LED) scenario gives for the Global South in 2050. If we use the value given for the Global North at 53 GJ/ae as minimum energy requirements, which still requires strong demand side measures, then inequality would need to be zero in the LED scenario and cut down by factors xx to yyy in all other scenarios. 

# Conclusions

Estimates of energy and carbon footprint inequality are increasingly being used to assign responsibility for climate change. At a global, regional, and within-country level, energy use and carbon emissions are often highly unequal [@piketty_carbon_2015 @kartha_carbon_2020 @gore_extreme_2015 @hubacek_global_2017 @ivanova_unequal_2020 @gore_t._confronting_2020 @wiedenhofer_unequal_2017 @golley_income_2012 @steenolsen_carbon_2016 @weber_quantifying_2008 @hardadi_implications_2020 @oswald_large_2020]. The proposed solution is often a call to reduce the carbon or energy inequality by reducing over-consumption, especially by the richest at the top of the economic distribution, which would then also reduce the energy and carbon footprints, everything else held equal. Complicating this picture, however, is the fact that energy and carbon intensities of consumption usually differ between economic groups. This is due to different consumption baskets and different access to technology. That lower-income groups tend to have higher energy and carbon intensities is an important finding from the environmental Kuznet's curve literature [@berthe_mechanisms_2015 @scruggs_political_1998]. This finding has not yet been well integrated with the current carbon and energy footprint inequality literature, that focuses more on assigning responsibility based on aggregate energy and carbon footprint inequality. 

In this study, we have found that, for Europe as a whole, lower-economic groups have higher energy and carbon intensities of consumption (although this is not necessarily true within each European country) [@sommer_carbon_2017 @kerkhof_determinants_2009]. These higher intensities come almost entirely from domestic electricity production and heating/cooling for housing, in a handful of Central and Eastern European countries. This is already an important focus of European climate policy, but reducing these intensities should be a major priority for investment fund allocation going forward, especially within a framework such as the EU's European Green Deal [@bianco_understanding_2019]. Efforts to break consumer lock-in to these high intensities must be occurring alongside policies that seek to reduce aggregate consumption and intensities higher up in the economic distribution [@royston_invisible_2018]. Bringing intensities of consumption for all economic groups in line with those of higher-economic groups in Europe with access to the cleanest and most efficient available technologies, would substantially reduce the European household energy and carbon footprints, everything else held equal. The unequal intensity structure hinders clear conclusions on footprint inequality. We have shown that in an important sector such as housing, lower-economic groups have almost the same level of footprint as higher-economic groups despite a fraction of the expenditure, because of their higher intensities. This can then be misleading in terms of assigning responsibility for climate change. Bringing energy and carbon intensities of all economic groups in line with the top group, and thus removing the inequality in intensity structure, would reduce the footprint, all else held equal, but *increase* energy and carbon inequality. The reduction of energy and carbon inequality is not a meaningful goal by itself. 

Current consumption inequality, however, is a barrier to achieving both scenario targets *and* providing minimum energy use (and minimum carbon in the short-term) for decent living to every European. At a global level, there is concern that achieving extreme poverty eradication may put global climate targets at risk [@hubacek_poverty_2017 @woodward_incrementum_2015 @alfredsson_why_2018]. In the European context, although less unequal than the globe as a whole, if/as lower-consumption groups increase their income and consumption, energy use and carbon emissions will increase if more efficient and cleaner technology is not adopted at a fast enough rate [@gough_recomposing_2017]. Achieving an average per capita/adult equivalent energy and carbon footprint in Europe, in scenarios that reach the Paris agreement goals, means either doing so at current consumption inequality levels and keeping lower-economic groups near or below minimum energy use levels for decent living, or reducing consumption inequality. 

<!-- The following line inserts a page break  -->

\newpage

# Associated Content

## Supplementary Information

Supplementary materials and methods, including extended discussion on limitations, and supplementary results.

# Author Information

## Corresponding Author

*email: 

## Author Contributions

I.S.J., P-P.P., and H.W. designed research; I.S.J. and P-P.P. performed research; I.S.J., P-P.P., J.T., and H.W.
interpreted results; and I.S.J., P-P.P., J.T., and H.W. wrote the paper.

## Notes

The authors declare no competing financial interest.

# Acknowledgements

Jakob Napiontek

\newpage

# References

<!-- The following line ensures the references appear here for the MS Word or HTML output files, rather than right at the end of the document (this will not work for PDF files):  -->

::: {#refs}
:::

\newpage

### Colophon

This report was generated on `r Sys.time()` using the following computational environment and dependencies:

```{r colophon, cache = FALSE}
# which R packages and versions?
if ("devtools" %in% installed.packages()) devtools::session_info()
```

The current Git commit details are:

```{r}
# what commit is this file at? 
if ("git2r" %in% installed.packages() & git2r::in_repository(path = ".")) git2r::repository(here::here())  
```
