Kolejna środa, a z nią kolejna dawka przykładów. Dzisiejszy temat to rozwinięcie parowych testów statystycznych z poprzedniego przykładu o większą ilość zbiorów danych. Oprócz tego znajdzie się też parę innych nowości. Gorąco zachęcam do lektury.
Eksperymenty na wielu zbiorach danych
Dzisiaj nauczymy się wykonywać eksperymenty nie tylko na wielu klasyfikatorach, ale również na wielu zbiorach danych. W szybki sposób bez rozpisywania się nad tym co już było przeprowadzimy całą procedurę krok po kroku.
Modele klasyfikacji
Tym razem zaczniemy od klasyfikatorów. Najpierw przygotowujemy słownik zawierający wybrane modele klasyfikacji, które zostaną użyte do dzisiejszych eksperymentów. Wybieramy trzy:
- Naiwny klasyfikator bayesowski - GNB
- Maszynę wektorów nośnych - SVM
- Klasyfikator k najbliższych sąsiadów - kNN
- Drzewo klasyfikacyjne i regresyjne - CART
from sklearn.naive_bayes import GaussianNB
from sklearn.neighbors import KNeighborsClassifier
from sklearn.tree import DecisionTreeClassifier
from sklearn.svm import SVC
clfs = {
'GNB': GaussianNB(),
'SVM': SVC(),
'kNN': KNeighborsClassifier(),
'CART': DecisionTreeClassifier(random_state=1234),
}
Dane
Teraz wybieramy dane na których chcemy przeprowadzić eksperymenty. W tym wypadku wybierzemy 20 różnych zbiorów danych rzeczywistych. Zestawy te można pobrać w postaci plików CSV z repozytorium benchmark_datasets.
datasets = ['australian', 'balance', 'breastcan', 'cryotherapy', 'diabetes',
'digit', 'ecoli4', 'german', 'glass2', 'heart', 'ionosphere',
'liver', 'monkthree', 'shuttle-c0-vs-c4', 'sonar', 'soybean',
'vowel0', 'waveform', 'wisconsin', 'yeast3']
Eksperyment
Mając już dane i wybrane metody, które chcemy przetestować. Na początku musimy się do eksperymentów trochę przygotować. Odczytujemy liczbę wybranych zbiorów za pomocą funkcji len(). Wybieramy liczbę podziałów oraz powtórzeń. Tworzymy bardzo dobrze już znany obiekt ze stratyfikowaną wielokrotną walidacją krzyżową rskf. Następnie tworzymy zerową macierz na wyniki scores. Macierz ta posiada aż trzy wymiary.
import numpy as np
from sklearn.model_selection import RepeatedStratifiedKFold
n_datasets = len(datasets)
n_splits = 5
n_repeats = 2
rskf = RepeatedStratifiedKFold(n_splits=n_splits, n_repeats=n_repeats, random_state=1234)
scores = np.zeros((len(clfs), n_datasets, n_splits * n_repeats))
Kolejny krok to już przeprowadzanie eksperymentów na wcześniej wybranych zbiorach danych. Najpierw w pętli wczytujemy te dane. Pamiętajmy o tym, żeby zgadzały się ścieżki do plików ze zbiorami. Następnie w kolejnej pętli wczytany aktualnie zbiór dzielimy na foldy. Klonujemy modele za pomocą funkcji clone i każdy z tych modeli uczymy oraz dokonujemy predykcji na odpowiednich podzbiorach danych. Na podstawie uzyskanych wyników obliczamy dokładność klasyfikacji i zapisujemy w macierzy scores. W ostatnim kroku zapisujemy wszystkie wyniki do pliku.
from sklearn.base import clone
from sklearn.metrics import accuracy_score
for data_id, dataset in enumerate(datasets):
dataset = np.genfromtxt("datasets/%s.csv" % (dataset), delimiter=",")
X = dataset[:, :-1]
y = dataset[:, -1].astype(int)
for fold_id, (train, test) in enumerate(rskf.split(X, y)):
for clf_id, clf_name in enumerate(clfs):
clf = clone(clfs[clf_name])
clf.fit(X[train], y[train])
y_pred = clf.predict(X[test])
scores[clf_id, data_id, fold_id] = accuracy_score(y[test], y_pred)
np.save('results', scores)
Analiza wyników
Mając już wyniki możemy przystąpić do ich analizy. Na samym początku musimy wczytać wcześniej zapisane wyniki. Warto zwrócić uwagę na wymiary macierzy, która została wczytana. Można to zrobić za pomocą odpowiedniego pola shape w obiekcie typu ndarray z biblioteki numpy.
scores = np.load('results.npy')
print("\nScores:\n", scores.shape)
Wypisane zostaną wymiary macierzy:
Scores:
(4, 20, 10)
Oznacza to, że w naszej macierzy znajdują się zebrane wszystkie wyniki dla 4 klasyfikatorów, 20 baz i 10 foldów. Teraz musimy te wyniki uśrednić po foldach. Pozwoli nam to zrobić funkcja mean z ustawieniem parametru axis=2 co oznacza, że uśredniony zostanie 3 wymiar (liczymy od 0, 1, 2 - trzeci wymiar). Następnie dokonamy transpozycji tej macierzy, tak aby w kolumnach (drugi wymiar) naszej macierzy znajdowały się metody, a w wierszach zbiory danych.
mean_scores = np.mean(scores, axis=2).T
print("\nMean scores:\n", mean_scores)
Zostanie wypisana taka tabela:
Mean scores:
[[0.78188406 0.65869565 0.70289855 0.81449275]
[0.904 0.9048 0.8496 0.7744 ]
[0.96267175 0.96999249 0.97290146 0.94511593]
[0.85 0.62222222 0.71666667 0.85555556]
[0.75004244 0.69529751 0.68164417 0.69337917]
[1. 0.99440491 0.99213064 1. ]
[0.83489903 0.98511853 0.98658911 0.96870061]
[0.7305 0.71 0.6895 0.691 ]
[0.47895903 0.92059801 0.90869324 0.87131783]
[0.84074074 0.65 0.65925926 0.74074074]
[0.89030181 0.94160966 0.83468813 0.88738431]
[0.56086957 0.69565217 0.65652174 0.62463768]
[0.9196724 0.96390663 0.98375102 0.9729484 ]
[0.99644734 0.99672056 0.99945355 1. ]
[0.67102207 0.79088269 0.7858885 0.73554007]
[1. 0.97777778 0.94666667 1. ]
[0.93725837 0.96155463 0.99290366 0.98481003]
[0.81 0.8668 0.8214 0.7521 ]
[0.95420349 0.96566804 0.95706064 0.93846865]
[0.30521203 0.95148444 0.94609155 0.93059082]]
Rangi
W kolejnym kroku przystąpimy do wyliczenia rang na podstawie uśrednionych wyników. Pamiętajmy, że im wyższa ranga tym metoda jest lepsza. W tym wypadku rangi przypisujemy od 1 do 4, gdzie 4 to jest liczba testowanych modeli. W przypadku remisów dokonujemy uśrednienia. Pomoże nam w tym funkcja rankdata() z biblioteki numpy. Na koniec uzyskane rangi zamieniamy z listy na macierz numpy i wypisujemy.
from scipy.stats import rankdata
ranks = []
for ms in mean_scores:
ranks.append(rankdata(ms).tolist())
ranks = np.array(ranks)
print("\nRanks:\n", ranks)
Uzyskane rangi:
Ranks:
[[3. 1. 2. 4. ]
[3. 4. 2. 1. ]
[2. 3. 4. 1. ]
[3. 1. 2. 4. ]
[4. 3. 1. 2. ]
[3.5 2. 1. 3.5]
[1. 3. 4. 2. ]
[4. 3. 1. 2. ]
[1. 4. 3. 2. ]
[4. 1. 2. 3. ]
[3. 4. 1. 2. ]
[1. 4. 3. 2. ]
[1. 2. 4. 3. ]
[1. 2. 3. 4. ]
[1. 4. 3. 2. ]
[3.5 2. 1. 3.5]
[1. 2. 4. 3. ]
[2. 4. 3. 1. ]
[2. 4. 3. 1. ]
[1. 4. 3. 2. ]]
Nadal taki sposób jest mało czytelny i trudno jest jednoznacznie określić, która metoda jest nalepsza. Dlatego teraz dokonujemy uśrednienia uzyskanych rang i wypisujemy uzyskane wyniki.
mean_ranks = np.mean(ranks, axis=0)
print("\nMean ranks:\n", mean_ranks)
Uśrednione rangi:
Models: ['GNB', 'SVM', 'kNN', 'CART']
Mean ranks: [2.25 2.85 2.5 2.4 ]
Na podstawie tych wyników można stwierdzić, że metody są na dość zbliżonym poziomie, ale nadal nie wiemy o tym, czy te wyniki są statystycznie istotne.
Testy parowe
Mając już rangi możemy przystąpić do testów parowych, które pierwszy raz pojawiły się na poprzednich przykładach. Dzisiaj zamiast testu T Studenta zostanie użyty test Wilcoxona. Użycie tego testu oznacza zaimportowanie oraz wywołanie funkcji ranksums z biblioteki scipy. Ustawiamy wartość alpha, która będzie przyjętym przez nas progiem ufności podczas sprawdzania istotności.
from scipy.stats import ranksums
alfa = .05
w_statistic = np.zeros((len(clfs), len(clfs)))
p_value = np.zeros((len(clfs), len(clfs)))
for i in range(len(clfs)):
for j in range(len(clfs)):
w_statistic[i, j], p_value[i, j] = ranksums(ranks.T[i], ranks.T[j])
Teraz utworzymy dużo czytelniejsze tabele za pomocą biblioteki tabulate. Odbywa to się podobnie jak w ostatnim przykładzie. Jedyna różnica to, że ustawimy dynamicznie nazwy nagłówków i nazwy kolumn z wykorzystaniem kluczy słownika modeli klasyfikacji.
from tabulate import tabulate
headers = list(clfs.keys())
names_column = np.expand_dims(np.array(list(clfs.keys())), axis=1)
w_statistic_table = np.concatenate((names_column, w_statistic), axis=1)
w_statistic_table = tabulate(w_statistic_table, headers, floatfmt=".2f")
p_value_table = np.concatenate((names_column, p_value), axis=1)
p_value_table = tabulate(p_value_table, headers, floatfmt=".2f")
print("\nw-statistic:\n", w_statistic_table, "\n\np-value:\n", p_value_table)
Po wypisaniu widzimy uzyskane statystyki i wartości p-wartości:
w-statistic:
GNB SVM kNN CART
---- ----- ----- ----- ------
GNB 0.00 -1.61 -0.62 -0.50
SVM 1.61 0.00 0.99 1.30
kNN 0.62 -0.99 0.00 0.28
CART 0.50 -1.30 -0.28 0.00
p-value:
GNB SVM kNN CART
---- ----- ----- ----- ------
GNB 1.00 0.11 0.53 0.62
SVM 0.11 1.00 0.32 0.19
kNN 0.53 0.32 1.00 0.78
CART 0.62 0.19 0.78 1.00
Następnie musimy wyznaczyć macierz przewag. Procedura przebiega w dość analogiczny sposób jak w ostatnich przykładach. Jedynie teraz używamy zmiennej w_statistic, do której zostały zapisane wyniki z parowego testu Wilcoxona.
advantage = np.zeros((len(clfs), len(clfs)))
advantage[w_statistic > 0] = 1
advantage_table = tabulate(np.concatenate(
(names_column, advantage), axis=1), headers)
print("\nAdvantage:\n", advantage_table)
Poniżej wyniki:
Advantage:
GNB SVM kNN CART
---- ----- ----- ----- ------
GNB 0 0 0 0
SVM 1 0 1 1
kNN 1 0 0 1
CART 1 0 0 0
Pamiętajmy jednak, że po pierwsze posiadamy wyznaczone rangi, a po drugie nie wiemy czy te wyniki są istotne statystycznie. Teraz możemy przejść do kluczowego momentu, który pozwoli nam to określić. Robimy to w dosłownie identyczny sposób jak było to pokazane w ostatnim przykładzie.
significance = np.zeros((len(clfs), len(clfs)))
significance[p_value <= alfa] = 1
significance_table = tabulate(np.concatenate(
(names_column, significance), axis=1), headers)
print("\nStatistical significance (alpha = 0.05):\n", significance_table)
Wypisujemy wyniki i sprawdzamy:
Statistical significance (alpha = 0.05):
GNB SVM kNN CART
---- ----- ----- ----- ------
GNB 0 0 0 0
SVM 0 0 0 0
kNN 0 0 0 0
CART 0 0 0 0
Niestety, w żadnym testowanym zestawieniu nie można odrzucić hipotezy zerowej. Oznacza to, że wyniki jakie ostatecznie uzyskaliśmy są nieistotne statystycznie po przeprowadzeniu odpowiednich testów parowych. Mimo różnic jakie występują pomiędzy przetestowanymi modelami klasyfikacji eksperymenty przeprowadzone na wybranych przez nas zbiorach danych nie pozwalają stwierdzić, która z tych metod jest lepsza z istotnością statystyczną. Jest to dla nas bardzo ważna lekcja mówiąca o tym, że lepszy nie zawsze musi oznaczać lepszy statystycznie.
