Load the necessary libraries
#library(tidyverse) #for data wrangling, includes dplyr and ggplot2
library(dplyr)
library(ggplot2)
library(pspearman) #for spearman correlation
library(tidybayes) #for data wrangling
library(performance) #for VIF estimation
library(GGally) #for plotting graphs
library(ggpubr) #for combining multiple plots
library(DHARMa) #for residuals and diagnostics
library(rstan) #for interfacing with STAN
library(brms) #for fitting models in STAN
library(ggmcmc) #for MCMC diagnostics
library(mice) #for imputing missing data
library(ppcor) #for estimating partial correlations
library(ggnewscale) #for manual plotting of Bayesian results
library(MCMCvis) #for automatic plotting of Bayesian results
library(tidybayes) #for tidymcmc Bayesian results
library(ggeffects) #for model testing
library(ggridges) #for plotting
library(stringr) #for pars variables
The aim of the analysis is to investigate whether a comprehensive operationalization of actor-specific adaptive capacity could usefully explain variations in observed responses among tourism operators, while controlling for contextual conditions related to the disturbance characteristics and governance settings where the tourism operators were based.
As described in greater detail in our manuscript, we conducted surveys with representatives of reef tourism companies (operators) in the APAC Region. The surveys requested information on the actions that each operator took in response to a specific climate disturbance and indicators of the operator’s adaptive capacity before the disturbance. Our study design included a priori treatment and control groups of tourism operators based on whether their locations had been directly affected by the climate disturbance. We used a combination of descriptive and multivariate statistical approaches to explore the relationships between adaptive capacity, adaptive responses, and contextual variables that were relevant to the relationship between adaptive capacity and adaptive responses.
data = read.csv('HB_adaptive_capacity_raw.csv')
glimpse(data)
## Rows: 231
## Columns: 38
## $ operator <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, …
## $ location <chr> "id", "id", "id", "id", "mi", "id", "id", "fp", "id…
## $ gov_effect <dbl> -0.04, -0.04, -0.04, -0.04, 0.44, -0.04, -0.04, 1.4…
## $ resp_region <chr> "eur", "eur", "other", "nam", "nam", "nam", "nam", …
## $ resp_age <chr> "45 - 54", "35 - 44", "35 - 44", "55 - 64", "55 - 6…
## $ resp_gender <int> 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, …
## $ dist_type <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ dist_severity <dbl> 0.00, 0.75, 0.50, 0.00, 0.50, 0.25, 0.25, 0.00, 0.0…
## $ divbl <int> 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 1, 0, …
## $ chgsites <int> 0, 1, 0, 0, 1, 1, 1, 0, 0, 0, 0, 1, 0, 0, 0, 1, 1, …
## $ chgact <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, …
## $ chgopmode <int> 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 1, …
## $ nrm <int> 1, 1, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, …
## $ insurance <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ monitor <int> 1, 1, 0, 0, 1, 1, 1, 0, 1, 0, 1, 0, 0, 0, 1, 0, 1, …
## $ relief <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ support <int> 0, 1, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ co2 <int> 0, 1, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, 0, 0, 1, 0, 1, …
## $ educate <int> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ none <int> 0, 0, 1, 1, 0, 0, 0, 1, 0, 0, 0, 0, 1, 1, 0, 0, 0, …
## $ topresponse <chr> "monitor", "monitor", "none", "none", "nrm", "chgsi…
## $ exclaccess <int> 1, 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
## $ particmgmt_tourism <int> 1, 2, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 2, 2, 2, 2, …
## $ particmgmt_coral <int> 1, 2, 1, 0, 0, 1, 1, 1, 1, 1, 1, 0, 0, 2, 1, 2, 2, …
## $ psgseats <int> 1, 1, 1, 5, 1, 2, 2, 1, 2, 2, 1, 3, 1, 1, 2, 1, 2, …
## $ savings <int> 2, 3, 2, 4, 1, 2, 3, 1, 3, 4, 1, 3, 2, 2, 3, 2, 4, …
## $ distmarket <int> 90, 90, 60, 90, 17, 135, 135, 5, 50, 18, 7, 46, 4, …
## $ othservices <int> 2, 3, 0, 0, 0, 2, 3, 0, 2, 2, 1, 7, 0, 0, 0, 2, 0, …
## $ accsites <int> 4, 5, 3, 5, 2, 5, 5, 2, 3, 5, 5, 5, 1, 4, 4, 5, 2, …
## $ edu <int> 3, 4, 3, 4, 1, 3, 3, 3, 3, 3, 3, 4, 4, 3, 4, 3, 1, …
## $ experience <int> 3, 3, 1, 4, 4, 3, 3, 3, 2, 3, 4, 3, 2, 4, 3, 4, 3, …
## $ employ_involve <int> 2, 3, 3, 2, 3, 3, 3, 0, 2, 2, 1, 2, 2, 0, 2, 2, 2, …
## $ ind_memberships <int> 1, 1, 0, 1, 0, 0, 1, 2, 0, 1, 0, 0, 1, 2, 0, 2, 2, …
## $ research_rely <int> 1, 2, 3, 0, 2, 0, 0, 0, 2, 1, 1, 2, 1, 0, 0, 0, 2, …
## $ govt_ties <int> 2, 1, 1, 1, 2, 1, 1, 2, 1, 1, 1, 1, 1, 2, 1, 1, 2, …
## $ reef_attach <int> 1, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 1, 0, 2, 2, 1, 2, …
## $ comp_concerns <int> 1, 0, 0, 3, 0, 1, 1, 0, 2, 2, 1, 2, 0, 0, 1, 1, 0, …
## $ adapt_conf <int> 2, 2, 2, 3, 0, 3, 3, 2, 2, 2, 1, 2, 0, 2, 0, 1, 2, …
We first evaluated the paired correlation between our predictors to see whether there are any highly correlated variables. The highest correlation was found, as we expected between our two indicators of management participation (in reef tourism mgmt., and in coral reef mgmt.). These predictors had a correlation of 74% using Spearman’s rank test. All other correlations were below 40%. We merged our two management participation indicators into a combined predictor that took the maximum of any of the two values.
For participation in management (x-scale): 0 for none, 1 for passive attendance (attended meetings but did not usually speak up), 2 for active participation (spoke up in meetings), and 3 for power (had a formal position / power over decision making)
cor.test(data$particmgmt_tourism, data$particmgmt_coral,
alternative = c("two.sided"),
method = c("spearman"))
## Warning in cor.test.default(data$particmgmt_tourism, data$particmgmt_coral, :
## Cannot compute exact p-value with ties
##
## Spearman's rank correlation rho
##
## data: data$particmgmt_tourism and data$particmgmt_coral
## S = 536050, p-value < 2.2e-16
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
## rho
## 0.7390672
data <- data %>%
mutate(particmgmt = pmax(particmgmt_coral,particmgmt_tourism)
)
ggplot(data, aes(x=particmgmt)) +
geom_bar()
We used the ‘mice’ package (Buuren and Groothuis-Oudshoorn, 2011; Zhang, 2016) in R modeling software to impute missing values in our dataset. For the 230 reef tourism operator surveys in our dataset, we had missing data for savings (15), competing concerns (7), adaptation confidence (5), education (1), accessible dive / snorkel sites (1), and reef attachment (1). Without imputing missing data, we would have had to exclude 20 operators from our sample and this would have led to the exclusion of other relevant data for these operators. We imputed the missing data based on the rest of the dataset but we did not use the disturbance type, severity, topresponse as predictors for missing values. We only used the combined mgmt. participation indicator.
imp0 <- mice(data, maxit = 0) #First, do a set-up run of mice(), without any iterations (maxit = 0). This is to get all the objects that we may need to adjust in order to customize the imputation to our data.
## Warning: Number of logged events: 4
pred <- imp0$predictorMatrix #in the code below we tell the model which variables not to use for predicting missing values
pred[c("savings", "accsites", "edu", "reef_attach", "comp_concerns", "adapt_conf"), "dist_type"] <- 0
pred[c("savings", "accsites", "edu", "reef_attach", "comp_concerns", "adapt_conf"), "dist_severity"] <- 0
pred[c("savings", "accsites", "edu", "reef_attach", "comp_concerns", "adapt_conf"), "operator"] <- 0
pred[c("savings", "accsites", "edu", "reef_attach", "comp_concerns", "adapt_conf"), "topresponse"] <- 0
pred[c("savings", "accsites", "edu", "reef_attach", "comp_concerns", "adapt_conf"), "particmgmt_coral"] <- 0
pred[c("savings", "accsites", "edu", "reef_attach", "comp_concerns", "adapt_conf"), "particmgmt_tourism"] <- 0
imp <- mice(data, predictorMatrix = pred, m = 1, print = FALSE,meth='pmm', seed=123) #impute data using predictive mean matching. We set a seed value to offset the random number generator (so that our imputed results have reproducibility, without seed there would be variations in the predictions each time we run the model)
densityplot(imp) #compare the distribution of the imputed values against the distribution of the observed values.
imp$imp$savings #check the imputed values for 'savings' indicator
imp$imp$comp_concerns
imp$imp$adapt_conf
imp$imp$edu
imp$imp$accsites
imp$imp$reef_attach
imputed_data <- complete(imp,1) #The missing values have been replaced with the imputed values and we now have a complete dataset to work with in our consequent analysis.
We collected data on the indicators in our surveys using multiple-choice categories to provide a consistent and directly comparable level of detail in the answers. Based on the distribution of the collected data, we transformed each indicator into a binary variable. While transforming some of our indicators into continuous or ordinal variables would likely have provided slightly more nuance, the nature of our research warranted the comparison of indicators at similar levels of complexity. We asked respondents about each indicator of adaptive capacity as they judged it to be before they were affected by the particular climate disturbance. Below we detail how we decided to cut-off ordinal and continuous into binary variables. These graphs do not include missing data point, which will be handled later.
For government effectiveness (‘gov_effect’): -0.37 for Fiji, -0.04 for Indonesia, 0.44 for Mariana Islands, 1.42 for the Hawaiian Islands, 1.42 for French Polynesia, 1.61 for Australia, and 1.81 for JP. For exclusivity of reef access (‘exclaccess’): 0 for open access (accessed same reef sites as everyone else), 1 for limited access (accessed sites that were shared with only a few other operators), and 2 for exclusive access (accessed sites where no one else was able to go). For participation in management (‘particmgmt_coral’): 0 for none, 1 for passive attendance (attended meetings but did not usually speak up), 2 for active participation (spoke up in meetings), and 3 for power (had a formal position / power over decision making). For savings, 1 for 0 – 3 months, 2 for 3 – 6 months, 3 for 6 – 12 months, and 4 for >12 months.
For boats (‘psgseats’): 1 for 0 – 10 seats, 2 for 10 – 20 seats, 3 for 20 – 50 seats, 4 for 50 – 100 seats, 5 for 100 – 200 seats, 6 for 200 – 300 seats, 7 for 300 – 400 seats, 8 for 400 – 500 seats, and 9 for >500 seats. For other services (‘othservices’): The number reflects the other types of products / services the operator was selling besides tourism activities (e.g. accommodation, first-aid training). For accessible dive / snorkel sites (‘accsites’): 1 for 0 – 2 sites, 2 for 2 – 5 sites, 3 for 5 – 10 sites, 4 for 10 – 20 sites, and 5 for >20 sites. For (GM) education (‘edu_ordinal’): 1 for secondary school, 2 for trade certification, 3 for graduate degree, and 4 for postgraduate degree.
For (GM) experience (‘experience’): 1 for 0 – 2 years, 2 for 2 – 5 years, 3 for 5 – 10 years, and 4 for >10 years. For employee involvement (‘employ_involve’): 0 for not involved, 1 for passively involved (they were involved prior to change), 2 for moderately involved (they were asked for, and voiced their opinions, and these were taken into account), and 3 for highly involved (they actively participated in the decision making process). For government ties (‘govt_ties’): 1 for weak ties (did not interact much with government stakeholder), 2 for moderate ties (had some interaction with government stakeholders, mainly through formal events), and 3 for strong ties (management or staff had informal ties with key government stakeholders). For reliance on research in decision making (‘research_rely’): 0 for not at all, 1 for low reliance, 2 for moderate reliance, and 3 for high reliance.
For industry memberships (‘ind_memberships’): The number reflects the number of tourism industry associations the operator was a member of. For reef attachment (‘reef_attach’): 0 for coral reef sites were not at all important in company’s identity, 2 for coral reef sites were somewhat important in company’s identity, and 3 for coral reef sites were major part of company’s identity. For competing concerns (‘comp_concerns’): 0 for no competing concerns, 1 for 1 – 2 competing concerns, 2 for 3 – 4 competing concerns, and 3 for >5 competing concerns. For adaptation confidence in ability to adapt to climate disturbances (‘adapt_conf’): 0 for not at all confident, 1 for low confidence, 2 for moderately confident, and 3 for highly confident.
We now apply these modifications to our data using the code below. We also create z-scored variables for the continuous predictors in our models: disturbance severity and distance to markets
data_modified <- imputed_data %>%
mutate(operator = factor(operator),
resp_gender = factor(resp_gender),
resp_region = factor (resp_region),
resp_age = factor(resp_age),
location = factor(location),
gov_effect_binary = ifelse(gov_effect > 0.5, 1, 0),
exclaccess_binary = ifelse(exclaccess > 0, 1, 0),
particmgmt_binary = ifelse(particmgmt > 1, 1, 0),
psgseats_binary = ifelse(psgseats > 2, 1, 0),
savings_binary = ifelse(savings > 2, 1, 0),
othservices_binary = ifelse(othservices > 0, 1, 0),
accsites_binary = ifelse(accsites > 3, 1, 0),
edu_binary = ifelse(edu > 2, 1, 0),
experience_binary = ifelse(experience > 3, 1, 0),
employ_involve_binary = ifelse(employ_involve > 1, 1, 0),
ind_memberships_binary = ifelse(ind_memberships > 0, 1, 0),
research_rely_binary = ifelse(research_rely > 0, 1, 0),
govt_ties_binary = ifelse(govt_ties > 1, 1, 0),
reef_attach_binary = ifelse(reef_attach > 1, 1, 0),
comp_concerns_binary = ifelse(comp_concerns > 0, 1, 0),
adapt_conf_binary = ifelse(adapt_conf > 1, 1, 0),
z.dist_severity = (dist_severity-mean(dist_severity))/(2*sd(dist_severity)),
z.distmarket = (distmarket-mean(distmarket))/(2*sd(distmarket))
)
We now calculate the partial correlations between the adaptive responses that operators adopted in response to climate disturbances. These partial correlations reflect whether particular responses were more frequently implemented together than others. We used Spearman’s Rank correlation because our responses are measures on a binary scale.
response_cor <- data_modified[, c("divbl","relief","support","insurance","chgopmode","chgact", "educate","chgsites","monitor","nrm","co2")]
pcor(response_cor,method="spearman")
## $estimate
## divbl relief support insurance chgopmode
## divbl 1.000000000 -0.020357860 0.07702295 0.01024531 0.20864800
## relief -0.020357860 1.000000000 0.19952106 -0.02632262 0.05180185
## support 0.077022948 0.199521055 1.00000000 -0.01375645 0.02830850
## insurance 0.010245306 -0.026322622 -0.01375645 1.00000000 -0.07616094
## chgopmode 0.208648005 0.051801846 0.02830850 -0.07616094 1.00000000
## chgact 0.037796189 0.003162198 0.06836009 -0.08821421 0.35515039
## educate -0.127592314 -0.114291373 -0.03133336 0.20047507 0.13208329
## chgsites 0.042165317 0.171738350 -0.06275669 0.18695861 0.18979091
## monitor 0.031665550 -0.107723152 0.07054858 -0.10563299 -0.03138290
## nrm -0.000335073 0.046177719 0.28126118 -0.00538762 -0.10650199
## co2 -0.012193354 0.054519949 -0.03944470 0.05810986 0.10597687
## chgact educate chgsites monitor nrm
## divbl 0.037796189 -0.12759231 0.04216532 0.03166555 -0.000335073
## relief 0.003162198 -0.11429137 0.17173835 -0.10772315 0.046177719
## support 0.068360085 -0.03133336 -0.06275669 0.07054858 0.281261181
## insurance -0.088214211 0.20047507 0.18695861 -0.10563299 -0.005387620
## chgopmode 0.355150391 0.13208329 0.18979091 -0.03138290 -0.106501985
## chgact 1.000000000 0.09800120 0.24784670 -0.07271320 0.106511089
## educate 0.098001196 1.00000000 -0.02280847 0.19694840 -0.012888401
## chgsites 0.247846702 -0.02280847 1.00000000 0.19071162 0.099265614
## monitor -0.072713196 0.19694840 0.19071162 1.00000000 0.253184790
## nrm 0.106511089 -0.01288840 0.09926561 0.25318479 1.000000000
## co2 0.006726194 -0.09624868 -0.05649930 0.24329687 0.326590348
## co2
## divbl -0.012193354
## relief 0.054519949
## support -0.039444704
## insurance 0.058109855
## chgopmode 0.105976873
## chgact 0.006726194
## educate -0.096248682
## chgsites -0.056499300
## monitor 0.243296870
## nrm 0.326590348
## co2 1.000000000
##
## $p.value
## divbl relief support insurance chgopmode
## divbl 0.000000000 0.762923308 2.531067e-01 0.879349540 1.773768e-03
## relief 0.762923308 0.000000000 2.825634e-03 0.696498264 4.424939e-01
## support 0.253106676 0.002825634 0.000000e+00 0.838494957 6.748579e-01
## insurance 0.879349540 0.696498264 8.384950e-01 0.000000000 2.584727e-01
## chgopmode 0.001773768 0.442493875 6.748579e-01 0.258472688 0.000000e+00
## chgact 0.575362391 0.962632918 3.105928e-01 0.190366218 5.323994e-08
## educate 0.057680843 0.089345380 6.424078e-01 0.002693854 4.935562e-02
## chgsites 0.531983880 0.010363087 3.520125e-01 0.005196749 4.544065e-03
## monitor 0.638883591 0.109456883 2.953207e-01 0.116553592 6.418817e-01
## nrm 0.996039085 0.493650869 2.107227e-05 0.936379511 1.135607e-01
## co2 0.856636039 0.418892263 5.588015e-01 0.388875769 1.153620e-01
## chgact educate chgsites monitor nrm
## divbl 5.753624e-01 0.057680843 0.5319838798 0.6388835912 9.960391e-01
## relief 9.626329e-01 0.089345380 0.0103630867 0.1094568834 4.936509e-01
## support 3.105928e-01 0.642407764 0.3520125277 0.2953207303 2.107227e-05
## insurance 1.903662e-01 0.002693854 0.0051967492 0.1165535916 9.363795e-01
## chgopmode 5.323994e-08 0.049355623 0.0045440655 0.6418816509 1.135607e-01
## chgact 0.000000e+00 0.145545398 0.0001912091 0.2807149326 1.135297e-01
## educate 1.455454e-01 0.000000000 0.7353902986 0.0032106414 8.485595e-01
## chgsites 1.912091e-04 0.735390299 0.0000000000 0.0043483617 1.403958e-01
## monitor 2.807149e-01 0.003210641 0.0043483617 0.0000000000 1.370727e-04
## nrm 1.135297e-01 0.848559511 0.1403958382 0.0001370727 0.000000e+00
## co2 9.206195e-01 0.152918672 0.4021780239 0.0002524891 6.500265e-07
## co2
## divbl 8.566360e-01
## relief 4.188923e-01
## support 5.588015e-01
## insurance 3.888758e-01
## chgopmode 1.153620e-01
## chgact 9.206195e-01
## educate 1.529187e-01
## chgsites 4.021780e-01
## monitor 2.524891e-04
## nrm 6.500265e-07
## co2 0.000000e+00
##
## $statistic
## divbl relief support insurance chgopmode
## divbl 0.000000000 -0.30201846 1.1458389 0.15197042 3.1643958
## relief -0.302018459 0.00000000 3.0200990 -0.39056291 0.7693785
## support 1.145838866 3.02009900 0.0000000 -0.20406038 0.4200513
## insurance 0.151970418 -0.39056291 -0.2040604 0.00000000 -1.1329399
## chgopmode 3.164395843 0.76937853 0.4200513 -1.13293988 0.0000000
## chgact 0.561008940 0.04690321 1.0163214 -1.31354903 5.6350876
## educate -1.908095283 -1.70639656 -0.4649771 3.03514292 1.9764280
## chgsites 0.625969422 2.58570829 -0.9326706 2.82281677 2.8671662
## monitor 0.469911653 -1.60714464 1.0490183 -1.57560570 -0.4657130
## nrm -0.004969936 0.68565769 4.3472714 -0.07991248 -1.5887156
## co2 -0.180870115 0.80986605 -0.5855152 0.86336736 1.5807931
## chgact educate chgsites monitor nrm co2
## divbl 0.56100894 -1.9080953 0.6259694 0.4699117 -0.004969936 -0.18087012
## relief 0.04690321 -1.7063966 2.5857083 -1.6071446 0.685657692 0.80986605
## support 1.01632139 -0.4649771 -0.9326706 1.0490183 4.347271371 -0.58551518
## insurance -1.31354903 3.0351429 2.8228168 -1.5756057 -0.079912476 0.86336736
## chgopmode 5.63508758 1.9764280 2.8671662 -0.4657130 -1.588715557 1.58079314
## chgact 0.00000000 1.4606236 3.7945535 -1.0813735 1.588852923 0.09976783
## educate 1.46062365 0.0000000 -0.3383923 2.9795752 -0.191181764 -1.43425746
## chgsites 3.79455355 -0.3383923 0.0000000 2.8815991 1.479655056 -0.83936081
## monitor -1.08137349 2.9795752 2.8815991 0.0000000 3.881814847 3.72046889
## nrm 1.58885292 -0.1911818 1.4796551 3.8818148 0.000000000 5.12515025
## co2 0.09976783 -1.4342575 -0.8393608 3.7204689 5.125150248 0.00000000
##
## $n
## [1] 231
##
## $gp
## [1] 9
##
## $method
## [1] "spearman"
data_modified <- data_modified %>%
mutate(top_mn = case_when(
topresponse == 'none' | topresponse == 'relief' ~ 'cope',
topresponse == 'chgsites' | topresponse == 'chgact' | topresponse == 'chgopmode' | topresponse == 'divbl' | topresponse == 'insurance' | topresponse == 'monitor' | topresponse == 'educate' ~ 'adapt',
topresponse == 'nrm' | topresponse == 'support'| topresponse == 'co2' ~ 'transform'),
top_mn = factor(top_mn,levels = c('cope', 'adapt', 'transform'))
)
ggplot(data_modified, aes(x=top_mn)) +
geom_bar()
table(data_modified$top_mn)
##
## cope adapt transform
## 63 92 76
The relationship between adaptive capacity and the primary responses to climate disturbances were explored using multinomial logistic regression in a Bayesian framework (Gelman et al., 2013). The Bayesian approach for binary regression extends to multinomial models. We used a baseline-category logit model in which the primary response of ‘coping measures (e.g. doing nothing or looking for relief)’ was used as the baseline category. The likelihood of other response clusters (‘adaptive’ and ‘transformative’ measures) were then estimated against the baseline category. As such with sample size ni and αi representing the probability of an operator either adopting a coping (αi1), adaptive (αi2), or transformative (αi3) primary response to a climate disturbance. The multinomial generalized linear model is parameterized in terms of the logarithm of the ratio of the probability of each category relative to that of our baseline category (coping: αi1), where β(j) is a vector of parameters for each category (Gelman et al., 2013, p. 426). We used weakly-informative normal priors on all population effects. Although we had some indications of the effect of our predictor variables on adaptive responses (Bartelet et al., 2022a), we decided to not use informative priors in our Bayesian models because most of the prior knowledge providing the foundation for our study is based on a different context, i.e., smallholder farming in Asia and Africa.
\[ \begin{align} y_i &\sim{} \mathcal{Multin}(n_i;α_i1, α_i2, α_i3)\\ with\ \alpha_i1 + \alpha_i2 + \alpha_i3 = 1\\ \log(\alpha_i2/\alpha_i1) = X_i \beta^j \\ \log(\alpha_i3/\alpha_i1) = X_i \beta^j \\ \beta_0 &\sim{} \mathcal{N}(0,3)\\ \beta_{1-19} &\sim{} \mathcal{N}(0,3)\\ \end{align} \] ## Test for collinearity
We evaluated whether there was potential for collinearity issues in our model by running a frequentist statistical model (logistic regression), in which we predict the likelihood of transformative action as primary responses versus doing something else. We found that all variance inflation factors (VIFs) are below the critical threshold of ‘3’.
data_modified <- data_modified %>%
mutate(top_transform = ifelse(top_mn == 'transform', 1, 0)
)
model.glm <- glm(top_transform ~ gov_effect_binary + z.distmarket +
dist_type + z.dist_severity +
exclaccess_binary + particmgmt_binary +
psgseats_binary + savings_binary +
othservices_binary + accsites_binary +
edu_binary + experience_binary + employ_involve_binary +
ind_memberships_binary + research_rely_binary + govt_ties_binary +
reef_attach_binary + comp_concerns_binary + adapt_conf_binary,
data=data_modified, family=binomial(link='logit'))
check_collinearity(model.glm)
We now fit the Bayesian model using the ‘categorical’ data distribution family reflecting the multinomial nature of our response data. We set adapt_delta to 0.99 (from the default of 0.95) to reduce the step size in our models which will lead to lower likelihood of divergences. This will usually slow down the sampling but lead to more robust samples.
Setting one prior on all effects makes the computation much easier. We use iteratively-determined priors by running the model first only with priors and seeing whether the outcomes are non-informative and covering the full ribbon (0 to 100%). We use zero slope coefficients to not prime the model on either positive or negative effects. We also use a reasonably wide variance (3) that is reasonable for logistic regression.
priors <- prior(normal(0, 3), class = 'Intercept', dpar = 'muadapt') +
prior(normal(0, 3), class = 'Intercept', dpar = 'mutransform') +
prior(normal(0,3), class = 'b', dpar = 'muadapt') +
prior(normal(0,3), class = 'b', dpar = 'mutransform')
We first run the model with priors and find that our selected priors are uninformative.
mod_top_multinomial = brms::brm(
top_mn ~ gov_effect_binary + z.distmarket +
dist_type + z.dist_severity +
exclaccess_binary + particmgmt_binary +
psgseats_binary + savings_binary +
othservices_binary + accsites_binary +
edu_binary + experience_binary + employ_involve_binary +
ind_memberships_binary + research_rely_binary + govt_ties_binary +
reef_attach_binary + comp_concerns_binary + adapt_conf_binary,
data = data_modified,
prior = priors,
sample_prior='only',
family = 'categorical',
iter = 5000, warmup = 1000,
chains = 3, cores = 3,
thin = 5,
refresh = 0,
control = list(adapt_delta = 0.99),
seed = 123)
## Compiling Stan program...
## Start sampling
mod_top_multinomial %>% ggpredict() %>% plot(add.data=TRUE, jitter=FALSE)
## $gov_effect_binary
##
## $z.distmarket
##
## $dist_type
##
## $z.dist_severity
##
## $exclaccess_binary
##
## $particmgmt_binary
##
## $psgseats_binary
##
## $savings_binary
##
## $othservices_binary
##
## $accsites_binary
##
## $edu_binary
##
## $experience_binary
##
## $employ_involve_binary
##
## $ind_memberships_binary
##
## $research_rely_binary
##
## $govt_ties_binary
##
## $reef_attach_binary
##
## $comp_concerns_binary
##
## $adapt_conf_binary
After checking the priors, run the Bayesian model using both prior and data information.
mod_top_multinomial <- mod_top_multinomial %>%
update(sample_prior='yes', refresh = 0)
## The desired updates require recompiling the model
## Compiling Stan program...
## Start sampling
preds <- posterior_predict(mod_top_multinomial, ndraws=250,summary=FALSE)
fert.resids <- createDHARMa(simulatedResponse=t(preds),
observedResponse = as.numeric(as.factor(data_modified$top_mn)),
fittedPredictedResponse=apply(preds,2,median),
integerResponse=FALSE)
testResiduals (fert.resids)
## $uniformity
##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: simulationOutput$scaledResiduals
## D = 0.040303, p-value = 0.8472
## alternative hypothesis: two-sided
##
##
## $dispersion
##
## DHARMa nonparametric dispersion test via sd of residuals fitted vs.
## simulated
##
## data: simulationOutput
## dispersion = 1.0138, p-value = 0.872
## alternative hypothesis: two.sided
##
##
## $outliers
##
## DHARMa outlier test based on exact binomial test with approximate
## expectations
##
## data: simulationOutput
## outliers at both margin(s) = 1, observations = 231, p-value = 1
## alternative hypothesis: true probability of success is not equal to 0.007968127
## 95 percent confidence interval:
## 0.0001095949 0.0238823445
## sample estimates:
## frequency of outliers (expected: 0.00796812749003984 )
## 0.004329004
## $uniformity
##
## Asymptotic one-sample Kolmogorov-Smirnov test
##
## data: simulationOutput$scaledResiduals
## D = 0.040303, p-value = 0.8472
## alternative hypothesis: two-sided
##
##
## $dispersion
##
## DHARMa nonparametric dispersion test via sd of residuals fitted vs.
## simulated
##
## data: simulationOutput
## dispersion = 1.0138, p-value = 0.872
## alternative hypothesis: two.sided
##
##
## $outliers
##
## DHARMa outlier test based on exact binomial test with approximate
## expectations
##
## data: simulationOutput
## outliers at both margin(s) = 1, observations = 231, p-value = 1
## alternative hypothesis: true probability of success is not equal to 0.007968127
## 95 percent confidence interval:
## 0.0001095949 0.0238823445
## sample estimates:
## frequency of outliers (expected: 0.00796812749003984 )
## 0.004329004
plot(fert.resids)
We tested model convergence and diagnostics (traceplots, autocorrelation, convergence, and effective sample size). Tests reveal no problems with our model convergence.
pars <- mod_top_multinomial$fit %>% get_variables()
pars_adapt <- str_extract(pars, '^b_muadapt_.*') %>% na.omit()
pars_transform <- str_extract(pars, '^b_mutransform_.*') %>% na.omit()
mod_top_multinomial$fit %>% stan_trace(pars=pars_adapt, inc_warmup=TRUE)
mod_top_multinomial$fit %>% stan_trace(pars=pars_transform, inc_warmup=TRUE)
mod_top_multinomial$fit %>% stan_ac(pars=pars_adapt) #autocorrelation test (you only want autocorrelation at beginning of plot. If evidence of autocorrelation, use more thinning. If lag 2 is greater than 0.2 you might have autocorrelation, adjust tinning to get beyond 0.2)
mod_top_multinomial$fit %>% stan_ac(pars=pars_transform)
mod_top_multinomial$fit %>% stan_rhat() #measure of convergence between chains. We want values to be less than 1.03 or 1.01. Values greater than 1.05 shows that chains were not consistent, they did not converge. There is a bar for every parameter. If problem with converge, run it longer and potentially narrow the priors.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
mod_top_multinomial$fit %>% stan_ess() #effective sample size. If chains don’t converge a lot of samples will be thrown away. This tests the efficiency of the sampler. If too many samples are rejected you don’t them contributing. Above 50% effective sample size is high enough. Less than 50% means sampler has drifted into unsupported areas. You need to tighten the priors.
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
test whether prior and posterior match for different predictors and responses. For example:
#mod_top_multinomial %>% get_variables()
mod_top_multinomial %>% hypothesis('muadapt_savings_binary=0') %>% plot()
mod_top_multinomial %>% hypothesis('mutransform_comp_concerns_binary=0') %>% plot()
We use median highest probability density intervals and present credibility intervals at both 80% (weak evidence) and 95% (strong evidence).
Strong evidence: log odds scale
mod_top_multinomial %>% tidy_draws() %>% dplyr::select(starts_with("b_")) %>%
summarise_draws(median,~HDInterval::hdi(.x, credMass = 0.95),rhat,ess_bulk,ess_tail) %>% knitr::kable()
variable | median | lower | upper | rhat | ess_bulk | ess_tail |
---|---|---|---|---|---|---|
b_muadapt_Intercept | 0.7896507 | -0.5707917 | 2.0941099 | 1.0008616 | 2571.411 | 2388.281 |
b_mutransform_Intercept | 0.0670525 | -1.2767953 | 1.5644713 | 1.0002720 | 2411.670 | 1700.096 |
b_muadapt_gov_effect_binary | -0.5225526 | -1.3643373 | 0.3846161 | 0.9996307 | 2421.821 | 2402.461 |
b_muadapt_z.distmarket | -0.6829090 | -1.5891268 | 0.1761063 | 1.0013428 | 2207.214 | 2446.568 |
b_muadapt_dist_type | -0.8526734 | -2.0980241 | 0.2723958 | 1.0001378 | 2283.985 | 2403.939 |
b_muadapt_z.dist_severity | 1.8699248 | 1.0372663 | 2.8064329 | 1.0002079 | 1998.903 | 2232.692 |
b_muadapt_exclaccess_binary | 0.1784635 | -0.7582171 | 1.0884664 | 1.0000032 | 2526.078 | 2491.499 |
b_muadapt_particmgmt_binary | 0.1266520 | -0.7636110 | 1.0242901 | 1.0024361 | 2573.241 | 2028.750 |
b_muadapt_psgseats_binary | -0.3566502 | -1.1862244 | 0.5644566 | 0.9998594 | 2515.490 | 2361.012 |
b_muadapt_savings_binary | 0.2700554 | -0.4824477 | 1.0701633 | 0.9999330 | 2304.356 | 2028.181 |
b_muadapt_othservices_binary | 0.3527455 | -0.4330590 | 1.1631295 | 0.9995843 | 2376.174 | 2321.149 |
b_muadapt_accsites_binary | -0.3448160 | -1.1143163 | 0.4192817 | 0.9997169 | 2161.249 | 2021.134 |
b_muadapt_edu_binary | -0.1493597 | -0.9888886 | 0.6771764 | 1.0003925 | 2498.604 | 2489.500 |
b_muadapt_experience_binary | -0.6657694 | -1.4866936 | 0.1902605 | 0.9997352 | 2341.851 | 2209.791 |
b_muadapt_employ_involve_binary | -0.1325657 | -0.8930583 | 0.7040542 | 1.0002235 | 2288.784 | 2282.263 |
b_muadapt_ind_memberships_binary | 0.8142634 | -0.0118199 | 1.6387235 | 1.0002858 | 2289.751 | 2251.340 |
b_muadapt_research_rely_binary | 0.5983230 | -0.2101715 | 1.5074950 | 0.9997970 | 2449.534 | 2403.939 |
b_muadapt_govt_ties_binary | -0.3869888 | -1.3596592 | 0.5281761 | 1.0002337 | 2317.363 | 2301.522 |
b_muadapt_reef_attach_binary | -0.2293992 | -0.9913047 | 0.5620005 | 1.0002022 | 2290.691 | 2289.754 |
b_muadapt_comp_concerns_binary | 0.0224119 | -0.7269714 | 0.8204324 | 1.0006043 | 2255.706 | 2281.905 |
b_muadapt_adapt_conf_binary | 0.0550231 | -0.8848561 | 0.9381379 | 0.9996399 | 2478.925 | 2026.598 |
b_mutransform_gov_effect_binary | -1.1892413 | -2.2207860 | -0.2254530 | 1.0004771 | 2474.323 | 2185.953 |
b_mutransform_z.distmarket | 0.4316454 | -0.4663618 | 1.3228894 | 0.9996506 | 2395.812 | 2397.392 |
b_mutransform_dist_type | -1.1525783 | -2.4497045 | 0.2416389 | 0.9998151 | 2378.624 | 2245.774 |
b_mutransform_z.dist_severity | 0.9272731 | 0.0018739 | 1.8896830 | 1.0002159 | 2255.098 | 2304.976 |
b_mutransform_exclaccess_binary | -0.2385664 | -1.2819340 | 0.7744782 | 0.9999840 | 2483.451 | 2125.298 |
b_mutransform_particmgmt_binary | 0.5052244 | -0.4361672 | 1.4344391 | 1.0004089 | 2445.127 | 2283.747 |
b_mutransform_psgseats_binary | 0.6958288 | -0.1800360 | 1.6172336 | 0.9999440 | 2400.006 | 2208.597 |
b_mutransform_savings_binary | 0.0103069 | -0.8387588 | 0.8880051 | 0.9996355 | 2234.431 | 2242.833 |
b_mutransform_othservices_binary | -0.0692409 | -0.8909545 | 0.7944451 | 0.9997495 | 2473.133 | 2471.211 |
b_mutransform_accsites_binary | 0.0485982 | -0.8703754 | 0.9321499 | 1.0001868 | 2092.162 | 2051.294 |
b_mutransform_edu_binary | -0.3284547 | -1.1369170 | 0.6146089 | 1.0007582 | 2557.956 | 2248.297 |
b_mutransform_experience_binary | -0.2103759 | -1.1678801 | 0.6649130 | 1.0011946 | 2248.228 | 2020.691 |
b_mutransform_employ_involve_binary | 0.7677624 | -0.0591760 | 1.7498694 | 1.0001015 | 2339.526 | 2263.377 |
b_mutransform_ind_memberships_binary | 1.0782039 | 0.2372243 | 1.9834708 | 0.9998260 | 2314.764 | 2362.630 |
b_mutransform_research_rely_binary | 1.0827900 | 0.1683884 | 1.9822058 | 0.9999343 | 2484.064 | 2211.670 |
b_mutransform_govt_ties_binary | 0.7790557 | -0.1886954 | 1.7402248 | 0.9996178 | 2452.987 | 2282.263 |
b_mutransform_reef_attach_binary | 0.1680752 | -0.6719244 | 1.0598517 | 0.9997563 | 2251.871 | 2099.048 |
b_mutransform_comp_concerns_binary | -0.6072639 | -1.5106822 | 0.2136596 | 1.0001265 | 2327.044 | 2234.775 |
b_mutransform_adapt_conf_binary | -1.2904208 | -2.2719004 | -0.3514563 | 0.9997993 | 2388.340 | 2362.428 |
Strong evidence: odds scale
mod_top_multinomial %>% tidy_draws() %>% exp() %>% dplyr::select(starts_with("b_")) %>%
summarise_draws(median,~HDInterval::hdi(.x, credMass = 0.95),rhat,ess_bulk,ess_tail) %>% knitr::kable()
variable | median | lower | upper | rhat | ess_bulk | ess_tail |
---|---|---|---|---|---|---|
b_muadapt_Intercept | 2.2026269 | 0.3015927 | 7.0769799 | 1.0013365 | 2571.411 | 2388.281 |
b_mutransform_Intercept | 1.0693516 | 0.1339271 | 3.4688070 | 1.0003494 | 2411.670 | 1700.096 |
b_muadapt_gov_effect_binary | 0.5930049 | 0.1543910 | 1.2597059 | 0.9996492 | 2421.821 | 2402.461 |
b_muadapt_z.distmarket | 0.5051455 | 0.1547902 | 1.0848722 | 1.0010931 | 2207.214 | 2446.568 |
b_muadapt_dist_type | 0.4262739 | 0.0629299 | 1.1863878 | 1.0005549 | 2283.985 | 2403.939 |
b_muadapt_z.dist_severity | 6.4878151 | 1.8870845 | 14.1827026 | 1.0002079 | 1998.903 | 2232.692 |
b_muadapt_exclaccess_binary | 1.1953792 | 0.3277001 | 2.5570613 | 1.0000032 | 2526.078 | 2491.499 |
b_muadapt_particmgmt_binary | 1.1350221 | 0.3403974 | 2.4776616 | 1.0023723 | 2573.241 | 2028.750 |
b_muadapt_psgseats_binary | 0.7000173 | 0.2142925 | 1.4396566 | 0.9998594 | 2515.490 | 2361.012 |
b_muadapt_savings_binary | 1.3100375 | 0.4859544 | 2.6177932 | 0.9999330 | 2304.356 | 2028.181 |
b_muadapt_othservices_binary | 1.4229689 | 0.5265716 | 2.8821483 | 0.9995842 | 2376.174 | 2321.149 |
b_muadapt_accsites_binary | 0.7083507 | 0.2606377 | 1.4033812 | 0.9997169 | 2161.249 | 2021.134 |
b_muadapt_edu_binary | 0.8612592 | 0.3092182 | 1.8009614 | 1.0002441 | 2498.604 | 2489.500 |
b_muadapt_experience_binary | 0.5138780 | 0.1774303 | 1.0255894 | 0.9997352 | 2341.851 | 2209.791 |
b_muadapt_employ_involve_binary | 0.8758454 | 0.3236181 | 1.8358319 | 1.0002879 | 2288.784 | 2282.263 |
b_muadapt_ind_memberships_binary | 2.2575123 | 0.7707048 | 4.6246852 | 1.0003575 | 2289.751 | 2251.340 |
b_muadapt_research_rely_binary | 1.8190658 | 0.6053523 | 3.8275378 | 0.9997970 | 2449.534 | 2403.939 |
b_muadapt_govt_ties_binary | 0.6790987 | 0.1658744 | 1.5293200 | 1.0002576 | 2317.363 | 2301.522 |
b_muadapt_reef_attach_binary | 0.7950112 | 0.2853164 | 1.5526568 | 1.0002867 | 2290.691 | 2289.754 |
b_muadapt_comp_concerns_binary | 1.0226649 | 0.4306755 | 2.0995293 | 1.0006932 | 2255.706 | 2281.905 |
b_muadapt_adapt_conf_binary | 1.0565652 | 0.2554295 | 2.2907905 | 0.9996399 | 2478.925 | 2026.598 |
b_mutransform_gov_effect_binary | 0.3044522 | 0.0746546 | 0.7009630 | 1.0005021 | 2474.323 | 2185.953 |
b_mutransform_z.distmarket | 1.5397890 | 0.4553948 | 3.2913945 | 0.9996506 | 2395.812 | 2397.392 |
b_mutransform_dist_type | 0.3158214 | 0.0533050 | 0.9960727 | 0.9996653 | 2378.624 | 2245.774 |
b_mutransform_z.dist_severity | 2.5276075 | 0.7277176 | 5.7413034 | 1.0002159 | 2255.098 | 2304.976 |
b_mutransform_exclaccess_binary | 0.7877564 | 0.2160790 | 1.9336998 | 0.9999840 | 2483.451 | 2125.298 |
b_mutransform_particmgmt_binary | 1.6573575 | 0.4917394 | 3.7234518 | 1.0004440 | 2445.127 | 2283.747 |
b_mutransform_psgseats_binary | 2.0053704 | 0.5681610 | 4.2707586 | 0.9999606 | 2400.006 | 2208.597 |
b_mutransform_savings_binary | 1.0103602 | 0.3254636 | 2.0607464 | 0.9996147 | 2234.431 | 2242.833 |
b_mutransform_othservices_binary | 0.9331019 | 0.3145282 | 1.9634308 | 0.9997495 | 2473.133 | 2471.211 |
b_mutransform_accsites_binary | 1.0497986 | 0.3207284 | 2.3042775 | 1.0000160 | 2092.162 | 2051.294 |
b_mutransform_edu_binary | 0.7200356 | 0.1866646 | 1.5573561 | 1.0007582 | 2557.956 | 2248.297 |
b_mutransform_experience_binary | 0.8102796 | 0.2202343 | 1.7301554 | 1.0012981 | 2248.228 | 2020.691 |
b_mutransform_employ_involve_binary | 2.1549389 | 0.6031091 | 4.7117546 | 1.0001244 | 2339.526 | 2263.377 |
b_mutransform_ind_memberships_binary | 2.9393961 | 0.8296821 | 6.3490269 | 0.9998909 | 2314.764 | 2362.630 |
b_mutransform_research_rely_binary | 2.9529067 | 0.9164544 | 6.4130129 | 0.9997787 | 2484.064 | 2211.670 |
b_mutransform_govt_ties_binary | 2.1794135 | 0.6798297 | 5.2131803 | 0.9996115 | 2452.987 | 2282.263 |
b_mutransform_reef_attach_binary | 1.1830256 | 0.3851240 | 2.5246507 | 0.9997563 | 2251.871 | 2099.048 |
b_mutransform_comp_concerns_binary | 0.5448396 | 0.1892597 | 1.1427649 | 1.0000649 | 2327.044 | 2234.775 |
b_mutransform_adapt_conf_binary | 0.2751550 | 0.0687799 | 0.6035765 | 0.9997993 | 2388.340 | 2362.428 |
Weak evidence: log odds scale
mod_top_multinomial %>% tidy_draws() %>% dplyr::select(starts_with("b_")) %>%
summarise_draws(median,~HDInterval::hdi(.x, credMass = 0.80),rhat,ess_bulk,ess_tail) %>% knitr::kable()
variable | median | lower | upper | rhat | ess_bulk | ess_tail |
---|---|---|---|---|---|---|
b_muadapt_Intercept | 0.7896507 | -0.0380937 | 1.7356287 | 1.0008616 | 2571.411 | 2388.281 |
b_mutransform_Intercept | 0.0670525 | -0.8337771 | 0.9951409 | 1.0002720 | 2411.670 | 1700.096 |
b_muadapt_gov_effect_binary | -0.5225526 | -1.1205329 | 0.0350963 | 0.9996307 | 2421.821 | 2402.461 |
b_muadapt_z.distmarket | -0.6829090 | -1.1981734 | -0.0546751 | 1.0013428 | 2207.214 | 2446.568 |
b_muadapt_dist_type | -0.8526734 | -1.6025678 | -0.0187438 | 1.0001378 | 2283.985 | 2403.939 |
b_muadapt_z.dist_severity | 1.8699248 | 1.2915806 | 2.4347703 | 1.0002079 | 1998.903 | 2232.692 |
b_muadapt_exclaccess_binary | 0.1784635 | -0.4386495 | 0.7553614 | 1.0000032 | 2526.078 | 2491.499 |
b_muadapt_particmgmt_binary | 0.1266520 | -0.4600502 | 0.6838187 | 1.0024361 | 2573.241 | 2028.750 |
b_muadapt_psgseats_binary | -0.3566502 | -0.8711731 | 0.2383807 | 0.9998594 | 2515.490 | 2361.012 |
b_muadapt_savings_binary | 0.2700554 | -0.2215155 | 0.7723556 | 0.9999330 | 2304.356 | 2028.181 |
b_muadapt_othservices_binary | 0.3527455 | -0.1768710 | 0.8669532 | 0.9995843 | 2376.174 | 2321.149 |
b_muadapt_accsites_binary | -0.3448160 | -0.8640243 | 0.1637284 | 0.9997169 | 2161.249 | 2021.134 |
b_muadapt_edu_binary | -0.1493597 | -0.6583141 | 0.4133722 | 1.0003925 | 2498.604 | 2489.500 |
b_muadapt_experience_binary | -0.6657694 | -1.1474395 | -0.0813802 | 0.9997352 | 2341.851 | 2209.791 |
b_muadapt_employ_involve_binary | -0.1325657 | -0.6602729 | 0.3924227 | 1.0002235 | 2288.784 | 2282.263 |
b_muadapt_ind_memberships_binary | 0.8142634 | 0.2848200 | 1.3499991 | 1.0002858 | 2289.751 | 2251.340 |
b_muadapt_research_rely_binary | 0.5983230 | 0.0317533 | 1.1710166 | 0.9997970 | 2449.534 | 2403.939 |
b_muadapt_govt_ties_binary | -0.3869888 | -0.9452191 | 0.2472671 | 1.0002337 | 2317.363 | 2301.522 |
b_muadapt_reef_attach_binary | -0.2293992 | -0.7550111 | 0.2277267 | 1.0002022 | 2290.691 | 2289.754 |
b_muadapt_comp_concerns_binary | 0.0224119 | -0.5096420 | 0.5048475 | 1.0006043 | 2255.706 | 2281.905 |
b_muadapt_adapt_conf_binary | 0.0550231 | -0.4734947 | 0.7062899 | 0.9996399 | 2478.925 | 2026.598 |
b_mutransform_gov_effect_binary | -1.1892413 | -1.7976446 | -0.4952732 | 1.0004771 | 2474.323 | 2185.953 |
b_mutransform_z.distmarket | 0.4316454 | -0.1568487 | 1.0036763 | 0.9996506 | 2395.812 | 2397.392 |
b_mutransform_dist_type | -1.1525783 | -2.0418327 | -0.3024858 | 0.9998151 | 2378.624 | 2245.774 |
b_mutransform_z.dist_severity | 0.9272731 | 0.3217585 | 1.5276140 | 1.0002159 | 2255.098 | 2304.976 |
b_mutransform_exclaccess_binary | -0.2385664 | -0.8939594 | 0.4510119 | 0.9999840 | 2483.451 | 2125.298 |
b_mutransform_particmgmt_binary | 0.5052244 | -0.1062662 | 1.1126471 | 1.0004089 | 2445.127 | 2283.747 |
b_mutransform_psgseats_binary | 0.6958288 | 0.1060918 | 1.2609721 | 0.9999440 | 2400.006 | 2208.597 |
b_mutransform_savings_binary | 0.0103069 | -0.5187061 | 0.6019658 | 0.9996355 | 2234.431 | 2242.833 |
b_mutransform_othservices_binary | -0.0692409 | -0.5846318 | 0.5337172 | 0.9997495 | 2473.133 | 2471.211 |
b_mutransform_accsites_binary | 0.0485982 | -0.5921316 | 0.6011409 | 1.0001868 | 2092.162 | 2051.294 |
b_mutransform_edu_binary | -0.3284547 | -0.9333058 | 0.2339310 | 1.0007582 | 2557.956 | 2248.297 |
b_mutransform_experience_binary | -0.2103759 | -0.8034390 | 0.3830542 | 1.0011946 | 2248.228 | 2020.691 |
b_mutransform_employ_involve_binary | 0.7677624 | 0.1737649 | 1.3233784 | 1.0001015 | 2339.526 | 2263.377 |
b_mutransform_ind_memberships_binary | 1.0782039 | 0.4694498 | 1.6282849 | 0.9998260 | 2314.764 | 2362.630 |
b_mutransform_research_rely_binary | 1.0827900 | 0.4311505 | 1.5979607 | 0.9999343 | 2484.064 | 2211.670 |
b_mutransform_govt_ties_binary | 0.7790557 | 0.1154874 | 1.4110606 | 0.9996178 | 2452.987 | 2282.263 |
b_mutransform_reef_attach_binary | 0.1680752 | -0.4217511 | 0.7080049 | 0.9997563 | 2251.871 | 2099.048 |
b_mutransform_comp_concerns_binary | -0.6072639 | -1.1572139 | -0.0521454 | 1.0001265 | 2327.044 | 2234.775 |
b_mutransform_adapt_conf_binary | -1.2904208 | -1.9138298 | -0.6752319 | 0.9997993 | 2388.340 | 2362.428 |
Conditional effects
mod_top_multinomial %>% conditional_effects(categorical = TRUE) %>% plot(points=TRUE)