remove magrittr dependency

DESCRIPTION

@@ -24,7 +24,7 @@ Description: Efficient methods for Bayesian inference of state space models
     models and discretised diffusion models are supported. 
     See Helske and Vihola (2021, <doi:10.32614/RJ-2021-103>) for details.
 License: GPL (>= 2)
-Depends: R (>= 3.5.0)
+Depends: R (>= 4.1.0)
 Suggests: 
     covr,
     ggplot2 (>= 2.0.0),
@@ -37,7 +37,6 @@ Suggests:
     sitmo,
     testthat
 Imports: 
-    magrittr,
     bayesplot,
     checkmate,
     coda (>= 0.18-1),

NAMESPACE

@@ -101,7 +101,6 @@ importFrom(dplyr,as_tibble)
 importFrom(dplyr,group_by)
 importFrom(dplyr,summarise)
 importFrom(dplyr,ungroup)
-importFrom(magrittr,"%>%")
 importFrom(posterior,as_draws)
 importFrom(posterior,as_draws_df)
 importFrom(posterior,default_convergence_measures)

NEWS.md

@@ -9,6 +9,8 @@ bssm 2.0.2 (Release date: 2023-09-04)
    * Added a default plot method for the `run_mcmc` output.
    * Fixed the aliases of the main help page to accomodate changes in roxygen2.
    * Removed explicit C++ version requirement as required by new CRAN policies.
+   * Removed `magrittr` dependency and switched to native pipe, leading to 
+     requirement for R 4.1.0+.
 
 bssm 2.0.1 (Release date: 2022-05-02)
 ==============

R/fitted.R

@@ -5,7 +5,6 @@
 #' 
 #' @export
 #' @importFrom stats fitted
-#' @importFrom magrittr %>%
 #' @importFrom dplyr group_by ungroup summarise as_tibble
 #' @importFrom diagis weighted_quantile weighted_var weighted_mean weighted_se
 #' @name fitted.mcmc_output
@@ -110,13 +109,13 @@ fitted.mcmc_output <- function(object, model,
     Variable = variables,
     Time = rep(time(model$y), each = nrow(pred)))
 
-  d %>% dplyr::group_by(.data$Variable, .data$Time) %>%
+  d  |>  dplyr::group_by(.data$Variable, .data$Time) |>
     dplyr::summarise(
       Mean = weighted_mean(.data$value, w),
       SD = sqrt(weighted_var(.data$value, w)),
       dplyr::as_tibble(as.list(weighted_quantile(.data$value, w, 
         probs = probs))),
-      "SE(Mean)" = as.numeric(sqrt(asymptotic_var(.data$value, w)))) %>% 
+      "SE(Mean)" = as.numeric(sqrt(asymptotic_var(.data$value, w)))) |> 
     dplyr::ungroup()
 }
 

R/models.R

@@ -898,8 +898,8 @@ bsm_lg <- function(y, sd_y, sd_level, sd_slope, sd_seasonal,
 #'
 #' # Traceplot using as.data.frame method for MCMC output
 #' library("dplyr")
-#' as.data.frame(mcmc_out) %>% 
-#'   filter(variable == "sd_level") %>% 
+#' as.data.frame(mcmc_out) |> 
+#'   filter(variable == "sd_level") |> 
 #'   ggplot(aes(y = value, x = iter)) + geom_line()
 #'   
 #' }

R/post_correction.R

@@ -210,34 +210,34 @@ suggest_N <- function(model, theta,
 #' # latent state
 #' library("dplyr")
 #' library("ggplot2")
-#' state_approx <- as.data.frame(out_approx, variable = "states") %>% 
-#'   group_by(time) %>%
+#' state_approx <- as.data.frame(out_approx, variable = "states") |> 
+#'   group_by(time) |>
 #'   summarise(mean = mean(value))
 #'   
-#' state_exact <- as.data.frame(out_is2, variable = "states") %>% 
-#'   group_by(time) %>%
+#' state_exact <- as.data.frame(out_is2, variable = "states") |> 
+#'   group_by(time) |>
 #'   summarise(mean = weighted.mean(value, weight))
 #' 
 #' dplyr::bind_rows(approx = state_approx, 
-#'   exact = state_exact, .id = "method") %>%
-#'   filter(time > 200) %>%
+#'   exact = state_exact, .id = "method") |>
+#'   filter(time > 200) |>
 #' ggplot(aes(time, mean, colour = method)) + 
 #'   geom_line() + 
 #'   theme_bw()
 #'
 #' # posterior means
 #' p_approx <- predict(out_approx, model, type = "mean", 
-#'   nsim = 1000, future = FALSE) %>% 
-#'   group_by(time) %>%
+#'   nsim = 1000, future = FALSE) |> 
+#'   group_by(time) |>
 #'   summarise(mean = mean(value))
 #' p_exact <- predict(out_is2, model, type = "mean", 
-#'   nsim = 1000, future = FALSE) %>% 
-#'   group_by(time) %>%
+#'   nsim = 1000, future = FALSE) |> 
+#'   group_by(time) |>
 #'   summarise(mean = mean(value))
 #' 
 #' dplyr::bind_rows(approx = p_approx, 
-#'   exact = p_exact, .id = "method") %>%
-#'   filter(time > 200) %>%
+#'   exact = p_exact, .id = "method") |>
+#'   filter(time > 200) |>
 #' ggplot(aes(time, mean, colour = method)) + 
 #'   geom_line() + 
 #'   theme_bw() 

R/predict.R

@@ -52,35 +52,35 @@
 #'   nsim = 1000)
 #' 
 #' library("dplyr")
-#' sumr_fit <- as.data.frame(mcmc_results, variable = "states") %>%
-#'   group_by(time, iter) %>% 
+#' sumr_fit <- as.data.frame(mcmc_results, variable = "states") |>
+#'   group_by(time, iter) |> 
 #'   mutate(signal = 
 #'       value[variable == "level"] + 
-#'       value[variable == "seasonal_1"]) %>%
-#'   group_by(time) %>%
+#'       value[variable == "seasonal_1"]) |>
+#'   group_by(time) |>
 #'   summarise(mean = mean(signal), 
 #'     lwr = quantile(signal, 0.025), 
 #'     upr = quantile(signal, 0.975))
 #' 
-#' sumr_pred <- pred %>% 
-#'   group_by(time, sample) %>%
+#' sumr_pred <- pred |> 
+#'   group_by(time, sample) |>
 #'   mutate(signal = 
 #'       value[variable == "level"] + 
-#'       value[variable == "seasonal_1"]) %>%
-#'   group_by(time) %>%
+#'       value[variable == "seasonal_1"]) |>
+#'   group_by(time) |>
 #'   summarise(mean = mean(signal),
 #'     lwr = quantile(signal, 0.025), 
 #'     upr = quantile(signal, 0.975)) 
 #'     
 #' # If we used type = "mean", we could do
-#' # sumr_pred <- pred %>% 
-#' #   group_by(time) %>%
+#' # sumr_pred <- pred |> 
+#' #   group_by(time) |>
 #' #   summarise(mean = mean(value),
 #' #     lwr = quantile(value, 0.025), 
 #' #     upr = quantile(value, 0.975)) 
 #'     
 #' library("ggplot2")
-#' rbind(sumr_fit, sumr_pred) %>% 
+#' rbind(sumr_fit, sumr_pred) |> 
 #'   ggplot(aes(x = time, y = mean)) + 
 #'   geom_ribbon(aes(ymin = lwr, ymax = upr), 
 #'    fill = "#92f0a8", alpha = 0.25) +
@@ -96,23 +96,23 @@
 #' meanrep <- predict(mcmc_results, model = model, type = "mean", 
 #'   future = FALSE, nsim = 1000)
 #'   
-#' sumr_yrep <- yrep %>% 
-#'   group_by(time) %>%
+#' sumr_yrep <- yrep |> 
+#'   group_by(time) |>
 #'   summarise(earnings = mean(value),
 #'     lwr = quantile(value, 0.025), 
-#'     upr = quantile(value, 0.975)) %>%
+#'     upr = quantile(value, 0.975)) |>
 #'   mutate(interval = "Observations")
 #'
-#' sumr_meanrep <- meanrep %>% 
-#'   group_by(time) %>%
+#' sumr_meanrep <- meanrep |> 
+#'   group_by(time) |>
 #'   summarise(earnings = mean(value),
 #'     lwr = quantile(value, 0.025), 
-#'     upr = quantile(value, 0.975)) %>%
+#'     upr = quantile(value, 0.975)) |>
 #'   mutate(interval = "Mean")
 #'     
-#' rbind(sumr_meanrep, sumr_yrep) %>% 
+#' rbind(sumr_meanrep, sumr_yrep) |> 
 #'   mutate(interval = 
-#'     factor(interval, levels = c("Observations", "Mean"))) %>%
+#'     factor(interval, levels = c("Observations", "Mean"))) |>
 #'   ggplot(aes(x = time, y = earnings)) + 
 #'   geom_ribbon(aes(ymin = lwr, ymax = upr, fill = interval), 
 #'    alpha = 0.75) +

R/run_mcmc.R

@@ -325,17 +325,17 @@ run_mcmc.lineargaussian <- function(model, iter, output_type = "full",
 #' library("ggplot2")
 #'
 #'  # compute summary statistics
-#' level_sumr <- d_states %>%
-#'   filter(variable == "level") %>%
-#'   group_by(time) %>%
+#' level_sumr <- d_states |>
+#'   filter(variable == "level") |>
+#'   group_by(time) |>
 #'   summarise(mean = diagis::weighted_mean(value, weight),
 #'     lwr = diagis::weighted_quantile(value, weight,
 #'       0.025),
 #'     upr = diagis::weighted_quantile(value, weight,
 #'       0.975))
 #'
 #' # visualize
-#' level_sumr %>% ggplot(aes(x = time, y = mean)) +
+#' level_sumr |> ggplot(aes(x = time, y = mean)) +
 #'   geom_line() +
 #'   geom_line(aes(y = lwr), linetype = "dashed", na.rm = TRUE) +
 #'   geom_line(aes(y = upr), linetype = "dashed", na.rm = TRUE) +
@@ -376,8 +376,8 @@ run_mcmc.lineargaussian <- function(model, iter, output_type = "full",
 #' # Note small number of iterations for CRAN checks
 #' out <- run_mcmc(model, iter = 4000, mcmc_type = "approx")
 #'
-#' sumr <- as.data.frame(out, variable = "states") %>%
-#'   group_by(time) %>% mutate(value = exp(value)) %>%
+#' sumr <- as.data.frame(out, variable = "states") |>
+#'   group_by(time) |> mutate(value = exp(value)) |>
 #'   summarise(mean = mean(value),
 #'     ymin = quantile(value, 0.05), ymax = quantile(value, 0.95))
 #' ggplot(sumr, aes(time, mean)) +

R/summary.R

@@ -87,9 +87,9 @@ summary.mcmc_output <- function(object, return_se = FALSE, variable = "theta",
   }
   if (variable %in% c("theta", "both")) {
     sumr_theta <- 
-      as.data.frame(object, variable = "theta", expand = TRUE) %>%
-      group_by(.data$variable) %>% 
-      summarise(as_tibble(as.list(summary_f(.data$value, .data$weight)))) %>% 
+      as.data.frame(object, variable = "theta", expand = TRUE) |>
+      group_by(.data$variable) |> 
+      summarise(as_tibble(as.list(summary_f(.data$value, .data$weight)))) |> 
       as.data.frame()
     if (variable == "theta") return(sumr_theta)
   }
@@ -117,9 +117,9 @@ summary.mcmc_output <- function(object, return_se = FALSE, variable = "theta",
 
     sumr_states <- 
       as.data.frame(object, variable = "states", expand = TRUE, 
-        times = times, states = states, use_times = use_times) %>%
-      group_by(.data$variable, .data$time) %>% 
-      summarise(as_tibble(as.list(summary_f(.data$value, .data$weight)))) %>% 
+        times = times, states = states, use_times = use_times) |>
+      group_by(.data$variable, .data$time) |> 
+      summarise(as_tibble(as.list(summary_f(.data$value, .data$weight)))) |> 
       as.data.frame()
     if (variable == "states") return(sumr_states)
   }

README.Rmd

@@ -137,7 +137,7 @@ set.seed(1)
 data("airquality", package = "datasets")
 
 # Covariates as matrix. For complex cases, check out as_bssm function
-xreg <- airquality %>% select(Wind, Temp) %>% as.matrix()
+xreg <- airquality |> select(Wind, Temp) |> as.matrix()
 
 model <- bsm_lg(airquality$Ozone,
   xreg = xreg,  
@@ -152,10 +152,10 @@ fit <- run_mcmc(model, iter = 20000, burnin = 5000)
 fit
 
 obs <- data.frame(Time = 1:nrow(airquality),
-  Ozone = airquality$Ozone) %>% filter(!is.na(Ozone))
+  Ozone = airquality$Ozone) |> filter(!is.na(Ozone))
 
 pred <- fitted(fit, model)
-pred %>%
+pred |>
   ggplot(aes(x = Time, y = Mean)) + 
   geom_ribbon(aes(ymin = `2.5%`, ymax = `97.5%`), 
     alpha = 0.5, fill = "steelblue") + 
@@ -185,7 +185,7 @@ Comparison:
 ```{r compare}
 pred2 <- fitted(fit2, model2)
 
-bind_rows(list(Gaussian = pred, Gamma = pred2), .id = "Model") %>%
+bind_rows(list(Gaussian = pred, Gamma = pred2), .id = "Model") |>
   ggplot(aes(x = Time, y = Mean)) + 
   geom_ribbon(aes(ymin = `2.5%`, ymax = `97.5%`, fill = Model), 
     alpha = 0.25) + 
@@ -214,7 +214,7 @@ We combine these two models as a bivariate Gaussian model with `ssm_mlg`:
 
 ```{r missing-values}
 # predictors (not including solar radiation) for ozone
-xreg <- airquality %>% select(Wind, Temp) %>% as.matrix()
+xreg <- airquality |> select(Wind, Temp) |> as.matrix()
 
 # Function which outputs new model components given the parameter vector theta
 update_fn <- function(theta) {
@@ -252,9 +252,9 @@ Draw predictions:
 pred <- fitted(fit, model)
 
 obs <- data.frame(Time = 1:nrow(airquality),
-  Solar = airquality$Solar.R) %>% filter(!is.na(Solar))
+  Solar = airquality$Solar.R) |> filter(!is.na(Solar))
 
-pred %>% filter(Variable == "Solar") %>%
+pred |> filter(Variable == "Solar") |>
   ggplot(aes(x = Time, y = Mean)) + 
   geom_ribbon(aes(ymin = `2.5%`, ymax = `97.5%`), 
     alpha = 0.5, fill = "steelblue") + 
@@ -265,9 +265,9 @@ pred %>% filter(Variable == "Solar") %>%
 
 
 obs <- data.frame(Time = 1:nrow(airquality),
-  Ozone = airquality$Ozone) %>% filter(!is.na(Ozone))
+  Ozone = airquality$Ozone) |> filter(!is.na(Ozone))
 
-pred %>% filter(Variable == "Ozone") %>%
+pred |> filter(Variable == "Ozone") |>
   ggplot(aes(x = Time, y = Mean)) + 
   geom_ribbon(aes(ymin = `2.5%`, ymax = `97.5%`), 
     alpha = 0.5, fill = "steelblue") +