<jupyter_start><jupyter_text>Rain Prediction, Model Evaluation & GeneralizationThis dataset contains about 10 years of daily weather observations from many locations across Australia. RainTomorrow is the target variable to predict.Objectives:- Developing an ML model that suits best for the datasets predicting the rain on the next day.- How would you evaluate your model and why do you consider the model as a good fit?- What are the steps you will prefer to do for improving your model for further analysis?  Import libraries and dataset<jupyter_code>import numpy as np # linear algebra
import pandas as pd # data processing, CSV file I/O (e.g. pd.read_csv)
import matplotlib.pyplot as plt # data visualization
import seaborn as sns # statistical data visualization
%matplotlib inline

import warnings
warnings.filterwarnings('ignore')
data = '../input/weather-dataset-rattle-package/weatherAUS.csv'
df = pd.read_csv(data)
# view dimensions of dataset

df.shape
# preview the dataset

df.head()
col_names = df.columns

col_names
# view summary of dataset

df.info()<jupyter_output><class 'pandas.core.frame.DataFrame'>
RangeIndex: 145460 entries, 0 to 145459
Data columns (total 23 columns):
 #   Column         Non-Null Count   Dtype  
---  ------         --------------   -----  
 0   Date           145460 non-null  object 
 1   Location       145460 non-null  object 
 2   MinTemp        143975 non-null  float64
 3   MaxTemp        144199 non-null  float64
 4   Rainfall       142199 non-null  float64
 5   Evaporation    82670 non-null   float64
 6   Sunshine       75625 non-null   float64
 7   WindGustDir    135134 non-null  object 
 8   WindGustSpeed  135197 non-null  float64
 9   WindDir9am     134894 non-null  object 
 10  WindDir3pm     141232 non-null  object 
 11  WindSpeed9am   143693 non-null  float64
 12  WindSpeed3pm   142398 non-null  float64
 13  Humidity9am    142806 non-null  float64
 14  Humidity3pm    140953 non-null  float64
 15  Pressure9am    130395 non-null  float64
 16  Pressure3pm    130432 non-null  float64
 17  Cloud9am       89572 non-null[...]<jupyter_text>Developing an ML model that suits best for the datasets predicting the rain on the next day Feature Selection<jupyter_code># Declare feature vector and target variable

X = df.drop(['RainTomorrow'], axis=1)
y = df['RainTomorrow']
# split X and y into training and testing sets

from sklearn.model_selection import train_test_split

X_train, X_test, y_train, y_test = train_test_split(X, y, test_size = 0.2, random_state = 0)
# check the shape of X_train and X_test

X_train.shape, X_test.shape<jupyter_output><empty_output><jupyter_text>Feature Engineering<jupyter_code># check data types in X_train

X_train.dtypes
# display categorical variables

categorical = [col for col in X_train.columns if X_train[col].dtypes == 'O']
categorical
# display numerical variables

numerical = [col for col in X_train.columns if X_train[col].dtypes != 'O']
numerical<jupyter_output><empty_output><jupyter_text>Handling missing values in numerical variables<jupyter_code># check missing values in numerical variables in X_train

X_train[numerical].isnull().sum()
# check missing values in numerical variables in X_test

X_test[numerical].isnull().sum()
# print percentage of missing values in the numerical variables in training set

for col in numerical:
    if X_train[col].isnull().mean()>0:
        print(col, round(X_train[col].isnull().mean(),4))
# impute missing values in X_train and X_test with respective column median in X_train

for df1 in [X_train, X_test]:
    for col in numerical:
        col_median=X_train[col].median()
        df1[col].fillna(col_median, inplace=True)
# check again missing values in numerical variables in X_train

X_train[numerical].isnull().sum()
# check missing values in numerical variables in X_test

X_test[numerical].isnull().sum()<jupyter_output><empty_output><jupyter_text>Now, we can see that there are no missing values in the numerical columns of training and test set.  Handling missing values in categorical variables<jupyter_code># print percentage of missing values in the categorical variables in training set

X_train[categorical].isnull().mean()
# print categorical variables with missing data

for col in categorical:
    if X_train[col].isnull().mean()>0:
        print(col, (X_train[col].isnull().mean()))
# impute missing categorical variables with most frequent value

for df2 in [X_train, X_test]:
    df2['WindGustDir'].fillna(X_train['WindGustDir'].mode()[0], inplace=True)
    df2['WindDir9am'].fillna(X_train['WindDir9am'].mode()[0], inplace=True)
    df2['WindDir3pm'].fillna(X_train['WindDir3pm'].mode()[0], inplace=True)
    df2['RainToday'].fillna(X_train['RainToday'].mode()[0], inplace=True)
# check missing values in categorical variables in X_train

X_train[categorical].isnull().sum()
# check missing values in categorical variables in X_test

X_test[categorical].isnull().sum()<jupyter_output><empty_output><jupyter_text>As a final check, I will check for missing values in X_train and X_test.<jupyter_code># check missing values in X_train

X_train.isnull().sum()
# check missing values in X_test

X_test.isnull().sum()<jupyter_output><empty_output><jupyter_text>We can see that there are no missing values in X_train and X_test.  The outliers in numerical variables<jupyter_code>def max_value(df3, variable, top):
    return np.where(df3[variable]>top, top, df3[variable])

for df3 in [X_train, X_test]:
    df3['Rainfall'] = max_value(df3, 'Rainfall', 3.2)
    df3['Evaporation'] = max_value(df3, 'Evaporation', 21.8)
    df3['WindSpeed9am'] = max_value(df3, 'WindSpeed9am', 55)
    df3['WindSpeed3pm'] = max_value(df3, 'WindSpeed3pm', 57)
X_train.Rainfall.max(), X_test.Rainfall.max()
X_train.Evaporation.max(), X_test.Evaporation.max()
X_train.WindSpeed9am.max(), X_test.WindSpeed9am.max()
X_train.WindSpeed3pm.max(), X_test.WindSpeed3pm.max()
X_train[numerical].describe()<jupyter_output><empty_output><jupyter_text>Encode categorical variables<jupyter_code>categorical
X_train[categorical].head()
# encode RainToday variable

import category_encoders as ce
encoder = ce.BinaryEncoder(cols=['RainToday'])
X_train = encoder.fit_transform(X_train)
X_test = encoder.transform(X_test)
X_train.head()
X_train = pd.concat([X_train[numerical], X_train[['RainToday_0', 'RainToday_1']],
                     pd.get_dummies(X_train.Location), 
                     pd.get_dummies(X_train.WindGustDir),
                     pd.get_dummies(X_train.WindDir9am),
                     pd.get_dummies(X_train.WindDir3pm)], axis=1)
X_train.head()
X_test = pd.concat([X_test[numerical], X_test[['RainToday_0', 'RainToday_1']],
                     pd.get_dummies(X_test.Location), 
                     pd.get_dummies(X_test.WindGustDir),
                     pd.get_dummies(X_test.WindDir9am),
                     pd.get_dummies(X_test.WindDir3pm)], axis=1)
X_test.head()<jupyter_output><empty_output><jupyter_text>The Feature Scaling<jupyter_code>X_train.describe()
cols = X_train.columns
from sklearn.preprocessing import MinMaxScaler
scaler = MinMaxScaler()
X_train = scaler.fit_transform(X_train)
X_test = scaler.transform(X_test)
X_train = pd.DataFrame(X_train, columns=[cols])
X_test = pd.DataFrame(X_test, columns=[cols])
X_train.describe()<jupyter_output><empty_output><jupyter_text>Model training<jupyter_code># train a logistic regression model on the training set
from sklearn.linear_model import LogisticRegression


# instantiate the model
logreg = LogisticRegression(solver='liblinear', random_state=0)

y_train.fillna(y_train.mode()[0], inplace=True)
y_test.fillna(y_train.mode()[0], inplace=True)

# fit the model
logreg.fit(X_train, y_train)<jupyter_output><empty_output><jupyter_text>Predict results<jupyter_code>y_pred_test = logreg.predict(X_test)

y_pred_test
from sklearn.metrics import accuracy_score

print('Model accuracy score: {0:0.4f}'. format(accuracy_score(y_test, y_pred_test)))<jupyter_output>Model accuracy score: 0.8487<jupyter_text>Check for overfitting and underfitting by the train-set and test-set accuracy<jupyter_code>y_pred_train = logreg.predict(X_train)
y_pred_train
# print the scores on training and test set

print('Training set score: {:.4f}'.format(logreg.score(X_train, y_train)))
print('Test set score: {:.4f}'.format(logreg.score(X_test, y_test)))<jupyter_output>Training set score: 0.8487
Test set score: 0.8487<jupyter_text>These two values are same in score. So, there is no question of overfitting.  In Logistic Regression, we use default value of C = 1. It provides good performance with approximately 85% accuracy on both the training and the test set. But the model performance on both the training and test set are very comparable. It is likely the case of underfitting. I will increase C and fit a more flexible model.<jupyter_code># fit the Logsitic Regression model with C=100

# instantiate the model
logreg100 = LogisticRegression(C=100, solver='liblinear', random_state=0)

# fit the model
logreg100.fit(X_train, y_train)
# print the scores on training and test set

print('Training set score: {:.4f}'.format(logreg100.score(X_train, y_train)))
print('Test set score: {:.4f}'.format(logreg100.score(X_test, y_test)))<jupyter_output>Training set score: 0.8487
Test set score: 0.8485<jupyter_text>We can see that, C=100 results in higher test set accuracy and also a slightly increased training set accuracy. Lets see, what happens if we use more regularized model than the default value of C=1, by setting C=0.01.<jupyter_code># fit the Logsitic Regression model with C=001

# instantiate the model
logreg001 = LogisticRegression(C=0.01, solver='liblinear', random_state=0)

# fit the model
logreg001.fit(X_train, y_train)
# print the scores on training and test set

print('Training set score: {:.4f}'.format(logreg001.score(X_train, y_train)))
print('Test set score: {:.4f}'.format(logreg001.score(X_test, y_test)))<jupyter_output>Training set score: 0.8427
Test set score: 0.8415<jupyter_text>So, if we use more regularized model by setting C=0.01, then both the training and test set accuracy decrease relative to the default parameters.  Compare model accuracy with null accuracy<jupyter_code># check class distribution in test set

y_test.value_counts()<jupyter_output><empty_output><jupyter_text>We can see that the occurences of most frequent class is 22726. So, we can calculate null accuracy by dividing 22726 by total number of occurences.<jupyter_code># check null accuracy score

null_accuracy = (22726/(22726+6366))
print('Null accuracy score: {0:0.4f}'. format(null_accuracy))<jupyter_output>Null accuracy score: 0.7812<jupyter_text>We can see that our model accuracy score is 0.8485 but null accuracy score is 0.7812. So, we can conclude that our Logistic Regression model is doing a very good job in predicting the class labels.  How would you evaluate your model and why do you consider the model as a good fit? Model Evaluation Confusion matrix<jupyter_code># Print the Confusion Matrix and slice it into four pieces

from sklearn.metrics import confusion_matrix

cm = confusion_matrix(y_test, y_pred_test)
print('Confusion matrix\n\n', cm)
print('\nTrue Positives(TP) = ', cm[0,0])
print('\nTrue Negatives(TN) = ', cm[1,1])
print('\nFalse Positives(FP) = ', cm[0,1])
print('\nFalse Negatives(FN) = ', cm[1,0])
# visualize confusion matrix with seaborn heatmap

cm_matrix = pd.DataFrame(data=cm, columns=['Actual Positive:1', 'Actual Negative:0'], 
                                 index=['Predict Positive:1', 'Predict Negative:0'])
sns.heatmap(cm_matrix, annot=True, fmt='d', cmap='YlGnBu')<jupyter_output><empty_output><jupyter_text>Classification metrices<jupyter_code>from sklearn.metrics import classification_report

print(classification_report(y_test, y_pred_test))
# classification accuracy
TP = cm[0,0]
TN = cm[1,1]
FP = cm[0,1]
FN = cm[1,0]
# print classification accuracy

classification_accuracy = (TP + TN) / float(TP + TN + FP + FN)
print('Classification accuracy : {0:0.4f}'.format(classification_accuracy))
# print classification error

classification_error = (FP + FN) / float(TP + TN + FP + FN)
print('Classification error : {0:0.4f}'.format(classification_error))
# print precision score

precision = TP / float(TP + FP)
print('Precision : {0:0.4f}'.format(precision))
recall = TP / float(TP + FN)
print('Recall or Sensitivity : {0:0.4f}'.format(recall))
# True Positive Rate is synonymous with Recall.

true_positive_rate = TP / float(TP + FN)
print('True Positive Rate : {0:0.4f}'.format(true_positive_rate))
# False Positive Rate

false_positive_rate = FP / float(FP + TN)
print('False Positive Rate : {0:0.4f}'.format(false_positive_rate))
# Specificity

specificity = TN / (TN + FP)
print('Specificity : {0:0.4f}'.format(specificity))<jupyter_output>Specificity : 0.7277<jupyter_text>Adjusting the threshold level<jupyter_code># print the first 10 predicted probabilities of two classes- 0 and 1

y_pred_prob = logreg.predict_proba(X_test)[0:10]
y_pred_prob
# store the probabilities in dataframe

y_pred_prob_df = pd.DataFrame(data=y_pred_prob, columns=['Prob of - No rain tomorrow (0)', 'Prob of - Rain tomorrow (1)'])
y_pred_prob_df
# print the first 10 predicted probabilities for class 1 - Probability of rain

logreg.predict_proba(X_test)[0:10, 1]
# store the predicted probabilities for class 1 - Probability of rain

y_pred1 = logreg.predict_proba(X_test)[:, 1]
# plot histogram of predicted probabilities

# adjust the font size 
plt.rcParams['font.size'] = 12

# plot histogram with 10 bins
plt.hist(y_pred1, bins = 10)

# set the title of predicted probabilities
plt.title('Histogram of predicted probabilities of rain')

# set the x-axis limit
plt.xlim(0,1)

# set the title
plt.xlabel('Predicted probabilities of rain')
plt.ylabel('Frequency')<jupyter_output><empty_output><jupyter_text>Lower the threshold<jupyter_code>from sklearn.preprocessing import binarize

for i in range(1,5):
    cm1=0
    y_pred1 = logreg.predict_proba(X_test)[:,1]
    y_pred1 = y_pred1.reshape(-1,1)
    y_pred2 = binarize(y_pred1, i/10)
    y_pred2 = np.where(y_pred2 == 1, 'Yes', 'No')
    cm1 = confusion_matrix(y_test, y_pred2)    
    print ('With',i/10,'threshold the Confusion Matrix is ','\n\n',cm1,'\n\n',
            'with',cm1[0,0]+cm1[1,1],'correct predictions, ', '\n\n', 
            cm1[0,1],'Type I errors( False Positives), ','\n\n',
            cm1[1,0],'Type II errors( False Negatives), ','\n\n',
           'Accuracy score: ', (accuracy_score(y_test, y_pred2)), '\n\n',
           'Sensitivity: ',cm1[1,1]/(float(cm1[1,1]+cm1[1,0])), '\n\n',
           'Specificity: ',cm1[0,0]/(float(cm1[0,0]+cm1[0,1])),'\n\n',
            '====================================================', '\n\n')<jupyter_output>With 0.1 threshold the Confusion Matrix is  

 [[13243  9483]
 [  578  5788]] 

 with 19031 correct predictions,  

 9483 Type I errors( False Positives),  

 578 Type II errors( False Negatives),  

 Accuracy score:  0.6541660937714836 

 Sensitivity:  0.9092051523719761 

 Specificity:  0.5827246325794244 

 ==================================================== 


With 0.2 threshold the Confusion Matrix is  

 [[17759  4967]
 [ 1355  5011]] 

 with 22770 correct predictions,  

 4967 Type I errors( False Positives),  

 1355 Type II errors( False Negatives),  

 Accuracy score:  0.782689399147532 

 Sensitivity:  0.7871504869619855 

 Specificity:  0.781439760626595 

 ==================================================== 


With 0.3 threshold the Confusion Matrix is  

 [[19732  2994]
 [ 2043  4323]] 

 with 24055 correct predictions,  

 2994 Type I errors( False Positives),  

 2043 Type II errors( False Negatives),  

 Accuracy score:  0.8268596177643338 

 Sensitivity:  0.67907634[...]<jupyter_text>ROC - AUC<jupyter_code># plot ROC Curve

from sklearn.metrics import roc_curve

fpr, tpr, thresholds = roc_curve(y_test, y_pred1, pos_label = 'Yes')
plt.figure(figsize=(6,4))
plt.plot(fpr, tpr, linewidth=2)
plt.plot([0,1], [0,1], 'k--' )
plt.rcParams['font.size'] = 12
plt.title('ROC curve for RainTomorrow classifier')
plt.xlabel('False Positive Rate (1 - Specificity)')
plt.ylabel('True Positive Rate (Sensitivity)')

plt.show()<jupyter_output><empty_output><jupyter_text>ROC-AUC<jupyter_code># compute ROC AUC

from sklearn.metrics import roc_auc_score

ROC_AUC = roc_auc_score(y_test, y_pred1)
print('ROC AUC : {:.4f}'.format(ROC_AUC))
# calculate cross-validated ROC AUC 

from sklearn.model_selection import cross_val_score

Cross_validated_ROC_AUC = cross_val_score(logreg, X_train, y_train, cv=5, scoring='roc_auc').mean()
print('Cross validated ROC AUC : {:.4f}'.format(Cross_validated_ROC_AUC))<jupyter_output>Cross validated ROC AUC : 0.8674<jupyter_text>What are the steps you will prefer to do for improving your model for further analysis? Model Generalization k-Fold Cross Validation<jupyter_code># Applying 5-Fold Cross Validation

from sklearn.model_selection import cross_val_score

scores = cross_val_score(logreg, X_train, y_train, cv = 5, scoring='accuracy')
print('Cross-validation scores:{}'.format(scores))
# compute Average cross-validation score

print('Average cross-validation score: {:.4f}'.format(scores.mean()))<jupyter_output>Average cross-validation score: 0.8483<jupyter_text>The Average cross-validation score is 0.8483 and our model original score was 0.8485. So, it's not an improvement.  Hyperparameter Optimization using GridSearch CV<jupyter_code>from sklearn.model_selection import GridSearchCV

parameters = [{'penalty':['l1','l2']}, 
              {'C':[1, 10, 100, 1000]}]

grid_search = GridSearchCV(estimator = logreg,  
                           param_grid = parameters,
                           scoring = 'accuracy',
                           cv = 5,
                           verbose=0)

grid_search.fit(X_train, y_train)
# examine the best model

# best score achieved during the GridSearchCV
print('GridSearch CV best score : {:.4f}\n\n'.format(grid_search.best_score_))
# print parameters that give the best results
print('Parameters that give the best results :','\n\n', (grid_search.best_params_))
# print estimator that was chosen by the GridSearch
print('\n\nEstimator that was chosen by the search :','\n\n', (grid_search.best_estimator_))
# calculate GridSearch CV score on test set

print('GridSearch CV score on test set: {0:0.4f}'.format(grid_search.score(X_test, y_test)))<jupyter_output>GridSearch CV score on test set: 0.8485
{'fraction_non_alphanumeric': 0.1163205710577547, 'fraction_numerical': 0.0427211767250703, 'mean_word_length': 4.831914222642699, 'pattern_counts': {'":': 0, '<': 87, '<?xml version=': 0, '>': 90, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '24', 'dir': 'jupyter-structured-clean-dedup', 'id': '137208', 'license': 'MIT', 'n_tokens_mistral': 6755, 'n_tokens_neox': 6071, 'n_words': 1700, 'path': 'rain-prediction-model-evaluation-generalization.ipynb', 'repo_name': 'saadbinmanjur/Rain-Prediction-Model-Evaluation-and-Generalization', 'size': '18492'}