RPI-Trends/iv_plotcoverage.R at main · BirdsCanada/RPI-Trends · 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
## Station coverage plots

# directory to post plots
in.dir <- paste("./Data/", max.yr, "/", sep = "")
plot.dir <- paste("./Plots/", max.yr, "/", sep = "")

site <- as.character(anal.param[t, "SiteCode"])
seas <- as.character(anal.param[t,"seas"])
name<-as.character(anal.param[t,"site"])


## Including Plots
#Plot using raw data the station coverage plots


  ############################################################################
  # READ IN DATA AND MANIPULATE - SITE AND SEASON SPECIFIC

  in.data <- read.csv(paste(in.dir, site, ".", seas, ".RawData.csv", sep = ""))
  #in.data <- subset(in.data, !is.na(datetime))
  #in.data$datetime <- as.POSIXct(in.data$datetime)

  tmp.data<-NULL

  tmp.data<-in.data %>% format_dates()
  tmp.data <- tmp.data %>% select(SiteCode, project_id, YearCollected, MonthCollected, DayCollected, doy, TimeCollected, DurationInHours, species_id, ObservationCount)
  tmp.data<-tmp.data %>% left_join(sp.names, by="species_id")

  # if no reason to drop any years at beginning (in anal.param), use min year in dataframe
   min.yr.filt <- ifelse(is.na(min.yr.filt), min(tmp.data$YearCollected), min.yr.filt)
   max.yr.filt <- ifelse(is.na(max.yr.filt), max.yr, max.yr.filt)

  ## Filter to min and max year

  tmp.data <- filter(tmp.data, YearCollected >= min.yr.filt & YearCollected <= max.yr.filt)

  #print("year range:"); print(range(tmp.data$YearCollected))

  ## Assign date and season variables
  tmp.data <- tmp.data %>%
    mutate(date = ymd(paste(YearCollected, MonthCollected, DayCollected, sep = "/")),
           doy = yday(date),
           season = if_else(doy < 180, "spring", "fall"))

#  tmp.data<-tmp.data %>% drop_na(species_code)

  #truncate spring and fall data
  if(seas =="fall"){
  tmp.data<-tmp.data %>% filter(doy>=213)
  }

    if(seas =="spring"){
    tmp.data<-tmp.data %>% filter(doy>=32 & doy <=181)
  }

  if(site == "484"){
    tmp.data <- tmp.data %>%  filter (!(species_code == "BLVU"),
                                      !(species_code == "TUVU")) %>%
      droplevels()
  }

  obsDays <- unique(subset(tmp.data, select = c("YearCollected", "doy")))
  obsDays <- summaryBy(doy ~ YearCollected, data = obsDays, FUN = c(length, min, max))


  # calculate number of observation hours/day

 obsHours <- tmp.data %>%
    #filter(ObservationCount > 0) %>%
    group_by(SiteCode, YearCollected, MonthCollected, DayCollected, doy, TimeCollected) %>%
    slice_max(DurationInHours) %>%
    select(SiteCode, YearCollected, MonthCollected, DayCollected, doy, TimeCollected, DurationInHours) %>%
    distinct() %>%
    ungroup() %>%
    group_by(SiteCode, YearCollected, MonthCollected, DayCollected, doy) %>% summarize(DurationInHours=sum(DurationInHours)) %>%
    ungroup() %>%
    as.data.frame() %>% drop_na()

  obsHours <- unique(subset(obsHours, select = c("YearCollected", "DurationInHours")))
  obsHours <- summaryBy(DurationInHours ~ YearCollected, data = obsHours, FUN = c(mean))

  maxhours<-max(obsHours$DurationInHours.mean)

  pdf(paste(plot.dir, site, "_", seas, "_", name, "_SamplingCoverPlot.pdf", sep=""),
      height = 10, width = 8, paper = "letter")
  par(mfrow = c(3, 1))
  plot(doy.length ~ YearCollected, data = obsDays,
       ylab = "Number Days Sampled", xlab = "Year",
       col = "black", pch = 20, cex = 1)
  plot(doy.min ~ YearCollected, data = obsDays,
       ylab = "Range of Dates Sampled",
       xlab = "Year", ylim = c(min(obsDays$doy.min), max(obsDays$doy.max)),
       col = c("black"), pch = c(20), cex = 1)
  points(doy.max ~ YearCollected, data = obsDays,
         col = "grey50", pch = 1, cex = 1)
  plot(DurationInHours.mean ~ YearCollected, data = obsHours,
       ylab = "Mean # hours sampled/day", xlab = "Year",
       col = "black", pch = 20, cex = 1, ylim=c(0, maxhours))

  while(!is.null(dev.list())) dev.off()