Блог

Мне трудно понять, почему моя подгонка по Гауссу к набору данных ( ydata ) не работает должным образом, если я сдвигаю интервал значений x, соответствующих этим данным ( xdata1 на xdata2 ). Гауссово записывается как:

где A — это просто коэффициент амплитуды. Изменяя некоторые значения данных, легко заставить его работать в обоих случаях, но также можно легко найти случаи, в которых он плохо работает xdata1 , а также в которых ковариация параметров не оценивается. Я использую scipy.optimize.curve_fit в Spyder с Python 3.7.1 в Windows 7.

 import numpy as np
from scipy.optimize import curve_fit
import matplotlib.pyplot as plt

xdata1 = np.linspace(-9,4,20, endpoint=True) # works fine
xdata2 = xdata1 2
ydata = np.array([8,9,15,12,14,20,24,40,54,94,160,290,400,420,300,130,40,10,8,4])

def gaussian(x, amp, mean, sigma):
    return amp*np.exp(-(((x-mean)**2)/(2*sigma**2)))/(sigma*np.sqrt(2*np.pi))

popt1, pcov1 = curve_fit(gaussian, xdata1, ydata)
popt2, pcov2 = curve_fit(gaussian, xdata2, ydata)

fig, ([ax1, ax2]) = plt.subplots(nrows=1, ncols=2,figsize=(9, 4))

ax1.plot(xdata1, ydata, 'b :', label='xdata1')
ax1.plot(xdata1, gaussian(xdata1, *popt1), 'r-', label='fit')
ax1.legend()
ax2.plot(xdata2, ydata, 'b :', label='xdata2')
ax2.plot(xdata2, gaussian(xdata2, *popt2), 'r-', label='fit')
ax2.legend()
  

Проблема в том, что ваша вторая попытка подгонки по гауссу застревает в локальном минимуме при поиске в пространстве параметров: curve_fit — это оболочка для least_squares, которая использует градиентный спуск для минимизации функции затрат, и это может застрять в локальных минимумах.

Вы должны попытаться предоставить разумные начальные параметры (используя p0 аргумент curve_fit), чтобы избежать этого:

  #...  your code

 y_max = np.max(y_data)
 max_pos = ydata[ydata==y_max][0]
 initial_guess = [y_max, max_pos, 1] # amplitude, mean, std

 popt2, pcov2 = curve_fit(gaussian, xdata2, ydata, p0=initial_guess)
  

Что, как вы можете видеть, обеспечивает разумную подгонку:

Вы должны написать функцию, которая может обеспечить разумные оценки начальных параметров. Здесь я просто нашел максимальное значение y и использовал это для определения начальных параметров. Я обнаружил, что это хорошо работает для подходящих нормальных распределений, но вы могли бы рассмотреть другие методы.

Редактировать:

Вы также можете решить проблему путем масштабирования амплитуды: амплитуда настолько велика, что пространство параметров искажается, а градиентный спуск просто следует направлению наибольшего изменения амплитуды и эффективно игнорирует сигму. Рассмотрим следующий график в пространстве параметров (цвет представляет собой сумму квадратов остатков подгонки по заданным параметрам, а белый крест показывает оптимальное решение):

Обязательно обратите внимание на разные масштабы для осей x и y.

Нужно сделать большое количество шагов размером с единицу в y (амплитуда), чтобы добраться до минимума из точки x, y = (0,0), где, поскольку вам нужно всего лишь меньше одного шага размером с единицу, чтобы добраться до минимума в x (сигма). Алгоритм просто выполняет шаги по амплитуде, поскольку это самый крутой градиент. Когда он достигает амплитуды, которая минимизирует функцию затрат, он просто останавливает алгоритм, поскольку он, по-видимому, сошелся, и практически не вносит изменений в параметр sigma.

Один из способов исправить это — масштабировать ваши ydata, чтобы не искажать пространство параметров: разделите ваши ydata на 100, и вы увидите, что ваша подгонка работает без предоставления каких-либо начальных параметров!