R/plotTsAnom.R

Defines functions plotTsAnom

Documented in plotTsAnom

#' Anomaly plot of time series
#' 
#' Series are illustrated by vertical lines extending from individual data
#' values to the long-term mean. The axes are not scaled in any way. Anomaly
#' plots are useful for visualizing shifts in time series levels.
#' 
#' Options are passed to the underlying \code{facet_wrap} function in
#' \pkg{ggplot2}. The main ones of interest are \code{ncol} for setting the
#' number of plotting columns and \code{scales = "free_y"} for allowing the y
#' scales of the different plots to be independent.
#' @author 
#' Alan Jassby, James Cloern
#' @param x matrix or vector time series
#' @param xlab optional x-axis label
#' @param ylab optional y-axis label
#' @param strip.labels labels for individual time series plots
#' @param ...  additional options
#' @return A plot and corresponding object of class \dQuote{ggplot}.
#' @seealso \code{\link{plotTs}}
#' @keywords Graphics ts
#' @export
#' @examples
#' 
#' # Spring bloom size for 6 stations in SF Bay
#' bloom <- aggregate(sfbayChla[, 1:6], 1, meanSub, sub=3:5)
#' plotTsAnom(bloom, ylab = 'Chl-a', strip.labels = paste('Station',
#'   substring(colnames(bloom), 2, 3)), ncol = 2, scales = "free_y")
#' 
plotTsAnom <-
function(x, xlab = NULL, ylab = NULL,
         strip.labels = colnames(x), ...) {

  # Validate arguments
  if (!is.ts(x)) stop("x must be of class 'ts'")
  if (missing(xlab)) xlab <- ""
  if (missing(ylab)) ylab <- ""

  if (is.matrix(x)) {  # a matrix time series

    # fill possible spaces in column names so melt+merge will work
    strip.labels <- strip.labels
    colnames(x) <- gsub(' ', '.', colnames(x))

    # Create data frame
    x.mean <- apply(x, 2, mean, na.rm = TRUE)
    x.mean.df <- data.frame(variable = factor(names(x.mean)), x.mean)
    d <- data.frame(time=as.Date(time(x)), x)
    d1 <- melt(d, id = 'time')
    d2 <- merge(d1, x.mean.df)
    d3 <- within(d2, variable <- factor(variable, levels = levels(variable),
                                        labels = strip.labels))
    d3 <- na.omit(d3)
    d3$ymin. <- with(d3, ifelse(value >= x.mean, x.mean, value))
    d3$ymax. <- with(d3, ifelse(value >= x.mean, value, x.mean))
    d3$colour. <- with(d3, value >= x.mean)
    # Plot
    ggplot(d3, aes_string(x="time", y="value", ymin="ymin.", ymax="ymax.",
                          colour="colour.")) +
      geom_linerange() +
      geom_hline(aes(yintercept = x.mean), size = 0.25) +
      labs(x = xlab, y = ylab) +
      facet_wrap(~ variable, ...) +
      theme(legend.position='none', panel.grid.minor = element_blank(),
        axis.text.x = element_text(angle=45, colour="grey50"))

  } else {  # a vector time series

    # Create data frame
    x.mean <- mean(x, na.rm = TRUE)
    d1 <- data.frame(time = as.Date(time(x)), x = as.numeric(x), x.mean)
    d1 <- na.omit(d1)
    d1$ymin. <- with(d1, ifelse(x >= x.mean, x.mean, x))
    d1$ymax. <- with(d1, ifelse(x >= x.mean, x, x.mean))
    d1$colour. <- with(d1, x >= x.mean)
    # Plot
    ggplot(d1, aes_string(x="time", y="x", ymin="ymin.", ymax="ymax.",
                          colour="colour.")) +
      geom_linerange() +
      geom_hline(aes(yintercept = x.mean), size = 0.25) +
      labs(x = xlab, y = ylab) +
      theme(legend.position='none', panel.grid.minor = element_blank(),
            axis.text.x = element_text(angle=45, colour="grey50"))
  }
}

Try the wql package in your browser

Any scripts or data that you put into this service are public.

wql documentation built on Aug. 10, 2022, 5:06 p.m.