Supervised, unsupervised and reinforcement learning; scikit-learn pipelines.

Three ways of learning

The level uses three families of methods. They differ in what they receive to learn from, not in the algorithm:

flowchart LR
  ML[Learning] --> S[Supervised]
  ML --> U[Unsupervised]
  ML --> R[Reinforcement]
  S --> SR[Regression]
  S --> SC[Classification]
  U --> UD[Dimensionality reduction]
  U --> UG[Clustering]
  R --> RP[Policy]

• Supervised. Each example comes with its features $x$ and an observed label $y$, and the model learns $f(x) \approx y$. If $y$ is a number (next week's return, tomorrow's volatility) it is regression; if it is a class (up or down, which barrier is hit first), classification.
• Unsupervised. There is only $x$. The model looks for structure: a few directions that summarize many columns (PCA), groups of assets that move together, market states that repeat.
• Reinforcement. An agent makes decisions in sequence, the environment responds and returns a reward; the agent learns a policy (what to do in each state) that maximizes the cumulative reward. There is no label for “the correct action”, only consequences.

Most of the level is supervised learning, and so are most of the mistakes: the label is information from the future by definition, and every step in the process that touches it in the wrong place leaks it into the past. The rest of the lesson builds the scikit-learn chain that avoids those mistakes by design.

The scikit-learn API: what each step learns

scikit-learn has a uniform API. An estimator learns from the data with fit and stores what it learns in attributes ending in an underscore (mean_, coef_). A transformer also has transform, which applies what it learned to new data; a predictor, predict. The rule that orders everything: fit is called only with training data, transform and predict with any data.

The data is the panel5 panel from lesson 5.1.1: 21 weekly features (without the revision_analistas leak), the sector and, as the label, the following week's return minus its cross-sectional mean (what matters for ranking stocks). Training up to 2018 and test from 2019, with a one-week gap:

import numpy as np
import pandas as pd

f = pd.read_parquet("panel5_features.parquet")
et = pd.read_parquet("panel5_etiquetas.parquet")
d = (f.merge(et[["fecha", "ticker", "ret_sig_semana"]],
             on=["fecha", "ticker"])
      .dropna(subset=["ret_sig_semana"]).reset_index(drop=True))
d["y"] = d["ret_sig_semana"] - d.groupby("fecha")[
    "ret_sig_semana"].transform("mean")
num = [c for c in f.columns
       if c not in ("fecha", "ticker", "sector", "revision_analistas")]
cols = ["fecha", "sector"] + num
ent = d[d["fecha"] < "2019-01-01"].reset_index(drop=True)
pru = d[d["fecha"] >= "2019-01-08"].reset_index(drop=True)
print(len(ent), "training rows;", len(pru), "test rows")
print(ent["fecha"].max().date(), "|", pru["fecha"].min().date())

145600 training rows; 145600 test rows
2018-12-28 | 2019-01-11

A StandardScaler learns the mean and standard deviation of each column in fit and uses them in transform. Applied to the test set, it uses the training ones: that is how a live system would see the data, because on Monday it does not know the means of the next seven years.

from sklearn.preprocessing import StandardScaler

esc = StandardScaler().fit(ent[["log_capitalizacion", "vol_63d"]])
print("learned means:", esc.mean_.round(3))
z = esc.transform(pru[["log_capitalizacion", "vol_63d"]])
print("test mean after transforming:", z.mean(axis=0).round(3))

learned means: [8.932 0.29 ]
test mean after transforming: [0.131 0.197]

The transformed test features do not have zero mean: the stocks' market cap grew between the two periods (0.13 standard deviations above the training mean) and volatility in 2019-2025 was higher (0.20 standard deviations). That is correct. If the mean came out exactly zero, the scaling would have been learned with the test data.

Pipeline: the chain as a single estimator

A Pipeline chains transformers and ends in a predictor, and it behaves like a single estimator: fit fits each step with the output of the previous one and predict applies the whole chain. The advantage is not convenience: it is that any scikit-learn validation tool refits all the steps with the training data of each split. For a stateful step not to see the test, it is enough for it to be inside.

Before the pipeline you need the metric and the first custom transformer. The metric is the mean weekly IC (lesson 5.1.3); the puntuacion_ic function has the signature of a scikit-learn scorer, (estimator, X, y), and takes the dates from X itself:

def ic_semanal(fecha, pred, y):
    """Rank correlation week by week."""
    s = pd.DataFrame({"fecha": np.asarray(fecha),
                      "p": np.asarray(pred), "y": np.asarray(y)})
    r = s.groupby("fecha")[["p", "y"]].rank()
    r = r - r.groupby(s["fecha"]).transform("mean")
    num_ = (r["p"] * r["y"]).groupby(s["fecha"]).sum()
    den = ((r["p"] ** 2).groupby(s["fecha"]).sum()
           * (r["y"] ** 2).groupby(s["fecha"]).sum()) ** 0.5
    return num_ / den


def puntuacion_ic(estimador, X, y):
    return ic_semanal(X["fecha"], estimador.predict(X), y).mean()

The custom transformer is a winsorizer: it clips each column to its 1% and 99% quantiles. It is stateful: the quantiles are learned in fit and applied in transform. A scikit-learn transformer inherits from BaseEstimator (which gives it get_params and set_params, needed for hyperparameter search) and from TransformerMixin (which gives it fit_transform); the parameters are stored in __init__ untouched, what is learned in attributes ending in _, and get_feature_names_out says what the output columns are called:

from sklearn.base import BaseEstimator, TransformerMixin


class Winsorizador(BaseEstimator, TransformerMixin):
    """Clips each column to the q and 1 - q quantiles of the fit."""

    def __init__(self, q=0.01):
        self.q = q

    def fit(self, X, y=None):
        X = pd.DataFrame(X)
        self.bajo_ = X.quantile(self.q)
        self.alto_ = X.quantile(1 - self.q)
        self.feature_names_in_ = np.asarray(X.columns, dtype=object)
        return self

    def transform(self, X):
        return pd.DataFrame(X).clip(self.bajo_, self.alto_, axis=1)

    def get_feature_names_out(self, input_features=None):
        return self.feature_names_in_


w = Winsorizador(0.01).fit(ent[num])
print(pd.DataFrame({"low": w.bajo_, "high": w.alto_})
      .loc[["ret_1s", "sue", "sentimiento"]].round(3))

               low   high
ret_1s      -0.112  0.116
sue         -2.990  3.662
sentimiento -2.344  2.307

With it, a complete linear pipeline: winsorize, impute the gaps in sentimiento (32% of the rows have no coverage) with the median and a missing-value indicator column, standardize and fit a ridge:

from sklearn.impute import SimpleImputer
from sklearn.linear_model import Ridge
from sklearn.pipeline import make_pipeline

lineal = make_pipeline(Winsorizador(0.01),
                       SimpleImputer(strategy="median",
                                     add_indicator=True),
                       StandardScaler(), Ridge(alpha=1000.0))
lineal.set_output(transform="pandas")
lineal.fit(ent[num], ent["y"])
print(list(lineal.named_steps))
print(lineal[:-1].get_feature_names_out()[-3:])
ic = ic_semanal(pru["fecha"], lineal.predict(pru[num]), pru["y"])
print(f"Mean IC 2019-2025: {ic.mean():.4f}")

['winsorizador', 'simpleimputer', 'standardscaler', 'ridge']
['ruido_4' 'ruido_5' 'missingindicator_sentimiento']
Mean IC 2019-2025: 0.0181

lineal[:-1] is the pipeline without the last step: it shows which columns the ridge receives, including the missingindicator_sentimiento indicator added by the imputer. The 2019-2025 IC is 0.018, close to the ridge in lesson 5.1.1 (0.019) and below the oracle's ceiling in that period (0.025).

ColumnTransformer and a stateless transformer

Not all features are treated the same way. The sector is a category and is encoded into binary columns; the numeric ones are turned into a cross-sectional rank: each stock's percentile among the 400 of its week. The rank removes extreme values and each era's level in one go (a 30% volatility is high in 2017 and normal in 2020). Lesson 5.2.4 compares this and other scalings; here what matters is its form: it is a stateless transformer, because each week is transformed only with that week's data. Its fit learns nothing from the features, so it cannot leak anything from the future.

ColumnTransformer applies a different transformer to each group of columns and concatenates the result:

flowchart TB
  X[fecha, sector, features] --> RF[RangoPorFecha]
  RF --> CT[ColumnTransformer]
  CT --> IM[Numeric: impute]
  CT --> OH[Sector: one-hot]
  IM --> M[Ridge]
  OH --> M

from sklearn.compose import ColumnTransformer
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import OneHotEncoder


class RangoPorFecha(BaseEstimator, TransformerMixin):
    """Percentile of each column within its date, minus 0.5."""

    def __init__(self, columnas):
        self.columnas = columnas

    def fit(self, X, y=None):
        self.resto_ = [c for c in X.columns
                       if c not in self.columnas + ["fecha"]]
        return self

    def transform(self, X):
        g = X.groupby("fecha")[self.columnas]
        return (g.rank(pct=True) - 0.5).join(X[self.resto_])

    def get_feature_names_out(self, input_features=None):
        return np.asarray(self.columnas + self.resto_, dtype=object)


prep = ColumnTransformer(
    [("num", SimpleImputer(strategy="constant", fill_value=0.0), num),
     ("sec", OneHotEncoder(handle_unknown="ignore",
                           sparse_output=False), ["sector"])],
    verbose_feature_names_out=False)
modelo = Pipeline([("rango", RangoPorFecha(num)), ("prep", prep),
                   ("ridge", Ridge())])
modelo.set_output(transform="pandas")
Xt = modelo[:-1].fit_transform(ent[cols])
print(Xt.shape)
Xt.iloc[:2, [0, 3, 20, 21, 22]].round(3)

(145600, 31)
   ret_1s  mom_12_1  ruido_5  sector_consumo_basico  sector_consumo_ciclico
0  -0.305    -0.350   -0.440                    0.0                     0.0
1  -0.285     0.047   -0.107                    0.0                     0.0

The rank leaves each feature between −0.5 and 0.5; the imputer sets to 0 (the median of their week) the stocks with no sentiment data, which the rank leaves empty. The one-hot encoding adds ten columns, one per sector. handle_unknown="ignore" is a production decision: if a new sector appears, it is encoded as zeros instead of breaking the system.

Splitting by date: validating a panel

To choose hyperparameters you need to validate within the training period. Standard cross-validation (KFold) shuffles the rows and, in financial series, promises results that do not exist; lesson 5.4.1 measures it in detail. scikit-learn's alternative for series, TimeSeriesSplit, trains on the past and tests on the future, but it splits by rows. In a panel, each date has 400 rows:

from sklearn.model_selection import TimeSeriesSplit

for ent_i, pru_i in TimeSeriesSplit(n_splits=5).split(ent):
    comunes = set(ent.loc[ent_i, "fecha"]) & set(
        ent.loc[pru_i, "fecha"])
    print(f"{len(ent_i):6d} training rows; "
          f"shared weeks: {len(comunes)}")

 24270 training rows; shared weeks: 1
 48536 training rows; shared weeks: 1
 72802 training rows; shared weeks: 1
 97068 training rows; shared weeks: 1
121334 training rows; shared weeks: 1

In all five splits the boundary falls in the middle of a week: some of that week's stocks are used for training and others for testing, with the same market and the same common return. In addition, the label of the last training week ends one Friday later, inside the test period. Both problems are fixed with a date-based splitter with a gap of at least the label horizon (one week for ret_sig_semana; for ret_sig_21d, five). A scikit-learn splitter is any object with split and get_n_splits:

class SplitPorFechas:
    """Expanding window by date: train on the past, leave
    `hueco` dates unused and test on the next block."""

    def __init__(self, n_splits=5, hueco=1):
        self.n_splits = n_splits
        self.hueco = hueco

    def get_n_splits(self, X=None, y=None, groups=None):
        return self.n_splits

    def split(self, X, y=None, groups=None):
        fechas = np.asarray(X["fecha"])
        unicas = np.unique(fechas)
        bloques = np.array_split(unicas, self.n_splits + 1)
        for b in bloques[1:]:
            i = np.searchsorted(unicas, b[0]) - 1 - self.hueco
            yield (np.flatnonzero(fechas <= unicas[i]),
                   np.flatnonzero(np.isin(fechas, b)))


cv = SplitPorFechas(n_splits=5, hueco=1)
for i_ent, i_pru in cv.split(ent):
    fe, fp = ent.loc[i_ent, "fecha"], ent.loc[i_pru, "fecha"]
    print(fe.min().date(), "to", fe.max().date(), "| test",
          fp.min().date(), "to", fp.max().date())

2012-01-13 to 2013-03-01 | test 2013-03-15 to 2014-05-09
2012-01-13 to 2014-05-02 | test 2014-05-16 to 2015-07-10
2012-01-13 to 2015-07-02 | test 2015-07-17 to 2016-09-09
2012-01-13 to 2016-09-02 | test 2016-09-16 to 2017-11-03
2012-01-13 to 2017-10-27 | test 2017-11-10 to 2018-12-28

Each split trains on everything before it (expanding window), leaves one week unused and tests on the next block of about 60 weeks. The figure draws the five splits and, separately, the final 2019-2025 test period, which the search does not touch. It is a simple walk-forward; lessons 5.4.2 and 5.4.3 add purging, embargo and rolling windows.

Hyperparameter search with a trading metric

GridSearchCV tries each combination in a grid with the splitter it is given and scores it with the scorer. The hyperparameters of a pipeline step are named step__parameter. Here the search covers the ridge's regularization, α, across nine orders of magnitude:

from sklearn.model_selection import GridSearchCV

rejilla = {"ridge__alpha": [1, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7,
                            1e8]}
busqueda = GridSearchCV(modelo, rejilla, scoring=puntuacion_ic,
                        cv=cv, n_jobs=1)
busqueda.fit(ent[cols], ent["y"])
res = pd.DataFrame(busqueda.cv_results_)
res = res.set_index("param_ridge__alpha")[
    ["mean_test_score", "std_test_score"]]
print(res.round(4))
print("chosen:", busqueda.best_params_)

                    mean_test_score  std_test_score
param_ridge__alpha
1.0                          0.0043          0.0135
10.0                         0.0043          0.0135
100.0                        0.0045          0.0134
1000.0                       0.0052          0.0128
10000.0                      0.0063          0.0115
100000.0                     0.0068          0.0104
1000000.0                    0.0070          0.0102
10000000.0                   0.0070          0.0102
100000000.0                  0.0070          0.0102
chosen: {'ridge__alpha': 10000000.0}

Two things that are not obvious. First: the validation IC rises with α and levels off at huge values. When α dominates $X^\top X$, the ridge gives $\hat\beta \approx X^\top y/\alpha$: the coefficients tend to zero, but their direction stops changing, and the IC only looks at the order of the predictions. With so little signal, heavy shrinkage is what generalizes best: the ridge almost becomes an average of the features weighted by their covariance with the label. Second: the standard deviation across splits (0.010 for the best one) is larger than the differences between values of α. Choosing between $\alpha = 10^6$ and $\alpha = 10^8$ with this data is choosing between noises; lesson 5.5.4 measures how much the validation itself can be overfitted when the grid is large.

GridSearchCV retrains the best model on the whole training set (refit=True, the default). It is evaluated once on 2019-2025 and compared with an untuned baseline: the sum of the ranks of four features with the sign given by the literature (momentum, value and SUE positive; one-week reversal negative).

mejor = busqueda.best_estimator_
ic_mod = ic_semanal(pru["fecha"], mejor.predict(pru[cols]),
                    pru["y"])
rk = pru.groupby("fecha")[["mom_12_1", "ret_1s", "valor_bp",
                           "sue"]].rank(pct=True)
base = rk["mom_12_1"] - rk["ret_1s"] + rk["valor_bp"] + rk["sue"]
ic_base = ic_semanal(pru["fecha"], base, pru["y"])
for nombre, ic in [("chosen ridge", ic_mod), ("baseline", ic_base)]:
    t = ic.mean() / ic.std() * np.sqrt(len(ic))
    print(f"{nombre:14s} IC {ic.mean():.4f}  t {t:.1f}")
coef = pd.Series(mejor["ridge"].coef_,
                 index=mejor[:-1].get_feature_names_out())
(coef / coef.abs().max()).sort_values().iloc[[0, 1, 2, -3, -2, -1]
                                            ].round(2)

chosen ridge   IC 0.0189  t 3.1
baseline       IC 0.0178  t 4.1
ret_1s                      -1.00
rsi_14                      -0.71
sector_telecomunicaciones   -0.62
sector_energia               0.53
mom_12_1                     0.55
sector_industria             0.65
dtype: float64

The chosen ridge gives an IC of 0.019 in 2019-2025 and the four-feature baseline 0.018, with a higher t-statistic (4.1 versus 3.1: it is more stable from one week to the next). The difference in means is smaller than the error of either one. This is the comparison you must always make: if a pipeline with 21 features, imputation, encoding and a search does not clearly beat a sum of four ranks, its complexity adds nothing.

The largest coefficients (relative to the largest) tell what it learned. The one-week reversal with a negative sign and momentum with a positive sign are effects of the generator; the negative RSI inherits the reversal. But three of the six are sectors: telecommunications negative, energy and industrials positive. The generator puts no sector effect in the expected return; those coefficients are the sectors that, by chance, returned more or less than the mean in 2012-2018. With such strong regularization, the ridge bets in proportion to each column's covariance with the label, and that of a binary sector column is simply its past return. A model can improve in validation by learning a period's noise; removing the sector or neutralizing the portfolio by sector (lesson 5.2.4) is a decision made on economic reasoning, not on what the search says.

The leak inside the validation: selecting with all the data

The most expensive pipeline mistake is not scaling with the test data, which here would change little. It is selecting with it. The experiment: 500 columns of pure noise, with no relationship at all with the label, over two years of data (2012-2013, 41,200 rows). The wrong way picks the 10 most correlated with the label using all the rows and then validates a ridge with them. The right way puts the selection inside the pipeline, so in each split the choice is made only with the training rows. As a control, the 10 columns chosen the wrong way are also evaluated on data the selection did not see (2014-2018):

from sklearn.feature_selection import SelectKBest, f_regression
from sklearn.model_selection import cross_val_score


class SinFecha(BaseEstimator, TransformerMixin):
    def fit(self, X, y=None):
        return self

    def transform(self, X):
        return X.drop(columns="fecha")


rng = np.random.default_rng(0)
ruido = [f"r{i:03d}" for i in range(500)]
datos = d[d["fecha"] < "2019-01-01"].reset_index(drop=True)
Z = pd.DataFrame(rng.standard_normal((len(datos), 500)),
                 columns=ruido).assign(fecha=datos["fecha"])
dos = (datos["fecha"] < "2014-01-01").to_numpy()
Z2, y2 = Z[dos].reset_index(drop=True), datos.loc[dos, "y"]

elegidas = SelectKBest(f_regression, k=10).fit(
    Z2[ruido], y2).get_feature_names_out()
mal = cross_val_score(make_pipeline(SinFecha(), Ridge(1.0)),
                      Z2[list(elegidas) + ["fecha"]], y2, cv=cv,
                      scoring=puntuacion_ic)
bien = cross_val_score(
    make_pipeline(SinFecha(), SelectKBest(f_regression, k=10),
                  Ridge(1.0)),
    Z2, y2, cv=cv, scoring=puntuacion_ic)
despues = Ridge(1.0).fit(Z2[elegidas], y2).predict(Z[~dos][elegidas])
ic_desp = ic_semanal(Z.loc[~dos, "fecha"], despues,
                     datos.loc[~dos, "y"]).mean()
for nombre, v in [("outside the pipeline", mal),
                  ("inside the pipeline", bien)]:
    print(f"{nombre:20s} IC {v.mean():.3f} ± "
          f"{v.std(ddof=1) / np.sqrt(len(v)):.3f}", v.round(3))
print(f"selected, 2014-2018  IC {ic_desp:.3f}")

outside the pipeline IC 0.028 ± 0.009 [0.035 0.02  0.022 0.006 0.058]
inside the pipeline  IC 0.008 ± 0.004 [ 0.004  0.018 -0.002  0.     0.019]
selected, 2014-2018  IC -0.001

The figure at the top summarizes the cell. The wrong way gives a mean validation IC of 0.028 with 500 noise columns: the same level as the oracle's ceiling in the panel. There is no signal at all; the selection had already seen the labels of each split's test weeks and chose the columns that, by chance, resembled them. With the selection inside the pipeline, the IC is 0.008 with an error of 0.004, consistent with zero. And the “chosen” columns predict nothing in the following five years. Hastie, Tibshirani and Friedman call it “the wrong way to do cross-validation” (section 7.10.2 of their book); in finance it is more dangerous than in other fields because noise is almost everything and any selection using the label finds something.

There are leaks the pipeline cannot prevent, because they are not in the process but in the data: a feature the vendor dates before it is known (revision_analistas, lesson 5.2.2) or labels that overlap with the test (lessons 5.3.4 and 5.4.2). The pipeline protects against procedural mistakes; data mistakes are found by auditing dates.

Exercises

1. Replace the pipeline's ridge with a scikit-learn HistGradientBoostingRegressor (max_iter=100, learning_rate=0.05, max_leaf_nodes=8, min_samples_leaf=2000, random_state=0) and evaluate it with cross_val_score and the same splitter. Does it beat the ridge in validation?

from sklearn.ensemble import HistGradientBoostingRegressor

hgb = Pipeline([
    ("rango", RangoPorFecha(num)), ("prep", prep),
    ("hgb", HistGradientBoostingRegressor(
        max_iter=100, learning_rate=0.05, max_leaf_nodes=8,
        min_samples_leaf=2000, random_state=0))])
hgb.set_output(transform="pandas")
v = cross_val_score(hgb, ent[cols], ent["y"], cv=cv,
                    scoring=puntuacion_ic)
print("boosting:", v.round(4), round(v.mean(), 4))
print("chosen ridge:", round(busqueda.best_score_, 4))

boosting: [ 0.0053 -0.0085  0.0005  0.0099 -0.0118] -0.0009
chosen ridge: 0.007

2. Check how much the ridge's 2019-2025 IC changes if the linear pipeline's StandardScaler is fitted on the whole 2012-2025 period (the “small” leak) instead of only on the training set.

todo = pd.concat([ent, pru])
w = Winsorizador(0.01).fit(ent[num])
imp = SimpleImputer(strategy="median", add_indicator=True)
imp.set_output(transform="pandas")
imp.fit(w.transform(ent[num]))
fuga = StandardScaler().fit(imp.transform(w.transform(todo[num])))
r = Ridge(1000.0).fit(fuga.transform(
    imp.transform(w.transform(ent[num]))), ent["y"])
p = r.predict(fuga.transform(imp.transform(w.transform(pru[num]))))
ic_fuga = ic_semanal(pru["fecha"], p, pru["y"]).mean()
ic_ok = ic_semanal(pru["fecha"], lineal.predict(pru[num]),
                   pru["y"]).mean()
print(f"no leak {ic_ok:.4f} | with leak {ic_fuga:.4f}")

no leak 0.0181 | with leak 0.0182

3. You are going to predict ret_sig_21d (the 21-session return, which overlaps between consecutive weeks) with the same pipeline. What would you change in the splitter and why?