[PYTHON] Analyse der Messdaten (2) -Hydrobacter und Anpassung, lmfit Empfehlung-

Vorwort

Eines Tages Treffen, Kyo "T-Kun, kannst du diese Grafik nicht ein bisschen besser anpassen?" "Können wir den Anpassungsbereich nicht ändern oder ihn beschweren?" T "... ich werde es versuchen"

Er möchte die Armatur modifizieren. Es ist oft erforderlich, das Analyseverfahren aus einer Laune heraus zu ändern, und es ist schwierig, damit umzugehen. Ich möchte es vor Ort tun können.

Dieses Mal werde ich anstelle von scipy.optimize.curve_fit eine neue Bibliothek namens lmfit einführen. Die Anpassung unter Verwendung von Scipy wurde in Analyse der Messdaten ①-Memorandum of Scipy Fitting- durchgeführt. Die Codierung erfolgt auf Jupyter & Jupyter empfohlen.

Ich hoffe ich kann lmfit gut vorstellen.

GithHub hat auch Notebook- und Beispieldaten. Von hier

Über lmfit

lmfit wurde als aufwärtskompatibel mit scipy.optimize entwickelt. Sie können Parameter und Gewichtsfehler beim Anpassen einschränken. Es scheint so, als würde man dem Objekt verschiedene Informationen hinzufügen und schließlich den Befehl fit ausführen, um es zu optimieren. Wenn Sie sich daran gewöhnen, ist es bequemer als scipy.

Dieses Experiment

Pulsdaten wurden vom Messgerät erhalten. Erstellen Sie zunächst ein Histogramm mit dem Spitzenwert des Impulses. Danach wird eine Anpassung durchgeführt und die Schärfe des Spektrums bewertet. Darüber hinaus werden wir versuchen, die Schärfe mithilfe der lmfit-Gewichtungsoption zu verbessern.

Daten lesen

Lesen Sie die vom Messgerät mit Pandas überstrichenen Pulsdaten ab. Verwenden Sie erneut tkinter, um die Dateipfade zu verbinden.

.ipynb

import tkinter 
from tkinter import filedialog
import numpy as np
import pandas as pd

root = tkinter.Tk()
root.withdraw()
file = filedialog.askopenfilename(
title='file path_',
filetypes=[("csv","csv")])
if file:
    pulse = pd.read_csv(file,engine='python',header=0,index_col=0)

pulse.iloc[:,[0,1,2]].plot()

Zeichnen wir nur 3 Daten. Es wurde herausgefegt, so dass der maximale Spitzenwert nahe dem 50. Punkt erhalten werden konnte.

Histogramm erstellen

Erstellen Sie ein Histogramm mit dem maximalen Spitzenwert. Beim Erstellen eines Histogramms müssen Sie die Schrittgröße auf der horizontalen Achse bestimmen. Dieses Mal wird gemäß der Sturges-Formel die Anzahl der Fächer auf k = 1 + Log2 (n) eingestellt.

* Np.histogram gibt ** width ** so viele wie die Anzahl zurück.

.ipynb

bins = int(1+np.log2(pulse.shape[1]))#Sturges offiziell
count, level = np.histogram(pulse.max(),bins=bins)
level = np.array(list(map(lambda x:(level[x]+level[x+1])/2,range(level.size-1))))
#Frequenz und Größe ausrichten
level = np.array(list(map(lambda x:(level[x]+level[x+1])/2,range(level.size-1))))
import matplotlib.pyplot as plt
plt.plot(level,count,marker='.',linestyle='None')
plt.show()

Auf der rechten Seite war ein kleiner Berg. Lassen Sie uns die Daten nach 50 trennen und erneut eine Verteilung für den Berg auf der rechten Seite vornehmen.

.ipynb

mask = pulse.max()>50#Erstellen Sie ein Array von Bool-Werten
pulse = pulse.loc[:,mask]#Pulsdaten maskieren

bins = int(1+np.log2(pulse.shape[1]))#Sturges offiziell
count, level = np.histogram(pulse.max(),bins=bins)
level = np.array(list(map(lambda x:(level[x]+level[x+1])/2,range(level.size-1))))
plt.plot(level,count,marker='o',linestyle='None')
plt.show()

Es ist so etwas wie eine Normalverteilung aufgetreten.

Passend zu lmfit

Wie Sie aus import ~ sehen können, verwendet lmfit einige Klassenobjekte zum Anpassen. Die Anpassungsfunktion ist eine Gaußsche Verteilung und lautet wie folgt. Das Zentrum des Peaks wird durch $ \ mu $ angegeben, und die Variation der Verteilung wird durch $ \ sigma $ angegeben. amp ist der Amplitudenbegriff. $f_{\left(x\right)}=\frac{amp}{\sqrt{2\pi\sigma^2}} \exp\left(-\frac{\left(x-\mu \right)^2}{2\sigma^2}\right)$ Geben Sie zunächst eine Funktion an, die zur Modellklasse passt. Erstellen Sie als Nächstes eine Parameterklasse und geben Sie ihr Anfangswerte. Die Model-Klasse verfügt über eine .make_params () -Methode, die nach dem zweiten Argument der angegebenen Funktion automatisch ein Parameters-Objekt erstellt.

.ipynb

from lmfit import Model, Parameters, Parameter
def gaussian(x,amp,mu,sigma):
    out = amp / (sigma * np.sqrt(2 * np.pi)) * np.exp(-(x - mu)**2 /(2 * sigma**2))
    return out
model = Model(gaussian)#Anpassungsfunktion zuweisen, Modellobjekt erstellen
pars = model.make_params()#Parameter Objekterzeugung
pars

Das erstellte Parameterobjekt besteht aus einer Reihe von Parameterobjekten. Auf jeden Parameter kann wie auf ein Wörterbuch zugegriffen werden. Von dort mit der .set-Methode ・ Anfangswert ・ Ober- und Untergrenze der Parameter, max, min: -np.inf ~ np.inf

• Ob der Wert vor und nach dem Anpassen geändert werden soll (Einschränkung), variieren = Richtig oder Falsch Sie können mit jedem von ihnen spielen.

.ipynb

pars['mu'].set(value=70)
pars['amp'].set(value=1)
pars['sigma'].set(value=1)
pars

Es ist sehr rau, aber ich habe den Anfangswert angegeben. Die Parameter können detailliert eingestellt werden und werden in einer DataFrame-ähnlichen Tabelle angezeigt. Es ist cool ...

Montage durchführen

.ipynb

result = model.fit(data=count,params=pars,x=level)
result

Das Ergebnis wird folgendermaßen angezeigt. Es scheint, dass ein neues ModelResult-Objekt erstellt wird, wenn .fit ausgeführt wird. Das generierte Objekt hatte eine .plot_fit () -Methode. Visualisieren Sie mit dieser Funktion.

Optimierte Parameter können mithilfe der Methode .best_values in einem Wörterbuch abgerufen werden.

.ipynb

opt = result.best_values
opt['mu']/opt['sigma']

Die Schärfe wurde mit $ \ frac {\ mu} {\ sigma} $ bewertet. Es war 23,67.

Ich habe versucht, die Schrittgröße auf der horizontalen Achse zu ändern

・ Scotts Beamter $h = \frac{3.49\sigma}{ n^\frac{1}{3} }$ Ich habe versucht, die Schrittgröße mit zu definieren. Da es eine große Sache ist, werde ich die verwenden, die ich in der Anpassung für $ \ sigma $ gefunden habe.

.ipynb

h = 3.49 * opt['sigma'] / (pulse.shape[1] **(1/3))#Scotts Beamter
r = pulse.max().max()-pulse.max().min()
bins = int(r/h)

count, level = np.histogram(pulse.max(),bins=bins)
level = np.array(list(map(lambda x:(level[x]+level[x+1])/2,range(level.size-1))))
plt.plot(level,count,marker='o',linestyle='None')
plt.show()

Die Schrittgröße ist enger geworden und der Maximalwert der Fraktion ist kleiner geworden, aber die Form der Verteilung scheint in Ordnung zu sein. Passen Sie an, indem Sie nur die Werte der zu ersetzenden Daten und x austauschen. Die Schärfe beträgt ... 24,44. Verbessert um 0,77.

Ich habe versucht, den Fehler zu gewichten

Dies ist ab hier der Fall.

Er wollte den Bereich von $ \ mu \ pm \ sigma $ stärken und Gewichte von weniger als 63 bis 0 anpassen. Es scheint, dass die Gewichtung durch Angabe der Option weight in .fit erfolgen kann.

.ipynb

s = (level>opt['mu']-opt['sigma']) * (level<opt['mu']+opt['sigma'])
weight = (np.ones(level.size) + 4*s)*(level>63)
result3 = model.fit(data=count,params=pars,x=level,weights=weight)

↓ Gewicht & Ergebnis Es sieht nicht viel anders aus, aber wenn man genau hinschaut, fühlt es sich an, als würde es ein wenig heruntergezogen. Die Schärfe beträgt ... 25,50. Nur 0,06 wurden schärfer.

Post-Fit-Punkte können mit der Methode .best_fit als Array aufgerufen werden. Es war praktisch, wenn man CSV macht. Lösen Sie das Problem sicher ...

Zusammenfassung

Dieses Mal wurden die Parameter automatisch mit Model.make_params () generiert. Es gibt jedoch auch eine Methode, um eine leere Klasse mit Parameters () direkt zu generieren und mit der Methode .add Variablen hinzuzufügen.

Ich war zuerst verwirrt, weil ich mehrere Objekte zum Anpassen verwendet habe. Es war leicht zu verwechseln, wo und welche Methode entspricht. Wenn Sie sich jedoch daran gewöhnen, werden Sie Parameter schätzen.

Die Tatsache, dass es nur wenige Websites gab, die lmfit auf Japanisch erklärten, war auch der Grund, warum ich anfing, diesen Artikel zu schreiben. Gibt es eine bessere Bibliothek ...? Ich würde es begrüßen, wenn Sie mich unterrichten könnten. Junioren werden in Zukunft keine Probleme mehr haben. Es war mein erstes Mal, dass ich beim Lesen der offiziellen Referenz Code schrieb (lacht). Aufgrund der Forschungsergebnisse bleibt keine andere Wahl, als dies zu tun.

In Zukunft möchte ich etwas mehr Technologie lernen, die mit meinem Leben verbunden zu sein scheint, und ich möchte auch in der Lage sein, die Technologie, die ich in verschiedenen Bereichen habe, zu nutzen. Ich möchte kreativ sein.

Das ist alles für Python-bezogene Fortschritte. Ich konnte nicht anders, als es dem Professor zu melden (niemand macht Python), deshalb möchte ich diese Gelegenheit nutzen, um es zu melden. Ich bin dem Professor zutiefst dankbar, dass er geduldig auf meine Datenanalyse gewartet hat.