9 Stochastisches Gradientenverfahren und Maschinelles Lernen

Dieses Kapitel ist kopiert (und, teilweise von mir und teilweise automatisch, übersetzt) aus dem Wikipedia Artikel (Stochastistic Gradient, Wikipedia contributors 2022). Dieses Kapitel steht unter der WP:CC BY-SA Lizenz.

Der stochastische Gradientenabstieg (oft als SGD abgekürzt) ist ein iteratives Verfahren zur Optimierung einer Zielfunktion mit geeigneten Glattheitseigenschaften (z. B. differenzierbar oder subdifferenzierbar). Sie kann als stochastische Näherung der Gradientenabstiegsoptimierung angesehen werden, da sie die tatsächliche Steigung (berechnet aus dem gesamten Datensatz) durch eine Schätzung davon ersetzt (berechnet aus einer zufällig ausgewählten Teilmenge der Daten). Insbesondere bei hochdimensionalen Optimierungsproblemen reduziert dies den sehr hohen Rechenaufwand, wodurch schnellere Iterationen, allerdings im Ausgleich für eine niedrigere Konvergenzrate, erreicht werden.

Während die Grundidee der stochastischen Approximation auf den Robbins-Monro-Algorithmus der 1950er Jahre zurückgeht, hat sich der stochastische Gradientenabstieg zu einer wichtigen Optimierungsmethode im maschinellen Lernen entwickelt.

9.1 Hintergrund

Sowohl in der Berechnung statistischer Schätzer als auch im Maschinellen Lernen spielt die Minimierung von Zielfunktionalen in Summenform \[\begin{equation*} Q(w) = \frac{1}{N}\sum_{i=1}^N Q_i(w) \end{equation*}\] eine Rolle, wobei der Parameter \(w\in \mathbb R^n\), der \(Q\) minimiert, gefunden oder geschätzt werden soll. Jede der Summandenfunktionen \(Q_i\) ist typischerweise assoziiert mit einem \(i\)-ten Datenpunkt (einer Beobachtung) beispielsweise aus einer Menge von Trainingsdaten.

Um obige Funktion zu minimieren, würde ein sogenannter Gradientenabstiegsverfahren den folgenden Minimierungsschritt

\[\begin{equation*} w^{k+1} := w^{k} - \eta \nabla Q(w^k) = w^k - \eta \frac{1}{N} \sum_{i=1}^N \nabla Q_i(w^k), \end{equation*}\] iterativ anwenden, wobei \(\eta\) die Schrittweite ist, die besonders in der ML community oft auch learning rate genannt wird.

Die Berechnung der Abstiegsrichtung erfordert hier also in jedem Schritt die Bestimmung von \(N\) Gradienten \(\nabla Q_i(w^k)\) der Summandenfunktionen. Wenn \(N\) groß ist, also beispielsweise viele Datenpunkte in einer Regression beachtet werden sollen, dann ist die Berechnung entsprechend aufwändig.

Andererseits entspricht die Abstiegsrichtung \[\begin{equation*} \frac{1}{N} \sum_{i=1}^N \nabla Q_i(w^k) \end{equation*}\] dem Mittelwert der Gradienten aller \(Q_i\)s am Punkt \(w_k\), der durch ein kleineres Sample

\[\begin{equation*} \frac{1}{N} \sum_{i=1}^N \nabla Q_i(w^k) \approx \frac{1}{|\mathcal J|} \sum_{j\in \mathcal J} \nabla Q_j(w^k), \end{equation*}\]

wobei \(\mathcal J \subset \{1, \dotsc, N\}\) eine Indexmenge ist, die den batch der zur Approximation gewählten \(Q_i\)s beschreibt.

9.2 Stochastisches Abstiegsverfahren

Beim stochastischen (oder “Online”) Gradientenabstieg wird der wahre Gradient von \(Q(w^k)\) durch einen Gradienten bei einer einzelnen Probe angenähert: \[\begin{equation*} w^{k+1} = w^k-\eta \nabla Q_j(w^k), \end{equation*}\] mit \(j\in \{1,\dotsc, N\}\) zufällig gewählt (ohne zurücklegen).

Während der Algorithmus den Trainingssatz durchläuft, führt er die obige Aktualisierung für jede Trainingsprobe durch. Es können mehrere Durchgänge über den Trainingssatz gemacht werden, bis der Algorithmus konvergiert. Wenn dies getan wird, können die Daten für jeden Durchlauf gemischt werden, um Zyklen zu vermeiden. Typische Implementierungen können eine adaptive Lernrate verwenden, damit der Algorithmus konvergiert.

Die wesentlichen Schritte als Algorithmus sehen wie folgt aus:

###################################################
# The basic steps of a stochastic gradient method #
###################################################

w = ...  # initialize the weight vector
eta = ... # choose the learning rate
I = [1, 2, ..., N]  # the full index set

for k in range(maxiter):
    J = shuffle(I)  # shuffle the indices
    for j in J:
        # compute the gradient of Qj at current w
        gradjk = nabla(Q(j, w))  
        # update the w vector
        w = w - eta*gradjk
    if convergence_criterion:
       break

###################################################

Die Konvergenz des stochastischen Gradientenabstiegsverfahren als Kombination von stochastischer Approximation und numerischer Optimierung ist gut verstanden. Allgemein und unter bestimmten Voraussetzung lässt sich sagen, dass das stochastische Verfahren ähnlich konvergiert wie das exakte Verfahren mit der Einschränkung, dass die Konvergenz fast sicher stattfindet.

In der Praxis hat sich der Kompromiss etabliert, der anstelle des Gradienten eines einzelnen Punktes \(\nabla Q_j(w_k)\), den Abstieg aus dem Mittelwert über mehrere Samples berechnet, also (wie oben beschrieben) \[\begin{equation*} \frac{1}{N} \sum_{i=1}^N \nabla Q_i(w^k) \approx \frac{1}{|\mathcal J|} \sum_{j\in \mathcal J} \nabla Q_j(w^k). \end{equation*}\] Im Algorithmus wird dann anstelle der zufälligen Indices \(j \in \{1, \dotsc, N\}\), über zufällig zusammengestellte Indexmengen \(\mathcal J \subset \{1, \dotsc, N\}\) iteriert.

Da die einzelnen Gradienten \(\nabla Q_j(w^K)\) unabhängig voneinander berechnet werden können, kann so ein batch Verfahren effizient auf Computern mit mehreren Prozessoren realisiert werden. Die Konvergenztheorie ist nicht wesentlich verschieden vom eigentlichen stochastischen Gradientenabstiegsverfahren, allerdings erscheint die beobachte Konvergenz weniger erratisch, da der Mittelwert statistische Ausreißer ausmitteln kann.

9.3 Aufgabe

Schreiben Sie ein Programm, das mit Hilfe eines neuronalen Netzes (NN) mit einer hidden layer \[\begin{equation*} \eta_i = NN(x_i):=\tanh \bigl (A_2 \tanh (A_1 x_i + b_1) + b_2\bigr ) \end{equation*}\]

für einen Datenpunkt \(x_i \in \mathbb R^{n_0}\), Gewichten \(A_1 \in \mathbb R^{n_1 \times n_0}\), \(b_1 \in \mathbb R^{n_1}\), \(A_2 \in \mathbb R^{1, n_1}\), \(b_2 \in \mathbb R^{1}\) und dem Ergebnisvektor \(\eta_i\in \mathbb R^{1}\), anhand der gemessenen Daten \(x_i\) die bekannte Pinguin Population penguin-data.json in zwei Gruppen aufteilt, so dass in der ersten Gruppe eine Spezies enthalten ist und in der anderen die beiden anderen Spezies.

Dazu kann eine Funktion \(\ell \colon X \mapsto \{-1, 1\}\) definiert werden, die die bekannten Pinguine \(x_i\) aus dem Datensatz \(X\) ihrer Gruppe zuordnet. Dann können die Koeffizienten des \(NN\) über das Optimierungsproblem \[\begin{equation*} \frac{1}{|X|}\sum_{x_i \in X} \|\ell(x_i)-NN(x_i)\|^2 \to \min_{A_1, b_1, A_2, b_2} \end{equation*}\] mittels des stochastischen (batch) Gradientenabstiegs bestimmt werden.

Hinweis: Für eine Funktion \(f \colon w \mapsto f(w) \in \mathbb R\), können sie den Gradienten \(\nabla_w f(w^*)\) an der Stelle \(w^*\) numerisch mit der Funktion scipy.optimize.approx_fprime bestimmen lassen.

Führen Sie die Optimierung auf einem Teil (z.B. 90%) der Pinguin Daten durch und testen sie wie gut ihr Netz funktioniert auf den verbleibenden Datenpunkten.

Der Beginn könnte also wie folgt aussehen (wobei die Größe der hidden layer zu \(n_1=2\) gesetzt ist).

import json
import numpy as np
from scipy.optimize import approx_fprime

with open('penguin-data.json', 'r') as f:
    datadict = json.load(f)

data = np.array(datadict['data'])
# centering the data
data = data - data.mean(axis=0)
lbls = np.array(datadict['target'])

# a dictionary that maps the labels(=targets) of the data into labels {1, -1}
# that will use for distinction of two groups
mplbldict = {0: np.array([-1]),
             1: np.array([-1]),
             2: np.array([1])}

# sizes of the layers
sxz, sxo, sxt = data.shape[1], 2, mplbldict[0].size
# define also the sizes of the weightmatrices

# parameters for the training -- these worked fine for me
batchsize = 30  # how many samples for the stochastic gradients
lr = 0.25  # learning rate
iterations = 200  # how many gradient steps

# the data
traindataratio = .9  # the ratio of training data vs. test data
ndata = data.shape[0]
trnds = np.int(ndata*traindataratio)
allidx = np.arange(ndata)
trnidx = np.random.choice(allidx, trnds, replace=False)
tstidx = np.setdiff1d(allidx, trnidx)

Und die Funktion, die das neuronale Netz realisiert, so:

def fwdnn(xzero, Aone=None, bone=None, Atwo=None, btwo=None):
    ''' a neural networks of two layers

    '''
    xone = np.tanh(Aone @ xzero + bone)
    xtwo = np.tanh(Atwo @ xone + btwo)
    return xtwo

Für die Umsetzung ist es hilfreich, alle Koeffizienten in einen Vektor \(w\) vorzuhalten. Damit kann dann optimiert werden und bei Bedarf die Matrizen wieder hergestellt werden:

def wvec_to_wmats(wvec):
    ''' helper to turn the vector of weights into the system matrices

    '''
    Aone = wvec[:sxz*sxo].reshape((sxo, sxz))
    cidx = sxz*sxo
    bone = wvec[cidx:cidx+sxo]
    cidx = cidx + sxo
    Atwo = wvec[cidx:cidx+sxo*sxt].reshape((sxt, sxo))
    cidx = cidx + sxo*sxt
    btwo = wvec[cidx:]
    if Aone.size + bone.size + Atwo.size + btwo.size == wvec.size:
        return Aone, bone, Atwo, btwo
    else:
        raise UserWarning('mismatch weightsvector/matrices')

Damit fehlen zum Programm insbesondere noch

die Zielfunktion der Optimierung in einem Datenpunkt \(x_i\)
eine Funktion, die den stochastischen Gradienten über einem kleinen batch von Daten ausrechnet
eine Iteration, die das Gradientenabstiegsverfahren zusammen mit der learning rate durchführt
ein Loop, der testet, wie das optimierte \(NN\), die verbliebenen Daten klassifiziert

Hinweis: Auch hier ist wieder sehr viel Zufall involviert. Im Zweifel lieber zweimal testen.

Hinweis: Der letzte Testloop könnte so aussehen:

print('***** testing the classification *****')
faillst = []  # list to collect the failures for later examination
for cti in tstidx:  # iteration over the test data points
    itrgt = data[cti, :]
    ilbl = mplbldict[lbls[cti]]
    # the prediction of the neural network -- Aonex, ... are the optimized weights
    nnlbl = fwdnn(itrgt, Aone=Aonex, bone=bonex, Atwo=Atwox, btwo=btwox)
    sccs = np.sign(ilbl) == np.sign(nnlbl)
    print(f'label: {ilbl.item()} -- nn: {nnlbl.item():.4f} -- success: {sccs}')
    if not sccs:
        faillst.append((cti, ilbl.item(), nnlbl.item(),
                        datadict['target_names'][lbls[cti]]))
    else:
        pass

print('\n***** Results *****')
print(f'{100-len(faillst)/tstidx.size*100:.0f}% was classified correctly')
print('***** Misses *****')
if len(faillst) == 0:
    print('None')
else:
    for cfl in faillst:
        cid, lbl, nnlbl, name = cfl
        print(f'ID: {cid} ({name} pinguin) was missclassified ' +
              f'with score {nnlbl:.4f} vs. {lbl}')

Referenzen

Wikipedia contributors: Stochastic gradient descent — Wikipedia, the free encyclopedia, https://en.wikipedia.org/w/index.php?title=Stochastic_gradient_descent&oldid=1098148439, (2022)