GeeseBehavior-Weather/16_plotting_weather_variables.R at main · s-cunningham/GeeseBehavior-Weather · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118

## 2022-08-31

library(tidyverse)
library(patchwork)

options(scipen=999, digits=3)

theme_set(theme_classic())

#### Read in bird data ####

# Read in weather covariates
wdat <- read_csv("files_for_models/weather_covars.csv")
range(wdat$prcp, na.rm=TRUE)
wdat$prcp[wdat$prcp==-99999] <- NA

# Read in ACC data
dat <- read_csv("files_for_models/daily_odba_behavior.csv")

# Join weather and ACC data
dat <- left_join(dat, wdat, by=c("animal_id", "date"))
dat <- as.data.frame(dat)

# read in migration dates
mdates <- read_csv("files_for_models/migration_dates.csv")

# save maximum duration of migration period
dur <- max(mdates$duration)

# subset to migration dates
dat$RelDay <- NA
dat$RevRelDay <- NA
dat$migration <- NA
un.id <- unique(dat$animal_id)
for (i in 1:length(un.id)) {

  # subset migration dates
  md <- mdates[mdates$animal_id==un.id[i],]

  # Create vector of days
  mdays <- md$start:md$end

  # Add column, indicate when is migration
  dat$RelDay[dat$animal_id==un.id[i] & (dat$julian %in% mdays)] <- 1:length(mdays)
  dat$RevRelDay[dat$animal_id==un.id[i] & (dat$julian %in% mdays)] <- length(mdays):1
  dat$migration[dat$animal_id==un.id[i] & (dat$julian %in% mdays)] <- "yes"

}

# Subset to just migration
dat <- dat[!is.na(dat$migration), ]

# log odba
dat$lnODBAmedian <- log(dat$median.odba)

# Reorganize columns
dat <- dat[,c(1,3,4,23,24,5:8,10,11,16,18,19,21)]

# Weather summary states
# dat %>% group_by(pop, variable) %>% summarize(range(values))

#### Plot weather variables ####

dat <- pivot_longer(dat, 14:15, names_to="variable", values_to="values") %>% as.data.frame()

prcp <- dat %>% filter(variable=="prcp")
temp <- dat %>% filter(variable=="mintemp")

# Plot precipitation
prcp_plot <- ggplot() +
              geom_line(data=prcp, aes(x=julian, y=values, group=animal_id, color=pop), alpha=0.5, size=1) +
              scale_x_continuous(breaks=c(60,91,121,152),
                                 labels=c("1-Mar","1-Apr","1-May","1-Jun")) +
              scale_color_manual(values=c("#2166ac","#b2182b"), labels=c("Greenland", "Midcontinent")) +
              facet_grid(.~year, scales="free_x") +
              xlab("Date") + ylab("Precipitation (mm)") +
              guides(color=guide_legend(title="Population")) +
              annotate("segment", x=-Inf, xend=Inf, y=-Inf, yend=-Inf)+
              annotate("segment", x=-Inf, xend=-Inf, y=-Inf, yend=Inf) +
              theme(legend.position='bottom',
                    legend.title=element_text(size=12),
                    legend.text=element_text(size=12),
                    legend.background=element_rect(fill=NA),
                    panel.border=element_blank(),
                    axis.line=element_line(),
                    axis.text.y=element_text(size=12),
                    axis.text.x=element_blank(),
                    axis.title.y=element_text(size=12),
                    axis.title.x=element_blank(),
                    strip.text=element_text(size=12),
                    strip.background=element_rect(fill="white", color=NA))

temp_plot <- ggplot() +
  geom_line(data=temp, aes(x=julian, y=values, group=animal_id, color=pop), alpha=0.5, size=1) +
  scale_x_continuous(breaks=c(60,91,121,152),
                     labels=c("1-Mar","1-Apr","1-May","1-Jun")) +
  scale_color_manual(values=c("#2166ac","#b2182b"), labels=c("Greenland", "Midcontinent")) +
  facet_grid(.~year, scales="free_x") +
  xlab("Date") + ylab("Temperature (C)") +
  guides(color=guide_legend(title="Population")) +
  annotate("segment", x=-Inf, xend=Inf, y=-Inf, yend=-Inf)+
  annotate("segment", x=-Inf, xend=-Inf, y=-Inf, yend=Inf) +
  theme(legend.position='bottom',
        legend.title=element_text(size=12),
        legend.text=element_text(size=12),
        legend.background=element_rect(fill=NA),
        panel.border=element_blank(),
        axis.line=element_line(),
        axis.text=element_text(size=12),
        axis.title=element_text(size=12),
        strip.text=element_blank(),
        strip.background=element_rect(fill="white"))


(prcp_plot / temp_plot + plot_layout(guides="collect")) +
  plot_annotation(tag_levels="a", tag_prefix="(", tag_suffix=")") &
  theme(legend.position = 'bottom')