[PYTHON] Introduction à OPTIMIZER ~ De la régression linéaire à Adam à Eve

À propos de cet article

Présentation de l'optimiseur utilisé dans l'apprentissage automatique.

Introduction Commençons par l'exemple de la régression linéaire, que vous aimez tous. La plupart du temps, vous mettez des données dans un package d'apprentissage automatique et vous vous retrouvez avec. Vous ne remarquerez peut-être même pas que l'optimiseur se tortille dans le package. Mais ils sont là.

Un peu de théorie

Maintenant, disons que la variable objective est $ y $, la variable explicative est $ X $ et le coefficient de régression est $ w $, bien que ce soit un peu brusque. En régression linéaire, nous voulons approcher la variable objective avec une combinaison linéaire de variables explicatives. C'est,

y \sim Xw \tag{1.1}

Nous cherchons à déterminer le coefficient $ w $ pour qu'il soit le plus «bon» possible. Ici, $ y $ et $ w $ sont traités comme des vecteurs verticaux, et $ X $ est traité comme une matrice. Si le nombre d'échantillons de données est $ N $ et le nombre de variables explicatives est $ M $, alors $ y $ est la dimension $ N $, $ w $ est la dimension $ M $ et $ X $ est la dimension $ N \ fois M $. est. Veuillez noter que cet article n'inclut pas de termes constants dans un souci de simplicité. Vous devez définir la fonction de perte pour déterminer $ w $ «bon». Erreur carrée

L(w) := (y-Xw)^2 \tag{1.2}

Est typique, donc je vais l'utiliser. Pour trouver la valeur de $ w $ qui minimise la fonction de perte $ L $, il est possible de différencier et de trouver le point où "= 0", alors différencions.

\frac{\partial L }{\partial w} = -2(X^Ty-X^TXw)=0 \tag{1.3}

Que,

w=\big(X^TX\big)^{-1}X^Ty \tag{1.4}

Il peut être obtenu. $ X ^ T $ est une matrice transposée de $ X $. Cela résout tout. Après cela, si vous pensez qu'il est correct de lancer mécaniquement $ X $ et $ y $ dans (1.4), c'est une grosse erreur. En premier lieu, (1.4) n'a pas de sens si l'inverse de $ X ^ TX $ n'existe pas. Prenons un exemple concret. Essayez de régler comme suit.

X=\begin{pmatrix}
1 & 1 & 1 \\
1 & 2 & 3
\end{pmatrix},\;
y=\begin{pmatrix}
0\\1
\end{pmatrix},\;
w=\begin{pmatrix}
a\\b\\c
\end{pmatrix} \tag{1.5}

Dans ce cas, si vous essayez de calculer concrètement $ X ^ TX $

X^TX=
\begin{pmatrix}
1 & 1 \\
1 & 2 \\
1 & 3 \\
\end{pmatrix}
\begin{pmatrix}
1 & 1 & 1 \\
1 & 2 & 3
\end{pmatrix}
=
\begin{pmatrix}
2 & 3 & 4 \\
3 & 5 & 7 \\
4 & 7 & 10 
\end{pmatrix} \tag{1.6}

Il est devenu. Si vous demandez également la formule matricielle

\text{det}\big(X^TX\big) = 0 \tag{1.7}

Il s'avère qu'il n'y a pas de matrice inverse. Alors, la solution de (1.3) n'existe-t-elle pas? Ce n'est pas vrai. Au contraire, il y en a d'innombrables. Sous le réglage de (1.5), (1.1) est une équation simultanée

a+b+c=0, \;\; a+2b+3c=1 \tag{1.8}

Est équivalent à Cependant, comme il existe deux équations pour trois variables, la solution ne peut pas être déterminée de manière unique. En fait, si vous résolvez $ a $ et $ b $ pour $ c $, un ensemble de solutions avec $ c $ comme paramètre.

a=c-1, \;\; b=1-2c \tag{1.9}

Il peut être obtenu.

Essayez d'utiliser le package

Maintenant, avez-vous des doutes ici? Qu'est-ce qui sera retourné si vous le jetez dans le package d'apprentissage automatique avec le paramètre (1.5)? Faisons le. Utilisez python.

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
from mpl_toolkits.mplot3d import Axes3D
sns.set_context('poster')

Créez les données correspondant à (1.5).

X = np.array([[1, 1, 1],[1, 2, 3]])
y = np.array([[0],[1]])

Je vais le jeter dans scikit-learn.

from sklearn.linear_model import LinearRegression
np.random.seed(0)
model = LinearRegression(fit_intercept=False)
model.fit(X,y)

Affiche le coefficient de régression après l'apprentissage.

model.coef_

Les résultats sont $ a = -0,5 $, $ b = 0 $, $ c = 0,5 $. Il rencontre certainement (1,8), donc c'est certainement l'un des (1,9). Cependant, comme il n'y a pas de matrice inverse de $ X ^ TX $, il est certain qu'elle n'a pas été calculée à partir de (1.4). Alors, qu'est-ce-qu'il s'est passé?

Gradient Decent Parlons un peu plus général, pas seulement la régression linéaire. Supposons maintenant que vous vouliez connaître la valeur minimale (ou maximale) d'une fonction $ f (x) $. L'argument $ x $ est un vecteur. À un certain moment, fixez un $ x_0 $ et considérez une petite fluctuation $ x_0 $ → $ x_0 + \ Delta x $ autour de $ x_0 $. Expansion en laissant jusqu'au terme de première commande de $ \ Delta x $

f(x_0+\Delta x)\simeq f(x_0)+\Delta x \cdot \frac{\partial f}{\partial x} (x_0) \tag{1.10}

Peut être écrit. "$ \ Cdot $" représente le produit interne des vecteurs. À propos, lorsque la direction de $ \ Delta x $ est modifiée de diverses manières, si le montant de la variation de $ f (x) $ se poursuit dans la direction du plus petit (ou plus grand), la valeur minimale (ou la valeur maximale) est bien recherchée. Je pense que je peux le faire. À cette fin, le produit intérieur doit être minimisé (ou maximisé).

\Delta x \propto \mp \frac{\partial f}{\partial x} (x_0) \tag{1.11}

Ça devrait être. En d'autres termes, si vous voulez trouver la valeur minimale, après avoir choisi la position initiale de manière appropriée, définissez l'indice du pas de temps sur $ t $.

x_{t+1} = x_{t} -\eta \frac{\partial f}{\partial x}(x_t) \tag{1.12}

Il vous suffit de répéter la recherche séquentiellement comme ceci. $ \ eta $ est une constante appelée taux d'apprentissage que vous définissez manuellement une "petite" valeur positive. Il s'agit de la technique appelée Gradient Decent, qui est l'optimiseur le plus basique. Dans le cas de la régression linéaire (1.1), si vous notez la formule de mise à jour de Gradient Decent

\begin{align}
w_{t+1} &= w_t - \eta \frac{\partial L}{\partial w}(w_t) \\
        &= w_t - \eta \big[-2(X^Ty-X^TXw_t) \big] \tag{1.13}
\end{align}

Et l'inverse de $ X ^ TX $ n'est pas nécessaire. De plus, l'existence d'innombrables solutions dans (1.9) est considérée comme correspondant à la valeur initiale de la recherche.

Implémentons une régression linéaire basée sur (1.13).

class LR_GradientDecent(object):
    
    def __init__(self, lr=0.01, n_iter=100, seed=0):
        self.lr = lr # learning rate
        self.n_iter = n_iter # number of iterations
        self.seed = seed # seed for initializing coefficients
        
    def fit(self, X, y):
        np.random.seed(self.seed)
        self.w = np.random.uniform(low=-1.0, high=1.0, size=(X.shape[1], 1)) # initializing coefficients
        self.loss_arr = np.array([])
        for _ in range(self.n_iter):
            self.loss_arr = np.append(self.loss_arr, (y-X.dot(self.w)).T.dot(y-X.dot(self.w)).flatten(), axis=0)
            dw = -2*(X.T.dot(y)-X.T.dot(X).dot(self.w)) # derivatives of the loss function
            self.w -= self.lr * dw # up-date coefficients
            
    def pred(self, X):
        return X.dot(self.w)

J'essaierai de l'utiliser.

model = LR_GradientDecent(n_iter=1000)
model.fit(X,y)

Tracons le changement de la valeur de la fonction de perte pendant l'entraînement.

plt.plot(model.loss_arr)
plt.xlabel('Number of Iterations')
plt.ylabel('Error')

Regardons également le coefficient de régression.

model.w.flatten()

coefs = []
err = []
for seed in range(200):
    model = LR_GradientDecent(n_iter=1000, seed=seed)
    model.fit(X, y)
    coefs.append(model.w.flatten())
    err.append(model.loss_arr[-1])
df_sol = pd.DataFrame(coefs)
df_sol.columns = ['a', 'b', 'c']
df_sol['error'] = err

df_sol

On peut voir que la valeur de convergence du coefficient de régression est complètement différente selon la façon dont la valeur initiale est prise. C'est la soi-disant dépendance de valeur initiale. Ensuite, tracez le coefficient de régression obtenu $ (a, b, c) $ sur une carte tridimensionnelle.

fig = plt.figure()
ax = Axes3D(fig)
# solutions found by gradient decent
ax.plot(df_sol['a'], df_sol['b'], df_sol['c'], 'o', ms=4, mew=0.5)
# exact solutions
c = np.linspace(0, 1, 100)
a = np.array([i - 1 for i in c])
b = np.array([1-2*i for i in c])
ax.plot(a, b, c, color='yellow', ms=4, mew=0.5, alpha=0.4)
#ax.set_xlim(-1,0)
#ax.set_ylim(-1,1)
#ax.set_zlim(0,1)
#ax.set_xlabel('a')
#ax.set_ylabel('b')
#ax.set_zlabel('c')
ax.view_init(elev=30, azim=45)
plt.show()

Dans la figure ci-dessus, le cercle bleu montre le coefficient de régression obtenu en balançant la valeur initiale, et la ligne jaune clair montre la solution exacte (1,9). En secouant la valeur initiale, vous pouvez voir comment les cercles bleus sont dispersés sur la ligne jaune.

Une dernière chose. En fait, LinearRegression de scikit-learn n'utilise pas la mise à jour Gradient Decent (1.13). Je ne suis pas familier avec les circonstances internes de scikit-learn, mais quand j'ai regardé le code source, j'ai senti qu'un algorithme spécialisé pour résoudre des équations linéaires appelé LSQR était utilisé. Cet algorithme recherche également la valeur minimale en répétant la recherche itérative, mais le coefficient de régression obtenu en touchant la graine n'a pas changé, probablement parce que la valeur initiale a été donnée. En pensant de cette façon, lorsque je place les données dans le package de régression linéaire et que le résultat est renvoyé, j'ai l'impression d'avoir terminé un travail, mais je me demande si j'obtiens une réponse vraiment significative. Dans cet article, je voudrais présenter les successeurs de Gradient Decent d'un point de vue plus général qu'un algorithme spécifique au modèle.

Stochastic Gradient Decent (SGD) Eh bien, dans le chapitre précédent, nous avons présenté Gradient Decent, qui était un traitement par lots. En d'autres termes, toutes les données sont traitées ensemble pour chaque étape. Dans l'exemple précédent, il n'y avait que deux échantillons de données, il n'y avait donc pas de problème, mais à mesure que la quantité de données augmente, le traitement par lots devient difficile. En outre, il peut ne pas être possible d'utiliser toutes les données ensemble, comme dans un flux de données en ligne. Par conséquent, nous avons besoin d'une méthode pour traiter le gradient décent de manière séquentielle, et celles-ci sont appelées Gradient décent stochastique (SGD). En fonction de la taille du morceau à traiter, il peut être traité pour chaque 1 échantillon, ou il peut être traité dans une certaine unité (mini-lot). En général, la fonction de perte $ L $ est la somme de chaque échantillon de données.

L=\sum_{n=1}^NL_n \tag{2.1}

Par conséquent, lors du traitement d'un échantillon par échantillon, l'algorithme de (1.12) est appliqué à chaque $ L_n $ en séquence. Dans le cas du traitement par mini-batch, l'algorithme de mise à jour est utilisé pour la somme $ \ sum_ {n \ in \ text {mini-batch}} L_n $ des échantillons inclus dans le mini-batch. À propos, SGD présente également certains avantages que Gradient Decent n'a pas. Dans Gradient Decent, il y a de fortes chances que vous restiez coincé dans la solution locale et que vous ne puissiez pas en sortir, mais dans SGD, vous serez plus susceptible de sortir de la solution locale en échangeant des données et en secouant le chemin de recherche. Cependant, il semble que le traitement par mini-lot soit préféré car le traitement à 1 échantillon tremblera trop et la convergence du calcul se détériorera.

Le contenu suivant est une histoire générale qui ne dépend pas de savoir s'il faut traiter toutes les données à la fois, 1 échantillon par échantillon ou mini-lot, j'ai donc décidé d'écrire la fonction de perte comme seulement $ L $. Je vais.

Successeur de SGD

SGD est un algorithme très simple, mais les performances en souffrent à mesure que la fonction de perte devient plus complexe. Par conséquent, divers algorithmes ont été proposés qui sont des améliorations de SGD. À partir de là, je voudrais vous présenter une méthode typique.

Momentum SGD Par exemple, considérons la fonction de perte suivante. Tout d'abord, veuillez voir la figure à gauche.

SGD	Momentum SGD

Supposons que la ligne noire représente la ligne de contour de la fonction de perte et que l'étoile violette représente la valeur minimale. Selon la formule de mise à jour de Gradient Decent (1.12), le mouvement de zigzag commence car l'inclinaison est poussée dans la direction du minimum et il est difficile d'atteindre la marque d'étoile cible. Cependant, si vous regardez de près, il semble que si vous prenez la "moyenne" de la partie qui va et vient, le mouvement en zigzag sera décalé et vous obtiendrez une force de propulsion vers l'avant (flèche bleue sur la figure de droite). C'est l'idée de Momentum SGD. Lorsqu'il est écrit dans une formule,

\begin{align}
& m_t = \gamma m_{t-1} + \eta \frac{\partial L}{\partial w}(w_t)　\tag{3.1}\\
& w_{t+1} = w_t - m_t \tag{3.2}
\end{align}

Ce sera. Définissez la valeur initiale de $ m_t $ sur $ m_0 = 0 $. $ M $ in (3.1) sert à compenser le mouvement en zigzag dans un terme appelé Momentum. Si vous mettez $ g_t ^ {\ eta}: = \ eta \ frac {\ partial L} {\ partial w} (w_t) $, (3.1) est

m_t=g_t^{\eta} + \gamma g_{t-1}^{\eta} + \gamma^2 g_{t-2}^{\eta} + \gamma^3 g_{t-3}^{\eta} + \cdots \tag{3.3}

Puisqu'il peut être réécrit comme, vous pouvez voir que les informations précédentes sur les garadients sont accumulées au taux d'atténuation $ \ gamma $. Le remplacement de $ m_t $ par $ g_t ^ \ eta $ dans (3.2) reviendra à SGD. Il semble que $ \ gamma $ soit souvent fixé à environ 0,9 $.

Nesterov accelerated gradient（NAG） Nesterov Accelerated Gradient (NAG) est une version légèrement modifiée de Momentum SGD. Dans la figure précédente, l'étoile violette était la destination. Même si je supprime le mouvement en zigzag, je crains qu'il s'arrête correctement à destination. Juste avant votre arrivée, ce serait bien si vous pouviez prédire votre destination et appliquer les freins. Donc, la formule de mise à jour

\begin{align}
& m_t = \gamma m_{t-1} + \eta \frac{\partial L}{\partial w}(w_t - m_{t-1})　\tag{3.4}\\
& w_{t+1} = w_t - m_t \tag{3.5}
\end{align}

Je vais essayer. (3.4) Notez que l'évaluation du gradient sur le côté droit n'est pas le point courant $ w_t $ mais le point de mouvement attendu $ w_t --m_ {t-1} $. Le pas de temps de $ m $ est $ t-1 $ car $ m_t $ lui-même ne peut pas être utilisé pour composer $ m_t $, c'est juste un point de déplacement "prédit".

AdaGrad AdaGrad est basé sur SGD, mais conçoit un taux d'apprentissage différent pour chaque composant du paramètre d'apprentissage $ w $ et un taux d'apprentissage élevé pour les composants moins mis à jour. Ci-dessous, par souci de simplification des symboles

g^t_i:=\frac{\partial L}{\partial w_i}(w^t) \tag{3.6}

Au fait, j'écrirai l'indice de pas de temps $ t $ dans le coin supérieur droit. $ i $ est un index qui spécifie les composants du vecteur $ w $. L'algorithme de mise à jour d'AdaGrad est donné ci-dessous

\begin{align}
& v^t_i = v^{t-1}_i + (g^t_i)^2 \tag{3.7} \\
& w^{t+1}_i = w^t_i - \frac{\eta}{\sqrt{v^t_i} + \epsilon} g^t_i　\tag{3.8}
\end{align}

Définissez la valeur initiale de $ v $ sur $ v_i ^ 0 = 0 $. $ \ eta $ est une constante donnée manuellement, également appelée taux d'apprentissage. $ \ Epsilon $ définit une petite valeur pour empêcher la division par $ 0 $. De plus, $ v ^ t_i $ provient de (3.7)

v^t_i = (g^t_i)^2 + (g^{t-2}_i)^2 + (g^{t-3}_i)^2 + \cdots \tag{3.9}

Par conséquent, la valeur de $ v ^ t_i $ est petite pour les paramètres d'apprentissage qui ne sont pas mis à jour très souvent, et le taux d'apprentissage réel $ \ eta / (\ sqrt {v ^ t_i} + \ epsilon) $ est grand. Devenir. Si vous n'ajoutez pas la racine carrée "$ \ sqrt {\ cdot} $" attachée à $ v ^ t_i $, les performances se détérioreront.

RMSprop Dans AdaGrad, comme indiqué dans (3.9), le taux d'apprentissage réel diminue à mesure que la mise à jour est répétée et la mise à jour de la valeur s'arrête. Le but de RMSprop est d'empêcher cela. Plutôt que d'accumuler toutes les informations passées, nous accumulerons les accumulations passées et les dernières informations dans un rapport de $ \ gamma: 1- \ gamma $. Autrement dit, l'algorithme de mise à jour ressemble à ceci:

\begin{align}
& v^t_i = \gamma v^{t-1}_i + (1-\gamma)(g^t_i)^2 \tag{3.10} \\
& w^{t+1}_i = w^t_i - \frac{\eta}{\sqrt{v^t_i}+\epsilon} g^t_i　\tag{3.11}
\end{align}

Définissez la valeur initiale de $ v $ sur $ v_i ^ 0 = 0 $. Si vous réécrivez (3.10)

v^t_i=\gamma^t v^0_i + (1-\gamma)\sum_{l=1}^t \gamma^{t-l}(g^l_i)^2 \tag{3.12}

Par conséquent, il peut être confirmé que les informations passées sont accumulées tout en diminuant de manière exponentielle, et il s'agit d'une sorte de méthode de moyenne appelée moyenne d'index mobile.

AdaDelta AdaDelta semble avoir été proposé indépendamment de RMSprop, mais il s'agit fonctionnellement d'une autre modification de RMSprop. Le but est d'augmenter la distance de recherche lorsque le paramètre d'apprentissage entre dans la zone non recherchée. L'algorithme de mise à jour ressemble à ceci:

\begin{align}
& v^t_i = \gamma v^{t-1}_i + (1-\gamma)(g^t_i)^2 \tag{3.13} \\
& s^t_i = \gamma s^{t-1}_i + (1-\gamma)(\Delta w^{t-1}_i)^2 \tag{3.14} \\
& \Delta w^t_i = -\frac{\sqrt{s^t_i + \epsilon}}{\sqrt{v^t_i + \epsilon}} g^t_i \tag{3.15} \\
& w^{t+1}_i = w^t_i + \Delta w^t_i　\tag{3.16}
\end{align}

Définissez les valeurs initiales comme $ v_i ^ 0 = 0 $ et $ s_i ^ 0 = 0 $. En regardant (3.15), $ s ^ t $ est inclus au lieu de $ \ eta $. Si $ w ^ t_i $ est mis à jour de manière significative, il entrera dans une zone non recherchée, mais (3.14) met également à jour $ s ^ t_i $ de manière significative, il est donc conçu pour augmenter la distance de recherche. Je vois parfois des expressions qui incluent le paramètre de taux d'apprentissage.

Adam Adam ressemble à une combinaison de RMSprop et de la méthode Momentum. L'algorithme de mise à jour est donné comme suit

\begin{align}
& m^t_i = \beta_1 m^{t-1}_i + (1-\beta_1)g^t_i \tag{3.17} \\
& v^t_i = \beta_2 v^{t-1}_i + (1-\beta_2)(g^t_i)^2 \tag{3.18} \\
& \hat{m}^t_i = \frac{m^t_i}{1-\beta_1^t} \tag{3.19} \\
& \hat{v}^t_i = \frac{v^t_i}{1-\beta_2^t} \tag{3.20} \\
& w^{t+1}_i = w^t_i - \frac{\eta}{\sqrt{\hat{v}^t_i}+\epsilon} \hat{m}^t_i　\tag{3.21}
\end{align}

Définissez les valeurs initiales comme $ v_i ^ 0 = 0 $ et $ m_i ^ 0 = 0 $. (3,17)

m^t_i=\beta_1^tm_0 + (1-\beta_1)\sum_{k=1}^t \beta_1^{t-k}g^k_i \tag{3.22}

Après la réécriture, remplacez la valeur initiale $ m_0 = 0 $, et comparez les commandes des deux côtés avec l'ordre de $ g ^ k_i $ comme $ g ^ k_i \ sim O [g] $.

O[m] \sim O[g](1-\beta_1) \sum_{k=1}^t \beta_1^{t-k} = O[g](1-\beta_1^t) \tag{3.23}

Ce sera. En premier lieu, $ m ^ t $ a été introduit avec l'espoir que ce soit la "moyenne" de $ g ^ t $, il semble donc désagréable qu'un écart de $ 1- \ beta_1 ^ t $ se produise. est. En particulier, quand il vaut $ \ beta_1 ^ t \ sim 0 $, il devient $ O [m] \ sim 0 $. (3.19) est la correction pour cela. La correction dans (3.20) est pour la même raison.

Eve Eve est une méthode proposée comme compatibilité ascendante avec Adam vers novembre 2016. Il a été souligné que lorsque la fonction de perte devient compliquée, en particulier dans l'apprentissage profond, il est plus probable que les zones plates appelées selles et plateau gênent l'apprentissage plutôt que de rester coincées dans la valeur minimale locale. Je vais. Cependant, étant donné que SGD n'utilise que les informations sur la différenciation de premier ordre de la fonction de perte, certains aspects ne peuvent pas être traités. Par conséquent, Eve utilise directement la fluctuation de la fonction de perte qui accompagne la mise à jour des paramètres d'apprentissage. Lorsque la quantité de fluctuation est faible, il est fort probable que vous soyez coincé dans le point de selle / plateau et la distance de recherche des paramètres d'apprentissage est étendue. Les autres étapes sont les mêmes que celles d'Adam et l'algorithme de mise à jour est le suivant.

\begin{align}
& m^t_i = \beta_1 m^{t-1}_i + (1-\beta_1)g^t_i \tag{3.24} \\
& v^t_i = \beta_2 v^{t-1}_i + (1-\beta_2)(g^t_i)^2 \tag{3.25} \\
& \hat{m}^t_i = \frac{m^t_i}{1-\beta_1^t} \tag{3.26} \\
& \hat{v}^t_i = \frac{v^t_i}{1-\beta_2^t} \tag{3.27} \\
& r^t = (\text{feed-back from the loss function}) \tag{3.28} \\
& d^t = \beta_3 d_{t-1} + (1-\beta_3) r^t \tag{3.29} \\
& w^{t+1}_i = w^t_i - \frac{\eta}{d^t \sqrt{\hat{v}^t_i}+\epsilon} \hat{m}^t_i　\tag{3.30}
\end{align}

Le point est de savoir comment configurer (3.28). (3,29) est la moyenne mobile de l'indice du retour d'expérience passé et du dernier retour d'information, et est ajoutée au dénominateur de (3,30). Définissez la valeur initiale de $ d $ sur $ d ^ 0 = 1 $. Maintenant, passons à l'explication de (3.28). Tout d'abord, définissez $ L ^ t: = L (w ^ t) $ pour la fonction de perte $ L (w) $. Fondamentalement, le montant de la variation non négative de $ L $ $ \ frac {L ^ {t} -L ^ {t-1}} {L ^ {t-1}} $ ($ L ^ t> L ^ {t -1} $), $ \ frac {L ^ {t-1} -L ^ {t}} {L ^ {t}} $ ($ L ^ {t-1}> L ^ t $ ) Est utilisé comme $ r ^ t $, mais le calcul peut ne pas être stable en raison du saut de valeurs, alors définissez un ensemble approprié de constantes $ 0 <k <K $ et $ k <r ^ Il sera arrondi pour que t <K $. Cependant, selon l'article, cette méthode simple ne semble pas fonctionner et Eve prend les mesures suivantes: Premièrement, pour $ L ^ 0 = L (w ^ 0) $, $ [\ frac {L ^ 0} {K + 1}, \ frac {L ^ 0} {k + 1}] $ (dans la figure) ②), $ [(k + 1) L ^ 0, (K + 1) L ^ 0] $ (⑤ sur la figure). De plus, coupez la section avec $ L ^ 0 $ comme limite comme indiqué sur la figure.

\frac{\partial f}{\partial x}(x_0+\Delta x)\simeq \frac{\partial f}{\partial x}(x_0)+ \frac{\partial f}{\partial x \partial x}(x_0) \Delta x \tag{4.1}

\Delta x = - \Big[\frac{\partial f}{\partial x \partial x}(x_0)\Big]^{-1}\frac{\partial f}{\partial x}(x_0) \tag{4.2}

x_{t+1} = x_t -\Big[\frac{\partial f}{\partial x \partial x}(x_t)\Big]^{-1}\frac{\partial f}{\partial x}(x_t) \tag{4.3}