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

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


statistics curve-fitting data-fitting centroid normal-distribution
person Sigh_at_psi    schedule 16.10.2018    source источник
comment
У меня есть пример Python для подгонки двух лоренцевских пиков к спектроскопии комбинационного рассеяния углерода углеродных нанотрубок с использованием генетического алгоритма для автоматической оценки исходных параметров здесь, замените данные и уравнения своими собственными, и это может работать довольно легко: bitbucket.org/zunzuncode/RamanSpectroscopyFit   -  person James Phillips    schedule 16.10.2018
comment
Для нахождения частоты в БПФ существует метод, использующий оконные функции Гаусса. В логарифмическом масштабе вы можете вычислить центр из двух последовательных интервалов. Это, однако, работает только для одного пика и лучше, если оно близко к максимальному. Вы могли бы попробовать больше от пика. Но если два гауссиана сильно перемешиваются, трудно понять, как избежать подгонки. Но я не вижу вычислительной сложности.   -  person mikuszefski    schedule 16.10.2018
comment
Это смесь моделей Гаусса. Обычный способ найти параметры такого типа модели — это так называемый алгоритм максимизации ожидания, который представляет собой просто обобщение уже имеющейся у вас процедуры вырезания и подгонки. Веб-поиск смеси гауссовых ожиданий и максимизации должен найти много ресурсов. Во многих языках программирования есть пакеты, реализующие этот алгоритм.   -  person Robert Dodier    schedule 16.10.2018

Ответы (1)


arrow_upward
0
arrow_downward

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

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

#!/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

гауссианцы

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

person mikuszefski    schedule 23.10.2018