Analisis de PCA en Python
La dimensionalidad en los datasets, suele ser un problema muy frecuente en los modelos que utilizan técnicas de machine learning o Deep learning. En este nuevo artículo, vamos a tratar de una forma lo mas simplificada posible, un método para reducir la dimensionalidad, en este caso utilizaremos el análisis de componentes principales (PCA) en Python, mediante la librería de sklearn.
Para poner un poco en contexto, vamos a dejar unos enlaces, con una explicacion matematica sobre lo que estamos haciendo, se puede consultar esta informacion aqui:
https://www.cienciadedatos.net/documentos/35_principal_component_analysis
https://towardsdatascience.com/pca-using-python-scikit-learn-e653f8989e60
Mediante el análisis de los componentes principales (PCA) vamos a transformar un dataset con demasiadas features, algunas relevantes, pero otras altamente correlacionadas, en un nuevo dataset, donde aplicamos PCA y reducimos su dimensionalidad, hasta un punto, que podamos afirmar, que no perdemos demasiada información.

Los modelos de Deep learning aplicados para trading,que podemos programar en Python, de una forma realmente sencilla funcionan muchísimo mejor, con datasets con pocas features, y muy relevantes, además de modelos de redes neuronales, con pocas capas y pocas neuronas por capa, además de nuestros preciados dropouts, que nos eliminan información de forma aleatoria, añadiendo dificultad adicional al modelo, y evitando la sobre optimización o el overfit de dicho modelo.
Este articulo será utilizado de base, para los futuros artículos, donde implementaremos técnicas de machine learning y Deep learning en Python para poder predecir la dirección del precio.
ATENCION! El proceso de formación de los precios es un proceso estocastico, un movimiento browniano, en el cual la aleatoriedad está MUY PRESENTE, y que incluso la propia predicción del precio altera el propio precio, al disponer de predicciones, es información que añadimos a nuestras decisiones o modelos, y, en consecuencia, se opera en función a ello, alterando la formación del precio.
En mi opinión es erróneo intentar predecir el precio, pues tal que retornos pasados no garantizan retornos futuros, precios pasados, no influyen de forma determinante en la formación de precios futuros.
Programando un analisis PCA en Python

Para empezar, lo primero que vamos a hacer, es importar todas las librerías necesarias, como viene siendo clásico en Python…
import pandas as pd
import numpy as np
import datetime
import yfinance as yf
import quandl
import matplotlib.pyplot as plt
import seaborn as sns
import ta
import warnings
from sklearn.preprocessing import MinMaxScaler, StandardScaler
from sklearn import metrics
from sklearn.decomposition import PCA
from sklearn.utils import class_weight
from sklearn.utils.class_weight import compute_class_weight
from sklearn.metrics import classification_report
import matplotlib.dates as mdates
import quantstats as qt
import warnings
%matplotlib inline
warnings.filterwarnings(action='ignore')
plt.style.use('ggplot')
plt.rcParams["figure.figsize"] = [20,9]
Una vez tenemos todas las librerías descargadas, vamos a crear nuestro dataset, en este caso he creado una combinación de features con diferentes datos, procesándolos en distintos datasets, para su posterior manipulación. Se han utilizado los futuros continuos del Nasdaq, el S&P500 el VIX y los bonos a 10 y 30 años americanos.
Una vez tenemos todas las librerías descargadas, vamos a crear nuestro dataset, en este caso he creado una combinación de features con diferentes datos, procesándolos en distintos datasets, para su posterior manipulación. Se han utilizado los futuros continuos del Nasdaq, el S&P500 el VIX y los bonos a 10 y 30 años americanos. Además, añadimos unos datos alternativos, tales como los nuevos máximos y mínimos, los COT o el índice VXST, que es un VIX a 9 días.
#Datos RAW de mercado
#indces
nq = pd.read_csv('NQ___CCS.csv',parse_dates=True,index_col=['Date'])
es = pd.read_csv('ES___CCS.csv',parse_dates=True,index_col=['Date'])
vx = pd.read_csv('VX___CCS.csv',parse_dates=True,index_col=['Date'])
#bonos
y10 = pd.read_csv('TY2__CCS.csv',parse_dates=True,index_col=['Date'])
y30 = pd.read_csv('US2__CCS.csv',parse_dates=True,index_col=['Date'])
#alternativos
newhighs = quandl.get("URC/NASDAQ_52W_HI", authtoken="xPLyKnHL4wCYWyyziuSs")
newlows = quandl.get("URC/NASDAQ_52W_LO", authtoken="xPLyKnHL4wCYWyyziuSs")
nqhighs =quandl.get("URC/NASDAQ_ADV", authtoken="xPLyKnHL4wCYWyyziuSs")
nqlows =quandl.get("URC/NASDAQ_DEC", authtoken="xPLyKnHL4wCYWyyziuSs")
cot_30 = quandl.get("CFTC/020601_FO_ALL", authtoken="xPLyKnHL4wCYWyyziuSs")
cot_30 = cot_30[['Dealer Longs','Dealer Shorts','Asset Manager Longs','Asset Manager Shorts','Leveraged Funds Longs','Leveraged Funds Shorts',
'Total Reportable Longs','Total Reportable Shorts','Non Reportable Longs','Non Reportable Shorts']]
cot_vx = quandl.get("CFTC/1170E1_FO_ALL", authtoken="xPLyKnHL4wCYWyyziuSs")
cot_vx = cot_vx[['Dealer Longs','Dealer Shorts','Asset Manager Longs','Asset Manager Shorts','Leveraged Funds Longs','Leveraged Funds Shorts',
'Total Reportable Longs','Total Reportable Shorts','Non Reportable Longs','Non Reportable Shorts']]
#vx9d
vx9 = pd.read_csv('https://www.cboe.com/publish/scheduledtask/mktdata/datahouse/vix9ddailyprices.csv',skiprows=3,parse_dates=True)
vx9=vx9[['Date','Close']]
vx9['Date']=vx9['Date'].str.strip('*')
vx9['Date'] = pd.to_datetime(vx9['Date'], infer_datetime_format=True)
vx9 = vx9.set_index('Date')
Una vez que ya tenemos toda la información raw (cruda, sin preprocesar) descargada, vamos a hacer lo propio con dichos datos, vamos a obtener features desde ellos, al final estamos añadiendo la misma información, desde diferentes puntos de vista. Creamos un nuevo dataframe llamado features, conde contendremos toda esta información.
features = pd.DataFrame()
#Features ES
features['es'] = es['Close']
features['es_vol'] = es['Volume']
features['es_oi'] = es['Open Interest']
features['es_delta_vol'] = features['es_vol'] - features['es_vol'].shift(1)
features['es_delta_vol2'] = features['es_vol'] - features['es_vol'].shift(10)
features['es_delta_vol3'] = features['es_vol'] - features['es_vol'].shift(25)
features['es_range_oh'] = es['Open'] - es['High'].shift(1)
features['es_range_oc'] = es['Open'] - es['Close'].shift(1)
features['es_range_ol'] = es['Open'] - es['Low'].shift(1)
features['es_range_oo'] = es['Open'] - es['Open'].shift(1)
features['es_vol_dolar'] = es['Volume'] * es['Close'] * 50
features['es_dist_max'] = es['Close'].rolling(251).max() - es['Close']
features['es_dist_min'] = es['Close'].rolling(251).min() - es['Close']
features['es_autocorr'] = es['Close'].rolling(20).corr() - es['Close'].shift(1).rolling(20).corr()
features['es_autocorr2'] = es['Close'].rolling(20).corr() - es['Close'].shift(5).rolling(20).corr()
#Features US
features['us'] = y30['Close']
features['us_vol'] = y30['Volume']
features['us_vol2'] = y10['Volume']
features['us_oi'] = y30['Open Interest']
features['us_delta_vol'] = features['us_vol'] - features['us_vol'].shift(1)
features['us_delta_vol2'] = features['us_vol'] - features['us_vol'].shift(10)
features['us_delta_vol3'] = features['us_vol'] - features['us_vol'].shift(25)
features['us_range_oh'] = y30['Open'] - y30['High'].shift(1)
features['us_range_oc'] = y30['Open'] - y30['Close'].shift(1)
features['us_range_ol'] = y30['Open'] - y30['Low'].shift(1)
features['us_range_oo'] = y30['Open'] - y30['Open'].shift(1)
features['us_vol_dolar'] = y30['Volume'] * y30['Close'] * 31.25
features['us_dist_max'] = y30['Close'].rolling(251).max() - y30['Close']
features['us_dist_min'] = y30['Close'].rolling(251).min() - y30['Close']
features['us_autocorr'] = y30['Close'].rolling(20).corr() - y30['Close'].shift(1).rolling(20).corr()
features['us_autocorr2'] = y30['Close'].rolling(20).corr() - y30['Close'].shift(5).rolling(20).corr()
#Features VIX
features['vx'] = vx['Close']
features['vx_vol'] = vx['Volume']
features['vx_oi'] = vx['Open Interest']
features['vx_delta_vol'] = features['vx_vol'] - features['vx_vol'].shift(1)
features['vx_delta_vol2'] = features['vx_vol'] - features['vx_vol'].shift(10)
features['vx_delta_vol3'] = features['vx_vol'] - features['vx_vol'].shift(25)
features['vx_range_oh'] = vx['Open'] - vx['High'].shift(1)
features['vx_range_oc'] = vx['Open'] - vx['Close'].shift(1)
features['vx_range_ol'] = vx['Open'] - vx['Low'].shift(1)
features['vx_range_oo'] = vx['Open'] - vx['Open'].shift(1)
features['vx_vol_dolar'] = vx['Volume'] * vx['Close'] * 1000
features['vx_dist_max'] = vx['Close'].rolling(251).max() - vx['Close']
features['vx_dist_min'] = vx['Close'].rolling(251).min() - vx['Close']
features['vx_autocorr'] = vx['Close'].rolling(20).corr() - vx['Close'].shift(1).rolling(20).corr()
features['vx_autocorr2'] = vx['Close'].rolling(20).corr() - vx['Close'].shift(5).rolling(20).corr()
#Alternative Features
features['yield_curve_30_10'] = (y30['Close'].pct_change() - y10['Close'].pct_change())
features['y30v'] = np.where(features['us_vol'] > features['us_vol2'],1,0)
features['inter_corr1'] = round(features['es'].rolling(20).corr(features['us']) - features['es'].shift(2).rolling(20).corr(features['us']),4)
features['nyse_amp'] = newhighs['Numbers of Stocks'] - newlows['Numbers of Stocks']
features['nq_amp'] = nqhighs['Numbers of Stocks'] - nqlows['Numbers of Stocks']
#Alternative CBOE
features['vxst'] = vx9['Close']
features['vxst_vx'] = features['vxst'].rolling(9).corr(features['vx'])
features['vxst_vx2'] = features['vxst'].rolling(20).corr(features['vx'])
features['delta_vxst_vx'] = features['vxst_vx'] - features['vxst_vx2']
features = features.dropna()
#Alternative COT
features['us_dl'] = cot_30['Dealer Longs']
features['us_ds'] = cot_30['Dealer Shorts']
features['us_aml'] = cot_30['Asset Manager Longs']
features['us_ams'] = cot_30['Asset Manager Shorts']
features['us_lfl'] = cot_30['Leveraged Funds Longs']
features['us_lfl'] = cot_30['Leveraged Funds Shorts']
features['us_longs'] = cot_30['Total Reportable Longs'] + cot_30['Non Reportable Longs']
features['us_shorts'] = cot_30['Total Reportable Shorts'] + cot_30['Non Reportable Shorts']
features['vx_dl'] = cot_vx['Dealer Longs']
features['vx_ds'] = cot_vx['Dealer Shorts']
features['vx_aml'] = cot_vx['Asset Manager Longs']
features['vx_ams'] = cot_vx['Asset Manager Shorts']
features['vx_lfl'] = cot_vx['Leveraged Funds Longs']
features['vx_lfl'] = cot_vx['Leveraged Funds Shorts']
features['vx_longs'] = cot_vx['Total Reportable Longs'] + cot_vx['Non Reportable Longs']
features['vx_shorts'] = cot_vx['Total Reportable Shorts'] + cot_vx['Non Reportable Shorts']
features.fillna(method='ffill', inplace=True)
Al final hemos creado un dataset con 69 columnas, algo totalmente inviable, repito lo mencionado anteriormente, “Keep it simple”, Pocas features, y muy relevantes, pocas capas y pocas neuronas.
Al tener cada feature(columna) unos valores máximos y mínimos distintos, debemos normalizar los datos, es una condición fundamental, para ello utilizaremos, la librería sklearn, muy conocida en Python y su método MinMaxScaler()
features.shape
X = features
X2 = X
X = MinMaxScaler().fit_transform(X)
X = pd.DataFrame(X,columns=X2.columns).set_index(X2.index)
X = X.dropna()
X
Una vez que tenemos los datos homogeneizados, el dataset esta listo para enfrentarse a un proceso de PCA que le reducirá la dimensionalidad. Para ello utilizaremos una vez más, la librería sklearn.
Empezaremos, intentando reducir el dataset con 69 columnas, hasta 3 columnas.
Hemos conseguido un 55% de explicar el dataset con tan solo 3 columnas, es un gran dato, pero no es suficiente, En mi opinión, un proceso de reducción de dimensionalidad mediante PCA necesitamos un 85% de explicación, en el menor número de columnas posibles.
pca = PCA(n_components=3)
principalComponents = pca.fit_transform(X)
principalDf = pd.DataFrame(data = principalComponents
, columns = ['pca1', 'pca2','pca3']).set_index(X.index)
print("shape of X_pca", principalDf.shape)
expl = pca.explained_variance_ratio_
print(expl)
print('suma:',sum(expl[0:16]))
principalDf
pca = PCA(n_components=4)
principalComponents = pca.fit_transform(X)
principalDf = pd.DataFrame(data = principalComponents
, columns = ['pca1', 'pca2','pca3','pca4']).set_index(X.index)
print("shape of X_pca", principalDf.shape)
expl = pca.explained_variance_ratio_
print(expl)
print('suma:',sum(expl[0:16]))
pca = PCA(n_components=5)
principalComponents = pca.fit_transform(X)
principalDf = pd.DataFrame(data = principalComponents
, columns = ['pca1', 'pca2','pca3','pca4','pca5']).set_index(X.index)
print("shape of X_pca", principalDf.shape)
expl = pca.explained_variance_ratio_
print(expl)
print('suma:',sum(expl[0:16]))
pca = PCA(n_components=8)
principalComponents = pca.fit_transform(X)
principalDf = pd.DataFrame(data = principalComponents
, columns = ['pca1', 'pca2','pca3','pca4','pca5','pca6','pca7','pca8']).set_index(X.index)
print("shape of X_pca", principalDf.shape)
expl = pca.explained_variance_ratio_
print(expl)
print('suma:',sum(expl[0:16]))
pca = PCA(n_components=10)
principalComponents = pca.fit_transform(X)
principalDf = pd.DataFrame(data = principalComponents
, columns = ['pca1', 'pca2','pca3','pca4','pca5','pca6','pca7','pca8','pca9','pca10']).set_index(X.index)
print("shape of X_pca", principalDf.shape)
expl = pca.explained_variance_ratio_
print(expl)
print('suma:',sum(expl[0:16]))
Vemos como aumentando el numero de outputs que nos devuelve el proceso de PCA, va aumentando la explicación del dataset, al final debemos encontrar el equilibro entre la mayor explicación posible (el numero mas cercano a 100) con el menor numero de columnas (features) posibles.
Conclusiones
Y para concluir, recordamos, que la mejor combinacion es, la que tiene menos columnas, y aporta un mayor grado de explicacion del dataset

Como podemos apreciar en la imagen,a partir de las 6 dimensiones de nuestro PCA, deja de aportar en exceso, y es mas, añades dimensionalidad al dataset y no estas ganando demasiado poder explicativo, por consecuencia, es contraproducente.
Si te ha gustado este articulo,te recompendamos tambien puedes revisar los articulos anterioress seguir aprendiendo! Portfolios1 y Portfolios2