Hessematrix durch finite Differenzen approximieren

Hallo,

ich möchte die Hessematrix einer Funktion f(x) mit $f : R^n \rightarrow R$ mit finiten Differenzen approximieren, wobei f(x) mindestens zweimal differenzierbar ist. Dazu entwickele ich den Gradienten im Punkt x wie folgt

$\nabla f(x + tp) = \nabla f(x) + t \nabla^2 f(x) p + O(t^2)$

bzw. nach der Hessematrix umgeformt

$\nabla^2 f(x) p = \frac{\nabla f(x + tp) - \nabla f(x)}{t}$

Nun setze ich $t = \sqrt{EPS}$ und $p = e_i$ als i-ten Basisvektor und erhalte damit die i-te Spalte der Hessematrix

$\frac{\nabla f(x + t e\_i) - \nabla f(x)}{t} = \nabla^2 f(x) e\_i = \frac{\partial f}{\partial x\_j \partial x\_i}, \quad j=1,...,n$

Das Problem ist nun: der Gradient wurde auch durch finite Differenzen berechnet; sollte aber trotzdem gehen. Allerdings weicht die approximierte Hessematrix stark von der analytisch berechneten Loesung ab und ist teilweise sehr unsymmetrisch.

Andererseits: wenn ich in der obigen Formel die tatsaechlichen Gradienten anstatt der approximierten Gradienten nehme, dann passt alles und der Approximationsfehler der Hessematrix (Norm) ist < 1e-7.
Die approximierten Gradienten werden allerdings auch bis auf < 1e-7 genau berechnet, sind also auch richtig.

Wenn man beide Approximationen, die fuer sich genommen funktionieren, kombiniert, dann kommt nur Murks raus. Warum? Ich habe auch schon zentrale Differenzenquotienten verwendet, fuehrt aber zum gleichen Ergebnis.

otze

1. Die Ableitung einer verrauschten Variable ist extrem instabil zu berechnen.
Das ist einfach zu sehen:wenn du 2 verrauschte variablen voneinander abziehst, verdoppelt sich das Rauschen der Differenz im schlimmsten Fall und dann teilst du das rauschen durch den winzigen Abstand -> dein Rauschen explodiert.

angenommen du kommst bis auf ε ran:
$\frac{(\nabla f(x + tp) \pm \epsilon) - (\nabla f(x)\pm \epsilon)}{t}=\nabla^2 f(x) \pm 2\frac{\epsilon} t$

2.benutze symmetrische differenzen für weniger rauschen.
3.anstatt den gradienten zuerst zu approximieren,approximier die hessematrix lieber direkt mithilfe eines symmetrischen dreipunktverfahrens.

otze schrieb:

angenommen du kommst bis auf ε ran:
$\frac{(\nabla f(x + tp) \pm \epsilon) - (\nabla f(x)\pm \epsilon)}{t}=\nabla^2 f(x) \pm 2\frac{\epsilon} t$

Danke, so ist es unmittelbar einsichtig.

Die Rechnung suggeriert allerdings auch, dass man den Fehler klein machen kann, wenn man $t$ gross gegenueber $\epsilon$ waehlt. Beispiel: $\epsilon = 10^{-7}, \, t = 10^{-4} \, \rightarrow 2 \frac{\epsilon}{t} = 2 \cdot 10^{-3}$ . Das funktioniert bei mir auch; allerdings entfernt man sich damit immer weiter vom eigentlichen Differenzial, so dass (3.) wohl doch die beste Option ist.

krümelkacker

Darf ich fragen, das das für eine Funktion ist? Was gibt sie zurück? Wenn das so eine Art Summe von ggf gewichteten quadratischen Abweichungen, also so etwas wie $f(x) = \frac{1}{2} \left\| F(x) \right\|_2^2$ mit $F: \mathbb{R}^n \mapsto \mathbb{R}^m$ , dann bekommst du die Hesse-Matrix von f über die Jacobi-Matrix von F. Und die exakte Jacobi-Matrix für so etwas lässt sich per automatischer Differentiation 1. Ordnung ausrechnen. Das ist gar nicht sooo schwer. Man findet hier und da schon fertige AD-Bibliotheken.

Das schöne an AD ist, dass du dabei keine Schrittweite wählen musst und deswegen auch die Probleme, die du mit zu hoher Schrittweite (Diskretisierungsfehler) aber auch zu niedriger Schrittweite (Auslöschung) bekommen würdest, komplett umgehen kannst. Eine sehr geile Sache also!