[R] How to identify the two largest peaks in a trimodal distribution
jim holtman
jholtman at gmail.com
Sat Oct 13 19:39:21 CEST 2007
You can use 'rle' to find where the direction changes and use that to
determine the peaks:
> # create test data
> x <- seq(-2,7, length=2000)
> y <- dnorm(x) + dnorm(x,3)
> plot(y,type='l')
> # find where direction changes from plus to minus
> z <- rle(diff(y) > 0) # find where breaks are
> z
Run Length Encoding
lengths: int [1:4] 452 325 325 897
values : logi [1:4] TRUE FALSE TRUE FALSE
> change.index <- cumsum(c(1, z$lengths)) # get index where change occurs
> change.index
[1] 1 453 778 1103 2000
> # print out index of max
> which.one <- change.index[which(z$values ==FALSE)]
> which.one
[1] 453 1103
> y[which.one]
[1] 0.4036175 0.4036175
>
On 10/13/07, Rob Knell <R.Knell at qmul.ac.uk> wrote:
> Hello all
>
> I'm trying to do a simulation that involves identifying the minimum
> point between two peaks of a (usually) bimodal distribution. I can do
> this easily if there are only two peaks:
>
> CnBdens<-density(Ys/Xs) #probability density function for ratio of Ys
> to Xs
>
> for(p in 1:512) ifelse(CnBdens$y[p]>CnBdens$y[p-1],peak1<-p,break)
> #identifies first peak in probability distribution
>
> for(p in 1:512) ifelse(CnBdens$y[512-p]>CnBdens$y[512-p
> +1],peak2<-512-p,break) #identifies second peak in probability
> distribution
>
> but the simulation sometimes produces a small third peak at one end
> of the distribution. Is there any simple way to identify the two
> highest maxima in a trimodal distribution?
>
> Thanks for any help
>
> Rob Knell
>
>
> School of Biological and Chemical Sciences
> Queen Mary, University of London
>
> 'Phone +44 (0)20 7882 7720
> Skype Rob Knell
>
> Research: http://www.qmw.ac.uk/~ugbt794
> Giant edible caterpillars: http://www.mopane.org
> Invertebrate macro photography: http://web.mac.com/rknell/iWeb/Bugsite
>
> "The truth is that they have no clue why the beetles had horns, it's
> the researchers who have sex on the brain and everything has to have
> a sexual explanation. And this is reasearch?!" Correspondent known as
> FairOpinion on Neo-Con website discussing my research.
>
>
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> 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.
>
--
Jim Holtman
Cincinnati, OH
+1 513 646 9390
What is the problem you are trying to solve?
More information about the R-help
mailing list