Overview

Vamos a investigar el rendimiento del combustible en coches, medido en millas por galon, para ello vamos a usar el conocido dataset publico mtcars. Pero nos vamos a centrar en averiguarlo para coches de cambio automatico comparado con coches de cambio manual. Tambien vamos a evaluar el impacto de introducir todas o algunas covariables en el ajuste lineal para acercarnos holisticamente a la realidad del conjunto de datos. Las conclusiones las podremos ir descubriendo a lo largo de este informe, pero al final habrá un executive summary.

Exploratory Analysis

datos <- mtcars
str(datos)
## 'data.frame':    32 obs. of  11 variables:
##  $ mpg : num  21 21 22.8 21.4 18.7 18.1 14.3 24.4 22.8 19.2 ...
##  $ cyl : num  6 6 4 6 8 6 8 4 4 6 ...
##  $ disp: num  160 160 108 258 360 ...
##  $ hp  : num  110 110 93 110 175 105 245 62 95 123 ...
##  $ drat: num  3.9 3.9 3.85 3.08 3.15 2.76 3.21 3.69 3.92 3.92 ...
##  $ wt  : num  2.62 2.88 2.32 3.21 3.44 ...
##  $ qsec: num  16.5 17 18.6 19.4 17 ...
##  $ vs  : num  0 0 1 1 0 1 0 1 1 1 ...
##  $ am  : num  1 1 1 0 0 0 0 0 0 0 ...
##  $ gear: num  4 4 4 3 3 3 3 4 4 4 ...
##  $ carb: num  4 4 1 1 2 1 4 2 2 4 ...
head(datos)
##                    mpg cyl disp  hp drat    wt  qsec vs am gear carb
## Mazda RX4         21.0   6  160 110 3.90 2.620 16.46  0  1    4    4
## Mazda RX4 Wag     21.0   6  160 110 3.90 2.875 17.02  0  1    4    4
## Datsun 710        22.8   4  108  93 3.85 2.320 18.61  1  1    4    1
## Hornet 4 Drive    21.4   6  258 110 3.08 3.215 19.44  1  0    3    1
## Hornet Sportabout 18.7   8  360 175 3.15 3.440 17.02  0  0    3    2
## Valiant           18.1   6  225 105 2.76 3.460 20.22  1  0    3    1

En un vistazo rapido podemos observar que cada observación pertenece a un modelo de coche, no repitiendose mas observaciones para el mismo modelo.

Segun la documentación del dataset mtcars podemos observar que las covariables que nos interesan para este estudio, inicialmente, son:

Procedo a convertir la covariable am en factor, y cambio la descripcion de sus levels para mejor entendimiento.

datos$am <- factor(datos$am)
levels(datos$am) <- c("automatic","manual")
head(table(datos$mpg, datos$am))
##       
##        automatic manual
##   10.4         2      0
##   13.3         1      0
##   14.3         1      0
##   14.7         1      0
##   15           0      1
##   15.2         2      0

Los valores de millas por galon mas bajos vemos que pertenecen a los coches automaticos.

tail(table(datos$mpg, datos$am))
##       
##        automatic manual
##   24.4         1      0
##   26           0      1
##   27.3         0      1
##   30.4         0      2
##   32.4         0      1
##   33.9         0      1

Los valores de millas por galon mas altos vemos, a su vez, que pertenecen a los coches de cambio manual.

Si observamos el Anexo - Grafico 1, podemos comprobar claramente que los coches de cambio manual obtienen mejor rendimiento en distancia recorrida por unidad de combustible que los coches automaticos, in average.

Hypothesis Tests

Vamos a corroborrar esta información formulando un Hypothesis Test comparando los datos de los coches automaticos con los de los coches manuales.

La Null Hypothesis va a ser que la diferencia entre tipos de transmission es cero.

Y la Hypothesis alternativa va a ser que la diferencia entre tipos de trasmission es distinta de cero, o lo que es lo mismo, que uno de los tipos de coche es mas efectivo que el otro en cuanto a la distancia recorrida por unidad de combustible.

datos.automatic <- subset(datos, as.character(am) == "automatic")
datos.manual <- subset(datos, as.character(am) == "manual")
datos.ht <- t.test(datos.manual$mpg, datos.automatic$mpg, paired = FALSE, var.equal = FALSE)
datos.ht$conf.int
## [1]  3.209684 11.280194
## attr(,"conf.level")
## [1] 0.95
datos.ht$p.value
## [1] 0.001373638

Viendo que el p-value es menor que 0.05 (5%) podemos considerar rechazar la null hypothesis (aceptando la Hypotesis alternativa).

Como ya sabiamos podemos asegurar que el rendimiento en cuanto a distancia recorrida para coches manuales es diferente al rendimiento de los coches automaticos. Y hasta nos podemos aventurar a decir que los coches manuales recorren mucha mas distancia por unidad de combustible que los coches automaticos.

Obviamente existen otras covariables que pueden intervenir variando y seguramente equilibrando esta diferencia, como el peso, los cilindros, los caballos de potencia, etc…

Simple Linear Regression

Vamos a usar solo una covariable am para el ajuste de regresión lineal.

fit.slm <- lm(mpg ~ am, data = datos)
summary(fit.slm)
## 
## Call:
## lm(formula = mpg ~ am, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -9.3923 -3.0923 -0.2974  3.2439  9.5077 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   17.147      1.125  15.247 1.13e-15 ***
## ammanual       7.245      1.764   4.106 0.000285 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.902 on 30 degrees of freedom
## Multiple R-squared:  0.3598, Adjusted R-squared:  0.3385 
## F-statistic: 16.86 on 1 and 30 DF,  p-value: 0.000285

Podemos sacar 2 conclusiones interesantes:

Multivariate Linear Regression

Vamos a incluir todas las covariables de la muestra de datos, para poder estudiar el impacto que producen en el modelo de regresión lineal. Algunas de ellas procedo a convertirlas en factores previamente.

datos$cyl<-relevel(as.factor(datos$cyl), ref='8')
datos$vs<-relevel(as.factor(datos$vs), ref='0')
datos$gear<-relevel(as.factor(datos$gear), ref='3')
datos$carb<-relevel(as.factor(datos$carb), ref='2')
fit.mlm <- lm(mpg ~ ., data = datos)
summary(fit.mlm)
## 
## Call:
## lm(formula = mpg ~ ., data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.5087 -1.3584 -0.0948  0.7745  4.6251 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept) 22.56362   17.74315   1.272   0.2229  
## cyl4         0.33616    7.15954   0.047   0.9632  
## cyl6        -2.31253    5.49944  -0.421   0.6801  
## disp         0.03555    0.03190   1.114   0.2827  
## hp          -0.07051    0.03943  -1.788   0.0939 .
## drat         1.18283    2.48348   0.476   0.6407  
## wt          -4.52978    2.53875  -1.784   0.0946 .
## qsec         0.36784    0.93540   0.393   0.6997  
## vs1          1.93085    2.87126   0.672   0.5115  
## ammanual     1.21212    3.21355   0.377   0.7113  
## gear4        1.11435    3.79952   0.293   0.7733  
## gear5        2.52840    3.73636   0.677   0.5089  
## carb1        0.97935    2.31797   0.423   0.6787  
## carb3        3.97899    3.85764   1.031   0.3187  
## carb4        2.07078    3.49316   0.593   0.5621  
## carb6        5.45692    5.82857   0.936   0.3640  
## carb8        8.22977    7.84852   1.049   0.3110  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.833 on 15 degrees of freedom
## Multiple R-squared:  0.8931, Adjusted R-squared:  0.779 
## F-statistic:  7.83 on 16 and 15 DF,  p-value: 0.000124
datos.origen <- mtcars
cor(datos.origen)[1,]
##        mpg        cyl       disp         hp       drat         wt 
##  1.0000000 -0.8521620 -0.8475514 -0.7761684  0.6811719 -0.8676594 
##       qsec         vs         am       gear       carb 
##  0.4186840  0.6640389  0.5998324  0.4802848 -0.5509251

Como sospechabamos incluir todas las covariables no es buena idea, podemos observar como el 7.245 anterior se convierte en 1.212, un valor que no es representativo, ni correctamente inferido. Lo que si hemos conseguido es ampliar el espectro de variabilidad de los datos, llegando a un 89%.

Respecto a la correlación, podemos observar una fuerte relación con las covariables: cyl, disp, hp y wt. Estas 4 covariables junto con am serian las apropiadas para construir nuestro propio modelo de regression lineal.

Our Own Multivariable Linear Regression

Con la información y las averiguaciones de la que disponemos en este momento, voy a construir nuestro modelo de regresion lineal.

fit.jclm <- lm(mpg ~ am + cyl + disp + hp + wt, data = datos)
summary(fit.jclm)
## 
## Call:
## lm(formula = mpg ~ am + cyl + disp + hp + wt, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.9374 -1.3347 -0.3903  1.1910  5.0757 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 31.146495   4.144394   7.515  7.2e-08 ***
## ammanual     1.806099   1.421079   1.271   0.2155    
## cyl4         2.717781   2.898149   0.938   0.3573    
## cyl6        -0.418285   2.243037  -0.186   0.8536    
## disp         0.004088   0.012767   0.320   0.7515    
## hp          -0.032480   0.013983  -2.323   0.0286 *  
## wt          -2.738695   1.175978  -2.329   0.0282 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.453 on 25 degrees of freedom
## Multiple R-squared:  0.8664, Adjusted R-squared:  0.8344 
## F-statistic: 27.03 on 6 and 25 DF,  p-value: 8.861e-10

Mi conclusión principal es que puedo perfectamente prescindir de la covariable disp al tener una variación, por aumento de una unidad, de +0.004 mpg, cantidad casi imperceptible.

Mi apuesta final es incluir las covariables cyl, hp y wt (a parte de am) en el modelo de regresion lineal. Ahora vamos a ver si esto coincide con los calculos automaticos de la funcion step de R.

Optimized Multivariable Linear Regression

Llegado este punto, vamos a hacer un calculo del mejor ajuste lineal posible.

fit.blm <- step(lm(mpg ~ ., data = datos), trace = 0)
summary(fit.blm)
## 
## Call:
## lm(formula = mpg ~ cyl + hp + wt + am, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.9387 -1.2560 -0.4013  1.1253  5.0513 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 31.54465    3.88461   8.120 1.34e-08 ***
## cyl4         2.16368    2.28425   0.947  0.35225    
## cyl6        -0.86767    1.71921  -0.505  0.61803    
## hp          -0.03211    0.01369  -2.345  0.02693 *  
## wt          -2.49683    0.88559  -2.819  0.00908 ** 
## ammanual     1.80921    1.39630   1.296  0.20646    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.41 on 26 degrees of freedom
## Multiple R-squared:  0.8659, Adjusted R-squared:  0.8401 
## F-statistic: 33.57 on 5 and 26 DF,  p-value: 1.506e-10

Incluyendo las covariables cyl (num of cylinders), hp (horse power) y wt (weight) el modelo se ajusta mucho mejor. Podemos ver que el 87% de la variabilidad esta explicada en el calculo, cantidad mas que aceptable.

Podemos ver que los coches de cambio manual tienen un 1.81 mpg mas comparados con los de cambio automatico, dejando las demas covariables constantes.

Podemos sacar otras conclusiones, cuando aumentamos en 1000 lb el weight, vemos que disminuye la distancia recorrida en 2.5 mpg, manteniendo las otras covariables constantes.

Asimismo, los coches con 4 cilindros tienen un aumento de distancia recorrida por galon de 2.2 mpg que los de 8 cilindros, manteniendo las demas covariables constantes.

Solo añadir que nuestras conclusiones y las calculadas por el comando step de R coinciden. Con lo cual solo queda comparar los modelos generados. Para mayor claridad en la comparación solo vamos a hacerlo con el mas simple y con el obtimizado. Cuanto mas diferentes sean, podemos asumir que el modelo optimizado es mucho mejor.

Simple Linear Regression vs Multivariate Linear Regression

Vamos a demostrar que el modelo optimizado multivariable es muy diferente al Modelo Simple de Regresión lineal. Con esto asumimos que esta mejor ajustado y es mas proximo a la realidad subyacente de los datos analizados.

anova(fit.slm, fit.blm)
## Analysis of Variance Table
## 
## Model 1: mpg ~ am
## Model 2: mpg ~ cyl + hp + wt + am
##   Res.Df    RSS Df Sum of Sq      F    Pr(>F)    
## 1     30 720.90                                  
## 2     26 151.03  4    569.87 24.527 1.688e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Si consideramos, en este contexto, que la null hypothesis es que los 2 modelos son iguales, y la alternativa que son diferentes. Podemos rechazar la null hypothesis porque el p-value resultante (1.688e-08 ~ 0) es extremadamente inferior al 5% (0.05). Con lo cual concluimos que estos dos modelos son muy diferentes.

En el Anexo - Grafico 2 podemos observar sus residuos, que a todas luces parecen no presentar ninguna anomalia.

Executive Summary

Podemos inferir finalmente que bajo cualquier circunstancia con los datos de la muestra, los coches de cambio manual obtienen mayor rendimiento del combustible utilizado que los coches de cambio automatico, es decir recorren mas millas por gallon consumido.

Añadiendo las covariables: cyl (num of cylinders), hp (horse power) y wt (weight) en el modelo simple de regresion lineal (mpg ~ am), podemos inferir que los coches de cambio manual recorren 1.81 millas por gallon (de media) mas que los coches de cambio automatico.

Este valor hay que tomarlo en cuenta pensando en la media de la muestra de los datos analizados.

Anexo

Grafico 1 - Distancia recorrida en mpg comparando el tipo de cambio del coche, automatico vs manual.

g <- ggplot(datos, aes(x = am, y = mpg))
g <- g + geom_boxplot(aes(fill = am))
g <- g + labs(title = "Distancia recorrida (mpg) coches automaticos vs manuales")
g <- g + labs(x = "Transmission", y = "Distancia (mpg)")
g

NOTA: de una forma aislada (tomando solo en cuenta 2 covariables) podemos ver que los coches de cambio manual hacen mejor uso del combustible (gallons) que los coches de cambio automatico, lo que es lo mismo, recorren mas distancia (millas).

Grafico 2 - Representación de los Residuos

par(mfrow = c(2,2))
plot(fit.blm)

NOTA: Los 4 graficos paracen correctos, no se ve nada fuera de lo normal. The QQ-plot se ajusta a la linea casi perfectamente y en los demas no se observa nada anomalo. No veo necesario hacer mas diagnosticos.