Попытка подобрать распределение Гаусса: ошибка scipy/optimize/minpack.py", строка 765, в curve_fit поднять ValueError ("`sigma` имеет неправильную форму.")

У меня довольно известная проблема, но решить ее пока невозможно. Это о функции curve_fit: я получаю сообщение об ошибке:

Error scipy/optimize/minpack.py", line 765, in curve_fit raise ValueError("`sigma` has incorrect shape.")

Вот код, не обращайте внимания на цикл, просто мне нужно 5 разных гистограмм:

for i in range(5):
  mean_o[i] = np.mean(y3[:,i])
  sigma_o[i] = np.std(y3[:,i])

## Histograms
# Number of bins
Nbins=100
binwidth = np.zeros(5)

# Fitting curves
def gaussian(x, a, mean, sigma):
  return a * np.exp(-((x - mean)**2 / (2 * sigma**2)))

for i in range(5):

  binwidth[i] = (max(y3[:,i]) - min(y3[:,i]))/Nbins
  bins_plot = np.arange(min(y3[:,i]), max(y3[:,i]) + binwidth[i], binwidth[i])
  plt.title('Distribution of O observable for redshift bin = '+str(z_ph_fid[i]))
  plt.hist(y3[:,i], bins=bins_plot, label='bin '+str(z_ph_fid[i]))
  plt.legend(loc='upper right')
  # Fitting and plot
  range_fit = np.linspace(min(y3[:,i]), max(y3[:,i]), len(y3[:,i]))
  popt, pcov = curve_fit(gaussian, range_fit, y3[:,i], mean_o[i], sigma_o[i])
  plt.plot(range_fit, gaussian(range_fit, *popt))
  # Save figure
  plt.savefig('chi2_gaussian_bin_'+str(i+1)+'.png')
  plt.close()

Первая гистограмма i=0 выглядит так:

Попытка подобрать распределение Гаусса: ошибка scipy/optimize/minpack.py", строка 765, в curve_fit поднять ValueError ("`sigma` имеет неправильную форму.")

Я хотел бы построить красную гауссову аппроксимацию гистограммы.

Если бы кто-нибудь мог помочь мне с этой ошибкой...

python scipy curve-fitting gaussian standard-deviation

08.05.2022 17:00

Почему в Python есть оператор "pass"?

Оператор pass в Python - это простая концепция, которую могут быстро освоить даже новички без опыта программирования.

Некоторые методы, о которых вы не знали, что они существуют в Python

Python - самый известный и самый простой в изучении язык в наши дни. Имея широкий спектр применения в области машинного обучения, Data Science,...

Основы Python Часть I

Вы когда-нибудь задумывались, почему в программах на Python вы видите приведенный ниже код?

LeetCode - 1579. Удаление максимального числа ребер для сохранения полной проходимости графа

Алиса и Боб имеют неориентированный граф из n узлов и трех типов ребер:

Оптимизация кода с помощью тернарного оператора Python

И последнее, что мы хотели бы показать вам, прежде чем двигаться дальше, это

Советы по эффективной веб-разработке с помощью Python

Как веб-разработчик, Python может стать мощным инструментом для создания эффективных и масштабируемых веб-приложений.

Перейти к ответу Данный вопрос помечен как решенный

Ответы 1

Ответ принят как подходящий

Есть две проблемы с ОП.

Первая проблема заключается в том, что код пытается подогнать случайные выборки к нормальному распределению. Это не правильно. Однако можно подогнать вывод гистограммы. Это показано в коде ниже. Лучше было бы использовать scipy.stats.norm.fit(), что позволяет подогнать случайные выборки. Это тоже показано.

Вторая проблема - сигма-форма. Здесь curve_fit на самом деле ожидает ошибок в y-данных, что, естественно, требует формы y-данных. Что нужно было сделать, так это: предоставить начальные значения для подгонки. Это также показано ниже.

Код выглядит так:

import matplotlib.pyplot as plt
import numpy as np

from scipy.stats import norm
from scipy.optimize import curve_fit

mean_o = list()
sigma_o = list()
y3 = list()

### generate some data
for i in range( 5 ):
    y3.append( norm.rvs( size=150000 ) )
y3 = np.transpose( y3 )

for i in range(5):
    mean_o.append( np.mean( y3[ :, i ] ) )
    sigma_o.append(  np.std( y3[ :, i ] ) )

## Histograms
# Number of bins
Nbins=100
binwidth = np.zeros(5)

# Fitting curves
def gaussian( x, a , mean, sigma ):
  return a * np.exp( -( ( x - mean )**2 / (2 * sigma**2 ) ) )

fig = plt.figure()
ax = { i : fig.add_subplot( 2, 3, i+1) for i in range( 5 ) }


for i in range(5):
    ymin = min(y3[:,i])
    ymax = max(y3[:,i])
    binwidth[i] = ( ymax - ymin) / Nbins
    bins_plot = np.arange( ymin, ymax + binwidth[i], binwidth[i])
    histdata = ax[i].hist(
        y3[:,i],
        bins=bins_plot,
        label='bin '+str(i)
    )

    range_fit = np.linspace( ymin, ymax, 250)
    # Fitting and plot version 1
    popt, pcov = curve_fit(
        gaussian,
        0.5 * ( histdata[1][:-1] + histdata[1][1:] ),
        histdata[0],
        p0=( max(histdata[0]), mean_o[i], sigma_o[i] ) )
    ax[i].plot(range_fit, gaussian( range_fit, *popt ) )
    ax[i].axvline( x=mean_o[i], ls=':', c='r' )

    # Fitting and plot version 2
    params = norm.fit( y3[ ::, i ], loc=mean_o[i], scale=sigma_o[i] )
    nth = gaussian(
        range_fit,
        len( y3[::, i]) * binwidth[i] / np.sqrt( 2 * np.pi ),
        *params
    )
    ax[i].plot(range_fit, nth, ls = "--" )
plt.tight_layout()
plt.show()

Большое спасибо, это работает! просто деталь, как поместить наилучшее значение на оси ax[i]? Я имею в виду, как добавить mean_o[i] по оси x каждого графика? С наилучшими пожеланиями

— 09.05.2022 13:47

@stefan Я сделаю обновление .... также забыл упомянуть важный момент.

— 10.05.2022 08:51

Спасибо, я не вижу модификаций для красной вертикальной линии (наиболее подходит), это то, что вы говорите? вы сделаете обновление? С уважением

— 10.05.2022 15:44

@стефан... извини. добавлено: ax[i].axvline( x=mean_o[i], ls=':', c='r' )

— 11.05.2022 08:55

09.05.2022 09:54