Brood No Round 1
Emily Runnion
2023-01-26
1 Input Pollen Data
Start by imputing pollen data and creating a new data frame with average pollen consumption on a per-colony basis
### Figure out average pollen consumption by treatment
pollen <- read_csv("pollen1.csv", col_types = cols(round = col_factor(levels = c("1",
"2")), treatment = col_factor(levels = c("1",
"2", "3", "4", "5")), replicate = col_factor(levels = c("1",
"2", "3", "4", "5", "6", "7", "9", "11",
"12")), start_date = col_date(format = "%m/%d/%Y"),
start_time = col_character(), end_date = col_date(format = "%m/%d/%Y"),
end_time = col_character()))
pollen$colony <- as.factor(pollen$colony)
pollen$pid <- as.factor(pollen$pid)
pollen$count <- as.factor(pollen$count)
pollen <- subset(pollen, pollen$round != 1)
pollen## # A tibble: 933 × 15
## round colony treatment replicate count start_date start_…¹ start…² end_date
## <fct> <fct> <fct> <fct> <fct> <date> <chr> <dbl> <date>
## 1 2 1.11R2 1 11 2 2022-08-22 10:30 1.18 2022-08-24
## 2 2 1.11R2 1 11 3 2022-08-24 12:00 1.15 2022-08-26
## 3 2 1.11R2 1 11 4 2022-08-26 10:30 1.09 2022-08-28
## 4 2 1.11R2 1 11 5 2022-08-28 1:00 1.26 2022-08-30
## 5 2 1.11R2 1 11 5 2022-08-30 11:30 1.14 2022-09-01
## 6 2 1.11R2 1 11 21 2022-09-01 10:30 1.07 2022-09-03
## 7 2 1.11R2 1 11 8 2022-09-03 12:20 0.844 2022-09-05
## 8 2 1.11R2 1 11 9 2022-09-05 12:40 1.30 2022-09-07
## 9 2 1.11R2 1 11 10 2022-09-07 10:30 1.22 2022-09-09
## 10 2 1.11R2 1 11 11 2022-09-09 11:30 1.26 2022-09-11
## # … with 923 more rows, 6 more variables: end_time <chr>, end_weight <dbl>,
## # whole_dif <chr>, difference <dbl>, pid <fct>, bees_alive <dbl>, and
## # abbreviated variable names ¹start_time, ²start_weightpollen <- na.omit(pollen)
range(pollen$difference)## [1] -0.30646 1.56542# get rid of negative numbers
pollen$difference[pollen$difference < 0] <- NA
pollen <- na.omit(pollen)
range(pollen$difference)## [1] 0.002715 1.5654202 Average pollen consumption per colony
pollen$whole_dif <- as.numeric(pollen$difference)
wd.1 <- lm(difference ~ treatment, data = pollen)
summary(wd.1)##
## Call:
## lm(formula = difference ~ treatment, data = pollen)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.4853 -0.2416 -0.1504 0.1576 1.0638
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.429541 0.024170 17.772 <2e-16 ***
## treatment2 0.072099 0.034886 2.067 0.0390 *
## treatment3 0.078078 0.034406 2.269 0.0235 *
## treatment4 0.058490 0.034988 1.672 0.0949 .
## treatment5 0.004982 0.035040 0.142 0.8870
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.3375 on 915 degrees of freedom
## Multiple R-squared: 0.009927, Adjusted R-squared: 0.005599
## F-statistic: 2.294 on 4 and 915 DF, p-value: 0.05773wd.emm <- emmeans(wd.1, "treatment")
summary(wd.emm)## treatment emmean SE df lower.CL upper.CL
## 1 0.430 0.0242 915 0.382 0.477
## 2 0.502 0.0252 915 0.452 0.551
## 3 0.508 0.0245 915 0.460 0.556
## 4 0.488 0.0253 915 0.438 0.538
## 5 0.435 0.0254 915 0.385 0.484
##
## Confidence level used: 0.95wd.summary <- pollen %>%
group_by(colony) %>%
summarize(whole.mean = mean(difference))
as.data.frame(wd.summary) # this is the data frame I will merge with subsequent data frames to contain average pollen consumption per colony as a new column ## colony whole.mean
## 1 1.11R2 0.2829509
## 2 1.12R2 0.1697964
## 3 1.1R2 0.5213458
## 4 1.2R2 0.3374200
## 5 1.3R2 0.4512891
## 6 1.4R2 0.6063016
## 7 1.5R2 0.7079516
## 8 1.7R2 0.7400381
## 9 1.9R2 0.2240081
## 10 2.11R2 0.4178270
## 11 2.12R2 0.4035568
## 12 2.1R2 0.6101589
## 13 2.2R2 0.5109234
## 14 2.3R2 0.3280036
## 15 2.4R2 0.3881394
## 16 2.5R2 0.7194915
## 17 2.7R2 0.5299685
## 18 2.9R2 0.5857152
## 19 3.11R2 0.4216410
## 20 3.12R2 0.3390993
## 21 3.1R2 0.3711948
## 22 3.2R2 0.3461010
## 23 3.3R2 0.8465806
## 24 3.4R2 0.4120433
## 25 3.5R2 0.8906211
## 26 3.7R2 0.6266680
## 27 3.9R2 0.4409511
## 28 4.11R2 0.3456011
## 29 4.12R2 0.2496295
## 30 4.1R2 0.7074755
## 31 4.2R2 0.3871742
## 32 4.3R2 0.5800074
## 33 4.4R2 0.8356247
## 34 4.5R2 0.2335967
## 35 4.7R2 0.6237400
## 36 4.9R2 0.5322950
## 37 5.11R2 0.4113960
## 38 5.12R2 0.2077741
## 39 5.1R2 0.3782286
## 40 5.2R2 0.4912468
## 41 5.3R2 0.2128778
## 42 5.4R2 0.4563045
## 43 5.5R2 0.7095479
## 44 5.7R2 0.6113189
## 45 5.9R2 0.4188290brood <- read_csv("brood.csv", col_types = cols(round = col_factor(levels = c("1",
"2")), treatment = col_factor(levels = c("1",
"2", "3", "4", "5")), replicate = col_factor(levels = c("1",
"2", "3", "4", "5", "7", "9", "11", "12"))))
brood$colony <- as.factor(brood$colony)
brood <- subset(brood, brood$round != 1)
brood <- merge(wd.summary, brood, by.x = "colony") # this is good because I calculated average pollen consumption two different ways (avg_pollen was calculated a couple months ago manually) and it's the same numbers as the group_by function that created whole.mean
mortality <- read_csv("surviving workers per colony.csv")
mortality$colony <- as.factor(mortality$colony)
brood <- merge(brood, mortality, by.x = "colony")3 Brood cells
plot(brood$treatment, brood$brood_cells)
range(brood$brood_cells)## [1] 0 96Total brood count not affected by treatment, Anova = 0.343328
b.gn <- glm(brood_cells ~ treatment + whole.mean + alive , data = brood, family = "poisson")
summary(b.gn) # We do have overdispersion here now that we took out round 1##
## Call:
## glm(formula = brood_cells ~ treatment + whole.mean + alive, family = "poisson",
## data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -4.7900 -0.9928 -0.0846 0.9184 3.2283
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 2.058742 0.119312 17.255 < 2e-16 ***
## treatment2 -0.005613 0.080289 -0.070 0.944267
## treatment3 0.045938 0.077836 0.590 0.555062
## treatment4 -0.086204 0.080672 -1.069 0.285260
## treatment5 -0.118587 0.086335 -1.374 0.169577
## whole.mean 2.347534 0.140916 16.659 < 2e-16 ***
## alive 0.078524 0.023289 3.372 0.000747 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 526.03 on 44 degrees of freedom
## Residual deviance: 142.41 on 38 degrees of freedom
## AIC: 388.25
##
## Number of Fisher Scoring iterations: 5Anova(b.gn)## Analysis of Deviance Table (Type II tests)
##
## Response: brood_cells
## LR Chisq Df Pr(>Chisq)
## treatment 5.817 4 0.2132312
## whole.mean 275.581 1 < 2.2e-16 ***
## alive 12.035 1 0.0005221 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1b.gnb <- glm.nb(brood_cells ~ treatment + whole.mean + alive , data = brood)
AIC(b.gn, b.gnb)## df AIC
## b.gn 7 388.2476
## b.gnb 8 359.7677Anova(b.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: brood_cells
## LR Chisq Df Pr(>Chisq)
## treatment 2.090 4 0.719301
## whole.mean 93.057 1 < 2.2e-16 ***
## alive 7.697 1 0.005531 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(b.gn)
plot(b.gnb)
brood.emm <- emmeans(b.gnb, "treatment")
pairs(brood.emm)## contrast estimate SE df z.ratio p.value
## treatment1 - treatment2 -0.0415 0.149 Inf -0.279 0.9987
## treatment1 - treatment3 -0.0677 0.150 Inf -0.451 0.9915
## treatment1 - treatment4 0.0454 0.150 Inf 0.302 0.9982
## treatment1 - treatment5 0.1241 0.155 Inf 0.799 0.9310
## treatment2 - treatment3 -0.0262 0.143 Inf -0.183 0.9997
## treatment2 - treatment4 0.0868 0.145 Inf 0.597 0.9756
## treatment2 - treatment5 0.1655 0.149 Inf 1.112 0.8002
## treatment3 - treatment4 0.1131 0.145 Inf 0.780 0.9366
## treatment3 - treatment5 0.1918 0.149 Inf 1.291 0.6967
## treatment4 - treatment5 0.0787 0.152 Inf 0.518 0.9856
##
## Results are given on the log (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimatessummary(brood.emm)## treatment emmean SE df asymp.LCL asymp.UCL
## 1 3.47 0.108 Inf 3.26 3.68
## 2 3.51 0.102 Inf 3.31 3.71
## 3 3.54 0.103 Inf 3.34 3.74
## 4 3.43 0.105 Inf 3.22 3.63
## 5 3.35 0.109 Inf 3.13 3.56
##
## Results are given on the log (not the response) scale.
## Confidence level used: 0.95sum_brood <- brood %>%
group_by(treatment) %>%
summarise( mean.b = mean(brood_cells),
n.b = length(brood_cells),
sd.b = sd(brood_cells)) %>%
mutate(seb = sd.b / sqrt(n.b))
sum_brood## # A tibble: 5 × 5
## treatment mean.b n.b sd.b seb
## <fct> <dbl> <int> <dbl> <dbl>
## 1 1 33.8 9 22.6 7.53
## 2 2 36.9 9 11.2 3.74
## 3 3 45.6 9 26.2 8.73
## 4 4 36.7 9 18.3 6.09
## 5 5 29.6 9 17.5 5.82ggplot(data = sum_brood, aes(x = treatment, y = mean.b, fill = treatment)) +
geom_col(position = "dodge", col = "black")+
coord_cartesian(ylim=c(20,60)) +
scale_fill_manual(values = c("peachpuff3", "darkseagreen", "lightseagreen", "darkolivegreen3", "darkolivegreen4"),
name = "Pristine Level",
labels = c("Treatment 1 (control)", "Treatment 2",
"Treatment 3", "Treatment 4", "Treatment 5")) +
geom_errorbar(aes(ymin = mean.b - seb,
ymax = mean.b + seb),
position = position_dodge(0.9), width = 0.4) +
labs(y = "Mean Brood Cellls",) +
ggtitle("Average of All Brood Cells per Treatment") +
scale_x_discrete(name = "Treatment",
labels = c("0 PPB", "150 PPB", "1,500 PPB", "15,000 PPB", "150,000 PPB")) +
theme_cowplot() +
theme_classic(base_size = 12)
3.0.1 honey pots
Treatment does not* (Anova sig but post hoc not) impact count of honey pots
Biggest difference is between treatments two and one
hp.gn <- glm(honey_pot ~ treatment + whole.mean + alive , data = brood, family = "poisson")
summary(hp.gn) # overdispersion not horrible but not great ##
## Call:
## glm(formula = honey_pot ~ treatment + whole.mean + alive, family = "poisson",
## data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.5571 -1.3070 0.1144 0.7176 2.7640
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.27439 0.30038 0.913 0.361000
## treatment2 0.49968 0.20306 2.461 0.013864 *
## treatment3 0.05395 0.21908 0.246 0.805488
## treatment4 0.36456 0.20731 1.759 0.078657 .
## treatment5 0.31604 0.21508 1.469 0.141731
## whole.mean 1.34884 0.36079 3.739 0.000185 ***
## alive 0.13772 0.05760 2.391 0.016810 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 112.51 on 44 degrees of freedom
## Residual deviance: 71.33 on 38 degrees of freedom
## AIC: 237.54
##
## Number of Fisher Scoring iterations: 5hp.gnb <- glm.nb(honey_pot ~ treatment + whole.mean + alive , data = brood)
summary(hp.gnb)##
## Call:
## glm.nb(formula = honey_pot ~ treatment + whole.mean + alive,
## data = brood, init.theta = 16.77011274, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.27540 -1.10841 0.07595 0.63373 2.40056
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.25184 0.33570 0.750 0.4531
## treatment2 0.48969 0.23602 2.075 0.0380 *
## treatment3 0.04259 0.25111 0.170 0.8653
## treatment4 0.36698 0.23976 1.531 0.1259
## treatment5 0.28991 0.24831 1.168 0.2430
## whole.mean 1.36528 0.42739 3.194 0.0014 **
## alive 0.14311 0.06465 2.213 0.0269 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(16.7701) family taken to be 1)
##
## Null deviance: 86.855 on 44 degrees of freedom
## Residual deviance: 55.915 on 38 degrees of freedom
## AIC: 237.78
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 16.8
## Std. Err.: 15.5
##
## 2 x log-likelihood: -221.782AIC(hp.gn, hp.gnb)## df AIC
## hp.gn 7 237.5396
## hp.gnb 8 237.7822plot(hp.gn)
plot(hp.gnb) # stick with glm.nb
Anova(hp.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: honey_pot
## LR Chisq Df Pr(>Chisq)
## treatment 6.8343 4 0.144906
## whole.mean 10.0691 1 0.001508 **
## alive 5.1175 1 0.023686 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Anova(hp.gn)## Analysis of Deviance Table (Type II tests)
##
## Response: honey_pot
## LR Chisq Df Pr(>Chisq)
## treatment 9.6044 4 0.0476460 *
## whole.mean 13.7448 1 0.0002094 ***
## alive 6.3582 1 0.0116841 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$honey_pot)
hp.emm <- emmeans(hp.gnb, "treatment")
summary(hp.emm)## treatment emmean SE df asymp.LCL asymp.UCL
## 1 1.50 0.184 Inf 1.14 1.86
## 2 1.99 0.147 Inf 1.70 2.28
## 3 1.54 0.170 Inf 1.21 1.87
## 4 1.87 0.155 Inf 1.56 2.17
## 5 1.79 0.163 Inf 1.47 2.11
##
## Results are given on the log (not the response) scale.
## Confidence level used: 0.95pairs(hp.emm)## contrast estimate SE df z.ratio p.value
## treatment1 - treatment2 -0.4897 0.236 Inf -2.075 0.2310
## treatment1 - treatment3 -0.0426 0.251 Inf -0.170 0.9998
## treatment1 - treatment4 -0.3670 0.240 Inf -1.531 0.5424
## treatment1 - treatment5 -0.2899 0.248 Inf -1.168 0.7700
## treatment2 - treatment3 0.4471 0.220 Inf 2.030 0.2517
## treatment2 - treatment4 0.1227 0.211 Inf 0.582 0.9777
## treatment2 - treatment5 0.1998 0.217 Inf 0.922 0.8885
## treatment3 - treatment4 -0.3244 0.225 Inf -1.443 0.6000
## treatment3 - treatment5 -0.2473 0.232 Inf -1.065 0.8245
## treatment4 - treatment5 0.0771 0.225 Inf 0.343 0.9970
##
## Results are given on the log (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimatessum_hp <- brood %>%
group_by(treatment) %>%
summarise( mean.h = mean(honey_pot),
n.h = length(honey_pot),
sd.h = sd(honey_pot)) %>%
mutate(sehp = sd.h / sqrt(n.h))
sum_hp## # A tibble: 5 × 5
## treatment mean.h n.h sd.h sehp
## <fct> <dbl> <int> <dbl> <dbl>
## 1 1 4.22 9 3.27 1.09
## 2 2 7.89 9 3.59 1.20
## 3 3 5.56 9 2.96 0.988
## 4 4 7 9 3.61 1.20
## 5 5 6.22 9 4.35 1.45ggplot(data = sum_hp, aes(x = treatment, y = mean.h, fill = treatment)) +
geom_col(position = "dodge", col = "black")+
coord_cartesian(ylim=c(2,10)) +
scale_fill_manual(values = c("peachpuff3", "darkseagreen", "lightseagreen", "darkolivegreen3", "darkolivegreen4"),
name = "Pristine Level",
labels = c("Treatment 1 (control)", "Treatment 2",
"Treatment 3", "Treatment 4", "Treatment 5")) +
geom_errorbar(aes(ymin = mean.h - sehp,
ymax = mean.h + sehp),
position = position_dodge(0.9), width = 0.4) +
labs(y = "Mean Honey Pots",) +
ggtitle("Average of Honey Pots per Treatment") +
scale_x_discrete(name = "Treatment",
labels = c("0 PPB", "150 PPB", "1,500 PPB", "15,000 PPB", "150,000 PPB")) +
theme_cowplot() +
theme_classic(base_size = 12)
3.0.2 Larvae
3.0.2.1 Live Larvae
ll.gn <- glm(live_larvae ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(ll.gn)##
## Call:
## glm(formula = live_larvae ~ treatment + whole.mean + alive +
## qro, family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -5.1687 -2.9099 -0.2689 1.6383 4.9706
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.398776 0.165575 8.448 < 2e-16 ***
## treatment2 -0.686302 0.112454 -6.103 1.04e-09 ***
## treatment3 -0.013370 0.094056 -0.142 0.886964
## treatment4 -0.689669 0.110931 -6.217 5.06e-10 ***
## treatment5 -0.369537 0.114038 -3.240 0.001193 **
## whole.mean 3.463393 0.226518 15.290 < 2e-16 ***
## alive 0.006448 0.033078 0.195 0.845444
## qroB3 -0.437644 0.130090 -3.364 0.000768 ***
## qroB4 -0.129719 0.116249 -1.116 0.264475
## qroB5 0.509452 0.090192 5.649 1.62e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 807.57 on 44 degrees of freedom
## Residual deviance: 321.59 on 35 degrees of freedom
## AIC: 519.61
##
## Number of Fisher Scoring iterations: 6ll.gnb <- glm.nb(live_larvae ~ treatment + whole.mean + alive + qro, data = brood)
summary(ll.gnb)##
## Call:
## glm.nb(formula = live_larvae ~ treatment + whole.mean + alive +
## qro, data = brood, init.theta = 1.837659942, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.5314 -1.2212 -0.1145 0.3618 1.8550
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.40055 0.55254 -0.725 0.4685
## treatment2 -0.60123 0.38524 -1.561 0.1186
## treatment3 -0.08114 0.38668 -0.210 0.8338
## treatment4 -0.38510 0.39074 -0.986 0.3244
## treatment5 -0.25743 0.40396 -0.637 0.5240
## whole.mean 5.26383 0.84791 6.208 5.37e-10 ***
## alive 0.17130 0.11553 1.483 0.1381
## qroB3 -0.68875 0.43232 -1.593 0.1111
## qroB4 -0.43316 0.47293 -0.916 0.3597
## qroB5 0.98401 0.32670 3.012 0.0026 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(1.8377) family taken to be 1)
##
## Null deviance: 111.620 on 44 degrees of freedom
## Residual deviance: 56.208 on 35 degrees of freedom
## AIC: 350.84
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 1.838
## Std. Err.: 0.505
##
## 2 x log-likelihood: -328.841AIC(ll.gn, ll.gnb) # poisson overdispersed, nb AIC better ## df AIC
## ll.gn 10 519.6073
## ll.gnb 11 350.8407plot(ll.gnb)
Anova(ll.gn)## Analysis of Deviance Table (Type II tests)
##
## Response: live_larvae
## LR Chisq Df Pr(>Chisq)
## treatment 82.652 4 < 2.2e-16 ***
## whole.mean 249.094 1 < 2.2e-16 ***
## alive 0.038 1 0.8453
## qro 64.784 3 5.579e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Anova(ll.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: live_larvae
## LR Chisq Df Pr(>Chisq)
## treatment 3.1025 4 0.540814
## whole.mean 31.2720 1 2.243e-08 ***
## alive 2.3196 1 0.127756
## qro 14.8069 3 0.001989 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$live_larvae)
ll.sum <- brood %>%
group_by(treatment) %>%
summarise(ll.m= mean(live_larvae),
ll.n = length(live_larvae),
sd.ll = sd(live_larvae))
ll.sum## # A tibble: 5 × 4
## treatment ll.m ll.n sd.ll
## <fct> <dbl> <int> <dbl>
## 1 1 26.7 9 26.6
## 2 2 13.8 9 9.72
## 3 3 31.2 9 23.1
## 4 4 17.6 9 14.3
## 5 5 16.4 9 13.43.0.2.2 Dead Larvae
dl.gn <- glm(dead_larvae ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(dl.gn)##
## Call:
## glm(formula = dead_larvae ~ treatment + whole.mean + alive +
## qro, family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.7801 -1.3938 -0.4008 0.4811 4.9399
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.3105 0.5940 -5.573 2.50e-08 ***
## treatment2 0.4410 0.3131 1.408 0.159014
## treatment3 0.6943 0.2912 2.384 0.017117 *
## treatment4 0.7590 0.3133 2.423 0.015407 *
## treatment5 -0.4820 0.3839 -1.255 0.209383
## whole.mean 3.1090 0.5636 5.516 3.47e-08 ***
## alive 0.4371 0.1271 3.438 0.000586 ***
## qroB3 0.1172 0.2904 0.403 0.686635
## qroB4 0.9654 0.3199 3.018 0.002546 **
## qroB5 1.3210 0.2130 6.201 5.61e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 359.50 on 44 degrees of freedom
## Residual deviance: 152.68 on 35 degrees of freedom
## AIC: 265.66
##
## Number of Fisher Scoring iterations: 6dl.gnb <- glm.nb(dead_larvae ~ treatment + whole.mean + alive + qro, data = brood)
summary(dl.gnb)##
## Call:
## glm.nb(formula = dead_larvae ~ treatment + whole.mean + alive +
## qro, data = brood, init.theta = 0.9784953413, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.2278 -1.1300 -0.3852 0.4047 1.7093
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.23848 1.00310 -2.232 0.0256 *
## treatment2 0.19356 0.60544 0.320 0.7492
## treatment3 0.82111 0.58943 1.393 0.1636
## treatment4 0.57303 0.59637 0.961 0.3366
## treatment5 -0.53367 0.67023 -0.796 0.4259
## whole.mean 1.31984 1.25295 1.053 0.2922
## alive 0.44211 0.21018 2.104 0.0354 *
## qroB3 -0.07176 0.65004 -0.110 0.9121
## qroB4 1.23503 0.69660 1.773 0.0762 .
## qroB5 1.18181 0.50830 2.325 0.0201 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(0.9785) family taken to be 1)
##
## Null deviance: 78.290 on 44 degrees of freedom
## Residual deviance: 47.138 on 35 degrees of freedom
## AIC: 211.39
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 0.978
## Std. Err.: 0.318
## Warning while fitting theta: alternation limit reached
##
## 2 x log-likelihood: -189.389AIC(dl.gn, dl.gnb) #poisson overdispersed and nb AIC better ## df AIC
## dl.gn 10 265.6575
## dl.gnb 11 211.3894plot(dl.gnb)
Anova(dl.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: dead_larvae
## LR Chisq Df Pr(>Chisq)
## treatment 5.3531 4 0.25295
## whole.mean 1.1475 1 0.28407
## alive 4.2331 1 0.03964 *
## qro 6.0359 3 0.10988
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$dead_larvae)
dl.sum <- brood %>%
group_by(treatment) %>%
summarise(dl.m= mean(dead_larvae),
dl.n = length(dead_larvae),
sd.dl = sd(dead_larvae))
dl.sum## # A tibble: 5 × 4
## treatment dl.m dl.n sd.dl
## <fct> <dbl> <int> <dbl>
## 1 1 1.78 9 1.99
## 2 2 3.22 9 5.26
## 3 3 5.78 9 6.22
## 4 4 6.67 9 14.8
## 5 5 1.56 9 1.883.0.2.3 All larvae
brood$all_larvae <- brood$dead_larvae + brood$live_larvae
al.gn <- glm(all_larvae ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(al.gn)##
## Call:
## glm(formula = all_larvae ~ treatment + whole.mean + alive + qro,
## family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -4.2637 -1.7945 0.0573 1.0396 4.8286
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.23488 0.15781 7.825 5.08e-15 ***
## treatment2 -0.55537 0.10388 -5.346 8.97e-08 ***
## treatment3 0.04734 0.08842 0.535 0.592377
## treatment4 -0.48166 0.10097 -4.770 1.84e-06 ***
## treatment5 -0.39058 0.10910 -3.580 0.000344 ***
## whole.mean 3.49686 0.20905 16.727 < 2e-16 ***
## alive 0.04738 0.03182 1.489 0.136471
## qroB3 -0.36443 0.11863 -3.072 0.002127 **
## qroB4 -0.02008 0.10792 -0.186 0.852394
## qroB5 0.62986 0.08258 7.627 2.40e-14 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 837.19 on 44 degrees of freedom
## Residual deviance: 217.65 on 35 degrees of freedom
## AIC: 437.64
##
## Number of Fisher Scoring iterations: 5al.gnb <- glm.nb(all_larvae ~ treatment + whole.mean + alive + qro, data = brood)
summary(al.gnb)##
## Call:
## glm.nb(formula = all_larvae ~ treatment + whole.mean + alive +
## qro, data = brood, init.theta = 4.160257125, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.3826 -0.9829 -0.1342 0.4692 2.1456
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.28749 0.38933 0.738 0.460266
## treatment2 -0.49318 0.26786 -1.841 0.065589 .
## treatment3 0.04778 0.26623 0.179 0.857575
## treatment4 -0.32787 0.27042 -1.212 0.225347
## treatment5 -0.30842 0.28191 -1.094 0.273936
## whole.mean 4.22456 0.58319 7.244 4.36e-13 ***
## alive 0.15878 0.08140 1.951 0.051095 .
## qroB3 -0.57345 0.30310 -1.892 0.058500 .
## qroB4 -0.10539 0.32233 -0.327 0.743709
## qroB5 0.83665 0.22594 3.703 0.000213 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(4.1603) family taken to be 1)
##
## Null deviance: 156.752 on 44 degrees of freedom
## Residual deviance: 52.139 on 35 degrees of freedom
## AIC: 349.22
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 4.16
## Std. Err.: 1.19
##
## 2 x log-likelihood: -327.217AIC(al.gn, al.gnb)## df AIC
## al.gn 10 437.6448
## al.gnb 11 349.2167Anova(al.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: all_larvae
## LR Chisq Df Pr(>Chisq)
## treatment 6.261 4 0.18047
## whole.mean 48.846 1 2.768e-12 ***
## alive 4.051 1 0.04415 *
## qro 22.086 3 6.259e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$all_larvae)
al.sum <- brood %>%
group_by(treatment) %>%
summarise(al.m= mean(all_larvae),
al.n = length(all_larvae),
sd.al = sd(all_larvae))
al.sum## # A tibble: 5 × 4
## treatment al.m al.n sd.al
## <fct> <dbl> <int> <dbl>
## 1 1 28.4 9 27.0
## 2 2 17 9 12.5
## 3 3 37 9 26.3
## 4 4 24.2 9 23.8
## 5 5 18 9 13.63.0.2.4 Larvae Outliers
out <- boxplot.stats(brood$dead_larvae)$out
out## [1] 15 12 14 15 46out_dl <- which(brood$dead_larvae %in% c(out))
out_dl## [1] 16 20 23 25 33brood[out_dl, ]## colony whole.mean round dose treatment replicate brood_cells honey_pot eggs
## 16 2.5R2 0.7194915 2 150 2 5 44 8 14
## 20 3.12R2 0.3390993 2 1500 3 12 29 9 0
## 23 3.3R2 0.8465806 2 1500 3 3 64 9 15
## 25 3.5R2 0.8906211 2 1500 3 5 96 9 16
## 33 4.4R2 0.8356247 2 15000 4 4 53 10 13
## dead_larvae live_larvae dead_pupae live_pupae dead_drones live_drones drones
## 16 15 29 8 12 0 0 11
## 20 12 0 9 0 0 0 1
## 23 14 49 24 1 0 0 10
## 25 15 71 12 18 1 2 15
## 33 46 32 6 0 0 0 27
## avg_pollen qro alive dead all_larvae
## 16 0.7194915 B4 4 1 44
## 20 0.3390993 B1 5 0 12
## 23 0.5044652 B3 5 0 63
## 25 0.4120433 B4 4 1 86
## 33 0.5800074 B5 5 0 78brood_dl_out <- brood[-c(22, 32, 34, 39, 45), ]
dl.gn <- glm(dead_larvae ~ treatment + whole.mean + alive + qro, data = brood_dl_out, family = "poisson")
summary(dl.gn)##
## Call:
## glm(formula = dead_larvae ~ treatment + whole.mean + alive +
## qro, family = "poisson", data = brood_dl_out)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.1749 -1.7276 -0.5128 1.0317 3.8002
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -4.0970 0.7610 -5.384 7.28e-08 ***
## treatment2 0.4082 0.3128 1.305 0.191876
## treatment3 0.9862 0.2983 3.306 0.000946 ***
## treatment4 0.8172 0.3480 2.348 0.018866 *
## treatment5 -0.6800 0.3946 -1.724 0.084796 .
## whole.mean 1.6703 0.6961 2.399 0.016419 *
## alive 0.7379 0.1733 4.259 2.05e-05 ***
## qroB3 0.3575 0.3349 1.067 0.285811
## qroB4 1.5396 0.3937 3.911 9.20e-05 ***
## qroB5 1.8263 0.2671 6.837 8.09e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 332.16 on 39 degrees of freedom
## Residual deviance: 127.07 on 30 degrees of freedom
## AIC: 235.05
##
## Number of Fisher Scoring iterations: 6dl.gnb <- glm.nb(dead_larvae ~ treatment + whole.mean + alive + qro, data = brood_dl_out)
summary(dl.gnb)##
## Call:
## glm.nb(formula = dead_larvae ~ treatment + whole.mean + alive +
## qro, data = brood_dl_out, init.theta = 1.14525421, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.1120 -1.0999 -0.4191 0.4628 1.7988
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.3027 1.0897 -2.113 0.03458 *
## treatment2 0.2029 0.5760 0.352 0.72464
## treatment3 1.0035 0.5702 1.760 0.07843 .
## treatment4 0.6289 0.6247 1.007 0.31405
## treatment5 -0.3689 0.6729 -0.548 0.58355
## whole.mean 0.8030 1.3254 0.606 0.54462
## alive 0.4968 0.2247 2.211 0.02703 *
## qroB3 0.0116 0.7081 0.016 0.98693
## qroB4 1.4292 0.7445 1.920 0.05489 .
## qroB5 1.4073 0.5438 2.588 0.00967 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(1.1453) family taken to be 1)
##
## Null deviance: 75.152 on 39 degrees of freedom
## Residual deviance: 42.407 on 30 degrees of freedom
## AIC: 195.58
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 1.145
## Std. Err.: 0.401
## Warning while fitting theta: alternation limit reached
##
## 2 x log-likelihood: -173.577Anova(dl.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: dead_larvae
## LR Chisq Df Pr(>Chisq)
## treatment 5.6898 4 0.22354
## whole.mean 0.4053 1 0.52438
## alive 5.4470 1 0.01960 *
## qro 8.0591 3 0.04481 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood_dl_out$treatment, brood_dl_out$dead_larvae)
3.0.2.5 cbind Larvae
#### cbind Larvae
larvae.mortality <- glm(cbind(live_larvae, dead_larvae) ~ treatment + whole.mean + qro, data = brood, family = binomial("logit"))
summary(larvae.mortality)##
## Call:
## glm(formula = cbind(live_larvae, dead_larvae) ~ treatment + whole.mean +
## qro, family = binomial("logit"), data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -7.603 -1.646 0.000 1.690 4.054
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 3.5754 0.4250 8.412 < 2e-16 ***
## treatment2 -1.3544 0.3359 -4.032 5.52e-05 ***
## treatment3 -1.0149 0.3083 -3.292 0.000994 ***
## treatment4 -1.7935 0.3205 -5.596 2.19e-08 ***
## treatment5 -0.4505 0.3853 -1.169 0.242319
## whole.mean -0.7257 0.6306 -1.151 0.249779
## qroB3 -0.4609 0.3392 -1.359 0.174302
## qroB4 -0.5556 0.3188 -1.743 0.081415 .
## qroB5 -0.6609 0.2122 -3.114 0.001843 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 329.89 on 42 degrees of freedom
## Residual deviance: 263.65 on 34 degrees of freedom
## AIC: 352.64
##
## Number of Fisher Scoring iterations: 5Anova(larvae.mortality)## Analysis of Deviance Table (Type II tests)
##
## Response: cbind(live_larvae, dead_larvae)
## LR Chisq Df Pr(>Chisq)
## treatment 43.321 4 8.874e-09 ***
## whole.mean 1.347 1 0.24585
## qro 10.093 3 0.01779 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(larvae.mortality)
lm.emm <- emmeans(larvae.mortality, "treatment")
pairs(lm.emm)## contrast estimate SE df z.ratio p.value
## treatment1 - treatment2 1.354 0.336 Inf 4.032 0.0005
## treatment1 - treatment3 1.015 0.308 Inf 3.292 0.0088
## treatment1 - treatment4 1.793 0.320 Inf 5.596 <.0001
## treatment1 - treatment5 0.450 0.385 Inf 1.169 0.7690
## treatment2 - treatment3 -0.339 0.270 Inf -1.258 0.7170
## treatment2 - treatment4 0.439 0.295 Inf 1.488 0.5705
## treatment2 - treatment5 -0.904 0.351 Inf -2.574 0.0753
## treatment3 - treatment4 0.779 0.244 Inf 3.187 0.0125
## treatment3 - treatment5 -0.564 0.335 Inf -1.684 0.4436
## treatment4 - treatment5 -1.343 0.356 Inf -3.777 0.0015
##
## Results are averaged over the levels of: qro
## Results are given on the log odds ratio (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimatessummary(lm.emm)## treatment emmean SE df asymp.LCL asymp.UCL
## 1 2.81 0.281 Inf 2.256 3.36
## 2 1.45 0.230 Inf 1.002 1.90
## 3 1.79 0.193 Inf 1.414 2.17
## 4 1.01 0.255 Inf 0.514 1.51
## 5 2.36 0.294 Inf 1.780 2.93
##
## Results are averaged over the levels of: qro
## Results are given on the logit (not the response) scale.
## Confidence level used: 0.95plot(brood$treatment, cbind(brood$dead_larvae, brood$live_larvae))
emmeans(larvae.mortality, pairwise ~ treatment, type = "response")## $emmeans
## treatment prob SE df asymp.LCL asymp.UCL
## 1 0.943 0.0151 Inf 0.905 0.966
## 2 0.810 0.0353 Inf 0.732 0.870
## 3 0.857 0.0236 Inf 0.804 0.898
## 4 0.734 0.0498 Inf 0.626 0.820
## 5 0.913 0.0232 Inf 0.856 0.949
##
## Results are averaged over the levels of: qro
## Confidence level used: 0.95
## Intervals are back-transformed from the logit scale
##
## $contrasts
## contrast odds.ratio SE df null z.ratio p.value
## treatment1 / treatment2 3.874 1.3014 Inf 1 4.032 0.0005
## treatment1 / treatment3 2.759 0.8506 Inf 1 3.292 0.0088
## treatment1 / treatment4 6.010 1.9263 Inf 1 5.596 <.0001
## treatment1 / treatment5 1.569 0.6045 Inf 1 1.169 0.7690
## treatment2 / treatment3 0.712 0.1922 Inf 1 -1.258 0.7170
## treatment2 / treatment4 1.551 0.4579 Inf 1 1.488 0.5705
## treatment2 / treatment5 0.405 0.1422 Inf 1 -2.574 0.0753
## treatment3 / treatment4 2.178 0.5320 Inf 1 3.187 0.0125
## treatment3 / treatment5 0.569 0.1906 Inf 1 -1.684 0.4436
## treatment4 / treatment5 0.261 0.0928 Inf 1 -3.777 0.0015
##
## Results are averaged over the levels of: qro
## P value adjustment: tukey method for comparing a family of 5 estimates
## Tests are performed on the log odds ratio scaleemm1 <- emmeans(larvae.mortality, ~ treatment, type = "response")
emm1## treatment prob SE df asymp.LCL asymp.UCL
## 1 0.943 0.0151 Inf 0.905 0.966
## 2 0.810 0.0353 Inf 0.732 0.870
## 3 0.857 0.0236 Inf 0.804 0.898
## 4 0.734 0.0498 Inf 0.626 0.820
## 5 0.913 0.0232 Inf 0.856 0.949
##
## Results are averaged over the levels of: qro
## Confidence level used: 0.95
## Intervals are back-transformed from the logit scalepairs(emm1)## contrast odds.ratio SE df null z.ratio p.value
## treatment1 / treatment2 3.874 1.3014 Inf 1 4.032 0.0005
## treatment1 / treatment3 2.759 0.8506 Inf 1 3.292 0.0088
## treatment1 / treatment4 6.010 1.9263 Inf 1 5.596 <.0001
## treatment1 / treatment5 1.569 0.6045 Inf 1 1.169 0.7690
## treatment2 / treatment3 0.712 0.1922 Inf 1 -1.258 0.7170
## treatment2 / treatment4 1.551 0.4579 Inf 1 1.488 0.5705
## treatment2 / treatment5 0.405 0.1422 Inf 1 -2.574 0.0753
## treatment3 / treatment4 2.178 0.5320 Inf 1 3.187 0.0125
## treatment3 / treatment5 0.569 0.1906 Inf 1 -1.684 0.4436
## treatment4 / treatment5 0.261 0.0928 Inf 1 -3.777 0.0015
##
## Results are averaged over the levels of: qro
## P value adjustment: tukey method for comparing a family of 5 estimates
## Tests are performed on the log odds ratio scaleemmdf <- as.data.frame(emm1)
p <- ggplot(data = emmdf, aes(x=treatment, y=prob, fill=treatment)) +
geom_col(position = "dodge", color = "black") +
coord_cartesian(ylim = c(0,1))
p + geom_errorbar(position=position_dodge(.9),width=.25, aes(ymax=asymp.UCL, ymin=asymp.LCL),alpha=.6)+
labs(x="Treatment", y="Probability of Survival", fill = "treatment") + theme_bw() +
scale_fill_manual(values = c("peachpuff3", "darkseagreen", "lightseagreen", "darkolivegreen3", "darkolivegreen4")) +
labs(y = "Probability of Survival",) +
ggtitle("Probability of larvae being alive upon colony dissection") +
scale_x_discrete(name = "Treatment",
labels = c("0 PPB", "150 PPB", "1,500 PPB", "15,000 PPB", "150,000 PPB")) +
theme_cowplot() +
theme_classic(base_size = 12)
3.0.3 Pupae
3.0.3.1 Live Pupae
lp.gn <- glm(live_pupae ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(lp.gn)##
## Call:
## glm(formula = live_pupae ~ treatment + whole.mean + alive + qro,
## family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.8137 -1.7257 -0.0934 0.7192 3.1344
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.30791 0.34731 0.887 0.3753
## treatment2 0.18431 0.21445 0.859 0.3901
## treatment3 -0.22210 0.23855 -0.931 0.3518
## treatment4 -0.54854 0.25672 -2.137 0.0326 *
## treatment5 -0.21140 0.27058 -0.781 0.4346
## whole.mean 3.71756 0.56044 6.633 3.28e-11 ***
## alive -0.14146 0.07282 -1.942 0.0521 .
## qroB3 -0.80480 0.33227 -2.422 0.0154 *
## qroB4 -0.33642 0.25730 -1.308 0.1910
## qroB5 -0.17451 0.22851 -0.764 0.4450
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 212.80 on 44 degrees of freedom
## Residual deviance: 136.34 on 35 degrees of freedom
## AIC: 263.48
##
## Number of Fisher Scoring iterations: 5lp.gnb <- glm.nb(live_pupae ~ treatment + whole.mean + alive + qro, data = brood)
summary(lp.gnb)##
## Call:
## glm.nb(formula = live_pupae ~ treatment + whole.mean + alive +
## qro, data = brood, init.theta = 1.353246607, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.5996 -1.2273 -0.1603 0.4228 1.4946
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.60720 0.69502 -0.874 0.38231
## treatment2 0.09051 0.48097 0.188 0.85073
## treatment3 -0.35289 0.50348 -0.701 0.48337
## treatment4 -0.46536 0.50757 -0.917 0.35923
## treatment5 -0.38140 0.52620 -0.725 0.46856
## whole.mean 4.03180 1.09439 3.684 0.00023 ***
## alive 0.02460 0.14556 0.169 0.86578
## qroB3 -0.51626 0.55342 -0.933 0.35089
## qroB4 -0.41144 0.59601 -0.690 0.48999
## qroB5 0.29408 0.42599 0.690 0.48998
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(1.3532) family taken to be 1)
##
## Null deviance: 71.812 on 44 degrees of freedom
## Residual deviance: 51.897 on 35 degrees of freedom
## AIC: 231.55
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 1.353
## Std. Err.: 0.473
##
## 2 x log-likelihood: -209.550AIC(lp.gn, lp.gnb)## df AIC
## lp.gn 10 263.4770
## lp.gnb 11 231.5498Anova(lp.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: live_pupae
## LR Chisq Df Pr(>Chisq)
## treatment 1.8933 4 0.7553832
## whole.mean 11.5823 1 0.0006658 ***
## alive 0.0277 1 0.8678000
## qro 2.1240 3 0.5470680
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$live_pupae)
lp.sum <- brood %>%
group_by(treatment) %>%
summarise(lp.m= mean(live_pupae),
sd.lp = sd(live_pupae))
lp.sum## # A tibble: 5 × 3
## treatment lp.m sd.lp
## <fct> <dbl> <dbl>
## 1 1 4.78 4.60
## 2 2 5.56 5.00
## 3 3 4.11 5.75
## 4 4 3 2.87
## 5 5 2.89 2.983.0.3.2 Dead Pupae
dp.gn <- glm(dead_pupae ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(dp.gn)##
## Call:
## glm(formula = dead_pupae ~ treatment + whole.mean + alive + qro,
## family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.1879 -1.4754 -0.6729 0.6997 6.6750
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.9980 0.6050 -3.302 0.000959 ***
## treatment2 1.2132 0.2670 4.544 5.53e-06 ***
## treatment3 1.1173 0.2664 4.194 2.74e-05 ***
## treatment4 0.0369 0.3145 0.117 0.906577
## treatment5 0.6690 0.3027 2.210 0.027108 *
## whole.mean 3.5946 0.5173 6.949 3.68e-12 ***
## alive 0.2490 0.1202 2.072 0.038227 *
## qroB3 -0.3758 0.2346 -1.602 0.109214
## qroB4 -0.5532 0.2704 -2.046 0.040784 *
## qroB5 -0.4223 0.2411 -1.751 0.079872 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 309.89 on 44 degrees of freedom
## Residual deviance: 144.89 on 35 degrees of freedom
## AIC: 283.4
##
## Number of Fisher Scoring iterations: 6dp.gnb <- glm.nb(dead_pupae ~ treatment + whole.mean + alive + qro, data = brood) # stick with this model overall
summary(dp.gnb)##
## Call:
## glm.nb(formula = dead_pupae ~ treatment + whole.mean + alive +
## qro, data = brood, init.theta = 2.166705709, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.2273 -1.0564 -0.3197 0.5204 2.5944
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.67904 0.76359 -2.199 0.027886 *
## treatment2 0.97034 0.44264 2.192 0.028365 *
## treatment3 1.08226 0.44245 2.446 0.014441 *
## treatment4 -0.04256 0.48451 -0.088 0.930001
## treatment5 0.47513 0.47846 0.993 0.320696
## whole.mean 3.31190 0.95322 3.474 0.000512 ***
## alive 0.23403 0.15536 1.506 0.131971
## qroB3 -0.55456 0.47270 -1.173 0.240726
## qroB4 -0.27993 0.52078 -0.538 0.590901
## qroB5 -0.38462 0.39122 -0.983 0.325538
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(2.1667) family taken to be 1)
##
## Null deviance: 99.941 on 44 degrees of freedom
## Residual deviance: 49.966 on 35 degrees of freedom
## AIC: 233.47
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 2.167
## Std. Err.: 0.764
##
## 2 x log-likelihood: -211.466AIC(dp.gn, dp.gnb)## df AIC
## dp.gn 10 283.3951
## dp.gnb 11 233.4660Anova(dp.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: dead_pupae
## LR Chisq Df Pr(>Chisq)
## treatment 12.0698 4 0.0168400 *
## whole.mean 12.5049 1 0.0004059 ***
## alive 2.1621 1 0.1414558
## qro 1.7479 3 0.6263385
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$dead_pupae)
em.dp <- emmeans(dp.gnb, "treatment")
summary(em.dp)## treatment emmean SE df asymp.LCL asymp.UCL
## 1 0.575 0.371 Inf -0.151 1.30
## 2 1.545 0.299 Inf 0.959 2.13
## 3 1.657 0.291 Inf 1.086 2.23
## 4 0.532 0.377 Inf -0.206 1.27
## 5 1.050 0.330 Inf 0.403 1.70
##
## Results are averaged over the levels of: qro
## Results are given on the log (not the response) scale.
## Confidence level used: 0.95pairs(em.dp)## contrast estimate SE df z.ratio p.value
## treatment1 - treatment2 -0.9703 0.443 Inf -2.192 0.1825
## treatment1 - treatment3 -1.0823 0.442 Inf -2.446 0.1033
## treatment1 - treatment4 0.0426 0.485 Inf 0.088 1.0000
## treatment1 - treatment5 -0.4751 0.478 Inf -0.993 0.8586
## treatment2 - treatment3 -0.1119 0.376 Inf -0.298 0.9983
## treatment2 - treatment4 1.0129 0.431 Inf 2.353 0.1285
## treatment2 - treatment5 0.4952 0.411 Inf 1.205 0.7486
## treatment3 - treatment4 1.1248 0.424 Inf 2.652 0.0614
## treatment3 - treatment5 0.6071 0.405 Inf 1.499 0.5629
## treatment4 - treatment5 -0.5177 0.469 Inf -1.104 0.8045
##
## Results are averaged over the levels of: qro
## Results are given on the log (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimatesdp.gnb1 <- glm.nb(dead_pupae ~ treatment + whole.mean + qro, data = brood)
dp.gnb2 <- glm.nb(dead_pupae ~ treatment + whole.mean + qro + alive, data = brood)
dp.gnb3 <- glm.nb(dead_pupae ~ treatment + whole.mean + qro , data = brood)
AIC(dp.gnb, dp.gnb1, dp.gnb2, dp.gnb3)## df AIC
## dp.gnb 11 233.4660
## dp.gnb1 10 233.6074
## dp.gnb2 11 233.4660
## dp.gnb3 10 233.6074plot(dp.gnb1)
plot(dp.gnb)
plot(dp.gnb3)
em.dp3 <- emmeans(dp.gnb, "treatment")
summary(em.dp3)## treatment emmean SE df asymp.LCL asymp.UCL
## 1 0.575 0.371 Inf -0.151 1.30
## 2 1.545 0.299 Inf 0.959 2.13
## 3 1.657 0.291 Inf 1.086 2.23
## 4 0.532 0.377 Inf -0.206 1.27
## 5 1.050 0.330 Inf 0.403 1.70
##
## Results are averaged over the levels of: qro
## Results are given on the log (not the response) scale.
## Confidence level used: 0.95pairs(em.dp3)## contrast estimate SE df z.ratio p.value
## treatment1 - treatment2 -0.9703 0.443 Inf -2.192 0.1825
## treatment1 - treatment3 -1.0823 0.442 Inf -2.446 0.1033
## treatment1 - treatment4 0.0426 0.485 Inf 0.088 1.0000
## treatment1 - treatment5 -0.4751 0.478 Inf -0.993 0.8586
## treatment2 - treatment3 -0.1119 0.376 Inf -0.298 0.9983
## treatment2 - treatment4 1.0129 0.431 Inf 2.353 0.1285
## treatment2 - treatment5 0.4952 0.411 Inf 1.205 0.7486
## treatment3 - treatment4 1.1248 0.424 Inf 2.652 0.0614
## treatment3 - treatment5 0.6071 0.405 Inf 1.499 0.5629
## treatment4 - treatment5 -0.5177 0.469 Inf -1.104 0.8045
##
## Results are averaged over the levels of: qro
## Results are given on the log (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimatesdp.sum <- brood %>%
group_by(treatment) %>%
summarise(dp.m= mean(dead_pupae),
sd.dp = sd(dead_pupae),
n.dp = length(dead_pupae)) %>%
mutate(sedp = sd.dp/ sqrt(n.dp))
dp.sum## # A tibble: 5 × 5
## treatment dp.m sd.dp n.dp sedp
## <fct> <dbl> <dbl> <int> <dbl>
## 1 1 2 2.06 9 0.687
## 2 2 7.89 12.1 9 4.04
## 3 3 8.89 7.27 9 2.42
## 4 4 2.89 3.02 9 1.01
## 5 5 3.89 4.28 9 1.43ggplot(data = dp.sum, aes(x = treatment, y = dp.m, fill = treatment)) +
geom_col(position = "dodge", col = "black")+
coord_cartesian(ylim=c(0,14)) +
scale_fill_manual(values = c("peachpuff3", "darkseagreen", "lightseagreen", "darkolivegreen3", "darkolivegreen4"),
name = "Pristine Level",
labels = c("Treatment 1 (control)", "Treatment 2",
"Treatment 3", "Treatment 4", "Treatment 5")) +
geom_errorbar(aes(ymin = dp.m - sedp,
ymax = dp.m + sedp),
position = position_dodge(0.9), width = 0.4) +
labs(y = "Mean Dead Pupae",) +
ggtitle("Average of Dead Pupae per Treatment") +
scale_x_discrete(name = "Treatment",
labels = c("0 PPB", "150 PPB", "1,500 PPB", "15,000 PPB", "150,000 PPB")) +
theme_cowplot() +
theme_classic(base_size = 12)
3.0.3.3 All Pupae
brood$all_pupae <- brood$dead_pupae + brood$live_pupae
ap.gn <- glm(all_pupae ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(ap.gn)##
## Call:
## glm(formula = all_pupae ~ treatment + whole.mean + alive + qro,
## family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.8914 -1.7980 -0.4441 0.7966 6.0279
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.26544 0.28789 0.922 0.35652
## treatment2 0.63796 0.16090 3.965 7.34e-05 ***
## treatment3 0.43271 0.16589 2.608 0.00910 **
## treatment4 -0.35467 0.19435 -1.825 0.06802 .
## treatment5 0.16720 0.19295 0.867 0.38620
## whole.mean 3.71191 0.37817 9.816 < 2e-16 ***
## alive -0.01685 0.05845 -0.288 0.77310
## qroB3 -0.52312 0.18711 -2.796 0.00518 **
## qroB4 -0.55219 0.18409 -2.999 0.00270 **
## qroB5 -0.33766 0.16290 -2.073 0.03819 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 350.04 on 44 degrees of freedom
## Residual deviance: 156.99 on 35 degrees of freedom
## AIC: 332.87
##
## Number of Fisher Scoring iterations: 5ap.gnb <- glm.nb(all_pupae ~ treatment + whole.mean + alive + qro, data = brood)
summary(ap.gnb)##
## Call:
## glm.nb(formula = all_pupae ~ treatment + whole.mean + alive +
## qro, data = brood, init.theta = 3.191596136, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.4545 -0.8815 -0.1592 0.4591 2.2075
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.15077 0.48607 -0.310 0.7564
## treatment2 0.54703 0.32439 1.686 0.0917 .
## treatment3 0.36441 0.33290 1.095 0.2737
## treatment4 -0.33364 0.34991 -0.954 0.3403
## treatment5 0.02044 0.35340 0.058 0.9539
## whole.mean 3.77931 0.72737 5.196 2.04e-07 ***
## alive 0.06842 0.10020 0.683 0.4947
## qroB3 -0.50376 0.36253 -1.390 0.1647
## qroB4 -0.43693 0.39774 -1.099 0.2720
## qroB5 -0.07309 0.28858 -0.253 0.8001
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(3.1916) family taken to be 1)
##
## Null deviance: 104.085 on 44 degrees of freedom
## Residual deviance: 51.164 on 35 degrees of freedom
## AIC: 281.25
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 3.19
## Std. Err.: 1.03
##
## 2 x log-likelihood: -259.249AIC(ap.gn, ap.gnb)## df AIC
## ap.gn 10 332.8670
## ap.gnb 11 281.2485Anova(ap.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: all_pupae
## LR Chisq Df Pr(>Chisq)
## treatment 8.9545 4 0.06225 .
## whole.mean 27.8628 1 1.302e-07 ***
## alive 0.4477 1 0.50343
## qro 2.6980 3 0.44056
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$all_pupae)
ap.sum <- brood %>%
group_by(treatment) %>%
summarise(al.m= mean(all_pupae),
al.n = length(all_pupae),
sd.al = sd(all_pupae)) %>%
mutate(sep = sd.al/sqrt(al.n))
ap.sum## # A tibble: 5 × 5
## treatment al.m al.n sd.al sep
## <fct> <dbl> <int> <dbl> <dbl>
## 1 1 6.78 9 6.22 2.07
## 2 2 13.4 9 14.1 4.68
## 3 3 13 9 9.5 3.17
## 4 4 5.89 9 5.16 1.72
## 5 5 6.78 9 5.89 1.96ggplot(data = ap.sum, aes(x = treatment, y = al.m, fill = treatment)) +
geom_col(position = "dodge", col = "black")+
coord_cartesian(ylim=c(0,20)) +
scale_fill_manual(values = c("peachpuff3", "darkseagreen", "lightseagreen", "darkolivegreen3", "darkolivegreen4"),
name = "Pristine Level",
labels = c("Treatment 1 (control)", "Treatment 2",
"Treatment 3", "Treatment 4", "Treatment 5")) +
geom_errorbar(aes(ymin = al.m - sep,
ymax = al.m + sep),
position = position_dodge(0.9), width = 0.4) +
labs(y = "Mean of all Pupae",) +
ggtitle("Average of All (dead and alive) Pupae per Treatment") +
scale_x_discrete(name = "Treatment",
labels = c("0 PPB", "150 PPB", "1,500 PPB", "15,000 PPB", "150,000 PPB")) +
theme_cowplot() +
theme_classic(base_size = 12)
##cbind pupae
pupae.mortality <- glm(cbind(live_pupae, dead_pupae) ~ treatment + whole.mean + qro, data = brood, family = binomial("logit"))
summary(pupae.mortality)##
## Call:
## glm(formula = cbind(live_pupae, dead_pupae) ~ treatment + whole.mean +
## qro, family = binomial("logit"), data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.507 -1.388 0.374 1.075 3.360
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.1605 0.5381 2.157 0.031025 *
## treatment2 -1.3652 0.3537 -3.859 0.000114 ***
## treatment3 -1.7369 0.3605 -4.818 1.45e-06 ***
## treatment4 -0.8751 0.4075 -2.148 0.031750 *
## treatment5 -1.6576 0.4124 -4.019 5.84e-05 ***
## whole.mean -0.6907 0.8159 -0.847 0.397251
## qroB3 -0.5694 0.4330 -1.315 0.188433
## qroB4 0.9133 0.3583 2.549 0.010805 *
## qroB5 0.8723 0.3042 2.868 0.004135 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 172.65 on 40 degrees of freedom
## Residual deviance: 127.41 on 32 degrees of freedom
## AIC: 215.17
##
## Number of Fisher Scoring iterations: 4Anova(pupae.mortality)## Analysis of Deviance Table (Type II tests)
##
## Response: cbind(live_pupae, dead_pupae)
## LR Chisq Df Pr(>Chisq)
## treatment 29.8204 4 5.324e-06 ***
## whole.mean 0.7204 1 0.396014
## qro 18.8401 3 0.000295 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(pupae.mortality)
lm.emm <- emmeans(pupae.mortality, "treatment")
pairs(lm.emm)## contrast estimate SE df z.ratio p.value
## treatment1 - treatment2 1.3652 0.354 Inf 3.859 0.0011
## treatment1 - treatment3 1.7369 0.361 Inf 4.818 <.0001
## treatment1 - treatment4 0.8751 0.407 Inf 2.148 0.2001
## treatment1 - treatment5 1.6576 0.412 Inf 4.019 0.0006
## treatment2 - treatment3 0.3717 0.289 Inf 1.286 0.6999
## treatment2 - treatment4 -0.4900 0.346 Inf -1.417 0.6165
## treatment2 - treatment5 0.2924 0.340 Inf 0.860 0.9113
## treatment3 - treatment4 -0.8618 0.360 Inf -2.393 0.1172
## treatment3 - treatment5 -0.0794 0.352 Inf -0.225 0.9994
## treatment4 - treatment5 0.7824 0.408 Inf 1.916 0.3085
##
## Results are averaged over the levels of: qro
## Results are given on the log odds ratio (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimatessummary(lm.emm)## treatment emmean SE df asymp.LCL asymp.UCL
## 1 1.133 0.319 Inf 0.507 1.7587
## 2 -0.233 0.243 Inf -0.708 0.2432
## 3 -0.604 0.248 Inf -1.091 -0.1176
## 4 0.258 0.341 Inf -0.411 0.9259
## 5 -0.525 0.293 Inf -1.099 0.0495
##
## Results are averaged over the levels of: qro
## Results are given on the logit (not the response) scale.
## Confidence level used: 0.95plot(brood$treatment, cbind(brood$dead_pupae, brood$live_pupae))
emmeans(pupae.mortality, pairwise ~ treatment, type = "response")## $emmeans
## treatment prob SE df asymp.LCL asymp.UCL
## 1 0.756 0.0589 Inf 0.624 0.853
## 2 0.442 0.0599 Inf 0.330 0.561
## 3 0.353 0.0567 Inf 0.251 0.471
## 4 0.564 0.0839 Inf 0.399 0.716
## 5 0.372 0.0685 Inf 0.250 0.512
##
## Results are averaged over the levels of: qro
## Confidence level used: 0.95
## Intervals are back-transformed from the logit scale
##
## $contrasts
## contrast odds.ratio SE df null z.ratio p.value
## treatment1 / treatment2 3.916 1.385 Inf 1 3.859 0.0011
## treatment1 / treatment3 5.680 2.048 Inf 1 4.818 <.0001
## treatment1 / treatment4 2.399 0.978 Inf 1 2.148 0.2001
## treatment1 / treatment5 5.246 2.164 Inf 1 4.019 0.0006
## treatment2 / treatment3 1.450 0.419 Inf 1 1.286 0.6999
## treatment2 / treatment4 0.613 0.212 Inf 1 -1.417 0.6165
## treatment2 / treatment5 1.340 0.455 Inf 1 0.860 0.9113
## treatment3 / treatment4 0.422 0.152 Inf 1 -2.393 0.1172
## treatment3 / treatment5 0.924 0.325 Inf 1 -0.225 0.9994
## treatment4 / treatment5 2.187 0.893 Inf 1 1.916 0.3085
##
## Results are averaged over the levels of: qro
## P value adjustment: tukey method for comparing a family of 5 estimates
## Tests are performed on the log odds ratio scaleemm2 <- emmeans(pupae.mortality, ~ treatment, type = "response")
emm2## treatment prob SE df asymp.LCL asymp.UCL
## 1 0.756 0.0589 Inf 0.624 0.853
## 2 0.442 0.0599 Inf 0.330 0.561
## 3 0.353 0.0567 Inf 0.251 0.471
## 4 0.564 0.0839 Inf 0.399 0.716
## 5 0.372 0.0685 Inf 0.250 0.512
##
## Results are averaged over the levels of: qro
## Confidence level used: 0.95
## Intervals are back-transformed from the logit scalepairs(emm2)## contrast odds.ratio SE df null z.ratio p.value
## treatment1 / treatment2 3.916 1.385 Inf 1 3.859 0.0011
## treatment1 / treatment3 5.680 2.048 Inf 1 4.818 <.0001
## treatment1 / treatment4 2.399 0.978 Inf 1 2.148 0.2001
## treatment1 / treatment5 5.246 2.164 Inf 1 4.019 0.0006
## treatment2 / treatment3 1.450 0.419 Inf 1 1.286 0.6999
## treatment2 / treatment4 0.613 0.212 Inf 1 -1.417 0.6165
## treatment2 / treatment5 1.340 0.455 Inf 1 0.860 0.9113
## treatment3 / treatment4 0.422 0.152 Inf 1 -2.393 0.1172
## treatment3 / treatment5 0.924 0.325 Inf 1 -0.225 0.9994
## treatment4 / treatment5 2.187 0.893 Inf 1 1.916 0.3085
##
## Results are averaged over the levels of: qro
## P value adjustment: tukey method for comparing a family of 5 estimates
## Tests are performed on the log odds ratio scaleemmdf2 <- as.data.frame(emm2)
p <- ggplot(data = emmdf2, aes(x=treatment, y=prob, fill=treatment)) +
geom_col(position = "dodge", color = "black") +
coord_cartesian(ylim = c(0,1))
p + geom_errorbar(position=position_dodge(.9),width=.25, aes(ymax=asymp.UCL, ymin=asymp.LCL),alpha=.6)+
labs(x="Treatment", y="Probability of Survival", fill = "treatment") + theme_bw() +
scale_fill_manual(values = c("peachpuff3", "darkseagreen", "lightseagreen", "darkolivegreen3", "darkolivegreen4")) +
labs(y = "Probability of Survival",) +
ggtitle("Probability of pupae being alive upon colony dissection") +
scale_x_discrete(name = "Treatment",
labels = c("0 PPB", "150 PPB", "1,500 PPB", "15,000 PPB", "150,000 PPB")) +
theme_cowplot() +
theme_classic(base_size = 12)
3.0.4 Dead Larvae + Dead Pupae
brood$all_d_lp <- brood$dead_larvae + brood$dead_pupae
dlp.gn <- glm(all_d_lp ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(dlp.gn)##
## Call:
## glm(formula = all_d_lp ~ treatment + whole.mean + alive + qro,
## family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -3.8308 -1.6570 -0.5231 1.1627 6.2517
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.07730 0.43090 -4.821 1.43e-06 ***
## treatment2 0.91289 0.20032 4.557 5.18e-06 ***
## treatment3 0.92349 0.19610 4.709 2.49e-06 ***
## treatment4 0.43548 0.21561 2.020 0.04341 *
## treatment5 0.21182 0.23447 0.903 0.36631
## whole.mean 3.55740 0.38049 9.350 < 2e-16 ***
## alive 0.36562 0.08801 4.154 3.27e-05 ***
## qroB3 -0.24110 0.18108 -1.331 0.18304
## qroB4 0.02634 0.20285 0.130 0.89667
## qroB5 0.45668 0.14752 3.096 0.00196 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 484.94 on 44 degrees of freedom
## Residual deviance: 194.21 on 35 degrees of freedom
## AIC: 362.89
##
## Number of Fisher Scoring iterations: 5dlp.gnb <- glm.nb(all_d_lp ~ treatment + whole.mean + alive + qro, data = brood)
summary(dlp.gnb)##
## Call:
## glm.nb(formula = all_d_lp ~ treatment + whole.mean + alive +
## qro, data = brood, init.theta = 2.195213809, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.3710 -0.8591 -0.2551 0.5116 2.2069
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.35803 0.65774 -2.065 0.03895 *
## treatment2 0.65736 0.39373 1.670 0.09500 .
## treatment3 0.83071 0.39310 2.113 0.03458 *
## treatment4 0.28839 0.40559 0.711 0.47707
## treatment5 -0.05105 0.42989 -0.119 0.90547
## whole.mean 2.62933 0.83989 3.131 0.00174 **
## alive 0.35279 0.13590 2.596 0.00943 **
## qroB3 -0.43310 0.42841 -1.011 0.31205
## qroB4 0.34250 0.46690 0.734 0.46322
## qroB5 0.39828 0.33935 1.174 0.24053
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(2.1952) family taken to be 1)
##
## Null deviance: 102.58 on 44 degrees of freedom
## Residual deviance: 47.63 on 35 degrees of freedom
## AIC: 276.98
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 2.195
## Std. Err.: 0.632
##
## 2 x log-likelihood: -254.977AIC(dlp.gn, dlp.gnb)## df AIC
## dlp.gn 10 362.8904
## dlp.gnb 11 276.9769Anova(dlp.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: all_d_lp
## LR Chisq Df Pr(>Chisq)
## treatment 7.8007 4 0.099158 .
## whole.mean 10.4570 1 0.001222 **
## alive 6.2967 1 0.012096 *
## qro 3.1277 3 0.372345
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$all_d_lp)
dlp.sum <- brood %>%
group_by(treatment) %>%
summarise(al.m= mean(all_d_lp),
al.n = length(all_d_lp),
sd.al = sd(all_d_lp))
dlp.sum## # A tibble: 5 × 4
## treatment al.m al.n sd.al
## <fct> <dbl> <int> <dbl>
## 1 1 3.78 9 2.99
## 2 2 11.1 9 13.2
## 3 3 14.7 9 12.6
## 4 4 9.56 9 16.1
## 5 5 5.44 9 5.773.0.5 Live Larvae + Live Pupae
brood$all_a_lp <- brood$live_larvae + brood$live_pupae
alp.gn <- glm(all_a_lp ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(alp.gn)##
## Call:
## glm(formula = all_a_lp ~ treatment + whole.mean + alive + qro,
## family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -5.5458 -2.4285 -0.3596 1.8114 5.0777
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.66892 0.14930 11.178 < 2e-16 ***
## treatment2 -0.50456 0.09818 -5.139 2.76e-07 ***
## treatment3 -0.03876 0.08735 -0.444 0.65725
## treatment4 -0.66432 0.10162 -6.537 6.27e-11 ***
## treatment5 -0.34495 0.10504 -3.284 0.00102 **
## whole.mean 3.48866 0.20954 16.649 < 2e-16 ***
## alive -0.01605 0.02999 -0.535 0.59253
## qroB3 -0.49024 0.12070 -4.062 4.87e-05 ***
## qroB4 -0.16064 0.10583 -1.518 0.12903
## qroB5 0.40953 0.08357 4.901 9.55e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 864.12 on 44 degrees of freedom
## Residual deviance: 333.87 on 35 degrees of freedom
## AIC: 544.85
##
## Number of Fisher Scoring iterations: 6alp.gnb <- glm.nb(all_a_lp ~ treatment + whole.mean + alive + qro, data = brood)
summary(alp.gnb)##
## Call:
## glm.nb(formula = all_a_lp ~ treatment + whole.mean + alive +
## qro, data = brood, init.theta = 2.183176069, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.7635 -1.0289 -0.1945 0.3398 2.0526
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.07191 0.49797 -0.144 0.88517
## treatment2 -0.44112 0.35015 -1.260 0.20775
## treatment3 -0.12109 0.35397 -0.342 0.73228
## treatment4 -0.41114 0.35704 -1.152 0.24952
## treatment5 -0.30374 0.36875 -0.824 0.41012
## whole.mean 5.12239 0.77379 6.620 3.59e-11 ***
## alive 0.15896 0.10423 1.525 0.12726
## qroB3 -0.62251 0.39020 -1.595 0.11063
## qroB4 -0.44898 0.43180 -1.040 0.29844
## qroB5 0.92239 0.29856 3.089 0.00201 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(2.1832) family taken to be 1)
##
## Null deviance: 115.765 on 44 degrees of freedom
## Residual deviance: 55.625 on 35 degrees of freedom
## AIC: 365.41
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 2.183
## Std. Err.: 0.583
##
## 2 x log-likelihood: -343.406AIC(alp.gn, alp.gnb)## df AIC
## alp.gn 10 544.8540
## alp.gnb 11 365.4059Anova(alp.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: all_a_lp
## LR Chisq Df Pr(>Chisq)
## treatment 2.389 4 0.664636
## whole.mean 35.794 1 2.194e-09 ***
## alive 2.395 1 0.121739
## qro 15.495 3 0.001439 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$all_a_lp)
alp.sum <- brood %>%
group_by(treatment) %>%
summarise(al.m= mean(all_a_lp),
al.n = length(all_a_lp),
sd.al = sd(all_a_lp))
alp.sum## # A tibble: 5 × 4
## treatment al.m al.n sd.al
## <fct> <dbl> <int> <dbl>
## 1 1 31.4 9 29.9
## 2 2 19.3 9 13.9
## 3 3 35.3 9 26.6
## 4 4 20.6 9 16.2
## 5 5 19.3 9 14.33.0.6 All Larvae + All Pupae
brood$all_lp <- brood$all_larvae + brood$all_pupae
lp.gn <- glm(all_lp ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(lp.gn)##
## Call:
## glm(formula = all_lp ~ treatment + whole.mean + alive + qro,
## family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -4.1545 -1.9008 -0.0325 1.0822 4.4637
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.55086 0.13812 11.228 < 2e-16 ***
## treatment2 -0.19105 0.08411 -2.271 0.02312 *
## treatment3 0.13157 0.07769 1.694 0.09035 .
## treatment4 -0.45485 0.08908 -5.106 3.29e-07 ***
## treatment5 -0.26132 0.09454 -2.764 0.00571 **
## whole.mean 3.55940 0.18318 19.431 < 2e-16 ***
## alive 0.03663 0.02789 1.313 0.18919
## qroB3 -0.41650 0.10002 -4.164 3.12e-05 ***
## qroB4 -0.16855 0.09280 -1.816 0.06932 .
## qroB5 0.40371 0.07262 5.559 2.71e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 937.80 on 44 degrees of freedom
## Residual deviance: 218.37 on 35 degrees of freedom
## AIC: 456.94
##
## Number of Fisher Scoring iterations: 5lp.gnb <- glm.nb(all_lp ~ treatment + whole.mean + alive + qro, data = brood)
summary(lp.gnb)##
## Call:
## glm.nb(formula = all_lp ~ treatment + whole.mean + alive + qro,
## data = brood, init.theta = 6.311776013, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.5242 -1.0381 -0.0703 0.4058 2.2768
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.72579 0.31822 2.281 0.022560 *
## treatment2 -0.17525 0.21634 -0.810 0.417890
## treatment3 0.09406 0.21768 0.432 0.665674
## treatment4 -0.37079 0.22259 -1.666 0.095748 .
## treatment5 -0.30270 0.23087 -1.311 0.189811
## whole.mean 4.09843 0.47651 8.601 < 2e-16 ***
## alive 0.15125 0.06627 2.282 0.022461 *
## qroB3 -0.56704 0.24475 -2.317 0.020514 *
## qroB4 -0.19264 0.26377 -0.730 0.465202
## qroB5 0.63013 0.18526 3.401 0.000671 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(6.3118) family taken to be 1)
##
## Null deviance: 188.11 on 44 degrees of freedom
## Residual deviance: 52.45 on 35 degrees of freedom
## AIC: 366.68
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 6.31
## Std. Err.: 1.79
##
## 2 x log-likelihood: -344.676AIC(lp.gn, lp.gnb)## df AIC
## lp.gn 10 456.9434
## lp.gnb 11 366.6761Anova(lp.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: all_lp
## LR Chisq Df Pr(>Chisq)
## treatment 6.693 4 0.15302
## whole.mean 71.683 1 < 2.2e-16 ***
## alive 5.435 1 0.01973 *
## qro 22.357 3 5.497e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(brood$treatment, brood$all_lp)
lp.sum <- brood %>%
group_by(treatment) %>%
summarise(lp.m= mean(all_lp),
lp.n = length(all_lp),
lp.al = sd(all_lp))
lp.sum## # A tibble: 5 × 4
## treatment lp.m lp.n lp.al
## <fct> <dbl> <int> <dbl>
## 1 1 35.2 9 30.6
## 2 2 30.4 9 19.4
## 3 3 50 9 34.6
## 4 4 30.1 9 26.8
## 5 5 24.8 9 17.1em.lp <- emmeans(lp.gnb, "treatment")
pairs(em.lp)## contrast estimate SE df z.ratio p.value
## treatment1 - treatment2 0.1753 0.216 Inf 0.810 0.9276
## treatment1 - treatment3 -0.0941 0.218 Inf -0.432 0.9928
## treatment1 - treatment4 0.3708 0.223 Inf 1.666 0.4553
## treatment1 - treatment5 0.3027 0.231 Inf 1.311 0.6843
## treatment2 - treatment3 -0.2693 0.209 Inf -1.291 0.6969
## treatment2 - treatment4 0.1955 0.217 Inf 0.901 0.8967
## treatment2 - treatment5 0.1274 0.221 Inf 0.576 0.9785
## treatment3 - treatment4 0.4648 0.215 Inf 2.164 0.1935
## treatment3 - treatment5 0.3968 0.216 Inf 1.836 0.3528
## treatment4 - treatment5 -0.0681 0.231 Inf -0.295 0.9983
##
## Results are averaged over the levels of: qro
## Results are given on the log (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimates####cbind all larvae and all pupae
pl_bind <- glm(cbind(all_a_lp, all_d_lp) ~ treatment + whole.mean + qro, data = brood, family = binomial("logit"))
summary(pl_bind )##
## Call:
## glm(formula = cbind(all_a_lp, all_d_lp) ~ treatment + whole.mean +
## qro, family = binomial("logit"), data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -7.653 -1.318 0.000 2.193 4.821
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 2.5891 0.3105 8.337 < 2e-16 ***
## treatment2 -1.6031 0.2240 -7.157 8.22e-13 ***
## treatment3 -1.1485 0.2136 -5.377 7.58e-08 ***
## treatment4 -1.2431 0.2302 -5.401 6.62e-08 ***
## treatment5 -0.9699 0.2474 -3.921 8.82e-05 ***
## whole.mean -0.9760 0.4614 -2.115 0.0344 *
## qroB3 -0.1974 0.2323 -0.850 0.3954
## qroB4 0.3607 0.2132 1.692 0.0906 .
## qroB5 0.1809 0.1510 1.198 0.2310
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 411.27 on 42 degrees of freedom
## Residual deviance: 331.16 on 34 degrees of freedom
## AIC: 470.26
##
## Number of Fisher Scoring iterations: 5Anova(pl_bind )## Analysis of Deviance Table (Type II tests)
##
## Response: cbind(all_a_lp, all_d_lp)
## LR Chisq Df Pr(>Chisq)
## treatment 63.817 4 4.567e-13 ***
## whole.mean 4.553 1 0.03286 *
## qro 6.895 3 0.07531 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1plot(pl_bind )
lm.emm <- emmeans(pl_bind , "treatment")
pairs(lm.emm)## contrast estimate SE df z.ratio p.value
## treatment1 - treatment2 1.6031 0.224 Inf 7.157 <.0001
## treatment1 - treatment3 1.1485 0.214 Inf 5.377 <.0001
## treatment1 - treatment4 1.2431 0.230 Inf 5.401 <.0001
## treatment1 - treatment5 0.9699 0.247 Inf 3.921 0.0008
## treatment2 - treatment3 -0.4546 0.172 Inf -2.644 0.0626
## treatment2 - treatment4 -0.3600 0.197 Inf -1.832 0.3552
## treatment2 - treatment5 -0.6333 0.208 Inf -3.048 0.0195
## treatment3 - treatment4 0.0946 0.181 Inf 0.524 0.9850
## treatment3 - treatment5 -0.1787 0.205 Inf -0.870 0.9078
## treatment4 - treatment5 -0.2733 0.228 Inf -1.198 0.7526
##
## Results are averaged over the levels of: qro
## Results are given on the log odds ratio (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimatessummary(lm.emm)## treatment emmean SE df asymp.LCL asymp.UCL
## 1 2.206 0.201 Inf 1.812 2.601
## 2 0.603 0.151 Inf 0.307 0.899
## 3 1.058 0.137 Inf 0.790 1.326
## 4 0.963 0.190 Inf 0.591 1.335
## 5 1.236 0.175 Inf 0.894 1.579
##
## Results are averaged over the levels of: qro
## Results are given on the logit (not the response) scale.
## Confidence level used: 0.95plot(brood$treatment, cbind(brood$dead_pupae, brood$live_pupae))
emmeans(pl_bind , pairwise ~ treatment, type = "response")## $emmeans
## treatment prob SE df asymp.LCL asymp.UCL
## 1 0.901 0.0180 Inf 0.860 0.931
## 2 0.646 0.0345 Inf 0.576 0.711
## 3 0.742 0.0262 Inf 0.688 0.790
## 4 0.724 0.0379 Inf 0.644 0.792
## 5 0.775 0.0305 Inf 0.710 0.829
##
## Results are averaged over the levels of: qro
## Confidence level used: 0.95
## Intervals are back-transformed from the logit scale
##
## $contrasts
## contrast odds.ratio SE df null z.ratio p.value
## treatment1 / treatment2 4.969 1.113 Inf 1 7.157 <.0001
## treatment1 / treatment3 3.154 0.674 Inf 1 5.377 <.0001
## treatment1 / treatment4 3.466 0.798 Inf 1 5.401 <.0001
## treatment1 / treatment5 2.638 0.652 Inf 1 3.921 0.0008
## treatment2 / treatment3 0.635 0.109 Inf 1 -2.644 0.0626
## treatment2 / treatment4 0.698 0.137 Inf 1 -1.832 0.3552
## treatment2 / treatment5 0.531 0.110 Inf 1 -3.048 0.0195
## treatment3 / treatment4 1.099 0.199 Inf 1 0.524 0.9850
## treatment3 / treatment5 0.836 0.172 Inf 1 -0.870 0.9078
## treatment4 / treatment5 0.761 0.174 Inf 1 -1.198 0.7526
##
## Results are averaged over the levels of: qro
## P value adjustment: tukey method for comparing a family of 5 estimates
## Tests are performed on the log odds ratio scaleemm3 <- emmeans(pl_bind , ~ treatment, type = "response")
emm3## treatment prob SE df asymp.LCL asymp.UCL
## 1 0.901 0.0180 Inf 0.860 0.931
## 2 0.646 0.0345 Inf 0.576 0.711
## 3 0.742 0.0262 Inf 0.688 0.790
## 4 0.724 0.0379 Inf 0.644 0.792
## 5 0.775 0.0305 Inf 0.710 0.829
##
## Results are averaged over the levels of: qro
## Confidence level used: 0.95
## Intervals are back-transformed from the logit scalepairs(emm3)## contrast odds.ratio SE df null z.ratio p.value
## treatment1 / treatment2 4.969 1.113 Inf 1 7.157 <.0001
## treatment1 / treatment3 3.154 0.674 Inf 1 5.377 <.0001
## treatment1 / treatment4 3.466 0.798 Inf 1 5.401 <.0001
## treatment1 / treatment5 2.638 0.652 Inf 1 3.921 0.0008
## treatment2 / treatment3 0.635 0.109 Inf 1 -2.644 0.0626
## treatment2 / treatment4 0.698 0.137 Inf 1 -1.832 0.3552
## treatment2 / treatment5 0.531 0.110 Inf 1 -3.048 0.0195
## treatment3 / treatment4 1.099 0.199 Inf 1 0.524 0.9850
## treatment3 / treatment5 0.836 0.172 Inf 1 -0.870 0.9078
## treatment4 / treatment5 0.761 0.174 Inf 1 -1.198 0.7526
##
## Results are averaged over the levels of: qro
## P value adjustment: tukey method for comparing a family of 5 estimates
## Tests are performed on the log odds ratio scaleemmdf3 <- as.data.frame(emm3)
p <- ggplot(data = emmdf3, aes(x=treatment, y=prob, fill=treatment)) +
geom_col(position = "dodge", color = "black") +
coord_cartesian(ylim = c(0,1))
p + geom_errorbar(position=position_dodge(.9),width=.25, aes(ymax=asymp.UCL, ymin=asymp.LCL),alpha=.6)+
labs(x="Treatment", y="Probability of Survival", fill = "treatment") + theme_bw() +
scale_fill_manual(values = c("peachpuff3", "darkseagreen", "lightseagreen", "darkolivegreen3", "darkolivegreen4")) +
labs(y = "Probability of Survival",) +
ggtitle("Combined Larvae and Pupae probability of being alive upon colony dissection") +
scale_x_discrete(name = "Treatment",
labels = c("0 PPB", "150 PPB", "1,500 PPB", "15,000 PPB", "150,000 PPB")) +
theme_cowplot() +
theme_classic(base_size = 12)
3.0.7 Eggs
eggs.gn <- glm(eggs ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(eggs.gn)##
## Call:
## glm(formula = eggs ~ treatment + whole.mean + alive + qro, family = "poisson",
## data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -4.836 -2.527 -1.114 1.250 6.211
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.84525 0.23951 3.529 0.000417 ***
## treatment2 -0.46167 0.14170 -3.258 0.001121 **
## treatment3 -1.25533 0.17510 -7.169 7.55e-13 ***
## treatment4 -0.77201 0.16135 -4.785 1.71e-06 ***
## treatment5 -1.09533 0.18912 -5.792 6.97e-09 ***
## whole.mean 2.56703 0.38105 6.737 1.62e-11 ***
## alive -0.05171 0.04621 -1.119 0.263104
## qroB3 0.98714 0.18240 5.412 6.24e-08 ***
## qroB4 1.43150 0.16957 8.442 < 2e-16 ***
## qroB5 0.95463 0.16708 5.714 1.11e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 624.61 on 44 degrees of freedom
## Residual deviance: 290.35 on 35 degrees of freedom
## AIC: 430.74
##
## Number of Fisher Scoring iterations: 6#eggs.gnb <- glm.nb(eggs ~ treatment + whole.mean + alive + qro, data = brood)
#summary(eggs.gnb)
#AIC(eggs.gn, eggs.gnb) #### I think we have too many zeros in eggs, I will switch to a zero-inflated negative binomial model
library(glmmTMB)
egg.zero <- glmmTMB(eggs ~ treatment + whole.mean + qro, family = "nbinom2", brood) # had to take out alive, because it made model not converge -- Anova from glm shows it's not impactful anywas
summary(egg.zero)## Family: nbinom2 ( log )
## Formula: eggs ~ treatment + whole.mean + qro
## Data: brood
##
## AIC BIC logLik deviance df.resid
## 269.7 287.8 -124.9 249.7 35
##
##
## Dispersion parameter for nbinom2 family (): 0.855
##
## Conditional model:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.2479 0.6562 0.378 0.70563
## treatment2 -0.4207 0.5942 -0.708 0.47889
## treatment3 -0.9663 0.5852 -1.651 0.09867 .
## treatment4 -0.7371 0.5806 -1.270 0.20428
## treatment5 -1.0130 0.5952 -1.702 0.08874 .
## whole.mean 3.3460 1.2335 2.713 0.00667 **
## qroB3 0.6962 0.5810 1.198 0.23083
## qroB4 0.9635 0.6375 1.511 0.13073
## qroB5 1.1071 0.4767 2.322 0.02021 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Anova(egg.zero)## Analysis of Deviance Table (Type II Wald chisquare tests)
##
## Response: eggs
## Chisq Df Pr(>Chisq)
## treatment 4.1574 4 0.385128
## whole.mean 7.3585 1 0.006675 **
## qro 6.7957 3 0.078702 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1AIC(eggs.gn, egg.zero)## df AIC
## eggs.gn 10 430.7424
## egg.zero 10 269.7371plot(brood$treatment, brood$eggs)
eggs.sum <- brood %>%
group_by(treatment) %>%
summarise(eggs.m= mean(eggs),
eggs.n = length(eggs),
eggs.al = sd(eggs))
eggs.sum## # A tibble: 5 × 4
## treatment eggs.m eggs.n eggs.al
## <fct> <dbl> <int> <dbl>
## 1 1 14.8 9 27.7
## 2 2 9.11 9 11.7
## 3 3 5.56 9 6.56
## 4 4 6.56 9 5.90
## 5 5 4.33 9 4.39em.egg <- emmeans(egg.zero, "treatment")
pairs(em.egg)## contrast estimate SE df t.ratio p.value
## treatment1 - treatment2 0.4207 0.594 35 0.708 0.9533
## treatment1 - treatment3 0.9663 0.585 35 1.651 0.4761
## treatment1 - treatment4 0.7371 0.581 35 1.269 0.7111
## treatment1 - treatment5 1.0130 0.595 35 1.702 0.4458
## treatment2 - treatment3 0.5456 0.596 35 0.916 0.8888
## treatment2 - treatment4 0.3163 0.594 35 0.532 0.9834
## treatment2 - treatment5 0.5923 0.571 35 1.038 0.8361
## treatment3 - treatment4 -0.2293 0.576 35 -0.398 0.9945
## treatment3 - treatment5 0.0467 0.592 35 0.079 1.0000
## treatment4 - treatment5 0.2760 0.588 35 0.470 0.9896
##
## Results are averaged over the levels of: qro
## Results are given on the log (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimatessummary(em.egg)## treatment emmean SE df lower.CL upper.CL
## 1 2.55 0.404 35 1.727 3.37
## 2 2.13 0.419 35 1.275 2.98
## 3 1.58 0.448 35 0.671 2.49
## 4 1.81 0.439 35 0.920 2.70
## 5 1.53 0.450 35 0.621 2.45
##
## Results are averaged over the levels of: qro
## Results are given on the log (not the response) scale.
## Confidence level used: 0.953.0.8 Count of drones
drones.gn <- glm(drones ~ treatment + whole.mean + alive + qro, data = brood, family = "poisson")
summary(drones.gn)##
## Call:
## glm(formula = drones ~ treatment + whole.mean + alive + qro,
## family = "poisson", data = brood)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.8161 -1.4434 -0.4646 0.9617 3.4167
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.02449 0.27207 -0.090 0.928273
## treatment2 -0.20393 0.16081 -1.268 0.204736
## treatment3 -0.57922 0.17020 -3.403 0.000666 ***
## treatment4 -0.07774 0.15616 -0.498 0.618598
## treatment5 -0.03768 0.16528 -0.228 0.819662
## whole.mean 2.72682 0.35191 7.749 9.28e-15 ***
## alive 0.18765 0.05783 3.245 0.001176 **
## qroB3 0.08375 0.17195 0.487 0.626205
## qroB4 0.08165 0.18022 0.453 0.650487
## qroB5 0.54761 0.13238 4.137 3.52e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for poisson family taken to be 1)
##
## Null deviance: 259.723 on 44 degrees of freedom
## Residual deviance: 97.101 on 35 degrees of freedom
## AIC: 275.05
##
## Number of Fisher Scoring iterations: 5Anova(drones.gn)## Analysis of Deviance Table (Type II tests)
##
## Response: drones
## LR Chisq Df Pr(>Chisq)
## treatment 16.240 4 0.0027131 **
## whole.mean 61.803 1 3.796e-15 ***
## alive 11.775 1 0.0006003 ***
## qro 18.232 3 0.0003939 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1drones.gnb <- glm.nb(drones ~ treatment + whole.mean + alive + qro, data = brood)
summary(drones.gnb)##
## Call:
## glm.nb(formula = drones ~ treatment + whole.mean + alive + qro,
## data = brood, init.theta = 7.713652119, link = log)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -2.1727 -1.1565 -0.1630 0.5915 2.1725
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.8279193 0.4265853 -1.941 0.052282 .
## treatment2 -0.2481215 0.2468556 -1.005 0.314835
## treatment3 -0.6606815 0.2568098 -2.573 0.010092 *
## treatment4 -0.0552440 0.2476179 -0.223 0.823456
## treatment5 -0.0648081 0.2571487 -0.252 0.801021
## whole.mean 3.3691460 0.5439340 6.194 5.86e-10 ***
## alive 0.2876680 0.0865995 3.322 0.000894 ***
## qroB3 0.0709023 0.2636990 0.269 0.788025
## qroB4 -0.0002418 0.2948032 -0.001 0.999346
## qroB5 0.7840522 0.2112366 3.712 0.000206 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for Negative Binomial(7.7137) family taken to be 1)
##
## Null deviance: 139.043 on 44 degrees of freedom
## Residual deviance: 49.747 on 35 degrees of freedom
## AIC: 261.88
##
## Number of Fisher Scoring iterations: 1
##
##
## Theta: 7.71
## Std. Err.: 3.04
## Warning while fitting theta: alternation limit reached
##
## 2 x log-likelihood: -239.877AIC(drones.gn, drones.gnb)## df AIC
## drones.gn 10 275.0533
## drones.gnb 11 261.8765Anova(drones.gnb)## Analysis of Deviance Table (Type II tests)
##
## Response: drones
## LR Chisq Df Pr(>Chisq)
## treatment 9.111 4 0.0583955 .
## whole.mean 37.301 1 1.012e-09 ***
## alive 11.625 1 0.0006507 ***
## qro 15.416 3 0.0014933 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1drone.em <- emmeans(drones.gnb, "treatment")
pairs(drone.em)## contrast estimate SE df z.ratio p.value
## treatment1 - treatment2 0.24812 0.247 Inf 1.005 0.8531
## treatment1 - treatment3 0.66068 0.257 Inf 2.573 0.0755
## treatment1 - treatment4 0.05524 0.248 Inf 0.223 0.9995
## treatment1 - treatment5 0.06481 0.257 Inf 0.252 0.9991
## treatment2 - treatment3 0.41256 0.252 Inf 1.638 0.4730
## treatment2 - treatment4 -0.19288 0.244 Inf -0.789 0.9339
## treatment2 - treatment5 -0.18331 0.247 Inf -0.742 0.9466
## treatment3 - treatment4 -0.60544 0.250 Inf -2.423 0.1093
## treatment3 - treatment5 -0.59587 0.257 Inf -2.318 0.1391
## treatment4 - treatment5 0.00956 0.255 Inf 0.038 1.0000
##
## Results are averaged over the levels of: qro
## Results are given on the log (not the response) scale.
## P value adjustment: tukey method for comparing a family of 5 estimatesdrones.sum <- brood %>%
group_by(treatment) %>%
summarise(m.d = mean(drones),
sd.d = sd(drones),
n.d = length(drones)) %>%
mutate(sed = sd.d / sqrt(n.d))
drones.sum## # A tibble: 5 × 5
## treatment m.d sd.d n.d sed
## <fct> <dbl> <dbl> <int> <dbl>
## 1 1 9.22 7.38 9 2.46
## 2 2 8.44 5.34 9 1.78
## 3 3 7.22 6.42 9 2.14
## 4 4 11.9 8.95 9 2.98
## 5 5 9.56 6.17 9 2.06ggplot(data = drones.sum, aes(x = treatment, y = m.d, fill = treatment)) +
geom_col(position = "dodge", col = "black")+
coord_cartesian(ylim=c(0,16)) +
scale_fill_manual(values = c("peachpuff3", "darkseagreen", "lightseagreen", "darkolivegreen3", "darkolivegreen4"),
name = "Pristine Level",
labels = c("Treatment 1 (control)", "Treatment 2",
"Treatment 3", "Treatment 4", "Treatment 5")) +
geom_errorbar(aes(ymin = m.d - sed,
ymax = m.d + sed),
position = position_dodge(0.9), width = 0.4) +
labs(y = "Mean Number of Drones",) +
ggtitle("Average of Drones Count per Treatment") +
scale_x_discrete(name = "Treatment",
labels = c("0 PPB", "150 PPB", "1,500 PPB", "15,000 PPB", "150,000 PPB")) +
theme_cowplot() +
theme_classic(base_size = 12)
LS0tDQp0aXRsZTogIkJyb29kIE5vIFJvdW5kIDEiDQphdXRob3I6ICJFbWlseSBSdW5uaW9uIg0KZGF0ZTogIjIwMjMtMDEtMjYiDQpvdXRwdXQ6DQogIGh0bWxfZG9jdW1lbnQ6DQogICAgdG9jOiB0cnVlDQogICAgdG9jX2RlcHRoOiA0DQogICAgbnVtYmVyX3NlY3Rpb25zOiB0cnVlDQogICAgdG9jX2Zsb2F0OiB0cnVlDQogICAgdGhlbWU6IGpvdXJuYWwNCiAgICBjb2RlX2Rvd25sb2FkOiB0cnVlDQotLS0NCg0KYGBge3Igc2V0dXAsIGluY2x1ZGU9RkFMU0V9DQprbml0cjo6b3B0c19jaHVuayRzZXQod2FybmluZyA9IEZBTFNFLCBtZXNzYWdlID0gRkFMU0UpDQpgYGANCg0KYGBge3IgbG9hZCBsaWJyYXJpZXMsIGluY2x1ZGU9RkFMU0V9DQpsaWJyYXJ5KHJlYWRyKQ0KbGlicmFyeShrYWJsZUV4dHJhKQ0KbGlicmFyeShzdGF0cykNCmxpYnJhcnkoZ2dwbG90MikNCmxpYnJhcnkoY2FyKQ0KbGlicmFyeShlbW1lYW5zKQ0KbGlicmFyeShNQVNTKQ0KbGlicmFyeShsbWU0KQ0KbGlicmFyeShibG1lY28pDQpsaWJyYXJ5KHRpZHl2ZXJzZSkNCmxpYnJhcnkoZHBseXIpDQpsaWJyYXJ5KGNvd3Bsb3QpDQpsaWJyYXJ5KGJlc3ROb3JtYWxpemUpDQpsaWJyYXJ5KHBsb3RseSkNCmxpYnJhcnkoYWdyaWNvbGFlKSANCmxpYnJhcnkoZ2dwdWJyKQ0KbGlicmFyeShnbHVlKQ0KbGlicmFyeShtdWx0Y29tcFZpZXcpDQpgYGANCg0KDQojIElucHV0IFBvbGxlbiBEYXRhIA0KDQpTdGFydCBieSBpbXB1dGluZyBwb2xsZW4gZGF0YSBhbmQgY3JlYXRpbmcgYSBuZXcgZGF0YSBmcmFtZSB3aXRoIGF2ZXJhZ2UgcG9sbGVuIGNvbnN1bXB0aW9uIG9uIGEgcGVyLWNvbG9ueSBiYXNpcyANCg0KDQpgYGB7cn0NCiMjIyBGaWd1cmUgb3V0IGF2ZXJhZ2UgcG9sbGVuIGNvbnN1bXB0aW9uIGJ5IHRyZWF0bWVudCANCg0KDQpwb2xsZW4gPC0gcmVhZF9jc3YoInBvbGxlbjEuY3N2IiwgY29sX3R5cGVzID0gY29scyhyb3VuZCA9IGNvbF9mYWN0b3IobGV2ZWxzID0gYygiMSIsIA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIjIiKSksIHRyZWF0bWVudCA9IGNvbF9mYWN0b3IobGV2ZWxzID0gYygiMSIsIA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAiMiIsICIzIiwgIjQiLCAiNSIpKSwgcmVwbGljYXRlID0gY29sX2ZhY3RvcihsZXZlbHMgPSBjKCIxIiwgDQogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIjIiLCAiMyIsICI0IiwgIjUiLCAiNiIsICI3IiwgIjkiLCAiMTEiLCANCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAiMTIiKSksIHN0YXJ0X2RhdGUgPSBjb2xfZGF0ZShmb3JtYXQgPSAiJW0vJWQvJVkiKSwgDQogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICBzdGFydF90aW1lID0gY29sX2NoYXJhY3RlcigpLCBlbmRfZGF0ZSA9IGNvbF9kYXRlKGZvcm1hdCA9ICIlbS8lZC8lWSIpLCANCiAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIGVuZF90aW1lID0gY29sX2NoYXJhY3RlcigpKSkNCg0KDQpwb2xsZW4kY29sb255IDwtIGFzLmZhY3Rvcihwb2xsZW4kY29sb255KQ0KcG9sbGVuJHBpZCA8LSBhcy5mYWN0b3IocG9sbGVuJHBpZCkNCnBvbGxlbiRjb3VudCA8LSBhcy5mYWN0b3IocG9sbGVuJGNvdW50KQ0KDQpwb2xsZW4gPC0gc3Vic2V0KHBvbGxlbiwgcG9sbGVuJHJvdW5kICE9IDEpDQpwb2xsZW4NCg0KDQpwb2xsZW4gPC0gbmEub21pdChwb2xsZW4pDQoNCnJhbmdlKHBvbGxlbiRkaWZmZXJlbmNlKQ0KDQojIGdldCByaWQgb2YgbmVnYXRpdmUgbnVtYmVycw0KcG9sbGVuJGRpZmZlcmVuY2VbcG9sbGVuJGRpZmZlcmVuY2UgPCAwXSA8LSBOQQ0KcG9sbGVuIDwtIG5hLm9taXQocG9sbGVuKQ0KcmFuZ2UocG9sbGVuJGRpZmZlcmVuY2UpDQoNCmBgYA0KDQoNCg0KDQojIEF2ZXJhZ2UgcG9sbGVuIGNvbnN1bXB0aW9uIHBlciBjb2xvbnkNCg0KYGBge3J9DQpwb2xsZW4kd2hvbGVfZGlmIDwtIGFzLm51bWVyaWMocG9sbGVuJGRpZmZlcmVuY2UpDQp3ZC4xIDwtIGxtKGRpZmZlcmVuY2UgfiB0cmVhdG1lbnQsIGRhdGEgPSBwb2xsZW4pDQpzdW1tYXJ5KHdkLjEpDQp3ZC5lbW0gPC0gZW1tZWFucyh3ZC4xLCAidHJlYXRtZW50IikNCnN1bW1hcnkod2QuZW1tKQ0KDQp3ZC5zdW1tYXJ5IDwtIHBvbGxlbiAlPiUgDQogIGdyb3VwX2J5KGNvbG9ueSkgJT4lDQogIHN1bW1hcml6ZSh3aG9sZS5tZWFuID0gbWVhbihkaWZmZXJlbmNlKSkNCg0KYXMuZGF0YS5mcmFtZSh3ZC5zdW1tYXJ5KSAgIyB0aGlzIGlzIHRoZSBkYXRhIGZyYW1lIEkgd2lsbCBtZXJnZSB3aXRoIHN1YnNlcXVlbnQgZGF0YSBmcmFtZXMgdG8gY29udGFpbiBhdmVyYWdlIHBvbGxlbiBjb25zdW1wdGlvbiBwZXIgY29sb255IGFzIGEgbmV3IGNvbHVtbiAgDQoNCmBgYA0KDQoNCmBgYHtyfQ0KDQpicm9vZCA8LSByZWFkX2NzdigiYnJvb2QuY3N2IiwgY29sX3R5cGVzID0gY29scyhyb3VuZCA9IGNvbF9mYWN0b3IobGV2ZWxzID0gYygiMSIsIA0KICAgICIyIikpLCB0cmVhdG1lbnQgPSBjb2xfZmFjdG9yKGxldmVscyA9IGMoIjEiLCANCiAgICAiMiIsICIzIiwgIjQiLCAiNSIpKSwgcmVwbGljYXRlID0gY29sX2ZhY3RvcihsZXZlbHMgPSBjKCIxIiwgDQogICAgIjIiLCAiMyIsICI0IiwgIjUiLCAiNyIsICI5IiwgIjExIiwgIjEyIikpKSkNCg0KYnJvb2QkY29sb255IDwtIGFzLmZhY3Rvcihicm9vZCRjb2xvbnkpDQoNCmJyb29kIDwtIHN1YnNldChicm9vZCwgYnJvb2Qkcm91bmQgIT0gMSkNCg0KDQpicm9vZCA8LSBtZXJnZSh3ZC5zdW1tYXJ5LCBicm9vZCwgYnkueCA9ICJjb2xvbnkiKSAgICAjIHRoaXMgaXMgZ29vZCBiZWNhdXNlIEkgY2FsY3VsYXRlZCBhdmVyYWdlIHBvbGxlbiBjb25zdW1wdGlvbiB0d28gZGlmZmVyZW50IHdheXMgKGF2Z19wb2xsZW4gd2FzIGNhbGN1bGF0ZWQgYSBjb3VwbGUgbW9udGhzIGFnbyBtYW51YWxseSkgYW5kIGl0J3MgdGhlIHNhbWUgbnVtYmVycyBhcyB0aGUgZ3JvdXBfYnkgZnVuY3Rpb24gdGhhdCBjcmVhdGVkIHdob2xlLm1lYW4NCg0KbW9ydGFsaXR5ICA8LSByZWFkX2Nzdigic3Vydml2aW5nIHdvcmtlcnMgcGVyIGNvbG9ueS5jc3YiKQ0KDQptb3J0YWxpdHkkY29sb255IDwtIGFzLmZhY3Rvcihtb3J0YWxpdHkkY29sb255KQ0KDQpicm9vZCA8LSBtZXJnZShicm9vZCwgbW9ydGFsaXR5LCBieS54ID0gImNvbG9ueSIpDQoNCg0KYGBgDQoNCg0KIyBCcm9vZCBjZWxscyANCg0KYGBge3J9DQoNCnBsb3QoYnJvb2QkdHJlYXRtZW50LCBicm9vZCRicm9vZF9jZWxscykNCg0KcmFuZ2UoYnJvb2QkYnJvb2RfY2VsbHMpDQoNCmBgYA0KDQoNClRvdGFsIGJyb29kIGNvdW50IG5vdCBhZmZlY3RlZCBieSB0cmVhdG1lbnQsIEFub3ZhID0gMC4zNDMzMjgNCg0KYGBge3J9DQoNCmIuZ24gPC0gZ2xtKGJyb29kX2NlbGxzIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICwgZGF0YSA9IGJyb29kLCBmYW1pbHkgPSAicG9pc3NvbiIpDQpzdW1tYXJ5KGIuZ24pICAjIFdlIGRvIGhhdmUgb3ZlcmRpc3BlcnNpb24gaGVyZSBub3cgdGhhdCB3ZSB0b29rIG91dCByb3VuZCAxDQoNCkFub3ZhKGIuZ24pDQoNCmIuZ25iIDwtIGdsbS5uYihicm9vZF9jZWxscyB+IHRyZWF0bWVudCArIHdob2xlLm1lYW4gKyBhbGl2ZSAsIGRhdGEgPSBicm9vZCkNCg0KQUlDKGIuZ24sIGIuZ25iKQ0KDQpBbm92YShiLmduYikNCg0KcGxvdChiLmduKQ0KcGxvdChiLmduYikNCg0KYnJvb2QuZW1tIDwtIGVtbWVhbnMoYi5nbmIsICJ0cmVhdG1lbnQiKQ0KcGFpcnMoYnJvb2QuZW1tKQ0Kc3VtbWFyeShicm9vZC5lbW0pDQoNCnN1bV9icm9vZCA8LSBicm9vZCAlPiUNCiAgZ3JvdXBfYnkodHJlYXRtZW50KSAlPiUNCiAgc3VtbWFyaXNlKCBtZWFuLmIgPSBtZWFuKGJyb29kX2NlbGxzKSwgDQogICAgICAgICAgICAgbi5iID0gbGVuZ3RoKGJyb29kX2NlbGxzKSwgDQogICAgICAgICAgICAgc2QuYiA9IHNkKGJyb29kX2NlbGxzKSkgJT4lDQogIG11dGF0ZShzZWIgPSBzZC5iIC8gc3FydChuLmIpKQ0KDQpzdW1fYnJvb2QNCg0KYGBgDQoNCmBgYHtyfQ0KZ2dwbG90KGRhdGEgPSBzdW1fYnJvb2QsIGFlcyh4ID0gdHJlYXRtZW50LCB5ID0gbWVhbi5iLCBmaWxsID0gdHJlYXRtZW50KSkgKw0KICBnZW9tX2NvbChwb3NpdGlvbiA9ICJkb2RnZSIsIGNvbCA9ICJibGFjayIpKw0KICBjb29yZF9jYXJ0ZXNpYW4oeWxpbT1jKDIwLDYwKSkgKw0KICBzY2FsZV9maWxsX21hbnVhbCh2YWx1ZXMgPSBjKCJwZWFjaHB1ZmYzIiwgImRhcmtzZWFncmVlbiIsICJsaWdodHNlYWdyZWVuIiwgImRhcmtvbGl2ZWdyZWVuMyIsICJkYXJrb2xpdmVncmVlbjQiKSwNCiAgICAgICAgICAgICAgICAgICAgbmFtZSA9ICJQcmlzdGluZSBMZXZlbCIsDQogICAgICAgICAgICAgICAgICAgIGxhYmVscyA9IGMoIlRyZWF0bWVudCAxIChjb250cm9sKSIsICJUcmVhdG1lbnQgMiIsIA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICJUcmVhdG1lbnQgMyIsICJUcmVhdG1lbnQgNCIsICJUcmVhdG1lbnQgNSIpKSArIA0KICBnZW9tX2Vycm9yYmFyKGFlcyh5bWluID0gbWVhbi5iIC0gc2ViLCANCiAgICAgICAgICAgICAgICAgICAgeW1heCA9IG1lYW4uYiArIHNlYiksDQogICAgICAgICAgICAgICAgcG9zaXRpb24gPSBwb3NpdGlvbl9kb2RnZSgwLjkpLCB3aWR0aCA9IDAuNCkgKw0KICBsYWJzKHkgPSAiTWVhbiBCcm9vZCBDZWxsbHMiLCkgKw0KICBnZ3RpdGxlKCJBdmVyYWdlIG9mIEFsbCBCcm9vZCBDZWxscyBwZXIgVHJlYXRtZW50IikgKw0KICBzY2FsZV94X2Rpc2NyZXRlKG5hbWUgPSAiVHJlYXRtZW50IiwgDQogICAgICAgICAgICAgICAgICAgbGFiZWxzID0gYygiMCBQUEIiLCAiMTUwIFBQQiIsICIxLDUwMCBQUEIiLCAiMTUsMDAwIFBQQiIsICIxNTAsMDAwIFBQQiIpKSArDQogIHRoZW1lX2Nvd3Bsb3QoKSArDQogIHRoZW1lX2NsYXNzaWMoYmFzZV9zaXplID0gMTIpDQpgYGANCg0KDQojIyMgaG9uZXkgcG90cw0KDQpUcmVhdG1lbnQgZG9lcyBub3QqIChBbm92YSBzaWcgYnV0IHBvc3QgaG9jIG5vdCkgaW1wYWN0IGNvdW50IG9mIGhvbmV5IHBvdHMgDQoNCkJpZ2dlc3QgZGlmZmVyZW5jZSBpcyBiZXR3ZWVuIHRyZWF0bWVudHMgdHdvIGFuZCBvbmUNCg0KYGBge3J9DQoNCmhwLmduIDwtIGdsbShob25leV9wb3QgfiB0cmVhdG1lbnQgKyB3aG9sZS5tZWFuICsgYWxpdmUgLCBkYXRhID0gYnJvb2QsIGZhbWlseSA9ICJwb2lzc29uIikNCnN1bW1hcnkoaHAuZ24pICMgb3ZlcmRpc3BlcnNpb24gbm90IGhvcnJpYmxlIGJ1dCBub3QgZ3JlYXQgDQpocC5nbmIgPC0gZ2xtLm5iKGhvbmV5X3BvdCB+IHRyZWF0bWVudCArIHdob2xlLm1lYW4gKyBhbGl2ZSAsIGRhdGEgPSBicm9vZCkNCnN1bW1hcnkoaHAuZ25iKQ0KQUlDKGhwLmduLCBocC5nbmIpDQoNCnBsb3QoaHAuZ24pDQpwbG90KGhwLmduYikgIyBzdGljayB3aXRoIGdsbS5uYg0KDQpBbm92YShocC5nbmIpDQpBbm92YShocC5nbikNCg0KcGxvdChicm9vZCR0cmVhdG1lbnQsIGJyb29kJGhvbmV5X3BvdCkNCg0KaHAuZW1tIDwtIGVtbWVhbnMoaHAuZ25iLCAidHJlYXRtZW50IikNCnN1bW1hcnkoaHAuZW1tKQ0KcGFpcnMoaHAuZW1tKQ0KDQoNCnN1bV9ocCA8LSBicm9vZCAlPiUNCiAgZ3JvdXBfYnkodHJlYXRtZW50KSAlPiUNCiAgc3VtbWFyaXNlKCBtZWFuLmggPSBtZWFuKGhvbmV5X3BvdCksIA0KICAgICAgICAgICAgIG4uaCA9IGxlbmd0aChob25leV9wb3QpLCANCiAgICAgICAgICAgICBzZC5oID0gc2QoaG9uZXlfcG90KSkgJT4lDQogIG11dGF0ZShzZWhwID0gc2QuaCAvIHNxcnQobi5oKSkNCg0Kc3VtX2hwDQoNCmBgYA0KDQpgYGB7cn0NCmdncGxvdChkYXRhID0gc3VtX2hwLCBhZXMoeCA9IHRyZWF0bWVudCwgeSA9IG1lYW4uaCwgZmlsbCA9IHRyZWF0bWVudCkpICsNCiAgZ2VvbV9jb2wocG9zaXRpb24gPSAiZG9kZ2UiLCBjb2wgPSAiYmxhY2siKSsNCiAgY29vcmRfY2FydGVzaWFuKHlsaW09YygyLDEwKSkgKw0KICBzY2FsZV9maWxsX21hbnVhbCh2YWx1ZXMgPSBjKCJwZWFjaHB1ZmYzIiwgImRhcmtzZWFncmVlbiIsICJsaWdodHNlYWdyZWVuIiwgImRhcmtvbGl2ZWdyZWVuMyIsICJkYXJrb2xpdmVncmVlbjQiKSwNCiAgICAgICAgICAgICAgICAgICAgbmFtZSA9ICJQcmlzdGluZSBMZXZlbCIsDQogICAgICAgICAgICAgICAgICAgIGxhYmVscyA9IGMoIlRyZWF0bWVudCAxIChjb250cm9sKSIsICJUcmVhdG1lbnQgMiIsIA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICJUcmVhdG1lbnQgMyIsICJUcmVhdG1lbnQgNCIsICJUcmVhdG1lbnQgNSIpKSArIA0KICBnZW9tX2Vycm9yYmFyKGFlcyh5bWluID0gbWVhbi5oIC0gc2VocCwgDQogICAgICAgICAgICAgICAgICAgIHltYXggPSBtZWFuLmggKyBzZWhwKSwNCiAgICAgICAgICAgICAgICBwb3NpdGlvbiA9IHBvc2l0aW9uX2RvZGdlKDAuOSksIHdpZHRoID0gMC40KSArDQogIGxhYnMoeSA9ICJNZWFuIEhvbmV5IFBvdHMiLCkgKw0KICBnZ3RpdGxlKCJBdmVyYWdlIG9mIEhvbmV5IFBvdHMgcGVyIFRyZWF0bWVudCIpICsNCiAgc2NhbGVfeF9kaXNjcmV0ZShuYW1lID0gIlRyZWF0bWVudCIsIA0KICAgICAgICAgICAgICAgICAgIGxhYmVscyA9IGMoIjAgUFBCIiwgIjE1MCBQUEIiLCAiMSw1MDAgUFBCIiwgIjE1LDAwMCBQUEIiLCAiMTUwLDAwMCBQUEIiKSkgKw0KICB0aGVtZV9jb3dwbG90KCkgKw0KICB0aGVtZV9jbGFzc2ljKGJhc2Vfc2l6ZSA9IDEyKQ0KYGBgDQoNCg0KIyMjIExhcnZhZQ0KDQojIyMjIExpdmUgTGFydmFlDQoNCmBgYHtyfQ0KDQpsbC5nbiA8LSBnbG0obGl2ZV9sYXJ2YWUgfiB0cmVhdG1lbnQgKyB3aG9sZS5tZWFuICsgYWxpdmUgICsgcXJvLCBkYXRhID0gYnJvb2QsIGZhbWlseSA9ICJwb2lzc29uIikNCnN1bW1hcnkobGwuZ24pDQpsbC5nbmIgPC0gZ2xtLm5iKGxpdmVfbGFydmFlIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kKQ0Kc3VtbWFyeShsbC5nbmIpDQpBSUMobGwuZ24sIGxsLmduYikgICAgICAjIHBvaXNzb24gb3ZlcmRpc3BlcnNlZCwgbmIgQUlDIGJldHRlciANCg0KcGxvdChsbC5nbmIpDQoNCkFub3ZhKGxsLmduKQ0KQW5vdmEobGwuZ25iKQ0KDQpwbG90KGJyb29kJHRyZWF0bWVudCwgYnJvb2QkbGl2ZV9sYXJ2YWUpDQoNCmxsLnN1bSA8LSBicm9vZCAlPiUNCiAgZ3JvdXBfYnkodHJlYXRtZW50KSAlPiUNCiAgc3VtbWFyaXNlKGxsLm09IG1lYW4obGl2ZV9sYXJ2YWUpLCANCiAgICAgICAgICAgIGxsLm4gPSBsZW5ndGgobGl2ZV9sYXJ2YWUpLCANCiAgICAgICAgICAgIHNkLmxsID0gc2QobGl2ZV9sYXJ2YWUpKQ0KDQpsbC5zdW0NCg0KYGBgDQojIyMjIERlYWQgTGFydmFlDQoNCmBgYHtyfQ0KDQpkbC5nbiA8LSBnbG0oZGVhZF9sYXJ2YWUgfiB0cmVhdG1lbnQgKyB3aG9sZS5tZWFuICsgYWxpdmUgICsgcXJvLCBkYXRhID0gYnJvb2QsIGZhbWlseSA9ICJwb2lzc29uIikNCnN1bW1hcnkoZGwuZ24pDQpkbC5nbmIgPC0gZ2xtLm5iKGRlYWRfbGFydmFlIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kKQ0Kc3VtbWFyeShkbC5nbmIpDQpBSUMoZGwuZ24sIGRsLmduYikgICAgICNwb2lzc29uIG92ZXJkaXNwZXJzZWQgYW5kIG5iIEFJQyBiZXR0ZXIgDQpwbG90KGRsLmduYikNCg0KQW5vdmEoZGwuZ25iKQ0KDQpwbG90KGJyb29kJHRyZWF0bWVudCwgYnJvb2QkZGVhZF9sYXJ2YWUpDQoNCmRsLnN1bSA8LSBicm9vZCAlPiUNCiAgZ3JvdXBfYnkodHJlYXRtZW50KSAlPiUNCiAgc3VtbWFyaXNlKGRsLm09IG1lYW4oZGVhZF9sYXJ2YWUpLCANCiAgICAgICAgICAgIGRsLm4gPSBsZW5ndGgoZGVhZF9sYXJ2YWUpLCANCiAgICAgICAgICAgIHNkLmRsID0gc2QoZGVhZF9sYXJ2YWUpKQ0KDQpkbC5zdW0NCg0KYGBgDQoNCiMjIyMgQWxsIGxhcnZhZQ0KDQpgYGB7cn0NCg0KYnJvb2QkYWxsX2xhcnZhZSA8LSBicm9vZCRkZWFkX2xhcnZhZSArIGJyb29kJGxpdmVfbGFydmFlDQoNCmFsLmduIDwtIGdsbShhbGxfbGFydmFlIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kLCBmYW1pbHkgPSAicG9pc3NvbiIpDQpzdW1tYXJ5KGFsLmduKQ0KYWwuZ25iIDwtIGdsbS5uYihhbGxfbGFydmFlIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kKQ0Kc3VtbWFyeShhbC5nbmIpDQpBSUMoYWwuZ24sIGFsLmduYikNCg0KQW5vdmEoYWwuZ25iKQ0KDQpwbG90KGJyb29kJHRyZWF0bWVudCwgYnJvb2QkYWxsX2xhcnZhZSkNCg0KDQphbC5zdW0gPC0gYnJvb2QgJT4lDQogIGdyb3VwX2J5KHRyZWF0bWVudCkgJT4lDQogIHN1bW1hcmlzZShhbC5tPSBtZWFuKGFsbF9sYXJ2YWUpLCANCiAgICAgICAgICAgIGFsLm4gPSBsZW5ndGgoYWxsX2xhcnZhZSksIA0KICAgICAgICAgICAgc2QuYWwgPSBzZChhbGxfbGFydmFlKSkNCg0KYWwuc3VtDQpgYGANCg0KIyMjIyBMYXJ2YWUgT3V0bGllcnMNCg0KYGBge3J9DQoNCm91dCA8LSBib3hwbG90LnN0YXRzKGJyb29kJGRlYWRfbGFydmFlKSRvdXQgDQpvdXQNCm91dF9kbCA8LSB3aGljaChicm9vZCRkZWFkX2xhcnZhZSAlaW4lIGMob3V0KSkgDQpvdXRfZGwNCg0KYnJvb2Rbb3V0X2RsLCBdDQoNCmJyb29kX2RsX291dCA8LSBicm9vZFstYygyMiwgMzIsIDM0LCAzOSwgNDUpLCBdIA0KDQoNCmRsLmduIDwtIGdsbShkZWFkX2xhcnZhZSB+IHRyZWF0bWVudCArIHdob2xlLm1lYW4gKyBhbGl2ZSAgKyBxcm8sIGRhdGEgPSBicm9vZF9kbF9vdXQsIGZhbWlseSA9ICJwb2lzc29uIikNCnN1bW1hcnkoZGwuZ24pDQpkbC5nbmIgPC0gZ2xtLm5iKGRlYWRfbGFydmFlIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kX2RsX291dCkNCnN1bW1hcnkoZGwuZ25iKQ0KQW5vdmEoZGwuZ25iKQ0KDQpwbG90KGJyb29kX2RsX291dCR0cmVhdG1lbnQsIGJyb29kX2RsX291dCRkZWFkX2xhcnZhZSkNCg0KYGBgDQoNCg0KIyMjIyBjYmluZCBMYXJ2YWUNCg0KYGBge3J9DQojIyMjIGNiaW5kIExhcnZhZQ0KDQpsYXJ2YWUubW9ydGFsaXR5IDwtIGdsbShjYmluZChsaXZlX2xhcnZhZSwgZGVhZF9sYXJ2YWUpIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIHFybywgZGF0YSA9IGJyb29kLCBmYW1pbHkgPSBiaW5vbWlhbCgibG9naXQiKSkNCnN1bW1hcnkobGFydmFlLm1vcnRhbGl0eSkNCg0KQW5vdmEobGFydmFlLm1vcnRhbGl0eSkNCg0KcGxvdChsYXJ2YWUubW9ydGFsaXR5KQ0KDQpsbS5lbW0gPC0gZW1tZWFucyhsYXJ2YWUubW9ydGFsaXR5LCAidHJlYXRtZW50IikNCnBhaXJzKGxtLmVtbSkNCnN1bW1hcnkobG0uZW1tKQ0KDQpwbG90KGJyb29kJHRyZWF0bWVudCwgY2JpbmQoYnJvb2QkZGVhZF9sYXJ2YWUsIGJyb29kJGxpdmVfbGFydmFlKSkNCg0KZW1tZWFucyhsYXJ2YWUubW9ydGFsaXR5LCBwYWlyd2lzZSB+IHRyZWF0bWVudCwgdHlwZSA9ICJyZXNwb25zZSIpDQoNCmVtbTEgPC0gZW1tZWFucyhsYXJ2YWUubW9ydGFsaXR5LCB+IHRyZWF0bWVudCwgdHlwZSA9ICJyZXNwb25zZSIpDQplbW0xDQpwYWlycyhlbW0xKQ0KDQplbW1kZiA8LSBhcy5kYXRhLmZyYW1lKGVtbTEpDQoNCnAgPC0gZ2dwbG90KGRhdGEgPSBlbW1kZiwgYWVzKHg9dHJlYXRtZW50LCB5PXByb2IsIGZpbGw9dHJlYXRtZW50KSkgKyANCiAgZ2VvbV9jb2wocG9zaXRpb24gPSAiZG9kZ2UiLCBjb2xvciA9ICJibGFjayIpICsNCiAgY29vcmRfY2FydGVzaWFuKHlsaW0gPSBjKDAsMSkpDQoNCnAgKyBnZW9tX2Vycm9yYmFyKHBvc2l0aW9uPXBvc2l0aW9uX2RvZGdlKC45KSx3aWR0aD0uMjUsIGFlcyh5bWF4PWFzeW1wLlVDTCwgeW1pbj1hc3ltcC5MQ0wpLGFscGhhPS42KSsNCiAgbGFicyh4PSJUcmVhdG1lbnQiLCB5PSJQcm9iYWJpbGl0eSBvZiBTdXJ2aXZhbCIsIGZpbGwgPSAidHJlYXRtZW50IikgKyB0aGVtZV9idygpICsNCiAgc2NhbGVfZmlsbF9tYW51YWwodmFsdWVzID0gYygicGVhY2hwdWZmMyIsICJkYXJrc2VhZ3JlZW4iLCAibGlnaHRzZWFncmVlbiIsICJkYXJrb2xpdmVncmVlbjMiLCAiZGFya29saXZlZ3JlZW40IikpICsNCiAgIGxhYnMoeSA9ICJQcm9iYWJpbGl0eSBvZiBTdXJ2aXZhbCIsKSArDQogIGdndGl0bGUoIlByb2JhYmlsaXR5IG9mIGxhcnZhZSBiZWluZyBhbGl2ZSB1cG9uIGNvbG9ueSBkaXNzZWN0aW9uIikgKw0KICBzY2FsZV94X2Rpc2NyZXRlKG5hbWUgPSAiVHJlYXRtZW50IiwgDQogICAgICAgICAgICAgICAgICAgbGFiZWxzID0gYygiMCBQUEIiLCAiMTUwIFBQQiIsICIxLDUwMCBQUEIiLCAiMTUsMDAwIFBQQiIsICIxNTAsMDAwIFBQQiIpKSArDQogIHRoZW1lX2Nvd3Bsb3QoKSArDQogIHRoZW1lX2NsYXNzaWMoYmFzZV9zaXplID0gMTIpDQpgYGANCg0KIyMjIFB1cGFlDQoNCiMjIyMgTGl2ZSBQdXBhZQ0KDQpgYGB7cn0NCmxwLmduIDwtIGdsbShsaXZlX3B1cGFlIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kLCBmYW1pbHkgPSAicG9pc3NvbiIpDQpzdW1tYXJ5KGxwLmduKQ0KbHAuZ25iIDwtIGdsbS5uYihsaXZlX3B1cGFlIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kKQ0Kc3VtbWFyeShscC5nbmIpDQpBSUMobHAuZ24sIGxwLmduYikNCg0KQW5vdmEobHAuZ25iKQ0KDQpwbG90KGJyb29kJHRyZWF0bWVudCwgYnJvb2QkbGl2ZV9wdXBhZSkNCg0KDQpscC5zdW0gPC0gYnJvb2QgJT4lDQogIGdyb3VwX2J5KHRyZWF0bWVudCkgJT4lDQogIHN1bW1hcmlzZShscC5tPSBtZWFuKGxpdmVfcHVwYWUpLCANCiAgICAgICAgICAgIHNkLmxwID0gc2QobGl2ZV9wdXBhZSkpDQoNCmxwLnN1bQ0KYGBgDQoNCg0KIyMjIyBEZWFkIFB1cGFlDQoNCmBgYHtyfQ0KDQpkcC5nbiA8LSBnbG0oZGVhZF9wdXBhZSB+IHRyZWF0bWVudCArIHdob2xlLm1lYW4gKyBhbGl2ZSAgKyBxcm8sIGRhdGEgPSBicm9vZCwgZmFtaWx5ID0gInBvaXNzb24iKQ0Kc3VtbWFyeShkcC5nbikNCmRwLmduYiA8LSBnbG0ubmIoZGVhZF9wdXBhZSB+IHRyZWF0bWVudCArIHdob2xlLm1lYW4gKyBhbGl2ZSAgKyBxcm8sIGRhdGEgPSBicm9vZCkgICAgIyBzdGljayB3aXRoIHRoaXMgbW9kZWwgb3ZlcmFsbCANCnN1bW1hcnkoZHAuZ25iKQ0KQUlDKGRwLmduLCBkcC5nbmIpDQoNCkFub3ZhKGRwLmduYikNCnBsb3QoYnJvb2QkdHJlYXRtZW50LCBicm9vZCRkZWFkX3B1cGFlKQ0KDQplbS5kcCA8LSBlbW1lYW5zKGRwLmduYiwgInRyZWF0bWVudCIpDQpzdW1tYXJ5KGVtLmRwKQ0KcGFpcnMoZW0uZHApDQoNCmRwLmduYjEgPC0gZ2xtLm5iKGRlYWRfcHVwYWUgfiB0cmVhdG1lbnQgKyB3aG9sZS5tZWFuICsgcXJvLCBkYXRhID0gYnJvb2QpDQpkcC5nbmIyIDwtIGdsbS5uYihkZWFkX3B1cGFlIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIHFybyArIGFsaXZlLCBkYXRhID0gYnJvb2QpDQpkcC5nbmIzIDwtIGdsbS5uYihkZWFkX3B1cGFlIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIHFybyAsIGRhdGEgPSBicm9vZCkNCkFJQyhkcC5nbmIsIGRwLmduYjEsIGRwLmduYjIsIGRwLmduYjMpDQoNCnBsb3QoZHAuZ25iMSkNCnBsb3QoZHAuZ25iKQ0KcGxvdChkcC5nbmIzKQ0KDQoNCmVtLmRwMyA8LSBlbW1lYW5zKGRwLmduYiwgInRyZWF0bWVudCIpDQpzdW1tYXJ5KGVtLmRwMykNCnBhaXJzKGVtLmRwMykNCg0KDQoNCmRwLnN1bSA8LSBicm9vZCAlPiUNCiAgZ3JvdXBfYnkodHJlYXRtZW50KSAlPiUNCiAgc3VtbWFyaXNlKGRwLm09IG1lYW4oZGVhZF9wdXBhZSksIA0KICAgICAgICAgICAgc2QuZHAgPSBzZChkZWFkX3B1cGFlKSwNCiAgICAgICAgICAgIG4uZHAgPSBsZW5ndGgoZGVhZF9wdXBhZSkpICU+JQ0KICBtdXRhdGUoc2VkcCA9IHNkLmRwLyBzcXJ0KG4uZHApKQ0KDQpkcC5zdW0NCg0KYGBgDQoNCmBgYHtyfQ0KZ2dwbG90KGRhdGEgPSBkcC5zdW0sIGFlcyh4ID0gdHJlYXRtZW50LCB5ID0gZHAubSwgZmlsbCA9IHRyZWF0bWVudCkpICsNCiAgZ2VvbV9jb2wocG9zaXRpb24gPSAiZG9kZ2UiLCBjb2wgPSAiYmxhY2siKSsNCiAgY29vcmRfY2FydGVzaWFuKHlsaW09YygwLDE0KSkgKw0KICBzY2FsZV9maWxsX21hbnVhbCh2YWx1ZXMgPSBjKCJwZWFjaHB1ZmYzIiwgImRhcmtzZWFncmVlbiIsICJsaWdodHNlYWdyZWVuIiwgImRhcmtvbGl2ZWdyZWVuMyIsICJkYXJrb2xpdmVncmVlbjQiKSwNCiAgICAgICAgICAgICAgICAgICAgbmFtZSA9ICJQcmlzdGluZSBMZXZlbCIsDQogICAgICAgICAgICAgICAgICAgIGxhYmVscyA9IGMoIlRyZWF0bWVudCAxIChjb250cm9sKSIsICJUcmVhdG1lbnQgMiIsIA0KICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgICJUcmVhdG1lbnQgMyIsICJUcmVhdG1lbnQgNCIsICJUcmVhdG1lbnQgNSIpKSArIA0KICBnZW9tX2Vycm9yYmFyKGFlcyh5bWluID0gZHAubSAtIHNlZHAsIA0KICAgICAgICAgICAgICAgICAgICB5bWF4ID0gZHAubSArIHNlZHApLA0KICAgICAgICAgICAgICAgIHBvc2l0aW9uID0gcG9zaXRpb25fZG9kZ2UoMC45KSwgd2lkdGggPSAwLjQpICsNCiAgbGFicyh5ID0gIk1lYW4gRGVhZCBQdXBhZSIsKSArDQogIGdndGl0bGUoIkF2ZXJhZ2Ugb2YgRGVhZCBQdXBhZSBwZXIgVHJlYXRtZW50IikgKw0KICBzY2FsZV94X2Rpc2NyZXRlKG5hbWUgPSAiVHJlYXRtZW50IiwgDQogICAgICAgICAgICAgICAgICAgbGFiZWxzID0gYygiMCBQUEIiLCAiMTUwIFBQQiIsICIxLDUwMCBQUEIiLCAiMTUsMDAwIFBQQiIsICIxNTAsMDAwIFBQQiIpKSArDQogIHRoZW1lX2Nvd3Bsb3QoKSArDQogIHRoZW1lX2NsYXNzaWMoYmFzZV9zaXplID0gMTIpDQpgYGANCg0KDQoNCiMjIyMgQWxsIFB1cGFlDQoNCmBgYHtyfQ0KYnJvb2QkYWxsX3B1cGFlIDwtIGJyb29kJGRlYWRfcHVwYWUgKyBicm9vZCRsaXZlX3B1cGFlDQoNCmFwLmduIDwtIGdsbShhbGxfcHVwYWUgfiB0cmVhdG1lbnQgKyB3aG9sZS5tZWFuICsgYWxpdmUgICsgcXJvLCBkYXRhID0gYnJvb2QsIGZhbWlseSA9ICJwb2lzc29uIikNCnN1bW1hcnkoYXAuZ24pDQphcC5nbmIgPC0gZ2xtLm5iKGFsbF9wdXBhZSB+IHRyZWF0bWVudCArIHdob2xlLm1lYW4gKyBhbGl2ZSAgKyBxcm8sIGRhdGEgPSBicm9vZCkNCnN1bW1hcnkoYXAuZ25iKQ0KQUlDKGFwLmduLCBhcC5nbmIpDQoNCkFub3ZhKGFwLmduYikNCg0KcGxvdChicm9vZCR0cmVhdG1lbnQsIGJyb29kJGFsbF9wdXBhZSkNCg0KDQphcC5zdW0gPC0gYnJvb2QgJT4lDQogIGdyb3VwX2J5KHRyZWF0bWVudCkgJT4lDQogIHN1bW1hcmlzZShhbC5tPSBtZWFuKGFsbF9wdXBhZSksIA0KICAgICAgICAgICAgYWwubiA9IGxlbmd0aChhbGxfcHVwYWUpLCANCiAgICAgICAgICAgIHNkLmFsID0gc2QoYWxsX3B1cGFlKSkgJT4lDQogIG11dGF0ZShzZXAgPSBzZC5hbC9zcXJ0KGFsLm4pKQ0KDQphcC5zdW0NCg0KDQpgYGANCmBgYHtyfQ0KDQpnZ3Bsb3QoZGF0YSA9IGFwLnN1bSwgYWVzKHggPSB0cmVhdG1lbnQsIHkgPSBhbC5tLCBmaWxsID0gdHJlYXRtZW50KSkgKw0KICBnZW9tX2NvbChwb3NpdGlvbiA9ICJkb2RnZSIsIGNvbCA9ICJibGFjayIpKw0KICBjb29yZF9jYXJ0ZXNpYW4oeWxpbT1jKDAsMjApKSArDQogIHNjYWxlX2ZpbGxfbWFudWFsKHZhbHVlcyA9IGMoInBlYWNocHVmZjMiLCAiZGFya3NlYWdyZWVuIiwgImxpZ2h0c2VhZ3JlZW4iLCAiZGFya29saXZlZ3JlZW4zIiwgImRhcmtvbGl2ZWdyZWVuNCIpLA0KICAgICAgICAgICAgICAgICAgICBuYW1lID0gIlByaXN0aW5lIExldmVsIiwNCiAgICAgICAgICAgICAgICAgICAgbGFiZWxzID0gYygiVHJlYXRtZW50IDEgKGNvbnRyb2wpIiwgIlRyZWF0bWVudCAyIiwgDQogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIlRyZWF0bWVudCAzIiwgIlRyZWF0bWVudCA0IiwgIlRyZWF0bWVudCA1IikpICsgDQogIGdlb21fZXJyb3JiYXIoYWVzKHltaW4gPSBhbC5tIC0gc2VwLCANCiAgICAgICAgICAgICAgICAgICAgeW1heCA9IGFsLm0gKyBzZXApLA0KICAgICAgICAgICAgICAgIHBvc2l0aW9uID0gcG9zaXRpb25fZG9kZ2UoMC45KSwgd2lkdGggPSAwLjQpICsNCiAgbGFicyh5ID0gIk1lYW4gb2YgYWxsIFB1cGFlIiwpICsNCiAgZ2d0aXRsZSgiQXZlcmFnZSBvZiBBbGwgKGRlYWQgYW5kIGFsaXZlKSBQdXBhZSBwZXIgVHJlYXRtZW50IikgKw0KICBzY2FsZV94X2Rpc2NyZXRlKG5hbWUgPSAiVHJlYXRtZW50IiwgDQogICAgICAgICAgICAgICAgICAgbGFiZWxzID0gYygiMCBQUEIiLCAiMTUwIFBQQiIsICIxLDUwMCBQUEIiLCAiMTUsMDAwIFBQQiIsICIxNTAsMDAwIFBQQiIpKSArDQogIHRoZW1lX2Nvd3Bsb3QoKSArDQogIHRoZW1lX2NsYXNzaWMoYmFzZV9zaXplID0gMTIpDQoNCmBgYA0KDQoNCg0KIyNjYmluZCBwdXBhZSANCg0KYGBge3J9DQpwdXBhZS5tb3J0YWxpdHkgPC0gZ2xtKGNiaW5kKGxpdmVfcHVwYWUsIGRlYWRfcHVwYWUpIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiAgKyBxcm8sIGRhdGEgPSBicm9vZCwgZmFtaWx5ID0gYmlub21pYWwoImxvZ2l0IikpDQpzdW1tYXJ5KHB1cGFlLm1vcnRhbGl0eSkNCg0KQW5vdmEocHVwYWUubW9ydGFsaXR5KQ0KDQpwbG90KHB1cGFlLm1vcnRhbGl0eSkNCg0KbG0uZW1tIDwtIGVtbWVhbnMocHVwYWUubW9ydGFsaXR5LCAidHJlYXRtZW50IikNCnBhaXJzKGxtLmVtbSkNCnN1bW1hcnkobG0uZW1tKQ0KDQpwbG90KGJyb29kJHRyZWF0bWVudCwgY2JpbmQoYnJvb2QkZGVhZF9wdXBhZSwgYnJvb2QkbGl2ZV9wdXBhZSkpDQoNCmVtbWVhbnMocHVwYWUubW9ydGFsaXR5LCBwYWlyd2lzZSB+IHRyZWF0bWVudCwgdHlwZSA9ICJyZXNwb25zZSIpDQoNCmVtbTIgPC0gZW1tZWFucyhwdXBhZS5tb3J0YWxpdHksIH4gdHJlYXRtZW50LCB0eXBlID0gInJlc3BvbnNlIikNCmVtbTINCnBhaXJzKGVtbTIpDQoNCmVtbWRmMiA8LSBhcy5kYXRhLmZyYW1lKGVtbTIpDQoNCnAgPC0gZ2dwbG90KGRhdGEgPSBlbW1kZjIsIGFlcyh4PXRyZWF0bWVudCwgeT1wcm9iLCBmaWxsPXRyZWF0bWVudCkpICsgDQogIGdlb21fY29sKHBvc2l0aW9uID0gImRvZGdlIiwgY29sb3IgPSAiYmxhY2siKSArDQogIGNvb3JkX2NhcnRlc2lhbih5bGltID0gYygwLDEpKQ0KDQpwICsgZ2VvbV9lcnJvcmJhcihwb3NpdGlvbj1wb3NpdGlvbl9kb2RnZSguOSksd2lkdGg9LjI1LCBhZXMoeW1heD1hc3ltcC5VQ0wsIHltaW49YXN5bXAuTENMKSxhbHBoYT0uNikrDQogIGxhYnMoeD0iVHJlYXRtZW50IiwgeT0iUHJvYmFiaWxpdHkgb2YgU3Vydml2YWwiLCBmaWxsID0gInRyZWF0bWVudCIpICsgdGhlbWVfYncoKSArDQogIHNjYWxlX2ZpbGxfbWFudWFsKHZhbHVlcyA9IGMoInBlYWNocHVmZjMiLCAiZGFya3NlYWdyZWVuIiwgImxpZ2h0c2VhZ3JlZW4iLCAiZGFya29saXZlZ3JlZW4zIiwgImRhcmtvbGl2ZWdyZWVuNCIpKSArDQogICBsYWJzKHkgPSAiUHJvYmFiaWxpdHkgb2YgU3Vydml2YWwiLCkgKw0KICBnZ3RpdGxlKCJQcm9iYWJpbGl0eSBvZiBwdXBhZSBiZWluZyBhbGl2ZSB1cG9uIGNvbG9ueSBkaXNzZWN0aW9uIikgKw0KICBzY2FsZV94X2Rpc2NyZXRlKG5hbWUgPSAiVHJlYXRtZW50IiwgDQogICAgICAgICAgICAgICAgICAgbGFiZWxzID0gYygiMCBQUEIiLCAiMTUwIFBQQiIsICIxLDUwMCBQUEIiLCAiMTUsMDAwIFBQQiIsICIxNTAsMDAwIFBQQiIpKSArDQogIHRoZW1lX2Nvd3Bsb3QoKSArDQogIHRoZW1lX2NsYXNzaWMoYmFzZV9zaXplID0gMTIpDQpgYGANCg0KIyMjIERlYWQgTGFydmFlICsgRGVhZCBQdXBhZQ0KDQpgYGB7cn0NCmJyb29kJGFsbF9kX2xwIDwtIGJyb29kJGRlYWRfbGFydmFlICsgYnJvb2QkZGVhZF9wdXBhZQ0KDQpkbHAuZ24gPC0gZ2xtKGFsbF9kX2xwIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kLCBmYW1pbHkgPSAicG9pc3NvbiIpDQpzdW1tYXJ5KGRscC5nbikNCmRscC5nbmIgPC0gZ2xtLm5iKGFsbF9kX2xwIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kKQ0Kc3VtbWFyeShkbHAuZ25iKQ0KQUlDKGRscC5nbiwgZGxwLmduYikNCg0KQW5vdmEoZGxwLmduYikNCg0KcGxvdChicm9vZCR0cmVhdG1lbnQsIGJyb29kJGFsbF9kX2xwKQ0KDQoNCmRscC5zdW0gPC0gYnJvb2QgJT4lDQogIGdyb3VwX2J5KHRyZWF0bWVudCkgJT4lDQogIHN1bW1hcmlzZShhbC5tPSBtZWFuKGFsbF9kX2xwKSwgDQogICAgICAgICAgICBhbC5uID0gbGVuZ3RoKGFsbF9kX2xwKSwgDQogICAgICAgICAgICBzZC5hbCA9IHNkKGFsbF9kX2xwKSkNCg0KZGxwLnN1bQ0KDQpgYGANCg0KIyMjIExpdmUgTGFydmFlICsgTGl2ZSBQdXBhZQ0KDQpgYGB7cn0NCmJyb29kJGFsbF9hX2xwIDwtIGJyb29kJGxpdmVfbGFydmFlICsgYnJvb2QkbGl2ZV9wdXBhZQ0KDQphbHAuZ24gPC0gZ2xtKGFsbF9hX2xwIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kLCBmYW1pbHkgPSAicG9pc3NvbiIpDQpzdW1tYXJ5KGFscC5nbikNCmFscC5nbmIgPC0gZ2xtLm5iKGFsbF9hX2xwIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kKQ0Kc3VtbWFyeShhbHAuZ25iKQ0KQUlDKGFscC5nbiwgYWxwLmduYikNCg0KQW5vdmEoYWxwLmduYikNCg0KcGxvdChicm9vZCR0cmVhdG1lbnQsIGJyb29kJGFsbF9hX2xwKQ0KDQoNCmFscC5zdW0gPC0gYnJvb2QgJT4lDQogIGdyb3VwX2J5KHRyZWF0bWVudCkgJT4lDQogIHN1bW1hcmlzZShhbC5tPSBtZWFuKGFsbF9hX2xwKSwgDQogICAgICAgICAgICBhbC5uID0gbGVuZ3RoKGFsbF9hX2xwKSwgDQogICAgICAgICAgICBzZC5hbCA9IHNkKGFsbF9hX2xwKSkNCg0KYWxwLnN1bQ0KYGBgDQoNCg0KIyMjIEFsbCBMYXJ2YWUgKyBBbGwgUHVwYWUNCg0KYGBge3J9DQoNCmJyb29kJGFsbF9scCA8LSBicm9vZCRhbGxfbGFydmFlICsgYnJvb2QkYWxsX3B1cGFlDQoNCmxwLmduIDwtIGdsbShhbGxfbHAgfiB0cmVhdG1lbnQgKyB3aG9sZS5tZWFuICsgYWxpdmUgICsgcXJvLCBkYXRhID0gYnJvb2QsIGZhbWlseSA9ICJwb2lzc29uIikNCnN1bW1hcnkobHAuZ24pDQpscC5nbmIgPC0gZ2xtLm5iKGFsbF9scCB+IHRyZWF0bWVudCArIHdob2xlLm1lYW4gKyBhbGl2ZSAgKyBxcm8sIGRhdGEgPSBicm9vZCkNCnN1bW1hcnkobHAuZ25iKQ0KQUlDKGxwLmduLCBscC5nbmIpDQoNCkFub3ZhKGxwLmduYikNCg0KcGxvdChicm9vZCR0cmVhdG1lbnQsIGJyb29kJGFsbF9scCkNCg0KDQpscC5zdW0gPC0gYnJvb2QgJT4lDQogIGdyb3VwX2J5KHRyZWF0bWVudCkgJT4lDQogIHN1bW1hcmlzZShscC5tPSBtZWFuKGFsbF9scCksIA0KICAgICAgICAgICAgbHAubiA9IGxlbmd0aChhbGxfbHApLCANCiAgICAgICAgICAgIGxwLmFsID0gc2QoYWxsX2xwKSkNCg0KbHAuc3VtDQoNCmVtLmxwIDwtIGVtbWVhbnMobHAuZ25iLCAidHJlYXRtZW50IikNCnBhaXJzKGVtLmxwKQ0KYGBgDQoNCiMjIyNjYmluZCBhbGwgbGFydmFlIGFuZCBhbGwgcHVwYWUgDQoNCg0KYGBge3IsIGZpZy53aWR0aD0xMH0NCnBsX2JpbmQgPC0gZ2xtKGNiaW5kKGFsbF9hX2xwLCBhbGxfZF9scCkgfiB0cmVhdG1lbnQgKyB3aG9sZS5tZWFuICsgcXJvLCBkYXRhID0gYnJvb2QsIGZhbWlseSA9IGJpbm9taWFsKCJsb2dpdCIpKQ0Kc3VtbWFyeShwbF9iaW5kICkNCg0KQW5vdmEocGxfYmluZCApDQoNCnBsb3QocGxfYmluZCApDQoNCmxtLmVtbSA8LSBlbW1lYW5zKHBsX2JpbmQgLCAidHJlYXRtZW50IikNCnBhaXJzKGxtLmVtbSkNCnN1bW1hcnkobG0uZW1tKQ0KDQpwbG90KGJyb29kJHRyZWF0bWVudCwgY2JpbmQoYnJvb2QkZGVhZF9wdXBhZSwgYnJvb2QkbGl2ZV9wdXBhZSkpDQoNCmVtbWVhbnMocGxfYmluZCAsIHBhaXJ3aXNlIH4gdHJlYXRtZW50LCB0eXBlID0gInJlc3BvbnNlIikNCg0KZW1tMyA8LSBlbW1lYW5zKHBsX2JpbmQgLCB+IHRyZWF0bWVudCwgdHlwZSA9ICJyZXNwb25zZSIpDQplbW0zDQpwYWlycyhlbW0zKQ0KDQplbW1kZjMgPC0gYXMuZGF0YS5mcmFtZShlbW0zKQ0KDQpwIDwtIGdncGxvdChkYXRhID0gZW1tZGYzLCBhZXMoeD10cmVhdG1lbnQsIHk9cHJvYiwgZmlsbD10cmVhdG1lbnQpKSArIA0KICBnZW9tX2NvbChwb3NpdGlvbiA9ICJkb2RnZSIsIGNvbG9yID0gImJsYWNrIikgKw0KICBjb29yZF9jYXJ0ZXNpYW4oeWxpbSA9IGMoMCwxKSkNCg0KcCArIGdlb21fZXJyb3JiYXIocG9zaXRpb249cG9zaXRpb25fZG9kZ2UoLjkpLHdpZHRoPS4yNSwgYWVzKHltYXg9YXN5bXAuVUNMLCB5bWluPWFzeW1wLkxDTCksYWxwaGE9LjYpKw0KICBsYWJzKHg9IlRyZWF0bWVudCIsIHk9IlByb2JhYmlsaXR5IG9mIFN1cnZpdmFsIiwgZmlsbCA9ICJ0cmVhdG1lbnQiKSArIHRoZW1lX2J3KCkgKw0KICBzY2FsZV9maWxsX21hbnVhbCh2YWx1ZXMgPSBjKCJwZWFjaHB1ZmYzIiwgImRhcmtzZWFncmVlbiIsICJsaWdodHNlYWdyZWVuIiwgImRhcmtvbGl2ZWdyZWVuMyIsICJkYXJrb2xpdmVncmVlbjQiKSkgKw0KICAgbGFicyh5ID0gIlByb2JhYmlsaXR5IG9mIFN1cnZpdmFsIiwpICsNCiAgZ2d0aXRsZSgiQ29tYmluZWQgTGFydmFlIGFuZCBQdXBhZSBwcm9iYWJpbGl0eSBvZiBiZWluZyBhbGl2ZSB1cG9uIGNvbG9ueSBkaXNzZWN0aW9uIikgKw0KICBzY2FsZV94X2Rpc2NyZXRlKG5hbWUgPSAiVHJlYXRtZW50IiwgDQogICAgICAgICAgICAgICAgICAgbGFiZWxzID0gYygiMCBQUEIiLCAiMTUwIFBQQiIsICIxLDUwMCBQUEIiLCAiMTUsMDAwIFBQQiIsICIxNTAsMDAwIFBQQiIpKSArDQogIHRoZW1lX2Nvd3Bsb3QoKSArDQogIHRoZW1lX2NsYXNzaWMoYmFzZV9zaXplID0gMTIpDQpgYGANCg0KDQoNCiMjIyBFZ2dzDQoNCmBgYHtyfQ0KDQplZ2dzLmduIDwtIGdsbShlZ2dzIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kLCBmYW1pbHkgPSAicG9pc3NvbiIpDQpzdW1tYXJ5KGVnZ3MuZ24pDQojZWdncy5nbmIgPC0gZ2xtLm5iKGVnZ3MgfiB0cmVhdG1lbnQgKyB3aG9sZS5tZWFuICsgYWxpdmUgICsgcXJvLCBkYXRhID0gYnJvb2QpDQojc3VtbWFyeShlZ2dzLmduYikNCiNBSUMoZWdncy5nbiwgZWdncy5nbmIpICAjIyMjIEkgdGhpbmsgd2UgaGF2ZSB0b28gbWFueSB6ZXJvcyBpbiBlZ2dzLCBJIHdpbGwgc3dpdGNoIHRvIGEgemVyby1pbmZsYXRlZCBuZWdhdGl2ZSBiaW5vbWlhbCBtb2RlbA0KDQoNCmxpYnJhcnkoZ2xtbVRNQikNCmVnZy56ZXJvIDwtIGdsbW1UTUIoZWdncyB+IHRyZWF0bWVudCArIHdob2xlLm1lYW4gICsgcXJvLCBmYW1pbHkgPSAibmJpbm9tMiIsIGJyb29kKSAgIyBoYWQgdG8gdGFrZSBvdXQgYWxpdmUsIGJlY2F1c2UgaXQgbWFkZSBtb2RlbCBub3QgY29udmVyZ2UgLS0gQW5vdmEgZnJvbSBnbG0gc2hvd3MgaXQncyBub3QgaW1wYWN0ZnVsIGFueXdhcyANCnN1bW1hcnkoZWdnLnplcm8pDQpBbm92YShlZ2cuemVybykNCg0KQUlDKGVnZ3MuZ24sIGVnZy56ZXJvKQ0KDQoNCnBsb3QoYnJvb2QkdHJlYXRtZW50LCBicm9vZCRlZ2dzKQ0KDQoNCmVnZ3Muc3VtIDwtIGJyb29kICU+JQ0KICBncm91cF9ieSh0cmVhdG1lbnQpICU+JQ0KICBzdW1tYXJpc2UoZWdncy5tPSBtZWFuKGVnZ3MpLCANCiAgICAgICAgICAgIGVnZ3MubiA9IGxlbmd0aChlZ2dzKSwgDQogICAgICAgICAgICBlZ2dzLmFsID0gc2QoZWdncykpDQoNCmVnZ3Muc3VtDQoNCmVtLmVnZyA8LSBlbW1lYW5zKGVnZy56ZXJvLCAidHJlYXRtZW50IikNCnBhaXJzKGVtLmVnZykNCnN1bW1hcnkoZW0uZWdnKQ0KDQpgYGANCg0KIyMjIENvdW50IG9mIGRyb25lcw0KDQpgYGB7cn0NCg0KZHJvbmVzLmduIDwtIGdsbShkcm9uZXMgfiB0cmVhdG1lbnQgKyB3aG9sZS5tZWFuICsgYWxpdmUgICsgcXJvLCBkYXRhID0gYnJvb2QsIGZhbWlseSA9ICJwb2lzc29uIikNCnN1bW1hcnkoZHJvbmVzLmduKQ0KQW5vdmEoZHJvbmVzLmduKQ0KZHJvbmVzLmduYiA8LSBnbG0ubmIoZHJvbmVzIH4gdHJlYXRtZW50ICsgd2hvbGUubWVhbiArIGFsaXZlICArIHFybywgZGF0YSA9IGJyb29kKQ0Kc3VtbWFyeShkcm9uZXMuZ25iKQ0KQUlDKGRyb25lcy5nbiwgZHJvbmVzLmduYikNCg0KQW5vdmEoZHJvbmVzLmduYikNCg0KZHJvbmUuZW0gPC0gZW1tZWFucyhkcm9uZXMuZ25iLCAidHJlYXRtZW50IikNCnBhaXJzKGRyb25lLmVtKQ0KDQpkcm9uZXMuc3VtIDwtIGJyb29kICU+JQ0KICBncm91cF9ieSh0cmVhdG1lbnQpICU+JQ0Kc3VtbWFyaXNlKG0uZCA9IG1lYW4oZHJvbmVzKSwgDQogICAgICAgICAgc2QuZCA9IHNkKGRyb25lcyksDQogICAgICAgICAgbi5kID0gbGVuZ3RoKGRyb25lcykpICU+JQ0KICBtdXRhdGUoc2VkID0gc2QuZCAvIHNxcnQobi5kKSkNCg0KZHJvbmVzLnN1bQ0KDQpgYGANCmBgYHtyfQ0KDQpnZ3Bsb3QoZGF0YSA9IGRyb25lcy5zdW0sIGFlcyh4ID0gdHJlYXRtZW50LCB5ID0gbS5kLCBmaWxsID0gdHJlYXRtZW50KSkgKw0KICBnZW9tX2NvbChwb3NpdGlvbiA9ICJkb2RnZSIsIGNvbCA9ICJibGFjayIpKw0KICBjb29yZF9jYXJ0ZXNpYW4oeWxpbT1jKDAsMTYpKSArDQogIHNjYWxlX2ZpbGxfbWFudWFsKHZhbHVlcyA9IGMoInBlYWNocHVmZjMiLCAiZGFya3NlYWdyZWVuIiwgImxpZ2h0c2VhZ3JlZW4iLCAiZGFya29saXZlZ3JlZW4zIiwgImRhcmtvbGl2ZWdyZWVuNCIpLA0KICAgICAgICAgICAgICAgICAgICBuYW1lID0gIlByaXN0aW5lIExldmVsIiwNCiAgICAgICAgICAgICAgICAgICAgbGFiZWxzID0gYygiVHJlYXRtZW50IDEgKGNvbnRyb2wpIiwgIlRyZWF0bWVudCAyIiwgDQogICAgICAgICAgICAgICAgICAgICAgICAgICAgICAgIlRyZWF0bWVudCAzIiwgIlRyZWF0bWVudCA0IiwgIlRyZWF0bWVudCA1IikpICsgDQogIGdlb21fZXJyb3JiYXIoYWVzKHltaW4gPSBtLmQgLSBzZWQsIA0KICAgICAgICAgICAgICAgICAgICB5bWF4ID0gbS5kICsgc2VkKSwNCiAgICAgICAgICAgICAgICBwb3NpdGlvbiA9IHBvc2l0aW9uX2RvZGdlKDAuOSksIHdpZHRoID0gMC40KSArDQogIGxhYnMoeSA9ICJNZWFuIE51bWJlciBvZiBEcm9uZXMiLCkgKw0KICBnZ3RpdGxlKCJBdmVyYWdlIG9mIERyb25lcyBDb3VudCBwZXIgVHJlYXRtZW50IikgKw0KICBzY2FsZV94X2Rpc2NyZXRlKG5hbWUgPSAiVHJlYXRtZW50IiwgDQogICAgICAgICAgICAgICAgICAgbGFiZWxzID0gYygiMCBQUEIiLCAiMTUwIFBQQiIsICIxLDUwMCBQUEIiLCAiMTUsMDAwIFBQQiIsICIxNTAsMDAwIFBQQiIpKSArDQogIHRoZW1lX2Nvd3Bsb3QoKSArDQogIHRoZW1lX2NsYXNzaWMoYmFzZV9zaXplID0gMTIpDQoNCmBgYA0KDQo=