Asymptotyczna normalność statystyki rzędu rozkładów ciężkich

Tło: Mam próbkę, którą chcę modelować z rozkładem o dużym ogonie. Mam pewne ekstremalne wartości, takie, że zasięg obserwacji jest stosunkowo duży. Moim pomysłem było modelowanie tego z uogólnioną dystrybucją Pareto, i tak zrobiłem. Teraz kwantyl 0,975 moich danych empirycznych (około 100 punktów danych) jest niższy niż kwantyl 0,975 Uogólnionego rozkładu Pareto, który dopasowałem do moich danych. Pomyślałem, czy jest jakiś sposób, aby sprawdzić, czy ta różnica jest czymś, o co należy się martwić?

Wiemy, że asymptotyczny rozkład kwantyli podano jako:

Pomyślałem więc, że dobrym pomysłem byłoby rozbudzenie mojej ciekawości, próbując wykreślić 95% przedziały ufności wokół kwantyla 0,975 uogólnionego rozkładu Pareto z tymi samymi parametrami, jakie otrzymałem z dopasowania moich danych.

Jak widać, pracujemy tutaj z pewnymi ekstremalnymi wartościami. A ponieważ spread jest tak ogromny, funkcja gęstości ma niezwykle małe wartości, co powoduje, że pasma ufności idą do rzędu przy użyciu wariancji powyższej formuły asymptotycznej normalności:$\pm 10^{12}$

$\pm 1.96\frac{0.975*0.025}{n({f_{GPD}(q_{0.975})})^2}$

Nie ma to więc żadnego sensu. Mam rozkład z tylko pozytywnymi wynikami, a przedziały ufności zawierają wartości ujemne. Coś się tu dzieje. Jeśli obliczę pasma wokół kwantyla 0,5, pasma nie są tak duże, ale wciąż ogromne.

Widzę, jak to wygląda w przypadku innej dystrybucji, a mianowicie dystrybucji . Symuluj obserwacji z rozkładu i sprawdź, czy kwantyle mieszczą się w pasmach ufności. Robię to 10000 razy, aby zobaczyć proporcje kwantyli 0,975 / 0,5 symulowanych obserwacji, które mieszczą się w pasmach ufności.$\mathcal{N}(1,1)$$n=100$$\mathcal{N}(1,1)$

    ################################################
# Test at the 0.975 quantile
################################################

#normal(1,1)

#find 0.975 quantile
q_norm<-qnorm(0.975, mean=1, sd=1)
#find density value at 97.5 quantile:
f_norm<-dnorm(q_norm, mean=1, sd=1)
#confidence bands absolute value:
band=1.96*sqrt((0.975*0.025)/(100*(f_norm)^2))
u=q_norm+band
l=q_norm-band

hit<-1:10000
for(i in 1:10000){
  d<-rnorm(n=100, mean=1, sd=1)
  dq<-quantile(d, probs=0.975)

  if(dq[[1]]>=l & dq[[1]]<=u) {hit[i]=1} else {hit[i]=0} 

}
sum(hit)/10000

#################################################################3
# Test at the 0.5 quantile  
#################################################################
#using lower quantile:

#normal(1,1)

#find 0.7 quantile
q_norm<-qnorm(0.7, mean=1, sd=1)
#find density value at 0.7 quantile:
f_norm<-dnorm(q_norm, mean=1, sd=1)
#confidence bands absolute value:
band=1.96*sqrt((0.7*0.3)/(100*(f_norm)^2))
u=q_norm+band
l=q_norm-band

hit<-1:10000
for(i in 1:10000){
  d<-rnorm(n=100, mean=1, sd=1)
  dq<-quantile(d, probs=0.7)

  if(dq[[1]]>=l & dq[[1]]<=u) {hit[i]=1} else {hit[i]=0} 

} 
sum(hit)/10000

EDYCJA : Naprawiłem kod i oba kwantyle dają około 95% trafień przy n = 100 i przy . Jeśli zwiększę standardowe odchylenie do , wtedy bardzo niewiele trafień będzie w pasmach. Tak więc pytanie wciąż jest aktualne.$\sigma=1$$\sigma=2$

EDYCJA 2 : Cofam to, co twierdziłem w pierwszym EDYCJI powyżej, jak wskazano w komentarzach pomocnego dżentelmena. Wygląda na to, że te CI są dobre dla normalnego rozkładu.

Czy ta asymptotyczna normalność statystyki rzędu jest po prostu bardzo złym miernikiem, jeśli chce się sprawdzić, czy prawdopodobne jest zaobserwowanie kwantylu przy określonym rozkładzie kandydata?

Intuicyjnie wydaje mi się, że istnieje związek między wariancją rozkładu (który, jak się uważa, stworzył dane, lub w moim przykładzie R, który wiemy, że utworzył dane) a liczbą obserwacji. Jeśli masz 1000 obserwacji i ogromną wariancję, te pasma są złe. Jeśli ktoś ma 1000 obserwacji i małą wariancję, pasma te mogą mieć sens.

Czy ktoś chce mi to wyjaśnić?

Zakładam, że twoje pochodzenie pochodzi z czegoś takiego jak ten na tej stronie .

Mam rozkład z tylko pozytywnymi wynikami, a przedziały ufności zawierają wartości ujemne.

Cóż, biorąc pod uwagę normalne przybliżenie, które ma sens. Nic nie stoi na przeszkodzie, aby normalne przybliżenie dawało wartości ujemne, dlatego jest to złe przybliżenie wartości ograniczonej, gdy próbka jest mała i / lub wariancja jest duża. Jeśli zwiększysz rozmiar próbki, interwały będą się zmniejszać, ponieważ wielkość próbki jest w mianowniku wyrażenia dla szerokości interwału. Wariancja wchodzi w problem poprzez gęstość: dla tej samej średniej wyższa wariancja będzie miała inną gęstość, wyższą na marginesach i niższą w pobliżu środka. Niższa gęstość oznacza szerszy przedział ufności, ponieważ gęstość jest w mianowniku wyrażenia.

Trochę googlowania znalazło tę stronę , między innymi, która wykorzystuje normalne przybliżenie do rozkładu dwumianowego w celu skonstruowania granic ufności. Podstawową ideą jest to, że każda obserwacja spada poniżej kwantyla z prawdopodobieństwem q , tak że rozkład jest dwumianowy. Gdy wielkość próbki jest wystarczająco duża (to ważne), rozkład dwumianowy jest dobrze przybliżony przez rozkład normalny ze średnią i wariancją . Zatem dolny limit ufności będzie miał indeks , a górny limit ufności będzie miał indeks . Istnieje możliwość, że albo albo$nq$$nq(1-q)$$j = nq - 1.96 \sqrt{nq(1-q)}$$k = nq - 1.96 \sqrt{nq(1-q)}$$k > n$$j < 1$ podczas pracy z kwantylami w pobliżu krawędzi, a odniesienie, które znalazłem, milczy na ten temat. Postanowiłem traktować maksimum lub minimum jako odpowiednią wartość.

W poniższym przepisaniu kodu skonstruowałem limit ufności na danych empirycznych i przetestowałem, czy teoretyczny kwantyl się w to mieści. Ma to dla mnie większy sens, ponieważ kwantyl obserwowanego zestawu danych jest zmienną losową. Pokrycie dla n> 1000 wynosi ~ 0,95. Dla n = 100 jest gorzej przy 0,85, ale należy się tego spodziewać w przypadku kwantyli w pobliżu ogonów o małych próbkach.

#find 0.975 quantile
q <- 0.975
q_norm <- qnorm(q, mean=1, sd=1)

#confidence bands absolute value (note depends on sample size)
n <- 10000
band <- 1.96 * sqrt(n * q * (1 - q))

hit<-1:10000
for(i in 1:10000){
  d<-sort(rnorm(n, mean=1, sd=1))
  dq<-quantile(d, probs=q)
  u <- ceiling(n * q + band)
  l <- ceiling(n * q - band)
  if (u > n) u = n
  if (l < 1) l = 1
  if(q_norm>=d[l] & q_norm<=d[u]) {hit[i]=1} else {hit[i]=0} 

}
sum(hit)/10000

Jeśli chodzi o określenie, która wielkość próby jest „wystarczająco duża”, cóż, większe jest lepsze. To, czy jakaś konkretna próbka jest „wystarczająco duża”, zależy silnie od aktualnego problemu i tego, jak bardzo jesteś wybredny w takich kwestiach, jak pokrycie limitów pewności.