INTRO
Poniżej druga część analizy danych dotyczących kontuzji dolnych kończyn graczy NFL (1 część tutaj). W tym wpisie spróbuję sprawdzić czy dane udostępnione przez NFL pozwalają na zbudowanie użytecznego modelu przewidującego czy dany gracz będzie kontuzjowany czy nie. Więcej informacji o samych danych i szczegółach wyzwania możecie znaleźć na platformie Kaggle pod tym linkiem.
PRZYGOTOWANIE DANYCH
Najpierw wgrywamy biblioteki, które będą nam potrzebne w trakcie budowania modelu.
import numpy as np
import pandas as pd
import gc
import matplotlib.pyplot as plt
%matplotlib inline
#models
from sklearn.tree import DecisionTreeClassifier, DecisionTreeRegressor
from sklearn.linear_model import LogisticRegression
from sklearn.ensemble import RandomForestClassifier
import xgboost as xgb
#validation
from sklearn.model_selection import KFold
from sklearn.model_selection import StratifiedKFold
#succes metrics
import ml_metrics as metrics
from ml_metrics import auc
from sklearn.metrics import accuracy_score, confusion_matrix, f1_score, recall_score, precision_score, classification_report
#imbalanced dataset
from imblearn.under_sampling import RandomUnderSampler
from imblearn.over_sampling import RandomOverSampler
from imblearn.over_sampling import SMOTE
import collections
pd.set_option('display.max_rows', 500)
pd.set_option('display.max_columns', 500)
pd.set_option('display.width', 1000)
Wgrywamy dane. „InjuryRecord” to oryginalny plik udostępniony przez NFL, a dwa pozostałe są wynikiem operacji przeprowadzonych w pierwszym wpisie dotyczącym wizualizacji danych (LINK).
injury_df_ext = pd.read_csv('injury_df_ext.csv')
injury_df = pd.read_csv('InjuryRecord.csv')
play_list_df = pd.read_csv('play_list_df.csv')
W dataframe „injury_df” mamy sporo brakujących danych w kolumnie „PlayKey”. Uzupełniamy brakujące dane wykorzystując kolumnę „GameID” i dodając do niej string ’-1′.
injury_df['new']='-1' injury_df['PlayKey'].fillna(injury_df['GameID']+injury_df['new'], inplace=True) del injury_df['new']
Dane w kolumnie „Temperature” zawierają całkiem sporą ilość wartości „-999”, która wskazuje na brak danych.
play_list_df['Temperature'].value_counts()
Żeby uzupełnić te brakujące dane zastępuję wartości „-999” średnią wartością temperatury dla każdego rodzaju stadionu.
play_list_df['Temperature'] = play_list_df['Temperature'].replace(-999, np.NaN)
play_list_df['Temperature'] = play_list_df.groupby('StadiumType_cat')['Temperature'].transform(lambda x: x.fillna(x.mean()))
Sprawdzamy jak wygląda histogram kolumny „Temperature” po naszych zmianach. Rozkład jest w miarę zbliżony do normalnego.
play_list_df['Temperature'].hist()
Łączymy zbiory, Mamy 3 kolumny, które powtarzają się w obu zbiorach („PlayerKey”,”GameID” i „PlayKey”.
play_list_df.columns, injury_df.columns
W pierwszym zbiorze mamy 16 kolumn, w drugim 9, więc po połączeniu powinniśmy dostać 16+9-3=22 kolumny. Nowy zbiór nazywamy oddzielną nazwą „play_injury_df”
play_injury_df = play_list_df.merge(injury_df, how='left') play_injury_df.info()
W tym momencie nasz dataframe wygląda następująco:
play_injury_df.head()
Usuwamy niepotrzebne kolumny „BodyPart” i „Surface”. „Surface” opisuje rodzaj nawierzchni co mamy już podane w kolumnie „FiledType”. „BodyPart” wskazuje na rodzaj kontuzji (kolano, staw skokowy, pięta itp). Nas interesuje fakt czy kontuzja miała miejsce czy nie, bez wdawania się w szczegóły dotyczące rodzaju kontuzji.
if 'BodyPart' in play_injury_df: del play_injury_df['BodyPart'] if 'Surface' in play_injury_df: del play_injury_df['Surface']
Po połączeniu zbiorów w kolumnach „DM_M1”, „DM_M7”, „DM_M28” i „DM_M42” mamy dużo brakujących danych. Wartości NaN zastępujemy zerami (0).
injury_dm = ['DM_M1', 'DM_M7', 'DM_M28', 'DM_M42']
for injury in injury_dm:
play_injury_df[injury] = play_injury_df[injury].fillna(0).astype(int)
Tworzymy nową kolumnę „IsInjured”, która będzie naszą zmienną docelową, którą model będzie przewidywał. „IsInjured” przyjmuje 1, gdy doszło do kontuzji oraz 0, gdy nie doszło do kontuzji.
play_injury_df['IsInjured']=play_injury_df[['DM_M1','DM_M7','DM_M28','DM_M42']].sum(axis=1) play_injury_df['IsInjured']=play_injury_df['IsInjured'].map(lambda x: 1 if x>0 else 0)
Nasz zbiór zawiera wiele cech, ale część z nich nie będzie pomocna przy budowaniu modelu. Może wprowadzać tylko niepotrzebny szum.
play_injury_df.columns
Usuwamy niepotrzebne kolumny.
play_injury_df = play_injury_df.drop(columns=['GameID', 'PlayerKey', 'PlayKey', 'StadiumType', 'Weather','Position','DM_M1', 'DM_M7',
'DM_M28', 'DM_M42'])
Zostaje nam 11 kolumn, ale tylko kilka z nich ma wartości numeryczne (int64 lub float64), które są wartościowe dla przyszłego modelu.
play_injury_df.info()
Zmieniamy zmienne kategorialne (object) na wartości numeryczne. Użyjemy do tego funkcji „get_dummies”.
cat_features = play_injury_df.select_dtypes(include= np.object).columns play_injury_df_cat = pd.get_dummies(play_injury_df, columns=cat_features, drop_first=True) play_injury_df_cat = play_injury_df_cat.drop(columns=['PlayerDay', 'PlayerGame', 'Temperature', 'PlayerGamePlay', 'IsInjured']) play_injury_df = pd.concat([play_injury_df, play_injury_df_cat], axis=1)
Powstaje nam zbiór z 50 kolumnami.
play_injury_df.info()
Żeby nie powtarzać wszystkich powyższych zmian zapisujemy nasz dataframe do do pliku CSV.
play_injury_df.to_csv('play_injury_df.csv', index=False)
TRENOWANIE I WALIDACJA MODELU
Ustalamy cechy, które będziemy wykorzystywać przy budowaniu modelu predykcyjnego.
feats = play_injury_df.select_dtypes(include= np.number).columns #in feats we list only numerical black_list = ['IsInjured'] #zmienna docelowa feats = [feat for feat in feats if feat not in black_list]
Tworzymy dane do trenowania modelu: X_train jako dane treningowe, y_train jako zmienna docelowa.
X_df = play_injury_df[feats] X_train = X_df.values y_train = play_injury_df['IsInjured'].values.astype(np.int8)
Lista modeli, których wyniki będziemy porównywać.
models = [
xgb.XGBClassifier(random_state=8888,max_depth=5, n_estimators=100, objective='binary:logistic'),
LogisticRegression(random_state=8888,solver='lbfgs', max_iter=1000),
DecisionTreeClassifier(random_state=8888),
RandomForestClassifier(n_estimators=50, max_depth=10, min_samples_split=10, min_samples_leaf=5,random_state=2019),
]
Tworzymy funkcję model_train_predict, która łączy w sobie trenowanie modelu i jego cross walidację (walidację krzyżową) z różnymi metrykami sukcesu (accuracy, recall, precision i f1 score). Walidacja krzyżowa polega na sprawdzaniu wyników modelu przy K-różnych podziałach naszych danych na zbiór treningowy (train) i testowy (test). Więcej o walidacji krzyżowej możecie przeczytać tutaj lub tutaj. Zdecydowałem się na użycie walidacji StratifiedKFold ze względu na bardzo niezbilansowane klasy w moim zbiorze danych (ilość przypadków, w której doszło do kontuzji jest dramatycznie mniejsza od przypadków, gdy kontuzji nie było. O metrykach sukcesu dla modeli klasyfikacyjnych możecie przeczytać tutaj. Zwykła dokładność predykcji (accuracy) nie jest wystarczającym parametrem świadczącym o jakości modelu (więcej informacji poniżej).
def model_train_predict(X, y, model, folds=3, succes_metric1=accuracy_score, succes_metric2=recall_score, succes_metric3=precision_score, succes_metric4=f1_score):
cv = StratifiedKFold(n_splits=folds, shuffle=True, random_state=8888)
accu_scores = []
f1_scores = []
recall_scores=[]
precision_scores = []
for train_idx, test_idx in cv.split(X,y):
X_fold_train, X_fold_test = X[train_idx], X[test_idx]
y_fold_train, y_fold_test = y[train_idx], y[test_idx]
model.fit(X_fold_train, y_fold_train)
y_pred = model.predict(X_fold_test)
accu_score = succes_metric1(y[test_idx], y_pred)
accu_scores.append(accu_score)
recall_score = succes_metric2(y[test_idx], y_pred)
recall_scores.append(recall_score)
precision_score = succes_metric3(y[test_idx], y_pred)
precision_scores.append(precision_score)
f1_score = succes_metric4(y[test_idx], y_pred)
f1_scores.append(f1_score)
print('Mean accuracy score= ', np.mean(accu_scores), ', standard deviation (std)=' , np.std(accu_scores) )
print('Mean recall score= ', np.mean(recall_scores), ', standard deviation (std)=' , np.std(recall_scores) )
print('Mean precision score= ', np.mean(precision_scores), ', standard deviation (std)=' , np.std(precision_scores) )
print('Mean f1 score= ', np.mean(f1_scores), ', standard deviation (std)=' , np.std(f1_scores) )
Jak widać na poniższych screenach modele nie poradziły sobie z tak niezbilansowanymi danymi. Bardzo wysoka dokładność wynika z tego, że model za każdym razem przewiduje brak kontuzji i w ponad 99,96% ma rację, bo przypadków z kontuzjami jest tak mało, że „gubią się” one wśród przypadków bez kontuzji. Metryka f1 score, czyli średnia harmoniczna wartości recall i precision, jest o wiele bardziej wiarygodna, gdyż uwzględnia dokładność przewidywania dla obu klas. W naszych danych jest tak mało przypadków klasy 1 (kontuzja), że model nie jest w stanie policzyć f1_score i zwraca ostrzeżenie (’warning’).
for i in range (0,len(models)):
print ('Model: ', models[i])
model = models[i]
model_train_predict(X_train, y_train, model)
print('______________________________________________________________________________')
Jak widzimy poniżej różnica ilości danych dla obu klas jest ogromna. Musimy coś zrobić, aby poprawić tą sytuację.
play_injury_df['IsInjured'].value_counts()
W przypadku tak nierównomiernego rozkładu elementów w klasach mamy do wyboru dwie możliwości:
- Oversampling – czyli zwiększenie ilości elementów w klasie mniejszościowej (spróbujemy dwóch algorytmów RandomOverSampler i SMOTE)
- Undersampling – czyli zmniejszenie ilości elementów w klasie większościowej (algorytm RandomUnderSampler)
Zacznijmy od udersamplingu i prostego algorytmu RandomUnderSampler, który losowo usuwa elementy z klasy większościowej, aż wyrówna ilość elementów w obu klasach.
rus = RandomUnderSampler(random_state=0) X_resampled, y_resampled = rus.fit_resample(X_train, y_train)
Undersampling pomógł w tym, że modele mogły wyliczyć już wszystkie oczekiwanie przez nas metryki sukcesu. Otrzymane wartości accuracy i f1_score jak na tak małą ilość danych też są całkiem przyzwoite. Ale czy taki model ma jakąś wartość, gdy usuwamy aż tyle danych? Moim zdaniem nie.
for i in range (0,len(models)):
print ('Model: ', models[i])
model = models[i]
model_train_predict(X_resampled, y_resampled, model)
print('______________________________________________________________________________')
Zobaczmy jak wyglądają wyniki modeli po zastosowaniu oversamplingu. Najpierw RandomOverSampler, a potem SMOTE.
res = RandomOverSampler(random_state=8888) X_resampled, y_resampled = res.fit_resample(X_train, y_train)
Ilość elementów w klasie mniejszościowej zrównała się z klasą większościową.
np.count_nonzero(y_resampled)
Wyniki dla RandomOverSampler:
SMOTE – podobnie jak wcześniej otrzymujemy dwie równe klasy
X_resampled, y_resampled = SMOTE().fit_resample(X_train, y_train) print(sorted(collections.Counter(y_resampled).items()))
Wyniki modeli:
WOW!!! Wyniki od 94 do 99,9% i to bez żadnej optymalizacji parametrów modelu. Czyżby oversampling był aż tak skuteczny? Czy aby na pewno nasz model jest w stanie tak dokładnie przewidzieć kontuzje graczy opierając się na sztucznym uzupełnianiu danych w mniejszościowej klasie?
Chyba niestety nie…:) Problem polega na tym, że użyliśmy oversamplingu, a dopiero potem zwalidowaliśmy krzyżowo nasz model. Trochę więcej na ten temat możecie przeczytać tutaj. Poniżej zmodyfikowana funkcja model_train_oversample_predict, która uwzględnia oversampling, ale dopiero po podziale danych na foldy.
def model_train_oversample_predict(X, y, model, folds=3, succes_metric1=accuracy_score, succes_metric2=recall_score, succes_metric3=precision_score, succes_metric4=f1_score):
cv = StratifiedKFold(n_splits=folds, shuffle=True, random_state=2019)
res = RandomOverSampler(random_state=8888)
accu_scores = []
f1_scores = []
recall_scores=[]
precision_scores = []
for train_idx, test_idx in cv.split(X,y):
X_fold_train, X_fold_test = X[train_idx], X[test_idx]
y_fold_train, y_fold_test = y[train_idx], y[test_idx]
X_resampled, y_resampled = res.fit_resample(X_fold_train, y_fold_train)
model.fit(X_resampled, y_resampled)
y_pred = model.predict(X_fold_test)
accu_score = succes_metric1(y[test_idx], y_pred)
accu_scores.append(accu_score)
recall_score = succes_metric2(y[test_idx], y_pred)
recall_scores.append(recall_score)
precision_score = succes_metric3(y[test_idx], y_pred)
precision_scores.append(precision_score)
f1_score = succes_metric4(y[test_idx], y_pred)
f1_scores.append(f1_score)
print('Mean accuracy score= ', np.mean(accu_scores), ', standard deviation (std)=' , np.std(accu_scores) )
print('Mean recall score= ', np.mean(recall_scores), ', standard deviation (std)=' , np.std(recall_scores) )
print('Mean precision score= ', np.mean(precision_scores), ', standard deviation (std)=' , np.std(precision_scores) )
print('Mean f1 score= ', np.mean(f1_scores), ', standard deviation (std)=' , np.std(f1_scores) )
Wyniki modeli przy użyciu model_train_oversample_predict. Jak widać jednak nasz klasyfikator nie jest taki super:(:(:(
Sprawdźmy jeszcze, które cechy są wg naszych modeli najbardziej istotne.
def draw_feature_importances(model, features):
importances = model.feature_importances_
indices = np.argsort(importances)[::-1]
plt.figure(figsize=(10, 5))
plt.title("Feature importances")
plt.bar(range(X_resampled.shape[1]), model.feature_importances_[indices],
color="b", align="center")
plt.xticks(range(X_resampled.shape[1]), [ features[x-3] for x in indices])
plt.xticks(rotation=90)
plt.xlim([-1, X_resampled.shape[1]])
plt.show()
for i in [0,2,3]:
draw_feature_importances(models[i], feats)
XGBClassifier
DecisionTreeClassifier
RandomForestClassifier
