Вычисление центроидов двух наложенных гауссовских функций - PullRequest
Новая
       103

Вычисление центроидов двух наложенных гауссовских функций

0 голосов
/ 16 октября 2018

Я пытаюсь найти решение следующей проблемы.У меня есть набор точек, которые должны моделировать сумму 2 гауссовских функций с центром в разных точках.Мне нужно найти эти две точки.До сих пор мой подход состоял в том, чтобы найти центр тяжести всего набора и вырезать набор дат ниже и выше него.Затем я вычисляю центр тяжести каждой части, и это мои центры.Этот подход, однако, отсекает информацию, скажем, левого гауссиана, которая просачивается в правую половину данных.Это делает процедуру неудачной, когда гауссиане находятся близко друг к другу.Есть ли способ сделать это более разумно?Из-за сложности вычислений я бы предпочел, чтобы решение не включало подгонку кривой.

1 Ответ

0 голосов
/ 23 октября 2018

Поскольку ОП не показывает никаких данных, неясно, насколько шумны эти данные.Кроме того, не ясно, как здесь определяется «близко друг к другу».В следующем у меня есть простое приближение, которое работает с низким уровнем шума и предположением, что в данных левой части преобладает левый гауссов, а в правой части преобладает правый гауссов.Это дает некоторые ограничения на положение, высоту и особенно стандартное отклонение.

Это, безусловно, работает для одного пика, но вполне подходит для смешанных двойных пиков (в пределах вышеупомянутых ограничений) 

#!/usr/bin/python
import matplotlib.pyplot as plt
import numpy as np

def gaussian( x, x0, s, a):
    return a * np.exp( -( x - x0 )**2 / ( 2 * s**2 ) ) / np.sqrt( 2 * np.pi * s**2 )

def get_x0( x1, x2, x3, y1, y2, y3 ):
    l12= np.log( y1 / y2 )
    l13= np.log( y1 / y3 )
    return ( ( x2 + x1 )/2. - ( x3 + x1 )/2. * l12/l13 * ( x3 - x1 ) / ( x2 - x1 ) ) / ( 1 - l12 / l13 * (x3 - x1 ) / ( x2 - x1 ) )


fig = plt.figure( )
ax = fig.add_subplot( 2, 1, 1 )

xL = np.linspace(-8, 8, 150 )
yL = np.fromiter( ( gaussian( x,-2.1, 1.2, 8 ) for x in xL ), np.float )
marker=[10,15,20]

x1 = xL[ marker[0] ]
x2 = xL[ marker[1] ]
x3 = xL[ marker[2] ]
y1 = yL[ marker[0] ]
y2 = yL[ marker[1] ]
y3 = yL[ marker[2] ]

print get_x0( x1, x2, x3, y1, y2, y3 )

ax.plot( xL, yL )
ax.scatter( [ x1, x2, x3 ],[ y1, y2, y3 ])

bx = fig.add_subplot( 2, 1, 2 )
yL = np.fromiter( ( gaussian( x,-2.1, 1.2, 8) + gaussian( x,0.7, 1.4, 6 ) for x in xL ), np.float )
marker=[10,15,20]
x1 = xL[ marker[0] ]
x2 = xL[ marker[1] ]
x3 = xL[ marker[2] ]
y1 = yL[ marker[0] ]
y2 = yL[ marker[1] ]
y3 = yL[ marker[2] ]
bx.scatter( [ x1, x2, x3 ],[ y1, y2, y3 ])
print get_x0( x1, x2, x3, y1, y2, y3 )
marker=[-20,-25,-30]
x1 = xL[ marker[0] ]
x2 = xL[ marker[1] ]
x3 = xL[ marker[2] ]
y1 = yL[ marker[0] ]
y2 = yL[ marker[1] ]
y3 = yL[ marker[2] ]
bx.scatter( [ x1, x2, x3 ],[ y1, y2, y3 ])
print get_x0( x1, x2, x3, y1, y2, y3 )
bx.plot( xL, yL )

plt.show()

Показывает: 

#Single
-2.0999999999999455
#Double
-2.0951188129317813
0.6998760921436634

, что довольно близко к -2.1 и 0.7

gaussians

В случае шумаможет потребоваться некоторое усреднение.

...