Introducción

A lo largo de los siguientes artículos vamos a ver una introducción sobre cómo analizar series temporales utilizando Python.

Antes de empezar a diseñar una estrategia de trading resulta necesario analizar detenidamente el activo sobre el que vamos a trabajar. Algunas preguntas que podemos hacernos son; ¿Cuál es su volatilidad? ¿Qué nivel mínimo de volatilidad se alcanzó históricamente?¿Cuánto llegó a retroceder el precio en la peor sesión? O la pregunta por excelencia; ¿Tiende a tener un comportamiento tendencial o anti-tendencial? Al analizar su comportamiento histórico respondemos a estas preguntas, y esto nos permitirá enfocar el diseño de nuestra estrategia de una manera mucho más eficiente.

Aquí vamos a crear una herramienta en Python muy sencilla, pero también muy práctica, que nos va a permitir caracterizar el activo que queremos analizar. El propósito final será obtener información sobre si debemos abordar el diseño de la estrategia para ese activo con un enfoque tendencial o anti-tendencial. También podemos encontrarnos con que no hay un sesgo claramente definido, lo cual supondría que no hay una ventaja evidente para diseñar un sistema basándonos directamente en uno de estos dos enfoques.

He dividido el análisis completo de una serie temporal en las siguientes fases. Estas serán realizadas siguiendo el orden secuencial que se muestran a continuación.

  • Estadística Descriptiva: En este apartado analizaremos visualmente la distribución de retornos, y obtendremos un resumen estadístico de la serie.
  • Análisis de volatilidad: Aquí analizaremos la volatilidad anualizada de la serie y su evolución histórica. En esta sección también compararemos la amplitud del rango en los días alcistas y bajistas para determinar si existe asimetría.
  • Exponente de Hurst: Aplicaremos esta herramienta para medir la persistencia de la serie temporal, lo cual nos dará información sobre el carácter tendencial, anti-tendencial o neutral de la serie.
  • Análisis de autocorrelación: A través de este análisis mediremos la relación entre el último precio de cierre y los precios de cierre anteriores.  Esto nos permitirá saber qué impacto tiene tienen los cierres anteriores sobre los cierres futuros. Este análisis nos sirve también para conocer si la serie tiene un comportamiento anti-tendencial, tendencial o neutral.
  • Análisis de la estrategia HCLC: Programaremos una estrategia basada en un canal de máximos y mínimos, y haremos un análisis paramétrico para observar su comportamiento histórico. Comprobaremos así si hay patrones explotables en el corto, medio y/o largo plazo utilizando un enfoque tendencial y/o anti-tendencial. Finalmente registraremos las curvas de capital de la mejor estrategia tendencial y de la mejor estrategia anti-tendencial.
  • Resumen de resultados: Recopilaremos los resultados obtenidos en el análisis de la estrategia HCLC.
  • Simulación de MonteCarlo: Para finalizar haremos una simulación de Montecarlo del sistema que ha mostrado el mejor comportamiento (tendencial o anti-tendencial).

Herramientas utilizadas

Para diseñar lo que anteriormente hemos llamado nuestra “herramienta” vamos a utilizar un notebook de Jupyter Notebook. Tras completar este cuaderno podremos utilizarlo para analizar cualquier activo que deseemos. Podemos acceder al entorno interactivo de Jupyter Notebook  a través de la instalación de Anaconda. En el sitio oficial hay dos versiones de Python disponibles para descarga; 2.6 y 3.7. La versión 2.6 está obsoleta por lo que no es conveniente instalarla. Aquí vamos a utilizar la versión actual, 3.7.

Vamos a acceder también la API de Yahoo Finance para descargar datos históricos de fin de día. Debemos tener en cuenta algo importante; los datos de Yahoo no son de buena calidad y suele ser habitual encontrar errores en ellos. Podemos encontrar picos irreales en las cotizaciones, y/o fallos en el cálculo de los splits o dividendos. Sin embargo, aquí hemos utilizado este proveedor ya que también tiene grandes ventajas; entre ellas el acceso a la API es gratuito sin necesidad de registro, y también dispone de un amplio número de activos con largos históricos. Por estas razones, como punto de partida Yahoo me parece útil para hacer un análisis exploratorio rápido y de manera cómoda. Si necesitamos hacer un análisis exhaustivo podremos adaptar fácilmente este NoteBook, bien para descargar datos de otro proveedor o bien para leer datos históricos de archivos propios con extensión CSV.

Estadística Descriptiva

Vamos a comenzar con la creación del Notebook. Lo primero que haremos será añadir las librerías de Python que vamos a emplear para realizar el análisis:

import pandas as pd
import numpy as np
import pandas_datareader as dr
import matplotlib.pyplot as plt
# Librería para mejorar la visualización del gráfico
plt.style.use('seaborn-darkgrid')

Vamos a explicar para qué sirve cada una de esas librerías.

> Numpy: Acrónico de Numerical Python, es la librería utilizada por excelencia para computación científica y análisis numérico. Incluye un amplio número de funcionalidades para realizar operaciones vectoriales y matriciales.

> Pandas: Construida sobre Numpy, es una librería que sirve para manipular y analizar datos en Python. Ofrece estructuras de datos y operaciones para manipular tablas numéricas y series temporales.

> Pandas_datareader: Es una sub-librería que permite obtener datos históricos de diferentes proveedores on-line y almacenarnos en un DataFrame. Desde aquí se pueden ver los proveedores actualmente soportados. Esta librería debe ser instalada manualmente ya que no viene por defecto en Anaconda. Si lo necesitas aquí encontrarás un video-tutorial para instalarla.

> Matplotlib: Es una librería que sirve para generar gráficos a partir de listas o arrays en Python y su extensión matemática Numpy.

A continuación vamos a configurar los parámetros necesarios para descargar los datos históricos del activo que queremos analizar. Aquí necesitamos configurar el nombre del activo, y el rango de fechas que vamos a descargar. Vamos a seleccionar el ETF SPDR S&P 500 (SPY) desde Febrero de 1993 hasta la fecha actual. Como se puede observar a continuación, a través de la librería Pandas_Datareader descargaremos estos datos y los almacenaremos en un nuevo DataFrame llamado “df”:

# Seleccionamos el activo a analizar y las fechas de estudio
from datetime import date
Activo = 'SPY'
FechaInicio = '1993-02-01'
FechaFinal = date.today()
# Descargamos la serie y mostramos el Tail.
df = dr.data.get_data_yahoo(Activo,start=FechaInicio,end=FechaFinal)
df = df[~df.index.duplicated()] # Eliminamos los duplicados en el índice (Problemas con criptomonedas)
df.tail()

La función df.tail() nos permitirá observar las últimas 5 filas de este DataFrame:

Vamos a crear ahora una nueva columna en nuestro DataFrame que va a mostrar el retorno aritmético porcentual del precio de cierre ajustado (Adj Close). Para ello, simplemente utilizaremos la función pct_change() de Pandas tal y como se muestra a continuación:

# Creamos una nueva columna que recoge el cambio porecentual del precio de Cierre.
df["R. Aritmetico"] = df["Adj Close"].pct_change() 
df.tail()

¡Y listo! Nuestro DataFrame quedará de este modo:

Ahora ya tenemos algo que empezar a analizar. Vamos a comenzar obteniendo un histograma de frecuencias para conocer la distribución histórica que han tenido estos retornos. En este histograma también vamos a añadir la función de densidad de probabilidad (FDP) de una distribución normal. Esto nos permitirá saber si la distribución histórica de los retornos se ajusta, o no, a una distribución normal. Construir este gráfico a través de Python es muy sencillo. Vamos a crearlo añadiendo el bloque de código que se muestra a continuación:

import seaborn as sns
from scipy import stats
from scipy.stats import norm

# Eliminamos la primera fila del DF; Contiene un valor NaN en el cálculo del R. Aritmético.
df = df.iloc[1:] 

# Dibujamos el historgrama de frecuencias. 
plt.figure(figsize=(15,8))
sns.set(color_codes = True)
ax = sns.distplot(df['R. Aritmetico'], bins=100, kde=False, fit=stats.norm, color='green')

# Obtenemos los parámetros ajustados de la distribución normal utilizados por SNS
(mu, sigma) = stats.norm.fit(df['R. Aritmetico'])

# Configuramos el título del gráfico, legendas y etiquetas.
plt.title('Distribución Histórica de Retornos Diarios',fontsize = 16)
plt.ylabel('Frecuencia')
plt.legend(["Distribución normal. fit ($\mu=${0:.2g}, $\sigma=${1:.2f})".format(mu, sigma), 
            "Distribución R. Aritmeticos"])

Como se puede observar aquí, en primer lugar hemos importado Seaborn. Seaborn es una librería de visualización de datos basada en Matplotlib que proporciona una interfaz de alto nivel para dibujar gráficos estadísticos atractivos e informativos. En las dos líneas de código siguientes importamos el módulo scipy.stats que contiene un amplio número de funciones estadísticas, y seguidamente importamos norm. Esta última función nos va a permitir crear la función de densidad de probabilidad (FDP) para la distribución normal.  En las siguientes líneas de código tenemos las instrucciones necesarias para crear nuestro histograma de frecuencias (ver comentarios dentro del código). Si ejecutamos dicho bloque de código obtendremos un gráfico como el siguiente:

El histograma consta de un eje-x (horizontal) y de un eje-y (vertical).  En nuestro caso, las barras verdes del histograma muestran la frecuencia (eje-y) con la que los retornos aritméticos se sitúan dentro de un determinado rango de valores del eje-x. La línea negra muestra la FDP de la distribución normal teórica. Dicha distribución se computa utilizando el promedio y la desviación típica de nuestra serie de retornos histórica.

Aquí podemos ver que nuestra distribución de retornos tiene cierta similitud con una distribución normal, sin embargo, esta no llega a ajustarse perfectamente. En realidad, las series de retornos que nos encontramos en los mercados rara vez, si no nunca, se ajustan perfectamente a una distribución normal. En general tienden a mostrar valores extremos que se desvían de su media con una probabilidad más alta de lo que se espera en una distribución normal. Lo cual produce distribuciones con colas largas y por ende mayor probabilidad de sufrir riesgo de cola (Tail Risk). El coeficiente de asimetría y la curtosis nos sirven para saber cuánto se desvía nuestra distribución de una distribución normal teórica. Obtendremos estos dos valores más adelante.

En el análisis de las series temporales financieras, por simplificación, se parte de la asunción de que las series de precios se ajustan a una distribución normal.  El principal problema que surge con esta teoría es que subestima el riesgo de obtener valores extremos. Es decir, si asumimos que la distribución es normal y en realidad sus colas son más largas de lo normal,  las probabilidades de obtener valores extremos serán más altas de lo esperado. Por ello, cuando esto ocurre resulta necesario descartar la normalidad, e utilizar herramientas más eficientes para medir el riesgo. Entre las cuales tenemos el VaR o el CVaR aplicando el método de cálculo histórico.

A continuación vamos a añadir el siguiente bloque de código para obtener la estadística descriptiva de los retornos. Aquí hemos añadido los estadísticos que considero más importantes, si queremos podremos personalizar este apartado añadiendo aquellos otros que consideremos. Dentro del código encontraremos comentarios que nos indican para qué sirven las distintas instrucciones:

# Calculamos la Tasa de Crecimiento Anual Compuesto y el resultado de Comprar y mantener.
Años  = df["R. Aritmetico"].count() / 252
CAGR = (df['Adj Close'].iloc[-1]/ df['Adj Close'].iloc[0]) ** (1 / Años) - 1
print('> Tasa de Crecimiento Anual Compuesto:','%.6s' % (100 * CAGR) ,'%')
print('> Buy & Hold:', '%.6s' % (100*((df['Adj Close'].iloc[-1] - df['Adj Close'].iloc[0])/ df['Adj Close'].iloc[0])),'%')
# Calculamos el Máximo Drawdown.
Maximo_Anterior = df["Adj Close"].cummax()
drawdowns = 100*((df["Adj Close"] - Maximo_Anterior )/Maximo_Anterior )
DD = pd.DataFrame({"Adj Close": df["Adj Close"], 
                   "Previous Peak": Maximo_Anterior , 
                   "Drawdown": drawdowns})
print ('> Máximo Drawdown Histórico:', '%.6s' % np.min(DD['Drawdown']), '%') 

# Obtenemos el promedio, desviación típica, máximo y mínimo valor y número de datos analizados:
print ('> Media Diaria:', '%.6s' % (100 * df["R. Aritmetico"].mean()),'%')
print ('> Desviación Típica Diaria:', '%.6s' % (100 * df["R. Aritmetico"].std(ddof=1)),'%')
print ('> Máxima Pérdida Diaria:', '%.6s' % (100 * df["R. Aritmetico"].min()),'%')
print ('> Máximo Beneficio Diario:', '%.6s' % (100 * df["R. Aritmetico"].max()),'%')
print ('> Dias Analizados:', '%.6s' % df["R. Aritmetico"].count())
print('<-------------------------------------------------->')

# Coeficiente de asimetria y curtosis de la distribución.
print ('> Coeficiente de Asimetría:', '%.6s' % df["R. Aritmetico"].skew())
print ('> Curtosis:', '%.6s' % df["R. Aritmetico"].kurt())
print('<-------------------------------------------------->')

# VaR Teórico obtenido a través de la distribución normal al 95% y 99% de confianza.
print('> VaR Modelo Gaussiano NC-95% :' , '%.6s' % (100 * norm.ppf(0.05, mu, sigma)),'%')
print('> VaR Modelo Gaussiano NC-99% :' , '%.6s' % (100 * norm.ppf(0.01, mu, sigma)),'%')
print('> VaR Modelo Gaussiano NC-99.7% :' , '%.6s' % (100 * norm.ppf(0.003, mu, sigma)),'%')

# VaR histórico al 95% y 99% de confianza.
print('> VaR Modelo Historico NC-95% :', '%.6s' % (100 * np.percentile(df["R. Aritmetico"],5)),'%')
print('> VaR Modelo Historico NC-99% :', '%.6s' % (100 * np.percentile(df["R. Aritmetico"],1)),'%')
print('> VaR Modelo Historico NC-99.7% :', '%.6s' % (100 * np.percentile(df["R. Aritmetico"],.3)),'%')

Si ejecutamos este bloque de código obtendremos los siguientes resultados:

En el primer apartado tenemos los estadísticos básicos de la serie; tasa de crecimiento anual compuesta, resultado de la “estrategia” Buy & Hold, máximo Drawdown histórico, etc.

En el apartado inferior tenemos el coeficiente de asimetría y el coeficiente de curtosis. Ambos sirven para conocer cómo se distribuyen los retornos en nuestra distribución, y cuanto se desvían de una distribución normal.

El coeficiente de asimetría (skewness) mide la asimetría de la distribución de los retornos sobre su media. Este valor puede ser positivo, negativo o neutro. En una distribución normal el coeficiente de asimetría es 0 ya que ambas colas están simétricamente balanceadas. En nuestra distribución este valor es positivo lo cual indica que la cola derecha es más larga que la izquierda. Dicho de otro modo, en la cola derecha encontramos valores más alejados de la media que en la cola izquierda. Lo cual en este caso particular es positivo.

La curtosis (Kurtosis) define en qué  grado las colas de la distribución difieren de las colas de una distribución normal. En una distribución normal la curtosis es 3. En nuestra distribución este valor es superior, lo cual nos indica que las colas son más largas de lo que se espera en una distribución normal. Dicho de otro modo, nos indica que nuestra distribución contiene valores que exceden las 3 desviaciones típicas de la media, las cuales no son excedidas con una probabilidad del 99.7% en una distribución normal. Cuando la curtosis excede la normalidad se dice que la distribución es leptocúrtica.

Tras analizar ambos coeficientes vemos que la serie no se ajusta totalmente a la normalidad. Vamos a ver ahora cuál es la limitación de estimar el riesgo de la serie asumiendo que sigue una distribución normal.

En el último apartado de las estadísticas tenemos el VaR (Value at Risk) calculado de dos modos distintos: En el primer caso se utiliza el modelo gaussiano, el cual asume que los retornos están normalmente distribuidos. En el segundo caso se utiliza el modelo basado en datos históricos, en el que se utiliza la distribución empírica, y en el que no se hacen asunciones sobre el tipo de distribución de la serie. En ambos casos se obtiene cuál es la peor pérdida esperada con un nivel de confianza del 95%, 99% y 99.7%.

Aquí se observa que con un nivel de confianza del 95% ambos métodos obtienen un valor muy similar, sin embargo, al 99%  y al 99.7% las diferencias son notables. En ambos casos el resultado esperado en el modelo empírico es peor que en el modelo gaussiano, lo cual nos confirma nuevamente que la cola izquierda de nuestra distribución es más larga que la normal.

Análisis de volatilidad

Vamos a analizar ahora la volatilidad histórica de nuestra serie temporal. Principalmente nos interesa conocer la evolución de la volatilidad durante todo el periodo. El código inferior nos sirve para crear un gráfico que muestra la evolución de los siguientes tres valores:

  • Volatilidad anualizada: Desviación típica de los últimos 14 días multiplicada por la raíz cuadrada de 252.
  • Media aritmética sobre la volatilidad anualizada (SMA): Muestra el promedio de los últimos 126 valores de la volatilidad anualizada.
  • Precio de cierre de ajustado.

Dentro del siguiente bloque de código encontraremos diferentes comentarios que explican la utilidad de las distintas instrucciones:

# Registramos Matplotlib converters para adatpar la fecha -DateTime- al método de dibujo de Matplotlib
from pandas.plotting import register_matplotlib_converters
register_matplotlib_converters()

# Calculamos la volatilidad histórica de 14 días, la volatilidad anualizada y su SMA de 252 días.
df['Volatilidad_Historica_14_Dias'] = 100*df["R. Aritmetico"].rolling(14).std()
df['Volatilidad_14_Dias_Anualizada'] = df['Volatilidad_Historica_14_Dias']*(252**0.5)
df['SMA_126_Volatilidad_Anualizada'] =  df ['Volatilidad_14_Dias_Anualizada'].rolling(126).mean()

# Creamos un gráfico en el que mostramos el precio de cierre, la volatilidad anualizada y su SMA.
fig, ax1 = plt.subplots(figsize=(15, 8))
ax2 = ax1.twinx()
ax1.plot(df['Volatilidad_14_Dias_Anualizada'], 'orange', linestyle='--')
ax1.plot(df['SMA_126_Volatilidad_Anualizada'], 'Green', linestyle='-')
ax2.plot(df['Adj Close'], 'black')

# Asingamos un nombre para el gráfico y para los ejes x e y.
plt.title('Evolucion Histórica del Precio y la Volatilidad',fontsize = 16)
ax1.set_xlabel('Fecha')
ax1.set_ylabel('Volatilidad Anualizada', color='black')
ax2.set_ylabel('Precio de Cierre', color='black')


# Configuramos la leyenda, las líneas de cuadrícula y mostramos el gráfico. 
ax1.legend(loc='upper left', frameon=True, borderpad=1)
ax1.grid(True)
ax2.grid(False)
plt.show()

Este código produce el siguiente gráfico:

Sin duda una imagen vale más que mil palabras. Aquí podemos ver la evolución histórica de la volatilidad anualizada,  mostrando los distintos picos y valles producidos a lo largo del histórico, y su relación con los movimientos históricos del precio.

El siguiente bloque de instrucciones nos permitirá analizar esta volatilidad más a fondo. Aquí comenzaremos obtenido la volatilidad histórica. Seguidamente localizaremos cuál ha sido el valor máximo y mínimo alcanzado históricamente, y registraremos las fechas en las que se han producido ambos valores. Para completar el análisis también nos interesa saber si existe alguna relación entre la volatilidad y la dirección de la sesión. Para responder a esta pregunta calcularemos el rango medio porcentual de los días alcistas y bajistas de un modo independientemente. Y seguidamente construiremos un ratio que comparará ambos valores para conocer la asimetría de la volatilidad.

# Volatilidad anualizada.
VAM = (252**0.5)*(100* df["R. Aritmetico"].std())
print('> Volatilidad Anualizada:','%.6s' % VAM ,'%')

# Obtenemos el valor mínimo y máximo de la volatilidad anualizada así como las fechas
# en las que ambos valores son producidos.
Fecha_Minima_Volatilidad = df.Volatilidad_14_Dias_Anualizada[df.Volatilidad_14_Dias_Anualizada == 
                          df['Volatilidad_14_Dias_Anualizada'].min()].index.strftime('%Y-%m-%d').tolist()
Fecha_Maxima_Volatilidad = df.Volatilidad_14_Dias_Anualizada[df.Volatilidad_14_Dias_Anualizada ==  
                          df['Volatilidad_14_Dias_Anualizada'].max()].index.strftime('%Y-%m-%d').tolist()

print ('> La Mínima Volatilidad Anualizada fue de', '%.6s' % 
       (df['Volatilidad_14_Dias_Anualizada'].min()),
       '%', 'registrada el', Fecha_Minima_Volatilidad[0] )
print ('> La Máxima Volatilidad Anualizada fue de', '%.6s' % 
       (df['Volatilidad_14_Dias_Anualizada'].max()),
       '%', 'registrada el', Fecha_Maxima_Volatilidad[0])

# Obtenemos el promedio sobre el rango porcentual de los días negativos.
df['DiasNegativos'] = np.where(df['R. Aritmetico'] < 0, 100*(df['High']-df['Low'])/df['Low'],0)
df_dias_negativos = df.loc[df['DiasNegativos'] != 0]
DN = df_dias_negativos['DiasNegativos'].mean()
print('> Rango Medio días Negativos:','%.4s' % DN ,'%')

# Obtenemos el promedio del rango porcentgual de los días positivos.
df['DiasPositivos'] = np.where(df['R. Aritmetico'] > 0, 100*(df['High']-df['Low'])/df['Low'],0)
df_dias_positivos = df.loc[df['DiasPositivos'] != 0]
DP = df_dias_positivos['DiasPositivos'].mean()
print('> Rango Medio días Positivos:','%.4s' % DP ,'%')

# Calulamos el ratio del rango entre días positivos y negativos.
print('> Ratio RDN/RDP','%.4s'%(DN/DP),'%')

Tras ejecutar este grupo de instrucciones obtendremos las siguientes estadísticas:

Aquí podemos ver que la volatilidad anualizada del SPY fue del 18,16% (utilizando todo el historial). La volatilidad mínima se registró en agosto de 2017 con un valor del 2,45%. Y la volatilidad máxima en octubre de 2008 con un valor del 110,78%. También vemos que el rango medio porcentual de los días negativos fue de 1,42%, y el de los días positivos de 1,21%. El Ratio RDP/RDN nos indica que la volatilidad media de los días negativos fue concretamente un 17% superior a la volatilidad media de los días positivos. Una asimetría esperada en los índices de renta variable.

Llegados aquí damos por terminado el primer artículo de esta serie. Esperamos que te haya gustado. En breve continuaremos con las publicaciones. Si te gustaría estar informado de las novedades de QuantSpace puedes suscribirte a nuestro boletín de noticias.