Determinar la ecuación de regresión lineal mútiple y predecir valores
Se muestra cómo utilizar la función de regresión lineal múltiple y con ello se determinan las ecuaciones de regesión lineal múltple para distintos datos de varios ejercicios.
library(dplyr)
library(knitr)
library(readr)
library(ggplot2)
library(patchwork) # PAra varias gráficos en el mismo renglón
options(scipen = 999) # Aplicar notación No científica a los datos
Los ejercicios son tomados de libros y artículos de regresión.
Se sometió a prueba un grupo de camiones ligeros con motores que utilizan diesel como combustible para saber si la humedad, la temperatura del aire y la presión barométrica influyen en la cantidad de óxido nitroso que emiten (en ppm). Las emisiones se midieron en distintos momentos y en diversas condiciones experimentales. (Walpole et al., 2012)
humedad <- c(72.4, 41.6, 34.3, 35.1, 10.7, 12.9, 8.3, 20.1, 72.2, 24.0, 23.2, 47.4, 31.5, 10.6, 11.2, 73.3, 75.4, 96.6, 107.4, 54.9)
temperatura <- c(76.3, 70.3, 77.1, 68.0, 79.0, 67.4, 66.8, 76.9, 77.7, 67.7, 76.8, 86.6, 76.9, 86.3, 86.0, 76.3, 77.9,
78.7, 86.8, 70.9)
presion <- c(29.18, 29.35, 29.24, 29.27, 29.78, 29.39, 29.69, 29.48, 29.09, 29.60, 29.38, 29.35, 29.63, 29.56, 29.48, 29.40, 29.28, 29.29, 29.03, 29.37)
oxido.nitroso <- c(0.90, 0.91, 0.96, 0.89, 1.00, 1.10, 1.15, 1.03, 0.77, 1.07, 1.07, 0.94, 1.10, 1.10, 1.10, 0.91, 0.87, 0.78, 0.82, 0.95)
datos <- data.frame(oxido.nitroso, humedad, temperatura, presion)
kable(datos, caption = "Factores ambientales que influyen en la formación de óxido nitroso en motores disel en camiones")
| oxido.nitroso | humedad | temperatura | presion |
|---|---|---|---|
| 0.90 | 72.4 | 76.3 | 29.18 |
| 0.91 | 41.6 | 70.3 | 29.35 |
| 0.96 | 34.3 | 77.1 | 29.24 |
| 0.89 | 35.1 | 68.0 | 29.27 |
| 1.00 | 10.7 | 79.0 | 29.78 |
| 1.10 | 12.9 | 67.4 | 29.39 |
| 1.15 | 8.3 | 66.8 | 29.69 |
| 1.03 | 20.1 | 76.9 | 29.48 |
| 0.77 | 72.2 | 77.7 | 29.09 |
| 1.07 | 24.0 | 67.7 | 29.60 |
| 1.07 | 23.2 | 76.8 | 29.38 |
| 0.94 | 47.4 | 86.6 | 29.35 |
| 1.10 | 31.5 | 76.9 | 29.63 |
| 1.10 | 10.6 | 86.3 | 29.56 |
| 1.10 | 11.2 | 86.0 | 29.48 |
| 0.91 | 73.3 | 76.3 | 29.40 |
| 0.87 | 75.4 | 77.9 | 29.28 |
| 0.78 | 96.6 | 78.7 | 29.29 |
| 0.82 | 107.4 | 86.8 | 29.03 |
| 0.95 | 54.9 | 70.9 | 29.37 |
g1 <- ggplot(data = datos, mapping = aes(x = humedad, y = oxido.nitroso)) +
geom_point(color = "forestgreen", size = 2) +
labs(title = 'oxido.nitroso ~ humedad', x = 'humedad') +
geom_smooth(method = "lm", se = FALSE, color = "black") +
theme_bw() +
theme(plot.title = element_text(hjust = 0.5))
g2 <- ggplot(data = datos, mapping = aes(x = temperatura, y = oxido.nitroso)) +
geom_point(color = "orange", size = 2) +
labs(title = 'oxido.nitroso ~ temperatura', x = 'tempertura') +
geom_smooth(method = "lm", se = FALSE, color = "black") +
theme_bw() +
theme(plot.title = element_text(hjust = 0.5))
g3 <- ggplot(data = datos, mapping = aes(x = presion, y = oxido.nitroso)) +
geom_point(color = "darkblue", size = 2) +
labs(title = 'oxido.nitroso ~ presion', x = 'presion') +
geom_smooth(method = "lm", se = FALSE, color = "black") +
theme_bw() +
theme(plot.title = element_text(hjust = 0.5))
g1 + g2 + g3
modelo <- lm(formula = oxido.nitroso ~ ., data = datos)
summary(modelo)
##
## Call:
## lm(formula = oxido.nitroso ~ ., data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.11799 -0.02526 0.01345 0.04103 0.06523
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -3.5077781 3.0048641 -1.167 0.26017
## humedad -0.0026250 0.0006549 -4.008 0.00101 **
## temperatura 0.0007989 0.0020451 0.391 0.70121
## presion 0.1541550 0.1013675 1.521 0.14784
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.05617 on 16 degrees of freedom
## Multiple R-squared: 0.8005, Adjusted R-squared: 0.763
## F-statistic: 21.4 on 3 and 16 DF, p-value: 0.000007609
b0 = modelo$coefficients[1]
b1 = modelo$coefficients[2]
b2 = modelo$coefficients[3]
b3 = modelo$coefficients[4]
El resumen del modelo con la función summary(modelo) identifica que las variables temperatura y humedad no son estadísticamente significativas dado que presentan valores por encima de 0.5 en \(Pr(>|t|)\), sólo la humedad es estadísticamente significativa.
\(b0\)= -3.5077781
\(b1\)= -0.002625
\(b2\)= 0.0007989
\(b3\)= 0.154155
\[yˆ=β0+β1(x1)+β2(x2)+β3(x3)\]
\[yˆ=-3.5077781-0.002625(x1)+0.000798(x2)+0.154155(x3)\]
Para 50% de humedad, una temperatura de 76˚F y una presión barométrica de 29.30, ¿cuánto es la cantidad estimada de óxido nitroso emitido? \[yˆ=−3.507778−0.002625(50.0)+0.000799(76.0)+0.1541553(29.30)=0.9384ppm\]
nuevo.dato <- data.frame(humedad = 50, temperatura = 76, presion =29.30)
prediccion <- predict(modelo, newdata = nuevo.dato)
paste("La cantidad estimada de óxido nitroso emitido es:", round(prediccion, 2))
## [1] "La cantidad estimada de óxido nitroso emitido es: 0.94"
En este caso que es un conjunto de datos el cual demuestra la cantidad de oxido nitroso que se puede presentar dentro de los tanques de camiones con motores de diesel. Para este análisis nuestras variables dependientes es el óxido nitroso y las variables independientes son humedad, temperatura y presión. Este modelo de regresión lineal múltiple como podemos ver en la dispersión de los datos en algunas variables si se llega a presentar una relación lineal sin embargo algunas de ellas son de forma negativa por lo que en algunos casos las predicciones no serán acertadas o cercanas al verdadero valor. En cuanto al modelo podemos ver que el R-squared es de 0.8005 (80%), Adjusted R-squared es de 0.763 (76%) lo que quiere decir que la relación entre la variable dependiente y las variables independientes es muy buena, sin embargo el modelo aún necesita mucho reforzamiento ya que los coeficientes del modelo no tienen un valor estadístico importante.
Se cree que la energía eléctrica que una planta química consume cada mes se relaciona con:
y <- c(240, 236, 290, 274, 301, 316, 300, 296, 267, 276, 288, 261)
x1 <- c(25, 31, 45, 60, 65, 72, 80, 84, 75, 60, 50, 38)
x2 <- c(24, 21, 24, 25, 25, 26, 25, 25, 24, 25, 25, 23)
x3 <- c(91, 90, 88, 87, 91, 94, 87, 86, 88, 91, 90, 89)
x4 <- c(100, 95, 110, 88, 94, 99, 97, 96, 110, 105, 100, 98)
datos <- data.frame(y, x1, x2, x3, x4)
kable(datos, caption = "Aspectos que se relacionan con el consumo de energía eléctrica en una plata química")
| y | x1 | x2 | x3 | x4 |
|---|---|---|---|---|
| 240 | 25 | 24 | 91 | 100 |
| 236 | 31 | 21 | 90 | 95 |
| 290 | 45 | 24 | 88 | 110 |
| 274 | 60 | 25 | 87 | 88 |
| 301 | 65 | 25 | 91 | 94 |
| 316 | 72 | 26 | 94 | 99 |
| 300 | 80 | 25 | 87 | 97 |
| 296 | 84 | 25 | 86 | 96 |
| 267 | 75 | 24 | 88 | 110 |
| 276 | 60 | 25 | 91 | 105 |
| 288 | 50 | 25 | 90 | 100 |
| 261 | 38 | 23 | 89 | 98 |
g1 <- ggplot(data = datos, mapping = aes(x = x1, y = y)) +
geom_point(color = "orange", size = 2) +
labs(title = 'consumo ~ tempertura', x = 'temperatura') +
geom_smooth(method = "lm", se = FALSE, color = "black") +
theme_bw() +
theme(plot.title = element_text(hjust = 0.5))
g2 <- ggplot(data = datos, mapping = aes(x = x2, y = y)) +
geom_point(color = "forestgreen", size = 2) +
labs(title = 'consumo ~ dias de mes', x = 'dias de mes') +
geom_smooth(method = "lm", se = FALSE, color = "black") +
theme_bw() +
theme(plot.title = element_text(hjust = 0.5))
g3 <- ggplot(data = datos, mapping = aes(x = x3, y = y)) +
geom_point(color = "purple", size = 2) +
labs(title = 'consumo ~ purezas del producto', x = 'purezas del producto') +
geom_smooth(method = "lm", se = FALSE, color = "black") +
theme_bw() +
theme(plot.title = element_text(hjust = 0.5))
g4 <- ggplot(data = datos, mapping = aes(x = x4, y = y)) +
geom_point(color = "red", size = 2) +
labs(title = 'consumo ~ produccion toneladas', x = 'produccion toneldas') +
geom_smooth(method = "lm", se = FALSE, color = "black") +
theme_bw() +
theme(plot.title = element_text(hjust = 0.5))
g1 + g2 + g3 + g4
Generando el modelo de regresión lineal múltiple, el consumo de energía eléctrica \(y\) en función ~ de las cuatro variables \(x1\),\(x2\),\(x3\),\(x4\)
modelo <- lm(formula = y ~ ., data = datos)
summary(modelo)
##
## Call:
## lm(formula = y ~ ., data = datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -18.758 -9.952 3.350 6.627 23.311
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -102.71324 207.85885 -0.494 0.636
## x1 0.60537 0.36890 1.641 0.145
## x2 8.92364 5.30052 1.684 0.136
## x3 1.43746 2.39162 0.601 0.567
## x4 0.01361 0.73382 0.019 0.986
##
## Residual standard error: 15.58 on 7 degrees of freedom
## Multiple R-squared: 0.7447, Adjusted R-squared: 0.5989
## F-statistic: 5.106 on 4 and 7 DF, p-value: 0.0303
b0 = modelo$coefficients[1]
b1 = modelo$coefficients[2]
b2 = modelo$coefficients[3]
b3 = modelo$coefficients[4]
b4 = modelo$coefficients[5]
\[yˆ=β0+β1(x1)+β2(x2)+β3(x3)+β4(x4)\] \[yˆ=-102.7132364+0.6053705(x1)+8.9236442(x2)+1.4374567(x3)+0.0136093(x3)\]
Para un mes en que \(x1=75˚F, x2=24\) días, \(x3=90\) y \(x4=98\) toneladas. ¿Cúal es la predicción de consumo de energía eléctrica? \[yˆ=–102.7132+0.6054(75)+8.9236(24)+1.4374(90)+0.0136(98)=287.56\]
nuevo.dato <- data.frame(x1 = 75, x2 = 24, x3 = 90, x4 = 98)
prediccion <- predict(modelo, newdata = nuevo.dato)
paste("La predicción de consumo de energía eléctrica es:", round(prediccion, 2))
## [1] "La predicción de consumo de energía eléctrica es: 287.56"
En este caso el conjunto de datos es de consumo de energía eléctrica en donde la variable dependiente es el consumo de energía eléctrica y las variables independientes son la temperatura ambiental promedio, el número de los días del mes, la pureza promedio del producto y las toneladas fabricadas del producto. Como podemos ver en los diagramas de dispersión en casi todas las variables se puede llegar a apreciar una regresión lineal positiva que casi cruza por algunos de los puntos sin embargo los valores del modelo R-squared vale 0.7447 (74%), Adjusted R-squared vale 0.5989 (59%) este siendo el ultimo que se toma en cuenta para casos de la “vida real” es muy bajo lo que quiere indicar que no hay una buena relación entre la variable dependiente y las variables independientes e igualmente en los coeficientes se puede observar que no hay relación estadística por lo cual este modelo no es muy apto para este caso.