Examen Final
Mariangel
2025-05-24
Carga de los Datos
Datos <- read_excel("C:/Users/User/Downloads/Reg_1.xlsx")
View(Datos)Descriptivo
Antes que todo, procederemos a hacer un análisis descriptivo de los datos para poder tener un adelanto de como se comportarán. Por lo tanto, procedemos primero a realizar un histograma y un summary para ver la distribución y cuales son sus medidas de tendencia central.
\(y:\): Ventas
\(x_{1}:\) Publicidad en Televesión.
\(x_{2}:\) Publicidad en Radio.
\(x_{3}:\) Publicidad en Periódicos.
Primero, veamos la variable respuesta y:
hist(Datos$Ventas,
main = "Histograma de las Ventas",
xlab = "Ventas",
col = "skyblue",
border = "white")
De acuerdo al histograma obtenido, la distribución de la variable Ventas presenta, a simple vista, una ligera asimetría hacia la derecha, ya que la mayor parte de los valores se concentran en los intervalos bajos. Por otro lado, los intervalos con mayor frecuencia fueron los comprendidos entre 8-10 y 10-12, con una frecuencia de 5 observaciones cada uno. Mientras que, los intervalos con menor frecuencia fueron los comprendidos entre 16-18 y 20-22, con una frecuencia de 1 observación. Esto nos da a entender que, la mayoría de las ventas tienden a ubicarse en los valores más bajos del rango observado.
summary(Datos$Ventas)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 8.00 10.75 12.50 13.43 16.00 22.00IQR(Datos$Ventas)## [1] 5.25Las ventas oscilan entre las 8 y 22 unidades aproximadamente. Por lo observado en el histograma, los datos parecen presentar una ligera asímetria por lo que resulta apropiado mencionar a la mediana. Las ventas tienen una mediana de 12.5 unidades, es decir, el 50% de las observaciones se encuentran por debajo o igual a este valor, con un rango intercuartílico de 5.25 unidades, lo que indica que el 50% central de las ventas esta entre 10.75 y 16 unidades. Por último, el 25% de las ventas está por debajo de las 10.75 unidades, mientras que el 75% esta por debajo de las 16 unidades.
En segundo lugar, veamos la variable explicativa \(x_1\) (Publicidad en Televesión):
hist(Datos$Publ_TV,
main = "Histograma de la publicidad en televesión",
xlab = "Unidades",
col = "skyblue",
border = "white")
De acuerdo al histograma obtenido, la distribución de la variable Publ_TV no parece presentar, a simple vista, una asimetría clara. Por otro lado, el intervalo con mayor frecuencia fue el comprendido entre 0.8-1.0 , con una frecuencia de 6 observaciones. Mientras que, el intervalo con menor frecuencia fue el comprendido entre 1.0-1.2, con una frecuencia de 1 observación. En general, los valores se encuentran relativamente dispersos entre los distintos intervalos, lo que sugiere que no hay una fuerte concentración en torno a un intervalo especifíco.
summary(Datos$Publ_TV)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.800 1.000 1.425 1.413 1.762 2.000sd(Datos$Publ_TV)## [1] 0.413227La publicidad en televesión oscila entre las 0.8 y 2.0 unidades aproximadamente. Por lo observado en el histograma, los datos no parecen presentar asimetría por lo que resulta apropiado mencionar a la media. El valor promedio de la publicidad en televesión es de 1.413 unidades, lo que sugiere que, se invierte alrededor de esa cantidad, con una desviación estándar de 0.41 unidades, lo que indica que los valores suelen fluctuar entre 1.41 ± 0.41. Por último, el 25% de la publicidad en televesión está por debajo de 1 unidade, mientras que el 75% esta por debajo de las 1.7 unidades.
Luego, veamos la variable explicativa \(x_2\) (Publicidad en Radio):
hist(Datos$Publ_radio,
main = "Histograma de la publicidad en radio",
xlab = "Unidades",
col = "skyblue",
border = "white")
De acuerdo al histograma obtenido, la distribución de la variable Publ_radio presenta, a simple vista, una asimetría hacia la derecha, ya que la mayor parte de los valores se concentran en los intervalos bajos. Por otro lado, el intervalo con mayor frecuencia fue el comprendido entre 60-70, con una frecuencia de 8 observaciones. Mientras que, los intervalos con menor frecuencia fueron los comprendidos entre 90-100 y 100-110, con una frecuencia de 1 observación. Esto nos da a entender que, la mayoría de las ventas tienden a ubicarse en los valores más bajos del rango observado.
summary(Datos$Publ_radio)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 50.00 60.00 67.00 70.60 76.25 110.00IQR(Datos$Publ_radio)## [1] 16.25La publicidad en radio oscila entre las 50 y 110 unidades aproximadamente. Por lo observado en el histograma, los datos parecen presentar una ligera asímetria por lo que resulta apropiado mencionar a la mediana. La publicidad en radio tienen una mediana de 67 unidades, es decir, el 50% de las observaciones se encuentran por debajo o igual a este valor, con un rango intercuartílico de 16.25 unidades, lo que indica que el 50% central de las ventas esta entre 60 y 76.25 unidades. Por último, el 25% de la publicidad en televesión está por debajo de las 60 unidades, mientras que el 75% esta por debajo de las 76.25 unidades.
Por último, veamos la variable explicativa \(x_3\) (Publicidad en periódicos):
hist(Datos$Publ_periodicos,
main = "Histograma de la publicidad en periódicos",
xlab = "Unidades",
col = "skyblue",
border = "white")
De acuerdo al histograma obtenido, la distribución de la variable Publ_periodicos presenta, a simple vista, una asimetría hacia la derecha, ya que la mayor parte de los valores se concentran en los intervalos bajos. Por otro lado, el intervalo con mayor frecuencia fue el comprendido entre 0.4-0.5, con una frecuencia de 9 observaciones aproximadamente. Mientras que, los intervalos con menor frecuencia fueron los comprendidos entre 0.8-0.9 y 1.0-1.1 aproximadamente, con una frecuencia de 2 observaciones cada uno. Esto nos da a entender que, la mayoría de inversions en la publicidad en periódicos tienden a ubicarse en los valores más bajos del rango observado.
summary(Datos$Publ_periodicos)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 0.300 0.400 0.450 0.561 0.575 1.100IQR(Datos$Publ_periodicos)## [1] 0.175La publicidad en periódicos oscila entre las 0.3 y 1.10 unidades aproximadamente. Por lo observado en el histograma, los datos parecen presentar una ligera asímetria por lo que resulta apropiado mencionar a la mediana. La publicidad en periódicos tienen una mediana de 0.450 unidades, es decir, el 50% de las observaciones se encuentran por debajo o igual a este valor, con un rango intercuartílico de 0.17 unidades, lo que indica que el 50% central de la publicidad en televesión esta entre 0.40 y 0.57 unidades. Por último, el 25% de la publicidad en periódicos está por debajo de las 0.40 unidades, mientras que el 75% esta por debajo de las 0.57 unidades. En general, no parece que se invierta mucho en publicidad en periódicos, ya que el rango de valores observados es bastante reducido.
NORMALIDAD
Veamos si los datos siguen una distribución normal para saber que método aplicar:
shapiro.test(Datos$Ventas)##
## Shapiro-Wilk normality test
##
## data: Datos$Ventas
## W = 0.94076, p-value = 0.2478shapiro.test(Datos$Publ_TV)##
## Shapiro-Wilk normality test
##
## data: Datos$Publ_TV
## W = 0.92377, p-value = 0.1171shapiro.test(Datos$Publ_radio)##
## Shapiro-Wilk normality test
##
## data: Datos$Publ_radio
## W = 0.9034, p-value = 0.04779shapiro.test(Datos$Publ_periodicos)##
## Shapiro-Wilk normality test
##
## data: Datos$Publ_periodicos
## W = 0.73628, p-value = 0.0001118Como hay dos p-valores menores al nivel de significancia establecido de 0.05, rechazamos la hipótesis nula, es decir, los datos no siguen una distribución aproximadamente normal.
CORRELACIÓN
Veamos que tan correlacionadas están las variables entre ellas:
chart.Correlation(Datos, histogram = TRUE, method = "spearman")
mat_cor <- cor(Datos, method = "spearman") # Calcula matriz de correlación
mat_cor1<-cor(Datos, method = "spearman")
significancia1<- cor.mtest(Datos,
method = "spearman",
conf.level = .95)
corrplot(mat_cor1,
p.mat = significancia1$p, #llamado del p-valor para cada coeficiente r
sig.level = 0.05) #definición del nivel de significancia
En cuanto a la distribución, todas parecen seguir una distribución normal, excepto por la publicidad en televesión y la publicidad en peiódicos en la cual se observa una ligera asimetría a la derecha, lo que explica que no hayan pasado la prueba de normalidad. A continuación, se realiza un análisis exhaustivo de la relación de las variables explicativas en torno a la variable respuesta.
Ventas vs Publ_TV: De acuerdo al coeficiente de correlación obtenido (\(R = 0.69\)), remarcado con altos niveles de significancia, se puede afirmar que existe una correlación positiva moderada entre la variable publicidad en televesión y las ventas. Esto implica que, en general, a medida que aumenta la inversión en publicidad televisiva, también tienden a aumentar las ventas.
Ventas vs Publ_radio: De acuerdo al coeficiente de correlación obtenido (\(R = 0.43\)), se puede afirmar que existe una correlación positiva baja entre la variable publicidad en radio y las ventas. Esto implica que, en general, a medida que aumenta la inversión en publicidad radial, también tienden a aumentar las ventas levemente.
Ventas vs Publ_periodicos: De acuerdo al coeficiente de correlación obtenido (\(R = 0.55\)), se puede afirmar que existe una correlación positiva moderada entre la variable publicidad en periódicos y las ventas. Esto implica que, en general, a medida que aumenta la inversión en publicidad radial, también tienden a aumentar las ventas.
Por otro lado, es importante resaltar que las variables explicativas no presentan una correlación fuerte entre sí. En particular, la correlación entre publicidad en televisión y radio fue de \(r = 0.50\), lo que indica una correlación positiva moderada. Por su parte, la correlación entre publicidad en periódicos y televisión fue de \(r = 0.31\), y entre periódicos y radio fue de \(r = 0.35\), lo que en ambos casos sugiere una correlación positiva débil. Por lo tanto, estas correlaciones relativamente “bajas” favorecen nuestro análisis.
Selección de variables
Para esto, aplicaremos el metodo híbrido o both ya que es el más compacto y nos llevará a mejores resultados:
Método Forward
A continuación aplicaremos el método forward, utilizando de criterio el AIC:
Primero iniciemos creando el modelo nulo:
M_0 <- lm(Datos$Ventas ~ 1) #De este manera se le indica a R que es el modelo nulo.
summary(M_0)##
## Call:
## lm(formula = Datos$Ventas ~ 1)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.425 -2.675 -0.925 2.575 8.575
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 13.425 0.885 15.17 4.52e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.958 on 19 degrees of freedomAIC(M_0)## [1] 114.7604Procedamos a crear modelos con las tres variables explicativas:
#Creamos cada uno de los modelos.
P_1 <- lm(Datos$Ventas ~ Datos$Publ_TV)
P_2 <- lm(Datos$Ventas ~ Datos$Publ_radio)
P_3 <- lm(Datos$Ventas ~ Datos$Publ_periodicos)
#Hacemos summary y calculamos su AIC
summary(P_1)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_TV)
##
## Residuals:
## Min 1Q Median 3Q Max
## -6.1657 -1.3483 -0.2438 0.5682 8.6780
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 4.873 2.569 1.896 0.07406 .
## Datos$Publ_TV 6.055 1.749 3.461 0.00279 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.151 on 18 degrees of freedom
## Multiple R-squared: 0.3996, Adjusted R-squared: 0.3662
## F-statistic: 11.98 on 1 and 18 DF, p-value: 0.002786AIC(P_1)## [1] 106.5571summary(P_2)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_radio)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.9881 -2.6192 -0.2503 2.2935 6.6915
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 6.57076 4.12869 1.591 0.129
## Datos$Publ_radio 0.09709 0.05724 1.696 0.107
##
## Residual standard error: 3.776 on 18 degrees of freedom
## Multiple R-squared: 0.1378, Adjusted R-squared: 0.08988
## F-statistic: 2.876 on 1 and 18 DF, p-value: 0.1071AIC(P_2)## [1] 113.7955summary(P_3)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_periodicos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.5977 -1.2775 0.0968 1.3886 4.7584
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.817 1.214 4.794 0.000145 ***
## Datos$Publ_periodicos 13.562 1.987 6.826 2.17e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.147 on 18 degrees of freedom
## Multiple R-squared: 0.7213, Adjusted R-squared: 0.7059
## F-statistic: 46.6 on 1 and 18 DF, p-value: 2.171e-06AIC(P_3)## [1] 91.20466En este caso el modelo con el AIC más pequeño fue el presentado cuando se agrega la variable Publ_periodicos. Comparando con el modelo anterior (En este caso el modelo nulo representado por M_0) obtenemos:
\[AIC_{M_{0}} = 114.76 > 91.20 = AIC_{P_{3}}\]
Como el AIC disminuyo al agregar la variable explicativa Publ_periodicos, procedemos a crear un nuevo modelo agregando la nueva variable explicativa:
M_1 <- lm(Datos$Ventas ~ Datos$Publ_periodicos)
summary(M_1)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_periodicos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.5977 -1.2775 0.0968 1.3886 4.7584
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 5.817 1.214 4.794 0.000145 ***
## Datos$Publ_periodicos 13.562 1.987 6.826 2.17e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.147 on 18 degrees of freedom
## Multiple R-squared: 0.7213, Adjusted R-squared: 0.7059
## F-statistic: 46.6 on 1 and 18 DF, p-value: 2.171e-06AIC(M_1)## [1] 91.20466Volvemos a repetir el mismo de crear modelos temporales, con las variables explicativas restantes, en este caso, 2:
#Creamos cada uno de los modelos.
P_1 <- lm(Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV)
P_2 <- lm(Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_radio)
#Hacemos summary y calculamos su AIC
summary(P_1)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.1554 -0.2380 0.2326 0.7616 3.2403
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.177 1.456 1.495 0.1532
## Datos$Publ_periodicos 11.364 1.717 6.618 4.35e-06 ***
## Datos$Publ_TV 3.450 1.030 3.349 0.0038 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.714 on 17 degrees of freedom
## Multiple R-squared: 0.8321, Adjusted R-squared: 0.8124
## F-statistic: 42.13 on 2 and 17 DF, p-value: 2.586e-07AIC(P_1)## [1] 83.0706summary(P_2)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_radio)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.3068 -1.0718 0.0713 1.0138 4.7133
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 3.96968 2.39753 1.656 0.116
## Datos$Publ_periodicos 12.97484 2.10255 6.171 1.03e-05 ***
## Datos$Publ_radio 0.03083 0.03444 0.895 0.383
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.158 on 17 degrees of freedom
## Multiple R-squared: 0.7339, Adjusted R-squared: 0.7026
## F-statistic: 23.44 on 2 and 17 DF, p-value: 1.297e-05AIC(P_2)## [1] 92.28361En este caso el modelo con el AIC más pequeño fue el presentado cuando se agrega la variable Publ_Publ_TV. Comparando con el modelo anterior (En este caso el modelo representado por M_1) obtenemos:
\[AIC_{M_{1}} = 91.20 > 83.07 = AIC_{P_{1}}\]
Como el AIC disminuyo al agregar la variable explicativa Publ_TV, procedemos a crear un nuevo modelo agregando la nueva variable explicativa en el modelo actual M1, es decir, nos queda:
M_2 <- lm(Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV)
summary(M_2)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.1554 -0.2380 0.2326 0.7616 3.2403
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.177 1.456 1.495 0.1532
## Datos$Publ_periodicos 11.364 1.717 6.618 4.35e-06 ***
## Datos$Publ_TV 3.450 1.030 3.349 0.0038 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.714 on 17 degrees of freedom
## Multiple R-squared: 0.8321, Adjusted R-squared: 0.8124
## F-statistic: 42.13 on 2 and 17 DF, p-value: 2.586e-07AIC(M_2)## [1] 83.0706Volvemos a repetir el mismo de crear modelos temporale, con las variable explicativa restante, obtenemos:
#Creamos cada uno de los modelos.
P_1 <- lm(Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV + Datos$Publ_radio)
#Hacemos summary y calculamos su AIC
summary(P_1)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV +
## Datos$Publ_radio)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.1512 -0.2184 0.2400 0.7606 3.2460
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.107701 2.054898 1.026 0.32029
## Datos$Publ_periodicos 11.347059 1.801580 6.298 1.06e-05 ***
## Datos$Publ_TV 3.431741 1.121366 3.060 0.00748 **
## Datos$Publ_radio 0.001477 0.029781 0.050 0.96106
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.767 on 16 degrees of freedom
## Multiple R-squared: 0.8321, Adjusted R-squared: 0.8007
## F-statistic: 26.44 on 3 and 16 DF, p-value: 1.941e-06AIC(P_1)## [1] 85.06753En este caso, como hay un solo modelo procedamos a comparar el AIC con el modelo anterior (En este caso el modelo representado por M_2) obtenemos:
\[AIC_{M_{2}} = 83.07 < 85.06 = AIC_{P_{1}}\]
Como en esta iteración ningún modelo presenta un AIC menor que el del modelo anterior, se concluye el proceso de selección. Nuestro modelo final queda definido como:
M_F <- lm(Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV)
summary(M_F)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.1554 -0.2380 0.2326 0.7616 3.2403
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.177 1.456 1.495 0.1532
## Datos$Publ_periodicos 11.364 1.717 6.618 4.35e-06 ***
## Datos$Publ_TV 3.450 1.030 3.349 0.0038 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.714 on 17 degrees of freedom
## Multiple R-squared: 0.8321, Adjusted R-squared: 0.8124
## F-statistic: 42.13 on 2 and 17 DF, p-value: 2.586e-07\[y = 2.17 + 11.36x_{1} + 3.45x_{2}\] Todo este proceso se realizó de forma manual, sin embargo, existe una función que permite simplificar significativamente la selección del modelo:
Modelo_Nulo <- lm(Ventas ~ 1, data = Datos) #Utilizamos ~ 1, para decir que solo sea la variable respuesta con el intercepto.
Modelo_Completo <- lm(Ventas ~ ., data = Datos) #La expresión ~ ., representa el modelo con todas las variables explicativas.
modelo_forward <- stepAIC(Modelo_Nulo, scope = list(lower = formula(Modelo_Nulo), upper = formula(Modelo_Completo)), direction = "forward") #Utilizamos stepAIC para que su criterio de selección sea por el AIC.## Start: AIC=56
## Ventas ~ 1
##
## Df Sum of Sq RSS AIC
## + Publ_periodicos 1 214.700 82.938 32.447
## + Publ_TV 1 118.937 178.700 47.800
## + Publ_radio 1 41.009 256.629 55.038
## <none> 297.637 56.003
##
## Step: AIC=32.45
## Ventas ~ Publ_periodicos
##
## Df Sum of Sq RSS AIC
## + Publ_TV 1 32.969 49.968 24.313
## <none> 82.938 32.447
## + Publ_radio 1 3.733 79.205 33.526
##
## Step: AIC=24.31
## Ventas ~ Publ_periodicos + Publ_TV
##
## Df Sum of Sq RSS AIC
## <none> 49.968 24.313
## + Publ_radio 1 0.0076785 49.960 26.310modelo_forward #Lo llamamos.##
## Call:
## lm(formula = Ventas ~ Publ_periodicos + Publ_TV, data = Datos)
##
## Coefficients:
## (Intercept) Publ_periodicos Publ_TV
## 2.177 11.364 3.450Como vemos el resultado coincide con el proceso hecho manualmente. Procedamos a desglosar este modelo:
modelo_forward$anova## Stepwise Model Path
## Analysis of Deviance Table
##
## Initial Model:
## Ventas ~ 1
##
## Final Model:
## Ventas ~ Publ_periodicos + Publ_TV
##
##
## Step Df Deviance Resid. Df Resid. Dev AIC
## 1 19 297.63750 56.00288
## 2 + Publ_periodicos 1 214.69992 18 82.93758 32.44712
## 3 + Publ_TV 1 32.96945 17 49.96813 24.31306En este caso, vemos que al desglosar el modelo con ANOVA, podemos notar varias cosas:
En la primera fila, correspondiente al modelo nulo (sin variables explicativas), la deviancia residual fue de aproximadamente 297.63 y el AIC de 56.
En la segunda fila, al incorporar la variable explicativa Publ_Periódicos, tanto la deviancia como el AIC disminuyeron notablemente, lo que sugiere que esta variable mejoró significativamente el ajuste del modelo.
En la tercera fila, al agregar también la variable Publ_TV, se observa una nueva disminución en la deviancia y el AIC, lo que indica que esta variable también aporta significativamente a explicar la variabilidad de la variable respuesta (Ventas).
Método Both
Se procede a aplicar el método Forward con cada una de las variables:
m1<-lm(Datos$Ventas ~ Datos$Publ_TV)
AIC(m1)## [1] 106.5571m2<-lm(Datos$Ventas ~ Datos$Publ_radio)
AIC(m2)## [1] 113.7955m3<-lm(Datos$Ventas ~ Datos$Publ_periodicos)
AIC(m3)## [1] 91.20466El menor AIC es el de \(91.20\) del modelo 3 con la variable Publ_periodicos, por tanto lo elegimos y aplicamos backward.
m3<-lm(Datos$Ventas ~ Datos$Publ_periodicos) #El modelo actual
AIC(m3)## [1] 91.20466m_0 <- lm(Datos$Ventas ~ 1)
AIC(m_0)## [1] 114.7604Al eliminar a Publ_periodicos no mejora el AIC, por tanto se mantiene en el modelo. Se procede a aplicar el método forward otra vez con el modelo actual:
m1<-lm(Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV)
AIC(m1)## [1] 83.0706m2<-lm(Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_radio)
AIC(m2)## [1] 92.28361El mejor AIC es el de \(83.07\) que corresponde al modelo con las variables Publ_periodicos y Publ_TV. Se verifica con backward si eliminar alguna mejora el AIC.
m3<-lm(Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV) #El modelo actual
AIC(m3)## [1] 83.0706m_01 <- lm(Datos$Ventas ~ Datos$Publ_periodicos)
AIC(m_01)## [1] 91.20466m_02 <- lm(Datos$Ventas ~ Datos$Publ_TV)
AIC(m_02)## [1] 106.5571El AIC empeora en cualquiera de los dos casos, se mantienen las dos variables y se repite el proceso con el método forward.
m1<-lm(Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV + Datos$Publ_radio)
AIC(m1)## [1] 85.06753El AIC no disminuye, por tanto el proceso termina, pues no se puede añadir ninguna variable que reduzca el AIC ni se puede eliminar ninguna variable para mejorar el AIC. Por tanto el modelo final incluye a las variables Publ_TV y Publ_radio.
No obstante, existe una función que condensa todo este proceso, se puede realizar de la siguiente forma:
modelo_inicial<-lm(Datos$Ventas ~ 1)
scope = list(lower = ~1, upper = Datos$Ventas ~ Datos$Publ_TV + Datos$Publ_radio + Datos$Publ_periodicos)
modelo_both <- stepAIC(modelo_inicial, trace=TRUE, direction="both", scope=scope, k=2)## Start: AIC=56
## Datos$Ventas ~ 1
##
## Df Sum of Sq RSS AIC
## + Datos$Publ_periodicos 1 214.700 82.938 32.447
## + Datos$Publ_TV 1 118.937 178.700 47.800
## + Datos$Publ_radio 1 41.009 256.629 55.038
## <none> 297.637 56.003
##
## Step: AIC=32.45
## Datos$Ventas ~ Datos$Publ_periodicos
##
## Df Sum of Sq RSS AIC
## + Datos$Publ_TV 1 32.969 49.968 24.313
## <none> 82.938 32.447
## + Datos$Publ_radio 1 3.733 79.205 33.526
## - Datos$Publ_periodicos 1 214.700 297.637 56.003
##
## Step: AIC=24.31
## Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV
##
## Df Sum of Sq RSS AIC
## <none> 49.968 24.313
## + Datos$Publ_radio 1 0.008 49.960 26.310
## - Datos$Publ_TV 1 32.969 82.938 32.447
## - Datos$Publ_periodicos 1 128.732 178.700 47.800modelo_both$anova ## Stepwise Model Path
## Analysis of Deviance Table
##
## Initial Model:
## Datos$Ventas ~ 1
##
## Final Model:
## Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV
##
##
## Step Df Deviance Resid. Df Resid. Dev AIC
## 1 19 297.63750 56.00288
## 2 + Datos$Publ_periodicos 1 214.69992 18 82.93758 32.44712
## 3 + Datos$Publ_TV 1 32.96945 17 49.96813 24.31306summary(modelo_both)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.1554 -0.2380 0.2326 0.7616 3.2403
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2.177 1.456 1.495 0.1532
## Datos$Publ_periodicos 11.364 1.717 6.618 4.35e-06 ***
## Datos$Publ_TV 3.450 1.030 3.349 0.0038 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.714 on 17 degrees of freedom
## Multiple R-squared: 0.8321, Adjusted R-squared: 0.8124
## F-statistic: 42.13 on 2 and 17 DF, p-value: 2.586e-07Interpretación de los coeficientes
Intercepto (\(2.17\)): Cuando no se invierte en la publicidad en periódicos o en la de televesión se espera que las ventas sean, en promedio, de 2.17 unidades. No es significativo (\(p ≈ 0.1532\)), pues su p-valor es mayor al nivel de significancia de 0.05.
Publ_periodicos (\(11.36\)): Por cada unidad que se invierte para publicidad en periódicos, manteniendo la publicidad en televesión constate, las ventas aumentan aproximadamente 11.36 unidades. Muy significativo (\(p ≈ 4.35e-06\)).
Publ_TV (\(3.45\)): Por cada unidad que se invierte para publicidad en televesión, manteniendo la publicidad en periódicos constate, las ventas aumentan aproximadamente 3.45 unidades. Es significativo (\(p ≈ 0.0038\)).
Calidad del ajuste
\(R^{2} = 0.83\): El modelo explica aproximadamente el 83.21% de la variabilidad de las ventas, lo cual es un buen ajuste. El 83% de la variabilidad en las ventas se debe a la regresión de la publicidad en televesión y en periódicos, mientras que el 17% se debe a otros factores no incluidos en el estudio.
\(R^{2}_{ajustado} = 0.81\) : Considerando el número de predictores, sigue siendo sólido. Además de que no varia mucho entre el coeficiente de determinación, lo que nos indica que las variables explicativas si aportan significativamente al modelo.
Error estándar residual = 1.71: En promedio, los valores predichos se desvían en ±1.71 unidades de los valores reales de las ventas.
Significancia global del modelo F-statistic = 42.13, con \(p ≈ 2.586e-07\): El modelo como un todo es altamente significativo, lo que indica que al menos una variable predictora explica las ventas.
Multicolinealidad
Antes de verificar supuestos, veamos si hay multicolinealidad entre las variables explicativas, para eso aplicamos la siguiente función:
vifs <- vif(M_F)
print(vifs)## Datos$Publ_periodicos Datos$Publ_TV
## 1.171033 1.171033Todos los VIF están por debajo de 5, lo que indica que no hay problemas de multicolinealidad, todas las variables son independientes entre sí y aportan información única.
Ecucaciones normales - Estimación minímos cuadrados
Primero calculamos la base matricial:
y <- Datos$Ventas #Vector que representa la variable respuesta
X <- cbind(
1, #Columna llena de unos para que se pueda escribir el intercepto
Datos$Publ_TV, #Los valores de la publicidad en televesión.
Datos$Publ_periodicos #Los valores de la publicidad en periódicos.
)Luego, procedemos a estimar los coeficientes por el método de mínimos cuadrados:
#t(X) transpuesta de esa matriz, %*% multiplicación, supongo que solve es inversa.
XtX <- t(X) %*% X
Xty <- t(X) %*% y
beta_hat <- solve(XtX) %*% XtyLos mostramos de manera más ordenada:
colnames(beta_hat) <- "Estimación"
rownames(beta_hat) <- c("Intercepto", "Publ_TV", "Publ_periodicos")
beta_hat## Estimación
## Intercepto 2.177303
## Publ_TV 3.449649
## Publ_periodicos 11.363758Comparando con el modelo anterior:
lm(Datos$Ventas ~ Datos$Publ_TV + Datos$Publ_periodicos)##
## Call:
## lm(formula = Datos$Ventas ~ Datos$Publ_TV + Datos$Publ_periodicos)
##
## Coefficients:
## (Intercept) Datos$Publ_TV Datos$Publ_periodicos
## 2.177 3.450 11.364Vemos que dan exactamente los mismos coeficientes. Por lo tanto, el cálculo se hizo de la manera adecuada.
Significancia de los B’s
IC <- confint(M_F)
IC## 2.5 % 97.5 %
## (Intercept) -0.8949249 5.249530
## Datos$Publ_periodicos 7.7409524 14.986563
## Datos$Publ_TV 1.2765196 5.622779El coeficiente estimado para la publicidad en televesión es estadísticamente significativo, ya que su intervalo de confianza no incluye al 0. Una inversión en publicidad para televesión se asocia, en promedio, con un incremento de 1.27 a 5.62 unidades. Por otro lado, el coeficiente estimado para la publicidad en periódicos también es estadísticamente significativos, puesto que tampoco incluye al 0. Es decir, por cada unidad adicional (manteniendo la otra variable independiente constante), las ventas aumentan, en promedio, entre 7.74 y 14.98 unidades.
Verificación de supuestos
Primero verifiquemos los residuos:
Residuos <- M_F$residualsVeamos si son normales:
shapiro.test(Residuos)##
## Shapiro-Wilk normality test
##
## data: Residuos
## W = 0.89609, p-value = 0.03486Como el p-valor es menor a 0.05, los residuos no siguen una distribución aproximadamente normal.
plot(M_F)
Podemos destacar varias cosas de los gráficos generados:
EL gráfico de “Residuals Vs Fitted” no se observa una forma curva o patrón, lo que sugiere que cumple con la suposición de linealidad e independencia de los errores. Sin embargo, la mayoría de datos se concentran en los valores más bajos. Por otro lado, se observa que la mayoría de sus datos están relativamente cerca del 0, sin embargo, las observaciones 7, 14 y 15 aparecen marcados, por lo que podrían ser valores atípicos.
El gráfico de “Q- Q Residuals” se destaca que algunos de los datos se encuentran sobre la línea teórica de la normalidad, y otros están relativamente cerca. Con la prueba shapiro pudimos comprobar que los residuos no siguen una distribución normal, ya que su p-valor fue menor a 0.05 que es lo que se sospecha del gráfico. Sin embargo, otra vez las observaciones 20, 23 y 24 se marcan y se encuentran en los extremos, lo cual puede estar desviando a la normalidad.
El gráfico de “Scale-Location” no se observa ningún patrón a simple vista, lo cual sugiere que cumple con la suposición de linealidad. No obstante, otra vez se destacan las observaciones 7, 14 y 15.
El gráfico de “Leverage” se observa que los puntos 7, 12 y 15 son de alto leverage. Sin embargo, analicemos más de cerca.
Puntos de alto Leverage y puntos influyentes
Primero, hacemos la prueba Bonferroni para identificar las observaciones atípicas que se pueden presentar en el modelo, ya tenemos algunas sospechas de las observaciones 7, 14 y 15.
outlierTest(M_F, cutoff = Inf) # Usa Bonferroni por defecto## rstudent unadjusted p-value Bonferroni p
## 15 -3.1776731 0.0058456 0.11691
## 7 2.2781917 0.0367870 0.73575
## 14 -2.2750562 0.0370140 0.74028
## 12 1.0801942 0.2960700 NA
## 17 -1.0329348 0.3170000 NA
## 9 0.9575068 0.3525600 NA
## 2 0.6794906 0.5065400 NA
## 18 -0.6042644 0.5541400 NA
## 6 0.5809693 0.5693600 NA
## 19 0.4425471 0.6640200 NADe acuerdo a los resultados obtenidos podemos destacar varias cosas:
La observación 15 presenta un valor de rstudent = -3.18, cuyo valor absoluto es mayor que 3, lo que sugiere que es una observación potencialmente atípica. Su p-valor no ajustado es 0.0058, lo cual indica significancia estadística al nivel del 5%. Sin embargo, al aplicar la corrección de Bonferroni, el p-valor ajustado asciende a 0.1169, lo cual ya no es significativo bajo el mismo umbral. A pesar de esto, el valor de rstudent sigue siendo alto, por lo que esta observación debería ser investigada más a fondo, ya que podría influir de forma importante en el modelo.
La observación 7 tiene un rstudent = 2.28, lo cual está dentro del rango usualmente sospechoso (entre 2 y 3). Su p-valor no ajustado es 0.0368, pero después de la corrección de Bonferroni este se incrementa a 0.7357, indicando que no hay evidencia suficiente para considerarlo un outlier significativo. Por lo tanto, se rechaza la hipótesis nula de que esta observación es atípica.
De forma similar, la observación 14 tiene un rstudent = -2.28. En términos absolutos, este valor también se encuentra en el rango dudoso. Su p-valor no ajustado es 0.0370, pero el p-valor ajustado con Bonferroni es 0.7403, lo que también lleva a rechazar la hipótesis nula. En consecuencia, no se considera un outlier estadísticamente significativo tras la corrección.
Ahora, procederemos a calcular cada \(h_{ii}\) para ver el leverage de cada unio:
leverage <- hatvalues(M_F)
plot(leverage, type = "h", main = "Leverage de cada observación",
ylab = "Leverage", xlab = "Observación")
abline(h = 2 * mean(leverage), col = "red", lty = 2) # umbral común
Podemos notar que la observación 12 presenta un leverage alto, ya que sobrepasa la línea del umbral. Es importante destacar que, según la prueba de Bonferroni, la observación 12 no se considera un outlier, ya que su valor p ajustado es mayor a 0.05. Esto indica que, aunque no se desvía significativamente, sí está alejada del centroide del espacio formado por las variables predictoras (Alto Leverage). En consecuencia, podemos clasificarla como una observación con leverage alto pero sin evidencia de ser un outlier, es decir, un “no outlier con leverage alto”. Este tipo de observaciones no necesariamente son problemáticas, siempre que sigan la tendencia general del modelo.
Por otro lado, notamos que la observación 20 presenta un leverage moderadamente bajo, ya que aunque sobrepasa la línea del umbral, no es tan notable como la observación 12. Además, como destacamos anteriormente, la observación 20 tampoco fue considera un outlier, pues su p-valor fue mayor a 0.05. Por lo cual, lo podemos clasificar como una observación con leverage moderamente bajo pero sin evidencia de ser un atípico (según prueba Bonferroni), es decir, un “no outlier con leverage alto”.
Verifiquemos ahora la distancia de Cook:
cooks_d <- cooks.distance(M_F)
plot(cooks_d, type = "h", main = "Distancia de Cook",
ylab = "Cook's Distance", xlab = "Observación")
abline(h = 4 / nrow(Datos), col = "red", lty = 2) # Umbral sugerido
plot(M_F, which=4)
cooks.distance(M_F)## 1 2 3 4 5 6
## 4.518902e-04 1.317799e-02 1.654905e-03 6.555043e-04 1.219012e-04 1.454636e-02
## 7 8 9 10 11 12
## 2.299086e-01 9.301275e-04 2.046781e-02 4.747484e-03 6.989973e-05 4.967809e-01
## 13 14 15 16 17 18
## 7.198508e-04 1.600363e-01 2.623559e-01 1.171736e-04 5.213884e-02 2.682911e-02
## 19 20
## 1.694546e-02 1.409202e-02cooks.distance(M_F)[which.max(cooks.distance(M_F))]## 12
## 0.4967809De acuerdo al gráfico podemos destacar varias cosas:
Recordar que: Que un punto tenga algo Leverage, significa que tiene un valor extremo en las variables predictoras, esto es, su combinación de valores para las variables independientes se encuentra muy alejada del centroide (el punto promedio de todas las observaciones).
La observación 7 aparece destacada como influyente según el gráfico de la distancia de Cook, lo cual resulta particularmente interesante dado que no fue identificada ni como un valor atípico ni como una observación con alto leverage. Esto indica que, a pesar de no presentar un comportamiento extremo en términos de valores individuales o influencia estructural (leverage), su exclusión del modelo sí podría generar cambios en los coeficientes estimados, es decir, influye en la estabilidad del modelo de regresión. Su valor no es excesivamente alto, por lo que, si bien se debe tener en cuenta, no necesariamente implica una distorsión grave en el modelo.
La observación 12 se identifica como influyente, siendo la que presenta la mayor distancia de Cook (0.49), superando el umbral de 0.2. Además, esta observación ya había sido señalada por tener leverage elevado, lo que confirma su impacto en el modelo y sugiere que puede afectar significativamente los coeficientes estimados.
Con la observación 15 ocurre algo similar a la 7. Aunque no fue clasificada como atípica tras aplicar la corrección de Bonferroni, su rstudent mayor a 3 ya la señalaba como una observación que merecía atención. Ahora, al ser identificada como influyente según la distancia de Cook, se refuerza la sospecha de que puede generar cambios importantes en los coeficientes del modelo, por lo que debe ser analizada con más detalle.
Veamos la influencia de cada uno:
influencePlot(M_F, id.method = "identify", main = "Gráfico de influencia",
sub = "Tamaño del punto = Cook's Distance")
## StudRes Hat CookD
## 7 2.2781917 0.1421067 0.22990855
## 12 1.0801942 0.5632812 0.49678091
## 15 -3.1776731 0.1068709 0.26235586
## 20 0.2912994 0.3203721 0.01409202De acuerdo al gráfico podemos notar varias cosas ya previamente mencionadas:
La observación 12 se observa que tiene algo leverage ya que esta bastante alejada horizontalmente hacia la derecha. Tiene un residual moderado, por lo cual, no nos da sospecha de ser un outlier. Destaca por su gran tamaño y color intenso, indicando una alta distancia de Cook. Esto la convierte en un punto altamente influyente que puede estar afectando de manera importante los coeficientes del modelo.
La observación 7 no presenta un leverage elevado, pero sí un residual mayor a 2, lo que inicialmente podría hacer pensar que es un outlier. Sin embargo, la prueba de Bonferroni nos permitió rechazar la hipótesis nula, por lo que no se considera como tal. A pesar de que su tamaño y color en el gráfico no son tan marcados, sigue siendo una observación influyente, aunque en menor medida que la 12.
La observación 15 no tiene un leverage alto, pero presenta un residual mayor a 3, lo que sí la clasifica como un outlier. Aunque su círculo no es muy grande ni de color muy intenso, se considera influyente. Es más influyente que la observación 7, pero no tanto como la 12.
plot(M_F, which = 5)
par(mfrow = c(2, 2))
plot(M_F) # Diagnóstico clásico: residuos, leverage, Cook, QQ
par(mfrow = c(1, 1)) # Restaurar layoutDe acuerdo con el gráfico podemos destacar varias cosas:
La observación 12 presenta un leverage alto, ya que se encuentra bastante desplazada hacia la derecha del gráfico. Su residuo estandarizado no es muy extremo, por lo que no parece ser un outlier en términos de residuo. Sin embargo, al estar cerca o incluso dentro de la curva de Cook’s distance = 0.5, se considera un punto altamente influyente, lo cual sugiere que podría estar afectando de manera significativa los coeficientes del modelo.
La observación 7 tiene un leverage bajo, es decir, no está muy alejada horizontalmente. Sin embargo, muestra un residuo estandarizado mayor a 2, lo que podría levantar sospechas de ser un outlier en términos del valor ajustado. Aun así, según la prueba de Bonferroni, se rechaza la hipótesis nula, por lo que no se clasifica como un outlier significativo. Además, no se encuentra cerca de las curvas de Cook’s distance, por lo tanto, no es muy influyente en el modelo.
La observación 15 también tiene un leverage bajo, pero presenta un residuo estandarizado cercano a -3, lo cual indica que sí es un outlier en cuanto a su valor ajustado. Aunque no se encuentra cerca de las curvas de Cook’s distance, su posición vertical extrema le otorga cierta influencia, más que la observación 7, pero menos que la 12.
Analicemos más de cerca:
summary(influence.measures(M_F))## Potentially influential observations of
## lm(formula = Datos$Ventas ~ Datos$Publ_periodicos + Datos$Publ_TV) :
##
## dfb.1_ dfb.Dt$P_ dfb.D$P_T dffit cov.r cook.d hat
## 12 0.23 1.07_* -0.84 1.23 2.22_* 0.50 0.56_*
## 15 0.13 0.60 -0.72 -1.10 0.31_* 0.26 0.11
## 20 -0.13 0.14 0.05 0.20 1.74_* 0.01 0.32De acuerdo a los resultados obtenidos podemos observar que:
De la observación 12: De acuerdo al dfb.Publ_periodicos se observa que es muy influyente, es decir, que si eliminamos esa observación, el coeficiente podría cambiar significativamente, mostrando mucha influencia en ese coeficiente. Para publicidad en televisión no es influyente. En dffit el valor ajustado cambia considerablemente. En cov.r la varianza cambiaría mucho si se quita ese punto. En cook.d tiene influencia moderada (0.50), mayor al umbral 0.2 que indica alta influencia. Además, tiene un leverage muy alto (0.56_*) resaltado con asterisco. Por lo tanto, es un punto muy influyente.
De la observación 15: De acuerdo al dfb.Publ_periodicos y dfb.Publ_TV se observa que no es muy influyente, pues aunque tiene valores relativamente altos, no son lo suficientemente grandes para que muestren mucha influencia en esos coeficientes. En dffit el valor ajustado cambia de manera considerable. En cov.r la varianza cambiaría mucho si se quita ese punto. En cook.d tiene una influencia no muy alta, pues no es mayor que 0.5 ni 1, pero si supera el umbral 0.2 pero por muy poco, por lo que no es tan influyente. Por último, no tiene un Leverage alto. Por lo que puede ser un outlier, pero no muy influyente.
De la observación 15: De acuerdo al dfb.Publ_periodicos y dfb.Publ_TV se observa que no es muy influyente, pues no tiene valores altos, por lo que no muestran mucha influencia en esos coeficientes. En dffit el valor ajustado cambia muy poco, por lo que tiene poco impacto en la predicción. En cov.r la varianza cambiaría mucho si se quita ese punto. En cook.d tiene una influencia no muy alta, pues no es mayor que 0.5 ni 1, y tampoco supera el umbral 0.2, por lo que no es influyente. Por último, tiene un Leverage moderadamente alto.
En conclusión, con todas estas pruebas que se hicieron la observación 12 fue la más destacada en la distancia de Cook, por lo que es la más influyente y puede cambiar drásticamente los coeficientes del modelo si se considera eliminarla. Por lo que hay que analizarla con cuidado y verificar herramientas para observarla de manera más adecuada. La observación 7 y 15 también vale la pena revisarlas, pero parecen ser más outliers que puntos muy influyentes, pues aunquen aparecen marcadas en la distancia de Cook no parecen influir de manera notable como lo hace la observación 12.
VARIABLES CATÉGORICAS - MRLM
Carga de datos
Meses <- c(2, 6, 8, 3, 2, 7, 9, 8, 4, 6)
Tipo <- c("eléctrico", "mecánico", "eléctrico", "mecánico", "eléctrico", "eléctrico", "mecánico", "mecánico", "eléctrico", "eléctrico")
Horas <- c(2.9, 3.0, 4.8, 1.8, 2.9, 4.9, 4.2, 4.8, 4.4, 4.5)
Datos_Johnson <- data.frame(Meses, Tipo, Horas)
Datos_Johnson## Meses Tipo Horas
## 1 2 eléctrico 2.9
## 2 6 mecánico 3.0
## 3 8 eléctrico 4.8
## 4 3 mecánico 1.8
## 5 2 eléctrico 2.9
## 6 7 eléctrico 4.9
## 7 9 mecánico 4.2
## 8 8 mecánico 4.8
## 9 4 eléctrico 4.4
## 10 6 eléctrico 4.5Creaación variable dummy
Datos_Johnson$Tipo_dummy <- ifelse(Datos_Johnson$Tipo == "eléctrico", 1, 0)
Datos_Johnson## Meses Tipo Horas Tipo_dummy
## 1 2 eléctrico 2.9 1
## 2 6 mecánico 3.0 0
## 3 8 eléctrico 4.8 1
## 4 3 mecánico 1.8 0
## 5 2 eléctrico 2.9 1
## 6 7 eléctrico 4.9 1
## 7 9 mecánico 4.2 0
## 8 8 mecánico 4.8 0
## 9 4 eléctrico 4.4 1
## 10 6 eléctrico 4.5 1Y1 <- Datos_Johnson$Horas
X1 <- cbind(1,
Datos_Johnson$Meses,
Datos_Johnson$Tipo_dummy
)Ahora, hagamos los cálculos con la estimación de minimos cuadrados:
#t(X) transpuesta de esa matriz, %*% multiplicación, supongo que solve es inversa.
XtX1 <- t(X1) %*% X1
Xty1 <- t(X1) %*% Y1
beta_hat1 <- solve(XtX1) %*% Xty1Los mostramos de manera más ordenada:
colnames(beta_hat1) <- "Estimación"
rownames(beta_hat1) <- c("Intercepto", "Meses", "Tipo de reparación")
beta_hat1## Estimación
## Intercepto 0.9304954
## Meses 0.3876161
## Tipo de reparación 1.2626935Ahora, hagamos lo mismo con la función:
modelo_multiple <- lm(Horas ~ Meses + Tipo_dummy, data = Datos_Johnson)
summary(modelo_multiple)##
## Call:
## lm(formula = Horas ~ Meses + Tipo_dummy, data = Datos_Johnson)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.49412 -0.24690 -0.06842 -0.00960 0.76858
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.93050 0.46697 1.993 0.086558 .
## Meses 0.38762 0.06257 6.195 0.000447 ***
## Tipo_dummy 1.26269 0.31413 4.020 0.005062 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.459 on 7 degrees of freedom
## Multiple R-squared: 0.8592, Adjusted R-squared: 0.819
## F-statistic: 21.36 on 2 and 7 DF, p-value: 0.001048Como podemos observar, nos da lo mismo por lo cual, se hicieron los cálculos de manera correcta. El modelo nos queda de la siguiente manera:
\[y = 0.93 + 0.38x_1+1.26 x_2\] - Intercepto (0.93): Una solicitud de servicio mecánico realizada en el mes 0 (es decir, cuando Meses = 0), se espera que tome 0.93 horas en promedio para su reparación. No es significativo ya que tiene un p-valor mayor a 0.05.
x1 (0.38): Por cada mes adicional, se espera que el tiempo promedio para su reparación aumente aproximadamente 0.38 horas, manteniendo el resto de variables constantes. Es altamente significativo.
x2 (0.38): Una solicitud de servicio eléctrico toma, en promedio, 1.26 horas más que un servicio mécanico, manteniendo constante el número de meses. Es significativo.
Coeficientes: El 85% de la variabilidad en el tiempo de reparación (en horas) se explica por el número de meses transcurrido y el tipo de reparación, mientras que el 15% restante se debe a otras variables no incluidas en el estudio. El \(R^{2}_{ajustado} = 0.819\) indica que aproximadamente el 81.9% de la variabilidad en el tiempo de reparación se explica por el modelo, ajustando por el número de variables que incluye. Esto significa que el modelo tiene un buen poder explicativo y que solo queda un 18.1% de variabilidad debido a factores no considerados o error.
Verificar significancia
modelo_reducido <- lm(Horas ~ Meses, data = Datos_Johnson)
modelo_completo <- lm(Horas ~ Meses + Tipo_dummy, data = Datos_Johnson)
anova(modelo_reducido, modelo_completo)## Analysis of Variance Table
##
## Model 1: Horas ~ Meses
## Model 2: Horas ~ Meses + Tipo_dummy
## Res.Df RSS Df Sum of Sq F Pr(>F)
## 1 8 4.8800
## 2 7 1.4751 1 3.4049 16.158 0.005062 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Al comparar el modelo reducido (solo Meses) con el modelo completo (Meses + Tipo_dummy), el valor de p = 0.005 es menor que 0.05, lo que indica que agregar la variable Tipo_dummy mejora significativamente el modelo. Es decir, el tipo de reparación aporta información importante para explicar el tiempo de reparación.