from sklearn.metrics import roc_curve, auc,confusion_matrix,f1_score, precision_score, recall_score,accuracy_score
from sklearn.linear_model import LogisticRegression as LR
from sklearn.model_selection import train_test_split
from imblearn.over_sampling import SMOTE
from keras.layers.core import Dense
from keras.models import Sequential
from xgboost import XGBClassifier
import matplotlib.pyplot as plt
from keras import initializers
from keras.layers import LSTM
import tensorflow as tf
import seaborn as sns
import pandas as pd
import numpy as np
import keras as ks
import random
import time
import os

import warnings
warnings.filterwarnings('ignore')

start_time = time.time()
rd = 369

###Fixed Neural Network Random Seeds
def seed_tensorflow(seed):
    os.environ['PYTHONHASHSEED'] = str(seed)
    random.seed(seed)
    np.random.seed(seed)
    tf.random.set_seed(seed)
    os.environ['TF_DETERMINISTIC_OPS'] = '1'
seed_tensorflow(rd)

dim = [0,1]


plt.rcParams['font.sans-serif']=['SimHei']
plt.rcParams['axes.unicode_minus']=False
path = 'C:/users/youcai/Data'
data = pd.read_csv(path + 'default of credit card clients.csv',index_col=0)

###Plotting the AUC curve
def make_figure_auc_roc(l1,l2,l3,h_ls,auc_avg,name='All model'):
    plt.figure(figsize=(15, 9))
    for true_label,score,modelname,h in zip(l1,l2,l3,h_ls):
        fpr, tpr, thresholds = roc_curve(true_label, score, pos_label=0, drop_intermediate=False)

        roc_auc = auc(fpr, tpr)
        roc_result = pd.DataFrame({'fpr' : fpr,'tpr' : tpr, 'tf' : tpr - (1-fpr),'thresholds' : thresholds})
        optimal_idx = np.argmax(tpr - fpr)
        optimal_threshold = thresholds[optimal_idx]
        y = tpr[optimal_idx]
        x = fpr[optimal_idx]
        plt.plot(fpr, tpr,label='%s AUC=%.2f'%(modelname,roc_auc))
        plt.plot([-0.005, 1.005], [-0.005, 1.005], color='orange',linestyle='--')
        plt.xlabel('1-specificity',fontsize=15)
        plt.ylabel('sensitive',fontsize=15)
        plt.scatter(x, y, c='red',s=15)
        plt.tick_params(labelsize=12)
        plt.title('%s avg auc：%s'%(name,auc_avg),fontsize=18)
        plt.legend(loc='best',fontsize=15)
        #plt.text(x+0.05,y-h,s='('+str(round(1-x, 2))+', '+str(round(y, 2))+')',fontsize=10)
    plt.show()

###Histogram of y-number distribution
def bar_count(modeldata,msg):

    print('\033[1;34m{0:*^100}\033[0m'.format('Bar charting'))

    plt.figure(figsize=(15, 9))
    v = pd.Series(modeldata[msg]).value_counts()
    vv = pd.Series(modeldata[msg]).value_counts() / modeldata.shape[0]
    ax = sns.countplot(x=msg,data=modeldata,order=v.index.to_list())
    for x, y in enumerate(v):
        text = ax.text(x, y, round(list(vv)[x], 3))
        text.set_ha('center')
    plt.tick_params(labelsize=12)
    plt.show()


def radar_fig(modeldata):

    print(modeldata)
    labels = modeldata.index.to_list()
    data = list(modeldata)
    dataLenth = len(labels)
    angles = np.linspace(0, 2 * np.pi, dataLenth, endpoint=False)


    data = np.concatenate((data, [data[0]]))
    angles = np.concatenate((angles, [angles[0]]))
    labels = np.concatenate((labels, [labels[0]]))

    plt.figure(figsize=(15,9),facecolor='white')
    plt.subplot(111, polar=True)
    plt.plot(angles, data, 'bo-', color='b', linewidth=2)
    plt.fill(angles, data, facecolor='b', alpha=0.25)
    plt.thetagrids(angles * 180 / np.pi, labels)
    plt.figtext(0.52,0.95,'XGB model feature importance',ha='center',fontsize=15)
    plt.grid(True)
    plt.tick_params(labelsize=12)
    plt.show()

##Calculation of auc
def auc_roc(true_label,score,class_,):

    fpr, tpr, thresholds = roc_curve(true_label, score, pos_label=class_, drop_intermediate=False)
    roc_auc = auc(fpr, tpr)
    return roc_auc

###xgboost model
def my_xgb(xtrain,ytrain,xtest,ytest,gamma=0.1,max_depth=3,subsample=0.5,scale_pos_weight=1,colsample_bytree=0.5,lr=0.01):

    xgb = XGBClassifier(gamma=gamma, max_depth=max_depth, subsample=subsample,scale_pos_weight=scale_pos_weight,colsample_bytree=colsample_bytree,learning_rate=lr,random_state=rd)  # 建立Xgboost并忽略警告
    xgb.fit(xtrain, ytrain, early_stopping_rounds=10, eval_metric='auc', eval_set=[(xtest, ytest)], verbose=False)
    pre_xgb = xgb.predict(xtest)
    pro_xgb = xgb.predict_proba(xtest)
    pret_xgb = xgb.predict(xtrain)
    prot_xgb = xgb.predict_proba(xtrain)

    featureimportance = xgb.feature_importances_
    # featureimportance = featureimportance / featureimportance.max()###
    plt.figure(figsize=(15, 9))
    df = pd.DataFrame(featureimportance,index=list(xtrain),columns=['Importance of XGB'])

    return pre_xgb,ytest,pro_xgb,xgb.classes_[0],pret_xgb,ytrain,prot_xgb,df

#Establishment of logistic regression (LR)
def my_lr(xtrain,ytrain,xtest,ytest,C=10,max_iter=100,solver='lbfgs'):#

    lr = LR(C=C,max_iter=max_iter,solver=solver,random_state=rd)
    lr.fit(xtrain, ytrain)
    pre_lr = lr.predict(xtest)
    pro_lr = lr.predict_proba(xtest)
    pret_lr = lr.predict(xtrain)
    prot_lr = lr.predict_proba(xtrain)
    return pre_lr, ytest, pro_lr, lr.classes_[0],pret_lr,ytrain,prot_lr

###y labels one-hotised
def tran_onehot_vector(label_train, dim_name):
    target_train = np.zeros([len(label_train), len(dim_name)])
    for i in range(len(label_train)):
        target_train[i, label_train[i]] = 1
    return target_train

###DNN network model
def my_bpnn(xtrain,ytrain,xtest,ytest,ininlr=0.01):

    model = Sequential()
    BpNeuralNet_temp = tran_onehot_vector(list(ytrain), dim)
    layer1 = Dense(128, input_dim=xtrain.shape[1], activation='relu',kernel_initializer=initializers.random_normal(stddev=0.01), bias_initializer='zeros')
    layer2 = Dense(64, activation='relu', kernel_initializer=initializers.random_normal(stddev=0.01))
    layer3 = Dense(32, activation='relu', kernel_initializer=initializers.random_normal(stddev=0.01))
    layer4 = Dense(BpNeuralNet_temp.shape[1], activation='sigmoid',kernel_initializer=initializers.random_normal(stddev=0.01))

    model.add(layer1)
    model.add(ks.layers.Dropout(0.2))
    model.add(layer2)
    model.add(ks.layers.Dropout(0.2))
    model.add(layer3)
    model.add(ks.layers.Dropout(0.2))
    model.add(layer4)
    lr = ininlr
    adam = tf.keras.optimizers.Adam(learning_rate=lr)
    model.compile(loss='binary_crossentropy', optimizer=adam, metrics=['accuracy'])

    # Activation = ['tanh','relu','sigmoid','softmax']
    # optimizer = ['sgd', 'rmsprop', 'adam', 'adagrad']
    # loss = ['binary_crossentropy','mean_squared_error','categorical_crossentropy','mae']
    batchsize = 300
    epochs = 100
    model.fit(xtrain, BpNeuralNet_temp, epochs=epochs, batch_size=batchsize,verbose=0)
    pro_bpnn = model.predict(xtest)
    pre_bpnn = np.argmax(pro_bpnn, axis=1)
    prot_bpnn = model.predict(xtrain)
    pret_bpnn = np.argmax(prot_bpnn, axis=1)
    return pre_bpnn, ytest, pro_bpnn,0,pret_bpnn,ytrain,prot_bpnn

###LSTM network model
def my_lstm(xtrain,ytrain,xtest,ytest,ininlr=0.01):


    xtrain = np.array(xtrain)
    ytrain = list(ytrain)
    xtest = np.array(xtest)
    ytest = list(ytest)
    model = Sequential()
    nb_time_steps = 1

    def data_format(X,y):

        X_seq = []
        y_seq = []
        for i in range(nb_time_steps, len(X)):
            X_seq.append(X[i - nb_time_steps+1:i+1])
            y_seq.append(y[i])
        X_seq = np.array(X_seq)
        y_t = np.array(y_seq)
        y_seq = tran_onehot_vector(list(y_t), dim)
        return X_seq,y_seq,y_t

    xtrain = np.array(xtrain)
    xtest = np.array(xtest)
    dataset_train = data_format(xtrain,ytrain)
    dataset_test = data_format(xtest,ytest)

    nb_lstm_outputs1 = 128
    nb_lstm_outputs2 = 64
    nb_input_vector = xtrain.shape[1]
    nb_output_vector = dataset_train[1].shape[1]

    layer1 = LSTM(units=nb_lstm_outputs1, input_shape=(nb_time_steps, nb_input_vector), return_sequences=True)
    layer2 = LSTM(units=nb_lstm_outputs2, input_shape=(nb_time_steps, nb_input_vector))
    layer3 = Dense(32, activation='relu')
    layer4 = Dense(16, activation='relu')
    layer5 = Dense(nb_output_vector, activation='sigmoid')

    model.add(layer1)
    model.add(ks.layers.Dropout(0.2))
    model.add(layer2)
    model.add(ks.layers.Dropout(0.2))
    model.add(layer3)
    model.add(ks.layers.Dropout(0.2))
    model.add(layer4)
    model.add(ks.layers.Dropout(0.2))
    model.add(layer5)
    lr = ininlr
    adam = tf.keras.optimizers.Adam(learning_rate=lr)
    model.compile(loss='binary_crossentropy', optimizer=adam, metrics=['accuracy'])

    batchsize = 300
    epochs = 100
    model.fit(dataset_train[0], dataset_train[1], epochs=epochs, batch_size=batchsize, verbose=0)
    pro_lstm = model.predict(dataset_test[0])
    pre_lstm = np.argmax(pro_lstm, axis=1)
    prot_lstm = model.predict(dataset_train[0])
    pret_lstm = np.argmax(prot_lstm, axis=1)
    return pre_lstm, dataset_test[2], pro_lstm,0,pret_lstm,dataset_train[2],prot_lstm

##Randomised cross-validation function
def rand_cross_vaild(rdls,x,y,func,name):

    print('\033[1;38m{0:*^40}\033[0m'.format('Summary of stochastic cross-validation results for the %s model'%name))

    acc_ls = []
    f1_ls = []
    p_ls = []
    r_ls = []
    auctrain_ls = []
    auctest_ls = []

    auc_train1 = []
    auc_train2 = []
    auc_train3 = []
    auc_train4 = []
    auc_test1 = []
    auc_test2 = []
    auc_test3 = []
    auc_test4 = []

    DD = pd.DataFrame()
    for k in rdls:

        xtrain, xtest, ytrain, ytest = train_test_split(x, y, test_size=0.2, random_state=k)
        smote = SMOTE(random_state=rd,sampling_strategy=1)
        xtrain, ytrain = smote.fit_resample(xtrain, ytrain)

        res = func(xtrain,ytrain,xtest,ytest)
        acc = accuracy_score(res[1], res[0])
        f1 = f1_score(res[1], res[0], average='binary')
        p = precision_score(res[1], res[0], average='binary')
        r = recall_score(res[1], res[0], average='binary')
        auctest = auc_roc(res[1], res[2][:, 0], res[3])
        auctrain = auc_roc(res[5], res[6][:, 0], res[3])
        f1_ls.append(f1)
        p_ls.append(p)
        r_ls.append(r)
        acc_ls.append(acc)
        auctrain_ls.append(auctrain)
        auctest_ls.append(auctest)


        auc_train1.append(res[5])
        auc_train2.append(res[6][:, 0])
        auc_train3.append('random seed:%s'%k)
        auc_train4.append(0.05)

        auc_test1.append(res[1])
        auc_test2.append(res[2][:, 0])
        auc_test3.append('random seed:%s'%k)
        auc_test4.append(0.05)

        if name == 'XGB':
            dd = res[7]
            DD = pd.concat([DD,dd],axis=1)
            DD = DD.mean(axis=1)

    if name == 'XGB':
        radar_fig(DD)

    auctest_avg = pd.Series(auctest_ls).mean()
    auctrain_avg = pd.Series(auctrain_ls).mean()
    make_figure_auc_roc(auc_train1,auc_train2,auc_train3,auc_train4,round(auctrain_avg,2),name='train'+' '+name)
    make_figure_auc_roc(auc_test1,auc_test2,auc_test3,auc_test4,round(auctest_avg,2),name='test'+' '+name)

    print('\033[1;36mAverage test set accuracy：%.2f\033[0m'%(pd.Series(acc_ls).mean()))
    print('\033[1;36mAverage test set f1 value：%.2f\033[0m' % (pd.Series(f1_ls).mean()))
    print('\033[1;36mAverage test set recall value：%.2f\033[0m' % (pd.Series(r_ls).mean()))
    print('\033[1;36mAverage test set precision rate value：%.2f\033[0m' % (pd.Series(p_ls).mean()))
    print('\033[1;36mAverage test set AUC value：%.2f\033[0m' % (auctest_avg))
    return acc_ls,f1_ls,p_ls,r_ls,auctest_ls



###Read the x and y of the data
x = data.iloc[:,:-1]
# x = (x - x.mean()) / x.std()##
x = (x - x.min()) / (x.max() - x.min())
y = data.iloc[:,-1]


rand_ls = np.random.randint(0,100,10)
reslr = rand_cross_vaild(rand_ls,x,y,my_lr,'LR')
resxgb = rand_cross_vaild(rand_ls,x,y,my_xgb,'XGB')



resbpnn = rand_cross_vaild(rand_ls,x,y,my_bpnn,'BPNN')
reslstm = rand_cross_vaild(rand_ls,x,y,my_lstm,'LSTM')

#
# ###The following plots the discounted distribution of each technical indicator on the different models for stochastic cross-validation
plt.figure(figsize=(15, 9))
plt.plot(reslr[0],label='LR，std=%s'%(round(np.array(reslr[0]).std(),5)), marker='o')
plt.plot(resxgb[0],label='XGB，std=%s'%(round(np.array(resxgb[0]).std(),5)), marker='*')
plt.plot(resbpnn[0],label='BPNN，std=%s'%(round(np.array(resbpnn[0]).std(),5)), marker='+')
plt.plot(reslstm[0],label='LSTM，std=%s'%(round(np.array(reslstm[0]).std(),5)), marker='^')
plt.legend(loc='best', fontsize=15)
plt.tick_params(labelsize=15)
plt.xticks([j for j in range(len(rand_ls))], [str(j) for j in rand_ls])
plt.xlabel('random %s times' % len(rand_ls), fontsize=18)
plt.ylabel('accuracy', fontsize=18)
plt.show()

plt.figure(figsize=(15, 9))
plt.plot(reslr[1],label='LR，std=%s'%(round(np.array(reslr[1]).std(),5)), marker='o')
plt.plot(resxgb[1],label='XGB，std=%s'%(round(np.array(resxgb[1]).std(),5)), marker='*')
plt.plot(resbpnn[1],label='BPNN，std=%s'%(round(np.array(resbpnn[1]).std(),5)), marker='+')
plt.plot(reslstm[1],label='LSTM，std=%s'%(round(np.array(reslstm[1]).std(),5)), marker='^')
plt.legend(loc='best', fontsize=15)
plt.tick_params(labelsize=15)
plt.xticks([j for j in range(len(rand_ls))], [str(j) for j in rand_ls])
plt.xlabel('random %s times' % len(rand_ls), fontsize=18)
plt.ylabel('F1-score', fontsize=18)
plt.show()

plt.figure(figsize=(15, 9))
plt.plot(reslr[4],label='LR，std=%s'%(round(np.array(reslr[4]).std(),5)), marker='o')
plt.plot(resxgb[4],label='XGB，std=%s'%(round(np.array(resxgb[4]).std(),5)), marker='*')
plt.plot(resbpnn[4],label='BPNN，std=%s'%(round(np.array(resbpnn[4]).std(),5)), marker='+')
plt.plot(reslstm[4],label='LSTM，std=%s'%(round(np.array(reslstm[4]).std(),5)), marker='^')
plt.legend(loc='best', fontsize=15)
plt.tick_params(labelsize=15)
plt.xticks([j for j in range(len(rand_ls))], [str(j) for j in rand_ls])
plt.xlabel('random %s times' % len(rand_ls), fontsize=18)
plt.ylabel('AUC', fontsize=18)
plt.show()


stop_time = time.time()
print('\033[1;34mFinal time at end of complete run (seconds)：%s\033[0m' % (round(stop_time - start_time, 3)))