Todas las entradas de: rvaquerizo

Longitud de las frases del Quijote con #rstats

Siempre he querido hacer cosas con Rstats y el Quijote y ayer se me ocurrió medir la longitud de las frases del Quijote y crear un histograma que describa esta longitud. Aunque confieso que no me lo he leído, me he quedado en el capítulo 7 u 8 (no recuerdo) el caso es que me pareció hipnótico con sus ritmos, es musical. Además tengo muchas ganas de meter mano al proyecto Gutemberg porque esos ritmos, esa musicalidad, el uso de palabras esdrújulas,… me llama la atención.
Bueno, al lío, todo el código está subido al repositorio por si lo queréis, pero hay algunas funciones y algunas ideas que me parecen interesantes.

library(dplyr)
library(ggplot2)
library(plotly)

#Leemos el fichero desde proyecto Gutemberg
ubicacion <- "https://www.gutenberg.org/cache/epub/2000/pg2000.txt"
quijote <- read.table (ubicacion,sep="\r", encoding="UTF-8")
quijote <- data.frame(quijote)
names(quijote) <- 'linea'

#Transformaciones e identificar el inicio del libro.
quijote <- quijote %>%
  mutate(linea = toupper(linea),
         inicio = grepl("EN UN LUGAR DE LA MANCHA",linea)>0)

Leemos directamente un txt desde Gutemberg y prefiero transformarlo en data frame para usar dplyr. Todas las palabras las pongo en mayúsculas e identifico donde empieza el Quijote, para evitar prólogos y demás. Ya tengo unos datos con los que poder trabajar:

#Marcamos lo que vamos a leer
desde <- which(quijote$inicio)
hasta <- nrow(quijote)

#Texto de trabajo
texto <- quijote[desde:hasta,1]

#El texto lo transformamos en una lista separada por espacios
texto_split = strsplit(texto, split=" ")

#Deshacemos esa lista y tenemos el data.frame
texto_col = as.character(unlist(texto_split))
texto_col = data.frame(texto_col)
names(texto_col) = 'palabra'

En este caso los datos los quiero de tal forma que disponga de un data frame con una sola variable que sea cada palabra del Quijote. Ahora voy a medir las frases identificando donde hay puntos en esas palabras:

#Identificamos donde tenemos puntos y un autonumérico del registro
texto_col <- texto_col %>% filter(!is.na(palabra)) %>%
  mutate(punto = ifelse(grepl('.',palabra,fixed=T),"FIN","NO"),
         posicion = row_number())

¿Qué se me ha ocurrido? Trabajar con autonuméricos, tengo identificados los puntos, ahora tengo que fijar una posición inicial y una posición final:

#Si unimos las posiciones con puntos con lag podemos calcular la longitud
pos_puntos1 <- filter(texto_col,punto=="FIN") %>% 
  select(posicion) %>% mutate(id = row_number())

pos_puntos2 <- pos_puntos1 %>% mutate(id = id + 1) %>%
  rename(posicion_final = posicion)

pos_puntos <- left_join(pos_puntos1,pos_puntos2) %>%
  mutate(longitud = ifelse(is.na(posicion_final), posicion, posicion - posicion_final))

Como no soy un tipo muy brillante opto por una opción sencilla de cruzar una tabla consigo misma, como me ponen los productos cartesianos “con talento”. La idea es seleccionar solo los registros que marcan el final de la frase, un autonumérico me marca cual es cada frase, ahora si hago una left join por el id de la frase y el id + 1 de la frase creo una especie de lag. La longitud de la frase será donde está el punto menos donde estaba el final de la anterior frase. Creo que me he explicado de pena, pero si veis el data frame final lo entenderéis mejor. Ahora ya pinto un histograma:

#GRaficamos la longitud
plot_ly(data = pos_puntos, x = ~longitud, type = "histogram") %>%
  layout(title = "Longitud de las frases del Quijote",
         xaxis = list(title = "Longitud"), yaxis = list(title = ""))

Y queda una gamma perfecta, yo diría que hasta bonita. Ahora quedaría identificar los parámetros de esta gamma y compararlos con otros libros, e incluso comparar lenguas. Pero esas tareas se las dejo a los “buenos”.

Ya no nos interesa el GDPR

GDPR en Google Trends un pico en mayo de 2018 y desaparece el interés. Será porque todas las compañías se han adecuado a ese marco general legislativo y no es necesario trabajar más en ello o será porque seguimos en modo “preproducción”. En cualquier caso es probable que esta legislación haya nacido anticuada, es probable que no la hayamos entendido e incluso es probable que no sea necesaria. 

Y recuerda, legislar también nos cuesta dinero. 

Gráfico con eje secundario en ggplot2

Los gráficos con eje secundario o con dos ejes son un tema que ya he puesto en el blog en varias ocasiones, hay un ejemplo con R que tenía sus problemas y hay un ejemplo con Python y matplotlib que particularmente me gusta por elegancia y sencillez. En esta entrada vamos a repetir el ejercicio y vamos a realizar un gráfico de columnas y líneas con 2 ejes, primario y secundario pero con ggplot2. Este tipo de gráficos son muy utilizados por los actuarios para representar frecuencias o siniestralidades y exposición. Para ilustrar el ejercicio vamos a emplear los mismos datos que usamos en el ejemplo con matplotlib pero vemos paso a paso como realizaríamos el gráfico Seguir leyendo Gráfico con eje secundario en ggplot2

Trucos simples para #rstats

En mi cuenta de twitter suelo poner algunos trucos sencillos de R, cosas que me surgen cuando estoy trabajando y que no me cuesta compartir en 2 minutos, por si puedo ayudar a alguien. Me acabo de dar cuenta que de verdad son útiles y que tenerlos en twitter desperdigados es un problema, así que he pensado en recopilarlos en una entrada del blog para que sea más sencillo buscarlos (incluso para mi). Aquí van algunos de esos trucos Seguir leyendo Trucos simples para #rstats

Muestrear no es pecado

Hace unos días que nació la web de Jose Luis Cañadas y ya está añadida a los enlaces de AyD, desde aquí nos hacemos eco de esta nueva bitácora:

https://muestrear-no-es-pecado.netlify.com/

Seguiremos muy de cerca a JL porque es un tipo de esos que merece ser escuchado por como mezcla pragmatismo y rigor dos de las principales características que tienen que tener los científicos de datos y que últimamente están vilipendiadas.

Las mujeres son más inteligentes pero no están en puestos directivos

A raíz de una noticia sobre la reestructuración del consejo directivo de un gran banco en España donde sólo una mujer ha sido elegida entre los 12 puestos de dirección general me ha surgido la oportunidad para explicar que es una distribución de probabilidad, que es una distribución normal y que es la media y la desviación típica.

Aquí tenéis en python un código que simula el IC de los hombres y el IC de las mujeres, no me he complicado mucho la vida ni he buscado datos al respecto pero leyendo un poco creo que deben de parecerse mucho a esto:

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns

hombres = np.random.normal(loc=60, scale=32, size=10000)
mujeres = np.random.normal(loc=70, scale=25, size=10000)

p1=sns.kdeplot(hombres, shade=True, color="r")
p1=sns.kdeplot(mujeres, shade=True, color="b")
sns.plt.show()

Bonito gráfico de densidades con varias variables hecho con seaborn y kdeplot de sintaxis sencilla que pone de manifiesto con unos datos simulados (con cierto talento) que en media la mujer es un 15% más inteligente que el hombre, pero la dispersión de la inteligencia del hombre es mayor, o como le he dicho a la persona que le he explicado que es la distribución normal, hay hombres muy tontos y hombres muy listos, muchos más que mujeres muy tontas y mujeres muy listas; así es la biología.

Pero ya os digo yo que la relación de hombres/mujeres inteligentes no es 1/12, así que esa importante entidad bancaria no está tomando decisiones correctas. Ellos verán lo que hacen pero no ha sido una medida inteligente, probablemente porque la ha tomado un hombre.

Mapa de México rápido (y sucio) y estático con #rstats

No sabía como mostraros el funcionamiento de getData del paquete raster para tener que evitaros ir a GADM y descargar los correspondientes mapas. Bueno, pues se me ha ocurrido hacer una entrada que tenga el menor número de líneas posibles y que genere un mapa. No me lo tengáis mucho en cuenta:

library(ggplot2)
library(raster)
library(dplyr)

#Obtenemos el mapa de GADM
mex <- getData("GADM", country = "MX", level = 2)

#El dato que vamos a pintar
prov <- data.frame(region=unique(mex@data$NAME_2))
prov$aleatorio <- runif(nrow(prov),0,100)

#Creamos el objeto mapa al que le añadimos el dato que necesitamos pintar
mex@data$name = mex@data$NAME_2
mex <- map_data(mex)
mex <- left_join(mex,prov)

#Pintamos el mapa
ggplot(data = mex, aes(x = long, y = lat, group = group)) +
  geom_polygon(aes(fill = aleatorio)) +
  scale_fill_continuous(low="white",high="blue") +
  labs(title = "Quick and dirty") + 
  theme_void()

Ahí lo tenéis getData se conecta a GADM donde vía ISO 3 le decimos que mapa queremos y el nivel que queremos y pintamos un mapa de México con ggplot2 en un pis pas. Comentad si no entendéis algo, hay miles de entradas que hacen lo mismo de forma más detallada pero en menos líneas ninguna. Saludos.

Porque no vamos a cobrar pensiones. Animación con R y pirámides de población

Estoy creando material para un módulo de un máster que voy a impartir y escribiendo sobre seguros de ahorro he llegado a crear esta animación:

Se trata de una animación con las pirámides de población de España desde 1975 hasta 2018 de 5 en 5 años. El sistema de pensiones español se basa en 5 principios:
1. principio de proporcionalidad
2. principio de universalidad
3. principio de gestión pública
4. principio de suficiencia
5. principio de reparto

La animación va directa contra el principio de reparto. En el sistema español nadie ha cotizado para garantizarse su pensión, los actuales trabajadores pagan las prestaciones de aquellos trabajadores jubilados. Si tras leer estas dos frases y mirar la animación sigues recelando de la migración de personas a España espero que tengas un buen plan de ahorro privado.

Esta animación está hecha con R, los datos están descargados del INE pero están ligeramente cocinados (no al estilo Tezanos). En https://github.com/analisisydecision/wordpress tenéis este Excel con el formato adecuado, el código empleado para realizar la animación está en https://github.com/analisisydecision/wordpress/blob/master/Piramide_poblacional.R Es un buen ejemplo de uso de plotrix y pyramid.plot espero que el código no tenga algún gazapo…

Data management con dplyr

Dos años con pandas y sckitlearn y ahora vuelvo a R. Y en mi regreso me propuse comenzar a trabajar con dplyr y mi productividad se está incrementando exponencialmente, creo que dplyr es LA HERRAMIENTA para el manejo de data frame con R, ni me imagino como puede funcionar sparlyr… Para aquellos que estéis iniciando vuestra andadura con R o para los que no estéis acostumbrados a dplyr he hecho una recopilación de las tareas más habituales que hago con esta librería. Se pueden resumir:

• Seleccionar columnas
• Seleccionar registros
• Crear nuevas variables
• Sumarizar datos
• Ordenar datos
• Uniones de datos

Como es habitual trabajamos con ejemplos data(iris); library(dplyr):

Seleccionar columnas select():

two.columns <- iris %>%
select(Sepal.Length,Sepal.Width)

columns = c(“Sepal.Length”,”Sepal.Width”)
two.columns <- iris %>%
select(columns)

Seleccionar registros filter():

setosa <- iris %>%
filter(Species==”setosa”)

species_to_select = c(“setosa”,”virginica”)
species <- iris %>%
filter(Species %in% species_to_select)
table(species$Species)

Crear nuevas variables mutate():

iris2 <- iris %>%
mutate(Sepal.Length.6 = ifelse(Sepal.Length >=6, “GE 6”, “LT 6”)) %>%
mutate(Sepal.Length.rela = Sepal.Length/mean(Sepal.Length))

Sumarizar group_by() summarize():

iris %>% group_by(Species) %>%
summarize(mean.Sepal.Length = mean(Sepal.Length),
sd.Sepal.Length = sd(Sepal.Length),
rows = n())

Ordenar datos arrange():

order1 <- iris %>%
arrange(Sepal.Length)

order2 <- iris %>%
arrange(desc(Sepal.Length))

iris %>% group_by(Species) %>%
summarize(mean.Sepal.Length = mean(Sepal.Length),
sd.Sepal.Length = sd(Sepal.Length),
rows = n()) %>%
arrange(mean.Sepal.Length)

Uniones de datos:

Inner_join():

iris2 <- iris %>%
mutate(id = row_number())

iris3 <- iris2 %>%
filter(Species==”setosa”) %>%
mutate(Sepal.Length.6 = ifelse(Sepal.Length >=6, “GE 6”, “LT 6”)) %>%
mutate(Sepal.Length.rela = Sepal.Length/mean(Sepal.Length)) %>%
select(id,Sepal.Length.6,Sepal.Length.rela)

iris4 <- iris2 %>% inner_join(iris3, by=c(“id”))

Left_join():

iris5 <- iris2 %>% left_join(iris3, by=c(“id”))

anti_join():

iris6 <- iris2 %>% anti_join(iris3)

Aquí tenéis una muestra de las posibilidades de dplyr y como se pueden combinar entre ellas. Creo que la sintaxis es bastante sencilla y se aprende con facilidad, si a mi no me esta costando mucho…

Crear una RESTful API con R con plumber

Podéis buscar info en la web acerca de lo que es una REST y una RESTful pero el objetivo de este trabajo es la creación de una API para “escorear” unos datos a partir de un modelo que hemos creado en R. Vamos a hacer lo más sencillo, un modelo de regresión lineal creado por R será guardado y una API con datos podrá llamar a este modelo mediante un cliente RESTful para obtener una predicción. Esta será la primera de una serie de entradas que le voy a dedicar a Carlos, un antiguo compañero mío y que me ha enseñado a desaprender y el primer guiño a Carlos será abandonar mi subversion local para conectar mi RStudio con GitHub, todo el trabajo que voy desarrollando lo tenéis en https://github.com/analisisydecision/Modelo1. Si echáis un vistazo al repositorio ya os podéis imaginar hacia donde irán encaminadas esta serie de entradas.

Bien, lo primero será crear y guardar el modelo con R:

#Programa de creación del modelo
Altura <-c(175,180,162,157,180,173,171,168,165,165)
Peso <-c(80,82,57,63,78,65,66,67,62,58)

modelo1 <- lm(Peso~Altura)
summary(modelo1)
save(modelo1, file = "modelo1/modelo1.rda")
#rm(modelo1)

Modelo de regresión lineal simple de alturas y pesos que guarda en la carpeta modelo1 el objeto con el modelo. Ahora quiero crear una API que, dada una altura, me estime el peso. Para ello creo un nuevo programa en R que debería llamarse despliegue pero que llamo depliegue_modelo1.R debido a que es bastante tarde. Este programa es una función para realizar una predicción que tiene el siguiente contenido:

library(jsonlite)

load("modelo1/modelo1.rda")

#* @post /prediccion
predict.peso <- function(Altura) {
data <- list(
Altura=Altura
)
prediccion <- predict.lm(modelo1, data )
return(prediccion)
}

Este código es el core de nuestra API a la que llamamos prediccion y que recibirá un json con ‘{“Altura”:XXX}’ y retornará la predicción del peso para esa altura. Y ahora viene plumber que es el “fontanero” que nos permite canalizar las llamadas a nuestra API prediccion. Esta llamada la hacemos con la función plumb:

library(plumber)
r <- plumb("depliegue_modelo1.R")
r$run(port=8000)

Ejecutado este código nuestra API esta funcionando en el puerto 8000 y sólo nos queda probarla y para ello yo recomiendo añadir una extensión de RESTClient a nuestro navegador habitual, en mi caso concreto es Chrome y he añadido una extensión Cliente de servicios Web RESTful y al ejecutarla tengo que modificar:

Como cuerpo de la solicitud pasamos un json con la altura y le damos a enviar, como respuesta debemos obtener [73.2291]. Lo que hemos hecho es una solicitud curl con la sintaxis:

curl -i -H Accept:application/json -X POST http://127.0.0.1:8000/prediccion -H Content-Type: application/json -d '{"Altura":175} '

Si ponéis esto en el terminal debe funcionar. ¿Qué os parece si empezamos a poner en producción modelos de R?