Oczekiwana wartość fałszywej korelacji

12

Rysujemy próbek, każda o rozmiarze , niezależnie od rozkładu Normal .Nn(μ,σ2)

Z próbek wybieramy następnie 2 próbki, które mają najwyższą (absolutną) korelację Pearsona ze sobą.N

Jaka jest oczekiwana wartość tej korelacji?

Dzięki [PS To nie zadanie domowe]

P Sellaz
źródło
2
(+1) To byłoby dość trudne zadanie domowe :-). Czy potrzebujesz ogólnej odpowiedzi, czy mógłbyś (być może) skoncentrować się na określonych wartościach lub ? Na przykład, możliwe jest opracowanie dobrych przybliżeń, gdy jest znacznie większe niż ; w innych przypadkach potrzebne byłyby różne przybliżenia. NnnN
whuber
1
Miałem nadzieję na ogólną odpowiedź, ale taką, w której założenie byłoby OK! Dla określonych wartości i nie byłoby to tak interesujące, ponieważ mogę spojrzeć na takie konkretne przypadki za pomocą symulacji (właśnie to robię w tej chwili), ale wciąż może być interesujące. n>>NNn
P Sellaz,
1
Myślę, że ogólne rozwiązanie jakiejkolwiek rzeczywistej użyteczności jest prawdopodobnie mało prawdopodobne, choć mogę się mylić. Jest to dość ściśle związane z niektórymi otwartymi problemami na styku geometrii i algebry liniowej. W aplikacjach potrzeba informacji o takich ilościach pojawia się, na przykład, w czujnikach skompresowanych.
kardynał
1
FWIW, jest to wynik właśnie przeprowadzonej symulacji: przy użyciu Normalnej (0,1) stwierdziłem, że średnia korelacja, (ponad 1000 symulacji) i liczba próbek są w przybliżeniu powiązane przez dla i przy użyciu modelu regresji liniowej. Dopasowanie modelu i zwykła diagnostyka były całkiem dobre. Odkryłem również, że średnia korelacja była w przybliżeniu normalnie rozłożona (choć nieco pochylona w prawo). ρN
ρ=0.025+0.113ln(N)0.008ln(N)2
n=1004Nn
P Sellaz,

Odpowiedzi:

9

Znalazłem następujący artykuł, który dotyczy tego problemu: Jiang, Tiefeng (2004). Rozkłady asymptotyczne największych wpisów macierzy korelacji próbek. The Annals of Applied Prawdopodobieństwo, 14 (2), 865-880

Jiang pokazuje asymptotyczny rozkład statystyki, gdzie jest korelacją między tym i tym losowym wektorem o długości (z ), wynosiLn=max1i<jN|ρij|ρijijnij

limnPr[nLn24logn+log(log(n))y]=exp(1a28πexp(y/2)),
gdzie zakłada się, że istnieje w dokumencie, a jest funkcją .a=limnn/NNn

Najwyraźniej ten wynik dotyczy wszelkich rozkładów dystrybucji z wystarczającą liczbą skończonych momentów ( edycja: patrz komentarz @ kardynała poniżej). Jiang zwraca uwagę, że jest to ekstremalny rozkład wartości typu I. Lokalizacja i skala są

σ=2,μ=2log(1a28π).

Oczekiwana wartość rozkładu EV typu I to , gdzie oznacza stałą Eulera. Jednakże, jak zauważono w komentarzach, konwergencja w dystrybucji sama w sobie nie gwarantuje konwergencji środków do dystrybucji ograniczającej.μ+σγγ

Jeśli moglibyśmy pokazać taki wynik w tym przypadku, to asymptotyczna oczekiwana wartośćwynosiłabynLn24logn+log(log(n))

limnE[nLn24logn+log(log(n))]=2log(a28π)+2γ.

Należy zauważyć, że dałoby to asymptotyczną wartość oczekiwaną największej korelacji do kwadratu, podczas gdy pytanie dotyczyło oczekiwanej wartości największej korelacji bezwzględnej. Więc nie 100% tam, ale blisko.

Zrobiłem kilka krótkich symulacji, które doprowadziły mnie do wniosku, że albo 1) jest problem z moją symulacją (prawdopodobnie), 2) jest problem z moją transkrypcją / algebrą (prawdopodobnie również), lub 3) przybliżenie nie jest poprawne dla zastosowane wartości i I. Być może PO może zważyć niektóre wyniki symulacji przy użyciu tego przybliżenia?nN

jmtroos
źródło
2
A na marginesie: bardzo podobało mi się to pytanie - zastanawiałem się wcześniej nad tym pytaniem. Byłem zaskoczony połączeniem z dystrybucją typu I - uznałem, że to całkiem fajne. Chciałbym tylko zrozumieć matematykę, która do tego doprowadziła ...
jmtroos,
1
(+1) Niezłe znalezisko !! Myślę, że możemy założyć, że dodatni pierwiastek kwadratowy tego jest równoważny z oczekiwaną wartością największej korelacji absolutnej? W twoim odczuciu oczekiwania nie możemy po prostu wyjąć wszystkich części obejmujących aby uzyskać: ? W każdym razie porównałem to do moich symulacji i wygląda to całkiem blisko! Mój kod R jest bardzo zaniedbany, więc postaram się uporządkować go później dzisiaj lub jutro i umieścić go ...Lnn
E[Ln2]=1n{2log(N2n28π)+2γ+4lognlog(log(n))}
P Sellaz
BTW, gazeta jest dostępna bezpośrednio tutaj projecteuclid.org/DPubS/Repository/1.0/…
P Sellaz
3
(+1) To jest bardzo fajny papier, a ja tylko go przejrzałem, ale musimy tu być trochę ostrożni . Kilka uwag: ( 1 ) Wyniki dotyczą reżimu , więc wymiar wektorów musi rosnąć mniej więcej proporcjonalnie do liczby wektorów branych pod uwagę dla tych wyników trzymać. ( 2 ) Nawet w tym przypadku wyniki nie dotyczą „żadnego” rozkładu; Rzeczywiście, warunki w artykule wymagają, aby zmienne losowe były „prawie wykładniczo ograniczone” w tym sensie, że zasadniczo potrzebujemy 30. momentu, aby być skończonym! (cd.)n/pγ(0,)
kardynał
3
(cd.) ( 3 ) Konwergencja w dystrybucji nie gwarantuje konwergencji środków do dystrybucji ograniczającej. W tym celu zwykle używamy czegoś podobnego do jednolitej integralności zbioru . Nie zostało to wykazane w artykule, a ponieważ zajmowanie się rozkładami wartości ekstremalnych może nie być prawdą. Jednym z moich ulubionych przykładów tego zjawiska jest sekwencja zmiennych losowych, które zbiegają się w rozkładzie do , ale można uzyskać środki do zbieżności z dowolną wybraną stałą dodatnią. {Ln}χ12
kardynał
2

Oprócz odpowiedzi udzielonej przez @jmtroos, poniżej znajdują się szczegóły mojej symulacji i porównanie z wyprowadzeniem @ jmtroos oczekiwań od Jiang (2004) , to znaczy:

E[Ln2]=1n{2log(N2n28π)+2γ+4lognlog(log(n))}

Wartości tego oczekiwania wydają się przekraczać symulowane wartości dla małego i poniżej dla dużego i wydają się nieco różnić w miarę wzrostuJednak różnice zmniejszają się wraz ze wzrostem , jak można się spodziewać, ponieważ artykuł twierdzi, że rozkład jest asymptotyczny. Próbowałem różnych . Poniższa symulacja wykorzystuje . Jestem całkiem nowy w R, więc wszelkie wskazówki i sugestie dotyczące ulepszenia mojego kodu będą mile widziane.NNNnn[100,500]n=200

set.seed(1)

ns <- 500
# number of simulations for each N

n <- 200
# length of each vector

mu <- 0
sigma <- 1
# parameters for the distribution we simulate from

par(mfrow=c(5,5))
x<-trunc(seq(from=5,to=n, length=20))
#vector of Ns

y<-vector(mode = "numeric")
#vector to store the mean correlations

k<- 1
#index for y

for (N in x) {
# loop over a range of N

    dt <- matrix(nrow=n,ncol=N)

    J <- vector(mode = "numeric")
    # vector to store the simulated largest absolute 
    # correlations for each N

    for (j in 1:ns) {
    # for each N, simulated ns times    

      for (i in 1:N) {
        dt[,i] <- rnorm(n,mu,sigma)
      }
      # perform the simulation

      M<-matrix(cor(dt),nrow=N,ncol=N)
      m <- M
      diag(m) <- NA
      J[j] <- max(abs(m), na.rm=TRUE)   
      # obtain the largest absolute correlation
      # these 3 lines came from stackoverflow
  }

    hist(J,main=paste("N=",N, " n=",n, " N(0,1)", "\nmean=",round(J[j],4))) 
    y[k]<-mean(J)
    k=k+1
}

lm1 <- lm(y~log(x))
summary(lm1)

logx_sq=log(x)^2
lm2<-lm(y~log(x)+logx_sq)
summary(lm2)
# linear models for these simulations

# Jiang 2004 paper, computation:

gamma = 0.5772
yy <- vector(mode = "numeric")
yy <- sqrt((2*log((x^2)/(sqrt(8*pi)*n^2)) + 2*gamma-(-4*log(n)+log(log(n))))/n)


plot(x,yy)
# plot the simulated correlations
points(x,y,col='red')
# add the points using the expectation
P Sellaz
źródło
Zobacz moje komentarze do drugiej odpowiedzi, która może (ale nie musi) pomóc wyjaśnić niektóre zauważone rozbieżności.
kardynał