Análisis descriptivo

A continuación, se hará un breve análisis univariado de las variables precio_millon, Area_construida y Tipo, con el fin de observar el comportamiento de cada variable por separado, para, posteriormente hacer un análisis bivariado que permita estudiar las relaciones entre ellas y a partir de ahí, ajustar el mejor modelo.

Para iniciar el procesamiento, se procede a cargar la base de datos y se eliminan las filas con datos faltantes:

Datos <- na.omit(read_excel("Datos_Vivienda.xlsx"))
attach(Datos)
head(Datos)
## # A tibble: 6 x 12
##   Zona       piso  Estrato precio_millon Area_contruida parqueaderos Banos
##   <chr>      <chr>   <dbl>         <dbl>          <dbl> <chr>        <dbl>
## 1 Zona Sur   2           6           880            237 2                5
## 2 Zona Oeste 2           4          1200            800 3                6
## 3 Zona Sur   3           5           250             86 NA               2
## 4 Zona Sur   NA          6          1280            346 4                6
## 5 Zona Sur   2           6          1300            600 4                7
## 6 Zona Sur   3           6           513            160 2                4
## # ... with 5 more variables: Habitaciones <dbl>, Tipo <chr>, Barrio <chr>,
## #   cordenada_longitud <dbl>, Cordenada_latitud <dbl>

Análisis univariado

Para el análisis univariado, tenemos dos variables numéricas y una nominal. Entonces, procedemos a visualizar las variables precio_millon, Area_construida, por medio de un gráfico boxplot.

par(mfrow=c(1,2))
boxplot(precio_millon, xlab = "Precio de vivienda (Millones)", ylab="Millones",
        col="green")
boxplot(Area_contruida, xlab="Area construida", ylab="Metros cuadrados",
        col="purple")

En el gráfico se puede observar que, ambas variables tienen un comportamiento asimétrico positivo (o a la derecha), lo que a primera vista permite ver que su mejor estimador es la Mediana (Me). Sin embargo, lo coprobaremos revisando sus estadísticos descriptivos:

summary(precio_millon)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    58.0   220.0   330.0   433.9   540.0  1999.0
summary(Area_contruida)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    30.0    80.0   123.0   174.9   229.0  1745.0

Los estadísticos descriptivos permiten ver, que efectivamente el promedio se vio afectado por los valores extremos tanto en el precio como en el área cosntruida de las viviendad, pues, vemos que 50% de las viviendas tiene un costo de 330 millones y un área contruida de 123m2, sin embargo sus promedios toman valores superiores y alejados de éstos, lo que corrobora la afirmación, de que la media no sería el estimador más adecuado para evaluar el comportamiento de los datos. Adicionalmente, la vivienda de menos precio cuesta 53 millones y la más costosa, 1999 millines y su area construida es de 30 y 1745 m2, respectivamente.

Sobre la variable Tipo, vemos que es dicotómica y nominal, en ella se guardan los datos referentes al tipo de vivienda y se distribuye como se muestra a continuación:

par(cex.axis= 0.7, cex.lab= 0.5, cex.sub= 0.7, las=1)
tip=round((prop.table(table(Tipo))*100),2)
tab=barplot(tip,col=c("mediumorchid4","mediumvioletred"),main="Tipos de viviendas",
        ylab="Frecuencia", ylim=c(0,100), cex.main=0.8)
text(tab, tip+3,format(tip))

Como se ve, el 61,31% de las viviendas, son apartamentos y el restante 38,69% de ellas, son casas.

Análisis bivariado.

A continución se muestra un gráfico de dispersión que permite ver la interacción entre las variables precio_millon y Area_construida.

plot(Area_contruida,precio_millon,col="mediumpurple4",xlab=expression("Área construida" (m^{2})),ylab="Precio (millones)")

En el gráfico se observa una aparente relación lineal positiva fuerte entre las variables sobre todo entre las viviendas con un área construida de 1500 m2 y cuyo costo es menor o igual a los 1500 millones. No obstante, para revisar si realmente las variables están relacionadas, se realiza un test de correlación por el método de Pearson.

cor.test(precio_millon, Area_contruida)
## 
##  Pearson's product-moment correlation
## 
## data:  precio_millon and Area_contruida
## t = 86.304, df = 8317, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.6758453 0.6985236
## sample estimates:
##      cor 
## 0.687352

Dado que el coficiente de correlación es 0.68, se puede determinar que efectivamente las variables precio_millon y Area_construida están linealmente correlacionadas. Por tanto, se puede realizar un modelo de regresión que muestre en qué medida el área construida de una vivienda, puede determinar su precio en millones.

Modelo de regresión lineal

El modelo que se propondrá busca modelar la relación entre las variables precio_millon (variable de respuesta), Area_construida y el Tipo (variables libres). Lo que se quiere estudiar es, la manera en que las covariables logran explicar el precio de las viviendas.

Como la variable Tipo es una variable nominal, es necesario recodificarla, para poder incluirla en el modelo. Esto se hará usando la función dummy_cols de la librería fastDummies, que creará una columna con valores binarios para Apartamentos y Casas por separado. Así pues, la data quedaría con éstas columas.

library(fastDummies)
Datos=dummy_cols(Datos,select_columns = "Tipo")
names(Datos)
##  [1] "Zona"               "piso"               "Estrato"           
##  [4] "precio_millon"      "Area_contruida"     "parqueaderos"      
##  [7] "Banos"              "Habitaciones"       "Tipo"              
## [10] "Barrio"             "cordenada_longitud" "Cordenada_latitud" 
## [13] "Tipo_Apartamento"   "Tipo_Casa"

Donde la variable auxiliar Tipo_Apartamento tendrá esta forma: \[ Tipo\_Apartamento=\left\{ \begin{array}{ll} 1 \ \ si \ \ Tipo=Apartamento\\ 0 \ \ si \ \ Tipo=Casa \end{array} \right.\]

Ahora, se ajusta el modelo por el método de Mínimos Cuadrados Orinarios (MCO), con la variable auxiliar Tipo_Apartamento.

mod=lm(precio_millon~Area_contruida+Tipo_Apartamento,data=Datos)
coefficients(mod)
##      (Intercept)   Area_contruida Tipo_Apartamento 
##         49.19021          1.79514        115.29391

El modelo ajustado es de la forma: \[Precio=49,19021+1,79514(Area\_Construida)+115.2939(Tipo\_Apartamento)\] Del modelo, se observa que el estimador de \(\hat{\beta_0}\) toma un valor bastante alejado de cero, por tanto, éste no se interpreta y se toma sólo como un valor necesario para ajustar correctamente el modelo. Por otro lado, el \(\hat{\beta_1}\) nos indica que por cada m2, el precio de la vivienda aumentará 1.8 millones aproximadamente, mientras que el \(\hat{\beta_2}\) permite ver que el el costo de las viviendas puede aumentar alrededor de 115.3 millones si la vivienda es un apartamento.

De lo anterior cabe mencionar que, el hecho de que el precio de la vivienda cueste 115.3 millones adicionales sólo por ser apartamento, es una cifra bastante grande, que puede estarse dando por posibles incumplimientos en los supuestos del modelo, pero éstos se evaluaran más adelante.

Intervalo de confianza del 95%

Se calcula un intervalo de confianza que nos ayude a determinar bajo qué rango de valores se pueden ubicar los verdaderos parámetros del modelo (\(\beta_0\), \(\beta_1\) y \(\beta_2\)).

confint(mod)
##                       2.5 %     97.5 %
## (Intercept)       35.135015  63.245410
## Area_contruida     1.753095   1.837186
## Tipo_Apartamento 102.953133 127.634687

Con un 95% de confianza, se tiene que: el \(\beta_0\) seguirá siendo un parámetro sin interpretación, pero se esperaría que su valor real se encuentre en un intervalo que va desde 35.2 a 63.3; el parámetro \(\beta_1\), se esperaría que tomase valores entre 1.75 y 1.84 (nótese que éste intervalo tiene una amplitud de 0.09); y que el parámetro \(\beta_2\), tome valores que pertenezcan al rango de valores que va desde 102.95 a 127.63.

#Rojo casa , Negro apartamento
mod1=lm(precio_millon~Area_contruida+Tipo_Apartamento+Area_contruida*Tipo_Apartamento,data=Datos)


plot(Datos$precio_millon ~ Datos$Area_contruida ,data=Datos,col=Datos$Tipo_Casa+1,ylab='Precio millon vivienda',xlab='Area construida',main='Casa-Apartamento')
abline(a=mod1$coefficients[1],b=mod1$coefficients[2],lwd=2)
abline(a=mod1$coefficients[1]+mod1$coefficients[3],b=mod1$coefficients[2]+mod1$coefficients[4],col=2,lwd=2)
legend("bottomright", legend = c("Casa", "Apartamento"),
       lwd = 3, col = c("red", "black"))

Donde se observa que la pendiente de las casas es mas pronunciada esto afectado al precio de la vivienda por lo que se observa una diferencia entre ambas aun cuando tienen el mismo area construida, por lo tanto el precio de las vivientas de tipo casa es mayor que las de apartamento.

Prueba de hipótesis para los estimadores \(\beta\)

Se tiene la siguiente hipótesis: \[H_0: \beta_i = 0 \ \ vs \ \ H_1:\beta_i \neq 0 \ \ \ \ | \ \ i=0,1,2 \] Para contrastar la hipótesis, la función summary() de R, calcula los respectivos valores críticos de la tabla t y los Valor-P correspondientes.

summary(mod)
## 
## Call:
## lm(formula = precio_millon ~ Area_contruida + Tipo_Apartamento, 
##     data = Datos)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2926.71  -122.58   -49.12    71.12  1276.36 
## 
## Coefficients:
##                   Estimate Std. Error t value Pr(>|t|)    
## (Intercept)       49.19021    7.17011    6.86 7.36e-12 ***
## Area_contruida     1.79514    0.02145   83.69  < 2e-16 ***
## Tipo_Apartamento 115.29391    6.29551   18.31  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 234.1 on 8316 degrees of freedom
## Multiple R-squared:  0.4929, Adjusted R-squared:  0.4928 
## F-statistic:  4042 on 2 and 8316 DF,  p-value: < 2.2e-16

Como se puede ver, el Valor-P de los diferentes parámetros es menor que la significancia \(\alpha=0.05\), de lo que se sigue que, la hipótesis de que el parámetro \(\beta_i=0\) se debe rechazar.

Indicador de bondad y ajuste R2

Del resúmen anterior, se obtuvo que la bondad de ajuste del modelo es del 0.49 lo que indica que el modelo sólo está explicando un \(49\%\) de la varianza de los datos. Por tanto, este no es un modelo confiable.

Precio promedio estimado para un apartamento de 110 metros cuadrados

Ahora, se desea estimar el precio promedio de los apartamentos con un área construida de 110 m2. El resultado muestra la media y el intervalo de las 10 primeras observaciones.

x.nuevo=data.frame(Area_contruida = 110)
predict(mod2, x.nuevo, interval = 'prediction') 
##        fit      lwr      upr
## 1 331.2958 -136.712 799.3036

fit = media, lwr y upr = intervalo de confiaza del 95% para el un apartamento de 110 metros cuadrados

Con un valor de 331.29 millones de pesos en promedio por un apartamento de 110 metros cuadrados, en comparación de un apartamento en la misma zona con 110 metros cuadrados que tiene un precio de 200 millones esa vivienda tiene un menor precio que el promedio de los apartamentos con estas características por lo tanto sería una buena oferta para comprar, teniendo en cuenta que no solo el área debe importar en el precio final de la vivienda, dado que pueden tener número de habitaciones diferentes u otros agregados de los apartamentos.

Validación de supuestos

A continuación, se procede a realizar la validación de los supuestos del modelo.

res.stud.MOD1 = studres(mod)
mod.fit.MOD1 = mod$fitted.values
par(mfrow=c(1,2))
plot(mod.fit.MOD1,res.stud.MOD1, ylab='Residuos estudentizados',
     xlab='Valores ajustados',main='(Grafico residuos estudentizados )')
abline(h=0,lty=2)
lines(lowess(res.stud.MOD1~mod.fit.MOD1), col = 2)
plot(mod.fit.MOD1,abs(res.stud.MOD1), 
     ylab='|Residuos estudentizados|',
     xlab='Valores ajustados',main='(Grafico residuos estudentizados absolutos)')
lines(lowess(abs(res.stud.MOD1)~mod.fit.MOD1), col = 2)

Del gráfico anterior, podemos observar que no parece haber un patrón determinado en los residuos del modelo, cuando los valores ajustados del modelo son menores a 1500; para los valores ajustados mayores a 1500, se observa que los residuos dejan de estar alrededor de cero, lo que muestra un indicio de heterocesdasticidad.

Prueba de homocedasticidad

Para motrar si efectivamente hay heterocedasticidad en el modelo, se realiza la prueba de Breush-Pagan en R.

library(lmtest)
bptest(mod, ~Area_contruida +I(Area_contruida^2) + Datos$Tipo_Apartamento+I(Datos$Tipo_Apartamento^{2}), data=Datos)
## 
##  studentized Breusch-Pagan test
## 
## data:  mod
## BP = 2551.9, df = 3, p-value < 2.2e-16

Teniendo en cuenta que la prueba de Breush-Pagan se realiza bajo la hipótesis nula de igual de varianza entre los errores y con una confianza del 95%, a partir del \(Valor-P= 2.2e^{-16}\), la \(H_0\) se rechaza, por tanto se demuestra que efectivamente el modelo incumple el supuesto de homocedasticidad.

Prueba de normalidad

Para evaluar el supuesto de normalidad, se realiza el test de Anderson-Darling, desde R.

\[H_0:La \ \ distribución \ \ de \ \ X \ \ es \ \ normal.\]

\[H_1:La \ \ distribución \ \ de \ \ X \ \ no \ \ es \ \ normal.\]

qqPlot(res.stud.MOD1)

## [1] 3324 1017
ad.test(res.stud.MOD1)
## 
##  Anderson-Darling normality test
## 
## data:  res.stud.MOD1
## A = 355.85, p-value < 2.2e-16

Del gráfico se puede observar que se la ubicación de cada punto quantil-quantil, no se ajusta a la línea de regresión normal que se pinta en azul, lo que muestra la falta de normalidad en la distribución de los errores. La aplicación de la prueba de Anderson-Darling, confirma el incumplimiento del supuesto de normalidad, pues, para una confianza del 95%, el \(Valor-P\) es casi cero.

Prueba de linealidad

Ahora, se usa un gráfico de la librería car de R, para visulizar si se cumple el supuesto de linealidad.

library(car)
par(mfrow=c(1,2))
crPlots(mod,'Area_contruida',xlab='Area construida')
crPlots(mod,'Tipo_Apartamento',xlab='Apartamentos')

En el gráfico anterior se puede ver que, hay linealidad entre el precio y el área construida de las viviendas, pues, aunque las líneas de regresión no paramétrica (línea morada) y la de regresión lineal simple (línea azul punteada), se separan un poco, no es lo suficientemente significativo. Sobre la relación entre el precio y el tipo de vivienda (puntualmente apartamentos), vemos que amabas líneas están unidas. También es importante tener en cuenta que la manera en que se ven los residuos corresponde a la cualidad binaria de la variable.

Transformación para mejorar el ajuste y supuestos del modelo

A continuación, se realiza la transfomación del modelo usando el método de Box-Cox. Entonces, primero se halla el valor que máximiza la verosimilitud del modelo; luego, se realiza la transformación.

boxcoxPrecio= MASS::boxcox(mod,lambda=seq(-3,3,length.out = 1000),
                               ylab='log-verosimilitud')

boxcoxPrecio$x[boxcoxPrecio$y ==max(boxcoxPrecio$y)] # valor que maximiza la log-verosimilitud
## [1] 0.03303303

Como el \(\hat{\lambda}\) es diferente de cero, se realiza la siguiente transformación.

mod1trans = lm((precio_millon^{0.03}-1)/(0.03*(precio_millon^{0.03}-1))~Area_contruida+Tipo_Apartamento,data=Datos)
resmod1trans = studres(mod1trans)

Luego de realizar la transformación del modelo, se procede a evaluar si se corrigieron los supuestos.

Prueba de homocedasticidad

Al realizar la prueba de homocedasticidad de Breush-Pagan, se puede ver que la transformación corrigió la heterocedasticidad que se presentó en el primer modelo, pues, a partir del \(Valor-P=0,6164\) sabemos que no se debe rechazar la hipótesis de igualdad de varianza entre los errores del modelo transformado.

bptest(mod1trans,~Area_contruida +I(Area_contruida^2) + Datos$Tipo_Apartamento+I(Datos$Tipo_Apartamento^{2}), data=Datos)
## 
##  studentized Breusch-Pagan test
## 
## data:  mod1trans
## BP = 1.7933, df = 3, p-value = 0.6164

Prueba de normalidad

Al realizar la prueba de normalidad de Anderson-Darling, vemos que el supuesto sigue sin cumplirse, por tanto se realizará otra transformación al modelo inicial, para ver si se logra corregir la falta a los dos supuestos.

res.stud.MOD2 = studres(mod1trans)
ad.test(res.stud.MOD2)
## 
##  Anderson-Darling normality test
## 
## data:  res.stud.MOD2
## A = 3213.3, p-value < 2.2e-16

Se realizará la transformación por Mínimos Cuadrados Ponderados (MCP).

mod.stdev.Precio = lm(abs(res.stud.MOD1)~Area_contruida+Tipo_Apartamento,data = Datos)

w = 1/mod.stdev.Precio$fitted.values^2
mod.Precio.mcp = lm(precio_millon~Area_contruida+Tipo_Apartamento,
                  data=Datos,weights = w)
summary(mod.Precio.mcp)
## 
## Call:
## lm(formula = precio_millon ~ Area_contruida + Tipo_Apartamento, 
##     data = Datos, weights = w)
## 
## Weighted Residuals:
##     Min      1Q  Median      3Q     Max 
## -2329.2  -188.9   -41.2   146.2  4143.4 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      29.18492    2.29726   12.70   <2e-16 ***
## Area_contruida    2.15770    0.02086  103.42   <2e-16 ***
## Tipo_Apartamento 47.11572    2.56329   18.38   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 301.5 on 8316 degrees of freedom
## Multiple R-squared:  0.5703, Adjusted R-squared:  0.5702 
## F-statistic:  5518 on 2 and 8316 DF,  p-value: < 2.2e-16

Al realizar las pruebas de normalidad y homocedasticidad, se tiene que para la transformación por MCP, los supuestos no se corrigen.

residuos=studres(mod.Precio.mcp)
ad.test(residuos)
## 
##  Anderson-Darling normality test
## 
## data:  residuos
## A = 115.04, p-value < 2.2e-16
bptest(mod.Precio.mcp,~Area_contruida +I(Area_contruida^2) + Datos$Tipo_Apartamento+I(Datos$Tipo_Apartamento^{2}), data=Datos)
## 
##  studentized Breusch-Pagan test
## 
## data:  mod.Precio.mcp
## BP = 2551.9, df = 3, p-value < 2.2e-16

Comparación de modelo inicial vs modelo transformado

A continuación se presenta la comparación de los ajustes y las pruebas de normalidad y homocedasticidad, entre el modelo inicial (ajuste por Mínimos Cuadrados Ordinarios), el transformado con el método de Box-Cox y el transformado por MCP.

Modelos Inicial Trans. Box-Cox Trans. MCP
R cuadrado 0,4928 0,4999 0,5702
Supuesto normalidad se rechaza Se rechaza Se rechaza
Supuesto homocedasticidad se rechaza No se rechaza Se rechaza

En la tabla se puede ver que, con las transformaciones se mejoró el ajuste del modelo, aunque el más significativo se logró con el método de MCP.

A modo de comentario, puede que el supuesto de normalidad no se corrija, debido a la existencia de la variable dummy Tipo de apartamento.