Apache Commons Math: Eigene Fitnessfunktion für Curve Fitting?

[Senior Sanchez]

Hallo miteinander,

ich möchte die Apache Commons Math-Bibliothek gerne zum Fitten eines Polynoms des Grad n an eine Messreihe verwenden. Als Ausgangsbasis habe ich dieses Code-Stück unter "12.6 Curve Fitting" verwendet (http://commons.apache.org/proper/commons-math/userguide/optimization.html).

Allerdings nutzt diese Funktion standardmäßig RMS als Fitnessfunktion, minimiert also die Summe der quadratischen Abstände. Allerdings habe ich eine eigene Fitnessfunktion, die ich dafür verwenden möchte. Nur wie kann ich das mit Apache Commons bewerkstelligen? Ich steige da irgendwie nicht durch. Natürlich könnte ich jetzt im Quelltext nachschauen, aber ich möchte ungerne da irgendwelche Modifikationen drin vornehmen nur um meine Fitnessfunktion darein zu quetschen.

Hat jemand eine Idee?

Danke!

[Milchkanne]

Rein aus interesse, was für eine Fittnessfunktion willst du denn verwenden? L1?

[Senior Sanchez]

Rein aus interesse, was für eine Fittnessfunktion willst du denn verwenden? L1?

Nope, das ist wesentlich ausgefeilter. Diese Funktion kombiniert mehrere Verfahren und gewichtet diese entsprechend.

[Milchkanne]

Also die CurveFitter Klasse scheint keine möglichkeit dafür zu bieten:
http://commons.apache.org/proper/commons-math/apidocs/org/apache/commons/math3/fitting/CurveFitter.html
Ich denke du musst das einfach als Optimierungsproblem formulieren und dann mit einem allgemeinen Optimierer lösen, wie in 12.5 General Case beschrieben.
Wenn du konkrete Fragen hast, kann ich dir auch weiterhelfen...
Deine Fittnessfunktion ist vermutlich auch nicht stetig differenzierbar oder? Das macht die Sache immer etwas schwieriger... aber wie gesagt, erklär das Problem dann kann man dir helfen.

[Senior Sanchez]

Das dachte ich mir schon, dass man das als Optimierungsproblem definieren muss, aber so richtig blicke ich da nicht durch, welche Softwarekomponenten zusammenspielen.

Ich gehe nicht davon aus, dass die stetig differenzierbar ist. Ich habe es aber auch nie untersucht, weil sie sehr aufwändig ist.

Das Problem habe ich ja schon erklärt. :) Ich möchte einfach ein Polynom gegen eine Messreihe fitten und dafür meine, eventuell nicht stetig differenzierbare Fitnessfunktion verwenden.

[Gast]

Ist deine Fitnessfunktion wenigstens konvex? Für konvexe nicht glatte Optimierung gibt es schon ein paar Sachen...

[Senior Sanchez]

Das habe ich nicht überprüft. Ich wüsste auch nicht, wie ich das jetzt auf die schnelle überprüfen kann.

[Milchkanne]

Dass sie konvex ist, würde ich ja mal stark annehmen. Zumindest kann ich mir keine sinnvolle nicht-konvexe Fittness vorstellen. Am einfachsten wäre ja dur würdest die hier einfach mal posten, oder hast du angst wir fangen dann an über die Sinnhaftigkeit der Funktion zu diskutieren? ;-)
Nutzlose Nebenfrage: Bist du auf Apache & Java festgelegt?

[Senior Sanchez]

Über die Sinnhaftigkeit werdet ihr nicht anfangen zu diskutieren, keine Bange. ;) Die Funktion hat schon ihre Berechtigung, nur wie gesagt, das ist nicht zwischen Tür und Angel erklärbar. Im wesentlichen bezieht die Funktion neben einer Abstandsmetrik auch solche Sachen wie Phasenverschiebung, Flächenintegrale und eine Korrelationsanalyse mit ein.

Auf Apache bin ich nicht festgelegt, Java würde ich aber gerne behalten, da meine ganze Toolchain darauf aufsetzt.

[Milchkanne]

Aber geschrieben hast du die Routine für die Fittness schon ja? Oder ist das bisher nur nicht diskretisiert auf dem Papier?

Und die Gefahr dass du mit der Diskussion über die Sinnhaftigkeit anfängst bestand wirklich nie... ;-)

[Gast]

Eine einfache Fitnessfunktion die nicht konvex ist wär z.b. eine abgeschnittene Betragsfunktion (bzw. abgeschnitten quadratisch wenns glatt sein soll)

|
-------\ | /--------
\ | /
\|/
--------------+------------
0

Wenn deine Funktion aus mehreren Teilen besteht kannst du die Teile einzeln analysieren und mit diesen Regeln feststellen ob die Gesamtfunktion konvex ist:
http://en.wikipedia.org/wiki/Convex_function#Convex_function_calculus

Wie gesagt, gerade wenn du eine kompliziertere Funktion hast ist es wichtig erstmal festzustellen ob konvex oder nicht. Für konvexe nicht glatte Funktionen wäre es dann machbar, für nicht konvex und nicht glatt wirds schwierig...

[Milchkanne]

Schönes Beispiel, danke. Manchmal kommt man auf die einfachsten Beispiele nicht, insbesondere wenn man noch nie lange drüber nachgedacht hat.

[Senior Sanchez]

Aber geschrieben hast du die Routine für die Fittness schon ja? Oder ist das bisher nur nicht diskretisiert auf dem Papier?

Und die Gefahr dass du mit der Diskussion über die Sinnhaftigkeit anfängst bestand wirklich nie... ;-)

Ich habe die geschrieben, ja. ;) Sie funktioniert auch in anderen Anwendungen und leistet dort genau was sie soll.

Hehe, das ist mir auch klar, dass ich mit der Sinnhaftigkeitsdiskussion nicht anfangen werde. Ich wollte nur darauf hinaus, dass die Funktion, so wie sie ist, sehr sinnvoll ist. Die Funktion selbst ist eine aufgebohrte Variante einer Grundfunktion, wie sie z.B. in der Airbagentwicklung mit genutzt wird und für diese Ursprungsvariante wird auch bald ein offizieller Standard verabschiedet werden.

Bezüglich der Diskussion ob konvex oder nicht. Sobald ich Zeit habe, werde ich mir das mal näher anschauen. Danke auf jeden Fall schon mal für die Erläuterung!

[pest]

einfach in Java programmieren und einen random-search algorithmus wie DDS verwenden, das sind keine 50 Zeilen. Ich fitte so sehr komplexe Modelle.

Generell ist das Auffinden eines (lokalen) Minimums nicht das Problem. Schwierig wird dann die Unsicherheitsanalyse der Paramater, was imo das Wichtigste ist.
Die meisten NonLinearFit-Methoden bassieren auf approximativen Konfidenzintervallen (Normalverteilung des Likelihoods).

Da du nicht parametrisch arbeitest bleiben dir am Ende nur Bootstrap-Methoden die aber zumindest bei einfachen Modellen gut durchführbar sind.

[AlecWhite]

Das dachte ich mir schon, dass man das als Optimierungsproblem definieren muss, aber so richtig blicke ich da nicht durch, welche Softwarekomponenten zusammenspielen.

Ich gehe nicht davon aus, dass die stetig differenzierbar ist. Ich habe es aber auch nie untersucht, weil sie sehr aufwändig ist.

Das Problem habe ich ja schon erklärt. :) Ich möchte einfach ein Polynom gegen eine Messreihe fitten und dafür meine, eventuell nicht stetig differenzierbare Fitnessfunktion verwenden.
Formuliere das bloss nicht als Optimierungs-Problem - da droht dir förmlich ein Methodenfehler.

Am Ende bekommst du zwar eine schöne Funktion raus - die ist aber jedoch quasi "zu gut" - das Rauschen in der Messreihe wird dann in das Polynom aufgenommen.

http://de.wikipedia.org/wiki/Overfitting

[Senior Sanchez]

Formuliere das bloss nicht als Optimierungs-Problem - da droht dir förmlich ein Methodenfehler.

Am Ende bekommst du zwar eine schöne Funktion raus - die ist aber jedoch quasi "zu gut" - das Rauschen in der Messreihe wird dann in das Polynom aufgenommen.

http://de.wikipedia.org/wiki/Overfitting

Das musst du mir jetzt mal erklären, denn aus meiner Sicht "lerne" ich kein Modell, dass ich dann auf unbekannte Daten anwenden möchte. Ich brauche also keine Generalisierungsfähigkeit, sondern möchte mit einer wesentlich einfacheren Funktion die Messreihe nur möglichst gut annähern.

[Senior Sanchez]

einfach in Java programmieren und einen random-search algorithmus wie DDS verwenden, das sind keine 50 Zeilen. Ich fitte so sehr komplexe Modelle.

Generell ist das Auffinden eines (lokalen) Minimums nicht das Problem. Schwierig wird dann die Unsicherheitsanalyse der Paramater, was imo das Wichtigste ist.
Die meisten NonLinearFit-Methoden bassieren auf approximativen Konfidenzintervallen (Normalverteilung des Likelihoods).

Da du nicht parametrisch arbeitest bleiben dir am Ende nur Bootstrap-Methoden die aber zumindest bei einfachen Modellen gut durchführbar sind.

Wofür steht DDS? Was meinst du denn genau mit der Unsicherheitsanalyse der Parameter?

[pest]

Wofür steht DDS?

Dynamically Dimensioned Search (http://www.cs.ubc.ca/~hutter/earg/papers09/2005WR004723.pdf)

Obwohl der Algorithmus sehr trivial ist, ist er gut für hohe Dimensionen geeignet. Vor allem wenn das Problem teilweise separierbar ist.
Ist dir ein möglichst genaues Ergebnis wichtig sind andere (ableitungsfreie) Algorithmen wie Differential Evolution (http://en.wikipedia.org/wiki/Differential_evolution), CMA-ES (http://en.wikipedia.org/wiki/CMA-ES) oder Nelder-Mead (http://de.wikipedia.org/wiki/Downhill-Simplex-Verfahren) geeignet.

Bis auf CMA-ES sind die auch recht einfach zu implementieren.

Ist das Problem glatt, diff'bar und unimodal gehen auch klassiche Quasi-Newton (http://de.wikipedia.org/wiki/Quasi-Newton-Verfahren) Verfahren


Was meinst du denn genau mit der Unsicherheitsanalyse der Parameter?

Das kommt auf den Kontext an :smile:
In einem frequentistischen Kontext geht man davon aus, dass die Parameter unbekannt aber fest sind.

Nehmen wir ein Beispiel: Wir wollen den Erwartungswert einer Stichprobe unbekannter Verteilung schätzen.

Ein konsistenter und erwartungsfreier Schätzer hierfür ist das arithmetische Mittel. Dieser Schätzer ist natürlich mit einem Fehler behaftet.

Ein 95%-Konfidenzintervall für den Schätzer beantwortet jetzt die Frage: Wenn ich 100 Stichproben derselben Größe habe gib mir einen Intervall an, dass den Erwartungswert (wahren Wert) in 95 Fällen beinhaltet.

Die Beantwortung der Frage beinhaltet das Wissen über die Verteilung des Schätzers, die man aber häufig nicht kennt.

Das arithmetische Mittel ist aber asymptotisch Normalverteilt (bei einer i.i.d Stichprobe),
das heißt ein 95%-Konfidenzintervall lässt sich (bei bekanntem Standardfehler) angeben durch [Schätzer-1.96*Standardfehler,Schätzer+1.96*Standardfehler]
Wenn man den Standardfehler schätzt nimmt man die entsprechenden Quantile der t-Verteilung.

Beim nichtparametrischen Bootstrap geht man anders vor. Man resampled die Stichprobe um neue Stichproben zu erhalten und damit eine empirische Verteilung für den Schätzer.

Beispiel: Schätzung des Erwartungswertes einer Stichprobe aus einer LogNormalverteilung (mu=0, sigma=1) (Rot: approximativ, Histogramm: Bootstrap mit 10000 Samples)
http://abload.de/img/logboot1g2bn4.png

http://abload.de/img/logboot2pzlfi.png

Bei der Regression ist es das gleiche. Es gibt unterschiedliche Quellen der Unsicherheit: z.B. Messfehler und Strukturelle Fehler des Modells

Deine Fitnesslandschaft könnte ja sehr flach sein, oder multimodal usw.

Im klassischen (additiven) Fall nichtlinearer Regression Daten=Modell+Fehler nimmt man an, dass die Fehler Normalverteilt sind. Dann ergibt sich als Maximum-Likelihood-Schätzung der Parameter die kleinste-Quadrate-Schätzung.

Der Maximum-Likelihood-Schätzer ist (unabhängig von der Fehlerannahme) asymptotisch normalverteilt und man kann Konfidenzintervalle angeben.

Beispiel: Schätzung des Parameters lambda einer Exponentialverteilung mit lambda=0.5 und Normalgestörten Daten (sigma=0.01)
http://abload.de/img/logboot36bzxt.png
Die "Estimated Variance" entspricht ungefähr der Varianz unseres eingebauten Fehlers (0.01²=0.0001), unserer Modell war schließlich perfekt.

eine Erhöhung der Streuung des Fehlers auf sigma=0.03 verbreitert das Konfidenzintervall
http://abload.de/img/logboot4jklzd.png

In deinem Fall (nicht-parametrisch) kann man immernoch Bootstrapping (http://en.wikipedia.org/wiki/Bootstrapping_%28statistics%29#Regression) verwenden.

Ich mache das mal für das obige Beispiel und nehme das Residuen-Resampling.

Das heißt, wir machen einen Fit. Merken uns die Residuen und die Antwort des Modells. Die Residuen werden dann resampled (1000 mal) und auf die Antwort addiert. Anschließend wird ein erneuter Fit durchgeführt.
Im Vergleich zum Beispiel mit dem Erwartungswert nehme ich hier die t-Verteilung weil der Standardfehler geschätzt wird.
Man sieht, dass das Ergebniss gleich dem von Mathematica ist.

http://abload.de/img/fit35eg7.png

[Senior Sanchez]

Danke schon mal für die Erläuterung, pest! Ich muss mir das heute Abend noch mal anschauen, da der Firmenproxy deine grafischen Darstellungen wegblockt.

Ich probiere übrigens auch noch andere Regressionsmodelle aus, also nicht nur Polynome, und da habe ich in einem Fall CMA-ES für mich entdeckt. Das liefert allgemein schon ziemlich gute Ergebnisse, wobei ich aber gerade auch sehe, dass CMA-ES, auch wenn es selten ist, mal in lokalen Minima hängen bleiben kann.

[Pinoccio]

Kurzer Hinweis: Ein Polynomialer Fit, der die quadratischen Abstände minimiert, geht algorithmitisch sehr einfach (Gauß und so). Einfach so eine andere Fitness-Funktion zu nehmen geht da also nicht.

mfg

[Senior Sanchez]

Kurzer Hinweis: Ein Polynomialer Fit, der die quadratischen Abstände minimiert, geht algorithmitisch sehr einfach (Gauß und so). Einfach so eine andere Fitness-Funktion zu nehmen geht da also nicht.

mfg

Das die Minimierung der quadratischen Abstände sehr einfach ist, weiß ich. :) Und das ich nicht einfach so die Fitness-Funktion ersetzen kann, dem bin ich mir auch bewusst, denn im Endeffekt kommt es eben auf den Optimierer an, der eben bestimmte Anforderungen an die Fitness-Funktion stellen.

[Senior Sanchez]

Pest, danke für deine Erläuterungen!

Ich bin aber gerade am Überlegen, an welcher Stelle ich das jetzt brauche. Kannst du mir das erklären?

[pest]

Du hast einen Punktschätzer für das Minimum einer Fitnessfunktion.

1. Deine Fitnesslandschaft ist extrem flach. Die Unsicherheit in den geschätzen Parametern ist groß

2. Deine Fitnesslandschaft ist ein peak. Die Unsicherheit in den geschätzten Parametern ist klein

3. du hast mehrere Minima: Es gibt mehrere optimale Parameter, d.h. Uneindeutigkeit

4. Parameter sind direkt abhängig und du hast unendlich viele Lösungen

Im Schnitt solltest du versuchen zumindest festzustellen ob solche Extremfälle auftreten können. Die ersteren zwei kannst du mit einer Sensitivitätsanalyse beurteilen. Die letzteren beiden mit mehrmaligem Ausführen der Optimierung.

Beim Durchstöbern auf Wikipedia bin ich immernoch unsicher, was der Unterschied zwischen Sensitivity und Uncertainty im Fall von Regression ist :freak:

[Senior Sanchez]

Ahh, okay, das hilft weiter! Da schaue ich mal.