Spis treści

1 Instrukcje
2 Kod
3 Kross-walidacja

Instrukcje

W tym ćwiczeniu zbudujemy klasyfikator bazujący na regresji logistycznej. Jego zadaniem będzie określanie prawdopodobieństwa przyjęcia kandydata na studia na podstawie wyników z dwóch egzaminów maturalnych (każdy przeskalowany na zakres 0-100%): z matematyki i z biologii.

Hipotetyczne dane znajdują się w pliku Plik:Reg log data1.txt.

Dla przypomnienia:

Hipoteza

w regresji logistycznej ma postać: [math]h_\theta(x) = \frac{1}{1+\exp(-\theta^Tx)}[/math].
ze względu na stabilność numeryczną obliczeń dobrze jest ograniczyć zakres zmienności exp np do zakresu [1e-8, 1e+8]:

f = exp(x)
if f < 1e-8:
        f = 1e-8
if f>1e8:
        f = 1e8

Parametry

parametry regresji znajdujemy przez maksymalizację funkcji log-wiarygodności:

[math]l(\theta) = \log L(\theta) = \sum_{j=1}^m y^{(j)} \log h(x^{(j)}) + (1 - y^{(j)}) \log (1 - h(x^{(j)}))[/math]

W tym ćwiczeniu zrobimy to za pomocą funkcji optymalizacyjnych z modułu scipy.optimize[1]. Wynikają z tego dwie konsekwencje:

Funkcje te są przystosowane do szukania minimów funkcji celu. Musimy więc podawać im jako argumenty funkcję minus log-wiarygodności
Niektóre algorytmy mogą działać szybciej jeśli zaimplementujemy jawnie postać pochodnej:

[math] \begin{array}{lcl} \frac{\partial}{\partial \theta_i} l(\theta) = \dots =\sum_{j=1}^m (y^{(j)}-h_\theta(x^{(j)}))x_i^{(j)} \end{array} [/math]

Wyniki

Wyniki regresji logistycznej możemy odbierać na dwa sposoby:

obliczyć wartość hipotezy dla badanego wejścia i dopasowanych parametrów: miara ta ma interpretację prawdopodobieństwa przynależności wejścia do klasy 1
dopisać funkcję wykonującą klasyfikację, tzn. porównanie wartości hipotezy z 1/2. Dla wartości hipotezy > 1/2 klasyfikacja zwraca 1, w przeciwnym razie 0.

Kod

Ten kod pomoże Ci zapoznać się z regresją logistyczną. W tym celu będziesz misiał uzupełnić kody funkcji:

sigmoida
funkcjaKosztu
predykcja
funkcjaKosztuReg

Proponuję robić to stopniowo w miarę jak będę sie pojawiać wywołania tych funkcji.

# -*- coding: utf-8 -*-
#
#importujemy potrzebne moduły i klasy
import matplotlib
matplotlib.use('TkAgg')
import numpy as np
import pylab as py
import scipy.optimize as so

#===========================================
# definicje funkcji, które trzeba uzupełnić:
#===========================================
def hipoteza(x, theta):
    '''ta funkcja zwraca wartość hipotezy dla danego wejścia x i parametrów theta'''
    # zaimplementuj hipotezę dla regresji logistycznej zgodnie ze wzorami z wykładu:
    #https://brain.fuw.edu.pl/edu/Uczenie_maszynowe_i_sztuczne_sieci_neuronowe/Wykład_6#Hipoteza
   
   # return h
   pass
 
def funkcjaLogWiarygodnosci(theta, X, y):
    '''Ta funkcja oblicza wartość funkcji log-wiarygodności  dla regresji logistycznej
    używając theta jako parametrów oraz X i y jako zbioru uczącego'''
    # zaimlementuj funkcję l(theta) z wykładu:
    #  https://brain.fuw.edu.pl/edu/Uczenie_maszynowe_i_sztuczne_sieci_neuronowe/Wykład_6#Funkcja_wiarygodno.C5.9Bci
    '''l=0.0
    for j in range(len(y)):
        h = hipoteza(X[j,:],theta)
        l +=  ... 
    return l'''
    pass

def minusFunkcjaLogWiarygodnosci(theta, X, y):
    return (-1.)*funkcjaLogWiarygodnosci(theta, X, y)
    
def pochodnaLogWiarygodnosci(theta, X, y):
    '''ta funkcja oblicza wartość pochodnej funkcji log-wiarygodności
    dla podaanych wartości theta, X i y'''
    # zaimplementuj wzory na dl/dtheta z
    # https://brain.fuw.edu.pl/edu/Uczenie_maszynowe_i_sztuczne_sieci_neuronowe/Wykład_6#Funkcja_wiarygodno.C5.9Bci
    '''
    dl_dtheta = np.zeros(len(theta))
    for i in range(len(theta)):
        for j in range(len(y)):    
            dl_dtheta[i] += ....
    return dl_dtheta
    '''
   pass

def  minusPochodnaLogWiarygodnosci(theta, X, y):
    return (-1)*pochodnaLogWiarygodnosci(theta, X, y)

def klasyfikacja(testX, theta):
    ''' Ta funkcja zwraca wynik klasyfikacji przykładu testX przy parametrach theta.
    Po obliczeniu hipotezy, jeśli otrzymane prawdopodobieństwo jest większe niż 0.5 to zwraca 1 w przeciwnym wypadku zwraca 0'''

    pass



#####################################################################
# definicje funkcji pomocniczych, których nie musisz modyfikować
#======================================================

def rysujDaneGrup(X, y, marker, xlabel, ylabel,legend_list):
    '''X - macierz wartości wejściowych zorganizowana tak, że kolejne przykłady są
    w wierszach, kolumny to kolejne wynmiary  wejścia,
    y - wektor określający przynależność do grupy, indeksy tego wektora odpowiadają wireszom macierzy X,
    marker - zestaw markerów do oznaczania elementów grup, markerów powinno być tyle ile jest grup'''
    p=[]
    for g in np.unique(y):
        g = int(g)
        tmp =py.plot(X[np.where(y==g),0],X[np.where(y==g),1],marker[g])
        p.append(tmp[0])
    py.legend(p,legend_list)
    # Dodajemy napisy
    py.xlabel(xlabel)
    py.ylabel(ylabel)
    
def rysujGranice(theta, X, y, marker, xlabel, ylabel,legend_list):
    rysujDaneGrup(X, y, marker, xlabel, ylabel,legend_list)
    # granica między obszarami dana jest równaniem
    # theta^T x = 0
    x_plot = np.array([np.min(X[:,1]), np.max(X[:,1])])
    y_plot = -1./theta[2]*(theta[1]*x_plot + theta[0])
    py.plot(x_plot,y_plot,'b')
#===================================================





######################### Program #################################
# Wczytanie danych
# Pierwsze dwie kolumny zawierają wyniki egzaminów,
# trzecia kolumna zawiera etykietę (przynależność do grupy)

data = np.loadtxt('reg_log_data1.txt',delimiter=',')
X = data[:, [0, 1]]
y = data[:, 2]


## ==================== Część1:  Rysunki ====================
#  Zawsze dobrze jest pooglądać dane aby nabrać wyczucia do problemu

print 'rysujemy dane: "+" oznacza przykłady z gdzie y = 1 zaś "o" te z y = 0'

rysujDaneGrup(X, y,marker = ('bo','r+'),xlabel='wynik z matematyki', ylabel='wynik z biologii' ,legend_list=(u'Nie przyjęty',u'Przyjęty'));

py.show()

## ============ Część 2: Obliczamy funkcję kosztu i gradient ============
#  W tej części maksymalizujemy funkcję log-wiarygodności dla regresji logistycznej
# Aby to zrobić musisz uzupełnić kody w funkcjach:
# hipoteza
# funkcjaLogWiarygodnosci
# pochodnaLogWiarygodnosci
# Szkielety tych funkcji znajdują sie na początku tego pliku


#  Teraz trzeba zmodyfikować macierz X i dodać jej z lewej strony kolumnę jedynek odopwiadającą
# parametrow theta[0]
N = len(y) # ilość przykładów w zbiorze uczącym
XX = np.concatenate((np.ones((N,1)), X),axis = 1)
# rozmiar wejścia rozszerzonego o jedynki
xDim = XX.shape[1]

# IInicjalizujemy parametry:
theta0 = np.zeros((xDim, 1));

# Oblicz funkcje log-wiarygodności i jej pochodną dla danych początkowych
logWiar = funkcjaLogWiarygodnosci(theta0, XX, y)
pochLogWiar = pochodnaLogWiarygodnosci(theta0, XX, y)

print 'wartość log-wiarygodności dla początkowej thety: '+ str(logWiar)
print 'pochodna log-wiarygodnosci dla poczatkowej thety: '+ str(pochLogWiar)




## ============= Część 3: Optymalizacja  =============
#  Funkcje optymalizujące zaczerpniemy z modułu scipy.optimize
# ponieważ funkcje te są zaimplementowane do mnimalizowania
# zamiast maksymalizować funkcję lowWiarygodności będziemy
# minimalizować tą funkcje przemnożoną przez -1
# czyli minusFunkcjaLogWiarygodnosci
#fprime=minusPochodnaLogWiarygodnosci,

theta_opt = so.fmin_bfgs(minusFunkcjaLogWiarygodnosci, theta0, fprime=minusPochodnaLogWiarygodnosci, args=(XX,y), disp= True)



# wypiszmy theta 
print 'Wartość log wiarygodnosci  dla optymalnych parametrów: '+str(funkcjaLogWiarygodnosci(theta_opt, XX, y))
print 'theta: '+str(theta_opt)


# narysujmy uzyskany podział
rysujGranice(theta_opt, X, y,marker = ('bo','r+'),xlabel='wynik z matematyki', ylabel='wynik z biologii' ,legend_list=(u'Nie przyjęty',u'Przyjęty'));

py.show()

## ============== Część 4: Przewidywanie ==============
#  PO dopasowaniu parametrów nadszedł czas aby zrobić predykcję.
# obliczmy jakie prawdopodobieństwo przyjęcia ma kandydat z wynikami
#  20 z matematyki
#  80 z biologii
#  Do przewidywania wykorzystujemy funkcję hipoteza, bo zgodnie z naszą interpretacją daje ona
# prawdopodobieństwo przyjęcia

prob = hipoteza([1, 20, 80], theta_opt)
print 'dla kandydata z wymnikami 20 z matematyki i 80 z biologii prawdopodobieństwo przyjęcia wynosi: ' +str(prob)

Kross-walidacja

Porównywanie klasyfikatorów na podstawie błędów generalizacji.

Najlepiej byłoby mieć możliwie mały błąd generalizacji — jak go oszacować? Można zastosować następujące techniki:

wiele zbiorów testowych: trzeba mieć dużo danych, żeby wystarczyło na rozsądny zbiór treningowy i kilka testowych
kross-walidacja: najprostsza wersja (leave-one-out):
- wybierz przypadek do odrzucenia
- trenuj klasyfikator na wzystkich przypadkach oprócz tego jednego — na tym jednym oblicz błąd generalizacji
- powtarzaj to dla każdego przypadku
  Zaleta — mażna efektywnie użyć całego zbioru danych do uczenia i testowania,
  Cena — wielokrotne uczenie sieci
Bootstrapowanie — wielokrotnie losuje się z powtórzeniami z pełnego zbioru dwie próby:
- do nauki klasyfikatora
- do testowania.

Krzywe charakterystyki roboczej odbiorcy (ROC)

Rozważmy nasz przykład klasyfikowania kandydatów na przyjętych i nie przyjętych. Możemy podjąć następujące decyzje:

		Stan faktyczny		Razem
		przyjęty	nie przyjęty	Razem
wynik klasyfikacji	przyjęty	TP	FP (błąd: wynik fałszywie pozytywny)	P'
wynik klasyfikacji	nie przyjęty	FN (błąd: wynik fałszywie negatywny)	TN	N'
Razem		P	N

Anglojęzyczna nomenklatura często stosowana do opisu tych możliwości to:

TP: true positive, hit
TN: true negative, correct rejection
FP: false positive, false alarm, Type I error
FN: false negative, with miss, Type II error

Miary oparte na tych definicjach to:

czułość (inne spotykane nazwy: sensitivity, true positive rate TPR, hit rate, recall)

[math]TPR = TP / P = TP / (TP+FN)[/math]

częstość fałszywych alarmów ( false positive rate (FPR), fall-out)

[math]FPR = FP / N = FP / (FP + TN)[/math]

dokładność ( accuracy (ACC))

[math]ACC = (TP + TN) / (P + N)[/math]

specyficzność: (specificity (SPC), True Negative Rate)

[math]SPC = TN / N = TN / (FP + TN) = 1 - FPR [/math]

precyzja: (positive predictive value (PPV), precision)

[math]PPV = TP / (TP + FP)[/math]

negative predictive value (NPV)

[math]NPV = TN / (TN + FN)[/math]

częstość fałszywych odkryć: (false discovery rate (FDR))

[math]FDR = FP / (FP + TP)[/math]

Prawdopodobieństwa podejmowania każdego rodzaju decyzji będą się zmieniały wraz z przesuwaniem progu podejmowania decyzji, czyli wartości hipotezy przy której zaliczamy przypadek do klasy 1.

Krzywa ROC: każdy punkt na tej krzywej otrzymywany jest dla ustalonej wartości progu i ma współrzędne (1−specyficzność, czułość).

Krzywa ROC przydaje się do porównywania różnych klasyfikatorów oraz do wyboru punktu pracy (progu)

Zastosowanie w naszym przykładzie

W części piątej naszego programu dodamy kross-walidację typu leave-one-out. Po kolei odłożymy po jednym przykładzie ze zbioru uczącego, na takim zredukownaym zbiorze nauczymy regresję, a następnie sprawdzimy która z poniższych możliwych sytuacji zachodzi:

TP: stan faktyczny jest pozytywny (y=1) i klasyfikator się nie myli (wynik = 1)
TN: stan faktyczny jest negatywny (y=0) i klasyfikator się nie myli (wynik = 0)
FP: wynik fałszywie pozytywny (fałszywy alarm): stan faktyczny jest negatywny (y=0) ale klasyfikator się myli (wynik = 1)
FN: przegapiony alarm: stan faktyczny jest pozytywny (y=1) i klasyfikator się myli (wynik = 0)

## =============Część 5: walidacja ==============
# zastosujemy kross-walidację typu leave-one-out
# przygotowujemy liczniki:
TP = 0 
TN = 0
FP = 0 
FN = 0

for v in range(len(y)):
    print v ,'/', len(y)
    # odkładamy przykład v do testowania
    testX = XX[v]
    testY = y[v]
    # robimy zredukowany zbiór uczący przez usunięcie przykładu v
    tenX = np.delete(XX,v,axis=0)
    tenY = np.delete(y,v)
    # uczymy regresję
    theta_opt = so.fmin_bfgs(minusFunkcjaLogWiarygodnosci, theta0, fprime=minusPochodnaLogWiarygodnosci, args=(tenX,tenY), disp= False)
   # klasyfikujemy odłożony przykład : proszę uzupełnić funkcję klasyfikacja na początku pliku
    wynik = klasyfikacja(testX, theta_opt)
    # aktualizujemy liczniki; proszę uzupełnić kod:
    if testY == 1:
        if wynik == 1:
            ...
        else:
            ...           
    else:
        if wynik == 1:
            ...
        else:
            ...
print 'TP: ', TP
print 'FP: ', FP
print 'TN: ', TN
print 'FN: ', FN

Dla naszego zbioru uczącego powinniśmy uzyskać:

TP: 55 FP: 6 TN: 34 FN: 5

Krzywa ROC

Aby wykreślić krzywą ROC należy przeprowadzić klasyfikację dla wielu możliwych wartości progu dla hipotezy, powyżej którego uznajemy przypadek za należący do klasy 1.

Modyfikując funkcję klasyfikacja, uzyskujemy następującą funkcje klasyfikującą zależną od progu:

def klasyfikacjaProg(testX, theta_opt,prog):
    prob = hipoteza(testX, theta_opt)
    if prob >prog:
        return 1
    else:
        return 0

Funkcję tą możemy wykorzystać do obliczenia liczebności poszczególnych przypadków klasyfikacji w zależności od progu:

def liczROC(XX,y,progi):
    '''funkcja oblicza FPR i TPR dla zadanych progów'''
    TP = np.zeros(len(progi))
    TN = np.zeros(len(progi))
    FP = np.zeros(len(progi))
    FN = np.zeros(len(progi))
    
    for v in range(len(y)):
        print v ,'/', len(y)
        testX = XX[v]
        testY = y[v]
        tenX = np.delete(XX,v,axis=0)
        tenY = np.delete(y,v)
        theta_opt = so.fmin_bfgs(minusFunkcjaLogWiarygodnosci, theta0, fprime=minusPochodnaLogWiarygodnosci, args=(tenX,tenY), disp= False)
        for ind, prog in enumerate(progi):
            wynik = klasyfikacjaProg(testX, theta_opt,prog)
           #==========================
           #      tu wstaw odpowiedni kawałek kodu
           #==========================
    TPR = TP/(TP+FN)
    FPR = FP/(FP+TN) 
    return (FPR,TPR)

Do wykreślenia krzywej ROC możesz użyć następującego kodu. Zaznaczamy w nim na wykresie wartości progów dla których osiągnięto konkretne wartości FPR i TPR.

progi = np.arange(0.0,1.1,0.1)
FPR,TPR= liczROC(XX,y,progi)
py.plot(FPR,TPR,'o')
py.plot(FPR,TPR)
for ind,pr in enumerate(progi):
    py.text(FPR[ind],TPR[ind],str(pr))
py.xlabel('FPR')
py.ylabel('TPR')
py.xlim((0,1))
py.ylim((0,1))
py.show()

Dopisz powyższe funkcje do głównego modułu.

Anonimowy

Szukaj

Uczenie maszynowe i sztuczne sieci neuronowe/Ćwiczenia 8

Przestrzenie nazw

Więcej

Działania na stronie

Spis treści

Instrukcje

Hipoteza

Parametry

Wyniki

Kod

Kross-walidacja

Porównywanie klasyfikatorów na podstawie błędów generalizacji.

Krzywe charakterystyki roboczej odbiorcy (ROC)

Zastosowanie w naszym przykładzie

Krzywa ROC

Nawigacja

Nawigacja

Narzędzia Wiki

Narzędzia Wiki

Anonimowy

Szukaj

Uczenie maszynowe i sztuczne sieci neuronowe/Ćwiczenia 8

Spis treści

Instrukcje

Hipoteza

Parametry

Wyniki

Kod

Kross-walidacja

Porównywanie klasyfikatorów na podstawie błędów generalizacji.

Krzywe charakterystyki roboczej odbiorcy (ROC)

Zastosowanie w naszym przykładzie

Krzywa ROC

Nawigacja

Narzędzia Wiki

Narzędzia dla stron