# Library imports here in first cell
# Importing the matplotlib.pyplot module and aliasing it as 'plt'
# This is a library for creating visualizations like plots, histograms, and charts.
import matplotlib.pyplot as plt

# Importing the numpy library and aliasing it as 'np'
# Numpy is used for numerical computations, such as handling arrays, mathematical operations, etc.
import numpy as np

# Importing the pandas library and aliasing it as 'pd'
# pandas is a library for data manipulation and analysis.
# It is used to read, process, clean, and analyze structured data like tables (e.g., CSV, Excel, or SQL files).
import pandas as pd

#For stansardisation and scaling of numerical values, utilise StandardScaler from sklearn.preprocessing library
from sklearn.preprocessing import StandardScaler

# Importing the LogisticRegression class from scikit-learn's linear_model module
# LogisticRegression is a machine learning model used for classification tasks and forms part 1 of this assignment
from sklearn.linear_model import LogisticRegression

# Importing evaluation metrics from scikit-learn
# classification_report: Provides metrics like precision, recall, F1-score for each class.
# confusion_matrix: Outputs a matrix showing the performance of a classification model.
from sklearn.metrics import classification_report, confusion_matrix, ConfusionMatrixDisplay

# Importing the train_test_split function from scikit-learn's model_selection module
# train_test_split is used to split a dataset into training and testing subsets.
# It is commonly used in machine learning workflows to evaluate model performance.
# Typical use: dividing the data into X_train, X_test, y_train, and y_test (features and target).
from sklearn.model_selection import train_test_split

# Imports evaluation metrics from the scikit-learn library.
#   accuracy_score: Computes the proportion of correct predictions (overall model accuracy).
#   precision_score: Measures the percentage of correctly identified positive predictions.
#   recall_score: Measures the ability to capture all actual positive samples (True Positives).
#   f1_score: Combines precision and recall into a single score (harmonic mean).
#   roc_auc_score: Computes the Area Under the Receiver Operating Characteristic (ROC) Curve.
#   roc_curve: Calculates points for the ROC curve (False Positive Rate, True Positive Rate) 
#              based on prediction probabilities.
from sklearn.metrics import accuracy_score, precision_score, recall_score, f1_score, roc_auc_score, roc_curve

# Imports utilities for splitting datasets and searching for model hyperparameters.
#   GridSearchCV: A tool for hyperparameter tuning which automates the search for the best 
#                 combination of hyperparameters for a model.
from sklearn.model_selection import GridSearchCV

# Import the libraries necessary to implement the SVM ML model
from sklearn.svm import SVC

# Read in the Telco dataset
telcoData = pd.read_csv('WA_Fn-UseC_-Telco-Customer-Churn.csv')

# Display the first few rows to verify it loaded correctly
print(telcoData.head())

# Display concise summary of the DataFrame, including non-null counts and data types
print(telcoData.info())

# Start of data validation in the dataset
# Check for missing values in each column
print(telcoData.isnull().sum())

# Count unique values in each column
unique_counts = telcoData.nunique()
print(unique_counts)

#It can be seen that there are 7 binary columns which will be mapped later. 
#However, out of 7k entries there are 73 values for tenure, below check min and max values and frequency
# Minimum value of tenure
min_value = telcoData['tenure'].min()
print(f"Minimum value: {min_value}")

# Maximum value of tenure
max_value = telcoData['tenure'].max()
print(f"Maximum value: {max_value}")

# Frequency of each unique value of tenure
value_counts = telcoData['tenure'].value_counts()
print(value_counts)

#We need to convert Churn from strings to boolean values of 0 and 1 - this will not be required to be converted in the binary list later
telcoData['Churn'] = telcoData['Churn'].map({'Yes': 1, 'No': 0})

#check churn only contains 0 and 1. 
print(telcoData['Churn'].unique())

#Plot a graph of tenure against churn to determine whether there is a link between the two.
#this affects whether to use numerical bin indices or one-hot encoding on the column

churn_by_tenure = telcoData.groupby('tenure')['Churn'].mean()
print(churn_by_tenure)

churn_by_tenure.plot(kind='bar', color='skyblue', figsize=(8, 5))
plt.title('Churn Rate by Tenure Group')
plt.xlabel('Tenure Group')
plt.ylabel('Churn Rate')
plt.xticks(rotation=45)
plt.show()

# The above graph shows that there is a strong correlation with tenure and churn - the longer tenure the less likely to churn
# therefore numerical bin indices will be implemented
# Below line creates a tenure group that only has numerical values called numerical bin indices as required.
telcoData['tenure_group'] = pd.cut(
    telcoData['tenure'],
    bins=[0, 12, 24, 36, 48, 60, telcoData['tenure'].max() + 1],
    labels=False,  # Numerical bin indices
    right=False
)

# Senior Citizen is already an int64 column which seems to contain 0 or 1, check for unique values
print(telcoData['SeniorCitizen'].unique())
is_binary = telcoData['SeniorCitizen'].isin([0,1]).all()
print(f"Column contains only 0s and 1s: {is_binary}")

# For all of the binary columns remaining from the original set found in the nunique check, there were 7 - 
# Five remain (Churn has already been converted, and SeniorCitizen is int64 with values of 0&1 already),
# binary columns should be made numerical 0 or 1 and checked for unique values

# define binary columns from the above results to ultimately make them numerical of 1 or 0.
binary_columns = ['gender', 'Partner', 'Dependents', 'PhoneService', 'PaperlessBilling']

# Loop through the binary columns set and print the unique values first to make sure the unique values are consistent
for col in binary_columns:        
    print(f"{col} unique values: {telcoData[col].unique()}")    # print the unique values to make sure that they are consistent.

# The values seem consistent and so we can run through the columns again and convert to 0&1, 
# note one column uses 'Male' & 'Female'.

# run through binary columns and map Male to 0 and female to 1 for gender and 
# Yes to 1 and No to 0 for all other columns, now that we know those are the discreet values.
for col in binary_columns:
    if col == 'gender':  # Special mapping for gender
        telcoData[col] = telcoData[col].map({'Male': 0, 'Female': 1})
    else:  # Standard mapping for Yes/No columns
        telcoData[col] = telcoData[col].map({'Yes': 1, 'No': 0})

#check that this has mapped correctly with head check
print(telcoData.head())

#Check telcoData columns
print(telcoData.columns)

#From the original nunique check there are 9 columns with 3 values
#now group the columns with 3 or 4 unique values
#first assign the three value columns to an array
three_value_columns = ['MultipleLines', 'InternetService', 'OnlineSecurity', 
                       'OnlineBackup', 'DeviceProtection', 'TechSupport', 
                       'StreamingTV', 'StreamingMovies', 'Contract']
#Then check for the unique values in  each for commonality so that group values can be assigned.
for col in three_value_columns:
    print(f"{col} unique values: {telcoData[col].unique()}")

#Use the pd.get_dummies function to convert the three value columns to two boolean columns.
telcoData = pd.get_dummies(telcoData, columns=three_value_columns, drop_first=True)

#There is one column (Payment method) that has four unique values, the same code can do this for that column 
telcoData = pd.get_dummies(telcoData, columns=['PaymentMethod'], drop_first=True)

#check that those columns are now converted
# Check the updated column names
print("Columns list:")
print(telcoData.columns)

# Verify a sample of the encoded columns
print("Sample of converted columns:")
print(telcoData[['MultipleLines_No phone service', 'MultipleLines_Yes']].head())

# Verify current column datatypes
print("Current column datatypes:")
print(telcoData.dtypes)

# Convert boolean columns to integers (0 and 1) as int64 type
bool_columns = telcoData.select_dtypes(include='bool').columns
telcoData[bool_columns] = telcoData[bool_columns].astype(int)

# Verify all data types are now numeric
print("Current column datatypes:")
print(telcoData.dtypes)

# There are only two columns left that require checking
# A) that they have complete numerical values and B) whether any other preparation is required.

# Check for non-numeric values in 'TotalCharges'
non_numeric_total = telcoData['TotalCharges'][pd.to_numeric(telcoData['TotalCharges'], errors='coerce').isna()]
print("TotalCharges:")
print(non_numeric_total)
print("TC End\n")
# Check for non-numeric values in 'MonthlyCharges'
non_numeric_monthly = telcoData['MonthlyCharges'][pd.to_numeric(telcoData['MonthlyCharges'], errors='coerce').isna()]
print("MonthlyCharges:")
print(non_numeric_monthly)

#There are rows in TotalCharges with non-numeric values, as detailed above.
#View their contents here to understand how to rectify
print(telcoData.loc[non_numeric_total.index, :])

#TotalCharges has "" or " " as its contents on those rows. This is fixed below:
telcoData['TotalCharges'] = telcoData['TotalCharges'].replace(' ', pd.NA)

#Then convert the column to numeric
telcoData['TotalCharges'] = pd.to_numeric(telcoData['TotalCharges'], errors='coerce')

#Now that TotalCharges and MonthlyCharges have been converted check datatypes across dataset
print (telcoData.dtypes)

#There are two ways to impute the NaN values, use the Mean or the Median, 
#If TotalCharges has a skewed distribution it is best to use median, but if the distribution is symmetrical mean is better.

#The following plots a histogram showing the distribution
plt.hist(telcoData['TotalCharges'], bins=30, edgecolor='k')
plt.title('Distribution of TotalCharges')
plt.xlabel('TotalCharges')
plt.ylabel('Frequency')
plt.show()

#The above shows a long tail and a skew to the left, which indiactes that the median would be more effective.
#The following line imputes the median.

telcoData['TotalCharges'] = telcoData['TotalCharges'].fillna(telcoData['TotalCharges'].median())

#Check the above has worked

print("Total of na entries in column TotalCharges:")
print(telcoData['TotalCharges'].isna().sum())  # Should be 0
print("\nDescription of Total Charges Column\n")
print(telcoData['TotalCharges'].describe())

#There are two columns that have a range of values which neeed to be scaled for encoding

# Initialize the scaler variable as an instance of StandardScaler
scaler = StandardScaler()

# Scale the numerical columns
telcoData[['MonthlyCharges', 'TotalCharges']] = scaler.fit_transform(telcoData[['MonthlyCharges', 'TotalCharges']])

#Check that the columns have been properly scaled.
print(telcoData[['MonthlyCharges', 'TotalCharges']].head())

#Check there are no NaN values left in the dataset
print(f"Min tenure_group: {telcoData['tenure_group'].min()}, Max tenure: {telcoData['tenure_group'].max()}")
print(telcoData['tenure_group'].isna().sum())

# Define the features (X) by dropping the target column customerID as it is unique for each chromosome
# tenure is dropped as it has been replaced by the tenure_group column
# and the Churn column as this is the target
X = telcoData.drop(columns=['customerID', 'tenure', 'Churn'])  
# Define the target (y) by selecting the target column Churn as that is the column to be predicted.
y = telcoData['Churn']

# Split into training and testing sets (80% training, 20% testing)
X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42)

# Check the shape of the splits and that there are no NaN entries
print(X_train.shape, X_test.shape, y_train.shape, y_test.shape)
print(X_train.isna().sum())

#Train the Logistic Regression Model
model = LogisticRegression()
model.fit(X_train, y_train)

# Evaluate the model
print(model.score(X_test, y_test))

#Display metrics for straight LR predictions
# Predictions on test set
y_pred = model.predict(X_test)
y_pred_proba = model.predict_proba(X_test)[:, 1]  # For AUC-ROC, use probabilities

# Calculate metrics
accuracy = accuracy_score(y_test, y_pred)
precision = precision_score(y_test, y_pred, average='weighted')
recall = recall_score(y_test, y_pred, average='weighted')
f1 = f1_score(y_test, y_pred, average='weighted')
auc = roc_auc_score(y_test, y_pred_proba)

# Print results
print(f"Accuracy: {accuracy:.2f}")
print(f"Precision: {precision:.2f}")
print(f"Recall: {recall:.2f}")
print(f"F1 Score: {f1:.2f}")
print(f"AUC-ROC: {auc:.2f}")

fpr, tpr, thresholds = roc_curve(y_test, y_pred_proba)
plt.plot(fpr, tpr, label="ROC Curve (AUC = {:.2f})".format(auc))
plt.xlabel("False Positive Rate")
plt.ylabel("True Positive Rate")
plt.title("ROC Curve")
plt.legend(loc="best")
plt.show()

#To display a confusion matrix of the results
# Generate confusion matrix
cm = confusion_matrix(y_test, y_pred)

# Display the confusion matrix
disp = ConfusionMatrixDisplay(confusion_matrix=cm, display_labels=model.classes_)
disp.plot(cmap="Blues")  # You can change the colormap if you like
plt.title("Confusion Matrix")
plt.show()

# Define the parameter grid
param_grid = {
    'C': [0.001, 0.01, 0.1, 1, 10],  # Regularization strength
    'penalty': ['l1', 'l2'],       # Regularization type
    #'class_weight': [None, 'balanced']
    'solver': ['liblinear', 'saga']  # Solvers compatible with l1 and l2
}

# Create a GridSearchCV object
grid_search = GridSearchCV(
    estimator=LogisticRegression(max_iter=5000),  # Base model
    param_grid=param_grid,                       # Parameter grid
    scoring='f1_weighted',                       # Optimize for weighted F1 score
    cv=5,                                       # 5-fold cross-validation
    verbose=1                                   # Print progress
)

# Perform the grid search on the training data
grid_search.fit(X_train, y_train)

# Best parameters and score
print("Best Hyperparameters:", grid_search.best_params_)
print("Best Cross-Validation F1 Score:", grid_search.best_score_)

# Train the final model with the best parameters
best_model = grid_search.best_estimator_

# Evaluate the final model on the test data
y_pred = best_model.predict(X_test)
y_pred_proba = best_model.predict_proba(X_test)[:, 1]

# Gather metrics and display
accuracy = accuracy_score(y_test, y_pred)
precision = precision_score(y_test, y_pred, average='weighted')
recall = recall_score(y_test, y_pred, average='weighted')
f1 = f1_score(y_test, y_pred, average='weighted')
auc = roc_auc_score(y_test, y_pred_proba)

print(f"Final Model Metrics:")
print(f"Accuracy: {accuracy:.2f}")
print(f"Precision: {precision:.2f}")
print(f"Recall: {recall:.2f}")
print(f"F1 Score: {f1:.2f}")
print(f"AUC-ROC: {auc:.2f}")

#Print new AUC-ROC curve
fpr, tpr, thresholds = roc_curve(y_test, y_pred_proba)
plt.plot(fpr, tpr, label="ROC Curve (AUC = {:.2f})".format(auc))
plt.xlabel("False Positive Rate")
plt.ylabel("True Positive Rate")
plt.title("ROC Curve")
plt.legend(loc="best")
plt.show()

#To display a confusion matrix of the results
# Generate confusion matrix
cm = confusion_matrix(y_test, y_pred)

# Display the confusion matrix
disp = ConfusionMatrixDisplay(confusion_matrix=cm, display_labels=model.classes_)
disp.plot(cmap="Blues")  # You can change the colormap if you like
plt.title("Confusion Matrix")
plt.show()

# Custom decision threshold
threshold = 0.52  # Adjust based on your needs
y_pred_custom = (y_pred_proba >= threshold).astype(int)

# Evaluate confusion matrix
cm_custom = confusion_matrix(y_test, y_pred_custom)
print("Confusion Matrix with Threshold 0.6:\n", cm_custom)

# Display the confusion matrix
disp = ConfusionMatrixDisplay(confusion_matrix=cm_custom, display_labels=model.classes_)
disp.plot(cmap="Blues")  # You can change the colormap if you like
plt.title("Confusion Matrix")
plt.show()

# Compute FPR, TPR, and thresholds
fpr, tpr, thresholds = roc_curve(y_test, y_pred_proba)

# Calculate metrics for each threshold
for i, threshold in enumerate(thresholds):
    print(f"Threshold: {threshold:.2f} | TPR: {tpr[i]:.2f} | FPR: {fpr[i]:.2f}")

# Plot the ROC curve
plt.plot(fpr, tpr, label="ROC Curve (AUC = {:.2f})".format(auc))
plt.axhline(y=0.5, color='r', linestyle='--', label="Random Guess")
plt.axvline(x=0.2, color='g', linestyle='--', label="Example Threshold")
plt.xlabel("False Positive Rate")
plt.ylabel("True Positive Rate")
plt.title("ROC Curve")
plt.legend(loc="best")
plt.show()

#The below defines a narrow range of thresholds and provides the metrics as they change.
thresholds = [0.50, 0.52, 0.54, 0.56, 0.58]
for threshold in thresholds:
    y_pred_custom = (y_pred_proba >= threshold).astype(int)
    precision = precision_score(y_test, y_pred_custom)
    recall = recall_score(y_test, y_pred_custom)
    f1 = f1_score(y_test, y_pred_custom)
    print(f"Threshold: {threshold:.2f} | Precision: {precision:.2f} | Recall: {recall:.2f} | F1 Score: {f1:.2f}")

# Now implement the SVM model and run that across the dataset.
# Train the SVM model
svm_model = SVC(probability=True)  # Enable probability estimates for AUC-ROC
svm_model.fit(X_train, y_train)

# Predictions
y_pred = svm_model.predict(X_test)
y_pred_proba = svm_model.predict_proba(X_test)[:, 1]  # For AUC-ROC

# Calculate metrics
accuracy = accuracy_score(y_test, y_pred)
precision = precision_score(y_test, y_pred, average='weighted')
recall = recall_score(y_test, y_pred, average='weighted')
f1 = f1_score(y_test, y_pred, average='weighted')
auc = roc_auc_score(y_test, y_pred_proba)

# Print metrics
print(f"Accuracy: {accuracy:.2f}")
print(f"Precision: {precision:.2f}")
print(f"Recall: {recall:.2f}")
print(f"F1 Score: {f1:.2f}")
print(f"AUC-ROC: {auc:.2f}")

# ROC Curve
fpr, tpr, thresholds = roc_curve(y_test, y_pred_proba)
plt.plot(fpr, tpr, label="ROC Curve (AUC = {:.2f})".format(auc))
plt.xlabel("False Positive Rate")
plt.ylabel("True Positive Rate")
plt.title("ROC Curve")
plt.legend(loc="best")
plt.show()

# Confusion Matrix
cm = confusion_matrix(y_test, y_pred)
print("Confusion Matrix:\n", cm)
# Display the confusion matrix
disp = ConfusionMatrixDisplay(confusion_matrix=cm, display_labels=model.classes_)
disp.plot(cmap="Blues")  # You can change the colormap if you like
plt.title("Confusion Matrix")
plt.show()

# Define the parameter grid
param_grid = {
    'C': [0.1, 1, 10, 100],  # Regularization parameter
    'kernel': ['linear', 'rbf', 'poly'],  # Kernels to try
    'gamma': ['scale', 'auto'],  # Kernel coefficient
    'degree': [2, 3]  # Only relevant for 'poly' kernel
}

# Create a GridSearchCV object
grid_search = GridSearchCV(
    estimator=SVC(probability=True),
    param_grid=param_grid,
    scoring='f1_weighted',  # Optimize for F1 score
    cv=5,  # 5-fold cross-validation
    verbose=1
)

# Perform the grid search
grid_search.fit(X_train, y_train)

# Best parameters and score
print("Best Hyperparameters:", grid_search.best_params_)
print("Best Cross-Validation F1 Score:", grid_search.best_score_)

# Use the best model for evaluation
best_svm_model = grid_search.best_estimator_

# Evaluate the final model
y_pred = best_svm_model.predict(X_test)
y_pred_proba = best_svm_model.predict_proba(X_test)[:, 1]

# Print out the metrics
accuracy = accuracy_score(y_test, y_pred)
precision = precision_score(y_test, y_pred, average='weighted')
recall = recall_score(y_test, y_pred, average='weighted')
f1 = f1_score(y_test, y_pred, average='weighted')
auc = roc_auc_score(y_test, y_pred_proba)

# Print metrics
print(f"Accuracy: {accuracy:.2f}")
print(f"Precision: {precision:.2f}")
print(f"Recall: {recall:.2f}")
print(f"F1 Score: {f1:.2f}")
print(f"AUC-ROC: {auc:.2f}")

# ROC Curve
fpr, tpr, thresholds = roc_curve(y_test, y_pred_proba)
plt.plot(fpr, tpr, label="ROC Curve (AUC = {:.2f})".format(auc))
plt.xlabel("False Positive Rate")
plt.ylabel("True Positive Rate")
plt.title("ROC Curve")
plt.legend(loc="best")
plt.show()

# Confusion Matrix
cm = confusion_matrix(y_test, y_pred)
print("Confusion Matrix:\n", cm)

# Display the confusion matrix
disp = ConfusionMatrixDisplay(confusion_matrix=cm, display_labels=best_svm_model.classes_)
disp.plot(cmap="Blues")  # You can change the colormap if you like
plt.title("Confusion Matrix")
plt.show()