Temat: Minimalna korelacja idealnych predyktorów...

Witam

Dziś lekko prowokacyjnie...

Jakoś tak cicho na tym forum, więc pomyślałem, że może uda się sprowokować jakąś ciekawą dyskusję...

Proponuję rozważenie następującego problemu:

1. mamy dwie zmienne ciągłe A i B, oraz jedną zmienną objaśnianą T dyskretną o wartościach 0 i 1

2. każda z naszych 2 zmiennych jest idealnym predyktorem, czyli umożliwia rozdzielenie wartości targetu, np.:

istnieje c takie, że:
A>=c ---> T=0
A< c ---> T=1

3. PYTANIE: jaka jest minimalna korelacja (liniowa lub rangowa) zmiennych A i B?

Na koniec parę słów o tym skąd takie pomysły na podobne rozważania... W wielu metodach wskazane jest używanie zmiennych o możliwie małej korelacji. Powodów jest wiele, chociażby stabilność numeryczna. Ale jednocześnie chcemy mieć zmienne jak najlepiej opisujące badane zjawisko... Te dwa cele mogą stać ze sobą w sprzeczności... pytnie jak bardzo...?

pozdrawiam
Michał M.

Temat: Minimalna korelacja idealnych predyktorów...

Bardzo fajny temat :)

1. Faktycznie, czy pójdziemy drogą MNK, czy ML, w przypadku modelu liniowego mamy te same estymatory i podobne problemy.

Jeśli bierzemy na tapetę metody oparte o "dzierganie" po macierzy ko(rela/warian)cji (gramian), czyli np. MNK, to - zależnie od siły tej współliniowości najpierw polecą oszacowania błędów, a w przypadku układu "obliczeniowo osobliwego" (czy dokładnie osobliwego) oczywiście metoda poleci na 1/wyznacznik (a więc m.in. odwracanie macierzy).

Dodatkowo, o czym nie mówi się na wykładach zbyt często, jeśli macierz danych przypomina macierz Lauchliego (co może "sprawiać wrażenie" współliniowości), to wtedy gramian (XtX) też będzie "obliczeniowo" osobliwy, co paskudnie widać na przykładzie PCA (i LDA), gdzie wartości własne dodatnio określonego gramianu (choć winny być tylko dodatnie) mogą być wskutek błędów num. ujemne.

Czyli sypnie się zarówno liniowy model prawdopodobieństwa (MNK), jak i regresja logistyczna (ML).

2. Jeśli dwie zmienne idealnie dyskryminują zbiór, to czy faktycznie model z nimi obiema będzie tak dobry? Może jednak wyrzucić którąś zmienną, skoro inna idealnie ją "naśladuje"? IMHO obowiązkowo powinno się przejść procedurę selekcji zmiennych, jak chociażby metoda korelacyjna Hellwiga, regresja krokowa (kryteria AIC, BIC), drzewa (CART, CHAID), PCA, LDA ( |Sw| / |Sw + Sb| ) czy Naiwny Bayes i wiele innych. Brzytwa Ochama to jednak dobra rzecz. A może obie zmienne są powiązane znaczeniowo i można z nich utworzyć intepretowalny syntetyczny wskaźnik?

Zawsze warto sprawdzić współczynnik podbicia wariancji (VIF) i współczynnik wskaźnik uwarunkowania dla i-tej zmiennej. Tutaj była taka dyskusja. A zatem patrzeć na te współczynniki, a nie na korelacje jako takie. Na te wskaźniki są "reguły kciuka".

Swoją drogą, ja stosowałem dotychczas regułę kciuka 70%, ale to jest tak na wyczucie. Opinie są różne: http://michaelguth.com/economist/kmemo.htm

3. Co jeszcze można zrobić? Ano można próbować:
- regularyzacji -> regresji grzbietowej (np. w przypadku regularyzacji Tichonowa - regresji SVM),
- uogólnionych równań estymujących (Artykuł z Journal of Statistical Software),

--------------------------------------
Przeprowadziłem tak na szybko prościutką symulację, jak to wygląda w praktyce......

Tu były wyniki mówiące, że już przy korelacji <0.9 dla modelu Y~X1+X2 SE oszacowań wsp. znacznie spada. Zostawiam tylko wykres, jak to wyglądało (R-korelacja predyktorów, A - stosunek SE(b)/b):




Wyciąłem wyniki symulacji, bo poległem na tym zadaniu. Dziś ją powtórzyłem na różnych innych modelach, na danych o innej skali, na różnej ich liczbie - i wszystko się rozjechało, daleko od przewidywalności w kwestii wielkości błędu i współczynnika uwarunkowania (który potrafi znacznie wzrosnąć dla korelacji 0.01). Wpływa na to liczba zmiennych, ich skala, wzajemna korelacja, korelacja ze zmienną zależną, rozrzut danych, postać modelu. Nie takie to proste... Pomyślę znów, jak to ugryźć. Pozostaje symulacja na szerokim zbiorze modeli i obserwacja ich zachowań.

Jedno jest pewne - ze spadkiem współzależności błędy znacznie, ale we zorze na SE mamy dwie korelacje, które się "ścigają".



Adrian Olszewski edytował(a) ten post dnia 10.11.11 o godzinie 16:04

Temat: Minimalna korelacja idealnych predyktorów...

Drugie podejście do współliniowości.

Na początek klasyczna regresja liniowa (logistyczna następnym razem), dwa predyktory. Model Z~X+Y+e.

Przeprowadziłem symulację, w trakcie której wielokrotnie (1000x) losowałem dane dla zmiennych X i Y z dwuwymiarowego rozkładu normalnego (generowany "ręcznie" za pomocą dekomp. Cholesky'ego) o zadanym, stopniowo rosnącym współczynniku korelacji i zadanej liczności. Następnie konstruowałem z nich zmienną zależną, podstawiałem zmienne do procedury liczącej LM, a następnie agregowałem informacje o estymatach parametrów modelu. Dzięki temu mogłem oszacować zachowanie estymat w funkcji korelacji predyktorów.

Interesowały mnie następujące parametry, które w jakiś sposób mówią o stabilności modelu:
1. rozrzut wartości estymat modelu (bety). Rozrzut ten szacuję odchyleniem międzykwartylnym (IQR)
2. wielkość błędu standardowego estymat (SE(bi)). Wielkość tę szacuję medianą
3. procent przypadków, gdy zmienna była istotna statystycznie (na poziomie istotności a=0.05) - SIG(nificance)
4. procent przypadków, gdy estymata miała wartość dodatnią, czyli pośrednio ile razy nastąpiła zmiana znaku (niestabilność) - SIGN
5. wskaźnik uwarunkowania gramianu XtX (dla zmiennych wystandaryzowanych jest to macierz kowariancji, w przeciwnym razie jest do niej proporcjonalny). To jest chyba najlepszy wskaźnik "problemów" - od niego zależy stabilność wielu metod numerycznych.

Sprawdzałem także na bieżąco korelację predyktorów ze zmienną zależną (Z jest "jitterowany", więc lepiej to kontrolować).

Powstały trzy tabelki dla różnej liczności zmiennych: 50 (przeciętny rozmiar próby w statystyce), 100, 500.

W zastosowaniach data miningowych mogą to być wielkości rzędy tysięcy (dziesiątek i setek), więc problem współliniowości będzie coraz mniej znaczący. (Pamiętamy, że takie rozmiary próby nie są jednak dobrą rzeczy, bo pozwalają udowodnić "cokolwiek", nawet związek kichania z burzą na Marsie i prowadzą do przeuczenia). Z drugiej zaś strony - w zagadnieniach DM jest wiele zmiennych, o wiele więcej, niż w klasycznych zagadnieniach statystycznych i nie wyobrażam sobie, by nie była wpierw przeprowadzana ich selekcja, względnie tworzenie wskaźników, agregatów (np. PCA).

Próba 50 elementowa (PS: zapomniałem ustawić rozsądną liczbę miejsc dziesiętnych, a nie chce mi się tego teraz przerabiać)

         Rs       Rxz       Ryz     b0.IQR     b1.IQR     b2.IQR       b0se

1 0.7078631 0.9018512 0.9002969 0.05795857 0.07903395 0.07757124 0.04641147

2 0.8054600 0.9283203 0.9282821 0.05990625 0.08784878 0.08992475 0.04875497

3 0.9010027 0.9562166 0.9548024 0.06082693 0.13296108 0.13989134 0.05018609

4 0.9504632 0.9692699 0.9689196 0.05970506 0.20469154 0.21273968 0.05093004

5 0.9803186 0.9779061 0.9782096 0.06458648 0.32540662 0.31290473 0.05201412

6 0.9950661 0.9803131 0.9801835 0.06410221 0.64140797 0.65389609 0.05271146

7 0.9990324 0.9801946 0.9799811 0.06559838 1.55295606 1.57099801 0.05394606



        b1se       b2se b0sig b1sig b2sig b0sign b1sign b2sign        Kmed

1 0.06338004 0.06334546 0.057 0.993 0.995  0.477  0.999  1.000    6.523409

2 0.08058206 0.08023386 0.059 0.985 0.982  0.496  1.000  0.997   10.185944

3 0.11576263 0.11465734 0.048 0.939 0.918  0.516  1.000  0.995   20.135575

4 0.17140403 0.17317291 0.045 0.834 0.841  0.497  0.986  0.991   40.614847

5 0.26593284 0.26430598 0.048 0.730 0.721  0.531  0.960  0.967  102.242335

6 0.54225371 0.54404725 0.048 0.495 0.520  0.475  0.887  0.900  404.595106

7 1.28148735 1.27838595 0.047 0.265 0.267  0.494  0.792  0.788 2073.245757

Próba 100 elementowa

         Rs       Rxz       Ryz     b0.IQR     b1.IQR     b2.IQR       b0se

1 0.6998292 0.9119868 0.9105295 0.02222196 0.03400405 0.03381368 0.01215638

2 0.8031049 0.9356706 0.9361619 0.02909592 0.04761750 0.04667387 0.01289200

3 0.9020949 0.9629271 0.9621930 0.03308900 0.08405306 0.08280329 0.02352077

4 0.9498030 0.9775420 0.9775663 0.03847896 0.11046114 0.11036531 0.02387422

5 0.9801408 0.9858614 0.9850667 0.03889692 0.17651843 0.18503063 0.02444852

6 0.9950860 0.9886473 0.9894029 0.03696740 0.38064108 0.39138538 0.02597552

7 0.9990096 0.9907751 0.9904346 0.04070566 0.90391555 0.89207716 0.02521315



        b1se       b2se b0sig b1sig b2sig b0sign b1sign b2sign        Kmed

1 0.01776170 0.01764451 0.048 1.000 1.000  0.529  1.000  1.000    6.406659

2 0.02331215 0.02300605 0.057 1.000 1.000  0.505  1.000  1.000    9.978332

3 0.05411416 0.05431975 0.045 0.981 0.988  0.492  1.000  1.000   20.337840

4 0.07748203 0.07774292 0.050 0.908 0.883  0.514  0.998  0.998   39.876157

5 0.12836581 0.12757731 0.063 0.710 0.721  0.526  0.986  0.977  100.402670

6 0.27984689 0.28141369 0.051 0.586 0.593  0.502  0.908  0.918  407.800890

7 0.61236620 0.61099116 0.054 0.473 0.467  0.499  0.847  0.828 2018.026543

Próba 500 elementowa

         Rs       Rxz       Ryz     b0.IQR     b1.IQR     b2.IQR       b0se

1 0.7018872 0.7858608 0.7843768 0.06232952 0.08550884 0.08823685 0.05163115

2 0.8000217 0.8117189 0.8115205 0.06295941 0.11075034 0.11062715 0.05168079

3 0.9005334 0.8385289 0.8390197 0.06683937 0.16257750 0.16172947 0.05176165

4 0.9503726 0.8523748 0.8529823 0.06798753 0.21923415 0.21407655 0.05173124

5 0.9800688 0.8602205 0.8603374 0.06813147 0.34198388 0.33676651 0.05171223

6 0.9949961 0.8642801 0.8639614 0.06835973 0.68308731 0.69333288 0.05170436

7 0.9989998 0.8655995 0.8655017 0.06912292 1.61015252 1.60500192 0.05170868



        b1se       b2se b0sig b1sig b2sig b0sign b1sign b2sign        Kmed

1 0.07214550 0.07206863 0.049 1.000 1.000  0.493  1.000  1.000    6.415322

2 0.08601955 0.08605199 0.038 1.000 1.000  0.494  1.000  1.000    9.812146

3 0.11897471 0.11867102 0.046 1.000 1.000  0.521  1.000  1.000   20.035299

4 0.16599973 0.16621896 0.052 1.000 1.000  0.514  1.000  1.000   40.240888

5 0.26038431 0.26047932 0.042 0.981 0.977  0.521  1.000  1.000  100.336240

6 0.51776654 0.51769726 0.045 0.502 0.475  0.502  0.977  0.967  399.522201

7 1.15799599 1.15822682 0.052 0.139 0.140  0.496  0.801  0.790 1997.955763

Wyniki mówią same za siebie. Niezależnie od rozmiarów próby "magiczną granicą" dla korelacji między predyktorami jest 0,9. Wtedy to następuje nagły skok wskaźnika uwarunkowania gramianu.

Zmienną parametr b0 (wyraz wolny) zostawmy, bo w tym modelu jest nieistotny statystycznie (co pięknie widać z tabelek; X i Y miały niewielką skalę, blisko 0), a nie chciałem szacować bez wyrazu wolnego;

Spójrzmy na rozrzut wartości b1 dla R=0.7 i R=0.99 - widzimy 10 krotny skok. b1 i b2 oraz ich SE mają niemal identyczne wartości (co było do przewidzenia patrząc na postać modelu). Błędy jednak rosną szybko, też niemal 10x. Początkowo wolno, a potem szybciej, zgodnie ze wskaźnikiem uwarunkowania gramianu.

Ze wzrostem R zaczyna zwiększać się odsetek przypadków, gdy dana zmienna przestawała być istotna (pierwsza cały czas trzyma się na a=0.05), co jest oczywiste, skoro rośnie SE. Podobnie dzieje się z liczbą zmian znaku estymaty. Powyższej R=0.98 dzieją się już "cuda".

Widać, że im większa liczebność próby, tym procesy te zachodzą dla coraz wyższego R.

Tak to wygląda na wykresach:







W kolejce czekają symulacje dla różnej skali rozkładów brzegowych (czyli naszych zmiennych), dla 3 zmiennych i dla regresji logistycznej.

Myślę, że to moje R=0.7 było sporo na wyrost, ale za to bezpieczne. Grunt, by nie mierzyć wzrostu w calach i centymetrach :) No cóż, kłaniają się lekcje geometrii: jedna prosta może należeć do nieskończenie wielu płaszczyzn. Zaś gdy dwie proste leżą bardzo blisko siebie, niewielka zmiana geometrii "chmury" punktów, do których są one dopasowane sprawi, że przechodząca przez te proste płaszczyzna zmieni drastycznie swoje nachylenie, równanie i znaki współczynników kierunkowych.

A inna nauka, płynąca z problemu uwarunkowania, to by standaryzować zmienne o dużym rozrzucie skali (choć pozmieniają się wartości własne, trzeba uważać).

--------------------

EDIT:
Symulacja dla 3 zmiennych.

         Rs       Rxz       Ryz       Rwz    b0.IQR    b1.IQR    b2.IQR

1 0.7041814 0.8241896 0.8226277 0.8251333 0.1548157 0.2581060 0.2365834

2 0.8029350 0.8611274 0.8611930 0.8621581 0.1587929 0.2980249 0.2767118

3 0.9000867 0.8976252 0.8975096 0.8976769 0.1507765 0.4115465 0.4298054

4 0.9507164 0.9158223 0.9156049 0.9162122 0.1594110 0.5464605 0.6027881

5 0.9803503 0.9271584 0.9268540 0.9270759 0.1525939 0.8669610 0.9736675

6 0.9950590 0.9316241 0.9314906 0.9314752 0.1517298 1.8905234 1.8966831

7 0.9990233 0.9329371 0.9330406 0.9328843 0.1609568 4.3891092 4.1485540



     b3.IQR      b0se      b1se      b2se      b3se b0sig b1sig b2sig b3sig

1 0.2438041 0.1171134 0.1811349 0.1802130 0.1798125 0.041 1.000 1.000 1.000

2 0.2925279 0.1175289 0.2185929 0.2193366 0.2186055 0.044 0.992 0.997 0.992

3 0.4349527 0.1174768 0.3060802 0.3077589 0.3071780 0.048 0.886 0.887 0.880

4 0.6012223 0.1170671 0.4306806 0.4297945 0.4312728 0.053 0.620 0.627 0.629

5 0.9101361 0.1172758 0.6831696 0.6800903 0.6771571 0.049 0.307 0.324 0.307

6 1.7558299 0.1173959 1.3648906 1.3566076 1.3536303 0.053 0.115 0.113 0.120

7 3.9722461 0.1173109 3.0312434 3.0357926 3.0508419 0.043 0.062 0.066 0.053



  b0sign b1sign b2sign b3sign           K

1  0.502  1.000  1.000  1.000    7.726145

2  0.481  1.000  1.000  1.000   11.611557

3  0.505  0.998  0.999  1.000   22.964347

4  0.515  0.990  0.989  0.992   45.156813

5  0.533  0.935  0.912  0.929  112.953399

6  0.492  0.778  0.730  0.784  452.517930

7  0.494  0.598  0.643  0.635 2252.201331

1. Blisko 8 krotne skoki rozrzutu i SE dla estymat dla R=0.995. Dla R=0.999 to już 20 razy. No ale to praktycznie funkcyjna zależność.
2. Dla R=0.999 wszystkie zmienne przestają być istotne na poz. istn. a=0.05
3. Dla tych samych korelacji estymaty częściej zmieniają znak
4. Nadal graniczną wartością R jest 0.9. Powyżej tej wartości:
a) zaczynają się wyraźnie zmieniać wartości różnych parametrów,
b) wsp. uwarunkowania gramianu (sqrt(λmax/λmin)) przekracza 30 (wartość graniczna z "reguły kciuka")Adrian Olszewski edytował(a) ten post dnia 14.11.11 o godzinie 01:20

Temat: Minimalna korelacja idealnych predyktorów...

Ciekawe symulacje dot. współliniowości i ładna ilustracja reguły kciuka VIF < 4, chociaż wtedy wartość, którą Adrian wskazujesz, powinna być r=0,866, żeby wyprodukować tę łatwą do zapamiętania wartość VIF :) [ 1/(1-(0.866)^2) ~= 4]. Swoją drogą, zawsze traktowałem predyktor jako potencjalnie (ale niekoniecznie faktycznie) problematyczny, gdy VIF>2, co się ładnie przekłada na Twoje r~=0.7

Zwykle obok VIF-ów obserwuję po prostu zachowanie i wartości parametrów i ich błędów (choćby: przed/po włączeniu do modelu potencjalnie problematycznego predyktora). W problemach wielowymiarowych ciekawsze wydają mi się problemy związane ze współliniowością dla zmiennych pomocniczych reprezentujących zmienne kategorialne (np. dwa zestawy zmiennych pomocniczych definiujące dwa ortogonalne zestawy porównań planowanych mogą produkować skorelowane ze sobą wzajemnie zmienne pomocnicze; ich korelację też można minimalizować, ale to trochę inna sprawa). Tutaj wydaje mi się, że można dużo łatwiej przekombinować i uzyskać niestabilny model, choć niekoniecznie VIFy będą wysokie.

A odpowiadając bezpośrednio na pytanie Michała:

- minimalny kwadrat korelacji liniowej (zlikwidujmy sobie znak korelacji ujemnej) dla tego problemu wynosi... dowolnie mała wartość niezerowa!

przykład:
x1....x2....y
ɛ....K-ɛ....0
K-ɛ....ɛ....0
-ɛ...ɛ-K....1
K-1...-ɛ....1

gdzie ɛ = bardzo mała wartość większa od zero
a K = dowolna liczba większa od zero

x1 i x2: predyktory spełniające kryterium idealnej separacji zmiennej y i przy takim układzie danych niemal nieskorelowane ze sobą.

Wystarczy wygenerować sobie mniej więcej taki sprytny i bardzo sztuczny układ danych:





I jest i idealna separacja po obu predyktorach i korelacja między nimi tym bliższa 0, im punkty danych bliższe skrzyżowaniu wartości granicznych mieszczących się w średnich x1 i x2 (standaryzacja).

Z korelacją rangową ta sztuczka się jednak nie uda....

A co do sprzeczności:
Praktyczny problem rozumiem, ale mimo to sprzeczności pomiędzy precyzyjnym opisem badanego zjawiska a wymaganiami modelu statystycznego nie widzę. To model ma się nadawać do opisu zjawiska, nie odwrotnie.Mariusz T. edytował(a) ten post dnia 18.11.11 o godzinie 20:33