Unidad 1 -Informe aplicación regresión lineal

1. Predicción de los precios de las acciones. Analizar el comportamiento de los precios de las Acciones de Ecopetrol según la variación del precio del barril de petróleo WTI producido en Colombia. Se tienen los siguientes precios.

prec_ecop=c(1090, 1170, 1160, 1230, 1155, 1165, 1205, 1170, 1150, 1130, 1110, 1105, 
1085, 1060, 1035, 1015, 955, 961) 
  
prec_petr=c(35.62,36.31,37.35,34.95,34.53,35.81,36.14,37.50,37.80,36.81,37.87,37.04,
36.76,35.97,33.97,33.27,31.41,30.44)

data_ecopetrol=data.frame(prec_ecop,prec_petr)

head(data_ecopetrol)

##   prec_ecop prec_petr
## 1      1090     35.62
## 2      1170     36.31
## 3      1160     37.35
## 4      1230     34.95
## 5      1155     34.53
## 6      1165     35.81

cor(data_ecopetrol$prec_ecop,data_ecopetrol$prec_petr,)

## [1] 0.7074373

El coeficiente de correlación entre las dos variables estudiadas es de 0.70, lo que indica una relación positiva, que si una variable la otra también.

A.Proponga un modelo de regresión lineal simple que permita predecir el valor de las Acciones de Ecopetrol con base en el Precio del barril de petróleo en Colombia. Indique la ecuación de regresión y el valor del R cuadrado.

mod=lm(prec_ecop~prec_petr)

summary(mod)

## 
## Call:
## lm(formula = prec_ecop ~ prec_petr)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -59.90 -40.74 -15.94  33.40 136.82 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept)  177.768    232.828   0.764  0.45627   
## prec_petr     26.192      6.542   4.004  0.00102 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 57.13 on 16 degrees of freedom
## Multiple R-squared:  0.5005, Adjusted R-squared:  0.4692 
## F-statistic: 16.03 on 1 and 16 DF,  p-value: 0.001024

Utilizando la función summary identicamos que el R cuadrado que tiene un valor de 0.46, lo que significa que el modelo explica el 46 % de las predicciones de las acciones de Ecopetrol.

B.Pruebe la significancia del modelo propuesto en “a)” plantee las hipótesis respectivas y use el concepto de Valor _p para tomar la decisión sobre las hipótesis. Use α = 0.05

Significancia de Intercepto b0

H0:β0=0

H1:β0≠0

Para el intercerpto, el valor p es de 0.45627, por tal razón, se rechaza H0 ya que el valor p es mayor que el de alfa

Significancia de pendiente b1

H0:β1=0

H1:β1≠0

Para la pendiente, el valor p es de 0.00102, por tal razón, se rechaza H0 ya que el valor p es menor que el de alfa

C.Interprete los coeficientes del modelo propuesto en “a)”

mod$coefficients

## (Intercept)   prec_petr 
##   177.76779    26.19213

Con base en los coeficientes obtenidos se puede afirmar que: -con base en Bo, Si no existe un precio del barril de petroleo, entonces una acción de Ecopetrol tendra un costo de $177.76779 aproximadamente. -Con base en B1, por cada dolar que incremente un barril de petroleo, la acción de Ecopetrol aumentara en $26.19213 dolares.

D. Haga un análisis de los residuos. ¿Qué supuesto no se cumple?

Supuesto 1: Los errores del modelo tienen media cero en el modelo

 summary(mod$residuals)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  -59.90  -40.74  -15.94    0.00   33.40  136.82

t.test(mod$residuals, mu = 0)

## 
##  One Sample t-test
## 
## data:  mod$residuals
## t = -4.2309e-16, df = 17, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -27.56364  27.56364
## sample estimates:
##     mean of x 
## -5.527407e-15

El valor P obtenido es 1 y en consecuencia es mayor que α = 0.05, por lo tanto el supuesto se cumple.

Supuesto 2: Los errores tienen varianza constante

library(lmtest)

lmtest::gqtest(mod)

## 
##  Goldfeld-Quandt test
## 
## data:  mod
## GQ = 0.17924, df1 = 7, df2 = 7, p-value = 0.9813
## alternative hypothesis: variance increases from segment 1 to 2

El valor P obtenido es 0.9813 y en consecuencia es mayor que α = 0.05, por lo tanto el supuesto se cumple.

Supuesto 3: Los errores del modelo se distribuyen normal

shapiro.test(mod$residuals)

## 
##  Shapiro-Wilk normality test
## 
## data:  mod$residuals
## W = 0.89259, p-value = 0.04276

El valor P obtenido es 0.04276 y en consecuencia es menor que α = 0.05, por lo tanto el supuesto no se cumple.

Supuesto 4: Los errores del modelo son independientes.

lmtest::dwtest(mod)

## 
##  Durbin-Watson test
## 
## data:  mod
## DW = 0.74504, p-value = 0.0004666
## alternative hypothesis: true autocorrelation is greater than 0

El valor P obtenido es 0.0004666 y en consecuencia es menor que α = 0.05, por lo tanto el supuesto no se cumple.

E. Concluya sobre la validez del modelo propuesto en a)

Teniendo en cuenta las pasos anteriores, se considera que no es modelo valido ya que no su cumplen todos los supuestos, solamente dos se cumplen, y con base en el R cuadrado, este modelo solo explica el 46,92% de las estimaciones, lo que generá bastante confianza para la toma de decisiones.

2. Los siguientes datos corresponden a la INFLACION y al SALARIO MINIMO LEGAL MENSUAL (SMLM) desde el año 1999 para Colombia.

year=1999:2015

inflacion=c(9.23, 8.75, 7.65, 6.99, 6.49, 5.50, 4.85, 4.48, 5.69, 7.67, 2.00, 3.17, 
3.73, 2.44,  1.94, 3.66, 6.77)

salario=c(236460, 260100, 286000, 309000, 332000, 358000, 381500, 408000, 433700, 
461500, 496900, 515000, 535600, 566700, 589500, 616027, 644350)

datos_salariocol=data.frame(year,inflacion,salario)

head(datos_salariocol)

##   year inflacion salario
## 1 1999      9.23  236460
## 2 2000      8.75  260100
## 3 2001      7.65  286000
## 4 2002      6.99  309000
## 5 2003      6.49  332000
## 6 2004      5.50  358000

A. Escriba la ecuación del modelo de regresión lineal simple

mod2=lm(datos_salariocol$salario~datos_salariocol$inflacion)
summary(mod2)

## 
## Call:
## lm(formula = datos_salariocol$salario ~ datos_salariocol$inflacion)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -75463 -63456 -42854  17623 263207 
## 
## Coefficients:
##                            Estimate Std. Error t value Pr(>|t|)    
## (Intercept)                  648486      58947   11.00  1.4e-08 ***
## datos_salariocol$inflacion   -39489      10151   -3.89  0.00145 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 94130 on 15 degrees of freedom
## Multiple R-squared:  0.5022, Adjusted R-squared:  0.469 
## F-statistic: 15.13 on 1 and 15 DF,  p-value: 0.00145

Identificando los coeficientes, la ecuación que describe el modelo de regresión Y(SML)=648486-39489*X(inflación)

B. Plantee y valide las hipótesis correspondientes a la linealidad general del modelo propuesto en a)

Significancia de Intercepto b0

H0:β0=0

H1:β0≠0

Para el intercerpto, el valor p es de 1.4e-08, por tal razón, se rechaza H0 ya que el valor p es menor que el de alfa.

Significancia de pendiente b1

H0:β1=0

H1:β1≠0

Para la pendiente, el valor p es de 0.00145, por tal razón, se rechaza H0 ya que el valor p es menor que el de alfa.

C.Indique e interprete el coeficiente de correlación del modelo propuesto en a)

cor(datos_salariocol$inflacion,datos_salariocol$salario)

## [1] -0.7086581

El coeficiente de correlación de -0.7086581 indica que existe una relación negativa entra el salario y la inflación, lo que se traduce como que a mayor inflación disminuye el salario minimo.

D.Interprete cada uno de los coeficientes del modelo propuesto en a)

mod2$coefficients

##                (Intercept) datos_salariocol$inflacion 
##                  648485.93                  -39489.33

Con base en los coeficientes obtenidos se puede afirmar que: -con base en Bo, Si la inflación es 0 entonces el salario minimo sera de $ 648.485 aproximadamente. -Con base en B1, por cada incremento en una unidad de inflación, el valor del salario minimo disminuira en $39.489,33 aproximadamente.

E.Construya una gráfica de residuales y haga un análisis cualitativo de los supuestos del modelo propuesto en a)

par(mfrow = c(2, 2))
plot(mod2)

Supuesto 1: Los errores del modelo tienen media cero en el modelo

 summary(mod2$residuals)

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##  -75463  -63456  -42854       0   17623  263207

t.test(mod2$residuals, mu = 0)

## 
##  One Sample t-test
## 
## data:  mod2$residuals
## t = -6.7462e-17, df = 16, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -46862.45  46862.45
## sample estimates:
##     mean of x 
## -1.491304e-12

El valor P obtenido es 1 y en consecuencia es mayor que α = 0.05, por lo tanto el supuesto se cumple.

Supuesto 2: Los errores tienen varianza constante

library(lmtest)

lmtest::gqtest(mod2)

## 
##  Goldfeld-Quandt test
## 
## data:  mod2
## GQ = 140.68, df1 = 7, df2 = 6, p-value = 3.171e-06
## alternative hypothesis: variance increases from segment 1 to 2

El valor P obtenido es 3.171e-06 y en consecuencia es menor que α = 0.05, por lo tanto el supuesto no se cumple.

Supuesto 3: Los errores del modelo se distribuyen normal

shapiro.test(mod2$residuals)

## 
##  Shapiro-Wilk normality test
## 
## data:  mod2$residuals
## W = 0.78826, p-value = 0.001407

El valor P obtenido es 0.001407 y en consecuencia es menor que α = 0.05, por lo tanto el supuesto no se cumple.

Supuesto 4: Los errores del modelo son independientes.

lmtest::dwtest(mod2)

## 
##  Durbin-Watson test
## 
## data:  mod2
## DW = 0.68432, p-value = 0.0002714
## alternative hypothesis: true autocorrelation is greater than 0

El valor P obtenido es 0.0002714 y en consecuencia es menor que α = 0.05, por lo tanto el supuesto no se cumple.

F.Comente sobre la conveniencia de usar el modelo propuesto en a) para predecir el SMLM para Colombia.

No es conveniente usar el modelo propuesto ya que solamente se cumple un supuesto de los cuatro planteados y además el R cuadrado solamente explica el 47% de las estimaciones que se realizan para predecir el salario minimo.

3. Con base en los datos de precios de vivienda de la actividad en clase realizar un informe que contenga los siguientes puntos utilizando R y RMarkdown (publicar el informe final en Rpubs presentando código, resultados e interpretaciones).

library(readxl)
library(ggplot2)
require(plotly)

## Loading required package: plotly

## 
## Attaching package: 'plotly'

## The following object is masked from 'package:ggplot2':
## 
##     last_plot

## The following object is masked from 'package:stats':
## 
##     filter

## The following object is masked from 'package:graphics':
## 
##     layout

bd_vivienda <- read_excel("C:/Users/PACHO/Downloads/bd_vivienda.xlsx")
ID=1:dim(bd_vivienda)[1]
bd_vivienda=data.frame(ID,bd_vivienda)
head(bd_vivienda)

##   ID       Zona piso Estrato precio_millon Area_contruida parqueaderos Banos
## 1  1   Zona Sur    2       6           880            237            2     5
## 2  2 Zona Oeste    2       4          1200            800            3     6
## 3  3   Zona Sur    3       5           250             86           NA     2
## 4  4   Zona Sur   NA       6          1280            346            4     6
## 5  5   Zona Sur    2       6          1300            600            4     7
## 6  6   Zona Sur    3       6           513            160            2     4
##   Habitaciones        Tipo                Barrio cordenada_longitud
## 1            4        Casa                 pance          -76.46300
## 2            7        Casa            miraflores          -76.46400
## 3            3 Apartamento           multicentro          -76.46400
## 4            5 Apartamento ciudad jardv<U+2260>n          -76.46400
## 5            5        Casa                 pance          -76.46438
## 6            4        Casa                 pance          -76.46438
##   Cordenada_latitud
## 1           3.43000
## 2           3.42800
## 3           3.42900
## 4           3.43300
## 5           3.43463
## 6           3.43463

A.Realice un filtro a la base de datos e incluya solo las ofertas de apartamentos, de la zona norte de la ciudad con precios inferiores a los 500 millones de pesos y áreas menores a 300 mt2. Presente los primeros 3 registros de la base y algunas tablas que comprueben la consulta. (Adicional un mapa con los puntos de la base, discutir si todos los puntos se ubican en la zona norte o se presentan valores en otras zonas, por que?).

pos=which(bd_vivienda$Tipo=="Apartamento" & bd_vivienda$Zona=="Zona Norte" & bd_vivienda$precio_millon < 500 & bd_vivienda$Area_contruida < 300)

datos_sub=bd_vivienda[pos,]
datos_sub <- transform(datos_sub, parqueaderos = as.numeric(parqueaderos))

## Warning in eval(substitute(list(...)), `_data`, parent.frame()): NAs
## introducidos por coerción

head(datos_sub,3)

##    ID       Zona piso Estrato precio_millon Area_contruida parqueaderos Banos
## 31 31 Zona Norte    2       3           135             56            1     1
## 71 71 Zona Norte   NA       3            78             54            2     1
## 89 89 Zona Norte   NA       5           340            106            2     2
##    Habitaciones        Tipo             Barrio cordenada_longitud
## 31            3 Apartamento torres de comfandi          -76.46745
## 71            3 Apartamento        chiminangos          -76.47820
## 89            3 Apartamento           la flora          -76.48200
##    Cordenada_latitud
## 31           3.40763
## 71           3.44898
## 89           3.43500

require(leaflet)

leaflet() %>% addCircleMarkers(lng = datos_sub$cordenada_longitud,lat = datos_sub$Cordenada_latitud,radius = 0.3,color = "black",label = datos_sub$ID) %>% addTiles()

En el mapa se observa la ubicación de las viviendas en la zona norte, pero llama la atención que se muestran unas viviendas hacia la zona sur , lo que puede implicar que las coordenadas de longitud y latitud no se encuentran correctas.

B. Realice un análisis exploratorio de datos enfocado en la correlación entre la variable respuesta (precio del apartamento) en función del área construida, estrato y si tiene parqueadero. Use gráficos interactivos con plotly e interprete los resultados.

g1=ggplot(data = datos_sub,aes(y=precio_millon,x=Area_contruida)) + geom_point() + geom_smooth()

ggplotly(g1)

## `geom_smooth()` using method = 'gam' and formula 'y ~ s(x, bs = "cs")'

De acuerdo con el gráfico de dispersión, se puede observar que existe una relación directa entre las variables precio y área construida.

grafico=ggplot(datos_sub,aes(y=precio_millon,x=Estrato)) + geom_point(colour = "black", size = 2) + geom_smooth()
grafico

## `geom_smooth()` using method = 'gam' and formula 'y ~ s(x, bs = "cs")'

## Warning: Computation failed in `stat_smooth()`:
## x has insufficient unique values to support 10 knots: reduce k.

De acuerdo con la gráfica, entre mayor sea el estrato mayor será el precio de la vivienda.

grafico=ggplot(datos_sub,aes(y=precio_millon,x=parqueaderos)) + geom_point(colour = "red", size = 2) + geom_smooth()
grafico

## `geom_smooth()` using method = 'gam' and formula 'y ~ s(x, bs = "cs")'

## Warning: Removed 322 rows containing non-finite values (stat_smooth).

## Warning: Computation failed in `stat_smooth()`:
## x has insufficient unique values to support 10 knots: reduce k.

## Warning: Removed 322 rows containing missing values (geom_point).

Según el gráfico no se puede definir una relación existente entre el numero de parqueaderos vs el precio de la vivienda.

C. Estime un modelo de regresión lineal múltiple con las variables del punto anterior e interprete los coeficientes si son estadísticamente significativos. Las interpretaciones deber están contextualizadas y discutir si los resultados son lógicos. Adicionalmente interprete el coeficiente R2 y discuta el ajuste del modelo e implicaciones (que podrían hacer para mejorarlo).

mod3 = lm(precio_millon ~ Area_contruida + Estrato +  parqueaderos, data = datos_sub )
summary(mod3)

## 
## Call:
## lm(formula = precio_millon ~ Area_contruida + Estrato + parqueaderos, 
##     data = datos_sub)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -184.798  -37.957    0.529   31.854  231.856 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    -172.21312   10.87566 -15.835  < 2e-16 ***
## Area_contruida    0.76566    0.08161   9.382  < 2e-16 ***
## Estrato          69.77587    2.89836  24.074  < 2e-16 ***
## parqueaderos     42.94832    5.24414   8.190 1.12e-15 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 56.62 on 751 degrees of freedom
##   (322 observations deleted due to missingness)
## Multiple R-squared:  0.7096, Adjusted R-squared:  0.7084 
## F-statistic: 611.6 on 3 and 751 DF,  p-value: < 2.2e-16

Con los coeficientes obtenidos la ecucación del modelo de regresión seria Y (Preciomillon) = −172.21312 + 0.76566(Area_contruida) + 69.77587(Estrato) + 42.94832(parqueaderos)

D. Realice la validación de supuestos del modelo e interprete los resultados (no es necesario corregir en caso de presentar problemas solo realizar sugerencias de que se podría hacer).

par(mfrow = c(2, 2))
plot(mod3)

Supuesto 1: Los errores del modelo tienen media cero en el modelo

 summary(mod3$residuals)

##      Min.   1st Qu.    Median      Mean   3rd Qu.      Max. 
## -184.7981  -37.9567    0.5289    0.0000   31.8539  231.8560

t.test(mod3$residuals, mu = 0)

## 
##  One Sample t-test
## 
## data:  mod3$residuals
## t = 6.0668e-16, df = 754, p-value = 1
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -4.037333  4.037333
## sample estimates:
##    mean of x 
## 1.247694e-15

El valor P obtenido es 1 y en consecuencia es mayor que α = 0.05, por lo tanto el supuesto se cumple.

Supuesto 2: Los errores tienen varianza constante

library(lmtest)

lmtest::gqtest(mod3)

## 
##  Goldfeld-Quandt test
## 
## data:  mod3
## GQ = 1.3017, df1 = 374, df2 = 373, p-value = 0.005503
## alternative hypothesis: variance increases from segment 1 to 2

El valor P obtenido es 0.005503 y en consecuencia es menor que α = 0.05, por lo tanto el supuesto no se cumple.

Supuesto 3: Los errores del modelo se distribuyen normal

shapiro.test(mod3$residuals)

## 
##  Shapiro-Wilk normality test
## 
## data:  mod3$residuals
## W = 0.99366, p-value = 0.002832

El valor P obtenido es 0.002832 y en consecuencia es menor que α = 0.05, por lo tanto el supuesto no se cumple.

Supuesto 4: Los errores del modelo son independientes.

lmtest::dwtest(mod3)

## 
##  Durbin-Watson test
## 
## data:  mod3
## DW = 1.7757, p-value = 0.0009303
## alternative hypothesis: true autocorrelation is greater than 0

El valor P obtenido es 0.0009303 y en consecuencia es menor que α = 0.05, por lo tanto el supuesto no se cumple.

E.Con el modelo identificado predecir el precio de un apartamento con 100 mt2, de estrato 4 y con parqueadero. ¿Si este apartamento lo están ofreciendo en 450 millones cual seria su opinión con base en el resultado del modelo considera que es una buena oferta?

#con un parqueadero:
predict(mod3,list(Area_contruida=100, Estrato = 4, parqueaderos = 1), interval = "confidence" )

##        fit      lwr      upr
## 1 226.4048 220.4613 232.3483

#con dos parqueaderos:
predict(mod3,list(Area_contruida=100, Estrato = 4, parqueaderos = 2), interval = "confidence" )

##        fit      lwr      upr
## 1 269.3531 260.3174 278.3888

#con tres parqueaderos:
predict(mod3,list(Area_contruida=100, Estrato = 4, parqueaderos = 3), interval = "confidence" )

##        fit      lwr      upr
## 1 312.3014 293.8641 330.7387

#con cuatro parqueaderos:
predict(mod3,list(Area_contruida=100, Estrato = 4, parqueaderos = 4), interval = "confidence" )

##        fit      lwr      upr
## 1 355.2497 326.7858 383.7136

F. Con las predicciones del modelo sugiera potenciales ofertas para una persona interesada en un apartamento en la zona norte con mas de 100 mt2 de área, de estrato 4, que tenga parqueadero y tenga encuentra que la persona tiene un crédito preaprobado de máximo 400 millones de pesos. Realice un análisis y presente en un mapa al menos 5 ofertas potenciales que debe discutir.

ofertas = which(bd_vivienda$Tipo=="Apartamento"& bd_vivienda$Zona=="Zona Norte"& bd_vivienda$precio_millon<=400 & bd_vivienda$Area_contruida>100 & bd_vivienda$Estrato == 4 & bd_vivienda$parqueaderos >= 1)
ofertas_2=bd_vivienda[ofertas,]
ofertas_final = head(ofertas_2,5)
ofertas_final

##        ID       Zona piso Estrato precio_millon Area_contruida parqueaderos
## 1264 1264 Zona Norte    4       4           380            123            1
## 1612 1612 Zona Norte   NA       4           370            117           NA
## 1668 1668 Zona Norte   NA       4           310            102           NA
## 2186 2186 Zona Norte    1       4           250            160           NA
## 2606 2606 Zona Norte   NA       4           350            130            1
##      Banos Habitaciones        Tipo           Barrio cordenada_longitud
## 1264     3            3 Apartamento         la flora          -76.51437
## 1612     4            4 Apartamento            acopi          -76.51687
## 1668     2            3 Apartamento             Cali          -76.51700
## 2186     2            4 Apartamento prados del norte          -76.51900
## 2606     2            3 Apartamento         la flora          -76.52100
##      Cordenada_latitud
## 1264           3.48618
## 1612           3.45319
## 1668           3.36971
## 2186           3.42158
## 2606           3.49000

oferta1 = predict(mod3,list(Area_contruida=123,Estrato=4,parqueaderos=1))
oferta2 = predict(mod3,list(Area_contruida=130,Estrato=4,parqueaderos=1))
oferta3 = predict(mod3,list(Area_contruida=108,Estrato=4,parqueaderos=1))
oferta4 = predict(mod3,list(Area_contruida=104,Estrato=4,parqueaderos=1))
oferta5 = predict(mod3,list(Area_contruida=125,Estrato=4,parqueaderos=2))
data.frame(oferta1,oferta2,oferta3,oferta4,oferta5)

##   oferta1  oferta2 oferta3  oferta4  oferta5
## 1 244.015 249.3746  232.53 229.4674 288.4946

Teniendo en cuenta el crédito preaprobado de 400 millones, con este dinero puede adquirir alguna de las 5 ofertas que se presentan.

library(leaflet)
leaflet() %>% addCircleMarkers(lng = ofertas_final$cordenada_longitud,lat = ofertas_final$Cordenada_latitud,radius = 5,color = "red") %>% addTiles()