Régression linéaire en ligne efficace

53

J'analyse des données pour lesquelles je souhaite effectuer une régression linéaire ordinaire. Toutefois, cela n’est pas possible car je traite d’un paramètre en ligne avec un flux continu de données d’entrée (qui deviendra rapidement trop volumineux pour la mémoire). pour mettre à jour les estimations de paramètres pendant la consommation. c'est-à-dire que je ne peux pas simplement tout charger en mémoire et effectuer une régression linéaire sur l'ensemble des données.

Je suppose un modèle de régression linéaire multivarié simple, à savoir

y=Ax+b+e

Quel est le meilleur algorithme pour créer une estimation constamment mise à jour des paramètres de régression linéaire et ?bAb

Idéalement:

  • Je voudrais un algorithme qui soit le plus complexité de l'espace et du temps par mise à jour, où est la dimensionnalité de la variable indépendante ( ) et la dimensionnalité de la variable dépendante ( ).N x M yO(NM)NxMy
  • J'aimerais pouvoir spécifier un paramètre pour déterminer le niveau de mise à jour des paramètres pour chaque nouvel échantillon, par exemple 0,000001 signifierait que le prochain échantillon fournirait un millionième de l'estimation du paramètre. Cela donnerait une sorte de décroissance exponentielle pour l’effet des échantillons du passé lointain.
Mikera
la source
2
Rechercher (1) Régression linéaire flexible, (2) Filtres de Kalman.
Jase

Réponses:

31

Maindonald décrit une méthode séquentielle basée sur les rotations de Givens . (Une rotation de Givens est une transformation orthogonale de deux vecteurs qui supprime une entrée donnée dans l'un des vecteurs.) A l'étape précédente, vous avez décomposé la matrice de conception en une matrice triangulaire via un transformation orthogonale pour que . (Il est rapide et facile d'obtenir les résultats de régression à partir d'une matrice triangulaire.) Si vous ajoutez une nouvelle ligne sous , vous prolongez d'une ligne non nulle. aussi, disT Q Q X = ( T , 0 ) v X ( T , 0 ) t T T t T QXTQQX=(T,0)vX(T,0)t. La tâche consiste à mettre à zéro cette ligne tout en conservant les entrées dans la position de diagonale. Une séquence de rotations de Givens fait ceci: la rotation avec la première ligne de zéros le premier élément de ; puis la rotation avec la deuxième ligne de le deuxième élément à zéro , et ainsi de suite. L’effet est de prémultiplier par une série de rotations, ce qui ne change pas son orthogonalité.TTtTQ

Lorsque la matrice de conception a colonnes (ce qui est le cas lors d'une régression sur variables plus une constante), le nombre de rotations nécessaires ne dépasse pas et chaque rotation modifie deux vecteurs . Le stockage nécessaire pour est . Ainsi, cet algorithme a un coût de calcul de dans le temps et dans l’espace.p p + 1 p + 1 T O ( ( p + 1 ) 2 ) O ( ( p + 1 ) 2 )p+1pp+1p+1TO((p+1)2)O((p+1)2)

Une approche similaire vous permet de déterminer l’effet de la suppression d’une ligne sur la régression. Maindonald donne des formules; donc faire Belsley, Kuh, et gallois . Ainsi, si vous recherchez une fenêtre mobile pour la régression, vous pouvez conserver les données de cette fenêtre dans un tampon circulaire, en joignant le nouveau datum et en supprimant l’ancien à chaque mise à jour. Cela double le temps de mise à jour et nécessite un stockage supplémentaire de pour une fenêtre de largeur . Il semble que serait l'analogue du paramètre d'influence.k 1 / kO(k(p+1))k1/k

Pour la décroissance exponentielle, je pense (de manière spéculative) que vous pouvez adapter cette approche aux moindres carrés pondérés, en attribuant à chaque nouvelle valeur un poids supérieur à 1. Il ne devrait pas être nécessaire de conserver une mémoire tampon des valeurs précédentes ou de supprimer des données anciennes.

Références

JH Maindonald, Calcul statistique. J. Wiley & Sons, 1984. Chapitre 4.

DA Belsley, E. Kuh, RE Welsch, Diagnostics de régression: identification des données d'influence et des sources de colinéarité. J. Wiley & Sons, 1980.

whuber
la source
1
La méthode décrite par Maindonald est-elle liée ou identique à l'algorithme de Gentleman? jstor.org/stable/2347147
onestop
6
Dans ce cas, voir également les extensions d’Alan Miller jstor.org/stable/2347583 . Une archive de son site logiciel Fortran se trouve maintenant sur jblevins.org/mirror/amiller
onestop le
5
Un algorithme explicite apparaît au bas de p. 4 de saba.kntu.ac.ir/eecd/people/aliyari/NN%20%20files/rls.pdf . Ceci peut être trouvé en googlant "les moindres carrés récursifs". Cela ne semble pas être une amélioration de l’approche Gentleman / Maindonald, mais au moins il est clairement et explicitement décrit.
whuber
2
Le dernier lien ressemble à la méthode que je voulais suggérer. L'identité matricielle qu'ils utilisent est connue ailleurs: identité de Sherman - Morrison - Woodbury. Son implémentation est également assez efficace numériquement, mais peut ne pas être aussi stable qu'une rotation de Givens.
cardinal
2
@suncoolsu Hmm ... Le livre de Maindonald a été publié récemment lorsque j'ai commencé à l'utiliser :-).
whuber
8

Je pense que la refonte de votre modèle de régression linéaire dans un modèle d’espace à états vous donnera ce que vous cherchez. Si vous utilisez R, vous voudrez peut-être utiliser le package dlm et consulter le livre de référence de Petris et al.

F. Tusell
la source
peut-être suis-je confus mais cela semble faire référence à un modèle de série chronologique? mon modèle est en réalité plus simple en ce sens que les échantillons ne sont pas une série chronologique (en fait, ce sont des échantillons indépendants (x-> y), ils sont simplement accumulés dans de grands volumes au fil du temps)
mikera
1
Oui, dans le cas général, il est utilisé pour les séries chronologiques avec des observations non indépendantes. mais vous pouvez toujours supposer une corrélation entre des observations successives, ce qui vous donne le cas particulier qui vous intéresse.
F. Tusell
7

Vous pouvez toujours effectuer une descente de gradient sur la somme des carrés coût WRT les paramètres de votre modèle . Il suffit d’en prendre le gradient mais ne pas opter pour la solution de formulaire fermé, mais uniquement pour le sens de la recherche.EW

Laissez le coût de l'échantillon de formation i'th compte tenu des paramètres . Votre mise à jour pour le ième paramètre est alorsE(i;W)W

WjWj+αE(i;W)Wj

où est un taux d'échelon, que vous devez choisir via une validation croisée ou une bonne mesure.α

C'est très efficace et la manière dont les réseaux de neurones sont formés. Vous pouvez traiter même de nombreux échantillons en parallèle (par exemple, environ 100) de manière efficace.

Bien entendu, des algorithmes d'optimisation plus sophistiqués (momentum, gradient conjugué, ...) peuvent être appliqués.

Bayerj
la source
Cela semble très semblable à cet article eprints.pascal-network.org/archive/00002147/01/… . Il a été implémenté dans un projet open source appelé jubatus.
saccharine
3

Surpris, personne d'autre n'en a encore parlé. La régression linéaire a une fonction objective quadratique. Ainsi, une étape newton Raphson à partir de n'importe quel point de départ vous mène directement aux optima. Supposons maintenant que vous ayez déjà effectué votre régression linéaire. La fonction objectif est:

L(β)=(yXβ)t(yXβ)
Le gradient devient et le hessian:
L(β)=2Xt(yXβ)
2L(β)=XtX

Maintenant, vous avez des données passées, fait une régression linéaire et vous associez à vos paramètres ( ). Le gradient à ce stade est zéro par définition. Le hessian est comme donné ci-dessus Un nouveau point de données ( ) arrive. Vous venez de calculer le gradient pour le nouveau point via:βxnew,ynew

Lnew(β)=2xnew(ynewxnewTβ)
et deviendra votre dégradé global (puisque le dégradé à partir des données existantes était nul) . Le hessian pour le nouveau point de données est:

2Lnew=xnewxnewT
.

Ajoutez ceci à la vieille Hesse donnée ci-dessus. Ensuite, il suffit de faire un pas de Newton Raphson.

βnew=βold+(2L)1Lnew

Et tu as fini.

ryu576
la source
1
J'aime l'idée pour sa simplicité, mais (a) pour ne pas dérouter les lecteurs, préférerais voir une définition explicite de " " et (b) croire que vous devez calculer le gradient correctement (ou afficher pourquoi être un facteur de 2 importe peu). Ce serait plus convaincant si vous pouviez donner un petit exemple démontrant que c'est correct. Pour plus grand il serait intéressant d'estimer l'effort de calcul. Inverser la Hesse ne prend-il pas temps? Lnewp,O(p3)
whuber
Merci, ajoutera plus de détails un peu plus tard aujourd'hui. Oui, inverser le hessian prend pour grand . Vous pouvez également essayer de maintenir l'inverse de hesse et de le mettre à jour directement à l'aide de la série de puissances ( ). Si vous avez un million de paramètres, la descente de gradient est plus ou moins votre seule option. p ( I - A )O(p3)p(IA)1=I+A+A2+
ryu576
2

L’ajustement des moindres carrés standard donne les coefficients de régression

β=(XTX)1XTY

β

XTXXTYM2+Mβ

Par exemple, si M = 1, le coefficient 1 est

β=i=1Nxiyii=1Nxi2

ainsi, chaque fois que vous obtenez un nouveau point de données, vous mettez à jour les deux sommes et calculez le ratio et vous obtenez le coefficient mis à jour.

XTXXTY(1λ)λ

Mark Higgins
la source
2
C'est bien de voir ce cas simple expliqué. Avez-vous remarqué, cependant, que la question concerne spécifiquement la régression multivariée ? Ce n'est pas si facile de mettre à jour le dénominateur de dans ce cas! β
whuber
Je pense que ma réponse fonctionne toujours: vous devez garder une trace de la matrice MxM et de la matrice Mx1 . Chaque élément de ces matrices est une somme comme dans l'exemple M = 1. Ou est-ce que je manque quelque chose? XTXXTY
Mark Higgins
6
Oui: en plus du calcul d'un produit de matrice et de l'application d'une matrice à un vecteur, vous devez maintenant inverser à chaque étape. C'est cher. Le but des algorithmes en ligne est de remplacer les étapes coûteuses en gros par des procédures de mise à jour moins coûteuses. XX
whuber
1
@whuber Soit dit en passant, Schraudolph, NN (2002) donne un moyen rapide et en ligne d'estimer pour une matrice et un vecteur changeants . Produits matriciels à courbure rapide pour la descente de gradient de second ordre. Vous laissez et comme . C x z t + 1 = z t + x - C z t z C - 1 x t C1xCxzt+1=zt+xCztzC1xt
Neil G
1

Le problème est plus facilement résolu lorsque vous réécrivez un peu les choses:

Y = y

X = [x, 1]

ensuite

Y = A * X

Une solution ponctuelle est trouvée en calculant

V = X '* X

et

C = X '* Y

notez que le V doit avoir une taille N-by-N et C une taille de N-M. Les paramètres que vous recherchez sont alors donnés par:

A = inv (V) * C

Puisque V et C sont calculés en faisant la somme de vos données, vous pouvez calculer A à chaque nouvel échantillon. Cela a cependant une complexité temporelle de O (N ^ 3).

Puisque V est carré et positif semi-défini, il existe une décomposition de LU qui rend l’inversion de V numériquement plus stable. Il existe des algorithmes pour effectuer des mises à jour de rang 1 à l'inverse d'une matrice. Trouvez-les et vous obtiendrez la mise en œuvre efficace que vous recherchez.

Les algorithmes de mise à jour de rang 1 se trouvent dans "Calculs matriciels" de Golub et van Loan. C'est un matériau difficile, mais il offre une vue d'ensemble complète de ces algorithmes.

Remarque: La méthode ci-dessus donne une estimation des moindres carrés à chaque étape. Vous pouvez facilement ajouter des poids aux mises à jour de X et Y. Lorsque les valeurs de X et Y deviennent trop grandes, elles peuvent être mises à l'échelle par un seul scalaire, sans affecter le résultat.

M. White
la source