[R] How to find local minimum between distributions using mixtools?
Luigi Marongiu
m@rong|u@|u|g| @end|ng |rom gm@||@com
Thu Oct 14 12:42:27 CEST 2021
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
More information about the R-help
mailing list