Regresión lineal simple
Gonzalo Pérez
28 de febrero de 2020
Ejemplo de motivación.
Una muestra aleatoria de 20 hombres de los Estados Unidos fue seleccionada y la información siguiente fue recolectada:
\(X=\) el número de gramos de grasa consumida por día.
\(Y=\) el total de colesterol en la sangre en miligramos.
Los datos son:
#Ingresaremos los datos
X=c(21,29, 43, 52, 56, 64, 77, 81, 84, 93, 98, 101, 107, 109, 113, 120, 127, 134, 148, 157)
Y=c(130, 163, 169, 136, 187, 193, 170, 115, 196, 237, 214, 239, 258, 283, 242, 289, 298, 271, 297, 316 )
Datos=as.data.frame(cbind(X,Y))
head(Datos)## X Y
## 1 21 130
## 2 29 163
## 3 43 169
## 4 52 136
## 5 56 187
## 6 64 193str(Datos)## 'data.frame': 20 obs. of 2 variables:
## $ X: num 21 29 43 52 56 64 77 81 84 93 ...
## $ Y: num 130 163 169 136 187 193 170 115 196 237 ...require(ggplot2)
require(ggiraph)
require(ggiraphExtra)
ggPoints(aes(x=X,y=Y),smooth=FALSE, data=Datos,interactive=TRUE)¿Qué se observa?
¿Se podría considerar que la relación es lineal?
Modelo de regresión lineal simple.
En el modelo de regresión lineal simple se asume lo siguiente:
\[y_i=\alpha+\beta x_i + \epsilon_i, i=1,..., n,\]
donde \(E(\epsilon_i)=0, \; V(\epsilon_i)=\sigma^2 \;\; \text{y} \;\; Cov(\epsilon_i, \epsilon_j)=0, \; i\neq j \; \; \forall i,j = 1,...,n.\)
En general, se realiza el supuesto adicional siguiente para realizar inferencias sobre los parámetros:
\[y_i \sim N(\mu_i, \sigma^2),\] con \(\mu_i=\alpha+\beta x_i\) y \(Cov(y_i, y_j)=0, \; i\neq j \; \; \forall i,j = 1,...,n.\)
Estimación
La estimación de los parámetros en el modelo de regresión lineal simple se puede realizar a partir de los estimadores obtenidos por el método de mínimos cuadrados o bien por máxima verosimilitud.
En el caso de mínimos cuadrados, se buscan los valores de \(\alpha\) y \(\beta\) que minimizan la suma de cuadrados de los errores definidos como \(e_i=y_i-\alpha-\beta x_i, \; i=1,...,n\). Es decir, minimizar \(\sum_{i=1}^n (y_i-\alpha-\beta x_i)^2\).
En el caso de máxima verosímilitud, se debe maximizar la log-verosimilitud:
\[l = log(L) =-\frac{n}{2}log(2\pi) - \frac{n}{2} log(\sigma^2)-\frac{1}{2 \sigma^2}\sum_{i=1}^n (y_i-\alpha-\beta x_i)^2. \] Los estimadores de \(\alpha\) y \(\beta\) obtenidos con ambos métodos coinciden y son:
\[\begin{align} \widehat{\alpha} = \overline{y} - \widehat{\beta}\overline{x} \end{align}\]
\[\begin{align} \widehat{\beta} = \frac{\sum_{i=1}^{n}(x_{i}-\overline{x})(y_{i}-\overline{y})}{\sum_{i=1}^{n}(x_{i}-\overline{x})^2} \end{align}\]
Nota: para encontrar los estimadores es necesario derivar con respecto a cada uno de los parámetros y resolver el sistema de ecuaciones que se obtiene al igualar las derivadas igual a cero
Ambos estimadores son insesgados, es decir, \(E(\widehat{\alpha})=\alpha\) y \(E(\widehat{\beta})=\beta\).
Para obtener las propiedades de los estimadores basta escribirlos como combinación lineal de los valores de \(y_i\)
Algunas estadísticas que se pueden calcular con los datos son:
\(\sum x_i=1,814\), \(\sum y_i=4,403\), \(\sum y_i^2=1,041,419\), \(\sum x_i^2=192,204\), \(\sum x_i y_i=438,848\),
\(\sum (x_i-\overline{x})^2=\sum x_i^2-n\overline{x}^2=27,674.20\),
\(\sum (y_i-\overline{y})^2=\sum y_i^2-n\overline{y}^2=72,098.55\), y
\(\sum (x_i-\overline{x})(y_i-\overline{y})=\sum x_iy_i-n\overline{x}\overline{y}=39,495.9\).
xbar=mean(Datos$X)
SSx=sum((Datos$X-xbar)^2)
ybar=mean(Datos$Y)
SSy=sum((Datos$Y-ybar)^2)
SSxy=sum((Datos$Y-ybar)*(Datos$X-xbar))
beta=SSxy/SSx
alpha=ybar-beta*xbar
options(digits=5)
print (c(alpha, beta))## [1] 90.7053 1.4272La linea azul corresponde a la recta ajustada
ggPoints(aes(x=X,y=Y),method="lm", data=Datos,interactive=TRUE)Las varianzas de los estimadores de \(\alpha\) y \(\beta\) son:
\[\begin{align} Var(\hat{\beta}) = \frac{\sigma^2}{\sum_{i=1}^{n}(x_{i}-\overline{x})^2} \end{align}\]
\[\begin{align} Var(\hat{\alpha}) = \sigma^2 [\frac{\sum_{i=1}^{n} x_i^2}{n\sum_{i=1}^{n}(x_{i}-\overline{x})^2}] \end{align}\]
Un estimador insesgado de \(\sigma^2\) es
\[\begin{align} \widehat{\sigma}^2= \dfrac{\sum_{i=1}^n (y_i-\widehat{\alpha}-\widehat{\beta} x_i)^2}{n-2} = \dfrac{\sum_{i=1}^n (y_i-\overline{y})^2-\widehat{\beta}^2\sum_{i=1}^n (x_i-\overline{x})^2}{n-2}\end{align} \]
n=length(Datos$X)
SSE=sum((Datos$Y-alpha-beta*Datos$X)^2)/(n-2)
print(SSE)## [1] 873.95SSEv2=(SSy-beta^2*SSx)/(n-2)
print(SSEv2)## [1] 873.95Usando esta estimación también se pueden obtener estimaciones de las varianzas de los parámetros \(\alpha\) y \(\beta\).
A partir de esto es posible observar que:
\[\dfrac{\hat{\beta} -\beta}{\sqrt{ \widehat{Var}(\hat{\beta})}} \sim t_{n-2}\]
Intervalos de confianza sobre \(\alpha\) y \(\beta\).
De donde es posible obtener un intervalo al \((1-\delta)\)% de confianza para \(\beta\).
\[[\hat{\beta} - t_{n-2, 1-\frac{\delta}{2}}\sqrt{ \widehat{Var}(\hat{\beta})}, \; \hat{\beta} + t_{n-2, 1-\frac{\delta}{2}}\sqrt{ \widehat{Var}(\hat{\beta})}] \]
De manera similar, un intervalo de confianza para \(\alpha\) es:
\[[\hat{\alpha} - t_{n-2, 1-\frac{\delta}{2}}\sqrt{ \widehat{Var}(\hat{\alpha})}, \; \hat{\alpha} + t_{n-2, 1-\frac{\delta}{2}}\sqrt{ \widehat{Var}(\hat{\alpha})}] \]
var.alpha=SSE*sum(Datos$X^2)/(n*SSx)
var.beta=SSE/SSx
print(sqrt(var.alpha))## [1] 17.421print(sqrt(var.beta))## [1] 0.17771#Intervalos al 95%
delta=.05
qt(1-delta/2, n-2,lower.tail = TRUE)## [1] 2.1009#alpha
print(c(alpha-qt(1-delta/2, n-2,lower.tail = TRUE)*sqrt(var.alpha),
alpha+qt(1-delta/2, n-2,lower.tail = TRUE)*sqrt(var.alpha)))## [1] 54.105 127.305#beta
print(c(beta-qt(1-delta/2, n-2,lower.tail = TRUE)*sqrt(var.beta),
beta+qt(1-delta/2, n-2,lower.tail = TRUE)*sqrt(var.beta)))## [1] 1.0538 1.8005También es posible realizar pruebas de hipótesis sobre combinaciones lineales de los parámetros.
Sea \(b_0\) una constante.
Pruebas de hipótesis sobre \(\alpha\) y \(\beta\).
a)
\[H_0: \beta = b_0 \] vs \[H_a: \beta \neq b_0 \]
b)
\[H_0: \beta \leq b_0 \] vs \[H_a: \beta > b_0 \]
En ambos casos, la estadística de prueba es
\[t^*= \dfrac{\hat{\beta} -b_0}{\sqrt{ \widehat{Var}(\hat{\beta})}}\]
La regla de decisión para una significancia \(\delta\) en el caso a) es:
Se rechaza \(H_0\) si \(|t^*|>t_{n-2, 1-\frac{\delta}{2}}\)
Para el caso b) es:
Se rechaza \(H_0\) si \(t^*>t_{n-2, 1-\delta}\)
**Nota. Cuando el objetivo es estudiar si hay una relación lineal entre la variable dependiente y la independiente se considera \(b_0=0\) con el caso a), mientras que si, por ejemplo, se quiere estudiar si la pendiente es positiva se considera \(b_0=0\) con el caso b).
Las pruebas de hipótesis asociadas al parámetro \(\alpha\) tienen las mimas reglas de decisión, sólo que en ese caso
\[t^*= \dfrac{\hat{\alpha} -a_0}{\sqrt{ \widehat{Var}(\hat{\alpha})}}\]
Ejercicio 1
En el ejemplo ¿se puede concluir que por cada gramo de grasa que se consume al día el aumento promedio en el total de colesterol en sangre supera un miligramo?
Escribir las hipótesis asociadas, estadística de prueba y regla de decisión. Usar \(\delta=.05\)
.
.
.
test=(beta-1)/sqrt(var.beta)
delta=.05
t= qt(1-delta, n-2,lower.tail = TRUE)
print(test)## [1] 2.4038print(t)## [1] 1.7341Con los resultados anteriores concluir.
.
.
.
Usando la función lm() de R.
fit1=lm(Y~X, data=Datos)
print(summary(fit1))##
## Call:
## lm(formula = Y ~ X, data = Datos)
##
## Residuals:
## Min 1Q Median 3Q Max
## -91.31 -11.86 6.74 16.51 36.73
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 90.705 17.421 5.21 5.9e-05 ***
## X 1.427 0.178 8.03 2.3e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 29.6 on 18 degrees of freedom
## Multiple R-squared: 0.782, Adjusted R-squared: 0.77
## F-statistic: 64.5 on 1 and 18 DF, p-value: 2.32e-07confint(fit1, level=.95)## 2.5 % 97.5 %
## (Intercept) 54.1053 127.3053
## X 1.0538 1.8005Ejercicio 2
Identificar lo que se calcula en la salida de R. Considerar la prueba de hipótesis tipo a) con \(a_0=0\) y \(b_0=0\).
Ejercicio 3
Estime el valor promedio del total de colesterol en la sangre para las personas cuya ingesta diaria de grasa es 50 gramos. Presente un intervalo de confianza al 95%.