Нахождение пересечения гауссовых распределений - PullRequest
0 голосов
/ 26 августа 2018

Как мне найти все точки, где распределения гауссовой смеси пересекаются в MATLAB?enter image description here

1 Ответ

0 голосов
/ 26 августа 2018

Общей категорией вашего вопроса является поиск пересечения двух кривых, что является управляемой, но нетривиальной задачей (самое сложное - убедиться, что вы перехватили все пересечения).

Но ваша проблема заключается в том, очень специфично: вы ищете пересечение двух гауссиан.Это очень хорошо: у нас обоих есть аналитическая формула для вашей функции, и гарантировано, что есть ровно два пересечения, если некоторые из ваших параметров не совпадают. *

Давайте предположим, что вашраспределения характеризуются с помощью mu1, mu2 и шкалы sigma1, sigma2.Тогда ваши гауссианы в положении x определяются функцией

1/sqrt(2*pi*sigma^2) * exp(-(x-mu)^2/2/sigma^2)

Оказывается, что мы можем полностью решить это уравнение для x на бумаге:

1/sqrt(2*pi*sigma1^2)*exp(-(x-mu1)^2/2/sigma1^2) == 1/sqrt(2*pi*sigma2^2)*exp(-(x-mu2)^2/2/sigma2^2)
sigma2/sigma1 == exp((x-mu1)^2/2/sigma1^2) * exp(-(x-mu2)^2/2/sigma2^2)
log(sigma2/sigma1) == (x-mu1)^2/2/sigma1^2) - (x-mu2)^2/2/sigma2^2

, что приводит кпараболическое уравнение ax^2 + bx + c == 0, где

a = 1/(2*sigma1^2) - 1/(2*sigma2^2);
b = mu2/(sigma2^2) - mu1/(sigma1^2);
c = mu1^2/(2*sigma1^2) - mu2^2/(2*sigma2^2) - log(sigma2/sigma1);

. Можно легко доказать, что дискриминант D = b^2 - 4 a c неотрицателен, поэтому в действительности уравнение имеет два действительных корня, когда параметры не вырождаются.Таким образом, две точки пересечения, с приведенными выше определениями,

D = b^2 - 4 * a * c;
x1 = (-b + sqrt(D))/(2*a);
x2 = (-b - sqrt(D))/(2*a);

Использование двух гауссианов с псевдослучайными параметрами:

% define parameters and Gaussian
mu1=1; sigma1=3; mu2=2; sigma2=4;

% intersections
a = 1/(2*sigma1^2) - 1/(2*sigma2^2);
b = mu2/(sigma2^2) - mu1/(sigma1^2);
c = mu1^2/(2*sigma1^2) - mu2^2/(2*sigma2^2) - log(sigma2/sigma1)
D = b^2 - 4 * a * c;
x1 = (-b + sqrt(D))/(2*a);
x2 = (-b - sqrt(D))/(2*a);

Просто чтобы доказать, что вышеуказанные точки пересечения являются правильными:

>> f = @(x,mu,sigma) 1/sqrt(2*pi*sigma^2) * exp(-(x-mu).^2/2/sigma^2);
>> f(x1,mu1,sigma1) - f(x1,mu2,sigma2)

ans =

   2.7756e-17

>> f(x2,mu1,sigma1) - f(x2,mu2,sigma2)

ans =

   1.0408e-17

Вышесказанное означает, что значения двух гауссианов в точках x1 и x2 равны друг другу в пределах точности машины, что так же хорошо, как любой числовой ответ.


* Первоначально я утверждал, что у нас всегда есть два пересечения, если гауссианы не совсем одинаковы.Ясно, что если и средние значения, и дисперсии одинаковы, мы имеем две вырожденные кривые и пересечения становятся бессмысленными.Но, как заметил Крис Луенго в комментарии , может случиться так, что есть только одно пересечение: когда дисперсии одинаковы, а средние значения различны (то есть у нас есть две кривые одинаковой формы, смещенные вдоль x).В этом случае a=0, следовательно, соответствующее уравнение равно b*x + c == 0, что дает нам x0 = -c/b для пересечения.Таким образом, более точный (но немного псевдокодовый) ответ (с учетом a, b и c) составляет

if a == 0 % or allow some tolerance... <=> sigma1 == sigma2
   if b == 0 % or allow some tolerance... <=> mu1 == mu2
      % degenerate curves: a == b == c == 0, f1(x)==f2(x) for all x
      disp('curves are degenerate...')
   else
      % single intersection: mu1 ~= mu2
      x1 = -c/b;
   end
else
   % two intersections; both parameters are different
   D = b^2 - 4 * a * c;
   x1 = (-b + sqrt(D))/(2*a);
   x2 = (-b - sqrt(D))/(2*a);
...