--- title: "Equitable Climate Change Actions, Transparent Planning, and Cooperation: a Path to Efficient and Sustainable Growth" author: "Ekaterina Melianova & Artem Volgin" output: html_document: toc: true toc_depth: 2 number_sections: yes fig_width: 8 fig_height: 6 theme: cosmo code_folding: hide --- ```{r message=FALSE, warning=FALSE, include=FALSE} # Libraries install.packages('deaR') library(foreign) library(dplyr) library(tidyr) library(stringr) library(tidyverse) library(quantmod) library(kableExtra) library(deaR) library(ggpubr) library(pheatmap) library(igraph) library(maditr) library(stm) library(cleanNLP) library(quanteda) # For nlp package_list <- c("cleanNLP") package_new <- package_list[!(package_list %in% installed.packages()[,"Package"])] if (length(package_new)) { install.packages(package_list) } if (reticulate::py_module_available("cleannlp") == FALSE) { reticulate::py_install("cleannlp") } if (reticulate::py_module_available("spacy") == FALSE) { reticulate::py_install("spacy") } # Some functions for the pre-processing source("/kaggle/input/cdpfunctions/cdp-functions.R") ``` ![Source: Renewable Energy Investment in South Africa](https://www.thirdway.co.za/wp-content/uploads/2017/05/Renewable-Energy-Investment-in-South-Africa.jpg) ```{r message=FALSE, warning=FALSE, include=FALSE} ####################################################### 2020 ####################################################### ###### Data setwd(paste0("../input/cdp-unlocking-climate-solutions/Cities/Cities Responses")) full_cities20 <- read.csv('2020_Full_Cities_Dataset.csv') %>% select('id' = Account.Number, 'org' = Organization, 'cntry' = Country, 'region' = CDP.Region, 'qstn' = Question.Number, 'qstname' = Question.Name, 'coln' = Column.Number, 'colname' = Column.Name, 'rown' = Row.Number, 'rowname' = Row.Name, 'resp' = Response.Answer) # names(full_cities20) subset <- c("3.2", "3.2a","3.3", "4.0", "4.15", "6.2a", "6.2a", "6.2","5.0", "5.0a", "5.0c", "5.0d", "5.0b", "5.2","5.3", "5.4", "5.5", "5.5a", "10.7a","0.5", "0.6", "1.0", '10.1','10.2', "1.0a", "6.0", "6.5", "8.0", "8.0a", "8.1", "8.4", "8.5", "8.5a", "2.0", "2.0b", "2.3a", "14.1", "12.4", "10.7", "13.0", "14.3", "14.4", "14.4a", '2.1', '3.0', '10.4', '12.3', "8.2", '0.4') cities20 <- full_cities20 %>% filter(qstn %in% subset) #length(unique(cities20$id)) ####################################################### General info 2020 ####################################################### ###### Q0.5: Population cities20_ <- cities20 %>% select(id, org, qstn, coln, 'popsize' = resp) %>% filter(qstn == 0.5 & coln == 1) cities20_$popsize <- as.numeric(as.character(cities20_$popsize)) # Fixing incorrect responses cities20_[cities20_$id == "841098",]$popsize <- 338000 cities20_[cities20_$id == "60318",]$popsize <- 494000 cities20_[cities20_$id == "54652",]$popsize <- 676000 cities20_[cities20_$id == "839650",]$popsize <- 124100 cities20_[cities20_$id == "826381",]$popsize <- 4412000 cities20_[cities20_$id == "841003",]$popsize <- 524700 cities20_[cities20_$id == "826211",]$popsize <- 399724 cities20_[cities20_$id == "42388",]$popsize <- 1383432 cities20_[cities20_$id == "845309",]$popsize <- 329675 ###### Q0.6: geography cities20_ <- cities20_ %>% full_join(cities20 %>% select(id, org, qstn, coln, 'landarea' = resp) %>% filter(qstn == 0.6 & coln == 1) %>% select(id, org, landarea), by = c('id', 'org')) ####################################################### Governance 2020 ####################################################### ###### Q1.0a # Selecting relevant columns Q1.0a <- cities20 %>% select(id, org, qstn, coln, resp) %>% filter(qstn == '1.0a' & coln == 1) %>% select(id, org, coln, resp) Q1.0a$resp <- as.character(Q1.0a$resp) # Number of sustainability goals Q1.0a_1 <- count_rows(Q1.0a, 1, sustain_goals) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q1.0a_1), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Risk assessment 2020 ####################################################### ###### Q2.0b # Selecting relevant columns Q2.0b <- cities20 %>% select(id, org, qstn, coln, resp) %>% filter(qstn == '2.0b' & coln %in% c(1,4,7,8)) %>% select(id, org, coln, resp) Q2.0b$resp <- as.character(Q2.0b$resp) ###### Boundary of assessment Q2.0b[Q2.0b$coln == 4,'resp'] <- boundary_question_recode(Q2.0b[Q2.0b$coln == 4,'resp']) #table(as.character(Q2.0b[Q2.0b$coln == 4,'resp'])) # Averaging for all assessments Q2.0b_4 <- mean_for_multiple_responses(Q2.0b, 4, boundary_assess) ###### Vulnerable populations # Recoding Q2.0b[Q2.0b$coln == 7,'resp'] <- ifelse(Q2.0b[Q2.0b$coln == 7,'resp'] == 'Yes', 1, ifelse(Q2.0b[Q2.0b$coln == 7,'resp'] == 'No', 0, ifelse(Q2.0b[Q2.0b$coln == 7,'resp'] == 'Question not applicable', 9999999, NA))) #table(as.character(Q2.0b[Q2.0b$coln == 7,'resp'])) # Averaging by all documents Q2.0b_7 <- mean_for_multiple_responses(Q2.0b, 7, vulnpop_assess) ###### Number of assessment publications and areas/sectors covered in total Q2.0b_1 <- count_rows(Q2.0b, 1, n_assess) ###### Number of areas/sectors covered in total Q2.0b_8 <- count_rows(Q2.0b, 8, areas_assess) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q2.0b_4, Q2.0b_1, Q2.0b_8, Q2.0b_7), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Adaptation targets 2020 ####################################################### ###### Q3.2a: plans that address climate change adaptation # Selecting relevant columns Q3.2a <- cities20 %>% select(id, org, qstn, coln, resp) %>% filter(qstn == '3.2a' & coln %in% c(1,3,4,6,8,9,10)) %>% select(id, org, coln, resp) Q3.2a$resp <- as.character(Q3.2a$resp) ###### Number of documents Q3.2a_1 <- count_rows(Q3.2a, 1, n_doc_adaptation) ##### Number of areas covered by plans Q3.2a_3_categories <- c('Energy', 'Transport (mobility)', 'Building and Infrastructure', 'Industry', 'ICT (Information and Communication Technology)', 'Spatial Planning', 'Agriculture and Forestry', 'Fishery', 'Water', 'Waste', 'Public Health and Safety', 'Business and Financial Service', 'Social Services') Q3.2a_3 <- count_rows(Q3.2a, 3, n_areas_adaptation_plan, select_from = T, categories = Q3.2a_3_categories) ##### Number of climate hazards factored into plan Q3.2a_4_categories <- as.data.frame(table(as.character(Q3.2a[Q3.2a$coln == 4,'resp']))) %>% select(Var1) %>% filter(grepl('>', Var1)) Q3.2a_4_categories <- as.character(Q3.2a_4_categories$Var1) Q3.2a_4 <- count_rows(Q3.2a, 4, n_hazards_plan, select_from = T, categories = Q3.2a_4_categories) ###### Boundary of plan Q3.2a[Q3.2a$coln == 6,'resp'] <- boundary_question_recode(Q3.2a[Q3.2a$coln == 6,'resp']) #table(as.character(Q3.2a[Q3.2a$coln == 6,'resp'])) # Averaging values of a question about boundaries for all plans Q3.2a_6 <- mean_for_multiple_responses(Q3.2a, 6, boundary_adaptation_plan) ##### Stage of implementation # Recoding Q3.2a[Q3.2a$coln == 8,'resp'] <- stage_question_recode(Q3.2a[Q3.2a$coln == 8,'resp']) #table(as.character(Q3.2a[Q3.2a$coln == 8,'resp'])) # checking the recoding results # Summing values of the question about a stage for all plans Q3.2a_8 <- mean_for_multiple_responses(Q3.2a, 8, stage_adaptation_plan) ###### Type of plan # Recoding to ordinal scale Q3.2a[Q3.2a$coln == 9,'resp'] <- ifelse(Q3.2a[Q3.2a$coln == 9,'resp'] == 'Standalone', 1, ifelse(Q3.2a[Q3.2a$coln == 9,'resp'] == 'Integrated mitigation / adaptation', 2, ifelse(Q3.2a[Q3.2a$coln == 9,'resp'] == 'Addressed in general city master plan', 3, ifelse(Q3.2a[Q3.2a$coln == 9,'resp'] == 'Addressed in city sector plan', 4, ifelse(Q3.2a[Q3.2a$coln == 9,'resp'] %in% as.character(str_subset(Q3.2a[Q3.2a$coln == 9,'resp'], "^Other, please specify")), 2, ifelse(Q3.2a[Q3.2a$coln == 9,'resp'] == 'Question not applicable', 9999999, NA)))))) #table(as.character(Q3.2a[Q3.2a$coln == 9,'resp'])) # checking the recoding results # Average type of a plan Q3.2a_9 <- mean_for_multiple_responses(Q3.2a, 9, adaptation_plan_type) ###### Assessment by the local government # Recoding Q3.2a[Q3.2a$coln == 10,'resp'] <- presence_or_absence_recode(Q3.2a[Q3.2a$coln == 10,'resp']) #table(as.character(Q3.2a[Q3.2a$coln == 10,'resp'])) # checking the recoding results # Averaging values of the question about local government's assessment for all plans Q3.2a_10 <- mean_for_multiple_responses(Q3.2a, 10, localgov_assess_adaptation_plan) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q3.2a_1, Q3.2a_3, Q3.2a_4, Q3.2a_6, Q3.2a_8, Q3.2a_9, Q3.2a_10), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Social equity 2020 ####################################################### ###### Q3.0 Selecting hazards that were addressed by actions Q3.0 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '3.0' & coln %in% c(1,2)) %>% select(-qstn) %>% arrange(org, coln, rown) Q3.0_2 <- Q3.0 %>% filter(coln %in% 2) %>% filter(!resp == "Question not applicable" & !resp == "" & !resp == 'No action currently taken') %>% left_join(Q3.0 %>% filter(coln %in% 1), by = c('id', 'org', 'rown')) %>% filter(!resp.y == '') %>% select(id, org, resp.y) %>% group_by(id, org) %>% mutate(resp = resp.y) %>% select(-resp.y) Q3.0_2 <- distinct(Q3.0_2) # Adding an indicator showing that a hazard was addressed Q3.0_2$selected = 1 # Selecting all hazards reported Q2.1 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '2.1' & coln %in% c(1,7)) %>% select(-qstn)%>% arrange(org, coln, rown) Q2.1$resp <- as.character(Q2.1$resp) # Merging the addressed hazards with all hazards reported Q2.1_Q3.0_2 <- Q2.1 %>% left_join(Q3.0_2, by = c('id', 'org', 'resp')) %>% filter(selected == 1, coln == 1) # Selecting reports with vulnerable populations and filtering only non-missing values Q2.1_7 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '2.1' & coln == 7) %>% select(-qstn) %>% arrange(org, coln, rown) %>% filter(!resp == '') # Merging the addressed hazards with their respective vulnerable populations Q2.1_Q3.0_2_Q2.1_7 <- Q2.1_Q3.0_2 %>% left_join(Q2.1_7, by = c('id', 'org', 'rown')) # Defining the number of vulnerable populations addressed df_vuln_addressed <- Q2.1_Q3.0_2_Q2.1_7 %>% group_by(id, org) %>% filter(!is.na(resp.y)) %>% # filtering those who did not report on vulnerable populations for a certain hazard summarise(vuln_addressed = n(), .groups = 'drop') # Defining the number of vulnerable populations in total df_vuln_total <- Q2.1_7 %>% group_by(id, org) %>% summarise(vuln_total = n(), .groups = 'drop') # Addressed per total df_vuln_addressed_per_total <- df_vuln_addressed %>% left_join(df_vuln_total, by = c('id', 'org')) %>% filter(!is.na(vuln_total)) %>% mutate(vuln_addressed_per_total = vuln_addressed/vuln_total) %>% select(id, org, vuln_addressed_per_total) #hist(df_vuln_addressed_per_total$vuln_addressed_per_total) ################################################################## #Proportion of actions addressed in hazards with social impact for vulnerable populations ################################################################## # Selecting all hazards reported Q2.1 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '2.1' & coln %in% c(1,5)) %>% select(-qstn)%>% arrange(org, coln, rown) Q2.1$resp <- as.character(Q2.1$resp) # Merging the addressed hazards with all hazards reported Q2.1_Q3.0_2 <- Q2.1 %>% left_join(Q3.0_2, by = c('id', 'org', 'resp')) %>% filter(selected == 1, coln == 1) # Selecting reports with vulnerable populations and filtering only non-missing values Q2.1_7 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '2.1' & coln == 5) %>% select(-qstn) %>% arrange(org, coln, rown) %>% filter(!resp == '') # Merging the addressed hazards with their respective vulnerable populations Q2.1_Q3.0_2_Q2.1_7 <- Q2.1_Q3.0_2 %>% left_join(Q2.1_7, by = c('id', 'org', 'rown')) # Defining the amount of social_impact addressed df_social_impact_addressed <- Q2.1_Q3.0_2_Q2.1_7 %>% group_by(id, org) %>% filter(!is.na(resp.y)) %>% # filtering those who did not report on social_impact for a certain hazard summarise(social_impact_addressed = n(), .groups = 'drop') # Defining the amount of social_impact in total df_social_impact_total <- Q2.1_7 %>% group_by(id, org) %>% summarise(social_impact_total = n(), .groups = 'drop') # Addressed per total df_social_impact_addressed_per_total <- df_social_impact_addressed %>% left_join(df_social_impact_total, by = c('id', 'org')) %>% filter(!is.na(social_impact_total)) %>% mutate(social_impact_addressed_per_total = social_impact_addressed/social_impact_total) %>% select(id, org, social_impact_addressed_per_total) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, df_vuln_addressed_per_total, df_social_impact_addressed_per_total), by = c('id', 'org'), type='full') ####################################################### Adaptation Actions 2020 ####################################################### ###### Q2.1 and Q3_1: Hazards in actions per hazards reported Q3.0 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '3.0' & coln %in% c(1,2)) %>% select(-qstn) %>% arrange(org, coln, rown) # Number of actions in total Q3.0_1 <- count_rows(Q3.0, 1, n_hazards) # Number of actions in total Q3.0_2 <- Q3.0 %>% filter(coln %in% 2) %>% filter(!resp == "Question not applicable" & !resp == "" & !resp == 'No action currently taken') %>% group_by(id, org) %>% summarise(n_adaptation_actions = n(), .groups = 'drop') # Merging Q3.0_12 <- Q3.0_1 %>% full_join(Q3.0_2, by = c('id', 'org')) # Number of UNIQUE hazards addressed Q3.0_12$n_hazards_addressed <- apply(Q3.0_12[, c('n_hazards', 'n_adaptation_actions')], 1, FUN = min) Q3.0_12$n_hazards_addressed <- ifelse(is.na(Q3.0_12$n_hazards) | is.na(Q3.0_12$n_adaptation_actions) | Q3.0_12$n_hazards == 9999999 | Q3.0_12$n_adaptation_actions == 9999999, 0, Q3.0_12$n_hazards_addressed) Q3.0_12[Q3.0_12 == 9999999 | is.na(Q3.0_12)] <- 0 # Actions per hazard Q3.0_12$n_actions_per_hazard <- Q3.0_12$n_adaptation_actions/Q3.0_12$n_hazards Q3.0_12[is.infinite(Q3.0_12$n_actions_per_hazard) | is.nan(Q3.0_12$n_actions_per_hazard), 'n_actions_per_hazard'] <- NA Q3.0_12 <- Q3.0_12 %>% select(id, org, n_hazards_addressed, n_actions_per_hazard) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q3.0_12), by = c('id', 'org'), type='full') ###### Q3.0 Q3.0 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '3.0' & coln %in% c(2,4,6,7)) %>% select(id, org, coln, rown, resp) %>% arrange(org, coln, rown) Q3.0$resp <- as.character(Q3.0$resp) ###### Number of actions Q3.0_2 <- count_rows(Q3.0, 2, n_adaptation_actions) ###### Status of actions (sum by all actions) Q3.0[Q3.0$coln == 4, 'resp'] <- ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Scoping', 1, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Pre-feasibility study/impact assessment', 2, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Pre-implementation', 3, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Implementation', 4, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Implementation complete but not in operation', 5, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Operation', 6, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == "Monitoring and reporting", 7, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] %in% as.character(str_subset(Q3.0[Q3.0$coln == 4, 'resp'], "^Other, please specify")), 3.5, NA)))))))) #table(as.character(Q3.0[Q3.0$coln == 4, 'resp'])) # Averaging values of a question about status for all actions Q3.0_4 <- mean_for_multiple_responses(Q3.0, 4, adaptation_action_status) ###### Number of co-benefits Q3.0_6 <- count_rows(Q3.0, 6, n_adaptation_cobenefits) ###### Number of Sectors/areas adaptation actions apply to Q3.0_7 <- count_rows(Q3.0, 7, n_adaptation_actions_areas) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q3.0_2, Q3.0_4, Q3.0_6, Q3.0_7), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Adaptation Goals 2020 ####################################################### ###### Q3.3 Q3.3 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '3.3' & coln %in% c(1,2)) %>% select(id, org, coln, rown, resp) %>% arrange(org, coln, rown) Q3.3$resp <- as.character(Q3.3$resp) ###### Number of goals Q3.3_1 <- count_rows(Q3.3, 1, n_goals) ###### Number of climate hazards that adaptation goals address Q3.3_2 <- count_rows(Q3.3, 2, n_hazards_goals) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q3.3_1, Q3.3_2), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Inventory 2020 ####################################################### ###### Q4.0: a city-wide emissions inventory Q4.0 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '4.0') %>% select(id, org, coln, rown, resp) %>% arrange(org, coln, rown) # Recoding into ordinal numeric Q4.0[,'resp'] <- presence_or_absence_recode(Q4.0[,'resp']) Q4.0 <- Q4.0 %>% select(id, org, 'inventory' = resp) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q4.0), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Emission Reduction 2020 ####################################################### ###### Q5.0: % of target achieved and reduction target in metric tonnes Q5.0a <- emission_reduction_targets(cities20, '5.0a', 10, 9) Q5.0b <- emission_reduction_targets(cities20, '5.0b', 8, 7) Q5.0c <- emission_reduction_targets(cities20, '5.0c', 12, 11) Q5.0d <- emission_reduction_targets(cities20, '5.0d', 10, 8) # cities with no target Q5_no_target <- cities20 %>% select(id, org, qstn, resp) %>% filter(qstn == '5.0' & resp == 'No target') %>% mutate(pcnt_target_achieved = 9999999, target_metric_tonnes = 9999999) %>% select(-c(qstn, resp)) %>% distinct(.keep_all = T) #length(unique(Q5_no_target$id)) # Combining all target types and summarising if a city fell into 2 categories Q5.0abcd <- rbind.data.frame(Q5.0a, Q5.0b, Q5.0c, Q5.0d, Q5_no_target) %>% group_by(id, org) %>% summarise_all(list(sum)) ##### Q5.4 Mitigation Targets Q5.4 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '5.4' & coln %in% c(1,4:7,9)) %>% select(-qstn) %>% arrange(org, coln, rown) Q5.4$resp <- as.character(Q5.4$resp) # Number of mitigation actions Q5.4_1 <- count_rows(Q5.4, 1, n_mitigation_actions) ###### Status of actions (sum by all actions) Q5.4[Q5.4$coln == 4, 'resp'] <- ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Scoping', 1, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Pre-feasibility study', 2, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Pre-implementation', 3, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Implementation', 4, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Implementation complete', 5, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Operation', 6, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == "Monitoring and reporting", 7, NA))))))) #table(as.character(Q5.4[Q5.4$coln == 4, 'resp'])) # Averaging values of a question about status for all actions Q5.4_4 <- mean_for_multiple_responses(Q5.4, 4, mitigation_action_status) ###### Number of co-benefits Q5.4_9 <- count_rows(Q5.4, 9, n_mitigation_cobenefits) ###### Estimated emissions reduction: sum for all actions per population size #Q5.4_5 <- sum_for_multiple_responses(Q5.4, 5, estimated_emissions_reduction) ###### Energy savings (MWh): sum for all actions per population size #Q5.4_6 <- sum_for_multiple_responses(Q5.4, 6, energy_savings) ###### Renewable energy production (MWh): sum for all actions per population size #Q5.4_7 <- sum_for_multiple_responses(Q5.4, 7, renewable_energy_production) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q5.0abcd, Q5.4_1, Q5.4_4, Q5.4_9), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ##### Q5.5a Mitigation Plan Q5.5a <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '5.5a' & coln %in% c(1,5,6,8,9)) %>% select(-qstn) %>% arrange(org, coln, rown) Q5.5a$resp <- as.character(Q5.5a$resp) ###### Number of documents Q5.5a_1 <- count_rows(Q5.5a, 1, n_doc_mitigation) ##### Number of areas covered by plan Q5.5a_5 <- count_rows(Q5.5a, 5, n_areas_mitigation_plan) ###### Boundary of a plan Q5.5a[Q5.5a$coln == 6,'resp'] <- boundary_question_recode(Q5.5a[Q5.5a$coln == 6,'resp']) #table(as.character(Q5.5a[Q5.5a$coln == 6,'resp'])) # Averaging values of a question about boundaries for all plans Q5.5a_6 <- mean_for_multiple_responses(Q5.5a, 6, boundary_mitigation_plan) ###### Stage of implementation # Recoding Q5.5a[Q5.5a$coln == 8,'resp'] <- stage_question_recode(Q5.5a[Q5.5a$coln == 8,'resp']) #table(as.character(Q5.5a[Q5.5a$coln == 8,'resp'])) # checking the recoding results # Averaging values of the question about a stage for all plans Q5.5a_8 <- mean_for_multiple_responses(Q5.5a, 8, stage_mitigation_plan) ###### Assessment by the local government # Recoding Q5.5a[Q5.5a$coln == 9,'resp'] <- presence_or_absence_recode(Q5.5a[Q5.5a$coln == 9,'resp']) table(as.character(Q5.5a[Q5.5a$coln == 9,'resp'])) # checking the recoding results # Averaging values of the question about local government's assessment for all plans Q5.5a_9 <- mean_for_multiple_responses(Q5.5a, 9, localgov_assess_mitigation_plan) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q5.5a_1, Q5.5a_5, Q5.5a_6, Q5.5a_8, Q5.5a_9), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Opportunities 2020 ####################################################### ###### Q6.0: Number of opportunities identified as a result of addressing climate change Q6.0 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '6.0' & coln == 1) %>% select(-qstn) %>% arrange(org, coln, rown) Q6.0$resp <- as.character(Q6.0$resp) # Number of opportunities Q6.0_1 <- count_rows(Q6.0, 1, n_opportunities) ###### Q6.2a: Collaboration with business Q6.2a <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '6.2a' & coln %in% c(1,2)) %>% select(-qstn) %>% arrange(org, coln, rown) Q6.2a$resp <- as.character(Q6.2a$resp) # Number of collaboration areas Q6.2a_1 <- count_rows(Q6.2a, 1, n_collab_areas) # Number of collaboration types Q6.2a_2 <- count_rows(Q6.2a, 2, n_collab_types) ###### Q6.5: Collaboration with business Q6.5 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '6.5' & coln %in% c(1,3,4)) %>% select(-qstn) %>% arrange(org, coln, rown) Q6.5$resp <- as.character(Q6.5$resp) # Number of unique project areas Q6.5_1 <- count_rows(Q6.5, 1, n_project_areas) ###### Stage of implementation # Recoding Q6.5[Q6.5$coln == 3,'resp'] <- ifelse(Q6.5[Q6.5$coln == 3,'resp'] == 'Scoping', 1, ifelse(Q6.5[Q6.5$coln == 3,'resp'] == 'Pre-feasibility study/impact assessment', 2, ifelse(Q6.5[Q6.5$coln == 3,'resp'] == 'Project feasibility', 3, ifelse(Q6.5[Q6.5$coln == 3,'resp'] == 'Project structuring', 4, ifelse(Q6.5[Q6.5$coln == 3,'resp'] == 'Transaction preparation', 5, ifelse(Q6.5[Q6.5$coln == 3,'resp'] == 'Implementation', 6, ifelse(Q6.5[Q6.5$coln == 3,'resp'] == "Post-implementation", 7, ifelse(Q6.5[Q6.5$coln == 3, 'resp'] == 'Question not applicable', 9999999, NA)))))))) #table(as.character(Q6.5[Q6.5$coln == 3,'resp'])) # checking the recoding results # Averaging values of the question about a stage for all projects Q6.5_3 <- mean_for_multiple_responses(Q6.5, 3, stage_project) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q6.0_1, Q6.2a_1, Q6.2a_2, Q6.5_1, Q6.5_3), by = c('id', 'org'), type='full') #length(unique(Q6.5_3$id)) ####################################################### Renewable Energy Target 2020 ####################################################### ###### Q8.0a Q8.0a <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '8.0a' & coln %in% c(9,10)) %>% select(-qstn) %>% arrange(org, coln, rown) Q8.0a$resp <- as.character(Q8.0a$resp) Q8.0a[Q8.0a$resp == 'Question not applicable', 'resp'] <- 9999999 # Average percentage of renewable energy / electricity of total energy or electricity in target year Q8.0a_9 <- mean_for_multiple_responses(Q8.0a, 9, percentage_renewable_energy) # Average percentage of targets achieved Q8.0a_10 <- mean_for_multiple_responses(Q8.0a, 10, renewable_energy_achieved) ###### Q8.1 the source mix of electricity consumed Q8.1 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '8.1' & coln %in% 5:9) %>% select(-qstn) %>% arrange(org, coln, rown) %>% group_by(id, org) %>% summarise(renewable_energy_source = ifelse(all(is.na(as.numeric(as.character(resp)))), NA, sum(as.numeric(as.character(resp)), na.rm = T)), .groups = 'drop') ####### Q8.5a city’s energy efficiency targets Q8.5a <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '8.5a' & coln %in% c(7,8)) %>% select(-qstn) %>% arrange(org, coln, rown) Q8.5a$resp <- as.character(Q8.5a$resp) Q8.5a[Q8.5a$resp == 'Question not applicable', 'resp'] <- 9999999 # Average percentage of energy efficiency improvement in target year compared to base year levels Q8.5a_7 <- mean_for_multiple_responses(Q8.5a, 7, percentage_energy_efficiency_improvement) # Average percentage of target achieved Q8.5a_8 <- mean_for_multiple_responses(Q8.5a, 8, energy_efficiency_achieved) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q8.0a_9, Q8.0a_10, Q8.1, Q8.5a_7, Q8.5a_8), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Transport 2020 ####################################################### ###### Q10.1 the mode share of each transport mode Q10.1 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '10.1' & coln %in% c(5,6)) %>% select(-qstn) %>% arrange(org, coln, rown) %>% group_by(id, org) %>% filter(!resp == "Question not applicable") %>% summarise(percentage_walking_cycling = ifelse(all(is.na(as.numeric(as.character(resp)))), NA, sum(as.numeric(as.character(resp)), na.rm = T)), .groups = 'drop') # Adding cities with 'not applicable" responses Q10.1 <- rbind.data.frame(Q10.1, cities20 %>% filter(qstn == '10.1' & coln %in% c(5,6)) %>% filter(resp == "Question not applicable") %>% select(id, org) %>% mutate(percentage_walking_cycling = 9999999)) %>% distinct(.keep_all = T) ###### Q10.4 Total number of vehicle types for private cars (Electric + Hybrid + Plug in hybrid + Hydrogen) Q10.4 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '10.4' & coln == 1 & rown %in% 2:5) %>% select(-qstn) %>% arrange(org, coln, rown)%>% group_by(id, org) %>% filter(!resp == "Question not applicable") %>% summarise(eco_private_cars = ifelse(all(is.na(as.numeric(as.character(resp)))), NA, sum(as.numeric(as.character(resp)), na.rm = T)), .groups = 'drop') # Adding cities with 'not applicable" responses Q10.4 <- rbind.data.frame(Q10.4, cities20 %>% filter(qstn == '10.4' & coln == 1 & rown %in% 2:5) %>% filter(resp == "Question not applicable") %>% select(id, org) %>% mutate(eco_private_cars = 9999999)) %>% distinct(.keep_all = T) ###### Q10.7a: zero emissions zone presence Q10.7 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '10.7') %>% select(-qstn) %>% arrange(org, coln, rown) %>% mutate(zero_emission_zone_presence = ifelse(resp == 'Yes', 1, ifelse(resp == 'No' | resp == 'Do not know', 0, NA))) %>% select(-c(coln, rown, resp)) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q10.1, Q10.4, Q10.7), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Food 2020 ####################################################### ###### Q12.3: presence of policies relating to food consumption Q12.3 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '12.3' & coln == 1) %>% select(-qstn) %>% arrange(org, coln, rown)%>% group_by(id, org) %>% mutate(food_policies_presence = ifelse(resp == 'Yes', 1, ifelse(resp == 'No' | resp == 'Do not know', 0, NA)))%>% select(-c(coln, rown, resp)) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q12.3), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Waste 2020 ####################################################### # A lot of missing values and variance is poor -> drop ###### Q13.0: Amount of solid waste generated (tonnes/year) per population size #Q13.0 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% # filter(qstn == '13.0' & coln == 1) %>% select(-qstn) %>% arrange(org, coln, rown) %>% # mutate(waste_amount = as.numeric(as.character(resp))) %>% # select(-c(coln, rown, resp)) # Merging with the selected variables so far #cities20_ <- plyr::join_all(list(cities20_, Q13.0), # by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ####################################################### Water Security 2020 ####################################################### ###### Q14.3 Q14.3 <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '14.3' & coln %in% 2:3) %>% select(-qstn) %>% arrange(org, coln, rown) Q14.3$resp <- as.character(Q14.3$resp) # Number of unique adaptation actions Q14.3_2 <- count_rows(Q14.3, 2, n_water_adaptation_actions) # Status Q14.3[Q14.3$coln == 3, 'resp'] <- ifelse(Q14.3[Q14.3$coln == 3, 'resp'] == 'Scoping', 1, ifelse(Q14.3[Q14.3$coln == 3, 'resp'] == 'Pre-feasibility study/impact assessment', 2, ifelse(Q14.3[Q14.3$coln == 3, 'resp'] == 'Pre-implementation', 3, ifelse(Q14.3[Q14.3$coln == 3, 'resp'] == 'Implementation', 4, ifelse(Q14.3[Q14.3$coln == 3, 'resp'] == 'Implementation complete but not in operation', 5, ifelse(Q14.3[Q14.3$coln == 3, 'resp'] == 'Operation', 6, ifelse(Q14.3[Q14.3$coln == 3, 'resp'] == "Monitoring and reporting", 7, ifelse(Q14.3[Q14.3$coln == 3, 'resp'] == 'Question not applicable', 9999999, NA)))))))) table(as.character(Q14.3[Q14.3$coln == 3, 'resp'])) # Averaging values for all water adaptation actions Q14.3_3 <- mean_for_multiple_responses(Q14.3, 3, water_adaptation_action_status) length(unique(Q14.3_3$id)) ###### Q14.4a: Amount of solid waste generated (tonnes/year) Q14.4a <- cities20 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '14.4a' & coln %in% c(1,5)) %>% select(-qstn) %>% arrange(org, coln, rown) Q14.4a$resp <- as.character(Q14.4a$resp) # Number of documents Q14.4a_1 <- count_rows(Q14.4a, 1, n_water_management) # Stage of implementation Q14.4a[Q14.4a$coln == 5, 'resp'] <- stage_question_recode(Q14.4a[Q14.4a$coln == 5, 'resp']) table(as.character(Q14.4a[Q14.4a$coln == 5,'resp'])) # checking the recoding results # Averaging values of the question about a stage for all plans Q14.4a_5 <- mean_for_multiple_responses(Q14.4a, 5, stage_water_management_plan) # Merging with the selected variables so far cities20_ <- plyr::join_all(list(cities20_, Q14.3_2, Q14.3_3, Q14.4a_1, Q14.4a_5), by = c('id', 'org'), type='full') #length(unique(cities20_$id)) ``` ```{r message=FALSE, warning=FALSE, include=FALSE} ####################################################### KPIs 2020 ####################################################### # Replacing 9999999 with 0 cities20_[cities20_ == 9999999] <- 0 # Norming by population size cities20_$target_metric_tonnes <- cities20_$target_metric_tonnes/cities20_$popsize #cities20_$estimated_emissions_reduction <- cities20_$estimated_emissions_reduction/cities20_$popsize #cities20_$energy_savings <- cities20_$energy_savings/cities20_$popsize #cities20_$renewable_energy_production <- cities20_$renewable_energy_production/cities20_$popsize #cities20_$waste_amount <- cities20_$waste_amount/cities20_$popsize # Exploring patterns of missing values cities20_items_only <- cities20_[,c(1,2,7:ncol(cities20_))] row_na <- cbind.data.frame(cities20_items_only$org, 'row_na_sum' = rowSums(is.na(cities20_items_only))) %>% arrange(row_na_sum) col_na <- cbind.data.frame(colnames(cities20_items_only), 'col_na_sum' = colSums(is.na(cities20_items_only))) %>% arrange(col_na_sum) rownames(col_na) <- NULL # Removing variables with >200 missing values cities20_items_only_ <- cities20_items_only %>% select(-(as.character(col_na[col_na$col_na_sum > 200, 'colnames(cities20_items_only)']))) # Recalculating missing patterns row_na <- cbind.data.frame(cities20_items_only_$org, 'row_na_sum' = rowSums(is.na(cities20_items_only_))) %>% arrange(row_na_sum) col_na <- cbind.data.frame(colnames(cities20_items_only_), 'col_na_sum' = colSums(is.na(cities20_items_only_))) %>% arrange(col_na_sum) rownames(col_na) <- NULL # Removing cities with >10 missing values cities20_items_only_2 <- cities20_items_only_ %>% filter(! org %in% as.character(row_na[row_na$row_na_sum > 10, 'cities20_items_only_$org'])) # Focusing only on cities with population size > 50000 df_cities_info <- rio::import("../input/df-cities-infoadditional/df_cities_info.additional.xlsx") %>% mutate(hazards_intensity = mean_Index_Current_impact*hazards_n, no2_per_population = no2_sum/population) df_cities_info <- df_cities_info %>% filter(population >= 50000, !is.na(gdp_per_capita)) %>% select(id, org, gdp_per_capita, hazards_intensity, population, mean_temp, sum_precip, no2_per_population) # Merging cities20_items_only_2 <- df_cities_info %>% left_join(cities20_items_only_2 %>% select(-org), by = 'id') #cities20_items_only_2 <- na.omit(cities20_items_only_2) # Imputing NA for 2020 ################################### library(Amelia) # Setting bounds for imputations max <- as.matrix(as.numeric(as.character(apply(cities20_items_only_2, 2 , max, na.rm=TRUE)))) max <- as.matrix(max[3:nrow(max),]) # Running imputation set.seed(123) amelia_res <- amelia(cities20_items_only_2, m = 10, idvars = c('id', 'org'), bounds = cbind(3:53, 0, max)) imputed_data2 <- amelia_res$imputations$imp10 #summary(imputed_data2) cities20_items_only_2 <- imputed_data2 ################################### ``` ```{r message=FALSE, warning=FALSE, include=FALSE} ####################################################### 2019 ####################################################### ###### Data setwd(paste0("../input/cdp-unlocking-climate-solutions/Cities/Cities Responses")) full_cities19 <- read.csv('2019_Full_Cities_Dataset.csv') %>% select('id' = Account.Number, 'org' = Organization, 'cntry' = Country, 'region' = CDP.Region, 'qstn' = Question.Number, 'qstname' = Question.Name, 'coln' = Column.Number, 'colname' = Column.Name, 'rown' = Row.Number, 'rowname' = Row.Name, 'resp' = Response.Answer) #names(full_cities19) subset <- c("3.2", "3.1a", "3.3", "4.0", "4.15", "6.1a", "6.2a", "6.2", "5.0", "5.0a", "5.0c", "5.0d", "5.0b", "5.2", "5.3", "5.4", "5.5", "5.5a", "10.7a", "0.5", "0.6", "1.0", '10.1','10.2', "1.0a", "6.0", "6.2", "8.0", "8.0a", "8.2", "8.4", "8.5", "8.6a", "2.0", "2.0b", "2.3a", "14.1", "10.7", "13.0", "14.4", "14.5a", '2.1', '3.0', '10.5', '12.4', "8.2", '0.4', '2.0') cities19 <- full_cities19 %>% filter(qstn %in% subset) # Fixing 'Not applicable' coding cities19$resp <- as.character(cities19$resp) cities19[cities19$resp %in% c('Not applicable.', 'Not applicable', 'not applicable') & is.na(cities19$resp) == F, 'resp'] <- 'Question not applicable' ####################################################### General info 2019 ####################################################### ###### Q0.5: Population cities19_ <- cities19 %>% select(id, org, qstn, coln, 'popsize' = resp) %>% filter(qstn == 0.5 & coln == 1) cities19_$popsize <- as.numeric(as.character(cities19_$popsize)) # Fixing incorrect responses cities19_[cities19_$id == "841098",]$popsize <- 338000 cities19_[cities19_$id == "60318",]$popsize <- 494000 cities19_[cities19_$id == "54652",]$popsize <- 676000 cities19_[cities19_$id == "839650",]$popsize <- 124100 #cities19_[cities19_$id == "826381",]$popsize <- 4412000 cities19_[cities19_$id == "841003",]$popsize <- 524700 #cities19_[cities19_$id == "826211",]$popsize <- 399724 #cities19_[cities19_$id == "42388",]$popsize <- 1383432 #cities19_[cities19_$id == "845309",]$popsize <- 329675 ###### Q0.6: geography cities19_ <- cities19_ %>% full_join(cities19 %>% select(id, org, qstn, coln, 'landarea' = resp) %>% filter(qstn == 0.6 & coln == 1) %>% select(id, org, landarea), by = c('id', 'org')) ####################################################### Governance 2019 ####################################################### ###### Q1.0a # Selecting relevant columns Q1.0a <- cities19 %>% select(id, org, qstn, coln, resp) %>% filter(qstn == '1.0a' & coln == 1) %>% select(id, org, coln, resp) Q1.0a$resp <- as.character(Q1.0a$resp) # Number of sustainability goals Q1.0a_1 <- count_rows(Q1.0a, 1, sustain_goals, year = 2019) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q1.0a_1), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Risk assessment 2019 ####################################################### ###### Q2.0b # Selecting relevant columns Q2.0b <- cities19 %>% select(id, org, qstn, coln, resp) %>% filter(qstn == '2.0b' & coln %in% c(1,4,6,8)) %>% select(id, org, coln, resp) Q2.0b$resp <- as.character(Q2.0b$resp) ###### Boundary of assessment Q2.0b[Q2.0b$coln == 4,'resp'] <- boundary_question_recode(Q2.0b[Q2.0b$coln == 4,'resp']) #table(as.character(Q2.0b[Q2.0b$coln == 4,'resp'])) # Averaging for all assessments Q2.0b_4 <- mean_for_multiple_responses(Q2.0b, 4, boundary_assess, year = 2019) ###### Vulnerable populations # Recoding Q2.0b[Q2.0b$coln == 8,'resp'] <- ifelse(Q2.0b[Q2.0b$coln == 8,'resp'] == 'Yes', 1, ifelse(Q2.0b[Q2.0b$coln == 8,'resp'] == 'No', 0, ifelse(Q2.0b[Q2.0b$coln == 8,'resp'] == 'Question not applicable', 9999999, NA))) #table(as.character(Q2.0b[Q2.0b$coln == 8,'resp'])) # Averaging by all documents Q2.0b_8 <- mean_for_multiple_responses(Q2.0b, 8, vulnpop_assess, year = 2019) ###### Number of assessment publications and areas/sectors covered in total Q2.0b_1 <- count_rows(Q2.0b, 1, n_assess, year = 2019) ###### Number of areas/sectors covered in total Q2.0b_6 <- count_rows(Q2.0b, 6, areas_assess, year = 2019) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q2.0b_4, Q2.0b_1, Q2.0b_8, Q2.0b_6), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Adaptation targets 2019 ####################################################### ###### Q3.1a: plans that address climate change adaptation # Selecting relevant columns Q3.1a <- cities19 %>% select(id, org, qstn, coln, resp) %>% filter(qstn == '3.1a' & coln %in% c(1,2,4,6,7,8)) %>% select(id, org, coln, resp) Q3.1a$resp <- as.character(Q3.1a$resp) ###### Number of documents Q3.1a_1 <- count_rows(Q3.1a, 1, n_doc_adaptation, year = 2019) ##### Number of areas covered by plans Q3.1a_2_categories <- c('Energy', 'Transport (mobility)', 'Building and Infrastructure', 'Industry', 'ICT (Information and Communication Technology)', 'Spatial Planning', 'Agriculture and Forestry', 'Fishery', 'Water', 'Waste', 'Public Health and Safety', 'Business and Financial Service', 'Social Services') Q3.1a_2 <- count_rows(Q3.1a, 2, n_areas_adaptation_plan, select_from = T, categories = Q3.1a_2_categories, year = 2019) ##### Number of climate hazards factored into plan #Q3.1a_4_categories <- as.data.frame(table(as.character(Q3.1a[Q3.1a$coln == 4,'resp']))) %>% select(Var1) %>% # filter(grepl('>', Var1)) #Q3.1a_4_categories <- as.character(Q3.1a_4_categories$Var1) #Q3.1a_4 <- count_rows(Q3.1a, 4, n_hazards_plan, select_from = T, categories = Q3.1a_4_categories) ###### Boundary of plan Q3.1a[Q3.1a$coln == 4,'resp'] <- boundary_question_recode(Q3.1a[Q3.1a$coln == 4,'resp']) #table(as.character(Q3.1a[Q3.1a$coln == 4,'resp'])) # Averaging values of a question about boundaries for all plans Q3.1a_4 <- mean_for_multiple_responses(Q3.1a, 4, boundary_adaptation_plan, year = 2019) ##### Stage of implementation # Recoding Q3.1a[Q3.1a$coln == 6,'resp'] <- stage_question_recode(Q3.1a[Q3.1a$coln == 6,'resp']) table(as.character(Q3.1a[Q3.1a$coln == 6,'resp'])) # checking the recoding results # Averaging values of the question about a stage for all plans Q3.1a_6 <- mean_for_multiple_responses(Q3.1a, 6, stage_adaptation_plan, year = 2019) ###### Type of plan # Recoding to ordinal scale Q3.1a[Q3.1a$coln == 7,'resp'] <- ifelse(Q3.1a[Q3.1a$coln == 7,'resp'] == 'Standalone', 1, ifelse(Q3.1a[Q3.1a$coln == 7,'resp'] == 'Integrated mitigation / adaptation', 2, ifelse(Q3.1a[Q3.1a$coln == 7,'resp'] == 'Addressed in general city master plan', 3, ifelse(Q3.1a[Q3.1a$coln == 7,'resp'] == 'Addressed in city sector plan', 4, ifelse(Q3.1a[Q3.1a$coln == 7,'resp'] %in% as.character(str_subset(Q3.1a[Q3.1a$coln == 7,'resp'], "^Other, please specify")), 2, ifelse(Q3.1a[Q3.1a$coln == 7,'resp'] == 'Question not applicable', 9999999, NA)))))) #table(as.character(Q3.1a[Q3.1a$coln == 7,'resp'])) # checking the recoding results # Average type of a plan Q3.1a_7 <- mean_for_multiple_responses(Q3.1a, 7, adaptation_plan_type, year = 2019) ###### Assessment by the local government # Recoding Q3.1a[Q3.1a$coln == 8,'resp'] <- presence_or_absence_recode(Q3.1a[Q3.1a$coln == 8,'resp']) #table(as.character(Q3.1a[Q3.1a$coln == 8,'resp'])) # checking the recoding results # Averaging values of the question about local government's assessment for all plans Q3.1a_8 <- mean_for_multiple_responses(Q3.1a, 8, localgov_assess_adaptation_plan, year = 2019) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q3.1a_1, Q3.1a_2, Q3.1a_4, Q3.1a_6, Q3.1a_7, Q3.1a_8), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Social equity 2019 ####################################################### ###### Q3.0 Selecting hazards that were addressed by actions Q3.0 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '3.0' & coln %in% c(1,2)) %>% select(-qstn) %>% arrange(org, coln, rown) Q3.0_2 <- Q3.0 %>% filter(coln %in% 2) %>% filter(!resp == "Question not applicable" & !resp == "" & !resp == 'No action currently taken') %>% left_join(Q3.0 %>% filter(coln %in% 1), by = c('id', 'org', 'rown')) %>% filter(!resp.y == '') %>% select(id, org, resp.y) %>% group_by(id, org) %>% mutate(resp = resp.y) %>% select(-resp.y) Q3.0_2 <- distinct(Q3.0_2) # Adding an indicator showing that a hazard was addressed Q3.0_2$selected = 1 # Selecting all hazards reported Q2.1 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '2.1' & coln %in% c(1,10)) %>% select(-qstn)%>% arrange(org, coln, rown) Q2.1$resp <- as.character(Q2.1$resp) # Merging the addressed hazards with all hazards reported Q2.1_Q3.0_2 <- Q2.1 %>% left_join(Q3.0_2, by = c('id', 'org', 'resp')) %>% filter(selected == 1, coln == 1) # Selecting reports with vulnerable populations and filtering only non-missing values Q2.1_7 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '2.1' & coln == 10) %>% select(-qstn) %>% arrange(org, coln, rown) %>% filter(!resp == '') # Merging the addressed hazards with their respective vulnerable populations Q2.1_Q3.0_2_Q2.1_7 <- Q2.1_Q3.0_2 %>% left_join(Q2.1_7, by = c('id', 'org', 'rown')) # Defining the number of vulnerable populations addressed df_vuln_addressed <- Q2.1_Q3.0_2_Q2.1_7 %>% group_by(id, org) %>% filter(!is.na(resp.y)) %>% # filtering those who did not report on vulnerable populations for a certain hazard summarise(vuln_addressed = n(), .groups = 'drop') # Defining the number of vulnerable populations in total df_vuln_total <- Q2.1_7 %>% group_by(id, org) %>% summarise(vuln_total = n(), .groups = 'drop') # Addressed per total df_vuln_addressed_per_total <- df_vuln_addressed %>% left_join(df_vuln_total, by = c('id', 'org')) %>% filter(!is.na(vuln_total)) %>% mutate(vuln_addressed_per_total = vuln_addressed/vuln_total) %>% select(id, org, vuln_addressed_per_total) ################################################################## #Proportion of actions addressed in hazards with social impact for vulnerable populations ################################################################## # Selecting all hazards reported Q2.1 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '2.1' & coln %in% c(1,5)) %>% select(-qstn)%>% arrange(org, coln, rown) Q2.1$resp <- as.character(Q2.1$resp) # Merging the addressed hazards with all hazards reported Q2.1_Q3.0_2 <- Q2.1 %>% left_join(Q3.0_2, by = c('id', 'org', 'resp')) %>% filter(selected == 1, coln == 1) # Selecting reports with vulnerable populations and filtering only non-missing values Q2.1_7 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '2.1' & coln == 5) %>% select(-qstn) %>% arrange(org, coln, rown) %>% filter(!resp == '') # Merging the addressed hazards with their respective vulnerable populations Q2.1_Q3.0_2_Q2.1_7 <- Q2.1_Q3.0_2 %>% left_join(Q2.1_7, by = c('id', 'org', 'rown')) # Defining the amount of social_impact addressed df_social_impact_addressed <- Q2.1_Q3.0_2_Q2.1_7 %>% group_by(id, org) %>% filter(!is.na(resp.y)) %>% # filtering those who did not report on social_impact for a certain hazard summarise(social_impact_addressed = n(), .groups = 'drop') # Defining the amount of social_impact in total df_social_impact_total <- Q2.1_7 %>% group_by(id, org) %>% summarise(social_impact_total = n(), .groups = 'drop') # Addressed per total df_social_impact_addressed_per_total <- df_social_impact_addressed %>% left_join(df_social_impact_total, by = c('id', 'org')) %>% filter(!is.na(social_impact_total)) %>% mutate(social_impact_addressed_per_total = social_impact_addressed/social_impact_total) %>% select(id, org, social_impact_addressed_per_total) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, df_vuln_addressed_per_total, df_social_impact_addressed_per_total), by = c('id', 'org'), type='full') ####################################################### Adaptation Actions 2019 ####################################################### ###### Q2.1 and Q3_1: Hazards in actions per hazards reported Q3.0 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '3.0' & coln %in% c(1,2)) %>% select(-qstn) %>% arrange(org, coln, rown) # Number of hazards in total Q3.0_1 <- count_rows(Q3.0, 1, n_hazards, year = 2019) # Number of actions in total Q3.0_2 <- Q3.0 %>% filter(coln %in% 2) %>% filter(!resp == "Question not applicable" & !resp == "" & !resp == 'No action currently taken') %>% group_by(id, org) %>% summarise(n_adaptation_actions = n(), .groups = 'drop') # Merging Q3.0_12 <- Q3.0_1 %>% full_join(Q3.0_2, by = c('id', 'org')) # Number of UNIQUE hazards addressed Q3.0_12$n_hazards_addressed <- apply(Q3.0_12[, c('n_hazards', 'n_adaptation_actions')], 1, FUN = min) Q3.0_12$n_hazards_addressed <- ifelse(is.na(Q3.0_12$n_hazards) | is.na(Q3.0_12$n_adaptation_actions) | Q3.0_12$n_hazards == 9999999 | Q3.0_12$n_adaptation_actions == 9999999, 0, Q3.0_12$n_hazards_addressed) Q3.0_12[Q3.0_12 == 9999999 | is.na(Q3.0_12)] <- 0 # Actions per hazard Q3.0_12$n_actions_per_hazard <- Q3.0_12$n_adaptation_actions/Q3.0_12$n_hazards Q3.0_12[is.infinite(Q3.0_12$n_actions_per_hazard) | is.nan(Q3.0_12$n_actions_per_hazard), 'n_actions_per_hazard'] <- NA Q3.0_12 <- Q3.0_12 %>% select(id, org, n_hazards_addressed, n_actions_per_hazard) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q3.0_12), by = c('id', 'org'), type='full') ###### Q3.0 Q3.0 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '3.0' & coln %in% c(2,4,5,7)) %>% select(id, org, coln, rown, resp) %>% arrange(org, coln, rown) Q3.0$resp <- as.character(Q3.0$resp) ###### Number of actions Q3.0_2 <- count_rows(Q3.0, 2, n_adaptation_actions, year = 2019) ###### Status of actions (sum by all actions) Q3.0[Q3.0$coln == 4, 'resp'] <- ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Scoping', 1, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Pre-feasibility study/impact assessment', 2, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Pre-implementation', 3, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Implementation', 4, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Implementation complete but not in operation', 5, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == 'Operation', 6, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] == "Monitoring and reporting", 7, ifelse(Q3.0[Q3.0$coln == 4, 'resp'] %in% as.character(str_subset(Q3.0[Q3.0$coln == 4, 'resp'], "^Other, please specify")), 3.5, NA)))))))) #table(as.character(Q3.0[Q3.0$coln == 4, 'resp'])) # Averaging values of a question about status for all actions Q3.0_4 <- mean_for_multiple_responses(Q3.0, 4, adaptation_action_status, year = 2019) ###### Number of co-benefits Q3.0_5 <- count_rows(Q3.0, 5, n_adaptation_cobenefits, year = 2019) ###### Number of Sectors/areas adaptation actions apply to #Q3.0_7 <- count_rows(Q3.0, 7, n_adaptation_actions_areas) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q3.0_2, Q3.0_4, Q3.0_5), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Adaptation Goals 2019 ####################################################### ###### Q3.3 Q3.2 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '3.2' & coln == 1) %>% select(id, org, coln, rown, resp) %>% arrange(org, coln, rown) Q3.2$resp <- as.character(Q3.2$resp) ###### Number of goals Q3.2_1 <- count_rows(Q3.2, 1, n_goals, year = 2019) ###### Number of climate hazards that adaptation goals address #Q3.2_2 <- count_rows(Q3.2, 2, n_hazards_goals) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q3.2_1), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Inventory 2019 ####################################################### ###### Q4.0: a city-wide emissions inventory Q4.0 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '4.0') %>% select(id, org, coln, rown, resp) %>% arrange(org, coln, rown) # Recoding into ordinal numeric Q4.0[,'resp'] <- presence_or_absence_recode(Q4.0[,'resp']) Q4.0 <- Q4.0 %>% select(id, org, 'inventory' = resp) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q4.0), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Emission Reduction 2019####################################################### ###### Q5.0: % of target achieved and reduction target in metric tonnes Q5.0a <- emission_reduction_targets(cities19, '5.0a', 10, 9) Q5.0b <- emission_reduction_targets(cities19, '5.0b', 8, 7) Q5.0c <- emission_reduction_targets(cities19, '5.0c', 12, 11) Q5.0d <- emission_reduction_targets(cities19, '5.0d', 10, 8) # cities with no target Q5_no_target <- cities19 %>% select(id, org, qstn, resp) %>% filter(qstn == '5.0' & resp == 'No target') %>% mutate(pcnt_target_achieved = 9999999, target_metric_tonnes = 9999999) %>% select(-c(qstn, resp)) %>% distinct(.keep_all = T) #length(unique(Q5_no_target$id)) # Combining all target types and summarising if a city fell into 2 categories Q5.0abcd <- rbind.data.frame(Q5.0a, Q5.0b, Q5.0c, Q5.0d, Q5_no_target) %>% group_by(id, org) %>% summarise_all(list(sum)) ##### Q5.4 Mitigation Targets Q5.4 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '5.4' & coln %in% c(1,4:7,9)) %>% select(-qstn) %>% arrange(org, coln, rown) Q5.4$resp <- as.character(Q5.4$resp) # Number of mitigation actions Q5.4_1 <- count_rows(Q5.4, 1, n_mitigation_actions, year = 2019) ###### Status of actions (sum by all actions) Q5.4[Q5.4$coln == 4, 'resp'] <- ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Scoping', 1, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Pre-feasibility study', 2, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Pre-implementation', 3, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Implementation', 4, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Implementation complete', 5, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == 'Operation', 6, ifelse(Q5.4[Q5.4$coln == 4, 'resp'] == "Monitoring and reporting", 7, NA))))))) #table(as.character(Q5.4[Q5.4$coln == 4, 'resp'])) # Averaging values of a question about status for all actions Q5.4_4 <- mean_for_multiple_responses(Q5.4, 4, mitigation_action_status, year = 2019) ###### Number of co-benefits Q5.4_9 <- count_rows(Q5.4, 9, n_mitigation_cobenefits, year = 2019) ###### Estimated emissions reduction: sum for all actions #Q5.4_5 <- sum_for_multiple_responses(Q5.4, 5, estimated_emissions_reduction, year = 2019) ###### Energy savings (MWh): sum for all actions #Q5.4_6 <- sum_for_multiple_responses(Q5.4, 6, energy_savings, year = 2019) ###### Renewable energy production (MWh): sum for all actions #Q5.4_7 <- sum_for_multiple_responses(Q5.4, 7, renewable_energy_production, year = 2019) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q5.0abcd, Q5.4_1, Q5.4_4, Q5.4_9), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ##### Q5.5a Mitigation Plan Q5.5a <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '5.5a' & coln %in% c(1,5,7,8)) %>% select(-qstn) %>% arrange(org, coln, rown) Q5.5a$resp <- as.character(Q5.5a$resp) ###### Number of documents Q5.5a_1 <- count_rows(Q5.5a, 1, n_doc_mitigation, year = 2019) ##### Number of areas covered by plan #Q5.5a_5 <- count_rows(Q5.5a, 5, n_areas_mitigation_plan) ###### Boundary of a plan Q5.5a[Q5.5a$coln == 5,'resp'] <- boundary_question_recode(Q5.5a[Q5.5a$coln == 5,'resp']) table(as.character(Q5.5a[Q5.5a$coln == 5,'resp'])) # Averaging values of a question about boundaries for all plans Q5.5a_5 <- mean_for_multiple_responses(Q5.5a, 5, boundary_mitigation_plan, year = 2019) ###### Stage of implementation # Recoding Q5.5a[Q5.5a$coln == 7,'resp'] <- stage_question_recode(Q5.5a[Q5.5a$coln == 7,'resp']) #table(as.character(Q5.5a[Q5.5a$coln == 7,'resp'])) # checking the recoding results # Averaging values of the question about a stage for all plans Q5.5a_7 <- mean_for_multiple_responses(Q5.5a, 7, stage_mitigation_plan, year = 2019) ###### Assessment by the local government # Recoding Q5.5a[Q5.5a$coln == 8,'resp'] <- presence_or_absence_recode(Q5.5a[Q5.5a$coln == 8,'resp']) #table(as.character(Q5.5a[Q5.5a$coln == 8,'resp'])) # checking the recoding results # Averaging values of the question about local government's assessment for all plans Q5.5a_8 <- mean_for_multiple_responses(Q5.5a, 8, localgov_assess_mitigation_plan, year = 2019) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q5.5a_1, Q5.5a_5, Q5.5a_7, Q5.5a_8), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Opportunities 2019 ####################################################### ###### Q6.0: Number of opportunities identified as a result of addressing climate change Q6.0 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '6.0' & coln == 1) %>% select(-qstn) %>% arrange(org, coln, rown) Q6.0$resp <- as.character(Q6.0$resp) # Number of opportunities Q6.0_1 <- count_rows(Q6.0, 1, n_opportunities, year = 2019) ###### Q6.1a: Collaboration with business Q6.1a <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '6.1a' & coln == 1) %>% select(-qstn) %>% arrange(org, coln, rown) Q6.1a$resp <- as.character(Q6.1a$resp) # Number of collaboration areas Q6.1a_1 <- count_rows(Q6.1a, 1, n_collab_areas, year = 2019) # Number of collaboration types #Q6.1a_2 <- count_rows(Q6.2a, 2, n_collab_types) ###### Q6.2: Collaboration with business Q6.2 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '6.2' & coln %in% c(1,3,4)) %>% select(-qstn) %>% arrange(org, coln, rown) Q6.2$resp <- as.character(Q6.2$resp) # Number of unique project areas Q6.2_1 <- count_rows(Q6.2, 1, n_project_areas, year = 2019) ###### Stage of implementation # Recoding Q6.2[Q6.2$coln == 3,'resp'] <- ifelse(Q6.2[Q6.2$coln == 3,'resp'] == 'Scoping', 1, ifelse(Q6.2[Q6.2$coln == 3,'resp'] == 'Pre-feasibility study/impact assessment', 2, ifelse(Q6.2[Q6.2$coln == 3,'resp'] == 'Project feasibility', 3, ifelse(Q6.2[Q6.2$coln == 3,'resp'] == 'Project structuring', 4, ifelse(Q6.2[Q6.2$coln == 3,'resp'] == 'Transaction preparation', 5, ifelse(Q6.2[Q6.2$coln == 3,'resp'] == 'Implementation', 6, ifelse(Q6.2[Q6.2$coln == 3,'resp'] == "Post-implementation", 7, ifelse(Q6.2[Q6.2$coln == 3, 'resp'] == 'Question not applicable', 9999999, NA)))))))) #table(as.character(Q6.2[Q6.2$coln == 3,'resp'])) # checking the recoding results # Averaging values of the question about a stage for all projects Q6.2_3 <- mean_for_multiple_responses(Q6.2, 3, stage_project, year = 2019) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q6.0_1, Q6.1a_1, Q6.2_1, Q6.2_3), by = c('id', 'org'), type='full') #length(unique(Q6.5_3$id)) ####################################################### Renewable Energy Target 2019 ####################################################### ###### Q8.0a Q8.0a <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '8.0a' & coln %in% c(8,9)) %>% select(-qstn) %>% arrange(org, coln, rown) Q8.0a$resp <- as.character(Q8.0a$resp) Q8.0a[Q8.0a$resp == 'Question not applicable', 'resp'] <- 9999999 # Average percentage of renewable energy / electricity of total energy or electricity in target year Q8.0a_8 <- mean_for_multiple_responses(Q8.0a, 8, percentage_renewable_energy, year = 2019) # Average percentage of targets achieved Q8.0a_9 <- mean_for_multiple_responses(Q8.0a, 9, renewable_energy_achieved, year = 2019) ###### Q8.2 the source mix of electricity consumed Q8.2 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '8.2' & coln %in% 5:9) %>% select(-qstn) %>% arrange(org, coln, rown) %>% group_by(id, org) %>% summarise(renewable_energy_source = ifelse(all(is.na(as.numeric(as.character(resp)))), NA, sum(as.numeric(as.character(resp)), na.rm = T)), .groups = 'drop') ####### Q8.6a city’s energy efficiency targets Q8.6a <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '8.6a' & coln %in% c(7,8)) %>% select(-qstn) %>% arrange(org, coln, rown) Q8.6a$resp <- as.character(Q8.6a$resp) Q8.6a[Q8.6a$resp == 'Question not applicable', 'resp'] <- 9999999 # Average percentage of energy efficiency improvement in target year compared to base year levels Q8.6a_7 <- mean_for_multiple_responses(Q8.6a, 7, percentage_energy_efficiency_improvement, year = 2019) # Average percentage of target achieved Q8.6a_8 <- mean_for_multiple_responses(Q8.6a, 8, energy_efficiency_achieved, year = 2019) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q8.0a_8, Q8.0a_9, Q8.2, Q8.6a_7, Q8.6a_8), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Transport 2019 ####################################################### ###### Q10.1 the mode share of each transport mode Q10.1 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '10.1' & coln %in% c(5,6)) %>% select(-qstn) %>% arrange(org, coln, rown) %>% group_by(id, org) %>% filter(!resp == "Question not applicable") %>% summarise(percentage_walking_cycling = ifelse(all(is.na(as.numeric(as.character(resp)))), NA, sum(as.numeric(as.character(resp)), na.rm = T)), .groups = 'drop') # Adding cities with 'not applicable" responses na_set_Q10.1 <- cities19 %>% filter(qstn == '10.1' & coln %in% c(5,6)) %>% filter(is.na(resp) | resp == "") %>% select(org) Q10.1 <- rbind.data.frame(Q10.1, full_cities19 %>% filter(!org %in% unique(Q10.1$org) & !org %in% na_set_Q10.1$org) %>% select(id, org) %>% mutate(percentage_walking_cycling = 9999999)) %>% distinct(.keep_all = T) ###### Q10.4 Total number of vehicle types for private cars (Electric + Hybrid + Plug in hybrid + Hydrogen) Q10.5 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '10.5' & coln == 1 & rown %in% 2:5) %>% select(-qstn) %>% arrange(org, coln, rown)%>% group_by(id, org) %>% filter(!resp == "Question not applicable") %>% summarise(eco_private_cars = ifelse(all(is.na(as.numeric(as.character(resp)))), NA, sum(as.numeric(as.character(resp)), na.rm = T)), .groups = 'drop') # Adding cities with 'not applicable" responses na_set_Q10.5 <- cities19 %>% filter(qstn == 'Q10.5' & coln == 1 & rown %in% 2:5) %>% filter(is.na(resp) | resp == "") %>% select(org) Q10.5 <- rbind.data.frame(Q10.5, full_cities19 %>% filter(!org %in% unique(Q10.5$org) & !org %in% na_set_Q10.5$org) %>% select(id, org) %>% mutate(eco_private_cars = 9999999)) %>% distinct(.keep_all = T) ###### Q10.7a: zero emissions zone presence Q10.7 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '10.7') %>% select(-qstn) %>% arrange(org, coln, rown) %>% mutate(zero_emission_zone_presence = ifelse(resp == 'Yes', 1, ifelse(resp == 'No' | resp == 'Do not know', 0, NA))) %>% select(-c(coln, rown, resp)) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q10.1, Q10.5, Q10.7), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Food 2019 ####################################################### ###### Q12.4: presence of policies relating to food consumption Q12.4 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '12.4' & coln == 1) %>% select(-qstn) %>% arrange(org, coln, rown)%>% group_by(id, org) %>% mutate(food_policies_presence = ifelse(resp == 'Yes', 1, ifelse(resp == 'No' | resp == 'Do not know', 0, NA)))%>% select(-c(coln, rown, resp)) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q12.4), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Waste 2019 ####################################################### ###### Q13.0: Amount of solid waste generated (tonnes/year) per population size #Q13.0 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% # filter(qstn == '13.0' & coln == 1 & rown == 1) %>% select(-qstn) %>% arrange(org, coln, rown) %>% # mutate(waste_amount = as.numeric(as.character(resp))) %>% # select(-c(coln, rown, resp)) # Merging with the selected variables so far #cities19_ <- plyr::join_all(list(cities19_, Q13.0), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ####################################################### Water Security 2019 ####################################################### ###### Q14.3 Q14.4 <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '14.4' & coln %in% 2:3) %>% select(-qstn) %>% arrange(org, coln, rown) Q14.4$resp <- as.character(Q14.4$resp) # Number of unique adaptation actions Q14.4_2 <- count_rows(Q14.4, 2, n_water_adaptation_actions, year = 2019) # Status Q14.4[Q14.4$coln == 3, 'resp'] <- ifelse(Q14.4[Q14.4$coln == 3, 'resp'] == 'Scoping', 1, ifelse(Q14.4[Q14.4$coln == 3, 'resp'] == 'Pre-feasibility study/impact assessment', 2, ifelse(Q14.4[Q14.4$coln == 3, 'resp'] == 'Pre-implementation', 3, ifelse(Q14.4[Q14.4$coln == 3, 'resp'] == 'Implementation', 4, ifelse(Q14.4[Q14.4$coln == 3, 'resp'] == 'Implementation complete but not in operation', 5, ifelse(Q14.4[Q14.4$coln == 3, 'resp'] == 'Operation', 6, ifelse(Q14.4[Q14.4$coln == 3, 'resp'] == "Monitoring and reporting", 7, ifelse(Q14.4[Q14.4$coln == 3, 'resp'] == 'Question not applicable', 9999999, NA)))))))) #table(as.character(Q14.4[Q14.4$coln == 3, 'resp'])) # Averaging values for all water adaptation actions Q14.4_3 <- mean_for_multiple_responses(Q14.4, 3, water_adaptation_action_status, year = 2019) ###### Q14.5a: Q14.5a <- cities19 %>% select(id, org, qstn, coln, rown, resp) %>% filter(qstn == '14.5a' & coln %in% c(1,5)) %>% select(-qstn) %>% arrange(org, coln, rown) Q14.5a$resp <- as.character(Q14.5a$resp) # Number of documents Q14.5a_1 <- count_rows(Q14.5a, 1, n_water_management, year = 2019) # Stage of implementation Q14.5a[Q14.5a$coln == 5, 'resp'] <- stage_question_recode(Q14.5a[Q14.5a$coln == 5, 'resp']) #table(as.character(Q14.5a[Q14.5a$coln == 5,'resp'])) # checking the recoding results # Averaging values of the question about a stage for all plans Q14.5a_5 <- mean_for_multiple_responses(Q14.5a, 5, stage_water_management_plan, year = 2019) # Merging with the selected variables so far cities19_ <- plyr::join_all(list(cities19_, Q14.4_2, Q14.4_3, Q14.5a_1, Q14.5a_5), by = c('id', 'org'), type='full') #length(unique(cities19_$id)) ``` ```{r message=FALSE, warning=FALSE, include=FALSE} ##################################################### KPIs 2020 2019 ####################################################### # Replacing 9999999 with 0 cities19_[cities19_ == 9999999] <- 0 # Filtering if popsize is missing (meaning these cities are not present in 2020 dataset) cities19_ <- cities19_ %>% filter(!is.na(popsize)) # Selecting only cities present in both 2020 and 2019 cities19_items_only <- cities19_ %>% filter(id %in% cities20_items_only_2$id) %>% select(-c(qstn, coln, popsize, landarea)) %>% select(-c(n_goals, # filtering items with low response rates in 2020 n_mitigation_cobenefits, stage_project, eco_private_cars)) cities20_items_only_2 <- cities20_items_only_2 %>% filter(id %in% cities19_items_only$id) # Adding contextual variables to 2019, the same as in 2020 cities19_items_only <- cities19_items_only %>% arrange(id) cities20_items_only_2 <- cities20_items_only_2 %>% arrange(id) cities19_items_only[, c('gdp_per_capita', "hazards_intensity", "population", "mean_temp", "sum_precip", "no2_per_population")] <- cities20_items_only_2[, c('gdp_per_capita', "hazards_intensity", "population", "mean_temp", "sum_precip", "no2_per_population")] # Norming by population size a range of items cities19_items_only$target_metric_tonnes <- cities19_items_only$target_metric_tonnes/cities19_items_only$population #cities19_items_only$estimated_emissions_reduction <- cities19_items_only$estimated_emissions_reduction/cities19_items_only$population #cities19_items_only$energy_savings <- cities19_items_only$energy_savings/cities19_items_only$population #cities19_items_only$renewable_energy_production <- cities19_items_only$renewable_energy_production/cities19_items_only$population #cities19_items_only$waste_amount <- cities19_items_only$waste_amount/cities19_items_only$population # Remove if NA >10 # Recalculating missing patterns #row_na <- cbind.data.frame(cities19_items_only$org, 'row_na_sum' = rowSums(is.na(cities19_items_only))) %>% arrange(row_na_sum) # Removing cities with >10 missing values #cities19_items_only <- cities19_items_only %>% filter(! org %in% as.character(row_na[row_na$row_na_sum > 10, 'cities19_items_only$org'])) # Imputing NA for 2019 ################################### # Setting bounds for imputations max <- as.matrix(as.numeric(as.character(apply(cities19_items_only, 2 , max, na.rm=TRUE)))) max <- as.matrix(max[3:nrow(max),]) # Running imputation set.seed(123) amelia_res <- amelia(cities19_items_only, m = 10, idvars = c('id', 'org'), bounds = cbind(3:50, 0, max)) data_imputed <- amelia_res$imputations$imp9 #summary(data_imputed) cities19_items_only <- data_imputed ################################### # Adding Year ID cities19_items_only$year <- 2019 cities20_items_only_2$year <- 2020 # Row bind all data to obtain a full variance for further scaling cities19_20 <- plyr::rbind.fill(cities20_items_only_2, cities19_items_only) # Scaling all items based on empirical data: max = max for 2 years, min = min for 2 years under focus # This is needed to ensure unification of a scale across years # 2019 cities19_norm <- cities19_20[cities19_20$year == 2019,] for (i in colnames(cities19_norm)[3:ncol(cities19_norm)]){ for (j in 1:nrow(cities19_norm)){ cities19_norm[j,i] <- (cities19_norm[j,i] - min(cities19_20[,i], na.rm = T)) / (max(cities19_20[,i], na.rm = T) - min(cities19_20[,i], na.rm = T)) } } # 2020 cities20_norm <- cities19_20[cities19_20$year == 2020,] for (i in colnames(cities20_norm)[3:ncol(cities20_norm)]){ for (j in 1:nrow(cities20_norm)){ cities20_norm[j,i] <- (cities20_norm[j,i] - min(cities19_20[,i], na.rm = T)) / (max(cities19_20[,i], na.rm = T) - min(cities19_20[,i], na.rm = T)) } } #summary(cities19_norm) #summary(cities20_norm) # Taking care of reverse scales cities19_norm$no2_per_population <- 1 - cities19_norm$no2_per_population cities20_norm$no2_per_population <- 1 - cities20_norm$no2_per_population ``` ```{r message=FALSE, warning=FALSE, include=FALSE} ############################################### KPI 2020 # Summing items by subgroups cities20_norm$Risk_Assessment <- rowMeans(cities20_norm[,c('boundary_assess', 'n_assess', 'areas_assess')]) cities20_norm$Adaptation_Plan <- rowMeans(cities20_norm[,c('n_doc_adaptation', 'n_areas_adaptation_plan', 'n_hazards_plan', 'boundary_adaptation_plan', 'stage_adaptation_plan', 'adaptation_plan_type', 'localgov_assess_adaptation_plan')]) cities20_norm$Adaptation_Actions <- rowMeans(cities20_norm[,c('n_hazards_addressed', 'n_actions_per_hazard', 'n_adaptation_actions', 'adaptation_action_status', 'n_adaptation_cobenefits')]) cities20_norm$Mitigation_Actions <- rowMeans(cities20_norm[,c('inventory', 'target_metric_tonnes', 'pcnt_target_achieved', 'n_mitigation_actions', 'mitigation_action_status', 'no2_per_population')]) cities20_norm$Mitigation_Plan <- rowMeans(cities20_norm[,c('n_doc_mitigation', 'n_areas_mitigation_plan', 'boundary_mitigation_plan', 'stage_mitigation_plan', 'localgov_assess_mitigation_plan')]) cities20_norm$Opportunities <- rowMeans(cities20_norm[,c('n_opportunities', 'n_collab_areas', 'n_collab_types', 'n_project_areas')]) cities20_norm$Renewable_Energy_Target <- rowMeans(cities20_norm[,c('percentage_renewable_energy', 'renewable_energy_achieved', 'renewable_energy_source', 'percentage_energy_efficiency_improvement', 'energy_efficiency_achieved')]) cities20_norm$Additional_Measures <- rowMeans(cities20_norm[,c('percentage_walking_cycling', 'zero_emission_zone_presence', 'food_policies_presence')]) cities20_norm$Water <- rowMeans(cities20_norm[,c('n_water_adaptation_actions', 'water_adaptation_action_status', 'n_water_management', 'stage_water_management_plan')]) cities20_norm$Social_Equity <- rowMeans(cities20_norm[,c('vulnpop_assess', 'vuln_addressed_per_total', 'social_impact_addressed_per_total')]) # TOTAL KPI cities20_KPIs <- cities20_norm[, c('id', 'org', "Risk_Assessment" , "Adaptation_Plan" , "Adaptation_Actions" , "Mitigation_Actions" , "Mitigation_Plan", "Opportunities", "Renewable_Energy_Target" , "Additional_Measures", "Water", "Social_Equity", 'hazards_intensity', 'population', 'mean_temp', 'sum_precip', 'gdp_per_capita')] # Computing the final index cities20_KPIs$KPI_total <- 100*rowMeans(cities20_KPIs[,3:(ncol(cities20_KPIs)-5)]) #hist(cities20_KPIs$KPI_total, breaks = 20) # summary(cities20_KPIs) ``` ```{r message=FALSE, warning=FALSE, include=FALSE} ############################################### KPI 2019 # Summing items by subgroups cities19_norm$Risk_Assessment <- rowMeans(cities19_norm[,c('boundary_assess', 'n_assess', 'areas_assess')]) cities19_norm$Adaptation_Plan <- rowMeans(cities19_norm[,c('n_doc_adaptation', 'n_areas_adaptation_plan', 'boundary_adaptation_plan', 'stage_adaptation_plan', 'adaptation_plan_type', 'localgov_assess_adaptation_plan')]) cities19_norm$Adaptation_Actions <- rowMeans(cities19_norm[,c('n_hazards_addressed', 'n_actions_per_hazard', 'adaptation_action_status', 'n_adaptation_cobenefits')]) cities19_norm$Mitigation_Actions <- rowMeans(cities19_norm[,c('inventory', 'target_metric_tonnes', 'pcnt_target_achieved', 'n_mitigation_actions', 'mitigation_action_status', 'no2_per_population')]) cities19_norm$Mitigation_Plan <- rowMeans(cities19_norm[,c('n_doc_mitigation', 'boundary_mitigation_plan', 'stage_mitigation_plan', 'localgov_assess_mitigation_plan')]) cities19_norm$Opportunities <- rowMeans(cities19_norm[,c('n_opportunities', 'n_collab_areas', 'n_project_areas')]) cities19_norm$Renewable_Energy_Target <- rowMeans(cities19_norm[,c('percentage_renewable_energy', 'renewable_energy_achieved', 'renewable_energy_source', 'percentage_energy_efficiency_improvement', 'energy_efficiency_achieved')]) cities19_norm$Additional_Measures <- rowMeans(cities19_norm[,c('percentage_walking_cycling', 'zero_emission_zone_presence', 'food_policies_presence')]) cities19_norm$Water <- rowMeans(cities19_norm[,c('n_water_adaptation_actions', 'water_adaptation_action_status', 'n_water_management', 'stage_water_management_plan')]) cities19_norm$Social_Equity <- rowMeans(cities19_norm[,c('vulnpop_assess', 'vuln_addressed_per_total', 'social_impact_addressed_per_total')]) # TOTAL KPI cities19_KPIs <- cities19_norm[, c('id', 'org', "Risk_Assessment" , "Adaptation_Plan" , "Adaptation_Actions" , "Mitigation_Actions" , "Mitigation_Plan", "Opportunities", "Renewable_Energy_Target" , "Additional_Measures", "Water", "Social_Equity", 'hazards_intensity', 'population', 'mean_temp', 'sum_precip', 'gdp_per_capita')] # Computing the final index cities19_KPIs$KPI_total <- 100*rowMeans(cities19_KPIs[,3:(ncol(cities19_KPIs)-5)]) #hist(cities19_KPIs$KPI_total, breaks = 20) #summary(cities19_KPIs) # Recode to UTF-8 2019 cities19_KPIs_org <- cities19_KPIs$org Encoding(cities19_KPIs_org) <- "UTF-8" cities19_KPIs$org <- cities19_KPIs_org ``` > The global disruption caused by the coronavirus pandemic has unleashed on human civilisation numerous environmental and climate consequences. Although the decline in human activity has temporarily led to the reduction in air and water pollutions and carbon/nitrogen oxides emissions (Le Quéré et al. 2020), enlarged investments in the structural transitions to sustainable energy resources, the outbreak triggered a bunch of unfavourable effects. The latter include the impediment to climate change diplomacy efforts, a worldwide recession capable of cutting down on green technologies subsidising (Newburger 2020), and others. More immediate threats of the pandemic brought down the environmental emergency, which has been insufficiently prioritised even before the coronavirus transmission even by leaders promoting sustainability over momentary profits. It is apparent that the problem of ample response to climate change is becoming increasingly acute, and different actors such as cities and corporations need to leverage practical up-to-date guidelines to ameliorate climate actions and build a brighter future for all populations. The current research attempts to contribute to these issues deploying mostly data collected by the Carbon Disclosure Project (CDP) and utilizing rigorous and diverse analytic instruments. > This notebook is organised as follows. First, we describe a performance assessment methodology and outline the reasoning behind its development. Second, we carry out an efficiency analysis, providing a scrupulous treatment of the social equity domain. This section suggests where concretely inefficient cities should channel their funds to catch up with their relatively successful peers and what strategies they may need to adapt to accelerate the introduction of efficient environmental measures especially at times of economic fallout. Third, with the help of Natural Language Processing some data-driven environmental areas, of which cities and corporations together may take advantage and fortify remunerative environmental outcomes, are put forward. Forth, we define some valuable points suitable for stimulating transparency in reporting practices of cities. The study concludes by summarising the key obtained findings and accentuating their climate policy relevance. # Methodology for the Key Performance Indicators Construction To bring in evidence of cities' environment-related activity, we elaborated an array of metrics, or Key Performance Indicators (KPIs), capturing how city leaders manage the climate change threats. To develop viable metrics, one needs to take on an approach that equitably addresses environmental concerns, aiming to promote a healthier and happier life when the most vulnerable populations are prioritised. The consideration of equity in our work has been made directly in the KPIs construction as well as in the employed methodology that estimates cities' effectiveness conditional on an empirical variation in cities' feasible resources that induce unequal chances to succeed in parrying environmental thrusts. We came up with the following distinct blocks of KPIs that broadly correspond to the general topics raised in the survey, mirror objectives delineated in Paris Agreements, UN’s Sustainable Development Goals and are strongly supported by cooperation movements and organisations such as C40 Cities or the Global Covenant of Mayors: 1. Social Equity: represents the accountability of vulnerable populations in risk assessments and actions. 2. Risk Assessment: reflects the extent of climate change risk and vulnerability assessment done. 3. Adaptation Plan: points on the development of a plan that addresses climate change adaptation. 4. Adaptation Actions: depicts the scope of adaptation actions undertaken to soften the risk to, and vulnerability of, a city’s services, infrastructure, citizens, and businesses from climate change. 5. Mitigation Plan: points on the development of a plan that addresses climate change mitigation. 6. Mitigation Actions: depicts the scope of mitigation actions undertaken to alleviate the risks caused by climate change. 7. Opportunities: indicates the expression of cities’ opportunities identification because of addressing climate change. 8. Renewable Energy Target: shows the readiness of goals connected with renewable energy or electricity. 9. Water: describes the achievements of cities in ensuring water security. 10. Additional Measures: reflects the scale of cities’ success in tailoring policies related to food consumption and transport. The table below demonstrates all these indicators and provides details on their measurement. To compute the metrics, we adopted a similar strategy used in the “Ease of Doing Business” (e.g., World Bank 2020) annual reports. Specifically, we applied a min-max normalisation to each item and averaged their values within the respective sub-metrics. To create the aggregated scores, all sub-indices were again averaged and scaled from 0 to 100. Scaling procedure accounted for the differences in data parameters in various measurement periods to make performance metrics dynamically comparable. In the present analysis, 2019 and 2020 datasets were exploited; therefore, cities participated only in both surveys were retained and observations with a substantial amount of missing values (> 10 unanswered questions of interest in 2020) were filtered, which left a sample size of 342. A relatively small fraction of missing responses within the resulting items dataset was imputed using the bootstrap EM algorithm for incomplete multivariate data. "Not Applicable" options were replaced with 0's since they manifest the absence of a response due to the preceding filtering in the survey. ```{r} # KPI description KPIs_description <- rio::import('../input/kpis-description/KPIs_Description.xlsx') KPIs_description[is.na(KPIs_description)] <- '' kbl(KPIs_description, align = "l") %>% kable_paper(full_width = F) %>% column_spec(1, bold = T) %>% scroll_box(width = "100%", height = "400px") ``` The graph below portrays the distribution of the total KPI scores by years. It is noticeable that the scores show prominently high and normal-like variability in both periods. Additionally, in 2020 the distribution of the overall KPI measure is lightly shifted to the right relative to the 2019 distribution, highlighting an overall improvement of cities' performance within a year (difference in averages for both years is statistically significant, around four percentage points). The examples of cities that display the highest positive drift from 2019 to 2020 in their total performance proxy (top 20) are shown in the following heatmap. ```{r message=FALSE, warning=FALSE, fig.height=4, fig.width=8} ####################################################### Reliability ####################################################### ###### Correlation # 2020 #cor.test(cities20_KPIs$KPI_total, cities20_KPIs$gdp_per_capita) KPI_top <- cities20_KPIs %>% arrange(desc(KPI_total)) #plot(cities20_KPIs$KPI_total, cities20_KPIs$gdp_per_capita) cities20_KPIs <- cities20_KPIs %>% filter(!id == 31090) # 2020 #cor.test(cities19_KPIs$KPI_total, cities19_KPIs$gdp_per_capita) KPI_top <- cities19_KPIs %>% arrange(desc(KPI_total)) #plot(cities19_KPIs$KPI_total, cities19_KPIs$gdp_per_capital) cities19_KPIs <- cities19_KPIs %>% filter(!id == 31090) ###### Densities # Combining data for 2019 and 2020 # Joint dataset cities20_KPIs_year <- cbind.data.frame(year = 2020, cities20_KPIs) cities19_KPIs_year <- cbind.data.frame(year = 2019, cities19_KPIs) cities19_20_KPIs_year <- rbind.data.frame(cities20_KPIs_year, cities19_KPIs_year) %>% arrange(year) ##### Difference in total KPI cities19_20_KPIs <- cities20_KPIs[,c('org', 'KPI_total')] %>% left_join(cities19_KPIs[,c('org', 'KPI_total')], by = 'org') %>% mutate(diff_KPI = KPI_total.x - KPI_total.y) %>% arrange(desc(diff_KPI)) colnames(cities19_20_KPIs) <- c('City', "Total KPI 2020", 'Total KPI 2019', 'Difference') # Selecting top 20 cities19_20_KPIs_top20 <- cities19_20_KPIs %>% slice(1:20) # Graphing mean_KPI_2019 <- cities19_20_KPIs_year %>% group_by(year) %>% summarise(mean_KPI=mean(KPI_total)) %>% filter(year==2019) mean_KPI_2020 <- cities19_20_KPIs_year %>% group_by(year) %>% summarise(mean_KPI=mean(KPI_total)) %>% filter(year==2020) cities19_20_KPIs_year %>% ggplot(aes(x = KPI_total, fill = factor(year))) + geom_density(alpha=.3) + labs(x = "Total KPI", y = "Density", fill="") + scale_fill_manual(values=c('darksalmon','deepskyblue2')) + ggtitle("KPI Distributions for 2019 and 2020") + geom_vline(aes(xintercept=mean(mean_KPI_2019$mean_KPI, na.rm=T)), # Ignore NA values for mean color="darksalmon", linetype="dashed", size=1) + geom_vline(aes(xintercept=mean(mean_KPI_2020$mean_KPI, na.rm=T)), # Ignore NA values for mean color="deepskyblue2", linetype="dashed", size=1) + theme_classic() + theme(plot.title = element_text(size=18, face="bold", hjust=0.5), legend.position="bottom") + ggplot2::annotate("text", label = paste0(as.character(round(mean(mean_KPI_2019$mean_KPI, na.rm=T), 1)), "\nMean KPI, 2019"), x = 25.5, y = 0.01, size = 4, colour = "darksalmon") + ggplot2::annotate("text", label = paste0(as.character(round(mean(mean_KPI_2020$mean_KPI, na.rm=T), 1)), "\nMean KPI, 2020"), x = 40.5, y = 0.01, size = 4, colour = "deepskyblue2") ``` ```{r message=FALSE, warning=FALSE} # Heatmap cities19_20_KPIs_top20[,c(1,3,2,4)] %>% kbl() %>% kable_styling(bootstrap_options = "striped", font_size = 12) %>% kable_paper(full_width = T) %>% # kable_styling(position = "center") %>% column_spec(4, color = "white", background = spec_color(cities19_20_KPIs_top20$Difference, end = 0.7), popover = paste("am:", cities19_20_KPIs_top20$Difference)) ``` However, KPIs along are not instrumental in determining the efficiency of cities as distinct decision-making units (DMUs): to evaluate cities’ relative productivity, one needs to make use of the differences in cities’ starting conditions, or inputs. The scatter plot below illustrates that inquiry by representing the relationship between the total KPI and cities attributes: GDP per capita, population size, and mean annual temperature. Data on GDP per capita were collected manually from multiple open sources: OECD, Wikipedia articles, and websites with local administrative statistics; the annual temperature was taken from the Google Earth Engine. It can be seen that on average higher GDP per capita, lower temperatures, and lower population sizes are associated with a larger performance indicator. This allows us to contend that cities pertaining to Global North seem to manage climate change more effectively than cities in the Global South. ```{r message=FALSE, warning=FALSE, fig.height=7, fig.width=7} options(scipen=10000) # Total KPI and cities attributes df_additional_info <- rio::import('/kaggle/input/df-cities-infoadditional/df_cities_info.additional.xlsx') cities20_KPIs_plot <- cities20_KPIs %>% dplyr::select(id, KPI_total) cities20_KPIs_plot <- cities20_KPIs_plot %>% left_join(df_additional_info) %>% mutate(population=population/1e6) ggplot(cities20_KPIs_plot, aes(x=gdp_per_capita, y=KPI_total)) + geom_point(aes(size=population, color = mean_temp)) + geom_smooth(method = lm, se = T, color = "darkgrey", size=0.5) + scale_color_gradient(low="deepskyblue2", high="darksalmon") + labs(x = "GDP per capita, $", y = "Total KPI", size="Population (millions)", col="Mean Annual Temperature, °C") + ggtitle("Total KPI and GDP per capita") + theme_classic() + theme(plot.title = element_text(size=18, face="bold", hjust=0.5), legend.position="bottom") ``` However, it might be unfair to assign the latter worse practices since the former cities possess vast structural opportunities to be effective: they are richer, less unequal, and carry the smaller intensity of hazards due to favourable climatic conditions and lower populations. Therefore, an equitable analysis would consider those baseline circumstances when gauging cities' achievement level in handling environmental threats. Computationally, this purpose is reachable through the Data Envelopment Analysis (DEA), which is a nonparametric analytical technique used in economics and operations research for the assessment of "the best-practice frontiers" (Sickles and Zelenyuk 2019). # Estimating the Efficiency of Cities’ Performance In brief, DEA looks at the decision surface formed by multiple linear equations connecting the given inputs (resources) of units with their respective outputs (performance indices). By using linear programming tools, DEA solves an objective function that obtains the optimal linear decision surface, or envelope, splitting the space such that all points are either on the envelope (efficient units) or beneath it (non-efficient units). Linear programming maximizes the following ratio: $$ \frac{virtual \space output}{virtual \space input} = \frac{u_oy_{1o} + ... + u_sy_{so}}{v_ox_{1o} + ... + v_mx_{mo}} $$ where s is the amount of outputs, m is the amount of inputs, (u,v) are the weights assigned to different inputs and outputs, $y_{jo}$ is the jth output of DMU observed $(j = 1, ..., s)$, $x_{io}$ is the ith input of DMU observed $(i = 1, ..., m)$. We employed GDP per capita, population size, and mean annual temperature as inputs (the last 2 were denoted as uncontrollable), whereas KPI indices – as outputs. An output-oriented model, assuming an output maximisation under the given inputs, with variable returns to scale was calculated. The histograms below reflect the amount of efficient and non-efficient cities as well as the distribution of the computed efficiency scores for non-efficient units for descriptive purposes. It can be noticed that the number of efficient cities slightly exceeds the number of non-efficient ones, and the efficiency score for many lagging units approaches 100%. Almost all cities exceed the efficiency of 50%, which implies that there are no dramatic imbalances in cities' efficient level of inputs used. ```{r, fig.height=5, fig.width=10} #################################### Data Envelopment Analysis ############################# ##### 2020 # Read data cities20_dear <- read_data(datadea = cities20_KPIs, dmus = 2, inputs = c('gdp_per_capita', 'population', 'mean_temp'), outputs = c("Risk_Assessment", "Adaptation_Plan", "Adaptation_Actions", "Mitigation_Actions", "Mitigation_Plan", "Opportunities", "Renewable_Energy_Target", "Additional_Measures", "Water", "Social_Equity"), nc_inputs = 2:3) # Run a basic model accounting for uncontrollable inputs deaR_model_20 <- model_basic(cities20_dear, orientation = 'oo', rts = 'vrs') #deaR_boot <- bootstrap_basic(cities20_dear, orientation = 'oo', rts = 'vrs') ##### Graphing # Barplot of efficient vs. non-efficient eff_20 <- data.frame(org = deaR_model_20$data$dmunames, eff_20 = 1/efficiencies(deaR_model_20)) eff_20 %>% mutate(iseff = ifelse(eff_20 == 1, 1 , 0)) -> eff_20 eff_20 <- eff_20 %>% left_join(cities20_KPIs[, c('org', 'gdp_per_capita')], by = 'org') eff_20$gdp_group <- ifelse(eff_20$gdp_per_capita >= mean(eff_20$gdp_per_capita), 1, 0) p1 <- eff_20 %>% ggplot(aes(x = factor(iseff), fill = factor(iseff))) + geom_bar(stat="count", colour="black", width=0.7) + ylab("Count") + xlab("") + scale_x_discrete(labels = c("Non-efficient", "Efficient")) + scale_fill_manual(values=c("coral3", "deepskyblue3")) + ggtitle("Efficient / Non-Efficient Cities") + geom_text(stat = 'count', aes(label = ..count..), vjust = -0.5) + theme_classic() + theme(plot.title = element_text(size=18, face="bold", hjust=0.5), legend.position="bottom") + guides(fill=FALSE) p2 <- eff_20 %>% filter(iseff == 0) %>% arrange(desc(eff_20)) %>% ggplot(aes(x = eff_20, fill = "coral3")) + geom_histogram(colour='black', breaks = c(seq( from = min(eff_20$eff_20), to = max(eff_20$eff_20[eff_20$iseff != 1]), length.out = 10 ))) + scale_fill_identity()+ ggtitle("Distribution of Efficency for\nNon-Efficient Cities") + ylab("Count") + xlab('Efficiency Score for Non-Efficient Cities') + theme_classic() + theme(plot.title = element_text(size=18, face="bold", hjust=0.5), legend.position="bottom") ggarrange(p1, p2, ncol = 2, nrow = 1) ``` DEA allows a researcher to calculate "slacks", which in our model specification mean the exact increase in outputs that can be done under the given inputs for less efficient units to achieve efficiency scores of successful units with similar input structure. The following heatmap visualises the resulting slacks by cities with above and below-average GDP levels; the higher intensity of the red colour represents greater performance improvement suggested by the model. The graph informs us that both groups of cities are quite efficient in terms of adaptation and mitigation actions, exploiting climate change opportunities, renewable energy performance, and regarding actions associated with water, food consumption, and transport. However, on average cities experience relatively pronounced problems with developing plans that consider climate change adaptation and mitigation. This signifies that generally, cities (irrespective of their affluence status) are capable of acting efficiently when addressing environmental hazards, but these practices seem to be poorly internalised by the whole policy system as the actions may not be carefully planned and outlined. Also, the graph shows that Social Equity issues are fulfilled in a substantially lower manner by less wealthy countries compared to wealthier countries even controlling for the discrepancy in environmental, demographic, and economic conditions in these two groups of cities. Policy leader will get ahead considerably from directing more efforts and investments into ensuring equitable treatment of vulnerable groups in less affluent cities. ```{r message=FALSE, warning=FALSE, fig.height=10, fig.width=10} # Slacks: heatmap by gdp # Preparing data slacks_20 <- as.data.frame(slacks(deaR_model_20)) slacks_20$org <- rownames(slacks_20) slacks_eff_20 <- slacks_20 %>% select(4:ncol(slacks_20)) %>% left_join(eff_20, by = 'org') colnames(slacks_eff_20) <- c("Risk_Assessment", "Adaptation_Plan", "Adaptation_Actions", "Mitigation_Actions", "Mitigation_Plan", "Opportunities", "Renewable_Energy_Target", "Additional_Measures", "Water", "Social_Equity",'org', 'eff_20', 'iseff', 'gdp_per_capita', 'gdp_group') rownames(slacks_eff_20) <- slacks_eff_20$org # Full high gdp slacks_eff_high_gdp <- slacks_eff_20 %>% filter(iseff == 0 & gdp_group == 1) %>% arrange(eff_20) %>% select(-c(org, eff_20, iseff, gdp_per_capita, gdp_group)) slacks_eff_high_gdp <- as.matrix(slacks_eff_high_gdp) # Full low gdp slacks_eff_low_gdp <- slacks_eff_20 %>% filter(iseff == 0 & gdp_group == 0) %>% arrange(eff_20) %>% select(-c(org, eff_20, iseff, gdp_per_capita, gdp_group)) slacks_eff_low_gdp <- as.matrix(slacks_eff_low_gdp) # Plotting colnames(slacks_eff_high_gdp) <- gsub("_", " ", colnames(slacks_eff_high_gdp)) colnames(slacks_eff_low_gdp) <- gsub("_", " ", colnames(slacks_eff_low_gdp)) pheatmap(slacks_eff_high_gdp, cluster_rows=T, cluster_cols=F, angle_col = "45", color = colorRampPalette(c("chartreuse3", "white", "brown2"))(100), border_color="white", display_numbers=F, legend=T, main="Slacks for Cities, GDP per capita above average", fontsize=14, fontsize_row=9, fontsize_col=12) pheatmap(slacks_eff_low_gdp, cluster_rows=T, cluster_cols=F, angle_col = "45", color = colorRampPalette(c("chartreuse3", "white", "brown2"))(100), border_color="white", display_numbers=F, legend=T, main="Slacks for Cities, GDP per capita below average", fontsize=14, fontsize_row=9, fontsize_col=12) ``` DEA output also provides peer efficient ultimate benchmark targets, i.e. for each non-efficient city an algorithm gauges a set of best-performing cities that have a similar scope of inputs. In other words, DEA will not assign a high-income peer to a low-income unit since the former do not possess the same volume of resources at their disposal. The directed network below portrays the connection between non-efficient cities and their more successful peers which the former should resemble. Red nodes depict cities with greater GDP per capita level, while blue nodes represent cities with smaller GDP per capita level. It is visible that Greater Sudbury, Boulder, Bloomington, Encinitas, are most frequently nominated (have the largest in-degree centrality weighted by lambdas) as efficient peers for richer cities, whereas Braga, Davao, Hermosillo, Pereira, Tagum, Providencia, Guimaraes seem to be leadership holders in a group of less affluent cities. ```{r fig.height=10, fig.width=10} ##### 2019 # Read data cities19_dear <- read_data(datadea = cities19_KPIs, dmus = 2, inputs = c('gdp_per_capita', 'population', 'mean_temp'), outputs = c("Risk_Assessment", "Adaptation_Plan", "Adaptation_Actions", "Mitigation_Actions", "Mitigation_Plan", "Opportunities", "Renewable_Energy_Target", "Additional_Measures", "Water", "Social_Equity"), nc_inputs = 2:3) # Run a basic model accounting for uncontrollable inputs deaR_model_19 <- model_basic(cities19_dear, orientation = 'oo', rts = 'vrs') ##### Malmquist index # Specifying read_data for the index data_malmquist <- read_malmquist(datadea = cities19_20_KPIs_year, percol = 'year', arrangement = "vertical", dmus = 3, inputs = c('gdp_per_capita', 'population', 'mean_temp'), outputs = c("Risk_Assessment", "Adaptation_Plan", "Adaptation_Actions", "Mitigation_Actions", "Mitigation_Plan", "Opportunities", "Renewable_Energy_Target", "Additional_Measures", "Water", "Social_Equity"), nc_inputs = 2:3) malmquist <- malmquist_index(data_malmquist, orientation = "oo", rts = 'vrs') mi <- data.frame(t(malmquist$mi)) colnames(mi) <- 'MI' tc <- as.data.frame(t(malmquist$tc)) colnames(tc) <- 'TC' pech <- as.data.frame(t(malmquist$pech)) colnames(pech) <- 'PTEC' # Malmquist output # Extracting efficiencies for 2019 eff_19 <- as.data.frame(100*1/efficiencies(deaR_model_19)) eff_19$org <- rownames(eff_19) rownames(eff_19) <- NULL colnames(eff_19) <- c('eff_19', 'org') # Multiplying efficiencies for 2020 by 100 eff_20$eff_20_normed <- 100*eff_20$eff_20 # Combining all pieces of the output Malmquist_output <- cbind.data.frame(mi, tc, pech) Malmquist_output$org <- rownames(Malmquist_output) Malmquist_output <- Malmquist_output[, c(4,1:3)] rownames(Malmquist_output) <- NULL Malmquist_output <- Malmquist_output %>% left_join(eff_20, by = 'org') %>% left_join(eff_19, by = 'org') %>% arrange(desc(gdp_per_capita, gdp_group)) %>% na.omit(Malmquist_output) # Selecting relevant columns Malmquist_output <- Malmquist_output[,c('org', 'MI', 'TC', 'PTEC', 'eff_19', 'eff_20_normed', 'gdp_group')] colnames(Malmquist_output) <- c('Cities', 'MI', 'TC', 'PTEC', 'Efficiency 2019', 'Efficiency 2020', 'gdp_group') # Sampling 20 cities for demonstration purposes set.seed(123) sample_id <- sample(Malmquist_output$Cities, 20) Malmquist_output_20 <- Malmquist_output[Malmquist_output$Cities %in% sample_id,] Malmquist_output_20 <- Malmquist_output_20 %>% dplyr::select(-gdp_group) %>% arrange(desc(`Efficiency 2019`)) Malmquist_output_20[,2:6] <- round(Malmquist_output_20[,2:6], 2) Malmquist_output <- Malmquist_output[!is.infinite(Malmquist_output$MI),] Malmquist_output$gdp_group <- NULL # Adding average by whole sample Averages <- Malmquist_output %>% summarise_if(is.numeric, mean) Malmquist_output_20 <- rbind.data.frame(Malmquist_output_20, cbind.data.frame('Cities' = 'Total average', Averages)) ``` ```{r message=FALSE, warning=FALSE, fig.height=10, fig.width=10} ##### References: peers # References: peers library(igraph) # Extracting lambdas lmbd <- lambdas(deaR_model_20) lmbd <- matrix(lmbd[complete.cases(lmbd), ], ncol = ncol(lmbd), dimnames = dimnames(lmbd)) # Make the matrix square in case dmu_eval != dmu_ref dmunames <- deaR_model_20$data$dmunames[deaR_model_20$dmu_eval] ref <- references(deaR_model_20) refnames <- refnames <- unique(unlist(lapply(ref, function (x) names(x)))) urefnames <- names(ref) effdmus <- dmunames[which(!dmunames %in% urefnames)] urefnames <- refnames[refnames %in% effdmus] allnames <- unique(unlist(dimnames(lmbd))) # An empty matrix for adjacency extendedlmbd <- matrix(0, nrow = length(allnames), ncol = length(allnames), dimnames = list(allnames, allnames)) extendedlmbd[dimnames(lmbd)[[1]],dimnames(lmbd)[[2]]] <- lmbd # Creating adjacency matrix adjmatrix <- extendedlmbd G <- graph_from_adjacency_matrix(adjmatrix, diag = FALSE, weighted = T, mode = 'directed') df_temp <- data.frame(names(V(G))) colnames(df_temp) <- "org" df_temp <- df_temp %>% left_join(df_additional_info %>% dplyr::select(org, gdp_per_capita)) G <- set_vertex_attr(G, 'gdp_per_capita', value = df_temp$gdp_per_capita) G <- delete_vertices(G, which(is.na(V(G)$gdp_per_capita))) # Adding weighted degree centrality as an attribute in_deg <- strength(G, mode = 'in', weights = get.edge.attribute(G)$weights) G <- set_vertex_attr(G, 'in_deg', value = in_deg) # Rem0ving efficient cities - isolated nodes Isolated = which(degree(G)==0) G2 = delete.vertices(G, Isolated) G2 <- readRDS("/kaggle/input/g2-graph/G2") ``` ```{r} # Network Plot with peers set.seed(25) par(mar=c(0,0,0,0)) # Set vertex size vertex_size = log(get.vertex.attribute(G2)$in_deg + 2)*3 vertex_size <- ifelse(vertex_size<1.5, 1.5, vertex_size) # Set vertex color gdp <- get.vertex.attribute(G2)$gdp_per_capita in_degree <- get.vertex.attribute(G2)$in_deg vertex_color <- ifelse( (gdp >= mean(gdp) & in_degree>=3), 'deepskyblue2', ifelse((gdp >= mean(gdp) & in_degree<3), adjustcolor('deepskyblue2', alpha=0.5), ifelse((gdp < mean(gdp) & in_degree>=3), 'darksalmon', adjustcolor('darksalmon', alpha=0.5) ))) # Set vertex label vertex_label <- get.vertex.attribute(G2)$name vertex_label_temp <- ifelse(in_degree>=3, vertex_label, NA) vertex_label <- dplyr::recode(vertex_label, "City of Boulder"="Boulder", "City of Greater Sudbury / Grand Sudbury"="Greater Sudbury", "Ayuntamiento de Hermosillo"="Hermosillo", "City of Bloomington"="Bloomington", "City Government of Davao"="Davao", "Alcaldia de Pereira"="Pereira", "Município de Guimarães"="Guimarães", "Município de Braga"="Braga", "Municipalidad de Providencia"="Providencia", "City of Encinitas, CA"="Encinitas", "Tagum City"="Tagum") vertex_label <- ifelse(in_degree>=3, vertex_label, NA) # Set vertex frame color vertex_frame_color <- ifelse(in_degree>=3, adjustcolor('black', alpha=0.5), "lightgrey") # Sed edge color edgelist <- data.frame(as_edgelist(G2)) edge_color <- ifelse(edgelist$X2 %in% c(na.omit(unique(vertex_label_temp, na.rm=T))), adjustcolor('darkgrey', alpha=0.5), adjustcolor('lightgrey', alpha=0.3)) # Set edge width edge_width <- E(G2)$weight set.seed(12) plot(G2, edge.arrow.size = 0.3, layout=layout.mds, # <<< OK2 vertex.label = vertex_label, vertex.size = vertex_size, vertex.color = vertex_color, vertex.label.color = "black", vertex.label.font = 2, vertex.label.family = "sans", vertex.label.dist=0.75, vertex.label.degree=-pi/2, vertex.frame.color=vertex_frame_color, edge.color= edge_color, # adjustcolor('lightgrey', alpha=0.5), edge.width=edge_width*3, edge.curved=T) ``` Each efficient peer has its weight or a lambda value, and a linear combination of 'weighted' peers constitutes the target position for a unit. The next table embraces names for all cities laying under the efficient frontier and prints their peers sorted by lambdas such that the first peer is the most important for a given unit, whereas the last one – the least. For instance, the following cities are the best role models for Chicago: San Francisco, Berlin, New York, Riga. City mayors can significantly simplify the delivery of solutions to climate change challenges by orienting towards such sets of peers when enacting government interventions. Peers analysis will replenish and bring in more credibility to the already existing C40 Network collaborations that were designed to connect city officials with peers, inspire innovation, advise city peers based on similar experiences, and impact international and national policy agendas. ```{r} # Full table with peers # Sorting by lambdas within each unit ref_sorted <- lapply(seq_along(ref), function(i) {data.frame('lambda' = ref[[i]]) %>% arrange(desc(lambda))}) # Extracting peers ref_dframe <- cbind.data.frame('Cities' = names(ref), 'Peers' = unlist(lapply(seq_along(ref_sorted), function(i){paste0(rownames(ref_sorted[[i]]), collapse = ', ')}))) # Plotting kbl(ref_dframe) %>% kable_paper() %>% scroll_box(width = "100%", height = "400px") ``` DEA enables a researcher to trace productivity changes over time. In a nutshell, this can happen due to changes in pure technical efficiencies (PTEC) of units (i.e., catching up with their frontier) and/or technological changes (TC) (i.e., shifts of a frontier over time). The Malmquist Index (MI) aggregates information about these two processes. The panel below shows the efficiency transformations of cities in 2020 compared to 2019 for 20 randomly sampled cities and their respective yardsticks for illustrative purposes. We can interpret a fraction of examples from this output. ```{r} # Printing kable(Malmquist_output_20, caption = 'Cities Efficiency over Time', booktabs = TRUE, "html", align = 'lccccc') %>% kable_styling("striped", full_width = F) ``` For instance, by 2020 the city of Helsinki has made a visible efficiency improvement (MI > 1) thanks to both its own better performance (PECH > 1) and shifts of the envelope (TC > 1). The city of Melbourne retained its overall productive level of 100; however, its overall productivity indicator slightly dropped down (MI = 0.98). This shift was composed of 2 processes: (1) an upward relocation of the efficient frontier (better performance of other cities) and (2) Melbourne's faintly diminished ability to make the maximal possible amount of outputs under the given level of inputs. In contrast, Kuala Lumpur city lost a sizeable part of its efficiency only due to its own worsened performance (PECH < 1), whereas the frontier has become even more favourable to this city (TC > 1). The total average change of efficiencies is relatively stagnant and remains at the level of 94%. Overall, such an analysis is profitable for detailing the nature of dynamic changes in the efficiency of the units. # Hazards – Action – Co-benefits: a Relational Approach ```{r message=FALSE, warning=FALSE, include=FALSE} library(visNetwork) library(igraph) library(rio) library(arules) library(arulesViz) ######################### ACTION-HAZARDS-BENIFITS NETWORK PREPROCESSING ############################# # --- Load the 2020 data setwd(paste0("../input/cdp-unlocking-climate-solutions/Cities/Cities Responses")) df_cities.20 <- read.csv('2020_Full_Cities_Dataset.csv', stringsAsFactors=FALSE) %>% dplyr::select('id' = Account.Number, 'org' = Organization, 'cntry' = Country, 'region' = CDP.Region, 'qstn' = Question.Number, 'qstname' = Question.Name, 'coln' = Column.Number, 'colname' = Column.Name, 'rown' = Row.Number, 'rowname' = Row.Name, 'resp' = Response.Answer) # Recode to UTF-8 temp_vec <- df_cities.20$resp Encoding(temp_vec) <- "UTF-8" df_cities.20$resp <- temp_vec temp_vec <- df_cities.20$colname Encoding(temp_vec) <- "UTF-8" df_cities.20$colname <- temp_vec temp_vec <- df_cities.20$org Encoding(temp_vec) <- "UTF-8" df_cities.20$org <- temp_vec # Recode 'Other' category df_cities.20$resp[grepl("^Other, please specify", as.character(df_cities.20$resp))] <- "Other" # --- Load the 2019 data df_cities.19 <- read.csv('2019_Full_Cities_Dataset.csv', stringsAsFactors=FALSE) %>% dplyr::select('id' = Account.Number, 'org' = Organization, 'cntry' = Country, 'region' = CDP.Region, 'qstn' = Question.Number, 'qstname' = Question.Name, 'coln' = Column.Number, 'colname' = Column.Name, 'rown' = Row.Number, 'rowname' = Row.Name, 'resp' = Response.Answer) # Recode to UTF-8 temp_vec <- df_cities.19$resp Encoding(temp_vec) <- "UTF-8" df_cities.19$resp <- temp_vec temp_vec <- df_cities.19$colname Encoding(temp_vec) <- "UTF-8" df_cities.19$colname <- temp_vec # Recode 'Other' category df_cities.19$resp[grepl("^Other:", as.character(df_cities.19$resp))] <- "Other" ############## Cities, Q3.0, 2020, Hazards, Action, Co-benifits: Association Rules # Loading df_q3_0 <- df_cities.20 %>% filter(qstn == "3.0", colname %in% c("Co-benefit area", "Climate hazards", "Action")) %>% mutate(resp=gsub("â€“", " - ", resp)) df_q3_0.hazards <- df_q3_0 %>% filter(colname %in% c("Climate hazards"), resp!="") %>% dplyr::select(id, rown, resp) %>% rename(hazard=resp) %>% mutate(hazard=str_replace(hazard, "\\> ", "")) %>% filter(hazard!="Question not applicable") df_q3_0.actions <- df_q3_0 %>% filter(colname %in% c("Action"), !(resp %in% c("", "Other"))) %>% rename(actions=resp) df_q3_0.benifits <- df_q3_0 %>% filter(colname %in% c("Co-benefit area"), !(resp %in% c("", "Other"))) %>% rename(benifits=resp) # Save freq tables for the future table.hazards.20 <- data.frame(table(df_q3_0.hazards$hazard)) table.actions.20 <- data.frame(table(df_q3_0.actions$actions)) table.benifits.20 <- data.frame(table(df_q3_0.benifits$benifits)) # Aggregate actions and co-benifits by id-rown df_q3_0.actions <- df_q3_0.actions %>% group_by(id, rown) %>% summarise_at(vars(actions), funs(paste(., collapse = ';'))) %>% ungroup() df_q3_0.benifits <- df_q3_0.benifits %>% group_by(id, rown) %>% summarise_at(vars(benifits), funs(paste(., collapse = ';'))) %>% ungroup() # Merge Hazards, Actions and Co-Benifits df_q3_0.wide.20 <- df_q3_0.hazards %>% full_join(df_q3_0.actions, by = c("id", "rown")) %>% full_join(df_q3_0.benifits, by = c("id", "rown")) %>% filter(rowSums(is.na(.)) == 0) %>% mutate(items=paste(hazard, actions, benifits, sep = ";")) %>% dplyr::select(items) ############## Cities, Q3.0, 2019, Hazards, Action, Co-benifits: Association Rules # Loading df_q3_0 <- df_cities.19 %>% filter(qstn == "3.0", colname %in% c("Co-benefit area", "Climate hazards", "Action")) %>% mutate(resp=gsub("â€“", " - ", resp)) df_q3_0.hazards <- df_q3_0 %>% filter(colname %in% c("Climate hazards"), resp!="") %>% dplyr::select(id, rown, resp) %>% rename(hazard=resp) %>% mutate(hazard=str_replace(hazard, "\\> ", "")) %>% filter(hazard!="Question not applicable") df_q3_0.actions <- df_q3_0 %>% filter(colname %in% c("Action"), !(resp %in% c("", "Other"))) %>% rename(actions=resp) df_q3_0.benifits <- df_q3_0 %>% filter(colname %in% c("Co-benefit area"), !(resp %in% c("", "Other"))) %>% rename(benifits=resp) # Save freq tables for the future table.hazards.19 <- data.frame(table(df_q3_0.hazards$hazard)) table.actions.19 <- data.frame(table(df_q3_0.actions$actions)) table.benifits.19 <- data.frame(table(df_q3_0.benifits$benifits)) # Aggregate actions and co-benifits by id-rown df_q3_0.actions <- df_q3_0.actions %>% group_by(id, rown) %>% summarise_at(vars(actions), funs(paste(., collapse = ';'))) %>% ungroup() df_q3_0.benifits <- df_q3_0.benifits %>% group_by(id, rown) %>% summarise_at(vars(benifits), funs(paste(., collapse = ';'))) %>% ungroup() # Merge Hazards, Actions and Co-Benifits df_q3_0.wide.19 <- df_q3_0.hazards %>% full_join(df_q3_0.actions, by = c("id", "rown")) %>% full_join(df_q3_0.benifits, by = c("id", "rown")) %>% filter(rowSums(is.na(.)) == 0) %>% mutate(items=paste(hazard, actions, benifits, sep = ";")) %>% dplyr::select(items) ############## 2019 and 2020 together # Combine main datasets df_q3_0.wide.rules <- rbind(df_q3_0.wide.19, df_q3_0.wide.20) df_q3_0.wide.rules <- df_q3_0.wide.rules %>% distinct(items, .keep_all = T) temp_vec <- df_q3_0.wide.rules$items Encoding(temp_vec) <- "UTF-8" df_q3_0.wide.rules$items <- temp_vec # Combine tables table.actions <- rbind(table.actions.19, table.actions.20) %>% mutate(Var1=as.character(Var1)) table.benifits <- rbind(table.benifits.19, table.benifits.20) %>% mutate(Var1=as.character(Var1)) table.hazards <- rbind(table.hazards.19, table.hazards.20) %>% mutate(Var1=as.character(Var1)) ``` ```{r message=FALSE, warning=FALSE, include=FALSE} ############## Association Rules # setwd(paste0("../output/kaggle/working")) export(df_q3_0.wide.rules, file = "df_q3_0.wide.rules.csv", "csv") tr <- read.transactions("df_q3_0.wide.rules.csv", format = 'basket', sep=';') ``` ```{r message=FALSE, warning=FALSE, include=FALSE} # Model association.rules <- apriori(tr, parameter = list(supp=0.001, conf=0.5, minlen=3, maxlen=3)) # Remove redundant rules. For example remove rule {A, B} -> {C}, if we have rule {A, B, D} -> {C} association.rules <- association.rules[!is.redundant(association.rules)] # Remove association rules with Hazards in RHS for (hazard in table.hazards$Var1){ association.rules <- subset(association.rules, subset = (!rhs %pin% hazard)) } # Select top-50 rules by Confidence and Lift topSubRules <- head(association.rules, n = 50, by = c("confidence", "lift")) # OK ``` ```{r message=FALSE, warning=FALSE, include=FALSE} # Save to igraph ig <- plot(topSubRules, method="graph", control=list(type="items") ) tf <- tempfile( ) saveAsGraph(topSubRules, file = tf, format = "dot" ) # clean up temp file if desired #unlink(tf) # let's bypass saveAsGraph and just use our igraph ig_df <- get.data.frame(ig, what = "both") ig_df$vertices$label <- gsub("â€“", "–", ig_df$vertices$label) # Make a variable with a node type ig_df$vertices$node_type <- ifelse(ig_df$vertices$label %in% table.actions$Var1, "Action", ifelse(ig_df$vertices$label %in% table.hazards$Var1, "Hazard", ifelse(ig_df$vertices$label %in% table.benifits$Var1, "Benifit", "Rule"))) # Dataframe with Nodes information nodes <- data.frame(group = ig_df$vertices$node_type, id = ig_df$vertices$name, value = ig_df$vertices$support, title = ifelse(ig_df$vertices$label == "", ig_df$vertices$name, ig_df$vertices$label), ig_df$vertices) # Add color for groups nodes$color <- ifelse(nodes$label %in% table.actions$Var1, "blue", ifelse(nodes$label %in% table.hazards$Var1, "darkred", ifelse(nodes$label %in% table.benifits$Var1, "green", "lightgray"))) ``` What specifically should cities do to adapt to climate hazards, maintain economic stability and acquire benefits from reacting to these hazards? To gain insights about this issue, we carried out an estimation allowing us to explore the relational structure between different classes of hazards, a variety of actions designed to deal with these hazards, and co-benefits areas associated with these actions. We leveraged question 3.0 from cities dataset to accomplish this task and devised a Road Map capable of navigating cities’ decision-making process. We merged responses from 2019 and 2020 cities databases and created a combined string characteristic embracing three distinct topics: climate hazards, actions, and co-benefits. As a result, each observation in the generated variable contained a path leading from climate hazards to actions undertaken and further to areas of co-benefit. The following example illustrates such a sequence: "chemical change saltwater intrusion – diversification of water supply – resource conservation – improved resource security - enhanced resilience". To analyse these series of records thoroughly and meaningfully, we utilized the Apriori algorithm – one the most popular methods in association rule mining (for details see e.g.: Zhang & Zhang 2003). Such techniques are usually employed to discover the connections between products’ sales in supermarkets. The main idea is to detect rules describing a regularity in the appearance of items within the analysed sequences. The metrics that were exploited in for this piece of our research are as follows: Let $X, Y$ be sets of items (e.g., certain hazards), $X \Rightarrow Y$ an association rule, $T$ a set of paths (e.g., hazard - action - co-benefit). 1. Support Support indicator reflects the frequency of item's appearance in the database $$ \mathrm{supp}(X) = \frac{|\{t \in T; X \subseteq t\}|}{|T|} $$ 2. Confindence Confidence indicator points on how frequent a rule has been established to be true. $$ {conf}(X \Rightarrow Y) = \mathrm{supp}(X \cup Y) / \mathrm{supp}(X) $$ 3. Lift Lift indicator potrays the ratio of the observed support to the expected support under the independence of $X$ and $Y$. $$ \mathrm{lift}(X\Rightarrow Y) = \frac{ \mathrm{supp}(X \cup Y)}{ \mathrm{supp}(X) \times \mathrm{supp}(Y) } $$ For the ease of interpretation, the analysis was restricted by the possibility of a rule to blend only 3 items, and rules where hazards appeared on the right side of it (i.e., when co-benefits or preventive actions were estimated to induce hazards) were removed. The interactive plot below displays the TOP-50 rules sorted by the confidence and lift constraining measures. One can observe that, for instance, flash/surface flood faced by tree planting and/or creation of green space can reduce GHG emission (hazards – actions – co-benefits type of bonds). Another path example of the same structure: storm and wind cyclone (hurricane/typhoon) tackled by the real-time risk monitoring results in disaster preparedness. We can visually discern other link types such as action-co-benefits-co-benefits: flood defences – development and operation & storage together with economic growth can enhance disaster preparedness. The picture is interactive, and an interested reader can play with it. Apparently, the information depicted in the plot is quite rich and can serve as additional guidance for city leaders at the policy formulation stage. ```{r} nodes$value <- ifelse(is.na(nodes$value), 2, 1) # VisNetwork visNetwork( nodes, edges = ig_df$edges ) %>% visEdges(arrows = "to") %>% # visNodes(color = node_colors) %>% visEvents(type = "once", startStabilizing = "function() { this.moveTo({scale:0.6})}") %>% visPhysics(stabilization = FALSE) %>% visOptions(highlightNearest = T) ``` # Cooperation between Cities and Corporations ```{r message=FALSE, warning=FALSE, include=FALSE} # --- CITIES # --- Load the 2020 data # Set the working directory df_cities.20 <- read.csv('/kaggle/input/cdp-unlocking-climate-solutions/Cities/Cities Responses/2020_Full_Cities_Dataset.csv', stringsAsFactors=FALSE) %>% dplyr::select('id' = Account.Number, 'org' = Organization, 'cntry' = Country, 'region' = CDP.Region, 'qstn' = Question.Number, 'qstname' = Question.Name, 'coln' = Column.Number, 'colname' = Column.Name, 'rown' = Row.Number, 'rowname' = Row.Name, 'resp' = Response.Answer) # Recode to UTF-8 temp_vec <- df_cities.20$resp Encoding(temp_vec) <- "UTF-8" df_cities.20$resp <- temp_vec temp_vec <- df_cities.20$colname Encoding(temp_vec) <- "UTF-8" df_cities.20$colname <- temp_vec # Recode 'Other' category df_cities.20$resp[grepl("^Other, please specify", as.character(df_cities.20$resp))] <- "Other" # --- CORPORATIONS # --- Load the 2020 data df_corps.20 <- read.csv('/kaggle/input/cdp-unlocking-climate-solutions/Corporations/Corporations Responses/Climate Change/2020_Full_Climate_Change_Dataset.csv', stringsAsFactors=FALSE) %>% dplyr::select('id' = account_number, 'org' = organization, 'qstn' = question_number, 'qstname' = question_unique_reference, 'coln' = column_number, 'colname' = column_name, 'rown' = row_number, 'rowname' = row_name, 'resp' = response_value) # Recode to UTF-8 temp_vec <- df_corps.20$resp Encoding(temp_vec) <- "UTF-8" df_corps.20$resp <- temp_vec temp_vec <- df_corps.20$colname Encoding(temp_vec) <- "UTF-8" df_corps.20$colname <- temp_vec # Recode 'Other' category df_corps.20$resp[grepl("^Other, please specify", as.character(df_corps.20$resp))] <- "Other" ``` ```{r message=FALSE, warning=FALSE, include=FALSE} ############## Q6.0 : NLP #################################################################### ### --- 2020 # df_q6_0 <- df_cities.20 %>% filter(qstn == "6.0") # df_q6_0.wide <- dcast(df_q6_0, id + rown ~ coln, value.var="resp") # df_q6_0.wide <- df_q6_0.wide %>% # rename("opportunity"="1", "text"="2") %>% # filter(text!="") # # # Translation # # Language of the city # city_lang <- df_cities.20 %>% # filter(qstname == "What language are you submitting your response in?") %>% # dplyr::select(id, resp) %>% mutate(resp=dplyr::recode(resp, "Chinese"='zh-CN', "English"="en", "French"="fr", # "Korean"="ko", "Other (only applies to free text): italiano"="it", # "Portuguese"="pt", "Spanish"="es")) %>% # rename(lang=resp) # # df_q6_0.wide <- df_q6_0.wide %>% left_join(city_lang) # df_q6_0.wide.en <- df_q6_0.wide %>% filter(lang=="en") # df_q6_0.wide.nonen <- df_q6_0.wide %>% filter(lang!="en") # # df_list <- list() # unique_langs <- unique(df_q6_0.wide.nonen$lang) # for (i in (1:length(unique_langs))){ # # df_temp <- df_q6_0.wide.nonen %>% filter(lang==unique_langs[i]) # df_temp <- translate(dataset = df_temp, # content.field = 'text', # google.api.key = 'AIzaSyCyADH7ABEw5_KB8e_pG2e9Qc6Aa720Jmg', # source.lang = unique_langs[i], # target.lang = 'en') # df_list[[i]] <- df_temp # print(unique_langs[i]) # # } # # df_q6_0.wide.nonen <- do.call(rbind, df_list) # df_q6_0.wide.nonen <- df_q6_0.wide.nonen %>% # dplyr::select(-c(text)) %>% # rename(text=translatedContent) # df_q6_0.wide <- rbind(df_q6_0.wide.nonen, df_q6_0.wide.en) # df_q6_0.wide <- df_q6_0.wide %>% dplyr::select(-lang) ############## Q2.4a: Opportunities ##################################################################### # Potenial connection to question Q3.2a in Cities. ### Preprocessing - 2020 # df_corp_q2_4a <- df_corps.20 %>% filter(qstn == "C2.4a") # df_corp_q2_4a.wide <- dcast(df_corp_q2_4a, id + rown ~ coln, value.var="resp") # df_corp_q2_4a.wide <- df_corp_q2_4a.wide %>% # rename("opp"="1", "value_chain"="2", "type"="3", "climate_driver"="4", # "finicial_impact"="5", "specific_descr"="6", "horizon"="7", "likehood"="8", # "magnitude_impact"="9", "provide_financial_figure"="10", # "financial_figure"="11", "min_financial_figure"="12", # "max_financial_figure"="13", "financial_figure_descr"="14", # "realize_cost"="15", "startegy_descr"="16", "comment"="17") %>% # mutate(total_descr=paste(specific_descr, financial_figure_descr, startegy_descr, comment)) %>% # filter(total_descr!=" ",) ############## Combined Prerpocessing ##################################################################### # Combine datasets # df_corp_q2_4a.wide.20 <- df_corp_q2_4a.wide %>% # dplyr::select(total_descr) %>% rename(text=total_descr) %>% # mutate(type="corporation", year="20") # df_q6_0.wide.20 <- df_q6_0.wide %>% # dplyr::select(text) %>% # mutate(type="city", year="20") # # # Combine in one dataset # df_corp_city_text <- rbind(df_corp_q2_4a.wide.20, # df_q6_0.wide.20) # # # NLP preprocessing # # Initilize spacy model # cnlp_download_spacy("en") # cnlp_init_spacy("en_core_web_sm") # # # Annotate the text # annotation <- cnlp_annotate(input = tolower(df_corp_city_text$text), verbose = 10) # df_corp_city_text$doc_id <- annotation$document$doc_id # # # Remove stopwords and non-alphabetical tokens # stopwords_vec <- stopwords::stopwords(language = "en",source = "smart") # stopwords_vec <- c(stopwords_vec, "e.g.", "la") # df_text_preprocessing <- annotation$token %>% # filter(!(lemma %in% stopwords_vec), # !(upos %in% c("DET", "PUNCT", "PROPN")), # (grepl("^[A-Za-z]+$", lemma, perl = T))) %>% # dplyr::select(doc_id, lemma) %>% # group_by(doc_id) %>% # summarise_at(vars(lemma), funs(paste(., collapse = ' '))) # df_corp_city_text <- df_corp_city_text %>% right_join(df_text_preprocessing) # # # saveRDS(df_corp_city_text, "df_corp_city_text") # LOAD PREPROCESSED FILE df_corp_city_text <- readRDS("/kaggle/input/nlpfiles/df_corp_city_text") ### Preprocessing with quanteda corpus_corp_city <- corpus(df_corp_city_text$lemma, docvars=df_corp_city_text %>% dplyr::select(-lemma)) dfmat_corp_city <- tokens(corpus_corp_city) %>% tokens_ngrams(n = c(1, 2)) %>% dfm() %>% dfm_trim(min_termfreq = 20) # Remove documents without vocabulary df_corp_city_text_temp <- df_corp_city_text[rowSums(dfmat_corp_city) != 0,] dfmat_corp_city <- dfmat_corp_city[rowSums(dfmat_corp_city) != 0,] # ~~ Model 2, K=20 K <- 20 stm_model.2 <- stm(documents = dfmat_corp_city, prevalence =~ type, data = dfmat_corp_city@docvars, K = K, max.em.its = 75, init.type = "Spectral") # Extract topic labels df_topiclabels <- data.frame(labelTopics(stm_model.2, n = 5)$frex) df_topiclabels <- df_topiclabels %>% mutate(topic=as.character(1:nrow(df_topiclabels)), topic_label=paste0("Topic ", topic, ": ", paste(str_to_title(X1), X2, X3, X4, X5, sep = ", "))) %>% dplyr::select(topic_label, topic) # Extract topics proportion proportion <- as.data.frame(colSums(stm_model.2$theta/nrow(stm_model.2$theta))) colnames(proportion) <- "topic_prop" df_topiclabels <- cbind(df_topiclabels, proportion) # Calculating the Effects prep.2 <- estimateEffect(1:K ~ type, stm_model.2, uncertainty = "Global", meta = dfmat_corp_city@docvars) # Extract values for plotting get_effects <- function(estimates, # estimates object variable, # variable for estimates type, # continuous or pointestimate ci = 0.95, # confidence interval moderator = NULL, # moderator for interaction modval = NULL, # cov values for difference cov_val1 = NULL, cov_val2 = NULL) { data <- plot.estimateEffect( x = estimates, covariate = variable, method = type, ci.level = ci, moderator = moderator, moderator.value = modval, cov.value1 = cov_val1, cov.value2 = cov_val2, omit.plot = TRUE ) # catching inconsistent stm naming conventions if (!'cis' %in% names(data)) { data$cis <- data$ci } if ('x' %in% names(data)) { data$uvals <- data$x } names(data$cis) <- data$topics names(data$means) <- data$topics tidy_stm <- data$topics %>% purrr::map(function(top) { top <- as.character(top) if (type == 'difference') { props <- tibble( difference = data$means[[top]], topic = top, lower = data$cis[[top]][[1]], upper = data$cis[[top]][[2]] ) } else { cis <- t(data$cis[[top]]) %>% as_tibble() %>% purrr::set_names(c('lower', 'upper')) props <- tibble( value = data$uvals, proportion = data$means[[top]], topic = top ) %>% bind_cols(cis) } }) tidy_stm <- tidy_stm %>% bind_rows() tidy_stm$topic <- as.factor(tidy_stm$topic) if (type == 'pointestimate') { tidy_stm$value <- as.factor(tidy_stm$value) } if (!is.null(moderator)) { tidy_stm$moderator <- modval } return(tidy_stm) } df_effects <- get_effects(prep.2, variable = "type", type = "difference", cov_val1 = "corporation", cov_val2 = "city") df_effects$topic <- as.character(sort(as.numeric(df_effects$topic), decreasing = T)) df_effects <- df_effects %>% left_join(df_topiclabels) ``` Our next task is to determine an array of practical and actionable points that can bring corporate and city ambition together. Cooperation between cities and organisations is vital since the mutual aid it delivers can lead to better attainments and savings in spending or efforts, which can be ploughed back to further facilitate environmental leadership. To disclose the overlap between the city and corporate aspirations we focused on 2020 responses about climate-related opportunities a city or a company has identified. Specifically, we selected a column titled "Describe how the city is maximizing this opportunity" from question 6.0 in cities' questionnaire. Besides, from the corporations' questionnaire, the following group of columns were considered: "Company-specific description", "Explanation of financial impact figure", "Strategy to realize opportunity and explanation of cost calculation" and "Comment". The latter text fields were consolidated into a single string feature called "Total description". Responses in non-English languages were translated with the "translateR" package in R, which uses Google API. We brought these two sources of data together to obtain an aggregated text corpus of environmental opportunities indicated by city and corporate leaders. To ascertain the intersecting topics of opportunities we leveraged the Structural Topic Modelling (STM) (Roberts, Stewart & Tingley 2019). Unlike the usual Topic Models, STM allows incorporating additional document-level variables into the estimation procedure. In our work, we plugged into the analysis a dichotomous categorical characteristic pointing on whether a respondent was a city or corporate representative. The final specification was made to embrace 20 different topics as according to several goodness-of-fit criteria (e.g., coherence and perplexity) this number enabled to keep an optimal match between the quality and quantity of the solution defined. For exploratory purposes, let us look at how the whole dataset of topics looks like. The histogram below exhibits the occurrence probability distribution for the topics revealed. It is noticeable that in general the most widespread topics are referred to urban aspects, energy efficiency or renewable energy issues, and financial points of the climate change. ```{r message=FALSE, warning=FALSE, fig.height=7, fig.width=15} # Plot topic proportions ggplot(df_topiclabels, aes(x = reorder(topic_label, topic_prop), y = topic_prop)) + geom_bar(stat = "identity", fill='deepskyblue3') + coord_flip() + theme_minimal() + scale_y_continuous(breaks = c(0, 0.05), limits = c(0, 0.16), expand = c(0, 0)) + #change breaks and limits as you need coord_flip() + geom_text(aes(label = scales::percent(topic_prop, 2)), #Scale in percent hjust = -0.25, size = 4, position = position_dodge(width = 1), inherit.aes = TRUE) + theme_classic() + theme(axis.title=element_text(size=14), legend.text=element_text(size=14), # legend.position = c(0.9, 0.1), legend.title=element_text(size=14), text = element_text(size=14), plot.title = element_text(size=22, hjust = 0.5)) + ggtitle("Proportion of Topics Related to Opportunities") + labs(x="", y="Topic Probability ") ``` To delineate the overlap between climate opportunities for the two entities under consideration, the respondent category (cities vs. corporations) effect on the probability of observing a topic in the document was evaluated. The results are visualised below. The x-axis of the graph represents the extent of "belonging" to corporations' responses (left side) or cities' responses (right side), y-axis marks the concrete topics. ```{r fig.height=10, fig.width=15} # Plot ggplot(df_effects, aes(x = reorder(topic_label, difference), y = difference)) + geom_point(data=df_effects, size=2) + ylim(-0.4, 0.4) + geom_errorbar(aes(ymin = lower, ymax = upper), width=0.25) + labs(x = "", y = "Corporation < ................... > City") + geom_hline(aes(yintercept=0), linetype="dashed", color = "gray40") + coord_flip() + theme_minimal() + theme(axis.title=element_text(size=14), legend.text=element_text(size=14), # legend.position = c(0.9, 0.1), legend.title=element_text(size=14), text = element_text(size=14), plot.title = element_text(size=22, hjust = 0.5)) + ggtitle("Difference in the Topics between Cities and Corporations") ``` There are 3 topics that cross the vertical line originating from zero, meaning that these topics are peculiar to cities and corporations in a relatively equal manner. Particularly, these intersecting responses include areas associated with the transportation system (topic 16), technological innovations (topic 3), and extreme climate events (topic 18). The topics that are potentially close to both cities and corporations relate to water consumption (topic 9) and energy efficiency issues (topic 4). Therefore, we can contend that the uncovered topics embody critical areas where cities’ and corporations’ intentions should harmonise above all to achieve bilateral benefits and foster mutual reciprocity in combating environmental threats. The next picture depicts the wordclouds for the three overlapping topics in greater detail. ```{r fig.height=8, fig.width=14} set.seed(42) par(mfrow=c(1,3)) stm::cloud(stm_model.2, topic = 18, scale = c(5,.5), color="coral3") stm::cloud(stm_model.2, topic = 16, scale = c(5,.5), color="deepskyblue3") stm::cloud(stm_model.2, topic = 3, scale = c(5,.5), color="forestgreen") mtext("Overlapping Opportunities between Cities and Corporations", side = 3, line = -9, outer = TRUE, cex=3) ``` # Towards Better Disclosure Practices ```{r message=FALSE, warning=FALSE, include=FALSE} ############## PARSING - Q2.0b - 2020 ######################################################################### # >>>>>>>>>> WE WILL UPLOAD PREPROCESSED FILES, # ALL THE REPORTS IN .pdf FORMAT ARE AVAILABLE BY THIS LINK: # https://drive.google.com/drive/folders/1zaD7mS1g421uAj3Y1oxM9P-C7n6qpTpn?usp=sharing ### --------- Q2.0b # df_q2_0b <- df_cities.20 %>% filter(qstn == "2.0b") # df_q2_0b.link <- df_q2_0b %>% # filter(colname == "Web link", # !(resp %in% c("Question not applicable", ""))) %>% # dplyr::select(id, resp) # # Extract the links # temp_vec <- sapply(df_q2_0b.link$resp, function(x) str_extract_all(x, "(?<=^|\\s)http[^\\s]+")) # temp_vec <- as.data.frame(sapply(temp_vec, function(x) paste(x, collapse = ' '))) # df_q2_0b.link$resp <- temp_vec[,1] # df_q2_0b.link <- df_q2_0b.link %>% filter(resp!="") # df_q2_0b.link <- df_q2_0b.link %>% separate(resp, paste0("link", as.character(1:9)), # sep=" ") # df_q2_0b.link <- df_q2_0b.link[,colSums(is.na(df_q2_0b.link)) != nrow(df_q2_0b.link)] # # Transform to long-format: id-link as observation # df_q2_0b.link.long <- melt(df_q2_0b.link, id.vars=c("id")) %>% # filter(!is.na(value)) %>% rename(link_num=variable, link=value) %>% # mutate(link_num=substr(link_num, 5,5), # id_link_num=paste(id, link_num, sep="_")) %>% # slice(-416) # Remove city with a large book as a report # # ### Download reports # setwd("WD") # <<< SET WORKING DIRECTORY # oldw <- getOption("warn") # options(warn = -1) # for (i in (1:nrow(df_q2_0b.link.long))) { # tryCatch(download.file(df_q2_0b.link.long$link[i], # destfile = paste0(df_q2_0b.link.long$id_link_num[i], ".pdf"), # mode = "wb", quiet = T), # error = function(e) print(paste(as.character(i), ' error ~~~~~'))) # print(i) # Sys.sleep(1) # } # options(warn = oldw) ############## PARSING - Q2.0b - 2019 ################################################################ ### --------- Q2.0b # df_q2_0b <- df_cities.19 %>% filter(qstn == "2.0b") # df_q2_0b.link <- df_q2_0b %>% # filter(colname == "Web link", # !(resp %in% c("Question not applicable", ""))) %>% # dplyr::select(id, resp) # # Extract the links # temp_vec <- sapply(df_q2_0b.link$resp, function(x) str_extract_all(x, "(?<=^|\\s)http[^\\s]+")) # temp_vec <- as.data.frame(sapply(temp_vec, function(x) paste(x, collapse = ' '))) # df_q2_0b.link$resp <- temp_vec[,1] # df_q2_0b.link <- df_q2_0b.link %>% filter(resp!="") # df_q2_0b.link <- df_q2_0b.link %>% separate(resp, paste0("link", as.character(1:9)), # sep=" ") # df_q2_0b.link <- df_q2_0b.link[,colSums(is.na(df_q2_0b.link)) != nrow(df_q2_0b.link)] # # Transform to long-format: id-link as observation # df_q2_0b.link.long <- melt(df_q2_0b.link, id.vars=c("id")) %>% # filter(!is.na(value)) %>% rename(link_num=variable, link=value) %>% # mutate(link_num=substr(link_num, 5,5), # id_link_num=paste(id, link_num, sep="_")) # # ## Download reports # setwd("WD") # <<< SET WORKING DIRECTORY # oldw <- getOption("warn") # options(warn = -1) # for (i in (1:nrow(df_q2_0b.link.long))) { # tryCatch(download.file(df_q2_0b.link.long$link[i], # destfile = paste0(df_q2_0b.link.long$id_link_num[i], ".pdf"), # mode = "wb", quiet = T), # error = function(e) print(paste(as.character(i), ' error ~~~~~'))) # print(i) # Sys.sleep(0.5) # } # options(warn = oldw) ############## >>>>> # LOAD Q2.0b 2019, 2020 REPORTS # ### Load the data from separate pdfs ------------- 2019 # setwd("WD") <-------------- FOLDER WITH 2019 REPORTS # pdf_files <- list.files(pattern = "pdf$") # oldw <- getOption("warn") # options(warn = -1) # pdf_texts <- list() # for (i in (1:length(pdf_files))) { # tryCatch(pdf_texts[[i]] <- pdf_text(pdf_files[i]), # error = function(e) print(paste(as.character(i), ' error ~~~~~'))) # print(i) # # Sys.sleep(1) # } # pdf_files <- pdf_files[!unlist(lapply(pdf_texts, is.null))] # pdf_texts <- pdf_texts[!unlist(lapply(pdf_texts, is.null))] # # df_texts <- as.data.frame(sapply(pdf_texts, function(x) paste(x, collapse = '---PAGE-BREAK---'))) # ---PAGE-BREAK--- # colnames(df_texts) <- "text" # df_texts$id_link_num <- substr(pdf_files, 1, (nchar(pdf_files)-4)) # # # Select only reports that are written in English # df_texts$lang <- detect_language(df_texts$text) # df_texts.en.2019 <- df_texts %>% filter(lang=="en") %>% dplyr::select(-c(lang)) # df_texts.en.2019$year <- "2019" # # ### Load the data from separate pdfs ------------- 2020 # setwd("WD") <-------------- FOLDER WITH 2019 REPORTS # pdf_files <- list.files(pattern = "pdf$") # oldw <- getOption("warn") # options(warn = -1) # pdf_texts <- list() # for (i in (1:length(pdf_files))) { # tryCatch(pdf_texts[[i]] <- pdf_text(pdf_files[i]), # error = function(e) print(paste(as.character(i), ' error ~~~~~'))) # print(i) # # Sys.sleep(1) # } # pdf_files <- pdf_files[!unlist(lapply(pdf_texts, is.null))] # pdf_texts <- pdf_texts[!unlist(lapply(pdf_texts, is.null))] # # # df_texts <- as.data.frame(sapply(pdf_texts, function(x) paste(x, collapse = '---PAGE-BREAK---'))) # ---PAGE-BREAK--- # colnames(df_texts) <- "text" # df_texts$id_link_num <- substr(pdf_files, 1, (nchar(pdf_files)-4)) # # # Select only reports that are written in English # df_texts$lang <- detect_language(df_texts$text) # df_texts.en.2020 <- df_texts %>% filter(lang=="en") %>% dplyr::select(-c(lang)) # df_texts.en.2020$year <- "2020" # # # # Merge together ------------- 2020 & 2019 # df_texts.en <- rbind(df_texts.en.2019, df_texts.en.2020) # df_texts.en <- df_texts.en %>% distinct(text, .keep_all = T) # df_texts.en <- df_texts.en %>% # mutate(id_link_num=substr(id_link_num, 1, (nchar(id_link_num)-2)), # id=paste(id_link_num, year, sep = "_")) %>% # dplyr::select(-c(id_link_num)) # # Save to RDS # setwd("C:/Users/Artem/YandexDisk/CDP/data/nlp_files") # # saveRDS(df_texts.en, "df_texts_en_2020_2019") # # # TO LONL # df_texts.en <- df_texts.en %>% separate(text, paste0("page", as.character(1:500)), sep="---PAGE-BREAK---") # # Transform to long-format: id-link as observation # df_texts.en.long <- melt(df_texts.en, id.vars=c("id")) %>% # filter(!is.na(value)) %>% rename(page=variable, text=value) %>% # mutate(page=as.character(page), # page=substr(page, 5, nchar(page)), # # text=tolower(text), # CHANGE FOR NER <<<<<<<<<<<<<<<<<<<<<<<<<<<<<< # ) # # ### NLP part # # Initilize spacy model # cnlp_init_spacy("en_core_web_sm") # # Annotate the text # annotation <- cnlp_annotate(input = df_texts.en.long, verbose = 100) # # Add City ID # df_temp <- annotation$token # df_temp <- df_temp %>% left_join(annotation$document %>% dplyr::select(id, doc_id)) # annotation$token <- df_temp # # Save the preprocessed file # setwd("WD") # saveRDS(annotation, "annotation_q2_0b_2019_2020") # annotation <- readRDS("annotation_q2_0b_2019_2020") # # # Remove stopwords and non-alphabetical tokens # stopwords_vec <- stopwords::stopwords(language = "en",source = "smart") # stopwords_vec <- c(stopwords_vec, "e.g.", "la") # df_text_preprocessing <- annotation$token %>% # filter(!(lemma %in% stopwords_vec), # !(upos %in% c("DET", "PUNCT", "PROPN")), # (grepl("^[A-Za-z]+$", lemma, perl = T))) %>% # dplyr::select(doc_id, id, lemma) %>% # # group_by(id) %>% # <<<<<<<<<<<<<<<<<<<<<<<<< BY CITY # group_by(doc_id) %>% # <<<<<<<<<<<<<<<<<<<<<<<<< BY PAGE # summarise_at(vars(lemma, id), funs(paste(., collapse = ' '))) # add id to vars # df_text_preprocessing$id <- unlist(lapply(df_text_preprocessing$id, function(x) str_split(x, " ")[[1]][1])) # # # saveRDS(df_text_preprocessing, "df_text_preprocessing_q2_0b_2019_2020_obspages") ### ---- WORD CLOUD # Load the preprocessed file df_text_preprocessing <- readRDS("/kaggle/input/nlpfiles/df_text_preprocessing_q2_0b_2019_2020_obspages") # Create Corpus and Document term matrix corpus_q2_0b <- corpus(df_text_preprocessing$lemma, docvars=df_text_preprocessing %>% dplyr::select(-c(lemma))) dfm_q2_0b <- tokens(corpus_q2_0b) %>% tokens_ngrams(n = c(1,2)) %>% dfm() %>% dfm_trim(min_termfreq = 30) ``` In our last section, we raise a problem about cities’ benefits from being transparent regarding the goal setting, actions planning and undertaking. Cities will gain substantially not only from participating in surveys promoting disclosure on climate change activity, but also from providing freely downloadable, understandable, and relevant reports to a greater audience. Particularly, we seek to figure out how clear and up-to-date cities are in their reports. We parsed cities' reports on climate change risk and vulnerability assessment based on the weblinks provided in question 2.0b in 2019 and 2020. Only publicly available .pdf documents were downloaded. In total, 95 reports were retrieved and further processed. These documents seem to entail an added value which may not be captured by the CDP questionnaire. First, let us inspect a word cloud with the most frequent words in the reports. ```{r message=FALSE, warning=FALSE, fig.height=5, fig.width=5} set.seed(42) # Plot World Cloud textplot_wordcloud(dfm_q2_0b) ``` ```{r message=FALSE, warning=FALSE, include=FALSE} ### ---- READABILITY # Load from RDS df_texts.en <- readRDS("/kaggle/input/nlpfiles/df_texts_en_2020_2019") # Load additional data df_id_name.20 <- df_cities.20 %>% dplyr::select(id, org) %>% distinct(id, .keep_all = T) %>% mutate(id=paste0(id, "_2020")) df_id_name.19 <- df_cities.19 %>% dplyr::select(id, org) %>% distinct(id, .keep_all = T) %>% mutate(id=paste0(id, "_2019")) df_id_name <- rbind(df_id_name.20, df_id_name.19) # Fix the document temp_string <- df_texts.en$text[81] temp_string <- strsplit(temp_string, "Appendix A:")[[1]][2] df_texts.en[81,]$text <- temp_string # Calculate the Readability score readability_index <- textstat_readability(df_texts.en$text) readability_index$id <- df_texts.en$id readability_index <- readability_index %>% left_join(df_id_name) readability_index <- readability_index %>% mutate(year=substr(id, (nchar(id)-3), nchar(id)), year=paste0("(", year, ")"), org=paste(org, year, sep = " ")) ``` The provision of the publicly available reports constitutes an important step towards environmental leadership; however, mayors need to make sure these reports are comprehensible or readable for the general population. The “readability” of a text is not an easily measurable concept, and there exist a lot of yardsticks to approach such a problem such as a Flesch-Kincaid readability test. The latter was developed to define how complicated a passage of English is to grasp. The formula for this reading-ease score is denoted by: $$206.835 - 1.015 \left( \frac{\text{total words}}{\text{total sentences}} \right) - 84.6 \left( \frac{\text{total syllables}}{\text{total words}} \right)$$ The distribution of cities' reports based on the Flesch–Kincaid test is displayed below. One can see that these metrics vary from 0 to almost 50 for the sample of cities under consideration. According to the conventional thresholds, texts scoring lower than 10 are extremely difficult to read. By comparison, Times magazine scores roughly 52, which corresponds to a fairly difficult-to-read text. In our analysis, better reporting practices were achieved by Park, Petaling Jaya, Miami, and Emeryville. The rest of the cities may need to catch up with the latter since improved readability of reports is a fruitful way to increase people's awareness about climate change hazards and boost more responsible behaviour on the part of both government and population. ```{r message=FALSE, warning=FALSE, fig.height=12, fig.width=10} # Plot ggplot(readability_index, aes(x = reorder(org, Flesch), y = Flesch)) + geom_point(size=2, color=adjustcolor("deepskyblue3", alpha.f = 0.8)) + coord_flip() + theme_minimal() + theme(axis.title=element_text(size=14), legend.text=element_text(size=14), # legend.position = c(0.9, 0.1), title=element_text(size=18), text = element_text(size=8)) + ggtitle("Readability of Reports") + labs(x="", "Flesch Index of Readability") ``` ```{r message=FALSE, warning=FALSE, include=FALSE} ### --- STM with pages as observations ---------------------------------------------------------------------- # Save the preprocessed file df_text_preprocessing <- readRDS("/kaggle/input/nlpfiles/df_text_preprocessing_q2_0b_2019_2020_obspages") # Add categorical variables to text # 2020 df_q2_0b.20 <- df_cities.20 %>% filter(qstn == "2.0b") df_q2_0b.vuln.20 <- df_q2_0b.20 %>% dplyr::filter(coln == 7, rown==1, resp != "") %>% dplyr::select(id, resp) %>% mutate(id=paste0(as.character(id), "_2020")) %>% rename(vulnerable_population=resp) df_q2_0b.year.20 <- df_q2_0b.20 %>% dplyr::filter(coln == 3, rown==1, resp != "") %>% dplyr::select(id, resp) %>% mutate(id=paste0(as.character(id), "_2020")) %>% rename(adoption_year=resp) df_q2_0b.author.20 <- df_q2_0b.20 %>% dplyr::filter(coln == 6, rown==1, resp != "") %>% dplyr::select(id, resp) %>% mutate(id=paste0(as.character(id), "_2020")) %>% rename(author=resp) # 2019 df_q2_0b.19 <- df_cities.19 %>% filter(qstn == "2.0b") df_q2_0b.vuln.19 <- df_q2_0b.19 %>% dplyr::filter(coln == 8, rown==1, resp != "") %>% dplyr::select(id, resp) %>% mutate(id=paste0(as.character(id), "_2019")) %>% rename(vulnerable_population=resp) df_q2_0b.year.19 <- df_q2_0b.19 %>% dplyr::filter(coln == 2, rown==1, resp != "") %>% dplyr::select(id, resp) %>% mutate(id=paste0(as.character(id), "_2019")) %>% rename(adoption_year=resp) df_q2_0b.author.19 <- df_q2_0b.19 %>% dplyr::filter(coln == 7, rown==1, resp != "") %>% dplyr::select(id, resp) %>% mutate(id=paste0(as.character(id), "_2019")) %>% rename(author=resp) # Merge df_q2_0b.year <- rbind(df_q2_0b.year.19, df_q2_0b.year.20) df_q2_0b.vuln <- rbind(df_q2_0b.vuln.19, df_q2_0b.vuln.20) df_q2_0b.author <- rbind(df_q2_0b.author.19, df_q2_0b.author.20) # Add to the dataset df_text_preprocessing <- df_text_preprocessing %>% data.frame() %>% left_join(df_q2_0b.year, by="id") %>% left_join(df_q2_0b.vuln, by="id") %>% left_join(df_q2_0b.author, by="id") # adoption_year preprocessing df_text_preprocessing <- df_text_preprocessing %>% filter(!is.na(adoption_year)) df_text_preprocessing <- df_text_preprocessing %>% mutate(adoption_year=ifelse(adoption_year %in% c("2016", "2017", "2018", "2019"), "2016-2019", "2015 and earlier")) # Remove pages with low number of characters df_text_preprocessing$len <- unlist(lapply(df_text_preprocessing$lemma, nchar)) df_text_preprocessing <- df_text_preprocessing %>% filter(len>100) # Create Corpus and Document term matrix corpus_q2_0b <- corpus(df_text_preprocessing$lemma, docvars=df_text_preprocessing %>% dplyr::select(-c(lemma))) dfm_q2_0b <- tokens(corpus_q2_0b) %>% tokens_ngrams(n = c(1,2)) %>% dfm() %>% dfm_trim(min_termfreq = 30) ncol(dfm_q2_0b) nrow(dfm_q2_0b) # Model 1. Base + adoption_year, K=20 K <- 20 # stm_model.1 <- stm(documents = dfm_q2_0b, # prevalence =~ adoption_year, # K = K, max.em.its = 75, init.type = "Spectral") # saveRDS(stm_model.1, "stm_model_pages_20") stm_model.1 <- readRDS("/kaggle/input/nlpfiles/stm_model_pages_20") # Estimate the effect of adoption_year on topics prep.1 <- estimateEffect(1:K ~ adoption_year, stm_model.1, uncertainty = "Global", meta = dfm_q2_0b@docvars) # Extract topic labels df_topiclabels <- data.frame(labelTopics(stm_model.1, n = 5)$frex) df_topiclabels <- df_topiclabels %>% mutate(topic=as.character(1:nrow(df_topiclabels)), topic_label=paste0("Topic ", topic, ": ", paste(str_to_title(X1), X2, X3, X4, X5, sep = ", "))) %>% dplyr::select(topic_label, topic) # Extract topics proportion proportion <- as.data.frame(colSums(stm_model.1$theta/nrow(stm_model.1$theta))) colnames(proportion) <- "topic_prop" df_topiclabels <- cbind(df_topiclabels, proportion) # Extract values for plotting get_effects <- function(estimates, # estimates object variable, # variable for estimates type, # continuous or pointestimate ci = 0.95, # confidence interval moderator = NULL, # moderator for interaction modval = NULL, # cov values for difference cov_val1 = NULL, cov_val2 = NULL) { data <- plot.estimateEffect( x = estimates, covariate = variable, method = type, ci.level = ci, moderator = moderator, moderator.value = modval, cov.value1 = cov_val1, cov.value2 = cov_val2, omit.plot = TRUE ) # catching inconsistent stm naming conventions if (!'cis' %in% names(data)) { data$cis <- data$ci } if ('x' %in% names(data)) { data$uvals <- data$x } names(data$cis) <- data$topics names(data$means) <- data$topics tidy_stm <- data$topics %>% purrr::map(function(top) { top <- as.character(top) if (type == 'difference') { props <- tibble( difference = data$means[[top]], topic = top, lower = data$cis[[top]][[1]], upper = data$cis[[top]][[2]] ) } else { cis <- t(data$cis[[top]]) %>% as_tibble() %>% purrr::set_names(c('lower', 'upper')) props <- tibble( value = data$uvals, proportion = data$means[[top]], topic = top ) %>% bind_cols(cis) } }) tidy_stm <- tidy_stm %>% bind_rows() tidy_stm$topic <- as.factor(tidy_stm$topic) if (type == 'pointestimate') { tidy_stm$value <- as.factor(tidy_stm$value) } if (!is.null(moderator)) { tidy_stm$moderator <- modval } return(tidy_stm) } df_effects <- get_effects(prep.1, variable = "adoption_year", type = "difference", cov_val1 = "2015 and earlier", cov_val2 = "2016-2019") df_effects$topic <- as.character(sort(as.numeric(df_effects$topic), decreasing = T)) df_effects <- df_effects %>% left_join(df_topiclabels) ``` Further, Structural Topic Modelling was implemented. Like previously, let us visually evaluate the occurrence probability distribution for the topics detected. Notably, the topics about funding and finance, risk and vulnerability assessment, temperature increase, and the effect of climate change are assessed as the most popular ones. ```{r message=FALSE, warning=FALSE, fig.height=7, fig.width=10} # Plot topic proportions ggplot(df_topiclabels, aes(x = reorder(topic_label, topic_prop), y = topic_prop)) + geom_bar(stat = "identity", fill='forestgreen') + coord_flip() + theme_minimal() + scale_y_continuous(breaks = c(0, 0.3), limits = c(0, 0.12), expand = c(0, 0)) + #change breaks and limits as you need coord_flip() + geom_text(aes(label = scales::percent(topic_prop, 2)), #Scale in percent hjust = -0.25, size = 4, position = position_dodge(width = 1), inherit.aes = TRUE) + theme_classic() + theme(axis.title=element_text(size=12), legend.text=element_text(size=12), # legend.position = c(0.9, 0.1), legend.title=element_text(size=12), text = element_text(size=12), plot.title = element_text(size=24, hjust = 5)) + ggtitle("Proportion of Topics in the Reports") + labs(x="", y="Topic Probability ") ``` Next, we added a categorical variable describing the reports’ year of publication or approval from the local government. The variable was binarized into two ranges: 2016-2019 and 2015 and earlier. Therefore, a structural topic model with 20 topics plus a year indicator was evaluated. ```{r message=FALSE, warning=FALSE, fig.height=7, fig.width=15} # Plot ggplot(df_effects, aes(x = reorder(topic_label, difference), y = difference)) + geom_point(data=df_effects, size=2) + ylim(-0.07, 0.07) + geom_errorbar(aes(ymin = lower, ymax = upper), width=0.25) + labs(x = "", y = "2015 and earlier < ................... > 2016-2019") + geom_hline(aes(yintercept=0), linetype="dashed", color = "gray40") + coord_flip() + theme_minimal() + theme_minimal() + theme(axis.title=element_text(size=12), legend.text=element_text(size=12), # legend.position = c(0.9, 0.1), legend.title=element_text(size=12), text = element_text(size=12), plot.title = element_text(size=20, hjust = 0.5)) + ggtitle("Difference in the Topics\nBased on the Reports' Year") ``` The plot above informs us that topics tending to 2016-2019 years may be more relevant and up-to-date and city leaders need to treat them properly in their environment-oriented activity to be able to effectively respond to the climatic threats of the modern world. These more recent topics outline a sea-level rise problem (topic 11), emphasise the necessity to account for a local level (topic 5), mention budget opportunities (topic 2), and refer to cities’ participation in global organisations (topic 4). Although a small sample size deprives us of determining far-reaching consequences coming from to such an analysis, we enabled to obtain a piece of evidence suggesting that cities vary in terms of how they keep up with the environmental trends, and this can be a potential future growth point. # Discussion and Conclusion To sum up, this piece of work was conducted to highlight practical points the can bright cities and corporations to environmental leadership and forward healthier life. We suggested a tactic to measure and assess cities' performance (KPIs) by preserving a socially equitable position in both the estimation algorithm and in the metrics deployed. We benchmarked cities based on their survey reports and external data about emissions abiding by equitable rules in the comparative algorithm, i.e. by accounting for cities' variation in their starting conditions (DGP per capita, population sizes, and temperature level). The efficiency exploration embraced the analysis of dynamic shifts that can indicate the underpinning drivers of such shifts: a city's performance or changes in the efficiency frontier (other cities' performance), or both. It was observed that generally, cities are lagging regarding indicators connected with assessment and planning, whereas in terms of actions they prove to be fairly efficient. Such a finding may symbolise poor internalisation of the environmental practices by cities’ policy systems: more liability in planning is needed because actions without careful arrangements may not be long-term or resilient and can effortlessly be stricken by economic slowdowns especially like in the current coronavirus era. Importantly, less wealthy cities are more vulnerable to social equity distortions in environmental responses than wealthier cities even controlling for the background characteristics; therefore, better treatment of vulnerable populations, undoubtfully, stands out as a priority for policy interventions and government investments. In this research, we also erected a network of cities’ peers which can be exploited as a practical tool for decision-makers. The network can instruct city officials in selecting real examples of more successful peers which should be followed while trying to accelerate the efficiency growth. Additionally, we showed what common environmental areas are estimated as essential by corporations and cities. The analysis provided evidence that opportunities in the transportation system, technological innovations, water consumption, and energy efficiency issues interest corporate and city leaders roughly to the same extent, meaning that these areas represent crucial points of prolific cooperation and should be approached through the concerted efforts. Finally, this examination unveiled possible ways to uplift transparency in cities' reporting strategies. We demonstrated that it may be necessary for cities to improve the readability of their reports to raise awareness of climate threats for broad masses of people and facilitate environment-friendly behaviour on the side of populations and corporations. Besides, it may be indispensable for cities to upgrade the dated reports with reasonable frequency to keep pace with the climate challenges of the 21st century. Overall, the current research seems to entail a practical value for government leaders proposing climate change policy initiatives as the analysis executed implements additional means for decision-making processes and can be used as a groundwork for a more scrutinised examination further. # References 1. Internet source: https://en.wikipedia.org/wiki/Flesch%E2%80%93Kincaid_readability_tests 2. Le Quéré, C., Jackson, R. B., Jones, M. W., Smith, A. J., Abernethy, S., Andrew, R. M., ... & Friedlingstein, P. (2020). Temporary reduction in daily global CO 2 emissions during the COVID-19 forced confinement. Nature Climate Change, 1-7. 3. Newburger, E. (2020). Coronavirus could weaken climate change action and hit clean energy investment, researchers warn. CNBC, March, 13. 4. Roberts, M. E., Stewart, B. M., & Tingley, D. (2019). Stm: An R package for structural topic models. Journal of Statistical Software, 91(1), 1-40. 5. Sickles, R. C., & Zelenyuk, V. (2019). Measurement of productivity and efficiency. Cambridge University Press. 6. World Bank. 2020. Doing Business 2020. Washington, DC: World Bank. DOI:10.1596/978-1-4648-1440-2. License: Creative Commons Attribution CC BY 3.0 IGO 7. Zhang, C., & Zhang, S. (2003). Association rule mining: models and algorithms (Vol. 2307). Springer.