diff --git a/scripts/ml/ml_functions/MultClfs_CVs.py b/scripts/ml/ml_functions/MultClfs_CVs.py new file mode 100755 index 0000000..2d0fb6e --- /dev/null +++ b/scripts/ml/ml_functions/MultClfs_CVs.py @@ -0,0 +1,453 @@ +#!/usr/bin/env python3 +# -*- coding: utf-8 -*- +""" +Created on Fri Mar 4 15:25:33 2022 + +@author: tanu +""" +#%% +import os, sys +import pandas as pd +import numpy as np +import pprint as pp +from copy import deepcopy +from sklearn import linear_model +from sklearn import datasets +from collections import Counter + +from sklearn.linear_model import LogisticRegression, LogisticRegressionCV +from sklearn.linear_model import RidgeClassifier, RidgeClassifierCV, SGDClassifier, PassiveAggressiveClassifier + +from sklearn.naive_bayes import BernoulliNB +from sklearn.neighbors import KNeighborsClassifier +from sklearn.svm import SVC +from sklearn.tree import DecisionTreeClassifier, ExtraTreeClassifier +from sklearn.ensemble import RandomForestClassifier, ExtraTreesClassifier, AdaBoostClassifier, GradientBoostingClassifier, BaggingClassifier +from sklearn.naive_bayes import GaussianNB +from sklearn.gaussian_process import GaussianProcessClassifier, kernels +from sklearn.gaussian_process.kernels import RBF, DotProduct, Matern, RationalQuadratic, WhiteKernel + +from sklearn.discriminant_analysis import LinearDiscriminantAnalysis, QuadraticDiscriminantAnalysis +from sklearn.neural_network import MLPClassifier + +from sklearn.svm import SVC +from xgboost import XGBClassifier +from sklearn.naive_bayes import MultinomialNB +from sklearn.preprocessing import StandardScaler, MinMaxScaler, OneHotEncoder + +from sklearn.compose import ColumnTransformer +from sklearn.compose import make_column_transformer + +from sklearn.metrics import make_scorer, confusion_matrix, accuracy_score, balanced_accuracy_score, precision_score, average_precision_score, recall_score +from sklearn.metrics import roc_auc_score, roc_curve, f1_score, matthews_corrcoef, jaccard_score, classification_report + +# added +from sklearn.model_selection import train_test_split, cross_validate, cross_val_score, LeaveOneOut, KFold, RepeatedKFold, cross_val_predict + +from sklearn.model_selection import train_test_split, cross_validate, cross_val_score +from sklearn.model_selection import StratifiedKFold,RepeatedStratifiedKFold, RepeatedKFold + +from sklearn.pipeline import Pipeline, make_pipeline + +from sklearn.feature_selection import RFE, RFECV + +import itertools +import seaborn as sns +import matplotlib.pyplot as plt + +from statistics import mean, stdev, median, mode + +from imblearn.over_sampling import RandomOverSampler +from imblearn.under_sampling import RandomUnderSampler +from imblearn.over_sampling import SMOTE +from sklearn.datasets import make_classification +from imblearn.combine import SMOTEENN +from imblearn.combine import SMOTETomek + +from imblearn.over_sampling import SMOTENC +from imblearn.under_sampling import EditedNearestNeighbours +from imblearn.under_sampling import RepeatedEditedNearestNeighbours + +from sklearn.model_selection import GridSearchCV +from sklearn.base import BaseEstimator +from sklearn.impute import KNNImputer as KNN +import json +import argparse +import re +import itertools +from sklearn.model_selection import LeaveOneGroupOut +from sklearn.decomposition import PCA +from sklearn.naive_bayes import ComplementNB +from sklearn.dummy import DummyClassifier + +#%% GLOBALS +#rs = {'random_state': 42} # INSIDE FUNCTION CALL NOW +#njobs = {'n_jobs': os.cpu_count() } # the number of jobs should equal the number of CPU cores + +scoring_fn = ({ 'mcc' : make_scorer(matthews_corrcoef) + , 'fscore' : make_scorer(f1_score) + , 'precision' : make_scorer(precision_score) + , 'recall' : make_scorer(recall_score) + , 'accuracy' : make_scorer(accuracy_score) + , 'roc_auc' : make_scorer(roc_auc_score) + , 'jcc' : make_scorer(jaccard_score) + }) +# for sel_cv INSIDE FUNCTION CALL NOW +#skf_cv = StratifiedKFold(n_splits = 10 +# #, shuffle = False, random_state= None) +# , shuffle = True, **rs) + +#rskf_cv = RepeatedStratifiedKFold(n_splits = 10 +# , n_repeats = 3 +# , **rs) + +mcc_score_fn = {'mcc': make_scorer(matthews_corrcoef)} +jacc_score_fn = {'jcc': make_scorer(jaccard_score)} + +############################################################################### +score_type_ordermapD = { 'mcc' : 1 + , 'fscore' : 2 + , 'jcc' : 3 + , 'precision' : 4 + , 'recall' : 5 + , 'accuracy' : 6 + , 'roc_auc' : 7 + , 'TN' : 8 + , 'FP' : 9 + , 'FN' : 10 + , 'TP' : 11 + , 'trainingY_neg': 12 + , 'trainingY_pos': 13 + , 'blindY_neg' : 14 + , 'blindY_pos' : 15 + , 'fit_time' : 16 + , 'score_time' : 17 + } + +scoreCV_mapD = {'test_mcc' : 'MCC' + , 'test_fscore' : 'F1' + , 'test_precision' : 'Precision' + , 'test_recall' : 'Recall' + , 'test_accuracy' : 'Accuracy' + , 'test_roc_auc' : 'ROC_AUC' + , 'test_jcc' : 'JCC' + } + +scoreBT_mapD = {'bts_mcc' : 'MCC' + , 'bts_fscore' : 'F1' + , 'bts_precision' : 'Precision' + , 'bts_recall' : 'Recall' + , 'bts_accuracy' : 'Accuracy' + , 'bts_roc_auc' : 'ROC_AUC' + , 'bts_jcc' : 'JCC' + } + +#gene_group = 'gene_name' +#%%############################################################################ +############################ +# MultModelsCl() +# Run Multiple Classifiers +############################ +# Multiple Classification - Model Pipeline +def MultModelsCl_CVs(input_df + , target + , tts_split_type + , resampling_type + #, group = None + , skf_cv_threshold = 10 #[None, 3, 5, 10] + + , add_cm = True # adds confusion matrix based on cross_val_predict + , add_yn = True # adds target var class numbers + , var_type = ['numerical', 'categorical','mixed'] + , scale_numeric = ['min_max', 'std', 'min_max_neg', 'none'] + + , return_formatted_output = True + + , random_state = 42 + , n_jobs = os.cpu_count() # the number of jobs should equal the number of CPU cores + ): + + ''' + @ param input_df: input features + @ type: df with input features WITHOUT the target variable + + @param target: target (or output) feature + @type: df or np.array or Series + + @param skv_cv: stratifiedK fold int or object to allow shuffle and random state to pass + @type: int or StratifiedKfold() + + @var_type: numerical, categorical and mixed to determine what col_transform to apply (MinMaxScalar and/or one-hot encoder) + @type: list + + returns + Dict containing multiple classification scores for each model and mean of each Stratified Kfold including training + ''' + +#%% Func globals + rs = {'random_state': random_state} + njobs = {'n_jobs': n_jobs} + + skf_cv = StratifiedKFold(n_splits = skf_cv_threshold + #, shuffle = False, random_state= None) + , shuffle = True,**rs) + + # rskf_cv = RepeatedStratifiedKFold(n_splits = skf_cv_threshold + # , n_repeats = 3 + # , **rs) + # logo = LeaveOneGroupOut() + + # select CV type: + # if group == None: + # sel_cv = skf_cv + # else: + # sel_cv = logo + #====================================================== + # Determine categorical and numerical features + #====================================================== + numerical_ix = input_df.select_dtypes(include=['int64', 'float64']).columns + numerical_ix + categorical_ix = input_df.select_dtypes(include=['object', 'bool']).columns + categorical_ix + + #====================================================== + # Determine preprocessing steps ~ var_type + #====================================================== + + if type(var_type) == list: + var_type = str(var_type[0]) + else: + var_type = var_type + + if var_type in ['numerical','mixed']: + if scale_numeric == ['none']: + t = [('cat', OneHotEncoder(), categorical_ix)] + if scale_numeric != ['none']: + if scale_numeric == ['min_max']: + scaler = MinMaxScaler() + if scale_numeric == ['min_max_neg']: + scaler = MinMaxScaler(feature_range=(-1, 1)) + if scale_numeric == ['std']: + scaler = StandardScaler() + + t = [('num', scaler, numerical_ix) + , ('cat', OneHotEncoder(), categorical_ix)] + + + if var_type == 'categorical': + t = [('cat', OneHotEncoder(), categorical_ix)] + + + col_transform = ColumnTransformer(transformers = t + , remainder='passthrough') + + + #====================================================== + # Specify multiple Classification Models + #====================================================== + models = [('AdaBoost Classifier' , AdaBoostClassifier(**rs) ) + # , ('Bagging Classifier' , BaggingClassifier(**rs, **njobs, bootstrap = True, oob_score = True, verbose = 3, n_estimators = 100) ) + # #, ('Bernoulli NB' , BernoulliNB() ) # pks Naive Bayes, CAUTION + # , ('Complement NB' , ComplementNB() ) + # , ('Decision Tree' , DecisionTreeClassifier(**rs) ) + # , ('Extra Tree' , ExtraTreeClassifier(**rs) ) + # , ('Extra Trees' , ExtraTreesClassifier(**rs) ) + # , ('Gradient Boosting' , GradientBoostingClassifier(**rs) ) + # , ('Gaussian NB' , GaussianNB() ) + # , ('Gaussian Process' , GaussianProcessClassifier(**rs) ) + # , ('K-Nearest Neighbors' , KNeighborsClassifier() ) + # , ('LDA' , LinearDiscriminantAnalysis() ) + # , ('Logistic Regression' , LogisticRegression(**rs) ) + # , ('Logistic RegressionCV' , LogisticRegressionCV(cv = 3, **rs)) + # , ('MLP' , MLPClassifier(max_iter = 500, **rs) ) + # , ('Multinomial NB' , MultinomialNB() ) + # , ('Passive Aggresive' , PassiveAggressiveClassifier(**rs, **njobs) ) + # , ('QDA' , QuadraticDiscriminantAnalysis() ) + # , ('Random Forest' , RandomForestClassifier(**rs, n_estimators = 1000, **njobs ) ) + # , ('Random Forest2' , RandomForestClassifier(min_samples_leaf = 5 + # , n_estimators = 1000 + # , bootstrap = True + # , oob_score = True + # , **njobs + # , **rs + # , max_features = 'auto') ) + # , ('Ridge Classifier' , RidgeClassifier(**rs) ) + # , ('Ridge ClassifierCV' , RidgeClassifierCV(cv = 3) ) + # , ('SVC' , SVC(**rs) ) + # , ('Stochastic GDescent' , SGDClassifier(**rs, **njobs) ) + , ('XGBoost' , XGBClassifier(**rs, verbosity = 0, use_label_encoder = False, **njobs) ) + , ('Dummy Classifier' , DummyClassifier(strategy = 'most_frequent') ) + ] + + mm_skf_scoresD = {} + + print('\n==============================================================\n' + , '\nRunning several classification models (n):', len(models) + ,'\nList of models:') + for m in models: + print(m) + print('\n================================================================\n') + + index = 1 + for model_name, model_fn in models: + print('\nRunning classifier:', index + , '\nModel_name:' , model_name + , '\nModel func:' , model_fn) + index = index+1 + + model_pipeline = Pipeline([ + ('prep' , col_transform) + , ('model' , model_fn)]) + + print('\nRunning model pipeline:', model_pipeline) + cv_modD = cross_validate(model_pipeline + , input_df + , target + , cv = skf_cv + #, groups = group + , scoring = scoring_fn + , return_train_score = True) + #============================== + # Extract mean values for CV + #============================== + mm_skf_scoresD[model_name] = {} + + for key, value in cv_modD.items(): + print('\nkey:', key, '\nvalue:', value) + print('\nmean value:', np.mean(value)) + mm_skf_scoresD[model_name][key] = round(np.mean(value),2) + + # ADD more info: meta data related to input df + mm_skf_scoresD[model_name]['resampling'] = resampling_type + mm_skf_scoresD[model_name]['n_training_size'] = len(input_df) + mm_skf_scoresD[model_name]['n_trainingY_ratio'] = round(Counter(target)[0]/Counter(target)[1], 2) + mm_skf_scoresD[model_name]['n_features'] = len(input_df.columns) + mm_skf_scoresD[model_name]['tts_split'] = tts_split_type + + ####################################################################### + #====================================================== + # Option: Add confusion matrix from cross_val_predict + # Understand and USE with caution + #====================================================== + if add_cm: + cmD = {} + + # Calculate cm + y_pred = cross_val_predict(model_pipeline + , input_df + , target + , cv = skf_cv + #, groups = group + , **njobs) + #_tn, _fp, _fn, _tp = confusion_matrix(y_pred, y).ravel() # internally + tn, fp, fn, tp = confusion_matrix(y_pred, target).ravel() + + # Build cm dict + cmD = {'TN' : tn + , 'FP': fp + , 'FN': fn + , 'TP': tp} + + # Update cv dict cmD + mm_skf_scoresD[model_name].update(cmD) + + #============================================= + # Option: Add targety numbers for data + #============================================= + if add_yn: + tnD = {} + + # Build tn numbers dict + tnD = {'n_trainingY_neg' : Counter(target)[0] + , 'n_trainingY_pos' : Counter(target)[1] } + + # Update cv dict with cmD and tnD + mm_skf_scoresD[model_name].update(tnD) + +#%% + #return(mm_skf_scoresD) + #============================ + # Process the dict to have WF + #============================ + if return_formatted_output: + CV_BT_metaDF = ProcessMultModelsCl(mm_skf_scoresD, cv_threshold_suffix = skf_cv_threshold) + return(CV_BT_metaDF) + else: + return(mm_skf_scoresD) + +#%% Process output function ################################################### +############################ +# ProcessMultModelsCl() +############################ +#Processes the dict from above if use_formatted_output = True + +def ProcessMultModelsCl(inputD = {} + , cv_threshold_suffix = 10 + #, blind_test_data = True + ): + + scoresDF = pd.DataFrame(inputD) + + #------------------------ + # Extracting split_name + #----------------------- + tts_split_nameL = [] + for k,v in inputD.items(): + tts_split_nameL = tts_split_nameL + [v['tts_split']] + + if len(set(tts_split_nameL)) == 1: + tts_split_name = str(list(set(tts_split_nameL))[0]) + print('\nExtracting tts_split_name:', tts_split_name) + + #---------------------- + # WF: CV results + #---------------------- + scoresDFT = scoresDF.T + + scoresDF_CV = scoresDFT.filter(regex='^test_.*$', axis = 1); scoresDF_CV.columns + # map colnames for consistency to allow concatenting + scoresDF_CV.columns = scoresDF_CV.columns.map(scoreCV_mapD); scoresDF_CV.columns + #scoresDF_CV['source_data'] = 'CV' + scoresDF_CV['source_data'] = 'CV_' + str(cv_threshold_suffix) + + + #---------------------- + # WF: Meta data + #---------------------- + metaDF = scoresDFT.filter(regex='^(?!test_.*$|bts_.*$|train_.*$).*'); metaDF.columns + + print('\nTotal cols in each df:' + , '\nCV df:', len(scoresDF_CV.columns) + , '\nmetaDF:', len(metaDF.columns)) + + #------------------------------------- + # Combine WF: CV + Metadata + #------------------------------------- + + combDF = pd.merge(scoresDF_CV, metaDF, left_index = True, right_index = True) + print('\nAdding column: Model_name') + combDF['Model_name'] = combDF.index + + #------------------------------------- + # Combine WF+Metadata: Final output + #------------------------------------- + + # if len(combDF.columns) == expected_ncols_out: + # print('\nPASS: Combined df has expected ncols') + # else: + # sys.exit('\nFAIL: Length mismatch for combined_df') + + # print('\nAdding column: Model_name') + # combDF['Model_name'] = combDF.index + + print('\n=========================================================' + , '\nSUCCESS: Ran multiple classifiers' + , '\n=======================================================') + + #resampling_methods_wf = combined_baseline_wf[['resampling']] + #resampling_methods_wf = resampling_methods_wf.drop_duplicates() + #, '\n', resampling_methods_wf) + + return combDF + +############################################################################### diff --git a/scripts/ml/ml_iterator_CVs.py b/scripts/ml/ml_iterator_CVs.py new file mode 100644 index 0000000..233f6fb --- /dev/null +++ b/scripts/ml/ml_iterator_CVs.py @@ -0,0 +1,132 @@ +#!/usr/bin/env python3 +# -*- coding: utf-8 -*- +""" +Created on Wed Jun 29 20:29:36 2022 + +@author: tanu +""" +import sys, os +import pandas as pd +import numpy as np +import re + +############################################################################### +homedir = os.path.expanduser("~") +sys.path.append(homedir + '/git/LSHTM_analysis/scripts/ml/ml_functions') +sys.path +############################################################################### +outdir = homedir + '/git/LSHTM_ML/output/genes/' + +#==================== +# Import ML functions +#==================== +from MultClfs_CVs import * +from GetMLData import * +from SplitTTS import * + +# skf_cv = StratifiedKFold(n_splits = 10 +# #, shuffle = False, random_state= None) +# , shuffle = True, random_state = 42) + +# #rskf_cv = RepeatedStratifiedKFold(n_splits = 10 +# # , n_repeats = 3 +# # , **rs) + +# param dict for getmldata() +gene_model_paramD = {'data_combined_model' : False + , 'use_or' : False + , 'omit_all_genomic_features': False + , 'write_maskfile' : False + , 'write_outfile' : False } +############################################################################### +#ml_genes = ["pncA", "embB", "katG", "rpoB", "gid"] + +ml_gene_drugD = { + 'pncA' : 'pyrazinamide' + #, 'embB' : 'ethambutol' + #, 'katG' : 'isoniazid' + #, 'rpoB' : 'rifampicin' + #, 'gid' : 'streptomycin' + } +gene_dataD={} + +split_types = ['none'] +split_data_types = ['complete'] + +for gene, drug in ml_gene_drugD.items(): + print ('\nGene:', gene + , '\nDrug:', drug) + gene_low = gene.lower() + gene_dataD[gene_low] = getmldata(gene, drug + , **gene_model_paramD) + + for split_type in split_types: + for data_type in split_data_types: + + out_filename = outdir + gene.lower() + '_' + split_type + '_' + data_type + '.csv' + + tempD = split_tts(gene_dataD[gene_low] + , data_type = data_type + , split_type = split_type + , oversampling = True + , dst_colname = 'dst' + , target_colname = 'dst_mode' + , include_gene_name = True + ) + paramD = { + 'baseline_paramD': { 'input_df' : tempD['X'] + , 'target' : tempD['y'] + , 'var_type' : 'mixed' + , 'resampling_type': 'none'} + + , 'smnc_paramD' : { 'input_df' : tempD['X_smnc'] + , 'target' : tempD['y_smnc'] + , 'var_type' : 'mixed' + , 'resampling_type' : 'smnc'} + + , 'ros_paramD' : { 'input_df' : tempD['X_ros'] + , 'target' : tempD['y_ros'] + , 'var_type' : 'mixed' + , 'resampling_type' : 'ros'} + + , 'rus_paramD' : { 'input_df' : tempD['X_rus'] + , 'target' : tempD['y_rus'] + , 'var_type' : 'mixed' + , 'resampling_type' : 'rus'} + + , 'rouC_paramD' : { 'input_df' : tempD['X_rouC'] + , 'target' : tempD['y_rouC'] + , 'var_type' : 'mixed' + , 'resampling_type' : 'rouC'} + } + + mmDD = {} + for k, v in paramD.items(): + print(k) + all_scoresDF = pd.DataFrame() + + # iterate over different cv thresholds + for skf_cv_threshold in [3,5,10]: + print('\nRunning CV threhhold:', skf_cv_threshold) + current_scoreDF = MultModelsCl_CVs(**paramD[k] + , skf_cv_threshold = skf_cv_threshold + , tts_split_type = split_type + , add_cm = True + , add_yn = True + , scale_numeric = ['min_max'] + , random_state = 42 + , n_jobs = os.cpu_count() + , return_formatted_output = True + ) + + all_scoresDF = pd.concat([all_scoresDF, current_scoreDF]) + mmDD[k] = all_scoresDF + + # Extracting the dfs from within the dict and concatenating to output as one df + for k, v in mmDD.items(): + out_wf = pd.concat(mmDD, ignore_index = True) + + #out_wf_f = out_wf.sort_values(by = ['resampling', 'source_data', 'MCC'], ascending = [True, True, False], inplace = False) + #out_wf_f.to_csv(out_filename, index = False) + out_wf.to_csv(out_filename, index = False) +