Zagadnienia

7. Klasyfikacja

Zadanie klasyfikacji polega na konstrukcji funkcji (klasyfikatora), która na podstawie zaobserwowanych cech będzie przydzielała obserwację do którejś z wcześniej zdefiniowanych grup. Do estymacji funkcji potrzebne są obserwacje, które już zostały sklasyfikowane, będziemy je nazywać próbą uczącą:

yi,Xi⁢⁢i∈1,…,n⁢ niezależne obserwacje, ⁢Xi∈Rp⁢⁢yi∈1,…,K.

Dane Xi oznaczają zaobserwowane cechy, yi grupę, do której obserwacja została zaklasyfikowana.

Cj oznaczaja zbiór tych indeksów i, że yi=j; C1,…,CK są rozłączne o licznościach odpowiednio n1,…,nK.

Funkcja wiarygodności dla opisanych danych wyraża się wzorem:

L⁢π1,…,πK,θ1,…,θK=∏i=1ng⁢yi,Xi=∏i=1nπyi⁢fθyi⁢Xi
pod warunkiem: ⁢∑k=1Kπk=1.

Logwiarygodność to logarytm funkcji wiarygodności:

l=log⁡L;
l=∑i=1nlog⁡πyi+∑i=1nlog⁡fθyi⁢Xi.

Zadanie klasyfikacji dzielimy na dwa kroki:

  1. Estymujemy parametry π1,…,πK oraz θ1,…,θK na podstawie zaobserwowanych par yi,Xi przy użyciu metody największej wiarygodności. Parametry πk możemy interpretować jako prawdopodobieństwa przynależności do danej grupy danych, a θk jako parametry rozkładu w danej grupie (na przykład dla wielowymiarowego rozkładu normalnego, byłyby to średnia μ i macierz kowariancji Σ).

  2. Obserwujemy nowe cechy Xn+1 i przyporządkowujemy im y^i+1 na podstawie zbudowanego przez nas klasyfikatora. Będziemy go także nazywać regułą decyzyjną.

Maksymalizujemy funkcję wiarygodności pod warunkiem ∑k=1Kπk=1 przy użyciu metody mnożników Lagrange'a:

maxπ1,…,πK,θ1,…,θK⁡F⁢π1,…,πK,θ1,…,θK=
maxπ1,…,πK,θ1,…,θK⁡n1⁢log⁡π1+…+nK⁢log⁡πK-λ⁢∑k=1Kπk-1+
+∑i∈C1log⁡fθ1⁢Xi+…⁢∑i∈CKlog⁡fθK⁢Xi.

Liczymy estymatory π^1,…,π^K:

∂⁡F∂⁡πk=nkπk-λ=0⁢⁢∀k=1⁢…,K; (7.1)
∂⁡F∂⁡λ=∑k=1Kπk-1=0; (7.2)

Z równań 7.1 otrzymujemy:

∀k=1⁢…,K⁢⁢nkλ=πk.

Sumujemy po k korzystając z równania 7.2 :

nλ=1⁢⇒⁢λ=n⁢⇒⁢π^k=nkn.

Estymację parametrów θk odłożymy do dalszej części wykładu.

7.1. Optymalna reguła decyzyjna

Zobaczmy teraz, jak można zdefiniować optymalny klasyfikator w zależności od funkcji straty karzącej za błędne sklasyfikowanie danych.

Definicja 7.1

Funkcja straty to funkcja przyporządkowująca nieujemną wielkość kary poprzez porównanie prawdy (założymy chwilowo, że ją znamy) do podjętej decyzji (wyliczonego estymatora):

L:K︸p⁢r⁢a⁢w⁢d⁢a×K︸d⁢e⁢c⁢y⁢z⁢j⁢a→R+︸k⁢a⁢r⁢a.

Przykładową funkcją kary dla ciągłego y jest L⁢y,y^=y-y^2.

Mając wyestymowaną regułę decyzyjną sensownym jest rozpatrywanie średniej straty dla naszego klasyfikatora:

Definicja 7.2

Ryzyko reguły decyzyjnej dla d:Rp→K:

Ryzyko=średnia strata reguły decyzyjnej ⁢d=
=R(d)=∑y∈K∫RpL(y,d(X))g(y,X)dX=

gdzie g⁢y,X jest gęstością łącznego rozkładu danych. Z twierdzenia o prawdopodobieństwie całkowitym:

=∑k=1K[∫RpL(k,d(X))f(X|k)dX]πk.
Definicja 7.3

Optymalna reguła decyzyjna d* to taka reguła decyzyjna, że

∀d⁢⁢R⁢d*≤R⁢d.
Definicja 7.4

Reguła bayesowska dB⁢X to reguła decyzyjna, która lokalnie dla danego X spełnia warunek:

dB⁢X=argmin1≤l≤K⁢Ey|X⁢L⁢y,l=
=argminl∑k=1KL(k,l)p(k|X)=

ze wzoru Bayesa:

=argminl[∑k=1KL(k,l)πkf(X|k)∑s=1Kπsf(X|s)]=
=argminl[∑k=1KL(k,l)πkf(X|k)].
Stwierdzenie 7.1
R⁢dB=R⁢d*.

Reguła bayesowska jest optymalną regułą decyzyjną.

Dla dowolnej reguły decyzyjnej d zachodzi:

R(d)=∑k=1K[∫RpL(k,d(X))f(X|k)dX]πk=
=∫Rp[∑k=1KL(k,d(X))πkf(X|k)]dX≥
≥∫Rp[min1≤l≤k∑k=1KL(k,l)πkf(X|k)]dX=
=∫Rp[∑k=1KL(k,dB(X))πkf(X|k)]dX=R(dB).
∎

7.2. Wielowymiarowy rozkład normalny

W dalszej części wykładu będziemy zakładać, że fθi mają rozkłady normalne. Dlatego przyjrzyjmy się bliżej własnościom wielowymiarowego rozkładu normalnego i estymacji jego parametrów metodą najwiękzej wiarygodności.

Definicja 7.5

Wektor losowy X=x1,…,xp ma rozkład wielowymiarowy normalny w Rp jeśli ∀ u∈Rp uT⁢X ma rozkład normalny w R. Oznaczmy ten rozkład poprzez N⁢μ,Σ, gdzie μ=E⁢x, Σ=Var⁢X.

Twierdzenie 7.1

Jeżeli X ma rokład normalny w Rp, to ∀ a∈Rk i macierzy A wymiaru k×p, A⁢X+a ma rozkład normalny w Rk.

∀ ⁢u∈Rk⁢⁢uT⁢A⁢X+a=uT⁢a⁢X+uT⁢a.
∎
Wniosek 7.1

Rozkłady brzegowe wielowymiarowego rozkładu normalnego są normalne w odpowiednich podprzestrzeniach Rp.

Twierdzenie 7.2

Fcja charakterystyczna zmiennej losowej X o rozkładzie normalnym w Rp jest postaci:

φX⁢t=ei⁢tT⁢μ-12⁢tT⁢Σ⁢t. (7.3)

Także na odwrót: jeżeli Σ jest symetryczną macierzą dodatnio określoną o wymiarach p×p, to φX określona w równaniu 7.3 jest funkcją charakterystyczną wektora losowego o rozkładzie normalnym w Rp.

Wniosek 7.2

Dowolna macierz symetryczna dodatnio określona o wymiarach p×p jest macierzą kowariancji wektora losowego o rokładzie normalnym w Rp.

Twierdzenie 7.3

Gęstość wielowymiarowego rozkładu normalnego N⁢μ,Σ:

f⁢X=12⁢πp2⁢det⁢Σ︸=Σ12⁢exp⁡-12⁢X-μT⁢Σ-1⁢X-μ.
Twierdzenie 7.4

Jeżeli X ma rozkład normalny w Rp: N⁢μ,Σ, to współrzędne wektora X są niezależne ⇔ Σ jest diagonalna. Dla rozkładu normalnego brak korelacji oznacza niezależność.

Twierdzenie 7.5

Jeżeli X∼N⁢μ,σ2⁢I, C jest macierzą ortonormalną o wymiarach p×p, to:

C⁢X∼N⁢C⁢μ,C⁢σ2⁢I⁢CT=N⁢C⁢μ,σ2⁢C⁢CT︸=I=N⁢C⁢μ,σ2⁢I.

7.2.1. Estymatory największej wiarygodności dla rozkładu normalnego N⁢μ,Σ

Niech X1,…,Xn będą niezależnymi wektorami losowymi z p-wymiarowego rozkładu N⁢μ,Σ. Znajdziemy estymatory dla parametrów μ i Σ. Łączna funkcja wiarygodności dla n wektorów losowych:

L⁢μ,Σ=∏i=1n12⁢πp2⁢Σ12⁢exp⁡-12⁢Xi-μT⁢Σ-1⁢Xi-μ.

Najpierw szukamy estymatora μ^; w tym celu opuszczamy wszystkie wyrazy nie zalezeżące od μ, które by się wyzerowały po policzeniu pochodnej. Dla prostoty obliczeń maksymalizujemy podwojoną logwiarygodność:

2⁢log⁡L⁢μ=2⁢l⁢μ=n⁢μT⁢Σ-1-2⁢n⁢μT⁢Σ-1⁢X¯,

gdzie X¯=1n⁢∑i=1nXi.

Przypomnijmy fakt:

Lemat 7.1

Oznaczmy: a,b – wektory tej samej długości p, A macierz o wymiarach p×p.

∂⁡aT⁢b∂⁡a=∂⁡bT⁢a∂⁡a=b.
∂⁡bT⁢A⁢b∂⁡b=A+AT⁢b︸=2⁢A⁢b⁢ jeśli A symetryczna.

Skorzystajmy z lematu 7.1 żeby obliczyć pochodną logwiarygodności:

1n⁢∂⁡2⁢l⁢μ∂⁡μ=2⁢Σ-1⁢μ-2⁢Σ-1⁢X¯=0,

stąd

μ^=X¯,

czyli estymatorem największej wiarygodności dla średniej rozkładu normalnego jest średnia arytmetyczna obserwacji.

Ponieważ optymalne μ^ nie zależy od Σ, przy obliczaniu Σ^ możemy wstawić X¯ za μ. Maksymalizujemy po Σ wyrażenie:

L(X¯,Σ)∝|Σ|-n2exp[-12∑i=1n(Xi-X¯)TΣ-1(Xi-X¯).]

Symbol ∝ oznacza proporcjonalność, możemy opuścić wszystkie stałe, które nie wpływają na wynik optymalizacji.

Ponieważ (Xi-X¯)TΣ-1(Xi-X¯) jest liczbą, a tr⁢liczba=liczba, oraz tr⁢A⁢B=tr⁢B⁢A, otrzymujemy:

L⁢X¯,Σ∝Σ-n2⁢exp⁡-12⁢∑i=1ntr⁢Xi-X¯⁢Xi-X¯T⁢Σ-1=

ślad macierzy jest funkcją liniową argumentu, więc zachodzi:

=|Σ|-n2exp[-12tr{Xi-X¯⁢Xi-X¯T︸=SΣ-1}]=

pomnóżmy i podzielmy przez Sn2

=S-n2⁢S⁢Σ-1n2⁢exp⁡-12⁢tr⁢S⁢Σ-1.

Ponieważ S-n2 nie zależy od Σ, możemy to wyrażenie opuścić. Podstawmy B=S⁢Σ-1:

L⁢X¯,B∝Bn2⁢exp⁡-12⁢tr⁢B.
Lemat 7.2

Dla macierzy kwadratowej A o wymiarach p×p zachodzi:

det⁢A=∏i=1pλi,
tr⁢A=∑i=1pλi,

gdzie λi to wartości własne macierzy.

Korzystając z lematu 7.2:

Bn2⁢exp⁡-12⁢tr⁢B=∏j=1pλjn2⁢e-12⁢λj

Zmaksymalizujmy to wyrażenie po każdej wartości własnej λj, co sprowadza się do maksymalizacji po λ funkcji:

F⁢λ=λn2⁢e-12⁢λ;
log⁡F⁢λ=n2⁢log⁡λ-12⁢λ;
∂⁡log⁡F⁢λ∂⁡λ=n2⁢λ-12=0;

skąd λ^=λ1^=…=λp^=n.

Macierzą o wszystkich wartościach własnych równych n jest n⁢I:

B=S⁢Σ-1=n⁢I,

skąd:

Σ^=1n⁢S=1n⁢Xi-X¯⁢Xi-X¯T,

czyli estymatorem największej wiarygodności dla macierzy kowariancji rozkładu normalnego jest obciążony estymator próbkowy macierzy kowariancji.

7.3. Klasyfikacja w modelu normalnym

Zrobimy dwa założenia dotyczące rozważanego wcześniej klasyfikatora:

  1. Funkcja straty jest postaci:

    L⁢k,l=1k≠l.
  2. W każdej z grup dane pochodzą z rozkładu normalnego, czyli fθk to gęstość rozkładu normalnego, θk=μk,Σk.

Dal zadanej funkcji straty optymalna (bayesowska) reguła decyzyjna będzie miała postać:

dB(X)=argminl[∑k=1KL(k,l)πkf(X|k)]=argminl[∑k=1K1k≠lπkf(X|k)]=
=argminl[∑k=1Kπkf(X|k)︸nie zależy od wyboru ⁢l-πlf(X|l)]=argmaxl[πlf(X|l)].

Znamy już postać szukanego klasyfikatora, potrzebujemy jeszcze estymatorów dla występujących w nim parametrów. Wiemy jak wyglądają estymatory π^k:

π^k=nkn,⁢nk=∑i=1n1yi=k.

Estymatory największej wiarygodności dla parametrów θk przy założeniu normalności rozkładów w grupach są postaci:

μk^=1nk⁢∑i∈CkXi=∑i=1nXi⁢1yi=k∑i=1n1yi=k;
Σk^=1nk⁢∑i∈CkXi-Xk¯⁢Xi-Xk¯T,

gdzie Xk¯ oznacza wektor średnich obserwacji dla Xi∈Ck.

Dla X niezależnego od próby uczącej: y1,X1,…,yn,Xn estymator reguły decyzyjnej ma postać:

d^⁢X=argmax1≤l≤K⁢πl^⁢fμ^l,Σ^l⁢X.

7.3.1. Kwadratowa (qda) i liniowa (lda) funkcja klasyfikacyjna

W zależności od założeń dotyczących parametrów, możemy otrzymać klasyfikator będący różną funkcją swojego argumentu X: albo kwadratową albo liniową.

Kwadratowa funkcja klasyfikacyjna (qda) nie wymaga dodatkowych założeń o parametrach:

d⁢X=argmaxl⁢πl⁢fμl,Σl=
=argmaxl[πl2⁢πp2⁢Σl12exp{-12(X-μl)TΣl-1(X-μl)}]=

po opuszczeniu wyrażeń niezależnych od l i zlogarytmowaniu:

=argmaxl⁢log⁡πl-12⁢log⁡Σl-12⁢X-μlT⁢Σl-1⁢X-μl,

czyli kwadratowa funkcja argumentu X.

Liniowa funkcja klasyfikacyjna (lda) wymaga założenia:

Σ1=…=ΣK=Σ.

Dzięki niemu mamy podwójny zysk obliczeniowy: o K-1 parametrów mniej do wyestymowania i liniową funkcję optymalizowaną:

d⁢X=argmaxl⁢log⁡πl-12⁢log⁡Σl-12⁢X-μlT⁢Σl-1⁢X-μl=

ponieważ log⁡Σ oraz XT⁢Σ-1⁢X nie zależy od l,

=argmaxl⁢log⁡πl+XT⁢Σ-1⁢μl-12⁢μlT⁢Σ-1⁢μl.

7.4. Metody porównywania klasyfikatorów

Chcemy znaleźć taką metodę porónywania, żeby każdą obserwację spośród y1,X1,…,yn,Xn wykorzystać do uczenia i testu, ale tak żeby testować tylko na tych obserwacjach, które nie były brane pod uwagę przy uczeniu klasyfikatorów.

Kroswalidacja m-krotna (walidacja krzyżowa) polega na podziale danych na m części (popularnymi wyborami są m=5, m=10): m-1 będzie tworzyć próbę uczącą, ostatnia będzie próbą testową. Estymujemy klasyfikatory na próbie uczącej, porównujemy metody na próbie testowej. Powtarzamy procedurę m razy tak, żeby każda z części była próbą testową. Dokładniej:

  1. Permutujemy obserwacje. Jeżeli dane mają jakąś strukturę, na przykład można je podzielić na klasy, permutujemy obserwacje w klasach.

  2. Dzielimy próbę na m części tak, żeby w każdej z grup było po tyle samo obserwacji z każdej klasy.

  3. Uczymy klasyfikatory na próbie uczącej - estymujemy parametry.

  4. Porównujemy metody na próbie testowej (np. poprzez estymację prawdopodobieństwa poprawnej predykcji)

Definicja 7.6

Prawdopodobieństwo poprawnej predykcji to dla danego klasyfikatora Pd⁢X=y. Np. jeżeli funkcja straty wyraża się wzorem L⁢k,l=1k≠l, możemy estymować prawdopodobieństwo poprawnej predykcji dla konkretnej próby treningowej i testowej następująco:

p⁢p⁢p^i=∑i∈próba testowa1d⁢Xi=yi∑i∈próba testowa1,

gdzie d jest klasyfikatorem wyestymowanym na podstawie próby uczącej. Uśrednione p⁢p⁢p^ jest dobrą metodą porównywania klasyfikatorów:

p⁢p⁢p^=∑i=1mp⁢p⁢p^im.

7.5. Przykłady w programie R

Klasyfikacja:

Treść automatycznie generowana z plików źródłowych LaTeXa za pomocą oprogramowania wykorzystującego LaTeXML.

Projekt współfinansowany przez Unię Europejską w ramach Europejskiego Funduszu Społecznego.

Projekt współfinansowany przez Ministerstwo Nauki i Szkolnictwa Wyższego i przez Uniwersytet Warszawski.