<jupyter_start><jupyter_text>Ridge Regression and Lasso RegressionTambién son llamados métodos de regularización. Se utilizan especialmente cuando el número de variables predictoras es mayor al número de observaciones. Con ellos:   1. Se logra que los coeficientes del modelo tiendan a cero.   2. Se minimiza el riesgo de _overfitting_.   3. Se reduce la varianza.   4. Se atenúa el efecto de la correlación entre las variables predictoras.   5. Se reduce la influencia en el modelo de las variables predictoras menos relevantes.**Ridge Regression:**- Se agrega una penalización equivalente al cuadrado de la magnitud de los coeficientes.- Se minimiza: $\displaystyle \sum_{i = 1}^{n} (y_{i} - \hat{y}_{i})^{2} + \lambda \sum_{j = 1}^{p} \beta_{j}^{2}$, donde$\displaystyle \sum_{i = 1}^{n} (y_{i} - \hat{y}_{i})^{2} = $ Suma del cuadrado de los residuales.$\displaystyle \sum_{j = 1}^{p} \beta_{j}^{2} = $ Se penaliza la suma de los coeficientes elevados al cuadrado.$\lambda = $ Determina el peso de la penalización._Notas:_- Si $\lambda = 0$ entonces se tiene el modelo de regresión lineal simple.- Si $\lambda \to \infty$ entonces la pendiente de la recta de los valores estimados tiende a cero. Las predicciones se hacen cada vez menos sensibles a las variables explicativas.- Se utiliza la validación cruzada (_Cross Validation_) para encontrar $\lambda$ que cumpla $0 < \lambda < \infty$. Con esto se determina cuál valor tiene la menor varianza.**Lasso Regression _(Least Absolute Shrinkage and Selection Operator)_:**- Se agrega una penalización equivalente al valor absoluto de la magnitud de los coeficientes.- Se minimiza: $\displaystyle \sum_{i = 1}^{n} (y_{i} - \hat{y}_{i})^{2} + \lambda \sum_{j = 1}^{p} |\beta_{j}|$_Nota:_- Si $\lambda \to \infty$ entonces la pendiente de la recta de los valores estimados es cero. Los pesos $\beta_{j}$ de las variables que no explican el modelo se hacen cero. Esto hace que la ecuación a minimizar sea más sencilla y más fácil de interpertar cuando se tienen muchas variables que no se utilizan.**Interpretación de los coeficientes:**Cada valor $\beta_{i}$ corresponde a una variable predictora. El valor indica qué tanto cambia **Y** por cada unidad que cambia **$X_{i}$** mientras las demás **$X_{j}$** (con $i \neq j$) permanecen constantes.Entre más grande sea el valor de $\beta_{i}$, más significativa es la _i-ésima_ variable predictora ($X_{i}$). La influencia puede ser positiva o negativa.Las variables con $\beta_{i} > 0$ y grandes, tienen mayor influencia en el modelo que variables con $\beta_{i}$ cercanos a cero.Las variables con $\beta_{i} < 0$ tienen un efecto negativo en las predicciones del modelo.**Ejemplos:**En este archivo se verán 2 ejemplos que muestran ambas regresiones.   - En _Ridge Regression_ se van a comparar 3 modelos        1. $\lambda = 5$        2. $\lambda = 0$        3. $\lambda$ elegida con _cross validation_   Para elegir el mejor de los 3 modelos vamos a calcular el error cuadrático medio: $ECM_{\lambda} = \dfrac{\displaystyle \sum_{i = 1}^{n} (real_{i} − estimado_{i})^{2}}{n}$.   Se elegirá el que tenga el menor $ECM_{\lambda}$      - En _LASSO Regression_ se va a definir el modelo utilizando $\lambda$ elegida con _cross validation_<jupyter_code># Se cargan las librerías que se van a utilizar en ambos ejemplos
import math
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import seaborn as sns
import sklearn

from sklearn.preprocessing import scale 
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error
from sklearn.metrics import r2_score
pd.set_option('display.max_columns', 25) # by default is 10, if change to None print ALL
## 1) EXTRAER DATOS
# Los datos pueden encontrarse en diferentes formatos, en nuestro caso están en formato csv.

# Se carga la base de datos
df = pd.read_csv('Hitters.csv') #Se encuentra en la misma carpeta que el jupyter notebook
print(df.shape) #Se tienen 322 renglones y 21 columnas
## 2) ANÁLISIS EXPLORATORIO
# Se realiza una descripción analítica de los datos.
df.head()# Se muestran los primeros 5 datos del data frame<jupyter_output><empty_output><jupyter_text>Los datos corresponden a la información de 322 jugadores de baseball en 1986 y 1987.Los términos utilizados se pueden encontrar en: https://www.mlb.com/glossary/La información que se tiene de cada uno de los jugadores es (21 columnas):   - **Player:** Nombre del jugador   - **AtBat:** Número de veces al bat en 1986   - **Hits:**  Número de hits en 1986 (cuando se le pega a la pelota y se llega al menos a la primera base)   - **HmRun:** Número de home runs en 1986   - **Runs:** Número de carreras en 1986   - **RBI:** Número de carreras impulsadas (carreras hechas por sus compañeros a causa suya) en 1986, (Runs Batted In)   - **Walks:** Número de bases por bolas en 1986   - **Years:** Número de años en ligas mayores   - **CAtBat:** Número de veces al bat durante su carrera   - **CHits:** Número de hits durante su carrera   - **CHmRun:** Número de home runs durante su carrera   - **CRuns:** Número de carreras durante su carrera   - **CRBI:** Número de carreras impulsadas durante su carrera   - **CWalks:** Número de bases por bolas durante su carrera   - **League:** Variable categórica que indica la liga del jugador al final de 1986 (A = American League -> Bateador designado batea en lugar del pitcher, N = National League -> Todos los jugadores batean)   - **Division:** Variable categórica que indica la división del jugador al final de 1986 (E = East, W = West)   - **PutOuts:** Número de outs en 1986   - **Assists:** Número de asistencias en 1986 (se acredita a un jugador defensivo que lanza o toca la pelota antes de un out)   - **Errors:** Número de errores en 1986   - **Salary:** Salario anual de 1987 en miles de dólares al día inaugural   - **NewLeague:** Factor con niveles A y N que indican la liga del jugador al inicio de 1987<jupyter_code># Se cuenta el número de NaN's por columna
df.isnull().sum() #Se tienen 59 NaN's en la columna del Salario
#En este caso vamos a eliminar los renglones de los jugadores con NaN en su salario.
#También vamos a eliminar la columna del nombre del jugador.
df = pd.read_csv('Hitters.csv').dropna().drop('Player', axis = 1)
print('Antes se tenían 322 renglones y 21 columnas.\nAhora se tienen ',df.shape[0],' renglones y ',df.shape[1],' columnas.')
print('Se eliminaron ',322-df.shape[0],' renglones.')
df.head()
# Se definen las varibles X (predictoras) y Y (dependiente)
Y = df.Salary
X = df.drop(['Salary'], axis = 1)#Se elimina la columna del salario
print(type(X))
X.head()
#Se definen las variables dummies para las variables categóricas
dummies = pd.get_dummies(X[['League', 'Division', 'NewLeague']])
dummies.head()
#Se eliminan las columnas 'League', 'Division', 'NewLeague' y se reemplazan por
#las variables dummies (una por cada columna).
#Recordemos los valores de cada variable:
# 'League' y 'NewLeague': A y N
# 'Division': E y W

X = X.drop(['League', 'Division', 'NewLeague'], axis = 1).astype('float64')
X = pd.concat([X, dummies[['League_N', 'Division_W', 'NewLeague_N']]], axis = 1)

print(X.shape[0],' renglones\n',X.shape[1],' columnas')
X.head()
# Se cuenta el número de NaN's por columna
df.isnull().sum() #No hay NaN's
#Se define un vector con las descripciones de cada columna en el data frame X:
descrip = ['Número de veces al bat en 1986',
           'Número de hits en 1986',
           'Número de home runs en 1986',
           'Número de carreras en 1986',
           'Número de carreras impulsadas en 1986',
           'Número de bases por bolas en 1986',
           'Número de años en ligas mayores',
           'Número de veces al bat durante su carrera',
           'Número de hits durante su carrera',
           'Número de home runs durante su carrera',
           'Número de carreras durante su carrera',
           'Número de carreras impulsadas durante su carrera',
           'Número de bases por bolas durante su carrera',
           'Número de outs en 1986',
           'Número de asistencias en 1986',
           'Número de errores en 1986',
           'Liga del jugador al final de 1986: 0 = American League, 1 = National League',
           'División del jugador al final de 1986: 0 = East, 1 = West',
           'Liga del jugador al inicio de 1987: 0 = American League, 1 = National League']
## 3) VISUALIZACIÓN DE LOS DATOS
# Para entender mejor los datos es necesario graficarlos.

#Histograma de la frecuencia relativa del salario con la densidad ajustada correspondiente.
sns.distplot(Y)
plt.title('Histograma de la frecuencia relativa del salario')
plt.xlabel('Salario')
plt.ylabel('Frecuencia relativa')
Y.describe()
#Histogramas de cada columna
for i in range(0,X.shape[1]): #[0,19)
    print('Mínimo = ',min(X.iloc[:,i]))
    print('Máximo = ',max(X.iloc[:,i]))
    plt.hist(X.iloc[:,i])
    plt.title('Histograma de la columna ' + X.columns[i])
    plt.xlabel(descrip[i])
    plt.ylabel('Número de jugadores')
    plt.show()
#Boxplots de cada columna
for i in range(0,(X.shape[1]-3)): #[0,16)
    print('Promedio = ',round(np.mean(X.iloc[:,i]),2))
    plt.title('Boxplot de la columna ' + X.columns[i])
    sns.boxplot(x=X.columns[i], data=X).set(xlabel=descrip[i])
    plt.show()
#No se muestran las gráficas de las columnas con variables dummies
#Se muestra la correlación entre las variables
sns.pairplot(df)
#Se muestra la correlación entre el salario y las variables predictivas
for i in [0,4,8,12]:
    sns.pairplot(data=df,
                  y_vars=['Salary'],
                  x_vars=['Salary',X.columns[i], X.columns[i+1], X.columns[i+2],X.columns[i+3]])
    plt.show()
#No se muestran las gráficas de las columnas con variables dummies
# Definimos el vector de lambdas, recordemos que el string 'lambda' está reservado para 
#las funciones con una sola instrucción.
#El vector tiene un gran rango de valores para poder cubrir varios escenarios.

print(np.linspace(10,-2,100)[:10])

lambdas = 10**np.linspace(10,-2,100)*0.5
'''
Con np.linspace(10,-2,100) se define un vector que va de 10 a -2, tiene 100 entradas .
Se eleva 10 a cada uno de esos valores.
Finalmente, cada una de esas entradas se multiplica por 0.5.
'''
print(lambdas[:10]) #Se muestran los primeros 10 valores del vector
print(lambdas[-10:]) #Se muestran los últimos 10 valores del vector
## 4) DIVIDIR LOS DATOS
# Se separan los datos en 2 grupos (usualmente 80% y 20%):
# i) Para entrenar al modelo (80%)
# ii) Para probar el modelo (20%)
X_train, X_test, Y_train, Y_test = train_test_split(X, Y, #Se indican los vectores que se van a dividir
                                                    test_size = 0.2, #Se indica el porcentaje de los datos para probar el modelo
                                                    random_state = 1) #Se fija la semilla

# Nota: Tomar la muestra aleatoria es muy importante porque en caso de que los datos estén
#ordenados el algoritmo no aprende adecuadamente. Por ejemplo si tenemos 80 sanos y 20 enfermos,
#sólo se tomarían los 80 sanos (por ser los primeros 80).<jupyter_output><empty_output><jupyter_text>Ridge Regression<jupyter_code># Se cargan las librerías que se van a utilizar
from sklearn.linear_model import Ridge, RidgeCV
#Se definen las variables que se van a utilizar
ridge = Ridge(normalize = True)#Para estandarizar las variables el parámetro "normalize" es verdadero
coefs = []

#Se llena la matriz con coeficientes asociadas a cada variable independiente y a cada valor de lambda
for k in lambdas:
    ridge.set_params(alpha = k)#Se toman los valores de lambda
    ridge.fit(X_train, Y_train)#Se ajusta el modelo
    coefs.append(ridge.coef_)#Se generan los coeficientes correspondientes
    
print(np.shape(coefs))
coefs[0]<jupyter_output><empty_output><jupyter_text>Esperamos que los coeficientes estimados sean cada vez más pequeños mientras lambda se hace mas grande.Se espera que los coeficientes estimados sean más grandes si lambda es pequeña.Si $\lambda \to \infty \,\, \Rightarrow \beta \to 0$ porque se quiere minimizar la función correspondiente.Entre menos pronunciada esté una pendiente implica que a cambios pequeños en **X**, se tienen cambios pequeños en **Y**.Al contrario si se tiene una pendiente pronunciada, entonces a cambios pequeños en **X**, se tienen cambios grandes en **Y**.<jupyter_code>#Grafiquemos
ax = plt.gca()
ax.plot(lambdas, coefs)
ax.set_xscale('log')
plt.axis('tight')
plt.xlabel('lambda')
plt.ylabel('beta')

#En esta gráfica vemos que entre más grande es lambda, más cercano a cero es beta<jupyter_output><empty_output><jupyter_text>$\lambda = 5$<jupyter_code>## 5) CONSTRUIR UN MODELO
mod_ridge5 = Ridge(alpha = 5, normalize = True)
mod_ridge5.fit(X_train, Y_train)
print(pd.Series(mod_ridge5.coef_, index = X.columns)) # Valor de los coeficientes
# Se grafican los valores de los coeficientes (betas).
eje_X = range(1,len(mod_ridge5.coef_)+1,1)
plt.scatter(eje_X, mod_ridge5.coef_, color = 'purple')
plt.title('Betas')
plt.xlabel('Índices')
plt.ylabel('Valor de coeficientes')
plt.show()<jupyter_output><empty_output><jupyter_text>Interpretación de los coeficientesCon la gráfica anterior podemos ver que la mayoría de los coeficientes están cercanos a cero.Esto nos indica que la mayoría de las variables no tienen una gran influencia en el modelo.Las variables con $\beta_{i} > 1$ son:|    Variable   | $\beta$ ||:-----------:|:-------:||    HmRun    |  1.40  ||    Years    |  2.32  ||   League_N  |  3.17  || NewLeague_N |  4.09  |Estas variables son las que más impactan en el sueldo del jugador.El coeficiente de la variable _Division_W_ es $-24.628$. Esta variable impacta de manera negativa en el modelo.Podemos concluir que jugar en la división del oeste implica tener un menor sueldo.<jupyter_code>## 6) PREDICCIONES
# Se hacen las predicciones con "X_test"
Y_pred = mod_ridge5.predict(X_test)
#print(X_test) #data frame
#print(Y_pred)
# Se grafican los resultados de la predicción.
plt.scatter(Y_test, Y_pred, color = 'blue')
plt.title('Predicciones')
plt.xlabel('Salario real')
plt.ylabel('Salario estimado')
plt.show()

# Nota: No estamos graficando contra ninguna variable explicativa.
#Los valores de las predicciones se graficaron contra el salario real.
## 7) EVALUACIÓN DEL MODELO CON lambda = 5
#Calculamos el error cuadrático medio: $\dfrac{\displaystyle \sum_{i = 1}^{n} (real_{i} − estimado_{i})^{2}}{n}$
ECM5 = round(mean_squared_error(Y_test, Y_pred),3) 
print(ECM5)<jupyter_output><empty_output><jupyter_text>$\lambda = 0$<jupyter_code>## 5) CONSTRUIR UN MODELO
mod_ridge0 = Ridge(alpha = 0, normalize = True)
mod_ridge0.fit(X_train, Y_train)
print(pd.Series(mod_ridge0.coef_, index = X.columns)) # Valor de los coeficientes
# Se grafican los valores de los coeficientes (betas).
eje_X = range(1,len(mod_ridge0.coef_)+1,1)
plt.scatter(eje_X, mod_ridge0.coef_, color = 'purple')
plt.title('Betas')
plt.xlabel('Índices')
plt.ylabel('Valor de coeficientes')
plt.show()<jupyter_output><empty_output><jupyter_text>Interpretación de los coeficientesCon la gráfica anterior podemos ver que la mayoría de los coeficientes están cercanos a cero.Esto nos indica que la mayoría de las variables no tienen una gran influencia en el modelo.Las variables con $\beta_{i} \notin [-1,1]$ son:|   Variable  | $\beta$ ||:-----------:|:-------:||    AtBat    |  -1.87  ||     Hits    |   6.58  ||    Walks    |   4.56  ||    Years    |   6.14  ||    CHmRun   |   1.59  ||    Errors   |   1.65  ||   League_N  |  51.95  ||  Division_W | -116.57 || NewLeague_N |  -3.25  |En este caso $\lambda = 0$ por lo que no importa el valor de los coeficientes ya que al multiplicarse por $\lambda$ sus valores se anulan.<jupyter_code>## 6) PREDICCIONES
# Se hacen las predicciones con "X_test"
Y_pred = mod_ridge0.predict(X_test)
# Se grafican los resultados de la predicción.
plt.scatter(Y_test, Y_pred, color = 'blue')
plt.title('Predicciones')
plt.xlabel('Salario real')
plt.ylabel('Salario estimado')
plt.show()

# Nota: No estamos graficando contra ninguna variable explicativa.
#Los valores de las predicciones se graficaron contra el salario real.
## 7) EVALUACIÓN DEL MODELO CON lambda = 0
#Calculamos el error cuadrático medio: $\dfrac{\displaystyle \sum_{i = 1}^{n} (real_{i} − estimado_{i})^{2}}{n}$
ECM0 = round(mean_squared_error(Y_test, Y_pred),3) 
print(ECM0)<jupyter_output><empty_output><jupyter_text>$\lambda$ definida con _cross validation_<jupyter_code>## 5) CONSTRUIR UN MODELO
#Le pasamos como parámetro el vector con los diferentes valores de lambda para que
#se elija la mejor con "cross validation"
ridgecv = RidgeCV(alphas = lambdas, scoring = 'neg_mean_squared_error', normalize = True)
ridgecv.fit(X_train, Y_train)
print('El valor de lambda encontrado con "cross validation" es: ',round(ridgecv.alpha_,3))
#Verificamos que el valor encontrado está en el vector de lambdas
ridgecv.alpha_ in lambdas
#Se imprimen los valores de los coeficientes
mod_ridgeCV = Ridge(alpha = ridgecv.alpha_, normalize = True)
mod_ridgeCV.fit(X_train, Y_train)
print(pd.Series(mod_ridgeCV.coef_, index = X.columns)) # Valor de los coeficientes
# Se grafican los valores de los coeficientes (betas).
eje_X = range(1,len(mod_ridgeCV.coef_)+1,1)
plt.scatter(eje_X, mod_ridgeCV.coef_, color = 'purple')
plt.title('Betas')
plt.xlabel('Índices')
plt.ylabel('Valor de coeficientes')
plt.show()<jupyter_output><empty_output><jupyter_text>Interpretación de los coeficientesCon la gráfica anterior podemos ver que la mayoría de los coeficientes están cercanos a cero.Esto nos indica que la mayoría de las variables no tienen una gran influencia en el modelo.Las variables con $\beta_{i} \notin [-1,1]$ son:|   Variable  | $\beta$ > 1 |  Variable  | $\beta$ < -1 ||:-----------:|:-----------:|:----------:|:------------:||     Hits    |     5.74    |    AtBat   |     -1.69    ||    Walks    |     4.20    |   HmRun   |     -2.40    ||    CHmRun   |     1.19    |    Years   |     -6.62    ||    Errors   |     1.70    | Division_W |    -117.51   ||   League_N  |    39.18    |      -     |       -      || NewLeague_N |     4.62    |      -     |       -      |Las variables con $\beta > 1$ son las que más aportan en el sueldo del jugador.El coeficiente de la variable _League_N_ es $39.18$. Podemos concluir que jugar en la liga nacional implica tener un mayor sueldo.Las variables con $\beta < -1$ son las que impactan de manera negativa en el sueldo del jugador.El coeficiente de la variable _Division_W_ es $-117.51$. Podemos concluir que jugar en la división del oeste implica tener un menor sueldo.<jupyter_code>## 6) PREDICCIONES
# Se hacen las predicciones con "X_test"
Y_pred = mod_ridgeCV.predict(X_test)
# Se grafican los resultados de la predicción.
plt.scatter(Y_test, Y_pred, color = 'blue')
plt.title('Predicciones')
plt.xlabel('Salario real')
plt.ylabel('Salario estimado')
plt.show()

# Nota: No estamos graficando contra ninguna variable explicativa.
#Los valores de las predicciones se graficaron contra el salario real.
## 7) EVALUACIÓN DEL MODELO CON lambda ELEGIDA CON "cross validation"
#Calculamos el error cuadrático medio: $\dfrac{\displaystyle \sum_{i = 1}^{n} (real_{i} − estimado_{i})^{2}}{n}$
ECMcv = round(mean_squared_error(Y_test, Y_pred),3) 
print(ECMcv)
# Comparamos el error cuadrático medio (ECM) de cada modelo:
print("ECM con lambda = 5 es: ",ECM5)
print("ECM con lambda = 0 es: ",ECM0)
print("ECM con lambda elegida es: ",ECMcv)<jupyter_output><empty_output><jupyter_text>Notamos que el error cuadrático medio menor es el del modelo en el que se elige $\lambda$ con _cross validation_.Elegimos ese modelo para seguir la evaluación. 7) EVALUACIÓN DEL MODELO ELEGIDOVeamos cómo se comporta el modelo:7.1 Calcular $R^{2}$ ajustada $ = 1 - \dfrac{(1 - R^{2}) (n-1)}{n - p - 1}$, donde$R^{2}:$ R cuadrada de los datos$n:$ Número de datos para entrenar al modelo$p:$ Número de variables independientes7.2 Calcular los errores absolutos $(real - estimado)$ y graficarlos7.3 Calcular los errores relativos $\left( \dfrac{\text{real - estimado}}{\text{real}} \right)$ y graficarlos7.4 Graficar valores estimados vs valores reales7.5 Calcular el error cuadrático: $(real − estimado)^{2}$<jupyter_code>#7.1 Calcular R^2 ajustada
r_cuadrada = r2_score(Y_test,Y_pred)
print('R^2 = ',round(r_cuadrada,3)) #Porcentaje de los datos explicados por el modelo

#R^2 ajustada
n = len(Y_train)
p = X_train.shape[1]
r_cuad_aj = 1 - (((1-r_cuadrada)*(n-1))/(n-p-1))

print('n = ',n)
print('p = ',p)
print('R^2_aj = ',round(r_cuad_aj,3))<jupyter_output><empty_output><jupyter_text>Interpretación de $R^{2}$ y $R_{aj}^{2}$Con $R^{2} =  0.44$ y $R_{aj}^{2} =  0.39$ vemos que el porcentaje de los datos explicados por el modelo es menor al $50\%$.Recordemos que en el análisis de los valores de los coeficientes vimos que la mayoría de las variables no eran significativas para el modelo.<jupyter_code>#7.2 Calcular los errores absolutos (real - estimado) y graficarlos
err_abs = Y_test-Y_pred

plt.scatter(Y_test, err_abs, color = 'blue')
plt.plot(Y_test, np.zeros(len(err_abs)), color = 'midnightblue') #Recta en Y = 0
plt.title('Errores absolutos (real - estimado)')
plt.xlabel('Salario real')
plt.ylabel('Errores absolutos')
plt.show()
#7.3 Calcular los errores relativos [(real - estimado)/real] y graficarlos
err_rel = err_abs/Y_test

plt.scatter(Y_test, err_rel, color = 'blue')
plt.plot(Y_test, np.zeros(len(err_abs)), color = 'midnightblue') #Recta en Y = 0
plt.title('Errores relativos [(real - estimado)/real]')
plt.xlabel('Salario real')
plt.ylabel('Errores relativos')
plt.show()

#Se tiene un mayor error en salarios bajos
#7.4 Graficar valores estimados vs valores reales
eje_X = range(1,len(Y_test)+1)
plt.plot(eje_X, Y_test, color = 'blue') #Recta de valores reales
plt.plot(eje_X, Y_pred, color = 'green') #Recta de valores estimados
plt.title('Valores estimados vs valores reales')
plt.xlabel('Índices')
plt.ylabel('Salario')
plt.show()<jupyter_output><empty_output><jupyter_text>En la siguiente gráfica vamos a ordenar de menor a mayor los vectores para poder apreciar de mejor manera el comportamiento que se tiene.<jupyter_code>#Los ordenamos de menor a mayor
eje_X = range(1,len(Y_test)+1)
plt.plot(eje_X, sorted(Y_test), color = 'blue') #Recta de valores reales
plt.plot(eje_X, sorted(Y_pred), color = 'green') #Recta de valores estimados
plt.title('Valores estimados vs valores reales')
plt.xlabel('Índices')
plt.ylabel('Salario')
plt.show()
#7.5 Calcular el error cuadrático = (real − estimado)^2
#print(np.around(err_abs,2))
err_cuad = pow(err_abs,2)#Vector

print("ECM con lambda = ",round(ridgecv.alpha_,3)," es: ",ECMcv)

#Graficamos los errores cuadráticos
Y= np.repeat(ECMcv, len(err_cuad))
plt.scatter(Y_test, err_cuad, color = 'blue')
plt.plot(Y_test,Y , color = 'lime') #Recta en Y = error cuadrático medio
plt.title('Errores cuadráticos: (real − estimado)^2')
plt.xlabel('Ganancia real')
plt.ylabel('Errores cuadráticos')
plt.show()<jupyter_output><empty_output><jupyter_text>Lasso Regression<jupyter_code># Se cargan las librerías que se van a utilizar
from sklearn.linear_model import Lasso, LassoCV
#Para estandarizar las variables el parámetro "normalize" es verdadero.
#Se define el máximo número de iteraciones en 10 mil
lasso = Lasso(max_iter = 10000, normalize = True)
coefs = []

for k in lambdas:
    lasso.set_params(alpha = k)
    lasso.fit(scale(X_train), Y_train)
    coefs.append(lasso.coef_)

np.shape(coefs)#Matriz con coeficientes asociadas a cada variable independiente y a cada valor de lambda
# Esperamos que los coeficientes estimados sean cada vez más pequeños mientras lambda se hace mas grande.
#Se espera que los coeficientes estimados sean más grandes si lambda es pequeña.

#Grafiquemos
ax = plt.gca()
ax.plot(lambdas, coefs)
ax.set_xscale('log')
plt.axis('tight')
plt.xlabel('lambda')
plt.ylabel('pesos')
## 5) CONSTRUIR UN MODELO
# Vamos a elegir lambda con "cross validation"
lassocv = LassoCV(alphas = None, cv = 10, max_iter = 100000, normalize = True)
lassocv.fit(X_train, Y_train)
#Se imprimen los valores de los coeficientes
print(pd.Series(lassocv.coef_, index = X.columns)) # Valor de los coeficientes
# Se grafican los valores de los coeficientes (betas).
eje_X = range(1,len(lassocv.coef_)+1,1)
plt.scatter(eje_X, lassocv.coef_, color = 'purple')
plt.title('Betas')
plt.xlabel('Índices')
plt.ylabel('Valor de coeficientes')
plt.show()<jupyter_output><empty_output><jupyter_text>Interpretación de los coeficientesAl imprimir los valores de los coeficientes vemos que 6 variables son cero. Las variables que no aportan nada al modelo son:   1. HmRun: Número de home runs en 1986   2. Runs: Número de carreras en 1986   3. RBI: Número de carreras impulsadas en 1986, (Runs Batted In)   4. CAtBat: Número de veces al bat durante su carrera   5. CHits: Número de hits durante su carrera   6. Assists: Número de asistencias en 1986Con la gráfica anterior podemos ver que la mayoría de los coeficientes son o están cercanos a cero. Esto nos indica que la mayoría de las variables no tienen una gran influencia en el modelo.Las variables con $\beta_{i} \notin [-1,1]$ son:|   Variable  | $\beta$ > 1 |  Variable  | $\beta$ < -1 ||:-----------:|:-----------:|:----------:|:------------:||     Hits    |     6.11    |    AtBat   |     -1.55    ||    Walks    |     3.66    |    Years   |     -5.71    ||    CHmRun   |    1.006    | Division_W |    -113.89   ||    Errors   |     1.23    |      -     |       -      ||   League_N  |    28.51    |      -     |       -      || NewLeague_N |    11.22    |      -     |       -      | **¿Qué pasa con _League_N, Division_W_ y _NewLeague_N_?** Las variables con $\beta > 1$ son las que más aportan en el sueldo del jugador.El coeficiente de la variable _League_N_ es $28.51$ y el de _NewLeague_N_ es $11.22$. Podemos concluir que jugar en la liga nacional implica tener un mayor sueldo.Las variables con $\beta < -1$ son las que impactan de manera negativa en el sueldo del jugador.El coeficiente de la variable _Division_W_ es $-113.89$. Podemos concluir que jugar en la división del oeste implica tener un menor sueldo.<jupyter_code>## 6) PREDICCIONES
# Se hacen las predicciones con "X_test"
lasso.set_params(alpha=lassocv.alpha_)
print('El valor de lambda encontrado con "cross validation" es: ',round(lassocv.alpha_,3))
lasso.fit(X_train, Y_train)
Y_pred = lasso.predict(X_test)

##Recordamos que el lambda encotrado en Ridge Regression fue de 0.012
# Se grafican los resultados de la predicción.
plt.scatter(Y_test, Y_pred, color = 'blue')
plt.title('Predicciones')
plt.xlabel('Salario real')
plt.ylabel('Salario estimado')
plt.show()

# Nota: No estamos graficando contra ninguna variable explicativa.
#Los valores de las predicciones se graficaron contra el salario real.<jupyter_output><empty_output><jupyter_text>7) EVALUACIÓN DEL MODELOVeamos cómo se comporta el modelo:7.1 Calcular $R^{2}$ ajustada $ = 1 - \dfrac{(1 - R^{2}) (n-1)}{n - p - 1}$, donde$R^{2}:$ R cuadrada de los datos$n:$ Número de datos para entrenar al modelo$p:$ Número de variables independientes7.2 Calcular los errores absolutos $(real - estimado)$ y graficarlos7.3 Calcular los errores relativos $\left( \dfrac{\text{real - estimado}}{\text{real}} \right)$ y graficarlos7.4 Graficar valores estimados vs valores reales7.5 Calcular el error cuadrático: $(real − estimado)^{2}$7.6 Calcular el error cuadrático medio: $\dfrac{\displaystyle \sum_{i = 1}^{n} (real_{i} − estimado_{i})^{2}}{n}$<jupyter_code>#7.1 Calcular R^2 ajustada
r_cuadrada = r2_score(Y_test,Y_pred)
print('R^2 = ',round(r_cuadrada,3)) #Porcentaje de los datos explicados por el modelo

#R^2 ajustada
n = len(Y_train)
p = X_train.shape[1]
r_cuad_aj = 1 - (((1-r_cuadrada)*(n-1))/(n-p-1))

print('n = ',n)
print('p = ',p)
print('R^2_aj = ',round(r_cuad_aj,3))
#7.1 Calcular R^2 ajustada
r_cuadrada = r2_score(Y_test,Y_pred)
print('R^2 = ',round(r_cuadrada,3)) #Porcentaje de los datos explicados por el modelo

#R^2 ajustada
n = len(Y_train)
p = 13 #Número de variables != 0
r_cuad_aj = 1 - (((1-r_cuadrada)*(n-1))/(n-p-1))

print('n = ',n)
print('p = ',p)
print('R^2_aj = ',round(r_cuad_aj,3))
len(Y_train)
X_train.shape[1]<jupyter_output><empty_output><jupyter_text>Interpretación de $R^{2}$ y $R_{aj}^{2}$Con $R^{2} =  0.45$ y $R_{aj}^{2} =  0.4$ vemos que el procentaje de los datos explicados por el modelo es menor al $50\%$.Notamos que el porcentaje es mayor comparado con los valores obtenidos en _Ridge Regression_.Recordemos que nuestros datos tienen 19 variables predictoras. En el análisis de los valores de los coeficientes vimos que, de dichas variables:   - 6 tienen coeficientes iguales a cero   - 6 variables tienen coeficientes mayores a 1   - 3 tienen coeficientes menores a -1   Las 4 restantes tienen coeficientes entre -1 y 1, por lo que decimos que no son significativas para el modelo.<jupyter_code>#7.2 Calcular los errores absolutos (real - estimado) y graficarlos
err_abs = Y_test-Y_pred

plt.scatter(Y_test, err_abs, color = 'blue')
plt.plot(Y_test, np.zeros(len(err_abs)), color = 'midnightblue') #Recta en Y = 0
plt.title('Errores absolutos (real - estimado)')
plt.xlabel('Salario real')
plt.ylabel('Errores absolutos')
plt.show()
#7.3 Calcular los errores relativos [(real - estimado)/real] y graficarlos
err_rel = err_abs/Y_test

plt.scatter(Y_test, err_rel, color = 'blue')
plt.plot(Y_test, np.zeros(len(err_abs)), color = 'midnightblue') #Recta en Y = 0
plt.title('Errores relativos [(real - estimado)/real]')
plt.xlabel('Salario real')
plt.ylabel('Errores relativos')
plt.show()

#Se tiene un mayor error en salarios bajos
#7.4 Graficar valores estimados vs valores reales
eje_X = range(1,len(Y_test)+1)
plt.plot(eje_X, Y_test, color = 'blue') #Recta de valores reales
plt.plot(eje_X, Y_pred, color = 'green') #Recta de valores estimados
plt.title('Valores estimados vs valores reales')
plt.xlabel('Índices')
plt.ylabel('Salario')
plt.show()
#Los ordenamos de menor a mayor
eje_X = range(1,len(Y_test)+1)
plt.plot(eje_X, sorted(Y_test), color = 'blue') #Recta de valores reales
plt.plot(eje_X, sorted(Y_pred), color = 'green') #Recta de valores estimados
plt.title('Valores estimados vs valores reales')
plt.xlabel('Índices')
plt.ylabel('Salario')
plt.show()
#7.5 Calcular el error cuadrático = (real − estimado)^2
#print(np.around(err_abs,2))
err_cuad = pow(err_abs,2)#Vector

#7.6 Calcular el error cuadrático medio = (1/n) * \sum (real − estimado)^2
'''
Indica qué tan cerca está la línea de la regresión lineal de los valores estimados.
i) Se elevan al cuadrado los errores absolutos.
ii) Se suman.
iii) Se divide el resultado entre el número de datos estimados.
'''
ECMlasso = round(mean_squared_error(Y_test, Y_pred),3)
print('Error cuadrático medio = ',round(ECMlasso,2))
print('Raíz cuadrada del error cuadrático medio = ',round(math.sqrt(ECMlasso),2))#Raíz cuadrada del error cuadrático medio

#Graficamos los errores cuadráticos
Y= np.repeat(ECMcv, len(err_cuad))
plt.scatter(Y_test, err_cuad, color = 'blue')
plt.plot(Y_test,Y , color = 'lime') #Recta en Y = error cuadrático medio
plt.title('Errores cuadráticos: (real − estimado)^2')
plt.xlabel('Ganancia real')
plt.ylabel('Errores cuadráticos')
plt.show()
# Comparamos el error cuadrático medio (ECM) de ambas regresiones:
# La comparación se hace tomando el modelo hecho con la lambda elegida con "cross validation"
if ECMlasso < ECMcv:
    print("ECM LASSO = ",ECMlasso," < ",ECMcv," = ECM Ridge \nElegir LASSO Regression")
elif ECMlasso > ECMcv:
    print("ECM LASSO = ",ECMlasso," > ",ECMcv," = ECM Ridge \nElegir Ridge Regression")
else:
    print("Ambas regresiones tienen el mismo ECM = ",ECMlasso,"\nElegir cualquier regresión")<jupyter_output><empty_output><jupyter_text>ConclusionesAl comparar ambos modelos con el error cuadrático medio, se elige _LASSO Regression_.Los coeficientes de cada modelo son muy similares. Podemos observar que _HmRun_ es la variable en la que más se difiere ya en en el caso de _Ridge Regression_ $\beta = -2.4$ y para _LASSO Regression_ $\beta = 0$.En _LASSO Regression_ se hacen cero las variables que tienen coeficientes casi cero en _Ridge Regression_._Nota:_ En Ciencia de Datos, lo que importa son las decisiones que se toman con respecto a los datos.  ¿Qué pasa si quitamos las variables que LASSO hace cero?<jupyter_code>#Definimos el nuevo data frame
X = X.drop(['HmRun','Runs','RBI','CAtBat','CHits','Assists'], axis = 1)
print(type(X))
print(X.shape)
X.head()
#División de los datos
X_train, X_test, Y_train, Y_test = train_test_split(X, Y,test_size = 0.2,random_state = 1)
#Se definen las variables que se van a utilizar
ridge = Ridge(normalize = True)#Para estandarizar las variables el parámetro "normalize" es verdadero
coefs = []

#Se llena la matriz con coeficientes asociadas a cada variable independiente y a cada valor de lambda
for k in lambdas:
    ridge.set_params(alpha = k)#Se toman los valores de lambda
    ridge.fit(X_train, Y_train)#Se ajusta el modelo
    coefs.append(ridge.coef_)#Se generan los coeficientes correspondientes
## 5) CONSTRUIR UN MODELO
#Le pasamos como parámetro el vector con los diferentes valores de lambda para que
#se elija la mejor con "cross validation"
ridgecv = RidgeCV(alphas = lambdas, scoring = 'neg_mean_squared_error', normalize = True)
ridgecv.fit(X_train, Y_train)
print('El valor de lambda encontrado con "cross validation" es: ',round(ridgecv.alpha_,4))
#Se imprimen los valores de los coeficientes
mod_ridgeCV = Ridge(alpha = ridgecv.alpha_, normalize = True)
mod_ridgeCV.fit(X_train, Y_train)
print(pd.Series(mod_ridgeCV.coef_, index = X.columns)) # Valor de los coeficientes
# Se grafican los valores de los coeficientes (betas).
eje_X = range(1,len(mod_ridgeCV.coef_)+1,1)
plt.scatter(eje_X, mod_ridgeCV.coef_, color = 'purple')
plt.title('Betas')
plt.xlabel('Índices')
plt.ylabel('Valor de coeficientes')
plt.show()
## 6) PREDICCIONES
# Se hacen las predicciones con "X_test"
Y_pred = mod_ridgeCV.predict(X_test)
# Se grafican los resultados de la predicción.
plt.scatter(Y_test, Y_pred, color = 'blue')
plt.title('Predicciones')
plt.xlabel('Salario real')
plt.ylabel('Salario estimado')
plt.show()
## 7) EVALUACIÓN DEL MODELO CON lambda ELEGIDA CON "cross validation"
#Calculamos el error cuadrático medio: $\dfrac{\displaystyle \sum_{i = 1}^{n} (real_{i} − estimado_{i})^{2}}{n}$
ECM = round(mean_squared_error(Y_test, Y_pred),3) 
print(ECM)

#7.1 Calcular R^2 ajustada
r_cuadrada = r2_score(Y_test,Y_pred)
print('R^2 = ',round(r_cuadrada,3)) #Porcentaje de los datos explicados por el modelo

#R^2 ajustada
n = len(Y_train)
p = X_train.shape[1]
r_cuad_aj = 1 - (((1-r_cuadrada)*(n-1))/(n-p-1))

print('n = ',n)
print('p = ',p)
print('R^2_aj = ',round(r_cuad_aj,3))<jupyter_output><empty_output><jupyter_text>Vemos que hay una mejora en la $R^{2}_{ajustada}$, comparada con el modelo con las 19 variables.  LASSO<jupyter_code>## 5) CONSTRUIR UN MODELO
# Vamos a elegir lambda con "cross validation"
lassocv = LassoCV(alphas = None, cv = 10, max_iter = 100000, normalize = True)
lassocv.fit(X_train, Y_train)

# Se grafican los valores de los coeficientes (betas).
eje_X = range(1,len(lassocv.coef_)+1,1)
plt.scatter(eje_X, lassocv.coef_, color = 'purple')
plt.title('Betas')
plt.xlabel('Índices')
plt.ylabel('Valor de coeficientes')
plt.show()<jupyter_output><empty_output><jupyter_text>En la función LassoCV(), se definen las $\alpha$ como _None_ (alphas = None), con este comando la función define automáticamente los valores de $\alpha$ que se van a utilizar. Ver [sklearn.linear_model.LassoCV](https://scikit-learn.org/stable/modules/generated/sklearn.linear_model.LassoCV.html).<jupyter_code>#Se imprimen los valores de los coeficientes
print(pd.Series(lassocv.coef_, index = X.columns)) # Valor de los coeficientes<jupyter_output><empty_output><jupyter_text>Ninguna variable tiene coeficiente igual a cero.<jupyter_code>## 6) PREDICCIONES
# Se hacen las predicciones con "X_test"
lasso.set_params(alpha=lassocv.alpha_)
print('El valor de lambda encontrado con "cross validation" es: ',round(lassocv.alpha_,3))
lasso.fit(X_train, Y_train)
Y_pred = lasso.predict(X_test)
# Se grafican los resultados de la predicción.
plt.scatter(Y_test, Y_pred, color = 'blue')
plt.title('Predicciones')
plt.xlabel('Salario real')
plt.ylabel('Salario estimado')
plt.show()
## 7) EVALUACIÓN DEL MODELO
#7.1 Calcular R^2 ajustada
r_cuadrada = r2_score(Y_test,Y_pred)
print('R^2 = ',round(r_cuadrada,3)) #Porcentaje de los datos explicados por el modelo

#R^2 ajustada
n = len(Y_train)
p = X_train.shape[1]
r_cuad_aj = 1 - (((1-r_cuadrada)*(n-1))/(n-p-1))

print('n = ',n)
print('p = ',p)
print('R^2_aj = ',round(r_cuad_aj,3))<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.11807504700656729, 'fraction_numerical': 0.020383137586178705, 'mean_word_length': 4.502774029089818, 'pattern_counts': {'":': 0, '<': 94, '<?xml version=': 0, '>': 95, 'https://': 2, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 7, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '21', 'dir': 'jupyter-structured-clean-dedup', 'id': '96286', 'license': 'MIT', 'n_tokens_mistral': 13909, 'n_tokens_neox': 13338, 'n_words': 4281, 'path': 'M2/2_4 Ridge and Lasso Regression.ipynb', 'repo_name': 'Pepoide45/Curso_Ciencia_de_Datos', 'size': '36697'}