Daremos un repaso rápido a algunos conceptos básicos de estadística para centrarnos en la detección y reemplazo de outliers con R. Por el camino aprenderemos a manejar fechas y horas, a convertir una variable numérica en categórica, los diagramas de caja y bigotes, y a crear nuestras propias funciones.
Índice de contenidos
- 1. Introducción
- 2. Leyendo el CSV de la actividad física del usuario
- 3. Manejando fechas y horas en R
- 4. Un poquito de estadística
- 5. Convertir una variable numérica en categórica
- 6. Valores atípicos (outliers)
- 7. Diagrama de caja y bigotes
- 8. Reemplazo de Outliers con R
- 9. Conclusiones
- Enlaces y referencias
1. Introducción
Continuando con la serie de tutoriales sobre R y RStudio, daremos un repaso a algunos conceptos básicos de estadística. Como ejemplo vehicular vamos a usar el proyecto Geria-TIC.
Se trata de un proyecto que mediante un wearable (una pulsera) registra una serie de parámetros de los ancianos como la actividad física, y los periodos de sueño y así poder establecer la relación con las caídas, problemas de insomnio e incontinencia urinaria.
El estudio sigue la filosofía de OpenData, en este caso promovida por en Red.es, y que se enmarca dentro de la Normativa RISP (Reutilización de la Información del Sector Público). En el siguiente enlace se puede consultar la información OpenData del proyecto Geria-TIC así como la documentación de sus ficheros y campos.
2. Leyendo el CSV de la actividad física del usuario
En un artículo anterior ya vimos cómo leer un archivo CSV. En este ejemplo nos vamos a centrar en el fichero que recoge la actividad física del usuario: user_activity.csv
Los campos que están en este fichero son:
- userId: idenfiticador del usuario
- registerStartDate: fecha de inicio del registro
- registerFinalDate: fecha final del registro
- steps: pasos caminados
- distance: distancia recorrida en metros
- runDistance: distancia recorrida en metros (corriendo)
- calories: calorias gastadas
- deepSleepTime: tiempo de sueño profundo (en minutos)
- shallowSleepTime: tiempo de sueño superficial (en minutos)
- wakeTime: tiempo para despertar (en minutos)
- sleepStartDate: hora a la que comienza el sueño
- sleepStopDate: hora a la que finaliza la fase de sueño
Como siempre, lo primero será definir nuestro directorio de trabajo:
#seteamos el directorio de trabajo
setwd("/lab/adictosaltrabajo/rstudio/04")
# Vamos a leer un CSV del proyecto GeriaTIC
user_activity <- read.csv("./user_activity.csv")
View(user_activity)
Vemos que nos ha recuperado las fechas como factores, pero eso no es lo que queremos.
3. Manejando fechas y horas en R
Lo primero que haremos será leer el CSV pero indicando que no queremos los strings como factores. Luego procedemos a la conversión de tipos a Date para fechas, y a POSIXlt para fecha-hora.
user_activity <- read.csv("./user_activity.csv", stringsAsFactors = F)
# para convertir string a fechas
user_activity$registerStartDate
<- as.Date(user_activity$registerStartDate, format="%Y/%m/%d")
user_activity$registerFinalDate
<- as.Date(user_activity$registerFinalDate, format="%Y/%m/%d")
# para convertir strings a fecha-hora
user_activity$sleepStartDate
<- strptime(user_activity$sleepStartDate, format="%Y/%m/%d %H:%M:%S", tz = "UTC")
user_activity$sleepStopDate
<- strptime(user_activity$sleepStopDate, format="%Y/%m/%d %H:%M:%S", tz = "UTC")
#nos quedamos con una version simplificada
activity <- user_activity[,c("userId","registerStartDate","steps","distance","calories")]
Fijémonos que con la función as.Date() cambiamos en tipo y con el parámetro format le indicamos el patrón del formato en que debe leer. La función strptime() funciona de forma similar, pero además tiene otro parámetro opcional tz que nos permite indicar el TimeZone.
Finalmente nos quedamos con un dataFrame que es una versión simplificada con sólo 5 columnas, pues nos vamos a centrar en la parte de la actividad física, y no tanto en los parámetros de sueño.
4. Un poquito de estadística
Es inevitable, y en algún momento teníamos que recordar ciertos conceptos elementales de estadística. Pero tranquilo, que empezaremos por los que estudiaste en el colegio:
- Media: también conocida como media aritmética o promedio. En una muestra es la suma de todos los valores obtenidos dividido entre el número de repeticiones.
- Mediana: si ordenamos todos los valores obtenidos, se trata del valor que deja a un lado la mitad de resultados y al otro la otra mitad.
- Moda: es el valor que más se repite en el muestreo.
- Cuartiles: Al igual que la mediana, los cuartiles dejan resultados a un lado y al otro. Pero dividiendo el muestre en cuatro partes. Así, el 25% de los resultados son menores que el primer cuartil (Q1), el 50% son menores que Q2 (es la mediana), el 75% menosres que el tercer cuartil (Q3)
- Percentiles: Igual pero con tanto por cieto. El percentil 95 el valor en que el 95% de los resultados son menores que él
Depende lo que se estudie, la media es un valor representativo de la población estudiada. Pero hay casos en que no.
Por ejempplo, los sueldos. Imaginemos un equipo de 5 personas y el jefe, donde los sueldos son (20.000, 22.000, 30.000, 30.000, 50.000, 100.000).
La media sale: 42.000 sin embargo sólo un empleado y el jefe superan esa cifra.
En este caso es más representativa la mediana: si ordenamos de menos a mayor y nos quedamos con la cifra de en medio, la mediana es 30:000.
Cuando el sueldo del que más gana en una empresa puede ser varias veces el del que menos gana, la media puede distorsionar la información.
5. Convertir una variable numérica en categórica
Volviendo a nuestro dataFrame de actividad de los ancianos, vamos a fijarnos en el número de pasos que dan. Ya vimos cómo obtener un histograma con R que nos aporte información de en qué rangos se mueve esta variable con un simple vistazo.
# Histograma numero de pasos hist(activity$steps, main="Frecuencia de número de pasos")
Vamos a definir una categorización, que a priori nos parece razonable:
- hasta 2000 pasos es bajo
- de 2000 – 4000 pasos sería medio
- de 4000 – 6000 pasos estaría alto
- mayor de 6000 es muy alto
# Convertir la variable numerica "pasos" en categorica
# para ello definimos los puntos de corte
breakPoints <- c(0, 2000, 4000, 6000, Inf)
categories <- c("Low", "Medium", "High", "Very High")
# y cortamos la variable número de pasos segun esta categorizacion
activity$steps.F <- cut(activity$steps, breaks = breakPoints, labels = categories)
En la imagen vemos que hay una nueva columna steps.F de tipo Factor, que puede tener 4 posibles valores. Pero si nos movemos un poco por el dataFrame vemos que se repite mucho el valor «Low». A lo mejor no hemos definido bien los puntos de corte. Vamos a obtener una serie de medidas estadísticas sobre la variable de pasos:
summary(activity$steps)
Y por consola nos salen los valores:
> summary(activity$steps)
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.0 353.0 878.5 1687.3 1947.5 33218.0
Vemos que el valor máximo es 33.218 pasos, mientras que la media son 1.687 pasos. Pero el 50% de las veces se han andado 878 pasos o menos. Esta mediana está tan distante de la media (casi el doble), que ya nos indica que pasa algo raro.
Tenemos unos pocos valores muy altos que nos están tergiversando el estudio.
6. Valores atípicos (outliers)
En el ejemplo que nos traemos entre manos, vemos que Q1 está en 353 pasos y que Q3 está en 1947 pasos. Entre Q1 y Q3 sabemos que están el 50% de los valores obtenidos en el estudio. A esta distancia se le llama rango intercuantílico (IQR: InterQuantile Range).
IQR = Q3 – Q1 = 1947.5 – 353.0 = 1594.5 pasos
Se define como valor atípico leve aquel que dista 1,5 veces el rango intercuantílico por debajo de Q1 o por encima de Q3
q < Q1 – 1,5 · IQR o bien q > Q3 + 1,5 · IQR
y valor atípico extremo aquel que dista 3 veces el rango intercantílico por debajo de Q1 o por encima de Q3
q < Q1 – 3 · IQR o bien q > Q2 + 3 · IQR
De hecho, si sospechamos que tenemos outliers por arriba, vamos a calcular cual sería el umbral:
Q1 + 1,5 · IQR # 1947.5 + 1,5 · 1594.5 = 4339.25
Todos los valores del muestro que superen 4339.25 son outliers. De ahí para arriba, está muy por encima y se consideran valores atípicos.
¿Y por abajo? ¿Tendremos outliers por debajo?
Vamos a ver cual es el umbral inferior:
Q1 - 1,5 · IQR # 353.0 - 1,5 · 1594.5 = -1241.5
El umbral inferior sería una cifra negativa, y el valor más bajo que tenemos es 0, el mínimo. Así que no tenemos ningún valor inferior al umbral. En este caso no hay valores atípicos por debajo.
Podemos quitar los outliers de una forma muy fácil
# quitamos los outliers superiores (inferiores no hay)
steps <- activity[activity$steps < 4339.25, c("steps")]
# y al hacer un summary de la variable vemos que los valores han cambiado
> summary(steps)
Min. 1st Qu. Median Mean 3rd Qu. Max.
0.0 327.0 778.5 1126.2 1682.0 4337.0
Es normal que la media, mediana y cuartiles hayan cambiado, pues los outliers tenían mucho peso y distorsionaban los datos. Aún así, este muestreo que es un subconjunto del original, tiene pinta que no está bien balanceado y tendrá sus propios outliers. No hay más que ver la distancia que hay por encima del Q3.
7. Diagrama de caja y bigotes
Una forma fácil de detectar outliers es con un diagrama de caja y bigotes. Para ilustrar este tipo de gráficos nos vamos a fijar en los pasos que da un anciano en concreto (userID = 10039)
steps_10039 <- activity[activity$userId == 10039, c("steps")]
Nos fijamos que empieza recorriendo pocos pasos, pero que va en aument.
> summary(activity_10039$steps)
Min. 1st Qu. Median Mean 3rd Qu. Max.
0 2057 2646 2624 3188 5842
Vamos a dibujar este diagrama de caja y bigotes.
boxplot(steps_10039)
En el diagrama vemos la caja, cuyo borde superior es Q3 y el inferior es Q1. Entre medias están el 50% de las ocurrencias. La altura de la caja es el rango intercuartílico. Y el bigote, la línea gruesa, es la mediana. Por encima y debajo se ven dos límites, que son los umbrales para los valores atípicos.
El superior es:
Q3 + 1,5 · IQR = 4884,5 pasos.
Y todo lo que esté por encima es un outlier.
Pero ¿y el límite inferior? Según lo que sabemos, debería ser: Q1 – 1,5 · IQR = 360,5 pasos. Sin embargo la función boxplot() nos pinta esa línea cercana a los 1.000 pasos. ¿Por qué? Muy sencillo. El día que menos andó fueron 1.045 pasos. De los 306 días registrados, 12 no andó nada (cero pasos). Así que el algoritmo que usa esta función considera esos días con cero pasos son valores atípicos, y dibuja la línea inferior en los 1.045 pasos y unos circulitos en el 0, para los outliers.
Podemos guardar los outliers en una variable a la vez que dibujamos el diagrama de caja y bigotes.
outliers_10039 <- boxplot(steps_10039)$out > outliers_10039 [1] 5344 4958 5017 5842 5139 0 0 0 0 0 0 0 0 0 0 0 0
Volviendo al resto de usuarios y ojeando el dataFrame, me doy cuenta, que casi todos los valores entre 20.000 y 30.000 pasos los da el usuario 10041. Este usuario me está distorsionando el estudio comparándolo con todos los demás ancianos. Es el andarín del grupo. Lo voy a sacar del dataFrame.
activity <- activity[activity$userId != 10041,]
Lo interesante sería poder dibujar este gráfico de caja y bigotes para cada usuario. Se puede hacer poniendo la variable «pasos» en función de la variable «usuario».
boxplot(activity$steps ~ activity$userId)
Cuando hacemos referencia a varias variables del mismo dataFrame es preferible indicarlo y nos ahorramos la sintasix del $.
# este boxplot sería equivalente al anterior
boxplot(steps ~ userId,
data = activity,
main = "Pasos por usuario")
8. Reemplazo de Outliers con R
Ya hemos visto como aprender a detectar outliers usando la función boxplot(), pero a lo mejor queremos definir nosotros nuestros propios umbrales de qué valores se quedan fuera, y reemplazarlos por la media o la mediana.
Nos definimos nuestra propia función que recibe como parámetros el conjunto de datos, el umbral inferior y el umbral superior. Y en la función decidimos reemplazar por la media aquellos valores atípicos que estén por debajo del umbral inferior, y por la mediana aquellos que estén por encima del umbral superior.
# reemplazo de outliers con R
outliersReplace <- function(data, lowLimit, highLimit){
data[data < lowLimit] <- mean(data)
data[data > highLimit] <- median(data)
data #devolvemos el dato
}
steps_10039_2 <- outliersReplace(steps_10039, 1045, 4884.5)
summary(steps_10039)
summary(steps_10039_2)
Y por consola veríamos:
> summary(steps_10039)
Min. 1st Qu. Median Mean 3rd Qu. Max.
0 2057 2646 2624 3188 5842
> summary(steps_10039_2)
Min. 1st Qu. Median Mean 3rd Qu. Max.
1045 2117 2645 2684 3147 4704
Podemos pintar un gráfico al lado del otro
par(mfrow = c(1,2)) boxplot(steps_10039, main = "Sin reemplazo de outliers con R") boxplot(steps_10039_2, main = "Con reemplazo de outliers con R")
Bueno, hemos definido nuestra primera función. Todo analista de datos tiene su buena colección de funciones a la que recurre habitualmente para tratar casuísticas habituales. R no deja de ser como cualquier otro lenguaje: cualquier tarea que se haga más de una vez es susceptible de ser programada.
Si volvemos al dataframe de activity, cada usuario tenía unos valores atípicos distintos. Nos interesará hacer una función que dado un dataFrame, particione los datos por un ID, y para cada ID encuentre los outliers de una de las columnas, reemplazando su valor por la mediana. Pero para no destruir los datos originales, podemos duplicar la columna con los datos a estudiar llamándola igual acabada en «.R» para indicar que tiene reemplazos.
Aunque tengamos en mente el dataFrame de activity y pensemos en los pasos por usuario, vamos a intentar hacer esta función genérica, que nos sirva para hacer el reemplazo de outliers con R sea cual sea el dataFrame, la columna por la que particionar los datos y la columna numérica a tratar.
# df es el dataFrame que recibimos (ej. activity)
# colNameData es la columna de los datos (ej. "steps")
# colNameBy es la columna por la que trocearemos (ej. "userId")
outliersReplace <- function(df, colNameData, colNameBy){
# creamos una nueva columna llamada igual que colNameData pero con .R
colNameData.R <- paste(colNameData, "R", sep=".")
df[colNameData.R] <- df[colNameData]
# obtenemos los IDs por los que partir el dataframe
IDs <- unique(df[,c(colNameBy)])
for (id in IDs){
data <- df[df[colNameBy] == id, c(colNameData) ]
Q <- quantile(data)
minimo <- Q[1] # valor minimo
Q1 <- Q[2] # primer cuartil
Me <- Q[3] # mediana
Q3 <- Q[4] # tercer cuartil
maximo <- Q[5] # valor maximo
IQR <- Q3 - Q1
lowLimit <- max(minimo, Q1 - 1.5*IQR)
highLimit <- min(maximo, Q3 + 1.5*IQR)
# todos los valores donde colNameBy es igual a id
# y el valor de colNameData es > Q3 + 1.5 * IQR
# lo reemplazamos por la mediana
df[df[colNameBy] == id & df[colNameData] > highLimit, c(colNameData.R)] <- Me
# lo mismo para el umbral inferior
df[df[colNameBy] == id & df[colNameData] < lowLimit, c(colNameData.R)] <- Me
cat(paste("El", colNameBy, id, "la mediana(", colNameData, ") ==", Me, "\n", sep=" " ))
}
df # devolvemos el valor del dataFrame
}
Y si ahora llamamos a la función, obtenemos por consola:
> activity <- outliersReplace(activity,"steps","userId") El userId 10036 la mediana( steps ) == 354 El userId 10037 la mediana( steps ) == 4345 El userId 10038 la mediana( steps ) == 0 El userId 10039 la mediana( steps ) == 2646.5 El userId 10043 la mediana( steps ) == 92 El userId 10044 la mediana( steps ) == 1141.5 El userId 10045 la mediana( steps ) == 450 El userId 10047 la mediana( steps ) == 155.5 El userId 10048 la mediana( steps ) == 1783 El userId 10056 la mediana( steps ) == 13.5 El userId 10057 la mediana( steps ) == 525 El userId 10058 la mediana( steps ) == 1135 El userId 10059 la mediana( steps ) == 633 El userId 10060 la mediana( steps ) == 939 El userId 10062 la mediana( steps ) == 1312
Y si vemos el dataFrame se observa la nueva columna que ha reemplazado los valores atípicos por la mediana de ese usuario.
Ahora podemos pintar los diagramas de caja y bigotes del número de pasos sin reemplazo y con reemplazo de outliers con R por la mediana de cada usuario.
par(mfrow = c(2,1)) # para ponerlos uno encima de otro boxplot(steps ~ userId, data = activity, main = "Sin reemplazo") boxplot(steps.R ~ userId, data = activity, main = "Con reemplazo")
9. Conclusiones
Lo principal que hemos aprendido en este tutorial es a detectar y a hacer reemplazo de outliers con R, pero también:
- hemos repasado algunos conceptos básicos de estadística
- aprendimos a manejar fechas y horas con R y a convertir tipos
- a convertir una variable numérica en factores (variables categóricas)
- vimos un nuevo tipo de gráfico: diagrama de caja y bigotes
- detectamos los valores atípicos o outliers
- y aprendimos a hacer el reemplazo de outliers con R por la media o mediana
- se vio cómo poner dos gráficos, uno al lado del otro, o uno encima de otro
- y creamos nuestras propias funciones para trabajar en bloque sobre el dataFrame
Enlaces y Referencias
- El código fuente de este tutorial en GitHub
https://github.com/eContento/rstudio - Wikipedia: medidas de tendencia central
- Método de detección temprana de Outliers (PDF). Trabajo de Grado de Juan Gabriel Moreno Castellanos
