TRES Color and Stress 2018-19
Conor Taff
Last Updated: September 10 2020 at 11:13
Overall goal
Is there an interraction between social environment experienced and the susceptibility to challenges? We did two different factorial experiments that coupled manipulations targeted at changing the social environment (signal color manipulation) with manipulations targeted at imposing realistic environmental challenges (flight handicap or simulated predations).
In 2018, manipulations were ordered so that birds were first dulled and then later in the season experienced one of the challenges.
In 2019, the order was reversed so that birds experienced the challenge (only predator, no handicaping in this year) and then subsequently received a dulling manipulation.
For each year, we can ask whether dulling or challenges influenced a number of female measures (mass loss, cort response, provisioning rate, visitation & trips) and nestling outcomes (number fledged, nestling mass, nestling fledge age). Note that in these years we also cross fostered nestlings so that effects on nestlings can be disentangled from those caused by maternal effects vs. experienced nestling environment. Not sure what level of detail we want to go into with nestlings.
Data that we have for these birds but that at the moment I am thinking would not be included in this paper: adult/nestling microbiome, full paternity assignment, nestling/adult feather quality, glucose, some diet.
In this first code chunk I'm loading packages and data to use in later plotting.
# Load packages ----
pacman::p_load(tidyverse, here, lme4, scales)
## not loading plyr because it conflicts with here and tidyverse
# Load data ----
# the data files include similar data from each year but experiments will mainly be analyzed separately
# data by female: one row per female with measures throughout the breeding season
d_fem <- read.delim(here::here("1_raw_data", "data_by_female.txt"))
# Remove records that were noted as exclude. These are mostly cases where either a nest failed so early
# that no treatment was delivered or something weird happened (e.g., treatments mixed up).
d_fem <- subset(d_fem, d_fem$exclude != "Yes")
# data by nestling: one row per nestling with some more individual level nestling info
d_nestling <- read.delim(here::here("1_raw_data", "data_by_nestling.txt"))
# provisioning data: one row per nest per day with provisiong info from rfid sensors
d_provision <- read.delim(here::here("1_raw_data", "daily_provision_data.txt"))
# visitation data: one row per nest per day with social visits and trips for each nest
d_social <- read.delim(here::here("1_raw_data", "nest_visit_data.txt"))
# rfid reads: total rfid reads at each box on each day (used to filter out malfunctioning rfid boards)
d_tot_rfid <- read.delim(here::here("1_raw_data", "total_rfid_reads.txt"))
# Set colors for different treatment groups ----
col_dull <- "slateblue"
col_sham <- "orange"
ch_pred <- "chartreuse3"
ch_con <- "orange"
ch_tape <- "coral3"
# Full treatments
col_CC <- "gray50" # 2018 & 2019, full control
col_CS <- "coral3" # 2018, sham color then taping challenge
col_CP <- "orange" # 2018, sham color then predator
col_DC <- "slateblue" # 2018, dulling folowed by sham coloring
col_DS <- "violet" # 2018, dulling folowed by taping challenge
col_DP <- "chartreuse3" # 2018, dulling followed by predator
col_PC <- "coral3" # 2019, predator first then sham color
col_PD <- "chartreuse3" # 2019, predator first then dulling
col_CD <- "slateblue" # 2019, control first, then dullingFemale Mass Loss
Females lose mass through the breeding season normally. The question here is whether the treatments we applied changed this mass loss trajectory.
## Female Mass Change ----
# Figures showing adult mass at each timepoint by treatment split by year (experiment)
long_mass <- d_fem %>%
pivot_longer(cols = c("fmass1", "fmass2", "fmass3"), names_to = "cap_num",
names_transform = list(
cap_num = ~ readr::parse_factor(.x, levels = c("fmass1", "fmass2", "fmass3"),
ordered = TRUE)),
values_to = "mass", values_drop_na = TRUE)
long_mass <- as.data.frame(long_mass)
long_mass$xpos <- as.numeric(long_mass$cap_num)
sum_mass <- d_fem %>%
pivot_longer(cols = c("fmass1", "fmass2", "fmass3"), names_to = "cap_num",
names_transform = list(
cap_num = ~ readr::parse_factor(.x, levels = c("fmass1", "fmass2", "fmass3"),
ordered = TRUE)),
values_to = "mass", values_drop_na = TRUE) %>%
group_by(as.factor(year), full_treatment, cap_num) %>%
summarise(n = n(), mu = mean(mass, na.rm = TRUE), se = sd(mass, na.rm = TRUE) / sqrt(n()))
sum_mass$x_pos <- as.numeric(sum_mass$cap_num)
for_plot <- data.frame(full_treatment = unique(sum_mass$full_treatment),
p_shape = c(21, 24, 22, 21, 24, 22, 21, 24, 24),
p_col = c(rep(col_sham, 3), rep(col_dull, 3), col_dull, col_sham, col_dull),
p_dodge = c(-.25, -.15, -.05, .05, .15, .25, .0833, -.0833, .25))
sum_mass2 <- as.data.frame(sum_mass)
colnames(sum_mass2)[1] <- "year"
sum_mass2 <- plyr::join(sum_mass2, for_plot, "full_treatment", "left", "first")
long_mass2 <- plyr::join(long_mass, for_plot, "full_treatment", "left", "first")
lm2_18 <- subset(long_mass2, long_mass2$year == "2018")
lm2_19 <- subset(long_mass2, long_mass2$year == "2019")
sm2_18 <- subset(sum_mass2, sum_mass2$year == "2018")
sm2_19 <- subset(sum_mass2, sum_mass2$year == "2019")
par(mfrow = c(1, 2))
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Dulling then Challenge", xlim = c(0.3, 3.7), ylim = c(15.5, 25),
xlab = "Capture Number", ylab = "Female Mass (g)")
axis(1, seq(0, 4, 1))
axis(2, seq(5, 30, 1), las = 2)
points(lm2_18$xpos + lm2_18$p_dodge, lm2_18$mass, pch = 16, col = alpha("black", 0.3))
for(i in 1:nrow(sm2_18)){
lines(rep(sm2_18$x_pos[i] + sm2_18$p_dodge[i], 2),
c(sm2_18$mu[i] - sm2_18$se[i], sm2_18$mu[i] + sm2_18$se[i]), lwd = 2)
}
points(sm2_18$x_pos + sm2_18$p_dodge, sm2_18$mu, pch = sm2_18$p_shape, cex = 1.4,
bg = as.character(sm2_18$p_col))
#legend("topright", c("No Color", "Dull Color", "Control", "Predator", "Handicap"),
# pch = c(21, 21, 21, 24, 22), pt.bg = c(col_sham, col_dull, rep("white", 3)), cex = 1.2, bty = "n")
rect(1.4, 15.7, 3.5, 16.2, col = alpha("chartreuse3", 0.3))
text(2.15, 15.95, "Dulling", pos = 4)
rect(2.2, 16.3, 3.5, 16.8, col = alpha("coral3", 0.2))
text(2.4, 16.55, "Challenge", pos = 4)
#2019
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Challenge then Dulling", xlim = c(0.3, 3.7), ylim = c(15.5, 25),
xlab = "Capture Number", ylab = "Female Mass (g)")
axis(1, seq(0, 4, 1))
axis(2, seq(5, 30, 1), las = 2)
points(lm2_19$xpos + lm2_19$p_dodge, lm2_19$mass, pch = 16, col = alpha("black", 0.3))
for(i in 1:nrow(sm2_19)){
lines(rep(sm2_19$x_pos[i] + sm2_19$p_dodge[i], 2),
c(sm2_19$mu[i] - sm2_19$se[i], sm2_19$mu[i] + sm2_19$se[i]), lwd = 2)
}
points(sm2_19$x_pos + sm2_19$p_dodge, sm2_19$mu, pch = sm2_19$p_shape, cex = 1.4,
bg = as.character(sm2_19$p_col))
legend("topright", c("Sham", "Dulled", "Control", "Predator", "Handicap"),
pch = c(21, 21, 21, 24, 22), pt.bg = c(col_sham, col_dull, rep("white", 3)), cex = 1, bty = "n")
rect(2.3, 15.7, 3.6, 16.2, col = alpha("chartreuse3", 0.3))
text(2.65, 15.95, "Dulling", pos = 4)
rect(1.4, 16.3, 2.5, 16.8, col = alpha("coral3", 0.2))
text(1.5, 16.55, "Challenge", pos = 4)
So the signature of mass loss across captures is very clear in all groups, but there is basically no suggestion that any of our treatments are resulting in a steaper mass loss than any others.
Female Cort Responses
## Female Cort
long_cort <- d_fem %>%
pivot_longer(cols = c("fbase1", "fstr1", "fdex1", "fbase2", "fbase3", "fstr3", "fdex3"), names_to = "cort_type",
names_transform = list(
cap_num = ~ readr::parse_factor(.x, levels = c("fbase1", "fstr1", "fdex1", "fbase2", "fbase3", "fstr3", "fdex3"),
ordered = TRUE)),
values_to = "cort", values_drop_na = TRUE)
long_cort <- as.data.frame(long_cort)
cort_pos <- data.frame(cort_type = c("fbase1", "fstr1", "fdex1", "fbase2", "fbase3", "fstr3", "fdex3"),
xpos = c(1, 2, 3, 5, 7, 8, 9))
long_cort <- plyr::join(long_cort, cort_pos, "cort_type", "left", "first")
sum_cort <- d_fem %>%
pivot_longer(cols = c("fbase1", "fstr1", "fdex1", "fbase2", "fbase3", "fstr3", "fdex3"), names_to = "cort_type",
names_transform = list(
cap_num = ~ readr::parse_factor(.x, levels = c("fbase1", "fstr1", "fdex1", "fbase2", "fbase3", "fstr3", "fdex3"),
ordered = TRUE)),
values_to = "cort", values_drop_na = TRUE) %>%
group_by(as.factor(year), full_treatment, cort_type) %>%
summarise(n = n(), mu = mean(cort, na.rm = TRUE), se = sd(cort, na.rm = TRUE) / sqrt(n()))
sum_cort <- as.data.frame(sum_cort)
sum_cort <- plyr::join(sum_cort, cort_pos, "cort_type", "left", "first")
for_plot <- data.frame(full_treatment = unique(sum_cort$full_treatment),
p_shape = c(21, 24, 22, 21, 24, 22, 21, 24, 24),
p_col = c(rep(col_sham, 3), rep(col_dull, 3), col_dull, col_sham, col_dull),
p_dodge = c(-.25, -.15, -.05, .05, .15, .25, .0833, -.0833, .25))
sum_cort <- plyr::join(sum_cort, for_plot, "full_treatment")
long_cort <- plyr::join(long_cort, for_plot, "full_treatment")
colnames(sum_cort)[2] <- "year"
lc2_18 <- subset(long_cort, long_cort$year == "2018")
lc2_19 <- subset(long_cort, long_cort$year == "2019")
sc2_18 <- subset(sum_cort, sum_cort$year == "2018")
sc2_19 <- subset(sum_cort, sum_cort$year == "2019")
#2018 plot
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Dulling then Challenge", xlim = c(0.3, 9.7), ylim = c(-20, 80),
xlab = "", ylab = "Corticosterone (ng/ul)")
axis(1, c(-1, 1, 2, 3, 5, 7, 8, 9, 15), c("", "B 1", "S 1", "D 1", "B 2", "B 3", "S 3", "D 3", ""))
axis(2, c(-50, seq(0, 200, 10)), las = 2)
points(lc2_18$xpos + lc2_18$p_dodge, lc2_18$cort, pch = 16, col = alpha("black", 0.3))
abline(h = 0)
for(i in 1:nrow(sc2_18)){
lines(rep(sc2_18$xpos[i] + sc2_18$p_dodge[i], 2),
c(sc2_18$mu[i] - sc2_18$se[i], sc2_18$mu[i] + sc2_18$se[i]), lwd = 2)
}
points(sc2_18$xpos + sc2_18$p_dodge, sc2_18$mu, pch = sc2_18$p_shape, cex = 1.4,
bg = as.character(sc2_18$p_col))
rect(4.4, -15, 9.6, -10, col = alpha("chartreuse3", 0.3))
text(6.5, -12.5, "Dulling", pos = 4)
rect(6.4, -9, 9.6, -4, col = alpha("coral3", 0.2))
text(7.3, -6.5, "Challenge", pos = 4)
legend("bottomleft", c("No Color", "Dull Color", "Control", "Predator", "Handicap"),
pch = c(21, 21, 21, 24, 22), pt.bg = c(col_sham, col_dull, rep("white", 3)), cex =0.9, bty = "n")
#2019 plot
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Challenge then Dulling", xlim = c(0.3, 9.7), ylim = c(-20, 80),
xlab = "", ylab = "Corticosterone (ng/ul)")
axis(1, c(-1, 1, 2, 3, 5, 7, 8, 9, 15), c("", "B 1", "S 1", "D 1", "B 2", "B 3", "S 3", "A 3", ""))
axis(2, c(-50, seq(0, 200, 10)), las = 2)
points(lc2_19$xpos + lc2_19$p_dodge, lc2_19$cort, pch = 16, col = alpha("black", 0.3))
abline(h = 0)
for(i in 1:nrow(sc2_19)){
lines(rep(sc2_19$xpos[i] + sc2_19$p_dodge[i], 2),
c(sc2_19$mu[i] - sc2_19$se[i], sc2_19$mu[i] + sc2_19$se[i]), lwd = 2)
}
points(sc2_19$xpos + sc2_19$p_dodge, sc2_19$mu, pch = sc2_19$p_shape, cex = 1.4,
bg = as.character(sc2_19$p_col))
rect(6, -15, 9.6, -10, col = alpha("chartreuse3", 0.3))
text(7, -12.5, "Dulling", pos = 4)
rect(3.8, -9, 6.6, -4, col = alpha("coral3", 0.2))
text(4.5, -6.5, "Challenge", pos = 4)
legend("bottomleft", c("No Color", "Dull Color", "Control", "Predator"),
pch = c(21, 21, 21, 24), pt.bg = c(col_sham, col_dull, rep("white", 2)), cex =0.9, bty = "n")
There doesn't seem to be much going on with cort differences between groups. Maybe in the first year the predator groups have higher stress response and dex than the non-predator groups, but that certainly isn't very strong. In any case, there doesn't seem to be any interaction between coloring treatment and challenge treatment.
Female Provisioning
## Provisioning ----
# Figure showing nestling provisioning
d_fem$uby <- paste(d_fem$unitbox, d_fem$year, sep = "_")
d_prov2 <- plyr::join(d_provision, d_fem, "uby", "left", "first")
d_prov3 <- subset(d_prov2, d_prov2$offset < 16 & d_prov2$doy < d_prov2$end_soc_date & d_prov2$f_reads > 50)
d_prov3 <- d_prov3[, c("offset", "f_feed", "m_feed", "year", "full_treatment", "experiment")]
d_prov3_18 <- subset(d_prov3, d_prov3$year == 2018)
ggplot(d_prov3_18, aes(x = offset, y = f_feed, col = full_treatment)) +
geom_point() + geom_smooth() + xlim(0, 15) + xlab("Days after hatching") +
ylab("Daily female provisioning trips") + ggtitle("Dulling then Challenge") +
theme_bw()
ggplot(d_prov3_18, aes(x = offset, y = m_feed, col = full_treatment)) +
geom_point() + geom_smooth() + xlim(0, 15) + xlab("Days after hatching") +
ylab("Daily male provisioning trips") + ggtitle("Dulling then Challenge") +
theme_bw()
d_prov3_19 <- subset(d_prov3, d_prov3$year == 2019)
ggplot(d_prov3_19, aes(x = offset, y = f_feed, col = full_treatment)) +
geom_point() + geom_smooth() + xlim(0.6, 14.5) + xlab("Days after hatching") +
ylab("Daily female provisioning trips") + ggtitle("Challenge then Dulling") +
theme_bw()
ggplot(d_prov3_19, aes(x = offset, y = m_feed, col = full_treatment)) +
geom_point() + geom_smooth() + xlim(0.6, 14.5) + xlab("Days after hatching") +
ylab("Daily male provisioning trips") + ggtitle("Challenge then Dulling") +
theme_bw()
There isn't anything super obvious going on with provisioning rate in any of the treatment groups. Maybe some hints of lower overall rates in predator groups, but not much. I could model these in a better way looking at hourly feeding and controlling for some other things like temperature and time of day, etc.
- I need to fix the colors/style on these but it's just using the ggplot defaults for now.
Nest Box Visitors
## Visits to box
soc_vis <- data.frame(uby = unique(d_fem$uby))
for(i in 1:nrow(soc_vis)){
sub <- subset(d_fem, d_fem$uby == soc_vis$uby[i])
if(sub$end_soc_date[1] - sub$cap1[1] + 1 > 0){
dats <- seq(from = sub$cap1[1] + 1, to = sub$end_soc_date[1], by = 1)
temp <- data.frame(uby = rep(sub$uby, length(dats)),
doy = dats)
if(i == 1){
soc_visits <- temp
}
if(i > 1){
soc_visits <- rbind(soc_visits, temp)
}
}
}
soc_visits <- plyr::join(soc_visits, d_fem, "uby", "left", "first")
#for(i in 1:nrow(soc_visits)){
# sub1 <- subset(d_social, as.character(d_social$ubox) == as.character(soc_visits$unitbox[i])
# & d_social$doy == soc_visits$doy[i])
# subf <- subset(sub1, sub1$sex == "Female")
# subm <- subset(sub1, sub1$sex == "Male")
# soc_visits$tot_f_vis[i] <- nrow(subf)
# soc_visits$tot_m_vis[i] <- nrow(subm)
# soc_visits$uni_f_vis[i] <- length(unique(subf$rfid))
# soc_visits$uni_m_vis[i] <- length(unique(subm$rfid))
# print(paste(i, " of ", nrow(soc_visits), sep = ""))
#}
for(i in 1:nrow(d_fem)){
sub <- subset(d_social, as.character(d_social$ubox) == as.character(d_fem$unitbox[i]))
subf <- subset(sub, sub$sex == "Female")
subm <- subset(sub, sub$sex == "Male")
subf1 <- subset(subf, subf$doy >= d_fem$cap1[i] & subf$doy < d_fem$cap2[i])
subf2 <- subset(subf, subf$doy >= d_fem$cap2[i])
subm1 <- subset(subm, subm$doy >= d_fem$cap1[i] & subm$doy < d_fem$cap2[i])
subm2 <- subset(subm, subm$doy >= d_fem$cap2[i])
d_fem$tot_f_vis[i] <- nrow(subf)
d_fem$tot_m_vis[i] <- nrow(subm)
d_fem$uni_f_vis[i] <- length(unique(subf$rfid))
d_fem$uni_m_vis[i] <- length(unique(subm$rfid))
d_fem$tot_f_vis1[i] <- nrow(subf1)
d_fem$tot_m_vis1[i] <- nrow(subm1)
d_fem$uni_f_vis1[i] <- length(unique(subf1$rfid))
d_fem$uni_m_vis1[i] <- length(unique(subm1$rfid))
d_fem$tot_f_vis2[i] <- nrow(subf2)
d_fem$tot_m_vis2[i] <- nrow(subm2)
d_fem$uni_f_vis2[i] <- length(unique(subf2$rfid))
d_fem$uni_m_vis2[i] <- length(unique(subm2$rfid))
}
##visitors to the box
long_vis <- d_fem %>%
pivot_longer(cols = c("uni_f_vis1", "uni_m_vis1", "uni_f_vis2", "uni_m_vis2"), names_to = "type",
names_transform = list(
type = ~ readr::parse_factor(.x, levels = c("uni_f_vis1", "uni_m_vis1", "uni_f_vis2", "uni_m_vis2"),
ordered = TRUE)),
values_to = "count", values_drop_na = TRUE)
long_vis <- as.data.frame(long_vis)
vis_pos <- data.frame(type = c("uni_f_vis1", "uni_m_vis1", "uni_f_vis2", "uni_m_vis2"),
xpos = seq(1, 4, 1))
long_vis <- plyr::join(long_vis, vis_pos, "type", "left", "first")
sum_vis <- d_fem %>%
pivot_longer(cols = c("uni_f_vis1", "uni_m_vis1", "uni_f_vis2", "uni_m_vis2"), names_to = "type",
names_transform = list(
type = ~ readr::parse_factor(.x, levels = c("uni_f_vis1", "uni_m_vis1", "uni_f_vis2", "uni_m_vis2"),
ordered = TRUE)),
values_to = "count", values_drop_na = TRUE) %>%
group_by(as.factor(year), full_treatment, type) %>%
summarise(n = n(), mu = mean(count, na.rm = TRUE), se = sd(count, na.rm = TRUE) / sqrt(n()))
sum_vis <- as.data.frame(sum_vis)
sum_vis <- plyr::join(sum_vis, vis_pos, "type", "left", "first")
colnames(sum_vis)[2] <- "year"
for_plot <- data.frame(full_treatment = unique(sum_vis$full_treatment),
p_shape = c(21, 24, 22, 21, 24, 22, 21, 24, 24),
p_col = c(rep(col_sham, 3), rep(col_dull, 3), col_dull, col_sham, col_dull),
p_dodge = c(-.25, -.15, -.05, .05, .15, .25, .0833, -.0833, .25))
sum_vis <- plyr::join(sum_vis, for_plot, "full_treatment")
long_vis <- plyr::join(long_vis, for_plot, "full_treatment")
vis2_18 <- subset(long_vis, long_vis$year == "2018")
vis2_19 <- subset(long_vis, long_vis$year == "2019")
svis2_18 <- subset(sum_vis, sum_vis$year == "2018")
svis2_19 <- subset(sum_vis, sum_vis$year == "2019")
## Visitors
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Dulling then Challenge", xlim = c(0.5, 4.5), ylim = c(-2, 25),
xlab = "", ylab = "Number of Unique Visitors")
axis(1, c(-1, 1, 2, 3, 4, 10), c("", "Females 1", "Males 1", "Females 2", "Males 2", ""))
axis(2, seq(-30, 200, 5), las = 2)
abline(h = 0)
points(vis2_18$xpos + vis2_18$p_dodge, vis2_18$count, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(svis2_18)){
lines(rep(svis2_18$xpos[i] + svis2_18$p_dodge[i], 2),
c(svis2_18$mu[i] - svis2_18$se[i], svis2_18$mu[i] + svis2_18$se[i]), lwd = 2)
}
points(svis2_18$xpos + svis2_18$p_dodge, svis2_18$mu, pch = svis2_18$p_shape, cex = 1.4,
bg = as.character(svis2_18$p_col))
legend(.5, 23, c("No Color", "Dull Color", "Control", "Predator", "Handicap"),
pch = c(21, 21, 21, 24, 22), pt.bg = c(col_sham, col_dull, rep("white", 3)), cex =0.9, bty = "n")
rect(0.5, 23.5, 4.5, 25, col = alpha("chartreuse3", 0.3))
text(2.3, 24.25, "Dulling", pos = 4)
rect(2.5, 21.5, 4.5, 23, col = alpha("coral3", 0.2))
text(3.2, 22.25, "Challenge", pos = 4)
lines(c(2.5, 2.5), c(-5, 23.5), lty = 2)
text(1.5, -1, "First to second capture")
text(3.5, -1, "After second capture")
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Challenge then Dulling", xlim = c(0.5, 4.5), ylim = c(-2, 25),
xlab = "", ylab = "Number of Unique Visitors")
axis(1, c(-1, 1, 2, 3, 4, 10), c("", "Females 1", "Males 1", "Females 2", "Males 2", ""))
axis(2, seq(-30, 200, 5), las = 2)
abline(h = 0)
points(vis2_19$xpos + vis2_19$p_dodge, vis2_19$count, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(svis2_19)){
lines(rep(svis2_19$xpos[i] + svis2_19$p_dodge[i], 2),
c(svis2_19$mu[i] - svis2_19$se[i], svis2_19$mu[i] + svis2_19$se[i]), lwd = 2)
}
points(svis2_19$xpos + svis2_19$p_dodge, svis2_19$mu, pch = svis2_19$p_shape, cex = 1.4,
bg = as.character(svis2_19$p_col))
legend(.5, 22, c("No Color", "Dull Color", "Control", "Predator"),
pch = c(21, 21, 21, 24), pt.bg = c(col_sham, col_dull, rep("white", 2)), cex =0.9, bty = "n")
rect(2.5, 23.5, 4.5, 25, col = alpha("chartreuse3", 0.3))
text(3.3, 24.25, "Dulling", pos = 4)
rect(0.5, 21.5, 2.5, 23, col = alpha("coral3", 0.2))
text(1.2, 22.25, "Challenge", pos = 4)
lines(c(2.5, 2.5), c(-5, 23.5), lty = 2)
text(1.5, -1, "First to second capture")
text(3.5, -1, "After second capture")
Hey! An actual difference that makes sense. It seems that dulling in these two years pretty consistently resulted in more unique visitors (both male and female) coming to the box. This is done a little bit of a different (simpler) way than in the previous year experiment. I'm only looking at days 6-15 and I haven't (yet) thought about initial coloration or interactions with different breeding stages. We did that in the preprint for 2017, although it gets really messy to explain and this is already complicated... For the 'challenge then dulling' experiment I really should split up at least into two periods since those birds actually were not dulled until after the second capture and that might be resulting in what looks like a slightly smaller effect size than in 2018. In 2019, there looks like maybe a trend for an interaction with the predation treatment too (and remember that predation treatment happened during incubation so wouldn't be too surprising if effect on visitors was a bit different).
One thing that I want to check on is that the predator difference (especially in the second year) might be driven in part by earlier failure in those predator nests so that they have fewer days to accrue unique visitors. In the previous paper the model we used was uniqe visitors each day as the response with multiple days of observation from the same bird controlled by a random effect. That probably makes sense to do again here but I just need to work on the code to get it into shape for that.
Female Trips to Other Boxes
### Trips to other boxes ----
## trips to other boxes
for(i in 1:nrow(d_fem)){
sub <- subset(d_social, as.character(d_social$rfid) == as.character(d_fem$fRFID[i]))
sub1 <- subset(sub, sub$doy >= d_fem$cap1[i] & sub$doy < d_fem$cap2[i])
sub2 <- subset(sub, sub$doy >= d_fem$cap2[i])
d_fem$tot_trip1[i] <- nrow(sub1)
d_fem$tot_trip2[i] <- nrow(sub2)
d_fem$uni_trip1[i] <- length(unique(sub1$ubox))
d_fem$uni_trip2[i] <- length(unique(sub2$ubox))
}
##trips to other boxes
long_vis <- d_fem %>%
pivot_longer(cols = c("uni_trip1", "uni_trip2"), names_to = "type",
names_transform = list(
type = ~ readr::parse_factor(.x, levels = c("uni_trip1", "uni_trip2"),
ordered = TRUE)),
values_to = "count", values_drop_na = TRUE)
long_vis <- as.data.frame(long_vis)
vis_pos <- data.frame(type = c("uni_trip1", "uni_trip2"),
xpos = seq(1, 4, 1))
long_vis <- plyr::join(long_vis, vis_pos, "type", "left", "first")
sum_vis <- d_fem %>%
pivot_longer(cols = c("uni_trip1", "uni_trip2"), names_to = "type",
names_transform = list(
type = ~ readr::parse_factor(.x, levels = c("uni_trip1", "uni_trip2"),
ordered = TRUE)),
values_to = "count", values_drop_na = TRUE) %>%
group_by(as.factor(year), full_treatment, type) %>%
summarise(n = n(), mu = mean(count, na.rm = TRUE), se = sd(count, na.rm = TRUE) / sqrt(n()))
sum_vis <- as.data.frame(sum_vis)
sum_vis <- plyr::join(sum_vis, vis_pos, "type", "left", "first")
colnames(sum_vis)[2] <- "year"
for_plot <- data.frame(full_treatment = unique(sum_vis$full_treatment),
p_shape = c(21, 24, 22, 21, 24, 22, 21, 24, 24),
p_col = c(rep(col_sham, 3), rep(col_dull, 3), col_dull, col_sham, col_dull),
p_dodge = c(-.25, -.15, -.05, .05, .15, .25, .0833, -.0833, .25))
sum_vis <- plyr::join(sum_vis, for_plot, "full_treatment")
long_vis <- plyr::join(long_vis, for_plot, "full_treatment")
vis2_18 <- subset(long_vis, long_vis$year == "2018")
vis2_19 <- subset(long_vis, long_vis$year == "2019")
svis2_18 <- subset(sum_vis, sum_vis$year == "2018")
svis2_19 <- subset(sum_vis, sum_vis$year == "2019")
## Trips
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Dulling then Challenge", xlim = c(0.5, 2.5), ylim = c(-.1, 25),
xlab = "", ylab = "Trips to Unique Boxes")
axis(1, c(-1, 1, 2, 10), c("", "Capture 1-2", "After Capture 2", ""))
axis(2, seq(-30, 200, 5), las = 2)
points(vis2_18$xpos + vis2_18$p_dodge, vis2_18$count, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(svis2_18)){
lines(rep(svis2_18$xpos[i] + svis2_18$p_dodge[i], 2),
c(svis2_18$mu[i] - svis2_18$se[i], svis2_18$mu[i] + svis2_18$se[i]), lwd = 2)
}
points(svis2_18$xpos + svis2_18$p_dodge, svis2_18$mu, pch = svis2_18$p_shape, cex = 1.4,
bg = as.character(svis2_18$p_col))
legend(.5, 23, c("No Color", "Dull Color", "Control", "Predator", "Handicap"),
pch = c(21, 21, 21, 24, 22), pt.bg = c(col_sham, col_dull, rep("white", 3)), cex =0.9, bty = "n")
rect(0.5, 23.5, 2.5, 25, col = alpha("chartreuse3", 0.3))
text(1.3, 24.25, "Dulling", pos = 4)
rect(1.5, 21.5, 2.5, 23, col = alpha("coral3", 0.2))
text(1.9, 22.25, "Challenge", pos = 4)
lines(c(1.5, 1.5), c(-5, 23.5), lty = 2)
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Challenge then Dulling", xlim = c(0.5, 2.5), ylim = c(-.1, 25),
xlab = "", ylab = "Trips to Unique Boxes")
axis(1, c(-1, 1, 2, 10), c("", "Capture 1-2", "After Capture 2", ""))
axis(2, seq(-30, 200, 5), las = 2)
points(vis2_19$xpos + vis2_19$p_dodge, vis2_19$count, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(svis2_19)){
lines(rep(svis2_19$xpos[i] + svis2_19$p_dodge[i], 2),
c(svis2_19$mu[i] - svis2_19$se[i], svis2_19$mu[i] + svis2_19$se[i]), lwd = 2)
}
points(svis2_19$xpos + svis2_19$p_dodge, svis2_19$mu, pch = svis2_19$p_shape, cex = 1.4,
bg = as.character(svis2_19$p_col))
legend(.5, 21.5, c("No Color", "Dull Color", "Control", "Predator"),
pch = c(21, 21, 21, 24), pt.bg = c(col_sham, col_dull, rep("white", 2)), cex =0.9, bty = "n")
rect(1.5, 23.5, 2.5, 25, col = alpha("chartreuse3", 0.3))
text(1.9, 24.25, "Dulling", pos = 4)
rect(0.5, 21.5, 1.5, 23, col = alpha("coral3", 0.2))
text(0.9, 22.25, "Challenge", pos = 4)
lines(c(1.5, 1.5), c(-5, 23.5), lty = 2)
So for the 2018 experiment it looks like dulled birds generally also made more trips to other boxes, but you can only really see this in the second time period when more trips are happening. In 2019 there isn't really much evidence for a difference, except that maybe the predator treatment birds made a bit fewer trips (although again those ones failed earlier so could be some effect of fewer days of observation, should model it like above).
I guess one interpretation for there being a trend towards social effects in 2019 but not as strong as in 2018 is just that the dulling wasn't applied until later and only applied twice instead of three times, so maybe it's expected to be a less obvious impact overall.
Reproductive Success
## Reproductive Success
long_rs <- d_fem %>%
pivot_longer(cols = c("clutch", "maxbrood", "numd6", "numband", "num_d15", "numfled"), names_to = "time_point",
names_transform = list(
time_point = ~ readr::parse_factor(.x, levels = c("clutch", "maxbrood", "numd6", "numband", "num_d15", "numfled"),
ordered = TRUE)),
values_to = "count", values_drop_na = TRUE)
long_rs <- as.data.frame(long_rs)
rs_pos <- data.frame(time_point = c("clutch", "maxbrood", "numd6", "numband", "num_d15", "numfled"),
xpos = seq(1, 6, 1))
long_rs <- plyr::join(long_rs, rs_pos, "time_point", "left", "first")
sum_rs <- d_fem %>%
pivot_longer(cols = c("clutch", "maxbrood", "numd6", "numband", "num_d15", "numfled"), names_to = "time_point",
names_transform = list(
time_point = ~ readr::parse_factor(.x, levels = c("clutch", "maxbrood", "numd6", "numband", "num_d15", "numfled"),
ordered = TRUE)),
values_to = "count", values_drop_na = TRUE) %>%
group_by(as.factor(year), full_treatment, time_point) %>%
summarise(n = n(), mu = mean(count, na.rm = TRUE), se = sd(count, na.rm = TRUE) / sqrt(n()))
sum_rs <- as.data.frame(sum_rs)
sum_rs <- plyr::join(sum_rs, rs_pos, "time_point", "left", "first")
colnames(sum_rs)[2] <- "year"
for_plot <- data.frame(full_treatment = unique(sum_rs$full_treatment),
p_shape = c(21, 24, 22, 21, 24, 22, 21, 24, 24),
p_col = c(rep(col_sham, 3), rep(col_dull, 3), col_dull, col_sham, col_dull),
p_dodge = c(-.25, -.15, -.05, .05, .15, .25, .0833, -.0833, .25))
sum_rs <- plyr::join(sum_rs, for_plot, "full_treatment")
long_rs <- plyr::join(long_rs, for_plot, "full_treatment")
rs2_18 <- subset(long_rs, long_rs$year == "2018")
rs2_19 <- subset(long_rs, long_rs$year == "2019")
srs2_18 <- subset(sum_rs, sum_rs$year == "2018")
srs2_19 <- subset(sum_rs, sum_rs$year == "2019")
#2018 plot
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Dulling then Challenge", xlim = c(0.5, 6.5), ylim = c(-0.1, 7.1),
xlab = "", ylab = "Count")
axis(1, c(-1, 1, 2, 3, 4, 5, 6, 10), c("", "Clutch", "Hatch", "Day 6", "Day 12", "Day 15", "Fledged", ""))
axis(2, seq(-1, 10, 1), las = 2)
points(rs2_18$xpos + rs2_18$p_dodge, rs2_18$count, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(srs2_18)){
lines(rep(srs2_18$xpos[i] + srs2_18$p_dodge[i], 2),
c(srs2_18$mu[i] - srs2_18$se[i], srs2_18$mu[i] + srs2_18$se[i]), lwd = 2)
}
points(srs2_18$xpos + srs2_18$p_dodge, srs2_18$mu, pch = srs2_18$p_shape, cex = 1.4,
bg = as.character(srs2_18$p_col))
abline(v = 1.5, lty = 2)
abline(v = 2.5, lty = 2)
abline(v = 3.5, lty = 2)
abline(v = 4.5, lty = 2)
abline(v = 5.5, lty = 2)
rect(0, 6.6, 3.5, 7, col = alpha("chartreuse3", 0.3))
text(1.7, 6.8, "Dulling", pos = 4)
rect(2, 6, 3, 6.4, col = alpha("coral3", 0.2))
text(2.1, 6.2, "Challenge", pos = 4)
legend("bottomleft", c("No Color", "Dull Color", "Control", "Predator", "Handicap"),
pch = c(21, 21, 21, 24, 22), pt.bg = c(col_sham, col_dull, rep("white", 3)), cex =0.9, bty = "n")
#2019 plot
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Challenge then Dulling", xlim = c(0.5, 6.5), ylim = c(-0.1, 7.1),
xlab = "", ylab = "Count")
axis(1, c(-1, 1, 2, 3, 4, 5, 6, 10), c("", "Clutch", "Hatch", "Day 6", "Day 12", "Day 15", "Fledged", ""))
axis(2, seq(-1, 10, 1), las = 2)
points(rs2_19$xpos + rs2_19$p_dodge, rs2_19$count, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(srs2_19)){
lines(rep(srs2_19$xpos[i] + srs2_19$p_dodge[i], 2),
c(srs2_19$mu[i] - srs2_19$se[i], srs2_19$mu[i] + srs2_19$se[i]), lwd = 2)
}
points(srs2_19$xpos + srs2_19$p_dodge, srs2_19$mu, pch = srs2_19$p_shape, cex = 1.4,
bg = as.character(srs2_19$p_col))
abline(v = 1.5, lty = 2)
abline(v = 2.5, lty = 2)
abline(v = 3.5, lty = 2)
abline(v = 4.5, lty = 2)
abline(v = 5.5, lty = 2)
rect(0, 6.6, 3.5, 7, col = alpha("chartreuse3", 0.3))
text(1.7, 6.8, "Dulling", pos = 4)
text(.9, 6.2, "Challenge", pos = 4)
arrows(0.5, 6.2, 0.9, 6.2, code = 1)
legend("bottomleft", c("No Color", "Dull Color", "Control", "Predator"),
pch = c(21, 21, 21, 24), pt.bg = c(col_sham, col_dull, rep("white", 2)), cex =0.9, bty = "n")
So there does seem to be a pretty clear effect of predator treatment on fledging success. In 2018 the pattern is very clear and appears right after the treatments were delivered (days 1-5 after hatching). It also looks like there could be a moderate decrease in fledging success from the handicap (wing taping) treatment, but not nearly as pronounced. In 2019, the same pattern is present though maybe a little bit weaker. Interestingly, in 2019 the challenges were applied during incubation but even though there was no difference in number hatched in the predator treatments, they had reduced numbers of nestlings later on (not sure how strong this evidence is).
Nestling Mass
## Nestling Mass ----
long_nm <- d_fem %>%
pivot_longer(cols = c("d6_av_mass", "d12_av_mass", "d15_av_mass"), names_to = "time_point",
names_transform = list(
time_point = ~ readr::parse_factor(.x, levels = c("d6_av_mass", "d12_av_mass", "d15_av_mass"),
ordered = TRUE)),
values_to = "mass", values_drop_na = TRUE)
long_nm <- as.data.frame(long_nm)
nm_pos <- data.frame(time_point = c("d6_av_mass", "d12_av_mass", "d15_av_mass"),
xpos = seq(1, 3, 1))
long_nm <- plyr::join(long_nm, nm_pos, "time_point", "left", "first")
sum_nm <- d_fem %>%
pivot_longer(cols = c("d6_av_mass", "d12_av_mass", "d15_av_mass"), names_to = "time_point",
names_transform = list(
time_point = ~ readr::parse_factor(.x, levels = c("d6_av_mass", "d12_av_mass", "d15_av_mass"),
ordered = TRUE)),
values_to = "mass", values_drop_na = TRUE) %>%
group_by(as.factor(year), full_treatment, time_point) %>%
summarise(n = n(), mu = mean(mass, na.rm = TRUE), se = sd(mass, na.rm = TRUE) / sqrt(n()))
sum_nm <- as.data.frame(sum_nm)
sum_nm <- plyr::join(sum_nm, nm_pos, "time_point", "left", "first")
colnames(sum_nm)[2] <- "year"
for_plot <- data.frame(full_treatment = unique(sum_nm$full_treatment),
p_shape = c(21, 24, 22, 21, 24, 22, 21, 24, 24),
p_col = c(rep(col_sham, 3), rep(col_dull, 3), col_dull, col_sham, col_dull),
p_dodge = c(-.25, -.15, -.05, .05, .15, .25, .0833, -.0833, .25))
sum_nm <- plyr::join(sum_nm, for_plot, "full_treatment")
long_nm <- plyr::join(long_nm, for_plot, "full_treatment")
nm2_18 <- subset(long_nm, long_nm$year == "2018")
nm2_19 <- subset(long_nm, long_nm$year == "2019")
snm2_18 <- subset(sum_nm, sum_nm$year == "2018")
snm2_19 <- subset(sum_nm, sum_nm$year == "2019")
##hbill and wing
long_hbw <- d_fem %>%
pivot_longer(cols = c("d12_av_hb", "d12_av_wing"), names_to = "measure",
names_transform = list(
measure = ~ readr::parse_factor(.x, levels = c("d12_av_hb", "d12_av_wing"),
ordered = TRUE)),
values_to = "length", values_drop_na = TRUE)
long_hbw <- as.data.frame(long_hbw)
hbw_pos <- data.frame(measure = c("d12_av_hb", "d12_av_wing"),
xpos = seq(1, 2, 1))
long_hbw <- plyr::join(long_hbw, hbw_pos, "measure", "left", "first")
sum_hbw <- d_fem %>%
pivot_longer(cols = c("d12_av_hb", "d12_av_wing"), names_to = "measure",
names_transform = list(
measure = ~ readr::parse_factor(.x, levels = c("d12_av_hb", "d12_av_wing"),
ordered = TRUE)),
values_to = "length", values_drop_na = TRUE) %>%
group_by(as.factor(year), full_treatment, measure) %>%
summarise(n = n(), mu = mean(length, na.rm = TRUE), se = sd(length, na.rm = TRUE) / sqrt(n()))
sum_hbw <- as.data.frame(sum_hbw)
sum_hbw <- plyr::join(sum_hbw, hbw_pos, "measure", "left", "first")
colnames(sum_hbw)[2] <- "year"
for_plot <- data.frame(full_treatment = unique(sum_hbw$full_treatment),
p_shape = c(21, 24, 22, 21, 24, 22, 21, 24, 24),
p_col = c(rep(col_sham, 3), rep(col_dull, 3), col_dull, col_sham, col_dull),
p_dodge = c(-.25, -.15, -.05, .05, .15, .25, .0833, -.0833, .25))
sum_hbw <- plyr::join(sum_hbw, for_plot, "full_treatment")
long_hbw <- plyr::join(long_hbw, for_plot, "full_treatment")
hbw2_18 <- subset(long_hbw, long_hbw$year == "2018")
hbw2_19 <- subset(long_hbw, long_hbw$year == "2019")
shbw2_18 <- subset(sum_hbw, sum_hbw$year == "2018")
shbw2_19 <- subset(sum_hbw, sum_hbw$year == "2019")
#2018 plot
par(mfrow = c(1, 1))
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Dulling then Challenge", xlim = c(0.5, 3.5), ylim = c(-0.1, 25),
xlab = "", ylab = "Average Nestling Mass (g)")
axis(1, c(-1, 1, 2, 3, 10), c("", "Day 6", "Day 12", "Day 15", ""))
axis(2, seq(-30, 30, 5), las = 2)
points(nm2_18$xpos + nm2_18$p_dodge, nm2_18$mass, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(snm2_18)){
lines(rep(snm2_18$xpos[i] + snm2_18$p_dodge[i], 2),
c(snm2_18$mu[i] - snm2_18$se[i], snm2_18$mu[i] + snm2_18$se[i]), lwd = 2)
}
points(snm2_18$xpos + snm2_18$p_dodge, snm2_18$mu, pch = snm2_18$p_shape, cex = 1.4,
bg = as.character(snm2_18$p_col))
abline(v = 1.5, lty = 2)
abline(v = 2.5, lty = 2)
legend("bottomright", c("No Color", "Dull Color", "Control", "Predator", "Handicap"),
pch = c(21, 21, 21, 24, 22), pt.bg = c(col_sham, col_dull, rep("white", 3)), cex =0.9, bty = "n")
## HB and wing
par(mfrow = c(1, 2))
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Dulling then Challenge", xlim = c(1.5, 2.5), ylim = c(15, 60),
xlab = "", ylab = "Length on Day 12 (mm)")
axis(1, c(-1, 1, 2, 10), c("", "Head + Bill", "Wing", ""))
axis(2, seq(-30, 200, 5), las = 2)
points(hbw2_18$xpos + hbw2_18$p_dodge, hbw2_18$length, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(shbw2_18)){
lines(rep(shbw2_18$xpos[i] + shbw2_18$p_dodge[i], 2),
c(shbw2_18$mu[i] - shbw2_18$se[i], shbw2_18$mu[i] + shbw2_18$se[i]), lwd = 2)
}
points(shbw2_18$xpos + shbw2_18$p_dodge, shbw2_18$mu, pch = shbw2_18$p_shape, cex = 1.4,
bg = as.character(shbw2_18$p_col))
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Dulling then Challenge", xlim = c(0.5, 1.5), ylim = c(19, 29),
xlab = "", ylab = "Length on Day 12 (mm)")
axis(1, c(-1, 1, 2, 10), c("", "Head + Bill", "Wing", ""))
axis(2, seq(-30, 200, 2), las = 2)
points(hbw2_18$xpos + hbw2_18$p_dodge, hbw2_18$length, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(shbw2_18)){
lines(rep(shbw2_18$xpos[i] + shbw2_18$p_dodge[i], 2),
c(shbw2_18$mu[i] - shbw2_18$se[i], shbw2_18$mu[i] + shbw2_18$se[i]), lwd = 2)
}
points(shbw2_18$xpos + shbw2_18$p_dodge, shbw2_18$mu, pch = shbw2_18$p_shape, cex = 1.4,
bg = as.character(shbw2_18$p_col))
legend("bottomright", c("No Color", "Dull Color", "Control", "Predator", "Handicap"),
pch = c(21, 21, 21, 24, 22), pt.bg = c(col_sham, col_dull, rep("white", 3)), cex =0.9, bty = "n")
## 2019 plot
par(mfrow = c(1, 1))
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Challenge then Dulling", xlim = c(0.5, 3.5), ylim = c(-0.1, 25),
xlab = "", ylab = "Average Nestling Mass (g)")
axis(1, c(-1, 1, 2, 3, 10), c("", "Day 6", "Day 12", "Day 15", ""))
axis(2, seq(-30, 30, 5), las = 2)
points(nm2_19$xpos + nm2_19$p_dodge, nm2_19$mass, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(snm2_19)){
lines(rep(snm2_19$xpos[i] + snm2_19$p_dodge[i], 2),
c(snm2_19$mu[i] - snm2_19$se[i], snm2_19$mu[i] + snm2_19$se[i]), lwd = 2)
}
points(snm2_19$xpos + snm2_19$p_dodge, snm2_19$mu, pch = snm2_19$p_shape, cex = 1.4,
bg = as.character(snm2_19$p_col))
abline(v = 1.5, lty = 2)
abline(v = 2.5, lty = 2)
legend("bottomright", c("No Color", "Dull Color", "Control", "Predator", "Handicap"),
pch = c(21, 21, 21, 24, 22), pt.bg = c(col_sham, col_dull, rep("white", 3)), cex =0.9, bty = "n")
## HB and wing
par(mfrow = c(1, 2))
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Challenge then Dulling", xlim = c(1.5, 2.5), ylim = c(15, 60),
xlab = "", ylab = "Length on Day 12 (mm)")
axis(1, c(-1, 1, 2, 10), c("", "Head + Bill", "Wing", ""))
axis(2, seq(-30, 200, 5), las = 2)
points(hbw2_19$xpos + hbw2_19$p_dodge, hbw2_19$length, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(shbw2_19)){
lines(rep(shbw2_19$xpos[i] + shbw2_19$p_dodge[i], 2),
c(shbw2_19$mu[i] - shbw2_19$se[i], shbw2_19$mu[i] + shbw2_19$se[i]), lwd = 2)
}
points(shbw2_19$xpos + shbw2_19$p_dodge, shbw2_19$mu, pch = shbw2_19$p_shape, cex = 1.4,
bg = as.character(shbw2_19$p_col))
plot(1, 1, type = "n", xaxs = "i", yaxs = "i", xaxt = "n", yaxt = "n", bty = "n",
main = "Challenge then Dulling", xlim = c(0.5, 1.5), ylim = c(19, 29),
xlab = "", ylab = "Length on Day 12 (mm)")
axis(1, c(-1, 1, 2, 10), c("", "Head + Bill", "Wing", ""))
axis(2, seq(-30, 200, 2), las = 2)
points(hbw2_19$xpos + hbw2_19$p_dodge, hbw2_19$length, pch = 16, col = alpha("black", 0.3))
#abline(h = 0)
for(i in 1:nrow(shbw2_19)){
lines(rep(shbw2_19$xpos[i] + shbw2_19$p_dodge[i], 2),
c(shbw2_19$mu[i] - shbw2_19$se[i], shbw2_19$mu[i] + shbw2_19$se[i]), lwd = 2)
}
points(shbw2_19$xpos + shbw2_19$p_dodge, shbw2_19$mu, pch = shbw2_19$p_shape, cex = 1.4,
bg = as.character(shbw2_19$p_col))
legend("bottomright", c("No Color", "Dull Color", "Control", "Predator"),
pch = c(21, 21, 21, 24), pt.bg = c(col_sham, col_dull, rep("white", 2)), cex =0.9, bty = "n")
There isn't really anything very obvious going on here. Maybe there is a hint of reduced growth in predation treatments. This is just plotting the average value for each nest, but I really probably should plot these as individual nestling points with nest as a random effect for the day 12 & 15 measures. One thing to consider is that because there seems to be higher nestling attrition in the predation treatments, it may just be that the lightest nestlings in those broods are already gone by the time day 12/15 measurements are taken.
Other Available Data
There are some other data types that we have on these birds and could look at or incorporate here in various ways (although there is already a lot in here). Those include:
- microbiome of adults and nestlings (save for different analyses?)
- paternity and cross fostering status of nestlings
- some diet data (save for different analyses)
- incubation data from HOBOS (would need to double check how many have data)
- glucose, BUTY, maybe some oxidative stress (no clear reason to include)
- nestling RFID data box entrance occupation, fledging time, visits (save for different paper)
- next year breeding data from adults in treatments?
- also these plots don't currently account for pre-treatment coloration in any way
Summary
I'm not sure exactly how this all gets put together into a paper. It seems like the main takeaways are:
- The color manipulation does seem to have some effect on social interactions, although maybe weaker in one year than the other.
- The predator treatment had a pretty clear impact on reproductive success that seems to be driven at least in part by changes in parental behavior (given that we see it even when the predator happens before hatching), but it isn't obvious what that change in behavior is since there isn't a clear change in provisioning rate (though provisioning should still fit in a more sophisticated model). There is maybe a hint that handicaping flight had some similar effects but not very strong.
- Despite the fact that we changed the social interactions and imposed a challenge, there really isn't any clear evidence or even suggestion that those two treatments interacted with each other, such that having one impacted the response to the other.
- There are lots of caveats and ways that other types of treatments could produce interactions between sociality and challenges, so I wouldn't say we have clear evidence to argue against that possibility, but these particular treatments don't seem to indicate an interaction and our evidence is most consistent with social behavior and environmental challenges operating indepednetly in this case.
Discussion Points
- Maybe our treatments just weren't severe enough to count as a real challenge?
- Maybe carryover effects in subsequent years?
- The biggest changes in social behavior seem to occur later in the provisioning period and maybe that is too late to have any effect on the challenges? Possible that manipulation of social environment during box settlement would be more salient but we don't have an easy way to manipulate then.
Analysis Info
sessionInfo()## R version 3.6.1 (2019-07-05)
## Platform: x86_64-apple-darwin15.6.0 (64-bit)
## Running under: OS X El Capitan 10.11.6
##
## Matrix products: default
## BLAS: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRblas.0.dylib
## LAPACK: /Library/Frameworks/R.framework/Versions/3.6/Resources/lib/libRlapack.dylib
##
## locale:
## [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] scales_1.1.1 lme4_1.1-23 Matrix_1.2-18 here_0.1
## [5] forcats_0.5.0 stringr_1.4.0 dplyr_1.0.1 purrr_0.3.4
## [9] readr_1.3.1 tidyr_1.1.1 tibble_3.0.3 ggplot2_3.3.2
## [13] tidyverse_1.3.0
##
## loaded via a namespace (and not attached):
## [1] Rcpp_1.0.5 lubridate_1.7.9 lattice_0.20-41 assertthat_0.2.1
## [5] rprojroot_1.3-2 digest_0.6.25 plyr_1.8.6 R6_2.4.1
## [9] cellranger_1.1.0 backports_1.1.8 reprex_0.3.0 evaluate_0.14
## [13] httr_1.4.2 pillar_1.4.6 rlang_0.4.7 readxl_1.3.1
## [17] rstudioapi_0.11 minqa_1.2.4 nloptr_1.2.2.2 blob_1.2.1
## [21] rmarkdown_2.3 labeling_0.3 splines_3.6.1 statmod_1.4.34
## [25] munsell_0.5.0 broom_0.7.0 compiler_3.6.1 modelr_0.1.8
## [29] xfun_0.16 pkgconfig_2.0.3 mgcv_1.8-31 htmltools_0.5.0
## [33] tidyselect_1.1.0 fansi_0.4.1 crayon_1.3.4 dbplyr_1.4.4
## [37] withr_2.2.0 MASS_7.3-51.6 grid_3.6.1 nlme_3.1-148
## [41] jsonlite_1.7.0 gtable_0.3.0 lifecycle_0.2.0 DBI_1.1.0
## [45] pacman_0.5.1 magrittr_1.5 cli_2.0.2 stringi_1.4.6
## [49] farver_2.0.3 fs_1.5.0 xml2_1.3.2 ellipsis_0.3.1
## [53] generics_0.0.2 vctrs_0.3.2 boot_1.3-25 tools_3.6.1
## [57] glue_1.4.1 hms_0.5.3 yaml_2.2.1 colorspace_1.4-1
## [61] rvest_0.3.6 knitr_1.29 haven_2.3.1