Climate Impact of Green and Non-Green Bonds at Country and Corporate Levels

Code
options(repos = c(CRAN = "https://cran.rstudio.com/"))

library(here) 
library(purrr)
library(readxl)
library(tidyr)
library(dplyr)
library(tidyverse)
Code
knitr::opts_knit$set(root.dir = "/Users/fangjun/Desktop/Sais/Year 2/Term 2/Sustainable Finance/Final Paper")
Code
annual_greenhouse_gas_air_emissions <- read.csv("emission by region.csv")
Code
library(here)
library(readxl)
# Construct the path
path_to_sheet <- here("green bond by industry.xlsx")
# Read the Excel file
green_bond_issuance_by_industry<- read_excel (path_to_sheet)
# View the data
green_bond_issuance_by_industry
# A tibble: 9 × 10
  Industry        `2014` `2015` `2016` `2017` `2018` `2019` `2020` `2021` `2022`
  <chr>            <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>
1 Energy            18.3   23.7   33.2   53     53.7   85    105.   201.   158. 
2 Buildings          7.5    8.3   17.9   46.3   46.9   81.4   82.6  165.   122. 
3 Transport          4.2    5.9   12.7   23.9   30.7   52.5   68.2   96.9   91.6
4 Water              2.7    4     10.7   19.8   18.1   24.9   18.1   36.7   33.2
5 Waste              1.1    1.5    4.4    6.1    7.6    9.7    7.8   24.6   24  
6 Land Use           1.5    0.6    1.7    5.3    7.3    8.9   15.5   30.9   26.5
7 Industry          NA      0.7    0.1    0.2    0.8    1.7    0.9    7.7    4.1
8 ICT               NA     NA     NA     NA      0.1    1.6    1.5    5.2    6.2
9 Unspecified A&R    1.2    0.9    3.9    3.5    3.5    2.6    2.4   11.4   18.6
Code
library(here)
library(readxl)
# Construct the path
path_to_sheet <- here("green bond by region.xlsx")
# Read the Excel file
green_bond_issuance_by_region<- read_excel (path_to_sheet)
# View the data
green_bond_issuance_by_region
# A tibble: 6 × 10
  Country       `2014` `2015`   `2016` `2017` `2018` `2019` `2020` `2021` `2022`
  <chr>          <dbl> <chr>     <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>
1 Africa           0.1 NA          0.2    0.3    0.2    0.9    1.2    0.4    0.3
2 Asia-Pacific     1.6 3.6        26.5   35.5   48.5   66.9   52.7  143.   133. 
3 Europe          17.9 19.6000…   25.8   60.6   67.1  120.   163.   294.   229. 
4 Latin America    0.2 1.10000…    1.6    4      2.4    5.2   10.9    9.1    3.1
5 North America    7.4 12.8       20.8   48.7   38.9   62.7   61.5  104.    76.6
6 Supranational    9.4 8.4        10.2    9.5   12.1   13.1   13.8   31.4   45.1
Code
library(dplyr)
green_bond_issuance_by_region <- green_bond_issuance_by_region %>%
  rename(Region = Country)
Code
library(dplyr)
#Select specific columns
selected_annual_greenhouse_gas_air_emissions <- select(annual_greenhouse_gas_air_emissions, 2, 11, 19:27)
Code
library(dplyr)

renamed_selected_annual_greenhouse_gas_air_emissions <- selected_annual_greenhouse_gas_air_emissions %>%
  mutate(
    Region = case_when(
      Country %in% c("Northern Africa", "Sub-Saharan Africa") ~ "Africa",
      Country %in% c("Central Asia", "Eastern Asia", "South-eastern Asia", "Southern Asia", "Western Asia") ~ "Asia-Pacific",
      Country %in% c("Eastern Europe", "Northern Europe", "Southern Europe", "Western Europe") ~ "Europe",
      Country == "Latin America and the Caribbean" ~ "Latin America",
      Country %in% c("Americas", "Northern America") ~ "North America",
      TRUE ~ "Other"
    )
  )
Code
library(dplyr)
#rearrange the columns to make "Region" the first column
renamed_selected_annual_greenhouse_gas_air_emissions <- renamed_selected_annual_greenhouse_gas_air_emissions %>%
  select(Region, everything())
Code
library(dplyr)
# Group by 'Region' and then summarise each numerical column by summing up
combined_annual_greenhouse_gas_air_emissions <- renamed_selected_annual_greenhouse_gas_air_emissions %>%
  group_by(Region) %>%
  summarise(across(where(is.numeric), sum, na.rm = TRUE))
Code
library(dplyr)
library(tidyr)

#| warning: false
# Convert all year columns to numeric
green_bonds_emissions <- combined_annual_greenhouse_gas_air_emissions %>%
  mutate(across(-Region, as.numeric, .names = "converted_{.col}"))

green_bonds_emissions_long <- green_bonds_emissions %>%
  pivot_longer(
    cols = starts_with("converted_"),  # Now selecting the converted columns
    names_to = "Year",
    names_prefix = "converted_X",  # Remove the 'converted_' prefix from the year names
    values_to = "GreenBondsEmissions"
  )

green_bonds_issuance <- green_bond_issuance_by_region %>%
  mutate(across(-Region, as.numeric, .names = "converted_{.col}"))

green_bonds_issuance_long <- green_bonds_issuance %>%
  pivot_longer(
    cols = starts_with("converted_"),  # Now selecting the converted columns
    names_to = "Year",
    names_prefix = "converted_",  # Remove the 'converted_' prefix from the year names
    values_to = "GreenBondsIssuance"
  )
Code
# Assuming emissions_long is also prepared
combined_green_bonds_emission_comparison <- left_join(green_bonds_issuance_long, green_bonds_emissions_long, by = c("Region", "Year"))
combined_green_bonds_emission_comparison
# A tibble: 54 × 22
   Region   `2014` `2015` `2016` `2017` `2018` `2019` `2020` `2021` `2022` Year 
   <chr>     <dbl> <chr>   <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl> <chr>
 1 Africa      0.1 NA        0.2    0.3    0.2    0.9    1.2    0.4    0.3 2014 
 2 Africa      0.1 NA        0.2    0.3    0.2    0.9    1.2    0.4    0.3 2015 
 3 Africa      0.1 NA        0.2    0.3    0.2    0.9    1.2    0.4    0.3 2016 
 4 Africa      0.1 NA        0.2    0.3    0.2    0.9    1.2    0.4    0.3 2017 
 5 Africa      0.1 NA        0.2    0.3    0.2    0.9    1.2    0.4    0.3 2018 
 6 Africa      0.1 NA        0.2    0.3    0.2    0.9    1.2    0.4    0.3 2019 
 7 Africa      0.1 NA        0.2    0.3    0.2    0.9    1.2    0.4    0.3 2020 
 8 Africa      0.1 NA        0.2    0.3    0.2    0.9    1.2    0.4    0.3 2021 
 9 Africa      0.1 NA        0.2    0.3    0.2    0.9    1.2    0.4    0.3 2022 
10 Asia-Pa…    1.6 3.6      26.5   35.5   48.5   66.9   52.7  143.   133.  2014 
# ℹ 44 more rows
# ℹ 11 more variables: GreenBondsIssuance <dbl>, X2014 <dbl>, X2015 <dbl>,
#   X2016 <dbl>, X2017 <dbl>, X2018 <dbl>, X2019 <dbl>, X2020 <dbl>,
#   X2021 <dbl>, X2022 <dbl>, GreenBondsEmissions <dbl>
Code
library(dplyr)

# Selecting specific columns
combined_green_bonds_emission_comparison2 <- combined_green_bonds_emission_comparison %>%
  select(Region, Year, GreenBondsIssuance, GreenBondsEmissions)
combined_green_bonds_emission_comparison2
# A tibble: 54 × 4
   Region       Year  GreenBondsIssuance GreenBondsEmissions
   <chr>        <chr>              <dbl>               <dbl>
 1 Africa       2014                 0.1              12576.
 2 Africa       2015                NA                12664.
 3 Africa       2016                 0.2              12835.
 4 Africa       2017                 0.3              13170.
 5 Africa       2018                 0.2              13471.
 6 Africa       2019                 0.9              13654.
 7 Africa       2020                 1.2              13116.
 8 Africa       2021                 0.4              13655.
 9 Africa       2022                 0.3              13815.
10 Asia-Pacific 2014                 1.6             108818.
# ℹ 44 more rows
Code
library(tidyr)
library(dplyr)

long_combined_green_bonds_emission_comparison2 <- combined_green_bonds_emission_comparison2 %>%
  pivot_longer(
    cols = c(GreenBondsIssuance, GreenBondsEmissions),
    names_to = "Type",
    values_to = "Value",
    values_drop_na = TRUE
  ) %>%
  mutate(Year = as.integer(Year)) # Convert Year to integer if it's not already
long_combined_green_bonds_emission_comparison2
# A tibble: 98 × 4
   Region  Year Type                  Value
   <chr>  <int> <chr>                 <dbl>
 1 Africa  2014 GreenBondsIssuance      0.1
 2 Africa  2014 GreenBondsEmissions 12576. 
 3 Africa  2015 GreenBondsEmissions 12664. 
 4 Africa  2016 GreenBondsIssuance      0.2
 5 Africa  2016 GreenBondsEmissions 12835. 
 6 Africa  2017 GreenBondsIssuance      0.3
 7 Africa  2017 GreenBondsEmissions 13170. 
 8 Africa  2018 GreenBondsIssuance      0.2
 9 Africa  2018 GreenBondsEmissions 13471. 
10 Africa  2019 GreenBondsIssuance      0.9
# ℹ 88 more rows
Code
library(dplyr)
library(ggplot2)

combined_green_bonds_emission_comparison2 %>%
 filter(!(Region %in% "Supranational")) %>%
 ggplot() +
 aes(x = Year, fill = Region, colour = Region, weight = GreenBondsIssuance) +
 geom_bar() +
 scale_fill_viridis_d(option = "viridis", direction = 1) +
 scale_color_viridis_d(option = "viridis", 
 direction = 1) +
 labs(x = "Year", y = "GHC Emission", title = "Regional Green Bonds Issuance Trends from 2014 to 2022", 
 subtitle = "Asia-Pacific Leads Increase in Issuancess as Europe Experiences Steady Growth Over Nine Years", 
 caption = "Source: IMF, Statistica") +
 theme_minimal() +
 theme(axis.text.x = element_text(angle = 45, hjust = 1), # Rotate x-axis labels
        axis.text.x.bottom = element_text(margin = margin(t = 5, r = 5))) + # Stagger x-axis labels)
 facet_wrap(vars(Region))

Code
library(dplyr)
library(ggplot2)

combined_green_bonds_emission_comparison2 %>%
 filter(!(Region %in% "Supranational")) %>%
 ggplot() +
 aes(x = Year, y = GreenBondsEmissions, fill = Region, 
 colour = Region) +
 geom_col() +
 scale_fill_viridis_d(option = "viridis", direction = 1) +
 scale_color_viridis_d(option = "viridis", 
 direction = 1) +
 labs(title = "\"Regional Greenhouse Gas Emission Trends from 2014 to 2022\"", subtitle = "Consistent High Emission Levels in Asia-Pacific and North America ", 
 caption = "Source: IMF") +
 theme_minimal() +
 theme(axis.text.x = element_text(angle = 45, hjust = 1), # Rotate x-axis labels
        axis.text.x.bottom = element_text(margin = margin(t = 5, r = 5))) + 
 facet_wrap(vars(Region))

Code
library(dplyr)

renamed_industry_selected_annual_greenhouse_gas_air_emissions <- selected_annual_greenhouse_gas_air_emissions %>%
  mutate(
    Industries = case_when(
      Industry == "Electricity, Gas, Steam and Air Conditioning Supply" ~ "Energy",
      Industry == "Construction" ~ "Buildings",
      Industry == "Transportation and Storage" ~ "Transport",
      Industry == "Water supply; sewerage, waste management and remediation activities" ~ "Water",
      Industry == "Total Households" ~ "Land Use",
      TRUE ~ "Other"
    )
  )%>%
 filter(Industries != "Other")
Code
library(dplyr)

renamed_green_bond_issuance_by_industry <- green_bond_issuance_by_industry%>%
rename(Industries = Industry) %>%
  filter(Industries %in% c("Energy", "Buildings", "Transport", "Water", "Land Use"))

renamed_green_bond_issuance_by_industry
# A tibble: 5 × 10
  Industries `2014` `2015` `2016` `2017` `2018` `2019` `2020` `2021` `2022`
  <chr>       <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>
1 Energy       18.3   23.7   33.2   53     53.7   85    105.   201.   158. 
2 Buildings     7.5    8.3   17.9   46.3   46.9   81.4   82.6  165.   122. 
3 Transport     4.2    5.9   12.7   23.9   30.7   52.5   68.2   96.9   91.6
4 Water         2.7    4     10.7   19.8   18.1   24.9   18.1   36.7   33.2
5 Land Use      1.5    0.6    1.7    5.3    7.3    8.9   15.5   30.9   26.5
Code
library(dplyr)
library(tidyr)

# Convert all year columns to numeric
green_bonds_emissions_by_industry <- renamed_industry_selected_annual_greenhouse_gas_air_emissions %>%
  mutate(across(-Industries, as.numeric, .names = "converted_{.col}"))
Warning: There were 2 warnings in `mutate()`.
The first warning was:
ℹ In argument: `across(-Industries, as.numeric, .names = "converted_{.col}")`.
Caused by warning:
! NAs introduced by coercion
ℹ Run `dplyr::last_dplyr_warnings()` to see the 1 remaining warning.
Code
green_bonds_emissions_by_industry_long <- green_bonds_emissions_by_industry %>%
  pivot_longer(
    cols = starts_with("converted_"),  # Now selecting the converted columns
    names_to = "Year",
    names_prefix = "converted_X",  # Remove the 'converted_' prefix from the year names
    values_to = "GreenBondsEmissions"
  )

converted_renamed_green_bond_issuance_by_industry <- renamed_green_bond_issuance_by_industry %>%
  mutate(across(-Industries, as.numeric, .names = "converted_{.col}"))

renamed_green_bond_issuance_by_industry_long <- converted_renamed_green_bond_issuance_by_industry %>%
  pivot_longer(
    cols = starts_with("converted_"),  # Now selecting the converted columns
    names_to = "Year",
    names_prefix = "converted_",  # Remove the 'converted_' prefix from the year names
    values_to = "GreenBondsIssuance"
  )
Code
combined_green_bonds_emission_industry_comparison <- left_join(renamed_green_bond_issuance_by_industry_long, green_bonds_emissions_by_industry_long, by = c("Industries", "Year"))
combined_green_bonds_emission_industry_comparison
# A tibble: 5,211 × 24
   Industries `2014` `2015` `2016` `2017` `2018` `2019` `2020` `2021` `2022`
   <chr>       <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>  <dbl>
 1 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
 2 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
 3 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
 4 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
 5 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
 6 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
 7 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
 8 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
 9 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
10 Energy       18.3   23.7   33.2     53   53.7     85   105.   201.   158.
# ℹ 5,201 more rows
# ℹ 14 more variables: Year <chr>, GreenBondsIssuance <dbl>, Country <chr>,
#   Industry <chr>, X2014 <dbl>, X2015 <dbl>, X2016 <dbl>, X2017 <dbl>,
#   X2018 <dbl>, X2019 <dbl>, X2020 <dbl>, X2021 <dbl>, X2022 <dbl>,
#   GreenBondsEmissions <dbl>
Code
library(dplyr)

# Selecting specific columns
combined_green_bonds_emission_industry_comparison2 <- combined_green_bonds_emission_industry_comparison %>%
  select(Industries, Year, GreenBondsIssuance, GreenBondsEmissions)
combined_green_bonds_emission_industry_comparison2
# A tibble: 5,211 × 4
   Industries Year  GreenBondsIssuance GreenBondsEmissions
   <chr>      <chr>              <dbl>               <dbl>
 1 Energy     2014                18.3            4246.   
 2 Energy     2014                18.3               2.43 
 3 Energy     2014                18.3            4394.   
 4 Energy     2014                18.3             101.   
 5 Energy     2014                18.3              44.6  
 6 Energy     2014                18.3             479.   
 7 Energy     2014                18.3             489.   
 8 Energy     2014                18.3               0.315
 9 Energy     2014                18.3               9.46 
10 Energy     2014                18.3            2614.   
# ℹ 5,201 more rows
Code
library(ggplot2)

ggplot(combined_green_bonds_emission_industry_comparison2) +
 aes(x = Year, y = GreenBondsIssuance, fill = Industries, 
 colour = Industries) +
 geom_col() +
 scale_fill_viridis_d(option = "viridis", direction = 1) +
 scale_color_viridis_d(option = "viridis", 
 direction = 1) +
 labs(title = "Dynamics of Green Bond Issuance Across Key Industries (2014-2022)", 
 subtitle = "Surge in Energy Sector Dominates as Transportation and Buildings Industries Show Notable Increases", 
 caption = "Source: Climate Bonds Initiatives") +
 theme_minimal() +
 theme(axis.text.x = element_text(angle = 45, hjust = 1), # Rotate x-axis labels
        axis.text.x.bottom = element_text(margin = margin(t = 5, r = 5)))+
 facet_wrap(vars(Industries))

Code
library(ggplot2)

ggplot(combined_green_bonds_emission_industry_comparison2) +
 aes(x = Year, y = GreenBondsEmissions, fill = Industries, 
 colour = Industries) +
 geom_col() +
 scale_fill_viridis_d(option = "viridis", direction = 1) +
 scale_color_viridis_d(option = "viridis", 
 direction = 1) +
 labs(title = "Trends in Green Bond Emission by Industry (2014-2022)", subtitle = "Steady Growth in Energy Sector Emission", 
 caption = "Source: IMF") +
 theme_minimal() +
 theme(axis.text.x = element_text(angle = 45, hjust = 1), # Rotate x-axis labels
        axis.text.x.bottom = element_text(margin = margin(t = 5, r = 5)))+
 facet_wrap(vars(Industries))

Code
library(here)
library(readxl)
# Construct the path
path_to_sheet <- here(" bond and company level data.xlsx")
# Read the Excel file
bond_and_company_level_data<- read_excel (path_to_sheet)
# View the data
bond_and_company_level_data
# A tibble: 18 × 10
   IssuerName           EUROutstanding\r\nAm…¹ Coupon MaturityDate Scope1 Scope2
   <chr>                <chr>                   <dbl> <chr>         <dbl>  <dbl>
 1 A2ASPA               400000000               1     16/07/2029      9.7    6.7
 2 EON1                 750000000               0.875 08/01/2025      3.6    6.6
 3 EON2                 850000000               1.25  19/10/2027      3.6    6.6
 4 ENGIESA1             480100000               0.875 27/03/2024      9.7    6.7
 5 ENGIESA2             750000000               2.12  30/03/2032      3.6    6.6
 6 FAURECIASE           400000000               2.38  15/06/2029      3.9    6.9
 7 IBERDROLA            750000000               3.12  22/11/2028      3.4    6.3
 8 KONINKLIJKEPHILIPSN… 750000000               0.5   22/05/2026      1.2    0.1
 9 KONINKLIJKEPHILIPSN… 650000000               2.12  05/11/2029      1.2    0.1
10 ORSTEDA/S            600000000               2.25  14/06/2028      3.4    6.3
11 POSTNLNV             300000000               0.625 23/09/2026      2.7    6.2
12 PROLOGIS             550000000               0.375 06/02/2028      2.8    6.2
13 RWEAG                1,000,000,000           2.12  24/05/2026      3.4    6.3
14 SEGROCAPITALSARL     650,000,000             1.25  23/03/2026      1      5.7
15 UNIBAIL-RODAMCO      643,748,000             2.5   26/02/2024      1      5.7
16 UPM-KYMMENEOYJ1      750,000.00c             0.125 19/11/2028     20.3    6.3
17 UPM-KYMMENEOYJ2      500,000,000             2.25  23/05/2029     98.6  103. 
18 VOLKSWAGENFINANCENV  1,250,000,000           0.875 22/09/2028      3.9    6.9
# ℹ abbreviated name: ¹​`EUROutstanding\r\nAmount`
# ℹ 4 more variables: `Scopes1+2` <dbl>, ITR <chr>, `Climate\r\nVaR` <chr>,
#   `SLB/SLGB\r\nOutstanding?` <chr>
Code
library(ggplot2)

ggplot(bond_and_company_level_data, aes(x = IssuerName, y = `EUROutstanding\r\nAmount`, fill = `EUROutstanding\r\nAmount`)) +
  geom_tile() +  # Removed the 'colour' aesthetic as it is typically not used with geom_tile()
  scale_fill_viridis_d(option = "viridis", direction = 1) +  # Fill scale for the tiles
  labs(title = "Green Bonds Value of Selected EU Companies", y = "EUR Outstanding Amount") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5, size = 8),  # Rotate x-axis labels
        axis.ticks.x = element_blank(),  # Remove x-axis ticks
        legend.position = "none")  # Remove the legend

Code
library(ggplot2)

ggplot(bond_and_company_level_data, aes(x = IssuerName, y = Scope1, fill = IssuerName)) +
  geom_col(show.legend = FALSE) +  # Remove the legend for the fill
  scale_fill_viridis_d(option = "viridis", direction = 1) +
  labs(title = "Financial Emission of Each Bond", subtitle = "UPM takes the lead", 
       caption = "Source: MSCI Research") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5), # Rotate x-axis labels
        axis.ticks.x = element_blank(), # Remove the x-axis ticks
        axis.title.x = element_blank(), # Remove the x-axis title if you wish
        legend.position = "none")  # Ensure the legend does not appear

Code
library(ggplot2)

ggplot(bond_and_company_level_data, aes(x = IssuerName, y = `Climate\r\nVaR`, fill = IssuerName)) +
  geom_tile() +
  scale_fill_viridis_d(option = "viridis", direction = 1) +
  labs(title = "Climate VaR for Each Green Bond", 
       y = "Climate Value at Risk",  # Add this if you want to label the y axis
       caption = "Source: MSCI Research") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 90, hjust = 1, vjust = 0.5), # Rotate x-axis labels
        axis.ticks.x = element_blank(), # Remove x-axis ticks
        legend.position = "none", # Remove the legend if it's redundant
        axis.title.x = element_blank()) # Optionally remove x-axis title for more space