From 8634b7fb07fc4d36eb59ea66784e3d1d8ae883d0 Mon Sep 17 00:00:00 2001 From: Ethan White Date: Sun, 20 Sep 2026 14:14:59 -0400 Subject: [PATCH] Add assignment covering R-time-series-modeling* material --- .gitignore | 4 + _quarto.yml | 3 + data/barco_oxygen_timeseries.csv | 93 ++++++++++ data/make_barco_timeseries.R | 57 ++++++ .../R-time-series-modeling-3/assignment.qmd | 136 ++++++++++++++ .../assignment_solution.qmd | 166 ++++++++++++++++++ lessons/R-time-series-modeling-3/index.qmd | 1 + 7 files changed, 460 insertions(+) create mode 100644 data/barco_oxygen_timeseries.csv create mode 100644 data/make_barco_timeseries.R create mode 100644 lessons/R-time-series-modeling-3/assignment.qmd create mode 100644 lessons/R-time-series-modeling-3/assignment_solution.qmd diff --git a/.gitignore b/.gitignore index 4b70f2b..5e17d5c 100644 --- a/.gitignore +++ b/.gitignore @@ -5,3 +5,7 @@ node_modules/ /.quarto/ **/*.quarto_ipynb .Rproj.user + +# Rendered answer key. Built by hand from assignment_solution.qmd; not published +# with the site and not worth committing as a 3 MB self-contained file. +lessons/*/assignment_solution.html diff --git a/_quarto.yml b/_quarto.yml index e95d794..78acbdc 100644 --- a/_quarto.yml +++ b/_quarto.yml @@ -17,6 +17,8 @@ project: # Raw R scratch files with no front matter; kept in the repo, not published. - "!lessons/R-state-space-models-1/r_tutorial_new.qmd" - "!lessons/R-state-space-models-1/r_tutorial_new_notes.qmd" + # Assignment answers don't render + - "!lessons/R-time-series-modeling-3/assignment_solution.qmd" resources: - data/** - media/** @@ -146,6 +148,7 @@ website: - lessons/R-time-series-modeling-3/material.qmd - lessons/R-time-series-modeling-3/r_tutorial.qmd - lessons/R-time-series-modeling-3/equation_review.qmd + - lessons/R-time-series-modeling-3/assignment.qmd - section: Intro to Ecological Forecasting href: lessons/intro-to-forecasting/index.qmd contents: diff --git a/data/barco_oxygen_timeseries.csv b/data/barco_oxygen_timeseries.csv new file mode 100644 index 0000000..4a99163 --- /dev/null +++ b/data/barco_oxygen_timeseries.csv @@ -0,0 +1,93 @@ +date,temperature,oxygen,chla +2018-11-01,20.514,6.045,3.818 +2018-12-01,17.485,8.073,1.465 +2019-01-01,17.063,8.568,0.952 +2019-02-01,18.853,9.088,0.313 +2019-03-01,20.739,8.13,0.702 +2019-04-01,24.223,7.87,0.427 +2019-05-01,28.802,6.873,0.326 +2019-06-01,30.683,6.041,1.053 +2019-07-01,31.618,5.215,1.343 +2019-08-01,31.21,5.753,1.933 +2019-09-01,29.889,5.944,1.985 +2019-10-01,27.373,5.974,1.516 +2019-11-01,21.554,7.296,1.583 +2019-12-01,18.169,8.692,1.776 +2020-01-01,17.889,8.306,1.896 +2020-02-01,17.643,8.598,2.11 +2020-03-01,22.232,8.268,1.834 +2020-04-01,25.234,7.421,1.963 +2020-05-01,26.835,7.166,1.287 +2020-06-01,29.708,6.846,2.856 +2020-07-01,31.421,6.105,2.615 +2020-08-01,31.424,6.069,1.782 +2020-09-01,29.677,6.694,2.222 +2020-10-01,26.943,7.319,2.663 +2020-11-01,23.219,7.944,3.103 +2020-12-01,17.062,8.936,1.767 +2021-01-01,15.935,9.931,1.603 +2021-02-01,17.782,8.951,0.755 +2021-03-01,21.266,8.249,2.182 +2021-04-01,23.866,7.896,2.803 +2021-05-01,27.235,7.307,2.384 +2021-06-01,29.275,6.722,2.856 +2021-07-01,30.214,6.569,2.763 +2021-08-01,31.061,6.431,5.04 +2021-09-01,29.795,6.226,6.85 +2021-10-01,27.068,6.384,5.002 +2021-11-01,20.028,7.736,1.861 +2021-12-01,19.067,8.371,1.541 +2022-01-01,16.626,8.683,1.076 +2022-02-01,17.129,9.446,1.177 +2022-03-01,21.452,8.336,1.184 +2022-04-01,24.654,7.715,1.295 +2022-05-01,28.374,7.357,1.301 +2022-06-01,31.134,7.136,1.239 +2022-07-01,31.668,6.914,1.177 +2022-08-01,31.726,6.693,1.115 +2022-09-01,29.816,6.471,1.052 +2022-10-01,24.809,7.276,1.077 +2022-11-01,22.007,7.85,1.164 +2022-12-01,18.61,8.861,1.499 +2023-01-01,16.637,10.069,1.33 +2023-02-01,19.889,9.46,1.381 +2023-03-01,22.851,8.084,1.652 +2023-04-01,24.992,7.731,1.462 +2023-05-01,27.254,7.274,1.656 +2023-06-01,29.074,7.533,2.169 +2023-07-01,31.946,6.798,2.464 +2023-08-01,32.384,5.759,4.632 +2023-09-01,30.506,5.692,4.566 +2023-10-01,25.356,6.643,3.425 +2023-11-01,21.345,7.603,2.764 +2023-12-01,17.428,9.219,2.382 +2024-01-01,15.846,10.074,2.385 +2024-02-01,16.977,9.081,1.159 +2024-03-01,21.789,8.363,1.257 +2024-04-01,24.408,7.734,1.539 +2024-05-01,28.731,6.965,1.759 +2024-06-01,30.753,6.698,2.762 +2024-07-01,31.659,6.174,3.564 +2024-08-01,31.272,6.329,3.799 +2024-09-01,29.907,5.624,5.108 +2024-10-01,25.864,6.847,4.115 +2024-11-01,23.521,6.957,3.69 +2024-12-01,17.272,8.28,1.544 +2025-01-01,14.236,6.04,1.662 +2025-02-01,18.754,8.868,0.635 +2025-03-01,19.989,8.627,1.672 +2025-04-01,25.806,7.763,1.152 +2025-05-01,28.934,6.954,1.281 +2025-06-01,31.018,6.529,1.509 +2025-07-01,31.654,6.635,1.912 +2025-08-01,32.337,6.224,3.655 +2025-09-01,29.175,6.52,3.687 +2025-10-01,25.67,6.858,3.069 +2025-11-01,20.663,7.983,2.045 +2025-12-01,17.942,8.757,2.105 +2026-01-01,16.273,9.586,1.481 +2026-02-01,15.79,9.838,1.366 +2026-03-01,21.961,8.543,1.223 +2026-04-01,24.691,8.079,1.271 +2026-05-01,28.03,7.481,1.741 +2026-06-01,31,7.168,2.571 diff --git a/data/make_barco_timeseries.R b/data/make_barco_timeseries.R new file mode 100644 index 0000000..c13caf3 --- /dev/null +++ b/data/make_barco_timeseries.R @@ -0,0 +1,57 @@ +# Build the monthly Lake Barco (BARC) limnology time series used by the +# "Time series modeling in R 3" assignment. +# +# Source: the aquatics target file from the NEON Ecological Forecasting +# Challenge, run by the Ecological Forecasting Initiative Research Coordination +# Network. https://projects.ecoforecast.org/neon4cast-docs/ +# +# The Challenge targets repackage two NEON data products for lake sites: +# DP1.20264.001 Temperature at specific depth in surface water -> temperature +# DP1.20288.001 Water quality -> oxygen, chla +# +# Thomas RQ, Boettiger C, Carey CC, et al. (2023) The NEON Ecological +# Forecasting Challenge. Front Ecol Environ 21(3): 112-113. doi:10.1002/fee.2616 +# +# Run from the repository root: Rscript data/make_barco_timeseries.R + +library(readr) +library(dplyr) +library(tidyr) +library(tsibble) +library(zoo) + +targets_url <- paste0( + "https://sdsc.osn.xsede.org/bio230014-bucket01/challenges/targets/", + "project_id=neon4cast/duration=P1D/aquatics-targets.csv.gz" +) + +# The window is fixed so that everyone works with the same 92 months. It starts +# in Nov 2018, the first month with all three sensors reporting, and ends in +# Jun 2026. Extend `end_month` to pick up newer data. +start_month <- yearmonth("2018 Nov") +end_month <- yearmonth("2026 Jun") + +barco <- read_csv(targets_url, show_col_types = FALSE) |> + filter(site_id == "BARC", variable %in% c("temperature", "oxygen", "chla")) |> + mutate(month = yearmonth(datetime)) |> + group_by(month, variable) |> + summarize(value = mean(observation, na.rm = TRUE), .groups = "drop") |> + pivot_wider(names_from = variable, values_from = value) |> + filter(month >= start_month, month <= end_month) |> + as_tsibble(index = month) |> + # A handful of months have no data at all because the sonde was out of the + # water. Interpolate them so that the series is gap free: fable's AR() errors + # out on missing values, and the assignment asks students to fit one. + fill_gaps() |> + mutate(across(c(temperature, oxygen, chla), \(x) round(na.approx(x), 3))) + +stopifnot( + nrow(barco) == 92, + !any(is.na(barco)) +) + +barco |> + as_tibble() |> + mutate(date = as.Date(month)) |> + select(date, temperature, oxygen, chla) |> + write_csv("data/barco_oxygen_timeseries.csv") diff --git a/lessons/R-time-series-modeling-3/assignment.qmd b/lessons/R-time-series-modeling-3/assignment.qmd new file mode 100644 index 0000000..cc54a03 --- /dev/null +++ b/lessons/R-time-series-modeling-3/assignment.qmd @@ -0,0 +1,136 @@ +--- +title: Assignment +order: 5 +--- + +## Overview + +This assignment covers the material from **Time series modeling in R 1, 2, and 3**. +You will fit the same family of models we fit in the tutorials -- white noise, autoregressive, +seasonal ARIMA, time-series linear models, and dynamic regression -- to a new ecological +time series, and evaluate each one the way we did in class. + +Work through the parts in order, which will expand the modeling approach +gradually fixing problems you diagnose as you go. + +Turn in a single R script (or Quarto document) on Canvas that contains your code, +generates your plots, +and includes written answers as comments (R script) or text block (Quarto document). + +::: {.callout-note} + +**Learning Objectives:** + +* Fit white noise, AR, ARIMA, TSLM, and dynamic regression models to a new dataset +* Use residual to check if a model has captured the structure in a time series +* Interpret model coefficients + +::: + +## Data + +We are exploring dissolved oxygen from **Lake Barco**, +the National Ecological Observatory Network (NEON) aquatic site at the Ordway-Swisher Biological Station in Florida. + +Dissolved oxygen is controlled by a combination of temperature and algae and serves as a good indicator of lake condition. + +The data is monthly and includes three variables: + +* `oxygen` -- dissolved oxygen concentration (mg/L), our response variable +* `temperature` -- water temperature (°C) +* `chla` -- chlorophyll-a (µg/L), a measure of how much algae is in the water + +### Where the data come from + +The data are NEON sensor measurements, obtained through the +[NEON Ecological Forecasting Challenge](https://projects.ecoforecast.org/neon4cast-docs/), +which is run by the Ecological Forecasting Initiative Research Coordination Network and +repackages the raw NEON products into tidy daily "target" files. + +Two NEON data products underlie the three columns: + +* [Temperature at specific depth in surface water (DP1.20264.001)](https://data.neonscience.org/data-products/DP1.20264.001) -- `temperature` +* [Water quality (DP1.20288.001)](https://data.neonscience.org/data-products/DP1.20288.001) -- `oxygen` and `chla` + +Both are collected at the [Lake Barco (BARC) field site](https://www.neonscience.org/field-sites/barc). + +If you want to see how the daily target file was reduced to the monthly series -- +or to rebuild it with more recent data -- the script is +[`make_barco_timeseries.R`](/data/make_barco_timeseries.R). + +> Thomas RQ, Boettiger C, Carey CC, et al. (2023) The NEON Ecological Forecasting Challenge. +> *Frontiers in Ecology and the Environment* 21(3): 112--113. +> [doi:10.1002/fee.2616](https://doi.org/10.1002/fee.2616) + +## Part 1: Look at the data + + +1. Download the [Lake Barco limnology data](/data/barco_oxygen_timeseries.csv) +2. Plot the oxygen time series along with its ACF and PACF using `gg_tsdisplay()`. +3. Briefly describe what you see. E.g., Is there a trend? Is there a seasonal signal? At what lags is the + autocorrelation strongest? + +## Part 2: White noise + +1. Fit a white noise model to `oxygen` using `MEAN()`. +2. Report the model with `report()`. +3. Plot the data with the fitted values on top using `augment()` and `autoplot()`. +4. Plot the residuals with `gg_tsresiduals()`. + +**Question:** What structure in the data is this model failing to capture and how do you now? + +## Part 3: Autoregressive models + +1. Based on the PACF from Part 1, pick the number of autoregressive terms and fit an `AR()` + model to `oxygen`. +2. Report the model. +3. Plot the fitted values over the data, and plot the residuals. + +**Question:** Describe how you chose the AR order. +Describe how this month's oxygen relates to last month's. +Compared to the white noise model, what did the AR model improve. +What structure is still left in the residuals? + +## Part 4: Seasonal ARIMA + +1. Fit an `ARIMA()` model to `oxygen`, letting `fable` select the structure. +2. Report the model. +3. Plot the fitted values over the data, and plot the residuals. + +**Question:** What model structure was selected +(How many AR terms, how many MA terms, how much differencing, and what seasonal terms). +How did the the residual ACF change compared to the AR model. + +## Part 5: Adding covariates with TSLM + +1. Fit a `TSLM()` model predicting `oxygen` from `temperature`. +2. Report the model. +3. Add `chla` as a second predictor and refit. +4. Plot the fitted values compared to the data, +and plot the residuals of the model with both `temperature and `chla` as drivers. + +**Question:** Describe the model coefficients including sign and significance. +The `TSLM()` model doesn't explicitly model seasonality, butit tracks the seasonal +ups and downs well. How is it doing that? Are there any issues with the residuals? + +## Part 6: Dynamic regression + +1. Fit an `ARIMA()` model with `temperature` as an external covariate. +2. Report the model. +3. Fit a second version that also includes `chla`. +4. Compare the two models using the AICc values (lower is better) in their reports and + plot the fitted values and residuals for your preferred model. + +**Question:** What error structure did `fable` choose for the two models? +Is the temperature coefficient similar to the one from the `TSLM()` model? +Compare the residuals of this model to the residuals from the `TSLM()` model in +Part 5. What did adding the ARIMA error structure change, and why does it matter for +interpreting the coefficients? + +## Part 7: Wrap-up + +Include a short paragraph the answers: + +* Which model would you choose to describe how dissolved oxygen behaves in Lake Barco, and why? +* One model suggests that oxygen is seasonal; another suggests why it is seasonal. + Which is which, and when might you prefer each one? diff --git a/lessons/R-time-series-modeling-3/assignment_solution.qmd b/lessons/R-time-series-modeling-3/assignment_solution.qmd new file mode 100644 index 0000000..cbf2914 --- /dev/null +++ b/lessons/R-time-series-modeling-3/assignment_solution.qmd @@ -0,0 +1,166 @@ +--- +title: "Assignment Solution: Time series modeling in R" +subtitle: Instructor answer key -- not published to the course website +format: + html: + toc: true + toc-depth: 2 + embed-resources: true + fig-width: 9 + fig-height: 6 +execute: + enabled: true + warning: false + message: false +--- + +::: {.callout-warning} + +This file is excluded from the site render in `_quarto.yml` so that it is never +published alongside the assignment. Render it by hand when you want to refresh it: + +```sh +quarto render lessons/R-time-series-modeling-3/assignment_solution.qmd +``` + +Unlike the rest of the site, this document **executes** its R code, so it needs +`tsibble`, `fable`, `feasts`, `ggtime`, and `dplyr` installed. + +::: + +```{r} +#| label: setup +library(tsibble) +library(fable) +library(feasts) +library(ggtime) +library(dplyr) +``` + +## Part 1: Look at the data + +```{r} +#| label: load +barco <- read.csv("../../data/barco_oxygen_timeseries.csv") |> + mutate(month = yearmonth(date)) |> + as_tsibble(index = month) +``` + +```{r} +#| label: tsdisplay +gg_tsdisplay(barco, oxygen) +``` + +## Part 2: White noise + +```{r} +#| label: mean-model +avg_model <- model(barco, MEAN(oxygen)) +report(avg_model) +``` + +```{r} +#| label: mean-fit +avg_aug <- augment(avg_model) +autoplot(avg_aug, oxygen) + + autolayer(avg_aug, .fitted, color = "red") +``` + +```{r} +#| label: mean-resid +gg_tsresiduals(avg_model) +``` + +## Part 3: Autoregressive models + +```{r} +#| label: ar-model +ar_model <- model(barco, AR(oxygen ~ order(2))) +report(ar_model) +``` + +```{r} +#| label: ar-fit +ar_aug <- augment(ar_model) +autoplot(ar_aug, oxygen) + + autolayer(ar_aug, .fitted, color = "orange") +``` + +```{r} +#| label: ar-resid +gg_tsresiduals(ar_model) +``` + +## Part 4: Seasonal ARIMA + +```{r} +#| label: arima-model +arima_model <- model(barco, ARIMA(oxygen)) +report(arima_model) +``` + +```{r} +#| label: arima-fit +arima_aug <- augment(arima_model) +autoplot(arima_aug, oxygen) + + autolayer(arima_aug, .fitted, color = "orange") +``` + +```{r} +#| label: arima-resid +gg_tsresiduals(arima_model) +``` + +## Part 5: Adding covariates with TSLM + +```{r} +#| label: tslm-temp +tslm_temp <- model(barco, TSLM(oxygen ~ temperature)) +report(tslm_temp) +``` + +```{r} +#| label: tslm-both +tslm_model <- model(barco, TSLM(oxygen ~ temperature + chla)) +report(tslm_model) +``` + +```{r} +#| label: tslm-fit +tslm_aug <- augment(tslm_model) +autoplot(tslm_aug, oxygen) + + autolayer(tslm_aug, .fitted, color = "orange") +``` + +```{r} +#| label: tslm-resid +gg_tsresiduals(tslm_model) +``` + +## Part 6: Dynamic regression + +```{r} +#| label: arimax-temp +arimax_temp <- model(barco, ARIMA(oxygen ~ temperature)) +report(arimax_temp) +``` + +```{r} +#| label: arimax-both +arimax_model <- model(barco, ARIMA(oxygen ~ temperature + chla)) +report(arimax_model) +``` + +```{r} +#| label: arimax-fit +arimax_aug <- augment(arimax_model) +autoplot(arimax_aug, oxygen) + + autolayer(arimax_aug, .fitted, color = "orange") +``` + +```{r} +#| label: arimax-resid +gg_tsresiduals(arimax_model) +``` + +## Part 7: Wrap-up diff --git a/lessons/R-time-series-modeling-3/index.qmd b/lessons/R-time-series-modeling-3/index.qmd index fcc3f3d..358d4b7 100644 --- a/lessons/R-time-series-modeling-3/index.qmd +++ b/lessons/R-time-series-modeling-3/index.qmd @@ -7,6 +7,7 @@ listing: - material.qmd - r_tutorial.qmd - equation_review.qmd + - assignment.qmd type: table fields: - title