Practic_1RMD.Rmd

---
title: "Tarea 1 - Algoritmos de clasificación"
author: "Jordi Vanrell Forteza"
date: "26/5/2021"
output: word_document
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = T, include = T, message = F, warning = F)
```

La tarea requiere escoger tres archivos de cada uno de los tres grupos provistos. Los grupos se corresponden con acciones de capitalización pequeña (nombrados _Stock20x.txt_), media (_Stock10x.txt_) y grande (_Stockx.txt_). Se descargan los ficheros y se meten en un mismo directorio. Se seleccionan aleatoriamente tres números entre 1 y 5 que marcarán la elección de las acciones en los tres grupos. Se leen con un bucle y se almacenan los nueve en una lista.

```{r}
# Se escogen 3 acciones de cada grupo (al azar)
set.seed(603) # semilla para completa reproducibilidad
idx <- as.character(sample(1:5, size = 3, replace = F))
idx <- c(paste0("", idx), paste0("10", idx), paste0("20", idx))

Stock <- vector(mode = "list", length = length(idx))
for (i in idx){ # Bucle de lectura de los archivos; se almacenan en lista
  Stock[[which(idx==i)]] <- read.table(paste0("Stock", i, ".txt"), header = T)
  rm(i)
}
```

Dicha lista consiste en 9 tablas de datos ordenadas según el índice de posición del vector `idx`: 5, 3, 1, 105, 103, 101, 205, 203 y 201.

# (a) Programario para la aplicación de los algoritmos de clasificación de transacciones.
```{r}
require(tidyverse) # Paquete para la manipulación de de data frames
```

Se procede a clasificar las transacciones según los algoritmos *Tick Rule*, *Quote Rule* y *Lee & Ready*. Para ello se arma un bucle que los implementará en cada una de las transacciones para cada activo. Debido a la propia idiosincasia de los algoritmos resulta más adecuado aplicarlos a través de bucles sobre objetos de tipo vector que sobre datos estructurados en *data frames*. Por ello se construyen tres bucles internos, uno por algoritmo, dentro de los cuales se codifican las reglas de clasificación mediante funciones condicionales. En cada caso se transforman los vectores en columnas de los *data frames* originales.

```{r}
for (d in 1:length(Stock)){
  Stock[[d]] <- Stock[[d]]  %>%
    arrange(day, time) %>% 
    # se asegura que las órdenes están correctamente ordenadas
    mutate(Q = round((ask+bid)/2, 4)) 
  # se obtiene el valor como la media de ask y bid
  
  day <- Stock[[d]]$day # se vectoriza day
  price <- Stock[[d]]$price # se vectoriza price
  
  # Tick Rule
  TR = c(); TR[1] = NA # se inicializa el vector con valor NA
  for (i in 2:length(price)){ # implementación del algoritmo
    TR[i] <- ifelse(day[i] == day[i-1] & price[i] > price[i-1], 1,
                    ifelse(day[i] == day[i-1] & price[i] < price[i-1], -1,
                           ifelse(day[i] == day[i-1] & price[i] == price[i-1],
                                  TR[i-1], NA)))
  }
  Stock[[d]]$TR <- TR; rm(TR) # se agrega el vector como columna al dataset
  
  # Quote Rule
  Q <- Stock[[d]]$Q # se vectoriza el valor
  QR = c() # se inicializa un vector vacío
  for (i in 1:length(price)){ # se implementa el algoritmo
    QR[i] <- ifelse(price[i] > Q[i], 1,
                    ifelse(price[i] < Q[i], -1, NA))
  }
  Stock[[d]]$QR <- QR; rm(QR) # se agrega el vector como columna al dataset
  
  # Lee & Ready
  LR = c(); LR[1] = ifelse(price[1] > Q[1], 1, # se inicializa el vector
                           ifelse(price[1] < Q[1], -1, NA))
  for (i in 2:length(price)){
    LR[i] <- ifelse(price[i] > Q[i], 1,
                    ifelse(price[i] < Q[i], -1, 
                           ifelse(price[i] == Q[i] & day[i] == day[i-1], 
                                  LR[i-1], NA)))
    rm(i)
  }
  Stock[[d]]$LR <- LR; rm(LR, day, price, Q, d)
  # se agrega el vector como columna al dataset
}
```

# (b) Precisión diaria por activo y algoritmo
```{r}
require(scales)
```

En primera instancia se requieren los datos de precisión de los algoritmos para cada activo y día del rango de los datos. Se define la precisión como la ratio de transacciones correctamente clasificadas entre el total de transacciones clasificadas (no se toman en cuenta las que el algoritmo no clasificada en cada caso). Para obtener esta información es necesario filtrar por día en cada una de las tablas de la lista; así se obtienen tres datos por tabla y día, uno por algoritmo. Estos tres se agregan y se almacenan en un *data frame*. El proceso se repite con un bucle hasta que el *data frame* contiene los datos para todos los días. La tabla completa, con filas y columnas convenientemente nombradas, se almacena en una lista y todo el proceso anterior se repite en bucle hasta que se ha efectuado para todos los activos. El resultado es una lista de nueve *data frames* con los datos de precisión requeridos, uno por cada activo.

```{r}
# se inicializa la lista de almacenamiento:
tabla_i_list <- vector(mode = "list", length = length(Stock))
for (s in 1:length(Stock)){
  subtabla_i <- data.frame() # se inicializa subtabla
  for (d in unique(Stock[[s]]$day)){
    stock_parse <- Stock[[s]] %>% filter(day == d) # filtro de día
    acc_f = function(R){ 
      # se define función interna de precisión (R es parámetro para algoritmo)
      mean(stock_parse[, R]==stock_parse[, "buysell"], na.rm = T)
    } # Aplicación de las funciones
    acc_TR = acc_f("TR"); acc_QR = acc_f("QR"); acc_LR = acc_f("LR")
    acc_df = data.frame(acc_TR, acc_QR, acc_LR) # se introducen en data.frame
    subtabla_i <- rbind(subtabla_i, acc_df) 
    # se agrega lo anterior a subtabla_i
    rm(d, acc_TR, acc_QR, acc_LR, acc_df, acc_f, stock_parse)
  }
  tabla_i_list[[s]] <- subtabla_i # se agrega la subtabla_i a tabla_i_list
  var <- colnames(tabla_i_list[[s]]) # se extraen nombres de columnas
  rownames(tabla_i_list[[s]]) <- paste0("day ", unique(Stock[[s]]$day)) 
  # se renombran las filas
  colnames(tabla_i_list[[s]]) <- paste0("Stock", idx[s], ", ", var) 
  # se añade nombre del stock a las columnas
  rm(subtabla_i, s, var)
}
```

A modo de ejemplo, en la Tabla I se disponen los resultados para el activo `Stock5`. Las filas contienen la dimensión día y en las columnas los algoritmos. En caso de querer explorar los resultados de los demás activos puede hacerse mediante la instrucción `tabla_i_list[[x]]`, sustituyendo la x por el índice del activo deseado. En estas los resultados están en tanto por uno.

```{r}
tabla_i <- tabla_i_list[[1]] 

var <- colnames(tabla_i) # Se disponen los datos como % con 2 dec.
tabla_i[, var] <- lapply(tabla_i[, var], percent_format(accuracy = .01))
rm(var)
```
```{r, eval=FALSE}
# Por cuestiones de compilación la tabla (esta y las siguientes) se ha generado de forma separada,de acuerdo con el código que sigue. Este es un código que NO puede compilarse en Rmd.
library(officer)
library(flextable)
library(magrittr)
# Para exportar tablas a formato Word
ft <- flextable(data = tabla_i %>% add_rownames()) %>% 
  theme_zebra %>% autofit
ft <- set_caption(ft, caption = "Tabla I", style = "Table Caption")
# Crea un archivo temp
tmp <- tempfile(fileext = ".docx")
# Crea un documento docx
read_docx() %>% body_add_flextable(ft) %>% print(target = tmp)
browseURL(tmp) # abre el documento
```

-AQUÍ VA TABLA I-

En segundo lugar se requiere la media y la desviación típica de las precisiones diarias por activo y algoritmo. Para ello se usa como fuente de partida la lista de tablas a partir de la cual se ha generado la Tabla I. Para cada una de ellas, y tras definir funciones para la media y la desviación típica, se calculan ambos indicadores para los tres algoritmos. Los resultados se disponen en la fila de un *data frame* donde van agregándose en bucle los de los nueve activos.

```{r}
tabla_ii <- data.frame() # Se inicializa tabla_ii como tabla vacía
for (s in 1:length(tabla_i_list)){
  var <- colnames(tabla_i_list[[s]]) 
  # Se toman los nombres de las columnas de la subtabla
  mean_f <- function(pos){ # Se define una función general para la media
    mean(tabla_i_list[[s]][, var[pos]], na.rm = T)
  }
  sd_f <- function(pos){ # Se define una función general para la desv típica
    sd(tabla_i_list[[s]][, var[pos]], na.rm = T)
  }
  mean_acc_TR <- mean_f(1); sd_acc_TR <- sd_f(1) # Aplicación de funciones a TR
  mean_acc_QR <- mean_f(2); sd_acc_QR <- sd_f(2) # Aplicación de funciones a QR
  mean_acc_LR <- mean_f(3); sd_acc_LR <- sd_f(3) # Aplicación de funciones a LR
  subtabla_ii <- data.frame(mean_acc_TR, sd_acc_TR, mean_acc_QR, 
                            sd_acc_QR, mean_acc_LR, sd_acc_LR)
  tabla_ii <- rbind(tabla_ii, subtabla_ii) 
  # Introducción en data frame y agregación
  rm(var, subtabla_ii, s, mean_f, sd_f, mean_acc_TR, sd_acc_TR, mean_acc_QR, 
     sd_acc_QR, mean_acc_LR, sd_acc_LR)
}
rownames(tabla_ii) <- paste0("Stock", idx) # Se renombran las filas
var <- colnames(tabla_ii) # Se disponen los datos como porcentajes con 2 dec
tabla_ii[, var] <- lapply(tabla_ii[, var], percent_format(accuracy = .01))
rm(var); tabla_ii
```

En la Tabla II se presentan los resultados. Cafa fila representa un activo y en las columnas se representan la media y la esviación típica en porcentaje.

```{r, eval=FALSE}
# Para exportar tablas a formato Word
ft <- flextable(data = tabla_ii %>% add_rownames()) %>% 
  theme_zebra %>% autofit
ft <- set_caption(ft, caption = "Tabla II", style = "Table Caption")
# Crea un archivo temp
tmp <- tempfile(fileext = ".docx")
# Crea un documento docx
read_docx() %>% body_add_flextable(ft) %>% print(target = tmp)
browseURL(tmp) # abre el documento
```

-AQUÍ VA LA TABLA II-

# (c) Precisión sobre el volumen por activo y algoritmo para las diferentes franjas temporales
```{r}
require(lubridate)
```

En este caso se requieren los datos de precisión de los algoritmos para cada activo y franja temporal de media hora medido sobre el volumen (y no sobre el número de transacciones). En este caso la precisión se mide como la suma de los mínimos de los volúmenes reales y estimados de compra y de venta ponderados sobre el volumen total intercambiado. Esta aproximación sobreestima la precisión real de los algoritmos en la medida en que los errores de identificación en las compras se compensan con los errores en las ventas. Sin embargo, vienen a ilustrar la capacidad de estos para predecir la proporción de compras y ventas en una franja temporal determinada. Del mismo modo que en (b), los resultados se almacenan en una lista de *data frames*, una tabla por activo.

Una vez más los resultados se obtienen a través de un triple bucle anidado. Dado un activo, para asignar cada operación a su franja temporal correspondiente se usan la variable time, los datos de fecha y la referencia de 34200 segundos correspondiente con la apertura del mercado a las 9:30. El formato fecha del paquete `lubridate` permite funciones de redondeo de fecha y hora a los 30 minutos. De esta información puede extraerse y agregarse la información sobre horas y minutos. El proceso únicamente requiere una corrección adicional: en ocasiones hay transacciones registradas a las 9:30:00 y 0 -ésimas, por lo que, por practicidad, se asignan al intervalo de las 10H 0M.

Los bucles interiores se hacen depender de los valores vectorizados de los algoritmos y las franjas temporales de media hora. Para cada combinación de los valores de dichos vectores se agregan los volúmenes en una matriz de confusión interna a partir de la cual se calcula la precisión tal y como se ha definido en este apartado. Los resultados se agregan primero por franja horaria para un algoritmo dado, a continuación se yuxtaponen las columnas correspondientes a los diferentes algoritmos para un mismo activo y por último se agregan las tablas, con filas y columnas renombradas, a la lista principal.

```{r}
alg <- c("TR", "QR", "LR")
for (d in 1:length(Stock)){
  for (a in alg){ # se reemplazan todos los NA por 0
    Stock[[d]][,a][is.na(Stock[[d]][,a])] <- 0
  }
}
```
```{r}
tabla_iii_list <- vector(mode = "list", length = length(Stock))
for (d in 1:length(Stock)){
  this_1 <- Stock[[d]] %>%
    mutate(time_stamp = ymd(paste0("2011-06-", day)) + hms("09:30:00") + time-34200,
           # 34200 son los segundos pasados la medianoche; 
           # time_stamp está redondeado a la baja en las -ésimas de segundo.
           ceil_t = ceiling_date(time_stamp, unit = "30 minutes"), 
           # se redondea al alza en intervalos de 30 minutos.
           ceil_halfh = hm(paste0(hour(ceil_t),":",minute(ceil_t)))) %>% 
    # se extren las horas y minutos
    select(-ceil_t)
  # Es posible que haya operaciones a las 9:30:00 y 0 -ésimas de s (no antes) 
  # así que, por practicidad, se asignan a las 10:00:00:
  this_1$ceil_halfh[this_1$ceil_halfh==hms("9H 30M 0S")] <- hms("10H 0M 0S")
  # Se construye un vector de texto con las franjas de tiempo:
  time <- unique(as.character(this_1$ceil_halfh))
  alg <- c("TR", "QR", "LR") # Se hace lo mismo con los algoritmos
  acc_at <- data.frame(aux = rep(0, length(time))) 
  # se inicializa el data.frame con una columna auxiliar
  for (a in alg){
    acc_t = data.frame() # se inicializa un data frame
    group_cols <- c("ceil_halfh", "buysell", a) 
    # criterio de agrupación dependiente de a
    for (t in time){ # bucle para la agregación de volúmenes según franja de t
      this_2 <- this_1 %>%
        group_by_at(group_cols) %>%
        summarise(vol = sum(vol)) %>%
        filter(ceil_halfh == t)
      this_2 <- this_2[this_2[,a]!=0,]
      # Precisión de un algoritmo dado para una franja de tiempo dada
      acc = (min(sum(this_2$vol[this_2[, "buysell"]==1]), 
                 sum(this_2$vol[this_2[, a]==1])) +
               min(sum(this_2$vol[this_2[, "buysell"]==-1]), 
                   sum(this_2$vol[this_2[, a]==-1])))/
        (sum(this_2$vol[this_2[, "buysell"]==1]) + 
           sum(this_2$vol[this_2[, "buysell"]==-1]))
      acc = data.frame(acc) # conversión de data.frame
      acc_t = rbind(acc_t, acc) 
      # agregación de resultados de t para un algoritmo dado
    }
    acc_at <- cbind(acc_at, acc_t) # agreg de algoritmos para un activo dado
  }
  acc_at <- acc_at[,2:4] # se desecha la columna auxiliar
  colnames(acc_at) <- paste0("Stock", idx[d], ", ", alg) 
  # se renombran las columnas con el nombre del activo
  rownames(acc_at) <- time # se renombran las filas con las franjas de tiempo
  tabla_iii_list[[d]] <- acc_at # se agrega la subtabla a la lista de tablas
  rm(a, acc, acc_at, acc_t, alg, d, group_cols, t, this_1, this_2, time)
}
```

A modo de ejemplo, en la Tabla III se disponen los resultados para el activo `Stock5`. Cada fila es un afranja de media hora y en las columnas se representan los tres algoritmos. En caso de querer explorar los resultados de los demás activos puede hacerse mediante la instrucción `tabla_iii_list[[x]]`, sustituyendo la x por el índice del activo deseado. En estas los resultados están en tanto por uno.

```{r}
tabla_iii <- data.frame() # Se inicializa tabla_ii como tabla vacía
for (s in 1:length(tabla_iii_list)){
  var <- colnames(tabla_iii_list[[s]]) # Se toman los nombres de las columnas de la subtabla
  mean_f <- function(pos){ # Se define una función general para la media
    mean(tabla_iii_list[[s]][, var[pos]], na.rm = T)
  }
  mean_acc_TR <- mean_f(1) # Aplicación de funciones a TR
  mean_acc_QR <- mean_f(2) # Aplicación de funciones a QR
  mean_acc_LR <- mean_f(3) # Aplicación de funciones a LR
  subtabla_iii <- data.frame(mean_acc_TR, mean_acc_QR, mean_acc_LR)
  tabla_iii <- rbind(tabla_iii, subtabla_iii) # Introducción en data frame y agregación
  rm(var, subtabla_iii, s, mean_f, mean_acc_TR, mean_acc_QR, mean_acc_LR)
}
rownames(tabla_iii) <- paste0("Stock", idx) # Se renombran las filas
var <- colnames(tabla_iii) # Se disponen los datos como porcentajes con 2 decimales
tabla_iii[, var] <- lapply(tabla_iii[, var], percent_format(accuracy = .01))
rm(var)
```
```{r, eval=FALSE}
# Para exportar tablas a formato Word
ft <- flextable(data = tabla_iii %>% add_rownames()) %>% 
  theme_zebra %>% autofit
ft <- set_caption(ft, caption = "Tabla III", style = "Table Caption")
# Crea un archivo temp
tmp <- tempfile(fileext = ".docx")
# Crea un documento docx
read_docx() %>% body_add_flextable(ft) %>% print(target = tmp)
browseURL(tmp) # abre el documento
```

-AQUÍ VA LA TABLA III-

# (d) Conclusiones

La precisión en la clasificación de las transacciones individuales varía entre activos y de un día a otro (Tabla I), pero puede decirse, a grandes rasgos, que *Quote Rule* clasifica mejor que *Tick Rule* cuando únicamente se tienen en cuenta las transacciones que consiguen clasificar (Tabla II). Sin embargo, la cantidad de transacciones que *Quote Rule* deja sin clasificar es considerable. Por ejemplo, para el activo `Stock5`, el `r round(100*mean(Stock[[1]]$QR==0),2)`% de las transacciones quedan sin clasificar por *Quote Rule*, pero solamente es el `r round(100*mean(Stock[[1]]$TR==0),2)`% por *Tick Rule*. El algoritmo de *Lee & Ready* intenta subsanar esta limitación aplicando *Tick Rule* donde *Quote Rule* no clasifica. El resultado es un empeoramiento de la precisión de *Quote Rule* a cambio de una proporción mayor de transacciones clasificadas. Comparadas *Tick Rule* y *Lee & Ready*, que clasifican una proporción similar de transacciones, este último supone una mejora media en la precisión del 1.4%. En cualquier caso, esta mejora va acompañada de un aumento de la desviación típica (`r round(mean(tabla_ii$sd_acc_TR),2)`% frente al `r round(mean(tabla_ii$sd_acc_LR),2)`%), con lo que resulta difícil aventurar cuál de los dos clasificaría mejor una serie cualquiera.

Por otro lado, se obtienen resultados de precisión sobre el volumen generalmente por encima del 90% y, en muchos casos, cercanos al 100% cuando se consideran volúmenes agregados y no transacción por transacción (Tabla III). Estos resultados dejan patente la capacidad de los tres algoritmos de predecir con un grado de acierto muy alto la proporción de compras y ventas producidas en un periodo de tiempo dado.