[R] which.minimums not which.min

[Philippe Grosjean]

Fred J. wrote:
> 
> Philippe Grosjean <phgrosjean at sciviews.org> wrote: What Fred is looking for is local minima/maxima, also known as turning 
> points, or pits/peaks in a series.  You can look at ?turnpoints in 
> pastecs library.
> 
>  > x <- c(1:4,0:5, 4, 11)
>  > x
>   [1]  1  2  3  4  0  1  2  3  4  5  4 11
>  > tp <- turnpoints(x)
>  > summary(tp)
> Turning points for: x
> 
> nbr observations  : 12
> nbr ex-aequos     : 0
> nbr turning points: 4 (first point is a peak)
> E(p) = 6.666667 Var(p) = 1.811111 (theoretical)
> 
>    point type       proba      info
> 1     4 peak 0.100000000 3.3219281
> 2     5  pit 0.002380952 8.7142455
> 3    10 peak 0.005952381 7.3923174
> 4    11  pit 0.666666667 0.5849625
>  > plot(tp) # Only useful for a longer and more complex series!
>  > # Get the position of peaks
>  > (1:length(x))[extract(tp, no.tp = FALSE, peak = TRUE, pit = FALSE)]
> [1]  4 10
> Warning message:
> arguments after the first two are ignored in: UseMethod("extract", e, n, 
> ...)
>  > (1:length(x))[extract(tp, no.tp = FALSE, peak = FALSE, pit = TRUE)]
> [1]  5 11
> Warning message:
> arguments after the first two are ignored in: UseMethod("extract", e, n, 
> ...)
>  > # By the way, there are warnings although it works well (I ask on R-Help)
> 
> Now, you can easily code your which.minima() function using turnpoints:
> 
> x <- c(1:4,0:5, 4, 11)
> x
> tp <- turnpoints(x)
> summary(tp)
> plot(tp) # Only useful for a longer and more complex series!
> # Get the position of peaks
> (1:length(x))[extract(tp, no.tp = FALSE, peak = TRUE, pit = FALSE)]
> (1:length(x))[extract(tp, no.tp = FALSE, peak = FALSE, pit = TRUE)]
> # By the way, there are warnings although it works well (I ask on R-Help)
> 
> which.minima <- function(x) {
>  if (!require(pastecs)) stop("pastecs library is required!")
>  x <- as.vector(x)
>  (1:length(x))[extract(turnpoints(x), no.tp = FALSE, peak = FALSE, pit = 
> TRUE)]
> }
> 
> which.minima(x)
> 
> Of course, you could optimize this code. This is just a rough solution!
> Best,
> 
> Philippe Grosjean
> 
> ..............................................<�}))><........
>   ) ) ) ) )
> ( ( ( ( (    Prof. Philippe Grosjean
>   ) ) ) ) )
> ( ( ( ( (    Numerical Ecology of Aquatic Systems
>   ) ) ) ) )   Mons-Hainaut University, Pentagone (3D08)
> ( ( ( ( (    Academie Universitaire Wallonie-Bruxelles
>   ) ) ) ) )   8, av du Champ de Mars, 7000 Mons, Belgium
> ( ( ( ( (
>   ) ) ) ) )   phone: + 32.65.37.34.97, fax: + 32.65.37.30.54
> ( ( ( ( (    email: Philippe.Grosjean at umh.ac.be
>   ) ) ) ) )
> ( ( ( ( (    web:   http://www.umh.ac.be/~econum
>   ) ) ) ) )          http://www.sciviews.org
> ( ( ( ( (
> ..............................................................
> 
>   Is there a limit as far as the data input to the function turnpoints? It seams to be unable to process a vector with 15000 records with number in the form of 734.2983, the data is known to have picks and pits
>    
>   The steps: 
>   Where pp is the vector with data of the above-mentioned style.
>   
>  >tp <- turnpoints(rev(pp))
>   > mx_id_rev <- (1:length(pp))[extract(tp, no.tp = FALSE, peak = TRUE, pit = FALSE)]
>   > tp
>   Turning points for: rev(pp) 
>    
>   nbr observations  : 11745 
>   nbr ex-aequos     : 10273 
>   nbr turning points: 1231 (first point is a pit)
>   E(p) = 7828.667 Var(p) = 2087.678 (theoretical)
>   > mx_id_rev
>   numeric(0)
>    
>   Regards,

There is certainly a limit, but not that low. My test gives (R under Win
XP with 1Gb RAM):

> x <- rnorm(15000)
> library(pastecs)
Loading required package: boot
> tp <-turnpoints(rev(x))
> tp
Turning points for: rev(x)

nbr observations  : 15000
nbr ex-aequos     : 0
nbr turning points: 10079 (first point is a peak)
E(p) = 9998.667 Var(p) = 2666.344 (theoretical)

> res <- (length(x):1)[extract(tp, no.tp = FALSE, peak = TRUE, pit = FALSE)]
> length(res)
[1] 6250
> res <- (length(x):1)[extract(tp, n= 15000, no.tp = FALSE, peak = TRUE, pit = FALSE)]
> length(res)

[1] 6250
[1] 5040
 > plot(x[1:500], type = "l")
 > points(res, x[res], col = 2)

That seems to work. So, there is a bug when you don't specify 'n'. I 
work on a patch for that and upload the new version of pastecs package 
soon. Thanks for pointing me this bug.
Best,

Philippe Grosjean