Генерация chi_squared с использованием наилучшего соответствия

Я генерирую 1000 итераций некоторого поддельного набора данных и использую curve_fit, чтобы найти наилучшую подходящую модель (гауссиан с определенным средним значением, смещением, усилителем...). Вместо того, чтобы угадывать среднее значение моей подходящей модели, используйте curve_fit, чтобы я мог имеют лучшие значения chi2. Но теперь вывод странный, и даже данные и модель не согласуются.Это точечные диаграммы моих данных и подходящей модели

Кроме того, если я построю гистограмму своих значений chi2, Значения chi2_fit не распространяются, как ожидалось

true_mean = 5 #Assume that we know the perfect fitting model is gaussian 
#with mean value 5. And use try_mean with different mean values to see how 
#does the chi2 behave
true_amp = 100
true_wid = 2
true_offset =10
x_values =  np.array([i for i in np.arange(0,10,0.4)])
exact_y_values = np.array([true_offset + true_amp*
                               np.exp(-((i-true_mean)/true_wid)**2)
                               for i in x_values])
    

    
try_mean = 4.8   # Notice the data is generated centered at 5;
# by comparing it to a 4.8 we expect disagreement.
try_amp = 100
try_wid = 2  
try_offset =10
try_these_y_values = np.array([try_offset + try_amp*
                                   np.exp(-((i-try_mean)/try_wid)**2)
                                   for i in x_values])

def func (x_values,offset,amp,mean,wid):
    return (offset + amp*np.exp(-((i-mean)/wid)**2))
#Excercise 2
#def func (x_values,offset,amp,mean,wid): return (offset + amp*np.exp(-((i-mean)/wid)**2))

chi2_fit=np.zeros(1000)
for i in range (1000):

    fake_exp_y_values = np.array([np.random.poisson(y)
                                      for y in exact_y_values])
    p01=[true_offset,true_amp,true_mean,true_wid]
    [popt,pcov]=opt.curve_fit(func,x_values, fake_exp_y_values,p01)
    y_values_fit=np.array([popt[0]+ popt[1]
                           *np.exp(
                               -((x_values-popt[2])/popt[3])**2)])
    residuals=fake_exp_y_values-y_values_fit
    y_err=np.clip(np.sqrt(fake_exp_y_values),1,9999)
    pulls=residuals/y_err
    chi2_fit[i]=np.sum(pulls**2)
    
plt.hist(chi2_fit)
plt.scatter(x_values,exact_y_values,color='k',ls='--',
            label='true model')
plt.scatter(x_values,y_values_fit,color='r')
plt.errorbar(x_values,y_values_fit,yerr=y_err)

Что следует изменить, чтобы получить разумный сюжет.что-то вроде этого И гистограмма должна быть такой


python numpy data-analysis chi-squared
person ZNNN    schedule 07.02.2021    source источник

Ответы (1)


arrow_upward
0
arrow_downward

Проблема заключается в определении вашей функции и использовании i:

def func (x_values,offset,amp,mean,wid):
    return (offset + amp*np.exp(-((i-mean)/wid)**2))

for i in range (1000):
    ...
    [popt,pcov]=opt.curve_fit(func,x_values, fake_exp_y_values,p01)
    ...

Исправленная версия:

def func (x_values,offset,amp,mean,wid):
    return (offset + amp*np.exp(-((x_values-mean)/wid)**2))

Также, пожалуйста, в следующий раз опубликуйте MRE. Этот не был минимальным и не воспроизводимым (мне пришлось добавить y_values_fit = y_values_fit.flatten() в конце, чтобы заставить его работать).

person rpoleski    schedule 07.02.2021
comment
О, эту функцию я определил для предыдущего вопроса, и все должно быть в порядке. Но я добавил flatted() в конце y_values_fit, и это сработало. Я немного смущен тем, как это решило проблему? Можете ли вы дать мне некоторые дополнительные объяснения, чтобы я мог лучше понять? В любом случае это решило мой вопрос, и я очень ценю вашу помощь! Спасибо - person ZNNN; 08.02.2021
comment
В опубликованной вами версии был аргумент x_values, но он не использовался. Вместо этого вы использовали переменную i, которая использовалась для других целей. - person rpoleski; 08.02.2021