punto di intersezione 2 curve normali

Question

Anche se penso che questa sia una domanda di base, non posso sembrare per scoprire come calcolare questo in R:

il punto di intersezione (ho bisogno del valore x) di 2 o più distribuzioni normali (montati su un istogramma), che hanno per esempio i seguenti parametri:

d=data.frame(mod=c(1,2),mean=c(14,16),sd=c(0.9,0.6),prop=c(0.6,0.4)) 

Con la media e la deviazione standard dei miei 2 curve, e puntellare le proporzioni di contributo di ciascun mod alla distribuzione.

Answer 1

È possibile utilizzare uniroot:

f <- function(x) dnorm(x, m=14, sd=0.9) * .6 - dnorm(x, m=16, sd=0.6) * .4 

uniroot(f, interval=c(12, 16)) 

$root 
[1] 15.19999 

$f.root 
[1] 2.557858e-06 

$iter 
[1] 5 

$estim.prec 
[1] 6.103516e-05 

---

ETA alcuni dell'esposizione:

uniroot è un cercatore di radice univariata, cioè data una funzione f di una variabile x, si trova il valore di x che risolve l'equazione f(x) = 0.

Per utilizzarlo, si forniscono la funzione f, insieme con un intervallo entro il quale il valore della soluzione si presume che giacere. In questo caso, f è solo la differenza tra i due densità; il punto in cui si intersecano sarà dove f è zero. Ho avuto l'intervallo (12, 16) in questo esempio facendo una trama e vedendo che intersecano attorno x = 15.

Answer 2

Ci dispiace, ma la risposta accettata non è buona. Vedi anche: intersection of two curves in matlab

È possibile ottenere entrambe le radici utilizzando una funzione come questa:

intersect <- function(m1, s1, m2, s2, prop1, prop2){ 

B <- (m1/s1^2 - m2/s2^2) 
A <- 0.5*(1/s2^2 - 1/s1^2) 
C <- 0.5*(m2^2/s2^2 - m1^2/s1^2) - log((s1/s2)*(prop2/prop1)) 

(-B + c(1,-1)*sqrt(B^2 - 4*A*C))/(2*A) 
} 

nel tuo caso:

> intersect(14,0.9,16,0.6,0.6,0.4) 
[1] 20.0 15.2 

Answer 3

Sono d'accordo con @Flounderer che si dovrebbe calcolare entrambe le radici, ma la funzione offerta è incompleta. Quando s1 = s2, A diventa zero e vengono prodotti Inf e NaN.

Una leggera modifica completa della funzione:

intersect <- function(m1, sd1, m2, sd2, p1, p2){ 

    B <- (m1/sd1^2 - m2/sd2^2) 
    A <- 0.5*(1/sd2^2 - 1/sd1^2) 
    C <- 0.5*(m2^2/sd2^2 - m1^2/sd1^2) - log((sd1/sd2)*(p2/p1)) 

    if (A!=0){ 
    (-B + c(1,-1)*sqrt(B^2 - 4*A*C))/(2*A) 
    } else {-C/B} 
} 
m1=0; sd1=2; m2=2.5; sd2=2; p1=.8; p2=.2 
(is=intersect(m1,sd1,m2,sd2,p1,p2)) 

xs = seq(-6, 6, by=.01) 
plot(xs, p1*dnorm(xs, m1, sd1), type= 'l') 
lines(xs, .2*dnorm(xs, m2,sd2)) 
abline(v=is) 

una soluzione generale si trova anche utilizzando uniroot.all:

library(rootSolve) 
f <- function(x, m1, sd1, m2, sd2, p1, p2){ 
    dnorm(x, m1, sd1) * p1 - dnorm(x, m2, sd2) * p2 } 
uniroot.all(f, lower=-6, upper=6, m1=m1, sd1=sd1, m2=m2, sd2=sd2, p1=p1, p2=p2)