Dobroć dopasowania do histogramów 2D

19

Mam dwa zestawy danych reprezentujących parametry gwiazd: obserwowany i modelowany. Za pomocą tych zestawów tworzę tak zwany schemat dwukolorowy (TCD). Próbkę można zobaczyć tutaj:

histogramy

Być obserwowane dane i B dane wydobyte z modelu (nieważne czarne linie, kropki reprezentują dane) Mam tylko jedno A schemat, ale może produkować tyle różnych B schematów, jak chcę, i to, co jest mi potrzebne aby utrzymać ten, który najlepiej pasuje A .

Potrzebuję więc niezawodnego sposobu sprawdzenia poprawności dopasowania schematu B (model) do schematu A (zaobserwowany).

Teraz tworzę histogram 2D ​​lub siatkę (tak to nazywam, może ma bardziej odpowiednią nazwę) dla każdego diagramu, dzieląc obie osie (100 pojemników na każdą) Następnie przechodzę przez każdą komórkę siatki i znajduję absolutną różnicę w liczbie między A i B dla tej konkretnej komórki. Po przejściu przez wszystkie komórki, to suma wartości dla każdej komórki, a więc kończy się z jednej pozytywnej parametru reprezentującego dobroć dopasowania ( ) pomiędzy A i B . Im najbliżej zera, tym lepsze dopasowanie. Zasadniczo tak wygląda ten parametr:solfa

a i j i j b i jsolfa=jajot|zajajot-bjajot|; gdzie jest liczba gwiazd na schemacie A dla danej komórki (oznaczonej przez ) oraz jest numerem B .zajajotjajotbjajot

Tak wyglądają te różnice w liczbie w każdej komórce w utworzonej przeze mnie siatce (zwróć uwagę, że nie używam wartości bezwzględnych w tym obrazie, ale zrobić z nich korzystać przy obliczaniu parametr):( a i j - b i j ) g f(zajajot-bjajot)(zajajot-bjajot)solfa

hess

Problem polega na tym, że doradzono mi, że może to nie być dobry estymator, głównie dlatego, że poza tym, że dopasowanie jest lepsze niż inne, ponieważ parametr jest niższy , naprawdę nie mogę nic więcej powiedzieć.


Ważne :

(dzięki @PeterEllis za podniesienie tego)

1- Punkty w B nie są związane jeden do jednego z punktów A . To ważna rzecz, aby pamiętać przy poszukiwaniu najlepszego dopasowania: liczba punktów A i B to nie koniecznie taka sama, a dobroć dopasowania testu należy również uwzględnić tej rozbieżności i spróbować go zminimalizować.

2- Liczba punktów w każdym zestawie danych B (model wyjściowy), który próbuję dopasować do A, nie jest stała.


W niektórych przypadkach widziałem test Chi-Squared :

O i E ija(Oja-mija)2)/mija ; gdzie to częstotliwość obserwowana (model), a to częstotliwość oczekiwana (obserwacja).Ojamija

ale problem polega na tym: co mam zrobić, jeśli wynosi zero? Jak widać na powyższym obrazku, jeśli utworzę siatkę tych diagramów w tym zakresie, będzie wiele komórek, w których wynosi zero.E imijamija

Przeczytałem również, że niektórzy ludzie zalecają test prawdopodobieństwa logarytmicznego Poissona do zastosowania w przypadkach takich jak ten, w których uczestniczą histogramy. Jeśli to prawda, naprawdę doceniłbym to, gdyby ktoś mógł poinstruować mnie, jak korzystać z tego testu w tym konkretnym przypadku (pamiętaj, że moja wiedza na temat statystyk jest dość słaba, więc proszę, zachowaj to tak proste, jak to możliwe :)

Gabriel
źródło
Czy punkty w B mają relację jeden-do-jednego z punktami w A (np. Każdy jest określoną gwiazdą), czy też jest bardziej abstrakcyjny?
Peter Ellis,
Hi @PeterEllis żadne punkty B nie są związane jeden do jednego z punktów A . W rzeczywistości to już inna ważną rzeczą, aby pamiętać przy poszukiwaniu najlepszego dopasowania: liczba punktów A i Bnie koniecznie równe.
Gabriel
Cześć - ciekawe pytanie, postaram się napisać prawidłową odpowiedź. Czy każda wersja B ma taką samą liczbę punktów, czy też się różnią?
Peter Ellis,
Różnią się również, tylko liczba punktów w A pozostaje stała. Nie masz pojęcia, jak bardzo byś mi pomógł, gdybyś pomógł mi to rozgryźć @PeterEllis.
Gabriel,
To pytanie bardzo przypomina ten temat: stats.stackexchange.com/questions/71036/... Gdzie udzieliłem odpowiedzi.
L Fischman

Odpowiedzi:

14

OK, gruntownie poprawiłem tę odpowiedź. Wydaje mi się, że zamiast binować dane i porównywać liczby w każdym bin, sugestia, którą zakopałem w mojej oryginalnej odpowiedzi na temat dopasowania oszacowania gęstości jądra 2d i porównania ich, jest znacznie lepszym pomysłem. Co więcej, w pakiecie ks Tarna Duonga dla R znajduje się funkcja kde.test (), która robi to łatwo.

Sprawdź dokumentację kde.test, aby uzyskać więcej informacji i argumentów, które możesz poprawić. Ale w zasadzie robi dokładnie to, czego chcesz. Zwracana wartość p to prawdopodobieństwo wygenerowania dwóch zestawów danych, które porównujesz w ramach hipotezy zerowej, że były one generowane z tego samego rozkładu. Więc im wyższa wartość p, tym lepsze dopasowanie między A i B. Zobacz mój przykład poniżej, gdzie łatwo zauważyć, że B1 i A są różne, ale że B2 i A są prawdopodobnie takie same (czyli w jaki sposób zostały wygenerowane) .

# generate some data that at least looks a bit similar
generate <- function(n, displ=1, perturb=1){
    BV <- rnorm(n, 1*displ, 0.4*perturb)
    UB <- -2*displ + BV + exp(rnorm(n,0,.3*perturb))
    data.frame(BV, UB)
}
set.seed(100)
A <- generate(300)
B1 <- generate(500, 0.9, 1.2)
B2 <- generate(100, 1, 1)
AandB <- rbind(A,B1, B2)
AandB$type <- rep(c("A", "B1", "B2"), c(300,500,100))

# plot
p <- ggplot(AandB, aes(x=BV, y=UB)) + facet_grid(~type) + 
    geom_smooth() +     scale_y_reverse() + theme_grey(9)
win.graph(7,3)
p +geom_point(size=.7)

wprowadź opis zdjęcia tutaj

> library(ks)
> kde.test(x1=as.matrix(A), x2=as.matrix(B1))$pvalue
[1] 2.213532e-05
> kde.test(x1=as.matrix(A), x2=as.matrix(B2))$pvalue
[1] 0.5769637

MOJA ORYGINALNA ODPOWIEDŹ PONIŻEJ JEST UTRZYMYWANA TYLKO PONIEWAŻ, ŻE SĄ TERAZ ŁĄCZA Z NIMI POZOSTAŁE, KTÓRE NIE UCZYNIĄ

Po pierwsze, mogą istnieć inne sposoby rozwiązania tego problemu.

Justel i wsp. Przedstawili wielowymiarowe rozszerzenie testu dobroci dopasowania Kołmogorowa-Smirnowa, które moim zdaniem można by zastosować w twoim przypadku, aby sprawdzić, jak dobrze każdy zestaw modelowanych danych pasuje do oryginału. Nie mogłem znaleźć implementacji tego (np. W R), ale może nie wyglądałem wystarczająco mocno.

Alternatywnie może istnieć sposób, aby to zrobić, dopasowując kopułę zarówno do oryginalnych danych, jak i do każdego zestawu modelowanych danych, a następnie porównując te modele. Istnieją implementacje tego podejścia w R i innych miejscach, ale nie jestem specjalnie z nimi zaznajomiony, więc nie próbowałem.

Ale aby odpowiedzieć bezpośrednio na twoje pytanie, przyjęte podejście jest rozsądne. Sugeruje się kilka punktów:

  • Chyba że twój zestaw danych jest większy niż się wydaje, myślę, że siatka 100 x 100 to zbyt wiele pojemników. Intuicyjnie mogę sobie wyobrazić, że wyciąganie wniosków z różnych zestawów danych jest bardziej odmienne niż tylko z powodu precyzji twoich pojemników, co oznacza, że ​​masz wiele pojemników z małą liczbą punktów, nawet gdy gęstość danych jest wysoka. Jednak ostatecznie jest to kwestia osądu. Z pewnością sprawdziłbym twoje wyniki za pomocą różnych podejść do binowania.

  • Po zakończeniu binowania i przekonwertowaniu danych na (w efekcie) tabelę zdarzeń z dwiema kolumnami i liczbą wierszy równą liczbie pojemników (10 000 w twoim przypadku), masz standardowy problem z porównaniem dwóch kolumn liczy. Sprawdziłby się albo test chi-kwadrat, albo dopasowanie jakiegoś modelu Poissona, ale jak mówisz, istnieje niezręczność z powodu dużej liczby zer. Każdy z tych modeli jest zwykle dopasowany, minimalizując sumę kwadratów różnicy, ważoną odwrotnością oczekiwanej liczby zliczeń; gdy zbliża się do zera, może powodować problemy.

Edytuj - reszty tej odpowiedzi nie uważam już za właściwe podejście.

Myślę, że dokładny test Fishera może nie być przydatny ani odpowiedni w tej sytuacji, w której krańcowe sumy wierszy w tabeli krzyżowej nie są ustalone. Daje wiarygodną odpowiedź, ale trudno mi pogodzić jej użycie z oryginalną pochodną z projektu eksperymentalnego. Zostawiam tutaj oryginalną odpowiedź, więc komentarze i pytania uzupełniające mają sens. Ponadto nadal może istnieć sposób odpowiedzi na pożądane podejście OP polegające na binowaniu danych i porównywaniu pojemników za pomocą jakiegoś testu opartego na średnich różnicach bezwzględnych lub kwadratowych. W takim podejściu nadal wykorzystano by , o której mowa poniżej, i sprawdzono niezależność, tj. Szukając wyniku, w którym kolumna A miałaby takie same proporcje jak kolumna B.nsol×2)

Podejrzewam, że rozwiązaniem powyższego problemu byłoby użycie dokładnego testu Fishera , stosując go do , gdzie jest całkowitą liczbą pojemników. Chociaż pełne obliczenia prawdopodobnie nie będą praktyczne ze względu na liczbę wierszy w tabeli, można uzyskać dobre oszacowania wartości p za pomocą symulacji Monte Carlo (implementacja R testu Fishera daje tę opcję jako opcję dla tabel, które są większy niż 2 x 2 i podejrzewam, że robią to inne pakiety). Te wartości p są prawdopodobieństwem, że drugi zestaw danych (z jednego z twoich modeli) ma taki sam rozkład w pojemnikach jak oryginał. Dlatego im wyższa wartość p, tym lepsze dopasowanie. nsol×2)nsol

Symulowałem niektóre dane, aby wyglądać trochę jak twoje, i odkryłem, że to podejście było dość skuteczne w identyfikacji, które z moich zestawów danych „B” zostały wygenerowane z tego samego procesu co „A”, a które nieco się różniły. Z pewnością bardziej skuteczny niż gołym okiem.

  • Przy takim podejściu do testowania niezależności zmiennych w tabeli kontyngencji nie ma znaczenia, że ​​liczba punktów w A jest różna od liczby w B (chociaż zauważ, że jest tonsol×2)problem, jeśli użyjesz tylko sumy różnic bezwzględnych lub różnic kwadratowych, jak pierwotnie proponujesz). Jednak nie ma znaczenia, że ​​każda z twoich wersji B ma inną liczbę punktów. Zasadniczo większe zestawy danych B będą miały tendencję do zwracania niższych wartości p. Mogę wymyślić kilka możliwych rozwiązań tego problemu. 1. Mógłbyś zredukować wszystkie swoje zestawy B danych do tego samego rozmiaru (rozmiar najmniejszego ze swoich zestawów B), pobierając losową próbkę tego rozmiaru ze wszystkich zestawów B, które są większe niż ten rozmiar. 2. Możesz najpierw dopasować dwuwymiarowe oszacowanie gęstości jądra do każdego ze swoich zestawów B, a następnie zasymulować dane z tego oszacowania o równych rozmiarach. 3. możesz użyć jakiejś symulacji, aby obliczyć stosunek wartości p do wielkości i użyć tego do „poprawienia” wartości p otrzymane z powyższej procedury, więc są one porównywalne. Prawdopodobnie istnieją też inne alternatywy. To, co zrobisz, będzie zależeć od tego, jak wygenerowano dane B, jak różne są rozmiary itp.

Mam nadzieję, że to pomaga.

Peter Ellis
źródło
Dokonałem kilku drobnych poprawek literówek; Mam nadzieję, że nie masz nic przeciwko. Może istnieć sposób formatowania rzeczy, aby nieco lepiej podkreślić główne pomysły, szczególnie w ostatnim punkcie. Ale ja też nie chciałem być nadgorliwy. Twoje zdrowie. :)
kardynał
nie ma problemu. Zmagałem się z dobrym sposobem na sformatowanie mojego ostatniego punktu wypunktowania - chciałem hierarchii numerowanej listy pod punktem wypunktowania. Ale nie mogłem znaleźć, jak to zrobić.
Peter Ellis,
Tak, ja też przez chwilę się tym bawiłem, bo tak właśnie zamierzałeś. Nie mogłem szybko wymyślić, jak to zrobić i bardzo wahałem się, czy wprowadzić hurtowe zmiany w układzie, więc pomyślałem, że po prostu skomentuję. :)
kardynał
1
Polecam książkę Wilcoxa jako ogólny tekst statystyk, który używa R (zobacz moją odpowiedź stats.stackexchange.com/questions/25632/... ). Chociaż nie ukrywa dokładnego tekstu Fishera, jest wystarczająco dużo szczegółów na temat konkretnego tekstu w Internecie, co będzie bardziej sensowne, jeśli masz doświadczenie w tej książce na podobnych testach. Może być lepszy tekst na temat tego rodzaju problemu z dobrością dopasowania, ale myślę, że książka Wilcoxa jest świetna jako ogólne wprowadzenie, które szybko przenosi cię na wyższy poziom.
Peter Ellis,
1
Łał. Odpowiedziałeś do cholery z tego. Gdyby istniała „najlepsza wymiana stosów”, byłoby to w niej.
Colin K