dvpart: Perform Streamflow Partitioning after Rutledge (1998)...

In wasquith-usgs/akqdecay: Many Distributional Statistics of USGS Daily Streamflow Hydrographs and Powerful Utilities

Description

Perform streamflow partitioning to base-flow using port of algorithms of the DVstats package and the function DVstats::part() of daily mean-streamflow values based on the daily value structure of dvget. The function fill_dvpartenv is a wrapper that can be used to fill an R environment with the output of the dvpart function. (Users interested in monthly streamflow partitioning should consult dvpart2mnpart.)

The algorithms implement the procedures of Rutledge (1998). For the akqdecay package, some new development was made on how to handle missing flows to help in data mining studies. First, instead of issuing stop() as done in DVstats::part(), the dvpart implementation returns NULL. Second, a more sophisticated method is optionally available (see Details).

Usage

dvpart(akdvtable, sdate="", edate="", cda=NA, site_no=NA,
                  fillgaps=FALSE, na.negflow=FALSE, ...)

Arguments

akdvtable	A USGS daily-value of streamflow table with some extensions unique to this package as returned by dvget;
sdate	Optional start date (default is earliest) and nomenclature matches that of the dataRetrieval package with string format of “YYYY-MM-DD” for year (YYYY), month (MM), and day (DD) respectively padded by zeros as needed. The akdvtable is consulted for the first day of record for the default is sdate = "";
edate	Optional ending date (default is earliest) and nomenclature matches that of the dataRetrieval package with string format of “YYYY-MM-DD” for year (YYYY), month (MM), and day (DD) respectively padded by zeros as needed. The akdvtable is consulted for the last day of record for the default is edate = "";
cda	Contributing drainage area of the streamgage in square miles. Only the first value is used without warning if a vector is given. The area can readily be determined from sites_to_SpatialPointsDataFrame, though a missing drainage area will default to unity inside dvpart;
site_no	Optional site number that will be taken if NA from the akdvtable as needed;
fillgaps	A logical to trigger the record gap infiller (see Details). This trigger when set to TRUE is an extension of the Rutledge (1998) algorithm as represented by DVstats::part();
na.negflow	A logical to trigger the conversion of FlowBase, FlowPart1, FlowPart2, and FlowPart3 to NA is the Flow (actual daily mean streamflow) is a negative value. The original algorithms truncate negative flows a “small” near $0^+ value. This trigger when set to TRUE is an extension of the Rutledge (1998) algorithm as represented by DVstats::part(); and
...	Other arguments to pass.

Details

Missing data and gaps in the record can actually be accommodated by a simple approach triggered by the fillgaps argument. It does require mentioning that there are two types of missing information. Real NAs could be delivered on the daily value retrieval, but also jumps or gaps in time too. The dvpart checks for both conditions and returns NULL if either is encountered—Except if the fillgaps is set.

The strategy has the following steps. First, construct an absolutely daily-continuous sequence of dates from sdate to edate. Second, flows that are less than or equal to zero are set to small number (see pmax use in the sources). Third, log-linearly interpolate for the daily-continuous sequence that streamflow across record gaps using the non-NA flows and respective dates in akdvtable. Fourth, streamflow separation is made on the daily-continuous streamflows instead of the non-NA original non-daily continuous data.

The fifth and last step is crucial. For only the dates in the akdvtable, extract the streamflow separation estimates and add them to the table. This means that the output table from this function still retains gaps in record.

Value

An R data.frame is returned with these expected columns.

agency_cd	The agency code for the data;
site_no	The streamgage identification number;
Date	The date, note that the capitalization is from the dataRetrieval package, elsewhere in akqdecay this will become lower case;
Flow	The streamflow in cubic feet per second, note that the capitalization is from the dataRetrieval package;
Flow_cd	A coding system for the streamflow (A, approved record; P, provisional, and W, working record);
site	A character representation of the site_no, which likely is already a character. The reasoning for having another site column is in case a user need some type of flexibility in later processing. The user could freely replace contents of this column;
year	The calendar year of the date;
decade	The decade of the daily value. The decade is assign by taking the year and the trailing digit has been stripped and replaced with zero. This is not a technique in which a “decade” is centered on an even step of 10—meaning, say that 1996–2005 is not the “2000 decade” but simply 01/01/2000–12/31/2009 is the “2000 decade;”
wyear	The water year;
month	The month;
FlowBase	The baseflow of the Flow (the column BaseQ from DVstats::part());
FlowPart1	The partition 1 baseflow of the Flow (the column Part1 from DVstats::part());
FlowPart2	The partition 2 baseflow of the Flow (the column Part2 from DVstats::part()); and
FlowPart3	The partition 3 baseflow of the Flow (the column Part3 from DVstats::part()).

Note

For the greater purposes of the akqdecay package, it was desired to have an independent subsystem for streamflow partitioning from external packages. The DVstats::part() as stated above is the basis with some code cleaning to the style of Asquith and Knight that includes the critical adjustment to processing on daily values in the canonical form of dvget. The part() implementation requires a short suite of utilities: smwrBase::eventNum(), smwrBase::na2miss(), and smwrBase::shiftData(). These are embedded within dvpart so that akqdecay is fullly independent of DVstats and smwrBase. Both na2miss() and shiftData() in original form have special handling for factor variables. This logic is not needed and has been stripped within the dvpart sources.

Rutledge (1998) is the basis for the algorithm. One oft question is why a drainage area is needed for the streamflow partitioning. Rutledge (1998, p. 10) states “The duration of time after a peak in streamflow during which the components of streamflow due to surface runoff and interflow are significant can be estimated [by empirical formula].” The formula is N = A^{0.2}, where N is the number of days after the peak and A is the drainage area in square miles. Algorithmically, the N is then cast as the integer \lceil N \rceil using the ceiling() function in R. Rutledge (1998, p. 3) also states that the lower limit is 1 square mile and suggests 500 square miles as an upper limit (2,000 square miles is also mentioned as an upper limit). No check against these limits is provided in dvpart and evidently not in DVstats::part() for that matter.

Author(s)

W.H. Asquith

References

De Cicco, L.A., Lorenz, D.L., 2017, smwrBase—Functions to import and manipulate hydrologic data: R package version 1.1.5, September 18, 2017, accessed April 28, 2019 at https://github.com/USGS-R/smwrBase.

Lorenz, D.L., De Cicco, L.A., 2017, DVstats—Functions to manipulate daily-values data: R package version 0.3.4, July 26, 2017, accessed April 28, 2019 at https://github.com/USGS-R/DVstats.

Rutledge, A.T., 1998, Computer programs for describing the recession of ground-water discharge and for estimating mean ground-water recharge and discharge from streamflow data—Update: U.S. Geological Survey Water-Resources Investigations Report 98–4148, 43 p., https://doi.org/10.3133/wri984148.

See Also

dvget, fill_dvpartenv, dvpart2mnpart

Examples

# USGS 14362000 Applegate River near Copper, Oregon
dv <- dvget("14362000", sdate="1989-10-01", edate="1999-09-30")
pdv <- dvpart(dv, cda=225) # The drainage area is 225 square miles.
plot(pdv$Date,  pdv$Flow, log="y", type="l", col=8)
lines(pdv$Date, pdv$FlowBase,  col=1)
lines(pdv$Date, pdv$FlowPart1, col=2)
lines(pdv$Date, pdv$FlowPart2, col=3)
lines(pdv$Date, pdv$FlowPart3, col=4) #


DV  <- dvget("07031740", edate="2017-12-31"); cda <- 788 # square miles
pDV <- dvpart(DV, cda=cda) # NULL will be return and here is the message
#         for 07031740 noncontinuous data between 1995-02-01 to 2017-11-28
#         gaps about: 1995-03-23[12days], 1995-06-25[53days], 1995-10-22[22days],
# 1996-02-12[18days], 1996-03-19[11days], 1996-04-20[26days], 1996-06-29[65days],
# 2017-01-21[97days], 2017-03-07[28days], 2017-04-05[17days], 2017-04-27[6days],
# 2017-06-16[3days] : total=358

pDV <- dvpart(DV, cda=cda, fillgaps=TRUE) # see Details
length(pDV$Date)  # 7,991. This means this many discharges NA or numeric
length(seq(pDV$Date[1],
       as.Date("2017-12-31"), by=1)) # 8370. This mean days desired/expected.
# Note 8370 - 7991 != 358 above because gaps at beginning or ending not counted.

DV  <- dvget("07031740", edate="2017-12-31", ignore.provisional=FALSE)
pDV <- dvpart(DV, cda=cda, fillgaps=FALSE)
#          for 07031740 noncontinuous data between 1995-02-01 to 2017-12-31
#          gaps about: 1995-03-23[12days], 1995-06-25[53days], 1995-10-22[22days],
# 1996-02-12[18days], 1996-03-19[11days], 1996-04-20[26days],
# 1996-06-29[65days] : total=207
# ** So we see improvement! Gaps only at beginning of gage. ** #

## Not run: 
# Let us now work on practical example with further implications.
# We want monthly baseflows for this gage that had the record gaps.
cda <- sites_to_SpatialPointsDataFrame("07031740")$CDA # 788 square miles
DV  <- dvget("07031740", edate="2017-12-31", ignore.provisional=FALSE)
pDV <- dvpart(DV, cda=cda, fillgaps=TRUE)
length(DV$Date)  # 8,170   See how no new days are created!
length(pDV$Date) # 8,170   This is a big deal (see Details).

MN <- dvpart2mnpart(pDV)
tail(MN); # shows July 2017 Flow as 351.3871, manually NWISWeb shows 351.4
          # shows July 2017 FlowBase as 275.7239
plot(MN$DateMean, MN$Flow, type="l", log="y",
                  xlab="Date", ylab="Streamflow, cfs")
lines(MN$DateMean, MN$FlowBase, col=2) #
## End(Not run)

## Not run: 
# This example demonstrates setting the baseflows to NA when flows are NA
cda <- sites_to_SpatialPointsDataFrame("07382500")$CDA # 715 square miles
DV  <- dvget("07382500", sdate="2000-01-01", edate="2012-12-31",
                         ignore.provisional=FALSE)
pDV <- dvpart(DV, cda=cda, fillgaps=TRUE, na.negflow=TRUE)
head(pDV[pDV$Flow < 0,]) # see how the baseflows are NAs
mDV <- dvpart2mnpart(pDV)
mDV[      mDV$Flow  < 0,] # October 2010 with baseflows as NAs to view all
mDV[is.na(mDV$FlowBase),] # months in which at least one daily streamflow
                          # was less than zero 
## End(Not run)

wasquith-usgs/akqdecay documentation built on Nov. 9, 2020, 1:13 p.m.