-
Notifications
You must be signed in to change notification settings - Fork 8
/
Copy pathindex.Rmd
86 lines (71 loc) · 3.08 KB
/
index.Rmd
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
---
title: "Portal Forecast"
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
library(dplyr)
library(portalr)
source("forecast_tools.R")
```
## Total Abundance Forecast
This is the forecast for next month's sampling of rodents at Portal.
```{r, echo=FALSE, message=FALSE, warning=FALSE}
obs_data = portalr::abundance("repo")
obs_data$total = rowSums(select(obs_data, -period))
new_moons = read.csv(text = RCurl::getURL("https://raw.githubusercontent.com/weecology/PortalData/master/Rodents/moon_dates.csv"))
obs_data_newmoon = inner_join(obs_data, new_moons, by = c("period" = "period"))
obs_data_newmoon$censusdate = as.Date(obs_data_newmoon$censusdate)
#for_data = read.csv("predictions/2017-02-24Allforecasts.csv")
for_data = compile_forecasts() %>%
dplyr::filter(level == 'All', model == 'Ensemble', date == max(date))
most_recent_forecast = max(for_data$date)
for_data = filter(for_data, date == most_recent_forecast)
for_data$model = "ensemble"
for_data$forecastdate = as.Date(paste(for_data$forecastyear,for_data$forecastmonth,'15',sep='-'))
forecast_viz(obs_data = obs_data_newmoon,
obs_date_col_name = "censusdate",
obs_val_col_name = "total",
for_data = for_data,
for_date_col_name = "forecastdate",
for_val_col_name = "estimate",
for_model_name = "ensemble",
for_lowerpi_col_name = "LowerPI",
for_upperpi_col_name = "UpperPI",
start_newmoon = 300,
ylabel = 'Total Abundance')
```
## Species-Level Forecasts
```{r speciesforecasts, echo=FALSE, message=FALSE, warning=FALSE}
source("forecast_tools.R")
data = compile_forecasts()
ensemble = dplyr::filter(data, level == 'All', model == 'Ensemble', date == max(date))
sp_predictions = get_sp_predicts(ensemble, 'All', lead_time = 1)
sp_predict = plot_species_forecast(sp_predictions)
plot(sp_predict)
```
```{r highabund_ts_forecasts, echo=FALSE, message=FALSE, warning=FALSE}
most_abund_sp = sp_predictions %>%
filter(species != 'total') %>%
arrange(desc(estimate)) %>%
head(3) %>%
select(species)
# load in rodent species table to get scientific names to display on plots
species_table = read.csv(
text = RCurl::getURL(
"https://raw.githubusercontent.com/weecology/PortalData/master/Rodents/Portal_rodent_species.csv"),stringsAsFactors = F,na.strings = '')
most_abund_sp_names = filter(species_table,speciescode %in% most_abund_sp$species) %>% select(speciescode,scientificname)
for (n in seq(dim(most_abund_sp_names)[1])) {
species_forecast = forecast_viz(obs_data = obs_data_newmoon,
obs_date_col_name = "censusdate",
obs_val_col_name = most_abund_sp_names$speciescode[n],
for_data = for_data,
for_date_col_name = "forecastdate",
for_val_col_name = "estimate",
for_model_name = "ensemble",
for_lowerpi_col_name = "LowerPI",
for_upperpi_col_name = "UpperPI",
start_newmoon = 300,
ylabel = most_abund_sp_names$scientificname[n])
plot(species_forecast)
}
```