Estadística y Machine Learning

Autor: Luis Fernando Apáez Álvarez
Clase 4: Regresión lineal

Modelo estadístico

Python

Comenzamos por cargar la librería a utilizar e importamos el conjunto de datos:

# Paquete estadistico
import statsmodels.api as sm

# Conjunto de datos
import seaborn as sns
df = sns.load_dataset('iris')
df.head()
##    sepal_length  sepal_width  petal_length  petal_width species
## 0           5.1          3.5           1.4          0.2  setosa
## 1           4.9          3.0           1.4          0.2  setosa
## 2           4.7          3.2           1.3          0.2  setosa
## 3           4.6          3.1           1.5          0.2  setosa
## 4           5.0          3.6           1.4          0.2  setosa

veamos la correlación de pearson entre las variables:

# Libreria para no ver advertencias en salida
import warnings
warnings.filterwarnings('ignore')

# Matriz de correlaciones
df.corr()
##               sepal_length  sepal_width  petal_length  petal_width
## sepal_length      1.000000    -0.117570      0.871754     0.817941
## sepal_width      -0.117570     1.000000     -0.428440    -0.366126
## petal_length      0.871754    -0.428440      1.000000     0.962865
## petal_width       0.817941    -0.366126      0.962865     1.000000

vemos una fuerte correlación entre el largo y ancho del pétalo de modo que:

# Definimos las variables
x = df.petal_width.values
y = df.petal_length.values

Graficamos un diagrama de dispersión:

import matplotlib.pyplot as plt

sns.scatterplot(x=x, y=y, hue=df.species)
plt.show()

Implementamos el modelo estadístico:

# A la matriz de predictores se le tiene que añadir una columna de 1's para el intercepto del modelo
X = sm.add_constant(x, prepend=True)

# Definicion del modelo
modelo = sm.OLS(endog=y, exog=X,)

# Ajustamos
modelo = modelo.fit()

# Vemos el resumen
print(modelo.summary())
##                             OLS Regression Results                            
## ==============================================================================
## Dep. Variable:                      y   R-squared:                       0.927
## Model:                            OLS   Adj. R-squared:                  0.927
## Method:                 Least Squares   F-statistic:                     1882.
## Date:                Fri, 05 May 2023   Prob (F-statistic):           4.68e-86
## Time:                        22:53:42   Log-Likelihood:                -101.18
## No. Observations:                 150   AIC:                             206.4
## Df Residuals:                     148   BIC:                             212.4
## Df Model:                           1                                         
## Covariance Type:            nonrobust                                         
## ==============================================================================
##                  coef    std err          t      P>|t|      [0.025      0.975]
## ------------------------------------------------------------------------------
## const          1.0836      0.073     14.850      0.000       0.939       1.228
## x1             2.2299      0.051     43.387      0.000       2.128       2.332
## ==============================================================================
## Omnibus:                        2.438   Durbin-Watson:                   1.430
## Prob(Omnibus):                  0.295   Jarque-Bera (JB):                1.966
## Skew:                           0.211   Prob(JB):                        0.374
## Kurtosis:                       3.369   Cond. No.                         3.70
## ==============================================================================
## 
## Notes:
## [1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

donde \(R^{2}=0.927\) lo que nos diría que el modelo es muy bueno. Continuando:

# coeficientes
modelo.params
## array([1.08355803, 2.2299405 ])
import numpy as np

# coeficientes
hbeta0 = modelo.params[0]
hbeta1 = modelo.params[1]

# Definimos las recta de regresion
r_reg = lambda x: hbeta0 + hbeta1 * x

# rango de graficacion
rango = np.arange(x.min(),x.max(),0.01)

# graficamos el scatter
sns.scatterplot(x=x, y=y, hue=df.species)

# graficamos la recta de regresion
plt.plot(rango, r_reg(rango))
plt.show()

R

Cargamos los datos y efectuamos la regresión

# Carga
df <- iris

# Definicion de las variables
x <- df$Petal.Width
y <- df$Petal.Length

# ajuste
fit <- lm(y~x)
summary(fit)
## 
## Call:
## lm(formula = y ~ x)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.33542 -0.30347 -0.02955  0.25776  1.39453 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.08356    0.07297   14.85   <2e-16 ***
## x            2.22994    0.05140   43.39   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4782 on 148 degrees of freedom
## Multiple R-squared:  0.9271, Adjusted R-squared:  0.9266 
## F-statistic:  1882 on 1 and 148 DF,  p-value: < 2.2e-16
plot(x, y, xlab='Ancho del pétalo', ylab='Largo del pétalo')
abline(fit)