UN EJEMPLO DETALLADO DE PRE-PROCESAMIENTO EN R - IMPUTACIÓN DE DATOS FALTANTES

Imputacion de datos faltantes

En este artículo, vamos a hacer todo el pre-procesamiento de una base de datos. Nos ocuparemos de las redundancias entre distintas variables, de los outliers y, como no, de los datos faltantes.

Este no es un proceso único. Distintas aproximaciones darían lugar a resultados diferentes. En todo caso, sirve de guía genérica para revisar los pasos que habría que dar.

Para la parte final de la imputación hemos usado un paquete de R llamado mice.

Carga de los datos y librerías necesarias

Librerías que necesitaremos:

#install.packages("mice")
#install.packages("DataExplorer")
#install.packages("mltools")
#install.packages("ggcorrplot")
library(DataExplorer) #plot missing esta aqui
library(MASS) #ordinal logistic regression
library(tidyverse) #glimpse, plot_missing
library(mice)
library(corrplot) #gráfica de la correlación
library(ggcorrplot)
library(mltools)
library(plotrix)
library(caret)
library(caTools)
library(nnet)
library(Hmisc)
library(rcompanion) #coeficiente V de crammer

Para ilustrar los distintos ejemplos de este artículo, usaremos los datos Melbourne Housing Dataset.

Vamos a leer los datos y les echaremos un primer vistazo:

melb_data <- read.csv("melb_data.csv")
glimpse(melb_data)
## Rows: 13,580
## Columns: 21
## $ Suburb        <chr> "Abbotsford", "Abbotsford", "Abbotsford", "Abbotsford", ~
## $ Address       <chr> "85 Turner St", "25 Bloomburg St", "5 Charles St", "40 F~
## $ Rooms         <int> 2, 2, 3, 3, 4, 2, 3, 2, 1, 2, 2, 3, 2, 2, 1, 2, 3, 3, 3,~
## $ Type          <chr> "h", "h", "h", "h", "h", "h", "h", "h", "u", "h", "u", "~
## $ Price         <dbl> 1480000, 1035000, 1465000, 850000, 1600000, 941000, 1876~
## $ Method        <chr> "S", "S", "SP", "PI", "VB", "S", "S", "S", "S", "S", "VB~
## $ SellerG       <chr> "Biggin", "Biggin", "Biggin", "Biggin", "Nelson", "Jelli~
## $ Date          <chr> "3/12/2016", "4/02/2016", "4/03/2017", "4/03/2017", "4/0~
## $ Distance      <dbl> 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2.5, 2~
## $ Postcode      <dbl> 3067, 3067, 3067, 3067, 3067, 3067, 3067, 3067, 3067, 30~
## $ Bedroom2      <dbl> 2, 2, 3, 3, 3, 2, 4, 2, 1, 3, 2, 3, 2, 2, 1, 2, 3, 2, 3,~
## $ Bathroom      <dbl> 1, 1, 2, 2, 1, 1, 2, 1, 1, 1, 2, 2, 2, 1, 1, 1, 2, 1, 1,~
## $ Car           <dbl> 1, 0, 0, 1, 2, 0, 0, 2, 1, 2, 1, 2, 1, 1, 1, 2, 1, 1, 0,~
## $ Landsize      <dbl> 202, 156, 134, 94, 120, 181, 245, 256, 0, 220, 0, 214, 0~
## $ BuildingArea  <dbl> NA, 79, 150, NA, 142, NA, 210, 107, NA, 75, NA, 190, 94,~
## $ YearBuilt     <dbl> NA, 1900, 1900, NA, 2014, NA, 1910, 1890, NA, 1900, NA, ~
## $ CouncilArea   <chr> "Yarra", "Yarra", "Yarra", "Yarra", "Yarra", "Yarra", "Y~
## $ Lattitude     <dbl> -37.7996, -37.8079, -37.8093, -37.7969, -37.8072, -37.80~
## $ Longitude     <dbl> 144.9984, 144.9934, 144.9944, 144.9969, 144.9941, 144.99~
## $ Regionname    <chr> "Northern Metropolitan", "Northern Metropolitan", "North~
## $ Propertycount <dbl> 4019, 4019, 4019, 4019, 4019, 4019, 4019, 4019, 4019, 40~

Puede verse que hay variables de todo tipo, tanto categóricas como numéricas, y también una importante colección de valores N/A con los que tendremos que lidiar.

Algunas de las variables - por el nombre que tienen - parecen redundantes. Tenemos, por ejemplo, información sobre el número de habitaciones, de baños y el área. Todas estas parecen decir, en esencia, lo mismo. Esto lo iremos examinando sobre la marcha.

Significado de las variables

Es importante, antes de continuar, conocer qué significa cada variable, puesto que nuestra intención es quitar columnas redundantes para dejar la base de datos lo más simplificada posible. Lo tenemos descrito aquí:

  • Rooms: Number of rooms

  • Price: Price in dollars

  • Method: S - property sold; SP - property sold prior; PI - property passed in;[…]

  • Type: br - bedroom(s); h - house,cottage,villa, […]

  • SellerG: Real Estate Agent

  • Date: Date sold

  • Distance: Distance from CBD

  • Regionname: General Region (West, North West, North, North east …etc)

  • Propertycount: Number of properties that exist in the suburb.

  • Bedroom2 : Scraped # of Bedrooms (from different source)

  • Bathroom: Number of Bathrooms

  • Car: Number of carspots

  • Landsize: Land Size

  • BuildingArea: Building Size

  • CouncilArea: Governing council for the area

Declarando variables categóricas como tal

siempre que tengamos variables categóricas lo primero es informar a R de que lo son para que las trate como factores. Si no lo hacemos, cuando le pidamos que nos muestre un summary de la variable las interpretará como numéricas y nos mostrará información sin sentido.

Todas estas lo son:

melb_data$Suburb <- factor(melb_data$Suburb)
melb_data$Type <- factor(melb_data$Type)
melb_data$SellerG <- factor(melb_data$SellerG)
melb_data$CouncilArea <- factor(melb_data$CouncilArea)
melb_data$Regionname<-factor(melb_data$Regionname)
melb_data$Method <- factor(melb_data$Method)

#el número de parkings también es categórica (ordinal), sólo toma valores 0,1,2
melb_data$Car <- factor(melb_data$Car)
summary(melb_data)
##             Suburb        Address              Rooms        Type    
##  Reservoir     :  359   Length:13580       Min.   : 1.000   h:9449  
##  Richmond      :  260   Class :character   1st Qu.: 2.000   t:1114  
##  Bentleigh East:  249   Mode  :character   Median : 3.000   u:3017  
##  Preston       :  239                      Mean   : 2.938           
##  Brunswick     :  222                      3rd Qu.: 3.000           
##  Essendon      :  220                      Max.   :10.000           
##  (Other)       :12031                                               
##      Price         Method             SellerG         Date          
##  Min.   :  85000   PI:1564   Nelson       :1565   Length:13580      
##  1st Qu.: 650000   S :9022   Jellis       :1316   Class :character  
##  Median : 903000   SA:  92   hockingstuart:1167   Mode  :character  
##  Mean   :1075684   SP:1703   Barry        :1011                     
##  3rd Qu.:1330000   VB:1199   Ray          : 701                     
##  Max.   :9000000             Marshall     : 659                     
##                              (Other)      :7161                     
##     Distance        Postcode       Bedroom2         Bathroom          Car      
##  Min.   : 0.00   Min.   :3000   Min.   : 0.000   Min.   :0.000   2      :5591  
##  1st Qu.: 6.10   1st Qu.:3044   1st Qu.: 2.000   1st Qu.:1.000   1      :5509  
##  Median : 9.20   Median :3084   Median : 3.000   Median :1.000   0      :1026  
##  Mean   :10.14   Mean   :3105   Mean   : 2.915   Mean   :1.534   3      : 748  
##  3rd Qu.:13.00   3rd Qu.:3148   3rd Qu.: 3.000   3rd Qu.:2.000   4      : 506  
##  Max.   :48.10   Max.   :3977   Max.   :20.000   Max.   :8.000   (Other): 138  
##                                                                  NA's   :  62  
##     Landsize         BuildingArea     YearBuilt           CouncilArea  
##  Min.   :     0.0   Min.   :    0   Min.   :1196                :1369  
##  1st Qu.:   177.0   1st Qu.:   93   1st Qu.:1940   Moreland     :1163  
##  Median :   440.0   Median :  126   Median :1970   Boroondara   :1160  
##  Mean   :   558.4   Mean   :  152   Mean   :1965   Moonee Valley: 997  
##  3rd Qu.:   651.0   3rd Qu.:  174   3rd Qu.:1999   Darebin      : 934  
##  Max.   :433014.0   Max.   :44515   Max.   :2018   Glen Eira    : 848  
##                     NA's   :6450    NA's   :5375   (Other)      :7109  
##    Lattitude        Longitude                          Regionname  
##  Min.   :-38.18   Min.   :144.4   Southern Metropolitan     :4695  
##  1st Qu.:-37.86   1st Qu.:144.9   Northern Metropolitan     :3890  
##  Median :-37.80   Median :145.0   Western Metropolitan      :2948  
##  Mean   :-37.81   Mean   :145.0   Eastern Metropolitan      :1471  
##  3rd Qu.:-37.76   3rd Qu.:145.1   South-Eastern Metropolitan: 450  
##  Max.   :-37.41   Max.   :145.5   Eastern Victoria          :  53  
##                                   (Other)                   :  73  
##  Propertycount  
##  Min.   :  249  
##  1st Qu.: 4380  
##  Median : 6555  
##  Mean   : 7454  
##  3rd Qu.:10331  
##  Max.   :21650  
## 

La variable DATE

El primer inconveniente que nos encontramos es la variable DATE. Con el formato fecha, no podemos usar ningún modelo de machine learning sobre ella, tenemos que transformarla para que nuestros algorítmos puedan digerirla.

DATE contiene la fecha (dd/mm/yyyy) en que la propiedad fue vendida. Vamos a extraer el año y el mes, ya que el día es una información demasiado volatil y no nos servirá de mucho.

Esta es una forma (no mala) de simplificar el modelo, agrupando por categorías, haciéndolo manejable a la par que conservamos gran parte de la información:

melb_data<-cbind(melb_data,year_sold=0,month_sold=0)

year_sold<-format(as.Date(melb_data$Date, format="%d/%m/%Y"),"%Y")
melb_data$year_sold<-as.numeric(year_sold)

month_sold<-format(as.Date(melb_data$Date, format="%d/%m/%Y"),"%m")
melb_data$month_sold<-as.numeric(month_sold)

# la variable original la podemos dejar o quitarla. Así el modelo será más simple:
melb_data$Date<-NULL

Un vistazo a la variable Price

La mayoría de modelos predictivos estarán interesados en deducir el precio de venta de una propiedad. Por lo tanto, aunque no sea nuestro objetivo, vamos a aprovechar para observar cómo es la variable price:

summary(melb_data$Price)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   85000  650000  903000 1075684 1330000 9000000
boxplot(melb_data$Price)

Vemos que la distribución de precios es muy sesgada a la derecha. De hecho es una lognormal, como puede verse a continuación:

hist(log(melb_data$Price),main="log(price)",xlab="",ylab="",col="tomato")

Tal y como decíamos, puede comprobarse en el gráfico la forma perfectamente normal de log(price).

UTILIDADES

Usaremos estas dos funciones más adelante.

Función para encontrar outliers con Q1-Q3

La función examina el contenido del vector x y devuelve todos los índices en los que hay outliers. Consideraremos outlier todo lo que esté fuera del rango [Q1 - 1.5 IQR, Q3 + 1.5 IQR].

# devuelve lista de índices con los outliers
find_outliers<-function(x)
{
        x<-as.double(x)
        Q1<-quantile(x,0.25,na.rm = TRUE)[[1]]
        Q3<-quantile(x,0.75,na.rm = TRUE)[[1]]
        IQR<-Q3-Q1
        
        outlier<-c(which(x<Q1-1.5*IQR),which(x>Q3+1.5*IQR))
        return(outlier) 
}

La función para calcular correlación extendida

Más adelante necesitaremos calcular la correlación por pares de un data frame, tanto para variables categóricas como numéricas. Esta función hace precisamente eso (nosotros la usaremos sólo para las categóricas).

Es necesario ponerla aquí para poder llamarla luego y que el código ejecute limpiamente, aunque su lugar ideal sería un apéndice al final de este texto.

require(tidyverse)
require(rcompanion)


# Calculate a pairwise association between all variables in a data-frame. In particular nominal vs nominal with Chi-square, numeric vs numeric with Pearson correlation, and nominal vs numeric with ANOVA.
# Adopted from https://stackoverflow.com/a/52557631/590437
mixed_assoc = function(df, cor_method="spearman", adjust_cramersv_bias=TRUE){
    df_comb = expand.grid(names(df), names(df),  stringsAsFactors = F) %>% set_names("X1", "X2")

    is_nominal = function(x) class(x) %in% c("factor", "character")
    # https://community.rstudio.com/t/why-is-purr-is-numeric-deprecated/3559
    # https://github.com/r-lib/rlang/issues/781
    is_numeric <- function(x) { is.integer(x) || is_double(x)}

    f = function(xName,yName) {
        x =  pull(df, xName)
        y =  pull(df, yName)

        result = if(is_nominal(x) && is_nominal(y)){
            # use bias corrected cramersV as described in https://rdrr.io/cran/rcompanion/man/cramerV.html
            cv = cramerV(as.character(x), as.character(y), bias.correct = adjust_cramersv_bias)
            data.frame(xName, yName, assoc=cv, type="cramersV")

        }else if(is_numeric(x) && is_numeric(y)){
            correlation = cor(x, y, method=cor_method, use="complete.obs")
            data.frame(xName, yName, assoc=correlation, type="correlation")

        }else if(is_numeric(x) && is_nominal(y)){
            # from https://stats.stackexchange.com/questions/119835/correlation-between-a-nominal-iv-and-a-continuous-dv-variable/124618#124618
            r_squared = summary(lm(x ~ y))$r.squared
            data.frame(xName, yName, assoc=sqrt(r_squared), type="anova")

        }else if(is_nominal(x) && is_numeric(y)){
            r_squared = summary(lm(y ~x))$r.squared
            data.frame(xName, yName, assoc=sqrt(r_squared), type="anova")

        }else {
            warning(paste("unmatched column type combination: ", class(x), class(y)))
        }

        # finally add complete obs number and ratio to table
        result %>% mutate(complete_obs_pairs=sum(!is.na(x) & !is.na(y)), complete_obs_ratio=complete_obs_pairs/length(x)) %>% rename(x=xName, y=yName)
    }

    # apply function to each variable combination
    map2_df(df_comb$X1, df_comb$X2, f)
}

Reducción del número de columnas basándonos en redundancias

Cuando tenemos una base de datos nueva, una de las primeras cosas que hay que mirar es las posibles redundancias entre los datos. No queremos columnas (variables) con información repetida, ya que esto complica nuestro análisis innecesariamente.

Variables categóricas

Veremos a continuación que hay variables categóricas dependientes entre sí. Primero, vamos a separar el data frame original en dos, uno que contenga las variables categóricas y otro que contenga las numéricas.

melb_data_categoric <- melb_data[, sapply(melb_data, is.factor)]
melb_data_numeric<-  melb_data[, sapply(melb_data, is.numeric)]

Calculemos la matriz de correlaciones, usando la función que hemos definido arriba, para las variables categóricas:

CORREL_CATEGORICAS<-mixed_assoc(melb_data_categoric)
CORREL_CATEGORICAS<-subset(CORREL_CATEGORICAS,CORREL_CATEGORICAS$x!=CORREL_CATEGORICAS$y)

#quitamos la diagonal y nos quedamos con las que tengan correlación alta
CORREL_CATEGORICAS<-subset(CORREL_CATEGORICAS,abs(CORREL_CATEGORICAS$assoc)>0.80)
CORREL_CATEGORICAS
##                         x           y  assoc     type complete_obs_pairs
## Cramer V...6  CouncilArea      Suburb 0.8740 cramersV              13580
## Cramer V...7   Regionname      Suburb 0.9887 cramersV              13580
## Cramer V...36      Suburb CouncilArea 0.8740 cramersV              13580
## Cramer V...43      Suburb  Regionname 0.9887 cramersV              13580
##               complete_obs_ratio
## Cramer V...6                   1
## Cramer V...7                   1
## Cramer V...36                  1
## Cramer V...43                  1

Es evidente, observando la matriz de correlaciones, que suburb es totalmente redundante. Nos quedaremos tan sólo con RegionName (también podemos hacerlo al revés). Podiamos deshacernos también de CouncilArea, pero mejor la dejamos ya que está un poco al límite de lo aceptable.

Para estas cosas, no tenemos una regla escrita. Pero si la tuviéramos seguro que una correlación de 0.9887 se consideraba suficiente. Con 0.8740 ya no es tan seguro.

melb_data$Suburb<-NULL

Variables numéricas

Vamos a calcular la matriz de correlación para las variables numéricas. Buscamos variables redundantes:

# correlación sólo para variables numéricas
A<-cor(melb_data_numeric,use="complete.obs") 
corrplot(A)

De la gráfica, es claro que podemos quitar una de Rooms o Bedroom2, dado que, tal y como se ve arriba, tienen una correlación muy alta:

cor(melb_data$Rooms,melb_data$Bedroom2,use="complete.obs")
## [1] 0.9441903

Por lo tanto, para simplificar el modelo:

melb_data$Bedroom2<-NULL

Variables categóricas que no sirven para nada

Tanto Address como SellerG tienen muchísimos valores diferentes. Las variables categóricas con tanta variedad hay que agruparlas. De otro modo, no aportan nada.

El problema es que no podemos agruparlas. La variable Address, basicamente, contiene información redundante con latitude y longitud. Y para sellerG no tenemos un criterio para agruparlo. No sabemos nada al respecto de los vendedores.

Por lo tanto, sólo nos queda eliminarlas:

melb_data$Address<-NULL
melb_data$SellerG<-NULL

Variables con información incorrecta - el caso de Landsize

La variable Landsize tiene un montón de registros con cero asignado. Este dato, dado que es el área en m2 de la construcción (o algo similar), es evidente que está mal. De hecho, el número de valores muy pequeños en Landsize - una casa construída en menos de 40m2 de superficie seria inconcevible - deforma la variable haciéndola inservible.

Primero veamos un resumen de la variable:

summary(melb_data$Landsize)
##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max. 
##      0.0    177.0    440.0    558.4    651.0 433014.0

Veamos cuántos registros mal informados - con valores demasiado pequeños - tenemos:

indices<-which(melb_data$Landsize<=40)
cat("landsize mal informado:",round(100*length(indices)/nrow(melb_data),0),"%\n")
## landsize mal informado: 14 %

14% es una cantidad moderada, así que la imputaremos a la mediana. Nótese que no tiene correlación alta con ninguna otra variable numérica de la que podamos extraer información sobre landsize con algún modelo de regresión.

No podemos imputar landsize a la media porque es una distribución muy sesgada.

melb_data$Landsize[indices]<-median(melb_data$Landsize,na.rm = TRUE)

Tratamiento de los outliers

Antes de nada, vamos a definir nuestra propia función para encontrar outliers:

¿Qué hacemos con los outliers?

Lo mejor para hacerse una idea es echarle un vistazo a los boxplot de las distintas variables. Estas son todas las variables.

colnames(melb_data) 
##  [1] "Rooms"         "Type"          "Price"         "Method"       
##  [5] "Distance"      "Postcode"      "Bathroom"      "Car"          
##  [9] "Landsize"      "BuildingArea"  "YearBuilt"     "CouncilArea"  
## [13] "Lattitude"     "Longitude"     "Regionname"    "Propertycount"
## [17] "year_sold"     "month_sold"

Tomaremos como ejemplo Distance.

x<-find_outliers(melb_data$Distance)
length(x)
## [1] 411

Tenemos más de 400 outliers en esta variable - y ésta es la situación con las demás. Por supuesto, no deben eliminarse. De hecho, si vemos el histograma nos daremos cuenta de que en realidad lo que tenemos es una distribución con mucha asimetría:

hist(melb_data$Distance,main="distance",xlab="",ylab="",col="green")

No vamos a tocar estas variables.

El super-outlier de landsize

De todos los valores que tiene landsize hay uno que llama poderosamente la atención. Es tan grande respecto de todos los demás que aparenta ser imposible:

boxplot(melb_data$Landsize,main="landsize",col = "green")

Vamos a examinar ese registro en concreto. Admitimos que landsize pueda ser una distribución muy asimétrica a derecha, pero cuesta creer este valor en concreto. Vamos a examinar a mano los distintos valores de ese registro. Queremos saber si puede tener una justificación un valor tan grande:

fila_maximo_landsize<-which(melb_data$Landsize==max(melb_data$Landsize))
head(melb_data[fila_maximo_landsize,])
##       Rooms Type   Price Method Distance Postcode Bathroom Car Landsize
## 11021     3    h 2700000     VB      2.1     3065        3   1   433014
##       BuildingArea YearBuilt CouncilArea Lattitude Longitude
## 11021           NA        NA       Yarra -37.79751  144.9831
##                  Regionname Propertycount year_sold month_sold
## 11021 Northern Metropolitan          5825      2017          8

Vemos que el precio de dicha propiedad está en la media y, aunque sí que es cierto que tiene tres cuartos de baño (no parece tanto), cuenta con una sola plaza de parking. Nada parece justificar semejante “landsize”.

El dato es totalmente anómalo, así que le imputaremos la mediana (también podría borrarse). En todo caso tomamos nota de dónde estaba y del valor, por si fuera importante.

melb_data$Landsize[fila_maximo_landsize]<-median(melb_data$Landsize,na.rm = TRUE)

No es que esto que hemos hecho sea perfecto pero, sin duda, los valores de landsize están mejor de cómo los teníamos antes de estas correcciones. Así ha quedado la variable:

boxplot(melb_data$Landsize, main="landsize tras correcciones varias")

Tratamiento de los datos faltantes

Como puede verse en el ejemplo anterior, borrar los registros donde tengamos un N/A no es una opción viable en muchas ocasiones. Si tuviéramos un 10 o un 15% de datos faltantes todavía sería admisible, pero en este caso fulminaríamos el 50% de la información que tenemos .

Tenemos que tener poderosas razones antes de borrar un registro que contenga N/A.

Localizando los datos faltantes

Una vez hemos simplificado un poco los datos, vamos a proceder a ocuparnos de los datos faltantes.

Lo primero que debemos hacer es localizar todas las celdas que tengan valor N/A. Y hay unas cuantas:

plot_missing(melb_data, missing_only = TRUE)

melb_missing <- melb_data[rowSums(is.na(melb_data)) > 0,]

Las variables BuildingArea y YearBuilt tienen una barbaridad de registros faltantes. Car tan sólo unos pocos, 62, y podremos imputarla sin problemas.

Borrado de las columnas con demasiados datos faltantes

Aquellas columnas que tengan demasiados datos faltantes - en este caso BuildingArea y YearBuilt son imposibles de imputar o de usar para sacar información, ya que la distribución subyacente es totalmente desconocida. Deben eliminarse.

El % a partir del cual se eliminan columnas no está escrito. Normalmente se acepta hasta un 20% o un 30%, pero en ocasiones podría admitirse un porcentaje mayor si la variable que vamos a borrar es de gran importancia en el modelo.

melb_data$YearBuilt<-NULL
melb_data$BuildingArea<-NULL
colnames(melb_data)
##  [1] "Rooms"         "Type"          "Price"         "Method"       
##  [5] "Distance"      "Postcode"      "Bathroom"      "Car"          
##  [9] "Landsize"      "CouncilArea"   "Lattitude"     "Longitude"    
## [13] "Regionname"    "Propertycount" "year_sold"     "month_sold"

Siempre hay que tener mucho cuidado cuando se borran registros o variables.

Columnas con datos faltantes

Después de nuestra criba de columnas inútiles, y de unas cuantas correcciones menores, tenemos que abordar el problema de los N/A, que ahora mismo están presentes sólo en la variable Car.

summary(melb_data$Car)
##    0    1    2    3    4    5    6    7    8    9   10 NA's 
## 1026 5509 5591  748  506   63   54    8    9    1    3   62

Como puee verse, tenemos 62 valores NA en car.

Imputando car

Car es una variable ordinal dado que el número de parkings de una casa tiene sentido ponerlo de menos a más (más parkings, más precio). Tendremos que usar regresión logística ordinal.

summary(melb_data$Car)
##    0    1    2    3    4    5    6    7    8    9   10 NA's 
## 1026 5509 5591  748  506   63   54    8    9    1    3   62

Imputando los datos faltantes con MICE

Para la imputación, podemos crear nosotros el modelo a mano o podemos usar una librería ya preparada. En este caso, aprovecharemos par intrudocir la librería MICE, de la que hablaremos más detalladamente en el futuro.

Imputar con mice es tan simple como esto:

library(mice)
ncol(melb_data)
## [1] 16
my_imp<-mice(melb_data,m=5,maxit=2,)   #maxit=2 es suficiente para nuestra prueba
## 
##  iter imp variable
##   1   1  Car
##   1   2  Car
##   1   3  Car
##   1   4  Car
##   1   5  Car
##   2   1  Car
##   2   2  Car
##   2   3  Car
##   2   4  Car
##   2   5  Car

Y ahora veamos cómo ha quedado la imputación:

my_imp$imp$Car
##        1 2 3 4  5
## 12222  0 1 1 2  2
## 12248  1 4 1 1  1
## 12260  1 4 1 1  1
## 12321  3 3 0 2  6
## 12363  3 2 2 8  2
## 12384  2 0 4 2  2
## 12398  2 1 1 1  1
## 12427  2 2 2 4  0
## 12463  1 1 1 2  1
## 12464  1 2 2 1  1
## 12465  1 1 3 2  1
## 12475  0 0 1 0  1
## 12490  2 1 1 2  2
## 12491  1 4 0 1  1
## 12498  0 2 2 1  0
## 12532  1 3 2 1  2
## 12607  1 2 1 2  0
## 12608  1 2 1 0  1
## 12609  0 0 0 0  0
## 12611  2 2 2 1  0
## 12612 10 1 5 1  0
## 12614  2 2 2 2  2
## 12634  2 0 1 1  0
## 12639  3 0 2 1  1
## 12724  2 2 2 2  2
## 12786  1 2 2 2  0
## 12806  1 1 2 1  2
## 12807  1 3 1 0  1
## 12816  1 2 1 0  1
## 12819  2 1 1 1  1
## 12820  0 2 1 4  0
## 12831 10 1 2 2  3
## 12842  2 4 2 2  2
## 12845  1 0 2 1  1
## 12885  2 1 3 2  0
## 12889  1 2 2 2  2
## 12930  1 1 1 1  0
## 12938  1 1 0 2  1
## 12976  1 4 2 0  0
## 12996  3 0 2 2  1
## 13049  2 1 1 0 10
## 13092  0 2 0 5  1
## 13098  1 2 3 1  2
## 13176  0 1 0 1  1
## 13178  2 1 2 2  2
## 13225  4 2 2 1  0
## 13226  0 1 1 1  1
## 13250  0 1 1 0  1
## 13254  2 1 1 1  0
## 13267  2 2 3 2  4
## 13290  1 1 1 1  2
## 13295  1 3 2 2  1
## 13298  2 1 1 1  2
## 13346  1 0 0 0  2
## 13373  1 0 1 1  1
## 13390  2 1 1 2  2
## 13461  0 1 3 1  5
## 13497  2 2 2 2  2
## 13509  1 1 0 2  6
## 13523  2 0 1 0  1
## 13525  0 1 2 2  1
## 13551  2 2 2 1  1

Una vez imputado todo, procedemos a cargar los valores:

melb_data<-complete(my_imp,5)

Con esto, ya tendríamos totalmente corregidos los datos y podemos afrontar cualquier modelo predictivo con unas garantías mínimas.

summary(melb_data)
##      Rooms        Type         Price         Method       Distance    
##  Min.   : 1.000   h:9449   Min.   :  85000   PI:1564   Min.   : 0.00  
##  1st Qu.: 2.000   t:1114   1st Qu.: 650000   S :9022   1st Qu.: 6.10  
##  Median : 3.000   u:3017   Median : 903000   SA:  92   Median : 9.20  
##  Mean   : 2.938            Mean   :1075684   SP:1703   Mean   :10.14  
##  3rd Qu.: 3.000            3rd Qu.:1330000   VB:1199   3rd Qu.:13.00  
##  Max.   :10.000            Max.   :9000000             Max.   :48.10  
##                                                                       
##     Postcode       Bathroom          Car          Landsize    
##  Min.   :3000   Min.   :0.000   2      :5608   Min.   :   41  
##  1st Qu.:3044   1st Qu.:1.000   1      :5534   1st Qu.:  306  
##  Median :3084   Median :1.000   0      :1040   Median :  440  
##  Mean   :3105   Mean   :1.534   3      : 749   Mean   :  590  
##  3rd Qu.:3148   3rd Qu.:2.000   4      : 507   3rd Qu.:  651  
##  Max.   :3977   Max.   :8.000   5      :  64   Max.   :76000  
##                                 (Other):  78                  
##         CouncilArea     Lattitude        Longitude    
##               :1369   Min.   :-38.18   Min.   :144.4  
##  Moreland     :1163   1st Qu.:-37.86   1st Qu.:144.9  
##  Boroondara   :1160   Median :-37.80   Median :145.0  
##  Moonee Valley: 997   Mean   :-37.81   Mean   :145.0  
##  Darebin      : 934   3rd Qu.:-37.76   3rd Qu.:145.1  
##  Glen Eira    : 848   Max.   :-37.41   Max.   :145.5  
##  (Other)      :7109                                   
##                       Regionname   Propertycount     year_sold   
##  Southern Metropolitan     :4695   Min.   :  249   Min.   :2016  
##  Northern Metropolitan     :3890   1st Qu.: 4380   1st Qu.:2016  
##  Western Metropolitan      :2948   Median : 6555   Median :2017  
##  Eastern Metropolitan      :1471   Mean   : 7454   Mean   :2017  
##  South-Eastern Metropolitan: 450   3rd Qu.:10331   3rd Qu.:2017  
##  Eastern Victoria          :  53   Max.   :21650   Max.   :2017  
##  (Other)                   :  73                                 
##    month_sold    
##  Min.   : 1.000  
##  1st Qu.: 5.000  
##  Median : 7.000  
##  Mean   : 7.052  
##  3rd Qu.: 9.000  
##  Max.   :12.000  
## 

NOTAS FINALES

En este ejemplo hemos borrado unas cuantas columnas por distintos motivos. Hay que tener mucho cuidado con esto. En nuestro caso, simplemente se trataba de dejar la base de datos preparada para un modelo predictivo, pero eliminar columnas es algo que siempre podremos hacer y, por lo tanto no es buena idea borrar datos sin saber qué haremos en el futuro con ellos.

También hemos dejado en la base de datos un montón de outliers que no hemos tocado. Sin tener información adicional no es tan fácil discernir cuál de ellos es un valor legítimo y cuál no. En nuestro caso, hemos visto que las distribuciones eran bastante sesgadas, por lo que tiene sentido que aparezca valores más extremos y, sin saber más, hemos optado por no tocar nada. Siempre es aconsejable el consejo de un experto en la materia antes de asumir que algo es un outlier.

Finalmente, hemos usado el paquete R mice para imputar la variable car. Esto no era en absoluto necesario, pero hemos querido aprovechar para introducir el uso de este paquete, del que hablaremos más en el futuro, con un ejemplo sencillo.

Comentarios

Entradas populares de este blog

UNA BREVE EXPLICACIÓN DE LA LEY DE LAS ESPERANZAS ITERADAS

Que se mueran los pobres - Un modelo de regresión logística para predecir la supervivencia en el Titanic