Trabajo Practico 1

1: Carga de datos

1.A: Lectura y analisis de estructura del archivo

Procedemos a cargar el archivo de entrada y verificar la cantidad de variables y casos. Para ello tenemos que descomprimir el archivo antes de cargalo, y tabmien especificar el encoding correcto (UTF-8) para la carga correcta de caracteres especiales.

Debido a que la carga de los ids cuesta mucho mas memoria que todo el resto y no proporciona informacion de utilidad para el analisis se descarta de entrada.

# carga libreria tidyverse
suppressPackageStartupMessages(library(tidyverse))

#carga de datos
ar_properties <- read.csv(unz('ar_properties.zip', 'ar_properties.csv'), row.names = NULL, stringsAsFactors = TRUE, encoding = "UTF-8") %>% select( -c('id') )

#verificaion de casos y variables
head(ar_properties)

nrow(ar_properties)

[1] 388891

1.B: Preprocesado y filtrado

Cargamos la libreria tidyverse para realizar los filtrados, y procedemos a extraer los casos de interes:

# Filtrado de los casos deseados
ar_properties <- ar_properties %>% filter(
    l1             == 'Argentina' &
    l2             == 'Capital Federal' &
    currency       == 'USD' &
    operation_type == 'Venta' &
    property_type  %in% c('Casa', 'Departamento', 'PH'))

1.C: Preprocesado y filtrado

Una vez filtrados los casos, seleccionamos las variables deseadas y reajustamos los niveles de las variables de categoria (descartamos las categorias que ya no figuran en nuestros datos). Verificamos nuevamente la cantidad de registros y variables que quedan

# Seleccion de las variables de interes
ar_properties <- ar_properties %>% select(
  l3,
  rooms,
  bedrooms,
  bathrooms,
  surface_total,
  surface_covered,
  price,
  property_type )

# eliminacion de categorias en desuso
ar_properties <- droplevels(ar_properties)

#verificacion de casos y variables
nrow(ar_properties)

[1] 61905

ncol(ar_properties)

[1] 8

2: Analisis exploratorio (I)

2.A: Analisis de valores distintos y faltantes

calculamos la cantidad de valores distintos y faltantes, y los analizamos:

# calculamos la cantidad de valores distintos
uniques.values <- apply(ar_properties, 2, n_distinct)
as.data.frame(uniques.values)

# calculamos la cantidad de valores faltantes
missing.values <- colSums(is.na(ar_properties))
as.data.frame(missing.values)

Curiosamente, notamos que para el campo l3 que corresponde al barrio porteño, hay 57 valores distintos (58 si tomamos en cuenta NA). Esto se contradice con los 48 barrios porteños existentes, por lo que analizamos los barrios inesperados:

barrios.portenos.officiales <- c("Agronomía", "Almagro", "Balvanera", "Barracas", "Belgrano", "Boedo", "Caballito", "Chacarita", "Coghlan", "Colegiales", "Constitución", "Flores", "Floresta", "Boca", "Paternal", "Liniers", "Mataderos", "Monserrat", "Monte Castro", "Pompeya", "Nuñez", "Palermo", "Parque Avellaneda", "Parque Chacabuco", "Parque Chas", "Parque Patricios", "Puerto Madero", "Recoleta", "Retiro", "Saavedra", "San Cristobal", "San Nicolás", "San Telmo", "Velez Sarsfield", "Versalles", "Villa Crespo", "Villa del Parque", "Villa Devoto", "Villa General Mitre", "Villa Lugano", "Villa Luro", "Villa Ortuzar", "Villa Pueyrredón", "Villa Real", "Villa Riachuelo", "Villa Santa Rita", "Villa Soldati", "Villa Urquiza")
sort(setdiff(unique(ar_properties$l3), barrios.portenos.officiales))

[1] "Abasto"               "Barrio Norte"         "Catalinas"            "Centro / Microcentro"
[5] "Congreso"             "Las Cañitas"          "Once"                 "Parque Centenario"   
[9] "Tribunales"          

Observamos como los barrios ‘extra’ corresponden a denominaciones no-oficiales de ciertas regions de la CABA.

2.B: Matriz de correlacion

Para el calculo de la matriz de correlacion, necesitamos poder extraer las variables del tipo numerico. Dado que esta funcion sera utilizada en pasos posteriores, será de utilidad crear una funcion que permita automatizar esta funcion; y ademas agregamos la opcion para seleccionar columnas especificas en lugar de todas las numericas. Luego procedemos a calcular la correlacion de las variables numericas.

# funcion para obtener las columnas especificadas, o las numericas en caso de no especificarse
get.numeric.columns <- function(df, columns.to.get=NULL) {
    if (is.null(columns.to.get))
        columns.to.get <- colnames(df)[unlist(lapply(df, is.numeric))]

    return(columns.to.get)
}

# calculamos la correlacion de las variables numericas, descartando los casos incompletos
cor(ar_properties[,get.numeric.columns(ar_properties)], use="complete.obs", method="pearson")

                     rooms   bedrooms  bathrooms surface_total surface_covered      price
rooms           1.00000000 0.92138719 0.61335026    0.06828238      0.07468335 0.48748747
bedrooms        0.92138719 1.00000000 0.61578024    0.06746895      0.07206826 0.43221753
bathrooms       0.61335026 0.61578024 1.00000000    0.06234262      0.06777010 0.59904254
surface_total   0.06828238 0.06746895 0.06234262    1.00000000      0.69656225 0.05095265
surface_covered 0.07468335 0.07206826 0.06777010    0.69656225      1.00000000 0.06257960
price           0.48748747 0.43221753 0.59904254    0.05095265      0.06257960 1.00000000

3: Preparacion de datos

3.A: Eliminacion de variables altamente correlacionadas

Dado que las variables rooms y bedrooms presentan alta correlacion, procedemos a eliminar una de ellas, la de mayor cantaidad de faltantes (bedrooms)

# eliminacion de la variable `bedrooms`
ar_properties <- ar_properties %>% select( -c('bedrooms') )

3.B: Filtrado de los casos incompletos

Filtramos los casos incompletos y validamos la cantidad de variables y casos

# eliminacion de la variable `bedrooms`
ar_properties <- ar_properties[complete.cases(ar_properties),]

#Validamos la cantidad de casos y variables
nrow(ar_properties)

[1] 51210

ncol(ar_properties)

[1] 7

4: Analisis exploratorio (II)

4.A: Obtencion de estadisticas generales

Obtenemos estadisticas de la variable precio. Generamos funciones que reutilizaremos mas tarde.

# Funcion para generar estadisticas con la media incluida
summary.with.mean <- function(data)
{
    summary <- summary(data)
    summary['mean'] = mean(data)
    return( summary )
}

# Caluclamos las estadisticas
summary.with.mean(ar_properties$price)

   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    mean 
   6000  119000  170000  251577  270000 6000000  251577 

# Y realizamos un histograma
ggplot(ar_properties, aes(x = price)) + geom_histogram()

Notamos que hay un gran rango de precios, y que la mayoria de los valores se concentran en los valores mas bajos. Por ese motivo, repetimos el histograma pero utiliando una escala logaritmica para poder observar mejor la distribucion.

ggplot(ar_properties, aes(x = price)) + geom_histogram(alpha=0.5, position="identity", aes(y = ..density..)) + scale_x_log10() + geom_density(alpha=0.5)

4.B: Obtencion de estadisticas por tipo de propiedad

Realizamos el mismo analisis anterior, pero segmentado por tipo de propiedad.

# Calculamos las estadisticas, discriminando por tipo de propiedad
tapply(ar_properties$price, ar_properties$property_type, summary.with.mean)

$Casa
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    mean 
  20000  235000  335000  434189  490000 5000000  434189 

$Departamento
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    mean 
   6000  115000  164000  246856  260000 6000000  246856 

$PH
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    mean 
  32000  137000  190000  218747  270000 1500000  218747 

# Realizamos el histograma, discriminando por tipo de propiedad
ggplot(ar_properties, aes(x=price, fill=property_type)) +
    geom_histogram(alpha=0.5, position="identity", aes(y = ..density..)) +
    geom_density(alpha=0.5) +
    scale_x_log10()

4.C: Boxplots

Graficamos los boxplots del precio por cada tipo de propiedad. Siguiendo la misma logica, escalamos el eje de precios logaritmicamente.

# Graficamos los boxplots
ggplot(ar_properties, aes(x=property_type, y=price, fill=property_type)) +
    geom_boxplot() +
    scale_y_log10()

Obervamos como los precios tienen una tendencia del tipo Casa > PH > depto (ya obervada tambien previamente en los histogramas)

4.D: Correlograma

suppressPackageStartupMessages(library(GGally))

# Graficamos el correlograma
ggcorr(ar_properties[,get.numeric.columns(ar_properties)], method = c("everything", "pearson")) 

5: Outliers

Analizamos como se componen los precios en los extremos superiores e inferiores de la distribucion

# Cortes de precios para los cuantiles en los extremos
quant.inf <- quantile(ar_properties$price, probs=seq(0, 0.005, length.out = 6) )
quant.sup <- quantile(ar_properties$price, probs=seq(.995, 1, length.out = 6) )

quant.inf

   0%  0.1%  0.2%  0.3%  0.4%  0.5% 
 6000 35000 45000 49900 53000 55000 

quant.sup

  99.5%   99.6%   99.7%   99.8%   99.9%    100% 
1979550 2300000 2500000 2879100 3500000 6000000 

Analizando los resultados de zonaProp, la propiedad mas barata en venta en capital federal esta por sobre 50.000 dolares. Combinando este dato con la informacion de los cuartiles obtenida previamente, podemos asegurar que descartar los valores inferiores a U$S 50.000 seria descartar menos del 0.5% de las muestras. En cuanto a valores superiores, se prosiguió a descartar el mismo porcentage de corte que para las muestras inferiores.

ar_properties.sin.outliers = ar_properties %>% filter(
    price > quant.inf['0.5%'] & price < quant.sup['99.5%'] )

Como caso de interes, tambien filtramos los outliers utilizando la tecnica de isolation tree, utilizando las variables de interes y filtrando una cantidad demejante de muestras, pero utilizando como criterio su puntaje de anomalia

# Cargamos la libreria necesaria
suppressPackageStartupMessages(library("solitude"))

# iniciamos el algoritmo utilizando solamente las variables que vamos a correlacionar
iso <- solitude::isolationForest$new()
iso$fit(ar_properties %>% select( c("rooms", "price", "surface_total") ) )

Building Isolation Forest ... done
Computing depth of terminal nodes ... 

done

ar_properties.sin.outliers.isolationtree <- ar_properties

# Usamos los resultados para filtras las muestras mas anomalas, y luego descartamos esta columna
ar_properties.sin.outliers.isolationtree$anomalyScore = iso$scores$anomaly_score
ar_properties.sin.outliers.isolationtree <- ar_properties.sin.outliers.isolationtree %>% filter( anomalyScore < 0.6 ) %>% select( -c("anomalyScore") )

6: Analisis exploratorios (III)

# Calculamos las estadisticas, discriminando por tipo de propiedad
tapply(ar_properties.sin.outliers$price, ar_properties.sin.outliers$property_type, summary.with.mean)

$Casa
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    mean 
  62000  235000  330000  397974  485000 1950000  397974 

$Departamento
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    mean 
  55400  115000  164000  234475  260000 1970000  234475 

$PH
   Min. 1st Qu.  Median    Mean 3rd Qu.    Max.    mean 
  58000  138000  195000  219302  272750 1500000  219302 

# Realizamos el histograma, discriminando por tipo de propiedad
ggplot(ar_properties.sin.outliers, aes(x=price, fill=property_type)) +
    geom_histogram(alpha=0.5, position="identity", aes(y = ..density..)) +
    geom_density(alpha=0.5) +
    scale_x_log10()

# Graficamos los boxplots
ggplot(ar_properties.sin.outliers, aes(x=property_type, y=price, fill=property_type)) +
    geom_boxplot() +
    scale_y_log10()

# Graficamos el correlograma
ggcorr(ar_properties.sin.outliers[,get.numeric.columns(ar_properties.sin.outliers)], method = c("everything", "pearson"))

Se puede observar (con especial detalle en los boxplots) como gran cantidad de los outliers desaparecieron. Cabe destacar que esta es una tecnica un poco ingenua para la eliminacion de los mismos, ya que no considera la posibilidad de outliers debido a la relacion de las variables entre sí, y solamente descarta los extremos en cada dimension (en este caso, la dimension precio)

7: Modelo lineal

7.A: Generacion del modelo lineal

#calculamos los modelos lineales indicados
linear.model.surface <- lm(price ~ surface_total, data = ar_properties.sin.outliers)
linear.model.rooms <- lm(price ~ rooms, data = ar_properties.sin.outliers)

7.B: Descripcion de los modelos lineales

superficie vs precio

summary(linear.model.surface)

Call:
lm(formula = price ~ surface_total, data = ar_properties.sin.outliers)

Residuals:
     Min       1Q   Median       3Q      Max 
-2103818  -119487   -69187    30879  1727990 

Coefficients:
               Estimate Std. Error t value Pr(>|t|)    
(Intercept)   237800.77     972.82  244.45   <2e-16 ***
surface_total     16.51       1.18   13.99   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 217300 on 50630 degrees of freedom
Multiple R-squared:  0.003851,  Adjusted R-squared:  0.003832 
F-statistic: 195.8 on 1 and 50630 DF,  p-value: < 2.2e-16

ambientes vs precio

summary(linear.model.rooms)

Call:
lm(formula = price ~ rooms, data = ar_properties.sin.outliers)

Residuals:
     Min       1Q   Median       3Q      Max 
-2239135   -84697   -26697    32289  1692289 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  -9332.7     1809.4  -5.158 2.51e-07 ***
rooms        89014.6      580.8 153.272  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 179900 on 50630 degrees of freedom
Multiple R-squared:  0.3169,    Adjusted R-squared:  0.3169 
F-statistic: 2.349e+04 on 1 and 50630 DF,  p-value: < 2.2e-16

7.B (bis): Descripcion de los modelos lineales - isolation tree

Verificamos que la cantidad de casos en ambos filtrados sea semejante

nrow(ar_properties.sin.outliers)

[1] 50632

nrow(ar_properties.sin.outliers.isolationtree)

[1] 50981

Analizamos los modelos:

superficie vs precio

summary(lm(price ~ surface_total, data = ar_properties.sin.outliers.isolationtree))

Call:
lm(formula = price ~ surface_total, data = ar_properties.sin.outliers.isolationtree)

Residuals:
     Min       1Q   Median       3Q      Max 
-2303976   -45505   -17164    19403  2630776 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)   33868.30    1252.81   27.03   <2e-16 ***
surface_total  2339.93      10.68  219.06   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 180900 on 50979 degrees of freedom
Multiple R-squared:  0.4849,    Adjusted R-squared:  0.4849 
F-statistic: 4.799e+04 on 1 and 50979 DF,  p-value: < 2.2e-16

ambientes vs precio

summary(lm(price ~ rooms, data = ar_properties.sin.outliers.isolationtree))

Call:
lm(formula = price ~ rooms, data = ar_properties.sin.outliers.isolationtree)

Residuals:
    Min      1Q  Median      3Q     Max 
-909837  -91525  -28961   38139 3030147 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) -38467.3     2166.3  -17.76   <2e-16 ***
rooms       101664.1      700.5  145.13   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 212000 on 50979 degrees of freedom
Multiple R-squared:  0.2924,    Adjusted R-squared:  0.2924 
F-statistic: 2.106e+04 on 1 and 50979 DF,  p-value: < 2.2e-16

Y realizamos el correlograma para esta version de filtrado

ggcorr(ar_properties.sin.outliers.isolationtree[,get.numeric.columns(ar_properties.sin.outliers.isolationtree)], method = c("everything", "pearson"))

7.C: Analisis / conclusiones

filtrado por extremos

Se puede observar en los valores de R-squared (mayor en el caso de ambiente-precio) que la correlacion es mas fuete para esta variable (ya se podia observar este efecto en el correlograma). Tambien observamos como los residuos son menores para este modelo lineal. Dado estos motivos, el modelo de precio-ambiente parece ser un mejor predictor de precio que el de sueprficie-ambiente.

filtrado por isolation tree

Cuando el filtrado se realiza utilizando otra tecnica, se observa como los resultados anteriores se invierten: el modelo de area-precio se vuelve un mejor predictor en este caso; pero ademas se potencia, dado que este modelo (filtrado por isolation tree, y modelando area vs precio) es la que expresa una mayor correlacion que todos los otros casos; producto de haber filtrado no solo los valores extremos, sino tambien usando tecnicas mas avanzadas para detectar outliers de valores intermedios, que habrian sesgado el modelo a uno menor performante que el de ambientes vs precio (esto se observa claramente en los dispersogramas a continuacion)

ggplot(ar_properties.sin.outliers, aes(x=rooms, y=price)) + geom_point(color="red") + geom_smooth(method='lm') + ggtitle("Precio vs Ambientes - filtrado ingenuo")

ggplot(ar_properties.sin.outliers, aes(x=surface_total, y=price)) + geom_point(color="red") + geom_smooth(method='lm') + ggtitle("Precio vs Area - filtrado ingenuo")

ggplot(ar_properties.sin.outliers.isolationtree, aes(x=rooms, y=price)) + geom_point(color="green") + geom_smooth(method='lm') + ggtitle("Precio vs Ambientes - filtrado por isolation tree")

ggplot(ar_properties.sin.outliers.isolationtree, aes(x=surface_total, y=price)) + geom_point(color="green") + geom_smooth(method='lm') + ggtitle("Precio vs Area - filtrado por isolation tree")