Sesion 2: Regresion lineal multiple
Sergio Diaz
2023-06-04
Regresion polinomial
En ocaciones la relacion entre una variable x e y no estan representados por una linea recta, es por ello que
podemos adaptar nuestro modelo incluyendo diversos grados de potencia de x.
Relaciones no lineales
df = [Link]("[Link]
plot(df$x,df$y)
14
12
df$y
10
8
6
4
5 10 15
df$x
1
En estos datos podemos ver como hay una evidente relacion entre las variables x e y, sin embargo si evaluamos
esto desde una perspectiva lineal encontraremos que no podra ser identificado por el coeficiente de correlacion
de pearson.
cor(df)
## x y
## x 1.00000000 0.07191547
## y 0.07191547 1.00000000
El resultado nos indica que hay una nula relacion entre ambas variables. Para mejorar esto vamos a darle a
la relacion de x e y el ser mas flexible por medio de sus potencias.
• Polinomio de grado 1:
yi = β 0 + β 1 xi
• Polinomio de grado 2:
yi = β0 + β1 xi + β2 x2i
• Polinomio de grado 3:
yi = β0 + β1 xi + β2 x2i + β3 x3i
• Polinomio de grado 4:
yi = β0 + β1 xi + β2 x2i + β3 x3i + β4 x4i
Creacion de valores cuadraticos, cubicos, . . .
Actualmente solo existen en df dos columnas x e y. El siguien codigo crea nuevas columnas x_2, x_3 y x_4
como las verciones polinomiales de x.
df$x_2 <- df$xˆ2
df$x_3 <- df$xˆ3
df$x_4 <- df$xˆ4
df$x_5 <- df$xˆ5
df$x_6 <- df$xˆ6
Ahora nuestro data frame cuenta con nuevas columnas para ser empleadas en el modelo.
head(df)
2
## x y x_2 x_3 x_4 x_5 x_6
## 1 0.60 15.52 0.3600 0.216000 0.1296000 0.0777600 0.0466560
## 2 0.82 15.16 0.6724 0.551368 0.4521218 0.3707398 0.3040067
## 3 1.12 14.13 1.2544 1.404928 1.5735194 1.7623417 1.9738227
## 4 1.87 10.87 3.4969 6.539203 12.2283096 22.8669390 42.7611759
## 5 1.90 11.84 3.6100 6.859000 13.0321000 24.7609900 47.0458810
## 6 2.11 10.62 4.4521 9.393931 19.8211944 41.8227202 88.2459396
Ahora estimemos los modelos necesarios.
pol1 <- lm(y ~ x,df)
pol2 <- lm(y ~ x + x_2,df)
pol3 <- lm(y ~ x + x_2 + x_3,df)
pol4 <- lm(y ~ x + x_2 + x_3 + x_4,df)
pol5 <- lm(y ~ x + x_2 + x_3 + x_4 + x_5,df)
pol6 <- lm(y ~ x + x_2 + x_3 + x_4 + x_5 + x_6,df)
Los objetos pol1, pol2, pol3 y pol4 tienen informacion de los modelos estimados y podemos obtener
mucha informacion de ellos. Ahora extraeremos los valores ajustados de la recta. por medio del vector
[Link] presente en estos modelos.
pol4$[Link]
## 1 2 3 4 5 6 7 8
## 17.116614 15.631749 13.828994 10.330995 10.218324 9.483132 7.367465 7.281249
## 9 10 11 12 13 14 15 16
## 6.679821 6.468449 6.226064 6.194958 6.208712 6.312381 6.338508 6.689787
## 17 18 19 20 21 22 23 24
## 6.978384 7.173423 7.650578 7.799127 7.809075 8.619047 8.844917 8.972454
## 25 26 27 28 29 30 31 32
## 9.395788 9.606804 9.911891 10.081381 10.132175 10.352988 10.393876 10.402162
## 33 34 35 36 37 38 39 40
## 10.444605 10.423397 10.419203 10.276909 10.192949 10.169896 10.081261 9.944090
## 41 42 43 44 45
## 9.935574 9.929372 9.936867 9.942827 10.021776
Antes de graficar si queremos conocer los parametros del modelo solo debemos buscar en el objeto
coefficients.
pol4$coefficients
## (Intercept) x x_2 x_3 x_4
## 21.931604158 -9.018120107 1.729673217 -0.125963151 0.003146755
Lo que debemos hacer es llevar estos numeros a nuestra ecuacion inicial.
yi = 21.93 + −9.02xi + 1.73x2i + −0.13x3i + 0.003x4i
Suavizamiento de la linea de tendencia
Ahora si veamos el cambio en la linea de tendencia
3
par(mfrow=c(1,2))
plot(df$x,df$y,main="Grado 1", xlab = "x", ylab = "y")
lines(df$x,pol1$[Link],col="red")
plot(df$x,df$y,main="Grado 2", xlab = "x", ylab = "y")
lines(df$x,pol2$[Link],col="red")
Grado 1 Grado 2
14
14
12
12
10
10
y
y
8
8
6
6
4
5 10 15 5 10 15
x x
par(mfrow=c(1,2))
plot(df$x,df$y,main="Grado 3", xlab = "x", ylab = "y")
lines(df$x,pol3$[Link],col="red")
plot(df$x,df$y,main="Grado 4", xlab = "x", ylab = "y")
lines(df$x,pol4$[Link],col="red")
4
Grado 3 Grado 4
14
14
12
12
10
10
y
y
8
8
6
6
4
5 10 15 4 5 10 15
x x
par(mfrow=c(1,2))
plot(df$x,df$y,main="Grado 5", xlab = "x", ylab = "y")
lines(df$x,pol5$[Link],col="red")
plot(df$x,df$y,main="Grado 6", xlab = "x", ylab = "y")
lines(df$x,pol6$[Link],col="red")
5
Grado 5 Grado 6
14
14
12
12
10
10
y
y
8
8
6
6
4
5 10 15 4 5 10 15
x x
Criterios de eleccion entre modelos
Como el analisis visual no resulta suficiente se requiere conocer algunas metricas objetivas sobre el ajuste.
Para ello se emplea la funcion summary() para obtener las principales metricas del modelo.
s_pol1 <- summary(pol1)
s_pol2 <- summary(pol2)
s_pol3 <- summary(pol3)
s_pol4 <- summary(pol4)
R cuadrado
El coeficiente de determinación, comúnmente conocido como R cuadrado(R2 ), es una medida estadística
que proporciona información sobre la proporción de la variabilidad de una variable dependiente que puede
ser explicada por una variable independiente o un conjunto de variables independientes en un modelo de
regresión.
print(paste("Polinomio 1: Rcuadrado=",s_pol1$[Link]))
## [1] "Polinomio 1: Rcuadrado= 0.00517183431182067"
6
print(paste("Polinomio 2: Rcuadrado=",s_pol2$[Link]))
## [1] "Polinomio 2: Rcuadrado= 0.280170159438744"
print(paste("Polinomio 3: Rcuadrado=",s_pol3$[Link]))
## [1] "Polinomio 3: Rcuadrado= 0.783911690866406"
print(paste("Polinomio 4: Rcuadrado=",s_pol4$[Link]))
## [1] "Polinomio 4: Rcuadrado= 0.868469940319047"
R cuadrado ajustado
El R cuadrado ajustado es una medida estadística relacionada con el coeficiente de determinación (R
cuadrado) que tiene en cuenta el número de variables independientes en un modelo de regresión. A diferencia
del R cuadrado, que tiende a aumentar al agregar más variables independientes al modelo, el R cuadrado
ajustado penaliza el aumento en el número de variables, evitando la sobreestimación del ajuste del modelo.
print(paste("Polinomio 1: Rcuadrado ajustado=",s_pol1$[Link]))
## [1] "Polinomio 1: Rcuadrado ajustado= -0.0179637044251137"
print(paste("Polinomio 2: Rcuadrado ajustado=",s_pol2$[Link]))
## [1] "Polinomio 2: Rcuadrado ajustado= 0.245892547983446"
print(paste("Polinomio 3: Rcuadrado ajustado=",s_pol3$[Link]))
## [1] "Polinomio 3: Rcuadrado ajustado= 0.768100351173703"
print(paste("Polinomio 4: Rcuadrado ajustado=",s_pol4$[Link]))
## [1] "Polinomio 4: Rcuadrado ajustado= 0.855316934350952"
Variables dummy
Otra manera de incluir una mas variables en un ecuacion que querramos ver en dos dimenciones es por
medio de variables dummy. La cual emplea variables categoricas para poder contribuir con brindarle mayor
informacion al modelo.
Caso de estudio
Una inmobiliaria acaba de finalizar la venta de un proyecto. Debido a las diversas promociones, campañas
y pagos anticipados de los clientes la venta de cada inmueble ha generado distintos ingresos por lo que se
requiere hacer un analisis que nos permita valorizar los atributos de los terrenos vendidos.
7
terrenos = [Link]("[Link]
head(terrenos)
## precio m2 tipo
## 1 40507 91.6 NORMAL
## 2 36811 92.9 CARRETERA
## 3 49010 93.0 PARQUE
## 4 50898 95.4 PARQUE
## 5 34838 87.6 CARRETERA
## 6 53411 95.9 PARQUE
Lo mas evidente en estos casos es observar la relacion entre el precio y el metro cuadrado.
plot(terrenos$m2,terrenos$precio)
50000
terrenos$precio
40000
30000
85 90 95 100 105
terrenos$m2
Efectos categoricos
La relacion es creciente sin embargo no parece haber un buen ajuste de los puntos. Continuamos analizando
los datos y observamos que cada lote es de un tipo distinto:
• Parque: Singnifica que el terreno esta frente a un parque por lo cual se pidio un precio mayor.
• Carretera: Esto singnifica que el lote esta cerca del ruido de una via principal por lo que al ser mas
dificil de vender perdio valor.
8
• Normal: Estos terrenos carecen de una ventaja o desventaja respecto de los otros grupos.
plot(terrenos$m2,terrenos$precio,col=factor(terrenos$tipo))
legend("topleft", legend = c("Parque", "Comun","Carretera"),
col=c("green", "red","black"),
lty = c("solid", "solid", "solid"))
Parque
Comun
Carretera
50000
terrenos$precio
40000
30000
85 90 95 100 105
terrenos$m2
Es evidente que el factor tipo influye en el precio del terreno vendido ya que los terrenos frente a parques
tienen una ventaja y los que estan cerca a la carretera una desventaja. Sin embargo el modelo de regresion
lineal simple se observa:
precioi = β0 + β1 m2i
Segun la grafica la pendiente beta parece ser la misma entre los grupos sin embargo es el intercepto el que
requiere una modificacion. Requeririamos lo siguiente:
• Un modelo para los lotes frente a parques
precioi = βparque + β1 m2i
• Un modelo para los lotes normales
precioi = β0 + β1 m2i
9
• Un modelo para los lotes frente a carreteras
precioi = βcarretera + β1 m2i
Incluso podriamos intuir:
βparque > β0 > βcarretera
Una posible salida seria estimar por separado las rectas de cada grupo, sin embargo correriamos el riesgo de
perjudicar la estimacion de β1 al reducir el tamaño de datos y realizar el proceso por separado.
Variable binaria
Para poder estimar los tres interceptos en un solo proceso de estimacion vamos a recurir ha representar la
columna tipo como 3 variables binarias, es decir que presenten valores de 0 y 1. De tal manera que se
creara una columna parque que tendra valor 1 tenga un parque en frente y 0 en caso contrario. El proceso
se repetira para las otras categorias y en estas nuevas columnas ninguna fila tendra mas de un 1 y tampoco
podra estar compuesta de puros 0.
terrenos$parque <- ifelse(terrenos$tipo == "PARQUE",1,0)
terrenos$normal <- ifelse(terrenos$tipo == "NORMAL",1,0)
terrenos$carretera <- ifelse(terrenos$tipo == "CARRETERA",1,0)
head(terrenos)
## precio m2 tipo parque normal carretera
## 1 40507 91.6 NORMAL 0 1 0
## 2 36811 92.9 CARRETERA 0 0 1
## 3 49010 93.0 PARQUE 1 0 0
## 4 50898 95.4 PARQUE 1 0 0
## 5 34838 87.6 CARRETERA 0 0 1
## 6 53411 95.9 PARQUE 1 0 0
La funcion ifelse() se emplea para evaluar un criterio logico y asignar un resultado especifico cuando
resulte ser verdadero o falso.
Modelando ecuacion con variables dummy
Habiamos indicado que requeriamos estimar βparque y βcarretera ya que para el grupo de lotes normales nos
servia β0 al no tener un atributo resaltante. Adicionalmente a esta condicion del caso de estudio resulta
inviable matematicamente incluir las tres variables binarias creadas reciente mente, es por ello que la variable
no incluida debe ser de preferencia un referente de las otras. El modelo final quedaria asi.
precioi = β0 + βparque parque + βcarretera carretera + β1 m2i
10
Grupo frente a un parque
De tal manera que cuando un terreno tenga un parque en frente βparque = 1 y βcarretera = 0, entonces:
precioi = β0 + βparque 1 + βcarretera 0 + β1 m2i
precioi = (β0 + βparque ) + β1 m2i
De esta manera el nuevo intercepto para el grupo frente a un parque seria (β0 + βparque ) y βparque seria la
estimacion del valor adicional que gana un terreno al estar frente a un parque.
Grupo normal
Cuando el terreno sea normal βparque = 0 y βcarretera = 0, entonces:
precioi = β0 + βparque 0 + βcarretera 0 + β1 m2i
precioi = β0 + β1 m2i
Este grupo referencia ingnoraria los otros parametros estimados.
Grupo cerca a una carretera
En el caso de cuando un terreno tenga una carretera cerca βparque = 0 y βcarretera = 1, entonces:
precioi = β0 + βparque 0 + βcarretera 1 + β1 m2i
precioi = (β0 + βcarretera ) + β1 m2i
De esta manera el nuevo intercepto para el grupo frente a un parque seria (β0 + βcarretera ) y βcarretera seria
la estimacion del valor adicional que gana un terreno al estar cerca de una carretera. Por la grafica observada
previamente podriamos intuir que este parametro es negativo.
Estimando el modelo
m_terrenos <- lm(precio ~ m2 + parque + carretera, terrenos)
m_terrenos$coefficients
## (Intercept) m2 parque carretera
## 4324.1299 404.6695 8240.4020 -6007.3261
b1 <- m_terrenos$coefficients["m2"]
b0 <- m_terrenos$coefficients["(Intercept)"]
bparque <- m_terrenos$coefficients["parque"]
bcarretera <- m_terrenos$coefficients["carretera"]
Interpretando:
• β0 : El valor fijo de los lotes normales es 4324.
• β1 : Un metro cuadrado adicional cuesta en promedio 404
11
• βparque : El valor agregado por tener un parque al frente es 8240
• βcarretera : El valor perjuicio economico por el ruido de la carretera es 6007.
Ahora nuestra ecuacion resulta:
precioi = 4324 + 8240parque − 6007carretera + 404m2i
Visualizando modelo
plot(terrenos$m2,terrenos$precio,col=factor(terrenos$tipo))
abline(a= b0+bparque, b=b1, col="green")
abline(a= b0,b=b1, col="red")
abline(a= b0+bcarretera, b=b1, col="black")
legend("topleft", legend = c("Parque", "Comun","Carretera"),
col=c("green", "red","black"),
lty = c("solid", "solid", "solid"))
Parque
Comun
Carretera
50000
terrenos$precio
40000
30000
85 90 95 100 105
terrenos$m2
12