## Libraries
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
from scipy.stats import spearmanr
import sklearn.linear_model as lm
from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split
from sklearn import metrics
import matplotlib.dates as mdate
import datetime as dt
%matplotlib inline
Cargar y normalizar datos meteorológicos¶
Fuente de datos climáticos:
NASA (POWER Data Access Viewer)
Fuente de datos consumo de gas:
ENARGAS
# Load data
filepath = './Regresion'
filename_climate = 'POWER_SinglePoint_Daily_19930101_20200712_34d59S_58d42W_cb0244f8.csv'
filename_gasm3 = 'Todos_010601a-Cuadro_I_5_1a.xlsx'
filename_gasN = 'Todos_010601b-Cuadro_I_5_1b.xlsx'
df = pd.read_csv(filepath + '/' + filename_climate, skiprows = 16)
df.replace(-999, value = 0)
df['T2M_MED'] = ((df['T2M_MAX'] + df['T2M_MIN'])/ 2).round(1)
df['GD'] = np.where(df['T2M_MED'] <= 20, (20 - df['T2M_MED']), 0).round(0)
df.columns = ['Lat', 'Lon', 'Year', 'Month', 'Day', 'TMax', 'TMin', 'RHum', 'Prec', 'Pres', 'Wind', 'Rad', 'TMed', 'GD']
df.head()
df_mes = pd.pivot_table(df, values = ['GD', 'Prec'], index = ['Year', 'Month'], aggfunc='sum')
df_mes1 = pd.pivot_table(df, values = ['TMed', 'RHum', 'Pres', 'Wind', 'Rad'], index = ['Year', 'Month'], aggfunc='mean')
df_mes1['Pres'] = df_mes1['Pres'] * 10
df_mes['TMed'] = df_mes1['TMed'].values.round(1)
df_mes['RHum'] = df_mes1['RHum'].values.round(1)
df_mes['Pres'] = df_mes1['Pres'].values.round(2)
df_mes['Wind'] = df_mes1['Wind'].values.round(2)
df_mes['Rad'] = df_mes1['Rad'].values.round(2)
df_mes = df_mes.iloc[:324]
df_mes
Cargar y normalizar datos de consumo de gas por usuario¶
dfm3 = pd.read_excel(filepath + '/' + filename_gasm3, header = 9, nrows = 324)
dfm3 = dfm3[['Mes', 'Capital Federal']]
dfusuarios = pd.read_excel(filepath + '/' + filename_gasN, header = 9, nrows = 324)
dfusuarios = dfusuarios[['Capital Federal']]
dfm3['usuarios'] = dfusuarios['Capital Federal'].values
dfm3['usuariomed'] = ((dfm3['Capital Federal'] * 1000) / dfm3['usuarios']).round(2)
dfm3
df_mes['Mes'] = dfm3['Mes'].values
df_mes['Gas'] = dfm3['usuariomed'].values.round(2)
df = df_mes[['Mes', 'TMed', 'GD', 'Wind', 'Prec', 'Rad', 'RHum', 'Pres', 'Gas']]
#df = df[df['Mes'] >= '2000-01-01']
df.reset_index(drop=True, inplace=True)
df
Graficar relación entre variables¶
Las gráficas bivariantes ayudan a comprender las relaciones que se establecen entre los datos. Las primeras gráficas abordarán la totalidad de las variables, en una matrix de gráficas de dispersión/histogramas, y una matriz de correlación. Se hará énfasis en mostrar las relaciones más estrechas entre las variables climáticas y el consumo de gas.
sns.pairplot(df)
plt.subplots(figsize=(6,6), dpi=100)
sns.set(font_scale=0.8)
cmap = sns.diverging_palette(220, 10, as_cmap=True)
sns.heatmap(df.corr(), annot = True, cmap = cmap)
Se observan en rojo las relaciones positivas y en azul las relaciones inversas. Las relaciones en blanco o más pálidas, son relaciones que no se puede afirmar que estén condicionadas. La relación que nos interesa es la de la última columna, que es la relación del consumo de gas del usuario medio con las variables. Allí podemos observar que:
- La temperatura media posee una correlación inversa del 0.92 con el consumo de gas.
- Los grados día de calefacción calculados poseen una correlación positiva del 0.96 con el consumo de gas.
- No existe una relación clara entre la velocidad del viento, las precipitaciones mensuales y la radiación media incidente en una superficie horizontal y el consumo de gas.
- Existe una relación positiva del 0.71 entre el consumo de gas y la humedad relativa, lo que tiene sentido dada su relación inversa con la temperatura (-0.8 de correlación).
- Existe una relación positiva del 0.78 entre el consumo de gas y la presión, lo que tiene sentido dada su relación inversa con la temperatura por el movimiento de las masas de aire (-0.86 de correlación).
plt.style.use('seaborn')
t = df['Mes']
data1 = df['TMed']
data2 = df['Gas']
fig, ax1 = plt.subplots(figsize=(12,4), dpi=100)
# x axis
color = 'tab:blue'
ax1.set_xlabel('Years')
# y1 axis
ax1.set_ylabel('Mean Temperature [ºC]', color=color)
ax1.plot(t, data1, color=color, lw=0.8, linestyle='--')
ax1.tick_params(axis='y', labelcolor=color)
ax1.set_ylim([25, 0])
# y2 axis
ax2 = ax1.twinx()
color = 'tab:red'
ax2.set_ylabel('Gas consumption[m3]', color=color)
ax2.plot(t, data2, color=color, lw=0.8)
ax2.tick_params(axis='y', labelcolor=color)
ax2.set_ylim([0, 250])
fig.tight_layout()
Para poder observar de manera más gráfica la relación inversa entre temperatura y consumo de gas, se ha invertido la escala de temperatura en el eje.
plt.style.use('seaborn')
t = df['Mes']
data1 = df['GD']
data2 = df['Gas']
fig, ax1 = plt.subplots(figsize=(12,4), dpi=100)
# x axis
color = 'tab:blue'
ax1.set_xlabel('Years')
# y1 axis
ax1.set_ylabel('Degree-Days [ºC]', color=color)
ax1.plot(t, data1, color=color, lw=0.8, linestyle='--')
ax1.tick_params(axis='y', labelcolor=color)
ax1.set_ylim([0, 500])
# y2 axis
ax2 = ax1.twinx()
color = 'tab:red'
ax2.set_ylabel('Gas consumption[m3]', color=color)
ax2.plot(t, data2, color=color, lw=0.8)
ax2.tick_params(axis='y', labelcolor=color)
ax2.set_ylim([0, 250])
fig.tight_layout()
sns.jointplot(x='TMed',y='Gas',data=df,kind='reg')
corr, _ = spearmanr(df['TMed'], df['Gas'])
print('Spearmans correlation: %.3f' % corr)
df1 = df[df['GD'] != 0]
sns.jointplot(x='GD',y='Gas',data=df1,kind='reg')
corr, _ = spearmanr(df1['GD'], df1['Gas'])
print('Spearmans correlation: %.3f' % corr)
Modelo de regresión lineal¶
La regresión lineal es un enfoque lineal para modelar la relación entre una respuesta escalar (o variable dependiente) y una o más variables explicativas (o variables independientes). En la regresión lineal, las relaciones se modelan utilizando funciones de predicción lineal cuyos parámetros de modelo desconocidos se estiman a partir de los datos. Tales modelos se llaman modelos lineales.
Entrenamiento del modelo¶
# X and Y arrays
X = df[['TMed', 'GD', 'Pres', 'RHum']]
y = df['Gas']
# train test split
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.3, random_state=101)
#creating and training the model
lm = LinearRegression()
lm.fit(X_train,y_train)
#model evaluation
print(lm.intercept_)
coeff_df = pd.DataFrame(lm.coef_, index = X.columns, columns = ['Coefficient'])
coeff_df
Interpretación de los coeficientes:
Manteniendo todas las variables fijas, un incremento de 1 en...
...Temperatura media (ºC) está asociado a un incremento de 2.29m3 de consumo de gas por usuario.
...Grados-día de calefacción (ºC) está asociado a un incremento de 0.52m3 de consumo de gas por usuario.
...Presión atmosférica media (hPa) está asociado a una reducción de 0.78m3 de consumo de gas por usuario.
...Humedad relativa media (%) está asociado a un incremento de 0.22m3 de consumo de gas por usuario.
Predicciones¶
predictions = lm.predict(X_test)
plt.scatter(y_test,predictions)
sns.distplot((y_test-predictions),bins=50);
MAE = metrics.mean_absolute_error(y_test, predictions)
MSE = metrics.mean_squared_error(y_test, predictions)
RMSE = np.sqrt(metrics.mean_squared_error(y_test, predictions))
print('MAE:', MAE)
print('MSE:', MSE)
print('RMSE:', RMSE)
print(min(df['Gas']))
print(max(df['Gas']))
print(max(df['Gas']) - min(df['Gas']))
print('MAE:',MAE/ (max(df['Gas']) - min(df['Gas'])))
print('RMSE:',RMSE/(max(df['Gas']) - min(df['Gas'])))
En un ejercicio simple para tratar de predecir el consumo de gas de acuerdo a las condiciones climáticas, el error predictivo del modelo es muy bajo por la estrecha correlación existente entre las variables explicativas (Temperatura, Grados-día de calefacción, Presión atmosférica, Humedad Relativa) y la respuesta escalar (Consumo de gas).