5.3.1 Python for Data Science and Machine Learning (Capstone Project)¶

Submitted by: Jose Regino Quiambao¶

Import Basic Libraries¶

In [1]:
#Import Pandas, NumPy, Matplotlib, Seaborn
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns

Import ML Libraries¶

In [15]:
#Import Scikit-Learn and XGBoost
from sklearn.model_selection import train_test_split, KFold, cross_val_score, GridSearchCV
from sklearn.metrics import roc_auc_score, accuracy_score, confusion_matrix, f1_score
from sklearn.preprocessing import MinMaxScaler, StandardScaler
from sklearn.ensemble import RandomForestClassifier
from sklearn.decomposition import PCA
from sklearn.svm import SVC
from sklearn.linear_model import LogisticRegression
from sklearn.neighbors import KNeighborsClassifier
from sklearn.cluster import KMeans
import xgboost as xgb

Exploratory Analysis¶

In [3]:
#Load the Dataset
df = pd.read_csv('SP901_CS_completedata.csv', sep=';', encoding='utf-8')
df
Out[3]:
PatientID Failure.binary Entropy_cooc.W.ADC GLNU_align.H.PET Min_hist.PET Max_hist.PET Mean_hist.PET Variance_hist.PET Standard_Deviation_hist.PET Skewness_hist.PET ... LZLGE.W.ADC LZHGE.W.ADC GLNU_area.W.ADC ZSNU.W.ADC ZSP.W.ADC GLNU_norm.W.ADC ZSNU_norm.W.ADC GLVAR_area.W.ADC ZSVAR.W.ADC Entropy_area.W.ADC
0 1 0 12.85352 46.256345 6.249117 17.825541 9.783773 6.814365 2.612479 0.688533 ... 0.006900 6201.93480 4.134000 239.289380 0.979180 0.018990 0.955860 1145.104960 0.025860 6.286320
1 2 1 12.21115 27.454540 11.005214 26.469077 15.426640 12.932074 3.598298 0.789526 ... 0.004230 16054.01263 8.376270 644.737020 0.956370 0.014610 0.932880 847.525370 0.041530 6.778530
2 3 0 12.75682 90.195696 2.777718 6.877486 4.295330 0.923425 0.962163 0.248637 ... 0.004530 6674.63840 13.116860 1165.702610 0.972680 0.025010 0.915370 1923.857050 0.071040 7.156850
3 4 1 13.46730 325.643330 6.296588 22.029843 10.334779 6.649795 2.580759 0.832011 ... 0.008880 17172.90951 23.847260 2760.412930 0.972030 0.010690 0.946580 1329.952900 0.038480 7.295210
4 5 0 12.63733 89.579042 3.583846 7.922501 4.454175 0.572094 0.757225 1.574845 ... 0.004050 13231.94294 8.144370 784.597290 0.964690 0.025260 0.937690 1116.386690 0.052230 7.051490
... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ... ...
192 193 0 11.95184 32.691265 12.213982 33.473818 16.208968 5.519894 3.300192 2.555792 ... -0.030152 38716.62434 29.548066 2709.805588 1.912158 -0.011360 1.864258 1972.882222 0.034874 13.936042
193 194 0 9.88702 60.481188 8.860044 21.524942 12.410224 4.226854 2.888088 1.822704 ... 0.028016 32577.07934 33.427936 4222.228292 1.842530 -0.014724 1.847794 3482.537034 0.055972 14.884836
194 195 0 12.84907 82.701566 11.543354 39.525156 18.867790 11.240758 4.716164 1.218402 ... -0.031432 24264.56392 163.442656 13608.298890 1.854554 -0.010420 1.769414 2216.464110 0.118514 14.217008
195 196 0 12.44606 72.223728 14.413852 49.234694 24.682116 25.364880 7.097038 1.554986 ... -0.028798 23482.25672 77.169252 7889.657780 1.910132 -0.011150 1.856080 1931.545086 0.037302 14.145064
196 197 0 12.13425 109.304666 11.545814 39.527616 18.870250 11.243218 4.718624 1.220862 ... -0.028972 24264.56638 163.445116 13608.301350 1.857014 -0.007960 1.771874 2216.466570 0.120974 14.219468

197 rows × 430 columns

In [4]:
#Checking if there's null
df['Failure.binary'].info()

<class 'pandas.core.series.Series'>
RangeIndex: 197 entries, 0 to 196
Series name: Failure.binary
Non-Null Count  Dtype
--------------  -----
197 non-null    int64
dtypes: int64(1)
memory usage: 1.7 KB
In [5]:
#Describe the target column
df['Failure.binary'].describe()
Out[5]:
count    197.000000
mean       0.340102
std        0.474950
min        0.000000
25%        0.000000
50%        0.000000
75%        1.000000
max        1.000000
Name: Failure.binary, dtype: float64
In [6]:
#Value counts of the target column
df['Failure.binary'].value_counts()
Out[6]:
0    130
1     67
Name: Failure.binary, dtype: int64

Splitting the data¶

In [7]:
# Split the data into features (X) and the target variable (y)
X = df.drop(columns=['PatientID', 'Failure.binary'])
y = df['Failure.binary']
In [8]:
print(X)
print(y)

     Entropy_cooc.W.ADC  GLNU_align.H.PET  Min_hist.PET  Max_hist.PET  \
0              12.85352         46.256345      6.249117     17.825541   
1              12.21115         27.454540     11.005214     26.469077   
2              12.75682         90.195696      2.777718      6.877486   
3              13.46730        325.643330      6.296588     22.029843   
4              12.63733         89.579042      3.583846      7.922501   
..                  ...               ...           ...           ...   
192            11.95184         32.691265     12.213982     33.473818   
193             9.88702         60.481188      8.860044     21.524942   
194            12.84907         82.701566     11.543354     39.525156   
195            12.44606         72.223728     14.413852     49.234694   
196            12.13425        109.304666     11.545814     39.527616   

     Mean_hist.PET  Variance_hist.PET  Standard_Deviation_hist.PET  \
0         9.783773           6.814365                     2.612479   
1        15.426640          12.932074                     3.598298   
2         4.295330           0.923425                     0.962163   
3        10.334779           6.649795                     2.580759   
4         4.454175           0.572094                     0.757225   
..             ...                ...                          ...   
192      16.208968           5.519894                     3.300192   
193      12.410224           4.226854                     2.888088   
194      18.867790          11.240758                     4.716164   
195      24.682116          25.364880                     7.097038   
196      18.870250          11.243218                     4.718624   

     Skewness_hist.PET  Kurtosis_hist.PET  Energy_hist.PET  ...  LZLGE.W.ADC  \
0             0.688533          -0.339727         0.005095  ...     0.006900   
1             0.789526          -0.319613         0.006297  ...     0.004230   
2             0.248637          -0.944246         0.005015  ...     0.004530   
3             0.832011           0.855861         0.003289  ...     0.008880   
4             1.574845           3.250288         0.008066  ...     0.004050   
..                 ...                ...              ...  ...          ...   
192           2.555792           3.463256        -0.030428  ...    -0.030152   
193           1.822704          -0.009474        -0.024528  ...     0.028016   
194           1.218402          -0.030052        -0.031068  ...    -0.031432   
195           1.554986           0.320610        -0.027516  ...    -0.028798   
196           1.220862          -0.027592        -0.028608  ...    -0.028972   

     LZHGE.W.ADC  GLNU_area.W.ADC    ZSNU.W.ADC  ZSP.W.ADC  GLNU_norm.W.ADC  \
0     6201.93480         4.134000    239.289380   0.979180         0.018990   
1    16054.01263         8.376270    644.737020   0.956370         0.014610   
2     6674.63840        13.116860   1165.702610   0.972680         0.025010   
3    17172.90951        23.847260   2760.412930   0.972030         0.010690   
4    13231.94294         8.144370    784.597290   0.964690         0.025260   
..           ...              ...           ...        ...              ...   
192  38716.62434        29.548066   2709.805588   1.912158        -0.011360   
193  32577.07934        33.427936   4222.228292   1.842530        -0.014724   
194  24264.56392       163.442656  13608.298890   1.854554        -0.010420   
195  23482.25672        77.169252   7889.657780   1.910132        -0.011150   
196  24264.56638       163.445116  13608.301350   1.857014        -0.007960   

     ZSNU_norm.W.ADC  GLVAR_area.W.ADC  ZSVAR.W.ADC  Entropy_area.W.ADC  
0           0.955860       1145.104960     0.025860            6.286320  
1           0.932880        847.525370     0.041530            6.778530  
2           0.915370       1923.857050     0.071040            7.156850  
3           0.946580       1329.952900     0.038480            7.295210  
4           0.937690       1116.386690     0.052230            7.051490  
..               ...               ...          ...                 ...  
192         1.864258       1972.882222     0.034874           13.936042  
193         1.847794       3482.537034     0.055972           14.884836  
194         1.769414       2216.464110     0.118514           14.217008  
195         1.856080       1931.545086     0.037302           14.145064  
196         1.771874       2216.466570     0.120974           14.219468  

[197 rows x 428 columns]
0      0
1      1
2      0
3      1
4      0
      ..
192    0
193    0
194    0
195    0
196    0
Name: Failure.binary, Length: 197, dtype: int64

Scaling the Data¶

In [9]:
# Scaling the data
scaler = MinMaxScaler()
scaled_data = scaler.fit_transform(X)

PCA¶

In [10]:
# Fit the PCA model
pca = PCA(n_components=2)  # Choose the number of principal components
pca_components = pca.fit_transform(scaled_data)

# Analyze the principal components
explained_variance = pca.explained_variance_ratio_
print(explained_variance)

principaldf = pd.DataFrame(data = pca_components, columns=['Principal Component 1','Principal Component 2'])
print("")
print(principaldf)

plt.scatter(pca_components[:, 0], pca_components[:, 1])
plt.xlabel('Principal Component 1')
plt.ylabel('Principal Component 2')
plt.title('Scatter Plot of Principal Components')
plt.show()

# Project the original data onto the principal components
pca_scores = pca.transform(scaled_data)

# Analyze the principal component scores
from sklearn.cluster import KMeans

kmeans = KMeans(n_clusters=3)
kmeans_labels = kmeans.fit_predict(pca_scores)

[0.67558736 0.08759435]

     Principal Component 1  Principal Component 2
0                -1.999969              -0.662641
1                -1.966676              -0.227635
2                -1.943384              -1.076428
3                -2.022240              -0.109832
4                -2.228796              -1.305999
..                     ...                    ...
192               5.367085               2.374854
193               5.766655               1.646730
194               5.709939               3.576730
195               5.703741               3.448861
196               5.758057               3.472670

[197 rows x 2 columns]

/Users/jr/anaconda3/lib/python3.11/site-packages/sklearn/cluster/_kmeans.py:1412: FutureWarning: The default value of `n_init` will change from 10 to 'auto' in 1.4. Set the value of `n_init` explicitly to suppress the warning
  super()._check_params_vs_input(X, default_n_init=10)
In [11]:
X = principaldf

X_train, X_temp, y_train, y_temp = train_test_split(X, y, test_size=0.4, random_state=42)
X_val, X_test, y_val, y_test = train_test_split(X_temp, y_temp, test_size=0.5, random_state=42)

# The data is first split into a training set (60%), a temporary set (40%).
# Then the temporary set is further split into validation (20%) and test (20%) sets.

print('Training set size:', len(X_train))
print('Validation set size:', len(X_val))
print('Testing set size:', len(X_test))

Training set size: 118
Validation set size: 39
Testing set size: 40

10-fold Cross Validation (Logistic Regression)¶

In [12]:
# Create the Logistic Regression model
model = LogisticRegression()

# Initialize the K-Fold cross-validation iterator with 10 folds
kf = KFold(n_splits=10, shuffle=True, random_state=42)

# Perform 10-fold cross-validation and calculate the accuracy for each fold
accuracies = cross_val_score(model, X, y, cv=kf, scoring='accuracy')

# Print the accuracy for each fold and the mean accuracy
for fold, accuracy in enumerate(accuracies, 1):
    print(f"Fold {fold}: Accuracy = {accuracy:.2f}")

mean_accuracy = np.mean(accuracies)
print(f"Mean Accuracy: {mean_accuracy:.2f}")

Fold 1: Accuracy = 0.60
Fold 2: Accuracy = 0.55
Fold 3: Accuracy = 0.65
Fold 4: Accuracy = 0.70
Fold 5: Accuracy = 0.75
Fold 6: Accuracy = 0.70
Fold 7: Accuracy = 0.70
Fold 8: Accuracy = 0.63
Fold 9: Accuracy = 0.74
Fold 10: Accuracy = 0.58
Mean Accuracy: 0.66

Grid Search¶

In [13]:
param_grid = {
    'C': [0.1, 1, 10],
    'kernel': ['linear', 'rbf'],
    'gamma': [0.001, 0.01, 0.1, 1],
}

# Creates an instance of the model to be used
svm = SVC()

# Creates a GridSearchCV object and specify the model, parameter grid, and the number of cross-validation folds
grid_search = GridSearchCV(svm, param_grid, cv=5)

# Fit the grid search object to your training data
grid_search.fit(X_train, y_train)

# Access the best parameters and the best estimator
best_params = grid_search.best_params_
best_estimator = grid_search.best_estimator_

# Evaluate the best model on your test data
accuracy = best_estimator.score(X_test, y_test)
print(f"Best model accuracy on the test set: {accuracy}")
print("Best parameters:", best_params)

Best model accuracy on the test set: 0.625
Best parameters: {'C': 0.1, 'gamma': 0.001, 'kernel': 'linear'}

Classification (Logistic Regression, K-Nearest Neighbors, SVM, Random Forest, XGBoost)¶

In [28]:
# Initialize and train the models
models = {
    'Logistic Regression': LogisticRegression(),
    'K-Nearest Neighbors': KNeighborsClassifier(n_neighbors=3),
    'SVM': SVC(),
    'Random Forest': RandomForestClassifier(),
    'XGBoost': xgb.XGBClassifier()
}

for model_name, model in models.items():
    model.fit(X_train, y_train)
    y_pred = model.predict(X_test)
    accuracy = accuracy_score(y_test, y_pred)
    print(f'{model_name} Accuracy: {accuracy:.2f}')
    
# Create a bar chart
models = list(accuracy_results.keys())
accuracy_values = list(accuracy_results.values())

plt.figure(figsize=(10, 6))
bars = plt.bar(models, accuracy_values, color='lightblue')
plt.ylim(0.1, 0.8)
plt.title('Accuracy of Different Classification Models')
plt.xlabel('Model')
plt.ylabel('Accuracy')
plt.xticks(rotation=45)
#Add data labels to the bars
for bar, accuracy in zip(bars, accuracy_values):
    plt.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.005, f'{accuracy:.2f}', ha='center', va='bottom')
plt.show()

Logistic Regression Accuracy: 0.62
K-Nearest Neighbors Accuracy: 0.61
SVM Accuracy: 0.62
Random Forest Accuracy: 0.58
XGBoost Accuracy: 0.52

Evaluation Metrics (ROC-AUC, Accuracy, F-1 Score)¶

In [38]:
# Initialize dictionaries to store the evaluation metrics
accuracy_results = {}
f1_results = {}
roc_auc_results = {}

# Train and evaluate the models
models = {
    'Logistic Regression': LogisticRegression(),
    'K-Nearest Neighbors': KNeighborsClassifier(n_neighbors=3),
    'SVM': SVC(probability=True),
    'Random Forest': RandomForestClassifier(),
    'XGBoost': xgb.XGBClassifier()
}

for model_name, model in models.items():
    model.fit(X_train, y_train)
    y_pred = model.predict(X_test)
    accuracy = accuracy_score(y_test, y_pred)
    f1 = f1_score(y_test, y_pred, average='weighted')
    if 'SVM' in model_name:
        y_score = model.predict_proba(X_test)[:, 1]
    else:
        y_score = model.decision_function(X_test) if hasattr(model, 'decision_function') else model.predict_proba(X_test)[:, 1]
    roc_auc = roc_auc_score(y_test, y_score)
    
    accuracy_results[model_name] = accuracy
    f1_results[model_name] = f1
    roc_auc_results[model_name] = roc_auc

# Create a bar chart for each evaluation metric
models = list(accuracy_results.keys())

plt.figure(figsize=(12, 6))

# Accuracy bar chart
plt.subplot(131)
accuracy_values = list(accuracy_results.values())
bars = plt.bar(models, accuracy_values, color='skyblue')
plt.ylim(0.1, 1.0)
plt.title('Accuracy')
plt.xticks(rotation=90)

for bar, accuracy in zip(bars, accuracy_values):
    plt.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.005, f'{accuracy:.2f}', ha='center', va='bottom')

# F1-score bar chart
plt.subplot(132)
f1_values = list(f1_results.values())
bars = plt.bar(models, f1_values, color='lightcoral')
plt.ylim(0.1, 1.0)
plt.title('F1-Score')
plt.xticks(rotation=90)

for bar, f1 in zip(bars, f1_values):
    plt.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.005, f'{f1:.2f}', ha='center', va='bottom')

# ROC-AUC bar chart
plt.subplot(133)
roc_auc_values = list(roc_auc_results.values())
bars = plt.bar(models, roc_auc_values, color='lightgreen')
plt.ylim(0.1, 1.0)
plt.title('ROC-AUC')
plt.xticks(rotation=90)

for bar, roc_auc in zip(bars, roc_auc_values):
    plt.text(bar.get_x() + bar.get_width() / 2, bar.get_height() + 0.005, f'{roc_auc:.2f}', ha='center', va='bottom')

plt.tight_layout()
plt.show()