Wydajna regresja liniowa online

53

Analizuję niektóre dane, w których chciałbym przeprowadzić zwykłą regresję liniową, jednak nie jest to możliwe, ponieważ mam do czynienia z ustawieniem on-line z ciągłym strumieniem danych wejściowych (które szybko stają się zbyt duże dla pamięci) i potrzebują zaktualizować oszacowania parametrów podczas ich zużycia. tzn. nie mogę po prostu załadować wszystkiego do pamięci i przeprowadzić regresji liniowej dla całego zestawu danych.

Zakładam prosty liniowy model regresji wielowymiarowej, tj

y = A x + b + e

$\mathbf y = \mathbf A\mathbf x + \mathbf b + \mathbf e$

Jaki jest najlepszy algorytm do tworzenia stale aktualizowanej oceny parametrów regresji liniowej i ? $\mathbf A$ $\mathbf b$

Idealnie:

Chciałbym algorytm, który jest najbardziej złożonością przestrzenną i czasową na aktualizację, gdzie jest wymiarowością zmiennej niezależnej ( ), a jest wymiarowością zmiennej zależnej ( ). $\mathcal O(N\cdot M)$ $N$ $\mathbf x$ $M$ $\mathbf y$
Chciałbym móc określić jakiś parametr, aby określić, jak bardzo parametry są aktualizowane przez każdą nową próbkę, np. 0,000001 oznaczałoby, że następna próbka dostarczy jedną milionową szacunku parametru. Dałoby to pewnego rodzaju wykładniczy rozkład dla efektu próbek w odległej przeszłości.

— mikera
źródło

2

Wyszukaj (1) Elastyczna regresja liniowa, (2) Filtry Kalmana.

— Jase,

31

Maindonald opisuje sekwencyjną metodę opartą na obrotach Givensa . (Obrót Givensa jest ortogonalną transformacją dwóch wektorów, która zeruje dany wpis w jednym z wektorów.) W poprzednim etapie dokonałeś dekompozycji macierzy projektowej na macierz trójkątną za pomocą transformacji ortogonalnej tak że . (Uzyskanie wyników regresji z macierzy trójkątnej jest szybkie i łatwe). Po dołączeniu nowego wiersza poniżej skutecznie przedłużasz $\mathbf{X}$ $\mathbf{T}$ $\mathbf{Q}$ $\mathbf{Q}\mathbf{X} = (\mathbf{T}, \mathbf{0})'$ $v$ $\mathbf{X}$ Również przez niezerowy rząd, powiedz . Zadanie polega na wyzerowaniu tego rzędu przy jednoczesnym zachowaniu wpisów w pozycjiprzekątnejWykonuje to sekwencja rotacji Givens: rotacja z pierwszym rzędem zeruje pierwszy element ; następnie obrót z drugim rzędem zeruje drugi element i tak dalej. Efektem jest przedwczesne przez serię obrotów, co nie zmienia jego ortogonalności. $(\mathbf{T}, \mathbf{0})'$ $t$ $\mathbf{T}$ $\mathbf{T}$ $t$ $\mathbf{T}$ $\mathbf{Q}$

Gdy macierz projektowa ma kolumny (co ma miejsce w przypadku regresji na zmiennych plus stała), liczba potrzebnych obrotów nie przekracza a każdy obrót zmienia dwa wektory . Pamięć potrzebna dla to . Zatem ten algorytm ma koszt obliczeniowy zarówno w czasie, jak i przestrzeni. $p+1$ $p$ $p+1$ $p+1$ $\mathbf{T}$ $O((p+1)^2)$ $O((p+1)^2)$

Podobne podejście pozwala określić wpływ regresji na usunięcie wiersza. Maindonald podaje formuły; podobnie jak Belsley, Kuh i Welsh . Dlatego jeśli szukasz ruchomego okna do regresji, możesz zachować dane dla okna w okrągłym buforze, przylegając do nowego układu odniesienia i upuszczając stary z każdą aktualizacją. Podwaja to czas aktualizacji i wymaga dodatkowej pamięci dla okna o szerokości . Wydaje się, że będzie analogiem parametru wpływu. $O(k (p+1))$ $k$ $1/k$

W przypadku rozkładu wykładniczego uważam (spekulacyjnie), że można dostosować to podejście do ważonych najmniejszych kwadratów, nadając każdej nowej wartości wagę większą niż 1. Nie powinno być potrzeby utrzymywania bufora poprzednich wartości ani usuwania starych danych.

Bibliografia

JH Maindonald, Obliczenia statystyczne. J. Wiley & Sons, 1984. Rozdział 4.

DA Belsley, E. Kuh, RE Welsch, Diagnostyka regresji: identyfikacja wpływowych danych i źródeł kolinearności. J. Wiley & Sons, 1980.

— Whuber
źródło

1

Czy metoda opisana przez Maindonald jest powiązana z algorytmem Dżentelmena lub taka sama? jstor.org/stable/2347147

— onestop

6

W takim przypadku zobacz także rozszerzenia Alana Millera jstor.org/stable/2347583 . Archiwum swojej stronie Fortran oprogramowania jest teraz w jblevins.org/mirror/amiller

— onestop

5

Wyraźny algorytm pojawia się na dole strony p. 4 z saba.kntu.ac.ir/eecd/people/aliyari/NN%20%20files/rls.pdf . Można to znaleźć w Google „rekurencyjnych najmniejszych kwadratów”. Nie wygląda to na ulepszenie podejścia Gentleman / Maindonald, ale przynajmniej jest jasno i wyraźnie opisane.

— whuber

2

Ostatni link wygląda jak metoda, którą zamierzałem zaproponować. Używana przez nich tożsamość matrycy znana jest w innych miejscach jako tożsamość Sherman - Morrison - Woodbury. Jest także dość efektywny numerycznie, ale może nie być tak stabilny jak obrót Givens.

— kardynał

2

@ suncoolsu Hmm ... Książka Maindonalda została niedawno wydana, kiedy zacząłem z niej korzystać :-).

— whuber

8

Myślę, że przekształcenie modelu regresji liniowej w model przestrzeni stanów da ci to, czego szukasz. Jeśli używasz R, możesz użyć pakietu dlm i zajrzeć do książki towarzyszącej autorstwa Petrisa i in.

— F. Tusell
źródło

może jestem zdezorientowany, ale wydaje się, że odnosi się to do modelu szeregów czasowych? mój model jest rzeczywiście prostszy tym, że próbki nie są szeregi czasowe (w praktyce są one niezależne (X-> Y) próbek, są one po prostu gromadzą się w dużych ilościach w czasie)

— mikera

1

Tak, w ogólnym przypadku jest to stosowane do szeregów czasowych z nie niezależnymi obserwacjami; ale zawsze możesz założyć związek między kolejnymi obserwacjami, co daje ci szczególny interesujący cię przypadek.

— F. Tusell,

7

$E$ $W$

$E(i; W)$ $W$

W_{j} \leftarrow W_{j} + α \frac{\partial E (i; W)}{\partial W_{j}}

$W_{j} \leftarrow W_j + \alpha \frac{\partial{E(i; W)}}{\partial{W_j}}$

$\alpha$

Jest to bardzo wydajne i sposób, w jaki sieci neuronowe są zwykle szkolone. Możesz efektywnie przetwarzać nawet wiele próbek (powiedzmy 100 lub więcej).

Oczywiście można zastosować bardziej wyrafinowane algorytmy optymalizacji (pęd, gradient sprzężony, ...).

— bayerj
źródło

Wydaje się bardzo podobny do tego papieru eprints.pascal-network.org/archive/00002147/01/... . Został zaimplementowany w projekcie open source o nazwie jubatus.

— sacharyna

3

Zaskoczony, jak dotąd nikt inny tego nie dotykał. Regresja liniowa pełni kwadratową funkcję celu. Tak więc krok Newtona Raphsona z dowolnego punktu początkowego prowadzi prosto do optymów. Powiedzmy, że już wykonałeś regresję liniową. Funkcja celu to:

L (β) = (y - X β)^{t} (y - X β)

$L(\beta) = (y - X \beta)^t (y - X \beta)$

\nabla L (β) = - 2 X^{t} (y - X β)

$\nabla L (\beta) = -2 X^t (y - X \beta)$

\nabla^{2} L (β) = X^{t} X

$\nabla^2 L (\beta) = X^t X$

$\beta$ $x_{new}, y_{new}$

\nabla L_{n e w} (β) = - 2 x_{n e w} (y_{n e w} - x_{n e w}^{T} β)

$\nabla L_{new}(\beta) = -2 x_{new} (y_{new}-x_{new}^T \beta)$

\nabla^{2} L_{n e w} = x_{n e w} x_{n e w}^{T}

$\nabla^2 L_{new} = x_{new}x_{new}^T$

Dodaj to do starego hessianu podanego powyżej. Następnie zrób krok Newtona Raphsona.

β_{n e w} = β_{o l d} + (\nabla^{2} L)^{- 1} \nabla L_{n e w}

$\beta_{new} = \beta_{old} + (\nabla^2L)^{-1} \nabla L_{new}$

I jesteś skończony.

— ryu576
źródło

1

\nabla L_{n e w}

$\nabla L_{new}$

p,

$p,$

O (p^{3})

$O(p^3)$

O (p^{3})

$O(p^3)$

p

$p$

(I - A)^{- 1} = I + A + A^{2} + \dots

$(I-A)^{-1}=I+A+A^2+ \dots$

2

Standardowe dopasowanie najmniejszych kwadratów daje współczynniki regresji

$\beta = ( X^T X )^{-1} X^T Y$

$\beta$

$X^T X$ $X^T Y$ $M^2+M$ $\beta$

Na przykład, jeśli M = 1, wówczas jeden współczynnik wynosi

$\beta = \frac{\sum_{i=1}^N{x_i y_i}}{\sum_{i=1}^N{x_i^2}}$

więc za każdym razem, gdy otrzymujesz nowy punkt danych, aktualizujesz obie sumy i obliczasz współczynnik i otrzymujesz zaktualizowany współczynnik.

$X^T X$ $X^T Y$ $(1-\lambda)$ $\lambda$

— Mark Higgins
źródło

2

β

$\beta$

X^{T} X

$X^T X$

X^{T} Y

$X^T Y$

6

X^{'} X

$X'X$

1

C^{- 1} x

$C^{-1}x$

C

$C$

x

$x$

z_{t + 1} = z_{t} + x - C z_{t}

$z_{t+1}=z_t + x - Cz_t$

z \to C^{- 1} x

$z\to C^{-1}x$

t \to \infty

$t\to\infty$

1

Problem jest łatwiejszy do rozwiązania, gdy przepisujesz trochę rzeczy:

Y = y

X = [x, 1]

następnie

Y = A * X

Obliczenie pozwala znaleźć rozwiązanie jednorazowe

V = X '* X

i

C = X '* Y

zauważ, że V powinien mieć rozmiar N-by-N, a C rozmiar N-by-M. Parametry, których szukasz, są następnie podawane przez:

A = inv (V) * C

Ponieważ zarówno V, jak i C są obliczane poprzez zsumowanie danych, możesz obliczyć A dla każdej nowej próbki. Ma to jednak złożoność czasową O (N ^ 3).

Ponieważ V jest kwadratowe i półokreślone dodatnio, istnieje rozkład LU, co sprawia, że odwracanie V jest bardziej stabilne numerycznie. Istnieją algorytmy do przeprowadzania aktualizacji rangi 1 odwrotności macierzy. Znajdź je, a uzyskasz efektywne wdrożenie, którego szukasz.

Algorytmy aktualizacji rangi 1 można znaleźć w „Obliczeniach macierzy” Goluba i van Loana. Jest to trudny materiał, ale ma kompleksowy przegląd takich algorytmów.

Uwaga: Powyższa metoda daje oszacowanie metodą najmniejszych kwadratów na każdym etapie. Możesz łatwo dodawać wagi do aktualizacji X i Y. Gdy wartości X i Y stają się zbyt duże, można je skalować za pomocą jednego skalara, bez wpływu na wynik.

— Mr. White
źródło