Skip to content
/ tvsews Public

Time vs. Space Early Warning Statistics

License

Unknown, MIT licenses found

Licenses found

Unknown
LICENSE
MIT
LICENSE.md
Notifications You must be signed in to change notification settings

cbuelo/tvsews

Repository files navigation

tvsews

tvsews is an R package that provides code and data for a whole-lake fertilization experiment that caused an algal bloom in order to compared the performance of early warning statistics (EWS) in temporal and spatial data.

Installation

You can install tvsews from GitHub with:

# install.packages("devtools")
devtools::install_github("cbuelo/tvsews")

To install the package and build a vignette that carries out analyses and displays results for an accepted manuscript (note this will take a few minutes to run):

devtools::install_github("cbuelo/tvsews", build_vignettes = T)
library(tvsews)
vignette("TvS-AlgalBloom-EWS")

Example

Example results from most recent experiment:

Set up and visualize time series data

library(tvsews)
library(dplyr)
library(ggplot2)
ex_data <- ts_data %>%
  filter(Year >= 2018) %>%
  select(Lake, Year, DOY, BGA_HYLB)

ex_bloom_fert_dates <- bloom_fert_dates %>%
  filter(Year >= 2018 & Lake == "R")

plot_state_vars(
  time_series = ex_data,
  bloom_fert_df = ex_bloom_fert_dates,
  var_rename_vec = c(`Phycocyanin (cells/mL)` = "BGA_HYLB"),
  legend_location = c(0.18, 0.7)
)

Calculate rolling window statistics and quickest detection alarms

ex_rw_stats <- calc_rolling_stats(
  data = ex_data,
  var_cols = "BGA_HYLB",
  widths = c(21)
) # warnings are missing data, okay
ex_qd_all <- run_qd(
  rolling_window_stats = ex_rw_stats,
  var_cols = "BGA_HYLB",
  widths = c(21),
  stats_to_qd = c("SD", "Ar1"),
  exp_lakes = c("R"),
  ref_lake = "L",
  ar1_alarm_rho = 0.95
)
#> [1] "Current variable QD-ing: BGA_HYLB"
ex_qd_formatted <- format_qd(
  qd_stats = ex_qd_all,
  var_cols = c("BGA_HYLB")
)
# plot the rolling window statistics and alarms
time_series_plot <- plot_temporal_EWS(
  rolling_window_stats = ex_rw_stats,
  qd_alarms = ex_qd_formatted,
  bloom_fert_df = ex_bloom_fert_dates,
  var_rename_vec = c(`Phyco` = "BGA_HYLB"),
  legend_location = c(0.2, 0.7)
)

Calculate and plot true and false alarm rates

alarm_rates <- calc_alarm_rates(
  qd_alarms = ex_qd_formatted,
  bloom_fert_df = bloom_fert_dates
)

plot_alarm_rates(
  qd_alarm_rates = alarm_rates,
  var_rename_vec = c("Phycocyanin" = "BGA_HYLB"),
  legend_title = "Positive Alarm Rate"
)

Plot example spatial data

plot_spatial_data(
  spatial_data = flame_data,
  samples_to_plot = data.frame(Year = 2019, DOY = c(165, 210, 228), stringsAsFactors = FALSE),
  var_cols = c("BGApc_ugL_tau")
)

Calculate spatial statistics

Just dates in pre-bloom fertilization period to limit computation time.

ex_spatial_data_pbfp <- flame_data %>%
  filter(Year == 2019 & between(DOY, 161, 201))

spatial_stats_pbfp <- calc_spatial_stats(
  spatial_data = ex_spatial_data_pbfp,
  var_cols = c("BGApc_ugL_tau")
) # ignore estimated time as not using all of default data

Plot spatial statistics

spatial_results_plot <- plot_spatial_EWS(spatial_stats_pbfp)

(Note to self: devtools::build_readme()to render)

About

Time vs. Space Early Warning Statistics

Resources

License

Unknown, MIT licenses found

Licenses found

Unknown
LICENSE
MIT
LICENSE.md

Stars

Watchers

Forks

Packages

No packages published

Languages