[R] error in seq.POSIXt?

[Don MacQueen]

Brandon,
Here is the explanation.

Think in terms of two different definitions of "day":
   (1) an interval of 24 hours, and
   (2) the interval from one day to the same time the next day.
These are identical, except for two days each year in time zones 
where an adjustment is made for the summer months. Then, (2) is 23 
hours when the interval includes the transition to summer time, and 
25 when it includes the transition back.

by='day' in seq.POSIXt() uses definition (1)
by='DSTday' uses definition (2)  (only available in R >= version 1.5.0).

temp.begin is in standard time. In 1974 and 1975, temp.end is in 
summer time, and so is one hour short of the 90th day, causing 
seq.POSIXt to return only 89 days. The other years, temp.end is also 
in standard time, so temp.end and temp.begin are an integer multiple 
of 24 hours apart.

I suggest using
   seq(temp.begin, length.out=90, by="DSTday")
instead of
   seq(temp.begin, length.out=90, by="day")
study the output of the code below to see why.


for(year in 1973:1976) {
   temp.begin <- ISOdatetime(year-1,12,1,0,0,0)
   temp.end <- ISOdatetime(year,2,28,0,0,0)
   ts1 <- seq(temp.begin, temp.end, by="day")
   ts2 <- seq(temp.begin, length.out=90, by="day")
   ts3 <- seq(temp.begin, length.out=90, by="DSTday")
   cat(
       format(temp.end,'%Y-%m-%d %H:%M %Z'),
       format(max(ts1),'%Y-%m-%d %H:%M %Z'),
       format(max(ts2),'%Y-%m-%d %H:%M %Z'),
       format(max(ts3),'%Y-%m-%d %H:%M %Z'),
       length(ts1),
       length(ts2),
       length(ts3),
       difftime(max(ts1),ts1[1],units='h'),
       difftime(max(ts2),ts2[1],units='h'),
       difftime(max(ts3),ts3[1],units='h'),
       (as.numeric(temp.end)-as.numeric(temp.begin))/(24*60*60),
       "\n")
}


In case anyone reading Brandon's first message thought he had found a 
bug, I have to say, no this is not a bug. The transitions to and from 
daylight savings time (to use the American misnomer) are like a 
minefield waiting to trip the unwary--and I have been caught by them 
many times.

-Don

At 3:27 PM -0600 5/31/02, Brandon Whitcher wrote:
>I am trying to extract only the winters (defined to be 01-Dec through
>28-Feb) of daily data from 1948-2002.  There are 90 days in each winter
>season.  I wrote the following code to gather the winter dates into a
>single vector:
>
>DJF <- NULL
>for(year in 1949:1999) {
>   temp.begin <- strptime(paste("01/12", year-1, sep="/"), "%d/%m/%Y")
>   temp.end <- strptime(paste("28/02", year, sep="/"), "%d/%m/%Y")
>   ## Gives wrong answer:
>   temp.seq <- seq(temp.begin, temp.end, by="day")
>   if(is.null(DJF))
>     DJF <- temp.seq
>   else
>     DJF <- c(DJF, temp.seq)
>   cat("Winter:", year, "\t", length(temp.seq), "\n")
>}
>
>with the following output:
>
>Winter: 1949	 90
>Winter: 1950	 90
>Winter: 1951	 90
>Winter: 1952	 90
>Winter: 1953	 90
>Winter: 1954	 90
>Winter: 1955	 90
>Winter: 1956	 90
>Winter: 1957	 90
>Winter: 1958	 90
>Winter: 1959	 90
>Winter: 1960	 90
>Winter: 1961	 90
>Winter: 1962	 90
>Winter: 1963	 90
>Winter: 1964	 90
>Winter: 1965	 90
>Winter: 1966	 90
>Winter: 1967	 90
>Winter: 1968	 90
>Winter: 1969	 90
>Winter: 1970	 90
>Winter: 1971	 90
>Winter: 1972	 90
>Winter: 1973	 90
>Winter: 1974	 89
>Winter: 1975	 89
>Winter: 1976	 90
>Winter: 1977	 90
>Winter: 1978	 90
>Winter: 1979	 90
>Winter: 1980	 90
>Winter: 1981	 90
>Winter: 1982	 90
>Winter: 1983	 90
>Winter: 1984	 90
>Winter: 1985	 90
>Winter: 1986	 90
>Winter: 1987	 90
>Winter: 1988	 90
>Winter: 1989	 90
>Winter: 1990	 90
>Winter: 1991	 90
>Winter: 1992	 90
>Winter: 1993	 90
>Winter: 1994	 90
>Winter: 1995	 90
>Winter: 1996	 90
>Winter: 1997	 90
>Winter: 1998	 90
>Winter: 1999	 90
>Winter: 2000	 90
>Winter: 2001	 90
>Winter: 2002	 90
>
>There are two years, 1974 and 1975, that produce only 89 days per winter.
>I could not figure out what the cause of this was, but a workaround is two
>rewrite the call to seq.POSIXt as:
>
>temp.seq <- seq(temp.begin, length.out=90, by="day")
>
>I am unsure why one works and one does not.  Any ideas?
>
>Brandon
>
>----------------------------------------------------------------------------
>   Geophysical Statistics Project
>   National Center for Atmospheric Research          +1 303 497 1709  voice
>   P.O. Box 3000, Boulder, CO  80307-3000            +1 303 497 1333  fax
>
>   whitcher at ucar.edu                            www.cgd.ucar.edu/~whitcher/
>----------------------------------------------------------------------------
>
>
>
>
>
>-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
>r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
>Send "info", "help", or "[un]subscribe"
>(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
>_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._


-- 
--------------------------------------
Don MacQueen
Environmental Protection Department
Lawrence Livermore National Laboratory
Livermore, CA, USA
--------------------------------------
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._