In CentreForHydrology/MESHr: pre- and post processing for MESH

knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)
library(MESHr)
library(ggplot2)
MESH_streamflows <- MESH_streamflows

simpleHydrograph

The simpleHydrograph function is intended to plot hydrographs from MESH output alone, i.e. without additional data being supplied from another source, such as WSC data. The function produces a ggplot2 object, which can then be manipulated to customise the plot.

The function requires a data frame of MESH flows. The first column must be an R date (named DATE) or a POSIXct datetime (named DATETIME). The easiest way to produce this data frame is by reading in the data using the function read_MESH_OutputTimeseries_csv. Functions will be added to read the output directly from the netCDF output.

Note that you can subset the data to restrict the range plotted.

The simplest plot shows all of the data in a single un-facetted graph. The stations are denoted by colour, and the type of flow (i.e. measured or simulated), by line type. By default, the stations are labelled by their MESH numbers.

simpleHydrograph(MESH_streamflows)

You can change the station names by specifying them in a vector.

station_names <- c("Station 1", "Station 2")
simpleHydrograph(MESH_streamflows, stationNames = station_names)

The function encourages the use of facetting to produce more readable plots. When facetting by station, it's a good idea to turn off the coloring of plots by station, by setting the option byStation = FALSE. You can change the x-axis scale and/or labels. Note that this example forces the facets to be on separate rows, so that the x-axes are as long as possible. Because this requires other functions in ggplot2, you should load it before proceeding.

Because we will be modifying the plot, it needs to be saved to a variable, in this case p. The value of the variable is updated with any number of ggplot2 commands. To plot the graph at any time, just type the name of the variable to show its value.

p <- simpleHydrograph(MESH_streamflows, stationNames = station_names, byStation = FALSE)
library(ggplot2)
p <- p + facet_wrap(~station, nrow = 2) + scale_x_date(date_labels = "%Y")
p

Because the magnitudes of the flows are so different for the stations, you can use separate y scales:

p <- p + facet_wrap(~station, nrow = 2, scales = "free_y") 
p

Alternatively, you can facet the plots by the year. Note that this option adds a variable to the plot called YEAR. The uppercase variable name is in keeping with the other MESH variables.

p <- simpleHydrograph(MESH_streamflows, stationNames = station_names, byYear = TRUE)
p <- p + facet_wrap(~YEAR, scales = "free_y")
p

You can also create gridded facets by year and station.

p <- simpleHydrograph(MESH_streamflows, stationNames = station_names, byStation = FALSE, byYear = TRUE)
p <- p + facet_grid(YEAR~station, scales = "free_y")
p

The default colors are a bit bland - you can easily change them. In this example, the y-axis labels are rotated, so they show up better.

plotcols <- c("red", "blue")
p <- p + scale_colour_manual(values = plotcols) +
  theme(strip.text.y = element_text(angle = 0))
p

You can plot only the measured or the simulated values. This this case, it's a good idea to suppress the flow type legend.

p <- simpleHydrograph(MESH_streamflows, stationNames = station_names, 
                      byStation = FALSE, byYear = TRUE, meas = FALSE) +
  facet_grid(YEAR~station, scales = "free_y") +
  scale_colour_manual(values = plotcols) +
  theme(strip.text.y = element_text(angle = 0)) +
  theme(legend.position = "none")

p

You can specify that the modelled data consist of both caibration and verification data, by specifying the beginning and/or the end dates of the calibration series. This adds the variable period to the plot, which can then be used for facetting. The value of the variable will be "Calibration" during the calibration period and "Validation" during the rest.

calEnd <- "2010-01-01"
p <- simpleHydrograph(MESH_streamflows, stationNames = station_names, 
                       byStation = FALSE, calEnd = calEnd)
p <- p + facet_grid(station~period) + scale_x_date(date_labels = "%Y")
p

You can easily change the scales (x, y, and colour) to make them more attractive. Note that this example modifies the plot just created, and stores it in a new varaiable.

plotcols <- c("red", "blue")
p1 <- p + facet_grid(station~period, scales = "free") + 
  scale_x_date(date_labels = "%Y") +
  scale_colour_manual(values = plotcols)
p1

hydroStats

You can use the function hydroStats to calculate the goodness of fit of your modelled hydrograph. This function is a wrapper for the function gof in the package hydroGOF, making it easier to use with MESH data.

stats <- hydroStats(MESH_streamflows, stationNames = station_names)
stats

Once the stats have been calculated, they can be added to a hydrograph. You can add the station names using the ggplot2 function annotate. To select the stat you want, specify its variable name. You can combine values together, or add station names.

NSE <- paste(stats$station, "NSE =", stats$NSE, sep = " ")
NSE

To add both values as a single string (which is required to annotate a single graph), you can paste them together

all_NSE <- paste(NSE, collapse = ", ")
all_NSE

Then use annotate to add the string to the graph. Note that you have to specify the location of the string, relative to the axes. Because the x-axis uses dates, you have to specify the location as a date.

p2 <- simpleHydrograph(MESH_streamflows, stationNames = station_names, byStation = TRUE) +
  scale_x_date(date_labels = "%Y")

dateLoc <- as.Date("2005-01-01", format = "%Y-%m-%d")

p3 <- p2 + annotate("text", x = dateLoc, y = 1500, label = all_NSE, hjust = "left")
p3

Since the text doesn't fit very well, you can wrap it before adding it to the graph.

wrapped_NSE <- paste(NSE, collapse = "\n")
p4 <- p2 + annotate("text", x = dateLoc, y = 1500, label = wrapped_NSE, hjust = "left")
p4

The hydroStats function also works with separate calibration and validation periods.

stats2 <- hydroStats(MESH_streamflows, stationNames = station_names, calEnd = calEnd)
stats2

CentreForHydrology/MESHr documentation built on Jan. 11, 2021, 8:34 p.m.