introduction

Cet article est l'article du 21e jour du Calendrier de l'Avent 2019 de ACCESS Co., Ltd. est.

Dans un certain cas de notre entreprise, nous devrions faire un ajustement carré minimum non linéaire, donc pour la préparation, nous examinerons le raccord carré minimum non linéaire par lmfit que nous avons utilisé auparavant et le résumerons ici. ..

lmfit

lmfit est un minimum de deux non-linéaires comme le sous-titre officiel dit «Minimisation non linéaire des moindres carrés et ajustement de courbe pour Python». Une bibliothèque pour l'ajustement de modèle utilisant la multiplication, étendue basée sur de nombreuses méthodes d'optimisation dans scipy.optimize. , A été développé.

traits

Les fonctionnalités suivantes sont répertoriées.

Paramètre Introduction des objets de classe. Cela facilite la gestion des paramètres du modèle.
Facilité de changer l'algorithme d'ajustement du modèle. --Amélioration de l'estimation de l'intervalle de confiance des paramètres. Comparé à celui de scipy.optimize.leastsq, l'intervalle de confiance peut être estimé de manière plus concise. ing.
Model L'ajustement de la courbe est amélioré en utilisant l'objet de classe. Cet objet de classe Model étend les fonctionnalités de scipy.optimize.curve_fit.
De nombreux [modèles intégrés] linéaires courants (https://lmfit.github.io/lmfit-py/builtin_models.html) sont disponibles et prêts à l'emploi.

Où je le trouve utile (subjectivement)

Facilité de manipulation des modèles et des paramètres de modèle

Vous pouvez facilement introduire la fonction que vous souhaitez utiliser dans l'ajustement de modèle que vous avez défini.
L'accès aux paramètres du modèle est facile avec la fonction de dictino de Python.
Les limites de paramètres / réglages fixes peuvent être effectués facilement. Il est facile de définir la plage pour rechercher la valeur optimale du paramètre lors de l'ajustement et de fixer la valeur du paramètre (l'éliminer en tant que paramètre libre).
Les paramètres peuvent être modifiés sans changer la fonction objectif.
Les paramètres peuvent être contraints algébriquement. La relation avec d'autres paramètres peut être restreinte.

Montage du modèle avec lmfit

Les données ont été créées pour cet article et sont tracées comme suit:

Comme vous pouvez le voir en un coup d'œil, ces données ont deux structures, et autour de 6,0 et 7,5, une structure comme une distribution normale peut être vue. Adaptez ces données.

modèle

Comme mentionné ci-dessus, il existe deux structures de type distribution normale, il est donc nécessaire d'introduire un modèle qui les représente. De plus, comme il semble qu'il existe un décalage constant non nul, il est nécessaire d'incorporer un composant qui le reproduit dans le modèle. La formule est la suivante.

f\left( x \right) = G_1 \left( x \right) +  G_2 \left( x \right) + C

Cela peut être exprimé en utilisant l'objet de classe lmfit.models comme suit.

import lmfit as lf
import lmfit.models as lfm

#Définition du modèle
model = lfm.ConstantModel() +  lfm.GaussianModel(prefix='gauss1_') + lfm.GaussianModel(prefix='gauss2_')

En spécifiant prefix, c'est pratique car les paramètres peuvent être distingués même si la même fonction de modèle est utilisée.

Lors de la présentation de votre propre modèle

Il peut être introduit d'un seul coup en passant la fonction définie dans l'objet de classe Model.


def constant(x, const):
    return const

model = lf.Model(constant) +  lfm.GaussianModel(prefix='gauss1_') + lfm.GaussianModel(prefix='gauss2_')

Paramètres

Même si vous déclarez un objet de classe Model, l'objet de classe Parameters pour cet objet de classe Model n'est pas déclaré, alors commencez par cette déclaration.

#Définissez la valeur initiale de chaque paramètre
par_val = {
    'c' : 20,
    'const' : 20,
    'gauss1_center' : 6.0,
    'gauss1_sigma' : 1.0,
    'gauss1_amplitude' : 400,
    'gauss2_center' : 7.5,
    'gauss2_sigma' : 1.0,
    'gauss2_amplitude' : 150
}

par_min = { 
    'c' : 0,
    'const' : 0,
    'gauss1_center' : 0,
    'gauss1_sigma' : 0,
    'gauss1_amplitude' : 0,
    'gauss2_center' : 0,
    'gauss2_sigma' : 0,
    'gauss2_amplitude' : 0
}

par_max = { 
    'c' : 100,
    'const' : 100,
    'gauss1_center' : 10,
    'gauss1_sigma' : 10,
    'gauss1_amplitude' : 1000,
    'gauss2_center' : 10,
    'gauss2_sigma' : 10,
    'gauss2_amplitude' : 1000
}

par_vary = { 
    'c' : True,
    'const' : True,
    'gauss1_center' : True,
    'gauss1_sigma' : True,
    'gauss1_amplitude' : True,
    'gauss2_center' : True,
    'gauss2_sigma' : True,
    'gauss2_amplitude' : True
}

#Objet de classe Déclaration de paramètres pour le modèle défini
params = model.make_params()
for name in model.param_names: 
    params[name].set(
        value=par_val[name], #valeur initiale
        min=par_min[name], #limite inférieure
        max=par_max[name], #limite supérieure
        vary=par_vary[name] #Déplacer le paramètre
    )

La relation entre les paramètres peut être contrainte algébriquement comme suit.

params['gauss2_center'].set(expr='gauss1_center*1.25')

raccord

result = model.fit(x=df.x, data=df.y, weights=df.dy**(-1.0), params=params, method='leastsq')

Remarque: Les données sont données dans l'objet DataFrame de pandas.

Confirmation du résultat de l'ajustement

summary
```
print(result.fit_report())
```

Vous pouvez voir le résumé des résultats de l'ajustement comme suit.

```
[[Model]]
    ((Model(constant) + Model(gaussian, prefix='gauss1_')) + Model(gaussian, prefix='gauss2_'))
[[Fit Statistics]]
    # fitting method   = leastsq
    # function evals   = 383
    # data points      = 50
    # variables        = 6
    chi-square         = 29.2390839
    reduced chi-square = 0.66452463
    Akaike info crit   = -14.8258350
    Bayesian info crit = -3.35369698
[[Variables]]
    const:             21.4200227 +/- 1.47331527 (6.88%) (init = 20)
    gauss1_amplitude:  88.4566498 +/- 1.39004141 (1.57%) (init = 400)
    gauss1_center:     5.99651706 +/- 0.00136123 (0.02%) (init = 6)
    gauss1_sigma:      0.09691089 +/- 0.00153157 (1.58%) (init = 1)
    gauss2_amplitude:  31.6555368 +/- 1.39743186 (4.41%) (init = 150)
    gauss2_center:     7.49564632 +/- 0.00170154 (0.02%) == 'gauss1_center*1.25'
    gauss2_sigma:      0.09911853 +/- 0.00446584 (4.51%) (init = 1)
    gauss1_fwhm:       0.22820770 +/- 0.00360656 (1.58%) == '2.3548200*gauss1_sigma'
    gauss1_height:     364.139673 +/- 4.89553113 (1.34%) == '0.3989423*gauss1_amplitude/max(2.220446049250313e-16, gauss1_sigma)'
    gauss2_fwhm:       0.23340629 +/- 0.01051624 (4.51%) == '2.3548200*gauss2_sigma'
    gauss2_height:     127.410415 +/- 4.77590041 (3.75%) == '0.3989423*gauss2_amplitude/max(2.220446049250313e-16, gauss2_sigma)'
[[Correlations]](unreported correlations are < 0.100)
    C(gauss2_amplitude, gauss2_sigma)     =  0.647
    C(gauss1_amplitude, gauss1_sigma)     =  0.636
    C(const, gauss2_amplitude)            = -0.555
    C(const, gauss1_amplitude)            = -0.552
    C(const, gauss1_sigma)                = -0.367
    C(const, gauss2_sigma)                = -0.359
    C(gauss1_amplitude, gauss2_amplitude) =  0.306
    C(gauss1_sigma, gauss2_amplitude)     =  0.203
    C(gauss1_amplitude, gauss2_sigma)     =  0.198
    C(gauss1_sigma, gauss2_sigma)         =  0.131
```

confidence interval

print('[[Confidence Intervals]]')
print(result.ci_report())

[[Confidence Intervals]]
                    99.73%    95.45%    68.27%    _BEST_    68.27%    95.45%    99.73%
const           :  -4.72534  -3.05035  -1.49522  21.42002  +1.48820  +3.02111  +4.65531
gauss1_amplitude:  -4.37979  -2.84434  -1.40305  88.45665  +1.41243  +2.88513  +4.47133
gauss1_center   :  -0.00435  -0.00278  -0.00131   5.99652  +0.00131  +0.00278  +0.00436
gauss1_sigma    :  -0.00479  -0.00312  -0.00154   0.09691  +0.00156  +0.00321  +0.00499
gauss2_amplitude:  -4.31612  -2.82491  -1.40294  31.65554  +1.43259  +2.94767  +4.61057
gauss2_sigma    :  -0.01344  -0.00889  -0.00446   0.09912  +0.00466  +0.00971  +0.01540

--Tracé des résultats de l'ajustement

```python
fig, gridspec = result.plot(data_kws={'markersize': 5})
```

Le graphique suivant est généré automatiquement.

-Résultat d'ajustement Les propriétés suivantes ont la valeur optimale, etc., alors n'hésitez pas à faire le reste.

```python
result.chisqr #Valeur du chi carré
result.nfree  #Degré de liberté
result.redchi #Valeur chi carré convertie
result.aic    # aic
result.bic    # bic    
result.best_fit  # best-ajustement de la valeur du modèle(Valeur du graphique orange du résultat de l'ajustement)
result.residual  # (best-fit model)-(data)La valeur du(Valeurs dans le panneau supérieur des résultats d'ajustement)
result.eval_components(x=df.x)  # best-Valeur pour chaque composant du modèle d'ajustement
result.best_values #Valeur optimale pour chaque paramètre:type de dict
result.result.params['const'].value  #Valeur optimale du paramètre
result.result.params['const'].stderr #Erreur standard des paramètres
```

c'est tout. Veuillez continuer à profiter de l'article de @ rheza_h le 22e jour.

[PYTHON] Montage du modèle avec lmfit