Anwendung von Lmfit

Wenn du dir nicht sicher bist, in welchem der anderen Foren du die Frage stellen sollst, dann bist du hier im Forum für allgemeine Fragen sicher richtig.
Antworten
Micha_uni
User
Beiträge: 15
Registriert: Dienstag 23. Mai 2017, 12:14

Hallo!
Ich habe versucht gemäß
https://lmfit.github.io/lmfit-py/model.html
eine Kurve zu fitten. Es klappt leider nicht.
Ich erhalte folgende Fehlermeldung:
gmodel = Model(gaussian)
NameError: name 'gaussian' is not defined
Obwohl Gauß im lmfit drin sein müsste.

Code: Alles auswählen

import numpy as np
import random
import matplotlib.pyplot as plt
from lmfit import  Model

eqstep = 500 ; mcstep = 500 ; gitter1 = 10
T = np.linspace(1.,3.5,50)

def initialzustand(zeile,spalte):
    return np.random.choice([1, -1], size=(zeile, spalte))

def nachbar_spin(konfig,gitter,x,y):
    links = (x, y-1)
    rechts = (x, (y+1) % gitter)
    oben  = (x-1, y)
    unten = ((x+1) % gitter, y)
    return[konfig[links[0], links[1]],
           konfig[rechts[0], rechts[1]],
           konfig[oben[0], oben[1]],
           konfig[unten[0],unten[1]]]

def energie(konfig,gitter,x,y):
    return 2 * konfig[x,y] * sum(nachbar_spin(konfig,gitter,x,y))

def metropolis():
    magnet = np.array([]) ; magnet2 = np.array([]) ; magnet4 = np.array([]) ;
    for temp in T:
        konfig = initialzustand(gitter1,gitter1)
        mag    = np.zeros(eqstep + mcstep)
        for sweep in range(eqstep + mcstep):
            for i in range(gitter1):
                for j in range(gitter1):
                    e = energie(konfig,gitter1,i,j)
                    if e <= 0:
                        konfig[i,j] *= -1
                    elif np.exp((-1.0 * e) / temp) > random.random():
                        konfig[i,j] *= -1
            mag[sweep] = (sum(sum(konfig))) / gitter1**2
        magnet   =   np.append(magnet,((sum(mag[eqstep:]) / mcstep)))
        magnet2  =   np.append(magnet2,sum(abs(mag[eqstep:])) / mcstep)
        magnet4  =   np.append(magnet4,sum(mag[eqstep:]**2) / mcstep)
        EW       =   magnet4 / magnet2**2
        Sus      =   magnet4 - magnet **2

    return Sus

EW = metropolis()

gmodel = Model(gaussian)
result = gmodel.fit(EW,T=T)
plt.plot(T,EW,'bo')
plt.plot(result.best_fit,'r--')
plt.show()
Weiß jemand, wo mein Fehler ist.
Danke euch!
LG
Micha
__deets__
User
Beiträge: 14529
Registriert: Mittwoch 14. Oktober 2015, 14:29

Die Fehlermeldung ist doch sehr klar. Es gibt den Namen gaussian nicht. Wenn ich auf deinen Link klicke, dann sehe ich aber, das da eine Funktion mit dem Namen existiert. Die wirst du auch anlegen muessen.
Micha_uni
User
Beiträge: 15
Registriert: Dienstag 23. Mai 2017, 12:14

Ich dachte es geht um Gauß.
Also wie kann ich dann meinen Code mit lmfit fitten.
Parameter habe ich nicht.
EW sind lediglich Punkte, die mir die Funktion metropolis rausgibt. Nun möchte ich diese Punkte gegen T fitten. Weiß aber nicht, wie es geht.
Weiß du, wie ich da vorgehen sollte?
LG
Micha
__deets__
User
Beiträge: 14529
Registriert: Mittwoch 14. Oktober 2015, 14:29

Nochmal: schau dir deinen Link an, und such nach dem Vorkommen von "gaussian". Das steht sogar im ersten Stueck Code den die ueberhaupt zeigen.
Micha_uni
User
Beiträge: 15
Registriert: Dienstag 23. Mai 2017, 12:14

Das ist ein Beispiel um die Abweichungen der Werte zu demonstrieren.
Wenn ich anstatt gaussian meine Funktion rein tue, kriege ich wieder Fehlermeldungen.
Bei gmodel = Model(EW) oder gmodel = Model(metropolis()) erhalte ich: TypeError: unsupported callable
Is auch klar, weil die EWs Werte der Funktion sind und nicht die Funktion selbst.
Bei gmodel=Model(metropolis) erhalte ich: IndexError: list index out of range
__deets__
User
Beiträge: 14529
Registriert: Mittwoch 14. Oktober 2015, 14:29

Die Funktion die du an Model übergibst muss eine Verteilung sein. So wie ich dein metropolis verstehe ist das eine Simulation aus der Daten fallen, die einer Verteilung entsprechen. Was ist das für eine? Gauß? Oder welche sonst?
Micha_uni
User
Beiträge: 15
Registriert: Dienstag 23. Mai 2017, 12:14

Die Funktion gibt mir ein paar verschiedene Sachen raus.
Sus gegen T, was ich gerade plotten möchte ist die Suszeptibilität. Der Graph ist ähnlich wie eine Gaußglocke nur nach links verzerrt.
magnet gegen T ist die Magnetisierung. Sie ist bis T=2.26 fast konstant 1 oder -1 dann fallen beide Kurven von -1 und 1 zu Null zusammen.
EW ist die Erwartungswert <M^2>/<|M|>^2. Er ist bis T=2.26 fast konstant 1 dann steigt schnell.

wenn mann die Graphen als isolierte Punkte, Also Z.B. " marker='o' " plotet, sieht man die Verläufe.
Das reicht mir aber nicht. Ich muss sie fitten. Weiß aber nicht, wie.
Bin blutiger Anfänger im Programmieren.

LG
Micha
__deets__
User
Beiträge: 14529
Registriert: Mittwoch 14. Oktober 2015, 14:29

Das fitting so wie ich es verstehe benoetigt eine zu fittende Verteilung (im Beispiel gaussian), und danach fittet es uebergebene Daten.

Code: Alles auswählen

model = Model(verteilung)
model.fit(daten, x=xwerte)
Was du aber beschreibst ist, das du versuchst

Code: Alles auswählen

model = Model(daten)
auszufuehren. Das kann ja nirgendwo hinfuehren, denn deine Datenpunkte sind ja keine Verteilungsfunktion.

Und du kannst natuerlich auch nur eine deiner Datensaetze fitten - die verschiedenen Groessen, die du beschreibst, gehorchen ja jeweils anderen Verteilungen. Du musst also pro Groesse ein Modell mit einer entsprechenden Verteilung fitten.
Micha_uni
User
Beiträge: 15
Registriert: Dienstag 23. Mai 2017, 12:14

model = Model(verteilung)
model.fit(daten, x=xwerte)

Für die Suszeptibilität hätte ich dann: daten= sus_aufgenommene Punkte, x = T
Woher weiß ich, welche Verteilung in der Aufnahme Herrscht?

Danke dir...
LG
__deets__
User
Beiträge: 14529
Registriert: Mittwoch 14. Oktober 2015, 14:29

Keine Ahnung. Anschauen? Eine Gauss-Glocke erkennt man ja ganz gut, und andere beliebte Klassiker wie Exponentialverteilung kann man ja auch sehen.

Ich glaube nicht, dass man das gross automatisiert bekommt, bzw. bestenfalls ueber statistische Tests ala chi-quadrat, nachdem man das mal mit mehreren ausprobiert hat. Da bin ich aber auch wirklich kein Experte.
Antworten