Luis Eduardo Pino Villarreal MD, MSc, MBA
Abordaremos hoy el primer modelo que todo científico de datos debe(ría) conocer, me refiero a la regresión lineal (RL).
La RL es uno de los modelos que intenta identificar una relación entre dos o más variables y utilizar dicha relación para predecir valores de una variable con base en los valores dados de otra variable. Como su nombre lo indica, en RL se asume que las relaciones entre variables pueden modelarse a través de la ecuación lineal (y = a + bx) y existe en ella una variable que se usa en la predicción a la cual llamaremos independiente, explicatoria o regresora y una variable predicha a la que llamaremos dependiente, objetivo o de respuesta. En RL entonces asumimos que las variables independientes están relacionadas linealmente con la variable dependiente según la siguiente ecuación:
𝑦 = 𝑐 + 𝑚𝑥
En ML esta ecuación la representamos mejor mediante:
𝑦 = 𝑤0 + 𝑤1𝑥
En donde: w0 es el intercepto en el eje Y, w1 es la pendiente de la línea, x es una variable independiente y y la variable independiente u objetivo.
Por ejemplo, si tenemos las cifras de tensión arterial (TA) de un paciente junto con datos de su edad, podemos usar esta información para predecir los niveles de TA de otro paciente de acuerdo con su edad, este es un problema modelado como de regresión en el cual la edad es una variable independiente y la TA la variable objetivo, entendiendo su relación así:
TA = w0 + w1 x (edad)
En el mundo real no existen problemas con una sola variable independiente y por tanto tendremos que utilizar otras variantes de RL como la regresión múltiple, pero para fines de explicación básica estos ejemplos funcionan bien.
El mejor gráfico relacional entre estos tipos de variables para el mundo de la RL es el scatterplot que nos permite visualizar la posible fuerza de relación entre dos variables (ver nuestra columna sobre visualización de datos), pero esencialmente nos sirve para intuir si nuestro problema podrá o no resolverse a través de un algoritmo lineal. Es importante entonces conocer cuáles son las medidas que nos permiten establecer la fortaleza de asociación entre variables, es decir las medidas de correlación.
Medidas de correlación y Cálculo de Línea
La más común de todas es el coeficiente de correlación que es una medida para determinar la fuerza de asociación entre dos variables. Dicho coeficiente va de -1 (fuerte asociación negativa) a 1 (fuerte asociación positiva), pasando por el cero (0) que indica ausencia de relación entre variables. En ML las matrices de correlación nos pueden dar una idea inicial de este coeficiente:
Ejemplo
En esta matriz podemos ver correlaciones positivas altas como la de las variables Power y Engine y otras negativas como Engine y Mileage (ejemplo tomado de un dataset de carros usados disponible en Kaggle).
Ahora bien, si tenemos un conjunto de datos tabulares con dos variables posiblemente correlacionadas, por ejemplo:
| Horas de estudio por día | Calificación final (/10) |
| 6 | 7 |
| 5 | 5 |
| 4 | 2 |
| 7 | 9 |
| 8 | 6 |
| 2 | 3 |
Tendremos la siguiente gráfica en un plano cartesiano (en este caso la variable dependiente es la calificación final (X) y la variable independiente es horas de estudio por día (Y))
Si observamos la gráfica podríamos pensar en trazar una línea que uniera estas variables estableciendo su relación lineal, la pregunta es cómo puede calcularse la “mejor línea” para cubrir dicha correlación. Regresemos entonces a nuestra ecuación base:
y= w0 + w1 * x
Recordemos que W0 es el intercepto y X1 es la pendiente. Si trazamos una línea cualquiera es posible que algunos de los puntos del plano cartesiano NO sean cubiertos por ella y excepto que la línea pase directamente por ese punto, habra siempre un error o residual que será la “distancia” entre el punto y la línea, este puede medirse así:
La predicción (yˆi) para cualquier x (xi en este caso) es:
yˆi = w0 + w1*xi
Entonces el valor de los residuales (e) para cualquier punto i es la diferencia:
ei = yi – yˆi
Teniendo en cuenta lo que es un residual tenemos entonces que la mejor línea será aquella que minimice la suma de los errores cuadrados (es decir los ei2)
Ei2 = ( y1 – yˆ1)2 + ( y2 – yˆ2)2 + (yn – yˆn)2
Los cálculos secuenciales para encontrar esta línea se denominan relaciones lineasles cuadráticas y pueden verse así:
Los dos métodos mas utilizados para calcular esta línea son el de mínimos cuadrados y el de gradiente descendente (atención a este último que será muy importante en modelos posteriores de ML).
Dado que no es objeto de la columna profundizar en estos métodos remito al lector a las referencias recomendadas al final de nuestra serie. Debo mencionar que los métodos a ejecutar difieren en su contexto matemático y en sus medidas de desempeño, en el caso de los mínimos cuadrados aplicaremos el R2 y el R2 ajustado mientras en el gradiente descendente utilizaremos funciones de costo. Los dos modelos incluyen cálculos de derivadas y por tanto es importante antes de profundizar en ellos desempolvar un libro de cálculo básico.
Cómo mencionamos previamente, en medicina no existen problemas que se ajusten a una regresión lineal simple, usualmente tenemos múltiples variables que influencia un desenlace. En este caso nos asistimos de los modelos de regresión lineal múltiple cuya ecuación es:
y= w0 + w1*x1 + w2*x2 + w3*x3…+wn*xn
En este caso y es la variable dependiente y x las diversas variables dependientes. Los métodos para calcular los coeficientes y la mejor línea son los descritos previamente.
Medidas de desempeño
Un elemento crucial en el desarrollo de los modelos de ML son las métricas de desempeño. Para el caso de las regresiones las medidas de desempeño que se utilizan son:
- R2 y R2 ajustado
El R2 es esencialmente la diferencia en la propagación o el “spread” entre los residuales iniciales y los encontrados post regresión. En ML el R2 es una medida interna de la muestra (in sample).
Puede expresarse como
R2 = 1 – ∑ ei2 / ∑ (y1 – yˆ)2
En general R2 expresa la fracción de la varianza en y que es explicada por la regresión y se expresa de 0 a 1, en donde 1 expresa una correlación máxima.
El R2 ajustado lo utilizamos especialmente en regresiones múltiples dado que su valor se aumenta con cada variable añadida.
- Error absoluto promedio (MAE)
Es la media de la diferencia entre los variables actuales (y) y predichos (yˆ) de la variable dependiente.
- Raiz cuadrada del error cuadrático promedio (RMCE)
Esta ecuación cuadrática se utiliza para estimar la desviación estandar de los residuales y se usa la potencia para disminuir el sesgo que produce el n-1, la ecuación se expresa:
Ventajas y Limitaciones de la Regresión Lineal
| Ventajas | Limitaciones |
| – Es un modelo de amplio uso y muy antiguo . – Simple y elegante – Computacionalmente muy eficiente. – Facil interpretación de sus salidas (outputs) | – Muy simple para la complejidad del mundo real – Asume relaciones lineales que son infrecuentes – Los outliers pueden afectar su desempeño – Asume independencia entre atributos, es dificl que pueda gestionar la multicolinealidad |
Acompañanos a ver un ejemplo de regresión lineal (Anexo 1): file:///C:/Users/57310/Desktop/Ejemplo_RegLineal.html
En nuestra próxima columna hablaremos sobre regresión logística, aquella en la cual la variable dependiente no es un número sino una o varias categorías y por tanto difiere de la regresión lineal simple y múltiple en algunas características y métricas.
Anexo 1. Ejemplo de regresión líneal.
Regresión Lineal by AIpocrates
A continuación un ejemplo básico de Regresión Lineal en una base de datos de carros
Primero importamos las librerías con las que trabajaremos
In [1]:
import numpy as np # importamos Numpy para funciones matemáticas from sklearn.linear_model import LinearRegression # importanos el modelo de regresión desde Sklearn import pandas as pd #Importamos Pandas para manejo de dataframes y math import matplotlib.pyplot as plt # importamos nuestro visor de gráficas %matplotlib inline # Este código nos permite trabajar las gráficas en el mismo archivo import seaborn as sns # nos sirve como apoyo en la visualización de datos estadísticos from sklearn.model_selection import train_test_split # Libreria para fraccionar los datos en entrenamiento y prueba
C:\Users\Rajet singh\Anaconda3\lib\site-packages\statsmodels\tools\_testing.py:19: FutureWarning: pandas.util.testing is deprecated. Use the functions in the public API at pandas.testing instead. import pandas.util.testing as tm
Cargando y revisando la base de datos
In [2]:
cData = pd.read_csv("auto-mpg.csv") # importamos la base de datos en formato CSV
cData.shape # vemos las características de la base (en este caso tiene 398 registros y 9 variables)
Out[2]:
(398, 9)
In [3]:
# Las 9 variables: # Numero # MPG (miles per gallon), # cylinders, # engine displacement (cu. inches), # horsepower, # vehicle weight (lbs.), # time to accelerate from O to 60 mph (sec.), # model year (modulo 100), and # origin of car (1. American, 2. European,3. Japanese). # veamos las primeras cinco columnas con este código: cData.head()
Out[3]:
| mpg | cylinders | displacement | horsepower | weight | acceleration | model year | origin | car name | |
|---|---|---|---|---|---|---|---|---|---|
| 0 | 18.0 | 8 | 307.0 | 130 | 3504 | 12.0 | 70 | 1 | chevrolet chevelle malibu |
| 1 | 15.0 | 8 | 350.0 | 165 | 3693 | 11.5 | 70 | 1 | buick skylark 320 |
| 2 | 18.0 | 8 | 318.0 | 150 | 3436 | 11.0 | 70 | 1 | plymouth satellite |
| 3 | 16.0 | 8 | 304.0 | 150 | 3433 | 12.0 | 70 | 1 | amc rebel sst |
| 4 | 17.0 | 8 | 302.0 | 140 | 3449 | 10.5 | 70 | 1 | ford torino |
La idea de este ejercicio será evaluar qué variable o variables se relacionan con el consumo de millas por galón (mpg) de este grupo de carros.
In [4]:
# vamos a quitar la variable car name ya que no contribuye al ejercicio
cData = cData.drop('car name', axis=1)
# Reemplazemos la variable origin con los nombres respectivos
cData['origin'] = cData['origin'].replace({1: 'america', 2: 'europe', 3: 'asia'})
cData.head()
Out[4]:
| mpg | cylinders | displacement | horsepower | weight | acceleration | model year | origin | |
|---|---|---|---|---|---|---|---|---|
| 0 | 18.0 | 8 | 307.0 | 130 | 3504 | 12.0 | 70 | america |
| 1 | 15.0 | 8 | 350.0 | 165 | 3693 | 11.5 | 70 | america |
| 2 | 18.0 | 8 | 318.0 | 150 | 3436 | 11.0 | 70 | america |
| 3 | 16.0 | 8 | 304.0 | 150 | 3433 | 12.0 | 70 | america |
| 4 | 17.0 | 8 | 302.0 | 140 | 3449 | 10.5 | 70 | america |
Creación de variables Dummy
Creamos variables binarias para facilitar el análisis de los datos
In [5]:
cData = pd.get_dummies(cData, columns=['origin']) cData.head()
Out[5]:
| mpg | cylinders | displacement | horsepower | weight | acceleration | model year | origin_america | origin_asia | origin_europe | |
|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 18.0 | 8 | 307.0 | 130 | 3504 | 12.0 | 70 | 1 | 0 | 0 |
| 1 | 15.0 | 8 | 350.0 | 165 | 3693 | 11.5 | 70 | 1 | 0 | 0 |
| 2 | 18.0 | 8 | 318.0 | 150 | 3436 | 11.0 | 70 | 1 | 0 | 0 |
| 3 | 16.0 | 8 | 304.0 | 150 | 3433 | 12.0 | 70 | 1 | 0 | 0 |
| 4 | 17.0 | 8 | 302.0 | 140 | 3449 | 10.5 | 70 | 1 | 0 | 0 |
Tratamiento de los datos perdidos
In [6]:
#Un breve resumen de las columnas: cData.describe()
Out[6]:
| mpg | cylinders | displacement | weight | acceleration | model year | origin_america | origin_asia | origin_europe | |
|---|---|---|---|---|---|---|---|---|---|
| count | 398.000000 | 398.000000 | 398.000000 | 398.000000 | 398.000000 | 398.000000 | 398.000000 | 398.000000 | 398.000000 |
| mean | 23.514573 | 5.454774 | 193.425879 | 2970.424623 | 15.568090 | 76.010050 | 0.625628 | 0.198492 | 0.175879 |
| std | 7.815984 | 1.701004 | 104.269838 | 846.841774 | 2.757689 | 3.697627 | 0.484569 | 0.399367 | 0.381197 |
| min | 9.000000 | 3.000000 | 68.000000 | 1613.000000 | 8.000000 | 70.000000 | 0.000000 | 0.000000 | 0.000000 |
| 25% | 17.500000 | 4.000000 | 104.250000 | 2223.750000 | 13.825000 | 73.000000 | 0.000000 | 0.000000 | 0.000000 |
| 50% | 23.000000 | 4.000000 | 148.500000 | 2803.500000 | 15.500000 | 76.000000 | 1.000000 | 0.000000 | 0.000000 |
| 75% | 29.000000 | 8.000000 | 262.000000 | 3608.000000 | 17.175000 | 79.000000 | 1.000000 | 0.000000 | 0.000000 |
| max | 46.600000 | 8.000000 | 455.000000 | 5140.000000 | 24.800000 | 82.000000 | 1.000000 | 1.000000 | 1.000000 |
In [9]:
# Reemplazamos los valores perdidos con NaN (not a number)
cData = cData.replace('?', np.nan)
cData[hpIsDigit['horsepower'] == False]
Out[9]:
| mpg | cylinders | displacement | horsepower | weight | acceleration | model year | origin_america | origin_asia | origin_europe | |
|---|---|---|---|---|---|---|---|---|---|---|
| 32 | 25.0 | 4 | 98.0 | NaN | 2046 | 19.0 | 71 | 1 | 0 | 0 |
| 126 | 21.0 | 6 | 200.0 | NaN | 2875 | 17.0 | 74 | 1 | 0 | 0 |
| 330 | 40.9 | 4 | 85.0 | NaN | 1835 | 17.3 | 80 | 0 | 0 | 1 |
| 336 | 23.6 | 4 | 140.0 | NaN | 2905 | 14.3 | 80 | 1 | 0 | 0 |
| 354 | 34.5 | 4 | 100.0 | NaN | 2320 | 15.8 | 81 | 0 | 0 | 1 |
| 374 | 23.0 | 4 | 151.0 | NaN | 3035 | 20.5 | 82 | 1 | 0 | 0 |
In [10]:
# Los datos perdidos los vamos a tratar reemplazandolos por una medida de tendencia central (mediana) que es para cada variable: cData.median()
Out[10]:
mpg 23.0 cylinders 4.0 displacement 148.5 horsepower 93.5 weight 2803.5 acceleration 15.5 model year 76.0 origin_america 1.0 origin_asia 0.0 origin_europe 0.0 dtype: float64
In [11]:
# Para reemplazarlos creamos una función lambda:
medianFiller = lambda x: x.fillna(x.median())
cData = cData.apply(medianFiller,axis=0)
cData['horsepower'] = cData['horsepower'].astype('float64') # converting the hp column from object / string type to float
Análisis Bivariado
Vamos a ver el análisis gráfico bivariado usando pairplot
In [12]:
cData_attr = cData.iloc[:, 0:7] sns.pairplot(cData_attr, diag_kind='kde') # Esto nos permite ver densidades y las posibles relaciones entre dos variables
Out[12]:
<seaborn.axisgrid.PairGrid at 0x2aba2117108>
Nuestra variable objetivo es mpg la cual parece tener algun grado de correlación con horsepower y weight. Vamos a construir el modelo para verificarlo.
Hacemos split de los datos antes de ajustar el modelo de RL
In [13]:
## Primero las independientes en donde retiramos a mpg y a origin_europe (esta no la queremos en el modelo) X = cData.drop(['mpg','origin_europe'], axis=1) # y en Y nuestra variable objetivo o dependiente y = cData[['mpg']]
In [14]:
# Rompemos los datos en sets de entrenamiento y de prueba en relación 70:30 con estado random 1 para mantener esta estructura de azar en las iteraciones X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.30, random_state=1)
Ajustando el modelo de RL
In [15]:
# elegimos el modelo LR desde la sklearn y lo ajustaremos por el intercepto (mejor linea) sin normalizacion regression_model = LinearRegression() regression_model.fit(X_train, y_train)
Out[15]:
LinearRegression(copy_X=True, fit_intercept=True, n_jobs=None, normalize=False)
Ahora calcularemos los coeficientes de regresion de acuerdo con el intercepto:
In [16]:
for idx, col_name in enumerate(X_train.columns):
print("The coefficient for {} is {}".format(col_name, regression_model.coef_[0][idx]))
The coefficient for cylinders is -0.3948079661648244 The coefficient for displacement is 0.028945510765487174 The coefficient for horsepower is -0.02175220772354676 The coefficient for weight is -0.007352032065147351 The coefficient for acceleration is 0.061919366007618326 The coefficient for model year is 0.8369338917644993 The coefficient for origin_america is -3.0012830009185123 The coefficient for origin_asia is -0.6060179643247369
El año del modelo tiene una fuerte correlación positiva con mpg, La mayoria de variables restantes tienen una correlación negativa siendo el origen en América la de mayor correlación negativa. Ahora calcularemos el intercepto de todo el modelo
In [17]:
intercept = regression_model.intercept_[0]
print("The intercept for our model is {}".format(intercept))
The intercept for our model is -18.283451116372067
Ahora calcularemos una métrica de desempeño que es el R2:
In [18]:
regression_model.score(X_train, y_train)
Out[18]:
0.8141025501610559
el R2 es de 0.81 que es elevado lo cual quiere decir que el modelo acertará en la correlación 82% de las veces, ahora calcularemos este R2 por fuera de la muestra:
In [19]:
#out of sample score (R^2) regression_model.score(X_test, y_test)
Out[19]:
0.8433135132808833
Al ser positivo el R2 fuera de muestra nuestro modelo tiene menor error cuadrático promedio lo que quiere decir que es reproducible.
Adicionemos términos de interacción
Esto nos permitirá hacer ajustes polinomiales para optimizar variables y mejorar (posiblemente) el rendimiento del modelo
In [20]:
from sklearn.preprocessing import PolynomialFeatures from sklearn import linear_model poly = PolynomialFeatures(degree=2, interaction_only=True) X_train2 = poly.fit_transform(X_train) X_test2 = poly.fit_transform(X_test) poly_clf = linear_model.LinearRegression() poly_clf.fit(X_train2, y_train) y_pred = poly_clf.predict(X_test2) # veamos ahora nuestro nuevo R2 en muestra print(poly_clf.score(X_train2, y_train))
0.9015975294607972
In [21]:
# y el R2 fuera de muestra print(poly_clf.score(X_test2, y_test))
0.8647441061208313
Nuestros R2 mejoraron ampliamente, pero esto ha sido a costa de aumentar variables (de 8 a 37):
In [22]:
print(X_train.shape) print(X_train2.shape)
(278, 8) (278, 37)
Este ha sido un pequeño ejemplo de como puede trabajarse un modelo de RL en python versión 3.3. Gracias por su lectura
Luis Eduardo Pino Villarreal MD, MSc. MBA
Fundador de AIpocrates
CEO OxLER
luispino@oxler.me http://www.aipocrates.org http://www.oxler.info
aipocrates.org 