CityWaterBalance provides a reproducible workflow for studying urban water systems. Any system may be modeled with preassembled data, but data for US cities can be gathered via web services using this package and dependencies geoknife and dataRetrieval.

Install

The latest release of CityWaterBalance is available from CRAN. To install:

install.packages("CityWaterBalance")

Development version available on Github:

devtools::install_github("USEPA/CityWaterBalance", 
                         build_vignettes = TRUE)
library("CityWaterBalance")

Usage overview

CityWaterBalance is based on a model of the urban water system, shown in the diagram below. This diagram specifies the network of water flows along with a mathematical solution for the changes in water storages (i.e., inflows - outflows) within the system.



The function, CityWaterBalance(), evaluates the urban system model according to the numbered flows described in the table below.

Usage examples

There are two ways to apply the CityWaterBalance workflow.

Option 1: Input preassembled data

See CityWaterBalance() inputs for details on data and parameter inputs. Data must be in self-consistent units. Each row of data must represent the same time period.

The following example is based on the included dataset (cwb_data) on flows of water in the Chicago metropolitan area.

# Specify parameters
p <- list("interc" = 0, "et_mult" = 1, "flow_mult" = 1, "open_wat" = 0.02, "run_mult" = 1, "run_css" = 0.35, "bf_mult" = 1, "nonrev" = 0.08,"ind_evap" = 0.012, "wast_gen" = 0.85, "pot_atm" = 0.13, "npot_infilt" = 0.5, "slud_evap" = 0, "leak_css" = 0.05, "dgw" = 0.5,"dgw_rep" = 0.5)

# Run model
m <- CityWaterBalance(cwb_data,p, print=FALSE)

Output from the model includes a list of 5 zoo series with values at each timestep for: 1) global flows, 2) internal flows, 3) state variables (storages, producers and consumers), 4) internal and 5) global water balances.

# Visualize output
gf <- m$global_flows
plotWaterBalance(gf, yl = "Flux (mm/month)")

Option 2: Input data gathered from web services

CityWaterBalance has other functions that assemble data for the model. At this time, these functions access US-based web services.

Specify spatial and temporal boundaries

Define an area of interest (AOI) and upload that geometry to the USGS Geo Data Portal (GDP). The GDP will give the geometry a name, which may start with ‘upload:’. Here we use a geometry that is already available to the GDP in order to automate the example.

geometry <- 'sample:Counties'
attribute <- 'STATE'
value <- 'RI'
area <- 2707
start <- "2010-01-01"
end <- "2010-12-31"

Get atmospheric data

latitude <- 41.5801
atm <- getAtmoFlows(start, end, geometry, attribute, value, latitude)

Get streamflow data

Choose streamgages to evaluate total inflow and outflow for the AOI.
NWIS mapper may be useful here.

ingages <- c("01112500")
outgages <-c("01113895","01114000","01117000","01118500")

inflows <- getStreamflow(start,end,ingages)
outflows<- getStreamflow(start,end,outgages)

Look at the streamflow data with the function plotStreamflow(). Gaps can be filled using gapfillStreamflow(). When time series for all gages are complete, use combineStreamflow() to aggegate total in/outflow.

Get water use data

Identify the states and counties in your AOI. County is the finest spatial scale for USGS water use data served by NWIS.

states <- c("RI")
counties <- list(c("Providence","Kent","Bristol","Newport","Washington"))
wu_raw <- getWaterUse(states,counties,years="ALL")

The above retrieves USGS water use data for the specified counties for all available years and withdrawal categories. To reorganize the output for use by CityWaterBalance() use:

wu <- combineWaterUse(start,end,wu_raw)

Get other data

Not all input data is currently available via web services. These inputs (i.e., wastewater effluent, sewer overflows, runoff, baseflow and deep groundwater recharge) must be gathered from other sources and converted to xts objects with the same temporal range and resolution as the other flows. Units must correspond with those specified for inputs to mergeData().

Finally, merge the data into a single xts of fluxes (i.e., flow/AOI) for input to CityWaterBalance().

# combine streamflow among gauges
inflow <- combineStreamflow(inflows, c(1))
outflow <- combineStreamflow(outflows, c(1,1,1,1))
model_data <- mergeData(area,atm,inflow,outflow,wu)

Solve

CityWaterBalance() solves for the changes in storage of system components. Acceptable solutions can be determined using ancillary observations (e.g., groundwater levels). Owing to uncertainty in input data, parameters, and ancillary observations of system storages, acceptable solutions will vary within a solution space.

Solutions within this space may be discovered, and the sensitivity of modeled flow explored, with the function getSolutions(). getSolutions() tests for parameter sets that satisfy criteria for the changes in all system storages. It will continue testing additional parameter sets until the mean flow solution differs by less than a user-specified tolerance (tol). Although the CityWaterBalance() model runs quickly, running it many times to achieve a stable solution may take quite a while.

The following example shows how to search for model solutions that satisfy a criteria for storage in the combined sewer system again using cwb_data. In this dataset, the combined sewer overflow (cso) do not span the complete analysis period. To compensate, the cso flow is set to zero and the acceptable range for the css balance is set to the mean and median of the original data. Other criteria are left at default values. The tol is set high so the example will run in a reasonable amount of time.

The boxplot call illustrates one way that results may be visualized.

d <- cwb_data
d$cso <- 0

params <- list("interc" = 0, "et_mult" = 1, "flow_mult" = 1, "open_wat" = 0.02, "run_mult" = 1, "run_css" = 0.35, "bf_mult" = 1, "nonrev" = 0.08,"ind_evap" = 0.012, "wast_gen" = 0.85, "pot_atm" = 0.13, "npot_infilt" = 0.5, "slud_evap" = 0, "leak_css" = 0.05, "dgw" = 0.5,"dgw_rep" = 0.5)

low <- median(cwb_data$cso,na.rm=TRUE)*nrow(cwb_data)
high <- mean(cwb_data$cso,na.rm=TRUE)*nrow(cwb_data)
css_crit <- c(low,high)

out <- getSolutions(data = d, p = params, n = 100, tol = 100, css_bal = css_crit)
boxplot(out[,order(colMeans(out),decreasing=TRUE)]/sum(d$prcp),las = 2, 
ylab = "flow/precipitation", ylim = c(0,1),col = c("lightblue1"), 
border = c("royalblue3"))