Por que model.evaluate() vs. cálculo de perda manual com model.predict() em tf.keras não batem?

Utilizo Keras e Tensorflow para treinar um Perceptron Multicamadas (MLP) "simples" para uma tarefa de regressão, onde utilizo o erro quadrático médio (MSE) como função de perda. Denomino meus dados de treinamento como x_train, y_traine meus dados de teste como x_test, y_test. Reconheci o seguinte: Para Ae Bdefinido da seguinte forma:

1. A = model.evaluate(x_test, y_test)e
2. B = loss(pred_test, y_test), onde pred_test = model.predict(x_test)estão as previsões fora da amostra obtidas do meu modelo,

os valores para Ae Bsão (ligeiramente) diferentes. Minha pergunta é de onde vem a diferença e o que posso fazer para que os valores coincidam. Abaixo, dou um exemplo mínimo reproduzível no qual tentei encontrar a resposta sozinho (sem sucesso). Minha primeira suspeita foi que isso é causado pelo cálculo em lotes; após algumas experiências com os tamanhos de lote, esse não parece ser o caso. Há perguntas relacionadas neste site, mas a resposta a esta pergunta sobre o mesmo (?) problema parece ser específica para CNNs. A discussão neste post afirma que a diferença é causada pela avaliação em lotes em model.evaluate, mas 1.) Eu realmente não vejo como a escolha do tamanho do lote deve afetar o resultado, já que no final a média é construída de qualquer maneira e 2.) mesmo se definir o tamanho do lote para o número de amostras, os resultados ainda são diferentes . Este é até o caso na resposta ao post mencionado . Por fim, há este tópico , onde o problema parece ser causado pela propriedade da métrica de ser, na verdade, uma variante em relação aos tamanhos de lote. No entanto, este não é o caso do MSE!

Aqui está o exemplo mínimo onde treino uma função de regressão em simulações:

import tensorflow as tf
import keras
import numpy as np
import random as random # for sims and seed setting

random.seed(10)

x = np.random.normal([0, 1, 2], [2,1,4], (200, 3))
y = x[:,0] + 0.01 * np.power(x[:,1], 2) + np.sqrt(np.abs(x[:,2] - 3)) + np.random.normal(0, 1, (200))
y = y[:,np.newaxis]

x_train = x[0:100,:]
y_train = y[0:100,:]
x_test = x[101:200,:]
y_test = y[101:200,:]

# MSE
def MSE(a,b):
    return tf.reduce_mean(tf.pow(a - b, 2))

# layers
Inputs_MLP = tf.keras.Input(batch_shape = (100,3), dtype = tf.float32)
Layer1_MLP = tf.keras.layers.Dense(16)(Inputs_MLP)
Outputs_MLP = tf.keras.layers.Dense(1)(Layer1_MLP)

# keras model
model_MLP = tf.keras.Model(Inputs_MLP, Outputs_MLP)
model_MLP.compile(loss = MSE)
history = model_MLP.fit(x = x_train, y = y_train, epochs=5, batch_size = 25)

# evaluation

# out-of-sample
model_MLP.evaluate(x_test, y_test, 100)
# 5.561294078826904
pred_MLP_test = model_MLP.predict(x_test, batch_size = 100)
MSE(pred_MLP_test, y_test)
# <tf.Tensor: shape=(), dtype=float64, numpy=5.561294010797092>

# in-sample
model_MLP.evaluate(x_train, y_train, 100)
# 5.460160732269287
pred_MLP_train = model_MLP.predict(x_train, batch_size = 100)
MSE(pred_MLP_train, y_train)
# <tf.Tensor: shape=(), dtype=float64, numpy=5.46016054713104>

A avaliação fora da amostra resulta em 5,561294078826904 uma vez e, por outro lado, em 5,561294010797092. Para este exemplo, a diferença é pequena, mas ainda me incomoda. Além disso, para outro exemplo (mais longo e complexo), a diferença é maior. Agradeço qualquer ajuda!