In questo articolo descriviamo come effettuare l’analisi predittiva delle serie temporali in python per verificare se esiste un modo semplice per migliorare le capacità e le abilità di previsione in modo da produrre risultati più accurati. Quando analizziamo i dati storici dei prezzi per la maggioranza degli asset finanziari, in vari intervalli di tempo (alcuni più lunghi, altri più brevi) molti set di dati che possono apparire completamente casuali. Sembrano abbastanza casuali da non permetterci di avere facilmente una probabile previsione dei valori e degli andamenti futuri. Per ulteriori approfondimenti sul tema delle Time Series è possibile consultare i tutorial presenti su DataTrading.info.
Sebbene trovare una soluzione per l’analisi predittiva delle serie temporali è molto difficile, ci sono alcuni strumenti che potrebbero aiutarci, come la scomposizione delle serie temporali.
Le componenti delle serie temporali
La scomposizione delle serie temporali è una tecnica che permette di suddividere una serie storica nelle sue singole “componenti” elementari. Una serie temporale è formata da al massimo 4 tipi di componenti diversi:
- Componente di trend T.
- Componente della stagionalità S.
- Componente ciclica C.
- Componente residua o erratica R.
Non è necessario che una serie temporale contenga tutte e 4 le componenti, potrebbe mancare le componente stagionale o di trend.
Prevedere correttamente una serie temporale è spesso piuttosto difficile. Se possiamo scomporre la serie in singole componenti e analizzarle separatamente, possiamo migliorare la previsione complessiva. Ogni componente ha proprietà e comportamenti specifici quindi possiamo usare i metodi ed approcci più adatti e performanti per ogni specifica componente.
Da notare che la scomposizione è principalmente usate per aiutare ad analizzare e comprendere le serie storiche, ma può anche rivelarsi utile quando effettuiamo un’analisi predittiva.
Come accennato in precedenza, una serie può essere scomposta in un alcune di queste componenti: Trend, Stagionale, Ciclica e Residua. È utile pensare alle componenti come elementi combinati in modo additivo o moltiplicativo.
Differiscono come segue:
- Un modello additivo suggerisce che i componenti siano sommati come segue: \(y_{t} = T_{t} + S_{t} + C_{t} + R_{t}\) mentre un modello moltiplicativo suggerisce che i componenti siano moltiplicati insieme come segue: \(y_{t} = T_{t} \times S_{t} \times C_{t} \times N_{t}\)
- Un modello additivo è lineare quindi le variazioni nel tempo sono costantemente apportate nella stessa quantità, mentre un modello moltiplicativo è non lineare, può essere quadratico o esponenziale, e le variazioni aumentano o diminuiscono nel tempo.
- Un trend lineare è una linea retta, mentre un trend non lineare è una linea curva.
- Una stagionalità lineare ha la stessa frequenza (periodicità dei cicli) e ampiezza (altezza dei cicli), mentre una stagionalità non lineare ha una frequenza e/o un’ampiezza crescente o decrescente nel tempo.
Sottolineiamo che le precedenti sono regole generali. Nella realtà non è possibile scomporre una serie temporale in modo pulito o perfetto in modello additivo o moltiplicativo. Possono esistere combinazioni dei due modelli, i trend possono cambiare direzione e i cicli di natura non ripetitiva possono mescolarsi con le componenti ricorrenti della stagionalità. Fondamentalmente i dati reali possono essere disordinati e non sempre rispettano le regole! Sono comunque utili nella definizione di un semplice framework da cui analizzare i nostri dati.
Decomposizione classica
Esistono vari metodi di decomposizione, con il metodo “base” noto come “decomposizione classica”. E’ una procedura relativamente semplice che costituisce anche il punto di partenza per la maggior parte degli altri metodi di decomposizione delle serie temporali. Questo metodo esiste da anni, essendo nato negli anni ’20, quindi non si tratta esattamente di qualcosa di “nuovo”.
Esistono due forme di decomposizione classica, una per ciascuno dei due modelli (additivo e moltiplicativo).
Per la decomposizione additiva il processo (assumendo un periodo stagionale M) prevede i seguenti passaggi:
- Calcolare la componente “trend-cycle” \(\hat{T}_{t}\) tramite una \(2 \times m \textrm-MA\) se M è un numero pari, oppure tramite \(m \textrm-MA\) se M è un numero dispari.
- Calcolare la serie detrendizzata: \(y_{t} – \hat{T}_{t}\).
- Calcolare la componente stagionale per ogni periodo M come la media dei valori detrendizzati di quel periodo. In altre parole, se usiamo dati mensili, la componente stagionale per un mese corrisponde alla media di tutti i valori detrendizzati di quel mede di dati. Successivamente si modifica i valori mensili per garantire che la loro somma sia pari a zero, quindi si unisco i valori di ogni periodo e si replica la sequenza per ogni anno di dati. Si ottiene l’output \(\hat{S}_{t}\)
- Calcolare la componente residua sottraendo le componenti stimante della “stagionalità” e “trend-cycle”: \(\hat{R}_{t} = y_{t} – \hat{S}_{t}\)
La decomposizione moltiplicativa classica è simile. In questo caso usiamo le divisioni al posto delle sottrazioni.
Sebbene la decomposizione classica sia ancora ampiamente utilizzata, non è raccomandata, poiché è soggetta ad alcune criticità. Ad esempio l’assenza di stime del ciclo di trend per le prime e ultime osservazioni (ad esempio se m=12 non ci sarebbe alcuna stima del ciclo di trend per il primo e gli ultimi 6 periodi di tempo, e quindi nemmeno la componente residua). Inoltre ci sono criticità relative alla ipotesi che la componente stagionale si ripete di anno in anno. Sotto questa ipotesi si tende a smussare eccessivamente i repentini aumenti e diminuzioni dei dati durante il calcolo del ciclo di trend. Infine il metodo non è robusto per “valori insoliti” o “outlier” nei dati. Ad ogni modo, oggi abbiamo a disposizione metodi migliori come le X11 Decomposition, SEATS Decomposition e STL Decomposition.
Decomposizione STL
La Decomposizione STL, cioè la “Seasonal and Trend decomposition using Loess”. è un metodo per stimare le relazioni non lineari. Questo metodo presenta alcuni vantaggi rispetto ai metodi di decomposizione classici, SEATS e X11:
- STL gestisce qualsiasi tipo di stagionalità, non solo mensile o trimestrale (a differenza delle SEATS e X11).
- L’uniformità del trend-cycle può essere controllata dall’utente.
- La componente stagionale può cambiare nel tempo, consentendo nuovamente all’utente di controllare questa variabile.
- È generalmente più robusta per i dati “insoliti” o anomali. Questi dati non influenzano le stime delle componenti trend-cycle e stagionali, tuttavia influenzano la componente residua.
Questa decomposizione presenta anche alcuni svantaggi, ad esempio può essere usata solo con la scomposizione additiva e non gestisce automaticamente la variazione del giorno di negoziazione/calendario.
Vediamo come implementare il codice python per effettuare la decomposizione tramite la libreria statsmodels. In questo esempio usiamo la serie temporale dei prezzi del mercato EURUSD. I dati storici usate possono essere scaricati al seguente link:
Dopo aver preparato i dati storici passiamo ad installare le librerie e moduli python di cui abbiamo bisogno. Prima di tutto abbiamo dobbiamo installare la classe STL presente nel modulo statsmodels in un ambiente virtuale. Quindi dobbiamo seguire i seguenti passaggi:
- Creare e attivare un ambiente virtuale.
- Installare Cython nell’ambiente virtuale
pip install cython - Installare la libreria statsmodels con il comando
pip install statsmodels
Successivamente andiamo direttamente nel codice principale ed importiamo i moduli necessari. Quindi leggiamo i dati EURUSD utilizzando Pandas, estraiamo i prezzi che ci interessano e applichiamo su di essi un filtro Hodrick-Prescott e una decomposizione SYL. Sottolineiamo che il file EURUSD fornito è “separato da tabulazioni” anziché il più comune “separato da virgole”, quindi impostiamo il parametro sep="\t" nella funzione di Pandas.
import pandas as pd
import statsmodels.api as sm
from statsmodels.tsa.seasonal import STL
from sklearn.metrics import mean_squared_error
from statsmodels.tsa.ar_model import AR
import matplotlib.pyplot as plt
import seaborn as sns
plt.rcParams['figure.figsize'] = [12, 9]
df = pd.read_csv('EURUSD.csv',sep='\t', index_col='Date')
df.index = pd.to_datetime(df.index)
df.sort_index(inplace=True)
df = df.resample('W').last()
series = df['Price']
Il filtro HP è una tecnica comunemente usata con le serie macroeconomiche che hanno un trend (movimenti a lungo termine), un ciclo economico e parti irregolari (fluttuazioni a breve termine). Costruisce la componente trend risolvendo un problema di ottimizzazione. Permette di ottenere una stima della componente di trend più uniforme, riducendo al minimo la distanza euclidea dalla serie originale. In altre parole, trova un l’equilibrio tra l’uniformità del trend e la sua vicinanza all’originale.
Per approfondire i concetti teorici si può consultare la documentazione sul filtro di Hodrick–Prescott.
Implementiamo quindi la scomposizione tramite il filtro HP.
cycle, trend = sm.tsa.filters.hpfilter(series, 50)
fig, ax = plt.subplots(3,1)
ax[0].plot(series)
ax[0].set_title('Price')
ax[1].plot(trend)
ax[1].set_title('Trend')
ax[2].plot(cycle)
ax[2].set_title('Cycle')
plt.show()
Quindi implementiamo la decomposizione SYL e visualizziamo le componenti ottenute.
result = STL(series).fit()
chart = result.plot()
plt.show()
A partire dai dati dei prezzi settimanali abbiamo ottenuto le serie delle componenti di Trend, di Stagionalità e del Residuo. Per iniziare l’analisi predittiva delle serie temporali possiamo provare a prevedere i movimenti dei prezzi EURUSD. In questo caso possiamo iniziare ad applicare un semplice “modello di persistenza”. Questa tipologia di modello consiste nell’assegnare l’ultimo valore osservato come previsione per il valore successivo, cioè si “persiste” in avanti l’ultimo prezzo. Non è un modello molto sofisticato, ma è un inizio e fornisce un’idea di base delle prestazioni che possiamo usare per confrontare il modello di regressione automatica. Valutiamo le prestazioni del modello tramite il Root Mean Squared Error (RMSE).
Analisi predittiva con il modello di persistenza
Nel seguente codice creiamo una serie di previsioni “1-period forward”, semplicemente spostando in avanti di una settimana l’ultimo prezzo e confrontando quel valore con il prezzo effettivo che è stato registrato in quel momento. Usiamo quindi la funzione mean_squared_error di scikit-learn per calcolare l’MSE, di cui prendiamo solo la radice quadrata per produrre l’RMSE.
Vediamo che l’RMSE è calcolato solo per il periodo corrispondente a un “dataset di test”, che in questo caso è l’ultimo 30% dei dati. In questo modo possiamo confrontare le prestazioni con modelli più sofisticati che richiedono un “addestramento” su un set di dati di “training”, prima di usarli per la previsione. Dobbiamo quindi confrontare le prestazioni delle previsioni sul set di dati di “test”.
In questo modo possiamo valutare le prestazioni di un modello quando è usato per l’analisi predittiva delle serie temporali solamente con i dati che non sono stati usati per l’ottimizzazione dei parametri.
Calcoliamo l’errore quadratico medio, o MSE, come media dei valori quadratici dell’errore di previsione. I valori dell’errore sono in unità quadrate dei valori previsti e un errore quadratico medio pari a zero indica capacità di previsione perfette o “nessun errore”.
L’MSE può essere ritrasformato nelle unità originali delle previsioni calcolando la sua radice quadrata, ottenendo quindi un RMSE. Ovviamente anche un RMSE pari a zero che indica ancora una volta capacità di previsione perfette.
Infine, tracciamo solo gli ultimi punti dati, permettendoci di vedere più chiaramente come i valori previsti si relazionano ai valori effettivamente osservati.
predictions = series.shift(1).dropna()
test_score = np.sqrt(mean_squared_error(series[int(len(series) * 0.7)+1:], predictions.iloc[int(len(series) * 0.7):]))
print('Test RMSE: %.5f' % test_score)
plt.plot(series.iloc[-25:], label='Price')
plt.plot(predictions[-25:], color='red', label='Prediction')
plt.legend()
plt.show()
Creiamo quindi un grafico a dispersione relativo alla variazione percentuale settimanale prevista rispetto a quella effettiva per verificare se esiste una relazione evidente. Inoltre calcoliamo l’errore medio assoluto (MAE) delle previsioni, cioè la differenza media tra i valori di previsione il valore osservato.
price_pred = pd.concat([series.iloc[-int(len(series) * 0.3):].pct_change(),
predictions.iloc[-int(len(series) * 0.3):].pct_change()], axis=1)
price_pred.dropna(inplace=True)
price_pred.columns = ['Price', 'preds']
fig, ax = plt.subplots()
ax = sns.regplot(data=price_pred, x=price_pred['Price'], y=price_pred['preds'])
plt.xlabel('Observations')
plt.ylabel('Predictions')
plt.title('EURUSD Observed vs Predicted Values')
ax.grid(True, which='both')
ax.axhline(y=0, color='#888888')
ax.axvline(x=0, color='#888888')
sns.despine(ax=ax, offset=0)
plt.xlim(-0.05, 0.05)
plt.ylim(-0.05, 0.05)
plt.show()
mae = round(abs(price_pred['Price'] - price_pred['preds']).mean(),4)
print(f'The MAE is {mae}')
Guardando il grafico sembra che non ci sia alcuna relazione tra le previsioni e i valori osservati quando consideriamo la variazione percentuale settimanale del prezzo EURUSD. Non è un evento eccezionale!
Se calcoliamo il “tasso di successo” delle probabilità di prevedere correttamente la direzione del movimento dell’EURUSD per prossima settimana, otteniamo un valore molto vicino al 50%, cioè come il lancio di una moneta.
price_pred['hit'] = np.where(np.sign(price_pred['Price']) == np.sign(price_pred['preds']), 1, 0)
print(f"Hit rate: {round((price_pred['hit'].sum() / price_pred['hit'].count()) * 100,2)}%")
Hit rate: 48.49%
Possiamo ora confrontare il punteggio del modello “di base” con i successivi modelli più complessi. Proviamo a calcolare le previsioni usando un modello autoregressivo (AR), cioè un un modello di regressione lineare che utilizza i valori ritardati delle variabili come input delle variabili successive.
Analisi predittiva con il modello Autoregressivo
Per l’analisi predittiva delle serie temporali con il modello autoregressivo dobbiamo definire il valore del ritardo (lag) da usare per la variabile di input. L’approccio più semplice è usare la classe AutoReg della libreria statsmodels. Questa classe procede all’addestramento di un modello di regressione lineare.
Nel seguente codice suddividiamo la serie dei prezzi EURUSD in set di dati di “addestramento” e “test” e creiamo una lista vuota per memorizzare le previsioni. Dobbiamo usare una lista perchè in questo approccio calcoliamo le previsioni settimanali una ad una, in modo “walk-forward”.
Eseguiamo l’iterazione attraverso un “ciclo for” con un numero di iterazioni pari alla lunghezza dei dati di test. Ad ogni iterazioni creiamo il modello AR e lo alimentiamo con il set di dati di addestramento per il training e il fitting del modello. A questo punto generiamo una previsione di 1 periodo in avanti (che assume la forma di un singolo valore) e aggiungiamo il risultato alla lista delle previsioni. Inoltre dai dati di test estraiamo il valore nella posizione corrispondente al numero dell’iterazione del ciclo e aggiungiamo il valore alla fine dei dati di addestramento.
Nell’esecuzione successiva del ciclo, addestriamo e adattiamo nuovamente il modello, ma questa volta abbiamo un’osservazione in più nei dati di addestramento. Questa osservazione è il valore che abbiamo aggiunto durante la precedente iterazione del ultimo ciclo for.
Con questo approccio i dati di addestramento aumentano man mano dato che aggiungiamo una nuova osservazione ad ogni ciclo. In questo modo il modello imita il passaggio nel tempo della vita reale dove ogni giorno abbiamo a disposizione un dato in più rispetto al giorno precedente. Niente di magico! Garantiamo semplicemente che non stiamo esponendo il modello a bias “previsionali” e addestrandolo usando dati che non sarebbero stati disponibili in quel momento.
Come prima, calcoliamo l’MSE e tracciamo gli ultimi 25 punti dati: i prezzi osservati rispetto alle previsioni fatte dal modello.
historic = series.iloc[:int(len(series) * 0.7)].to_list()
test = series.iloc[int(len(series) * 0.7):]
predictions = []
for i in range(len(test)):
model = AutoReg(historic, lags=10)
model_fit = model.fit()
pred = model_fit.predict(start=len(historic), end=len(historic), dynamic=False)
predictions.append(pred[0])
historic.append(test[i])
predictions = pd.Series(predictions, index=test.index)
test_score = np.sqrt(mean_squared_error(test, predictions))
print('Test MSE: %.5f' % test_score)
# Grafico dei risultati
plt.plot(test.iloc[-25:], label='Prices')
plt.plot(predictions.iloc[-25:], color='red', label='Prediction')
plt.legend()
plt.show()
Anche in questo caso otteniamo un MSE molto vicino a 0,013, in realtà è leggermente più alto pari a 0,1318. I due modelli producono un MSE così simile perchè nessuno dei due modelli sembra essere in grado di fare previsioni. Se calcoliamo il movimento settimanale medio del prezzo dell’EURUSD nel periodo coperto dai dati di test, ho il vago sospetto che otteniamo un valori molto vicino a… avete indovinato – 0,013.
È come dire che le previsioni del modello non sono migliori, in media, al supporre che il prezzo non si muoverà, che è ovviamente una previsione inutile.
Tracciare il grafico a dispersione per l’ultimo modello conferma i nostri sospetti. Non vediamo nessuna relazione evidente tra le previsioni dei movimenti percentuali settimanali dei prezzi e i movimenti effettivamente osservati, inoltre il MAE è praticamente lo stesso del primo modello.
price_pred = pd.concat([series.iloc[-int(len(series) * 0.3):].pct_change(),
predictions.iloc[-int(len(series) * 0.3):].pct_change()], axis=1)
price_pred.dropna(inplace=True)
price_pred.columns = ['Price', 'preds']
fig, ax = plt.subplots()
ax = sns.regplot(data=price_pred, x=price_pred['Price'], y=price_pred['preds'])
plt.xlabel('Observations')
plt.ylabel('Predictions')
plt.title('EURUSD Observed vs Predicted Values')
ax.grid(True, which='both')
ax.axhline(y=0, color='#888888')
ax.axvline(x=0, color='#888888')
sns.despine(ax=ax, offset=0)
plt.xlim(-0.05, 0.05)
plt.ylim(-0.05, 0.05)
plt.show()
mae = round(abs(test.pct_change() - predictions.pct_change()).mean(),10)
print(f'The MAE is {mae}')
Dal grafico sembra, anche in questo caso,- che non ci sia alcuna relazione tra le previsioni e i valori osservati quando si considera la variazione percentuale settimanale del prezzo EURUSD.
Se calcoliamo il “tasso di successo” di prevedere correttamente la direzione del movimento della prossima settimana per l’EURUSD, ancora una volta otteniamo un valore leggermente inferiore al 50%.
price_pred['hit'] = np.where(np.sign(price_pred['Price']) == np.sign(price_pred['preds']), 1, 0)
print(f"Hit rate: {round((price_pred['hit'].sum() / price_pred['hit'].count()) * 100,2)}%")
Hit rate: 48.16%
Analisi predittiva delle componenti delle serie
Cosa possiamo fare per migliorare l’analisi predittiva delle serie temporali? Come possiamo usare i metodi di decomposizione per migliorare le nostre capacità predittive? Come accennato in precedenza, ogni componente ha caratteristiche diverse quindi devono essere trattate individualmente in modo più dettagliato. Se possiamo suddividere una serie temporale nelle sue componenti e usare i modelli previsionali sulle singole componenti, allora per ottenere una previsione completa è sufficiente ricombinare le previsioni delle singole componenti. In altre parole possiamo ottenere una previsione più accurata se sommiamo insieme le singole previsioni.
Iniziamo scomponendo la serie dei prezzi dell’EURUSD tramite il filtro Hodrick-Prescott e memorizziamo ogni componente. Tracciamo quindi il grafico dell’andamento di ogni singola componenti.
cycle, trend = sm.tsa.filters.hpfilter(series, 50)
fig, ax = plt.subplots(3,1)
ax[0].plot(series)
ax[0].set_title('Price')
ax[1].plot(trend)
ax[1].set_title('Trend')
ax[2].plot(cycle)
ax[2].set_title('Cycle')
plt.show()
Dobbiamo quindi eseguire il nostro modello AR per le singole componenti della serie e tracciamo i grafici per confrontare le previsioni con i valori osservati
component_dict = {'cycle': cycle, 'trend': trend}
prediction_results = []
for component in ['trend', 'cycle']:
historic = component_dict[component].iloc[:int(len(series) * 0.7)].to_list()
test = component_dict[component].iloc[int(len(series) * 0.7):]
predictions = []
for i in range(len(test)):
model = AutoReg(historic, lags=10)
model_fit = model.fit()
pred = model_fit.predict(start=len(historic), end=len(historic), dynamic=False)
predictions.append(pred[0])
historic.append(test[i])
predictions = pd.Series(predictions, index=test.index, name=component)
prediction_results.append(predictions)
test_score = np.sqrt(mean_squared_error(test, predictions))
print(f'Test for {component} MSE: {test_score}')
# grafico risultati
plt.plot(test.iloc[:], label='Observed '+component)
plt.plot(predictions.iloc[:], color='red', label='Predicted '+component)
plt.legend()
plt.show()
Sembra che sia stata tracciata una sola linea per la componente di trend ma non è così perchè i valori previsti sono così vicini ai valori effettivi che una linea è nascosta dall’altra.
Successivamente “ricomponiamo” i nostri dati in modo additivo e calcoliamo l’RMSE confrontando le prestazioni delle previsione combinate con i valori osservati.
recomposed_preds = pd.concat(prediction_results,axis=1).sum(axis=1)
recomposed_preds.name = 'recomposed_preds'
plt.plot(series.iloc[int(len(series) * 0.7):], label='Observed')
plt.plot(recomposed_preds, color='red', label='Predicted')
plt.legend()
plt.show()
test_score = np.sqrt(mean_squared_error(series.iloc[int(len(series) * 0.7):], recomposed_preds))
print(f'RMSE: {test_score}')
Otteniamo un RMSE intorno a 0,0091, minore rispetto ai due tentativi precedenti che producevano valori più vicini a 0,013. Abbiamo una diminuzione del 33%. Quindi possiamo fare alcuni approfondimenti relativi a questo approccio di scomposizione/ricomposizione.
Tracciamo il grafico a dispersione e calcoliamo il MAE per questo modello.
price_pred = pd.concat([series.iloc[-int(len(series) * 0.3):].pct_change(),
recomposed_preds.iloc[-int(len(series) * 0.3):].pct_change()], axis=1)
price_pred.dropna(inplace=True)
price_pred.columns = ['Price', 'recomposed_preds']
fig, ax = plt.subplots()
ax = sns.regplot(data=price_pred, x=price_pred['Price'], y=price_pred['recomposed_preds'])
plt.xlabel('Observations')
plt.ylabel('Predictions')
plt.title('EURUSD Observed vs Predicted Values')
ax.grid(True, which='both')
ax.axhline(y=0, color='#888888')
ax.axvline(x=0, color='#888888')
sns.despine(ax=ax, offset=0)
plt.xlim(-0.05, 0.05)
plt.ylim(-0.05, 0.05)
plt.show()
mae = round(abs(series.iloc[-int(len(series) * 0.3):].pct_change() -
recomposed_preds.iloc[-int(len(series) * 0.3):].pct_change()).mean(),10)
print(f'The MAE is {mae}')
Anche il MAE è diminuito, da una cifra precedente di circa 0,013 a 0,0088 – anche in questo caso una differenza di circa il 33%. Non solo, ma la distribuzione dei punti dati sul grafico a dispersione ha assunto una forma leggermente diversa: la linea di regressione mostra una relazione più positiva tra le osservazioni e le nostre previsioni.
Come mostrato di seguito, ora siamo in grado di prevedere la direzione del movimento del prezzo per la settimana successiva il 60% delle volte, anziché un 50% casuale come negli ultimi due tentativi precedenti.
price_pred['hit'] = np.where(np.sign(price_pred['Price']) == np.sign(price_pred['recomposed_preds']), 1, 0)
print(f"Hit rate: {round((price_pred['hit'].sum() / price_pred['hit'].count()) * 100,2)}%")
Hit rate: 60.2%
Nonostante il RMSE, il MAE e l’Hit Rate siano migliorati in modo significativo con l’approccio decompose/recompose, senza effettuare l’ottimizzare dei parametri del modello, non è opportuno cantare vittoria troppo in fretta.
Abbiamo descritto l’applicazione del metodo di decomposizione del filtro HP per l’analisi predittiva delle serie temporali. Bisogna stare attenti quando si applicano nuovi metodi se non comprendiamo appieno come funziona la logica “dietro le quinte”, possiamo essere colti di sorpresa.
C’è un piccolo problema con il metodo del filtro HP. L’algoritmo di decomposizione fa uso di osservazioni che precedono e seguono la stima corrente. Questo approccio introduce nell’analisi un bias di “look-ahead”, che è uno dei bias più comuni che i trader introducono nei modelli di trading/predittivi. In questo caso, evitare il bias è facile se conosciamo esattamente come si calcola il filtro HP. E’ invece impossibile se non studiamo e lavoriamo per familiarizzare con il metodo di calcolo.
Sottolineiamo che il backtesting di modelli di trading algoritmico/sistematico è un’area che è ASSOLUTAMENTE sensibile e piena di vari tipi di bias, con il bias “look-ahead” in cima alla lista.
Quanto sopra descritto implica che il filtro HP produce serie delle componenti che contengono informazioni del “futuro”. Queste rende generalmente le previsioni più accurate e di conseguenza anche le previsioni ricomposte sono più accurate in modo fuorviante.
Vediamo cosa succede quando torniamo a usare un metodo di scomposizione che non “imbroglia” e non si avvale di “informazioni future” durante la scomposizione di una serie temporale. Consideriamo il metodo di STL Decomposition. Di seguito vediamo una scomposizione e una visualizzazione dei risultati.
result = STL(series).fit()
result.plot()
plt.show()
result tramite le paroli chiave trend, seasonal o resid. Ad esempio, stampiamo le prime 5 righe della componente della stagionalità.
print(result.seasonal.head())
Date
2000-01-09 0.060011
2000-01-16 0.047314
2000-01-23 0.041816
2000-01-30 0.012704
2000-02-06 0.024431
Freq: W-SUN, Name: season, dtype: float64
Ora dobbiamo scrivere il codice per eseguire il modello AR, come in precedenza.
component_dict = {'seasonal': result.seasonal, 'trend': result.trend, 'residual': result.resid}
prediction_results = []
for component in ['seasonal', 'trend', 'residual']:
historic = component_dict[component].iloc[:int(len(series) * 0.7)].to_list()
test = component_dict[component].iloc[int(len(series) * 0.7):]
predictions = []
for i in range(len(test)):
model = AutoReg(historic, lags=10)
model_fit = model.fit()
pred = model_fit.predict(start=len(historic), end=len(historic), dynamic=False)
predictions.append(pred[0])
historic.append(test[i])
predictions = pd.Series(predictions, index=test.index, name=component)
prediction_results.append(predictions)
test_score = np.sqrt(mean_squared_error(test, predictions))
print(f'Test for {component} MSE: {test_score}')
# Grafico risultati
plt.plot(test.iloc[:], label='Observed '+component)
plt.plot(predictions.iloc[:], color='red', label='Predicted '+component)
plt.legend()
plt.show()
Successivamente “ricomponiamo” i dati in modo additivo e calcoliamo l’RMSE confrontando le prestazioni delle previsione combinate risultanti con i valori osservati, come abbiamo fatto per il modello con il filtro HP.
recomposed_preds = pd.concat(prediction_results,axis=1).sum(axis=1)
plt.plot(series.iloc[int(len(series) * 0.7):], label='Observed')
plt.plot(recomposed_preds, color='red', label='Predicted')
plt.legend()
plt.show()
test_score = np.sqrt(mean_squared_error(series.iloc[int(len(series) * 0.7):], recomposed_preds))
print(f'RMSE: {test_score}')
Test for residual MSE: 0.009840348913376536
Il punteggio RMSE sembra essersi spostato nella direzione sbagliata (0,0098). E’ più vicino ai valori che abbiamo ottenuto nei due primi modelli, il modello “Persistenza” e il modello AR applicati ai non scomposti, cioè usando la serie dei prezzi grezzi.
Questo non è di buon auspicio. Vediamo ora il grafico a dispersione delle variazioni di prezzo percentuali settimanali (i valori osservati rispetto a quelli previsti) e calcoliamo il MAE.
price_pred = pd.concat([series.iloc[-int(len(series) * 0.3):].pct_change(), recomposed_preds.iloc[-int(len(series) * 0.3):].pct_change()], axis=1)
price_pred.dropna(inplace=True)
price_pred.columns = ['Price', 'recomposed_preds']
fig, ax = plt.subplots()
ax = sns.regplot(data=price_pred, x=price_pred['Price'], y=price_pred['recomposed_preds'])
plt.xlabel('Observations')
plt.ylabel('Predictions')
plt.title('EURUSD Observed vs Predicted Values')
ax.grid(True, which='both')
ax.axhline(y=0, color='#888888')
ax.axvline(x=0, color='#888888')
sns.despine(ax=ax, offset=0)
plt.xlim(-0.05, 0.05)
plt.ylim(-0.05, 0.05)
plt.show()
mae = round(abs(series.iloc[-int(len(series) * 0.3):].pct_change() -
recomposed_preds.iloc[-int(len(series) * 0.3):].pct_change()).mean(),10)
print(f'The MAE is {mae}')
The MAE is 0.0120908732
Anche il MAE è ritornato intorno ai valori visti in precedenza, e giusto per confermare ciò che ora ci dovrebbe essere chiaro, l’Hit Ratio è tornato al livello del 50%. Questo significa che “il modello non è migliore per prevedere la direzione della variazione del prezzo nella prossima settimana rispetto al lancio di una moneta”.
price_pred = pd.concat([series.iloc[-int(len(series) * 0.3):].pct_change(),
recomposed_preds.iloc[-int(len(series) * 0.3):].pct_change()], axis=1)
price_pred.dropna(inplace=True)
price_pred.columns = ['Price', 'preds']
price_pred['hit'] = np.where(np.sign(price_pred['Price']) == np.sign(price_pred['preds']), 1, 0)
print(f"Hit rate: {round((price_pred['hit'].sum() / price_pred['hit'].count()) * 100,2)}%")
Hit rate: 51.17%
Abbiamo visto che l’esito del metodo di decomposizione utilizzato influisce enormemente sulle previsioni finali. I parametri del modello AR non sono stati stimati usando le osservazioni future, né tantomeno per la decomposizione STL, quindi non è sbagliato che non siamo riusciti a trovare una soluzione alla possibilità di effettuare una previsione.
Sottolineiamo, come accennato in precedenza, che dobbiamo avere molta attenzione quando si usano determinati metodi e strumenti come la decomposizione. Anche se il modello AR usato per prevedere i valori delle serie scomposte emesse dal filtro HP non soffre di bias di previsione, il modo in cui le serie sono separate durante la decomposizione stessa usa i movimenti futuri sconosciuti al momento della stima e quindi contamina l’output predittivo del modello.
A mio parere, un euro risparmiato evitando trappole ed inutili perdite causate da un’analisi fuorviante o semplicemente errata vale tanto quanto un euro guadagnato a seguito di un’analisi riuscita. Quindi non è stata una perdita di tempo.
Codice completo
In questo articolo abbiamo descritto come effettuare l’analisi predittiva delle serie temporali in Python. Per il codice completo riportato in questo articolo, si può consultare il seguente repository di github:
https://github.com/datatrading-info/AnalisiDatiFinanziari
