Ejercicio de regresión lineal

Con la intención de comparar el desempeño de dos clases de discos duros (0:SDD, 1: HDD). Este desempeño es medido a través de la variable Y: tiempo de respuesta del disco (segundos), la cual se relaciona, posiblemente bajo una dependencia no lineal, de X: la carga del sistema (Número de consultas por minuto)

i. Relación entre el tiempo de respuesta y la carga de trabajo diferenciando los dos tipos de disco duro.

Los datos que se utilizan en este análisis son de "DatosTaller. A continuación, subimos la base de datos.

library(readxl)
DatosTaller <- read_excel("C:/Users/home/Desktop/Sexto semestre/Electiva R/DatosTaller.xlsx")
View(DatosTaller)
attach(DatosTaller)

Primero vamos a comparar gráficamente el desempeño de los discos SDD y HDD por separado.

Disco 0:SDD

plot(x = SDD$Carga0, y = SDD$Tiempo0, main = "Desempeño disco SDD",
     xlab = "Carga", ylab= "Tiempo")

Disco 1:HDD

plot(x = HDD$Carga1, y = HDD$Tiempo1, main = "Desempeño disco HDD",
     xlab = "Carga", ylab= "Tiempo")

Ahora bien, asignamos como valor numérico nuestras variables para el análisis de correlación:

CARGA0 <- as.numeric(Carga0)
TIEMPO0 <- as.numeric(Tiempo0)
CARGA1 <- as.numeric(Carga1)
TIEMPO1 <- as.numeric(Tiempo1)

Realizamos para ambos casos los test de correlacion de Pearson

cor.test(x= CARGA0, y= TIEMPO0, alternative = "less", method = "pearson" )
## 
##  Pearson's product-moment correlation
## 
## data:  CARGA0 and TIEMPO0
## t = 28.334, df = 10, p-value = 1
## alternative hypothesis: true correlation is less than 0
## 95 percent confidence interval:
##  -1.0000000  0.9979347
## sample estimates:
##       cor 
## 0.9938293
cor.test(x= CARGA1, y= TIEMPO1, alternative = "less", method = "pearson" )
## 
##  Pearson's product-moment correlation
## 
## data:  CARGA1 and TIEMPO1
## t = 12.024, df = 11, p-value = 1
## alternative hypothesis: true correlation is less than 0
## 95 percent confidence interval:
##  -1.0000000  0.9871297
## sample estimates:
##       cor 
## 0.9640003

El análisis por separado de los dos tipos de discos y su relación entre la carga y el tiempo arroja una correlación muy cercana a 1, esto significa que existe una fuerte relacion entre ambos datos.Algo por destacar, es que existe una mayor correlacion entre el tiempo y la carga en el tipo de disco SDD de 0.9938 frente a 0.96400 del disco HDD. Por último,se aplica el método de PEARSON ya que se nota (por medio de las gráficas) una normalidad en los datos.

## The following objects are masked from DatosTaller (pos = 3):
## 
##     Carga, Conf, Tiempo
CARGA <- as.numeric(Carga)
TIEMPO <-as.numeric(Tiempo)
Tipodedisco <- as.character(Conf)
library(ggplot2)

ggplot(data = DatosTaller)+
  geom_line(aes(x=CARGA , y= TIEMPO, col= Tipodedisco))+
  geom_point(aes(x=CARGA , y=TIEMPO, col= Tipodedisco))+
  geom_text(data = NULL,x=8, y=1, label="CorrelaciónSSD = 0.99",col="gray")+
  geom_text(data = NULL, x=8, y=1.5, label="CorrelaciónHDD = 0.96",col="gray")

¿Se puede ver una evidencia lineal?

De acuerdo al análisis empírico, se puede evidenciar una tendencia lineal, es decir, a medida que aumenta la carga, el tiempo también,dicho comportamiento se distingue en HDD Y SDD. Algo por destacar es que SDD en la gráfica conjunta tiene una tendencia exponencial y HDD una especie de parábola hacia abajo. Además, la correlación que evalúa la tendencia entre dos variables, en este caso, entre el tiempo y la carga, arrojó valores muy cercanos a uno, indicando una fuerte relación entre los mismos.

ii. Modelo de regresión simple.

Para ajustar el modelo, se hace mediante una regresón lineal, en dónde nuestra variable regresada es el Tiempo en función de nuestra variable independiente Carga.

\[Tiempo= β_0+β_1*Carga+U_1\]

Modelo1 <- lm(TIEMPO~ CARGA, data= DatosTaller)
summary(Modelo1)
## 
## Call:
## lm(formula = TIEMPO ~ CARGA, data = DatosTaller)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.8902 -0.4350  0.1068  0.4260  1.0536 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.16756    0.29248   0.573    0.575    
## CARGA        0.45320    0.05139   8.819 2.54e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6135 on 15 degrees of freedom
## Multiple R-squared:  0.8383, Adjusted R-squared:  0.8275 
## F-statistic: 77.78 on 1 and 15 DF,  p-value: 2.541e-07

El modelo presenta una bondad de ajuste de 0.8275, lo cual significa que el modelo explica el 82,75% de las variaciones en la variable explicativa, en este caso, el tiempo. De manera general, el modelo presenta un p-value de 2.541e-07 que frente a un p-valor<0.05, se rechaza la Ho de que el modelo NO es significativo en su conjunto.

Algunos supuestos del modelo clásico de regresión lineal

Exogeneidad estricta

Este supuesto indica que no existe correlación entre el término de error y las variables explicativas, por lo que se espera que dicha correlación sea cero.

u=residuals(Modelo1)
cor(u, CARGA)
## [1] 4.595231e-17
library(ggplot2)

plot1 <- ggplot(data = DatosTaller, aes(TIEMPO, u)) +
  geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
  theme_bw()
plot2 <- ggplot(data = DatosTaller, aes(CARGA, u)) +
  geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
  theme_bw()
library(gridExtra)

grid.arrange(plot1, plot2)
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

homocedasticidad, varianza de los errores constante

## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## Loading required package: carData
ncvTest(Modelo1)
## Non-constant Variance Score Test 
## Variance formula: ~ fitted.values 
## Chisquare = 0.08056436, Df = 1, p = 0.77653
bptest(Modelo1)
## 
##  studentized Breusch-Pagan test
## 
## data:  Modelo1
## BP = 0.16549, df = 1, p-value = 0.6842
##Al realizar las pruebas, se encuentra que el test de Breusch-Pagan arroja un p-valor de 0.6842 que es > a 0.05, es decir, aceptamos la Ho de homocedasticidad en los errores.

Normalidad en los errores

library(normtest) 
library(nortest)
library(moments)
shapiro.test(u)
## 
##  Shapiro-Wilk normality test
## 
## data:  u
## W = 0.95141, p-value = 0.4792
pearson.test(u)
## 
##  Pearson chi-square normality test
## 
## data:  u
## P = 2.3529, p-value = 0.6711
##Las pruebas de Shapiro y pearson arrojan p-valores mayores a 0.05, por lo tanto se acepta la Ho de que hay normalidad en los errores.
qqnorm(Modelo1$residuals)
qqline(Modelo1$residuals)

iii. Modelo 2: Interacción entre Carga y Tipos de discos.

\[Tiempo= β_0+β_1*Carga+β_2*Tipodedisco+β_3 Carga*Tipodedisco+U_2\]

Conf <- as.numeric(DatosTaller$Conf)
Conf
##  [1] 1 0 1 0 1 1 0 0 0 1 0 1 1 1 1 0 0
Tipodedisco <- as.factor(Conf)
Modelo2 <- lm(TIEMPO ~ CARGA+ Tipodedisco + CARGA*Tipodedisco)
summary(Modelo2)
## 
## Call:
## lm(formula = TIEMPO ~ CARGA + Tipodedisco + CARGA * Tipodedisco)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.52162 -0.21362  0.03254  0.22885  0.43311 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        -1.25562    0.25399  -4.944 0.000269 ***
## CARGA               0.68923    0.04548  15.155 1.22e-09 ***
## Tipodedisco1        2.13235    0.31159   6.844 1.18e-05 ***
## CARGA:Tipodedisco1 -0.34169    0.05515  -6.195 3.24e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3068 on 13 degrees of freedom
## Multiple R-squared:  0.965,  Adjusted R-squared:  0.9569 
## F-statistic: 119.3 on 3 and 13 DF,  p-value: 1.036e-09
u2<- residuals(Modelo2)

El modelo presenta una bondad de ajuste de 0.9569, lo cual significa que el modelo explica el 95,69% de las variaciones en la variable explicativa, en este caso, el tiempo. Presenta un p-value de 1.036e-09 que frente a un p-valor<0.05, rechaza la Ho de que el modelo NO es significativo en su conjunto. Además, el ser disco 1: HDD aumenta el tiempo de respuesta del disco en 2.1323 frente a ser disco 0:SDD, dejando todo lo demás constante con un nivel de significancia de (p-valor < 0,001). Para el caso de la interacción entre el tipo de disco y la carga, al dejar todo lo demás constante, disminuye el tiempo de respuesta del disco en 0.34169 con un nivel de significancia estadística del 99.999%

Algunos supuestos para el modelo 2

homocedasticidad

ncvTest(Modelo2)
## Non-constant Variance Score Test 
## Variance formula: ~ fitted.values 
## Chisquare = 1.149731, Df = 1, p = 0.28361
##La prueba de non- constant variance arroja un p-valor de 0.28361 > 0.05, por lo tanto se acepta la Ho de que hay homocedasticidad en los errores.

spreadLevelPlot(Modelo2)

## 
## Suggested power transformation:  0.6469484
bptest(Modelo2)
## 
##  studentized Breusch-Pagan test
## 
## data:  Modelo2
## BP = 5.3269, df = 3, p-value = 0.1494
##La prueba de Breusch-Pagan arroja un p-valor de 0.1494 > 0.05, por lo tanto se acepta la Ho de que hay homocedasticidad en los errores.

Multicolinealidad

library(corrplot)
## corrplot 0.90 loaded
library(MASS)
vif(Modelo2)
##             CARGA       Tipodedisco CARGA:Tipodedisco 
##          3.132308          4.369211          6.235145
sqrt(vif(Modelo2)) > 2 
##             CARGA       Tipodedisco CARGA:Tipodedisco 
##             FALSE              TRUE              TRUE
car::vif(Modelo2)
##             CARGA       Tipodedisco CARGA:Tipodedisco 
##          3.132308          4.369211          6.235145
#Esta prueba nos dice que si los numeros son mayores a 2, existe multicolinealidad, en este caso sucede esto, sin embargo, hay que tener en cuenta la interacción entre Carga y el tipo de disco que se aplicó a este modelo.

Criterios de información: Akaike(AIC) y Bayesiano (BIC)

AIC (Modelo1, Modelo2)
##         df      AIC
## Modelo1  3 35.50385
## Modelo2  5 13.50758
BIC(Modelo1, Modelo2)
##         df      BIC
## Modelo1  3 38.00349
## Modelo2  5 17.67364
##Estos son criterios de comparación para escoger el mejor modelo. Se escoge el que arroje el valor más pequeño. En este caso, por criterios de Akaike y Bayesiano, se escoge el modelo 2 que es el de la interacción.

iv. Ajuste del modelo.

\[Tiempo= β_0+β_1*Carga+β_2*Tipodedisco+β_3 Carga*Tipodedisco+U_2\]

ggplot(data= DatosTaller) +
  geom_point(aes(x= CARGA, y= TIEMPO))+
  geom_line(aes(x= CARGA, y= predict(Modelo1)),
            col="yellow") + 
  geom_line(aes(x= CARGA, y= predict(Modelo2)),
            col="purple")

v. Test Anova

El análisis de varianza (ANOVA) prueba la hipótesis de que las medias de dos o más poblaciones son iguales. Los ANOVA evalúan la importancia de uno o más factores al comparar las medias de la variable de respuesta en los diferentes niveles de los factores.

Comparación entre modelos

anova(Modelo1, Modelo2)
## Analysis of Variance Table
## 
## Model 1: TIEMPO ~ CARGA
## Model 2: TIEMPO ~ CARGA + Tipodedisco + CARGA * Tipodedisco
##   Res.Df    RSS Df Sum of Sq      F    Pr(>F)    
## 1     15 5.6455                                  
## 2     13 1.2234  2    4.4221 23.494 4.822e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Como se puede observar, el test de anova arrojó un p-value de 4.822e-05 para el modelo 2 con una significancia estadística de 99.999% que frente a un p-value < 0.05 rechaza la Ho, por lo tanto, los modelos no presentan similar ajuste, el modelo 2 presenta mejor ajuste que el modelo 1. La inclusión de la variable cualitativa, mejora significativamente el ajuste del modelo.

vi. Supuestos

plot1.2 <- ggplot(data = DatosTaller, aes(Tipodedisco, u2)) +
  geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
  theme_bw()
plot2 <- ggplot(data = DatosTaller, aes(CARGA, u2)) +
  geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
  theme_bw()
plot3 <- ggplot(data = DatosTaller, aes(TIEMPO, u2)) +
  geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
  theme_bw()
grid.arrange(plot1.2, plot2, plot3)
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

verificando homocedasticidad, varianza de los errores constante

library("car")
ncvTest(Modelo2)
## Non-constant Variance Score Test 
## Variance formula: ~ fitted.values 
## Chisquare = 1.149731, Df = 1, p = 0.28361
library(lmtest)
bptest(Modelo2)
## 
##  studentized Breusch-Pagan test
## 
## data:  Modelo2
## BP = 5.3269, df = 3, p-value = 0.1494
##La prueba de Breusch-Pagan arroja un p-valor de 0.1494 > 0.05, por lo tanto se acepta la Ho de que hay homocedasticidad en los errores.

No Multicolinealidad

Define que las variables explicativas en una muestra no pueden ser constantes. Además, no deben existir relaciones lineales exactas entre variables explicativas (no multicolinealidad exacta).

library(corrplot)
vif(Modelo2)
##             CARGA       Tipodedisco CARGA:Tipodedisco 
##          3.132308          4.369211          6.235145
sqrt(vif(Modelo2)) > 2 
##             CARGA       Tipodedisco CARGA:Tipodedisco 
##             FALSE              TRUE              TRUE
car::vif(Modelo2)
##             CARGA       Tipodedisco CARGA:Tipodedisco 
##          3.132308          4.369211          6.235145

Normalidad en los errores:

La normalidad de los errores permite la estimación por intervalos de confianza no sólo para los coeficientes de regresión, sino también para la predicción. Cuando los errores no son normales, los intervalos y las pruebas de hipótesis no son exactas y pueden llegar a ser inválidas. Para la comprobación de la normalidad, se realiza la prueba formal mediante el Test de Pearson.

library(normtest) 
library(nortest)
library(moments)
pearson.test(u2)
## 
##  Pearson chi-square normality test
## 
## data:  u2
## P = 2.3529, p-value = 0.6711
## El test de Pearson arroja un p-valor de 0.6711 > 0.05, por lo tanto se acepta la Ho de que hay normalidad en los errores.

Distribución de los errores

median(u2)
## [1] 0.03253846
skewness.norm.test(u2)
## 
##  Skewness test for normality
## 
## data:  u2
## T = -0.25378, p-value = 0.632
kurtosis.norm.test(u2)
## 
##  Kurtosis test for normality
## 
## data:  u2
## T = 2.0677, p-value = 0.2505
##Decimos que una distribución es normal si tenemos una asimetría de 0 y una kurtosis de 3, los valores que arrojan los test no son exactos pero si muy cercanos a los anteriormente mencionados.
residplot <- function(Modelo2, nbreaks=10) {       
  z <- rstudent(Modelo2)            
  hist(z, breaks=nbreaks, freq=FALSE, xlab="Studentized Residual",
       main="Distribution of Errors")            
  rug(jitter(z), col="brown")            
  curve(dnorm(x, mean=mean(z), sd=sd(z)),
        add=TRUE, col="blue", lwd=2)     
  lines(density(z)$x, density(z)$y,      
        col="red", lwd=2, lty=2)         
  legend("topright",                     
         legend = c( "Normal Curve", "Kernel Density Curve"),         
         lty=1:2, col=c("blue","red"), cex=.7)        }
residplot(Modelo2)

###El gráfico anterior, nos muestra la distribución de los errores para el modelo 2, se puede ver que este presenta un ajuste similar a la curva azul, un poco corrida al lado derecho y con una media de 0.03253, la cuál es muy cercana a cero, por lo tanto, es pertinente decir que este se asemeja a una distribucion normal.

vii. Conclusiones.

Dada la información proporcionada en el presente trabajo, se encontró que la variación del tiempo de respuesta del disco guarda una estrecha relación con la carga del sistema del disco.

En primer lugar se compararon los rendimientos de los diferentes discos (SDD y HDD), para esto se utilizó la correlación que evalúa la tendencia entre dos variables, en este caso, entre el tiempo y la carga, ambos arrojaron valores muy altos, sin embargo, el disco SDD mostró una mayor relación entre carga y tiempo que el HDD. Esto quiere decir que existe un mayor rendimiento entre el tiempo de respuesta y la carga del sistema del disco SDD.

Para darle mayor validez a nuestro análisis, se utiliza la regresión que es un proceso estadístico para estimar las relaciones entre variables, primero se realizan dos modelos, uno simple y otro que incluya la variable cualitativa y su interacción con la carga del sistema. Estos ajustes del modelo arrojaron (gracias a los criterios de información y ANOVA) que el mejor modelo que explica la relación entre las variables Tiempo y Carga es el modelo que posee la interacción y distincion entre discos. Los resultados son significativos y poseen una buena bondad del ajuste; además, al momento de comprobarse algunos supuestos, en su totalidad cumplen, esto nos dice que este modelo es el mejor para explicar las variaciones que tiene el tiempo y la carga y que además, la interacción se da de manera positiva cuando el tipo de disco es SDD.

anova(Modelo2)
## Analysis of Variance Table
## 
## Response: TIEMPO
##                   Df  Sum Sq Mean Sq  F value    Pr(>F)    
## CARGA              1 29.2722 29.2722 311.0401 1.845e-10 ***
## Tipodedisco        1  0.8098  0.8098   8.6046   0.01164 *  
## CARGA:Tipodedisco  1  3.6123  3.6123  38.3833 3.243e-05 ***
## Residuals         13  1.2234  0.0941                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1