Una delle aree più problematiche del trading quantitativo è l’ottimizzazione di una strategia di previsione per migliorarne le prestazioni.

I trader quantistici esperti sono ben consapevoli che è fin troppo facile generare una strategia con capacità predittive stellari durante un backtest. Tuttavia, alcuni backtest possono mascherare il pericolo di un modello overfit , che può portare a una drastica sottoperformance quando viene implementata una strategia.

In questo articolo descriviamo un approccio per ridurre il problema dell’overfitting di un modello di machine learning, utilizzando una tecnica nota come cross-validation.

Per prima cosa introduciamo la definizione della cross-validation e poi descriviamo il funzionamento. In secondo luogo, costruiamo un modello di previsione utilizzando un indice azionario e quindi applichiamo due metodi di validazione incrociata a questo esempio: il validation set approach e k-fold cross-validation. Infine discuteremo il codice per le simulazioni utilizzando Python, Pandas , Matplotlib e Scikit-Learn .

Questo articolo è il “successore spirituale” di un precedente articolo scritto di recente sul compromesso tra bias e varianza . In quell’articolo abbiamo menzionato la convalida incrociata come un mezzo per risolvere alcuni dei problemi causati dal compromesso bias-varianza.

Il nostro obiettivo è infine creare una serie di strumenti statistici che possono essere utilizzati all’interno di un motore di backtesting per aiutarci a ridurre al minimo il problema dell’overfitting di un modello e quindi limitare le perdite future a causa di una strategia scarsamente performante.

Panoramica della Cross-Validation

Nel precedente articolo sul compromesso bias-varianza sono state introdotte le definizioni di errore di test e flessibilità:

  • Errore di test: l’errore medio, dove la media è calcolata tra molte osservazioni, associato alle prestazioni predittive di un particolare modello statistico quando è valutato su nuove osservazioni che non sono state utilizzate per addestrare il modello .
  • Flessibilità : i gradi di libertà a disposizione del modello per “adattarsi” ai dati di addestramento. Una regressione lineare è molto rigida (ha solo due gradi di libertà) mentre un polinomio di alto grado è molto flessibile (e come tale può avere molti gradi di libertà).

Con questi concetti in mente possiamo ora definire la cross-validation:

L’obiettivo della cross-validation è stimare l’errore di test associato a un modello statistico o selezionare il livello di flessibilità appropriato per un particolare metodo statistico.

Ancora una volta, possiamo ricordare dall’articolo sul compromesso bias-varianza che l’ errore di addestramento associato a un modello può sottovalutare notevolmente l’errore di test del modello. La convalida incrociata ci fornisce la capacità di stimare più accuratamente l’errore di test, che non conosceremo mai nella pratica.

La convalida incrociata consiste nell’escludere specifici sottoinsiemi dai dati di addestramento al fine di usarli come osservazioni di test. In questo articolo discuteremo i vari modi in cui tali sottoinsiemi vengono distribuiti e implementeremo i metodi usando Python su un modello di previsione di esempio basato su dati storici precedenti.

Esempio di previsione

Per rendere concreta la seguente discussione teorica prenderemo in considerazione lo sviluppo di una nuova strategia di trading basata sulla previsione dei livelli di prezzo di un indice azionario. Prenderemo in considerazione l’indicie S&P500, che contiene un raggruppamento ponderato delle cinquecento aziende più grandi società quotate (per capitalizzazione di mercato) a Wall Street. Allo stesso modo potremmo scegliere l’Euro Stoxx 50  o il DAX .

Per questa strategia considereremo semplicemente il prezzo di chiusura delle barre giornaliere storiche Open-High-Low-Close (OHLC) come predittori e il prezzo di chiusura del giorno successivo come risposta. Quindi stiamo tentando di prevedere il prezzo di domani utilizzando i prezzi storici giornalieri.

Un’osservazione sarà costituito da una coppia di vettori , \(X\) e \(y\), che contengono rispettivamente i valori predittori e il valore di risposta. Se consideriamo un ritardo giornaliero di \(p\) giorni, \(X\) ha \(p\) componenti. Ciascuno di questi componenti rappresenta il prezzo di chiusura di un giorno precedente. \(X_p\) rappresenta il prezzo di chiusura di oggi (noto), mentre \(X_{p-1}\) rappresenta il prezzo di chiusura di ieri, mentre \(X_1\) rappresenta il prezzo di \(p-1\) giorni fa.

\(Y\) contiene un solo valore, vale a dire il prezzo di chiusura di domani, ed è quindi uno scalare. Ogni osservazione è una tupla \((X, y)\). Considereremo una serie di \(n\) osservazioni corrispondenti a \(n\) giorni di informazioni storiche sui prezzi del SP500.

Il nostro obiettivo è trovare un modello statistico che tenti di prevedere il livello di prezzo del SP500 in base ai prezzi dei giorni precedenti. Se dovessimo ottenere una previsione accurata, potremmo usarla per generare segnali di trading di base. Questo articolo è principalmente interessato alla parte precedente del modello, quella della componente predittiva.

Useremo la convalida incrociata in due modi: in primo luogo per stimare l’errore di test di particolari metodi di apprendimento statistico (cioè le loro separate prestazioni predittive), e in secondo luogo per selezionare la flessibilità ottimale del metodo scelto al fine di minimizzare gli errori associati a bias e varianza.

Descriveremo ora i diversi modi di eseguire la convalida incrociata, iniziando con l’ approccio dell’insieme di convalida e poi infine con la convalida incrociata k-fold . In ogni caso useremo Pandas e Scikit-Learn per implementare questi metodi.

Validation Set Approach

L’approccio del set di convalida alla cross-validation è molto semplice da eseguire. Essenzialmente prendiamo l’insieme delle osservazioni (\(n\) giorni di dati) e le dividiamo casualmente in due metà. Una metà è nota come set di addestramento mentre la seconda metà è nota come set di convalida . Il modello si adatta utilizzando solo i dati nel set di addestramento, mentre il suo errore di test viene stimato utilizzando solo il set di convalida.

Questo è facilmente riconoscibile come una tecnica spesso utilizzata nel trading quantitativo come meccanismo per valutare le prestazioni predittive. Tuttavia, è più comune trovare due terzi dei dati utilizzati per il set di addestramento, mentre il terzo rimanente viene utilizzato per la convalida. Inoltre è più comune mantenere l’ordine delle serie temporali in modo tale che i primi due terzi rappresentino cronologicamente i primi due terzi dei dati storici.

Meno frequente è l’applicazione di questo metodo per la randomizzazione delle osservazioni in ciascuno dei due set. Ancora meno frequente è una discussione su quali sottili problemi possono sorgere quando si esegue questa randomizzazione.

In primo luogo, e soprattutto in situazioni con dati limitati, la procedura può portare ad un’elevata varianza per la stima dell’errore di test dovuto alla randomizzazione dei campioni. Questo è un tipico “trucco” quando si esegue l’approccio del set di convalida alla cross-validation. È fin troppo facile ottenere un basso errore di test semplicemente tramite un caso fortunato nel dividere appropriatamente attraverso la fortuna cieca nel ricevere una divisione del campione casuale appropriata. Quindi il vero errore di test (cioè il potere predittivo) può essere notevolmente sottostimato .

In secondo luogo, si noti che nella divisione 50-50 dei dati  addestramento/testing tralasciamo la metà di tutte le osservazioni. Quindi stiamo riducendo le informazioni che altrimenti sarebbero utilizzate per addestrare il modello. Quindi è probabile che abbia prestazioni peggiori rispetto a se avessimo usato tutte le osservazioni, comprese quelle nel set di convalida. Ciò porta a una situazione in cui potremmo effettivamente sovrastimare l’errore di test per l’intero set di dati.

Al fine di ridurre l’impatto di questi problemi, prenderemo in considerazione una suddivisione più sofisticata dei dati nota come convalida incrociata k-fold.

k-Fold Cross Validation

La convalida incrociata K-fold migliora validation set approach dividendo le \(n\) osservazioni in  \(k\) sottoinsiemi che si escludono a vicenda e di dimensioni approssimativamente uguali note come “fold”.

Il primo fold diventa un set di convalida, mentre i restanti \(k-1\) fold (aggregati insieme) diventano il set di addestramento. Il modello si adatta al set di addestramento e il suo errore di test viene stimato sul set di convalida. Questa procedura viene ripetuta \(k\) volte, con ciascuna ripetizione che offre un fold come set di convalida, mentre i restanti \(k-1\) vengono utilizzati per l’addestramento.

Ciò consente di calcolare una stima complessiva del test, \(\text{CV}_k\), che è una media di tutti i singoli errori quadratici medi, \(\text{MSE}_i\), per ogni fold:

\(\begin{eqnarray} \text {CV} _k = \frac {1} {k} \sum ^ {k} _ {i = 1} \text {MSE} _i \end{eqnarray}\)

La domanda ovvia che ci si pone in questa fase è come scegliere il valore di \(k\)? La risposta semplice (basata su studi empirici) è scegliere \(k = 5\) o \( k = 10 \). La risposta completa a questa domanda si riferisce sia alla spesa computazionale sia , ancora una volta, al compromesso bias-varianza.

Leave-One-Out Cross Validation

Possiamo effettivamente scegliere \(k = n\), cioè adattiamo il modello \(n\) volte, con una sola osservazione tralasciata per ogni adattamento. Questo è noto come validazione incrociata leave-one-out (LOOCV). Può essere molto costoso in termini di calcolo, soprattutto se \(n\) è grande e il modello ha una procedura di adattamento costosa.

Sebbene LOOCV sia vantaggioso per ridurre il bias , poiché quasi tutti i campioni vengono utilizzati per l’adattamento, in realtà soffre del problema dell’elevata varianza. Questo perché stiamo calcolando l’errore di test ogni volta su una singola risposta per ogni osservazione nel set di dati.

La convalida incrociata k-fold riduce la varianza a scapito dell’introduzione di qualche bias in più, dato che alcune delle osservazioni non sono utilizzate per l’addestramento. Con \(k = 5\) o \(k = 10\) il compromesso bias-varianza è generalmente ottimizzato.

Implementazione in Python

Siamo abbastanza fortunati quando lavoriamo con Python e il suo ecosistema di librerie, poiché gran parte del “lavoro pesante” è già stato implementato e così risparmiamo molto tempo e mal di testa! Utilizzando Pandas, Scikit-Learn e Matplotlib, possiamo creare rapidamente alcuni esempi per mostrare l’utilizzo e le problematiche relative alla convalida incrociata.

Se non hai ancora configurato un ambiente di ricerca Python, ti consiglio vivamente di scaricare il pacchetto Anaconda di Continuum Analytics che fornisce tutte le librerie che utilizzeremo in questo articolo e un IDE pronto per l’uso chiamato Spyder.

Ottenere i dati

Il primo compito è ottenere i dati e metterli in un formato che possiamo usare. In realtà abbiamo già eseguito questa procedura in un articolo precedente , ma vale la pena provare ad avere questi articoli il più autonomi possibile! Quindi, puoi utilizzare il seguente codice per ottenere dati storici di qualsiasi serie temporale finanziaria disponibile su Yahoo Finanza, nonché i valori di ritardo predittivi giornalieri associati:

import datetime
import numpy as np
import pandas as pd
import sklearn
import pandas_datareader as pdr


def create_lagged_series(symbol, start_date, end_date, lags=5):
    """
    Si crea un DataFrame pandas che memorizza i rendimenti percentuali dei
    prezzi di chiusura rettificata di un titolo azionario ottenuta da Yahoo
    Finance, insieme a una serie di rendimenti ritardati dai giorni di negoziazione
    precedenti (i valori predefiniti ritardano di 5 giorni).
    Sono inclusi anche il volume degli scambi, così come la direzione del giorno
    precedente.
    """

    # Ottieni informazioni sulle azioni da Yahoo Finance
    ts = pdr.get_data_yahoo(symbol, start_date-datetime.timedelta(days=365), end_date)

    # Crea un nuovo Dataframe per i ritardi
    tslag = pd.DataFrame(index=ts.index)
    tslag["Today"] = ts["Adj Close"]
    tslag["Volume"] = ts["Volume"]

    # Crea la serie traslata dei ritardi dei prezzi di chiusura del 
    # periodo (giorno) precedente
    for i in range(0,lags):
        tslag["Lag%s" % str(i+1)] = ts["Adj Close"].shift(i+1)

    # Crea il DataFrame dei ritorni
    tsret = pd.DataFrame(index=tslag.index)
    tsret["Volume"] = tslag["Volume"]
    tsret["Today"] = tslag["Today"].pct_change()*100.0

    # Se qualsiasi valore dei ritorni percentuali è uguale a zero, si imposta 
    # a un numero piccolo (per non avere problemi con il QDA in scikit-learn)
    for i,x in enumerate(tsret["Today"]):
        if (abs(x) < 0.0001):
            tsret["Today"][i] = 0.0001

    # Crea la serie dei ritorni precedenti percentuali
    for i in range(0,lags):
        tsret["Lag%s" % str(i+1)] = tslag["Lag%s" % str(i+1)].pct_change()*100.0

    # Crea la serie "Direction" (+1 o -1) che indica un giorno up/down
    tsret["Direction"] = np.sign(tsret["Today"])
    tsret = tsret[tsret.index >= start_date]

    return tsret

Da notare che non memorizziamo i valori del prezzo di chiusura direttamente nelle colonne “Today” o “Lags”. Invece, stiamo memorizzando il rendimento percentuale tra i prezzi di chiusura di un giorno e del giorno precedente.

Dobbiamo ottenere i dati per i prezzi giornalieri del SP500 per un periodo di tempo adeguato. Si considera dal 1 ° gennaio 2004 al 31 dicembre 2004. Tuttavia questa è una scelta arbitraria. Si può regolare l’intervallo di tempo come meglio si crede. Per ottenere i dati e inserirli in un Pandas DataFrame chiamato sp500_lags possiamo utilizzare il seguente codice:

if __name__ == "__main__":
    symbol = "^GSPC"
    start_date = datetime.datetime(2004, 1, 1)
    end_date = datetime.datetime(2004, 12, 31)
    sp500_lags = create_lagged_series(symbol, start_date, end_date, lags=5)

A questo punto abbiamo i dati necessari per iniziare a creare una serie di modelli statistici di machine learning..

Validation Set Approach

Ora che abbiamo i dati finanziari necessari per creare una serie di modelli di regressione predittiva, possiamo utilizzare i metodi di convalida incrociata sopra riportati per ottenere stime per l’errore di test.

Il primo compito è importare i modelli da Scikit-Learn. Scegliamo un modello di regressione lineare con caratteristiche polinomiali. Questo ci fornisce la possibilità di scegliere diversi gradi di flessibilità semplicemente aumentando il grado dell’ordine polinomiale delle features. Inizialmente si considera l’approccio del set di convalida per la cross-validation.

Scikit-Learn fornisce un approccio a set di convalida tramite il metodo train_test_split trovato nel modulo cross_validation. Successivamente dobbiamo importare il metodo KFold per la convalida incrociata k-fold, così come il modello di regressione lineare stesso. Dobbiamo importare il calcolo MSE così come Pipeline PolynomialFeatures. Gli ultimi due metodi ci consentono di creare facilmente un insieme di modelli di regressione lineare di feature polinomiali con una minima codifica aggiuntiva:

..
from sklearn.model_selection import train_test_split, KFold
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_squared_error
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import PolynomialFeatures
..

Una volta importati i moduli, possiamo creare un DataFrame SP500 che utilizza i rendimenti in ritardo dei cinque giorni precedenti come predittori. Possiamo quindi creare dieci separate suddivisioni casuali dei dati in un set di addestramento e convalida.

Infine, per gradi multipli delle features polinomiali della regressione lineare, possiamo calcolare l’errore di test. Questo ci fornisce dieci separate curve di errore di test, ogni valore delle quali mostra il test MSE per un grado diverso del kernel polinomiale:

..
..

def validation_set_poly(random_seeds, degrees, X, y):
    """
    Utilizza il metodo train_test_split per creare un set
    di addestramento e un set di convalida (50% per ciascuno)
    utilizzando separati campionamenti casuali "random_seeds"
    per modelli di regressione lineare di varia flessibilità
    """
    sample_dict = dict([("seed_%s" % i,[]) for i in range(1, random_seeds+1)])
    # Esegui un ciclo su ogni suddivisione casuale in una suddivisione train-test
    for i in range(1, random_seeds+1):
        print("Random: %s" % i)
        # Aumenta il grado di ordine polinomiale della regressione lineare
        for d in range(1, degrees+1):
            print("Degree: %s" % d)
            # Crea il modello, divide gli insiemi e li addestra
            polynomial_features = PolynomialFeatures(
                degree=d, include_bias=False
            )
            linear_regression = LinearRegression()
            model = Pipeline([
                ("polynomial_features", polynomial_features),
                ("linear_regression", linear_regression)
            ])
            X_train, X_test, y_train, y_test = train_test_split(
                X, y, test_size=0.5, random_state=i
            )
            model.fit(X_train, y_train)
            # Calcola il test MSE e lo aggiunge al
            # dizionario di tutte le curve di test
            y_pred = model.predict(X_test)
            test_mse = mean_squared_error(y_test, y_pred)
            sample_dict["seed_%s" % i].append(test_mse)
        # Converte queste liste in array numpy per calcolare la media
        sample_dict["seed_%s" % i] = np.array(sample_dict["seed_%s" % i])
    # Crea la serie delle "medie dei test MSE" colcolando la media
    # del test MSE per ogni grado dei modelli di regressione lineare,
    # attraverso tutti i campionamenti casuali
    sample_dict["avg"] = np.zeros(degrees)
    for i in range(1, random_seeds+1):
        sample_dict["avg"] += sample_dict["seed_%s" % i]
    sample_dict["avg"] /= float(random_seeds)
    return sample_dict

..
..

Possiamo usare Matplotlib per tracciare il grafico di questi dati. Dobbiamo importare pylab e quindi creare una funzione per tracciare le curve di errore di test:

..
import pylab as plt
..
..
def plot_test_error_curves_vs(sample_dict, random_seeds, degrees):
    fig, ax = plt.subplots()
    ds = range(1, degrees+1)
    for i in range(1, random_seeds+1):
        ax.plot(ds, sample_dict["seed_%s" % i], lw=2,
                label='Test MSE - Sample %s' % i)

    ax.plot(ds, sample_dict["avg"], linestyle='--', color="black",
                    lw=3, label='Avg Test MSE')
    ax.legend(loc=0)
    ax.set_xlabel('Degree of Polynomial Fit')
    ax.set_ylabel('Mean Squared Error')
    ax.set_ylim([0.0, 4.0])
    fig.set_facecolor('white')
    plt.show()
..
..

Abbiamo selezionato il grado delle nostre features polinomiali al variare tra \(d = 1\) e \(d = 3\), prevedendo così un ordine cubico nelle nostre features. La seguente Figura 1 mostra le dieci diverse suddivisioni casuali dei dati di addestramento e test, insieme alla media del test MSE (la linea tratteggiata nera):

trading-algoritmico-cross-val-fig1
Figura 1 - Le curve di test MSE per più divisioni di convalida dell'addestramento per una regressione lineare con features polinomiali di grado crescente.

È immediatamente evidente quanta variazione ci sia tra diverse suddivisioni casuali in un set di addestramento e convalida. Poiché non c’è una grande quantità di segnale predittivo nell’utilizzo dei prezzi di chiusura storici dei giorni precedenti del SP500, vediamo che all’aumentare del grado delle features polinomiali, effettivamente il test MSE aumenta.

Inoltre è chiaro che il set di convalida soffre di una elevata varianza. Il test MSE medio per l’approccio del set di validazione sul modello di grado \(d = 3\) è di circa 1,9.

Per ridurre al minimo questo problema, implementeremo ora la convalida incrociata k-fold sullo stesso set di dati dell’SP500.

k-Fold Cross Validation

Poiché ci siamo già occupati delle importazioni di cui sopra, ci limitiamo a delineare le nuove funzioni per eseguire la convalida incrociata k-fold. Sono quasi identiche alle funzioni utilizzate per la divisione del test di addestramento. Tuttavia, dobbiamo usare l’oggetto  KFold per iterare su \(k \) “fold”.

In particolare l’oggetto KFold fornisce un iteratore che ci permette di indicizzare correttamente i campioni nel data set e creare fold separati di training/test. In questo esempio abbiamo scelto \(k = 10\).

Come con l’approccio dell’insieme di convalida, creiamo una pipeline di trasformazione delle features polinomiali e applichiamo un modello di regressione lineare. Quindi calcoliamo il test MSE e costruiamo curve di test MSE separate per ogni fold. Infine, creiamo una curva MSE media tra le fold:

..
..
def k_fold_cross_val_poly(folds, degrees, X, y):
    n = len(X)
    kf = KFold(n, n_folds=folds)
    kf_dict = dict([("fold_%s" % i,[]) for i in range(1, folds+1)])
    fold = 0
    for train_index, test_index in kf:
        fold += 1
        print("Fold: %s" % fold)
        X_train, X_test = X.ix[train_index], X.ix[test_index]
        y_train, y_test = y.ix[train_index], y.ix[test_index]
        # Aumenta il grado di ordine polinomiale della regressione lineare
        for d in range(1, degrees+1):
            print("Degree: %s" % d)
            # Crea il modello e lo addestra
            polynomial_features = PolynomialFeatures(
                degree=d, include_bias=False
            )
            linear_regression = LinearRegression()
            model = Pipeline([
                ("polynomial_features", polynomial_features),
                ("linear_regression", linear_regression)
            ])
            model.fit(X_train, y_train)
            # Calcola il test MSE e lo aggiunge al 
            # dizionario di tutte le curve di test
            y_pred = model.predict(X_test)
            test_mse = mean_squared_error(y_test, y_pred)
            kf_dict["fold_%s" % fold].append(test_mse)
        # Converte queste liste in array numpy per calcolare la media
        kf_dict["fold_%s" % fold] = np.array(kf_dict["fold_%s" % fold])
    # Crea la serie dei "test MSE medi" calcolando la media dei 
    # test MSE per ogni grado del modello di regressione lineare,
    # tramite ogni k folds.
    kf_dict["avg"] = np.zeros(degrees)
    for i in range(1, folds+1):
        kf_dict["avg"] += kf_dict["fold_%s" % i]
    kf_dict["avg"] /= float(folds)
    return kf_dict
..
..

Possiamo tracciare queste curve con la seguente funzione:

..
..
def plot_test_error_curves_kf(kf_dict, folds, degrees):
    fig, ax = plt.subplots()
    ds = range(1, degrees+1)
    for i in range(1, folds+1):
        ax.plot(ds, kf_dict["fold_%s" % i], lw=2, label='Test MSE - Fold %s' % i)

    ax.plot(ds, kf_dict["avg"], linestyle='--', color="black", 
                                lw=3, label='Avg Test MSE')
    ax.legend(loc=0)
    ax.set_xlabel('Degree of Polynomial Fit')
    ax.set_ylabel('Mean Squared Error')
    ax.set_ylim([0.0, 4.0])
    fig.set_facecolor('white')
    plt.show()
..
..

L’output è riportato nella seguente Figura 2:

trading-algoritmico-cross-val-fig2
Figura 2 - Curve dei test MSE per più fold di convalida incrociata k-fold per una regressione lineare con features polinomiali di grado crescente.

Si noti che la variazione tra le curve di errore è molto inferiore rispetto al validation set approch. Questo è l’effetto desiderato dell’esecuzione della convalida incrociata. In particolare, per \(d = 3\) abbiamo un errore di test medio ridotto di circa 0,8.

La convalida incrociata fornisce generalmente una stima molto migliore del vero test MSE, a scapito di qualche lieve bias. Questo di solito è un compromesso accettabile nelle applicazioni di machine learning.

Negli articoli futuri prenderemo in considerazione approcci di ricampionamento alternativi , inclusi Bootstrap, Bootstrap Aggregation (“Bagging”) e Boosting. Si tratta di tecniche più sofisticate che ci aiuteranno a selezionare meglio i nostri modelli e (si spera) a ridurre ulteriormente i nostri errori.

Codice Python Completo

Di seguito il codice Python completo per il file cross_validation.py:

import datetime
import numpy as np
import pandas as pd
import sklearn
import pandas_datareader as pdr
import pylab as plt

from sklearn.model_selection import train_test_split, KFold
from sklearn.linear_model import LinearRegression
from sklearn.metrics import mean_squared_error
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import PolynomialFeatures



def create_lagged_series(symbol, start_date, end_date, lags=5):
    """
    Si crea un DataFrame pandas che memorizza i rendimenti percentuali dei
    prezzi di chiusura rettificata di un titolo azionario ottenuta da Yahoo
    Finance, insieme a una serie di rendimenti ritardati dai giorni di negoziazione
    precedenti (i valori predefiniti ritardano di 5 giorni).
    Sono inclusi anche il volume degli scambi, così come la direzione del giorno precedente.
    """

    # Ottieni informazioni sulle azioni da Yahoo Finance
    ts = pdr.get_data_yahoo(symbol, start_date-datetime.timedelta(days=365), end_date)

    # Crea un nuovo Dataframe per i ritardi
    tslag = pd.DataFrame(index=ts.index)
    tslag["Today"] = ts["Adj Close"]
    tslag["Volume"] = ts["Volume"]

    # Crea la serie traslata dei ritardi dei prezzi di chiusura del periodo (giorno) precedente
    for i in range(0,lags):
        tslag["Lag%s" % str(i+1)] = ts["Adj Close"].shift(i+1)

    # Crea il DataFrame dei ritorni
    tsret = pd.DataFrame(index=tslag.index)
    tsret["Volume"] = tslag["Volume"]
    tsret["Today"] = tslag["Today"].pct_change()*100.0

    # Se uno qualsiasi dei valori dei ritorni percentuali è uguale a zero, si impostano
    # a un numero piccolo (per non avere problemi con il modello QDA in scikit-learn)
    for i,x in enumerate(tsret["Today"]):
        if (abs(x) < 0.0001):
            tsret["Today"][i] = 0.0001

    # Crea la serie dei ritorni precedenti percentuali
    for i in range(0,lags):
        tsret["Lag%s" % str(i+1)] = tslag["Lag%s" % str(i+1)].pct_change()*100.0

    # Crea la serie "Direction" (+1 o -1) che indica un giorno up/down
    tsret["Direction"] = np.sign(tsret["Today"])
    tsret = tsret[tsret.index >= start_date]

    return tsret


def validation_set_poly(random_seeds, degrees, X, y):
    """
    Utilizza il metodo train_test_split per creare un set
    di addestramento e un set di convalida (50% per ciascuno)
    utilizzando separati campionamenti casuali "random_seeds"
    per modelli di regressione lineare di varia flessibilità
    """
    sample_dict = dict([("seed_%s" % i,[]) for i in range(1, random_seeds+1)])
    # Esegui un ciclo su ogni suddivisione casuale in una suddivisione train-test
    for i in range(1, random_seeds+1):
        print("Random: %s" % i)
        # Aumenta il grado di ordine polinomiale della regressione lineare
        for d in range(1, degrees+1):
            print("Degree: %s" % d)
            # Crea il modello, divide gli insiemi e li addestra
            polynomial_features = PolynomialFeatures(
                degree=d, include_bias=False
            )
            linear_regression = LinearRegression()
            model = Pipeline([
                ("polynomial_features", polynomial_features),
                ("linear_regression", linear_regression)
            ])
            X_train, X_test, y_train, y_test = train_test_split(
                X, y, test_size=0.5, random_state=i
            )
            model.fit(X_train, y_train)
            # Calcola il test MSE e lo aggiunge al
            # dizionario di tutte le curve di test
            y_pred = model.predict(X_test)
            test_mse = mean_squared_error(y_test, y_pred)
            sample_dict["seed_%s" % i].append(test_mse)
        # Converte queste liste in array numpy per calcolare la media
        sample_dict["seed_%s" % i] = np.array(sample_dict["seed_%s" % i])
    # Crea la serie delle "medie dei test MSE" colcolando la media
    # del test MSE per ogni grado dei modelli di regressione lineare,
    # attraverso tutti i campionamenti casuali
    sample_dict["avg"] = np.zeros(degrees)
    for i in range(1, random_seeds+1):
        sample_dict["avg"] += sample_dict["seed_%s" % i]
    sample_dict["avg"] /= float(random_seeds)
    return sample_dict


def plot_test_error_curves_vs(sample_dict, random_seeds, degrees):
    fig, ax = plt.subplots()
    ds = range(1, degrees+1)
    for i in range(1, random_seeds+1):
        ax.plot(ds, sample_dict["seed_%s" % i], lw=2, label='Test MSE - Sample %s' % i)
    ax.plot(ds, sample_dict["avg"], linestyle='--', color="black", lw=3, label='Avg Test MSE')
    ax.legend(loc=0)
    ax.set_xlabel('Degree of Polynomial Fit')
    ax.set_ylabel('Mean Squared Error')
    ax.set_ylim([0.0, 4.0])
    fig.set_facecolor('white')
    plt.show()


def k_fold_cross_val_poly(folds, degrees, X, y):
    n = len(X)
    kf = KFold(n, n_folds=folds)
    kf_dict = dict([("fold_%s" % i,[]) for i in range(1, folds+1)])
    fold = 0
    for train_index, test_index in kf:
        fold += 1
        print("Fold: %s" % fold)
        X_train, X_test = X.ix[train_index], X.ix[test_index]
        y_train, y_test = y.ix[train_index], y.ix[test_index]
        # Aumenta il grado di ordine polinomiale della regressione lineare
        for d in range(1, degrees+1):
            print("Degree: %s" % d)
            # Crea il modello e lo addestra
            polynomial_features = PolynomialFeatures(
                degree=d, include_bias=False
            )
            linear_regression = LinearRegression()
            model = Pipeline([
                ("polynomial_features", polynomial_features),
                ("linear_regression", linear_regression)
            ])
            model.fit(X_train, y_train)
            # Calcola il test MSE e lo aggiunge al
            # dizionario di tutte le curve di test
            y_pred = model.predict(X_test)
            test_mse = mean_squared_error(y_test, y_pred)
            kf_dict["fold_%s" % fold].append(test_mse)
        # Converte queste liste in array numpy per calcolare la media
        kf_dict["fold_%s" % fold] = np.array(kf_dict["fold_%s" % fold])
    # Crea la serie dei "test MSE medi" calcolando la media dei
    # test MSE per ogni grado del modello di regressione lineare,
    # tramite ogni k folds.
    kf_dict["avg"] = np.zeros(degrees)
    for i in range(1, folds+1):
        kf_dict["avg"] += kf_dict["fold_%s" % i]
    kf_dict["avg"] /= float(folds)
    return kf_dict


def plot_test_error_curves_kf(kf_dict, folds, degrees):
    fig, ax = plt.subplots()
    ds = range(1, degrees+1)
    for i in range(1, folds+1):
        ax.plot(ds, kf_dict["fold_%s" % i], lw=2, label='Test MSE - Fold %s' % i)
    ax.plot(ds, kf_dict["avg"], linestyle='--', color="black", lw=3, label='Avg Test MSE')
    ax.legend(loc=0)
    ax.set_xlabel('Degree of Polynomial Fit')
    ax.set_ylabel('Mean Squared Error')
    ax.set_ylim([0.0, 4.0])
    fig.set_facecolor('white')
    plt.show()



if __name__ == "__main__":
    symbol = "^FTSE"
    symbol = "^GSPC"
    start_date = datetime.datetime(2004, 1, 1)
    end_date = datetime.datetime(2004, 12, 31)
    sp500_lags = create_lagged_series(symbol, start_date, end_date, lags=5)

    # Uso tutti e venti i ritorni di 2 giorni precedenti come  
    # valori di predizione, con "Today" come risposta
    X = sp500_lags[[
        "Lag1", "Lag2", "Lag3", "Lag4", "Lag5",
        # "Lag6", "Lag7", "Lag8", "Lag9", "Lag10",
        # "Lag11", "Lag12", "Lag13", "Lag14", "Lag15",
        # "Lag16", "Lag17", "Lag18", "Lag19", "Lag20"
    ]]
    y = sp500_lags["Today"]
    degrees = 3

    # Visualizza le curve dell'errore di test per il set di validazione
    random_seeds = 10
    sample_dict_val = validation_set_poly(random_seeds, degrees, X, y)
    plot_test_error_curves_vs(sample_dict_val, random_seeds, degrees)

    # Visualizza le curve dell'errore di test per il set di k-fold CV
    folds = 10
    kf_dict = k_fold_cross_val_poly(folds, degrees, X, y)
    plot_test_error_curves_kf(kf_dict, folds, degrees)

Recommended Posts