In this analysis, I attempt to analyze whether there is a relationship between age and capture probability in Wood Thrushes.
migrant %>%
# Filter for valid calendar age and year values
filter(str_detect(calendar_age, "Y$"),
str_detect(year, "^20"),
!year %in% c("2012","2013")) %>%
# Make a season variable to specify banding season of captures
mutate(
month = as.numeric(month),
season_yr = case_when(
month < 7 ~ "Spring",
TRUE ~ "Fall") %>%
str_c(year, sep = " ")) %>%
# Do not keep captures in the same season
distinct(band_full,
calendar_age,
season_yr) %>%
# Rearrange calendar_age variable so that ages are in order
mutate(calendar_age = calendar_age %>%
fct_relevel("ASY",
"AHY",
"SY",
"HY"),
# Rearrange seasons so they are in order
season_yr = season_yr %>%
fct_relevel("Spring 2014", "Fall 2014", "Spring 2015", "Fall 2015", "Spring 2016", "Fall 2016", "Spring 2017", "Fall 2017", "Spring 2018", "Fall 2018", "Spring 2019", "Fall 2019", "Spring 2020", "Fall 2020", "Spring 2021", "Fall 2021", "Spring 2022", "Fall 2022")) %>%
# Plot captures over time, with fill as age of birds
ggplot(aes(x = season_yr, fill = calendar_age)) +
geom_bar() +
theme(axis.text.x = element_text(angle = 90, hjust = 1)) +
labs(title = "Wood Thrushes captured at TREES over time",
x = "Banding Season",
y = "Number of Individuals",
color = "Calendar Age") +
scale_y_continuous(expand = c(0,0),
breaks = seq(0, 1500, 100),
limits = c(0,1500)) +
scale_fill_brewer(palette = "Spectral")
As seen in the graph, the age of Wood Thrushes varies between banding seasons. In addition, the number of individual also varies between banding seasons, but the reasons are unknown (may be due to varying banding effort or disease).
I wanted to know how many seasons an individual was detected again at TREES.
mean_visits <-
migrant %>%
filter(species_code %in% "WOTH") %>%
mutate(
month = as.numeric(month),
season = case_when(
month < 7 ~ "Spring",
TRUE ~ "Fall")) %>%
# Do not keep captures in the same season
distinct(band_full, season, year) %>%
# Get the number of visits to TREES
summarize(visits = n(),
.by = band_full) %>%
# Find average no. of visits to TREES
summarize(mean_visits = mean(visits)) %>%
# Get it as a double
pull(mean_visits)
migrant %>%
filter(species_code %in% "WOTH") %>%
mutate(
month = as.numeric(month),
season = case_when(
month < 7 ~ "Spring",
TRUE ~ "Fall")) %>%
# Do not keep captures in the same season
distinct(band_full, season, year) %>%
# Get the number of visits to TREES
summarize(visits = n(),
.by = band_full) %>%
ggplot(aes(x = visits)) +
geom_bar() +
labs(title = "Number of seasons at TREES",
x = "No. of seasons",
y = "No. of individuals") +
scale_y_continuous(expand = c(0,0),
limits = c(0, 300)) +
geom_vline(xintercept = mean_visits,
color = "red") +
annotate("text",
x = mean_visits + 0.05,
y = 150,
label = "Mean visit for all individuals",
angle = 90,
color = "white")
The mean visit value is 1.226361, which is higher than one, and shows that Wood Thrushes return to TREES to spend the season.
I wanted to find out if younger or older birds tend to be detected at TREES. Below is a plot that proportions of the ages of each Wood Thrush when they were first banded.
migrant %>%
filter(
# Extract Wood Thrushes from data
species_code %in% "WOTH",
# Extract newly banded birds
code %in% "N",
# Extract valid ages
str_detect(calendar_age, "Y$")) %>%
summarize(count = n(),
.by = calendar_age) %>%
# Calculating proportion of individuals
mutate(prop = count / sum(count),
# Reordering calendar age
calendar_age = calendar_age %>%
fct_relevel("HY",
"SY",
"AHY",
"ASY")) %>%
arrange(calendar_age)
## # A tibble: 4 × 3
## calendar_age count prop
## <fct> <int> <dbl>
## 1 HY 163 0.509
## 2 SY 106 0.331
## 3 AHY 26 0.0812
## 4 ASY 25 0.0781
# Graph the results
migrant %>%
filter(
# Extract Wood Thrushes from data
species_code %in% "WOTH",
# Extract newly banded birds
code %in% "N",
# Extract valid ages
str_detect(calendar_age, "Y$")) %>%
summarize(count = n(),
.by = calendar_age) %>%
# Calculating proportion of individuals
mutate(prop = count / sum(count),
# Reordering calendar age
calendar_age = calendar_age %>%
fct_relevel("HY",
"SY",
"AHY",
"ASY")) %>%
arrange(calendar_age) %>%
ggplot(aes(calendar_age, prop)) +
geom_col() +
labs(title = "Age when first captured",
x = "Calendar age",
y = "Proportion of individuals") +
scale_y_continuous(expand = c(0,0),
limits = c(0,0.6),
breaks = seq(0,0.6,0.05))
Note that “HY” is a first-year bird, “SY” is a second-year bird, and “ASY” means after-second-year, which means that the bird is age 3 or higher. “AHY” signifies that the bird is not a first-year bird, but it is unclear whether it is an ASY or SY bird.
51% of Wood Thrushes captured are 1st year birds. 33% are second-year birds. The other 15% may be second-year or third-year and above.
I wanted to find out if younger or older birds tend to return to TREES. Below is a plot that includes only recaptured birds, and their ages when they were first banded.
# Generate a vector of band numbers for recaptured Wood Thrushes
recaps <-
migrant %>%
filter(species_code %in% "WOTH") %>%
mutate(
month = as.numeric(month),
season = case_when(
month < 7 ~ "Spring",
TRUE ~ "Fall")) %>%
# Do not keep captures in the same season
distinct(band_full, season, year) %>%
# Get the number of visits to TREES
summarize(visits = n(),
.by = band_full) %>%
# Get the individuals with more than 1 visit per season
filter(visits > 1) %>%
# Get band numbers of those individuals as a vector
pull(band_full)
# Do more younger birds return to TREES?
migrant %>%
filter(band_full %in% recaps,
code == "N",
str_detect(calendar_age, "Y$")) %>%
select(band_full,
calendar_age) %>%
drop_na(calendar_age) %>%
mutate(`Calendar age when first banded` = calendar_age %>%
fct_relevel("HY",
"SY",
"AHY",
"ASY")) %>%
summarize(count = n(),
.by = `Calendar age when first banded`) %>%
arrange(`Calendar age when first banded`)
## # A tibble: 4 × 2
## `Calendar age when first banded` count
## <fct> <int>
## 1 HY 19
## 2 SY 17
## 3 AHY 4
## 4 ASY 4
# Generate a bar plot of above results
migrant %>%
filter(band_full %in% recaps,
code == "N",
str_detect(calendar_age, "Y$")) %>%
select(band_full,
calendar_age) %>%
drop_na(calendar_age) %>%
mutate(calendar_age = calendar_age %>%
fct_relevel("HY",
"SY",
"AHY",
"ASY")) %>%
ggplot(aes(calendar_age)) +
geom_bar() +
scale_y_continuous(expand = c(0,0),
limits = c(0,20)) +
labs(title = "Initial age of recaptures",
x = "Calendar Age")
It seems that more first year and second year birds return to TREES than third year or above, but the trend seems similar to the graph above.
This may be a function of time or survival and detection differences.
To find out whether age when first banded influences detection, I ran several Cormack-Jolly-Seber models.
I firstly generated a capture history for each Wood Thrush individual. Each 0 or 1 signifies whether the individual was detected at each banding season. I then included initial_age as a column which states the age at first capture. The value 1 means HY, 2 means SY, 3 means AHY and 4 means ASY.
initial_age <-
migrant %>%
filter(species_code %in% "WOTH",
code == "N",
str_detect(calendar_age, "Y$"),
rowid != 3594) %>%
select(band_full,
initial_age = calendar_age) %>%
drop_na(initial_age) %>%
mutate(initial_age = initial_age %>%
fct_recode("1" = "HY",
"2" = "SY",
"3" = "AHY",
"4" = "ASY") %>%
as.numeric())
ch_woth <-
migrant %>%
filter(species_code %in% "WOTH",
str_detect(year, "^20")) %>%
mutate(
month = as.numeric(month),
season =
case_when(month < 7 ~ "S",
TRUE ~ "F"),
value = 1) %>%
distinct(band_full, season, year, value) %>%
unite(season, season, year, sep = "") %>%
mutate(season = season %>%
fct_relevel("F2014", "S2015", "F2015", "S2016", "F2016", "S2017", "F2017", "S2018", "F2018", "S2019", "F2019", "S2020", "F2020", "S2021", "F2021", "S2022", "F2022")) %>%
arrange(season) %>%
pivot_wider(id_cols = band_full,
names_from = season,
values_from = value,
values_fill = 0) %>%
left_join(initial_age, by = "band_full") %>%
unite(ch, F2014:F2022, sep = "") %>%
process.data(model = "CJS",
groups = "initial_age")
head(ch_woth$data)
## band_full ch initial_age group
## 1 1202-15853 10100100010000000 3 3
## 2 1202-15854 10000000000000000 <NA> <NA>
## 3 1202-15856 10000000000000000 3 3
## 4 2561-30093 10000000000000000 3 3
## 5 2561-30094 10000000000000000 3 3
## 6 2561-30095 10000000000000000 3 3
Note that only the first 6 rows of the data are shown above. There are 330 rows in total. I then ran several models that varied in whether survival or detection was dependent on time or initial age. The following code chunks are the CJS models being built. Please scroll down until you see a table of values!
# When set to null, the probability remains constant
null = list(formula=~1)
# When set to initial_age, the probability varies with initial_age
initial_age = list(formula=~initial_age)
# When set to time, the probability varies with time
time = list(formula=~time)
# Intercept Only
p0.phi0 <-
mark(data = ch_woth,
model = "CJS",
model.parameters = list(Phi = null, p = null),
brief = TRUE)
# Survival dependent on initial age
p0.phiinitial_age <-
mark(data = ch_woth,
model = "CJS",
model.parameters = list(Phi = initial_age, p = null),
brief = TRUE)
# Detection dependent on initial age
pinitial_age.phi0 <-
mark(data = ch_woth,
model = "CJS",
model.parameters = list(Phi = null, p = initial_age),
brief = TRUE)
# Detection and survival dependent on initial age
pinitial_age.phiinitial_age <-
mark(data = ch_woth,
model = "CJS",
model.parameters = list(Phi = initial_age,
p = initial_age),
brief = TRUE)
# Detection dependent on time
ptime.phi0 <-
mark(data = ch_woth,
model = "CJS",
model.parameters = list(Phi = null,
p = time),
brief = TRUE)
##
## Note: only 14 parameters counted of 17 specified parameters
## AICc and parameter count have been adjusted upward
# Survival dependent on time
p0.phitime <-
mark(data = ch_woth,
model = "CJS",
model.parameters = list(Phi = time,
p = null),
brief = TRUE)
##
## Note: only 13 parameters counted of 17 specified parameters
##
## AICc and parameter count have been adjusted upward
# Detection and survival dependent on time
ptime.phitime <-
mark(data = ch_woth,
model = "CJS",
model.parameters = list(Phi = time,
p = time),
brief = TRUE)
##
## Note: only 30 parameters counted of 32 specified parameters
##
## AICc and parameter count have been adjusted upward
# Detection dependent on time + Survival dependent on initial age
ptime.phiinitial_age <-
mark(data = ch_woth,
model = "CJS",
model.parameters = list(Phi = initial_age,
p = time),
brief = TRUE)
##
## Note: only 17 parameters counted of 20 specified parameters
##
## AICc and parameter count have been adjusted upward
# Detection dependent on initial age + Survival dependent on time
pinitial_age.phitime <-
mark(data = ch_woth,
model = "CJS",
model.parameters = list(Phi = initial_age,
p = time),
brief = TRUE)
##
## Note: only 17 parameters counted of 20 specified parameters
##
## AICc and parameter count have been adjusted upward
The lowest AICc values signify a model with the best fit.
collect.models()
## model npar AICc DeltaAICc weight
## 7 Phi(~1)p(~time) 17 512.7987 0.000000 7.156194e-01
## 6 Phi(~initial_age)p(~time) 20 516.0309 3.232197 1.421737e-01
## 8 Phi(~initial_age)p(~time) 20 516.0309 3.232197 1.421737e-01
## 3 Phi(~time)p(~1) 17 534.7682 21.969570 1.213530e-05
## 9 Phi(~time)p(~time) 32 534.7992 22.000559 1.194872e-05
## 4 Phi(~1)p(~initial_age) 5 536.8559 24.057178 4.272994e-06
## 1 Phi(~1)p(~1) 2 537.1976 24.398922 3.601832e-06
## 2 Phi(~initial_age)p(~1) 5 540.1132 27.314558 8.383033e-07
## 5 Phi(~initial_age)p(~initial_age) 8 541.8592 29.060562 3.501558e-07
## Deviance
## 7 236.9894
## 6 233.4866
## 8 233.4866
## 3 258.9590
## 9 223.9946
## 4 286.7729
## 1 293.2592
## 2 290.0303
## 5 285.5202
As seen from the table above, the model with best fit suggests that detection varies with time, and survival is constant. The next best fit model suggests that survival is dependent on age when first banded, and detection is dependent on time.
Below I plot the detection probabilities across time.
custom_labels <- c("Fall 2014", "Spring 2015", "Fall 2015", "Spring 2016", "Fall 2016", "Spring 2017", "Fall 2017", "Spring 2018", "Fall 2018", "Spring 2019", "Fall 2019", "Spring 2020", "Fall 2020", "Spring 2021", "Fall 2021", "Spring 2022", "Fall 2022")
tibble(season = 0:16,
p = ptime.phi0$results$real$estimate) %>%
ggplot(aes(x = season, y = p)) +
geom_line() +
scale_x_continuous(breaks = 0:16, labels = custom_labels) +
scale_y_continuous(breaks = seq(0,0.7,0.05),
expand = c(0,0),
limits = c(0,0.7)) +
labs(title = "Probability of detecting Wood Thrushes across banding seasons",
y = "Detection Probability",
x = "Banding Season") +
theme(axis.text.x = element_text(angle = 90, hjust = 1))