---
title: "Regresión lineal simple"
author: "Guillermo Ayala Gallego"  
date: "`r Sys.Date()`"  
address: Departamento de Estadística e I.O. Universidad de Valencia  
format: 
   revealjs:  
      scrollable: true
      echo: true
      embed-resources: true
---


# Recta de mínimos cuadrados

## $\ldots$
  
-  Nuestros datos son $(x_i,y_i)$ con $i=1,\ldots,n$.
- Un modelo sencillo podría ser

$$
    y_i  =  \beta_0 + \beta_1 x_i, \ i = 1, \ldots,n
$$

- No es posible.
- Una  **buena**  solución aproximada.
- Nos conformamos con
  
$$
y_i  \approx   \beta_0 + \beta_1 x_i, \ i = 1, \ldots,n
$$


## $\ldots$
   
- Una posibilidad para determinar unos **buenos** valores para
     $\beta_0$ y $\beta_1$ es considerar y encontrar los valores de los
	 coeficientes que hacen mínima esta expresión.

$$
S_t = \sum_{i=1}^n (y_i - \beta_0 - \beta_1 x_i)^2.
$$
  
- Haciendo algún cálculo se obtiene
  
$$
\hat{\beta}_1 = \frac{\sum_{i=1}^n x_i y_i - (\sum_{i=1}^nx_i \sum_{i=1}^n
y_i)/n}{\sum_{i=1}^n x_i^2 - (\sum_{i=1}^nx_i)^2/n}.
$$
	
 $$
\hat{\beta}_0 = \bar{y} - \hat{\beta}_1 \bar{x}.
$$

## Recta de mínimos cuadrados$\ldots$  
  
  - Si denotamos
		
  $$
	  S_{xy} = \sum_{i=1}^n (x_i-\bar{x}) (y_i - \bar{y}),
 $$

y 
		
$$
S_{xx} = \sum_{i=1}^n (x_i-\bar{x})^2, \ \ 
S_{yy} = \sum_{i=1}^n (y_i-\bar{y})^2,
$$
		
entonces 

$$
\hat{\beta}_1 = \frac{S_{xy}}{S_{xx}}.
$$
			
## Predicción y residuo  
  
-  Nos planteamos encontrar unas buenas aproximaciones para 
    $y_i$ utilizando   $x_i$.
- Lo natural sería tomar el siguiente valor la  **predicción**
    o  **valor ajustado** para $y_i$

$$
    \hat{y}_i = \hat{\beta}_0 + \hat{\beta}_1 x_i.
$$
- El  **residuo** es 

$$
e_i = y_i - \hat{y}_i.
$$ 

## Coeficiente de determinación: $R^2$  
  
- Es una medida de la calidad del ajuste. 

$$
    R^2 = \frac{\sum_{i=1}^n (\hat{y}_i - \bar{y})^2}{\sum_{i=1}^n 
      (y_i - \bar{y})^2} = \frac{SS(regression)}{SS(total)}.
$$

-  Se tiene 

$$
R^2 = \frac{S_{xy}^2}{S_{xx}S_{yy}},
$$

que no es más que el cuadrado del coeficiente de correlación entre la
predictora y la respuesta.

# Modelo  

## Regresión lineal simple  


- $Y_i = \beta_0 + \beta_1 x_i + \epsilon_i$ para $i=1,\ldots,n$.
- $\epsilon_i \sim N(0,\sigma^2)$.
- Los errores $\epsilon_i$ son independientes. 

## Regresión lineal simple$\ldots$  
  
- Estamos asumiendo que la distribución **condicionada**
    de $Y_i$ al predictor $x_i$ es normal con media $\beta_0 + \beta_1 x_i$
    y varianza $\sigma^2$, 

$$
Y_i | x_i \sim N(\beta_0 + \beta_1 x_i,\sigma^2),
$$

y que los distintos $Y_i$ son independientes entre sí. 

  

##  Representación matricial  

   
- Tenemos  el **vector de respuestas aleatorias**,
     $\mathbf Y$; la **matriz de modelo**,
     $\mathbf X$; el **vector de coeficientes**, $\mathbf \beta$ y el 
     **vector de errores aleatorios**, $\mathbf \epsilon$ dados por
	 
$$
\mathbf Y = 
\begin{bmatrix}
   Y_1 \\
   \vdots \\
    Y_n
\end{bmatrix};
\ \ 
\mathbf X = 
\begin{bmatrix}
  1 & x_1 \\
  \vdots&  \vdots \\
    1 & x_n
\end{bmatrix};
\ \ 
\mathbf \beta = 
\begin{bmatrix}
 \beta_0 \\
  \beta_1
\end{bmatrix};
\ \ 
\mathbf \epsilon = 
\begin{bmatrix}
   \epsilon_1 \\
   \vdots \\
    \epsilon_n
\end{bmatrix}
$$
- El modelo de regresión lineal simple se puede formular como
$$
  \mathbf Y = \mathbf X \mathbf \beta + \mathbf \epsilon.
  $$
   - El vector de errores aleatorios $\mathbf \epsilon$ sigue una
     distribución normal multivariante 
$$
  \mathbf \epsilon \sim N_n(\mathbf 0_n,\sigma^2 \mathbf I_n).
$$
   


## Verosimilitud  
  
  - La función de verosimilitud viene dada por
    \begin{multline} 
      L(\beta_0,\beta_1,\sigma^2;y_1,\ldots,y_n) =
      \prod_{i=1}^n \frac{1}{\sqrt{2\pi} \sigma}
      e^{-\frac{1}{2 \sigma^2} (y_i - \beta_0 - \beta_1 x_i)^2} = \\
      \frac{1}{(2 \pi)^{n/2} \sigma^n}
      e^{-  \frac{1}{2 \sigma^2} \sum_{i=1}^n (y_i - \beta_0 - \beta_1 x_i)^2} 
    \end{multline}
  -  Los estimadores máximo verosímiles de los coeficientes corresponden
    con los estimadores mínimo cuadráticos.
  

## Paquetes  

- Cargamos los paquetes que necesitamos.  

```{r}
#| label: rls1
pacman::p_load(Biobase,ggplot2)
```

## Datos  
- Leemos los datos.

- La primera es la más ortodoxa. 
```{r}
#| label: rls2
data(gse25171,package="tamidata2")
sel0 = which("261892_at"==fData(gse25171)[,"PROBEID"])
df0 = data.frame(pData(gse25171)[,c("time","Pi")],
                 expression=exprs(gse25171)[sel0,])
```
```{r}
#| label: rls4
#| fig-cap: Expression vs time
#| output-location: slide
(p = ggplot(df0,aes(x=time,y=expression)) + geom_point())
```
- Podemos querer guardar el gráfico.   

```{r}
#| label: rls5
#| eval: false  
ggsave("gse25171_261892_at.png",p)
```
## Recta de mínimos cuadrados  

- Ajustamos. 

```{r}
#| label: rls6
fit = lm(expression ~ time, data = df0)
```  

- ¿Qué clase de objeto tenemos?
```{r}
#| label: rls7
class(fit)
```
-  Vemos lo que el método `print` nos muestra de un objeto de clase
`lm`.
```{r}
#| label: rls8
fit
```
- Los coeficientes de la recta de mínimos cuadrados los podemos obtener
con el método  `coef` o  `coefficients`.

```{r}
#| label: rls9
coef(fit)
coefficients(fit)
```
- Añadimos la recta de mínimos cuadrados al diagrama de puntos.  
```{r}
#| label: rls10
(p = ggplot(df0,aes(x=time,y=expression))+geom_point() +
    geom_smooth(method='lm', se = FALSE))

```
- El método `predict` aplicado al objeto de clase `lm` nos proporciona
las predicciones para las observaciones. 

```{r}
#| label: rls11
predict(fit)
```
- El método `residuals` o `resid` nos proporciona los residuos. 

```{r}
#| label: rls12
residuals(fit)
resid(fit)
```
- El método `summary` aplicado a un objeto `lm`.  
```{r}
#| label: rls13
(fit.s = summary(fit))
```
-  Tenemos un objeto de clase

```{r}
#| label: rls14
class(fit.s)
```
- Este resumen tiene muchas cosas como podemos ver
```{r}
#| label: rls15
str(fit.s)
```
- Una versión simplificada puede ser
```{r}
#| label: rls16
attributes(fit.s)
```
- En particular el coeficiente de determinación lo tenemos con 
```{r}
#| label: rls17
fit.s$r.squared
```
- La matriz de modelo se obtiene con 
```{r}
#| label: rls18
model.matrix(fit)
``` 

## Contrastes e intervalos para la pendiente  

- En el `summary` tenemos los contrastes.

```{r}
#| label: rls19
summary(fit)
```  
- Los intervalos de confianza con nivel 0.95 vienen con
```{r}
#| label: rls20
confint(fit)
```

