Do analizy danych z eksperymentu biofizyki próbuję obecnie dopasować krzywą za pomocą wysoce nieliniowego modelu. Funkcja modelu wygląda następująco:
Tutaj szczególnie duże znaczenie ma wartość .
Wykres dla tej funkcji:
(Zauważ, że funkcja modelu opiera się na dokładnym matematycznym opisie systemu i wydaje się działać bardzo dobrze --- po prostu zautomatyzowane dopasowania są trudne).
Oczywiście funkcja modelu jest problematyczna: strategie dopasowania, które wypróbowałem do tej pory, zawodzą z powodu ostrej asymptoty przy , szczególnie przy zaszumionych danych.
Rozumiem ten problem tutaj: proste dopasowanie najmniejszych kwadratów (grałem z regresją zarówno liniową, jak i nieliniową w MATLAB; głównie Levenberg-Marquardt) jest bardzo wrażliwe na asymptotę pionową, ponieważ małe błędy w x są bardzo wzmocnione .
Czy ktoś mógłby wskazać mi odpowiednią strategię, która mogłaby obejść ten problem?
Mam podstawową wiedzę na temat statystyki, ale wciąż jest dość ograniczona. Chciałbym się uczyć, gdybym tylko wiedział, od czego zacząć :)
Bardzo dziękuję za radę!
Edytuj Błagam o wybaczenie za zapomnienie o wspominaniu błędów. Jedyny znaczący hałas występuje w i jest on addytywny.
Edytuj 2 Dodatkowe informacje o tle tego pytania. Powyższy wykres modeluje zachowanie rozciągania polimeru. Jak zauważył @whuber w komentarzach, potrzebujesz aby uzyskać wykres podobny do powyższego.
Co do tego, jak ludzie dopasowywali tę krzywą do tego momentu: wydaje się, że ludzie zazwyczaj odcinają pionową asymptotę, dopóki nie znajdą dobrego dopasowania. Wybór odcięcia jest jednak nadal arbitralny, co powoduje, że procedura dopasowania jest niewiarygodna i odtwarzalna.
Edycja 3 i 4 Naprawiono wykres.
źródło
Odpowiedzi:
Metody, których użylibyśmy, aby dopasować to ręcznie (to znaczy eksploracyjną analizę danych), mogą zadziałać wyjątkowo dobrze z takimi danymi.
Chciałbym nieco sparametryzować model , aby jego parametry były dodatnie:
Dla danego , załóżmy, że istnieje unikalny rzeczywisty spełniający to równanie; nazwij to lub, dla zwięzłości, gdy zrozumiemy .x f ( y ; a , b ) f ( y ) ( a , b )y x fa( y; a , b ) fa( y) ( a , b )
Obserwujemy zbiór uporządkowanych par gdzie odbiegają od przez niezależne losowe zmienne ze średnimi zerowymi. W tej dyskusji założę, że wszystkie mają wspólną wariancję, ale rozszerzenie tych wyników (przy użyciu ważonych najmniejszych kwadratów) jest możliwe, oczywiste i łatwe do wdrożenia. Oto symulowany przykład takiej kolekcji wartości wartości , i wspólnej wariancji .x i f ( y i ; a , b ) 100 a = 0,0001 b = 0,1 σ 2 = 4( xja, yja) xja fa( yja; a , b ) 100 a = 0,0001 b = 0,1 σ2)= 4
Jest to (celowo) trudny przykład, co można docenić przez niefizyczne (ujemne) wartości i ich niezwykły rozkład (który zwykle wynosi jednostki poziome , ale może wynosić do lub na osi ). Jeśli uda nam się uzyskać rozsądne dopasowanie do tych danych, które są bliskie oszacowania zastosowanych wartości , i , naprawdę zrobimy to dobrze.± 2 5 6 x a b σ 2x ± 2 5 6 x za b σ2)
Wyposażenie eksploracyjne jest iteracyjne. Każdy etap składa się z dwóch etapów: oszacowanie (w oparciu o dane i wcześniejszych oszacowań i z i , z których poprzednie wartości przewidywanych można uzyskać za ), a następnie oszacuj . Ponieważ błędy są w x , dopasowania szacują na podstawie , a nie na odwrót. Aby najpierw uporządkować błędy w , gdy jest wystarczająco duże,b a b x i x i b x I ( r i ) x xza za^ b^ za b x^ja xja b xja ( yja) x x
W związku z tym, możemy aktualizować montując ten model z najmniejszych kwadratów (zawiadomienie to ma tylko jeden parametr - zboczu, --and nie przechwycenia) oraz biorąc odwrotność współczynnika jako zaktualizowanego oszacowania .za^ za za
Następnie, gdy jest wystarczająco małe, dominuje odwrotny kwadratowy wyraz i stwierdzamy (ponownie w celu pierwszego uporządkowania błędów), żex
Po ponownym użyciem metody najmniejszych kwadratów (z tylko zbocze termin ) otrzymujemy zaktualizowanej estymaty B przez pierwiastek kwadratowy z zamontowanym zbocza.b b^
Aby zrozumieć, dlaczego to działa, surowy rozpoznawcza przybliżenie tego ataku można uzyskać przez wykreślenie przed 1 / r 2 i dla mniejszych X í . Jeszcze lepiej, ponieważ x i mierzy z błędem, a Y i zmienia się monotonicznie z x i , powinniśmy skupić się na danych z większych wartości 1 / r 2 í . Oto przykład z naszego symulowanego zbiorze pokazując największą połowę y Ixja 1 / y2)ja xja xja yja xja 1 / y2)ja yja na czerwono, najmniejsza połowa na niebiesko, a linia przechodząca przez początek pasuje do czerwonych punktów.
Punkty w przybliżeniu się wyrównują, chociaż przy małych wartościach i y występuje trochę krzywizny . (Zwróć uwagę na wybór osi: ponieważ x jest pomiarem, zwykle wykreśla się go na osi pionowej .) Skupiając dopasowanie na czerwonych punktach, w których krzywizna powinna być minimalna, powinniśmy uzyskać rozsądne oszacowanie b . Podana w tytule wartość 0,096 jest pierwiastkiem kwadratowym nachylenia tej linii: to tylko 4 % mniej niż prawdziwa wartość!x y x b 0,096 4
W tym momencie prognozowane wartości można zaktualizować za pomocą
Iteruj, aż oszacowania ustabilizują się (co nie jest gwarantowane) lub przejdą przez małe przedziały wartości (których nadal nie można zagwarantować).
Okazuje się, że jest trudne do oszacowania, chyba że mamy dobry zestaw bardzo dużych wartości x , ale to b - co determinuje pionową asymptotę na oryginalnym wykresie (w pytaniu) i jest w centrum pytania - można przypiąć dość dokładnie, pod warunkiem, że w pionowej asymptocie znajdują się pewne dane. W naszym przykładzie systemem iteracje są zbieżne dla a = 0.000196 (który jest prawie dwukrotnie prawidłowa wartość 0,0001 ) oraz b = 0,1073 (który jest zbliżony do właściwej wartości 0,1za x b za^= 0,000196 0,0001 b^= 0,1073 0,1 ). Ten wykres pokazuje dane jeszcze raz, na które nakładają się (a) prawdziwa krzywa w kolorze szarym (przerywanym) i (b) oszacowana krzywa w kolorze czerwonym (ciągłym):
To dopasowanie jest tak dobre, że trudno jest odróżnić prawdziwą krzywą od dopasowanej krzywej: pokrywają się prawie wszędzie. Nawiasem mówiąc, szacowana wariancja błędu jest bardzo zbliżona do prawdziwej wartości 4 .3,73 4
Z tym podejściem wiążą się pewne problemy:
Szacunki są tendencyjne. Odchylenie staje się widoczne, gdy zestaw danych jest mały i stosunkowo niewiele wartości znajduje się w pobliżu osi x. Dopasowanie jest systematycznie trochę niskie.
Procedura szacowania wymaga metody powiedzieć „duże” z „małych” wartości . Mógłbym zaproponować eksploracyjne sposoby identyfikacji optymalnych definicji, ale w praktyce można je pozostawić jako stałe „dostrajające” i zmieniać je, aby sprawdzić czułość wyników. Ustawiłem je arbitralnie, dzieląc dane na trzy równe grupy zgodnie z wartością yi i używając dwóch zewnętrznych grup.yja yja
Procedura nie będzie działać dla wszystkich możliwych kombinacji i B lub wszystkich możliwych zakresów danych. Powinien on jednak działać dobrze, ilekroć w zestawie danych jest wystarczająca liczba krzywych, aby odzwierciedlić obie asymptoty: pionową na jednym końcu i ukośną na drugim końcu.za b
Kod
W Mathematica napisano poniżej .
Zastosuj to do danych (podanych przez wektory równoległea = b = 0
x
iy
uformowanych w macierz dwukolumnowądata = {x,y}
) aż do zbieżności, zaczynając od oszacowań :źródło
Zobacz ważne pytania opublikowane przez @probabilityislogic
-
Edytuj, aby uwzględnić dodatkowe informacje:
Mamy teraz, że błędy są w x i addytywnie. Nadal nie wiemy, czy wariancja jest stała w tej skali.
Nie jestem pewien, czy to poprawi rzeczy! Sądzę, że istnieją metody tego typu, ale tak naprawdę to nie jest mój obszar.
Wspomniałem w komentarzach, że możesz chcieć spojrzeć na regresję odwrotną, ale szczególna forma twojej funkcji może uniemożliwić daleko do tego.
Możesz nawet utknąć w próbach dość solidnych metod na błędy w x w tej liniowej formie.
-
źródło
Po kilku tygodniach eksperymentów w tym konkretnym przypadku najlepsza wydaje się inna technika: dopasowanie Total Least Squares . Jest to wariant zwykłego (nieliniowego) dopasowania najmniejszych kwadratów, ale zamiast mierzyć błędy dopasowania wzdłuż jednej z osi (co powoduje problemy w wysoce nieliniowych przypadkach, takich jak ta), bierze pod uwagę obie osie.
Istnieje mnóstwo artykułów, samouczków i książek na ten temat, chociaż nieliniowy przypadek jest bardziej nieuchwytny. Jest nawet dostępny kod MATLAB .
źródło