Uczenie maszynowe i sztuczne sieci neuronowe/Ćwiczenia 8
Z Brain-wiki
Spis treści
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:
#http://brain.fuw.edu.pl/edu/index.php/Uczenie_maszynowe_i_sztuczne_sieci_neuronowe/Wyk%C5%82ad_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:
# http://brain.fuw.edu.pl/edu/index.php/Uczenie_maszynowe_i_sztuczne_sieci_neuronowe/Wyk%C5%82ad_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
# http://brain.fuw.edu.pl/edu/index.php/Uczenie_maszynowe_i_sztuczne_sieci_neuronowe/Wyk%C5%82ad_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)
Walidacja
Teoria do tej części znajduje się tu: wykład
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.
