Análisis de los Modelos Lineales mediante la librería de Machine Learning (ML) de Scikit-Learn de Python
Boston Housing Data
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
import sklearn
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split
import sklearn.metrics as metrics
data = pd.read_csv('/work/Data/housing.data',header=None,sep='\s+')
columns = ['CRIM', 'ZN', 'INDUS', 'CHAS', 'NOX', 'RM','AGE','DIS','RAD','TAX','PTRATIO','B','LSTAT','MEDV']
data.columns = columns
data.head()
print(data.shape)
print("El dataset de precios de casas en Boston cuenta con 506 filas y 14 variables (columnas).")
print(data.size)
print(f'El dataset de precios de casas en Boston tiene en total {data.size} datos.')
data.isnull().sum()
print('El dataset de precios de casas en Boston no presenta valores nulos.')
Análisis estadístico del dataset
data.describe()
Con el objetivo de ser específico en el análisis, tomaré 5 variables del total que posiblemente expliquen el comportamiento de la variable Valor medio de las viviendas ocupadas por sus propietarios en miles de dólares. Las 5 variables que tomaré en cuenta son: CRIM, ZN, INDUS, DIS y TAX.
new_columns = ['CRIM', 'ZN', 'INDUS', 'DIS', 'TAX', 'MEDV']
sns.set_theme(style='ticks')
sns.pairplot(data[new_columns], height=1)
plt.show()
coef_matrix = np.corrcoef(data[new_columns].values.T)
fig, ax = plt.subplots(figsize=(10,8),dpi=80)
ax = sns.heatmap(coef_matrix,cbar=True,annot=True, yticklabels=new_columns, xticklabels=new_columns)
plt.title('Análisis de correlaciones entre las variables predictoras',size=15)
plt.show()
fig, ax = plt.subplots(figsize=(12,8),dpi=60)
ax = sns.histplot(data = data, x = data['MEDV'],kde=True, color='orange')
plt.title('Histograma del valor medio de las viviendas ocupadas por sus propietarios en miles de dólares \n', size=14)
for i in ['bottom', 'left']:
ax.spines[i].set_color('black')
ax.spines[i].set_linewidth(1.5)
right_side = ax.spines["right"]
right_side.set_visible(False)
top_side = ax.spines["top"]
top_side.set_visible(False)
ax.set_axisbelow(True)
ax.grid(color='gray', linewidth=1, axis='y', alpha=0.4)
plt.show()
fig, ax = plt.subplots(figsize=(9,5), dpi = 80)
ax = sns.boxplot(data=data, x=data['MEDV'], color='#038554')
plt.title('Diagrama de caja de la variable valor medio de las viviendas ocupadas \npor sus propietarios en miles de dólares \n', size=14)
for i in ['bottom', 'left']:
ax.spines[i].set_color('black')
ax.spines[i].set_linewidth(1.5)
right_side = ax.spines["right"]
right_side.set_visible(False)
top_side = ax.spines["top"]
top_side.set_visible(False)
ax.set_axisbelow(True)
ax.grid(color='gray', linewidth=1, axis='y', alpha=0.4)
plt.show()
fig, axs = plt.subplots(5,1, figsize=(10,14), dpi=60)
sns.boxplot(ax = axs[0], x=data['CRIM'], color='red').set(title='CRIM')
sns.boxplot(ax = axs[1], x=data['ZN'], color='green').set(title='ZN')
sns.boxplot(ax = axs[2], x=data['INDUS'], color='purple').set(title='INDUS')
sns.boxplot(ax = axs[3], x=data['DIS'], color='brown').set(title='DIS')
sns.boxplot(ax = axs[4], x=data['TAX'],color='coral').set(title='TAX')
plt.show()
Estandarización de los datos y Modelo de Regresión Lineal (MCO)
new_columns_x = ['CRIM', 'ZN', 'INDUS', 'DIS', 'TAX']
new_columns_y = ['MEDV']
x = data[new_columns_x].values
y = data[new_columns_y].values
X_train, X_test, y_train, y_test = train_test_split(x,y, test_size=0.20, random_state=42)
standardscaler_x = StandardScaler().fit(x)
standardscaler_y = StandardScaler().fit(y)
X_train = standardscaler_x.transform(X_train)
X_test = standardscaler_x.transform(X_test)
y_train = standardscaler_y.transform(y_train)
y_test = standardscaler_y.transform(y_test)
model = LinearRegression()
model.fit(X_train, y_train)
y_prediccion = model.predict(X_test)
y_prediccion.shape
print(X_train.mean(axis=0))
print('La media de los datos estandarizados fluctúa alrededor del 0.')
print(X_train.std(axis=0))
print('Mientras que la varianza fluctúa alrededor del 1.')
print(f'Los coeficientes del modelo de regresión lineal: {model.coef_}')
print(f'El intercepto del modelo de regresión lineal: {model.intercept_}')
from sklearn.feature_selection import f_regression
f_stat_p_value = f_regression(X_train,np.ravel(y_train), force_finite=True)
f_stat_p_value
print("="*15 + "F Estadístico y p-value de los coeficientes del modelo" + "="*15)
print(f'CRIM: {f_stat_p_value[0][0]} -- {f_stat_p_value[1][0]}' )
print(f'ZN: {f_stat_p_value[0][1]} -- {f_stat_p_value[1][1]}' )
print(f'INDUS: {f_stat_p_value[0][2]} -- {f_stat_p_value[1][2]}' )
print(f'DIS: {f_stat_p_value[0][3]} -- {f_stat_p_value[1][3]}' )
print(f'TAX: {f_stat_p_value[0][4]} -- {f_stat_p_value[1][4]}' )
print("="*84)
from sklearn.metrics import r2_score
from sklearn.metrics import mean_squared_error
from sklearn.metrics import mean_absolute_error
from sklearn.metrics import explained_variance_score #Best possible score is 1.0, lower values are worse.
r2 = r2_score(y_test,y_prediccion)
mse = mean_squared_error(y_test, y_prediccion)
mae = mean_absolute_error(y_test, y_prediccion)
evs = explained_variance_score(y_test, y_prediccion)
print("="*30 + "Métricas" + "="*30)
print(f'El coeficiente de determinación del modelo es: {round(r2,2)}')
print("="*69)
print(f'El error cuadrático medio del modelo es: {round(mse,2)}')
print("="*69)
print(f'El error absoluto medio del modelo es: {round(mae,2)}')
print("="*69)
print(f'La función de puntuación de regresión de varianza explicada es: {round(evs,2)}*')
print("*La mejor puntuación posible es 1.")
print("="*69)
fig, ax = plt.subplots(figsize=(10,8),dpi=65)
residuals = np.subtract(y_test,y_prediccion)
plt.scatter(y_prediccion,residuals)
plt.xlabel('y_prediccion')
plt.ylabel('residuales')
plt.axhline(y=0, color='red')
plt.title('Gráfico de dispersión de los valores predichos contra los residuos', size=14)
for i in ['bottom', 'left']:
ax.spines[i].set_color('black')
ax.spines[i].set_linewidth(1.5)
right_side = ax.spines["right"]
right_side.set_visible(False)
top_side = ax.spines["top"]
top_side.set_visible(False)
ax.set_axisbelow(True)
ax.grid(color='gray', linewidth=1, axis='y', alpha=0.4)
plt.show()
fig, ax = plt.subplots(figsize=(10,8),dpi=70)
plt.plot(y_prediccion, label='y_predicción')
plt.plot(y_test, label='y_test')
plt.ylim(-5,5)
plt.title('Ajuste del modelo de regresión lineal contra la data de test')
plt.legend()
for i in ['bottom', 'left']:
ax.spines[i].set_color('black')
ax.spines[i].set_linewidth(1.5)
right_side = ax.spines["right"]
right_side.set_visible(False)
top_side = ax.spines["top"]
top_side.set_visible(False)
ax.set_axisbelow(True)
ax.grid(color='gray', linewidth=1, axis='y', alpha=0.4)
plt.show()
Modelo lineal mediante el método Ridge
from sklearn import linear_model
n_alphas = 500
alphas = np.logspace(-10,10,n_alphas)
coeficients = []
for a in alphas:
ridge = linear_model.Ridge(alpha=a, fit_intercept=False)
ridge.fit(X_train,np.ravel(y_train))
coeficients.append(ridge.coef_)
ax = plt.gca()
ax.plot(alphas, coeficients)
ax.set_xscale("log")
ax.set_xlim(ax.get_xlim()[::-1])
plt.xlabel("alpha")
plt.ylabel("pesos")
plt.title("Coeficientes de Ridge en función de la regularización")
plt.axis("tight")
plt.show()
from sklearn.linear_model import Ridge
model_Ridge = Ridge(alpha=0.01).fit(X_train,y_train)
y_prediccion_Ridge = model_Ridge.predict(X_test)
y_prediccion_Ridge
r2_Ridge = r2_score(y_test,y_prediccion_Ridge)
mse_Ridge = mean_squared_error(y_test, y_prediccion_Ridge)
mae_Ridge = mean_absolute_error(y_test, y_prediccion_Ridge)
evs_Ridge = explained_variance_score(y_test, y_prediccion_Ridge)
print("="*30 + "Métricas" + "="*30)
print(f'El coeficiente de determinación del modelo es: {round(r2_Ridge,2)}')
print("="*69)
print(f'El error cuadrático medio del modelo es: {round(mse_Ridge,2)}')
print("="*69)
print(f'El error absoluto medio del modelo es: {round(mae_Ridge,2)}')
print("="*69)
print(f'La función de puntuación de regresión de varianza explicada es: {round(evs_Ridge,2)}*')
print("*La mejor puntuación posible es 1.")
print("="*69)
fig, ax = plt.subplots(figsize=(10,8),dpi=70)
plt.plot(y_prediccion_Ridge, label='y_predicción_Ridge')
plt.plot(y_test, label='y_test')
plt.ylim(-5,5)
plt.title('Ajuste del modelo de regresión lineal contra la data de test')
plt.legend()
for i in ['bottom', 'left']:
ax.spines[i].set_color('black')
ax.spines[i].set_linewidth(1.5)
right_side = ax.spines["right"]
right_side.set_visible(False)
top_side = ax.spines["top"]
top_side.set_visible(False)
ax.set_axisbelow(True)
ax.grid(color='gray', linewidth=1, axis='y', alpha=0.4)
plt.show()
Modelo lineal mediante el método Lasso
# Cálculo optimizado del hiperparámetro alpha para el método Lasso que permite reducir el MSE.
# En función del método de Lasso por descenso coordinado
from sklearn.linear_model import LassoCV
from sklearn.pipeline import make_pipeline
import time
start_time = time.time()
model = make_pipeline(StandardScaler(), LassoCV(cv=10)).fit(X_train, np.ravel(y_train))
fit_time = time.time() - start_time
ymin, ymax = 0, 1
lasso = model[-1]
plt.semilogx(lasso.alphas_, lasso.mse_path_, linestyle=":")
plt.plot(
lasso.alphas_,
lasso.mse_path_.mean(axis=-1),
color="black",
label="Promedio de los pliegues",
linewidth=2,
)
plt.axvline(lasso.alpha_, linestyle="--", color="black", label="alpha: CV estimado")
plt.ylim(ymin, ymax)
plt.xlabel(r"$\alpha$")
plt.ylabel("Error cuadrático medio")
plt.legend()
_ = plt.title(
f"Error cuadrático medio por cada pliegue (Tiempo de entrenamiento: {fit_time:.2f}s)"
)
lasso.alpha_
from sklearn.linear_model import Lasso
model_Lasso = Lasso(alpha=0.0036108728151931844).fit(X_train,y_train)
y_prediccion_Lasso = model_Lasso.predict(X_test)
y_prediccion_Lasso
r2_Lasso = r2_score(y_test,y_prediccion_Lasso)
mse_Lasso = mean_squared_error(y_test, y_prediccion_Lasso)
mae_Lasso = mean_absolute_error(y_test, y_prediccion_Lasso)
evs_Lasso = explained_variance_score(y_test, y_prediccion_Lasso)
print("="*30 + "Métricas" + "="*30)
print(f'El coeficiente de determinación del modelo es: {round(r2_Lasso,2)}')
print("="*69)
print(f'El error cuadrático medio del modelo es: {round(mse_Lasso,2)}')
print("="*69)
print(f'El error absoluto medio del modelo es: {round(mae_Lasso,2)}')
print("="*69)
print(f'La función de puntuación de regresión de varianza explicada es: {round(evs_Lasso,2)}*')
print("*La mejor puntuación posible es 1.")
print("="*69)
fig, ax = plt.subplots(figsize=(10,8),dpi=70)
plt.plot(y_prediccion_Lasso, label='y_predicción_Lasso')
plt.plot(y_test, label='y_test')
plt.ylim(-5,5)
plt.title('Ajuste del modelo de regresión lineal contra la data de test')
plt.legend()
for i in ['bottom', 'left']:
ax.spines[i].set_color('black')
ax.spines[i].set_linewidth(1.5)
right_side = ax.spines["right"]
right_side.set_visible(False)
top_side = ax.spines["top"]
top_side.set_visible(False)
ax.set_axisbelow(True)
ax.grid(color='gray', linewidth=1, axis='y', alpha=0.4)
plt.show()
Conclusiones
Al analizar la tabla de métricas de los tres métodos lineales desarrollados en el presente análisis técnico de Machine Learning, puedo concluir que tanto el método de Regresión Lineal (Mínimos Cuadrados Ordinarios) y el método Ridge minimizan la varianza de los datos (error) para explicar con mayor precisión el impacto de las variables exógenas o features sobre la variable endógena.
Cabe resaltar que mediante la métrica del Coeficiente de determinación se busca maximizar su valor, puesto que representa el ajuste de los datos predichos contra los datos reales. Además, la Puntuación de Varianza Explicada se encuentra en el rango de 0 y 1 donde se procura que el modelo se encuentre cerca del valor máximo. Por el contrario, debe establecerse como objetivo que el modelo desarrollado cuente con valores mínimos en el Error Cuadrático Medio y el Error Cuadrático Absoluto.