Ethekwini Electricity Usage Analysis¶
Objectives¶
- To visualise the trend of Ethekwini Munipality Electricity usage from 1996-2019.
- To Develop a model to predict electricity demand based on usage and access data.
In [288]:
#Import important libriries
import pandas as pd
pd.set_option('display.max_columns', 500)
pd.set_option('display.max_rows', 500)
from pandas import read_csv
In [289]:
cd
In [290]:
df1 = pd.read_excel(r"C:\Users\Student 95\Downloads\electricity-use-ethekwini.xlsx")
df1.head()
Out[290]:
In [409]:
# droping variable that are uneeded to the contribution of the gaol
In [291]:
df1 = df1.drop(columns=['Percent energy growth (%)','Percent growth (%)','Percent loss (%)',
'Power factor at system peak (%)','Average monthly load factor (%)',
'Percent growth','Maximum kVA','Energy (kWh) purchased','Number of customers'], axis=1)
In [292]:
df1.shape
Out[292]:
In [293]:
df1.head()
Out[293]:
In [294]:
df2 = pd.read_excel(r"C:\Users\Student 95\Downloads\electricity-customers-and-revenue.xlsx")
df2.head()
Out[294]:
In [295]:
df2.shape
Out[295]:
In [410]:
# Renaming the variable
In [296]:
df2 = df2.rename(columns={'Unnamed: 0': 'Year'})
In [411]:
# merging Datasets
In [297]:
df = pd.merge(df1,df2, on = "Year", how = "inner")
df.head()
Out[297]:
In [301]:
power.columns
Out[301]:
In [303]:
power = power[['A: Total number of customers', 'B: Total use','C: Total revenue: total','D: Cents per total',
'E: Use per total','F: Rands per total','Year', 'Energy (kWh) sold']]
power.head()
Out[303]:
In [304]:
columns_rearranged = ['Year', 'A: Total number of customers','C: Total revenue: total','D: Cents per total',
'E: Use per total','F: Rands per total', 'Energy (kWh) sold']
In [412]:
# arraging the datasets
In [305]:
power = power.reindex(columns=columns_rearranged)
power.head()
Out[305]:
In [306]:
#We understand the nature of our dataset and the type of attributes
power.info()
In [307]:
power.head()
Out[307]:
In [308]:
#Data visualization libraries
import seaborn as sns
import matplotlib.pyplot as plt
%matplotlib inline
In [309]:
#Check the Null values for the entire DataFrame
power.isnull().sum().sum()
Out[309]:
In [310]:
power.head()
Out[310]:
In [311]:
power.columns[power.isna().any()].tolist()
Out[311]:
In [312]:
power.head()
Out[312]:
In [313]:
#Drop any column with NaN values
power = power.dropna(axis=1, how='all')
EXPLORATORY DATA ANALYSIS¶
In [392]:
#We plot the histogram to check the distribution of the datasets in each attributes
power.hist(column=columns_rearranged, rwidth=0.95, figsize=(15,8))
Out[392]:
In [389]:
#We do the count plot of the Energy (kWh) sold in different years range to check if it correspond with the dataset.
from matplotlib.pyplot import figure
figure(figsize=(16,8))
sns.countplot(data=power, x='Year', hue='Energy (kWh) sold')
Out[389]:
In [400]:
# Plotting the Histogram barplot to compare the frequency distribution of D: Cents per total variable.
# We also create a range to compare and analyse different distribution.
bins=[40,50,60,70,80,90,100,110,120,130,140,150,160,170,180,190,200]
plt.hist(power['D: Cents per total'], rwidth=0.75, bins=bins, color='g')
plt.xlabel('Cents range')
plt.ylabel('Total no of cents')
plt.title('Cents Per Total')
Out[400]:
Feature-Feature Relationships¶
- We explore the relationship between the attributes.
In [314]:
# The scatter charts between the pair of variables to understand the relationship.
sns.pairplot(power)
Out[314]:
In [316]:
from pandas.plotting import scatter_matrix
scatter_matrix(power, alpha=0.2, figsize=(16, 10), diagonal='kde')
Out[316]:
Correlation Analysis and Feature Selection¶
In [ ]:
power=power.astype('int64', )
In [324]:
power.columns
Out[324]:
In [326]:
power['C: Total revenue: total'] = power ['C: Total revenue: total'].astype('int64')
In [328]:
power['D: Cents per total'] = power ['D: Cents per total'].astype('int64')
In [329]:
power['E: Use per total'] = power ['E: Use per total'].astype('int64')
In [330]:
power['F: Rands per total'] = power ['F: Rands per total'].astype('int64')
In [331]:
power['Energy (kWh) sold'] = power ['Energy (kWh) sold'].astype('int64')
In [332]:
power.info()
In [333]:
#Important libraries
from sklearn.feature_selection import SelectKBest
from sklearn.feature_selection import chi2
X = power.iloc[:,1:6]
y = power.iloc[:,-1]
In [334]:
#Apply SeleckKBest class to extarct top 6 best features
bestfeatures = SelectKBest(score_func=chi2, k=5)
fit = bestfeatures.fit(X,y)
In [335]:
dfscores = pd.DataFrame(fit.scores_)
dfcolumns = pd.DataFrame(X.columns)
In [336]:
#Concat two DataFrames for better visualization
featureScores = pd.concat([dfcolumns, dfscores], axis=1)
featureScores.columns = ['col','Score']
In [337]:
featureScores
Out[337]:
In [338]:
# Print the 6 best features
print(featureScores.nlargest(4,'Score'))
In [340]:
sns.countplot(data=power, x='Energy (kWh) sold', palette='hls')
plt.show()
plt.savefig('count_plot')
Correlation Matrix with Heatmap¶
In [339]:
#The correlation of each feauture in the dataset
corrmat = power.corr()
top_corr_features = corrmat.index
plt.figure(figsize=(16,8))
#We plot the heatmap
sns.heatmap(power[top_corr_features].corr(), annot=True, cmap='RdYlGn')
Out[339]:
BUILDING AND DEPLOYING MODELS¶
Split the dataset¶
To understand the model perfomance we split the datasets into training set and test set.
In [342]:
# Split X and y into training and testing set.
from sklearn.model_selection import train_test_split
In [343]:
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.25, random_state=0)
The 75% of the data will be used for model training and 25% for model testing.
In [355]:
# import the class
from sklearn.linear_model import LogisticRegression
In [356]:
#instantiate the model
logreg = LogisticRegression()
#fit the model with train data
logreg.fit(X_train, y_train)
#fit the model with predictor test data(X_test)
y_pred = logreg.predict(X_test)
In [357]:
y_pred
Out[357]:
In [358]:
# import the metrics
from sklearn import metrics
cf_matrix = metrics.confusion_matrix(y_test, y_pred)
cf_matrix
Out[358]:
In [348]:
# Import the libraries (pipeline and models)
from sklearn.pipeline import Pipeline
from sklearn.ensemble import RandomForestClassifier
from sklearn.neighbors import KNeighborsClassifier
from sklearn.tree import DecisionTreeClassifier
from sklearn.linear_model import LogisticRegression
In [344]:
pipeline_dt=Pipeline([('dt_classifier',DecisionTreeClassifier(random_state=0))])
In [345]:
pipeline_rf=Pipeline([('rf_classifier',RandomForestClassifier())])
In [346]:
pipeline_knn=Pipeline([('kn_classifier',KNeighborsClassifier())])
In [349]:
pipeline_lr=Pipeline([('lr_classifier',LogisticRegression())])
In [350]:
#Make the list of pipelines
pipelines = [pipeline_dt,pipeline_rf,pipeline_knn,pipeline_lr]
In [351]:
best_accuracy=0.0
best_classifier=0
best_pipeline=""
In [352]:
#Dictionery of pipelines and classifier type for ease of reference
pipe_dict = {0: 'Decision Tree', 1: 'RandomForest', 2: 'KNeighbors', 3:'Logistic Regression'}
#Fit the pipilines
for pipe in pipelines:
pipe.fit(X_train, y_train)
In [353]:
for i,model in enumerate(pipelines):
print("{} Test Accuracy: {}".format(pipe_dict[i],model.score(X_test,y_test)))
In [354]:
for i,model in enumerate(pipelines):
if model.score(X_test,y_test)>best_accuracy:
best_accuracy=model.score(X_test,y_test)
best_pipeline=model
best_classifier=i
print('Classifier with best accuracy:{}'.format(pipe_dict[best_classifier]))
The Predicted Results for each Agorithm¶
Decision Tree Predicted Results¶
In [359]:
y_pred_0 = pipeline_dt.predict(X_test)
In [360]:
y_pred_0
Out[360]:
Random Forest Predicted Results¶
In [364]:
y_pred_1 = pipeline_rf.predict(X_test)
In [365]:
y_pred_1
Out[365]:
KNeighborsClassifier Predicted Results¶
In [362]:
y_pred_3 = pipeline_knn.predict(X_test)
In [363]:
y_pred_3
Out[363]:
LogisticRegression Predicted Results¶
In [366]:
y_pred_4=pipeline_lr.predict(X_test)
In [367]:
y_pred_4
Out[367]:
Testing the Model Performance of the Decision Tree¶
In [401]:
# Use the forest's predict method on the test data
predictions = pipeline_dt.predict(X_test)
In [402]:
# Calculate the absolute errors
errors = abs(predictions - y_test)
In [404]:
import numpy as np
In [405]:
# Print out the mean absolute error (mae)
print('Mean Absolute Error:', round(np.mean(errors), 2), 'degrees.')
In [406]:
# Calculate mean absolute percentage error (MAPE)
mape = 100 * (errors / y_test)
In [407]:
# Calculate and display accuracy
accuracy = 100 - np.mean(mape)
In [408]:
print('Accuracy:', round(accuracy, 2), '%.')
Our model has learned to predict using the train and test dataset with an accuracy of 97.9%¶
Ensemble Learning (The performance of boosting models)¶
Ensemble techniques combine individual models together to improve the stabilty and predictive power of the model. One of the ensemble methods are Boosting and Gradient Boosting model.
Boosting technique reduce bias by adding all the weak learners sequantially, each classifier trying to correct its predecessor.
Gradient Boosting minimizes the loss function(MSE)of a models by adding weak learners using a gradient descesnt procedure.
Building Models¶
- Gradient BoostingRegressor
- XGBoost
In [371]:
#From sklearn we load Gradient Boosting
from sklearn.ensemble import GradientBoostingRegressor
In [372]:
%%time
GB_regressor = GradientBoostingRegressor()
GB_regressor.fit(X_train,y_train)
Out[372]:
In [373]:
y_predictions = GB_regressor.predict(X_train)
In [374]:
y_predictions
Out[374]:
In [376]:
from sklearn.metrics import r2_score, mean_squared_error
In [377]:
print('Score of Model:', r2_score(y_train,y_predictions))
print('Mean square error:', mean_squared_error(y_train, y_predictions))
Test predictions
In [378]:
y_predictions = GB_regressor.predict(X_test)
In [379]:
y_predictions
Out[379]:
In [380]:
print('Score of Model:', r2_score(y_test,y_predictions))
print('Mean square error:', mean_squared_error(y_test,y_predictions))
This is pretty! good our model is doing well with the train dataset prediction of 97% accuracy.¶
XGBoost Modelling¶
In [381]:
from xgboost import XGBRegressor
In [382]:
%%time
xgb_regressor = XGBRegressor()
xgb_regressor.fit(X_train, y_train)
Out[382]:
In [383]:
y_prediction2 = xgb_regressor.predict(X_train)
In [384]:
y_prediction2
Out[384]:
In [385]:
print('Score of Model :', r2_score(y_train,y_prediction2))
print('Mean square error :', mean_squared_error(y_train,y_prediction2))
Test predictions
In [386]:
y_prediction3 = xgb_regressor.predict(X_test)
In [387]:
y_prediction3
Out[387]:
In [388]:
print('Score of Model :', r2_score(y_test,y_prediction3))
print('Mean square error :', mean_squared_error(y_test,y_prediction3))
Observations¶
- As we apply above the ensemble models we come to observe the performance of these model in productionization Gradient Boost is the best model with 99% accuracy compared to Xtreme Gradient Boost with 85% accuracy.
- We also note that Gradient Boost is the slowest model with 979 ms compared to Xtreme GB which ran 945 ms.
- We can improve our ensembe learning models by doing more feature engineering:
- Removing outliers
- Dropping columns with low variance of about 0.01 threshold value
- Checking some correlation of the indexes and remove those columns with same correlation as it denotes duplicates
- Implementing the models using hyperparameter tuning.
In [ ]: