2013-06-07 8 views
7

Chociaż myślę, że jest to podstawowe pytanie, nie mogę wydawać się, aby dowiedzieć się, jak obliczyć to w R:punkt przecięcia 2 normalnych krzywych

punktem przecięcia (muszę X-wartości) 2 lub więcej rozkładu normalnego (montowane na histogramie), które mają, na przykład następujące parametry:

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

ze średnią i odchylenie standardowe z 2 moich krzywych i prop proporcji udziału każdego mod rozkładu.

Odpowiedz

12

Można użyć 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 niektóre ekspozycja:

uniroot jest jednowymiarowy korzeń finder, czyli dana funkcja f jednej zmiennej x, stwierdzi wartość x że rozwiązuje równanie f(x) = 0.

Aby z niego skorzystać, należy podać funkcję f wraz z przedziałem, w którym zakłada się wartość rozwiązania. W tym przypadku, f jest tylko różnicą między dwiema gęstościami; punkt, w którym się przecinają, to miejsce, gdzie f jest równe zero. Dostałem interwał (12, 16) w tym przykładzie, wykonując fabułę i widząc, że przecinały się wokół x = 15.

+0

+1, ale można dodać trochę wyjaśnienie co to robi/jak to działa? Dzięki –

+0

Dzięki za tekst. To jest świetne! –

+0

dzięki, działa idealnie !! – Wave

0

Przepraszamy, ale zaakceptowana odpowiedź nie jest dobra. Zobacz także: intersection of two curves in matlab

można uzyskać zarówno korzenie przy użyciu funkcji takich jak to:

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) 
} 

w przypadku:

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

Jak wybrać najbardziej znaczący punkt? – Spidfire

0

zgadzam się z @Flounderer że należy obliczyć zarówno korzenie, ale oferowana funkcja jest niekompletna. Gdy s1 = s2, A staje się zerem, a Inf i NaN są wytwarzane.

Nieznaczna modyfikacja dopełnia funkcję:

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) 

ogólne rozwiązanie jest również znaleźć za pomocą 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) 
Powiązane problemy