Jestem bardzo nowy w częściowych najmniejszych kwadratach (PLS) i staram się zrozumieć wynik funkcji R plsr()
w pls
pakiecie. Symulujmy dane i uruchom PLS:
library(pls)
n <- 50
x1 <- rnorm(n); xx1 <- scale(x1)
x2 <- rnorm(n); xx2 <- scale(x2)
y <- x1 + x2 + rnorm(n,0,0.1); yy <- scale(y)
p <- plsr(yy ~ xx1+xx2, ncomp=1)
Spodziewałem się, że następujące liczby i
> ( w <- loading.weights(p) )
Loadings:
Comp 1
xx1 0.723
xx2 0.690
Comp 1
SS loadings 1.0
Proportion Var 0.5
> a <- w["xx1",]
> b <- w["xx2",]
> a^2+b^2
[1] 1
są obliczane w celu maksymalizacji
> cor(y, a*xx1+b*xx2)
[,1]
[1,] 0.9981291
ale nie jest tak dokładnie:
> f <- function(ab){
+ a <- ab[1]; b <- ab[2]
+ cor(y, a*xx1+b*xx2)
+ }
> optim(c(0.7,0.6), f, control=list(fnscale=-1))
$par
[1] 0.7128259 0.6672870
$value
[1] 0.9981618
Czy to błąd numeryczny, czy też źle rozumiem naturę i ? b
Chciałbym również wiedzieć, jakie są te współczynniki:
> p$coef
, , 1 comps
yy
xx1 0.6672848
xx2 0.6368604
EDYCJA : Teraz widzę, co p$coef
to jest:
> x <- a*xx1+b*xx2
> coef(lm(yy~0+x))
x
0.9224208
> coef(lm(yy~0+x))*a
x
0.6672848
> coef(lm(yy~0+x))*b
x
0.6368604
Więc myślę, że mam rację co do natury i . b
EDYCJA: W świetle komentarzy udzielonych przez @chl uważam, że moje pytanie nie jest wystarczająco jasne, dlatego podaję więcej szczegółów. W moim przykładzie jest wektor odpowiedzi i dwukolumnowa macierz predyktorów, a ja używam znormalizowanej wersji z i znormalizowanej wersji z (wyśrodkowana i podzielona przez odchylenia standardowe). Definicja pierwszych PLS składnika jest z i wybiera się w celu uzyskania wartości maksymalnej produktu wewnętrznej .X ˜ Y Y ˜ X X t 1 t 1 = a ˜ X 1 + b ˜ X 2b ⟨ t 1 , ~ Y ⟩ T 1 TZatem jest to równoważne z maksymalizacją korelacji między i , prawda?
źródło
pls
tym dokumencie JSS znajduje się dobry przegląd pakietu i regresji PLS .?coef.mvr
Odpowiedzi:
Regresja PLS opiera się na algorytmach iteracyjnych (np. NIPALS, SIMPLS). Twój opis głównych pomysłów jest poprawny: szukamy jednego (PLS1, jednej zmiennej odpowiedzi / wielu predyktorów) lub dwóch (PLS2, z różnymi trybami, wielu zmiennych odpowiedzi / wielu predyktorów) wektora (-ów) wag, (i ) powiedzmy, aby utworzyć kombinację liniową pierwotnej zmiennej (zmiennych), tak że kowariancja między Xu i Y (Yv, dla PLS2) jest maksymalna. Skupmy się na wyodrębnieniu pierwszej pary wag powiązanych z pierwszym składnikiem. Formalnie kryterium optymalizacji czyta W twoim przypadku jest jednoznaczne, więc sprowadza się do maksymalizacjiu v
X=[x_1;x_2]
Powinienem podziękować Arthurowi Tenenhausowi, który wskazał mi właściwy kierunek.
Używanie wektorów wagi jednostkowej nie jest ograniczające, a niektóre pakiety (u winieta oferuje również dobrą dyskusję (str. 26–29). Szczególnie ważny jest również fakt, że większość procedur PLS (przynajmniej ta, którą znam w R) zakłada, że udostępniasz niestandardowe zmienne, ponieważ centrowanie i / lub skalowanie jest obsługiwane wewnętrznie (jest to szczególnie ważne na przykład podczas sprawdzania poprawności krzyżowej ).
pls. regression
w wersji plsgenomics , oparte na kodzie z wcześniejszego pakietu Wehrensapls.pcr
) zwracają niestandardowe wektory wagi (ale z ukrytymi składnikami wciąż o normie 1), jeśli są wymagane. Ale większość pakietów PLS zwróci znormalizowany , w tym ten, którego użyłeś, zwłaszcza te implementujące algorytm SIMPLS lub NIPALS; Znalazłem dobry przegląd obu podejść w prezentacji Barry'ego M. Wise'a, właściwości regresji częściowych najmniejszych kwadratów (PLS) i różnic między algorytmami , ale chemometriiBiorąc pod uwagę ograniczenie , wektor jestu′u=1 u
Za pomocą małej symulacji można uzyskać w następujący sposób:
Możesz porównać powyższe wyniki (
u=[0.5792043;0.8151824]
w szczególności) z tym, co dadzą pakiety R. Np. Używając NIPALS z pakietu chemometrii (inna implementacja, o której wiem, że jest dostępna w pakiecie mixOmics ), otrzymalibyśmy:Podobne wyniki można uzyskać za pomocą
plsr
domyślnego algorytmu PLS jądra:We wszystkich przypadkach możemy sprawdzić, czy ma długość 1.u
Pod warunkiem, że zmienisz funkcję, aby zoptymalizować na czytającą
a
u
następnie normalizować (u <- u/sqrt(crossprod(u))
), powinieneś być bliżej powyższego rozwiązania.Sidenote : Jako kryterium (1) jest równa znajduje się w lewym pojedynczej wektora z SVD odpowiadający największej wartości własnej:
W bardziej ogólnym przypadku (PLS2) sposobem na podsumowanie powyższego jest stwierdzenie, że pierwsze wektory kanoniczne PLS są najlepszym przybliżeniem macierzy kowariancji X i Y w obu kierunkach.
Bibliografia
źródło