Almost finished gradient descent

programming/lecture/programming.tex

@@ -119,7 +119,6 @@ wie eine einzelne Variable gel\"oscht wird.
 
 Der Datentyp bestimmt, wie die im Speicher abgelegten Bitmuster
 interpretiert werden. Die wichtigsten Datentpyen sind folgende:
-
 \begin{itemize}
 \item \codetterm{integer} - Ganze Zahlen. Hier gibt es mehrere
   Unterarten, die wir in \matlab{} (meist) ignorieren k\"onnen.
@@ -148,7 +147,7 @@ unterschiedlichem Speicherbedarf und Wertebreich.
 Daten gespeichert werden. Dennoch lohnt es sich, sich ein wenig mit
 den Datentypen auseinanderzusetzen (Box \ref{daqbox}). 
 
-\begin{ibox}[tp]{\label{daqbox}Digitalisierung von Messdaten}
+\begin{ibox}[t]{\label{daqbox}Digitalisierung von Messdaten}
   Szenario: Die elektrische Aktivit\"at (z.B. die Membranspannung)
   einer Nervenzelle wird gemessen. Die gemessene Spannung wird mittels
   Messkarte digitalisiert und auf dem Computer f\"ur weitere Analysen

regression/code/errorSurface.m

@@ -1,6 +1,3 @@
-clear
-close all
-
 load('lin_regression.mat');
 
 % compute mean squared error for a range of sloopes and intercepts:
@@ -23,7 +20,3 @@ zlabel('error')
 set(gca,'xtick', (-5:2.5:5))
 grid on
 view(3)
-
-set(gcf, 'paperunits', 'centimeters', 'papersize', [15, 15], ... 
-    'paperposition', [0., 0., 15, 15])
-saveas(gcf, 'error_surface', 'pdf')

regression/code/lsqGradient.m

@@ -1,7 +1,17 @@
-function gradient = lsqGradient(parameter, x, y)
+function gradient = lsqGradient(x, y, parameter)
-h = 1e-6;
+% The gradient of the least square error
+%
+% Arguments: x, the input values
+%            y, the corresponding measured output values
+%            parameter, vector containing slope and intercept 
+%                       as the 1st and 2nd element 
+%
+% Returns: the gradient as a vector with two elements
 
-partial_m = (lsqError([parameter(1)+h, parameter(2)],x,y) - lsqError(parameter,x,y))/ h;
+  h = 1e-6;   % stepsize for derivatives
-partial_n = (lsqError([parameter(1), parameter(2)+h],x,y) - lsqError(parameter,x,y))/ h;
+
+  partial_m = (lsqError(x, y, [parameter(1)+h, parameter(2)]) - lsqError(x, y, parameter))/ h;
+partial_n = (lsqError(x, y, [parameter(1), parameter(2)+h]) - lsqError(x, y, parameter))/ h;
 
   gradient = [partial_m, partial_n];
+end

regression/lecture/derivative.py

@@ -0,0 +1,53 @@
+import numpy as np
+import matplotlib.pyplot as plt
+
+#plt.xkcd()
+fig = plt.figure( figsize=(2.5,3.4) )
+ax = fig.add_subplot(1, 1, 1)
+
+# parabula:
+x1 = -0.2
+x2 = 0.9
+x = np.linspace(x1, x2, 200)
+y = x*x
+ax.set_xlim(x1, x2)
+ax.set_ylim(-0.2, 0.7)
+ax.plot(x, y, c='b', lw=4, zorder=0)
+# secant:
+x = np.asarray([0.1, 0.7])
+y = x*x
+ax.set_xticks(x)
+ax.set_yticks(y)
+ax.set_xticklabels(['$x$','$x+\Delta x$'])
+ax.set_yticklabels(['',''])
+ax.scatter(x, y, c='#CC0000', edgecolor='none', s=150, zorder=10)
+# function values:
+ax.plot([x[0], x[0], x1],[-0.2, y[0], y[0]], '--k', lw=1, zorder=6)
+ax.plot([x[1], x[1], x1],[-0.2, y[1], y[1]], '--k', lw=1, zorder=6)
+ax.text(x1+0.05, y[0]+0.05, '$f(x)$', zorder=6)
+ax.text(x1+0.05, y[1]+0.05, '$f(x+\Delta x)$', zorder=6)
+# slope triangle:
+ax.plot([x[0], x[1], x[1]],[y[0], y[0], y[1]], '-k', lw=2, zorder=7)
+ax.text(np.mean(x), y[0]-0.08, '$\Delta x$', ha='center', zorder=7)
+ax.text(x[1]+0.05, np.mean(y), '$f(x+\Delta x)-f(x)$', va='center', rotation='vertical', zorder=7)
+# secant line:
+m = np.diff(y)/np.diff(x)
+xl = [x1, x2]
+yl = m*(xl-x[0])+y[0]
+ax.plot(xl, yl, c='#CC0000', lw=3, zorder=7)
+
+# derivative:
+md = 2.0*x[0]
+yl = md*(xl-x[0])+y[0]
+ax.plot(xl, yl, c='#FFFF66', lw=3, zorder=5)
+
+# limit:
+for ml in np.linspace(md, m, 5)[1:] :
+    yl = ml*(xl-x[0])+y[0]
+    xs = 0.5*(ml+np.sqrt(ml*ml-4.0*(ml*x[0]-y[0])))
+    ax.scatter([xs], [xs*xs], c='#FFCC00', edgecolor='none', s=80, zorder=3)
+    ax.plot(xl, yl, c='#FFCC00', lw=2, zorder=3)
+
+plt.tight_layout();
+plt.savefig('derivative.pdf')
+#plt.show();

regression/lecture/error_surface.py

@@ -41,8 +41,7 @@ def plot_error_plane(ax, x, y, m, n):
                     linewidth=0, shade=True)
     # Minimum:
     mini = np.unravel_index(np.argmin(error_surf), error_surf.shape)
-    #ax.scatter([m], [n], [0.0], color='#cc0000')
+    ax.scatter(slopes[mini[1]], intercepts[mini[0]], [0.0], color='#cc0000', s=60)
-    ax.scatter(slopes[mini[1]], intercepts[mini[0]], [0.0], color='#cc0000')
 
 
 if __name__ == "__main__":

regression/lecture/gradient.py

@@ -0,0 +1,42 @@
+import numpy as np
+import matplotlib.pyplot as plt
+
+#plt.xkcd()
+fig = plt.figure( figsize=(3.5,3.7) )
+ax = fig.add_subplot(1, 1, 1)
+
+def gaussian(x, y) :
+    return np.exp(-0.5*(x*x+y*y))
+
+def gaussiangradient(x, y) :
+    return np.asarray([-x*np.exp(-0.5*(x*x+y*y)), -y*np.exp(-0.5*(x*x+y*y))])
+
+# Gauss contour lines:
+x1 = -2.0
+x2 = 2.0
+x = np.linspace(x1, x2, 200)
+y = np.linspace(x1, x2, 200)
+X, Y = np.meshgrid(x, y)
+Z = gaussian(X, Y)
+ax.set_xlabel('x')
+ax.set_ylabel('y', rotation='horizontal')
+ax.set_xticks(np.arange(x1, x2+0.1, 1.0))
+ax.set_yticks(np.arange(x1, x2+0.1, 1.0))
+ax.contour(X, Y, Z, zorder=0)
+
+# gradients:
+xxg = [-1.1, 1.4, -1.0]
+yyg = [-1.1, 0.6, 0.75]
+for xg, yg in zip(xxg, yyg) :
+    zg = gaussian(xg, yg)
+    g = gaussiangradient(xg, yg)
+    ax.scatter(xg, yg, c='k', s=100, zorder=10)
+    ax.quiver([xg, xg], [yg, yg], [g[0], 0.0], [0.0, g[1]], units='xy', angles='uv', scale_units='x', scale=0.5, zorder=10)
+    ax.quiver([xg], [yg], [g[0]], [g[1]], units='xy', angles='uv', scale_units='x', scale=0.5, zorder=10, lw=2)
+ax.text(-0.8, -1.5, '$\partial f/\partial x$', ha='center', zorder=20)
+ax.text(-1.55, -0.8, '$\partial f/\partial y$', va='center', rotation='vertical', zorder=20)
+ax.text(-0.4, -0.8, r'$\nabla f$', ha='left', zorder=20)
+
+plt.tight_layout();
+plt.savefig('gradient.pdf')
+#plt.show();

regression/lecture/regression.tex

@@ -201,83 +201,124 @@ bestimmt werden soll erh\"oht die Anzahl der n\"otigen
 Berechnungen. Wir suchen also ein Verfahren, dass das Minimum der
 Kostenfunktion mit m\"oglichst wenigen Berechnungen findet.
 
+\begin{ibox}[t]{\label{differentialquotientbox}Differenzenquotient und Ableitung}
+  \includegraphics[width=0.33\textwidth]{derivative}
+  \hfill
+  \begin{minipage}[b]{0.63\textwidth}
+    Der Differenzenquotient 
+    \[ m = \frac{f(x + \Delta x) - f(x)}{\Delta x} \] einer Funktion $y = f(x)$ ist
+    die Steigung der Sekante (rot) durch die beiden Punkte $(x,f(x))$ und
+    $(x+\Delta x,f(x+\Delta x))$ mit dem Abstand $\Delta x$.
+
+    Die Steigung einer Funktion $y=f(x)$ an einer Stelle $x$ (gelb) wird durch
+    die Ableitung $f'(x)$ der Funktion an dieser Stelle berechnet.  Die
+    Ableitung ist \"uber den Grenzwert (orange) des Differenzenquotienten f\"ur
+    unendlich kleine Abst\"ande $\Delta x$ definiert:
+    \[f'(x) = \frac{{\rm d} f(x)}{{\rm d}x} = \lim\limits_{\Delta x \to 0} \frac{f(x + \Delta x) - f(x)}{\Delta x} \]
+  \end{minipage}
+\end{ibox}
+
+\begin{ibox}[t]{\label{partialderivativebox}Partielle Ableitungen und Gradient}
+  Bei Funktionen
+  \[ z = f(x,y) \]
+  die von mehreren Variablen, z.B. $x$ und $y$ abh\"angen,
+  kann die Steigung in Richtung jeder dieser Variablen
+  mit den partiellen Ableitungen 
+  \[ \frac{\partial f(x,y)}{\partial x} = \lim\limits_{\Delta x \to 0} \frac{f(x + \Delta x,y) - f(x,y)}{\Delta x} \]
+  und
+  \[ \frac{\partial f(x,y)}{\partial y} = \lim\limits_{\Delta y \to 0} \frac{f(x, y + \Delta y) - f(x,y)}{\Delta y} \]
+  definiert \"uber den jeweiligen Differenzenquotienten
+  (Box~\ref{differentialquotientbox}) berechnet werden.  \vspace{1ex}
+
+  \begin{minipage}[t]{0.44\textwidth}
+    \mbox{}\\[-2ex]
+    \includegraphics[width=1\textwidth]{gradient}
+  \end{minipage}
+  \hfill
+  \begin{minipage}[t]{0.52\textwidth}
+    Z.B. lauten die partiellen Ableitungen von 
+    \[ f(x,y) = x^2+y^2 \]
+    \[ \frac{\partial f(x,y)}{\partial x} = 2x \; , \quad \frac{\partial f(x,y)}{\partial y} = 2y \; .\]
+
+    Der Gradient ist der aus den partiellen Ableitungen gebildete Vektor
+    \[ \nabla f(x,y) = \left( \begin{array}{c} \frac{\partial f(x,y)}{\partial x} \\[1ex] \frac{\partial f(x,y)}{\partial y} \end{array} \right) \]
+    und zeigt in Richtung des st\"arksten Anstiegs der Funktion $f(x,y)$.
+  \end{minipage}
+
+  \vspace{1ex} Die Abbildung zeigt die Kontourlinien einer bivariaten
+  Gau{\ss}glocke $f(x,y) = \exp(-(x^2+y^2)/2)$ und den Gradienten mit
+  seinen partiellen Ableitungen an drei verschiedenen Stellen.
+\end{ibox}
+
 
 \section{Gradient}
 
-Man kann sich den Optimierungsprozess veranschaulichen wenn man sich
+Wenn eine Kugel an einem beliebigen Startpunkt auf der Fehlerfl\"ache
-vorstellt, dass eine Parameterkombination einen Punkt auf der
+\figref{errorsurfacefig} losgelassen werden w\"urde, dann w\"urde sie
-Fehlerfl\"ache definiert. Wenn von diesem Punkt aus eine Kugel
+entlang des steilsten Gef\"alles zum Minimum der Fehlerfl\"ache rollen
-losgelassen w\"urde, dann w\"urde sie, dem steilsten Gef\"alle
+und dort zum Stehen kommen. Um den Weg der Kugel mit einem
-folgend, zum Minimum der Fehlerfl\"ache rollen und dort zum Stehen
+Computerprogramm nachvollziehen zu k\"onnen, ben\"otigen wir
-kommen. Um dem Computer zu sagen, in welche Richtung die Position
+Information \"uber die Richtung des steilsten Gef\"alles an der
-ver\"andert werden soll muss also die Steigung an der aktellen
+jeweils aktuellen Position.
-Position berechnet werden.
+
-
+Der Gradient (Box~\ref{partialderivativebox}) der Kostenfunktion
-Die Steigung einer Funktion an einer Stelle ist die Ableitung der
+\[ \nabla e(m,b) = \left( \begin{array}{c} \frac{\partial
-Funktion an dieser Stelle, die durch den Differenzquotienten f\"ur
+      e(m,b)}{\partial m} \\[1ex] \frac{\partial f(m,b)}{\partial
-unendlich kleine Schritte $h$ bestimmt wird.
+      b} \end{array} \right) \]
-\[f'(x) = \lim\limits_{h \rightarrow 0} \frac{f(x + h) - f(x)}{h} \]
+bzgl. der beiden Parameter $m$ und $b$ der Geradengleichung ist ein
-Bei unserem Fittingproblem h\"angt die Fehlerfunktion von zwei
+Vektor, der in Richtung des steilsten Anstiegs der Kostenfunktion
-Parametern ab und wir bestimmen die Steigung partiell f\"ur jeden
+$e(m,b)$ zeigt.  Die L\"ange des Gradienten gibt die St\"arke des
-Parameter einzeln. Die partielle Ableitung nach $m$ sieht so aus:
+Anstiegs an (\figref{gradientquiverfig})).  Da wir aber abw\"arts zum
-\[\frac{\partial g(m,n)}{\partial m} = \lim\limits_{h \rightarrow 0} \frac{g(m + h, n) - g(m,n)}{h}\]
+Minimum laufen wollen, m\"ussen wir die dem Gradienten
-Da wir die Ver\"anderung des Fehlers in Abh\"angigkeit der beiden
+entgegengesetzte Richtung einschlagen.
-Parameter bewerten, ist der Gradient an der Stelle $(m,n)$ ein
-zweielementigen Vektor, der aus den partiellen Ableitungen nach $m$ und
-nach $n$ besteht. Die Richtung des Gradienten zeigt die Richtung der
-gr\"o{\ss}ten Steigung an, seine L\"ange repr\"asentiert die St\"arke
-des Gef\"alles.
-
-\[\bigtriangledown g(m,n) = \left( \frac{\partial g(m,n)}{\partial m}, \frac{\partial g(m,n)}{\partial n}\right)\]
-Die partiellen Ableitungen m\"ussen nicht analytisch berechnet werden
-sondern wir n\"ahern sie numerisch durch einen sehr kleinen Schritt
-an.
-\[\frac{\partial g(m,n)}{\partial m} = \lim\limits_{h \rightarrow
-  0} \frac{g(m + h, n) - g(m,n)}{h} \approx \frac{g(m + h, n) -
-  g(m,n)}{h}\]
 
-Die graphische Dargestellung der Gradienten (Abbildung
+Die partiellen Ableitungen m\"ussen nicht analytisch berechnet werden
-\ref{gradientquiverfig}) zeigt, dass diese in die Richtung der
+sondern k\"onnen numerisch entsprechend dem Differenzenquotienten
-gr\"o{\ss}ten Steigung zeigen. Um zum Minimum der Fehlerfunktion zu
+(Box~\ref{differentialquotientbox}) mit kleinen Schrittweiten $\Delta
-gelangen sollte man also die entgegengesetzte Richtung einschlagen.
+m$ und $\Delta b$ angen\"ahert werden:
+\[\frac{\partial e(m,b)}{\partial m} = \lim\limits_{\Delta m \to
+  0} \frac{e(m + \Delta m, b) - e(m,b)}{\Delta m} \approx \frac{e(m + \Delta m, b) -
+  e(m,b)}{\Delta m} \; . \]
 
 \begin{figure}
-  \includegraphics[width=0.75\columnwidth]{figures/error_gradient}
+  \includegraphics[width=0.75\columnwidth]{error_gradient}
-  \caption{\textbf{Darstellung der Gradienten f\"ur jede
+  \titlecaption{Der Gradienten der Fehlerfl\"ache.} 
-      Parameterkombination.} Jeder Pfeil zeigt die Richtung und die
+  {Jeder Pfeil zeigt die Richtung und die
     Steigung f\"ur verschiedene Parameterkombination aus Steigung und
-    y-Achsenabschnitt an. Die Kontourlinien im Hintergrund
+    $y$-Achsenabschnitt an. Die Kontourlinien im Hintergrund
     illustrieren die Fehlerfl\"ache. Warme Farben stehen f\"ur
     gro{\ss}e Fehlerwerte, kalte Farben f\"ur kleine. Jede
     Kontourlinie steht f\"ur eine Linie gleichen
     Fehlers.}\label{gradientquiverfig}
 \end{figure}
 
-\begin{exercise}{lsqGradient.m}{}\label{gradientexercise}
+\begin{exercise}{lsqGradient.m}{}\label{gradientexercise}%
-  Implementiert eine Funktion \code{lsqGradient}, die den aktuellen
+  Implementiere eine Funktion \code{lsqGradient}, die die $x$- und
-  Parametersatz als 2-elementigen Vektor, die x-Werte und die y-Werte
+  $y$-Werte der Messdaten sowie den Parametersatz $(m, b)$ der
-  als Argumente entgegennimmt und den Gradienten zur\"uckgibt.
+  Geradengleichung als 2-elementigen Vektor als Argumente
+  entgegennimmt und den Gradienten an dieser Stelle zur\"uckgibt.
 \end{exercise}
 
 \begin{exercise}{errorGradient.m}{}
-  Benutzt die Funktion aus vorheriger \"Ubung (\ref{gradientexercise}),
+  Benutzt die Funktion aus der vorherigen \"Ubung (\ref{gradientexercise}),
-  um f\"ur die jede Parameterkombination aus der Fehlerfl\"ache
+  um f\"ur jede Parameterkombination aus der Fehlerfl\"ache
   (\"Ubung \ref{errorsurfaceexercise}) auch den Gradienten zu
   berechnen und darzustellen. Vektoren im Raum k\"onnen mithilfe der
   Funktion \code{quiver} geplottet werden.
 \end{exercise}
 
+
 \section{Gradientenabstieg}
 
-Zu guter Letzt muss ``nur'' noch der Gradientenabstieg implementiert
+Zu guter Letzt muss nur noch der Gradientenabstieg implementiert
-werden. Die daf\"ur ben\"otigten Zutaten sollten wir aus den
+werden. Die daf\"ur ben\"otigten Zutaten haben wir aus den
-vorangegangenen \"Ubungen haben. Wir brauchen: 1. Die Fehlerfunktion
+vorangegangenen \"Ubungen bereits vorbereitet. Wir brauchen: 1. Die Fehlerfunktion
 (\code{meanSquareError.m}), 2. die Zielfunktion (\code{lsqError.m})
 und 3. den Gradienten (\code{lsqGradient.m}).  Der Algorithmus
 f\"ur den Abstieg lautet:
 \begin{enumerate}
 \item Starte mit einer beliebigen Parameterkombination $p_0 = (m_0,
-  n_0)$.
+  b_0)$.
 \item \label{computegradient} Berechne den Gradienten an der akutellen Position $p_i$.
 \item Wenn die L\"ange des Gradienten einen bestimmten Wert
   unterschreitet, haben wir das Minum gefunden und k\"onnen die Suche
@@ -287,14 +328,14 @@ f\"ur den Abstieg lautet:
   wird (z.B. \code{norm(gradient) < 0.1}).
 \item \label{gradientstep} Gehe einen kleinen Schritt ($\epsilon =
   0.01$) in die entgegensetzte Richtung des Gradienten:
-  \[p_{i+1} = p_i - \epsilon \cdot \bigtriangledown g(m_i, n_i)\]
+  \[p_{i+1} = p_i - \epsilon \cdot \nabla e(m_i, b_i)\]
 \item Wiederhole die Schritte \ref{computegradient} -- \ref{gradientstep}.
 \end{enumerate}
 
 Abbildung \ref{gradientdescentfig} zeigt den Verlauf des
 Gradientenabstiegs. Von einer Startposition aus wird die Position
 solange ver\"andert, wie der Gradient eine bestimmte Gr\"o{\ss}e
-\"uberschreitet. An den Stellen, an denen der Gradient sehr stark ist
+\"uberschreitet. An den Stellen, an denen der Gradient sehr stark ist,
 ist auch die Ver\"anderung der Position gro{\ss} und der Abstand der
 Punkte in Abbildung \ref{gradientdescentfig} gro{\ss}.
 
@@ -318,3 +359,37 @@ Punkte in Abbildung \ref{gradientdescentfig} gro{\ss}.
   \item Erstelle einen Plot, der den besten Fit in die Daten plottet.
   \end{enumerate}
 \end{exercise}
+
+
+\section{Fazit}
+Mit dem Gradientenabstieg haben wir eine wichtige Methode zur
+Bestimmung eines globalen Minimums einer Kostenfunktion
+kennengelernt. F\"ur den Fall des Kurvenfits mit einer
+Geradengleichung zeigt der mittlere quadratische Abstand als
+Kostenfunktion in der Tat ein einziges klar definiertes Minimum.  
+
+Wie wir im n\"achsten Kapitel sehen werden, kann die Position des
+Minimums bei Geradengleichungen sogar analytisch bestimmt werden, der
+Gradientenabstieg w\"are also gar nicht n\"otig
+\matlabfun{polyfit}. Bei allen nichtlinearen Funktionen, wie z.B. der
+Exponentialfunktion $f(x;\lambda) = \exp(\lambda x)$ gibt es aber
+keine analytische L\"osung, und das Minimum der Kostenfunktion muss
+numerisch, z.B. mit dem Gradientenabstiegsverfhren bestimmt werden.
+
+Um noch schneller das Minimum zu finden, kann der Gradientenabstieg
+auf vielf\"altige Wei{\ss}e verfeinert werden. z.B. kann die
+Schrittweite an die St\"arke des Gradienten angepasst werden. Diese
+numerischen Tricks sind in bereits vorhandenen Funktionen
+implementiert.  Allgemeine Funktionen sind f\"ur beliebige
+Kostenfunktionen gemacht \matlabfun{fminsearch}, w\"ahrend spezielle
+Funktionen z.B. f\"ur die Minimierung des quadratischen Abstands bei
+einem Kurvenfit angeboten werden \matlabfun{lsqcurvefit}.
+
+\begin{important}
+  Das Finden des globalen Minimums ist leider nur selten so leicht wie
+  bei einem Geradenfit. Oft gibt es viele Nebenminima, in denen der
+  Gradientenabstieg enden kann, obwohl das gesuchte globale Minimum
+  noch weit entfernt ist. Darum ist es meist sehr wichtig, wirklich
+  gute Startwerte f\"ur die zu bestimmenden Parameter der
+  Kostenfunktion zu haben.
+\end{important}