added plot method for metabolism attribute, added metabolim and units…

… attribute to swmpr object

NAMESPACE

@@ -2,6 +2,7 @@
 
 S3method(comb,swmpr)
 S3method(ecometab,swmpr)
+S3method(plot_met,swmpr)
 S3method(setstep,swmpr)
 export(aggregate.swmpr)
 export(all_params)
@@ -22,6 +23,7 @@ export(na.approx.swmpr)
 export(param_names)
 export(parser)
 export(plot.swmpr)
+export(plot_met)
 export(plot_summary)
 export(plot_summary.swmpr)
 export(qaqc)

R/swmpr_analyze.R

@@ -414,8 +414,7 @@ decomp.swmpr <- function(swmpr_in, param, type = 'additive', frequency = 'daily'
 #' @param ... additional arguments passed to other methods, including \code{\link[wq]{decompTs}} 
 #' 
 #' @return  
-#' A \code{\link[ggplot2]{ggplot}} object if \code{vals_out = TRUE} (default), otherwise 
-#' a monthly time series matrix of class \code{\link[stats]{ts}}.
+#' A \code{\link[ggplot2]{ggplot}} object if \code{vals_out = TRUE} (default), otherwise a monthly time series matrix of class \code{\link[stats]{ts}}.
 #' 
 #' @details
 #' This function is a simple wrapper to the \code{\link[wq]{decompTs}} function in the wq package, also described in Cloern and Jassby (2010).  The function is similar to \code{\link{decomp.swmpr}} (which is a wrapper to \code{\link[stats]{decompose}}) with a few key differences.  The \code{\link{decomp.swmpr}} function decomposes the time series into a trend, seasonal, and random components, whereas the current function decomposes into the grandmean, annual, seasonal, and events components.  For both functions, the random or events components, respectively, can be considered anomalies that don't follow the trends in the remaining categories.  
@@ -861,7 +860,7 @@ plot_summary.swmpr <- function(swmpr_in, param, years = NULL, ...){
 #' 
 #' @param swmpr_in Input swmpr object which must include time series of dissolved oxygen, 
 #' @param depth_val numeric value for station depth if time series is not available
-#' @param units chr indicating units of output for oxygen, either as mmol or grams
+#' @param met_units chr indicating units of output for oxygen, either as mmol or grams
 #' @param trace logical indicating if progress is shown in the console
 #' @param ... arguments passed to other methods
 #' 
@@ -906,10 +905,9 @@ plot_summary.swmpr <- function(swmpr_in, param, years = NULL, ...){
 #' Thebault J, Schraga TS, Cloern JE, Dunlavey EG. 2008. Primary production and carrying capacity of former salt ponds after reconnection to San Francisco Bay. Wetlands. 28(3):841-851.
 #' 
 #' @seealso 
-#' \code{\link{comb}}
+#' \code{\link{comb}}, \code{\link{plot_met}}
 #' 
 #' @examples
-#'
 #' ## import water quality and weather data
 #' data(apadbwq)
 #' data(apaebmet)
@@ -927,7 +925,7 @@ plot_summary.swmpr <- function(swmpr_in, param, years = NULL, ...){
 #' res <- attr(res, 'metabolism')
 #' 
 #' ## output units in grams of oxygen
-#' res <- ecometab(dat, trace = TRUE, units = 'grams')
+#' res <- ecometab(dat, trace = TRUE, met_units = 'grams')
 #' res <- attr(res, 'metabolism')
 ecometab <- function(swmpr_in, ...) UseMethod('ecometab')
 
@@ -936,12 +934,12 @@ ecometab <- function(swmpr_in, ...) UseMethod('ecometab')
 #' @export
 #' 
 #' @method ecometab swmpr
-ecometab.swmpr <- function(swmpr_in, depth_val = NULL, units = 'mmol', trace = FALSE, ...){
+ecometab.swmpr <- function(swmpr_in, depth_val = NULL, met_units = 'mmol', trace = FALSE, ...){
 
   stat <- attr(swmpr_in, 'station')
 
   # stop if units not mmol or grams
-  if(any(!(grepl('mmol|grams', units))))
+  if(any(!(grepl('mmol|grams', met_units))))
     stop('Units must be mmol or grams')
 
   # stop if input data does not include wq and met
@@ -1092,7 +1090,7 @@ ecometab.swmpr <- function(swmpr_in, depth_val = NULL, units = 'mmol', trace = F
   row.names(out) <- 1:nrow(out)
 
   # change units to grams
-  if('grams' %in% units){
+  if('grams' %in% met_units){
 
     # convert metab data to g m^-2 d^-1
     # 1mmolO2 = 32 mg O2, 1000mg = 1g, multiply by 32/1000
@@ -1103,9 +1101,93 @@ ecometab.swmpr <- function(swmpr_in, depth_val = NULL, units = 'mmol', trace = F
 
   # append to metabolism attribute
   attr(swmpr_in, 'metabolism') <- out
+  attr(swmpr_in, 'met_units') <- met_units
 
   if(trace) tictoc::toc()
 
   return(swmpr_in)
 
+}
+
+######
+#' Plot ecosystem metabolism for a swmpr object
+#'
+#' Plot gross production, total respiration, and net ecosystem metabolism for a swmpr object. 
+#'
+#' @param swmpr_in input swmpr object
+#' @param pretty logical indicating use of predefined plot aesthetics
+#' @param ...
+#'
+#' @export
+#' 
+#' @import ggplot2
+#'
+#' @details 
+#' A plot will only be returned if the \code{metabolism} attribute for the \code{\link{swmpr}} object is not \code{NULL}.  The \code{\link{ecometab}} function is used to create metabolism estimates.  See the examples.
+#' 
+#' By default, \code{pretty = TRUE} will return a \code{\link[ggplot2]{ggplot}} object with predefined aesthetics.  Setting \code{pretty = FALSE} will return the plot with minimal modifications to the \code{\link[ggplot2]{object}}.  Use the latter approach for easier customization of the plot.  
+#' 
+#' @return 
+#' A \code{\link[ggplot2]{ggplot}} object which can be further modified.
+#' 
+#' @seealso 
+#' \code\link{ecometab}
+#' 
+#' @examples
+#' ## import water quality and weather data
+#' data(apadbwq)
+#' data(apaebmet)
+#' 
+#' ## qaqc, combine
+#' wq <- qaqc(apadbwq)
+#' met <- qaqc(apaebmet)
+#' dat <- comb(wq, met)
+#' 
+#' ## estimate metabolism
+#' res <- ecometab(dat, trace = TRUE)
+#' 
+#' ## plot
+#' plot_met(res)
+plot_met <- function(swmpr_in, ...) UseMethod('plot_met')
+
+#' @rdname plot_met
+#'
+#' @export
+#'
+#' @method plot_met swmpr
+plot_met.swmpr <- function(swmpr_in, pretty = TRUE, ...){
+
+  # get metabolism estimates
+  metabolism <- attr(swmpr_in, 'metabolism')
+  met_units <- attr(swmpr_in, 'met_units')
+
+  if(is.null(metabolism)) 
+    stop('No metabolism data, use the ecometab function')
+
+  # format for plotting
+  to_plo <- metabolism[, c('date', 'Pg', 'Rt', 'NEM')]
+  to_plo <- tidyr::gather(to_plo, Estimate, Value, -date)
+
+  # plot
+  p <- ggplot(to_plo, aes(x = date, y = Value, group = Estimate)) +
+    geom_line()
+
+  if(!pretty) 
+    return(p)
+
+  # ylabs
+  ylab <- expression(paste('mmol ', O [2], ' ', m^-2, d^-1))
+  if(met_units == 'grams')
+    ylab <- expression(paste('g ', O [2], ' ', m^-2, d^-1))
+
+  p <- p + 
+    geom_line(aes(colour = Estimate)) +
+    geom_point(aes(colour = Estimate)) +
+    theme_bw() +
+    theme(axis.title.x = element_blank())
+    scale_y_continuous(ylab)
+
+  return(p)
+
+
 }

R/swmpr_misc.R

@@ -11,7 +11,7 @@
 #' @return Returns a swmpr object to be used with S3 methods
 #' 
 #' @details 
-#' This function is a simple wrapper to \code{\link[base]{structure}} that is used internally within other functions to create a swmpr object.  The function does not have to be used explicitly.  Attributes of a swmpr object include \code{names}, \code{row.names}, \code{class}, \code{station}, \code{parameters}, \code{qaqc_cols}, \code{date_rng}, \code{timezone}, \code{stamp_class}, and \code{metabolism} (if present). 
+#' This function is a simple wrapper to \code{\link[base]{structure}} that is used internally within other functions to create a swmpr object.  The function does not have to be used explicitly.  Attributes of a swmpr object include \code{names}, \code{row.names}, \code{class}, \code{station}, \code{parameters}, \code{qaqc_cols}, \code{date_rng}, \code{timezone}, \code{stamp_class}, \code{metabolism} (if present), and \code{met_units} (if present). 
 #' 
 swmpr <- function(stat_in, meta_in){
 
@@ -44,7 +44,8 @@ swmpr <- function(stat_in, meta_in){
     date_rng = range(stat_in$datetimestamp),
     timezone = timezone, 
     stamp_class = class(stat_in$datetimestamp),
-    metabolism = NULL
+    metabolism = NULL, 
+    met_units = NULL
     )
 
 }

man/decomp_cj.Rd

@@ -19,8 +19,7 @@ decomp_cj(swmpr_in, ...)
 \item{vals_out}{logical indicating of numeric output is returned, default is \code{FALSE} to return a plot.}
 }
 \value{
-A \code{\link[ggplot2]{ggplot}} object if \code{vals_out = TRUE} (default), otherwise
-a monthly time series matrix of class \code{\link[stats]{ts}}.
+A \code{\link[ggplot2]{ggplot}} object if \code{vals_out = TRUE} (default), otherwise a monthly time series matrix of class \code{\link[stats]{ts}}.
 }
 \description{
 Decompose monthly SWMP time series into grandmean, annual, seasonal, and event series using \code{\link[wq]{decompTs}}, as described in Cloern and Jassby 2010.

man/ecometab.Rd

@@ -7,7 +7,7 @@
 \usage{
 ecometab(swmpr_in, ...)
 
-\method{ecometab}{swmpr}(swmpr_in, depth_val = NULL, units = "mmol",
+\method{ecometab}{swmpr}(swmpr_in, depth_val = NULL, met_units = "mmol",
   trace = FALSE, ...)
 }
 \arguments{
@@ -17,12 +17,12 @@ ecometab(swmpr_in, ...)
 
 \item{depth_val}{numeric value for station depth if time series is not available}
 
-\item{units}{chr indicating units of output for oxygen, either as mmol or grams}
+\item{met_units}{chr indicating units of output for oxygen, either as mmol or grams}
 
 \item{trace}{logical indicating if progress is shown in the console}
 }
 \value{
-A \code{\link[base]{data.frame}} of daily integrated metabolism estimates.
+The original \code{\link{swmpr}} object is returned that includes a metabolism attribute as a \code{\link[base]{data.frame}} of daily integrated metabolism estimates.  See the examples for retrieval.
 \describe{
  \item{\code{date}}{The metabolic day for each estimate, defined as the approximate 24 hour period between sunsets on adjacent calendar days.}
  \item{\code{ddo}}{Mean hourly rate of change for DO, mmol/m3/hr}
@@ -44,7 +44,7 @@ Estimate ecosystem metabolism using the Odum open-water method.  Estimates of da
 \details{
 Input data must be combined water quality and weather datasets.  This is typically done using the \code{\link{comb}} function after creating separate swmpr objects.
 
-The open-water method is used to infer net ecosystem metabolism using a mass balance equation that describes the change in dissolved oxygen (DO) over time as a function of photosynthetic rate, minus respiration rate, corrected for air-sea gas exchange. The air-sea gas exchange is estimated as the difference between the DO saturation concentration and observed DO concentration, multiplied by a volumetric reaeration coefficient (see the appendix in Thebault et al. 2008).  The diffusion-corrected DO flux estimates are averaged during day and night for each 24 hour period in the time series, where flux is an hourly rate of DO change. DO flux is averaged during night hours for respiration and averaged during day hours for net production. Respiration rates are assumed constant during day and night such that total daily rates are calculated as hourly respiration multiplied by 24. The metabolic day is considered the approximate 24 hour period between sunsets on two adjacent calendar days.  Respiration is subtracted from daily net production estimates to yield gross production.
+The open-water method is a common approach to quantify net ecosystem metabolism using a mass balance equation that describes the change in dissolved oxygen over time from the balance between photosynthetic and respiration rates, corrected for air-sea gas diffusion at the surface. The air-sea gas exchange is estimated as the difference between the DO saturation concentration and observed DO concentration, multiplied by a volumetric reaeration coefficient (see the appendix in Thebault et al. 2008).  The diffusion-corrected DO flux estimates are averaged during day and night for each 24 hour period in the time series, where flux is an hourly rate of DO change. DO flux is averaged during night hours for respiration and averaged during day hours for net production. Respiration rates are assumed constant during day and night such that total daily rates are calculated as hourly respiration multiplied by 24. The metabolic day is considered the approximate 24 hour period between sunsets on two adjacent calendar days.  Respiration is subtracted from daily net production estimates to yield gross production.
 
 Volumetric rates for gross production and total respiration are based on total depth of the water column, which is assumed to be mixed.  Water column depth is based on mean value for the depth variable across the time series in the \code{\link{swmpr}} object and is floored at 1 meter for very shallow stations.  Additionally, the volumetric reaeration coefficient requires an estimate of the anemometer height of the weather station, which is set as 10 meters by default.  The metadata should be consulted for exact height.  The value can be changed manually using a \code{height} argument, which is passed to \code{\link{calcKL}}.
 
@@ -69,9 +69,11 @@ res <- ecometab(dat, trace = TRUE)
 
 ## change height (m) of weather station anemometer
 res <- ecometab(dat, trace = TRUE, height = 5)
+res <- attr(res, 'metabolism')
 
 ## output units in grams of oxygen
-res <- ecometab(dat, trace = TRUE, units = 'grams')
+res <- ecometab(dat, trace = TRUE, met_units = 'grams')
+res <- attr(res, 'metabolism')
 }
 \references{
 Caffrey JM, Murrell MC, Amacker KS, Harper J, Phipps S, Woodrey M. 2013. Seasonal and inter-annual patterns in primary production, respiration and net ecosystem metabolism in 3 estuaries in the northeast Gulf of Mexico. Estuaries and Coasts. 37(1):222-241.
@@ -81,6 +83,6 @@ Odum HT. 1956. Primary production in flowing waters. Limnology and Oceanography.
 Thebault J, Schraga TS, Cloern JE, Dunlavey EG. 2008. Primary production and carrying capacity of former salt ponds after reconnection to San Francisco Bay. Wetlands. 28(3):841-851.
 }
 \seealso{
-\code{\link{comb}}
+\code{\link{comb}}, \code{\link{plot_met}}
 }
 

man/plot_met.Rd

@@ -0,0 +1,49 @@
+% Generated by roxygen2 (4.1.0): do not edit by hand
+% Please edit documentation in R/swmpr_analyze.R
+\name{plot_met}
+\alias{plot_met}
+\alias{plot_met.swmpr}
+\title{Plot ecosystem metabolism for a swmpr object}
+\usage{
+plot_met(swmpr_in, ...)
+
+\method{plot_met}{swmpr}(swmpr_in, pretty = TRUE, ...)
+}
+\arguments{
+\item{swmpr_in}{input swmpr object}
+
+\item{...}{}
+
+\item{pretty}{logical indicating use of predefined plot aesthetics}
+}
+\value{
+A \code{\link[ggplot2]{ggplot}} object which can be further modified.
+}
+\description{
+Plot gross production, total respiration, and net ecosystem metabolism for a swmpr object.
+}
+\details{
+A plot will only be returned if the \code{metabolism} attribute for the \code{\link{swmpr}} object is not \code{NULL}.  The \code{\link{ecometab}} function is used to create metabolism estimates.  See the examples.
+
+By default, \code{pretty = TRUE} will return a \code{\link[ggplot2]{ggplot}} object with predefined aesthetics.  Setting \code{pretty = FALSE} will return the plot with minimal modifications to the \code{\link[ggplot2]{object}}.  Use the latter approach for easier customization of the plot.
+}
+\examples{
+## import water quality and weather data
+data(apadbwq)
+data(apaebmet)
+
+## qaqc, combine
+wq <- qaqc(apadbwq)
+met <- qaqc(apaebmet)
+dat <- comb(wq, met)
+
+## estimate metabolism
+res <- ecometab(dat, trace = TRUE)
+
+## plot
+plot_met(res)
+}
+\seealso{
+\code\link{ecometab}
+}
+

man/swmpr.Rd

@@ -18,6 +18,6 @@ Returns a swmpr object to be used with S3 methods
 Wrapper for creating a swmpr object
 }
 \details{
-This function is a simple wrapper to \code{\link[base]{structure}} that is used internally within other functions to create a swmpr object.  The function does not have to be used explicitly.  Attributes of a swmpr object include \code{names}, \code{row.names}, \code{class}, \code{station}, \code{parameters}, \code{qaqc_cols}, \code{date_rng}, \code{timezone}, and \code{stamp_class}.
+This function is a simple wrapper to \code{\link[base]{structure}} that is used internally within other functions to create a swmpr object.  The function does not have to be used explicitly.  Attributes of a swmpr object include \code{names}, \code{row.names}, \code{class}, \code{station}, \code{parameters}, \code{qaqc_cols}, \code{date_rng}, \code{timezone}, \code{stamp_class}, \code{metabolism} (if present), and \code{met_units} (if present).
 }