Supplemental.Rmd

---
title: "Rethinking Stakeholder Engagement: A Multidimensional Framework Inspired by Meta-Analysis of Food, Energy, and Water Research"
author: "-"
date: "`r Sys.Date()`"
output:
  pdf_document:
    toc: yes

---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```

\newpage

# Supplemental Materials Summary

The benefits for involving stakeholders in research projects are well established.  We investigated recent papers related to food-energy-water systems to assess whether stakeholders were engaged, when they were engaged, which stakeholders were involved, the level of their involvement, whether the needs of stakeholders were considered, and whether stakeholder engagement was evaluated and monitored. Stakeholder engagement theory and research suggest that stakeholder engagement is effective only when it is considered as an integrated process, but our review suggests that this is not common practice.  We developed a framework to assist scientists in engaging stakeholders in an integrated and thoughtful manner. Our framework includes six dimensions:  1. When should stakeholders be engaged and for how long? 2. Who should be engaged (identification and diversity)? 3. What role should stakeholders play? 4. How can researchers meet stakeholder needs?   5. What methods of engagement should be used?  6. How to evaluate and monitor the engagement? 

We conducted two literature searches: an initial search in 2020 and a follow-up search in 2023 to capture literature published between 2020 and 2023. We conducted the searches in two different online databases, ScienceDirect and WorldCat, to ensure a comprehensive identification of relevant literature. We identified 177 publications from our initial search and 540 publications in our follow up search, resulting in a total of 717 publications. We then manually screened the papers and removed additional papers that did not meet the criteria above on reading the manuscript, resulting in 489 total papers for analysis.

This supplemental materials provides:

* A summary of the variables quantified from the literature reviews;
* Chi-square testing of whether a solution was proposed or not vs. 1) stakeholder engagement as well as 2) whether a computational model was used;
* Summary statistics graphs of the literature analysis; and
* Regression and odds ratio testing for a variety of factors


```{r echo=FALSE, warning=FALSE, message=FALSE}

library(ggplot2)
library(epiDisplay)

crcdata <- read.csv("./data/crc_data_ES2.csv")


#set YN to 01
require(dplyr)
crcdata <- crcdata %>%
  mutate(solution_proposed_YN = ifelse(solution_proposed_YN == "N",0,1))

require(dplyr)
crcdata <- crcdata %>%
  mutate(solution_implemented_YN = ifelse(solution_implemented_YN == "N",0,1))

require(dplyr)
crcdata <- crcdata %>%
  mutate(S_stakeholder_engagment_YN = ifelse(S_stakeholder_engagment_YN == "N",0,1))

#fix colnames
names(crcdata)[names(crcdata) == 'STE_IAP2_NA'] <- 'STE_IAP2_data_gathering'
names(crcdata)[names(crcdata) == 'G_na'] <- 'G_nogeography'

names(crcdata)[43] <- "SC_researcher"

names(crcdata)[51] <- "G_notdescribed"


```


\newpage

# Variable Summary

Below is a list of the categorical variables generated from the literature reviews

```{r echo=FALSE, warning=FALSE, message=FALSE}
library(knitr)
library(kableExtra)
df <- data.frame(Variables = c("Year", "Solution Proposed", "Solution Implemented", "Solution Type", "Computational Model Used", "Researcher Type", "Stakeholder Type", "Stakeholder Engagment Scale - Ghodsvali", "Stakeholder Engagment Scale - IAP2", "Stakeholder Engagement Scale - Local", "Geographical Type"), Description = linebreak(c("Year of citation", "Was a solution proposed?", "Was a solution implemented?","If a solution was proposed, what was the solution type? Groups include: Technology, Policy, Institutional, Social, Economic, Ecological, and Educational.", "Was a computational model used?", "What was the research type?  Groups include: NGO, English, Math, Computer Science, Physics, Engineering, Interdisciplinary, Social Science, Economics, Agriculture, and Other","What was the stakeholder type? Groups include: Farmers, Combined Government, Combined Coalition, Combined Industry, Migrants, Youth, Public, Univerity, and Experts", "If a stakeholder was engaged, categorization of the engagement using the Ghodsvali scale. Groups include: Nominal, Instrumental, Representation, and Transformative", "If a stakeholder was engaged, categorization of the engagement using the IAP2 scale. Groups include: Data Gathering, Inform, Consult, Involve, Collaborate, and Empower", "If a stakeholder was engaged, categorization of the engagement using a customized scale. Groups include: Researcher, Data Gathering, Inform, Perspectives, Planning, Identify, Envision, and Implement", "What the geography type? Groups include: Not Described, Local, Regional, National, Multinational, Global, and No Geography")))
kable(df, col.names = c("Variable Name", "Description"), escape = F, caption = "Table T1: Variable Descriptions") %>% column_spec(2, width = "2in") %>% column_spec(c(1:2), width = "2.5in") %>% kable_styling(full_width = FALSE,latex_options = c("hold_position", "scale_down"), font_size = 8) %>% collapse_rows(columns = 2) %>% kable_classic_2(full_width = F) %>% row_spec(0,bold=TRUE)

```


<!-- \newpage -->

<!-- # Variable Summary REFINED -->

<!-- Below is a list of the categorical variables generated from the literature reviews REFINED -->

<!-- ```{r echo=FALSE, warning=FALSE, message=FALSE} -->
<!-- library(knitr) -->
<!-- library(kableExtra) -->
<!-- df <- data.frame(Variables = c("Year", "Solution Proposed", "Solution Implemented",  "Stakeholder Type", "Stakeholder Engagment Scale - Ghodsvali", "Stakeholder Engagment Scale - IAP2", "Stakeholder Engagement Scale - Local"), Description = linebreak(c("Year of citation", "Was a solution proposed?", "Was a solution implemented?","What was the stakeholder type? Groups include: Farmers, Combined Government, Combined Coalition, Combined Industry, Migrants, Youth, Public, Univerity, and Experts", "If a stakeholder was engaged, categorization of the engagement using the Ghodsvali scale. Groups include: Nominal, Instrumental, Representation, and Transformative", "If a stakeholder was engaged, categorization of the engagement using the IAP2 scale. Groups include: Data Gathering, Inform, Consult, Involve, Collaborate, and Empower", "If a stakeholder was engaged, categorization of the engagement using a customized scale. Groups include: Researcher, Data Gathering, Inform, Perspectives, Planning, Identify, Envision, and Implement"))) -->
<!-- kable(df, col.names = c("Variable Name", "Description"), escape = F, caption = "Table T1: Variable Descriptions") %>% column_spec(2, width = "2in") %>% column_spec(c(1:2), width = "2.5in") %>% kable_styling(full_width = FALSE,latex_options = c("hold_position", "scale_down"), font_size = 8) %>% collapse_rows(columns = 2) %>% kable_classic_2(full_width = F) %>% row_spec(0,bold=TRUE)  -->

<!-- ``` -->


\newpage

# Chi-Square Testing
## Chi Square Testing: solution proposed or not vs. vs. stakeholder engagement

Chi Square and Fishers Exact Test on contingency table with Solution/No Solution as the explanatory variable, and engaged stakeholder/did not engage stakeholder as the response variable. 

**ChiSquare = 26: Fishers Exact Test Odds Ratio: 10: Not Independent**

Both chi square and fishers exact test were significant, with a chi square approximation of ~26, which is well above the critical value (3.84 with one degree of freedom).  Fishers Exact Test returned an odds ratio of ~10. The alternative hypothesis: true odds ratio is not equal to 1, therefore the null hypothesis is rejected - the groups are not independent. 

The Fishers Exact Test defaults to associating the odds ratio (which can represent effect size) with the first cell. In this instance "The odds of having a solution is 10 times that for an engaged stakeholder".  You could flip the response and explanatory variables, but the odds ratio would stay the same.

For more info on this topic see: Kim HY. Statistical notes for clinical researchers: Chi-squared test and Fisher's exact test. Restor Dent Endod. 2017 May;42(2):152-155. doi: 10.5395/rde.2017.42.2.152. Epub 2017 Mar 30. PMID: 28503482; PMCID: PMC5426219.


```{r, echo=FALSE}
library(Barnard)

crcdata$solution_proposed_YN <- as.factor(crcdata$solution_proposed_YN)
crcdata$S_stakeholder_engagment_YN <- as.factor(crcdata$S_stakeholder_engagment_YN)
crcdata$solution_implemented_YN <- as.factor(crcdata$solution_implemented_YN)

levels(crcdata$S_stakeholder_engagment_YN) = c("N", "Y")
levels(crcdata$S_stakeholder_engagment_YN) = c("N", "Y")
levels(crcdata$solution_proposed_YN) = c("N", "Y")
levels(crcdata$solution_implemented_YN) = c("N", "Y")

crcdata_table <- table(crcdata$solution_proposed_YN, crcdata$S_model_YN)
colnames(crcdata_table) <- c("M", "NM")

#custom add
crcdata_table <- rbind(crcdata_table[2,], crcdata_table[1,])
crcdata_table[1,] <- rev(crcdata_table[1,])

solution_stakeholder <- crcdata_table

 # solution_stakeholder <- 
 #   as.matrix(
 #     read.table(text=
 #                  "  E     NE           
 #       S             14    4      
 #       NS            76    400",
 #                header=TRUE,
 #                row.names=1)
 #   )

names(dimnames(solution_stakeholder)) <- c('solution','stakeholder')
solution_stakeholder

solution_stakeholder <- as.table(solution_stakeholder)

summary(solution_stakeholder)

fisher.test(solution_stakeholder)

btest <- c(solution_stakeholder[1,1:2], solution_stakeholder[2,1:2])

barnard.test(13,5,95,370)

```

\newpage

## Chi Square Testing: solution proposed or not vs. computational model used or not? 

Chi Square and Fishers Exact Test on the contingency table with solution proposed (solution/no solution) as the explanatory variable, and whether a computational model used (model/no model) as the response variable.

**ChiSquare = .57: Fishers Exact Test Odds Ratio = .668: Independent**

Both chi square and fishers exact test (FET) were insignificant/borderline, with a chi square approximation of ~.57, which is well below to the critical value (3.84 for one degree of freedom).  FET returned an odds ratio of under 1, therefore the null is accepted - the groups are independent.

The FET defaults to associating the odds ratio (which can represent effect size) with the first cell. In this instance "The odds of having a solution is .67 times that for having a model".  You could flip the response and explanatory variables, but the odds ratio would stay the same.

```{r, echo=FALSE}
library(Barnard)

crcdata$S_model_YN <- as.factor(crcdata$S_model_YN)

crcdata$S_model_YN <- factor(crcdata$S_model_YN, levels = c("Y", "N"))

crcdata_table <- table(crcdata$solution_proposed_YN, crcdata$S_model_YN)
colnames(crcdata_table) <- c("M", "NM")

crcdata_table <- rbind(crcdata_table[2,], crcdata_table[1,])
#crcdata_table[2,] <- rev(crcdata_table[2,])


solution_model <- 
 crcdata_table

 
 # solution_model <- 
 #   as.matrix(
 #     read.table(text=
 #                  "  M     NM      
 #       S             11    7      
 #       NS            377   94",
 #                header=TRUE,
 #                row.names=1)
 #   )

names(dimnames(solution_model)) <- c('solution','model')
solution_model

solution_model <- as.table(solution_model)

summary(solution_model)

fisher.test(solution_model)

barnard.test(13,5,370,95)
```

\newpage

# Summary Statistics Graphs

## Were solutions proposed in the set of all papers?

```{r echo=FALSE, warning=FALSE, message=FALSE}
ggplot(crcdata,aes(x=solution_proposed_YN, fill=year)) + ggtitle("Solution Proposed: Y=18") + labs(x="Was a solution proposed?") +
  geom_bar(stat = "count", position = "dodge", fill="#99CC99") + theme_minimal() + theme_minimal()+
  theme(text = element_text(),
        axis.text.x = element_text(angle = 0, hjust = 1)) 

```
\newpage

## Were solutions implemented in the set of all papers?

```{r echo=FALSE, warning=FALSE, message=FALSE}

ggplot(crcdata,aes(x=solution_implemented_YN, fill=year)) + ggtitle("Solution Implemented: Y=11") + labs(x="Was a solution implemented?") + geom_bar(stat = "count", position = "dodge", fill="#99CC99") + theme_minimal()+
  theme(text = element_text(),
        axis.text.x = element_text(angle = 0, hjust = 1)) 

```
\newpage

## What were the solution types?

```{r echo=FALSE, warning=FALSE, message=FALSE}

# par(mar = c(8,4,4,4))
# plot3 <- barplot(colSums(crcdata[,4:10]), las=2, main="Solution Types: N=18", names.arg = c("Technological", "Policy", "Institutional", "Social", "Economic", "Ecological", "Educational"))

plot3_data <- colSums(crcdata[,4:10])
p3a <- as.data.frame(plot3_data)
p3a$group <- c("Technological", "Policy", "Institutional", "Social", "Economic", "Ecological", "Educational")

# 
ggplot(p3a, aes(x=group, y=plot3_data)) + geom_bar(stat="identity",fill="#99CC99") + theme_bw() + theme(axis.text.x = element_text(angle = 90, vjust = 1, hjust=1)) +
  xlab("Solution Types") + ylab("Count")


```
\newpage

## Was a computational model used?

```{r echo=FALSE, warning=FALSE, message=FALSE}
ggplot(crcdata,aes(x=S_model_YN, fill=year)) + ggtitle("Computational Model Used: Y=366") + labs(x="Was a computational model used?") +
  geom_bar(stat = "count", position = "dodge", fill="#99CC99") + theme_minimal()
 
``` 
\newpage

## Were stakeholders engaged?

```{r echo=FALSE, warning=FALSE, message=FALSE}

# nrow(subset(crcdata, solution_proposed_YN == "1"))/nrow(subset(crcdata, solution_proposed_YN == "0"))
# 
# nrow(subset(crcdata, solution_implemented_YN == "1"))/nrow(subset(crcdata, solution_implemented_YN == "0"))


#stakeholder
#par(mar = c(4,4,4,4))
#barplot(table(crcdata$S_stakeholder_engagment_YN), main = "Stakeholder Engagment: Y=92")

ggplot(crcdata,aes(x=S_stakeholder_engagment_YN, fill=year)) + ggtitle("Stakeholder Engagement: YES=92") + labs(x="Were stakeholders engaged?") +
  geom_bar(stat = "count", position = "dodge", fill="#99CC99") + theme_minimal()

```
\newpage

## All FEWS papers by year

```{r echo=FALSE, warning=FALSE, message=FALSE}
#year
crcdata_remove2011 <- subset(crcdata, year != "2011")
ggplot(crcdata_remove2011, aes(x=year)) + 
    geom_bar(stat="count")+theme_bw() + scale_x_continuous(breaks = crcdata$year) +
  geom_bar(stat = "count", position = "dodge", fill="#99CC99") + theme_minimal()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1))


```


```{r echo=FALSE, warning=FALSE, message=FALSE, results=FALSE}
#year

png("./images/fews_year.png", res=300, width = 1500, height = 1500 )
ggplot(crcdata_remove2011, aes(x=year)) + 
    geom_bar(stat="count")+theme_bw() + scale_x_continuous(breaks = crcdata$year) +
  geom_bar(stat = "count", position = "dodge", fill="#99CC99") + theme_minimal() + theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1))
dev.off()


```


\newpage




## Level of stakeholder engagement by year - Ghodsvali scale

```{r echo=FALSE, warning=FALSE, message=FALSE}
#year
library(dplyr)

crcdata_new <- crcdata %>% 
  mutate(ghodsvali = sprintf("%d_%d_%d_%d", STE_G_nominal, STE_G_instrumental, STE_G_representation, STE_G_transformative))

crcdata_new <- crcdata_new[!(crcdata_new$ghodsvali %in% "0_0_0_0"),]

library(dplyr)
crcdata_new <- crcdata_new %>%
    mutate(ghodsvali = recode(ghodsvali, '0_0_0_0' = 'none', '1_0_0_0' = 'nominal', '0_1_0_0' =  'instrumental', '0_0_1_0' = 'representative',  '0_0_0_1' = 'transformative' ))

crcdata_new$ghodsvali <- as.factor(crcdata_new$ghodsvali)


crcdata_new$ghodsvali <- factor(crcdata_new$ghodsvali, levels = c("nominal", "instrumental", "representative", "transformative"))

my_colors <- c("#99CC99", "#669999", "#6666CC", "#CC9966")

ggplot(crcdata_new, aes(fill=ghodsvali, x=year)) + 
    geom_bar(position="stack", stat="count")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + scale_x_continuous(breaks = crcdata_new$year) + scale_fill_manual(values = my_colors) + scale_fill_discrete(name = "Ghodsvali Scale")

```


```{r echo=FALSE, warning=FALSE, message=FALSE, results=FALSE}

png("./images/fews_ghodsvali.png", res=300, width = 1500, height = 1500 )
ggplot(crcdata_new, aes(fill=ghodsvali, x=year)) + 
    geom_bar(position="stack", stat="count")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + scale_x_continuous(breaks = crcdata_new$year) + scale_fill_manual(values = my_colors) + scale_fill_discrete(name = "Ghodsvali Scale")
dev.off()

```

\newpage

## Level of stakeholder engagement by year - IAP2 scale

```{r echo=FALSE, warning=FALSE, message=FALSE}
#year
library(dplyr)

crcdata_new <- crcdata %>% 
  mutate(IAP2 = sprintf("%d_%d_%d_%d_%d_%d", STE_IAP2_data_gathering,STE_IAP2_inform,STE_IAP2_consult,STE_IAP2_involve,STE_IAP2_collab,STE_IAP2_empower))

crcdata_new <- crcdata_new[!(crcdata_new$IAP2 %in% "0_0_0_0_0_0"),]

library(dplyr)
crcdata_new <- crcdata_new %>%
    mutate(IAP2 = recode(IAP2, '1_0_0_0_0_0' = 'data gathering', '0_1_0_0_0_0' = 'inform', '0_0_1_0_0_0' =  'consult', '0_0_0_1_0_0' = 'involve',  '0_0_0_0_1_0' = 'collab', '0_0_0_0_0_1' =  'empower'))


crcdata_new$IAP2 <- as.factor(crcdata_new$IAP2)


crcdata_new$IAP2 <- factor(crcdata_new$IAP2, levels = c("data gathering", "inform", "consult", "involve", "collab", "empower"))

my_colors <- c("#99CC99", "#669999", "#6666CC", "#9966CC", "#CC9966", "#FFCC66")

ggplot(crcdata_new, aes(fill=IAP2, x=year)) + 
    geom_bar(position="stack", stat="count")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + scale_x_continuous(breaks = crcdata_new$year) + scale_fill_manual(values = my_colors) + scale_fill_discrete(name = "IAP2 Scale")


```

```{r echo=FALSE, warning=FALSE, message=FALSE, results=FALSE}

png("./images/fews_iap2.png", res=300, width = 1500, height = 1500 )
ggplot(crcdata_new, aes(fill=IAP2, x=year)) + 
    geom_bar(position="stack", stat="count")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + scale_x_continuous(breaks = crcdata_new$year) + scale_fill_manual(values = my_colors) + scale_fill_discrete(name = "IAP2 Scale")
dev.off()

```

\newpage

## Level of stakeholder engagement by year - Local scale

```{r echo=FALSE, warning=FALSE, message=FALSE}
#year
library(dplyr)

crcdata_new <- crcdata %>% 
  mutate(local = sprintf("%d_%d_%d_%d_%d_%d_%d_%d", SC_researcher,SC_datagathering,SC_inform,SC_perspectives,SC_plan,SC_identify,SC_envision,SC_implement))


crcdata_new <- crcdata_new[!(crcdata_new$local %in% "0_0_0_0_0_0_0_0"),]

library(dplyr)
crcdata_new <- crcdata_new %>%
    mutate(local = dplyr::recode(local, '1_0_0_0_0_0_0_0' ='researcher', '0_1_0_0_0_0_0_0'='data gathering', '0_0_1_0_0_0_0_0' = 'inform', '0_0_0_1_0_0_0_0' =  'perspectives', '0_0_0_0_1_0_0_0' = 'plan',  '0_0_0_0_0_1_0_0' = 'identify', '0_0_0_0_0_0_1_0' =  'envision',  '0_0_0_0_0_0_0_1' = 'implement',  '1_0_0_0_0_0_0_1' ='researcher', '1_1_0_0_0_0_0_0' ='researcher' ))

crcdata_new$local <- as.factor(crcdata_new$local)


crcdata_new$local <- factor(crcdata_new$local, levels = c("researcher", "data gathering", "inform", "perspectives", "plan", "identify", "envision", "implement"))


my_colors <- c("#66FF33", "#339900", "#99CC99", "#669999", "#6666CC", "#9966CC", "#CC9966", "#FFCC66")

ggplot(crcdata_new, aes(fill=local, x=year)) + 
    geom_bar(position="stack", stat="count")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + scale_x_continuous(breaks = crcdata_new$year) + scale_fill_manual(values = my_colors)+ scale_fill_discrete(name = "Local Scale")




```
\newpage

## Stakeholder engagement by year

```{r echo=FALSE, warning=FALSE, message=FALSE}
#year
library(dplyr)
crcdata_st <- read.csv("./data/crc_data_ES2.csv")



my_colors <- c("#99CC99",  "#CC9966")

crcdata_stremove2011 <- subset(crcdata_st, year != "2011")
colnames(crcdata_stremove2011)[22] <- "Stakeholder Engagement"
ggplot(crcdata_stremove2011, aes(fill=`Stakeholder Engagement`, x=year)) + 
    geom_bar(position="stack", stat="count")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + scale_x_continuous(breaks = crcdata_st$year) + guides(col=guide_legend("Stakeholder Engagement? Yes or No")) + scale_fill_manual(values = my_colors)



```

\newpage

<!-- ## Level of Stakeholder Engagement by Year - Local -->

```{r echo=FALSE, warning=FALSE, message=FALSE}
#year

# library(dplyr)
# 
# crcdata_new <- crcdata %>% 
#   mutate(IAP2 = sprintf("%d_%d_%d_%d_%d_%d", STE_IAP2_data_gathering,STE_IAP2_inform,STE_IAP2_consult,STE_IAP2_involve,STE_IAP2_collab,STE_IAP2_empower))
# 
# crcdata_new <- crcdata_new[!(crcdata_new$IAP2 %in% "0_0_0_0_0_0"),]
# 
# library(dplyr)
# crcdata_new <- crcdata_new %>%
#     mutate(IAP2 = recode(IAP2, '1_0_0_0_0_0' = 'data gathering', '0_1_0_0_0_0' = 'inform', '0_0_1_0_0_0' =  'consult', '0_0_0_1_0_0' = 'involve',  '0_0_0_0_1_0' = 'collab', '0_0_0_0_0_1' =  'empower'))
# 
# ggplot(crcdata_new, aes(fill=IAP2, x=year)) + 
#     geom_bar(position="stack", stat="count")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + scale_x_continuous(breaks = crcdata_new$year)


```
\newpage

## Researcher types

```{r echo=FALSE, warning=FALSE, message=FALSE}
#researcher
df <- t(t(colSums(crcdata[,12:17])))
df <- df[order(df[,1], decreasing = TRUE),]
par(mar = c(6,4,4,4))
#barplot(df, las=2, main = "Researcher Types")

df <- data.frame(df)
rownames(df) <- c("Engineering", "Physics", "Imterdisciplinary", "Computer Science", "Math", "NGO")
df$group <- row.names(df)
df$group <- factor(df$group, levels=c("Engineering", "Physics", "Imterdisciplinary", "Computer Science", "Math", "NGO"))

ggplot(df, aes(x=group, y=df)) +
    geom_bar(stat='identity', fill="#99CC99")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + guides(col=guide_legend("Researcher Types")) +labs(y= "Count", x = "Researcher Type")


```
\newpage

## Stakeholder types

```{r echo=FALSE, warning=FALSE, message=FALSE}

#stakeholder types
df <- t(t(colSums(crcdata[,23:32])))
df <- df[order(df[,1], decreasing = TRUE),]
par(mar = c(11,4,4,4))
#barplot(df, las=2, main = "Stakeholder Types")


df <- data.frame(df)
rownames(df) <- c("Combined Government", "Combined Industry", "Farmers", "University", "Combined Coalition", "Public", "Experts", "Tribal Nations", "Youth", "Migrants")
df$group <- row.names(df)
df$group <- factor(df$group, levels=c("Combined Government", "Combined Industry", "Farmers", "University", "Combined Coalition", "Public", "Experts", "Tribal Nations", "Youth", "Migrants"))

ggplot(df, aes(x=group, y=df)) +
    geom_bar(stat='identity', fill="#99CC99")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + guides(col=guide_legend("Stakeholder Types")) +labs(y= "Count", x = "Stakeholder Type")

```
\newpage

## Ghodsvali scale breakdown

```{r echo=FALSE, warning=FALSE, message=FALSE}

#ghodsvali scale
df <- t(t(colSums(crcdata[,33:36])))
df <- df[order(df[,1]),]
par(mar = c(11,4,4,4))

#barplot(c(df[3], df[4], df[2], df[1]), las=2, main="Ghodsvali Scale Breakdown")

df <- data.frame(df)
rownames(df) <- c("Transformative", "Representation", "Nominal", "Instrumental")
df$group <- row.names(df)
library(dplyr)

# create a vector with letters in the desired order
x <- c("Nominal", "Instrumental",        "Representation",   "Transformative" )

df <- df %>%
  slice(match(x, group))



#rownames(df) <- c("Transformative", "Representation", "Nominal", "Instrumental")
#df$group <- row.names(df)
#df <- c(df[3,], df[4,], df[2,], df[1,])
df$group <- factor(df$group, levels=c("Nominal", "Instrumental", "Representation", "Transformative"))

ggplot(df, aes(x=group, y=df)) +
    geom_bar(stat='identity', fill="#99CC99")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + guides(col=guide_legend("Ghodsvali Scale Breakdown")) +labs(y= "Count", x = "Ghodsvali Scale")

```
\newpage

## IAP2 scale breakdown

```{r echo=FALSE, warning=FALSE, message=FALSE}

#IAP2 scale
df <- t(t(colSums(crcdata[,37:42])))
df <- df[order(df[,1], decreasing = TRUE),]
par(mar = c(12,4,4,4))
#barplot(c(df[2], df[6], df[1], df[3], df[4], df[5]), las=2, main = "IAP2 Scale Breakdown")


df <- data.frame(df)
rownames(df) <- c("Consulting", "Data Gathering", "Involve", "Collaborative", "Empowered", "Informed")
df$group <- row.names(df)
library(dplyr)

# create a vector with letters in the desired order
x <- c("Data Gathering", "Informed",        "Consulting",   "Involve", "Collaborative", "Empowered" )

df <- df %>%
  slice(match(x, group))



#rownames(df) <- c("Transformative", "Representation", "Nominal", "Instrumental")
#df$group <- row.names(df)
#df <- c(df[3,], df[4,], df[2,], df[1,])
df$group <- factor(df$group, levels=c("Data Gathering", "Informed",        "Consulting",   "Involve", "Collaborative", "Empowered"))

ggplot(df, aes(x=group, y=df)) +
    geom_bar(stat='identity', fill="#99CC99")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + guides(col=guide_legend("Ghodsvali Scale Breakdown")) +labs(y= "Count", x = "IAP2 Scale")




```
\newpage

## Local scale breakdown

```{r echo=FALSE, warning=FALSE, message=FALSE}

#local scale
df <- t(t(colSums(crcdata[,43:50])))
#df <- df[order(df[,1], decreasing = TRUE),]
#par(mar = c(9,4,4,4))
#barplot(df, las=2, main="Local Scale Breakdown")


#barplot(c(df[1], df[3], df[6], df[8], df[4], df[2], df[5], df[7]), las=2, main="Local Scale")





df <- data.frame(df)
rownames(df) <- c("Researcher", "Data Gathering", "Inform", "Perspectives", "Plan", "Identify", "Envision", "Implement")
df$group <- row.names(df)
library(dplyr)

# create a vector with letters in the desired order
#x <- c("Researcher", "Data Gathering",        "Inform",   "Envision", "Collaborative", "Empowered" )

#df <- df %>%
#  slice(match(x, group))



#rownames(df) <- c("Transformative", "Representation", "Nominal", "Instrumental")
#df$group <- row.names(df)
#df <- c(df[3,], df[4,], df[2,], df[1,])
df$group <- factor(df$group, levels=c("Researcher", "Data Gathering", "Inform", "Perspectives", "Plan", "Identify", "Envision", "Implement"))

ggplot(df, aes(x=group, y=df)) +
    geom_bar(stat='identity', fill="#99CC99")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + guides(col=guide_legend("Local Scale Breakdown")) +labs(y= "Count", x = "Local Scale")



```
\newpage

## Geographic location breakdown

```{r echo=FALSE, warning=FALSE, message=FALSE}

#geography
df <- t(t(colSums(crcdata[,52:57])))
df <- df[order(df[,1], decreasing = TRUE),]
par(mar = c(8,4,4,4))
#barplot(df, las=2, main="Geographic Location Breakdown")




df <- data.frame(df)
rownames(df) <- c("Regional", "Local", "National", "No Geography", "Multinational", "Global")
df$group <- row.names(df)
library(dplyr)

# create a vector with letters in the desired order
#x <- c("Researcher", "Data Gathering",        "Inform",   "Envision", "Collaborative", "Empowered" )

#df <- df %>%
#  slice(match(x, group))



#rownames(df) <- c("Transformative", "Representation", "Nominal", "Instrumental")
#df$group <- row.names(df)
#df <- c(df[3,], df[4,], df[2,], df[1,])
df$group <- factor(df$group, levels=c("Regional", "Local", "National", "No Geography", "Multinational", "Global"))

ggplot(df, aes(x=group, y=df)) +
    geom_bar(stat='identity', fill="#99CC99")+theme_bw()+theme(axis.text.x = element_text(angle = 90, vjust = 0.5, hjust=1)) + guides(col=guide_legend("Geography")) +labs(y= "Count", x = "Geography")

```


```{r echo=FALSE, warning=FALSE, message=FALSE}

#Odds of stakeholder scale predicting whether a solution was proposed or not

ghodsvali.fit<-glm(solution_proposed_YN ~STE_G_nominal+STE_G_instrumental+STE_G_representation+STE_G_transformative,family=binomial,data=crcdata)
IAP2.fit<-glm(solution_proposed_YN ~STE_IAP2_data_gathering+STE_IAP2_inform+STE_IAP2_consult+STE_IAP2_involve+STE_IAP2_collab+STE_IAP2_empower,family=binomial,data=crcdata)

local.fit<-glm(solution_proposed_YN ~SC_researcher+SC_datagathering+SC_inform+SC_perspectives+SC_plan+SC_identify+SC_envision+SC_implement,family=binomial,data=crcdata)

```

\newpage
# Regression Testing

## Ghodsvali scale regression

Ghodsvali scale regression testing on whether a solution was proposed or not

```{r echo=FALSE, warning=FALSE, message=FALSE}

summary(ghodsvali.fit)

```
\newpage

## Ghodsvali scale odds

Odds of Ghodsvali scale predicting whether a solution was proposed or not

```{r echo=FALSE, warning=FALSE, message=FALSE}

logistic.display(ghodsvali.fit)


```

\newpage

## IAP2 scale regression 

IAP2 scale regression testing on whether a solution was proposed or not

```{r echo=FALSE, warning=FALSE, message=FALSE}

summary(IAP2.fit)

```

\newpage

## IAP2 scale odds

Odds of IAP2 scale predicting whether a solution was proposed or not

```{r echo=FALSE, warning=FALSE, message=FALSE}

logistic.display(IAP2.fit)

```

\newpage

## Local scale regression

Local scale regression predicting whether a solution was proposed or not

```{r echo=FALSE, warning=FALSE, message=FALSE}

summary(local.fit)

```

\newpage

## Local scale odds

Odds of Local scale predicting whether a solution was proposed or not

```{r echo=FALSE, warning=FALSE, message=FALSE}

logistic.display(local.fit)

```


\newpage

## Engagement vs. solution - regression testing

Regression of whether stakeholder engagement (Y/N) predicts if a solution was proposed (Y/N)

```{r echo=FALSE, warning=FALSE, message=FALSE}

#engagement of stakeholders vs. solutions

stakeholder_vs_solutions <-glm(solution_proposed_YN ~S_stakeholder_engagment_YN,family=binomial,data=crcdata)

summary(stakeholder_vs_solutions)
logistic.display(stakeholder_vs_solutions)

```

\newpage

## Engagement vs. solution - odds

Odds of whether stakeholder engagement (Y/N) predicts if a solution was proposed (Y/N)

```{r echo=FALSE, warning=FALSE, message=FALSE}

#engagement of stakeholders vs. solutions

stakeholder_vs_solutions <-glm(solution_proposed_YN ~S_stakeholder_engagment_YN,family=binomial,data=crcdata)

#summary(stakeholder_vs_solutions)
logistic.display(stakeholder_vs_solutions)

```


\newpage

## Use of computational model vs. solution - regression testing

Regression of whether a computational model was used predicts if a solution was proposed (Y/N)

```{r echo=FALSE, warning=FALSE, message=FALSE}

#Use of computational model vs. solutions

model_vs_solutions <-glm(solution_proposed_YN ~S_model_YN,family=binomial,data=crcdata)

summary(model_vs_solutions)
logistic.display(model_vs_solutions)

```


\newpage

## Use of computational model vs. solution - odds

Odds of whether a computational model was used predicts if a solution was proposed (Y/N)

```{r echo=FALSE, warning=FALSE, message=FALSE}

#Use of computational model vs. solutions

model_vs_solutions <-glm(solution_proposed_YN ~S_model_YN,family=binomial,data=crcdata)

#summary(model_vs_solutions)
logistic.display(model_vs_solutions)

```
\newpage

## Diversity of stakeholders vs solution

Regression of whether Diversity of stakeholders predicts if a solution was proposed (Y/N).  In order to represent diversity, we have used a simple ratio calcuation which sums the number of stakeholders involved divided by the total number of possible stakeholder options.  A ratio which is closer to 1 has a greater level of stakeholder diversity.

```{r echo=FALSE, warning=FALSE, message=FALSE}

#Diversity of stakeholders vs solutions (using ratio)

crcdata$ST_ratio <- rowSums(crcdata[,23:31])/9

stratio_vs_solutions <-glm(solution_proposed_YN ~ ST_ratio,family=binomial,data=crcdata)

summary(stratio_vs_solutions)
#logistic.display(stratio_vs_solutions)

```

\newpage

## Diversity of stakeholders vs solution - odds

Odds whether Diversity of stakeholders predicts if a solution was proposed (Y/N). In order to represent diversity, we have used a simple ratio calcuation which sums the number of stakeholders involved divided by the total number of possible stakeholder options.  A ratio which is closer to 1 has a greater level of stakeholder diversity.

```{r echo=FALSE, warning=FALSE, message=FALSE}

#Diversity of stakeholders vs solutions (using ratio)

crcdata$ST_ratio <- rowSums(crcdata[,23:31])/9

stratio_vs_solutions <-glm(solution_proposed_YN ~ ST_ratio,family=binomial,data=crcdata)

#summary(stratio_vs_solutions)
logistic.display(stratio_vs_solutions)

```

\newpage

## Diversity of stakeholders vs solution - full predictors


```{r echo=FALSE, warning=FALSE, message=FALSE}

#Diversity of stakeholders vs solutions (using ratio)

crcdata$ST_ratio <- rowSums(crcdata[,23:31])/10

stratio_vs_solutions <-glm(solution_proposed_YN ~ ST_farmers + ST_combined_gov  + ST_combined_coalition + ST_combined_industry + ST_migrants + ST_public + ST_university + ST_experts,family=binomial,data=crcdata)

summary(stratio_vs_solutions)
#logistic.display(stratio_vs_solutions)

```

\newpage
## Diversity of researchers vs solution - regression testing

Regression of whether Diversity of researchers predicts if a solution was proposed (Y/N). In order to represent diversity, we have used a simple ratio calculation which sums the number of researcher types involved, divided by the total number of possible researcher options.  A ratio which is closer to 1 has a greater level of researcher diversity.

```{r echo=FALSE, warning=FALSE, message=FALSE}

#Diversity of researchers vs solutions (using ratio)


crcdata$R_ratio <- rowSums(crcdata[,12:21])/10

rratio_vs_solutions <-glm(solution_proposed_YN ~ R_ratio,family=binomial,data=crcdata)

summary(rratio_vs_solutions)
#logistic.display(rratio_vs_solutions)

```


\newpage
## Diversity of researchers vs solution - odds

Odds of whether Diversity of researchers predicts if a solution was proposed (Y/N). In order to represent diversity, we have used a simple ratio calculation which sums the number of researcher types involved, divided by the total number of possible researcher options.  A ratio which is closer to 1 has a greater level of researcher diversity.

```{r echo=FALSE, warning=FALSE, message=FALSE}

#Diversity of researchers vs solutions (using ratio)


crcdata$R_ratio <- rowSums(crcdata[,12:21])/10

rratio_vs_solutions <-glm(solution_proposed_YN ~ R_ratio,family=binomial,data=crcdata)

#summary(rratio_vs_solutions)
logistic.display(rratio_vs_solutions)

```


<!-- ## Interdisciplinary Researchers vs Solution - Regression Testing -->

<!-- In order to accurately reflect diversity of researchers, we run a binomial logistic regression of the interdisciplinary researchers category to assess predictability on stakeholder solution (Y/N). -->

<!-- ```{r echo=FALSE, warning=FALSE, message=FALSE} -->

<!-- #Diversity of researchers vs solutions (using ratio) -->


<!-- crcdata$R_ratio <- rowSums(crcdata[,12:21])/10 -->

<!-- rratio_vs_solutions <-glm(solution_proposed_YN ~ R_interdis,family=binomial,data=crcdata) -->

<!-- summary(rratio_vs_solutions) -->
<!-- #logistic.display(rratio_vs_solutions) -->

<!-- ``` -->

\newpage

## Stakeholder type vs level of engagement (Ghodsvali) - regression testing

```{r echo=FALSE, warning=FALSE, message=FALSE}

#stakeholder type vs level of engagement ghodsvali

ST_engagement<-lm(cbind(ST_farmers,ST_combined_gov,ST_tribal,ST_combined_coalition,ST_combined_industry,ST_migrants,ST_youth,ST_public,ST_university,ST_experts)~STE_G_nominal+STE_G_instrumental+STE_G_representation+STE_G_transformative,data=crcdata)

summary(ST_engagement)


```
<!-- \newpage -->
<!-- ## Stakeholder type vs level of engagement (Ghodsvali) -->
<!-- ) -->
<!-- ```{r echo=FALSE, warning=FALSE, message=FALSE} -->

<!-- #stakeholder type vs level of engagement ghodsvali -->

<!-- ST_engagement_ghodsvali<-lm(cbind(ST_farmers,ST_combined_gov,ST_tribal,ST_combined_coalition,ST_combined_industry,ST_migrants,ST_youth,ST_public,ST_university,ST_experts)~STE_G_nominal+STE_G_instrumental+STE_G_representation+STE_G_transformative,data=crcdata) -->


<!-- ``` -->

\newpage

## Stakeholder type vs level of engagement (IAP2) - regression testing

```{r echo=FALSE, warning=FALSE, message=FALSE}

#stakeholder type vs level of engagement IAP2

ST_engagement_IAP2<-lm(cbind(ST_farmers,ST_combined_gov,ST_tribal,ST_combined_coalition,ST_combined_industry,ST_migrants,ST_youth,ST_public,ST_university,ST_experts)~STE_IAP2_data_gathering+STE_IAP2_inform+STE_IAP2_consult+STE_IAP2_involve+STE_IAP2_collab+STE_IAP2_empower,data=crcdata)

summary(ST_engagement_IAP2)

```

\newpage

## Stakeholder type vs level of engagement (local) - regression testing

```{r echo=FALSE, warning=FALSE, message=FALSE}

#stakeholder type vs level of engagement local

ST_engagement_IAP2<-lm(cbind(ST_farmers,ST_combined_gov,ST_tribal,ST_combined_coalition,ST_combined_industry,ST_migrants,ST_youth,ST_public,ST_university,ST_experts)~SC_researcher+SC_datagathering+SC_inform+SC_perspectives+SC_plan+SC_identify+SC_envision+SC_implement,data=crcdata)

summary(ST_engagement_IAP2)


```

\newpage

## Stakeholder type vs solution - regression testing

```{r echo=FALSE, warning=FALSE, message=FALSE}

#stakeholder type vs solution

ST_solution <- glm(solution_proposed_YN~ST_farmers+ST_combined_gov+ST_tribal+ST_combined_coalition+ST_combined_industry+ST_migrants+ST_youth+ST_public+ST_university+ST_experts, family=binomial, data=crcdata)

summary(ST_solution)

```

\newpage

## Geographic area vs solution - regression testing

```{r echo=FALSE, warning=FALSE, message=FALSE}

#geographic area vs solution

geog_solution <- glm(solution_proposed_YN~G_notdescribed+G_local+G_regional+G_national+G_multination+G_global, family=binomial, data=crcdata)
summary(geog_solution)
logistic.display(geog_solution)


```

\newpage

## Geographic area vs solution - odds

```{r echo=FALSE, warning=FALSE, message=FALSE}

#geographic area vs solution

logistic.display(geog_solution)


```

\newpage

## stakeholder type vs geographic area

```{r echo=FALSE, warning=FALSE, message=FALSE}

#geographic area vs solution

geog_stakeholdertype <- lm(cbind(ST_farmers,ST_combined_gov,ST_tribal,ST_combined_coalition,ST_combined_industry,ST_migrants,ST_youth,ST_public,ST_university,ST_experts)~G_notdescribed+G_local+G_regional+G_national+G_multination+G_global, data=crcdata)

summary(geog_stakeholdertype)

```

\newpage

## Geographic area vs engagment (Ghodsvali)

```{r echo=FALSE, warning=FALSE, message=FALSE}

#geographic area vs solution

geog_engagment_ghodsvali <- lm(cbind(G_notdescribed,G_local,G_regional,G_national,G_multination,G_global)~STE_G_nominal+STE_G_instrumental+STE_G_representation+STE_G_transformative, data=crcdata)

summary(geog_engagment_ghodsvali)

```
\newpage

## Geographic area vs engagment (IAP2)
```{r echo=FALSE, warning=FALSE, message=FALSE}


geog_engagment_IAP2 <- lm(cbind(G_notdescribed,G_local,G_regional,G_national,G_multination,G_global)~STE_IAP2_data_gathering+STE_IAP2_inform+STE_IAP2_consult+STE_IAP2_involve+STE_IAP2_collab+STE_IAP2_empower, data=crcdata)

summary(geog_engagment_IAP2)

```
\newpage

## Geographic area vs engagment (local)
```{r echo=FALSE, warning=FALSE, message=FALSE}


geog_engagment_local <- lm(cbind(G_notdescribed,G_local,G_regional,G_national,G_multination,G_global)~SC_researcher+SC_datagathering+SC_inform+SC_perspectives+SC_plan+SC_identify+SC_envision+SC_implement, data=crcdata)

summary(geog_engagment_local)

```


\newpage

## Geographic area vs engagment (local)

```{r echo=FALSE, warning=FALSE, message=FALSE}


geog_engagment_local <- lm(cbind(G_notdescribed,G_local,G_regional,G_national,G_multination,G_global)~SC_researcher+SC_datagathering+SC_inform+SC_perspectives+SC_plan+SC_identify+SC_envision+SC_implement, data=crcdata)

summary(geog_engagment_local)

```