[R] Turning points in a series

Thu Sep 17 12:35:51 CEST 2009

Hello,

I don't see what's wrong with turnpoints() from pastecs. It is easy to 
use, and provides additional information for each turnpoints, i.e., 
probability of occurrence against the null hypothesis that the series is 
purely random, and the number of bits of information associated with the 
points according to Kendall's information theory. See ?turnpoints.

Rewriting the function is a nice exercise, but a small explanation on 
how to use turnpoints() is much easier. So, nobody is able to tell that 
simply using:

 > library(pastecs)
 > turnpoints(dat$count)

does the job? Am I the only one interested by the extra information 
provided by turnpoints()?

Here is a more extensive example that shows also how to get date or 
counts associated to pits/peaks:

txt <- "
y  m   d  count
93 02 07 3974.6
93 02 08 3976.7
93 02 09 3955.2
93 02 10 3955.0
93 02 11 3971.8
93 02 12 3972.8
93 02 13 3961.0
93 02 14 3972.8
93 02 15 4008.0
93 02 16 4004.2
93 02 17 3981.2
93 02 18 3996.8
93 02 19 4028.2
93 02 20 4029.5
93 02 21 3953.4
93 02 22 3857.3
93 02 23 3848.3
93 02 24 3869.8
93 02 25 3898.1
93 02 26 3920.5
93 02 27 3936.7
93 02 28 3931.9
"
con <- textConnection(txt)
dat <- read.table(con, header = TRUE)
close(con)
dat$date <- as.Date(paste(dat$y, dat$m, dat$d), format = "%y %m %d")

library(pastecs)
tp <- turnpoints(dat$count)
tp
summary(tp)
# Indicate which turnpoints are significant (see ?turnpoints)
plot(tp, level = 0.05)
# Another plot
plot(dat$count, type = "l")
lines(tp)

# Get counts for all turnpoints
allcounts <- <- dat$count[extract(tp, no.tp = FALSE, peak = TRUE, pit = 
TRUE)]

# Get dates for all turnpoints
alldates <- dat$date[extract(tp, no.tp = FALSE, peak = TRUE, pit = TRUE)]
alldates
# Get dates for "informative" turnpoints (5%) only (see ?turnpoints)
alldates[tp$proba < 0.05]
# Get dates for peaks only
dat$date[extract(tp, no.tp = FALSE, peak = TRUE, pit = FALSE)]
# Etc...

Best,

Philippe Grosjean

..............................................<°}))><........
  ) ) ) ) )
( ( ( ( (    Prof. Philippe Grosjean
  ) ) ) ) )
( ( ( ( (    Numerical Ecology of Aquatic Systems
  ) ) ) ) )   Mons-Hainaut University, Belgium
( ( ( ( (
..............................................................

(Ted Harding) wrote:
> On 17-Sep-09 08:10:47, ogbos okike wrote:
>> Good morning once more. My problem of yesterday has been addressed.
>> Having learned a few tricks from that, I wish to ask another question
>> in connection with that. My data is a cosmic ray data consisting of
>> dates and counts.
>> When I plot a graph of counts versus dates, the resultant signal
>> shows a number of maximum and minimum points. These minimum points
>> (turning points) are of interest to me. Reading these dates and counts
>> off from the plot is difficult as I am dealing with a large data.
>> I have been looking at turnpoints function in pastecs library but have
>> not been able to figure out the appropriate commands that one can use
>> to find the minima/maxima (turning points) or pits/peaks in a series.
>> My data is of the form shown below where y stands for year, m month,
>> d day and finally count. Is there a way I could find these minima
>> together with the dates they occurred?
>> I would be indebted to those of you who will show me the way out of
>> these problem.
>> Thank you.
>> Best regards
>> Ogbos
>>
>>
>> y  m   d  count
>> 93 02 07 3974.6
>> 93 02 08 3976.7
>> 93 02 09 3955.2
>> 93 02 10 3955.0
>> 93 02 11 3971.8
>> 93 02 12 3972.8
>> 93 02 13 3961.0
>> 93 02 14 3972.8
>> 93 02 15 4008.0
>> 93 02 16 4004.2
>> 93 02 17 3981.2
>> 93 02 18 3996.8
>> 93 02 19 4028.2
>> 93 02 20 4029.5
>> 93 02 21 3953.4
>> 93 02 22 3857.3
>> 93 02 23 3848.3
>> 93 02 24 3869.8
>> 93 02 25 3898.1
>> 93 02 26 3920.5
>> 93 02 27 3936.7
>> 93 02 28 3931.9
> 
> The following simple function TP() (for "Turning Point") locates
> the positions i where x[i] is greater than both of its immediate
> neighbours (local maximum) or less than both of its neighbours
> (local minimum).
> 
>   TP <- function(x){
>     L <- length(x)
>     which( ((x[1:(L-2)]<x[2:(N-1)])&(x[2:(L-1)]>x[3:L]))
>           |((x[1:(L-2)]>x[2:(N-1)])&(x[2:(L-1)]<x[3:L])) ) + 1
>   }
> 
> Applied to your series "count" above:
> 
>   TP(count)
>   # [1]  2  4  6  7  9 11 14 17 21
> 
> If you assign these values to an index:
> 
>   ix <- TP(count)
>   rbind(d[ix],count[ix])
> 
>   # [1,]    8.0   10   12.0   13   15   17.0   20.0   23.0   27.0
>   # [2,] 3976.7 3955 3972.8 3961 4008 3981.2 4029.5 3848.3 3936.7
> 
> Of course, this is only a very simplistic view of "turning point",
> and will pick out everything which is a local minimum or maximum.
> 
> The above function can be extended (in a fairly obvious way) to
> identify each position i where x[i] is greater than its neighbours
> out to 2 on either side, or less than these neighbours; or more
> generally out to k on either side.
> 
> A lot depends on how you want to interpret "turning point". With
> your "count" series, it might be that you were only interested
> in identifying the relatively extreme turning points, such as
> i=4 (maybe), i=9 (maybe), i=14, i=17, i=21(maybe).
> 
> Hoping this helps,
> Ted.
> 
> --------------------------------------------------------------------
> E-Mail: (Ted Harding) <Ted.Harding at manchester.ac.uk>
> Fax-to-email: +44 (0)870 094 0861
> Date: 17-Sep-09                                       Time: 09:47:53
> ------------------------------ XFMail ------------------------------
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
>