6

Currently I am working on a river discharge data analysis. I have the daily discharge record from 1935 to now. I want to extract the annual maximum discharge for each hydrolocial year (start from 01/11 to next year 31/10). However, I found that the hydroTSM package can only deal with the natural year. I tried to use the "zoo" package, but I found it's difficult to compute, as each year have different days. Does anyone have some idea? Thanks.

the data looks like:

01-11-1935 663
02-11-1935 596
03-11-1935 450
04-11-1935 381
05-11-1935 354
06-11-1935 312

my code:

mydata<-read.table("discharge")
colnames(mydata) <- c("date","discharge")

library(zoo)
z<-zooreg(mydata[,2],start=as.Date("1935-11-1"))

mydta$date <- as.POSIXct(dat$date)

q.month<-daily2monthly(z,FUN=max,na.rm = TRUE,date.fmt = "%Y-%m-%d",out.fmt="numeric")
q.month.plain=coredata(q.month)

z.month<-zooreg(q.month.plain,start=1,frequency=12)

4 Answers 4

8

With dates stored in a vector of class Date, you can just use cut() and tapply(), like this:

## Example data
df <- data.frame(date = seq(as.Date("1935-01-01"), length = 100, by = "week"),
                 flow = (runif(n = 100, min = 0, max = 1000)))

## Use vector of November 1st dates to cut data into hydro-years
breaks <- seq(as.Date("1934-11-01"), length=4, by="year")
df$hydroYear <- cut(df$date, breaks, labels=1935:1937)

## Find the maximum flow in each hydro-year
with(df, tapply(flow, hydroYear, max))
#     1935     1936     1937 
# 984.7327 951.0440 727.4210 


## Note: whenever using `cut()`, I take care to double-check that 
## I've got the cuts exactly right
cut(as.Date(c("1935-10-31", "1935-11-01")), breaks, labels=1935:1937)
# [1] 1935 1936
# Levels: 1935 1936 1937
3

Here is a one-liner to do that.

First convert the dates to "yearmon" class. This class represents a year month as the sum of a year as the integer part and a month as the fractional part (Jan = 0, Feb = 1/12, etc.). Add 2/12 to shift November to January and then truncate to give just the years. Aggregate over those. Although the test data we used starts at the beginning of the hydro year this solution works even if the data does not start on the beginning of the hydro year.

# test data
library(zoo)
z <- zooreg(1:1000, as.Date("2000-11-01")) # test input

aggregate(z, as.integer(as.yearmon(time(z)) + 2/12), max)

This gives:

2001 2002 2003 
 365  730 1000 
2

Try the xts package, which works together with zoo:

require(zoo)    
require(xts)

dates = seq(Sys.Date(), by = 'day', length = 365 * 3)
y = cumsum(rnorm(365 * 3))    
serie = zoo(y, dates)

# if you need to specify `start` and `end`
# serie = window(serie, start = "2015-06-01")

# xts function
apply.yearly(serie, FUN = max)
6
  • many thanks!! Each hydrolocial year (start from 01/11 to next year 31/10)? how to choose the start of the computing year?
    – uared1776
    Commented Feb 27, 2014 at 16:31
  • You can use the window function to specify start and end, see edit.
    – Fernando
    Commented Feb 27, 2014 at 16:34
  • thank you very much for your swift reply. I tried, however it still calculate the nature year (using the "window" function makes the calculation start from the given date and end at 31 December of that year, the following year still calculate from 1 Jan to 31 December). My purpose is to calculate start from 1 Nove to 31 Oct second year, and coutinuously... e.g. the maxiums: 1/11/2002 to 31/10/2003, 1/11/2003 to 31/10/2004, 1/11/2004 to 31/10/2005.....
    – uared1776
    Commented Feb 27, 2014 at 16:59
  • @user3361298 -- any reason, then, that you accepted this answer over mine, which has the virtue of doing what you actually asked for? If you're wanting only solutions that use the xts package, you should add that to your question... Commented Feb 27, 2014 at 17:38
  • 1
    @Josh I know you didn't mean to offend, english is not my first language, so sorry if i appear rough too. Maybe the OP found it simpler to use xts...but i'm happy if he chooses your solution over mine.
    – Fernando
    Commented Feb 27, 2014 at 18:35
0

You can use the apply.seasonal function from lfstat package which operates over xts objects.

To solve your case in one line:

apply.seasonal(mydata, varying = "yearly", fun = max, origin = 11)

origin = 11 means that the hydrological year will start in november.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.