# R-Code for Replication of Results and Figures in Möhring et al., 2025 # The code was run on R version 4.2.1. # ************************************* Replication code starts from here ************************************* ###### ### i) Read in data #### library(ggplot2) library(tidyverse) library(diagis) library(plotrix) s <- function(x){plotrix::std.error(x,na.rm=T)} m <- function(x){mean(x,na.rm=T)} mw <- function(x,w){diagis::weighted_mean(x,w,na.rm = T)} # m <- function(x,w){wtd.mean(x,w,normwt=T,na.rm = T)} sw <- function(x,w){weighted_se(x,w,na.rm = T)} # s <- function(x,w){sqrt(wtd.var(x,w,normwt=T,na.rm = T))} qw<- function(x,w){wtd.quantile(x,w,normwt = T,na.rm=T, probs=c(0.5))} path <- "path" datsur <- read.csv2(c(path,"/Mohring_etal_data_replication.csv")) ### iii) Replicate Fig. 1 and Fig. 2 #### chall <- datsur %>% mutate(row.n = row_number()) %>% select(FS1 = Q4.FS.2.1.sufficient.quantity, FS2 = Q4.FS.2.2.healthy.food, FS3 = Q4.FS.2.3.safe.food.and.feed, FS4 = Q4.FS.2.4.resilience.of.food, ENV1 = Q5.EE.2.1.pollution.drinking.water, ENV2 = Q5.EE.2.2.pollution.soil, ENV3 = Q5.EE.2.3.pollution.marine.ecosystems, ENV4 = Q5.EE.2.4.pollution.freshwater.ecosystems, ENV5 = Q5.EE.2.5.loss.of.genes.biodiversity, HH1 = Q6.HH.2.1.effects.farm.workers, HH2 = Q6.HH.2.2.effects.bystanders, HH3 = Q6.HH.2.3.effects.on.consumers, HH4= Q6.HH.2.4.restore.trust.in.food, SOC1 = Q7.SSE.2.1.equal.distrib.local.actors, SOC2 = Q7.SSE.2.2.equal.distrib.global.actors, SOC3 = Q7.SSE.2.3.adapt.to.future.pest.pressure, SOC4 = Q7.SSE.2.4.access.to.education.tools.technology, SOC5 = Q7.SSE.2.5.safe.working.conditions, ECON1 = Q8.EC.2.1.efficient.affordable.PM, ECON2 = Q8.EC.2.2.short.term.livelihood.keeping, ECON3 = Q8.EC.2.3.long.term.livelihood.keeping, ECON4 = Q8.EC.2.4.economic.resilience.extreme.pest.events, ECON5 = Q8.EC.2.5.reduction.threat.of.pests.in.productivity, ECON6 = Q8.EC.2.6.reduction.of.indirect.costs.of.PM, overview_continent_all, Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, pesticide_use_ha, pesticide_use_outputvalue, Q3.ISPM.sustainable.pest.management.in.your.crops, row.n) %>% pivot_longer(-c(overview_continent_all,Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, pesticide_use_ha, pesticide_use_outputvalue, Q3.ISPM.sustainable.pest.management.in.your.crops, profession.short, row.n), names_to = "Key_challenge",values_to="Importance_for_sustainable_pest_management") %>% mutate(category_challenge = ifelse(str_detect(Key_challenge,"FS"),"Food Security", ifelse(str_detect(Key_challenge,"ENV"),"Environment", ifelse(str_detect(Key_challenge,"HH"),"Human Health", ifelse(str_detect(Key_challenge,"SOC"),"Social", ifelse(str_detect(Key_challenge,"ECON"),"Economy","Other")))))) imp_wi<- datsur %>% mutate(row.n = row_number()) %>% select("Food Security"= Q4.FS.1.importance.pest.management, Environment= Q5.EE.1.importance.pest.management, "Human Health"= Q6.HH.1.importance.pest.management, Social= Q7.SSE.1.importance.pest.management, Economy= Q8.EC.1.importance.pest.management, row.n) %>% pivot_longer(-c(row.n), names_to = "category_challenge",values_to="Importance_pest_management") exp_we<- datsur %>% mutate(row.n = row_number()) %>% select("Food Security"=Q2.1.PME.expertise.key.challenges.food.security,"Environment"=Q2.2.expertise.key.challenges.environmental.effects, "Human Health"=Q2.3.PME.expertise.key.challenges.health, "Social"=Q2.4.PME.expertise.key.challenges.social.security, "Economy"=Q2.5.PME.expertise.key.challenges.economy, row.n) %>% pivot_longer(-row.n,names_to = "category_challenge",values_to="Strength_Expertise") chall_wi <- left_join(chall, imp_wi, by=c("row.n","category_challenge")) %>% mutate(importance.weight = Importance_pest_management/10) chall_wi_e <- left_join(chall_wi, exp_we, by=c("row.n","category_challenge")) %>% mutate(experience.weight = Strength_Expertise/10) %>% mutate(exp.imp.weight = importance.weight*experience.weight, ag_intensity = pesticide_use_ha/pesticide_use_outputvalue) # ag_intensity is output_value/ha # overall chall_wi_e <- chall_wi_e %>% group_by(Key_challenge) %>% mutate(Importance_for_sustainable_pest_management_mean_ind = mean(Importance_for_sustainable_pest_management, na.rm=T), Importance_for_sustainable_pest_management_se_ind = s(Importance_for_sustainable_pest_management)) %>% ungroup() all_inds_mean_dat <- chall_wi_e %>% select(Key_challenge, Importance_for_sustainable_pest_management_mean_ind, Importance_for_sustainable_pest_management_se_ind) %>% distinct(Key_challenge, Importance_for_sustainable_pest_management_mean_ind, Importance_for_sustainable_pest_management_se_ind) qq_all_inds = quantile(chall_wi_e$Importance_for_sustainable_pest_management, probs = seq(0, 1, .25),na.rm = T) qq_all_inds_mean_025 = quantile(all_inds_mean_dat$Importance_for_sustainable_pest_management_mean_ind, probs = seq(0, 1, .25),na.rm = T) qq_all_inds_mean_033 = quantile(all_inds_mean_dat$Importance_for_sustainable_pest_management_mean_ind, probs = seq(0, 1, .3333),na.rm = T) all_inds_mean_dat2 <- all_inds_mean_dat %>%mutate(l_bound =Importance_for_sustainable_pest_management_mean_ind-2.639*Importance_for_sustainable_pest_management_se_ind, u_bound =Importance_for_sustainable_pest_management_mean_ind+2.639*Importance_for_sustainable_pest_management_se_ind ) %>% mutate(quartile = ifelse(Importance_for_sustainable_pest_management_mean_ind <= 4.711407, "Q1", ifelse(Importance_for_sustainable_pest_management_mean_ind <= 5.513583, "Q2", ifelse(Importance_for_sustainable_pest_management_mean_ind <= 6.419376, "Q3","Q4")))) %>% mutate(third = ifelse(Importance_for_sustainable_pest_management_mean_ind <=5.080189, "T1", ifelse(Importance_for_sustainable_pest_management_mean_ind <=6.120597, "T2","T3"))) general_mean <- mean(chall_wi_e$Importance_for_sustainable_pest_management, na.rm=T) # mean of responses across all indicators plot_dat_ind <- chall_wi_e %>% group_by(category_challenge, Key_challenge) %>% summarise(mean= mean(((Importance_for_sustainable_pest_management-general_mean)/general_mean)*100, na.rm = T), se = s(x=((Importance_for_sustainable_pest_management-general_mean)/general_mean)*100), n=n(), mean_import= mean(Importance_pest_management, na.rm = T), mean_impl = mean(Q3.ISPM.sustainable.pest.management.in.your.crops, na.rm = T)) %>% ungroup() # per region: between -67% : +93% glob <- chall_wi_e %>% filter(overview_continent_all %in% "Global") %>% group_by(Key_challenge) %>% summarize(Importance_for_sustainable_pest_management_mean_ind_glob =mean(Importance_for_sustainable_pest_management, na.rm =T)) %>% ungroup() chall_wi_e_reg <- left_join(chall_wi_e,glob,by="Key_challenge") plot_dat_eur <- chall_wi_e_reg %>% filter(overview_continent_all %in% "Europe") %>% group_by(category_challenge, Key_challenge) %>% summarise(mean= m(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), se = s(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), n=n(), mean_import= m(x=Importance_pest_management), mean_impl = mean(Q3.ISPM.sustainable.pest.management.in.your.crops,na.rm=T)) %>% ungroup() # plot_dat_eur %>% summarize(min = min(mean), max=max(mean)) # range -67%% : -2% plot_dat_na <- chall_wi_e_reg %>% filter(overview_continent_all %in% "North America") %>% group_by(category_challenge, Key_challenge) %>% summarise(mean= m(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), se = s(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), n=n(), mean_import= m(x=Importance_pest_management), mean_impl = mean(Q3.ISPM.sustainable.pest.management.in.your.crops,na.rm=T)) %>% ungroup() # plot_dat_na %>% summarize(min = min(mean), max=max(mean)) # range -24% : +8% plot_dat_sa <- chall_wi_e_reg %>% filter(overview_continent_all %in% "South America") %>% group_by(category_challenge, Key_challenge) %>% summarise(mean= m(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), se = s(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), n=n(), mean_import= m(x=Importance_pest_management), mean_impl = mean(Q3.ISPM.sustainable.pest.management.in.your.crops,na.rm=T)) %>% ungroup() # plot_dat_sa %>% summarize(min = min(mean), max=max(mean)) # range -0.5% : +86% plot_dat_as <- chall_wi_e_reg %>% filter(overview_continent_all %in% "Asia") %>% group_by(category_challenge, Key_challenge) %>% summarise(mean= m(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), se = s(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), n=n(), mean_import= m(x=Importance_pest_management), mean_impl = mean(Q3.ISPM.sustainable.pest.management.in.your.crops,na.rm=T)) %>% ungroup() # plot_dat_as %>% summarize(min = min(mean), max=max(mean)) # range -12% : +80% plot_dat_af <- chall_wi_e_reg %>% filter(overview_continent_all %in% "Africa") %>% group_by(category_challenge, Key_challenge) %>% summarise(mean= m(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), se = s(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), n=n(), mean_import= m(x=Importance_pest_management), mean_impl = mean(Q3.ISPM.sustainable.pest.management.in.your.crops,na.rm=T)) %>% ungroup() # plot_dat_af %>% summarize(min = min(mean), max=max(mean)) # range -5% : +93% plot_dat_oc <- chall_wi_e_reg %>% filter(overview_continent_all %in% "Oceania") %>% group_by(category_challenge, Key_challenge) %>% summarise(mean= m(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), se = s(((x=Importance_for_sustainable_pest_management- Importance_for_sustainable_pest_management_mean_ind_glob)/Importance_for_sustainable_pest_management_mean_ind_glob)*100), n=n(), mean_import= m(x=Importance_pest_management), mean_impl = mean(Q3.ISPM.sustainable.pest.management.in.your.crops,na.rm=T)) %>% ungroup() # plot_dat_oc %>% summarize(min = min(mean), max=max(mean)) # # range -173% : +17% # Comparison of CV across regions plot_dat_cv_reg <- chall_wi_e_reg %>% filter(!(overview_continent_all %in% "Global")) %>% group_by(overview_continent_all,category_challenge, Key_challenge) %>% summarise(cv = sd(Importance_for_sustainable_pest_management+10, na.rm=T)/m(Importance_for_sustainable_pest_management+10)) %>% ungroup() plot_dat_cv_reg2 <- plot_dat_cv_reg %>% select(-category_challenge) %>% pivot_wider(names_from = Key_challenge, values_from= cv) # write.csv(plot_dat_cv_reg2,"Q:/User/Möhring/Paper/Challenges for sustainable pest management/R/Revision paper/CV_across_regions.csv") # Add labels on top of each bar label_data <- plot_dat_ind %>% mutate(Tot=mean) %>% mutate(Tot_lab = ifelse(Tot<0,0,Tot)) angle= round(90 - 360 * (c(1:nrow(label_data))-0.5) /nrow(label_data),0)*-1 # I substract 0.5 because the letter must have the angle of the centre of the bars. Not extreme right(1) or extreme left (0) label_data <- label_data %>% mutate(angle = ifelse(category_challenge %in% c("Human Health","Social"), angle+180, angle)) # label with code and text label_data$lab2 <- "" label_data$lab2[label_data$Key_challenge =="ECON1"] <- "ECON1: Cost-efficient" label_data$lab2[label_data$Key_challenge =="ECON2"] <- "ECON2: Short-run Income" label_data$lab2[label_data$Key_challenge =="ECON3"] <- "ECON3: Long-run Income " label_data$lab2[label_data$Key_challenge =="ECON4"] <- "ECON4: Economic Resilience " label_data$lab2[label_data$Key_challenge =="ECON5"] <- "ECON5: Productivity Growth " label_data$lab2[label_data$Key_challenge =="ECON6"] <- "ECON6: Indirect Costs" label_data$lab2[label_data$Key_challenge =="ENV1"] <- "ENV1: Drinking Water" label_data$lab2[label_data$Key_challenge =="ENV2"] <- "ENV2: Soils" label_data$lab2[label_data$Key_challenge =="ENV3"] <- "ENV3: Marine Eco-Sys" label_data$lab2[label_data$Key_challenge =="ENV4"] <- "ENV4: Freshwater Eco-Sys " label_data$lab2[label_data$Key_challenge =="ENV5"] <- "ENV5: Biodiversity" label_data$lab2[label_data$Key_challenge =="FS1"] <- "FS1: Food Quantity" label_data$lab2[label_data$Key_challenge =="FS2"] <- "FS2: Healthy Diets" label_data$lab2[label_data$Key_challenge =="FS3"] <- "FS3: Safe Food" label_data$lab2[label_data$Key_challenge =="FS4"] <- "FS4: Resilience Hunger" label_data$lab2[label_data$Key_challenge =="HH1"] <- " HH1: Farm Workers" label_data$lab2[label_data$Key_challenge =="HH2"] <- " HH2: Residents" label_data$lab2[label_data$Key_challenge =="HH3"] <- " HH3: Consumers" label_data$lab2[label_data$Key_challenge =="HH4"] <- " HH4: Healthy Food Choice" label_data$lab2[label_data$Key_challenge =="SOC1"] <- "SOC1: Local Equality" label_data$lab2[label_data$Key_challenge =="SOC2"] <- "SOC2: Global Equality" label_data$lab2[label_data$Key_challenge =="SOC3"] <- " SOC3: Adaption Capacity" label_data$lab2[label_data$Key_challenge =="SOC4"] <- "SOC4: Equal Access PM" label_data$lab2[label_data$Key_challenge =="SOC5"] <- "SOC5: Safe Work" # Plot 1 ggplot(label_data, aes(Key_challenge,mean, fill=category_challenge)) + geom_bar(stat="identity") + geom_errorbar(aes(ymin=mean-1.95*se, ymax=mean+1.95*se), color= "black", width=.3, linewidth=1.0, alpha=0.8)+ # geom_point(size=5,alpha=0.5, color= "grey")+ # geom_errorbar(aes(ymin=mean-se*qt(p=.05/2, df=n-1, lower.tail=FALSE), # ymax= mean+se*qt(p=.05/2, df=n-1, lower.tail=FALSE)), # alpha=0.5, color= "grey", width=.3, size=1.2)+ coord_polar("x", start=0,direction = -1)+ scale_fill_manual(values = c("#332288", "#117733", "#44AA99", "#88CCEE", "#DDCC77"))+ ylab("Expert rating (global average)")+ theme(legend.position = "none", axis.ticks = element_blank(),axis.title.y = element_blank(), axis.title.x = element_blank(), axis.text = element_blank(), panel.background = element_blank()) + geom_hline(aes(yintercept=-40,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-30,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-20,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-10,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=0,alpha=0.1),color="black", linewidth=0.65)+ geom_hline(aes(yintercept=10, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=20, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=30, alpha=0.1),color="darkgray")+ geom_text(label="-40",x=0.4,y=-40, color="black", size=3, angle =-1, alpha=0.1)+ geom_text(label="-30",x=0.4,y=-30, color="black", size=3, angle =-1, alpha=0.1)+ geom_text(label="-20",x=0.4,y=-20, color="black", size=3, angle =-1, alpha=0.1)+ geom_text(label="-10",x=0.4,y=-10, color="black", size=3, angle =-1, alpha=0.1)+ geom_text(label=" 0",x=0.4,y=0, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label="+10",x=0.4,y=10, color="black", size=3, angle =-1, alpha=0.1)+ geom_text(label="+20",x=0.4,y=20, color="black", size=3, angle =-1)+ geom_text(label="+30",x=0.4,y=30, color="black", size=3, angle =-1)+ geom_text(aes(label=lab2, y=Tot_lab+30), color="black", size=3, angle= label_data$angle) # Plot 2 - Europe - no labels label_data_eur <- label_data %>% arrange(Key_challenge) plot_dat_eur <- plot_dat_eur %>% arrange(Key_challenge) label_data_eur$mean <- plot_dat_eur$mean ggplot(label_data_eur, aes(Key_challenge,mean, fill=category_challenge)) + geom_bar(stat="identity") + geom_errorbar(aes(ymin=mean-1.95*se, ymax=mean+1.95*se), color= "black", width=.3, linewidth=1, alpha=0.8)+ coord_polar("x", start=0,direction = -1)+ scale_fill_manual(values = c("#332288", "#117733", "#44AA99", "#88CCEE", "#DDCC77"))+ ylab("Expert rating (global average)")+ theme( legend.position = "none", axis.ticks = element_blank(),axis.title.y = element_blank(), axis.title.x = element_blank(), axis.text = element_blank(), panel.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing panel outline panel.grid.major = element_blank(), # get rid of major grid panel.grid.minor = element_blank(), # get rid of minor grid plot.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing plot outline legend.background = element_rect(fill = "transparent"), legend.box.background = element_rect(fill = "transparent"), legend.key = element_rect(fill = "transparent")) + geom_hline(aes(yintercept=-100,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-80,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-60,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-40,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-20,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=0,alpha=0.1),color="black", linewidth=1.5)+ geom_hline(aes(yintercept=20, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=40, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=60, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=80, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=100, alpha=0.1),color="darkgray")+ geom_text(label=" 0",x=0.5,y=10, color="black", size=5, angle =-1, alpha=0.3)+ geom_text(label="+100",x=0.5,y=100, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" +50",x=0.5,y=50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" -50",x=0.5,y=-50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label="-100",x=0.5,y=-95, color="black", size=4, angle =-1, alpha=0.3) # geom_text(label="-40",x=0.4,y=-40, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-30",x=0.4,y=-30, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-20",x=0.4,y=-20, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-10",x=0.4,y=-10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+10",x=0.4,y=10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+20",x=0.4,y=20, color="black", size=3, angle =-1)+ # geom_text(label="+30",x=0.4,y=30, color="black", size=3, angle =-1) # Plot 2 - Oceania - no labels label_data_oc <- label_data %>% arrange(Key_challenge) plot_dat_oc <- plot_dat_oc %>% arrange(Key_challenge) label_data_oc$mean <- plot_dat_oc$mean ggplot(label_data_oc, aes(Key_challenge,mean, fill=category_challenge)) + geom_bar(stat="identity") + geom_errorbar(aes(ymin=mean-1.95*se, ymax=mean+1.95*se), color= "black", width=.3, linewidth=1, alpha=0.8)+ coord_polar("x", start=0,direction = -1)+ scale_fill_manual(values = c("#332288", "#117733", "#44AA99", "#88CCEE", "#DDCC77"))+ ylab("Expert rating (global average)")+ theme( legend.position = "none", axis.ticks = element_blank(),axis.title.y = element_blank(), axis.title.x = element_blank(), axis.text = element_blank(), panel.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing panel outline panel.grid.major = element_blank(), # get rid of major grid panel.grid.minor = element_blank(), # get rid of minor grid plot.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing plot outline legend.background = element_rect(fill = "transparent"), legend.box.background = element_rect(fill = "transparent"), legend.key = element_rect(fill = "transparent")) + geom_hline(aes(yintercept=-100,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-80,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-60,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-40,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-20,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=0,alpha=0.1),color="black", linewidth=1.5)+ geom_hline(aes(yintercept=20, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=40, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=60, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=80, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=100, alpha=0.1),color="darkgray")+ geom_text(label=" 0",x=0.5,y=10, color="black", size=5, angle =-1, alpha=0.3)+ geom_text(label="+100",x=0.5,y=100, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" +50",x=0.5,y=50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" -50",x=0.5,y=-50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label="-100",x=0.5,y=-95, color="black", size=4, angle =-1, alpha=0.3) # geom_text(label="-40",x=0.4,y=-40, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-30",x=0.4,y=-30, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-20",x=0.4,y=-20, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-10",x=0.4,y=-10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+10",x=0.4,y=10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+20",x=0.4,y=20, color="black", size=3, angle =-1)+ # geom_text(label="+30",x=0.4,y=30, color="black", size=3, angle =-1) # Plot 2 - North America - no label label_data_na <- label_data %>% arrange(Key_challenge) plot_dat_na <- plot_dat_na %>% arrange(Key_challenge) label_data_na$mean <- plot_dat_na$mean ggplot(label_data_na, aes(Key_challenge,mean, fill=category_challenge)) + geom_bar(stat="identity") + geom_errorbar(aes(ymin=mean-1.95*se, ymax=mean+1.95*se), color= "black", width=.3, linewidth=1, alpha=0.8)+ coord_polar("x", start=0,direction = -1)+ scale_fill_manual(values = c("#332288", "#117733", "#44AA99", "#88CCEE", "#DDCC77"))+ ylab("Expert rating (global average)")+ theme(legend.position = "none", axis.ticks = element_blank(),axis.title.y = element_blank(), axis.title.x = element_blank(), axis.text = element_blank(), panel.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing panel outline panel.grid.major = element_blank(), # get rid of major grid panel.grid.minor = element_blank(), # get rid of minor grid plot.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing plot outline legend.background = element_rect(fill = "transparent"), legend.box.background = element_rect(fill = "transparent"), legend.key = element_rect(fill = "transparent")) + geom_hline(aes(yintercept=-100,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-80,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-60,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-40,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-20,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=0,alpha=0.1),color="black", linewidth=1.5)+ geom_hline(aes(yintercept=20, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=40, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=60, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=80, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=100, alpha=0.1),color="darkgray")+ geom_text(label=" 0",x=0.5,y=10, color="black", size=5, angle =-1, alpha=0.3)+ geom_text(label="+100",x=0.5,y=100, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" +50",x=0.5,y=50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" -50",x=0.5,y=-50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label="-100",x=0.5,y=-95, color="black", size=4, angle =-1, alpha=0.3) # geom_text(label="-40",x=0.4,y=-40, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-30",x=0.4,y=-30, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-20",x=0.4,y=-20, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-10",x=0.4,y=-10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+10",x=0.4,y=10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+20",x=0.4,y=20, color="black", size=3, angle =-1)+ # geom_text(label="+30",x=0.4,y=30, color="black", size=3, angle =-1) # Plot 2 - South America - no label label_data_sa <- label_data %>% arrange(Key_challenge) plot_dat_sa <- plot_dat_sa %>% arrange(Key_challenge) label_data_sa$mean <- plot_dat_sa$mean ggplot(label_data_sa, aes(Key_challenge,mean, fill=category_challenge)) + geom_bar(stat="identity") + geom_errorbar(aes(ymin=mean-1.95*se, ymax=mean+1.95*se), color= "black", width=.3, linewidth=1, alpha=0.8)+ coord_polar("x", start=0,direction = -1)+ scale_fill_manual(values = c("#332288", "#117733", "#44AA99", "#88CCEE", "#DDCC77"))+ ylab("Expert rating (global average)")+ theme(legend.position = "none", axis.ticks = element_blank(),axis.title.y = element_blank(), axis.title.x = element_blank(), axis.text = element_blank(), panel.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing panel outline panel.grid.major = element_blank(), # get rid of major grid panel.grid.minor = element_blank(), # get rid of minor grid plot.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing plot outline legend.background = element_rect(fill = "transparent"), legend.box.background = element_rect(fill = "transparent"), legend.key = element_rect(fill = "transparent")) + geom_hline(aes(yintercept=-100,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-80,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-60,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-40,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-20,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=0,alpha=0.1),color="black", linewidth=1.5)+ geom_hline(aes(yintercept=20, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=40, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=60, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=80, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=100, alpha=0.1),color="darkgray")+ geom_text(label=" 0",x=0.5,y=10, color="black", size=5, angle =-1, alpha=0.3)+ geom_text(label="+100",x=0.5,y=100, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" +50",x=0.5,y=50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" -50",x=0.5,y=-50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label="-100",x=0.5,y=-95, color="black", size=4, angle =-1, alpha=0.3) # geom_text(label="-40",x=0.4,y=-40, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-30",x=0.4,y=-30, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-20",x=0.4,y=-20, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-10",x=0.4,y=-10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+10",x=0.4,y=10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+20",x=0.4,y=20, color="black", size=3, angle =-1)+ # geom_text(label="+30",x=0.4,y=30, color="black", size=3, angle =-1) # Plot 2 - Asia - no label label_data_as <- label_data %>% arrange(Key_challenge) plot_dat_as <- plot_dat_as %>% arrange(Key_challenge) label_data_as$mean <- plot_dat_as$mean ggplot(label_data_as, aes(Key_challenge,mean, fill=category_challenge)) + geom_bar(stat="identity") + geom_errorbar(aes(ymin=mean-1.95*se, ymax=mean+1.95*se), color= "black", width=.3, linewidth=1, alpha=0.8)+ coord_polar("x", start=0,direction = -1)+ scale_fill_manual(values = c("#332288", "#117733", "#44AA99", "#88CCEE", "#DDCC77"))+ ylab("Expert rating (global average)")+ theme(legend.position = "none", axis.ticks = element_blank(),axis.title.y = element_blank(), axis.title.x = element_blank(), axis.text = element_blank(), panel.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing panel outline panel.grid.major = element_blank(), # get rid of major grid panel.grid.minor = element_blank(), # get rid of minor grid plot.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing plot outline legend.background = element_rect(fill = "transparent"), legend.box.background = element_rect(fill = "transparent"), legend.key = element_rect(fill = "transparent")) + geom_hline(aes(yintercept=-100,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-80,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-60,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-40,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-20,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=0,alpha=0.1),color="black", linewidth=1.5)+ geom_hline(aes(yintercept=20, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=40, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=60, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=80, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=100, alpha=0.1),color="darkgray")+ geom_text(label=" 0",x=0.5,y=10, color="black", size=5, angle =-1, alpha=0.3)+ geom_text(label="+100",x=0.5,y=100, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" +50",x=0.5,y=50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" -50",x=0.5,y=-50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label="-100",x=0.5,y=-95, color="black", size=4, angle =-1, alpha=0.3) # geom_text(label="-40",x=0.4,y=-40, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-30",x=0.4,y=-30, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-20",x=0.4,y=-20, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-10",x=0.4,y=-10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+10",x=0.4,y=10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+20",x=0.4,y=20, color="black", size=3, angle =-1)+ # geom_text(label="+30",x=0.4,y=30, color="black", size=3, angle =-1) # Plot 2 - Africa - no label label_data_af <- label_data %>% arrange(Key_challenge) plot_dat_af <- plot_dat_af %>% arrange(Key_challenge) label_data_af$mean <- plot_dat_af$mean ggplot(label_data_af, aes(Key_challenge,mean, fill=category_challenge)) + geom_bar(stat="identity") + geom_errorbar(aes(ymin=mean-1.95*se, ymax=mean+1.95*se), color= "black", width=.3, linewidth=1, alpha=0.8)+ coord_polar("x", start=0,direction = -1)+ scale_fill_manual(values = c("#332288", "#117733", "#44AA99", "#88CCEE", "#DDCC77"))+ ylab("Expert rating (global average)")+ theme(legend.position = "none", axis.ticks = element_blank(),axis.title.y = element_blank(), axis.title.x = element_blank(), axis.text = element_blank(), panel.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing panel outline panel.grid.major = element_blank(), # get rid of major grid panel.grid.minor = element_blank(), # get rid of minor grid plot.background = element_rect(fill = "transparent", colour = NA_character_), # necessary to avoid drawing plot outline legend.background = element_rect(fill = "transparent"), legend.box.background = element_rect(fill = "transparent"), legend.key = element_rect(fill = "transparent")) + geom_hline(aes(yintercept=-100,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-80,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-60,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-40,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=-20,alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=0,alpha=0.1),color="black", linewidth=1.5)+ geom_hline(aes(yintercept=20, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=40, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=60, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=80, alpha=0.1),color="darkgray")+ geom_hline(aes(yintercept=100, alpha=0.1),color="darkgray")+ geom_text(label=" 0",x=0.5,y=10, color="black", size=5, angle =-1, alpha=0.3)+ geom_text(label="+100",x=0.5,y=100, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" +50",x=0.5,y=50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label=" -50",x=0.5,y=-50, color="black", size=4, angle =-1, alpha=0.3)+ geom_text(label="-100",x=0.5,y=-95, color="black", size=4, angle =-1, alpha=0.3) # geom_text(label="-40",x=0.4,y=-40, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-30",x=0.4,y=-30, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-20",x=0.4,y=-20, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="-10",x=0.4,y=-10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+10",x=0.4,y=10, color="black", size=3, angle =-1, alpha=0.1)+ # geom_text(label="+20",x=0.4,y=20, color="black", size=3, angle =-1)+ # geom_text(label="+30",x=0.4,y=30, color="black", size=3, angle =-1) ### iv) Replicate analyses of correlations and p-vals #### plot_quick<- datsur %>% mutate(row.n = row_number()) %>% select(gdp_capita, gdp_share_agri,HDI, PHDI, gfsi22, risk_score, SDG_1.5_affected_disaster, SDG_2.1_undernourishment, SDG_2.3_agricult_prod, SDG_2A_agricult_invest, SDG_3.9_chem_poisoning, SDG_6.3_water_quality, SDG_8.4_material_consump, SDG_8.8_work_injuries, SDG_12.3_food_waste, SDG_12.4_Rotterdam, SDG_14.1_Ocean_pollution, SDG_15.1_Freshwater_prot, SDG_15.5_red_list, FS1 = Q4.FS.2.1.sufficient.quantity, FS2 = Q4.FS.2.2.healthy.food, FS3 = Q4.FS.2.3.safe.food.and.feed, FS4 = Q4.FS.2.4.resilience.of.food, ENV1 = Q5.EE.2.1.pollution.drinking.water, ENV2 = Q5.EE.2.2.pollution.soil, ENV3 = Q5.EE.2.3.pollution.marine.ecosystems, ENV4 = Q5.EE.2.4.pollution.freshwater.ecosystems, ENV5 = Q5.EE.2.5.loss.of.genes.biodiversity, HH1 = Q6.HH.2.1.effects.farm.workers, HH2 = Q6.HH.2.2.effects.bystanders, HH3 = Q6.HH.2.3.effects.on.consumers, HH4= Q6.HH.2.4.restore.trust.in.food, SOC1 = Q7.SSE.2.1.equal.distrib.local.actors, SOC2 = Q7.SSE.2.2.equal.distrib.global.actors, SOC3 = Q7.SSE.2.3.adapt.to.future.pest.pressure, SOC4 = Q7.SSE.2.4.access.to.education.tools.technology, SOC5 = Q7.SSE.2.5.safe.working.conditions, ECON1 = Q8.EC.2.1.efficient.affordable.PM, ECON2 = Q8.EC.2.2.short.term.livelihood.keeping, ECON3 = Q8.EC.2.3.long.term.livelihood.keeping, ECON4 = Q8.EC.2.4.economic.resilience.extreme.pest.events, ECON5 = Q8.EC.2.5.reduction.threat.of.pests.in.productivity, ECON6 = Q8.EC.2.6.reduction.of.indirect.costs.of.PM) name_plot.x <- names(plot_quick)[!(names(plot_quick)%in% c("SDG_1.5_affected_disaster", "SDG_2.1_undernourishment", "SDG_2.3_agricult_prod", "SDG_2A_agricult_invest", "SDG_3.9_chem_poisoning", "SDG_6.3_water_quality", "SDG_8.4_material_consump", "SDG_8.8_work_injuries", "SDG_12.3_food_waste", "SDG_12.4_Rotterdam", "SDG_14.1_Ocean_pollution", "SDG_15.1_Freshwater_prot", "SDG_15.5_red_list", "gdp_capita", "gdp_share_agri","HDI", "PHDI", "gfsi22", "risk_score"))] name_plot.y <- names(plot_quick)[(names(plot_quick) %in% c("SDG_1.5_affected_disaster", "SDG_2.1_undernourishment", "SDG_2.3_agricult_prod", "SDG_2A_agricult_invest", "SDG_3.9_chem_poisoning", "SDG_6.3_water_quality", "SDG_8.4_material_consump", "SDG_8.8_work_injuries", "SDG_12.3_food_waste", "SDG_12.4_Rotterdam", "SDG_14.1_Ocean_pollution", "SDG_15.1_Freshwater_prot", "SDG_15.5_red_list", "gdp_capita", "gdp_share_agri","HDI", "PHDI", "gfsi22", "risk_score"))] # kendall (rank) correlations of all cor.table = data.frame(matrix(nrow = length(name_plot.y), ncol = length(name_plot.x))) colnames(cor.table) = name_plot.x row.names(cor.table) = name_plot.y for( x in 1:length(name_plot.x)){ for( y in 1:length(name_plot.y)){ c <- cor.test(x = plot_quick[,name_plot.x[[x]]], y = plot_quick[,name_plot.y[[y]]], use= "complete.obs", method = "kendall") t <- paste0( format(round(c$estimate,2), nsmall = 2), " (",format(round(c$p.value,2), nsmall = 2), ")") cor.table[y,x] <- t } } # Turn those combinations of Pot SPM X SDG Levels NA, that are not part of Table 1 # 1.5 - ECON4 # 2.1 - FS1, FS2, FS3 # 2.3 - ECON5, ECON6, SOC3 # 2A - ECON1 # 3.9 - HH1, HH2, HH3 # 6.3 - ENV1 # 8.4 - ECON1, ECON2, ECON3, ECON5 # 8.8 - SOC5 # 12.3 - FS1 # 12.4 - ENV1, ENV2, ENV3, ENV4, ENV5, HH1, HH2, HH3, HH4, FS1, FS2, FS3, ECON6 # 14.1 - ENV3 # 15.1 - ENV2, ENV4 # 15.5 - ENV5 cor.table[name_plot.y[6+1],name_plot.x[! name_plot.x %in% c("ECON4")]] <- NA cor.table[name_plot.y[6+2],name_plot.x[! name_plot.x %in% c("FS1", "FS2", "FS3" )]] <- NA cor.table[name_plot.y[6+3],name_plot.x[! name_plot.x %in% c( "ECON5", "ECON6","SOC3")]] <- NA cor.table[name_plot.y[6+4],name_plot.x[! name_plot.x %in% c( "ECON1")]] <- NA cor.table[name_plot.y[6+5],name_plot.x[! name_plot.x %in% c( "HH1", "HH2", "HH3")]] <- NA cor.table[name_plot.y[6+6],name_plot.x[! name_plot.x %in% c( "ENV1")]] <- NA cor.table[name_plot.y[6+7],name_plot.x[! name_plot.x %in% c( "ECON1", "ECON2", "ECON3", "ECON5")]] <- NA cor.table[name_plot.y[6+8],name_plot.x[! name_plot.x %in% c( "SOC5")]] <- NA cor.table[name_plot.y[6+9],name_plot.x[! name_plot.x %in% c( "FS1")]] <- NA cor.table[name_plot.y[6+10],name_plot.x[! name_plot.x %in% c( "ENV1", "ENV2", "ENV3", "ENV4", "ENV5", "HH1", "HH2", "HH3", "HH4", "FS1", "FS2", "FS3", "ECON6")]] <- NA cor.table[name_plot.y[6+11],name_plot.x[! name_plot.x %in% c( "ENV3")]] <- NA cor.table[name_plot.y[6+12],name_plot.x[! name_plot.x %in% c( "ENV2", "ENV4")]] <- NA cor.table[name_plot.y[6+13],name_plot.x[! name_plot.x %in% c( "ENV5")]] <- NA ### v) Replicate Fig. 3 and underlying analyses #### library(coefplot) library(car) chall <- datsur %>% mutate(row.n = row_number()) %>% select(FS1 = Q4.FS.2.1.sufficient.quantity, FS2 = Q4.FS.2.2.healthy.food, FS3 = Q4.FS.2.3.safe.food.and.feed, FS4 = Q4.FS.2.4.resilience.of.food, ENV1 = Q5.EE.2.1.pollution.drinking.water, ENV2 = Q5.EE.2.2.pollution.soil, ENV3 = Q5.EE.2.3.pollution.marine.ecosystems, ENV4 = Q5.EE.2.4.pollution.freshwater.ecosystems, ENV5 = Q5.EE.2.5.loss.of.genes.biodiversity, HH1 = Q6.HH.2.1.effects.farm.workers, HH2 = Q6.HH.2.2.effects.bystanders, HH3 = Q6.HH.2.3.effects.on.consumers, HH4= Q6.HH.2.4.restore.trust.in.food, SOC1 = Q7.SSE.2.1.equal.distrib.local.actors, SOC2 = Q7.SSE.2.2.equal.distrib.global.actors, SOC3 = Q7.SSE.2.3.adapt.to.future.pest.pressure, SOC4 = Q7.SSE.2.4.access.to.education.tools.technology, SOC5 = Q7.SSE.2.5.safe.working.conditions, ECON1 = Q8.EC.2.1.efficient.affordable.PM, ECON2 = Q8.EC.2.2.short.term.livelihood.keeping, ECON3 = Q8.EC.2.3.long.term.livelihood.keeping, ECON4 = Q8.EC.2.4.economic.resilience.extreme.pest.events, ECON5 = Q8.EC.2.5.reduction.threat.of.pests.in.productivity, ECON6 = Q8.EC.2.6.reduction.of.indirect.costs.of.PM, overview_continent_all, pesticide_use_ha, pesticide_use_outputvalue, risk_score, outputvalue, cropland, schooling, SDG_2.3_agricult_prod, SDG_2A_agricult_invest,ai_count, ai_bans, organic_share, rel_attainable_yield, rel_yield_gap, Actual.loss,Efficiency.PM,Potential.loss, gfsi22, Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, Q3.ISPM.sustainable.pest.management.in.your.crops, row.n) %>% pivot_longer(-c(overview_continent_all, pesticide_use_ha, pesticide_use_outputvalue, risk_score, outputvalue, cropland, schooling, SDG_2.3_agricult_prod, SDG_2A_agricult_invest,ai_count, ai_bans, organic_share, rel_attainable_yield, rel_yield_gap, Actual.loss,Efficiency.PM,Potential.loss,gfsi22, Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, Q3.ISPM.sustainable.pest.management.in.your.crops, profession.short, row.n), names_to = "Key_challenge",values_to="Importance_for_sustainable_pest_management") %>% mutate(category_challenge = ifelse(str_detect(Key_challenge,"FS"),"Food Security", ifelse(str_detect(Key_challenge,"ENV"),"Environment", ifelse(str_detect(Key_challenge,"HH"),"Human Health", ifelse(str_detect(Key_challenge,"SOC"),"Social", ifelse(str_detect(Key_challenge,"ECON"),"Economy","Other")))))) imp_wi<- datsur %>% mutate(row.n = row_number()) %>% select("Food Security"= Q4.FS.1.importance.pest.management, Environment= Q5.EE.1.importance.pest.management, "Human Health"= Q6.HH.1.importance.pest.management, Social= Q7.SSE.1.importance.pest.management, Economy= Q8.EC.1.importance.pest.management, row.n) %>% pivot_longer(-c(row.n), names_to = "category_challenge",values_to="Importance_pest_management") exp_we<- datsur %>% mutate(row.n = row_number()) %>% select("Food Security"=Q2.1.PME.expertise.key.challenges.food.security,"Environment"=Q2.2.expertise.key.challenges.environmental.effects, "Human Health"=Q2.3.PME.expertise.key.challenges.health, "Social"=Q2.4.PME.expertise.key.challenges.social.security, "Economy"=Q2.5.PME.expertise.key.challenges.economy, row.n) %>% pivot_longer(-row.n,names_to = "category_challenge",values_to="Strength_Expertise") chall_wi <- left_join(chall, imp_wi, by=c("row.n","category_challenge")) %>% mutate(importance.weight = Importance_pest_management/10) chall_wi_e <- left_join(chall_wi, exp_we, by=c("row.n","category_challenge")) %>% mutate(experience.weight = Strength_Expertise/10) %>% mutate(exp.imp.weight = importance.weight*experience.weight) # regression data reg <- chall_wi_e %>% filter(overview_continent_all %in% c("Global","Africa","Asia","Europe", "North America", "South America", "Oceania")) %>% filter(Q1.4.PME.crops.expertise %in% c("Arable crops","General (all)","Horticulture")) %>% filter(field.expertise.short %in% c("PMS","Ecol","Toxi","Soc-Eco","BioTech", "Other")) %>% rename(continent = overview_continent_all, crop = Q1.4.PME.crops.expertise, field = field.expertise.short, scope = Q1.1.PME.specify.larger.vs.specific) %>% mutate(continent = relevel(as.factor(continent), ref="Global"), crop = relevel(as.factor(crop), ref = "General (all)"), field = relevel(as.factor(field), ref= "PMS"), Key_challenge = as.factor(Key_challenge), pesticide_use_ha_eff= pesticide_use_ha/Efficiency.PM, pesticide_use_ha_act =pesticide_use_ha/Actual.loss, implementation.level = Q3.ISPM.sustainable.pest.management.in.your.crops, importance.category= Importance_pest_management, expertise.category = Strength_Expertise, pot.pest.damages = Potential.loss, # those are just renamed ag_intensity = outputvalue/cropland, ag_output_work = SDG_2.3_agricult_prod/1000, ag_invest = SDG_2A_agricult_invest/1000, loss_potential = pot.pest.damages*rel_attainable_yield) # ag_intensity is output_value/ha model <- Importance_for_sustainable_pest_management ~ crop + scope + field + implementation.level + expertise.category + importance.category +Potential.loss + risk_score + rel_attainable_yield # Analyses and Fig. 3 a <- summary(lm(data = reg %>% filter(Key_challenge == "FS1"), formula = model )) b <- summary(lm(data = reg %>% filter(Key_challenge == "FS2"), formula = model )) c <- summary(lm(data = reg %>% filter(Key_challenge == "FS3"), formula = model )) d <- summary(lm(data = reg %>% filter(Key_challenge == "FS4"), formula = model )) e <- summary(lm(data = reg %>% filter(Key_challenge == "ECON1"), formula = model )) f <- summary(lm(data = reg %>% filter(Key_challenge == "ECON2"), formula = model )) g <- summary(lm(data = reg %>% filter(Key_challenge == "ECON3"), formula = model )) h <- summary(lm(data = reg %>% filter(Key_challenge == "ECON4"), formula = model )) i <- summary(lm(data = reg %>% filter(Key_challenge == "ECON5"), formula = model )) j <- summary(lm(data = reg %>% filter(Key_challenge == "ECON6"), formula = model )) k <- summary(lm(data = reg %>% filter(Key_challenge == "ENV1"), formula = model )) l <- summary(lm(data = reg %>% filter(Key_challenge == "ENV2"), formula = model )) m <- summary(lm(data = reg %>% filter(Key_challenge == "ENV3"), formula = model )) n <- summary(lm(data = reg %>% filter(Key_challenge == "ENV4"), formula = model )) o <- summary(lm(data = reg %>% filter(Key_challenge == "ENV5"), formula = model )) p <- summary(lm(data = reg %>% filter(Key_challenge == "HH1"), formula = model )) q <- summary(lm(data = reg %>% filter(Key_challenge == "HH2"), formula = model )) r <- summary(lm(data = reg %>% filter(Key_challenge == "HH3"), formula = model )) s <- summary(lm(data = reg %>% filter(Key_challenge == "HH4"), formula = model )) t <- summary(lm(data = reg %>% filter(Key_challenge == "SOC1"), formula = model )) u <- summary(lm(data = reg %>% filter(Key_challenge == "SOC2"), formula = model )) v <- summary(lm(data = reg %>% filter(Key_challenge == "SOC3"), formula = model )) w <- summary(lm(data = reg %>% filter(Key_challenge == "SOC4"), formula = model )) x <- summary(lm(data = reg %>% filter(Key_challenge == "SOC5"), formula = model )) results <- list(a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r, s, t, u, v, w, x) nvar <- length(rownames(results[[1]]$coefficients))-1 nresults <- length(results) coeff <- c() pval <- c() for ( i in 1:nresults){ coeff <- c(coeff, as.numeric(results[[i]]$coefficients[1:nvar+1,"Estimate"])) pval <- c(pval, as.numeric(results[[i]]$coefficients[1:nvar+1,"Pr(>|t|)"])) } test_real_a <- data.frame( key_challenge = c(rep("FS1",nvar), rep("FS2",nvar), rep("FS3",nvar), rep("FS4",nvar), rep("ECON1",nvar), rep("ECON2",nvar), rep("ECON3",nvar), rep("ECON4",nvar), rep("ECON5",nvar), rep("ECON6",nvar), rep("ENV1",nvar), rep("ENV2",nvar), rep("ENV3",nvar), rep("ENV4",nvar), rep("ENV5",nvar), rep("HH1",nvar), rep("HH2",nvar), rep("HH3",nvar), rep("HH4",nvar), rep("SOC1",nvar), rep("SOC2",nvar), rep("SOC3",nvar), rep("SOC4",nvar), rep("SOC5",nvar)), reg_var = as.factor(c(rep(rownames(results[[1]]$coefficients)[-1], nresults))), reg_coeff = coeff, reg_p= pval ) test_real <- data.frame( key_challenge = c(rep("FS1",nvar), rep("FS2",nvar), rep("FS3",nvar), rep("FS4",nvar), rep("ECON1",nvar), rep("ECON2",nvar), rep("ECON3",nvar), rep("ECON4",nvar), rep("ECON5",nvar), rep("ECON6",nvar), rep("ENV1",nvar), rep("ENV2",nvar), rep("ENV3",nvar), rep("ENV4",nvar), rep("ENV5",nvar), rep("HH1",nvar), rep("HH2",nvar), rep("HH3",nvar), rep("HH4",nvar), rep("SOC1",nvar), rep("SOC2",nvar), rep("SOC3",nvar), rep("SOC4",nvar), rep("SOC5",nvar)), reg_var = as.factor(c(rep(rownames(results[[1]]$coefficients)[-1], nresults))), reg_dir = as.factor (ifelse(coeff>0, "Positive","Negative")), reg_p= as.factor( ifelse( pval>0.05, "not sign.", ifelse(pval<0.05 & pval>0.01, "<0.05","<0.01"))) ) test_real$key_challenge <- factor(test_real$key_challenge, levels = unique(test_real$key_challenge[order(test_real$key_challenge, decreasing = T)])) test_real2 <- test_real %>% mutate(reg_var = factor(reg_var, levels = c("cropArable crops", "cropHorticulture", "Potential.loss", "rel_attainable_yield", "risk_score", "fieldBioTech", "fieldToxi", "fieldEcol","fieldSoc-Eco", "fieldOther", "scopespecific","expertise.category", "importance.category", "implementation.level"))) ggplot(test_real2, aes(reg_var, key_challenge)) + geom_raster(aes( fill= reg_dir, alpha = reg_p), hjust =0, vjust =0.00) + geom_vline(xintercept=seq(0, 14, by=1))+ geom_hline(yintercept=seq(0, 24, by=1))+ scale_fill_manual(values = c("Positive"= "#8BCF4F","Negative" = "#185D8C"))+ # original "lighter" colors: "Positive"= "#b2df8a","Negative" = "#1F78B4" scale_alpha_manual(values = c("not sign." = 0,"<0.05"=0.4,"<0.01"=1))+ theme(axis.title = element_blank(), axis.text.x = element_text(angle = 45, vjust = -1, hjust=0, face = "bold"), axis.text.y = element_text( vjust = 1, hjust=0, face = "bold"))+ labs(fill = "Coefficient\nSign", alpha = "Coefficient\nP-value") + guides(fill = guide_legend(order=1), alpha = guide_legend(order=2))+ scale_y_discrete( expand = c(0,0))+ scale_x_discrete(position = "top", expand = c(0,0), labels = c("Type: Arable", "Type: Horticulture", "Pot. Pest Pressure", "Attainable Yield", "Level Pesticide Risk", "Field: BioTech", "Field: Toxicology", "Field: Ecology","Field: Socio-Econ", "Field: Other", "Scope: Local expert","Expertise in Category", "Importance of PM", "Level Sustainable PM" )) # Analyses and Fig. with Bonferroni correction a <- summary(lm(data = reg %>% filter(Key_challenge == "FS1"), formula = model )) b <- summary(lm(data = reg %>% filter(Key_challenge == "FS2"), formula = model )) c <- summary(lm(data = reg %>% filter(Key_challenge == "FS3"), formula = model )) d <- summary(lm(data = reg %>% filter(Key_challenge == "FS4"), formula = model )) e <- summary(lm(data = reg %>% filter(Key_challenge == "ECON1"), formula = model )) f <- summary(lm(data = reg %>% filter(Key_challenge == "ECON2"), formula = model )) g <- summary(lm(data = reg %>% filter(Key_challenge == "ECON3"), formula = model )) h <- summary(lm(data = reg %>% filter(Key_challenge == "ECON4"), formula = model )) i <- summary(lm(data = reg %>% filter(Key_challenge == "ECON5"), formula = model )) j <- summary(lm(data = reg %>% filter(Key_challenge == "ECON6"), formula = model )) k <- summary(lm(data = reg %>% filter(Key_challenge == "ENV1"), formula = model )) l <- summary(lm(data = reg %>% filter(Key_challenge == "ENV2"), formula = model )) m <- summary(lm(data = reg %>% filter(Key_challenge == "ENV3"), formula = model )) n <- summary(lm(data = reg %>% filter(Key_challenge == "ENV4"), formula = model )) o <- summary(lm(data = reg %>% filter(Key_challenge == "ENV5"), formula = model )) p <- summary(lm(data = reg %>% filter(Key_challenge == "HH1"), formula = model )) q <- summary(lm(data = reg %>% filter(Key_challenge == "HH2"), formula = model )) r <- summary(lm(data = reg %>% filter(Key_challenge == "HH3"), formula = model )) s <- summary(lm(data = reg %>% filter(Key_challenge == "HH4"), formula = model )) t <- summary(lm(data = reg %>% filter(Key_challenge == "SOC1"), formula = model )) u <- summary(lm(data = reg %>% filter(Key_challenge == "SOC2"), formula = model )) v <- summary(lm(data = reg %>% filter(Key_challenge == "SOC3"), formula = model )) w <- summary(lm(data = reg %>% filter(Key_challenge == "SOC4"), formula = model )) x <- summary(lm(data = reg %>% filter(Key_challenge == "SOC5"), formula = model )) results <- list(a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r, s, t, u, v, w, x) nvar <- length(rownames(results[[1]]$coefficients))-1 nresults <- length(results) coeff <- c() pval <- c() for ( i in 1:nresults){ coeff <- c(coeff, as.numeric(results[[i]]$coefficients[1:nvar+1,"Estimate"])) pval <- c(pval, as.numeric(results[[i]]$coefficients[1:nvar+1,"Pr(>|t|)"])) } test_real_a <- data.frame( key_challenge = c(rep("FS1",nvar), rep("FS2",nvar), rep("FS3",nvar), rep("FS4",nvar), rep("ECON1",nvar), rep("ECON2",nvar), rep("ECON3",nvar), rep("ECON4",nvar), rep("ECON5",nvar), rep("ECON6",nvar), rep("ENV1",nvar), rep("ENV2",nvar), rep("ENV3",nvar), rep("ENV4",nvar), rep("ENV5",nvar), rep("HH1",nvar), rep("HH2",nvar), rep("HH3",nvar), rep("HH4",nvar), rep("SOC1",nvar), rep("SOC2",nvar), rep("SOC3",nvar), rep("SOC4",nvar), rep("SOC5",nvar)), reg_var = as.factor(c(rep(rownames(results[[1]]$coefficients)[-1], nresults))), reg_coeff = coeff, reg_p= pval ) test_real <- data.frame( key_challenge = c(rep("FS1",nvar), rep("FS2",nvar), rep("FS3",nvar), rep("FS4",nvar), rep("ECON1",nvar), rep("ECON2",nvar), rep("ECON3",nvar), rep("ECON4",nvar), rep("ECON5",nvar), rep("ECON6",nvar), rep("ENV1",nvar), rep("ENV2",nvar), rep("ENV3",nvar), rep("ENV4",nvar), rep("ENV5",nvar), rep("HH1",nvar), rep("HH2",nvar), rep("HH3",nvar), rep("HH4",nvar), rep("SOC1",nvar), rep("SOC2",nvar), rep("SOC3",nvar), rep("SOC4",nvar), rep("SOC5",nvar)), reg_var = as.factor(c(rep(rownames(results[[1]]$coefficients)[-1], nresults))), reg_dir = as.factor (ifelse(coeff>0, "Positive","Negative")), reg_p= as.factor( ifelse( pval>0.01, "not sign.", ifelse(pval<0.01 & pval>0.0021, "Bonf. A: 0.01","Bonf. B: 0.0021"))) ) test_real$key_challenge <- factor(test_real$key_challenge, levels = unique(test_real$key_challenge[order(test_real$key_challenge, decreasing = T)])) test_real2 <- test_real %>% mutate(reg_var = factor(reg_var, levels = c("cropArable crops", "cropHorticulture", "Potential.loss", "rel_attainable_yield", "risk_score", "fieldBioTech", "fieldToxi", "fieldEcol","fieldSoc-Eco", "fieldOther", "scopespecific","expertise.category", "importance.category", "implementation.level"))) ggplot(test_real2, aes(reg_var, key_challenge)) + geom_raster(aes( fill= reg_dir, alpha = reg_p), hjust =0, vjust =0.00) + geom_vline(xintercept=seq(0, 14, by=1))+ geom_hline(yintercept=seq(0, 24, by=1))+ scale_fill_manual(values = c("Positive"= "#b2df8a","Negative" = "#1f78b4"))+ scale_alpha_manual(values = c("not sign." = 0,"Bonf. A: 0.01"=0.5,"Bonf. B: 0.0021"=1))+ theme(axis.title = element_blank(), axis.text.x = element_text(angle = 45, vjust = -1, hjust=0, face = "bold"), axis.text.y = element_text( vjust = 1, hjust=0, face = "bold"))+ labs(fill = "Coefficient\nSign", alpha = "Coefficient\nP-value") + guides(fill = guide_legend(order=1), alpha = guide_legend(order=2))+ scale_y_discrete( expand = c(0,0))+ scale_x_discrete(position = "top", expand = c(0,0), labels = c("Type: Arable", "Type: Horticulture", "Pot. Pest Pressure", "Attainable Yield", "Level Pesticide Risk", "Field: BioTech", "Field: Toxicology", "Field: Ecology","Field: Socio-Econ", "Field: Other", "Scope: Local expert","Expertise in Category", "Importance of PM", "Level Sustainable PM" )) # Export regression tables library(stargazer) modela<- Importance_for_sustainable_pest_management ~ crop + Potential.loss + rel_attainable_yield + risk_score + field + scope + + expertise.category + importance.category +implementation.level a <- lm(data = reg %>% filter(Key_challenge == "FS1"), formula = modela ) b <- lm(data = reg %>% filter(Key_challenge == "FS2"), formula = modela ) c <- lm(data = reg %>% filter(Key_challenge == "FS3"), formula = modela ) d <- lm(data = reg %>% filter(Key_challenge == "FS4"), formula = modela ) e <- lm(data = reg %>% filter(Key_challenge == "ECON1"), formula = modela ) f <- lm(data = reg %>% filter(Key_challenge == "ECON2"), formula = modela ) g <- lm(data = reg %>% filter(Key_challenge == "ECON3"), formula = modela ) h <- lm(data = reg %>% filter(Key_challenge == "ECON4"), formula = modela ) i <- lm(data = reg %>% filter(Key_challenge == "ECON5"), formula = modela ) j <- lm(data = reg %>% filter(Key_challenge == "ECON6"), formula = modela ) k <- lm(data = reg %>% filter(Key_challenge == "ENV1"), formula = modela ) l <- lm(data = reg %>% filter(Key_challenge == "ENV2"), formula = modela ) m <- lm(data = reg %>% filter(Key_challenge == "ENV3"), formula = modela ) n <- lm(data = reg %>% filter(Key_challenge == "ENV4"), formula = modela ) o <- lm(data = reg %>% filter(Key_challenge == "ENV5"), formula = modela ) p <- lm(data = reg %>% filter(Key_challenge == "HH1"), formula = modela ) q <- lm(data = reg %>% filter(Key_challenge == "HH2"), formula = modela ) r <- lm(data = reg %>% filter(Key_challenge == "HH3"), formula = modela ) s <- lm(data = reg %>% filter(Key_challenge == "HH4"), formula = modela ) t <- lm(data = reg %>% filter(Key_challenge == "SOC1"), formula = modela ) u <- lm(data = reg %>% filter(Key_challenge == "SOC2"), formula = modela ) v <- lm(data = reg %>% filter(Key_challenge == "SOC3"), formula = modela ) w <- lm(data = reg %>% filter(Key_challenge == "SOC4"), formula = modela ) x <- lm(data = reg %>% filter(Key_challenge == "SOC5"), formula = modela ) results_s <- list(a, b, c, d, e, f, g, h, i, j, k, l, m, n, o, p, q, r, s, t, u, v, w, x) stargazer(results_s, type="html", star.cutoffs = c(0.05, 0.01, 0.001), single.row = T, covariate.labels = c("Type: Arable", "Type: Horticulture", "Pot. Pest Pressure", "Attainable Yield", "Level Pesticide Risk", "Field: BioTech", "Field: Ecology", "Field: Other","Field: Socio-Econ","Field: Toxicology", "Scope: Local expert","Expertise in Category", "Importance of PM", "Level Sustainable PM", "Intercept"), column.labels=c("FS1","FS2","FS3","FS4", "ECON1","ECON2","ECON3","ECON4","ECON5","ECON6", "ENV1","ENV2","ENV3","ENV4","ENV5", "HH1","HH2","HH3","HH4", "SOC1","SOC2","SOC3","SOC4","SOC5")) ### vi) Supplementary Information #### # Importance of pest management for categories: Supp. Fig. 1. imp<- datsur %>% mutate(row.n = row_number()) %>% select("Food Security"= Q4.FS.1.importance.pest.management, Environment= Q5.EE.1.importance.pest.management, "Human Health"= Q6.HH.1.importance.pest.management, Social= Q7.SSE.1.importance.pest.management, Economy= Q8.EC.1.importance.pest.management, overview_continent_all, Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, row.n) %>% pivot_longer(-c(overview_continent_all,Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, row.n), names_to = "Key_challenge",values_to="Importance_pest_management") exp_we<- datsur %>% mutate(row.n = row_number()) %>% select("Food Security"=Q2.1.PME.expertise.key.challenges.food.security,"Environment"=Q2.2.expertise.key.challenges.environmental.effects, "Human Health"=Q2.3.PME.expertise.key.challenges.health, "Social"=Q2.4.PME.expertise.key.challenges.social.security, "Economy"=Q2.5.PME.expertise.key.challenges.economy, row.n) %>% pivot_longer(-row.n,names_to = "Key_challenge",values_to="Strength_Expertise") imp_we <- left_join(imp, exp_we, by=c("row.n","Key_challenge")) %>% mutate(expertise.weight = Strength_Expertise/10) %>% group_by(Key_challenge) %>% mutate(mean= mw(x=Importance_pest_management, w=expertise.weight), se = sw(x=Importance_pest_management, w= expertise.weight), n=n()) ggplot(imp_we, aes(x= Key_challenge, y=Importance_pest_management, fill = Key_challenge))+ stat_summary(fun.y = mean, geom = "bar", fill= c("#332288", "#117733", "#44AA99", "#88CCEE", "#DDCC77")) + stat_summary(fun.data = mean_se, geom = "errorbar",alpha=1, color= "black", width=.3, size=1.2)+ geom_jitter( color = "gray",position = position_jitter(width = 0.1, height = 0.25),alpha=0.5)+ ylab("Importance of pest management per category")+ xlab("")+ ylim(0,10)+ theme(axis.text.x= element_text(size=12, angle = 15), axis.title = element_text(size=12, face= "bold"))+ theme(legend.position = "none") # Sample properties : Suppl. Fig. 2 dat <- read.csv("cords_local_experts.csv") spamcrop <- read.csv("spam2020V1r0_global_H_TA.csv") qq = quantile(spamcrop$MAIZ_A, probs = seq(0, 1, .2)) qq spamcrop <- spamcrop %>% filter(as.numeric(MAIZ_A)> 0) %>% mutate(Maize_ha = ifelse(as.numeric(MAIZ_A)< 7.6, "very low", ifelse(as.numeric(MAIZ_A)< 68.1, "low", ifelse(as.numeric(MAIZ_A)< 251.2, "medium", ifelse(as.numeric(MAIZ_A)< 818.7, "high","very high")) ))) %>% mutate(Maize_ha = factor(Maize_ha, levels=c("very low","low","medium","high","very high"))) library(ggplot2) ggplot() + borders("world", fill = "white", colour = "grey80")+ geom_point(data = spamcrop , aes(x = x, y = y , colour = Maize_ha) , size = 0.0003)+ scale_color_brewer(palette= "Greens")+ guides(colour = guide_legend(override.aes = list(size=10)))+ geom_point(data = dat , aes(x = lon, y = lat) , size = 2, color = "deepskyblue2") # Plots of individual responses per indicator: Suppl. Fig. 3-7 datsur <- datsur %>% mutate(survey_type_short = ifelse(survey_type == "Literature","Literature","Networks_organisations")) chall <- datsur %>% select(FS1 = Q4.FS.2.1.sufficient.quantity, FS2 = Q4.FS.2.2.healthy.food, FS3 = Q4.FS.2.3.safe.food.and.feed, FS4 = Q4.FS.2.4.resilience.of.food, ENV1 = Q5.EE.2.1.pollution.drinking.water, ENV2 = Q5.EE.2.2.pollution.soil, ENV3 = Q5.EE.2.3.pollution.marine.ecosystems, ENV4 = Q5.EE.2.4.pollution.freshwater.ecosystems, ENV5 = Q5.EE.2.5.loss.of.genes.biodiversity, HH1 = Q6.HH.2.1.effects.farm.workers, HH2 = Q6.HH.2.2.effects.bystanders, HH3 = Q6.HH.2.3.effects.on.consumers, HH4= Q6.HH.2.4.restore.trust.in.food, SOC1 = Q7.SSE.2.1.equal.distrib.local.actors, SOC2 = Q7.SSE.2.2.equal.distrib.global.actors, SOC3 = Q7.SSE.2.3.adapt.to.future.pest.pressure, SOC4 = Q7.SSE.2.4.access.to.education.tools.technology, SOC5 = Q7.SSE.2.5.safe.working.conditions, ECON1 = Q8.EC.2.1.efficient.affordable.PM, ECON2 = Q8.EC.2.2.short.term.livelihood.keeping, ECON3 = Q8.EC.2.3.long.term.livelihood.keeping, ECON4 = Q8.EC.2.4.economic.resilience.extreme.pest.events, ECON5 = Q8.EC.2.5.reduction.threat.of.pests.in.productivity, ECON6 = Q8.EC.2.6.reduction.of.indirect.costs.of.PM, overview_continent_all, Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, survey_type_short) %>% pivot_longer(-c(overview_continent_all,Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, survey_type_short), names_to = "Key_challenge",values_to="Importance_for_sustainable_pest_management") %>% mutate(category_challenge = ifelse(str_detect(Key_challenge,"FS"),"Food Security", ifelse(str_detect(Key_challenge,"ENV"),"Environment", ifelse(str_detect(Key_challenge,"HH"),"Human Health", ifelse(str_detect(Key_challenge,"SOC"),"Social", ifelse(str_detect(Key_challenge,"ECON"),"Economy","Other")))))) #3 ggplot(chall , aes(x= Key_challenge, y= Importance_for_sustainable_pest_management))+ ## add half-violin from {ggdist} package ggdist::stat_halfeye( fill="darkgrey", ## custom bandwidth adjust = .9, ## adjust if neccesary width = .95, ## move geom if necessary justification = 0, ## remove slab interval .width = 0, point_colour = NA) + # geom_jitter(position = position_jitter(width = 0.1, height = 0.25),alpha=0.5, color="grey")+ xlab("Indicators assessment framework")+ ylab("Expert estimation: expected effect of transformation")+ theme(axis.text.x= element_text( angle = 20), axis.title = element_text(size=12, face= "bold"))+ facet_wrap(~category_challenge, scales = "free_x") #4 s <- function(x){plotrix::std.error(x,na.rm=T)} m <- function(x){mean(x,na.rm=T)} ggplot(chall %>% group_by(category_challenge, Key_challenge) %>% summarise(mean= m(x=Importance_for_sustainable_pest_management), se = s(x=Importance_for_sustainable_pest_management), n=n()), aes(x= Key_challenge, y= mean))+ geom_errorbar(aes(ymin=mean-se*qt(p=.05/2, df=n-1, lower.tail=FALSE), ymax= mean+se*qt(p=.05/2, df=n-1, lower.tail=FALSE)), alpha=1, width=.2, position=position_dodge(0.2))+ geom_point(position=position_dodge(0.2))+ xlab("Indicators assessment framework")+ ylab("Expert estimation: expected effect of transformation (unweighted)")+ ylim(-10,10) + geom_hline(yintercept = 0, linetype= "dashed", color = "grey")+ theme(axis.text.x= element_text( angle = 15), axis.title = element_text(size=12, face= "bold"))+ facet_wrap(~category_challenge, scales = "free_x") #5 ggplot(chall %>% group_by(category_challenge, Key_challenge, field.expertise.short) %>% summarise(mean= m(x=Importance_for_sustainable_pest_management), se = s(x=Importance_for_sustainable_pest_management), n=n()), aes(x= Key_challenge, y= mean, color= field.expertise.short))+ geom_errorbar(aes(ymin=mean-se*qt(p=.05/2, df=n-1, lower.tail=FALSE), ymax= mean+se*qt(p=.05/2, df=n-1, lower.tail=FALSE)), alpha=1, width=.2, position=position_dodge(0.2))+ geom_point(position=position_dodge(0.2))+ xlab("Indicators assessment framework")+ ylab("Expert estimation: expected effect of transformation (unweighted)")+ ylim(-10,10) + geom_hline(yintercept = 0, linetype= "dashed", color = "grey")+ theme(axis.text.x= element_text( angle = 15), axis.title = element_text(size=12, face= "bold"))+ facet_wrap(~category_challenge, scales = "free_x")+ labs(color = "Research \nField") #6 ggplot(chall , aes(x= Key_challenge, y= Importance_for_sustainable_pest_management, color = field.expertise.short))+ geom_jitter(position = position_jitter(width = 0.3, height = 0.25),alpha=0.5)+ xlab("Indicators assessment framework")+ ylab("Expert estimation: expected effect of transformation")+ theme(axis.text.x= element_text( angle = 20), axis.title = element_text(size=12, face= "bold"))+ facet_wrap(~category_challenge, scales = "free_x")+ scale_color_discrete(name="Research\nField")+ theme(legend.justification=c(1,0), legend.position=c(1,0)) #7 ggplot(chall %>% group_by(category_challenge, Key_challenge, survey_type_short) %>% summarise(mean= m(x=Importance_for_sustainable_pest_management), se = s(x=Importance_for_sustainable_pest_management), n=n()), aes(x= Key_challenge, y= mean, color= survey_type_short))+ geom_errorbar(aes(ymin=mean-se*qt(p=.05/2, df=n-1, lower.tail=FALSE), ymax= mean+se*qt(p=.05/2, df=n-1, lower.tail=FALSE)), alpha=1, width=.2, position=position_dodge(0.2))+ geom_point(position=position_dodge(0.2))+ xlab("Indicators assessment framework")+ ylab("Expert estimation: expected effect of transformation (unweighted)")+ ylim(-10,10) + geom_hline(yintercept = 0, linetype= "dashed", color = "grey")+ theme(axis.text.x= element_text( angle = 15), axis.title = element_text(size=12, face= "bold"))+ facet_wrap(~category_challenge, scales = "free_x")+ labs(color = "Type\nSurvey\nInvitation") # the share of negative responses per indicator share_ind<- chall %>% group_by(Key_challenge) %>% summarise(pos = sum(Importance_for_sustainable_pest_management>-1, na.rm=T), neg = sum(Importance_for_sustainable_pest_management<0, na.rm=T), share= round(sum(Importance_for_sustainable_pest_management<0, na.rm=T)/sum(Importance_for_sustainable_pest_management>-11, na.rm=T),2)) # The share of respondents with at least one negative response datsur <- datsur %>% mutate(survey_type_short = ifelse(survey_type == "Literature","Literature","Networks_organisations")) chall2 <- datsur %>% select(FS1 = Q4.FS.2.1.sufficient.quantity, FS2 = Q4.FS.2.2.healthy.food, FS3 = Q4.FS.2.3.safe.food.and.feed, FS4 = Q4.FS.2.4.resilience.of.food, ENV1 = Q5.EE.2.1.pollution.drinking.water, ENV2 = Q5.EE.2.2.pollution.soil, ENV3 = Q5.EE.2.3.pollution.marine.ecosystems, ENV4 = Q5.EE.2.4.pollution.freshwater.ecosystems, ENV5 = Q5.EE.2.5.loss.of.genes.biodiversity, HH1 = Q6.HH.2.1.effects.farm.workers, HH2 = Q6.HH.2.2.effects.bystanders, HH3 = Q6.HH.2.3.effects.on.consumers, HH4= Q6.HH.2.4.restore.trust.in.food, SOC1 = Q7.SSE.2.1.equal.distrib.local.actors, SOC2 = Q7.SSE.2.2.equal.distrib.global.actors, SOC3 = Q7.SSE.2.3.adapt.to.future.pest.pressure, SOC4 = Q7.SSE.2.4.access.to.education.tools.technology, SOC5 = Q7.SSE.2.5.safe.working.conditions, ECON1 = Q8.EC.2.1.efficient.affordable.PM, ECON2 = Q8.EC.2.2.short.term.livelihood.keeping, ECON3 = Q8.EC.2.3.long.term.livelihood.keeping, ECON4 = Q8.EC.2.4.economic.resilience.extreme.pest.events, ECON5 = Q8.EC.2.5.reduction.threat.of.pests.in.productivity, ECON6 = Q8.EC.2.6.reduction.of.indirect.costs.of.PM, id= X,overview_continent_all, Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, survey_type_short) %>% pivot_longer(-c(id, overview_continent_all,Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, survey_type_short), names_to = "Key_challenge",values_to="Importance_for_sustainable_pest_management") %>% mutate(category_challenge = ifelse(str_detect(Key_challenge,"FS"),"Food Security", ifelse(str_detect(Key_challenge,"ENV"),"Environment", ifelse(str_detect(Key_challenge,"HH"),"Human Health", ifelse(str_detect(Key_challenge,"SOC"),"Social", ifelse(str_detect(Key_challenge,"ECON"),"Economy","Other")))))) share_id<- chall2 %>% group_by(id) %>% summarise(pos = sum(Importance_for_sustainable_pest_management>-1, na.rm=T), neg = sum(Importance_for_sustainable_pest_management<0, na.rm=T), share= sum(Importance_for_sustainable_pest_management<0, na.rm=T)/sum(Importance_for_sustainable_pest_management>-11, na.rm=T)) share_id %>% mutate(share_tot= sum(share>0, na.rm=T)/n(), n= n()) # Analyse responses solution part: Supplementary Fig. 8 library(dplyr) library(ggplot2) library(tidyr) library(stringr) library(Hmisc) library(diagis) library(zoo) library(plotrix) datsur <- datsur %>% mutate(Q1.4.PME.crops.expertise = replace(Q1.4.PME.crops.expertise,Q1.4.PME.crops.expertise== "Horticulture (vegetables, wine and fruits)","Horticulture")) datsur$Q1.6.PME.main.field.expertise[datsur$Q1.6.PME.main.field.expertise=="Bio/Tech (developing mechanical, synthetic, biological crop management tools)"] <- "Bio/Tech" datsur <- datsur %>% mutate(field.expertise.short = ifelse(Q1.6.PME.main.field.expertise %in% c("Agronomy", "Bio/Tech", "Entomology","Plant Pathology"),"PMS", ifelse(Q1.6.PME.main.field.expertise %in% c("(Agro)Ecology"),"Ecol", ifelse(Q1.6.PME.main.field.expertise %in% c("(Eco)Toxicity","Environmental Sciences","Human Health"),"Toxi", ifelse(Q1.6.PME.main.field.expertise %in% c("Socio-Economics"),"Soc-Eco","Other"))))) datsur <- datsur %>% mutate(profession.short = ifelse(Q1.5.PME.profession.description %in% c("Research"),"Research", ifelse(Q1.5.PME.profession.description %in% c("Policy"),"Policy","Production"))) solu <- datsur %>% mutate(row.n = row_number()) %>% select( Awareness = Q9.KS.1.1.increase.awareness.of.neg.impacts.of.PM, "Risk management" = Q9.KS.1.2.instruments.to.manage.higher.risk.of.sustainable.production, "Econ Support" = Q9.KS.1.3.provision.of.econ.support, Education = Q9.KS.1.4.provision.of.education, Legislation = Q9.KS.1.5.adjustment.legislation, Markets = Q9.KS.1.6.abolish.barriers.trade.market.access.for.sustainable.food, "Value-Chain" = Q9.KS.1.7.support.food.value.chain.actors, Substitutes = Q9.KS.1.8.pesticide.substitutes, overview_continent_all, pesticide_use_ha, pesticide_use_outputvalue, risk_score, HDI, outputvalue, cropland, schooling, SDG_2.3_agricult_prod, SDG_2A_agricult_invest,ai_count, ai_bans, organic_share, rel_attainable_yield, rel_yield_gap, Actual.loss,Efficiency.PM,Potential.loss, gfsi22, Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, Q3.ISPM.sustainable.pest.management.in.your.crops, row.n) %>% pivot_longer(-c(overview_continent_all, pesticide_use_ha, pesticide_use_outputvalue, risk_score, HDI, outputvalue, cropland, schooling, SDG_2.3_agricult_prod, SDG_2A_agricult_invest,ai_count, ai_bans, organic_share, rel_attainable_yield, rel_yield_gap, Actual.loss,Efficiency.PM,Potential.loss, gfsi22, Q1.4.PME.crops.expertise, Q1.1.PME.specify.larger.vs.specific, field.expertise.short, profession.short, Q3.ISPM.sustainable.pest.management.in.your.crops, row.n), names_to = "Solution",values_to="Importance_solution") %>% distinct() fac_s<- solu %>% group_by(Solution, overview_continent_all) %>% summarise(mean= m(x= Importance_solution), se = s(x= Importance_solution), n=n()) %>% filter (overview_continent_all =="Europe") %>% arrange(mean) ggplot(solu %>% group_by(Solution) %>% summarise(mean= m(x= Importance_solution), se = s(x= Importance_solution), n=n()) %>% mutate(Solution=factor(Solution, levels=fac_s$Solution)), aes(x= Solution, y= mean))+ geom_errorbar(aes(ymin=mean-se*qt(p=.05/2, df=n-1, lower.tail=FALSE), ymax= mean+se*qt(p=.05/2, df=n-1, lower.tail=FALSE)), alpha=1, width=.2, position=position_dodge(0.5))+ geom_point(position=position_dodge(0.5))+ xlab("Categories Potential Solutions")+ ylab("Weight solution (0-100 points)")+ theme(axis.text.x= element_text( angle = 15, size= 12), axis.title = element_text(size=12, face= "bold")) # Supplementary Fig. 9 ggplot(as.data.frame(table(datsur$overview_continent_all)) %>% arrange(Freq), aes(x= Var1, y=Freq))+ geom_bar(stat="identity", width=0.7)+ geom_text(aes(label= Freq), vjust=1.6, position = position_dodge(0.9), size=5)+ labs(fill = "Expertise: \ncontinent")+ ylab("Number of respondents")+ # xlab("Expertise: specific vs. larger (continental/global)")+ xlab("")+ theme_bw()+ theme(axis.text.x= element_text(size=11, angle = 90, vjust=-0.0001), legend.position='none') # Supplementary Fig. 10 ggplot(as.data.frame(table(datsur$Q1.4.PME.crops.expertise)), aes(x= Var1, y=Freq))+ geom_bar(stat="identity", width=0.7)+ geom_text(aes(label= Freq), vjust=1.6, position = position_dodge(0.9), size=5)+ labs(fill = "Expertise: \ncrop type")+ ylab("Number of respondents")+ # xlab("Expertise: specific vs. larger (continental/global)")+ xlab("")+ theme_bw()+ theme(axis.text.x= element_text(size=11, angle = 90, vjust=-0.0001), legend.position='none') # Supplementary Fig. 11 library(ggthemes) expert <- distinct(chall2 %>% select(id,field.expertise.short)) ggplot(as.data.frame(table(expert$field.expertise.short)) , aes(x= Var1, y=Freq, fill=Var1))+ scale_fill_colorblind()+ geom_bar(stat="identity", width=0.7)+ geom_text(aes(label= Freq), vjust=1.6, position = position_dodge(0.9), size=4)+ ylab("Number of respondents")+ labs(fill = "Research field")+ # xlab("Expertise: specific vs. larger (continental/global)")+ xlab("")+ theme_bw()+ theme(axis.text.x= element_blank()) # Supplementary Fig. 12 library(ggdist) exp<- datsur %>% select(Food_security=Q2.1.PME.expertise.key.challenges.food.security,Environment=Q2.2.expertise.key.challenges.environmental.effects, Human_Health=Q2.3.PME.expertise.key.challenges.health, Social=Q2.4.PME.expertise.key.challenges.social.security, Economy=Q2.5.PME.expertise.key.challenges.economy) %>% pivot_longer(everything(),names_to = "Field_Expertise",values_to="Strength_Expertise") ggplot(exp , aes(x= Field_Expertise, y=Strength_Expertise))+ ## add half-violin from {ggdist} package ggdist::stat_halfeye( ## custom bandwidth adjust = .5, ## adjust height width = .6, ## move geom to the right justification = -.2, ## remove slab interval .width = 0, point_colour = NA) + geom_jitter(position = position_jitter(width = 0.2, height = 0.1))+ xlab("Field of Expertise")+ ylab("Indicated strength of expertise")+ theme(axis.text.x= element_text(size=12, angle = 15), axis.title = element_text(size=12, face= "bold")) #14 chall2 <- chall %>% filter(!(overview_continent_all %in% c("NA", "Global"))) chall3 <- chall2[!is.na(chall2$overview_continent_all),] ggplot(chall3 , aes(x= Key_challenge, y= Importance_for_sustainable_pest_management, color = overview_continent_all))+ geom_jitter(position = position_jitter(width = 0.3, height = 0.25),alpha=0.5)+ xlab("Indicators assessment framework")+ ylab("Expert estimation: expected effect of transformation")+ theme(axis.text.x= element_text( angle = 20), axis.title = element_text(size=12, face= "bold"))+ facet_wrap(~category_challenge, scales = "free_x")+ scale_color_discrete(name="Region")+ theme(legend.justification=c(1,0), legend.position=c(1,0))