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
