2011-08-25 11 views
8

Pracuję nad implementacją funkcji gęstości prawdopodobieństwa wielowymiarowego Gaussa w C++ i utknąłem na tym, jak najlepiej radzić sobie z przypadkami, w których wymiar> 2.Implementacja wielowymiarowej funkcji gęstości prawdopodobieństwa gaussowskiej dla> 2 wymiarów w C++

pDF gaussowski można zapisać jako

multivariate gaussian pdf

gdzie (a) „i a” oznacza transpozycję z «matrycą» utworzonej przez odjęcie średniej z wszystkich elementów x. W tym równaniu k jest liczbą wymiarów, które posiadamy, a sigma reprezentuje macierz kowariancji, która jest macierzą kx k. Wreszcie, | X | oznacza wyznacznik macierzy X.

W przypadku jednowymiarowym wdrożenie pdf jest banalne. Nawet w przypadku biwariatu (k = 2) jest to banalne. Jednakże, gdy przekraczamy dwa wymiary, implementacja jest znacznie trudniejsza.

w dwuwymiarowym przypadku, to mamy

bivariate gaussian pdf

gdzie Rho jest korelacja pomiędzy X i Y, z korelacji równy

correlation between two random variables X and Y

W tym przypadku, może użyć Eigen::Matrix<double, Eigen::Dynamic, Eigen::Dynamic>, aby zaimplementować pierwsze równanie, lub po prostu obliczyć wszystko za pomocą drugiego równania, bez korzystania z uproszczonego interfejsu algebry liniowej Eigen.

Moje myśli na próbę wieloczynnikowej przypadku prawdopodobnie rozpocznie rozszerzając powyższe równania do wielowymiarowej przypadku

multivariate pdf

z

multivariate pdf

moje pytania są następujące:

  1. Czy byłoby właściwe/zaleca się użycie boost::multi_array dla macierzy n-wymiarowej , czy też powinienem spróbować wykorzystać Eigen?
  2. Czy powinienem mieć oddzielne funkcje dla przypadków jednowymiarowych/dwuwymiarowych, , czy też powinienem je po prostu skopiować do wielowymiarowego przypadku, używając boost :: multi_array (lub odpowiedniej alternatywy)?
+0

Oof! Cóż, co próbowaliście do tej pory? : D –

+1

Odpowiednią odpowiedzią jest oczywiście użycie biblioteki numerycznej obsługującej operacje macierzy. Czy nie zapewnia tego uBLAS/LaPack? W każdym razie użycie 'multi_array' (lub czegokolwiek self made) jest bardzo prawdopodobnie * nie * dobrym pomysłem. –

Odpowiedz

1

jestem trochę z mojego elementu tutaj, ale niektóre myśli:

pierwsze, z widokiem programowania, odpowiedź akcji jest „profil”. To znaczy, najpierw zakoduj go w jaśniejszy sposób. Następnie wyprofiluj swoje wykonanie, aby sprawdzić, czy optymalizacja jest opłacalna. IMHO prawdopodobnie lepiej używać biblioteki macierzowej, aby trzymać się bliżej oryginalnej matematyki.

Z widoku matematycznego: Jestem nieco wątpliwy co do formuły, którą podajesz w przypadku wielowymiarowym. To mi nie pasuje. Wyrażenie Z powinno być kwadratową postacią, a Twoja Z nie jest. Chyba że czegoś mi brakuje.

Oto opcja, o której nie wspomniałeś, ale może mieć sens. Szczególnie, jeśli zamierzasz wielokrotnie oceniać plik PDF dla pojedynczej dystrybucji. Zacznij od obliczenia głównej składowej twojego rozkładu. Oznacza to, że własność własna dla Σ. Podstawowe kierunki komponentów są ortogonalne. W podstawowej zasadzie składowej, kowariancje krzyżowe są zerowe, więc plik PDF ma prostą formę. Jeśli chcesz ocenić, zmień podstawę na podstawie danych wejściowych na podstawową podstawę komponentu, a następnie wykonaj prostsze obliczenia w formacie PDF.

Uważa się, że można obliczyć zmianę macierzy bazowej i podstawowych komponentów raz na pierwszym planie, a następnie dokonać tylko pojedynczego mnożenia macierzy (zmiany podstawy) na ocenę, zamiast dwóch macierzowych multiplikacji potrzebnych do oceń (x-μ)' Σ (x-μ) w standardowej bazie.

+0

Gdzie gdzie jest mój ukochany znak TeX? MathOverflow wspiera to ... – Managu

+0

Innymi słowy, przekształć kwadratową formę '(x-μ) 'Σ (x-μ)' w diagonalną formę (a la http://en.wikipedia.org/wiki/Quadratic_form#Real_quadratic_forms) i ocenić na właściwej podstawie. – Managu

0

Zasadniczo wdrożyłem exp -część równania dla trójwymiarowego przypadku w this question. Użyłem najpierw biblioteki komputerowej o nazwie OpenCV. Ale zauważyłem, że interfejs C++ był bardzo wolny. Potem wypróbowałem interfejs C, który był nieco szybszy. W końcu zdecydowałem się zignorować elastyczność i czytelność, więc wdrożyłem to bez żadnych bibliotek i było znacznie szybciej.

Próbuję powiedzieć, że: gdy wydajność jest ważna, należy rozważyć zastosowanie specjalnych przypadków dla najczęściej używanych liczb wymiarów przy jak najmniejszym obciążeniu, jak to możliwe. W przeciwnym razie wybierz łatwość konserwacji ponad prędkość.

Nota prawna: Nie wiem nic na temat szybkości Eigen lub boost::multi_array (prawdopodobnie to jest to, do czego naprawdę zmierza to pytanie?).

Powiązane problemy