Bezstronna ocena macierzy kowariancji dla wielokrotnie cenzurowanych danych

22

Analizy chemiczne próbek środowiskowych są często cenzurowane poniżej limitów sprawozdawczych lub różnych limitów wykrywalności / ilościowych. Te ostatnie mogą się różnić, zwykle proporcjonalnie do wartości innych zmiennych. Na przykład, próbka o wysokim stężeniu jednego związku może wymagać rozcieńczenia do analizy, co spowoduje proporcjonalne zawyżenie limitów cenzury dla wszystkich innych związków analizowanych jednocześnie w tej próbce. Jako inny przykład, czasami obecność związku może zmienić odpowiedź testu na inne związki („interferencja matrycy”); kiedy laboratorium wykryje to, odpowiednio zwiększy swoje limity raportowania.

Szukam praktycznego sposobu oszacowania całej macierzy wariancji-kowariancji dla takich zestawów danych, zwłaszcza gdy wiele związków doświadcza ponad 50% cenzury, co często ma miejsce. Konwencjonalny model dystrybucji polega na tym, że logarytmy (prawdziwych) stężeń są rozkładane wielonormalnie, co wydaje się dobrze pasować w praktyce, więc przydatne byłoby rozwiązanie tej sytuacji.

(Przez „praktyczny” rozumiem metodę, którą można niezawodnie zakodować w co najmniej jednym ogólnie dostępnym środowisku oprogramowania, takim jak R, Python, SAS itp., W sposób, który wykonuje się wystarczająco szybko, aby obsługiwać iteracyjne ponowne obliczenia, takie jak wielokrotne przypisywanie, i który jest dość stabilny [dlatego niechętnie badam implementację BŁĘDU, chociaż rozwiązania bayesowskie są ogólnie mile widziane])

Z góry dziękuję za przemyślenia na ten temat.

Whuber
źródło
Właśnie dlatego dobrze rozumiem problem cenzury: kiedy rozcieńczasz próbkę, stężenie związku spada tak nisko, że przyrząd testowy może nie wykryć jego obecności. Czy to dokładne ponowne sformułowanie problemu cenzury?
Tak, to prawda: rozcieńczenie o współczynnik D zwiększa również wszystkie granice wykrywalności o współczynnik D. (Problem interferencji macierzy jest trudniejszy do oszacowania, a ogólna sytuacja jest niezwykle złożona. Aby uprościć to, konwencjonalny model jest taki, że zestaw testów na jednej próbce daje wektor (x [1], ..., x [k ]), w którym x [I] albo są liczbami rzeczywistymi lub są przedziały liczb rzeczywistych, zwykle z lewego punktu końcowego po -nieskończoności; odstępie identyfikuje zestaw, w której rzeczywista wartość jest wielorakie).
whuber
Dlaczego limity wykrywania wzrosłyby? Czy nie są one cechą przyrządu testowego, a nie testowanej próbki?
Jako przykład załóżmy, że granica wykrywalności przyrządu wynosi 1 mikrogram na litr (ug / l). Próbka jest rozcieńczana 10: 1 (z wielką precyzją, więc nie martwimy się tutaj błędem), a instrument odczytuje „<1”; to jest niewykrywalne dla rozcieńczonej próbki. Laboratorium wnioskuje, że stężenie w próbce jest mniejsze niż 10 * 1 = 10 ug / L i podaje je jako takie; to znaczy jako „<10”.
whuber
1
@amoeba Widzę, że powinienem był wyjaśnić te rzeczy w samym pytaniu. Odpowiedzi są następujące: PCA; wymiarowość będzie zmieniać się od 3 do kilkuset; rozmiary próbek zawsze znacznie przekraczają wymiarowość, ale wskaźniki cenzury mogą być bardzo wysokie (wymagana jest zdolność do obsługi do 50%, a pożądane jest do 95%).
whuber

Odpowiedzi:

3

Nie w pełni zinternalizowałem problem interferencji macierzy, ale oto jedno podejście. Pozwolić:

Y jest wektorem, który reprezentuje stężenie wszystkich docelowych związków w nierozcieńczonej próbce.

Z oznacza odpowiedni wektor w rozcieńczonej próbce.

dre jest współczynnikiem rozcieńczenia, tzn. próbka jest rozcieńczana : 1.re

Nasz model to:

YN.(μ,Σ)

Z=Yre+ϵ

gdzie reprezentuje błąd wynikający z błędów rozcieńczania.ϵN.(0,σ2) ja)

Wynika stąd, że:

ZN.(μre,Σ+σ2) ja)

Oznacz powyższy rozkład przez .fZfaZ(.)

Niech będzie obserwowanymi stężeniami, a reprezentuje próg przyrządu testowego, poniżej którego nie może wykryć związku. Następnie dla związku mamy:τOτjath

Oja=Zjaja(Zja>τ)+0ja(Zjaτ)

Bez utraty ogólności niech pierwsze związków będzie takie, aby były poniżej progu. Następnie funkcję prawdopodobieństwa można zapisać jako:k

L.(O1,...Ok,Ok+1,...On|-)=[ja=1ja=kP.r(Zjaτ)][ja=k+1ja=nfa(Oja|-)]

gdzie

fa(Oja|-)=jotjafaZ(Oja|-)ja(Oja>τ)

Oszacowanie polega zatem na wykorzystaniu albo maksymalnego prawdopodobieństwa, albo pomysłów bayesowskich. Nie jestem pewien, na ile powyższe jest wykonalne, ale mam nadzieję, że dostarczy ci kilku pomysłów.


źródło
Dziękuję bardzo za tę myśl. Rzeczywiście jest to standardowe i dobrze udokumentowane podejście do wielokrotnej cenzury. Jedna trudność polega na jej trudności: te całki są niezwykle trudne do obliczenia. Czai się tu również problem modelowania: wartość d jest zwykle dodatnio skorelowana z Y , jak sugeruje pierwszy akapit mojego opisu.
whuber
2

Inną bardziej wydajną obliczeniowo opcją byłoby dopasowanie macierzy kowariancji poprzez dopasowanie momentu za pomocą modelu, który został nazwany „dychomizowanym gaussowskim”, tak naprawdę tylko modelem kopuły Gaussa.

Niedawny artykuł Macke i in. 2010 opisuje procedurę zamkniętej formy dopasowania tego modelu, która obejmuje tylko (ocenzurowaną) empiryczną macierz kowariancji i obliczenie niektórych dwuwymiarowych normalnych prawdopodobieństw. Ta sama grupa (laboratorium Bethge'a z MPI Tuebingen) opisała również hybrydowe dyskretne / ciągłe modele gaussowskie, które prawdopodobnie są tutaj potrzebne (tj. Ponieważ Gaussowskie RV nie są w pełni „dychotomizowane” - tylko te poniżej progu).

Krytycznie nie jest to oszacowanie ML i obawiam się, że nie wiem, jakie są jego właściwości uprzedzające.

jpillow
źródło
@jp Dziękuję: przyjrzę się temu. (Może to zająć trochę czasu ...)
whuber
1

Ile związków jest w twojej próbce? (Lub, jak duża jest omawiana macierz kowariancji?).

Alan Genz ma bardzo ładny kod w różnych językach (R, Matlab, Fortran; patrz tutaj ) do obliczania całek wielowymiarowych normalnych gęstości w hiperprostokątach (tj. Rodzajów całek, których potrzebujesz do oceny prawdopodobieństwa, jak zauważono przez użytkownik 28).

Użyłem tych funkcji („ADAPT” i „QSIMVN”) dla całek o wielkości do około 10-12 wymiarów, a kilka funkcji na tej stronie reklamuje całki (i powiązane pochodne, których możesz potrzebować) w przypadku problemów do wymiaru 100. Nie nie wiem, czy jest to wystarczająca liczba wymiarów dla twoich celów, ale jeśli tak, to prawdopodobnie pozwoli ci znaleźć szacunki maksymalnego prawdopodobieństwa na podstawie wzrostu gradientu.

jpillow
źródło
Och, przepraszam - jestem tu nowy i nie zauważyłem, jak dawno temu został opublikowany - prawdopodobnie za późno, aby być bardzo pomocnym!
jpillow
@jp Jest to wciąż ważny problem, więc upływ czasu między pytaniem a odpowiedzią nie ma większego znaczenia. Dziękuję za odpowiedź!
whuber