| Duración de Llamada | 23 | 29 | 49 | 64 | 74 | 87 | 96 | 97 | 109 | 119 | 149 | 145 | 154 | 166 |
| No Componentes | 1 | 2 | 3 | 4 | 4 | 5 | 6 | 6 | 7 | 8 | 9 | 9 | 10 | 10 |
Regresión Lineal Simple en Rstudio
Regresión Lineal Simple
La regresión lineal simple es un método estadístico utilizado para entender y predecir la relación entre dos variables, una independiente (predictora) X_{i} para i = 1,2,\dots, n y una dependiente (respuesta) Y_{i} con i = 1,2,\dots, n, a través del modelo.
Y_{i} = \beta_{0} + \beta_{1}X_{i} + \varepsilon_{i}
Donde \beta_{0} (intercepto) y \beta_{1} (pendiente) son constantes llamadas coeficientes del Modelo o parámetros del modelo, y \varepsilon es un error estocástico.
Supuestos del modelo
El modelo de regresión lineal tiene varios supuestos que deben cumplirse para que los resultados sean válidos y confiables.
X_{i} es constante y conocida
Los errores \varepsilon_{i} son aleatorios e independientes
No existen valores atípicos
Se asume que los errores del modelo de regresión se distribuyen normalmente con media 0 y varianza constante (Normalidad y Homocedasticidad)
\varepsilon_{i} \sim N (0, \sigma^{2})
de estos supuesto se deriva
\begin{aligned} E[Y_{i}|X_{i} = x_{i}] &= E[\beta_{0} + \beta_{1}x_{i} + \varepsilon]\\ &= \beta_{0} + \beta_{1}x_{i} + E[\varepsilon] \\ &= \beta_{0} + \beta_{1}x_{i} \end{aligned}
y
\begin{aligned} Var[Y_{i}|X_{i} = x_{i}] &= Var[\beta_{0} + \beta_{1}x_{i} + \varepsilon_{i}] \\ &= Var[\varepsilon_{i}] = \sigma^{2} \end{aligned}
El objetivo es encontrar la línea que mejor se ajuste a los datos observados, minimizando la diferencia entre los valores predichos por la línea y los valores reales observados.
Se puede demostrar que mediante una estimación por mínimos cuadrados Ordinarios (MCO), los coenficientes estimados son
\hat{\beta}_{1} = \frac{\sum_{i = 1}^{n}(y_{i} - \bar{y})(x_{i} - \bar{x})}{\sum_{i = 1}^{n}(x_{i} - \bar{x})^{2}} \quad;\quad \hat{\beta}_{0} = \bar{y} - \hat{\beta}_{1}\bar{x}
Del cual se obtiene el modelo estimado
\hat{Y} = \hat{\beta}_{0} + \hat{\beta}_{1}X
que será la recta que se ajuste a los datos.
Interpretación de los Coeficientes estimados
Interpretación de \beta_{0}
\beta_{0} representa el valor esperado en la respuesta, cuando la variable predictora X toma el valor de 0.
\begin{aligned} E[Y_{i}|X_{i} = 0] &= E[\beta_{0} + \beta_{1}(0) + \varepsilon]\\ &= \beta_{0} + 0+ E[\varepsilon] \\ &= \beta_{0} \end{aligned}
Interpretación de \beta_{1}
\beta_{1} representa el cambio esperado en la respuesta, cuando la variable predictora X aumenta en una unidad. Si se define al cambio como la diferencia en el valor esperado de la respuesta entonces
\begin{aligned} E[Y_{i}|X_{i} = x_{i} + 1] - E[Y_{i}|X_{i} = x_{i}] &= E[\beta_{0} + \beta_{1}(x_{i} + 1) + \varepsilon] - E[\beta_{0} + \beta_{1}(x_{i}) + \varepsilon]\\ &= (\beta_{0} + \beta_{1}(x_{i} + 1) + E[\varepsilon]) - (\beta_{0} + \beta_{1}(x_{i}) + E[\varepsilon]) \\ &= (\beta_{0} + \beta_{1}x_{i} + \beta_{1}) - (\beta_{0} + \beta_{1}x_{i}) \\ &= \beta_{1} \end{aligned}
Ejemplo
Una empresa que comercializa y repara pequeñas computadoras. Para estudiar la relación entre la duración de una llamada de servicio y el número de componentes electrónicos en la computadora que deben ser reparados o reemplazados, se tomó una muestra de registros de llamadas. Los datos consisten en la duración de las llamadas en minutos (la variable de respuesta) y el número de componentes reparados (la variable predictora)
la relación lineal se puede ver en el siguiente gráfico de dispersión
y <- c(23, 29, 49, 64, 74, 87, 96, 97, 109, 119, 149, 145, 154, 166)
x <- c(1, 2, 3, 4, 4, 5, 6, 6, 7, 8, 9, 9, 10, 10)
plot(x, y,
xlab = "Numero de Componentes",
ylab = "Duración de la Llamada", main = "Correlación = 0.9936")
Estimando \bar{x}, \bar{y}
(barx <- mean(x))[1] 6(bary <- mean(y))[1] 97.21429Calculamos S_{xy} = \sum_{i = 1}^{n}(y_{i} - \bar{y})(x_{i} - \bar{x}) \text{ y } S_{xx} = \sum_{i = 1}^{n}(x_{i} - \bar{x})^{2}
(Sxy <- sum((y - bary)*(x - barx)))[1] 1768(Sxx <- sum((x - barx)^{2}))[1] 114Calculamos la estimación \hat{\beta}_{1} = \frac{\sum_{i = 1}^{n}(y_{i} - \bar{y})(x_{i} - \bar{x})}{\sum_{i = 1}^{n}(x_{i} - \bar{x})^{2}} = \frac{S_{xy}}{S_{xx}}
(beta1 <- Sxy/Sxx)[1] 15.50877Luego \hat{\beta}_{0} = \bar{y} - \hat{\beta}_{1}\bar{x}
(beta0 <- bary - beta1*barx)[1] 4.161654Por lo tanto la recta estimada es
\hat{Y_{i}} = 4.1616 + 15.508X_{i}
Se puede interpretar que:
\hat{\beta}_{0}: Cuando el número de componentes que deben ser reparados es igual a 0, se espera que la duración de la llamada sea de 4.1616 minutos.
\hat{\beta}_{1}: Se espera que, por cada componente adicional que se deba repararm la duración de la llamada aumente en 15.508 minutos.
Ajustamos la recta estimada a los valores observados
plot(x, y,
xlab = "Numero de Componentes",
ylab = "Duración de la Llamada", main = "Correlación = 0.9936");curve(Vectorize(4.1616 + 15.508*x), from = 0, to = 10,add = TRUE, col = "red")
Ajuste del modelo con la función lm
La función lm() se utiliza para ajustar modelos de regresión lineal. Permite modelar la relación entre una variable dependiente (o respuesta) y una o más variables independientes (o predictoras) mediante una ecuación lineal. La sintaxis básica de la función lm() es modelo <- lm(formula, data) junto a la función summary() se proporciona un resumen de las estadísticas del modelo. Para el ejemplo anterior:
data <- data.frame(x = x, y = y)
modelo <- lm(y ~ x, data = data)
summary(modelo)
Call:
lm(formula = y ~ x, data = data)
Residuals:
Min 1Q Median 3Q Max
-9.2318 -3.3415 -0.7143 4.7769 7.8033
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 4.162 3.355 1.24 0.239
x 15.509 0.505 30.71 8.92e-13 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 5.392 on 12 degrees of freedom
Multiple R-squared: 0.9874, Adjusted R-squared: 0.9864
F-statistic: 943.2 on 1 and 12 DF, p-value: 8.916e-13Podemos ver en el resumen del modelo los coeficientes estimados, los errores estándar, los valores p, coeficiente R^{2}, entre otros.
Verificación de los Supuestos del Modelo
La verificación de los supuestos del modelo es crucial para garantizar la validez de las inferencias y predicciones realizadas. Observemos la siguiente gráfica de los residuos del modelo.
library(ggplot2)
diagnorm = function(fit.model){
# fit.model: objeto con el resultado del modelo normal
# lineal homocedástico obtenido a través de la función "lm"
if (!is(fit.model, "lm")) stop("fit.model debe ser un objeto retornado por la función lm.")
X = model.matrix(fit.model)
n = nrow(X)
p = ncol(X)
H = X%*%solve(t(X)%*%X)%*%t(X)
h = diag(H)
lms = summary(fit.model)
s = lms$sigma
r = resid(lms)
ts = r/(s*sqrt(1-h))
di = (1/p)*(h/(1-h))*(ts^2)
si = lm.influence(fit.model)$sigma
tsi = r/(si*sqrt(1-h))
a = max(tsi)
b = min(tsi)
#
dados2 = data.frame(x=1:length(tsi), tsi, fit=fitted(fit.model))
f1 = ggplot(dados2, aes(x=x, y=tsi)) + geom_point() +
geom_hline(yintercept=c(-2, 0, 2), linetype="dashed", color="red") +
labs(x="Índice", y="Residuos Studentizados") + theme_bw()
#
f2 = ggplot(dados2, aes(x=fit, y=tsi)) + geom_point() +
geom_hline(yintercept=c(-2, 0, 2), linetype="dashed", color="red") +
labs(x="Valores ajustados", y="Residuos Studentizados") + theme_bw()
#
f3 = ggplot(dados2, aes(x=tsi, y=after_stat(density))) + geom_histogram(bins=11, fill="grey", color="black") +
labs(x="Residuos Studentizados", y="Densidad") + theme_bw()
#
ident = diag(n)
epsilon = matrix(0,n,100)
e = matrix(0,n,100)
e1 = numeric(n)
e2 = numeric(n)
#
for(i in 1:100){
epsilon[,i] = rnorm(n,0,1)
e[,i] = (ident - H)%*%epsilon[,i]
u = diag(ident - H)
e[,i] = e[,i]/sqrt(u)
e[,i] = sort(e[,i]) }
#
for(i in 1:n){
eo = sort(e[i,])
e1[i] = (eo[2]+eo[3])/2
e2[i] = (eo[97]+eo[98])/2 }
#
med = apply(e,1,mean)
faixa = range(tsi,e1,e2)
#
dados3 = data.frame(e1, e2, med)
q1 = data.frame(qqnorm(tsi, plot.it=FALSE))
q2 = data.frame(qqnorm(e1, plot.it=FALSE))
q3 = data.frame(qqnorm(e2, plot.it=FALSE))
q4 = data.frame(qqnorm(med, plot.it=FALSE))
f4 = ggplot() + geom_line(data=q2, aes(x=x, y=y)) + ylim(faixa) +
geom_line(data=q3, aes(x=x, y=y)) + geom_line(data=q4, aes(x=x, y=y), color="red", linetype="dashed") +
geom_point(data=q1, aes(x=x, y=y), size=0.6) + labs(x="Percentiles de la N(0,1)", y="Residuos Studentizados") + theme_bw()
gridExtra::grid.arrange(f1, f2, f3, f4, nrow=2)
return(tsi)
} # Crea la función para realizar los gráficos en ggplot (SOLO COPIAR LA FUNCIÓN)
diagnorm(modelo) # Grafica los gráficos para los residuos 
1 2 3 4 5 6
0.71831325 -1.33183234 -0.32675848 -0.41591246 1.63421422 1.02613790
7 8 9 10 11 12
-0.22427554 -0.03949066 -0.70472596 -2.03462727 1.06408274 0.24327244
13 14
-1.10685845 1.47824066 Verificar Homcedasticidad:
El primer gráfico muestra los residuos (diferencias entre los valores observados y los valores predichos por el modelo) estudentizados en el eje vertical y los índices en el eje horizontal, por otra parte, el segundo gráfico muestra los residuos estudentizados vs los ajustados.
Se utiliza para evaluar el supuesto de homocedasticidad. La dispersión uniforme de los puntos alrededor de una línea horizontal sugiere homocedasticidad.
Verificar Normalidad:
El tercer gráfico muestra un histigrama de los residuos Studentizados, mientras que el cuerto gráfico compara los percentiles observados de los residuos Studentizados con los cuantiles esperados de una distribución normal.
Se utiliza para evaluar el supuesto de normalidad de los residuos. Si los puntos en el gráfico siguen aproximadamente una línea diagonal, indica que los residuos siguen una distribución normal. Desviaciones significativas de la línea diagonal pueden indicar violaciones del supuesto de normalidad.
Para los gráficos del ejemplo podemos Interpretar lo siguiente:
Los residuos no muestran patrones claros, lo que es bueno, pero hay una ligera curvatura, lo que podría indicar una relación no lineal no capturada por el modelo.
Los puntos siguen en gran medida la línea diagonal, pero hay desviaciones en los extremos, lo que sugiere que los residuos podrían no ser completamente normales.
Los puntos están dispersos aleatoriamente, lo que es positivo para la homocedasticidad, aunque hay una tendencia a una mayor dispersión para valores ajustados más grandes.
Gráficos similares se obtienen con la siguiente instrucción
par(mfrow = c(2,2))
plot(modelo)
par(mfrow = c(1,1))Predicciones
Una predicción puntual para un nuevo valor X_{i} se puede obtener como
Y_{pred} = \hat{\beta_{0}} + \hat{\beta_{1}}X_{nuevo} Donde Y_{pred} es la predicción para el nuevo valor X_{nuevo}.
La función predict() se utiliza para hacer predicciones utilizando un modelo de regresión lineal. La sintaxis básica es predict(modelo, newdata), donde modelo es el objeto lm ajustado y newdata es un dataframe que contiene los valores de las variables predictoras para los cuales se desean hacer las predicciones.
Por ejemplo, si deseamos predecir el tiempo que dura una llamada si los componentes a reparar o remplazar son 12, 13 o 14.
(predict(modelo, newdata = data.frame(x = c(12, 13, 14)))) 1 2 3
190.2669 205.7757 221.2845 Entonces
Se predice que cuando se deben reparar o reemplazar 12 componentes, se espera que la llamada dure aproximadamente 190.27 minutos.
Se predice que cuando se deben reparar o reemplazar 13 componentes, se espera que la llamada dure aproximadamente 205.78 minutos.
Se predice que cuando se deben reparar o reemplazar 14 componentes, se espera que la llamada dure aproximadamente 221.28 minutos.
Al usar las predicciones de un modelo de regresión lineal, es esencial tener en cuenta su contexto y las limitaciones del modelo. Aunque una predicción como una llamada de 221 minutos puede parecer poco plausible en condiciones normales, puede ser plausible bajo circunstancias excepcionales o con tareas complejas de reparación. Por lo tanto, es crucial combinar las predicciones con el juicio humano y estar preparado para ajustarlas según sea necesario, reconociendo que los modelos estadísticos son simplificaciones de la realidad y pueden no capturar todas las complejidades del proceso analizado.
Intervalos de Confianza para los Parámetros
Los intervalos de confianza son cruciales en la interpretación de un modelo de regresión, ya que indican la precisión de las estimaciones de los coeficientes. Si el intervalo de confianza de un coeficiente no incluye el valor cero, esto sugiere que la variable predictora asociada contribuye significativamente al modelo. Por el contrario, si el intervalo de confianza incluye el valor cero, la variable predictora podría no ser relevante. Estos intervalos proporcionan una guía valiosa para evaluar la importancia de cada variable en el modelo y tomar decisiones informadas sobre su inclusión o exclusión.
Para el coeficiente \beta_0 con un nivel de confianza (1 - \alpha)\times 100\% el IC esta dado por:
IC({\beta}_0, 1 - \alpha) = \hat{\beta}_0 \pm t_{(n-2,1-\alpha/2)}\sqrt{\frac{\hat{\sigma}^2}{n} \frac{\sum_{i = 1}^{n} x_{i}^{2}}{\sum_{i=1}^{n}(x_i - \bar{x})^2}}
Para el coeficiente \beta_1 con un nivel de confianza (1 - \alpha)\times 100\% el IC esta dado por:
IC(\beta_{1}, 1 - \alpha) = \hat{\beta}_{1} \pm t_{n-2,1-\alpha/2}\sqrt{\frac{\hat{\sigma}^2}{\sum_{i=1}^{n}(x_i - \bar{x})^2}}
La función confint en R calcula intervalos de confianza para los coeficientes estimados de un modelo de regresión lineal.
Para el ejemplo anterior
confint(modelo, #Objeto resultante de la función lm
level = 0.95 # Nivel de confianza deseado
) 2.5 % 97.5 %
(Intercept) -3.148482 11.47179
x 14.408512 16.60903