An MLP on a stock panel: regularize and compare.

The problem and the yardstick

A multilayer perceptron (MLP) is the simplest network: dense layers with a nonlinearity between them. With enough width it can approximate any continuous function, and that is why it is tempting for tabular data: it should find the interactions and nonlinear shapes that a linear regression does not see. In the level 5 panel those interactions exist (the one-week reversal weakens with volume, the SUE decays with the days, sentiment saturates), so there is something to find. The problem is scale: the oracle, the true alpha, explains around 0.09% of the weekly cross-sectional variance (data README). A network that can learn anything will learn the 99.9% of noise before the 0.1% of signal.

The lesson follows the design of 5.5.3 so that the numbers are comparable: the 21 features without the leak (revision_analistas), transformed into normal ranks by week (lesson 5.2.4), and as the label the normal rank of ret_sig_semana. The splits are by date with a one-week gap between them, which purges the overlap of the weekly labels (lesson 5.4.2):

flowchart TB
  E[Train 2012-2016] --> V[Validation 2017-2018]
  V --> R[Retrain 2012-2018]
  R --> P[Test 2019-2025]

Validation chooses the architecture, the regularization and the number of epochs; the test is looked at only once, at the end, with the chosen model retrained on train and validation. The metric is the IC: the rank correlation, week by week, between the prediction and the next return, averaged. It is what a long-short portfolio uses (lesson 5.1.3).

The networks are trained with the JAX backend, set in the first line before importing Keras: in lesson 5.6.1 the same training of a small model was several times faster than with PyTorch.

import os
os.environ["KERAS_BACKEND"] = "jax"      # before the first import of keras
import numpy as np
import pandas as pd
from scipy.stats import norm
from sklearn.linear_model import Ridge
import lightgbm as lgb
import keras
from keras import layers

f = pd.read_parquet("panel5_features.parquet")
e = pd.read_parquet("panel5_etiquetas.parquet")
d = (f.merge(e[["fecha", "ticker", "ret_sig_semana"]],
             on=["fecha", "ticker"])
      .dropna(subset=["ret_sig_semana"])
      .reset_index(drop=True))
FE = [c for c in f.columns[3:] if c != "revision_analistas"]

def g(s):
    """Normal rank within each week (lesson 5.2.4); missing, 0."""
    por = s.groupby(d["fecha"])
    u = (por.rank() - 0.5) / por.transform("count")
    return pd.Series(norm.ppf(u), index=s.index).fillna(0.0)

X = np.column_stack([g(d[k]) for k in FE]).astype("float32")
y = g(d["ret_sig_semana"]).to_numpy().astype("float32")
sem = d["fecha"].to_numpy()
fe = np.datetime64
ent = sem < fe("2016-12-26")                          # 2012-2016
val = (sem >= fe("2017-01-06")) & (sem < fe("2018-12-24"))
tr_todo = sem < fe("2018-12-24")                      # 2012-2018
te = sem >= fe("2019-01-01")                          # 2019-2025

def ic(score, mask):
    """Rank IC of each week between score and the label."""
    t = pd.DataFrame({"s": score, "f": sem[mask],
                      "y": d.loc[mask, "ret_sig_semana"].to_numpy()})
    t[["s", "y"]] = t.groupby("f")[["s", "y"]].rank()
    return t.groupby("f")[["s", "y"]].corr().xs("s", level=1)["y"]

print("backend:", keras.config.backend(), "| features:", len(FE))
print({k: (int(m.sum()), len(np.unique(sem[m]))) for k, m in
       [("train", ent), ("validation", val), ("test", te)]})

backend: jax | features: 21
{'train': (103600, 259), 'validation': (41200, 103), 'test': (146000, 365)}

There are 103,600 training rows (259 weeks of 400 stocks), 41,200 validation rows and 146,000 test rows.

The baselines

Before the network, the models it has to beat: Ridge ($\alpha = 10$) and the LightGBM of lesson 5.5.3 with 25 trees, the number its validation chose. Both are trained on 2012-2016 for validation and on 2012-2018 for the test. The validation squared error of Ridge and that of always predicting 0 are also measured:

def lgbm(n=25):
    return lgb.LGBMRegressor(
        n_estimators=n, learning_rate=0.02, max_depth=2, num_leaves=4,
        min_child_samples=3000, subsample=0.5, subsample_freq=1,
        colsample_bytree=0.8, reg_lambda=10, random_state=0, n_jobs=4,
        deterministic=True, force_row_wise=True, verbose=-1)

base = {}
for nombre, modelo in [("ridge", lambda: Ridge(alpha=10.0)),
                       ("lightgbm", lgbm)]:
    v_ = ic(modelo().fit(X[ent], y[ent]).predict(X[val]), val)
    p_ = ic(modelo().fit(X[tr_todo], y[tr_todo]).predict(X[te]), te)
    base[nombre] = p_
    print(f"{nombre:<9} validation IC {v_.mean():.4f}   test "
          f"{p_.mean():.4f} ± {p_.std() / np.sqrt(len(p_)):.4f}")
mse_cero = np.mean(y[val] ** 2)
mse_ridge = np.mean((Ridge(alpha=10.0).fit(X[ent], y[ent])
                     .predict(X[val]) - y[val]) ** 2)
print(f"Validation MSE: predict 0 {mse_cero:.4f}, "
      f"Ridge {mse_ridge:.4f}")

ridge     validation IC 0.0155   test 0.0180 ± 0.0045
lightgbm  validation IC 0.0077   test 0.0139 ± 0.0047
Validation MSE: predict 0 0.9968, Ridge 0.9969

Ridge has a validation IC of 0.0155 and a test IC of 0.0180 ± 0.0045; LightGBM, 0.0077 and 0.0139: the same figures as in lesson 5.5.3. The last line is the first warning: Ridge's validation MSE, 0.9969, is worse than that of predicting 0, 0.9968, and yet Ridge ranks the stocks with a positive IC. With an IC of 0.02, the most a prediction can lower the squared error of a label with variance 1 is of the order of $\text{IC}^2 \approx 0.0004$; a minimal mismatch in scale or level hides it. The loss is for training, not for choosing.

A textbook MLP and its overfitting

The recipe of any tutorial: two hidden layers of 128 neurons with ReLU, Adam with rate 0.001, batches of 512 and 20 epochs. A custom callback (lesson 5.6.4) measures the validation IC at the end of each epoch:

class ICValidacion(keras.callbacks.Callback):
    """Measures the validation IC at the end of each epoch."""
    def on_train_begin(self, logs=None):
        self.ic = []

    def on_epoch_end(self, epoca, logs=None):
        p = self.model.predict(X[val], batch_size=8192, verbose=0)
        self.ic.append(ic(p.ravel(), val).mean())

keras.utils.set_random_seed(0)
manual = keras.Sequential([
    keras.Input((len(FE),)),
    layers.Dense(128, activation="relu"),
    layers.Dense(128, activation="relu"),
    layers.Dense(1)])
manual.compile(keras.optimizers.Adam(1e-3), "mse")
cb = ICValidacion()
h = manual.fit(X[ent], y[ent], epochs=20, batch_size=512, verbose=0,
               validation_data=(X[val], y[val]), callbacks=[cb])
curva_manual = pd.DataFrame({"train loss": h.history["loss"],
                             "val. loss": h.history["val_loss"],
                             "val. IC": cb.ic},
                            index=range(1, 21)).rename_axis("epoch")
curva_manual.iloc[[0, 1, 3, 7, 11, 15, 19]].round(4)

       train loss  val. loss  val. IC
epoch
1          1.0134     1.0058   0.0116
2          0.9986     1.0040   0.0067
4          0.9917     1.0086   0.0090
8          0.9811     1.0084   0.0005
12         0.9689     1.0170  -0.0005
16         0.9550     1.0317   0.0061
20         0.9389     1.0498   0.0035

The figure at the top plots both losses epoch by epoch. The training loss falls without stopping, from 1.013 to 0.939: 7% less, when the whole signal is worth less than 0.1%. That is memory, not learning. The validation loss reaches its minimum in the third epoch (1.0032) without ever going below that of predicting 0, and then rises to 1.050. The validation IC jumps between −0.0027 and 0.0118 with no trend, and in the last epochs it falls below zero: below Ridge in every epoch.

With 19,457 parameters (21×128 + 128×128 + 128 weights and their biases) and 103,600 rows of a label that is almost pure noise, the result was predictable. But there is a second cause, less obvious, that affects even a linear model.

The scale of the step and of the initialization

In this panel the good coefficients are tiny: Ridge's are worth 0.0054 on average in absolute value. Glorot initialization, the default for Dense, gives random weights of several hundredths, and Adam with rate 0.001 moves each weight by something of the order of 0.001 in every batch, in the direction set by a gradient that is almost all noise. The cell trains a single Dense(1) layer, which is a linear regression like the one in 5.6.1, in two ways:

def lineal(init, lr, semilla=0, epocas=12):
    keras.utils.set_random_seed(semilla)
    m = keras.Sequential([keras.Input((len(FE),)),
                          layers.Dense(1, kernel_initializer=init)])
    m.compile(keras.optimizers.Adam(lr), "mse")
    m.fit(X[ent], y[ent], epochs=epocas, batch_size=1024, verbose=0)
    w = m.layers[-1].get_weights()[0].ravel()
    p = m.predict(X[val], batch_size=8192, verbose=0).ravel()
    return ic(p, val).mean(), np.abs(w).mean()

w_ridge = np.abs(Ridge(alpha=10.0).fit(X[ent], y[ent]).coef_).mean()
for init, lr in [("glorot_uniform", 1e-3), ("zeros", 3e-4)]:
    ic_v, w_m = lineal(init, lr)
    print(f"{init:<15}lr {lr:.0e}: validation IC {ic_v:.4f}, "
          f"mean |weight| {w_m:.4f}")
print(f"Ridge: mean |weight| {w_ridge:.4f}")

glorot_uniform lr 1e-03: validation IC 0.0013, mean |weight| 0.0308
zeros          lr 3e-04: validation IC 0.0150, mean |weight| 0.0049
Ridge: mean |weight| 0.0054

With the default configuration, the linear layer ends up with weights of 0.031 on average, six times Ridge's, and a validation IC of 0.0013: noise. With the output initialized to zero and a rate of 0.0003, the weights end at 0.0049 and the IC at 0.0150, on a par with Ridge (0.0155). The model is the same; what changes is that the optimizer does not have to undo a random starting point larger than the signal, nor add with each step a noise larger than the effect it is looking for.

The rule of thumb for networks with financial data: initialize the last layer to zero (the network starts by predicting the mean), small rates and large batches. The other layers do need random initial weights: with everything at zero, the neurons of a layer would be identical and would receive the same gradient.

Regularization: dropout, weight decay and width

On that basis, three regularization tools and a fourth, the simplest, which is to use a small network:

• Dropout: in each batch, each hidden neuron is switched off with probability $p$ (here 0.3) and the rest are scaled by $1/(1-p)$. It forces the network not to depend on specific combinations of neurons. It only acts during training.
• Weight decay: each step shrinks the weights toward zero, $w \leftarrow w - \eta\,\lambda\,w$, apart from the gradient. It is the decoupled version of the L2 penalty (AdamW, by Loshchilov and Hutter): with Adam it is not the same as adding $\lambda\lVert w\rVert^2$ to the loss, because Adam rescales each gradient.
• Early stopping: stop training when validation stops improving (next section).
• Capacity: 16 neurons in one layer instead of 64×2.

Each configuration is trained for 8 epochs with three seeds. With only one, the effect of the regularization could not be told apart from the randomness of the initialization, of dropout and of the order of the batches:

def mlp(ancho, capas, dropout=0.0, wd=0.0, lr=3e-4):
    m = keras.Sequential([keras.Input((len(FE),))])
    for _ in range(capas):
        m.add(layers.Dense(ancho, activation="relu"))
        if dropout:
            m.add(layers.Dropout(dropout))
    m.add(layers.Dense(1, kernel_initializer="zeros"))
    opt = (keras.optimizers.AdamW(lr, weight_decay=wd) if wd
           else keras.optimizers.Adam(lr))
    m.compile(opt, "mse")
    return m

CONFIG = {"64x2 unregularized": (64, 2, 0.0, 0.0),
          "64x2 dropout 0.3": (64, 2, 0.3, 0.0),
          "64x2 weight decay 0.01": (64, 2, 0.0, 1e-2),
          "16x1 dropout + decay": (16, 1, 0.3, 1e-2)}
EPOCAS = 8
pred_val, filas = {}, []
for nombre, cfg in CONFIG.items():
    ics = []
    for s in range(3):
        keras.utils.set_random_seed(s)
        m = mlp(*cfg)
        m.fit(X[ent], y[ent], epochs=EPOCAS, batch_size=1024,
              verbose=0)
        p = m.predict(X[val], batch_size=8192, verbose=0).ravel()
        pred_val[(nombre, s)] = p
        ics.append(ic(p, val).mean())
    filas.append({"model": nombre, "parameters": m.count_params(),
                  "mean IC": np.mean(ics), "min": min(ics),
                  "max": max(ics)})
tabla_reg = pd.DataFrame(filas).set_index("model")
tabla_reg.round(4)

                        parameters  mean IC     min     max
model
64x2 unregularized            5633   0.0130  0.0113  0.0147
64x2 dropout 0.3              5633   0.0116  0.0111  0.0122
64x2 weight decay 0.01        5633   0.0128  0.0113  0.0144
16x1 dropout + decay           369   0.0126  0.0107  0.0163

The figure puts the twelve marks on one scale. Four readings:

• No configuration beats Ridge consistently. The means range from 0.0116 to 0.0130, all below 0.0155. Only one seed, that of the 16x1 with dropout and decay, goes above Ridge (0.0163).
• The seed weighs as much as the configuration. Within the 16x1 the three seeds range from 0.0107 to 0.0163; between the means of the four configurations there is 0.0014. With these numbers you cannot claim that one regularization is better than another.
• Regularizing does not rescue the large network. With the right basis (output at zero, low rate, few epochs), the unregularized 64x2 no longer overfits as much as the textbook MLP, and dropout does not improve it: with 8 epochs and a low rate, the short training itself acts as regularization.
• The small model costs fifteen times fewer parameters (369 against 5,633) for the same result. It is the one taken to the test.

Early stopping

keras.callbacks.EarlyStopping watches a validation metric, stops when it does not improve in patience epochs and, with restore_best_weights=True, returns the weights of the best epoch. Watching val_loss is the default:

keras.utils.set_random_seed(0)
m = mlp(*CONFIG["16x1 dropout + decay"])
parada = keras.callbacks.EarlyStopping(
    monitor="val_loss", patience=3, restore_best_weights=True)
cb = ICValidacion()
h = m.fit(X[ent], y[ent], epochs=30, batch_size=1024, verbose=0,
          validation_data=(X[val], y[val]), callbacks=[cb, parada])
perd = np.array(h.history["val_loss"])
print(f"epochs trained: {len(perd)}; best validation loss "
      f"in epoch {perd.argmin() + 1}")
print(f"validation IC in that epoch: {cb.ic[perd.argmin()]:.4f}; "
      f"best IC: {max(cb.ic):.4f} (epoch {np.argmax(cb.ic) + 1})")
print("validation losses:", np.round(perd, 5).tolist())

epochs trained: 5; best validation loss in epoch 2
validation IC in that epoch: 0.0134; best IC: 0.0134 (epoch 2)
validation losses: [0.99664, 0.99664, 0.99665, 0.99671, 0.99673]

The validation loss barely moves: the five epochs lie between 0.99664 and 0.99673, differences in the fourth decimal place. The stop is decided by a noise of 0.00001. Here it coincides with the best IC (epoch 2, 0.0134), but that is luck: in the textbook MLP, the IC had no relation to the loss. The consistent thing is to stop with the metric that matters, the validation IC; lesson 5.6.4 writes that stopping callback. And any stopping criterion uses validation to choose: the number of epochs is one more hyperparameter, and validation stops being a clean estimate (the test is still clean).

The seed lottery

The three seeds of the chosen model give different ICs. Training several networks and averaging their predictions is the cheapest way to reduce that variance (an ensemble like those in lesson 5.10.3):

nombre = "16x1 dropout + decay"
por_semilla = [ic(pred_val[(nombre, s)], val).mean() for s in range(3)]
media = np.mean([pred_val[(nombre, s)] for s in range(3)], axis=0)
print("IC per seed:", np.round(por_semilla, 4).tolist())
print(f"IC of the average of the three: {ic(media, val).mean():.4f}")
corr = np.corrcoef([pred_val[(nombre, s)] for s in range(3)])
print("correlation between predictions:",
      np.round(corr[np.triu_indices(3, 1)], 3).tolist())

IC per seed: [0.0109, 0.0163, 0.0107]
IC of the average of the three: 0.0141
correlation between predictions: [0.771, 0.732, 0.783]

The predictions of the three networks have correlations between 0.73 and 0.78: they are similar, but they are not the same network. The average gives 0.0141, above two of the three seeds and without depending on having picked the good one. Choosing the 0.0163 seed would be the usual selection bias (lessons 2.5.1 and 5.5.4): the maximum of three noises.

The test, only once

The chosen model (16x1 with dropout 0.3 and weight decay 0.01, 8 epochs, average of three seeds) is retrained on 2012-2018 and measured on 2019-2025, together with the baselines and the oracle. The important column is the difference from Ridge week by week, with its standard error: the 365 weeks are the same for everyone, and the paired comparison removes the common noise.

preds = []
for s in range(3):
    keras.utils.set_random_seed(s)
    m = mlp(*CONFIG["16x1 dropout + decay"])
    m.fit(X[tr_todo], y[tr_todo], epochs=EPOCAS, batch_size=1024,
          verbose=0)
    preds.append(m.predict(X[te], batch_size=8192, verbose=0).ravel())
ic_mlp = ic(np.mean(preds, axis=0), te)

verdad = pd.read_parquet("panel5_verdad.parquet",
                         columns=["fecha", "ticker", "alfa_verdadero"])
alfa = d[["fecha", "ticker"]].merge(verdad, how="left")
ic_oraculo = ic(alfa.loc[te, "alfa_verdadero"].to_numpy(), te)

filas = []
for nombre, serie in [("oracle", ic_oraculo),
                      ("MLP (3 seeds)", ic_mlp),
                      ("ridge", base["ridge"]),
                      ("lightgbm", base["lightgbm"])]:
    dif = serie - base["ridge"]
    filas.append({"model": nombre, "IC": serie.mean(),
                  "s.e.": serie.std() / np.sqrt(len(serie)),
                  "minus ridge": dif.mean(),
                  "s.e. diff.": dif.std() / np.sqrt(len(dif))})
pd.DataFrame(filas).set_index("model").round(4)

                   IC    s.e.  minus ridge  s.e. diff.
model
oracle         0.0256  0.0052       0.0076      0.0043
MLP (3 seeds)  0.0172  0.0042      -0.0008      0.0017
ridge          0.0180  0.0045       0.0000      0.0000
lightgbm       0.0139  0.0047      -0.0041      0.0022

The honest result:

• The MLP gives an IC of 0.0172 ± 0.0042 and Ridge 0.0180 ± 0.0045. The difference, −0.0008 ± 0.0017, is half its standard error: there is no evidence that the MLP improves on Ridge, nor that it is worse.
• LightGBM stays below Ridge (−0.0041 ± 0.0022), as in lesson 5.5.3.
• The oracle, 0.0256, marks the ceiling: the three models capture between half and two thirds of it. What is missing lies in the interactions and in the effects that change with the regime, exactly what a network should find and what, with 14 years of weekly data, cannot be told apart from noise.

In a tie, the simple model wins: Ridge trains in a second, does not depend on a seed, its coefficients can be read and audited, and its out-of-sample behavior is more predictable. This agrees with the literature on tabular data (Grinsztajn, Oyallon and Varoquaux, 2022: with medium-sized datasets, trees are still ahead of networks). In the other direction, Gu, Kelly and Xiu (2020) found that networks and trees were the best methods for predicting US stock returns, thanks to nonlinear interactions, with many more features and decades of monthly data. The two are compatible: the MLP needs much more information than a signal this weak in a few features provides.

When an MLP is worth it

None of the above says that networks are useless in finance; it says when they add nothing. The cases in which a network has a real advantage are those that go beyond a small table of features:

• Many features and many interactions, with large samples (hundreds of characteristics, decades, thousands of assets), as in Gu, Kelly and Xiu.
• Structured inputs: sequences (lesson 5.6.3), text, order book images, graphs (lesson 5.7.4).
• Tailored losses: optimizing the portfolio or a utility function directly (lessons 5.6.4 and 5.7.5).
• Shared representations: an embedding of the asset or the sector learned together with the model (functional API, lesson 5.6.1).

For a tabular panel like this one, the sensible order is: Ridge as the baseline, boosting as the first nonlinear alternative (5.5.3) and a small, regularized and averaged MLP as the third, only if the purged validation justifies it.

Exercises

1. Compute the test IC of the MLP and of Ridge year by year. In how many of the seven years does the MLP win?

anual = pd.DataFrame({"mlp": ic_mlp, "ridge": base["ridge"]})
anual = anual.groupby(anual.index.year.rename("year")).mean()
anual["mlp wins"] = anual["mlp"] > anual["ridge"]
anual.round(4)

         mlp   ridge  mlp wins
year
2019 -0.0011  0.0021     False
2020  0.0054  0.0062     False
2021  0.0216  0.0110      True
2022  0.0313  0.0311      True
2023 -0.0106 -0.0041     False
2024  0.0439  0.0497     False
2025  0.0299  0.0302     False

2. Build a wide and deep model with the functional API: a linear branch (Dense(1) on the features) plus the branch of the 16x1 MLP, added with layers.Add. Train it like the previous ones and compare its validation IC with the one in the table.

ics = []
for s in range(3):
    keras.utils.set_random_seed(s)
    x_in = keras.Input((len(FE),))
    lin = layers.Dense(1, kernel_initializer="zeros")(x_in)
    h = layers.Dropout(0.3)(layers.Dense(16, activation="relu")(x_in))
    prof = layers.Dense(1, kernel_initializer="zeros")(h)
    ap = keras.Model(x_in, layers.Add()([lin, prof]))
    ap.compile(keras.optimizers.AdamW(3e-4, weight_decay=1e-2), "mse")
    ap.fit(X[ent], y[ent], epochs=EPOCAS, batch_size=1024, verbose=0)
    p = ap.predict(X[val], batch_size=8192, verbose=0).ravel()
    ics.append(ic(p, val).mean())
print("Validation IC per seed:", np.round(ics, 4).tolist())

Validation IC per seed: [0.0158, 0.0152, 0.0161]

3. The lesson chose the number of epochs (8) and the configuration by looking at a single validation split (2017-2018). Propose a more robust procedure with the tools of module 5.4 and state what it costs with networks.