8.2 Práctica Guiada
Para este ejercicio utilizaremos los datos provistos por Properati: https://www.properati.com.ar/data/
(Primero acondicionamos la base original, para quedarnos con una base más fácil de trabajar, y que contiene unicamente los datos interesantes. También eliminamos previamente los outliers)
## # A tibble: 10 x 9
## id created_on l3 rooms bathrooms surface_total surface_covered
## <chr> <date> <chr> <dbl> <dbl> <dbl> <dbl>
## 1 YnEW… 2019-02-07 Boedo 2 1 48 41
## 2 BOzp… 2019-04-01 Saav… 3 2 200 84
## 3 QC9J… 2019-05-31 Pale… 4 2 100 100
## 4 JbTZ… 2019-04-18 Pale… 3 2 105 105
## 5 B4pO… 2019-02-22 Puer… 3 3 165 155
## 6 7X+I… 2019-01-02 Alma… 4 2 165 165
## 7 bSPU… 2019-05-03 Flor… 4 2 158 149
## 8 iA37… 2019-04-25 Flor… 3 1 60 55
## 9 j9UV… 2019-03-24 Flor… 4 1 162 115
## 10 NdU4… 2019-03-23 Alma… 3 1 60 60
## # … with 2 more variables: price <dbl>, property_type <chr>8.2.1 Primer modelo
##
## Call:
## lm(formula = price ~ l3 + rooms + bathrooms + surface_total +
## property_type, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2152396 -72315 -4170 46095 5284708
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.597e+05 1.389e+04 -11.504 < 2e-16 ***
## l3Agronomía 2.023e+04 2.605e+04 0.776 0.437561
## l3Almagro -7.963e+03 1.296e+04 -0.615 0.538880
## l3Balvanera -2.617e+04 1.360e+04 -1.924 0.054392 .
## l3Barracas 5.374e+02 1.592e+04 0.034 0.973061
## l3Barrio Norte 6.417e+04 1.331e+04 4.820 1.44e-06 ***
## l3Belgrano 1.212e+05 1.290e+04 9.396 < 2e-16 ***
## l3Boca -4.937e+04 2.020e+04 -2.444 0.014546 *
## l3Boedo -1.422e+04 1.554e+04 -0.915 0.360008
## l3Caballito 6.372e+03 1.298e+04 0.491 0.623515
## l3Catalinas -2.557e+04 1.059e+05 -0.241 0.809215
## l3Centro / Microcentro -3.332e+04 1.970e+04 -1.692 0.090682 .
## l3Chacarita 2.842e+04 1.609e+04 1.766 0.077375 .
## l3Coghlan 5.894e+04 1.644e+04 3.585 0.000337 ***
## l3Colegiales 3.945e+04 1.456e+04 2.710 0.006734 **
## l3Congreso -2.853e+04 1.610e+04 -1.772 0.076373 .
## l3Constitución -2.955e+04 1.767e+04 -1.672 0.094507 .
## l3Flores -2.403e+04 1.364e+04 -1.762 0.078056 .
## l3Floresta -1.222e+04 1.517e+04 -0.806 0.420517
## l3Las Cañitas 1.194e+05 1.759e+04 6.787 1.16e-11 ***
## l3Liniers -2.030e+04 1.592e+04 -1.275 0.202324
## l3Mataderos -3.333e+04 1.612e+04 -2.067 0.038733 *
## l3Monserrat -9.575e+03 1.570e+04 -0.610 0.541989
## l3Monte Castro 1.873e+04 1.794e+04 1.044 0.296448
## l3Nuñez 9.193e+04 1.373e+04 6.694 2.20e-11 ***
## l3Once -2.204e+04 1.599e+04 -1.379 0.167936
## l3Palermo 1.275e+05 1.272e+04 10.023 < 2e-16 ***
## l3Parque Avellaneda -1.670e+04 2.200e+04 -0.759 0.447906
## l3Parque Centenario -3.829e+04 1.524e+04 -2.513 0.011984 *
## l3Parque Chacabuco -1.405e+03 1.569e+04 -0.090 0.928651
## l3Parque Chas 2.208e+04 2.268e+04 0.974 0.330203
## l3Parque Patricios -1.129e+04 1.769e+04 -0.638 0.523388
## l3Paternal -2.800e+03 1.550e+04 -0.181 0.856664
## l3Pompeya -6.163e+04 2.211e+04 -2.787 0.005321 **
## l3Puerto Madero 5.285e+05 1.457e+04 36.273 < 2e-16 ***
## l3Recoleta 1.293e+05 1.310e+04 9.871 < 2e-16 ***
## l3Retiro 7.510e+04 1.571e+04 4.779 1.76e-06 ***
## l3Saavedra 3.673e+04 1.486e+04 2.472 0.013448 *
## l3San Cristobal -1.199e+04 1.480e+04 -0.810 0.417905
## l3San Nicolás -4.628e+03 1.535e+04 -0.302 0.762992
## l3San Telmo 1.762e+04 1.461e+04 1.207 0.227615
## l3Tribunales -4.234e+04 2.556e+04 -1.656 0.097677 .
## l3Velez Sarsfield 1.644e+03 2.487e+04 0.066 0.947319
## l3Versalles 4.489e+03 1.989e+04 0.226 0.821391
## l3Villa Crespo 1.071e+04 1.304e+04 0.822 0.411190
## l3Villa del Parque 2.439e+04 1.471e+04 1.658 0.097237 .
## l3Villa Devoto 3.089e+04 1.440e+04 2.145 0.031955 *
## l3Villa General Mitre -2.568e+04 2.025e+04 -1.268 0.204694
## l3Villa Lugano -1.002e+05 1.749e+04 -5.728 1.02e-08 ***
## l3Villa Luro 7.194e+03 1.618e+04 0.445 0.656600
## l3Villa Ortuzar 2.825e+04 2.043e+04 1.383 0.166821
## l3Villa Pueyrredón 2.685e+04 1.591e+04 1.688 0.091475 .
## l3Villa Real 1.340e+04 2.593e+04 0.517 0.605265
## l3Villa Riachuelo -5.138e+04 4.990e+04 -1.030 0.303147
## l3Villa Santa Rita 6.239e+03 1.910e+04 0.327 0.743936
## l3Villa Soldati -9.215e+04 3.637e+04 -2.534 0.011282 *
## l3Villa Urquiza 4.076e+04 1.333e+04 3.057 0.002238 **
## rooms 5.200e+04 8.992e+02 57.836 < 2e-16 ***
## bathrooms 1.418e+05 1.462e+03 97.039 < 2e-16 ***
## surface_total 5.800e+00 1.141e+00 5.083 3.72e-07 ***
## property_typeDepartamento 4.770e+03 5.268e+03 0.906 0.365198
## property_typePH -4.786e+04 5.693e+03 -8.407 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 210300 on 52184 degrees of freedom
## Multiple R-squared: 0.4868, Adjusted R-squared: 0.4862
## F-statistic: 811.6 on 61 and 52184 DF, p-value: < 2.2e-16¿ Qué pasó con las variables no numéricas? ¿Son significativos los estimadores? ¿cuales? ¿Cómo se leen los valores de los estimadores?
\[ \Delta y = \beta_1 \Delta x \]
8.2.2 Feature engineering.
Se denomina ingeniería de features a la construcción de nuevas variables a partir de las originales.
Por ejemplo, dado que muchos de los barrios no explican significativamente los cambios en los precios, no esta bueno conservarlos todos. A su vez, no sabemos respecto a qué barrio se compara.
Una solución puede ser agrupar los barrios en tres categorías respecto a su efecto en el precio:
- Alto
- Medio
- Bajo
En particular, podemos notar de esta primera regresión que algunos barrios tienen un efecto significativo en subir el valor de la propiedad, como Belgrano o Recoleta.
Para construir la nueva variable, podemos ver el precio promedio del metro cuadrado por barrio
df_barrios <- df %>%
group_by(l3) %>%
summarise(precio_m2 = mean(price/surface_total))
ggplot(df_barrios,aes(precio_m2)) +
geom_dotplot(method = 'histodot', binwidth = 100)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 871.2 2031.8 2147.3 2345.8 2560.0 6061.1Con este gráfico vemos que que hay muchos barrios con un precio promedio cercano a 2500 dólares el \(m^2\).
Podemos dividr los tres grupos al rededor de los quartiles 1 y 3.
- <2000 bajo
- 2000-2500 medio
- >2500 alto
df_barrios <- df_barrios %>%
mutate(barrio= case_when(precio_m2<=2100 ~ 'bajo',
precio_m2>2100 & precio_m2<2400 ~ 'medio',
#between(precio_m2,2100,2400)~'medio', # Otra forma de escribirlo
precio_m2>=2400 ~ 'alto'))
df_barrios %>%
sample_n(10)## # A tibble: 10 x 3
## l3 precio_m2 barrio
## <chr> <dbl> <chr>
## 1 Liniers 2057. bajo
## 2 Parque Chas 2305. medio
## 3 Belgrano 3421. alto
## 4 Versalles 1881. bajo
## 5 Abasto 2635. alto
## 6 Caballito 2687. alto
## 7 Nuñez 3264. alto
## 8 Villa Pueyrredón 2292. medio
## 9 Monte Castro 2121. medio
## 10 Villa Devoto 2345. medioCon esta nueva variable podemos modificar la tabla original.
y volvemos a calcular el modelo
lm_fit <- lm(price~ barrio+ rooms + bathrooms + surface_total + property_type,data = df)
summary(lm_fit)##
## Call:
## lm(formula = price ~ barrio + rooms + bathrooms + surface_total +
## property_type, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2143881 -67303 -11936 39417 5319777
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.196e+05 6.403e+03 -18.677 < 2e-16 ***
## barriobajo -9.165e+04 3.304e+03 -27.742 < 2e-16 ***
## barriomedio -7.399e+04 2.640e+03 -28.029 < 2e-16 ***
## rooms 4.709e+04 9.366e+02 50.274 < 2e-16 ***
## bathrooms 1.633e+05 1.508e+03 108.275 < 2e-16 ***
## surface_total 5.683e+00 1.208e+00 4.706 2.54e-06 ***
## property_typeDepartamento 2.395e+04 5.410e+03 4.427 9.58e-06 ***
## property_typePH -4.085e+04 5.969e+03 -6.844 7.80e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 222700 on 52238 degrees of freedom
## Multiple R-squared: 0.4243, Adjusted R-squared: 0.4242
## F-statistic: 5499 on 7 and 52238 DF, p-value: < 2.2e-16Si queremos que compare contra ‘barrio medio’ podemos convertir la variable en factor y explicitar los niveles
df <- df %>%
mutate(barrio = factor(barrio, levels = c('medio', 'alto','bajo')))
lm_fit <- lm(price~ barrio+ rooms + bathrooms + surface_total + property_type,data = df)
summary(lm_fit)##
## Call:
## lm(formula = price ~ barrio + rooms + bathrooms + surface_total +
## property_type, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2143881 -67303 -11936 39417 5319777
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.936e+05 6.375e+03 -30.362 < 2e-16 ***
## barrioalto 7.399e+04 2.640e+03 28.029 < 2e-16 ***
## barriobajo -1.766e+04 3.789e+03 -4.662 3.14e-06 ***
## rooms 4.709e+04 9.366e+02 50.274 < 2e-16 ***
## bathrooms 1.633e+05 1.508e+03 108.275 < 2e-16 ***
## surface_total 5.683e+00 1.208e+00 4.706 2.54e-06 ***
## property_typeDepartamento 2.395e+04 5.410e+03 4.427 9.58e-06 ***
## property_typePH -4.085e+04 5.969e+03 -6.844 7.80e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 222700 on 52238 degrees of freedom
## Multiple R-squared: 0.4243, Adjusted R-squared: 0.4242
## F-statistic: 5499 on 7 and 52238 DF, p-value: < 2.2e-16Cuál es el mejor modelo? Con los barrios colapsados en tres niveles, o con todos los barrios?
8.2.3 Transformaciones log.
df %>%
sample_n(10000) %>%
ggplot(., aes(surface_total, price))+
geom_point(alpha=0.75)+
geom_smooth(method = 'lm')
Nuestros datos no tienen varianza constante, tienen heterocedasticidad. Que la varianza sea constante es uno de los supestos del modelo lineal.
df %>%
sample_n(10000) %>%
ggplot(., aes(surface_total, price))+
geom_point(alpha=0.75)+
geom_smooth(method = 'lm')+
scale_y_log10()+
scale_x_log10()
La transformación logarítmica parece ser una buena idea para nuestros datos
lm_fit <- lm(log(price)~ barrio+ rooms + bathrooms + log(surface_total) + property_type,data = df)
summary(lm_fit)##
## Call:
## lm(formula = log(price) ~ barrio + rooms + bathrooms + log(surface_total) +
## property_type, data = df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.2874 -0.1782 -0.0097 0.1591 3.6961
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 8.449601 0.016879 500.589 < 2e-16 ***
## barrioalto 0.279975 0.003759 74.486 < 2e-16 ***
## barriobajo -0.134716 0.005393 -24.979 < 2e-16 ***
## rooms 0.020981 0.001680 12.486 < 2e-16 ***
## bathrooms 0.171951 0.002311 74.400 < 2e-16 ***
## log(surface_total) 0.697528 0.003919 177.996 < 2e-16 ***
## property_typeDepartamento 0.238179 0.007904 30.132 < 2e-16 ***
## property_typePH 0.023308 0.008551 2.726 0.00642 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.317 on 52238 degrees of freedom
## Multiple R-squared: 0.7759, Adjusted R-squared: 0.7759
## F-statistic: 2.584e+04 on 7 and 52238 DF, p-value: < 2.2e-16El \(R^2\) pasó de 0.6977 a 0.7833!!
Cómo se interpretan ahora los nuevos parámetros?
| Modelo | Variable Dependiente | Variable Independiente | Interpretación \(\beta_1\) |
|---|---|---|---|
| Nivel-Nivel | \(y\) | \(x\) | \(\Delta y = \beta_1 \Delta x\) |
| Nivel-log | \(y\) | \(log(x)\) | \(\Delta y = (\beta_1/100)\% \Delta x\) |
| Log-nivel | \(log(y)\) | \(x\) | \(\% \Delta y = (100 \beta_1) \Delta x\) |
| Log-log | \(log(y)\) | \(log(x)\) | \(\% \Delta y = \beta_1 \% \Delta x\) |
8.2.4 Predicciones
Para predecir un nuevo caso, podemos construir un dataframe con las variables. Por ejemplo
caso_nuevo <- tibble(barrio='alto',
rooms=2,
bathrooms=2,
property_type='Departamento',
surface_total=78)
predict(lm_fit,newdata = caso_nuevo)## 1
## 12.39254Pero debemos recordar que este es el valor del logarítmo del percio. Por lo tanto tenemos que realizar el caminio inverso con la función exp
## 1
## 2409988.2.5 Para seguir practicando
Un problema de lo que vimos en esta práctica es que las salidas de summary(lm_fit) es una impresión en la consola. Es muy difícil seguir trabajando con esos resultados. Para resolver esto hay un par de librerías que incorporan el modelado lineal al flujo del tidyverse: