Czy w przypadku macierzy losowej SVD nie powinien niczego wyjaśniać? Co ja robię źle?

Gdybym zbudował matrycę 2-D złożoną wyłącznie z losowych danych, oczekiwałbym, że komponenty PCA i SVD w zasadzie niczego nie wyjaśnią.

Zamiast tego wydaje się, że pierwsza kolumna SVD wydaje się wyjaśniać 75% danych. Jak to możliwe? Co ja robię źle?

Oto fabuła:

wprowadź opis zdjęcia tutaj

Oto kod R:

set.seed(1)
rm(list=ls())
m <- matrix(runif(10000,min=0,max=25), nrow=100,ncol=100)
svd1 <- svd(m, LINPACK=T)
par(mfrow=c(1,4))
image(t(m)[,nrow(m):1])
plot(svd1$d,cex.lab=2, xlab="SVD Column",ylab="Singluar Value",pch=19)

percentVarianceExplained = svd1$d^2/sum(svd1$d^2) * 100
plot(percentVarianceExplained,ylim=c(0,100),cex.lab=2, xlab="SVD Column",ylab="Percent of variance explained",pch=19)

cumulativeVarianceExplained = cumsum(svd1$d^2/sum(svd1$d^2)) * 100
plot(cumulativeVarianceExplained,ylim=c(0,100),cex.lab=2, xlab="SVD column",ylab="Cumulative percent of variance explained",pch=19)

Aktualizacja

Dziękuję @Aaron. Jak zauważyłeś, poprawka polegała na dodaniu skalowania do macierzy, aby liczby były wyśrodkowane wokół 0 (tj. Średnia wynosi 0).

m <- scale(m, scale=FALSE)

Oto poprawiony obraz, pokazujący macierz z losowymi danymi, pierwsza kolumna SVD jest bliska 0, zgodnie z oczekiwaniami.

Poprawiony obraz

r pca svd Contango
źródło

Twoja macierz jest w przybliżeniu równomiernie rozłożona na kostce jednostkowej w . SVD oblicza momenty bezwładności na temat pochodzenia . W „całkowita wariancja” musi być razy większa od interwału jednostkowego, czyli . Łatwo obliczyć, że moment wzdłuż głównej osi sześcianu (emanujący z początku) wynosi a wszystkie pozostałe momenty - z racji symetrii - równe . Dlatego pierwsza wartość własna wynosi całości. Dla to

[0, 1]^{100}

$[0,1]^{100}$

R^{100}

$\mathbb{R}^{100}$

R^{n}

$\mathbb{R}^n$

n

$n$

1 / 3

$1/3$

n / 3 - (n - 1) / 12

$n/3-(n-1)/12$

1 / 12

$1/12$

(n / 3 - (n - 1) / 12) / (n / 3) = 3 / 4 + 1 / (4 n)

$(n/3-(n-1)/12)/(n/3)=3/4+1/(4n)$

n = 100

$n=100$

75.25

$75.25$ %, widoczne na trzecim wykresie.

whuber

Odpowiedzi:

Pierwszy komputer wyjaśnia, że zmienne nie są wyśrodkowane wokół zera. Skalowanie najpierw lub centrowanie zmiennych losowych wokół zera da oczekiwany wynik. Na przykład jeden z tych:

m <- matrix(runif(10000,min=0,max=25), nrow=100,ncol=100)
m <- scale(m, scale=FALSE)

m <- matrix(runif(10000,min=-25,max=25), nrow=100,ncol=100)

Aaron opuścił Stack Overflow
źródło

Podnosisz dobrą rację, ale myślę, że to tylko część historii. Rzeczywiście, zgaduję, że PO spróbuje każdego z nich i nadal będzie zaskoczony wynikiem. Faktem jest, że ponieważ liczby pojedyncze są uporządkowane w danych wyjściowych, nie pojawią się (a nawet nie będą) równomiernie rozmieszczone, jak można naiwnie oczekiwać na podstawie „losowych” danych. W tym przypadku rozkład Marchenko-Pastur reguluje ich zachowanie.

kardynał

@Aaron Dziękuję, miałeś absolutną rację. Dodałem wykres skorygowanego wyjścia powyżej, aby pokazać, jak piękny jest wynik.

Contango,

@cardinal Dziękuję za komentarz, masz absolutną rację (patrz wykres utworzony przez poprawiony kod powyżej). Wierzę, że wartości SVD stałyby się mniej równomiernie rozłożone w miarę zmniejszania się macierzy, ponieważ mniejsza matryca miałaby większe szanse na posiadanie wzorów, które nie byłyby zmiażdżone przez prawo wielkich liczb.

Contango,

Dodam bardziej wizualną odpowiedź do twojego pytania, poprzez porównanie z zerowym modelem. Procedura losowo tasuje dane w każdej kolumnie, aby zachować ogólną wariancję, podczas gdy kowariancja między zmiennymi (kolumnami) zostaje utracona. Odbywa się to kilka razy, a wynikowy rozkład liczby pojedynczej w losowej macierzy jest porównywany z wartościami pierwotnymi.

Używam prcompzamiast svddo rozkładu macierzy, ale wyniki są podobne:

set.seed(1)
m <- matrix(runif(10000,min=0,max=25), nrow=100,ncol=100)

S <- svd(scale(m, center = TRUE, scale=FALSE))
P <- prcomp(m, center = TRUE, scale=FALSE)
plot(S$d, P$sdev) # linearly related

Porównanie modeli zerowych jest wykonywane na macierzy centrowanej poniżej:

library(sinkr) # https://github.com/marchtaylor/sinkr

# centred data
Pnull <- prcompNull(m, center = TRUE, scale=FALSE, nperm = 100)
Pnull$n.sig
boxplot(Pnull$Lambda[,1:20], ylim=range(Pnull$Lambda[,1:20], Pnull$Lambda.orig[1:20]), outline=FALSE, col=8, border="grey50", log="y", main=paste("m (center=FALSE); n sig. =", Pnull$n.sig))
lines(apply(Pnull$Lambda, 2, FUN=quantile, probs=0.95))
points(Pnull$Lambda.orig[1:20], pch=16)

Poniżej znajduje się wykres pudełkowy permutowanej macierzy z 95% kwantylem każdej liczby pojedynczej pokazanej jako linia ciągła. Pierwotne wartości PCA mto kropki. wszystkie leżą poniżej linii 95% - dlatego ich amplituda jest nie do odróżnienia od szumu przypadkowego.

Tę samą procedurę można wykonać w niecentralnej wersji mz tym samym wynikiem - Brak znaczących pojedynczych wartości:

# centred data
Pnull <- prcompNull(m, center = FALSE, scale=FALSE, nperm = 100)
Pnull$n.sig
boxplot(Pnull$Lambda[,1:20], ylim=range(Pnull$Lambda[,1:20], Pnull$Lambda.orig[1:20]), outline=FALSE, col=8, border="grey50", log="y", main=paste("m (center=TRUE); n sig. =", Pnull$n.sig))
lines(apply(Pnull$Lambda, 2, FUN=quantile, probs=0.95))
points(Pnull$Lambda.orig[1:20], pch=16)

Dla porównania spójrzmy na zestaw danych z nieprzypadkowym zestawem danych: iris

# iris dataset example
m <- iris[,1:4]
Pnull <- prcompNull(m, center = TRUE, scale=FALSE, nperm = 100)
Pnull$n.sig
boxplot(Pnull$Lambda, ylim=range(Pnull$Lambda, Pnull$Lambda.orig), outline=FALSE, col=8, border="grey50", log="y", main=paste("m (center=FALSE); n sig. =", Pnull$n.sig))
lines(apply(Pnull$Lambda, 2, FUN=quantile, probs=0.95))
points(Pnull$Lambda.orig[1:20], pch=16)

Tutaj pierwsza wartość pojedyncza jest znacząca i wyjaśnia ponad 92% całkowitej wariancji:

P <- prcomp(m, center = TRUE)
P$sdev^2 / sum(P$sdev^2)
# [1] 0.924618723 0.053066483 0.017102610 0.005212184

Marc w pudełku
źródło

+1. Przykład zestawu danych Iris jest interesujący, ponieważ patrząc na pierwsze dwa komputery osobiste (jak np. We własnym poście tutaj stats.stackexchange.com/a/88092 ) jest całkiem jasne, że drugi niesie jakiś sygnał. Test permutacji (inaczej tasowanie) wskazuje, że tylko pierwszy jest „znaczący”. Oczywiste jest, że tasowanie ma tendencję do niedoceniania liczby komputerów: duża wariancja pierwszego prawdziwego komputera zostanie „rozłożona” na tasowane komputery i podniesie je wszystkie, zaczynając od drugiego. Można opracować bardziej czułe testy, które to uwzględniają, ale jest to rzadko wykonywane.

ameba mówi Przywróć Monikę

@amoeba - Doskonały komentarz. Od jakiegoś czasu zastanawiam się nad efektem „rozprzestrzeniania się”. Przypuszczam, że test weryfikacji krzyżowej może być jednym z bardziej wrażliwych, na który się powołujesz (np. Twoja odpowiedź tutaj )? Byłoby świetnie, gdybyś mógł podać przykład / odniesienie.

Marc w pudełku

Zwykle wolę korzystać z walidacji krzyżowej (na podstawie błędu rekonstrukcji, zgodnie z moją odpowiedzią tutaj ), ale tak naprawdę nie jestem pewien, czy nie cierpi z powodu tego samego rodzaju niewrażliwości, czy nie. Może warto wypróbować to w zestawie danych Iris. Jeśli chodzi o podejście oparte na tasowaniu, nie znam żadnego odniesienia do rozliczania tego „rozprzestrzeniania się”, po prostu znam ludzi, którzy ostatnio nad tym pracowali. Myślę, że wkrótce chcieli to napisać. Chodzi o to, aby wprowadzić pewne czynniki zmniejszania skali dla wariantów wyższych komputerów z tasowaniem.

ameba mówi Przywróć Monikę

@amoeba - Dzięki za ten link. Wiele mi to wyjaśnia. Za szczególnie interesujące uważam, że krzyżowa walidacja w PCA wykorzystuje metody, które mogą działać na zestawach danych z brakującymi wartościami. Podjąłem kilka prób zastosowania tego podejścia i (jak twierdzisz) podejście tasowania modelu zerowego rzeczywiście nie docenia liczby znaczących komputerów. Jednak w przypadku zestawu danych tęczówki konsekwentnie zwracam pojedynczy komputer z powodu błędu rekonstrukcji. Interesujące, biorąc pod uwagę to, co wspomniałeś o fabule. Możliwe, że gdybyśmy oceniali błąd w oparciu o prognozy dotyczące gatunków, wyniki mogą być inne.

Marc w pudełku

Z ciekawości wypróbowałem to na danych Iris. W rzeczywistości otrzymuję dwa znaczące komputery z metodą krzyżowej weryfikacji. Zaktualizowałem mój połączony post, zobacz tam.

ameba mówi Przywróć Monikę