Свертка случайно распределенных данных с гауссианами

Question

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

Я хочу заменить каждую точку данных y [x_i] с фиксированной ширинойГауссова и сложить их вместе.Он должен дать мне:

Мой код очень примитивен и медленен:

def gaussian(x, mu, sig):
    return 1/(sig*np.sqrt(2*np.pi))*np.exp(-np.power(x - mu, 2.) / (
        2 * np.power(sig, 2.)))

def gaussian_smoothing(x, y, sig=0.5, n=1000):

    x_new = np.linspace(x.min()-10*sig, x.max()+10*sig, n)
    y_new = np.zeros(x_new.shape)

    for _x, _y in zip(x, y):
        y_new += _y*gaussian(x_new, _x, sig)

    return x_new, y_new

Для больших наборов данных требуетсядолгое время выполнять такое сглаживание.Я смотрел на np.convolve.Однако он показывает, что он применим только к равномерно распределенным данным, и шаг x для данных и гауссиан должны быть одинаковыми.Какой бы самый быстрый способ выполнить такую ​​операцию.

Answer 1

yon попытайтесь оценить это как смесь Guassian с меньшим количеством компонентов (как EM-алгоритм), используя sklearn:

import matplotlib.pyplot as plt
from numpy.random import choice
from sklearn import mixture
import scipy.stats
import numpy

# generate some data
x = numpy.array([1.,1.1,1.6,2.,2.1,2.2,2.9,3.,8.,62.,62.2,63.,63.4,64.5,65.,67.,69.])
# generate weights to it
y = numpy.random.rand(x.shape[0])
# normalize weigth to 1
y /= y.sum()  

# resamlple to 5000 samples with equal weights according to original weights
x_rsmp = numpy.array([choice(x, p=y) for _ in range(5000)])
x_rsmp.sort()
x_rsmp = x_rsmp.reshape(-1,1)

# define number of components - this must be user seelcted or estimated
n_comp = 2
# fit the mixture
gmm = mixture.GaussianMixture(n_components=n_comp, covariance_type='full')
gmm.fit(x_rsmp)

# plot it
fig = plt.figure()
ax = fig.add_subplot(111)
x_gauss = numpy.linspace(-10,100,1000)

for n_c in range(n_comp):
    norm_pdf = scipy.stats.norm.pdf(x_gauss, gmm.means_[n_c,0], gmm.covariances_[n_c,0])
    ax.plot(x_gauss, norm_pdf, label='gauss %d' % (n_c+1))

ax.stem(x,y,'gray')
plt.legend()

Это дает n_c гауссовские компоненты со средним значением gmm.means_ и ковариациями gmm.covariances_.