Lineárna klasifikácia

Dva typické problémy pri učení strojov sú regresia a klasifikácia. Pri regresii sa snažíme predpovedať kontinuálnu hodnotu napr. príjem na základe IQ. Cieľom klasifikácie je predpovedať kategóriu napr. pohlavie alebo typ vzdelania. Minule som sa zaoberal regresiou. Klasifikácia prebieha analogicky. V tomto článku demonštrujem klasifikáciu pomocou maximalizácie vierohodnosti.

Aby bola prezentácia konkrétnejšia vygenerujme si demo dáta. Dané sú dva kontinuálne prediktory a dve kategórie. Povedzme, že zisťujeme pohlavie jedincov určitého druhu chobotnice a prvým prediktorom je dĺžka hektakotylového ramena zatiaľčo druhý prediktor popisuje hmotnosť jedinca. Dáta simulujem ako kombináciu dvoch gaussovských rozdelení.

import pylab as plt
import numpy as np
np.random.seed(6)
m=np.random.multivariate_normal([10, 150],[[2,5],[5,30]],100)
f1=np.random.multivariate_normal([6, 140],[[4,-5],[-5,30]],100)
plt.close()

plt.plot(f1[:,0],f1[:,1],'or')
plt.plot(m[:,0],m[:,1],'xb')
plt.xlabel('Prediktor 1')
plt.xlabel('Prediktor 2')
plt.xlim([0, 15])
plt.ylim([100,200])
plt.show()

Kategórie sa z časti prelínajú. To ale nevadí, aspoň bude úloha zaujímavejšia.

Dáta sú dané vo forme Nx1 vektorov y, x_1, x_2, pričom elementy y určujú kategóriu daného jedinca y_n \in 0,1. Našou úlohou je získať lineárny model s vhodnými parametrami.

x=np.concatenate((m,f1))
x=np.concatenate((np.ones((x.shape[0],1)),x),axis=1)
y=np.concatenate((np.zeros(100),np.ones(100)))

Pri regresii sme optimalizovali p(y|x,b_0,b_1,\sigma)=\prod_{n=1}^N p(y_n| x_n,b_0,b_1,\sigma)=\prod_{n=1}^N \mathcal{N}(y_n|b_0 + b_1 x_n,\sigma). Pri klasifikácii môžeme postupovať analogicky, akurát gausovské rozdelenie v poslednom kroku nevyhovuje. Namiesto neho potrebujeme rozdelenie, ktoré by modelovalo dichotomické premenné. Vhodným rozdelením je Bernouliho rozdelenie \mathcal{B}(y|q)= q^y (1-q)^(1-y). Ak y=1 tak \mathcal{B}(y|q)= q a ak y=0 tak \mathcal{B}(y|q)= 1-q. Zvolíme teda
p(y|b,x_1,x_2)=\prod_{n=1}^N p(y_n| x_1,x_2,b)=\prod_{n=1}^N \mathcal{B}(y_n|x_1,x_2,b)= q^y (1-q)^(1-y)
Ideálne by sme zvolili q= b_0 + b_1 x_{n1}+b_2 x_{n2}. Tým sa snažíme docieliť, aby v lineárnej závislosti od parametrov a dát q= 1 ak y= 1 a q= 0 ak y= 0. Problém je stále v tom, že b_0 + b_1 x_{n1}+b_2 x_{n2} poskytuje kontinuálne hodnoty. Preto musíme zvoliť funkciu, ktorá by mapovala kontinuálne hodnoty (-\infty, \infty) do rozmedzia (0,1). Sigmoidálna funkcia nám to umožní: \sigma(x)=\frac{1}{1+e^{-1}}

Sigmoidálna funkcia nie je jedinou možnou voľbou, avšak táto funkcia má viaceré atraktívne vlastnosti. Je monotónna na celom svojom intervale a jej derivácia má jednoduchú formu \frac{\delta \sigma(x)}{\delta x}=\sigma(x) (1-\sigma(x)). Platí teda q_n= \sigma(b_0 + b_1 x_{n1}+b_2 x_{n2}).

Týmto sme s formuláciou vierohodnosti modelu hotoví a môžeme nájsť optimálne parametre. Nájdeme deriváciu vierohodnosti a nájdeme parametre pre prípad, že derivácia je rovná nule. Problémom je, že pri takomto postupe sa nám žiaľ nepodarí parametre “vyslobodiť”. Parametre sú na sebe závislé. Musíme zvoliť iteratívnu optimalizáciu. Použijeme Newtonovu metódu. Animácia z Wikipédie ilustruje jej fungovanie:

Našim cieľom je nájsť x pre ktoré platí f'(x)=0. Na začiatku zvolíme náhodne x_1. Nájdeme deriváciu y=f(x) v bode x_1. Táto tvorí tangentu k f(x) v bode x_1 Následne nájdeme x_2 pre ktoré je tangenta rovná nule, t.j. pretína x-ovú osu. Pritom platí

f'(x)=\frac{dy}{dx}=\frac{f(x_t)-f(x_{t+1})}{x_t-x_{t+1}}= \frac{f(x_t)}{x_t-x_{t+1}}

Vyriešime pre x_{t+1}:

x_{t+1}=x_t - \frac{f(x_t)}{f'(x_t)}

V našom prípade chceme nájsť parametre w, pre ktoré gradient vierohodnosti \nabla E je rovný nule. Ak aplikujeme Newtonovu metódu na náš problém získame w_{t+1}=w_t -\nabla \nabla E^{-1} \nabla E. Z formulácie vierohodnosti určíme potrebné derivácie: \nabla E=\sum_n (q_n - y_n)x_n\nabla \nabla E=\sum_n q_n (1 - q_n)x_n x _n^T, pričom q_n= \sigma(w_t x_n). V našom konkrétnom prípade platí q_n= \sigma(b_0 + b_1 x_{n1}+b_2 x_{n2}), kedže w= (b_0, b_1, b_2).

Takéto iteratívne riešenie nie je tragédia. Vierohodnosť je konvexnou funkciou a Newtonova metóda zaručene nájde jej globálne maximum.

Sigmoidálnej funkcii sa nazýva aj logistická funkcia a preto sa odvodený lineárny klasifikátor nesie meno logistická regresia. Implementácia môže vyzerať nasledovne:

import sys
def sigmoid(x):
    return 1.0 / (1.0 + np.exp(-x))

def newton(f,x0,fprime,tol=1e-8,maxiter=100):
    """ Newton method for optimization of with vectors"""
    print 'Newton Optimization started'
    xold=np.copy(x0)
    for i in range(maxiter):
        fxold=f(xold)
        if np.linalg.norm(np.zeros(xold.size)-fxold)<tol:
            print '\nOptimization terminated successfully after %d iterations'%i
            return xold

        xnew=xold-np.array(np.matrix(fprime(xold)).I).dot(fxold)
        sys.stdout.write('.')
        xold=xnew
    print '\nOptimization failed to converge after %d iterations'%maxiter
    return xnew

class LogisticRegression():
    """
    A simple logistic regression model
    """

    def __init__(self,x,y,alpha=0.0):
        """ N - nr of samples, M - nr of features
            x - Nx(M+1) numpy array, the first column is a vector of ones
            y - Nx1 numpy array
            alpha - L2 regularization strength
        """
        # Set the data.
        self.x = x
        self.alpha=alpha
        self.y= y
        self.n = y.shape[0]
        # Initialize parameters to zero, for lack of a better choice.
        self.w = np.zeros(self.x.shape[1])

    def prediction(self,x=None,w=None):
        if x is None: x=self.x
        if w is None: w=self.w
        q = np.zeros(x.shape[0])
        for i in range(x.shape[0]):
            q[i] = sigmoid(np.dot(w, x[i]))
        return q

    def dE(self, w):
        """ Error gradient of the data under
            the current settings of parameters. """
        q = self.prediction(w=w)
        return self.x.T.dot(q-self.y)+w*self.alpha
    def ddE(self,w):
        """ Computes Hessian of the data under
            the current settings of parameters. """
        q = self.prediction(w=w)
        R=np.diag(q*(1-q))
        return self.x.transpose().dot(R.dot(self.x))+self.alpha

    def train(self):
        """ Find dE=0 with Newton method """
        self.w = newton(self.dE,self.w,fprime=self.ddE)

    def checkCorrectness(self,xTest,yTest,w=None):
        """ Gives performance in percent correct """
        if not w is None: self.w=w
        t = np.squeeze(yTest) > 0
        y = np.squeeze(self.prediction(x=xTest) > 0.5)
        return 100 * np.sum(t==y)/ np.size(yTest)

Nakoniec trénujeme klasifikátor na našich dátach a diagnostikujeme správnosť predpovedí modelu:

lr=LogisticRegression(x,y)
lr.train()
print lr.checkCorrectness(x,y)

p2s=np.arange(0,20,0.1)
p1s=np.arange(120,181,1)
D=np.nan*np.ones((p1s.size,p2s.size))
xx=np.ones((1,3))
for i in range(p1s.size):
    for j in range(p2s.size):
        xx[0,1]=p2s[j]
        xx[0,2]=p1s[i]
        D[i,j]=lr.prediction(x=xx)
plt.imshow(D,extent=(p2s[0],p2s[-1],p1s[-1],p1s[0]))
plt.colorbar()
plt.plot(f1[:,0],f1[:,1],'or')
plt.plot(m[:,0],m[:,1],'ob')
plt.xlabel('Prediktor 1')
plt.ylabel('Prediktor 2')

Predpovede modelu sú znázornené farebne pre všetky možné hodnoty prediktorov 1 a 2 v rozmedzí (0,20) a (120,180) respektívne. Predpovede spadajú vďaka sigmoidálnej funkcii do rozhrania (0,1). V prípade lineárneho klasifikátora tvorí hranicu línia. Hranica nie je skoková, ale ako vidieť na prechodnom farebnom odtieni kontinuálna. Model tým vyjadruje neistotu kategorizovania v oblasti kde sa množiny bodov prelínajú.

Lineárny klasifikátor trpí viacerými nedostatkami.

Klasifikácia je citlivá na extrémne hodnoty. Ak pridám jednu od veci hodnotu (, ktorá mohla vzniknúť napr. ako chyba pri prepise dát), tak táto môže silne pohnúť hranicou (extrémnu hodnotu v grafike nižšie nevidieť).

x=np.concatenate((m,f1,np.array([200,200],ndmin=2)))


Neprekvapujúco, ak je hranica medzi kategóriami nelineárna, klasifikátor pracuje mizerne.
Pridajme napríklad ďalšiu množinu bodov a zopakujme analýzy.

f2=np.random.multivariate_normal([11, 165],[[10,10],[10,20]],100)
x=np.concatenate((m,f1,f2))
y=np.concatenate((np.zeros(100),np.ones(200)))


Lineárny klasifikátor nedokáže rekonštruovať  ten oblúk, v ktorom červené bodky ohraničujú modrú množinu. Namiesto toho zvolí ako hranicu líniu, ktorá ide krížom cez množinu modrých bodov a klasifikátor dosiahne slabých 65% správnych predpovedí. V ďalšom článku sa pozrieme na nelineárne klasifikačné metódy.

Pridaj komentár

Zadajte svoje údaje, alebo kliknite na ikonu pre prihlásenie:

WordPress.com Logo

Na komentovanie používate váš WordPress.com účet. Log Out / Zmeniť )

Twitter picture

Na komentovanie používate váš Twitter účet. Log Out / Zmeniť )

Facebook photo

Na komentovanie používate váš Facebook účet. Log Out / Zmeniť )

Google+ photo

Na komentovanie používate váš Google+ účet. Log Out / Zmeniť )

Connecting to %s