diff --git a/analysis/preprocessing/full_code.Rmd b/analysis/preprocessing/full_code.Rmd
index edb5dfae3f06145abf26607a0972fac7b3bbb0ad..ed175d2f12d96797a5a5020023f3368119ce4384 100644
--- a/analysis/preprocessing/full_code.Rmd
+++ b/analysis/preprocessing/full_code.Rmd
@@ -1124,7 +1124,7 @@ data_dir_income_stratified_footprints = here("analysis", "preprocessing", "incom
##########################################################################################################################################################
# EUROSTAT Household Budget Survey (HBS)
-## load 'Mean consumption expenditure by income quintile (hbs_exp_t133)' table
+## load 'mean consumption expenditure by income quintile (hbs_exp_t133)' table
hbs_exp_t133 = read_csv(paste0(data_dir_income_stratified_footprints, "/hbs_exp_t133.csv"))
## rename and arrange by country
mean_expenditure_by_quintile = hbs_exp_t133 %>%
@@ -3216,34 +3216,39 @@ results_formatted = results_recombined %>%
write.csv(results_formatted, paste0(data_dir_income_stratified_footprints, "/results_formatted_method1_ixi.csv"))
write_rds(results_formatted, paste0(data_dir_income_stratified_footprints, "/results_formatted_method1_ixi.rds"))
-################################################### !!!! method 1 - PXP version - PPS HH NO RENT !!!! ####################################################
-##########################################################################################################################################################
-##########################################################################################################################################################
-
+################################################### main paper method, EXIOBASE product-by-product version ############################################### ##########################################################################################################################################################
-# Exiobase - pxp version
+# EXIOBASE product-by-product version
+## set study years
years_exb_pxp = c(2005,2010)
+## create NULL dataframe for dis-aggregated final demand
disaggregated_final_demand = NULL
+## create NULL dataframe for 'total intensity vectors'
TIVs = NULL
+## create NULL dataframe for domestic 'total intensity vectors'
domestic_TIVs = NULL
+## create NULL dataframe for rest-of-Europe 'total intensity vectors'
europe_TIVs = NULL
+## create NULL dataframe for national footprints
national_fp = NULL
+## create NULL dataframe for national territorial energy use/emissions
national_territorial = NULL
+## 'for' loop which populates the above NULL dataframes
for (i in years_exb_pxp){
year_current = i
+ ## EXIOBASE final demand table without labels (.csv file)
Exiobase_FD = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/FD_",year_current,"_pxp.csv"))[,-1]
- # select household final demand vectors for relevant countries - figure out how to soft code this
-
+ ## select household (hh) final demand vectors for study countries
AT = Exiobase_FD[,1]
BE = Exiobase_FD[,8]
BG = Exiobase_FD[,15]
@@ -3277,48 +3282,44 @@ for (i in years_exb_pxp){
Eurostat_countries = cbind(AT,BE,BG,CY,CZ,DE,DK,EE,EL,ES,FI,FR,HR,HU,IE,IT,LT,LU,LV,MT,NL,NO,PL,PT,RO,SE,SI,SK,TR,UK)
- # labels
-
- Exiobase_T_labels = read.csv(paste0(data_dir_income_stratified_footprints, "/data/Exiobase_T_labels_pxp_w_coicop_mapping_no_rent.csv")) %>%
+ ## read in EXIOBASE production sector labels written and saved in previous code chunk
+ Exiobase_T_labels = read.csv(paste0(data_dir_income_stratified_footprints, "/data/Exiobase_T_labels_pxp_w_coicop_mapping.csv")) %>%
mutate(V1 = dplyr::recode(V1,"GR" = "EL","GB" = "UK"))
-
- # hh fd with production sector labels
-
+ ## bind EXIOBASE hh final demand with production sector labels
hh_fd_with_production_sector_labels = cbind(Exiobase_T_labels,Eurostat_countries) %>% rename(geo = V1, sector = V2)
- # assumption of same purchase structure between quintiles of domestic and foreign final demand
-
- # replicate each cell of each country's hh final demand as many times as there are income groups in the HBS data - in this preliminary case:5
-
+ ## replicate each cell of each country's hh final demand as many times as there are income quantiles in the HBS data - in this case: 5
cells_repeat = data.frame(hh_fd_with_production_sector_labels %>% slice(rep(1:n(), each = 5)))
+ ## create quintile names
quintiles = data.frame(rep(c("QUINTILE1","QUINTILE2","QUINTILE3","QUINTILE4","QUINTILE5"),200)) %>% rename_at(1,~"quintile")
+ ## bind replicated hh final demand cells to quintile names
replicated = cbind(cells_repeat,quintiles) %>% rename(country_of_production = geo)
- # make fd data long
-
+ ## convert to long data
replicated_long = replicated %>% gather(geo, value,-sector,-coicop,-quintile,-five_sectors,-country_of_production)
+ ## create year column
year = as.character(rep(year_current,nrow(replicated_long)))
+ ## bind year with long hh fd data
replicated_long = cbind(year,replicated_long)
-
+ ## populate NULL 'dis-aggregated_final_demand' dataframe with hh fd from each study year
disaggregated_final_demand = rbind(disaggregated_final_demand, replicated_long)
-
- # TIVs
-
- # CO2 - combustion - air
-
+ # load total intensity vectors
+ ## CO2 - combustion - air
Exiobase_TIV_co2_combustion_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_co2_combustion_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_co2_combustion_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_combustion_air_", year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CO2_combustion_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_co2_combustion_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_combustion_air_", year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3341,14 +3342,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CO2_combustion_europe,TIV_CO2_combustion_not_europe)
- # CO2 - noncombustion - cement - air
-
+ ## CO2 - noncombustion - cement - air
Exiobase_TIV_co2_noncombustion_cement_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_co2_noncombustion_cement_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_co2_noncombustion_cement_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_noncombustion_cement_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CO2_noncombustion_cement_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_co2_noncombustion_cement_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_noncombustion_cement_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3371,14 +3373,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CO2_noncombustion_cement_europe,TIV_CO2_noncombustion_cement_not_europe)
- # CO2 - noncombustion - lime - air
-
+ ## CO2 - noncombustion - lime - air
Exiobase_TIV_co2_noncombustion_lime_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_co2_noncombustion_lime_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_co2_noncombustion_lime_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_noncombustion_lime_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CO2_noncombustion_lime_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_co2_noncombustion_lime_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_noncombustion_lime_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3401,14 +3404,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CO2_noncombustion_lime_europe,TIV_CO2_noncombustion_lime_not_europe)
- # CO2 - agriculture - peat decay - air
-
+ ## CO2 - agriculture - peat decay - air
Exiobase_TIV_co2_agriculture_peatdecay_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_co2_agriculture_peatdecay_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_co2_agriculture_peatdecay_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_agriculture_peatdecay_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CO2_agriculture_peatdecay_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_co2_agriculture_peatdecay_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_agriculture_peatdecay_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3431,14 +3435,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CO2_agriculture_peatdecay_europe,TIV_CO2_agriculture_peatdecay_not_europe)
- # CO2 - waste - biogenic - air
-
+ ## CO2 - waste - biogenic - air
Exiobase_TIV_co2_waste_biogenic_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_co2_biogenic_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_co2_waste_biogenic_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_biogenic_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CO2_waste_biogenic_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_co2_waste_biogenic_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_biogenic_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3461,14 +3466,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CO2_waste_biogenic_europe,TIV_CO2_waste_biogenic_not_europe)
- # CO2 - waste - fossil - air
-
+ ## CO2 - waste - fossil - air
Exiobase_TIV_co2_waste_fossil_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_co2_waste_fossil_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_co2_waste_fossil_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_waste_fossil_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CO2_waste_fossil_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_co2_waste_fossil_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_waste_fossil_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3491,16 +3497,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CO2_waste_fossil_europe,TIV_CO2_waste_fossil_not_europe)
-
-
- # CH4 - combustion -air
-
+ ## CH4 - combustion - air
Exiobase_TIV_ch4_combustion_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_combustion_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_combustion_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_combustion_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_combustion_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_combustion_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_combustion_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3523,14 +3528,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_combustion_europe,TIV_CH4_combustion_not_europe)
- # CH4 - noncombustion - gas - air
-
+ ## CH4 - noncombustion - gas - air
Exiobase_TIV_ch4_noncombustion_gas_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_noncombustion_gas_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_noncombustion_gas_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_gas_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_noncombustion_gas_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_noncombustion_gas_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_gas_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3553,14 +3559,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_noncombustion_gas_europe,TIV_CH4_noncombustion_gas_not_europe)
- # CH4 - noncombustion - oil - air
-
+ ## CH4 - noncombustion - oil - air
Exiobase_TIV_ch4_noncombustion_oil_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_noncombustion_oil_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_noncombustion_oil_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_oil_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_noncombustion_oil_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_noncombustion_oil_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_oil_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3583,14 +3590,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_noncombustion_oil_europe,TIV_CH4_noncombustion_oil_not_europe)
- # CH4 - noncombustion - anthracite - air
-
+ ## CH4 - noncombustion - anthracite - air
Exiobase_TIV_ch4_noncombustion_anthracite_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_noncombustion_anthracite_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_noncombustion_anthracite_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_anthracite_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_noncombustion_anthracite_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_noncombustion_anthracite_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_anthracite_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3613,14 +3621,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_noncombustion_anthracite_europe,TIV_CH4_noncombustion_anthracite_not_europe)
- # CH4 - noncombustion - bituminouscoal - air
-
+ ## CH4 - noncombustion - bituminouscoal - air
Exiobase_TIV_ch4_noncombustion_bituminouscoal_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_noncombustion_bituminouscoal_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_noncombustion_bituminouscoal_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_bituminouscoal_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_noncombustion_bituminouscoal_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_noncombustion_bituminouscoal_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_bituminouscoal_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3643,14 +3652,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_noncombustion_bituminouscoal_europe,TIV_CH4_noncombustion_bituminouscoal_not_europe)
- # CH4 - noncombustion - cokingcoal - air
-
+ ## CH4 - noncombustion - cokingcoal - air
Exiobase_TIV_ch4_noncombustion_cokingcoal_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_noncombustion_cokingcoal_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_noncombustion_cokingcoal_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_cokingcoal_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_noncombustion_cokingcoal_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_noncombustion_cokingcoal_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_cokingcoal_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3673,14 +3683,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_noncombustion_cokingcoal_europe,TIV_CH4_noncombustion_cokingcoal_not_europe)
- # CH4 - noncombustion - lignite - air
-
+ ## CH4 - noncombustion - lignite - air
Exiobase_TIV_ch4_noncombustion_lignite_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_noncombustion_lignite_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_noncombustion_lignite_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_lignite_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_noncombustion_lignite_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_noncombustion_lignite_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_lignite_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3703,14 +3714,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_noncombustion_lignite_europe,TIV_CH4_noncombustion_lignite_not_europe)
- # CH4 - noncombustion - subbituminouscoal - air
-
+ ## CH4 - noncombustion - subbituminouscoal - air
Exiobase_TIV_ch4_noncombustion_subbituminouscoal_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_noncombustion_subbituminouscoal_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_noncombustion_subbituminouscoal_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_subbituminouscoal_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_noncombustion_subbituminouscoal_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_noncombustion_subbituminouscoal_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_subbituminouscoal_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3733,14 +3745,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_noncombustion_subbituminouscoal_europe,TIV_CH4_noncombustion_subbituminouscoal_not_europe)
- # CH4 - noncombustion - oilrefinery - air
-
+ ## CH4 - noncombustion - oilrefinery - air
Exiobase_TIV_ch4_noncombustion_oilrefinery_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_noncombustion_oilrefinery_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_noncombustion_oilrefinery_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_oilrefinery_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_noncombustion_oilrefinery_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_noncombustion_oilrefinery_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_oilrefinery_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3763,14 +3776,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_noncombustion_oilrefinery_europe,TIV_CH4_noncombustion_oilrefinery_not_europe)
- # CH4 - agriculture - air
-
+ ## CH4 - agriculture - air
Exiobase_TIV_ch4_agriculture_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_agriculture_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_agriculture_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_agriculture_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_agriculture_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_agriculture_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_agriculture_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3793,14 +3807,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_agriculture_europe,TIV_CH4_agriculture_not_europe)
- # CH4 - waste - air
-
+ ## CH4 - waste - air
Exiobase_TIV_ch4_waste_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_ch4_CO2eq_waste_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_ch4_waste_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_waste_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_CH4_waste_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_ch4_waste_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_waste_air_", year_current, "_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3823,16 +3838,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_CH4_waste_europe,TIV_CH4_waste_not_europe)
-
-
- # N2O - combustion - air
-
+ ## N2O - combustion - air
Exiobase_TIV_n2o_combustion_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_n2o_CO2eq_combustion_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_n2o_combustion_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_n2o_CO2eq_combustion_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_N2O_combustion_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_n2o_combustion_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_n2o_CO2eq_combustion_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3855,14 +3869,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_N2O_combustion_europe,TIV_N2O_combustion_not_europe)
- # N2O - agriculture - air
-
+ ## N2O - agriculture - air
Exiobase_TIV_n2o_agriculture_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_n2o_CO2eq_agriculture_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_n2o_agriculture_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_n2o_CO2eq_agriculture_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_N2O_agriculture_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_n2o_agriculture_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_n2o_CO2eq_agriculture_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3885,14 +3900,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_N2O_agriculture_europe,TIV_N2O_agriculture_not_europe)
- # SF6 - air
-
+ ## SF6 - air
Exiobase_TIV_sf6_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_sf6_CO2eq_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_sf6_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_sf6_CO2eq_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_SF6_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_sf6_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_sf6_CO2eq_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3915,14 +3931,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_SF6_europe,TIV_SF6_not_europe)
- # HFC - air
-
+ ## HFC - air
Exiobase_TIV_hfc_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_hfc_CO2eq_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_hfc_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_hfc_CO2eq_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_HFC_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_hfc_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_hfc_CO2eq_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3945,14 +3962,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_HFC_europe,TIV_HFC_not_europe)
- # PFC - air
-
+ ## PFC - air
Exiobase_TIV_pfc_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_pfc_CO2eq_air_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_pfc_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_pfc_CO2eq_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_PFC_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_pfc_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_pfc_CO2eq_air_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -3975,14 +3993,15 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_PFC_europe,TIV_PFC_not_europe)
- # Energy use
-
+ ## energy use
Exiobase_TIV_energy_use_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_energy_carrier_use_",year_current,"_pxp.csv"))[,-1]
+ ### domestic part
Exiobase_TIV_country_breakdown_energy_use_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_energy_carrier_use_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
gather(country, TIV_energy_domestic)
+ ### non-domestic part
Exiobase_TIV_europe_breakdown_energy_use_bp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_energy_carrier_use_",year_current,"_pxp.csv"))[,-1] %>%
row_to_names(row_number = 1) %>%
mutate_at(vars(AT:ZA), funs(as.numeric(as.character(.)))) %>%
@@ -4005,9 +4024,7 @@ for (i in years_exb_pxp){
ZA) %>%
select(TIV_energy_europe,TIV_energy_not_europe)
-
- # join with labels
-
+ # bind production sector labels with total intensity vectors
TIV_with_labels = cbind(Exiobase_T_labels,
t(Exiobase_TIV_co2_combustion_bp),
t(Exiobase_TIV_co2_noncombustion_cement_bp),
@@ -4057,17 +4074,17 @@ for (i in years_exb_pxp){
TIV_energy = "t(Exiobase_TIV_energy_use_bp)") %>%
mutate(V1 = dplyr::recode(V1,"GR" = "EL","GB" = "UK"))
+ # create year column
year = as.character(rep(year_current,nrow(TIV_with_labels)))
- look = cbind(year,TIV_with_labels) %>%
+ # bind year with labeled total intensity vectors
+ TIVs_intermediate = cbind(year,TIV_with_labels) %>%
rename(country_of_production = V1, sector = V2)
+ # populate NULL total intensity vectors dataframe
+ TIVs = rbind(TIVs, TIVs_intermediate)
- TIVs = rbind(TIVs,look)
-
-
- # join domestic_TIVs with labels
-
+ # bind production sector labels with domestic total intensity vectors
domestic_TIV_with_labels = cbind(Exiobase_T_labels,
Exiobase_TIV_country_breakdown_co2_combustion_bp,
Exiobase_TIV_country_breakdown_co2_noncombustion_cement_bp %>% select(-country),
@@ -4095,9 +4112,11 @@ for (i in years_exb_pxp){
mutate(V1 = dplyr::recode(V1,"GR" = "EL","GB" = "UK"),
country = dplyr::recode(country, "GR" = "EL", "GB" = "UK"))
+ # create year column for domestic TIVs
year_domestic = as.character(rep(year_current,nrow(domestic_TIV_with_labels)))
- look_domestic = cbind(year_domestic,domestic_TIV_with_labels) %>%
+ # bind year with labeled domestic total intensity vectors
+ TIVs_domestic_intermediate = cbind(year_domestic,domestic_TIV_with_labels) %>%
rename(country_of_production = V1, sector = V2, geo = country, year = year_domestic) %>%
mutate(TIV_CO2_combustion_domestic = as.numeric(TIV_CO2_combustion_domestic),
TIV_CO2_noncombustion_cement_domestic = as.numeric(TIV_CO2_noncombustion_cement_domestic),
@@ -4123,10 +4142,10 @@ for (i in years_exb_pxp){
TIV_PFC_domestic = as.numeric(TIV_PFC_domestic),
TIV_energy_domestic = as.numeric(TIV_energy_domestic))
- domestic_TIVs = rbind(domestic_TIVs, look_domestic)
-
- # europe TIVs with labels
+ # populate NULL domestic TIVs dataframe
+ domestic_TIVs = rbind(domestic_TIVs, TIVs_domestic_intermediate)
+ # bind production sector labels with European TIVs
europe_TIV_with_labels = cbind(Exiobase_T_labels,
Exiobase_TIV_europe_breakdown_co2_combustion-bp,
Exiobase_TIV_europe_breakdown_co2_noncombustion_cement_bp,
@@ -4153,9 +4172,11 @@ for (i in years_exb_pxp){
Exiobase_TIV_europe_breakdown_energy_use_bp) %>%
mutate(V1 = dplyr::recode(V1,"GR" = "EL","GB" = "UK"))
+ # create year column for European TIVs
year_europe = as.character(rep(year_current,nrow(europe_TIV_with_labels)))
- look_europe = cbind(year_europe,europe_TIV_with_labels) %>%
+ # bind year with labeled European total intensity vectors
+ TIVs_europe_intermediate = cbind(year_europe,europe_TIV_with_labels) %>%
rename(country_of_production = V1, sector = V2, year = year_europe) %>%
mutate(TIV_CO2_combustion_europe = as.numeric(TIV_CO2_combustion_europe),
TIV_CO2_noncombustion_cement_europe = as.numeric(TIV_CO2_noncombustion_cement_europe),
@@ -4181,13 +4202,11 @@ for (i in years_exb_pxp){
TIV_PFC_europe = as.numeric(TIV_PFC_europe),
TIV_energy_europe = as.numeric(TIV_energy_europe))
- europe_TIVs = rbind(europe_TIVs, look_europe)
-
-
- # total national footprints
-
- # FD labels
+ # populate NULL Europe TIVs dataframe
+ europe_TIVs = rbind(europe_TIVs, TIVs_europe_intermediate)
+ # calculating total national footprints
+ ## read in EXIOBASE final demand labels and calculate national footprints for each study satellite extension
Exiobase_FD_labels = as.data.frame(t(read.csv(paste0(data_dir_exiobase, "/Exiobase_FD_labels_pxp.csv")))[-1,-3]) %>%
mutate(V1 = dplyr::recode(V1,"GR" = "EL","GB" = "UK"))
@@ -4237,9 +4256,7 @@ for (i in years_exb_pxp){
national_energy_footprints = Exiobase_FD * t(Exiobase_TIV_energy_use_bp)
-
- # together
-
+ ## bind FD labels and national footprints together
national_footprints_w_labels = cbind(Exiobase_FD_labels,
rowSums(t(national_CO2_combustion_footprints)),
rowSums(t(national_CO2_noncombustion_cement_footprints)),
@@ -4266,16 +4283,17 @@ for (i in years_exb_pxp){
rowSums(t(national_energy_footprints))) %>%
mutate(V1 = dplyr::recode(V1,"GR" = "EL","GB" = "UK"))
+ ## create year column for national footprints
year_national_fp = as.character(rep(year_current,nrow(national_footprints_w_labels)))
-
-
- # direct FD emissions
-
+
+ # satellite extensions direct from household final demand
+ ## load satellite extensions direct from household final demand table from 'data_dir_exiobase' and convert to numeric
direct_FD_extensions = read.csv(paste0(data_dir_exiobase, "/IOT_", year_current, "_pxp/satellite/F_hh.csv", sep = ""),row.names=NULL,as.is=TRUE)[3:1106,3:345]
direct_FD_extensions[is.na(direct_FD_extensions)]=0
direct_FD_extensions = mapply(direct_FD_extensions, FUN = as.numeric)
direct_FD_extensions = matrix(data=direct_FD_extensions,ncol=343,nrow=1104)
+ ## extract each study satellite extension
direct_FD_co2_combustion = direct_FD_extensions[24,]
direct_FD_co2_noncombustion_cement = direct_FD_extensions[93,]
direct_FD_co2_noncombustion_lime = direct_FD_extensions[94,]
@@ -4300,7 +4318,7 @@ for (i in years_exb_pxp){
direct_FD_pfc = direct_FD_extensions[426,]
direct_FD_energy = direct_FD_extensions[470,]
-
+ ## bind all extensions together
direct_FD_fp = data.frame(direct_FD_co2_combustion,
direct_FD_co2_noncombustion_cement,
direct_FD_co2_noncombustion_lime,
@@ -4325,7 +4343,8 @@ for (i in years_exb_pxp){
direct_FD_pfc,
direct_FD_energy)
- look_national_fp = as.data.frame(cbind(year_national_fp,
+ ## bind national footprints and direct extensions from hh final demand together (with year column)
+ national_fp_intermediate = as.data.frame(cbind(year_national_fp,
national_footprints_w_labels,
direct_FD_fp)) %>%
rename(year = year_national_fp,
@@ -4404,14 +4423,13 @@ for (i in years_exb_pxp){
energy,
direct_FD_energy)
-
- national_fp = rbind(national_fp, look_national_fp)
+ ## populate NULL national footprints dataframe
+ national_fp = rbind(national_fp, national_fp_intermediate)
# national territorial
-
+ ## load saved satellite extensions table (as .csv file) from 'data_dir_exiobase' and extract each study satellite extension
satellite = read.csv(paste0(data_dir_exiobase, "/IOT_", year_current, "_pxp/satellite/satellite_",year_current,"_pxp.csv"))[,-1]
-
CO2_combustion_air = satellite[24,]
CO2_noncombustion_cement_air = satellite[93,]
@@ -4472,7 +4490,7 @@ for (i in years_exb_pxp){
energy_carrier_use = satellite[470,]
-
+ ## bind together
territorial = data.frame(t(CO2_combustion_air),
t(CO2_noncombustion_cement_air),
t(CO2_noncombustion_lime_air),
@@ -4518,9 +4536,11 @@ for (i in years_exb_pxp){
SF6 = 20,
HFC = 21, PFC = 22, energy = 23)
+ ## create year column
year_territorial = as.character(rep(year_current,nrow(territorial)))
- look_territorial = as.data.frame(cbind(year_territorial,
+ ## bind year column, production sector labels and territorial satellite extensions
+ territorial_intermediate = as.data.frame(cbind(year_territorial,
Exiobase_T_labels,
territorial)) %>%
rename(year = year_territorial,
@@ -4528,33 +4548,33 @@ for (i in years_exb_pxp){
sector = V2) %>%
select(-coicop,-five_sectors)
- national_territorial = rbind(national_territorial, look_territorial)
-
+ ## populate NULL territorial dataframe
+ national_territorial = rbind(national_territorial, territorial_intermediate)
}
+# write national territorial extensions as .csv and .rds files (save to 'data_dir_income_stratified_footprints')
write.csv(national_territorial, paste0(data_dir_income_stratified_footprints, "/national_territorial_pxp.csv"))
write_rds(national_territorial, paste0(data_dir_income_stratified_footprints, "/national_territorial_pxp.rds"))
-
+# write national footprints as .csv and .rds files (save to 'data_dir_income_stratified_footprints')
write.csv(national_fp, paste0(data_dir_income_stratified_footprints, "/national_fp_pxp.csv"))
write_rds(national_fp, paste0(data_dir_income_stratified_footprints, "/national_fp_pxp.rds"))
-# calculate quintile shares within each sector
+# calculate quintile shares within each HBS coicop category
shares = join_expenditures %>%
group_by(coicop,geo,year) %>%
mutate(share = pps_coicop/sum(pps_coicop))
-# pre-processing
-
+# join disaggregated EXIOBASE hh final demand dataframe with HBS coicop shares and calculate disaggregated hh final demand using the shares
fd_exiobase = disaggregated_final_demand %>%
left_join(shares, by = c("year","geo","coicop","quintile")) %>%
mutate(disaggregated_fd = value*share) %>%
select(year,geo,quintile,country_of_production,sector,coicop,disaggregated_fd) %>%
spread(quintile,disaggregated_fd)
-# direct from FD - to go back to results without direct FD fp, do not run this next chunk and do not bind_rows with 'results'
-
+# direct extensions from households
+## load EUROSTAT energy use data (no Turkey)
env_ac_pefasu_no_TR = read_csv(paste0(data_dir_income_stratified_footprints, "/env_ac_pefasu_1_Data.csv")) %>%
filter(TIME == 2015) %>%
mutate(geo = dplyr::recode(GEO,"Austria" = "AT",
@@ -4596,15 +4616,17 @@ env_ac_pefasu_no_TR = read_csv(paste0(data_dir_income_stratified_footprints, "/e
HH_OTH = other_activities_by_households/total_activities_by_households) %>%
select(geo,HH_HEAT,HH_TRA,HH_OTH)
-
+## create column for Turkey using Bulgarian shares between HH_HEAT, HH_TRA, and HH_OTH
env_ac_pefasu_TR = env_ac_pefasu_no_TR %>%
filter(geo == "BG") %>%
mutate(geo = dplyr::recode(geo,
"BG" = "TR"))
+## bind EUROSTAT energy use data (no Turkey) with created Turkey column
env_ac_pefasu = rbind(env_ac_pefasu_no_TR,env_ac_pefasu_TR) %>%
gather(sector,share_of_total_energy,-geo)
+## load EUROSTAT emissions data
env_ac_ainah_r2 = read_csv(paste0(data_dir_income_stratified_footprints, "/env_ac_ainah_r2_1_Data.csv")) %>%
filter(TIME == 2015) %>%
mutate(geo = dplyr::recode(GEO,"Austria" = "AT",
@@ -4647,22 +4669,25 @@ env_ac_ainah_r2 = read_csv(paste0(data_dir_income_stratified_footprints, "/env_a
HH_OTH = other_activities_by_households/total_activities_by_households) %>%
select(geo,airpol,HH_HEAT,HH_TRA,HH_OTH)
-
+### create separate dataframe for CO2
env_ac_ainah_r2_co2 = env_ac_ainah_r2 %>%
filter(airpol == "Carbon dioxide") %>%
select(-airpol) %>%
gather(sector,share_of_total_co2,-geo)
+### for CH4
env_ac_ainah_r2_ch4 = env_ac_ainah_r2 %>%
filter(airpol == "Methane") %>%
select(-airpol) %>%
gather(sector,share_of_total_ch4,-geo)
+### for N2O
env_ac_ainah_r2_n2o = env_ac_ainah_r2 %>%
filter(airpol == "Nitrous oxide") %>%
select(-airpol) %>%
gather(sector,share_of_total_n2o,-geo)
+## filter saved EXIOBASE national footprints dataframe, select direct extensions and repeat each cell 3 times
direct_FD_fp_long = national_fp %>%
filter(fd_category == "Final consumption expenditure by households",
geo %in% c("AT",
@@ -4701,8 +4726,10 @@ direct_FD_fp_long = national_fp %>%
direct_FD_energy) %>%
slice(rep(1:n(), each = 3))
+## create column with direct household sector names (3 times each)
sector = rep(c("HH_HEAT","HH_TRA","HH_OTH"), nrow(direct_FD_fp_long)/3)
+## bind names with EXIOBASE direct extensions data, map hh sectors to COICOP categories, then bind with EUROSTAT energy and emissions data and dis-aggregate EXIOBASE direct extensions data using EUROSTAT energy and emissions shares between the three sectors/categories
direct_FD_fp_long_disagg = cbind(sector,direct_FD_fp_long) %>%
mutate(coicop = ifelse(sector == "HH_TRA","CP072",
ifelse(sector == "HH_HEAT","CP045","CP05")),
@@ -4744,6 +4771,7 @@ direct_FD_fp_long_disagg = cbind(sector,direct_FD_fp_long) %>%
disaggregated_direct_FD_n2o,
disaggregated_direct_FD_energy)
+## spread quintiles across columns for CO2
direct_FD_co2 = direct_FD_fp_long_disagg %>%
select(year,geo,sector,quintile,coicop,five_sectors,disaggregated_direct_FD_co2) %>%
spread(quintile,disaggregated_direct_FD_co2) %>%
@@ -4762,6 +4790,7 @@ direct_FD_co2 = direct_FD_fp_long_disagg %>%
q2_co2_domestic+q3_co2_domestic+
q4_co2_domestic+q5_co2_domestic)
+## spread quintiles across columns for CH4
direct_FD_ch4 = direct_FD_fp_long_disagg %>%
select(year,geo,sector,quintile,coicop,five_sectors,disaggregated_direct_FD_ch4) %>%
spread(quintile,disaggregated_direct_FD_ch4) %>%
@@ -4780,7 +4809,7 @@ direct_FD_ch4 = direct_FD_fp_long_disagg %>%
q2_ch4_domestic+q3_ch4_domestic+
q4_ch4_domestic+q5_ch4_domestic)
-
+## spread quintiles across columns for N2O
direct_FD_n2o = direct_FD_fp_long_disagg %>%
select(year,geo,sector,quintile,coicop,five_sectors,disaggregated_direct_FD_n2o) %>%
spread(quintile,disaggregated_direct_FD_n2o) %>%
@@ -4799,6 +4828,7 @@ direct_FD_n2o = direct_FD_fp_long_disagg %>%
q2_n2o_domestic+q3_n2o_domestic+
q4_n2o_domestic+q5_n2o_domestic)
+## spread quintiles across columns for energy use
direct_FD_energy = direct_FD_fp_long_disagg %>%
select(year,geo,sector,quintile,coicop,five_sectors,disaggregated_direct_FD_energy) %>%
spread(quintile,disaggregated_direct_FD_energy) %>%
@@ -4817,7 +4847,7 @@ direct_FD_energy = direct_FD_fp_long_disagg %>%
q2_energy_domestic+q3_energy_domestic+
q4_energy_domestic+q5_energy_domestic)
-
+# join spread direct from hh fd extensions to each other and create co2eq as well
direct_FD_fp_wide = direct_FD_co2 %>%
left_join(direct_FD_ch4, by = c("year","geo",
"sector","coicop",
@@ -4870,8 +4900,7 @@ direct_FD_fp_wide = direct_FD_co2 %>%
-q5_n2o_domestic,
-n2o_total_domestic)
-
-
+# create intermediate results dataframe (without direct extensions data) by joining dis-aggregated EXIOBASE fd with total intensity vectors (plus domestic, plus Europe) and calculating co2, co2eq and energy use (total, domestic, and Europe)
results = fd_exiobase %>%
left_join(TIVs, by = c("year", "country_of_production", "coicop", "sector")) %>%
left_join(europe_TIVs, by = c("year", "country_of_production", "coicop", "sector", "five_sectors")) %>%
@@ -5004,10 +5033,11 @@ results = fd_exiobase %>%
q5_energy_europe = QUINTILE5*(TIV_energy_europe - TIV_energy_domestic),
energy_total_europe = q1_energy_europe+q2_energy_europe+q3_energy_europe+q4_energy_europe+q5_energy_europe)
+# bind intermediate results dataframe with direct from hh fd extensions dataframe
results_with_direct_FD_fp = bind_rows(results,direct_FD_fp_wide)
-### create compressed results_pxp rds file
-
+# create compressed results rds file
+## clean names of results file
dat_all = results_with_direct_FD_fp %>%
clean_names()
@@ -5016,7 +5046,7 @@ sectors = dat_all %>%
distinct(sector) %>%
mutate(sector_id = row_number())
-# if interested in looking at a sectoral breakdown of the product-by-product version results, un-comment line below
+## write EXIOBASE production sector ids to derived data folder (as .csv file) for use in creating figures in the main paper (commented out atm - not used by current main paper figures)
#write_csv(sectors, paste0(here("/analysis/data/derived/si/sectors_method1_pxp.csv")))
# convert aggregated sector labels to IDs
@@ -5024,6 +5054,7 @@ sectors_agg = dat_all %>%
distinct(five_sectors) %>%
mutate(sector_agg_id = row_number())
+## write aggregated sector ids to derived data folder (as .csv file) for use in creating figures in the main paper. If interested in looking at a sectoral breakdown of the product-by-product version results, un-comment line below
#write_csv(sectors_agg, paste0(here("analysis/data/derived/si/sectors_agg_method1_pxp.csv")))
# convert COICOP labels to IDs
@@ -5031,6 +5062,7 @@ coicop = dat_all %>%
distinct(coicop) %>%
mutate(coicop_id = row_number())
+## write COICOP consumption category ids to derived data folder (as .csv file) for use in creating figures in the main paper (commented out atm - not used by current main paper figures)
#write_csv(coicop, paste0(here("analysis/data/derived/si/coicop_method1_pxp.csv")))
# replace sector text labels with numerical IDs (save space)
@@ -5052,7 +5084,7 @@ dat_results = dat_compressed %>%
group_by(year, geo, sector_id) %>%
summarise_if(is.numeric, sum, na.rm = TRUE)
-## extract final demand and pivot long
+# extract final demand and pivot long
cols_final_demand = c("quintile1", "quintile2", "quintile3", "quintile4", "quintile5")
tmp_fd = dat_results %>%
select(year, geo, sector_id, cols_final_demand) %>%
@@ -5062,7 +5094,7 @@ tmp_fd = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-## extract co2 and pivot long
+# extract co2 and pivot long
cols_co2 = c("q1_co2", "q2_co2", "q3_co2", "q4_co2", "q5_co2")
tmp_co2 = dat_results %>%
select(year, geo, sector_id, cols_co2) %>%
@@ -5072,7 +5104,7 @@ tmp_co2 = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-## extract co2 domestic and pivot long
+# extract co2 domestic and pivot long
cols_co2_domestic = c("q1_co2_domestic", "q2_co2_domestic", "q3_co2_domestic", "q4_co2_domestic", "q5_co2_domestic")
tmp_co2_domestic = dat_results %>%
select(year, geo, sector_id, cols_co2_domestic) %>%
@@ -5082,7 +5114,7 @@ tmp_co2_domestic = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-## extract co2 europe and pivot long
+# extract co2 europe and pivot long
cols_co2_europe = c("q1_co2_europe", "q2_co2_europe", "q3_co2_europe", "q4_co2_europe", "q5_co2_europe")
tmp_co2_europe = dat_results %>%
select(year, geo, sector_id, cols_co2_europe) %>%
@@ -5092,8 +5124,7 @@ tmp_co2_europe = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-
-## extract co2eq and pivot long
+# extract co2eq and pivot long
cols_co2eq = c("q1_co2eq", "q2_co2eq", "q3_co2eq", "q4_co2eq", "q5_co2eq")
tmp_co2eq = dat_results %>%
select(year, geo, sector_id, cols_co2eq) %>%
@@ -5103,7 +5134,7 @@ tmp_co2eq = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-## extract co2eq domestic and pivot long
+# extract co2eq domestic and pivot long
cols_co2eq_domestic = c("q1_co2eq_domestic", "q2_co2eq_domestic", "q3_co2eq_domestic", "q4_co2eq_domestic", "q5_co2eq_domestic")
tmp_co2eq_domestic = dat_results %>%
select(year, geo, sector_id, cols_co2eq_domestic) %>%
@@ -5113,7 +5144,7 @@ tmp_co2eq_domestic = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-## extract co2eq europe and pivot long
+# extract co2eq europe and pivot long
cols_co2eq_europe = c("q1_co2eq_europe", "q2_co2eq_europe", "q3_co2eq_europe", "q4_co2eq_europe", "q5_co2eq_europe")
tmp_co2eq_europe = dat_results %>%
select(year, geo, sector_id, cols_co2eq_europe) %>%
@@ -5123,7 +5154,7 @@ tmp_co2eq_europe = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-## extract energy use and pivot long
+# extract energy use and pivot long
cols_energy = c("q1_energy","q2_energy","q3_energy","q4_energy","q5_energy")
tmp_energy = dat_results %>%
select(year, geo, sector_id, cols_energy) %>%
@@ -5133,7 +5164,7 @@ tmp_energy = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-## extract energy domestic and pivot long
+# extract energy domestic and pivot long
cols_energy_domestic = c("q1_energy_domestic","q2_energy_domestic","q3_energy_domestic","q4_energy_domestic","q5_energy_domestic")
tmp_energy_domestic = dat_results %>%
select(year, geo, sector_id, cols_energy_domestic) %>%
@@ -5143,7 +5174,7 @@ tmp_energy_domestic = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-## extract energy europe and pivot long
+# extract energy europe and pivot long
cols_energy_europe = c("q1_energy_europe","q2_energy_europe","q3_energy_europe","q4_energy_europe","q5_energy_europe")
tmp_energy_europe = dat_results %>%
select(year, geo, sector_id, cols_energy_europe) %>%
@@ -5153,8 +5184,7 @@ tmp_energy_europe = dat_results %>%
mutate(quint = parse_number(quintile)) %>%
select(-quintile)
-### TODO: also convert to other indicators to this format (as blocks above)
-### TODO: left join all indicators back to "results_formated" like her with co2
+# recombine pivoted extensions
results_recombined = tmp_fd %>%
left_join(tmp_co2, by=c("year", "geo", "sector_id", "quint")) %>%
left_join(tmp_co2_domestic, by=c("year", "geo", "sector_id", "quint")) %>%
@@ -5166,14 +5196,13 @@ results_recombined = tmp_fd %>%
left_join(tmp_energy_domestic, by=c("year", "geo", "sector_id", "quint")) %>%
left_join(tmp_energy_europe, by = c("year", "geo", "sector_id", "quint"))
-
-
# finally re-join aggregated sector IDs
results_formatted = results_recombined %>%
left_join(sector_mapping, by="sector_id") %>%
ungroup() %>%
select(-coicop_id)
+# write formatted results files (.csv and .rds) to 'data_dir_income_stratified_footprints'
write.csv(results_formatted, paste0(data_dir_income_stratified_footprints, "/results_formatted_method1_pxp.csv"))
write_rds(results_formatted, paste0(data_dir_income_stratified_footprints, "/results_formatted_method1_pxp.rds"))