diff --git a/analysis/preprocessing/full_code.Rmd b/analysis/preprocessing/full_code.Rmd
index e5e5aee96722f7b5adba8652560401130db45926..9e8deca663e24a24f2f154c924573b44b77d80c3 100644
--- a/analysis/preprocessing/full_code.Rmd
+++ b/analysis/preprocessing/full_code.Rmd
@@ -55,83 +55,93 @@ Each year is available as a .zip file ('IOT_year_ixi' or 'IOT_year_pxp') from th
 
 ```{r exiobase, eval = FALSE}
 
-# EXIOBASE_cluster_ixi_version
-
-# data directories (on cluster)
-#data_dir_exiobase = paste("/",file.path("data","metab","Exiobase", fsep=.Platform$file.sep),sep="")
+# set data directory for EXIOBASE files
 data_dir_exiobase = here("analysis", "preprocessing", "EXIOBASE")
 
+##### EXIOBASE industry-by-industry version
+
+# set study years
 years_exiobase_ixi = c(2005,2010,2015)
 
+# 'for' loop which writes 'total intensity vectors' (and row-wise breakdowns) for all study satellite extensions and years to 'data_dir_exiobase' using downloaded EXIOBASE files
 for (i in years_exiobase_ixi){
   
 year_current = i
 
+# read in A table as table
 A = read.delim(paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/A.txt"),header = F)
 
+# write A table as .csv file
 write.csv(A, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/A.csv"))
 
+# read in A table as .csv file and extract the data only (no labels), and convert to numeric
 A = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/A.csv",sep = ""),row.names=NULL,as.is=TRUE)[4:7990,4:7990]
 A[is.na(A)]=0
 A = mapply(A, FUN = as.numeric)
 A = matrix(data = A, ncol = 7987, nrow = 7987)
 
-L = solve(diag(dim(A)[1])-A)  # this solves the Leontief inverse initially
+# solve the Leontief inverse
+L = solve(diag(dim(A)[1])-A)  
 L[is.na(L)]=0
 
-
-# final demand
+# read in final demand table (Y) as table
 FD = read.delim(paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/Y.txt"),header = F)
 
+# write final demand table (Y) as .csv file
 write.csv(FD, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/Y.csv"))
 
-
+# read in final demand table as .csv file and extract the production sector labels
 Exiobase_T_labels_ixi = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/Y.csv"))[4:7990,1:3]
 
+# write a .csv file with only the production sector labels
 write.csv(Exiobase_T_labels_ixi, paste0(data_dir_exiobase, "/Exiobase_T_labels_ixi.csv"))
 
+# read in final demand table as .csv file and extract the final demand category labels
 Exiobase_FD_labels_ixi = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/Y.csv"))[1:3,4:346]
 
+# write a .csv file with only the final demand category labels
 write.csv(Exiobase_FD_labels_ixi, paste0(data_dir_exiobase, "/Exiobase_FD_labels_ixi.csv"))
 
-
+# read in final demand table as .csv file and extract the data only (no labels), and convert to numeric
 FD = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/Y.csv",sep=""),row.names=NULL,as.is=TRUE)[4:7990,4:346]
 FD[is.na(FD)]=0
 FD = mapply(FD, FUN = as.numeric)
 FD = matrix(data=FD,ncol=343,nrow=7987)
 
+# write a .csv file with final demand data only (no labels)
 write.csv(FD, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/FD_",year_current,"_ixi.csv"))
 
-
-# total output
+# calculate total output
 total_output = L %*% rowSums(FD)
 
+# write total output as a .csv file
 write.csv(total_output, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/total_output_",year_current,"_ixi.csv"))
 
-
-# direct environmental vectors
+# read in satellite extensions table (F) as table
 satellite = read.delim(paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/satellite/F.txt"),header = F)
 
+# write satellite extensions table (F) as .csv file
 write.csv(satellite, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/satellite/F.csv"))
 
-
+# read in satellite extensions table (F) as .csv file and extract the data only (no labels), and convert to numeric
 satellite = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/satellite/F.csv",sep=""),row.names=NULL,as.is=TRUE)[3:1115,3:7989]
 satellite[is.na(satellite)]=0
 satellite = mapply(satellite, FUN = as.numeric)
 satellite = matrix(data=satellite,ncol=7987,nrow=1113)
 
-
+# write a .csv file with satellite extensions data only (no labels)
 write.csv(satellite, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/satellite/satellite_",year_current,"_ixi.csv"))
 
-
-# direct environmental vectors on final demand
-
+# read in satellite extensions on final demand table (F_hh) as table
 satellite_FD = read.delim(paste0(data_dir_exiobase, "/IOT_", year_current, "_ixi/satellite/F_hh.txt"),header = F)
 
+# write satellite extensions on final demand table (F_hh) as .csv file
 write.csv(satellite_FD, paste0(data_dir_exiobase, "/IOT_", year_current, "_ixi/satellite/F_hh.csv"))
 
+## extract the relevant satellite extensions from the satellite table, calculate the 'total intensity 
+## vectors' (and their row-wise breakdowns), and save to 'data_dir_exiobase'
 
-# CO2 combustion air
+# CO2 - combustion - air
 CO2_combustion_air = satellite[24,]
 DIV_co2_combustion_air = CO2_combustion_air/total_output
 DIV_co2_combustion_air[is.na(DIV_co2_combustion_air)]=0
@@ -149,8 +159,7 @@ TIV_country_breakdown_co2_combustion_air_w_labels = t(TIV_breakdown_co2_combusti
 
 write.csv(TIV_country_breakdown_co2_combustion_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_co2_combustion_air_",year_current,"_ixi.csv"))
 
-
-# CO2 non-combustion air
+# CO2 - non-combustion - air
 ## cement
 CO2_noncombustion_cement_air = satellite[93,]
 DIV_co2_noncombustion_cement_air = CO2_noncombustion_cement_air/total_output
@@ -187,8 +196,7 @@ TIV_country_breakdown_co2_noncombustion_lime_air_w_labels = t(TIV_breakdown_co2_
 
 write.csv(TIV_country_breakdown_co2_noncombustion_lime_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_co2_noncombustion_lime_air_",year_current,"_ixi.csv"))
 
-
-# CO2 agriculture peat decay air
+# CO2 - agriculture - peat decay - air
 CO2_agriculture_peatdecay_air = satellite[428,]
 DIV_co2_agriculture_peatdecay_air = CO2_agriculture_peatdecay_air/total_output
 DIV_co2_agriculture_peatdecay_air[is.na(DIV_co2_agriculture_peatdecay_air)]=0
@@ -206,7 +214,7 @@ TIV_country_breakdown_co2_agriculture_peatdecay_air_w_labels = t(TIV_breakdown_c
 
 write.csv(TIV_country_breakdown_co2_agriculture_peatdecay_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_co2_agriculture_peatdecay_air_",year_current,"_ixi.csv"))
 
-# CO2 waste air
+# CO2 - waste - air
 ## biogenic
 CO2_waste_biogenic_air = satellite[438,]
 DIV_co2_waste_biogenic_air = CO2_waste_biogenic_air/total_output
@@ -243,8 +251,7 @@ TIV_country_breakdown_co2_waste_fossil_air_w_labels = t(TIV_breakdown_co2_waste_
 
 write.csv(TIV_country_breakdown_co2_waste_fossil_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_co2_waste_fossil_air_",year_current,"_ixi.csv"))
 
-
-# CH4 combustion air
+# CH4 - combustion - air
 CH4_combustion_air = satellite[25,]
 CH4_combustion_air = CH4_combustion_air*28
 DIV_ch4_combustion_air = CH4_combustion_air/total_output
@@ -263,8 +270,7 @@ TIV_country_breakdown_ch4_combustion_air_w_labels = t(TIV_breakdown_ch4_combusti
 
 write.csv(TIV_country_breakdown_ch4_combustion_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_ch4_CO2eq_combustion_air_",year_current,"_ixi.csv"))
 
-
-# CH4 noncombustion air
+# CH4 - non-combustion - air
 ## gas
 CH4_noncombustion_gas_air = satellite[68,]
 CH4_noncombustion_gas_air = CH4_noncombustion_gas_air*28
@@ -417,8 +423,7 @@ TIV_country_breakdown_ch4_noncombustion_oilrefinery_air_w_labels = t(TIV_breakdo
 
 write.csv(TIV_country_breakdown_ch4_noncombustion_oilrefinery_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_ch4_CO2eq_noncombustion_oilrefinery_air_",year_current,"_ixi.csv"))
 
-
-# CH4 agriculture air
+# CH4 - agriculture - air
 CH4_agriculture_air = satellite[427,]
 CH4_agriculture_air = CH4_agriculture_air*28
 DIV_ch4_agriculture_air = CH4_agriculture_air/total_output
@@ -437,8 +442,7 @@ TIV_country_breakdown_ch4_agriculture_air_w_labels = t(TIV_breakdown_ch4_agricul
 
 write.csv(TIV_country_breakdown_ch4_agriculture_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_ch4_CO2eq_agriculture_air_",year_current,"_ixi.csv"))
 
-
-# CH4 waste air
+# CH4 - waste - air
 CH4_waste_air = satellite[436,]
 CH4_waste_air = CH4_waste_air*28
 DIV_ch4_waste_air = CH4_waste_air/total_output
@@ -457,8 +461,7 @@ TIV_country_breakdown_ch4_waste_air_w_labels = t(TIV_breakdown_ch4_waste_air_w_l
 
 write.csv(TIV_country_breakdown_ch4_waste_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_ch4_CO2eq_waste_air_",year_current,"_ixi.csv"))
 
-
-# N2O combustion air
+# N2O - combustion - air
 N2O_combustion_air = satellite[26,]
 N2O_combustion_air = N2O_combustion_air*265
 DIV_n2o_combustion_air = N2O_combustion_air/total_output
@@ -477,8 +480,7 @@ TIV_country_breakdown_n2o_combustion_air_w_labels = t(TIV_breakdown_n2o_combusti
 
 write.csv(TIV_country_breakdown_n2o_combustion_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_n2o_CO2eq_combustion_air_",year_current,"_ixi.csv"))
 
-
-# N2O agriculture air
+# N2O - agriculture - air
 N2O_agriculture_air = satellite[430,]
 N2O_agriculture_air = N2O_agriculture_air*265
 DIV_n2o_agriculture_air = N2O_agriculture_air/total_output
@@ -497,8 +499,7 @@ TIV_country_breakdown_n2o_agriculture_air_w_labels = t(TIV_breakdown_n2o_agricul
 
 write.csv(TIV_country_breakdown_n2o_agriculture_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_n2o_CO2eq_agriculture_air_",year_current,"_ixi.csv"))
 
-
-# SF6 air
+# SF6 - air
 SF6_air = satellite[424,]
 SF6_air = SF6_air*23500
 DIV_sf6_air = SF6_air/total_output
@@ -517,8 +518,7 @@ TIV_country_breakdown_sf6_air_w_labels = t(TIV_breakdown_sf6_air_w_labels %>%
 
 write.csv(TIV_country_breakdown_sf6_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_sf6_CO2eq_air_",year_current,"_ixi.csv"))
 
-
-# HFC air
+# HFC - air
 HFC_air = satellite[425,]
 DIV_hfc_air = HFC_air/total_output
 DIV_hfc_air[is.na(DIV_hfc_air)]=0
@@ -536,8 +536,7 @@ TIV_country_breakdown_hfc_air_w_labels = t(TIV_breakdown_hfc_air_w_labels %>%
 
 write.csv(TIV_country_breakdown_hfc_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_hfc_CO2eq_air_",year_current,"_ixi.csv"))
 
-
-# PFC air
+# PFC - air
 PFC_air = satellite[426,]
 DIV_pfc_air = PFC_air/total_output
 DIV_pfc_air[is.na(DIV_pfc_air)]=0
@@ -555,8 +554,7 @@ TIV_country_breakdown_pfc_air_w_labels = t(TIV_breakdown_pfc_air_w_labels %>%
 
 write.csv(TIV_country_breakdown_pfc_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/TIV_country_breakdown_pfc_CO2eq_air_",year_current,"_ixi.csv"))
 
-
-# energy carrier use
+# Energy carrier - use
 energy_carrier_use = satellite[470,]
 write.csv(energy_carrier_use, paste0(data_dir_exiobase, "/IOT_",year_current,"_ixi/satellite/energy_carrier_use_",year_current,"_ixi.csv"))
 DIV_e_u = energy_carrier_use/total_output
@@ -577,79 +575,90 @@ write.csv(TIV_country_breakdown_e_u_w_labels, paste0(data_dir_exiobase, "/IOT_",
 }
 
 
-# EXIOBASE_cluster_pxp_version
+##### EXIOBASE product-by-product version
 
+# set study years
 years_exiobase_pxp = c(2005,2010)
 
+# 'for' loop which writes 'total intensity vectors' (and row-wise breakdowns) for all study satellite extensions and years to 'data_dir_exiobase' using downloaded EXIOBASE files
 for (i in years_exiobase_pxp){
   
 year_current = i
 
+# read in A table as table
 A = read.delim(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/A.txt"),header = F)
 
+# write A table as .csv file
 write.csv(A, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/A.csv"))
 
+# read in A table as .csv file and extract the data only (no labels), and convert to numeric
 A = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/A.csv",sep = ""),row.names=NULL,as.is=TRUE)[4:9803,4:9803]
 A[is.na(A)]=0
 A = mapply(A, FUN = as.numeric)
 A = matrix(data = A, ncol = 9800, nrow = 9800)
 
-L = solve(diag(dim(A)[1])-A)  # this solves the Leontief inverse initially
+# solve the Leontief inverse
+L = solve(diag(dim(A)[1])-A)
 L[is.na(L)]=0
 
-
-# final demand
+# read in final demand table (Y) as table
 FD = read.delim(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/Y.txt"),header = F)
 
+# write final demand table (Y) as .csv file
 write.csv(FD, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/Y.csv"))
 
-
+# read in final demand table as .csv file and extract the production sector labels
 Exiobase_T_labels_pxp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/Y.csv"))[4:9803,1:3]
 
+# write a .csv file with only the production sector labels
 write.csv(Exiobase_T_labels_pxp, paste0(data_dir_exiobase, "/Exiobase_T_labels_pxp.csv"))
 
+# read in final demand table as .csv file and extract the final demand category labels
 Exiobase_FD_labels_pxp = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/Y.csv"))[1:3,4:346]
 
+# write a .csv file with only the final demand category labels
 write.csv(Exiobase_FD_labels_pxp, paste0(data_dir_exiobase, "/Exiobase_FD_labels_pxp.csv"))
 
-
+# read in final demand table as .csv file and extract the data only (no labels), and convert to numeric
 FD = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/Y.csv",sep=""),row.names=NULL,as.is=TRUE)[4:9803,4:346]
 FD[is.na(FD)]=0
 FD = mapply(FD, FUN = as.numeric)
 FD = matrix(data=FD,ncol=343,nrow=9800)
 
+# write a .csv file with final demand data only (no labels)
 write.csv(FD, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/FD_",year_current,"_pxp.csv"))
 
-
-# total output
+# calculate total output
 total_output = L %*% rowSums(FD)
 
+# write total output as a .csv file
 write.csv(total_output, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/total_output_",year_current,"_pxp.csv"))
 
-
-# direct environmental vectors
+# read in satellite extensions table (F) as table
 satellite = read.delim(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/satellite/F.txt"),header = F)
 
+# write satellite extensions table (F) as .csv file
 write.csv(satellite, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/satellite/F.csv"))
 
-
+# read in satellite extensions table (F) as .csv file and extract the data only (no labels), and convert to numeric
 satellite = read.csv(paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/satellite/F.csv",sep=""),row.names=NULL,as.is=TRUE)[3:1106,3:9802]
 satellite[is.na(satellite)]=0
 satellite = mapply(satellite, FUN = as.numeric)
 satellite = matrix(data=satellite,ncol=9800,nrow=1104)
 
-
+# write a .csv file with satellite extensions data only (no labels)
 write.csv(satellite, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/satellite/satellite_",year_current,"_pxp.csv"))
 
-
-# direct environmental vectors on final demand
-
+# read in satellite extensions on final demand table (F_hh) as table
 satellite_FD = read.delim(paste0(data_dir_exiobase, "/IOT_", year_current, "_pxp/satellite/F_hh.txt"),header = F)
 
+# write satellite extensions on final demand table (F_hh) as .csv file
 write.csv(satellite_FD, paste0(data_dir_exiobase, "/IOT_", year_current, "_pxp/satellite/F_hh.csv"))
 
+## extract the relevant satellite extensions from the satellite table, calculate the 'total intensity 
+## vectors' (and their row-wise breakdowns), and save to 'data_dir_exiobase'
 
-# CO2 combustion air
+# CO2 - combustion - air
 CO2_combustion_air = satellite[24,]
 DIV_co2_combustion_air = CO2_combustion_air/total_output
 DIV_co2_combustion_air[is.na(DIV_co2_combustion_air)]=0
@@ -667,8 +676,7 @@ TIV_country_breakdown_co2_combustion_air_w_labels = t(TIV_breakdown_co2_combusti
 
 write.csv(TIV_country_breakdown_co2_combustion_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_combustion_air_",year_current,"_pxp.csv"))
 
-
-# CO2 non-combustion air
+# CO2 - non-combustion - air
 ## cement
 CO2_noncombustion_cement_air = satellite[93,]
 DIV_co2_noncombustion_cement_air = CO2_noncombustion_cement_air/total_output
@@ -705,8 +713,7 @@ TIV_country_breakdown_co2_noncombustion_lime_air_w_labels = t(TIV_breakdown_co2_
 
 write.csv(TIV_country_breakdown_co2_noncombustion_lime_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_noncombustion_lime_air_",year_current,"_pxp.csv"))
 
-
-# CO2 agriculture peat decay air
+# CO2 - agriculture - peat decay - air
 CO2_agriculture_peatdecay_air = satellite[428,]
 DIV_co2_agriculture_peatdecay_air = CO2_agriculture_peatdecay_air/total_output
 DIV_co2_agriculture_peatdecay_air[is.na(DIV_co2_agriculture_peatdecay_air)]=0
@@ -724,8 +731,7 @@ TIV_country_breakdown_co2_agriculture_peatdecay_air_w_labels = t(TIV_breakdown_c
 
 write.csv(TIV_country_breakdown_co2_agriculture_peatdecay_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_agriculture_peatdecay_air_",year_current,"_pxp.csv"))
 
-
-# CO2 waste air
+# CO2 - waste - air
 ## biogenic
 CO2_waste_biogenic_air = satellite[438,]
 DIV_co2_waste_biogenic_air = CO2_waste_biogenic_air/total_output
@@ -762,8 +768,7 @@ TIV_country_breakdown_co2_waste_fossil_air_w_labels = t(TIV_breakdown_co2_waste_
 
 write.csv(TIV_country_breakdown_co2_waste_fossil_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_co2_waste_fossil_air_",year_current,"_pxp.csv"))
 
-
-# CH4 combustion air
+# CH4 - combustion - air
 CH4_combustion_air = satellite[25,]
 CH4_combustion_air = CH4_combustion_air*28
 DIV_ch4_combustion_air = CH4_combustion_air/total_output
@@ -782,8 +787,7 @@ TIV_country_breakdown_ch4_combustion_air_w_labels = t(TIV_breakdown_ch4_combusti
 
 write.csv(TIV_country_breakdown_ch4_combustion_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_combustion_air_",year_current,"_pxp.csv"))
 
-
-# CH4 non-combustion air
+# CH4 - non-combustion - air
 ## gas
 CH4_noncombustion_gas_air = satellite[68,]
 CH4_noncombustion_gas_air = CH4_noncombustion_gas_air*28
@@ -936,8 +940,7 @@ TIV_country_breakdown_ch4_noncombustion_oilrefinery_air_w_labels = t(TIV_breakdo
 
 write.csv(TIV_country_breakdown_ch4_noncombustion_oilrefinery_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_noncombustion_oilrefinery_air_",year_current,"_pxp.csv"))
 
-
-# CH4 agriculture air
+# CH4 - agriculture - air
 CH4_agriculture_air = satellite[427,]
 CH4_agriculture_air = CH4_agriculture_air*28
 DIV_ch4_agriculture_air = CH4_agriculture_air/total_output
@@ -956,8 +959,7 @@ TIV_country_breakdown_ch4_agriculture_air_w_labels = t(TIV_breakdown_ch4_agricul
 
 write.csv(TIV_country_breakdown_ch4_agriculture_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_agriculture_air_",year_current,"_pxp.csv"))
 
-
-# CH4 waste air
+# CH4 - waste - air
 CH4_waste_air = satellite[436,]
 CH4_waste_air = CH4_waste_air*28
 DIV_ch4_waste_air = CH4_waste_air/total_output
@@ -976,8 +978,7 @@ TIV_country_breakdown_ch4_waste_air_w_labels = t(TIV_breakdown_ch4_waste_air_w_l
 
 write.csv(TIV_country_breakdown_ch4_waste_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_ch4_CO2eq_waste_air_",year_current,"_pxp.csv"))
 
-
-# N2O combustion air
+# N2O - combustion - air
 N2O_combustion_air = satellite[26,]
 N2O_combustion_air = N2O_combustion_air*265
 DIV_n2o_combustion_air = N2O_combustion_air/total_output
@@ -996,8 +997,7 @@ TIV_country_breakdown_n2o_combustion_air_w_labels = t(TIV_breakdown_n2o_combusti
 
 write.csv(TIV_country_breakdown_n2o_combustion_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_n2o_CO2eq_combustion_air_",year_current,"_pxp.csv"))
 
-
-# N2O agriculture air
+# N2O - agriculture - air
 N2O_agriculture_air = satellite[430,]
 N2O_agriculture_air = N2O_agriculture_air*265
 DIV_n2o_agriculture_air = N2O_agriculture_air/total_output
@@ -1016,8 +1016,7 @@ TIV_country_breakdown_n2o_agriculture_air_w_labels = t(TIV_breakdown_n2o_agricul
 
 write.csv(TIV_country_breakdown_n2o_agriculture_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_n2o_CO2eq_agriculture_air_",year_current,"_pxp.csv"))
 
-
-# SF6 air
+# SF6 - air
 SF6_air = satellite[424,]
 SF6_air = SF6_air*23500
 DIV_sf6_air = SF6_air/total_output
@@ -1036,8 +1035,7 @@ TIV_country_breakdown_sf6_air_w_labels = t(TIV_breakdown_sf6_air_w_labels %>%
 
 write.csv(TIV_country_breakdown_sf6_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_sf6_CO2eq_air_",year_current,"_pxp.csv"))
 
-
-# HFC air
+# HFC - air
 HFC_air = satellite[425,]
 DIV_hfc_air = HFC_air/total_output
 DIV_hfc_air[is.na(DIV_hfc_air)]=0
@@ -1055,8 +1053,7 @@ TIV_country_breakdown_hfc_air_w_labels = t(TIV_breakdown_hfc_air_w_labels %>%
 
 write.csv(TIV_country_breakdown_hfc_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_hfc_CO2eq_air_",year_current,"_pxp.csv"))
 
-
-# PFC air
+# PFC - air
 PFC_air = satellite[426,]
 DIV_pfc_air = PFC_air/total_output
 DIV_pfc_air[is.na(DIV_pfc_air)]=0
@@ -1074,8 +1071,7 @@ TIV_country_breakdown_pfc_air_w_labels = t(TIV_breakdown_pfc_air_w_labels %>%
 
 write.csv(TIV_country_breakdown_pfc_air_w_labels, paste0(data_dir_exiobase, "/IOT_",year_current,"_pxp/TIV_country_breakdown_pfc_CO2eq_air_",year_current,"_pxp.csv"))
 
-
-# energy carrier use
+# Energy carrier - use
 energy_carrier_use = satellite[470,]
 DIV_e_u = energy_carrier_use/total_output
 DIV_e_u[is.na(DIV_e_u)]=0