2017-10-17 1 views
2

으로 가우스 맞추기를 제한하려면 파이썬으로 데이터를 평가해야합니다. 불행히도 아직 내 동료 학생들의 대본이 없기 때문에 프로그래밍에 익숙하지 않습니다.내 학사 논문의 프레임 워크에서 curve_fit

이 데이터 세트가 있고 scipy.optimize.curve_fit를 사용하여 가우스로 맞추려고합니다. 특히 축의 끝에서 사용할 수없는 카운트가 많이 있기 때문에, 나는 맞추어야 할 부분을 제한하고 싶습니다. 내가 bounds=([2400,-np.inf, -np.inf],[2600, np.inf, np.inf])를 사용하려고하면 docs.scipy.org에

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

x=np.arange(5120) 
y=array([ 0.81434599, 1.17054264, 0.85279188, ..., 1.  , 
    1.  , 13.56291391]) #most of the data isn't interesting 
#to me, part of interest see below 

def Gauss(x, a, x0, sigma): 
    return a * np.exp(-(x - x0)**2/(2 * sigma**2)) 

mean = sum(x * y)/sum(y) 
sigma = np.sqrt(sum(y * (x - mean)**2)/sum(y)) 

popt,pcov = curve_fit(Gauss, x, y, p0=[max(y), mean, sigma], 
maxfev=360000) 

plt.plot(x,y,label='data') 
plt.plot(x,Gauss(x, *popt), 'r-',label='fit') 

내가 curve_fit

에 대한 일반적인 설명을 발견했습니다

사진은 raw data

이것은 내가 지금까지 무엇을 가지고 , ValueError가 표시됩니다. x0은 실행할 수 없습니다. 여기에 어떤 문제가 있습니까? 이 질문에 단 댓글에서 제안

는 또한 popt,pcov = curve_fit(Gauss, x[2400:2600], y[2400:2600], p0=[max(y), mean, sigma], maxfev=360000) 으로 범위를 제한하려고 : "오류를 그래프 가우스 적합을받을 때"나는 단지 비록 직선을 얻을이 경우에 유래 에서.

사진 : Confinement with x[2400:2600],y[2400:2600] as arguments of curve_fit

난 정말 당신이 나를 도울 수 있기를 바랍니다. 난 단지 내 데이터의 작은 부분에 맞는 방법이 필요합니다. 미리 감사드립니다!

흥미로운 데이터 :

y=array([ 0.93396226, 1.00884956, 1.15457413, 1.07590759, 
0.88915094, 1.07142857, 1.10714286, 1.14171123, 1.06666667, 
0.84975369, 0.95480226, 0.99388379, 1.01675978, 0.83967391, 
0.9771987 , 1.02402402, 1.04531722, 1.07492795, 0.97135417, 
0.99714286, 1.0248139 , 1.26223776, 1.1533101 , 0.99099099, 
1.18867925, 1.15772871, 0.95076923, 1.03313253, 1.02278481, 
0.93265993, 1.06705539, 1.00265252, 1.02023121, 0.92076503, 
0.99728997, 1.03353659, 1.15116279, 1.04336043, 0.95076923, 
1.05515588, 0.92571429, 0.93448276, 1.02702703, 0.90056818, 
0.96068796, 1.08493151, 1.13584906, 1.1212938 , 1.0739645 , 
0.98972603, 0.94594595, 1.07913669, 0.98425197, 0.87762238, 
0.96811594, 1.02710843, 0.99392097, 0.91384615, 1.09809264, 
1.00630915, 0.93175074, 0.87572254, 1.00651466, 0.78772379, 
1.12244898, 1.2248062 , 0.97109827, 0.94607843, 0.97900262, 
0.97527473, 1.01212121, 1.16422287, 1.20634921, 0.97275204, 
1.01090909, 0.99404762, 1.00561798, 1.01146132, 1.08695652, 
0.97214485, 1.03525641, 0.99096386, 1.05135952, 1.16451613, 
0.90462428, 0.76876877, 0.47701149, 0.27607362, 0.21580547, 
0.20598007, 0.16766467, 0.15533981, 0.19745223, 0.15407855, 
0.18925831, 0.26997245, 0.47603834, 0.596875 , 0.85126582, 0.96  
, 1.06578947, 1.08761329, 0.89548023, 0.99705882, 1.07142857, 
0.95677233, 0.86119874, 1.02857143, 0.98250729, 0.94214876, 
1.04166667, 0.96024465, 1.07022472, 1.10344828, 1.04859335, 
0.96655518, 1.06424581, 1.01754386, 1.03492063, 1.18627451, 
0.91036415, 1.03355705, 1.09116809, 0.96083551, 1.01298701, 
1.03691275, 1.02923977, 1.11612903, 1.01457726, 1.06285714, 
0.98186528, 1.16470588, 0.86645963, 1.07317073, 1.09615385, 
1.21192053, 0.94385027, 0.94244604, 0.88390501, 0.95718654, 
0.9691358 , 1.01729107, 1.01119403, 1.20350877, 1.12890625, 
1.06940063, 0.90410959, 1.14662757, 0.97093023, 1.03021148, 
1.10629921, 0.97118156, 1.10693642, 1.07917889, 0.9484127 , 
1.07581227, 0.98006645, 0.98986486, 0.90066225, 0.90066225, 
0.86779661, 0.86779661, 0.96996997, 1.01438849, 0.91186441, 
0.91290323, 1.03745318, 1.0615942 , 0.97202797, 1.16608997, 
0.94182825, 1.08333333, 0.9076087 , 1.18181818, 1.20618557, 
1.01273885, 0.93606138, 0.87457627, 0.90575916, 1.09756098, 
0.99115044, 1.13380282, 1.04333333, 1.04026846, 1.0297619 , 
1.04334365, 1.03395062, 0.92553191, 0.98198198, 1.  , 
0.9439528 , 1.02684564, 1.1372549 , 0.96676737, 0.99649123, 
1.07051282, 1.10367893, 1.0866426 , 1.15384615, 0.99667774]) 
+1

시도가 추가을 추가 일정한. 데이터 의미가 0이 아니기 때문입니다. 또한 범위를 더 좁히지 않으려 고합니다. 약 2500에 해당 피크를 맞추려고합니다. – joaquinn

+0

x0은 첫 번째 경계 인덱스에 대해 유효하지 않습니다. 그리고 오류에 대해 이야기 할 때 전체 스택 추적을 추가하십시오. – sascha

답변

1

당신은 이것에 대한 유용한 lmfit 모듈 (https://lmfit.github.io/lmfit-py/를) 찾을 수 있습니다. 커브 피팅을 매우 쉽게 수행 할 수 있도록 설계되었으며 Gaussian과 같은 공통 피크에 대한 내장 모델을 포함하며 매개 변수 경계를 설정하는 것과 같은 많은 유용한 기능을 제공합니다. lmfit와 데이터에 맞추기는 다음과 같습니다

import numpy as np 
import matplotlib.pyplot as plt 

from lmfit.models import GaussianModel, ConstantModel 

y = np.array([.....]) # uses your shorter data range 
x = np.arange(len(y)) 

# make a model that is a Gaussian + a constant: 
model = GaussianModel(prefix='peak_') + ConstantModel() 

# make parameters with starting values: 
params = model.make_params(c=1.0, peak_center=90, 
          peak_sigma=5, peak_amplitude=-5) 

# it's not really needed for this data, but you can put bounds on 
# parameters like this (or set .vary=False to fix a parameter) 
params['peak_sigma'].min = 0   # sigma > 0 
params['peak_amplitude'].max = 0  # amplitude < 0 
params['peak_center'].min = 80 
params['peak_center'].max = 100 

# run fit 
result = model.fit(y, params, x=x) 

# print, plot results 
print(result.fit_report()) 
plt.plot(x, y) 
plt.plot(x, result.best_fit) 
plt.show() 

[[Model]] 
    (Model(gaussian, prefix='peak_') + Model(constant)) 
[[Fit Statistics]] 
    # function evals = 54 
    # data points  = 200 
    # variables  = 4 
    chi-square   = 1.616 
    reduced chi-square = 0.008 
    Akaike info crit = -955.625 
    Bayesian info crit = -942.432 
[[Variables]] 
    peak_sigma:  4.03660814 +/- 0.204240 (5.06%) (init= 5) 
    peak_center:  91.2246614 +/- 0.200267 (0.22%) (init= 90) 
    peak_amplitude: -9.79111362 +/- 0.445273 (4.55%) (init=-5) 
    c:    1.02138228 +/- 0.006796 (0.67%) (init= 1) 
    peak_fwhm:  9.50548558 +/- 0.480950 (5.06%) == '2.3548200*peak_sigma' 
    peak_height:  -0.96766623 +/- 0.041854 (4.33%) == '0.3989423*peak_amplitude/max(1.e-15, peak_sigma)' 
[[Correlations]] (unreported correlations are < 0.100) 
    C(peak_sigma, peak_amplitude) = -0.599 
    C(peak_amplitude, c)   = -0.328 
    C(peak_sigma, c)    = 0.196 

를 인쇄이 같은 플롯 할 것입니다 :

enter image description here

+0

피팅 결과 중 하나 (이 경우 peak_center)를 새 배열로 인쇄 할 수 있습니까? 'np.append (graph, peak_center) '와 같은 것? – Kike

+1

예, 가장 적합한 매개 변수는 Parameter 객체의 명칭 인'result.params'에 보관됩니다. 'result.params [ 'peak_center']. value'가 필요합니다. 각 매개 변수에는 더 많은 속성이 있습니다. 예를 들어 불확실성은'result.params [ 'peak_center']. stderr'에 있습니다. –