Fix incomplete symmetrize + document it

R/layer_cdc_flatline_quantiles.R

@@ -32,10 +32,20 @@
 #'   These samples are spaced evenly on the (0, 1) scale, F_X(x) resulting
 #'   in linear interpolation on the X scale. This is achieved with
 #'   [stats::quantile()] Type 7 (the default for that function).
+#' @param symmetrize Logical. If `TRUE`, does two things: (i) forces the
+#'   "empirical" CDF of residuals to be symmetric by pretending that for every
+#'   actually-observed residual X we also observed another residual -X, and (ii)
+#'   at each ahead, forces the median simulated value to be equal to the point
+#'   prediction by adding or subtracting the same amount to every simulated
+#'   value. Adjustments in (ii) take place before propagating forward and
+#'   simulating the next ahead. This forces any 1-ahead predictive intervals to
+#'   be symmetric about the point prediction, and encourages larger aheads to be
+#'   more symmetric.
 #' @param nonneg Logical. Force all predictive intervals be non-negative.
-#'   Because non-negativity is forced _before_ propagating forward, this
-#'   has slightly different behaviour than would occur if using
-#'   [layer_threshold()].
+#'   Because non-negativity is forced _before_ propagating forward, this has
+#'   slightly different behaviour than would occur if using [layer_threshold()].
+#'   Thresholding at each ahead takes place after any shifting from
+#'   `symmetrize`.
 #'
 #' @return an updated `frosting` postprocessor. Calling [predict()] will result
 #'   in an additional `<list-col>` named `.pred_distn_all` containing 2-column
@@ -232,6 +242,9 @@ slather.layer_cdc_flatline_quantiles <-
 propogate_samples <- function(
     r, p, quantile_levels, aheads, nsim, symmetrize, nonneg) {
   max_ahead <- max(aheads)
+  if (symmetrize) {
+    r <- c(r, -r)
+  }
   samp <- quantile(r, probs = c(0, seq_len(nsim - 1)) / (nsim - 1), na.rm = TRUE)
   res <- list()