Preparations

if(!require(lme4)) install.packages("lme4")       #for GLMMs
if(!require(glmmTMB)) install.packages("glmmTMB")    #for glmmTMB
if(!require(emmeans)) install.packages("emmeans")   #for marginal means etc
if(!require(MuMIn)) install.packages("MuMIn")      #for model comparison
if(!require(broom)) install.packages("broom")     #for tidying outputs
if(!require(DHARMa))  install.packages("DHARMa")   #for residual diagnostics
if(!require(ggeffects)) install.packages("ggeffects") #for partial plots
if(!require(mgcv))  install.packages("mgcv")     #for GAMs
if(!require(performance)) install.packages("performance") #for model diagnostics
if(!require(see)) install.packages("see")       #for model diagnostics
if(!require(GGally))  install.packages("GGally")   #for correlation matrix
if(!require(tidyverse)) install.packages("tidyverse") #for data wrangling etc
if(!require(cowplot)) install.packages("cowplot")   #for adding theme
if(!require(rstanarm))  install.packages("rstanarm") #for fitting models in STAN
if(!require(brms))  install.packages("brms")     #for fitting models in STAN
if(!require(coda))  install.packages("coda")     #for diagnostics
if(!require(bayesplot)) install.packages("bayesplot") #for diagnostics

Scenario

The persistence of Australian terrestrial and aquatic functional groups was evaluated by experts against three levels of threat intensity. This analysis investigates the overall confidence the expert’s had in their assessments.

Terrestrial

Read in the data

terrestrial = read_csv('../data/Data_Table_1_Persistence_Data_Terrestrial.csv', trim_ws=TRUE)
## 
## -- Column specification --------------------------------------------------------
## cols(
##   Expert_ID = col_double(),
##   Category = col_character(),
##   Group = col_character(),
##   Threat = col_character(),
##   Level = col_double(),
##   Estimate = col_character(),
##   Persistence = col_double(),
##   Confidence = col_double(),
##   Persistence_fit = col_double(),
##   Confidence_req = col_double()
## )
glimpse(terrestrial)
## Rows: 66,033
## Columns: 10
## $ Expert_ID       <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ Category        <chr> "Mammals", "Mammals", "Mammals", "Mammals", "Mammal...
## $ Group           <chr> "M01", "M01", "M01", "M01", "M01", "M01", "M01", "M...
## $ Threat          <chr> "Fire", "Fire", "Fire", "Fire", "Fire", "Fire", "Fi...
## $ Level           <dbl> 1, 1, 1, 2, 2, 2, 3, 3, 3, 1, 1, 1, 2, 2, 2, 3, 3, ...
## $ Estimate        <chr> "Upper", "Best", "Lower", "Upper", "Best", "Lower",...
## $ Persistence     <dbl> 1.0, 1.0, 1.0, 1.0, 1.0, 0.5, 1.0, 0.5, 0.0, 1.0, 1...
## $ Confidence      <dbl> 90, 90, 90, 60, 60, 60, 60, 60, 60, 90, 90, 90, 60,...
## $ Persistence_fit <dbl> 1.000, 1.000, 1.000, 1.000, 1.000, 0.333, 1.000, 0....
## $ Confidence_req  <dbl> 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80,...

Reformat variables

terrestrial = terrestrial %>% mutate(
                       Category=factor(Category),
                       Func_Group=factor(Group),
                       Threat=factor(Threat),
                       Estimate=factor(Estimate),
                       Level=factor(Level),
                       Pers_fit_trans=sqrt(1-Persistence_fit))
glimpse(terrestrial)
## Rows: 66,033
## Columns: 12
## $ Expert_ID       <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ Category        <fct> Mammals, Mammals, Mammals, Mammals, Mammals, Mammal...
## $ Group           <chr> "M01", "M01", "M01", "M01", "M01", "M01", "M01", "M...
## $ Threat          <fct> Fire, Fire, Fire, Fire, Fire, Fire, Fire, Fire, Fir...
## $ Level           <fct> 1, 1, 1, 2, 2, 2, 3, 3, 3, 1, 1, 1, 2, 2, 2, 3, 3, ...
## $ Estimate        <fct> Upper, Best, Lower, Upper, Best, Lower, Upper, Best...
## $ Persistence     <dbl> 1.0, 1.0, 1.0, 1.0, 1.0, 0.5, 1.0, 0.5, 0.0, 1.0, 1...
## $ Confidence      <dbl> 90, 90, 90, 60, 60, 60, 60, 60, 60, 90, 90, 90, 60,...
## $ Persistence_fit <dbl> 1.000, 1.000, 1.000, 1.000, 1.000, 0.333, 1.000, 0....
## $ Confidence_req  <dbl> 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80,...
## $ Func_Group      <fct> M01, M01, M01, M01, M01, M01, M01, M01, M01, M02, M...
## $ Pers_fit_trans  <dbl> 0.0000000, 0.0000000, 0.0000000, 0.0000000, 0.00000...

Subset threats by Fire, Grazing, Buffalo,Cane toad, Cat,Pig, Perennial grasses (e.g. gamba grass).

terrestrial_sub = terrestrial %>% filter(Threat == "Fire" | Threat == "Grazing" | Threat == "Buffalo"| Threat == "Cane toad"| Threat == "Cat"| Threat == "Pig"| Threat == "Perennial grasses") %>% droplevels()
unique(terrestrial_sub$Threat)
## [1] Fire              Grazing           Buffalo           Cane toad        
## [5] Cat               Perennial grasses Pig              
## Levels: Buffalo Cane toad Cat Fire Grazing Perennial grasses Pig

Exploratory Data Analysis

#ggpairs(terrestrial_sub,title="Correlogram with ggpairs", #columns=c("Category","Level","Persistence_fit","Threat","Expert_ID"),
#        ggplot2::aes(colour="Category"))+
#  scale_color_brewer(palette = "Set1")+
#  scale_fill_brewer(palette = "Set1")
ggplot(terrestrial_sub,aes(y=Persistence_fit,x=Level,fill=Threat))+
  geom_boxplot(outlier.shape = NA)+
  facet_grid(.~Category)+
  scale_fill_viridis_d()+
  theme_bw()

ggplot(data=terrestrial_sub)+
  geom_histogram(aes(x=Persistence_fit,fill=Threat))+
  facet_grid(.~Category)+
  scale_fill_viridis_d()+
  theme_bw()
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Persistence summary figure

terr_sum_fig_prelim<-ggplot(data=subset(terrestrial_sub,Estimate=="Best"&Category!="Vegetation"), aes(y=Persistence_fit,x=Level,fill=Category,colour=Category))+
  facet_wrap(.~Threat,ncol=2, strip.position="top")+
  geom_boxplot(outlier.shape = NA)+
  scale_x_discrete("Threat level")+
  scale_y_continuous("Probability of persistence", limits=c(0.0,1.05),expand=c(0,0))+
  theme_bw()+
  theme(strip.background = element_blank(),strip.placement = "outside")+
  scale_colour_viridis_d()+
  scale_fill_viridis_d()
 
terr_sum_fig_prelim

#save_plot("../docs/terr_sum_fig",terr_sum_fig,nrow = 7, ncol=1, base_asp = 1.3)

Create individual plots for terrestrial categories for the threats: Fire, Grazing, Buffalo, Cane toad, Cat, Perennial grasses, Pig

Function to increase legend size:

addBigLegend <- function(myPlot, pointSize = 16, textSize = 16, spaceLegend = 2) {
    myPlot +
        guides(shape = guide_legend(override.aes = list(size = pointSize)),
               color = guide_legend(override.aes = list(size = pointSize))) +
        theme(legend.title = element_text(size = textSize), 
              legend.text  = element_text(size = textSize),
              legend.key.size = unit(spaceLegend, "lines"))
}

Suitable colour tones can be found here.

terr_test<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Fire"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  #geom_point(aes(colour=Category),alpha=0.3,position=position_jitterdodge(0.2))+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  #scale_colour_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Test/Fire")+
  theme(plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_test

terr_test_large<-addBigLegend(terr_test)
terr_test_large

terr_fire<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Fire"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Fire")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_fire

terr_grazing<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Grazing"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Grazing")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_grazing

terr_buffalo<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Buffalo"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Buffalo")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_buffalo

terr_toad<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Cane toad"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Cane toad")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_toad

terr_cat<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Cat"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Cat")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_cat

terr_grass<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Perennial grasses"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Perennial grasses")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_grass

terr_pig<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Pig"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Pig")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_pig

Extract legend

legend_terr<-get_legend(terr_test_large+
                          theme(legend.direction = "vertical", legend.justification = "center"))

Summary figure terrestrial categories

terr_sum<-plot_grid(terr_grazing,legend_terr,terr_toad,terr_cat,terr_pig,terr_grass,terr_buffalo,terr_fire, labels = c("a","","b","c","d","e","f","g"), ncol = 2, align="hv",axis = "b",
                    rel_widths = c(1,1,1,1,1,1,1,1)) #
terr_sum

save_plot("../docs/terr_sum_earth.pdf",terr_sum,base_asp=1.3,ncol=2,nrow=4)

Figure including small mammals

Subset terrestrial data set, where small mammals func group M08, is included into category.

terrestrial_sub2<-terrestrial_sub %>% filter(Func_Group=="M08") %>% droplevels()
terrestrial_sub3 = terrestrial_sub2 %>% mutate(
                       Category=factor(Func_Group))
terrestrial_sub_final = bind_rows(terrestrial_sub,terrestrial_sub3)
terr_test<-ggplot(data=subset(terrestrial_sub_final,Estimate=="Best" & Category!="Vegetation" & Threat=="Fire"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category, levels = c("Amphibians","Reptiles","Birds","Mammals","M08"), labels = c("Amphibians","Reptiles","Birds","Mammals","Small ground mammals"))),alpha=0.9,outlier.shape = NA)+
  #geom_point(aes(colour=Category),alpha=0.3,position=position_jitterdodge(0.2))+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown","grey85"))+
  #scale_colour_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Test/Fire")+
  theme(plot.title = element_text(hjust = 0.5))+
  guides(fill=guide_legend(title="Category"))
  #scale_fill_viridis_d()
terr_test

terr_test_large<-addBigLegend(terr_test)
terr_test_large

terr_fire<-ggplot(data=subset(terrestrial_sub_final,Estimate=="Best" & Category!="Vegetation" & Threat=="Fire"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category, levels = c("Amphibians","Reptiles","Birds","Mammals","M08"), labels = c("Amphibians","Reptiles","Birds","Mammals","Small ground mammals"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown","grey85"))+
  ggtitle("Fire")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_fire

terr_grazing<-ggplot(data=subset(terrestrial_sub_final,Estimate=="Best" & Category!="Vegetation" & Threat=="Grazing"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category, levels = c("Amphibians","Reptiles","Birds","Mammals","M08"), labels = c("Amphibians","Reptiles","Birds","Mammals","Small ground mammals"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown","grey85"))+
  ggtitle("Grazing")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_grazing

terr_buffalo<-ggplot(data=subset(terrestrial_sub_final,Estimate=="Best" & Category!="Vegetation" & Threat=="Buffalo"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category, levels = c("Amphibians","Reptiles","Birds","Mammals","M08"), labels = c("Amphibians","Reptiles","Birds","Mammals","Small ground mammals"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown","grey85"))+
  ggtitle("Buffalo")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_buffalo

terr_toad<-ggplot(data=subset(terrestrial_sub_final,Estimate=="Best" & Category!="Vegetation" & Threat=="Cane toad"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category, levels = c("Amphibians","Reptiles","Birds","Mammals","M08"), labels = c("Amphibians","Reptiles","Birds","Mammals","Small ground mammals"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown","grey85"))+
  ggtitle("Cane toad")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_toad

terr_cat<-ggplot(data=subset(terrestrial_sub_final,Estimate=="Best" & Category!="Vegetation" & Threat=="Cat"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category, levels = c("Amphibians","Reptiles","Birds","Mammals","M08"), labels = c("Amphibians","Reptiles","Birds","Mammals","Small ground mammals"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown","grey85"))+
  ggtitle("Cat")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_cat

terr_grass<-ggplot(data=subset(terrestrial_sub_final,Estimate=="Best" & Category!="Vegetation" & Threat=="Perennial grasses"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category, levels = c("Amphibians","Reptiles","Birds","Mammals","M08"), labels = c("Amphibians","Reptiles","Birds","Mammals","Small ground mammals"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown","grey85"))+
  ggtitle("Perennial grasses")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_grass

terr_pig<-ggplot(data=subset(terrestrial_sub_final,Estimate=="Best" & Category!="Vegetation" & Threat=="Pig"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category, levels = c("Amphibians","Reptiles","Birds","Mammals","M08"), labels = c("Amphibians","Reptiles","Birds","Mammals","Small ground mammals"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown","grey85"))+
  ggtitle("Pig")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_pig

Extract legend

legend_terr<-get_legend(terr_test_large+
                          theme(legend.direction = "vertical", legend.justification = "center"))

Summary figure terrestrial categories

terr_sum_M08<-plot_grid(terr_grazing,legend_terr,terr_toad,terr_cat,terr_pig,terr_grass,terr_buffalo,terr_fire, labels = c("a","","b","c","d","e","f","g"), ncol = 2, align="hv",axis = "b",
                    rel_widths = c(1,1,1,1,1,1,1,1)) #
terr_sum_M08

save_plot("../docs/terr_sum_earth_small_mammals.pdf",terr_sum_M08,base_asp=1.3,ncol=2,nrow=4)

Confidence summary figure

terr_sum_fig2_prelim<-ggplot(data=subset(terrestrial_sub,Estimate=="Best"&Category!="Vegetation"), aes(y=Confidence,x=Level,fill=Category))+
  facet_wrap(.~Threat,ncol=2, strip.position="top")+
  geom_boxplot(outlier.shape = NA)+
  scale_x_discrete("Threat level")+
  scale_y_continuous("Confidence")+
  theme_bw()+
  theme(strip.background = element_blank(),strip.placement = "outside")+
  #scale_colour_viridis_d()+
  scale_fill_viridis_d()
 
terr_sum_fig2_prelim

#save_plot("../docs/terr_sum_fig",terr_sum_fig,nrow = 7, ncol=1, base_asp = 1.3)
terr_grazing2<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Grazing"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Confidence (%)", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Grazing")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
terr_grazing2

terr_toad2<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Cane toad"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Confidence (%)", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Cane toad")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_toad2

terr_cat2<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Cat"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Cat")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_cat2

terr_pig2<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Pig"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Confidence (%)", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Pig")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_pig2

terr_grass2<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Perennial grasses"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Perennial grasses")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_grass2

terr_buffalo2<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Buffalo"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="Confidence (%)", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Buffalo")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_buffalo2

terr_fire2<-ggplot(data=subset(terrestrial_sub,Estimate=="Best" & Category!="Vegetation" & Threat=="Fire"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("goldenrod", "burlywood","sienna","brown"))+
  ggtitle("Fire")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
terr_fire2

terr_sum_conf<-plot_grid(terr_grazing2,legend_terr,terr_toad2,terr_cat2,terr_pig2,terr_grass2,terr_buffalo2,terr_fire2, labels = c("a","","b","c","d","e","f","g"), ncol = 2, align="hv",axis = "b",
                    rel_widths = c(1,1,1,1,1,1,1,1)) #
terr_sum_conf

save_plot("../docs/terr_sum_conf_earth.pdf",terr_sum_conf,base_asp=1.3,ncol=2,nrow=4)

Statistical analysis (Bayesian)

Fit the model

Advice on the implementation of Zero-One-Inflated-Beta models can be found here.

Model validation

load(file="../data/Testfit.RData")
preds <-posterior_predict(fit,nsamples=250,summary=FALSE) #was 250 instead of 30

ZOIB.resids<-createDHARMa(simulatedResponse = t(preds),
                          observedResponse = terrestrial_sub$Persistence_fit,
                          fittedPredictedResponse = apply(preds,2,median))
plot(ZOIB.resids)

Predictions

Summary figure

7 Threats as panels, in each predict persistence_fit for each category’s levels (1,2,3)

Aquatic

Read in data

aquatic = read_csv('../data/Data_Table_2_Persistence_Data_Aquatic.csv', trim_ws=TRUE)
## 
## -- Column specification --------------------------------------------------------
## cols(
##   Expert_ID = col_double(),
##   Category = col_character(),
##   Group = col_character(),
##   Threat = col_character(),
##   Level = col_double(),
##   Estimate = col_character(),
##   Persistence = col_double(),
##   Confidence = col_double(),
##   Persistence_fit = col_double(),
##   Confidence_req = col_double()
## )
glimpse(aquatic)
## Rows: 4,914
## Columns: 10
## $ Expert_ID       <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ Category        <chr> "Fish", "Fish", "Fish", "Fish", "Fish", "Fish", "Fi...
## $ Group           <chr> "F01", "F01", "F01", "F02", "F02", "F02", "F03", "F...
## $ Threat          <chr> "Altered flow regime", "Altered flow regime", "Alte...
## $ Level           <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ Estimate        <chr> "Upper", "Best", "Lower", "Upper", "Best", "Lower",...
## $ Persistence     <dbl> 1.0, 0.8, 0.5, 1.0, 0.8, 0.5, 1.0, 0.8, 0.5, 1.0, 0...
## $ Confidence      <dbl> 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70,...
## $ Persistence_fit <dbl> 1.000, 0.800, 0.457, 1.000, 0.800, 0.457, 1.000, 0....
## $ Confidence_req  <dbl> 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80,...

Exploratory Data Analysis

aquatic = aquatic %>% mutate(
                       Category=factor(Category),
                       Func_Group=factor(Group),
                       Threat=factor(Threat),
                       Estimate=factor(Estimate),
                       Level=factor(Level),
                       Pers_fit_trans=sqrt(1-Persistence_fit))

glimpse(aquatic)
## Rows: 4,914
## Columns: 12
## $ Expert_ID       <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ Category        <fct> Fish, Fish, Fish, Fish, Fish, Fish, Fish, Fish, Fis...
## $ Group           <chr> "F01", "F01", "F01", "F02", "F02", "F02", "F03", "F...
## $ Threat          <fct> Altered flow regime, Altered flow regime, Altered f...
## $ Level           <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ Estimate        <fct> Upper, Best, Lower, Upper, Best, Lower, Upper, Best...
## $ Persistence     <dbl> 1.0, 0.8, 0.5, 1.0, 0.8, 0.5, 1.0, 0.8, 0.5, 1.0, 0...
## $ Confidence      <dbl> 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70,...
## $ Persistence_fit <dbl> 1.000, 0.800, 0.457, 1.000, 0.800, 0.457, 1.000, 0....
## $ Confidence_req  <dbl> 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80, 80,...
## $ Func_Group      <fct> F01, F01, F01, F02, F02, F02, F03, F03, F03, F04, F...
## $ Pers_fit_trans  <dbl> 0.0000000, 0.4472136, 0.7368853, 0.0000000, 0.44721...

Check aquatic threats:

unique(aquatic$Threat)
## [1] Altered flow regime   Aquatic weeds         Buffalos             
## [4] Cane toads            Grazing               Longitudinal barriers
## [7] Pigs                 
## 7 Levels: Altered flow regime Aquatic weeds Buffalos Cane toads ... Pigs

Persistence summary figure

aqua_test<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Cane toads"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Test/Cane toad")+
  theme(plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_test

aqua_test_large<-addBigLegend(aqua_test)
aqua_test_large

aqua_grazing<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Grazing"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Grazing")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_grazing

aqua_toads<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Cane toads"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Cane toad")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_toads

aqua_weeds<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Aquatic weeds"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Aquatic weeds")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_weeds

aqua_pigs<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Pigs"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Pig")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_pigs

aqua_flow<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Altered flow regime"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Altered flow regime")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_flow

aqua_buffalo<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Buffalos"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Buffalo")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_buffalo

aqua_barrier<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Longitudinal barriers"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Longitudinal barriers")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_barrier

Extract legend

legend_aqua<-get_legend(aqua_test_large+
                          theme(legend.direction = "vertical", legend.justification = "center"))

Summary figure terrestrial categories

aqua_sum<-plot_grid(aqua_grazing,legend_aqua,aqua_toads,aqua_weeds,aqua_pigs,aqua_flow,aqua_buffalo,aqua_barrier, labels = c("a","","b","c","d","e","f","g"), ncol = 2, align="hv",axis = "b",
                    rel_widths = c(1,1,1,1,1,1,1,1)) #
aqua_sum

save_plot("../docs/aqua_sum_water.pdf",aqua_sum,base_asp=1.3,ncol=2,nrow=4)

Figure including grunters

Subset aquatic data set, where grunter func group F05, is included into category.

aquatic2<-aquatic %>% filter(Func_Group=="F05") %>% droplevels()
aquatic3 = aquatic2 %>% mutate(
                       Category=factor(Func_Group))
aquatic_final = bind_rows(aquatic,aquatic3)
aqua_test<-ggplot(data=subset(aquatic_final,Estimate=="Best" & Threat=="Cane toads"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category,levels = c("Waterbirds","Turtles","Fish","F05"),labels = c("Waterbirds","Turtles","Fish","Small-bodied\nmigratory invertivores"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E","grey85"))+
  ggtitle("Test/Cane toad")+
  theme(plot.title = element_text(hjust = 0.5))+
  guides(fill=guide_legend(title="Category"))
  #scale_fill_viridis_d()
aqua_test

aqua_test_large<-addBigLegend(aqua_test)
aqua_test_large

aqua_grazing<-ggplot(data=subset(aquatic_final,Estimate=="Best" & Threat=="Grazing"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category,levels = c("Waterbirds","Turtles","Fish","F05"),labels = c("Waterbirds","Turtles","Fish","Grunters"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E","grey85"))+
  ggtitle("Grazing")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_grazing

aqua_toads<-ggplot(data=subset(aquatic_final,Estimate=="Best" & Threat=="Cane toads"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category,levels = c("Waterbirds","Turtles","Fish","F05"),labels = c("Waterbirds","Turtles","Fish","Grunters"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E","grey85"))+
  ggtitle("Cane toad")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_toads

aqua_weeds<-ggplot(data=subset(aquatic_final,Estimate=="Best" & Threat=="Aquatic weeds"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category,levels = c("Waterbirds","Turtles","Fish","F05"),labels = c("Waterbirds","Turtles","Fish","Grunters"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E","grey85"))+
  ggtitle("Aquatic weeds")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_weeds

aqua_pigs<-ggplot(data=subset(aquatic_final,Estimate=="Best" & Threat=="Pigs"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category,levels = c("Waterbirds","Turtles","Fish","F05"),labels = c("Waterbirds","Turtles","Fish","Grunters"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E","grey85"))+
  ggtitle("Pig")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_pigs

aqua_flow<-ggplot(data=subset(aquatic_final,Estimate=="Best" & Threat=="Altered flow regime"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category,levels = c("Waterbirds","Turtles","Fish","F05"),labels = c("Waterbirds","Turtles","Fish","Grunters"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E","grey85"))+
  ggtitle("Altered flow regime")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_flow

aqua_buffalo<-ggplot(data=subset(aquatic_final,Estimate=="Best" & Threat=="Buffalos"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category,levels = c("Waterbirds","Turtles","Fish","F05"),labels = c("Waterbirds","Turtles","Fish","Grunters"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="Probability of persistence", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E","grey85"))+
  ggtitle("Buffalo")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_buffalo

aqua_barrier<-ggplot(data=subset(aquatic_final,Estimate=="Best" & Threat=="Longitudinal barriers"), aes(y=Persistence_fit,x=Level))+
  geom_boxplot(aes(fill=factor(Category,levels = c("Waterbirds","Turtles","Fish","F05"),labels = c("Waterbirds","Turtles","Fish","Grunters"))),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(0.0,1.0))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E","grey85"))+
  ggtitle("Longitudinal barriers")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_barrier

Extract legend

legend_aqua<-get_legend(aqua_test_large+
                          theme(legend.direction = "vertical", legend.justification = "center"))

Summary figure terrestrial categories

aqua_sum_F05<-plot_grid(aqua_grazing,legend_aqua,aqua_toads,aqua_weeds,aqua_pigs,aqua_flow,aqua_buffalo,aqua_barrier, labels = c("a","","b","c","d","e","f","g"), ncol = 2, align="hv",axis = "b",
                    rel_widths = c(1,1,1,1,1,1,1,1)) #
aqua_sum_F05

save_plot("../docs/aqua_sum_water_grunters.pdf",aqua_sum_F05,base_asp=1.3,ncol=2,nrow=4)

Confidence summary figure

aqua_grazing2<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Grazing"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Confidence (%)", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Grazing")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_grazing2

aqua_toads2<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Cane toads"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Confidence (%)", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Cane toad")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_toads2

aqua_weeds2<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Aquatic weeds"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Aquatic weeds")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_weeds2

aqua_pigs2<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Pigs"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="Confidence (%)", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Pig")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_pigs2

aqua_flow2<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Altered flow regime"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="")+
  scale_y_continuous(name="", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Altered flow regime")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_flow2

aqua_buffalo2<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Buffalos"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="Confidence (%)", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Buffalo")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_buffalo2

aqua_barrier2<-ggplot(data=subset(aquatic,Estimate=="Best" & Threat=="Longitudinal barriers"), aes(y=Confidence,x=Level))+
  geom_boxplot(aes(fill=Category),alpha=0.9,outlier.shape = NA)+
  scale_x_discrete(name="Threat level")+
  scale_y_continuous(name="", limits=c(50,100))+
  theme_cowplot()+
  background_grid(major=c("y"))+
  scale_fill_manual(values = c("#CFE8F3", "#73BFE2", "#12719E"))+
  ggtitle("Longitudinal barriers")+
  theme(legend.position = "none",plot.title = element_text(hjust = 0.5))
  #scale_fill_viridis_d()
aqua_barrier2

aqua_sum_conf<-plot_grid(aqua_grazing2,legend_aqua,aqua_toads2,aqua_weeds2,aqua_pigs2,aqua_flow2,aqua_buffalo2,aqua_barrier2, labels = c("a","","b","c","d","e","f","g"), ncol = 2, align="hv",axis = "b",
                    rel_widths = c(1,1,1,1,1,1,1,1)) #
aqua_sum_conf

save_plot("../docs/aqua_sum_conf_water.pdf",aqua_sum_conf,base_asp=1.3,ncol=2,nrow=4)

Statistical analysis (Bayesian)

Advice on the implementation of Zero-One-Inflated Beta models can be found here.

Fit the model

Model validation

Predictions

Summary figure

Final figure