<jupyter_start><jupyter_text>Note to self, remove a hold out set of data to play around with at the end<jupyter_code>import pandas as pd
import numpy as np
import seaborn as sns
import matplotlib.pyplot as plt
from textwrap import wrap
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LogisticRegression 
from sklearn import metrics
from sklearn.linear_model import LassoCV
from sklearn.model_selection import KFold
from sklearn.model_selection import cross_val_score
from sklearn.preprocessing import binarize

#Read pre processed data into pandas df and drop the index column that comes in from the spreadsheet
df=pd.read_excel('Data/HorseColicData_DeathTarget.xlsx',index=False)
df.shape
df=df.drop('Unnamed: 0',axis=1)
df
SimpleBinary=df
SimpleBinary['Target'].value_counts()
#Final check to ensure no nan values
SimpleBinary.isna().sum()
#Create features for model consumption
X = SimpleBinary.drop(columns= ['outcome','Target'])
X.columns
#Create target variable for model consumption
y = SimpleBinary.Target 
y.value_counts()<jupyter_output><empty_output><jupyter_text>Split the data into train and test, with a test size of 25% of the total dataset. randomstate 40 is used for reproducability of exercise. instantiate the ligistic regression model and fit the training data to it.<jupyter_code>X_train, X_test, Y_train, Y_test = train_test_split(X, y, test_size = .25, random_state = 40)

logreg = LogisticRegression(penalty = 'none', max_iter= 10000,random_state = 40) ## instantiate

logreg.fit(X_train, Y_train) ## fit
#Check the accuracy scores of the classification model using the test data
score = logreg.score(X_test,Y_test)
score
predictions = logreg.predict(X_test)
probabilities= logreg.predict_proba(X_test)
problive=[]
probdie=[]
for i in range(len(probabilities)):
    problive.append(probabilities[i][0])
    probdie.append(probabilities[i][1])
ModelChart = {'live':problive,'die':probdie,'predictedoutcome':list(predictions)}
Chartdf = pd.DataFrame(ModelChart)
Chartdf['predictedoutcome'].replace(to_replace =0, value ="Lives",inplace=True)
Chartdf['predictedoutcome'].replace(to_replace =1, value ="Dies",inplace=True)
Chartdf
#Some styling assistance from matplotlib example docs
#https://matplotlib.org/gallery/lines_bars_and_markers/barchart.html#sphx-glr-gallery-lines-bars-and-markers-barchart-py

sns.set(rc={'axes.facecolor':'#333333', 'figure.facecolor':'#333333'})
plt.rcParams['font.family'] = 'sans-serif'
plt.rcParams['font.sans-serif'] = 'montserrat'
plt.rcParams['text.color'] = 'white'
plt.rcParams['axes.labelcolor']= 'white'
plt.rcParams['xtick.color'] = 'white'
plt.rcParams['ytick.color'] = 'white'

plt.rcParams['patch.edgecolor'] = '#333333'

labels = Chartdf[0:10]['predictedoutcome'].astype(str)
labels = [ '\n'.join(wrap(l, 12)) for l in labels ]

A_means = Chartdf[0:10]['live'].round(3)
B_means = Chartdf[0:10]['die'].round(3)

x = np.arange(len(labels))  # the label locations
width = 0.35  # the width of the bars

fig, ax = plt.subplots(figsize=(16, 8))

rects1 = ax.bar(x - width/2, A_means, width, label='Liklihood of Lives',color=['#EAA700'])
rects2 = ax.bar(x + width/2, B_means, width, label='Liklihood of Death',color=['#9E1B34'])

# Add some text for labels, title and custom x-axis tick labels, etc.
ax.set_ylabel('Log Liklihood of Possible Outcome for Each Prediction')
ax.set_xlabel('Predicted Outcome for Each Observation')
ax.set_title('Logistic Regression Classifier Visualization of Maximum Liklihood',size=14)
ax.set_xticks(x)
ax.set_xticklabels(labels)
ax.legend()


def autolabel(rects):
    """Attach a text label above each bar in *rects*, displaying its height."""
    for rect in rects:
        height = rect.get_height()
        ax.annotate('{}'.format(height),
                    xy=(rect.get_x() + rect.get_width() / 2, height),
                    xytext=(0, 3),  # 3 points vertical offset
                    textcoords="offset points",
                    ha='center', va='bottom')

ax.grid(False)

autolabel(rects1)
autolabel(rects2)

fig.tight_layout()

plt.show()
#Generate a confusion matrix for the test data and predictions for the above model. code from class examples
cm = metrics.confusion_matrix(Y_test, predictions)<jupyter_output><empty_output><jupyter_text>In the below confusion matrix we see an overall prediction accuracy of approximately 74%. False posative outcomes outweigh false negative outcomes at 15 vs 9.<jupyter_code>plt.figure(figsize=(9,9))
sns.heatmap(cm, annot=True, fmt=".3f", linewidths=.5, square = True, cmap = 'RdBu');
plt.ylabel('Actual label');
plt.xlabel('Predicted label');
all_sample_title = 'Accuracy Score: {0}'.format(score)
plt.title(all_sample_title, size = 15);
plt.savefig('toy_Digits_ConfusionSeabornCodementor.png')
plt.show()<jupyter_output>findfont: Font family ['sans-serif'] not found. Falling back to DejaVu Sans.<jupyter_text>We will down attempt to do some feature selection, removing potential noise form uncorrelated variables. In this first step those features with a correlation of less than .1 to the target feature are removed from the dataframe.<jupyter_code>#remove values with correlation less than .1
Corrs=SimpleBinary[SimpleBinary.columns].corr()['Target']
remove2=Corrs[np.sqrt(Corrs**2) < 0.1].index
remove2
SimpleBinary2=SimpleBinary.drop(list(remove2), axis=1)
SimpleBinary2<jupyter_output><empty_output><jupyter_text>We will now refit the model with these variables removed and reassess its performance.<jupyter_code>X = SimpleBinary2.drop(columns= ['outcome','Target'])
#X = pd.get_dummies(X)
X.columns
y = SimpleBinary2.Target ## creating target values
y.value_counts()

X_train, X_test, Y_train, Y_test = train_test_split(X, y, test_size = .25, random_state = 40)

logreg = LogisticRegression(penalty = 'none', max_iter= 10000,random_state = 40) ## instantiate

logreg.fit(X_train, Y_train) ## fit
## Accuracy Scores
trainscore = logreg.score(X_train,Y_train)
testscore = logreg.score(X_test,Y_test)

print(trainscore)
print(testscore)
predictions = logreg.predict(X_test)
predictions
cm = metrics.confusion_matrix(Y_test, predictions)<jupyter_output><empty_output><jupyter_text>We actually observe a slightly reduced test accuracy based on the removal of the above variables. One additional prediction becomes a false posative.<jupyter_code>plt.figure(figsize=(9,9))
sns.heatmap(cm, annot=True, fmt=".3f", linewidths=.5, square = True, cmap = 'RdBu');
plt.ylabel('Actual label');
plt.xlabel('Predicted label');
all_sample_title = 'Accuracy Score: {0}'.format(testscore)
plt.title(all_sample_title, size = 15);
plt.savefig('toy_Digits_ConfusionSeabornCodementor.png')
plt.show()<jupyter_output><empty_output><jupyter_text>Now we will attempt feature selection using the lasso meathod, as opposed to simply discarding the variables with low simple correlation to target.<jupyter_code>X = SimpleBinary.drop(columns= ['outcome','Target'])
X.columns
y = SimpleBinary.Target ## creating target values
y.value_counts()
X_train, X_test, Y_train, Y_test = train_test_split(X, y, test_size = .25, random_state = 40)
#Assistance from https://github.com/justmarkham/scikit-learn-videos
lasso = LassoCV().fit(X_train, Y_train)
importance = np.abs(lasso.coef_)
#feature_names = np.array(X_train.values)
plt.bar(height=importance, x=X_train.columns)
plt.xticks(rotation=90)
plt.title("Feature importances via coefficients")
plt.show()
importance
LassoOutcome = {'feature':X_train.columns,'importance':list(importance)}
LassoOutcome = pd.DataFrame(LassoOutcome)
LassoOutcome
#remove values with importance of 0
Lassoremove=LassoOutcome[LassoOutcome['importance'] == 0]
list(Lassoremove['feature'])
SimpleBinary3=SimpleBinary.drop(list(Lassoremove['feature']), axis=1)
SimpleBinary3<jupyter_output><empty_output><jupyter_text>Now lets retrain the classifier model using these features.<jupyter_code>X = SimpleBinary3.drop(columns= ['outcome','Target'])
#X = pd.get_dummies(X)
X.columns
y = SimpleBinary3.Target ## creating target values
y.value_counts()
X_train, X_test, Y_train, Y_test = train_test_split(X, y, test_size = .25, random_state = 40)

logreg = LogisticRegression(penalty = 'none', max_iter= 10000,random_state = 40) ## instantiate

logreg.fit(X_train, Y_train) ## fit
## Accuracy Scores
trainscore = logreg.score(X_train,Y_train)
testscore = logreg.score(X_test,Y_test)

print(trainscore)
print(testscore)
predictions = logreg.predict(X_test)
predictions
cm = metrics.confusion_matrix(Y_test, predictions)<jupyter_output><empty_output><jupyter_text>We now observe a test accuracy equal to that of the baseline with the removal of the above variables identified by the lasso regression. While no improvement was noted, excluding uneeded features makes for a better model and we will now use this dataset moving forward for parameter testing and final evaluation.<jupyter_code>plt.figure(figsize=(9,9))
sns.heatmap(cm, annot=True, fmt=".3f", linewidths=.5, square = True, cmap = 'RdBu');
plt.ylabel('Actual label');
plt.xlabel('Predicted label');
all_sample_title = 'Accuracy Score: {0}'.format(testscore)
plt.title(all_sample_title, size = 15);
plt.savefig('ConfusionMatrix_Lasso.png')
plt.show()<jupyter_output><empty_output><jupyter_text>Validating the model via Cross Validation, k folds used as dataset isnt huge orig run time * k times. this is a more accurate measurment of model accracySource used for assistance in coding:https://github.com/justmarkham/scikit-learn-videos/blob/master/07_cross_validation.ipynband https://github.com/justmarkham/scikit-learn-videos/blob/master/09_classification_metrics.ipynb<jupyter_code>#Below is a visualization of how the folds used in this method will be evaluated, splitting out the test and train data for each fold
kf = KFold(n_splits=5, shuffle=False).split(range(len(SimpleBinary3)))

# print the contents of each training and testing set
print('{} {:^61} {}'.format('Iteration', 'Training set observations', 'Testing set observations'))
for iteration, data in enumerate(kf, start=1):
    print('{:^9} {} {:^25}'.format(iteration, data[0], str(data[1])))
# 10-fold cross-validation 
scores = cross_val_score(logreg, X, y, cv=10, scoring='accuracy')
print(scores)
#The scores for these ten folds are now averaged to come up with the final score for the model
print(scores.mean())<jupyter_output>0.7715465465465465<jupyter_text>Now that we have a more accurate measure of the accuracy of our model, we can use the follwing functions to fine tune the parameters of the logistic regression classifier to obtain the best possible results from our model Cfloat, default=1.0Inverse of regularization strength; must be a positive float. Like in support vector machines, smaller values specify stronger regularization.<jupyter_code>Cvals=[0.001,0.01,0.1,1,10,100]
k_scores = []
for c in Cvals:
    logreg = LogisticRegression(max_iter= 10000,random_state = 40,C=c)

    scores = cross_val_score(logreg, X, y, cv=10, scoring='accuracy')
    k_scores.append(scores.mean())
print(k_scores)
#Select the C value with the best outcome
BestC=Cvals[k_scores.index(max(k_scores))]<jupyter_output><empty_output><jupyter_text>solver{‘newton-cg’, ‘lbfgs’, ‘liblinear’, ‘sag’, ‘saga’}, default=’lbfgs’Algorithm to use in the optimization problem.<jupyter_code>solve = ['newton-cg', 'lbfgs', 'liblinear', 'sag', 'saga']
k_scores = []
for s in solve:
    logreg = LogisticRegression(max_iter= 1000000,random_state = 40,solver=s,C=BestC) ## instantiate

    scores = cross_val_score(logreg, X, y, cv=10, scoring='accuracy')
    k_scores.append(scores.mean())
print(k_scores)
#Select the solver with the best outcome
BestS=solve[k_scores.index(max(k_scores))]<jupyter_output><empty_output><jupyter_text>We see an improved outcome from the original default parameters set by the model, from 74 to 77<jupyter_code># 10-fold cross-validation using the parameters which the above functions produced
logreg = LogisticRegression(max_iter= 1000000,random_state = 40,C=BestC,solver=BestS) ## instantiate
logreg.fit(X_train, Y_train) ## fit
scores = cross_val_score(logreg, X, y, cv=10, scoring='accuracy')
print(scores.mean())<jupyter_output>0.7716216216216216<jupyter_text>One final way to measure the true effectiveness of our model is to consider the distribution of the possible outcomes (life or death). Does out model predict simply based on the numerical split/distribution within the data, or is it more accurate than that?<jupyter_code>#Checking null accuracy
print(Y_test.value_counts().head(2) / len(Y_test))
y.value_counts().head(2) / len(y)<jupyter_output><empty_output><jupyter_text>A prediction of both outcomes does not simply follow the distribtuon percentage of the targets within the dataset, a posative outcome. The final manual parameter which we could adjust to achieve greater accuracy would be to alter the final prediction output of the model based on the liklihood output. In other words, a liklihood greater than the defualt .5 can be specified when making the final prediction off of the proability array. This is the final parameter which should be considered and only once the final model and features have been selected. It is often a matter of domain knowledge. For us we will see if weighting a "death" outcome can allow us to reduce the number of false negative results. erring on the side of caution for the model used by veteranarians so that they be more lean towards more aggressive medical intervention when death could be a result of the evaluation of the horse.<jupyter_code># make class predictions for the testing set
y_pred_class = logreg.predict(X_test)
print('True:', Y_test.values[0:25])
print('Pred:', y_pred_class[0:25])
#Checking Distribution
y_pred_prob = logreg.predict_proba(X_test)[:, 1]
plt.hist(y_pred_prob, bins=10)
plt.xlim(0, 1)
plt.title('Distribution of predicted probabilities')
plt.xlabel('Predicted probability of death')
plt.ylabel('Frequency')
PossibleWeights=np.arange(.4, .6, 0.01).tolist()
BestAccuracy=[]
for i in range(len(PossibleWeights)):
    BestAccuracy.append(metrics.accuracy_score(Y_test, binarize([y_pred_prob], PossibleWeights[i])[0]))
Bestweight=PossibleWeights[BestAccuracy.index(max(BestAccuracy))]
Bestweight
lowfalseneg=[]
for i in range(len(PossibleWeights)):
    lowfalseneg.append(metrics.confusion_matrix(Y_test,binarize([y_pred_prob], PossibleWeights[i])[0])[1][0])
LowestFalseNegative=PossibleWeights[lowfalseneg.index(min(lowfalseneg))]
LowestFalseNegative
print('Lowest false negative count possible: ' + str(min(lowfalseneg)))
y_pred_classAccuracy=binarize([y_pred_prob], Bestweight)
plt.figure(figsize=(9,9))
sns.heatmap((metrics.confusion_matrix(Y_test, y_pred_classAccuracy[0])), annot=True, fmt=".3f", linewidths=.5, square = True);
plt.ylabel('Actual label');
plt.xlabel('Predicted label');
all_sample_title = 'Accuracy Score: {0}'.format((metrics.accuracy_score(Y_test, y_pred_classAccuracy[0])))
plt.title(all_sample_title, size = 15);
plt.savefig('toy_Digits_ConfusionSeabornCodementor.png')
plt.show()
def evaluate_threshold(threshold):
    print('Sensitivity:', tpr[thresholds > threshold][-1])
    print('Specificity:', 1 - fpr[thresholds > threshold][-1])
evaluate_threshold(Bestweight)

y_pred_classlowfalseneg=binarize([y_pred_prob], LowestFalseNegative)
plt.figure(figsize=(9,9))
sns.heatmap((metrics.confusion_matrix(Y_test, y_pred_classlowfalseneg[0])), annot=True, fmt=".3f", linewidths=.5, square = True);
plt.ylabel('Actual label');
plt.xlabel('Predicted label');
all_sample_title = 'Accuracy Score: {0}'.format((metrics.accuracy_score(Y_test, y_pred_classlowfalseneg[0])))
plt.title(all_sample_title, size = 15);
plt.savefig('toy_Digits_ConfusionSeabornCodementor.png')
plt.show()
def evaluate_threshold(threshold):
    print('Sensitivity:', tpr[thresholds > threshold][-1])
    print('Specificity:', 1 - fpr[thresholds > threshold][-1])
evaluate_threshold(LowestFalseNegative)

def evaluate_threshold(threshold):
    print('Sensitivity:', tpr[thresholds > threshold][-1])
    print('Specificity:', 1 - fpr[thresholds > threshold][-1])
evaluate_threshold(.5)<jupyter_output>Sensitivity: 0.5925925925925926
Specificity: 0.8<jupyter_text>Not finding any major gains in accuracy and losing a good bit of accuracy to acheive the lowest possible outcome of false negatives, we will refrain from modfying the model in this way. The default split of .5 will be retained<jupyter_code>#ROC curve to measure
fpr, tpr, thresholds = metrics.roc_curve(Y_test, y_pred_prob)
plt.plot(fpr, tpr)
plt.xlim([0.0, 1.0])
plt.ylim([0.0, 1.0])
plt.title('ROC curve for diabetes classifier')
plt.xlabel('False Positive Rate (1 - Specificity)')
plt.ylabel('True Positive Rate (Sensitivity)')
plt.grid(True)
#Checking area under curve as alternative to classification accuracy 
print(metrics.roc_auc_score(Y_test, y_pred_prob))
cross_val_score(logreg, X, y, cv=10, scoring='roc_auc').mean()<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.11667356226727348, 'fraction_numerical': 0.020036645191796206, 'mean_word_length': 6.729556875285518, 'pattern_counts': {'":': 0, '<': 62, '<?xml version=': 0, '>': 67, 'https://': 4, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '16', 'dir': 'jupyter-structured-clean-dedup', 'id': '310022', 'license': 'BSD-2-Clause', 'n_tokens_mistral': 5912, 'n_tokens_neox': 5508, 'n_words': 1494, 'path': 'Notebooks/Logistic Classifier for Death Outcome.ipynb', 'repo_name': 'yeagercmbpd/Identifying-Serious-Cases-of-Colic-in-Equines', 'size': '16919'}