[R] How to find local minimum between distributions using mixtools?

[Richard O'Keefe]

I presumed there was some reason why the 'optimise' function did not suit you.

optimize(function (x) a*dnorm(x, x1, s1) + (1-a)*dnorm(x, x2, s2),
range(c(x1, x2)))

should do the trick.

On Thu, 14 Oct 2021 at 23:42, Luigi Marongiu <marongiu.luigi using gmail.com> wrote:
>
> Thank you, I hoped there might be an automatic method more than
> function analysis...
>
> On Thu, Oct 14, 2021 at 9:06 AM Richard O'Keefe <raoknz using gmail.com> wrote:
> >
> > Do you really want the minimum?
> > It sounds as though your model is a*N(x1,s1) + (1-a)*N(x2,s2) where
> > you use mixtools to estimate
> > the parameters.  Finding the derivative of that is fairly
> > straightforward calculus, and solving for the
> > derivative being zero gives you extrema (you want the one between x1 and x2).
> >
> > In practice I might start with a quick and dirty hack like
> > x <- seq(from min(x1,x2), to=max(x1,x2), length=399)
> > p <- a*dnorm(x, x1, s1) + (1-a)*dnorm(x, x2, s2)
> > m <- min(p)
> > x[x == m]
> >
> > On Wed, 13 Oct 2021 at 22:12, Luigi Marongiu <marongiu.luigi using gmail.com> wrote:
> > >
> > > Hello,
> > > I have two peaks distributions (and some noise values in between), and
> > > I used mixtools to characterize the two overlapping data. How can I
> > > find the minimum between the peak distributions?
> > > The sigma value of the second population is close to that value, but I
> > > am not sure if this is really the cut-off point between the two
> > > distributions (by eye, I would say the value is 0.125 against 0.1 of
> > > the sigma value). Is there a better approach?
> > > Thanks
> > >
> > > ```
> > > set.seed(10)
> > > negative_mr <- runif(50, 0, 0.099)
> > > negative_fcn <- runif(50, 1, 40)
> > > positive_mr <- c(runif(30, 0.2, 0.5), runif(20, 0.4, 0.5))
> > > positive_fcn <- c(runif(30, 25, 40), runif(20, 10, 25))
> > > uncertain_mr <- runif(10, 0.099, 0.2)
> > > uncertain_fcn <- runif(10, 2, 40)
> > > df <- data.frame(Y=c(negative_mr, uncertain_mr, positive_mr),
> > > X=c(negative_fcn, uncertain_fcn, positive_fcn),
> > > GROUP=c(rep("negative", length(negative_mr)),
> > > rep("uncertain", length(uncertain_mr)),
> > > rep("positive", length(positive_mr))))
> > > library(mixtools)
> > > model = normalmixEM((x = df$Y))
> > > summary(model)
> > > # plot
> > > plot(model, which=2)
> > > Cut_off <- model$sigma[2]
> > > Cut_offE <- 0.125
> > > abline(v=Cut_off, col="blue", lwd=2)
> > > abline(v=Cut_offE, col="blue", lwd=2, lty=2)
> > > plot(df$Y[df$GROUP=="negative"]~df$X[df$GROUP=="negative"],
> > > pch=1, ylim=c(0,0.5), xlim=c(0,41), cex=1.5,
> > > xlab=expression(bold("X")),
> > > ylab=expression(bold("Y")))
> > > points(df$Y[df$GROUP=="positive"]~df$X[df$GROUP=="positive"],
> > > pch=16, cex=1.5)
> > > points(df$Y[df$GROUP=="uncertain"]~df$X[df$GROUP=="uncertain"],
> > > pch=16, cex=1.5, col="grey")
> > > abline(h=Cut_off, col="blue", lwd=2)
> > > abline(h=Cut_offE, col="blue", lwd=2, lty=2)
> > > ```
> > >
> > >
> > > --
> > > Best regards,
> > > Luigi
> > >
> > > ______________________________________________
> > > R-help using r-project.org mailing list -- To UNSUBSCRIBE and more, see
> > > 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.
>
>
>
> --
> Best regards,
> Luigi