2

I calculate the density function for the following data:

> dput(mydat)  
c(-20, -13, 30, 4, -4, 34, 27, 19, 13.5, 15, 13, 18, 10, 12, 
21, -0.769999999999996, 2.5, -7, 0, -30.6, 6.39999999999999, 
-18.6, -0.199999999999989, -20.4, -19.9, 4.60000000000001, -19.4, 
4.5, -9, -15, 9, -1, -14, 8, 6, -17, 5, 7)  

> myden = density(mydat)  # default kernel and bandwidth  

which gives me this result:

enter image description here

I want to find the location of the two density peaks. I initially thought of using diff() on myden$y and then check for all locations where there is a sign change, using this as a condition to select the X-axis values. I tried it on a few test vectors but I was not getting the expected result and I suspect it is not that simplistic.

Is there a straightforward way to accomplish this? I want a solution that will be repeatable because I will do this as part of a random simulation study with ~ e+05 realisations, and it could happen that the number of peaks vary across the simulation.

avg
  • 793
  • 2
  • 12
  • 24
  • Why do you think the diff() approach isn't working as expected? I just saw the diff output and where the sign changes from +ve to -ve seems to be around a peak. Diff should approximate the differential sort of logic one would normally apply to get the local maxima. – TheComeOnMan Oct 02 '13 at 11:18
  • @Thomas thanks! yes, its a duplicate. it will have to be flagged as such and i dont know if i have the privileges to do that. – avg Oct 02 '13 at 11:20
  • @Codoremifa after looking at the link Thomas posted, i remembered from calculus class that one has to test this with the second differential.. – avg Oct 02 '13 at 11:24

2 Answers2

2

Use which.max:

myden$x[which.max(myden$y)]
# [1] 5.91428

You can test this visually:

plot(myden, col='red')
abline(v=myden$x[which.max(myden$y)])

enter image description here

Thomas
  • 43,637
  • 12
  • 109
  • 140
  • I want to find the value for both (or multiple, if there are more) peaks, not just the global maximum.. – avg Oct 02 '13 at 11:08
  • @AdvaitGodbole Ah, then see the question and answers I just posted on your original question. – Thomas Oct 02 '13 at 11:11
2

I often use pastecs::turnpoints to find local maxima and minima.

Carl Witthoft
  • 20,573
  • 9
  • 43
  • 73
  • E.g., http://stats.stackexchange.com/q/30750/11849 – Roland Oct 02 '13 at 11:36
  • @Roland where's that "I posted it" disclaimer? :-) //shuffles off to link to his own `turnpoints` posts *_O – Carl Witthoft Oct 02 '13 at 11:40
  • `turnpoints` works like a charm. However, I can imagine that for noisy data (similar to what is written on this post http://stackoverflow.com/questions/14319826/finding-local-maxima-minima-in-r/14320172#14320172), I might have to look for a work-around..But thats a different matter and for now, this should work. – avg Oct 02 '13 at 12:03
  • 1
    If your data is noisy, just smooth it. – Roland Oct 02 '13 at 12:28