Odchylenie estymatora momentu rozkładu logarytmicznego

25

Robię eksperyment liczbowy, który polega na próbkowaniu logarytmicznego rozkładu i próbuję oszacować momenty dwiema metodami:XLN(μ,σ)E[Xn]

  1. Patrząc na średnią próbnąXn
  2. Oszacowanie i przy użyciu przykładowych środków dla , a następnie wykorzystując fakt, że dla rozkładu logarytmicznego mamy .μσ2log(X),log2(X)E[Xn]=exp(nμ+(nσ)2/2)

Pytanie brzmi :

Odkryłem eksperymentalnie, że druga metoda działa znacznie lepiej niż pierwsza, gdy utrzymuję stałą liczbę próbek i zwiększam o jakiś czynnik T. Czy istnieje jakieś proste wytłumaczenie tego faktu?μ,σ2

Dołączam cyfrę, na której oś x to T, zaś oś y to wartości porównujące prawdziwe wartości (pomarańczowa linia), do wartości szacunkowych. metoda 1 - niebieskie kropki, metoda 2 - zielone kropki oś y jest w skali logarytmicznejE[X2]E[X2]=exp(2μ+2σ2)

Prawdziwe i szacunkowe wartości dla $ \ mathbb {E} [X ^ 2] $.  Niebieskie kropki to średnie próbki dla $ \ mathbb {E} [X ^ 2] $ (metoda 1), podczas gdy zielone kropki są szacowanymi wartościami przy użyciu metody 2. Pomarańczowa linia jest obliczana na podstawie znanego $ \ mu $, $ \ sigma $ według tego samego równania jak w metodzie 2. Oś y jest w skali logarytmicznej

EDYTOWAĆ:

Poniżej znajduje się minimalny kod Mathematica do wygenerowania wyników dla jednego T, z wynikiem:

   ClearAll[n,numIterations,sigma,mu,totalTime,data,rmomentFromMuSigma,rmomentSample,rmomentSample]
(* Define variables *)
n=2; numIterations = 10^4; sigma = 0.5; mu=0.1; totalTime = 200;
(* Create log normal data*)
data=RandomVariate[LogNormalDistribution[mu*totalTime,sigma*Sqrt[totalTime]],numIterations];

(* the moment by theory:*)
rmomentTheory = Exp[(n*mu+(n*sigma)^2/2)*totalTime];

(*Calculate directly: *)
rmomentSample = Mean[data^n];

(*Calculate through estimated mu and sigma *)
muNumerical = Mean[Log[data]]; (*numerical \[Mu] (gaussian mean) *)
sigmaSqrNumerical = Mean[Log[data]^2]-(muNumerical)^2; (* numerical gaussian variance *)
rmomentFromMuSigma = Exp[ muNumerical*n + (n ^2sigmaSqrNumerical)/2];

(*output*)
Log@{rmomentTheory, rmomentSample,rmomentFromMuSigma}

Wydajność:

(*Log of {analytic, sample mean of r^2, using mu and sigma} *)
{140., 91.8953, 137.519}

powyżej, drugi wynik to średnia próbki , która jest poniżej dwóch pozostałych wynikówr2

użytkownik29918
źródło
2
Bezstronny estymator nie oznacza, że ​​niebieskie kropki powinny znajdować się w pobliżu oczekiwanej wartości (krzywa pomarańczowa). Estymator może być bezstronny, jeśli ma wysokie prawdopodobieństwo, że jest zbyt niski, a małe (być może znikomo małe) prawdopodobieństwo, że jest zbyt wysoki. Tak się dzieje, gdy T rośnie, a wariancja staje się bzdura ogromna (patrz moja odpowiedź).
Matthew Gunn
Aby uzyskać obiektywne estymatory, zobacz stats.stackexchange.com/questions/105717 . Wartości średniej i wariancji UMVUE podano w odpowiedziach i komentarzach do nich.
whuber

Odpowiedzi:

22

W tych wynikach jest coś zagadkowego

  1. pierwsza metoda zapewnia obiektywny estymator , a mianowicie ma jako środek. Dlatego niebieskie kropki powinny znajdować się wokół oczekiwanej wartości (pomarańczowa krzywa);E[X2]
    1Ni=1NXi2
    E[X2]
  2. druga metoda zapewnia tendencyjny estymator , a mianowicie gdy i są obiektywnymi estymatorami odpowiednio i , a zatem dziwne jest, że zielone kropki są wyrównane z pomarańczową krzywą.E[X2]
    E[exp(nμ^+n2σ^2/2)]>exp(nμ+(nσ)2/2)
    μ^σ^²μσ²

ale wynikają one z problemu, a nie z obliczeń numerycznych: powtórzyłem eksperyment w R i otrzymałem następujący obraz z tym samym kodem koloru i tą samą sekwencją i , który reprezentuje każdy estymator podzielony przez prawdziwe oczekiwania:μTσT

Dwa empiryczne drugie momenty, oparte na symulacjach logarytmicznych 10 normal

Oto odpowiedni kod R:

moy1=moy2=rep(0,200)
mus=0.14*(1:200)
sigs=sqrt(0.13*(1:200))
tru=exp(2*mus+2*sigs^2)
for (t in 1:200){
x=rnorm(1e5)
moy1[t]=mean(exp(2*sigs[t]*x+2*mus[t]))
moy2[t]=exp(2*mean(sigs[t]*x+mus[t])+2*var(sigs[t]*x+mus[t]))}

plot(moy1/tru,col="blue",ylab="relative mean",xlab="T",cex=.4,pch=19)
abline(h=1,col="orange")
lines((moy2/tru),col="green",cex=.4,pch=19)

Stąd rzeczywiście zapada się drugi moment empiryczny, gdy rośnie i , co przypisałbym ogromnemu wzrostowi wariancji wspomnianego drugiego momentu empirycznego, gdy rośnie i .μσμσ

Moje wyjaśnienie tego dziwnego zjawiska jest takie, że chociaż oczywiście jest średnią z , nie jest to centralna wartość: w rzeczywistości mediana jest równa . Reprezentując zmienną losową jako gdzie , jasne jest, że gdy jest duża wystarczy, że zmienna losowa prawie nigdy nie ma wielkości . Innymi słowy, jeśli toE[X2]X2X2e2μX2exp{2μ+2σϵ}ϵN(0,1)σσϵσ2XLN(μ,σ)

P(X2>E[X2])=P(log{X2}>2μ+2σ2)=P(μ+σϵ>μ+σ2)=P(ϵ>σ)=1Φ(σ)
które mogą być dowolnie małe.
Xi'an
źródło
1
Jestem też zdziwiony.
Dodam
Dobrze. Dzięki! Po dodaniu liczb widzę teraz, że mój niewielki rozmiar próbki naprawdę nie był w stanie sprostać zadaniu!
user29918
2
@ user29918: Przepraszam, nie widzę problemu jako próby, ale raczej fakt, że log-normal staje się bardzo wypaczony, gdy staje się nieskończonością, co oznacza, że ​​staje się bezużyteczny. σ
Xi'an
2
@ Xi'an Dobre rzeczy! . To oddaje dokładnie w równaniach to, co (raczej nieprecyzyjnie) próbowałem wyrazić słowami, że wraz ze wzrostem staje się coraz bardziej prawdopodobne (a dla dużej prawie pewne), że obserwacja jest poniżej średniej. Rzeczywiście prawdopodobieństwo jest tak wysokie, że jest bardzo prawdopodobne, że cała próbka jest poniżej średniej! P(X2>E[X2])=1Φ(σ)σσ
Matthew Gunn
2
Ten typ asymptotyczny nie jest zbyt pomocny, ponieważ liczba symulacji potrzebnych do prawidłowego przybliżenia momentów rośnie wykładniczo szybko z . σ
Xi'an
13

Pomyślałem, że rzuciłem figi, które pokazują, że wykresy user29918 i Xi'an są spójne. Ryc. 1 pokazuje, co zrobił użytkownik29918, a ryc. 2 (w oparciu o te same dane) robi to, co Xi'an zrobił dla swojego wykresu. Ten sam wynik, inna prezentacja.

Co się dzieje, gdy T rośnie, wariancje stają się ogromne, a estymator staje się jak próba oszacowania średniej populacji Powerball Lotto poprzez zakup biletów Lotto! W dużej części czasu nie docenisz wypłaty (ponieważ żadna obserwacja próbki nie trafi w dziesiątkę), a niewielki procent czasu znacznie przeszacujesz wypłatę (ponieważ w próbie jest zwycięzca głównej wygranej). Średnia z próby jest obiektywnym oszacowaniem, ale nie jest sprecyzowana, nawet przy tysiącach losowań! W rzeczywistości, ponieważ coraz trudniej jest wygrać w lotto, Twoja średnia próbki będzie mniejsza niż średnia populacji w przeważającej części czasu.1nixi2

Dalsze komentarze:

  1. Bezstronny estymator nie oznacza, że ​​estymator powinien być blisko! Niebieskie kropki nie muszą być bliskie oczekiwaniom. Na przykład. pojedyncza obserwacja wybrana losowo daje obiektywne oszacowanie średniej populacji, ale nie można oczekiwać, że estymator będzie bliski.
  2. Problem pojawia się, gdy wariancja staje się absolutnie astronomiczna. Gdy wariancja idzie w gówno, szacunek dla pierwszej metody opiera się na kilku obserwacjach. Zaczynasz też mieć małe, maleńkie prawdopodobieństwo NIESAMOWITEJ, NIESAMOWITEJ, NIESAMOWITEJ dużej liczby ...
  3. To jest intuicyjne wyjaśnienie. Xi'an ma bardziej formalne pochodzenie. Jego wynik implikuje, że gdy staje się duża, bardzo mało prawdopodobne jest, aby kiedykolwiek narysować obserwację powyżej średniej, nawet przy tysiącach obserwacji . Mój język „wygrywania w lotto” odnosi się do wydarzenia, w którym . P(X2>E[X2])=1Φ(σ)σX2>E[X2]wprowadź opis zdjęcia tutaj

wprowadź opis zdjęcia tutaj

Matthew Gunn
źródło