binary_classifier.py

import os
import sys

import math
from itertools import product

import numpy
import pandas

import lightgbm

from sklearn.preprocessing import MinMaxScaler
from sklearn.model_selection import train_test_split
from sklearn.metrics import accuracy_score, f1_score, precision_score, recall_score, balanced_accuracy_score # evan got mad at me for balanced_accuracy :( but alik said its fine so take that FWIW

from scipy.stats import pearsonr

import shap

import gridcv

'''
Vars
'''
rootPath = ""
dataPath = ""
targetPath = ""
savePath = ""

artificialBalance = False # if classes are imbalanced, set true to artificially inflate minority class. lightgbm includes a way to penalize minority class heavily, so not needed necessarily and not implemented
dropna = False 
testSetSize = 0.2
validSetSize = 0.2

subject = ""
session = ""

foldAccuracy = pandas.DataFrame(columns = ['Fold 1', 'Fold 2', 'Fold 3', 'Fold 4', 'Fold 5'])
foldPrecision = pandas.DataFrame(columns = ['Fold 1', 'Fold 2', 'Fold 3', 'Fold 4', 'Fold 5'])
foldRecall = pandas.DataFrame(columns = ['Fold 1', 'Fold 2', 'Fold 3', 'Fold 4', 'Fold 5'])
foldF1 = pandas.DataFrame(columns = ['Fold 1', 'Fold 2', 'Fold 3', 'Fold 4', 'Fold 5'])

'''
Data Loading
'''
input = pandas.read_csv(dataPath)
targets = pandas.read_csv(targetPath)

'''
Data Preprocessing
'''

if dropna:
    input = input.dropna(0)
    print("dropped nan rows")

'''
Data Splitting 
'''
xTrain, xTest, yTrain, yTest = train_test_split(input, 
                                                    targets, 
                                                    test_size=testSetSize, 
                                                    shuffle=False) #shuffle set to false to prevent temporal leakage


'''
5 fold CV splitting
'''
# in Liu et al., 2024: for each subject and session 5 models trained and assessed for each fold with shap values extracted
#                      hyperparamter tuning is also performed for each fold on the train/valid

indicesPerWindow = math.ceil(xTrain.shape[0] / 9)
indices = [i * indicesPerWindow for i in range(9 + 1)]
indices[-1] = xTrain.shape[0]

fold1XTrain = xTrain.iloc[indices[0]:indices[4], :]
fold1YTrain = yTrain.iloc[indices[0]:indices[4], :]
fold1XTest = xTrain.iloc[indices[4]:indices[5], :]
fold1YTest = yTrain.iloc[indices[4]:indices[5], :]

fold2XTrain = xTrain.iloc[indices[1]:indices[5], :]
fold2YTrain = yTrain.iloc[indices[1]:indices[5], :]
fold2XTest = xTrain.iloc[indices[5]:indices[6], :]
fold2YTest = yTrain.iloc[indices[5]:indices[6], :]

fold3XTrain = xTrain.iloc[indices[2]:indices[6], :]
fold3YTrain = yTrain.iloc[indices[2]:indices[6], :]
fold3XTest = xTrain.iloc[indices[6]:indices[7], :]
fold3YTest = yTrain.iloc[indices[6]:indices[7], :]

fold4XTrain = xTrain.iloc[indices[3]:indices[7], :]
fold4YTrain = yTrain.iloc[indices[3]:indices[7], :]
fold4XTest = xTrain.iloc[indices[7]:indices[8], :]
fold4YTest = yTrain.iloc[indices[7]:indices[8], :]

fold5XTrain = xTrain.iloc[indices[4]:indices[8], :]
fold5YTrain = yTrain.iloc[indices[4]:indices[8], :]
fold5XTest = xTrain.iloc[indices[8]:indices[9], :]
fold5YTest = yTrain.iloc[indices[8]:indices[9], :]

'''
Model Training
'''
# There exists more robust hyperparameter tuning methods (bayesian etc.) Grid search is just fine for our purposes.
paramGrid = {
    'num_leaves': [3, 5, 13, 29],
    'bagging_fraction': [ 0.7, 0.8, 0.9, 1.0],
    'bagging_freq': [0, 4, 8, 12],
    'feature_fraction': [0.7, 0.8, 0.9, 1.0],
    'learning_rate': numpy.logspace(-3, -1, 5),
}

combs = list(product(
    paramGrid['num_leaves'],
    paramGrid['bagging_fraction'],
    paramGrid['bagging_freq'],
    paramGrid['feature_fraction'],
    paramGrid['learning_rate']
))
print("tuning hyperparameters")
bestParamsFold1, _ = gridcv.tuneLGBClassifierBinaryParameters(combs, fold1XTrain, fold1YTrain, fold1XTest, fold1YTest)
bestParamsFold2, _ = gridcv.tuneLGBClassifierBinaryParameters(combs, fold2XTrain, fold2YTrain, fold2XTest, fold2YTest)
bestParamsFold3, _ = gridcv.tuneLGBClassifierBinaryParameters(combs, fold3XTrain, fold3YTrain, fold3XTest, fold3YTest)
bestParamsFold4, _ = gridcv.tuneLGBClassifierBinaryParameters(combs, fold4XTrain, fold4YTrain, fold4XTest, fold4YTest)
bestParamsFold5, _ = gridcv.tuneLGBClassifierBinaryParameters(combs, fold5XTrain, fold5YTrain, fold5XTest, fold5YTest)

fold1GBMData = lightgbm.Dataset(fold1XTrain, label=fold1YTrain)
fold1Model = lightgbm.train(bestParamsFold1, fold1GBMData)
fold2GBMData = lightgbm.Dataset(fold2XTrain, label=fold2YTrain)
fold2Model = lightgbm.train(bestParamsFold2, fold2GBMData)
fold3GBMData = lightgbm.Dataset(fold3XTrain, label=fold3YTrain)
fold3Model = lightgbm.train(bestParamsFold3, fold3GBMData)
fold4GBMData = lightgbm.Dataset(fold4XTrain, label=fold4YTrain)
fold4Model = lightgbm.train(bestParamsFold4, fold4GBMData)
fold5GBMData = lightgbm.Dataset(fold5XTrain, label=fold5YTrain)
fold5Model = lightgbm.train(bestParamsFold5, fold5GBMData)

'''
Model testing

'''
fold1Preds = fold1Model.predict(fold1XTest)
fold2Preds = fold2Model.predict(fold2XTest)
fold3Preds = fold3Model.predict(fold3XTest)
fold4Preds = fold4Model.predict(fold4XTest)
fold5Preds = fold5Model.predict(fold5XTest)


'''
Model evaluation
'''
fold1PredLabels = (fold1Preds >= 0.5).astype(int)
fold2PredLabels = (fold2Preds >= 0.5).astype(int)
fold3PredLabels = (fold3Preds >= 0.5).astype(int)
fold4PredLabels = (fold4Preds >= 0.5).astype(int)
fold5PredLabels = (fold5Preds >= 0.5).astype(int)

foldAccuracy.at[0, 'Fold 1'] = accuracy_score(fold1YTest, fold1PredLabels)
foldAccuracy.at[0, 'Fold 2'] = accuracy_score(fold2YTest, fold2PredLabels)
foldAccuracy.at[0, 'Fold 3'] = accuracy_score(fold3YTest, fold3PredLabels)
foldAccuracy.at[0, 'Fold 4'] = accuracy_score(fold4YTest, fold4PredLabels)
foldAccuracy.at[0, 'Fold 5'] = accuracy_score(fold5YTest, fold5PredLabels)

foldPrecision.at[0, 'Fold 1'] = precision_score(fold1YTest, fold1PredLabels, zero_division=0)
foldPrecision.at[0, 'Fold 2'] = precision_score(fold2YTest, fold2PredLabels, zero_division=0)
foldPrecision.at[0, 'Fold 3'] = precision_score(fold3YTest, fold3PredLabels, zero_division=0)
foldPrecision.at[0, 'Fold 4'] = precision_score(fold4YTest, fold4PredLabels, zero_division=0)
foldPrecision.at[0, 'Fold 5'] = precision_score(fold5YTest, fold5PredLabels, zero_division=0)

foldRecall.at[0, 'Fold 1'] = recall_score(fold1YTest, fold1PredLabels, zero_division=0)
foldRecall.at[0, 'Fold 2'] = recall_score(fold2YTest, fold2PredLabels, zero_division=0)
foldRecall.at[0, 'Fold 3'] = recall_score(fold3YTest, fold3PredLabels, zero_division=0)
foldRecall.at[0, 'Fold 4'] = recall_score(fold4YTest, fold4PredLabels, zero_division=0)
foldRecall.at[0, 'Fold 5'] = recall_score(fold5YTest, fold5PredLabels, zero_division=0)

foldF1.at[0, 'Fold 1'] = f1_score(fold1YTest, fold1PredLabels, zero_division=0)
foldF1.at[0, 'Fold 2'] = f1_score(fold2YTest, fold2PredLabels, zero_division=0)
foldF1.at[0, 'Fold 3'] = f1_score(fold3YTest, fold3PredLabels, zero_division=0)
foldF1.at[0, 'Fold 4'] = f1_score(fold4YTest, fold4PredLabels, zero_division=0)
foldF1.at[0, 'Fold 5'] = f1_score(fold5YTest, fold5PredLabels, zero_division=0)



'''
Pull Shap values
'''

explainer1 = shap.TreeExplainer(fold1Model)
shapValues1 = explainer1(fold1XTest)
explainer2 = shap.TreeExplainer(fold2Model)
shapValues2 = explainer2(fold2XTest)
explainer3 = shap.TreeExplainer(fold3Model)
shapValues3 = explainer3(fold3XTest)
explainer4 = shap.TreeExplainer(fold4Model)
shapValues4 = explainer4(fold4XTest)
explainer5 = shap.TreeExplainer(fold5Model)
shapValues5 = explainer5(fold5XTest)

'''
IO
'''

shapValues1 = pandas.DataFrame(shapValues1.values, columns = fold1XTest.columns)
shapValues1.to_csv(os.path.join(savePath, f"{subject}_{session}_shapValuesFold1.csv"), index=False)
shapValues2 = pandas.DataFrame(shapValues2.values, columns = fold2XTest.columns)
shapValues2.to_csv(os.path.join(savePath, f"{subject}_{session}_shapValuesFold2.csv"), index=False)
shapValues3 = pandas.DataFrame(shapValues3.values, columns = fold3XTest.columns)
shapValues3.to_csv(os.path.join(savePath, f"{subject}_{session}_shapValuesFold3.csv"), index=False)
shapValues4 = pandas.DataFrame(shapValues4.values, columns = fold4XTest.columns)
shapValues4.to_csv(os.path.join(savePath, f"{subject}_{session}_shapValuesFold4.csv"), index=False)
shapValues5 = pandas.DataFrame(shapValues5.values, columns = fold5XTest.columns)
shapValues5.to_csv(os.path.join(savePath, f"{subject}_{session}_shapValuesFold5.csv"), index=False)

foldAccuracy.to_csv(os.path.join(savePath, f"{subject}_{session}_foldAccuracy.csv"), index=False)
foldPrecision.to_csv(os.path.join(savePath, f"{subject}_{session}_foldPrecision.csv"), index=False)
foldRecall.to_csv(os.path.join(savePath, f"{subject}_{session}_foldRecall.csv"), index=False)
foldF1.to_csv(os.path.join(savePath, f"{subject}_{session}_foldF1.csv"), index=False)

fold1Model.save_model(os.path.join(savePath, f"{subject}_{session}_fold1Model.txt"))
fold2Model.save_model(os.path.join(savePath, f"{subject}_{session}_fold2Model.txt"))
fold3Model.save_model(os.path.join(savePath, f"{subject}_{session}_fold3Model.txt"))
fold4Model.save_model(os.path.join(savePath, f"{subject}_{session}_fold4Model.txt"))
fold5Model.save_model(os.path.join(savePath, f"{subject}_{session}_fold5Model.txt"))

'''
Evaluate on holdout set

Have to look at fold results and pick best model params
'''

holdoutParams = bestParamsFold3 # change this
holdoutGBMData = lightgbm.Dataset(xTrain, label=yTrain)
holdoutModel = lightgbm.train(holdoutParams, holdoutGBMData)
holdoutPreds = holdoutModel.predict(xTest)

holdoutPredLabels = (holdoutPreds >= 0.5).astype(int)
holdoutAccuracy = accuracy_score(yTest, holdoutPredLabels)
holdoutPrecision = precision_score(yTest, holdoutPredLabels)
holdoutRecall = recall_score(yTest, holdoutPredLabels)
holdoutF1 = f1_score(yTest, holdoutPredLabels)