Implementação do algoritmo Levenberg-Marquardt para otimizar redes neurais no TensorFlow

Este é um tutorial da biblioteca TensorFlow. Considere um pouco mais profundo do que nos artigos sobre reconhecimento de números manuscritos. Este é um tutorial sobre métodos de otimização. Aqui você não pode prescindir da matemática. Tudo bem se você esqueceu completamente. Lembre-se. Não haverá evidência formal e conclusões complexas, apenas o mínimo necessário para a compreensão intuitiva. Para começar, um pequeno histórico de como esse algoritmo pode ser útil para otimizar uma rede neural.

Seis meses atrás, um amigo me pediu para mostrar como criar uma rede neural em Python. Sua empresa produz instrumentos para medições geofísicas. Várias sondas diferentes durante a perfuração medem um conjunto de sinais associados aos parâmetros do ambiente ao redor do poço. Em alguns casos complexos, calcule com precisão os parâmetros ambientais dos sinais por um longo tempo, mesmo em um computador potente, e é necessário interpretar os resultados das medições em campo. Havia uma idéia de contar centenas de milhares de casos em um cluster e treinar uma rede neural neles. Como a rede neural é muito rápida, ela pode ser usada para determinar parâmetros consistentes com os sinais medidos, exatamente no processo de perfuração. Os detalhes estão no artigo:

Kushnir, D., Velker, N., Bondarenko, A., Dyatlov, G., & Dashevsky, Y. (2018, 29 de outubro). Simulação em tempo real da ferramenta de resistividade azimutal profunda no modelo de falhas 2D usando redes neurais (russo). Sociedade de Engenheiros de Petróleo. doi: 10.2118 / 192573-RU

Uma noite, mostrei como as keras poderiam implementar uma rede neural simples, e um amigo no trabalho começou a treinar nos dados contados. Depois de alguns dias, discutimos o resultado. Do meu ponto de vista, ele parecia promissor, mas um amigo disse que precisava de cálculos com a precisão do dispositivo. E se o erro quadrático médio fosse em torno de 1, seria necessário 1e-3. 3 pedidos menos. Mil vezes.

Experimentos com arquitetura de rede neural, normalização de dados e abordagens de otimização renderam quase nada. Depois de algumas semanas, um amigo ligou e disse que instalou o MatLab e resolveu o problema pelo método Levenberg-Marquardt (a seguir, chamaremos LM ). Foi otimizado por um longo tempo (vários dias), não funcionou na GPU, mas o resultado foi o correto. Pareceu um desafio.

A pesquisa rápida de um otimizador LM pronto para keras ou TensorFlow falhou. Encontrei apenas a biblioteca pyrenn, mas sua funcionalidade me pareceu fraca. Eu decidi implementá-lo eu mesmo. À primeira vista, tudo parecia simples, e duas noites deveriam ter sido suficientes. Demorou mais tempo. Havia dois problemas:

1. TensorFlow. Um monte de artigos, mas quase todos os níveis ", mas vamos escrever um olá mundo reconhecimento de dígitos manuscritos".
2. Matemática Eu esqueci muito, e os autores de artigos matemáticos não se importam com pessoas como eu: fórmulas sólidas sem explicação, “obviamente!” e assim por diante.

Como resultado, ele escreveu um artigo para aqueles que esqueceram a matemática e querem entender o TensorFlow um pouco mais fundo, mas sem hardcore. O artigo tem muito texto e pouco código. A opção oposta, quando há pouco texto e muito código, é aqui o Jupyter Notebook Levenberg-Marquardt .

Conheça o recurso Rosenbrock

Geraremos dados de treinamento pela função Rosenbrock , que geralmente é usada como referência para algoritmos de otimização:

$$ $$f (x, y) = (a-x) ^ 2 + b (y-x ^ 2) ^ 2$$

Por que ela é boa?

• Horário bonito. É chamado Vale de Rosenbrock e a função banana intraduzível de Rosenbrock .
• O mínimo global está dentro de um vale longo, estreito e parabólico. Encontrar um vale é trivial e um mínimo global é difícil.
• Existe uma opção multidimensional. Não é tão fácil criar uma boa função para muitas variáveis.

Começaremos a escrever o código dele, conectando as bibliotecas necessárias para mais trabalhos:

import numpy as np import tensorflow as tf import math def rosenbrock(x, y, a, b): return (a - x)**2 + b*(y - x**2)**2 

Afirmamos o problema

Como estávamos falando de um dispositivo de medição, vamos continuar usando a analogia. Nosso dispositivo em um mundo fictício pode medir coordenadas $$$(x, y)$ e altura $$$z$ . Os físicos estudaram o mundo e disseram: " Sim, é Rosenbrock! Conhecendo as coordenadas, você pode calcular com precisão a altura, não é necessário medi-la ". Em outras palavras, os cientistas nos deram um modelo $$$z = rosenbrock (x, y, a, b)$ que depende dos parâmetros $$$(a, b)$ . Esses parâmetros, embora constantes em um mundo fictício, são desconhecidos. Eles precisam ser encontrados.

Realizamos uma série de experimentos que deram $$$m$ pontos $$$(x_1, y_1, z_1), (x_2, y_2, z_2), ..., (x_m, y_m, z_m)$ :

 # (2.5, 2.5) -   ,  ,     data_points = np.array([[x, y, rosenbrock(x, y, 2.5, 2.5)] for x in np.arange(-2, 2.1, 2) for y in np.arange(-2, 2.1, 2)]) m = data_points.shape[0] 

A primeira maneira de otimizar é tentar adivinhar os parâmetros. Usamos a biblioteca Numpy:

 x, y = data_points[:, 0], data_points[:, 1] z = data_points[:, 2] #   =5  b=5? a_guess, b_guess = 5., 5. #  -hat   ,   , #   ,   ,   .   #     ^   - # .    hat. z_hat = rosenbrock(x, y, a_guess, b_guess) 

Como entender como estamos errados? Contagem de resíduos - tamanhos de erro. $$$m$ pontos dão $$$m$ resíduos - você precisa de um indicador integral. Nós quadrado cada resíduo em um quadrado e calculamos a média:

$$ $$MSE (a, b) = \ frac {1} {m} \ sum_ {i = 1} ^ {m} (z_ {i} - \ widehat {z_ {i}}) ^ 2$$

Essa medida de proximidade é chamada de erro médio quadrático (a seguir denominado mse ):

 # r - residuals () r = z - z_hat # mse loss = np.mean(r**2) print(loss) 

 [Out]: 3868.2291666666665 

Ao minimizar o mse , resolvemos o problema dos mínimos quadrados ( minimização de quadrados não lineares ):

Pode-se ver que os parâmetros não foram adivinhados.

Formulamos o problema no TensorFlow

O modelo tem a forma $$$z = rosenbrock (x, y, a, b)$ . Trazemos para o formulário $$$y = f (x, p)$ (geralmente matemática escreve $$$\ beta$ em vez de $$$p$ mas os programadores não usam beta). Agora o modelo tem a forma $$$y = rosenbrock (x, p)$ onde $$$y$ - altura $$$x$ É o vetor de coordenadas de dois elementos (componente) e $$$p$ - vetor de parâmetros.

Os programadores costumam pensar em vetores como matrizes unidimensionais. Isso não está totalmente correto. Uma matriz de números é um meio de representar um vetor. Você pode representar um vetor como uma matriz de dimensão $$$N$ matriz bidimensional $$$1 \ vezes N$ e até uma matriz $$$N \ vezes 1$ nos casos em que o fato de o vetor ser um vetor de coluna (por exemplo, multiplicar uma matriz por ele) é importante:

$$ $$\ begin {bmatrix} x_1 \\\ vdots \\ x_N \ end {bmatrix}$$

O TensorFlow usa o conceito de tensor . Um tensor , como uma matriz, pode ser unidimensional (para representar um vetor ), bidimensional (para um vetor de matriz ou coluna ) e qualquer dimensão maior.

 #    ('placeholder' ,    #      ) x = tf.placeholder(tf.float64, shape=[m, 2]) y = tf.placeholder(tf.float64, shape=[m]) #   ('variable' ,    ) #     (5, 5) p = tf.Variable([5., 5.], dtype=tf.float64) #  y_hat = rosenbrock(x[:, 0], x[:, 1], p[0], p[1]) #  r = y - y_hat # mse (mean squared error) loss = tf.reduce_mean(r**2) 

O código TensorFlow não é diferente no código Numpy. O conteúdo é enorme. O código numpy calcula o valor mse. O código TensorFlow não realiza nenhum cálculo, forma um gráfico de fluxo de dados que mse pode calcular . Um momento muito tolerante ao cérebro é o trabalho da função rosenbrock . Usamos nos dois casos. Mas quando passamos as matrizes Numpy, ele executa os cálculos pela fórmula e retorna os números. E quando transferimos os tensores para o TensorFlow, ele forma um subgráfico do fluxo de dados e retorna sua borda na forma de um tensor. Milagres de polimorfismo, mas não os abuse:

Graças à presença de um gráfico de fluxo de dados, o TensorFlow, em particular, pode calcular derivadas automaticamente (usando a técnica de diferenciação automática no modo reverso ).

Um momento de matemática. Blocos "para quem esqueceu" serão escondidos em um spoiler.

Construímos um gráfico de fluxo de dados, vamos executar o cálculo do mse:

 #        #      placeholder (  ) feed_dict = {x: data_points[:,0:2], y: data_points[:,2]} #       TensorFlow session = tf.Session() #     session.run(tf.global_variables_initializer()) #   ()  loss (mse) current_loss = session.run(loss, feed_dict) print(current_loss) 

 [Out]: 3868.2291666666665 

O resultado é o mesmo que com o Numpy. Então eles não estavam enganados.

Comece a otimizar

Infelizmente, não foi possível adivinhar os parâmetros. Mas então nós:

1. Definimos o critério de otimização - o valor mínimo de mse.
2. Os parâmetros variáveis ​​foram determinados: vetor $$$p$ com componentes $$$a$ , $$$b$ Funções Rosenbrock.
3. Ainda não pensamos em limitações, mas elas não estão lá.

Na última etapa, construímos um gráfico de fluxo de dados com um tensor de perda finita ( função de perda ). O objetivo da otimização é encontrar o valor do vetor de parâmetros $$$p$ em que o valor da função de perda é mínimo. Tivemos sorte, o gráfico desta função é muito simples (côncavo e sem mínimos locais):

Introdução à otimização. Para começar, escrevemos um ciclo generalizado:

 # :   mse,   ,  #   mse,        placeholder def train(target_loss, max_steps, loss_tensor, train_step_op, inputs): step = 0 current_loss = session.run(loss_tensor, inputs) #           while current_loss > target_loss and step < max_steps: step += 1 #    1, 2, 4, 8, 16...  if math.log(step, 2).is_integer(): print(f'step: {step}, current loss: {current_loss}') #    session.run(train_step_op, inputs) current_loss = session.run(loss_tensor, inputs) print(f'ENDED ON STEP: {step}, FINAL LOSS: {current_loss}') 

Otimizamos pelo método da descida de gradiente mais rápida (SGD)

As ações desse método podem ser comparadas ao montar um esquiador ousado, que sempre desce a ladeira (na direção mais íngreme). Nesse caso, apenas a inclinação no ponto de localização é levada em consideração. E se a inclinação for forte, o esquiador voa uma longa distância antes da próxima mudança. Com uma inclinação fraca, ele se move em pequenos passos. Talvez como voar para longe em uma árvore (o algoritmo diverge ) e fica preso em um poço ( mínimo local ).

Você pode escrever da seguinte forma (alterar $$$\ boldsymbol {p}$ em $$$\ boldsymbol {p} - ...$ ):

$$ $$\ boldsymbol {p} \ rightarrow \ boldsymbol {p} - \ alpha [\ nabla_ {p} perda (\ boldsymbol {p})]$$

Gorduroso $$$\ boldsymbol {p}$ enfatiza que este é o ponto da localização real - o valor do vetor de parâmetro na etapa atual. No primeiro passo, esse é o nosso palpite (5, 5). Existem dois pontos interessantes na fórmula: $$$\ alpha$ - taxa de aprendizagem ( taxa de aprendizagem ), $$$\ nabla_ {p} perda$ - gradiente ( gradiente ) da função de perda pelo vetor de parâmetros.

 #        grad = tf.gradients(loss, p)[0] #   learning_rate = 0.0005 #   ,     apply_gradients - #         opt = tf.train.GradientDescentOptimizer(learning_rate=1) #           sgd = opt.apply_gradients([(learning_rate*grad, p)]) #   ,    40000  session.run(tf.global_variables_initializer()) train(1e-10, 40000, loss, sgd, feed_dict) print('PARAMETERS:', session.run(p)) 

 [Out]: step: 1, current loss: 3868.2291666666665 step: 2, current loss: 1381.5379689135807 [...] ENDED ON STEP: 582, FINAL LOSS: 9.698531012270816e-11 PARAMETERS: [2.50000205 2.49999959] 

Foram necessárias 582 etapas:

Otimizamos pelo método Adam

Não vamos avançar nos métodos de gradiente, mas há muitas variações. Você pode ler sobre eles no artigo Métodos de otimização de redes neurais . No TensorFlow, muitos otimizadores já estão implementados. Por exemplo, Adam:

 #       , #      adm = tf.train.AdamOptimizer(15).minimize(loss) #   ,    40000  session.run(tf.global_variables_initializer()) train(1e-10, 40000, loss, adm, feed_dict) print('PARAMETERS:', session.run(p)) 

 [Out]: step: 1, current loss: 3868.2291666666665 step: 2, current loss: 34205.72916492336 [...] ENDED ON STEP: 317, FINAL LOSS: 2.424142714263483e-12 PARAMETERS: [2.49999969 2.50000008] 

Gerenciado em 317 etapas. Muito mais rapido

Otimizamos pelo método de Newton

As ações dos métodos de segunda ordem podem ser comparadas à pilotagem de um snowboardista de freeride racional que pondera o próximo ponto de sua rota por um longo tempo e leva em consideração não apenas a inclinação do local, mas também a curvatura.

De fato, os métodos de descida de gradiente e os de segunda ordem tentam adivinhar ( aproximar ) a função no ponto atual. Os métodos de gradiente focam apenas a inclinação do gráfico da função no ponto - a primeira derivada. Os métodos de segunda ordem, além do viés, levam em consideração a curvatura , a segunda derivada: "se a curvatura persistir, onde será o mínimo?" Nós calculamos e vamos lá:

Para construir essa aproximação e calcular o ponto mínimo estimado, você pode usar a série de Taylor . Para o caso unidimensional, a aproximação por um polinômio de segunda ordem no ponto $$$a$ é assim:

$$ $$f (x) \ approx f (a) + \ frac {f '(a) (x-a)} {1!} + \ frac {f' '(a) (x-a) ^ 2} {2!}$$

O mínimo é atingido em $$$x = a - \ frac {f '(a)} {f' '(a)}$ . O caso multidimensional parece mais sério:

Aproximação de um polinômio de segunda ordem para uma função de um vetor por meio de um gradiente e uma matriz Hessiana em um ponto $$$a$ é assim:

$$ $$f (x) \ aproximadamente f (a) + (xa) ^ \ intercal [\ nabla_ {x} f (a)] + \ frac {1} {2!} (xa) ^ \ intercal [\ boldsymbol {H } f_ {x} (a)] (xa)$$

O mínimo é atingido em $$$x = a - [\ boldsymbol {H} f_ {x} (a)] ^ {- 1} [\ nabla_ {x} f (a)]$ . A forma praticamente coincide com o caso unidimensional: substituímos a primeira derivada por um gradiente, a segunda por uma matriz Hessiana e fizemos uma correção para trabalhar com vetores. É impossível dividir um vetor por uma matriz, portanto, a multiplicação pela matriz inversa é usada. T significa transpor . A fórmula implica que, por padrão, um vetor é uma coluna. Transpose transforma um vetor de coluna em um vetor de linha . Ao implementar no TensorFlow, isso deve ser levado em consideração, mas na direção oposta: por padrão, o vetor é uma sequência (tensor unidimensional). Apenas no caso: a transposição não é uma rotação de 90 graus, é a transformação de linhas em colunas na mesma ordem.

Portanto, a etapa do método Newton tem a seguinte forma:

$$ $$\ boldsymbol {p} \ rightarrow \ boldsymbol {p} - [\ boldsymbol {H} loss_ {p} (\ boldsymbol {p})] ^ {- 1} [\ nabla_ {p} perda (\ boldsymbol {p} )]$$

O TensorFlow tem tudo para implementar este método:

 #        hess = tf.hessians(loss, p)[0] #    - grad_col = tf.expand_dims(grad, -1) # ,      dp = tf.matmul(tf.linalg.inv(hess), grad_col) #  -  - dp = tf.squeeze(dp) #  p  dp    newton = opt.apply_gradients([(dp, p)]) #   ,    40000  session.run(tf.global_variables_initializer()) train(1e-10, 40000, loss, newton, feed_dict) print('PARAMETERS:', session.run(p)) 

 [Out]: step: 1, current loss: 3868.2291666666665 step: 2, current loss: 105.04357496954218 step: 4, current loss: 9.96663526704236 ENDED ON STEP: 6, FINAL LOSS: 5.882202372519996e-20 PARAMETERS: [2.5 2.5] 

Bastante 6 etapas:

Otimizado pelo algoritmo de Gauss-Newton

O método de Newton tem uma desvantagem - a matriz hessiana. Graças ao TensorFlow, podemos contá-lo em uma linha de código. Segundo o wiki, Johann Karl Friedrich Gauss fez a primeira menção de seu método em 1809. O cálculo da matriz de Hessian para vários parâmetros para o método dos mínimos quadrados pode levar muito tempo. Agora podemos assumir que o algoritmo de Gauss-Newton usa a aproximação da matriz de Hessian através da matriz de Jacobi para simplificar os cálculos. Mas, do ponto de vista da história, não é assim: Ludwig Otto Hesse (que desenvolveu a matriz com o nome dele) nasceu em 1811 - 2 anos após a primeira menção do algoritmo. E Carl Gustav Jacobi tinha 5 anos.

O algoritmo de Gauss-Newton não funciona com a função de perda. Funciona com a função residual $$$r (p)$ . Esta função aceita um vetor de entrada de parâmetros $$$p$ e retorna um vetor de resíduos . No nosso caso, o vetor $$$p$ consiste em 2 componentes (parâmetros $$$a$ e $$$b$ Funções de Rosenbrock) e o vetor residual de $$$m$ componente (de acordo com o número de experiências). A função vetorial do argumento vetorial é obtida. Sua derivada:

Para não nos confundirmos com a abundância de caracteres, assumimos que $$$\ boldsymbol {J} _ {r}$ - matriz de Jacobi de funções residuais no ponto atual $$$\ boldsymbol {p}$ . Então o algoritmo de Gauss-Newton pode ser escrito da seguinte maneira:

$$ $$\ boldsymbol {p} \ rightarrow \ boldsymbol {p} - [\ boldsymbol {J} _ {r} ^ \ intercal \ boldsymbol {J} _ {r}] ^ {- 1} \ boldsymbol {J} _ {r } ^ \ intercal r (\ boldsymbol {p})$$

A gravação na forma coincide completamente com a gravação do método de Newton. Somente em vez da matriz Hessiana é usada $$$\ boldsymbol {J} _ {r} ^ \ intercal \ boldsymbol {J} _ {r}$ em vez do gradiente $$$\ boldsymbol {J} _ {r} ^ \ intercal r (\ boldsymbol {p})$ . A seguir, veremos por que essa aproximação pode ser usada. Enquanto isso, vamos prosseguir com a implementação no TensorFlow:

 #  ,  TensorFlow     , #   ,        #  .  ,   : # 1)       tf.unstack(r) # 2)      tf.gradients(r_i, p) # 3)       tf.stack #      ,     #       j = tf.stack([tf.gradients(r_i, p)[0] for r_i in tf.unstack(r)]) jT = tf.transpose(j) #     - r_col = tf.expand_dims(r, -1) #      hess_approx = tf.matmul(jT, j) grad_approx = tf.matmul(jT, r_col) # ,      dp = tf.matmul(tf.linalg.inv(hess_approx), grad_approx) #  -  - dp = tf.squeeze(dp) #  p  dp    ng = opt.apply_gradients([(dp, p)]) #   ,    40000  session.run(tf.global_variables_initializer()) train(1e-10, 40000, loss, ng, feed_dict) 

 [Out]: step: 1, current loss: 3868.2291666666665 step: 2, current loss: 14.653025157673625 step: 4, current loss: 4.3918079172783016e-07 ENDED ON STEP: 4, FINAL LOSS: 3.374364957618591e-17 PARAMETERS: [2.5 2.5] 

Suficiente 4 etapas. Menor que para o método de Newton.

Como pode ser visto no código, a função de perda não é usada na otimização, apenas nos critérios de parada e registro. Como o algoritmo de otimização sabe qual função minimizar? A resposta é surpreendente: de jeito nenhum! Gauss-Newton minimiza apenas o erro quadrático médio .

Corrija a parte matemática do artigo

Repetimos toda a matemática que precisávamos. Vamos corrigi-lo um pouco para focar ainda mais a programação e o TensorFlow. Você pode precisar de um lápis para rastrear a sequência de ações matemáticas.

Existe um modelo $$$y = f (x, p)$ onde $$$x$ - vetor $$$p$ - vetor de parâmetros de dimensão $$$n$ e $$$y$ - escalar. Das experiências recebidas $$$m$ pontos $$$(x_ {1}, y_ {1}), ..., (x_ {m}, y_ {m})$ ( pares de dados ). A função residual do vetor depende apenas do vetor de parâmetro: $$$r (p) = (r_ {1} (p), ... r_ {m} (p))$ onde $$$r_ {k} (p) = y_ {k} - \ widehat {y_ {k}} = y_ {k} - f (x_ {k}, p)$ . , $$$p$ , $$$x_{k}, y_{k}$ ? , $$$x_{k}, y_{k}$ , .

$$$p$ , ( sum of squared error — sse residual sum-of-squares — rss ) . mse sse , $$$m$ . . :

$$ $$loss(p) = r_{1}^2(p) + \cdots + r_{m}^2(p) = \sum_{k=1}^{m} r_{k}^2(p)$$

$$$p$ $$$(p)$ .

, . — . — , $$$r^2$ $$$2r \frac{\partial r}{\partial p}$ . :

$$ $$\nabla_{p}loss = (\sum_{k=1}^{m}2r_{k}\frac{\partial r_{k}}{\partial p_{1}}, \cdots, \sum_{k=1}^{m}2r_{k}\frac{\partial r_{k}}{\partial p_{n}})$$

. :

$$ $$[\boldsymbol{H}loss_{p}]_{ij} = \frac{\partial^2 loss}{\partial p_{i} \partial p_{j}} = \sum_{k=1}^{m}(2\frac{\partial r_{k}}{\partial p_{i}}\frac{\partial r_{k}}{\partial p_{j}} + 2r_{k}\frac{\partial^2 r_{k}}{\partial p_{i}\partial p_{j}})$$

. , , $$${(uv)}'={u}'v+u{v}'$ .
Ótimo! .

, , , — $$$2r_{k}\frac{\partial^2 r_{k}}{\partial p_{i}\partial p_{j}}$ . , , $$$r_{k}$ , . — . , ? -.

:

$$ $$\boldsymbol{J}_{r} = \begin{pmatrix} \frac{\partial r_{1}}{\partial p_{1}} & \cdots & \frac{\partial r_{1}}{\partial p_{n}}\\ \vdots & \ddots & \vdots \\ \frac{\partial r_{m}}{\partial p_{1}} & \cdots & \frac{\partial p_{m}}{\partial p_{n}} \end{pmatrix}$$

, , . Note que:

$$ $$2\boldsymbol{J}_{r}^\intercal\boldsymbol{J}_{r} \approx \boldsymbol{H}loss_{p}$$

"" . ( ). , — $$$2r_{k}\frac{\partial^2 r_{k}}{\partial p_{i}\partial p_{j}}$, .
( ):

$$ $$2\boldsymbol{J}_{r}^\intercal r = \nabla_{p}loss$$

, , - — , mse .

. , , . $$$m$$$$(x_{1}, y_{1}), ..., (x_{m}, y_{m})$, $$$y = rosenbrock(x, p)$. $$$p$, .

, : " . - ! ". , , , ( supervised learning ). , . : ( training set ) — ; — ( prediction model ) ; — , .

( multi-layer perceptron neural network mlp ). , , :

• ( starting values ) . Xavier'a, .
• ( overfitting ). — . , . — .
• ( scaling of the input ). , .

9 . 500:

 #    def get_random_rosenbrock_data_points(m): result = np.zeros((m, 3)) result[:, 0] = np.random.uniform(-2, 2, m) result[:, 1] = np.random.uniform(-2, 2, m) result[:, 2] = rosenbrock(result[:, 0], result[:, 1], 2.5, 2.5) return result m = 500 data_points = get_random_rosenbrock_data_points(m) # overfitting   ,      validation_data_points = get_random_rosenbrock_data_points(m) 

500 . — ( learner ), ( outcome measurement ) ( features ) .

( network diagram ). MatLab:

( input ). $$$W$( weights ) 2x10, $$$b$( bias ) 10, ( activation ). () ( hidden layer ) 10 . , , ( output ).

, , ( $$$tanh$):

$$ $$\begin{matrix} h_{1} = tanh(xW_{1} + b_{1})\\ \widehat{y} = h_{1}W_{2} + b_{2} \end{matrix}$$

:

$$ $$h_1 = tanh(\begin{bmatrix} x_1 & x_2 \end{bmatrix}\begin{bmatrix} w^{(1)}_{1,1} & \cdots& w^{(1)}_{1,10} \\ w^{(1)}_{2,1} &\cdots& w^{(1)}_{2,10} \end{bmatrix} + \begin{bmatrix} b^{(1)}_1 & \cdots & b^{(1)}_{10} \end{bmatrix}) \\ \widehat{y} = \begin{bmatrix}h^{(1)}_1 & \cdots & h^{(1)}_{10}\end{bmatrix} \begin{bmatrix} w^{(2)}_{1,1} \\ \vdots \\ w^{(2)}_{1,10} \\ \end{bmatrix} + b_2$$

. $$$W_{1}$"" $$$h_{1}$, - $$$W_{2}$. 41 . , .

$$$m \times 2$, . - $$$\widehat{y}$de $$$m$:

 #     10 "" n_hidden = 10 #      Xavier'a initializer = tf.contrib.layers.xavier_initializer() #    x = tf.placeholder(tf.float64, shape=[m, 2]) y = tf.placeholder(tf.float64, shape=[m, 1]) #         W1 = tf.Variable(initializer([2, n_hidden], dtype=tf.float64)) b1 = tf.Variable(initializer([1, n_hidden], dtype=tf.float64)) #   ,  tanh   h1 = tf.nn.tanh(tf.matmul(x, W1) + b1) #        W2 = tf.Variable(initializer([n_hidden, 1], dtype=tf.float64)) b2 = tf.Variable(initializer([1], dtype=tf.float64)) #   y_hat = tf.matmul(h1, W2) + b2 #  r = y - y_hat #   mse     loss = tf.reduce_mean(tf.square(r)) #      placeholder feed_dict = {x: data_points[:,0:2], y: data_points[:,2:3]} validation_feed_dict = {x: validation_data_points[:,0:2], y: validation_data_points[:,2:3]} 

Adam

Adam $$$rosenbrock$. mse :

 #    adm = tf.train.AdamOptimizer(1e-2).minimize(loss) session.run(tf.global_variables_initializer()) #   ,    40000  train(1e-10, 40000, loss, adm, feed_dict) print('VALIDATION LOSS: '+str(session.run(loss, validation_feed_dict))) 

 [Out]: step: 1, current loss: 671.4242576535694 [...] ENDED ON STEP: 40000, FINAL LOSS: 0.22862158574440725 VALIDATION LOSS: 0.29000289644978866 

. : , , .

$$$rosenbrock$2 . :

• . 9 , 500. .
• . - $$$p$, .

:

 #      y   x def jacobian(y, x): loop_vars = [ tf.constant(0, tf.int32), tf.TensorArray(tf.float64, size=m), ] #  -   #      _, jacobian = tf.while_loop( lambda i, _: i < m, #           #   (-),   x     lambda i, res: (i+1, res.write(i, tf.reshape(tf.gradients(y[i], x), (-1,)))), loop_vars) #       return jacobian.stack() #       r_flat = tf.squeeze(r) #        #       parms = [W1, b1, W2, b2] parms_sizes = [tf.size(p) for p in parms] j = tf.concat([jacobian(r_flat, p) for p in parms], 1) jT = tf.transpose(j) #           hess_approx = tf.matmul(jT, j) grad_approx = tf.matmul(jT, r) 

$$$\boldsymbol{J}r_{p}$. , 4 $$$W_1, b_1, W_2, b_2$. 4 $$$\boldsymbol{J}r_{W_1}, \boldsymbol{J}r_{b_1}, \boldsymbol{J}r_{W_2}, \boldsymbol{J}r_{b_2}$tf.concat .

. tf.while_loop , $$$r_i$, , stack .

$$$r_i$$$$W_1$: $$$\begin{bmatrix} \frac{\partial r_i}{\partial w^{(1)}_{1,1}} & \cdots& \frac{\partial r_i}{\partial w^{(1)}_{1,10}} \\ \frac{\partial r_i}{\partial w^{(1)}_{2,1}} & \cdots& \frac{\partial r_i}{\partial w^{(1)}_{2,10}} \end{bmatrix}$. tf.reshape (-1,) $$$\begin{bmatrix} \frac{\partial r_i}{\partial w^{(1)}_{1,1}} & \cdots& \frac{\partial r_i}{\partial w^{(1)}_{1,10}} & \frac{\partial r_i}{\partial w^{(1)}_{2,1}} & \cdots& \frac{\partial r_i}{\partial w^{(1)}_{2,10}} \end{bmatrix}$.

. - . — TensorFlow . — - - $$$W_1, b_1, W_2, b_2$. -. Levenberg-Marquardt Jupyter Notebook rosenbrock_train.py . , TensorFlow . - , ( ) , , .

-

hess_approx grad_approx -. $$$rosenbrock$, . :

1. : $$$\Delta \boldsymbol{p} = \begin{bmatrix}\Delta w^{(1)}_{1,1} & \cdots & \Delta w^{(1)}_{2,10} & \Delta b^{(1)}_1 & \cdots & \Delta b^{(1)}_{10} & \Delta w^{(2)}_{1,1} & \cdots & \Delta w^{(2)}_{1,10} & \Delta b_2\end{bmatrix}$
2. :
   $$$\Delta W_{1} = \begin{bmatrix}\Delta w^{(1)}_{1,1} & \cdots & \Delta w^{(1)}_{2,10} \end{bmatrix}$, $$$\Delta b_{1} = \begin{bmatrix} \Delta b^{(1)}_1 & \cdots & \Delta b^{(1)}_{10} \end{bmatrix}$, $$$\Delta W_{2} = \begin{bmatrix} \Delta w^{(2)}_{1,1} & \cdots & \Delta w^{(2)}_{1,10} \end{bmatrix}$, $$$\Delta b_{2} = \begin{bmatrix} \Delta b_2\end{bmatrix}$.
3. , :
   $$$\Delta W_{1} = \begin{bmatrix} \Delta w^{(1)}_{1,1} & \cdots & \Delta w^{(1)}_{1,10} \\ \Delta w^{(1)}_{2,1} &\cdots & \Delta w^{(1)}_{2,10} \end{bmatrix}$, $$$\Delta W_{2} = \begin{bmatrix} \Delta w^{(2)}_{1,1} \\ \vdots \\ \Delta w^{(2)}_{1,10} \\ \end{bmatrix}$
4. .

 # 1.     dp_flat = tf.matmul(tf.linalg.inv(hess_approx), grad_approx) # 2.     dps = tf.split(dp_flat, parms_sizes, 0) # 3.     for i in range(len(dps)): dps[i] = tf.reshape(dps[i], parms[i].shape) # 4.  :       gn = opt.apply_gradients(zip(dps, parms)) #   session.run(tf.global_variables_initializer()) train(1e-10, 100, loss, gn, feed_dict) 

 [Out]: step: 1, current loss: 548.8468777701685 step: 2, current loss: 49648941.340197295 InvalidArgumentError: Input is not invertible. 

- . , . - , .

, .

-

. Matlab trainlm . . MathWorks.

- : $$$\boldsymbol{p} \rightarrow \boldsymbol{p}-[\boldsymbol{J}_{r}^\intercal\boldsymbol{J}_{r}]^{-1}\boldsymbol{J}_{r}^\intercal r(\boldsymbol{p})$. - :

$$ $$\boldsymbol{p} \rightarrow \boldsymbol{p}-[\boldsymbol{J}_{r}^\intercal\boldsymbol{J}_{r}+\mu \boldsymbol{I}]^{-1}\boldsymbol{J}_{r}^\intercal r(\boldsymbol{p})$$

$$$\ mu$$$$I$$$$n$( ). $$$\ mu$, -. , . , LM -.

:

 mu = tf.placeholder(tf.float64, shape=[1]) n = tf.add_n(parms_sizes) I = tf.eye(n, dtype=tf.float64) # 1.     dp_flat = tf.matmul(tf.linalg.inv(hess_approx + tf.multiply(mu, I)), grad_approx) # 2.     dps = tf.split(dp_flat, parms_sizes, 0) # 3.     for i in range(len(dps)): dps[i] = tf.reshape(dps[i], parms[i].shape) # 4.  :       lm = opt.apply_gradients(zip(dps, parms)) 

$$$\ mu$? LM - . , . , $$$\ mu$, . — , mse . , :

 #       store = [tf.Variable(tf.zeros(p.shape, dtype=tf.float64)) for p in parms] #  TensorFlow       save_parms = [tf.assign(s, p) for s, p in zip(store, parms)] restore_parms = [tf.assign(p, s) for s, p in zip(store, parms)] #   mu    3. feed_dict[mu] = np.array([3.]) step = 0 session.run(tf.global_variables_initializer()) #    mse current_loss = session.run(loss, feed_dict) #    100   while current_loss > 1e-10 and step < 100: step += 1 #  1, 2, 4...   if math.log(step, 2).is_integer(): print(f'step: {step}, mu: {feed_dict[mu][0]} current loss: {current_loss}') #    session.run(save_parms) # ,     mse while True: #    session.run(lm, feed_dict) new_loss = session.run(loss, feed_dict) if new_loss > current_loss: #  -  mu  10     feed_dict[mu] *= 10 session.run(restore_parms) else: #  -  mu  10     feed_dict[mu] /= 10 current_loss = new_loss break print(f'ENDED ON STEP: {step}, FINAL LOSS: {current_loss}') print('VALIDATION LOSS: '+str(session.run(loss, validation_feed_dict))) 

 [Out]: step: 1, mu: 3.0 current loss: 692.6211687622557 [...] ENDED ON STEP: 100, FINAL LOSS: 0.012346989371823602 VALIDATION LOSS: 0.01859463694102034 

100 LM mse 10 , 40 .

. , . , rosenbrock_train.py .

2D . . . , " " ( curse of dimentionality , Bellman, 1961). . .

:

$$ $$f(\boldsymbol{x}) = \sum_{i=1}^{N-1}\left [ 100(x_{i+1} - x_{i}^2)^2 + (1-x_{i})^2 \right ], \boldsymbol{x}=[x_1 \cdots x_{N}]\in \mathbb{R}^N$$

rosenbrock_train.py get_rand_rosenbrock_points .

-

- : " ! 4 , 300! ". , ( ) -. , , . - . . : ? , . . , - :

1. 10 000 6D .
2. 3 12, 10, 8 (311 ).
3. .
4. 3.5 .

. - 2 . LM . 20 .

rosenbrock_train.py . . , .

Conclusão

, . " ", , . , . , 273 . - , .

, :

1. .
2. ( ) -:
   [1] Petros Drineas, Ravi Kannan, and Michael W. Mahoney. 2006. Fast Monte Carlo Algorithms for Matrices I: Approximating Matrix Multiplication. SIAM J. Comput. 36, 1 (July 2006), 132-157. DOI= http://dx.doi.org/10.1137/S0097539704442684
   [2] Adelman, M., & Silberstein, M. (2018). Faster Neural Network Training with Approximate Tensor Operations. CoRR, abs/1805.08079.

, - . , . "".