Lab 9.pdf

(519 KB) Pobierz
Laboratorium 9
Wymagania teoretyczne
Metoda regresji logistycznej
Ocena jakości modelu – miary i sposoby ich estymacji
Na tym laboratorium zwrócimy uwagę na ocenę jakości modelu regresji logistycznej.
Pewien ogląd daje miara nazywana odchyleniem (deviance).
1
Odzwierciedla ona stopień
dopasowania modelu do danych. Kolejną miarą jest błąd klasyfikacji (wzór 6 z wykładu). W
celu jego oszacowania musimy się nauczyć najpierw klasyfikować obiekty w programie R.
1. Klasyfikacja – obliczanie błędu
Przypomnijmy sobie, że po zastosowaniu funkcji
glm
mamy dostęp do prawdopodobieństw
a
posteriori.
Duże wartości takiego prawdopodobieństwa przemawia za klasyfikacją obiektu do
klasy kodowanej jedynką. Można przyjąć umowę, że dokonujemy tej klasyfikacji jeśli
prawdopodobieństwo
a posteriori
jest większe od �½.
2
Jak to zrobić w R?
Dostęp do omawianych prawdopodobieństw jest następujący
Model$fitted.values
Zobaczmy co będzie wynikiem następującego warunku
Model$fitted.values > 0.5
# jest nim wektor wartości logicznych
Gdyby teraz TRUE zamienić na jedynki, a FALSE na zera mielibyśmy wektor z
klasyfikacjami, tzn. i-ta współrzędna tego wektora mówiłyby do jakiej klasy (,,zer” czy
,,jedynek”
3
) należy i-ty obiekt w tabeli danych. Może czytelnik myśli już jak to zrobić za
pomocą pętli lub warunku? OK, ale ponownie przekonamy się, że R jest ,,sprytny” czy
inaczej mówiąc – wygodny w obsłudze. Wystarczy komenda:
as.numeric(Model$fitted.values > 0.5)
Wynik ten zapiszmy jeszcze do obiektu
pred = as.numeric(Model$fitted.values > 0.5)
OK, wyświetlmy zatem klasyfikacje dla 10 pierwszych obiektów
pred[1:10]
A teraz wyświtlmy ich prawdziwą przynależność do klas (czyli wartość zmiennej Bankrupt z
tabeli danych)
dane[1:10,6]
Jak widać na ogół klasy się zgadzają, ale pierwszy obiekt (pierwsza firma) została błędnie
sklasyfikowana. Model popełnia tu błąd. Policzymy teraz ile jest tych błędów w stosunku do
liczby obiektów.
Zobaczmy co będzie wynikiem warunku
pred!=dane[,6]
# UWAGA: symbol != oznacza różny
1
Jej wartość jest wyświetlana komendą
summary
jak widzieliśmy na poprzednim laboratorium (chodzi o
Residual deviance).
2
Zob. też komentarz pod wzorem 2 z wykładu.
3
W naszym przykładzie nie bankrutów lub bankrutów.
Zliczenie wartości logicznych TRUE da nam liczbę błędnych klasyfikacji. Oczywiście nie
chcemy się ograniczać do analizy wzrokowej
4
. Istnieje sprytne rozwiązanie:
sum(pred!=dane[,6])
# UWAGA: argumentem funkcji
sum
jest tu wektor wartości
# logicznych, które będą przez nią potraktowane jak liczbowe zera
# i jedynki
Ostatecznie błąd klasyfikacji jest równy:
sum(pred!=dane[,6])/50
Możemy jeszcze automatycznie zliczać liczbę obiektów w zbiorze danych.
nrow(dane)
# zlicza liczbę wierszy
Zatem
sum(pred!=dane[,6]) / nrow(dane)
# błąd klasyfikacji
2. Błąd przewidywań – losowy podział zbioru danych
Pamiętajmy o głównym celu budowania modelu regresji logistycznej. Jest nim
przewidywanie zdarzeń w przyszłości, a więc klasyfikacja obiektów, dla których nie będzie
znana wartość zmiennej objaśnianej. Mamy na to prosty sposób – omówiony pod wzorem (6)
z wykładu.
Dokonamy zatem losowego podziału zbioru danych na zbiór uczący 2/3 jego liczebności oraz
testowy 1/3 jego liczebności. Wykorzystamy funkcję
sample,
która losuje określoną liczbę
liczb (ze zwracaniem lub bez zwracania) z wyspecyfikowanego zakresu.
Np.
sample(1:50,5)
# wylosuje 5 liczb z zakresu od 1 do 50
Wyjaśnij w ramach ćwiczenia:
co zrobi poniższy kod
5
oraz wyświetl wynik każdej
komendy w trakcie swego domowego ćwiczenia.
N = nrow(dane)
nr = sample(1:N, N/3)
zb.test = dane[nr, ]
Nt = nrow(zb.test)
zb.ucz = dane[-nr, ]
3. Błąd przewidywań – funkcja
predict
Mamy przygotowany podział losowy oryginalnego zbioru danych. Teraz należy zbudować
model wykorzystując tylko zbiór uczący:
M.ucz = glm(Bankrupt~. , family = binomial, data = zb.ucz)
Dokonamy teraz za pomocą tego modelu klasyfikacji obiektów ze zbioru testowego.
Nie możemy już wykorzystać dostępu M.ucz$fitted.values ! Dlaczego ?
6
4
5
Tu mamy 50 obiektów a w zbiorze
spam
kilka tysięcy !
Odpowiedź na końcu dokumentu :)
Prawdopodobieństwa a posteriori dla obiektów zbioru testowego obliczymy wykorzystując
funkcję
predict:
pred = predict (MODEL, NOWE_DANE, TYP_WYNIKU )
MODEL: określenie modelu, który ma być wykorzystany do przewidywań (wartości
zmiennej objaśnianej)
NOWE_DANE: określenie zbioru obiektów, dla którego chcemy przewidywać (wartości
zmiennej objaśnianej)
7
.
TYP_WYNIKU: jeżeli chcemy prawdopodobieństwa a posteriori to należy wybrać opcję
”response”
Kontynuując nasze ćwiczenie piszemy komendę:
pred = predict (object = M.ucz, newdata = zb.test, type =
"
response
"
)
round(pred,3)
# wyświetlenie z dokładnością do 3 miejsca po przecinku
UWAGA: zauważmy, że obiekt (spółka) o nr 39 budzi pewne wątpliwości. Na podstawie
prawdopodobieństwa
a posteriori
0,502 będzie sklasyfikowana do bankrutów, ale to niewiele
się różni od wartości granicznej 0,5! Można łatwo sobie wyświetlić dane dot. tej spółki
dane[39,]
i zobaczyć, że w rzeczywistości nie zbankrutowała.
Pozostaje obliczyć błąd klasyfikacji. Czy potrafisz? Sprawdź się:
………………………………………………………………………………………….
………………………………………………………………………………………….
Rozwiązanie znajdziesz na końcu dokumentu!
Zadanie domowe
Zbudować model regresji logistycznej dla zbioru
spam.csv.
Oszacować jego błąd klasyfikacji
na zbiorze testowym.
Zadanie domowe dla wytrwałych !
Obliczyć błąd klasyfikacji szacując go na zbiorze testowym 30 razy w pętli (oczywiście za
każdym razem dokonując losowania na zbiór testowy i uczący). Ostatecznie wyliczyć średnią
tych błędów testowych i błąd standardowy.
6
7
Odpowiedź na końcu dokumentu :)
UWAGA: gdybyśmy wstawili tu ten sam zbiór, na którym był budowany model (czyli w naszym ćwiczeniu
zbiór
zb.ucz),
to otrzymalibyśmy te same wyniki co w M.ucz$fitted.values
Odpowiedzi do ćwiczeń
Wyjaśnij
N = nrow(dane)
nr = sample(1:N, N/3)
zb.test = dane[nr, ]
Nt = nrow(zb.test)
zb.ucz = dane[-nr, ]
# zliczenie liczby obiektów w całym zbiorze
# wylosowanie 1/3 numerów z zakresu od 1 do N (u nas N=50)
# przypisanie obiektów o wylosowanych numerach do nowego
# zbioru, o nazwie zb.test, który będziemy traktowali jako testowy
# zliczenie liczby obiektów w zbiorze testowym i zapamiętanie
# pod nazwą Nt
# przypisanie pozostałych obiektów do zbioru o nazwie zb.ucz,
# który będzie traktowany jako zbiór uczący
Nie możemy już wykorzystać dostępu M.ucz$fitted.values ! Dlaczego?
Gdyż tam mamy prawdopodobieństwa a posteriori obiektów ze zbioru uczącego, a nie
testowego.

Plik z chomika:

Inne pliki z tego folderu:

Inne foldery tego chomika:

Zgłoś jeśli naruszono regulamin