Dopasowanie rozkładu log-normalnego w R vs. SciPy

Dopasowałem logarytmiczny model za pomocą R z zestawem danych. Wynikowymi parametrami były:

meanlog = 4.2991610 
sdlog = 0.5511349

Chciałbym przenieść ten model do Scipy, z którego nigdy wcześniej nie korzystałem. Korzystając z Scipy, udało mi się uzyskać kształt i skalę 1 oraz 3.1626716539637488e + 90 - bardzo różne liczby. Próbowałem również użyć exp z meanlog i sdlog, ale nadal otrzymuję dziwny wykres.

Przeczytałem każdy dokument, który mogę na scipy i nadal jestem zdezorientowany co do znaczenia parametrów kształtu i skali w tym przypadku. Czy sensownie byłoby samemu kodować tę funkcję? To wydaje się być podatne na błędy, ponieważ jestem nowy w scipy.

SCIPY Lognormal (NIEBIESKI) vs. R Lognormal (RED):

Masz jakieś przemyślenia w jakim kierunku? Nawiasem mówiąc, dane są bardzo dobrze dopasowane do modelu R, więc jeśli wygląda na coś innego w Pythonie, możesz je udostępnić.

Dziękuję Ci!

Aktualizacja:

Używam Scipy 0.11

Oto podzbiór danych. Rzeczywista próbka wynosi 38k +, ze średnią 81.53627:

Podzbiór:

x
[60, 170, 137, 138, 81, 140, 78, 46, 1, 168, 138, 148, 145, 35, 82, 126, 66, 147, 88, 106, 80, 54, 83, 13, 102, 54, 134, 34]
numpy.mean (x)
99,071428571428569

Alternatywnie:

Pracuję nad funkcją przechwytywania pliku pdf:

def lognoral(x, mu, sigma):
    a = 1 / (x * sigma * numpy.sqrt(2 * numpy.pi) )
    b = - (numpy.log(x) - mu) ^ 2 / (2 * sigma ^ 2)
    p = a * numpy.exp(b)
    return p

Daje mi to jednak następujące liczby (wypróbowałem kilka na wypadek pomieszania znaczenia sdlog i meanlog):

>>> lognormal(54,4.2991610, 0.5511349)
0.6994656085799437
 >>> lognormal(54,numpy.exp(4.2991610), 0.5511349)
0.9846125119455129
>>> lognormal(54,numpy.exp(4.2991610), numpy.exp(0.5511349))
0.9302407837304372

jakieś pomysły?

Aktualizacja:

ponowne uruchomienie z sugestią „UPQuark”:

kształt, loc, skala (1.0, 50.03445923295007, 19.074457156766517)

Kształt wykresu jest jednak bardzo podobny, a szczyt ma miejsce około 21.

Przeszedłem przez kod źródłowy, aby dojść do następującej interpretacji zwariowanej lognormalnej rutyny.

$\frac{x-\text{loc}}{\text{scale}} \sim \text{Lognormal}(\sigma)$

gdzie jest parametrem „shape”. $\sigma$

Równoważność między parametrami scipy a parametrem R jest następująca:

loc - Brak odpowiednika, to jest odejmowane od twoich danych, tak że 0 staje się minimalnym zakresem danych.

skala - , gdzie jest średnią logu wariacji. (Przy dopasowywaniu zwykle używasz przykładowej średniej z dziennika danych).$\exp{\mu}$$\mu$

kształt - standardowe odchylenie dziennika zmiennej.

Zadzwoniłem lognorm.pdf(x, 0.55, 0, numpy.exp(4.29))tam, gdzie są odpowiednio argumenty (x, kształt, loc, skala), i wygenerowałem następujące wartości:

x pdf

10 0,000106

20 0,002275

30 0,006552

40 0,009979

50 0,11557

60 0,113479

70 0,103327

80 0,008941

90 0,007494

100 0,006155

które wydają się całkiem dobrze pasować do twojej krzywej R.

Logarytmiczna dystrybucja w SciPy pasuje do ogólnej struktury dla wszystkich dystrybucji w SciPy. Wszystkie mają słowo kluczowe dotyczące skali i lokalizacji (domyślnie 0 i 1, jeśli nie podano inaczej). Pozwala to na przesunięcie i skalowanie wszystkich rozkładów w stosunku do znormalizowanej specyfikacji z wyraźnymi implikacjami dla statystyki rozkładu. Rozkłady zazwyczaj mają również jeden lub więcej parametrów „kształtu” (chociaż niektóre, podobnie jak rozkład normalny, nie wymagają żadnych dodatkowych parametrów).

Chociaż to ogólne podejście ładnie ujednolica wszystkie dystrybucje, dla lognormal może powodować pewne zamieszanie ze względu na sposób, w jaki inne pakiety definiują parametry. Nadal bardzo łatwo jest dopasować dowolny rozkład logarytmiczny, jeśli masz na myśli log (średnią rozkładu podstawowego) i sdlog (standardowe odchylenie rozkładu podstawowego).

Najpierw upewnij się, że ustawiłeś parametr położenia na 0. Następnie ustaw parametr kształtu na wartość sdlog. Na koniec ustaw parametr skali na math.exp (meanlog). Zatem rv = scipy.stats.lognorm (0.5511349, scale = math.exp (4.2991610)) utworzy obiekt dystrybucji, którego pdf dokładnie pasuje do twojej wygenerowanej krzywej R. Jako x = numpy.linspace (0,180,1000); plot (x, rv.pdf (x)) zweryfikuje.

Zasadniczo rozkład logarytmiczny SciPy jest uogólnieniem standardowego rozkładu logarytmicznego, który pasuje do standardu dokładnie przy ustawieniu parametru położenia na 0.

Podczas dopasowywania danych za pomocą metody .fit można również użyć słów kluczowych, f0..fn, floc i fshape, aby zachować stały dowolny parametr kształtu, położenia i / lub skali i pasować tylko do innych zmiennych. W przypadku rozkładu logarytmicznego jest to bardzo przydatne, ponieważ zwykle wiesz, że parametr położenia powinien być ustawiony na 0. Zatem scipy.stats.lognorm.fit (zestaw danych, floc = 0) zawsze zwróci parametr położenia jako 0 i zmienia tylko inne parametry kształtu i skali.

Scipy lognormal fit zwraca kształt, lokalizację i skalę. Właśnie uruchomiłem następujące na tablicy przykładowych danych cenowych:

shape, loc, scale = st.lognorm.fit(d_in["price"])

Daje mi to rozsądne szacunki 1,0, 0,09, 0,86, a kiedy go kreślisz, powinieneś wziąć pod uwagę wszystkie trzy parametry.

Parametr kształtu jest odchyleniem standardowym leżącego u podstaw rozkładu normalnego, a skala jest wykładniczą średniej średniej.

Mam nadzieję że to pomoże.

Wydaje się, że rozkład w Scipy dla logarytmu normalnego nie jest taki sam jak w R, lub ogólnie nie jest taki sam jak rozkład, który znam. John D Cook poruszył tę kwestię : http://www.johndcook.com/blog/2010/02/03/statystyczny-distribution-in-scipy/ http://www.johndcook.com/distribution_scipy.html

Jednak nie znalazłem niczego rozstrzygającego na temat używania funkcji gęstości logarytmicznej w Pythonie. Jeśli ktoś chciałby to dodać, nie krępuj się.

Moim dotychczasowym rozwiązaniem jest użycie logarytmicznego pliku pdf ocenianego w zakresie od 0 do 180 (wyłącznie) i używanego jako słownik w skrypcie python.