library(tidyverse) library(car) library(lme4) library(ggplot2) library(readxl) library(emmeans) tim <- read.csv("Data/17 sex ratios.csv") summary(tim) tim$Clutch<-as.factor(tim$Clutch) Sex<-cbind(tim$Successes,tim$Failures) Sex model<-glmer(Sex~Treatment+(1|Clutch), data = tim, family = "binomial") #dropping clutch as random effect# model<-glm(Sex~Treatment, data = tim, family = "binomial") model<-glmer(cbind(tim$Successes,tim$Failures)~Treatment+(1|Clutch), data = tim, family = "binomial") summary(model) drop1(model, test = "Chisq") Anova(model) t.emm <- emmeans(model, "Treatment", type="response") t.emm emmeans(model, list(pairwise ~ Treatment), adjust = "tukey") sexratio<-summary(t.emm, infer = TRUE, type = "response") sexratio sexratio %>% ggplot(aes(x=Treatment, y = prob)) + geom_bar(stat = "identity", color="black", width = .45, fill = "white") + geom_errorbar(aes(ymin=asymp.LCL, ymax=asymp.UCL), width=.2, position = position_dodge(.9)) + theme(panel.grid.major = element_line(linetype = "blank"), panel.grid.minor = element_line(linetype = "blank"), axis.text = element_text(colour = "black"), axis.text.x = element_text(angle = -40), panel.background = element_rect(fill = "white"), legend.position = "none") +labs(x = "Heat Wave Treatment", y = "Probability of female production ± 95% C.I.") + theme(axis.line = element_line(linetype = "solid"), axis.text = element_text(vjust = 0), axis.text.x = element_text(vjust = 0)) + annotate("text", x = "Day 10-24", y = .4, label = "a") + annotate("text", x = "Day 38-52", y = .4, label = "a") + annotate("text", x = "Day 17-31", y = .95, label = "b") + annotate("text", x = "Day 24-38", y = 1.05, label = "b") + annotate("text", x = "Day 31-45", y = 1.025, label = "b") + theme(axis.text = element_text(hjust = 0.25)) + theme(axis.title = element_text(family = "Times"), axis.text = element_text(family = "Times"), axis.text.x = element_text(family = "Times"), axis.text.y = element_text(family = "Times"), plot.title = element_text(family = "Times"), strip.text = element_text(family = "Times")) + theme(axis.title = element_text(size = 12), axis.text = element_text(size = 12)) + theme(axis.text.x = element_text(vjust = 0.5)) + theme(axis.title = element_text(vjust = 0), axis.text = element_text(hjust = 0)) + theme(axis.text.x = element_text(vjust = 1.25, angle = -45))