Различные ограничения для параметра соответствия в модели lmfit - PullRequest
Новая
       27

Различные ограничения для параметра соответствия в модели lmfit

0 голосов
/ 08 мая 2018

Я пытаюсь с помощью lmfit создать функцию подгонки пиковых чисел voigt / Gaussian / Lorentizan-peak.Поэтому я написал следующую функцию: 

def apply_fit_mix_multy(data,modelPeak,peakPos,amplitud,**kwargs):
peakPos=np.array(peakPos)
Start=kwargs.get('Start',data[0,0])
length_data=len(data)-1
End=kwargs.get('End',data[length_data,0])
StartPeak=kwargs.get('StartPeak',data[0,0])
EndPeak=kwargs.get('EndPeak',data[length_data,0])
BackFunc=kwargs.get('BackFunc',False)
BackCut=kwargs.get('BackCut',False)
dataN=data_intervall(data,Start,End)
y=dataN[:, 1]
x=dataN[:, 0]
amplitud=amplitud
center=peakPos

mod = None
for i in range(len(peakPos)):
    this_mod = make_model(i,amplitud,center,modelPeak)
    if mod is None:
        mod = this_mod
    else:
        mod = mod + this_mod

bgy=[list() for f in range(len(x))]
if(BackFunc==True):
    bg,bgx=BackFunc
    for i in range(len(x)):
        bgy[i]=bg.best_values.get('c')        

elif(BackCut!=False):
    slope,intercept=back_ground_cut(data,BackCut[0],BackCut[1])  
    for i in range(len(x)):
        bgy[i]=slope*x[i]+intercept       

if(BackCut!=False):
    print('Background substraction model is used! (Sign=Sign-backgr(linear between two points))')
    y=y-bgy
    out = mod.fit(y, x=x)
else:  
    print('Combination model is used! (offset+Gauss/Lor/Voig)')
    offset=ConstantModel()
    mod=mod+offset

out = mod.fit(y, x=x)#out is the fitted function

area=[list() for f in range(len(peakPos))]
comps=out.eval_components(x=x)
if(BackCut!=False):
    for i in range(len(peakPos)):
        area[i]=simps(comps['peak'+str(i)+'_'],x=x,even='avg')-simps(bgy,x=x,even='avg')
    fit_dict={'signal':y, 'convol':out.best_fit,'x':x,'peak_area':area,'backgr':bgy,'comps':comps}
else:
    for i in range(len(peakPos)):
        area[i]=simps(comps['peak'+str(i)+'_'],x=x,even='avg')
    fit_dict={'convol':out.best_fit,'x':x,'peak_area':area,'comps':comps} #comps is inf. of  sperate peaks

return fit_dict

Функция считывает в наборе данных modelPeak (например, GaussianModel) начальное предположение о положениях и амплитудах пиков (peakPos, амплитуды).

В первой части я инициализирую модель пиков (сколько пиков ...) 

    for i in range(len(peakPos)):
    this_mod = make_model(i,amplitud,center,modelPeak)
    if mod is None:
        mod = this_mod
    else:
        mod = mod + this_mod

С помощью функции make_model: 

def make_model(num,amplitud,center,mod):
pref = "peak{0}_".format(num)
model = mod(prefix = pref)
model.set_param_hint(pref+'amplitud', value=amplitud[num], min=0, max=5*amplitud[num])
model.set_param_hint(pref+'center', value=center[num], min=center[num]-0.5, max=center[num]+0.5)
if(num==0):
    model.set_param_hint(pref+'sigma', value=0.3, min=0.01, max=1)
else:
    model.set_param_hint(pref+'sigma', value=0.3, min=0.01, max=1)
#print('Jetzt',center[num],amplitud[num])
return model

вот теперь моя проблема: IIнапример, чтобы соответствовать 3 пикам. Я хочу, чтобы, например, сигма первого пика изменялась во время подгонки, в то время как сигмы других пиков зависят от сигмы первого пика!любая идея?THX математика

К вашему сведению, вот как выглядит подгонка: введите описание изображения здесь

1 Ответ

0 голосов
/ 08 мая 2018

Если я понимаю ваш длинный вопрос (было бы полезно удалить посторонние вещи - а их довольно много), вы захотите создать Model с несколькими пиками, позволяющими sigma с 1-го пика свободно изменяться и ограничивать sigma, чтобы другие пики зависели от этого.

Для этого вы можете либо использовать подсказки параметров (как вы используете в своей функции make_model()), либо задавать выражения для параметров после создания объекта «Параметры». Для первого подхода, что-то вроде этого 

def make_model(num,amplitud,center,mod):
    pref = "peak{0}_".format(num)
    model = mod(prefix = pref)
    model.set_param_hint(pref+'amplitud', value=amplitud[num], min=0, max=5*amplitud[num])
    model.set_param_hint(pref+'center', value=center[num], min=center[num]-0.5, max=center[num]+0.5)
    if(num==0):
        model.set_param_hint(pref+'sigma', value=0.3, min=0.01, max=1)
    else:
        ## instead of 
        # model.set_param_hint(pref+'sigma', value=0.3, min=0.01, max=1)
        ## set peakN_sigma == peak0_sigma
        model.set_param_hint(pref+'sigma', expr='peak0_sigma') 
        ## or maybe set peakN_sigma == N * peak0_sigma 
        model.set_param_hint(pref+'sigma', expr='%d*peak0_sigma' % num)
    return model

Вы также можете сделать полную модель (несколько упрощенную из вашего кода, но ту же идею): 

model = (VoigtModel(prefix='peak0_') + VoigtModel(prefix='peak1_') +
         VoigtModel(prefix='peak2_') + LinearModel(prefix='const_'))

# create parameters with default values
params = model.make_params(peak0_amplitude=10, peak0_sigma=2, ....)

# set constraints for `sigma` params:
params['peak1_sigma'].expr = 'peak0_sigma'
params['peak2_sigma'].expr = 'peak0_sigma'

# similarly, set bounds as needed:
params['peak1_sigma'].min = 0
params['peak1_amplitude'].min = 0

Надеюсь, это поможет ...

...