Comment corréler deux séries temporelles avec des lacunes et des bases de temps différentes?

10

J'ai posé cette question sur StackOverflow et on m'a recommandé de la poser ici.


J'ai deux séries chronologiques de données d'accéléromètre 3D qui ont des bases de temps différentes (horloges démarrées à des moments différents, avec un très léger fluage pendant la période d'échantillonnage), ainsi que contenant de nombreuses lacunes de taille différente (en raison des retards associés à l'écriture pour séparer périphériques flash).

Les accéléromètres que j'utilise sont les GCDC X250-2 bon marché . J'exécute les accéléromètres à leur gain le plus élevé, donc les données ont un bruit de fond important.

Les séries chronologiques ont chacune environ 2 millions de points de données (sur une heure à 512 échantillons / s) et contiennent environ 500 événements d'intérêt, où un événement typique s'étend sur 100 à 150 échantillons (200 à 300 ms chacun). Beaucoup de ces événements sont affectés par des interruptions de données lors des écritures flash.

Ainsi, les données ne sont pas vierges et ne sont même pas très jolies. Mais mon inspection du globe oculaire montre qu'il contient clairement les informations qui m'intéressent. (Je peux publier des tracés, si nécessaire.)

Les accéléromètres sont dans des environnements similaires mais ne sont que modérément couplés, ce qui signifie que je peux voir à l'œil les événements qui correspondent à chaque accéléromètre, mais je n'ai pas réussi jusqu'à présent à le faire dans le logiciel. En raison de limitations physiques, les appareils sont également montés dans différentes orientations, où les axes ne correspondent pas, mais ils sont aussi proches de l'orthogonal que je pourrais les faire. Ainsi, par exemple, pour les accéléromètres 3 axes A & B, + Ax correspond à -By (haut-bas), + Az correspond à -Bx (gauche-droite) et + Ay correspond à -Bz (recto-verso) .

Mon objectif initial est de corréler les événements de choc sur l'axe vertical, bien que j'aimerais finalement a) découvrir automatiquement la cartographie des axes, b) corréler l'activité sur les as cartographiés, et c) extraire les différences de comportement entre les deux accéléromètres (comme la torsion) ou flexion).

La nature des données de séries chronologiques rend numpy.correlate () de Python inutilisable. J'ai également examiné le package R's Zoo, mais je n'ai pas fait de progrès avec. J'ai cherché de l'aide dans différents domaines de l'analyse du signal, mais je n'ai fait aucun progrès.

Quelqu'un a-t-il des indices sur ce que je peux faire ou des approches que je devrais rechercher?

Mise à jour le 28 février 2011: Ajout de quelques parcelles ici montrant des exemples de données.

BobC
la source
1
@BobC, l'un des modérateurs peut peut-être faire migrer votre message vers ce site. Ce serait le plus raisonnable. Quant à vos questions techniques, tout d'abord, utilisez-vous la FFT pour faire la corrélation? Cela devrait être possible pour 2 millions de points de données sur un ordinateur à moitié décent. Votre rapport signal / bruit semble relativement élevé, vous devriez donc être en affaires. Une coupe rapide et sale consisterait à remplir les données manquantes soit avec le dernier échantillon disponible, soit avec des zéros. Le fluage des différences d'intervalles d'échantillonnage peut être la "caractéristique" la plus difficile à gérer de vos données.
Cardinal
@cardinal: J'ai effectivement essayé une FFT, seulement pour obtenir des ordures. Les caractéristiques «intéressantes» facilement visibles dans les données sont indiscernables du bruit dans la FFT. Cependant, je n'ai fait des FFT que sur l'ensemble des données: une FFT à fenêtre mobile fournirait peut-être de meilleurs résultats, mais je n'ai pas encore trouvé de moyen efficace sur le plan informatique de l'implémenter. Je soupçonne qu'une transformation en ondelettes pourrait aider, mais je ne la connais pas (mais j'apprends lentement à ce sujet).
BobC
1
@BobC, ce que je voulais dire, avez-vous envisagé une implémentation basée sur la FFT pour calculer la corrélation? La convolution directe est , mais une implémentation basée sur FFT réduirait cet , ce qui la rendrait possible. Quant à regarder la FFT elle-même, avec 2 millions de points de données, votre résolution en fréquence sera très élevée. Tout fluage d'échantillonnage et autres éléments sont tenus d'éliminer le signal par fréquence. Mais, vous devriez pouvoir agréger sur de nombreux bacs pour faire sortir le signal du bruit. Quelque chose comme une approche Welch ou peut-être une technique de fenêtrage personnalisée. O(n2)O(nlogn)
Cardinal
@BobC, du haut de ma tête, il semble qu'une variante d'un algorithme de chevauchement et d'ajout ou de chevauchement et d'enregistrement pourrait être utilisée pour faire une FFT à fenêtre coulissante. Faire glisser les échantillons dans une fenêtre équivaut simplement à un décalage de phase, donc tout ce que vous avez à faire est de compenser les échantillons qui «tombent» à l'extrémité gauche et ceux qui «entrent» à l'extrémité droite.
Cardinal
Salut, j'ai une question similaire. J'ai 2 séries temporelles, chacune représentée par une matrice avec sa première colonne correspondant aux valeurs et la deuxième colonne correspondant à la différence de temps (depuis la valeur précédente) Comment trouver la corrélation entre ces 2 matrices? J'ai essayé de faire xcorr2 () mais cela ne semble pas correct et faire xcorr calculerait probablement la corrélation avec seulement les valeurs à considérer, mais je veux aussi tenir compte du temps. Je suis vraiment confus ici, une FFT va-t-elle aider? Comment proposeriez-vous que je m'y prenne?

Réponses:

12

La question concerne le calcul de la corrélation entre deux séries temporelles échantillonnées irrégulièrement (processus stochastiques unidimensionnels) et son utilisation pour trouver le décalage temporel où elles sont corrélées au maximum (leur "différence de phase").

Ce problème n'est généralement pas traité dans l'analyse des séries chronologiques, car les données des séries chronologiques sont supposées être collectées systématiquement (à intervalles réguliers). C'est plutôt le domaine de la géostatistique , qui concerne les généralisations multidimensionnelles des séries chronologiques. L'ensemble de données géostatistiques archétypales consiste en des mesures d'échantillons géologiques à des emplacements irrégulièrement espacés.

Avec un espacement irrégulier, les distances entre paires d'emplacements varient: deux distances peuvent ne pas être identiques. La géostatistique surmonte cela avec le variogramme empirique . Ceci calcule une valeur "typique" (souvent la moyenne ou la médiane) de 2/2 - la "semi-variance" - où dénote une valeur mesurée au point et la distance entre et est contrainte de se situer dans un intervalle appelé "décalage". Si nous supposons que le processus est stationnaire et a une covariance, alors l'espérance de la semi-variance est égale à la covariance maximale (égale à pour tout ) moins la covariance entre(z(p)z(q))2/2z(p)ppqZVar(Z(p))pZ(p) et . Ce regroupement en décalages résout le problème d'espacement irrégulier.Z(q)

Lorsqu'une paire ordonnée de mesures est effectuée à chaque point, on peut calculer de façon similaire le variogramme croisé empirique entre les et les et ainsi estimer la covariance à n'importe quel décalage . Vous voulez la version unidimensionnelle du variogramme croisé. Les packages R gstat et sgeostat , entre autres, estimeront les variogrammes croisés. Ne vous inquiétez pas que vos données sont unidimensionnelles; si le logiciel ne fonctionne pas directement avec eux, introduisez simplement une deuxième coordonnée constante: cela les fera apparaître en deux dimensions.(z(p),w(p))zw

Avec deux millions de points, vous devriez pouvoir détecter de petits écarts par rapport à la stationnarité. Il est possible que la différence de phase entre les deux séries temporelles varie également dans le temps. Faire face à cela en calculant le variogramme croisé séparément pour différentes fenêtres espacées tout au long de la période.

@cardinal a déjà évoqué la plupart de ces points dans les commentaires. La principale contribution de cette réponse est de pointer vers l'utilisation de progiciels de statistiques spatiales pour faire votre travail pour vous et d'utiliser des techniques de géostatistique pour analyser ces données. En ce qui concerne l'efficacité du calcul, notez que la convolution complète (cross-variogramme) n'est pas nécessaire: vous n'avez besoin que de ses valeurs proches de la différence de phase. Cela fait l'effort , pas , où est le nombre de retards à calculer, donc cela pourrait être faisable même avec un logiciel prêt à l'emploi. Sinon, l'algorithme de convolution directe est facile à mettre en œuvre.O(nk)O(n2)k

whuber
la source
@whuber, bons commentaires et suggestions. Si je lis bien la question, je pense qu'une des principales préoccupations est l'incertitude du moment de l'échantillonnage. Cela peut être un peu différent du cadre géostatistique typique où, je crois, l'espacement est irrégulier, mais toujours supposé connu (au moins avec une grande précision). Je pense qu'un modèle approximatif est, si le ème point de la série un est au temps , pour fixe alors le ème point de la série 2 est à où est probablement de l'ordre de quelques millisecondes et est probablement minuscule. ntn=nttnτn=tn+α+βnαβ
Cardinal
@Cardinal Je n'ai pas obtenu cela de la question. Je ne peux pas penser à un moyen d'estimer la qui ne serait pas assez intensif en calcul. Peut-être en divisant la série chronologique en groupes où l'effet net de est négligeable? ββ
whuber
@whuber, @BobC, je fais une supposition semi-éclairée basée sur l'expérience passée en traitant des problèmes et des problèmes similaires. La plupart des approches que j'ai vues nécessitent beaucoup de calculs et ne parviennent pas à impressionner. Une tentative pourrait être via quelque chose comme la déformation temporelle dynamique ou ce que Ramsay et Silverman appellent l' enregistrement de courbe . Je ne sais pas si l'une ou l'autre de ces données serait réalisable avec cet ensemble de données de taille.
Cardinal
Il me faudra un peu pour que mon cerveau s'enroule autour de ça. Je vais commencer par les exemples dans les packages R que vous avez mentionnés.
BobC
@BobC, le modèle approximatif que j'ai donné pour la synchronisation asynchronicité est-il proche de ce que vous avez? Je pense que c'est "décalage initial aléatoire" + "erreur linéaire" où ce dernier est dû à une petite différence constante dans l'intervalle d'échantillonnage entre vos deux appareils. Ensuite, il y a une petite erreur aléatoire supplémentaire, par exemple en raison de l'interruption du traitement de deux uC différents en plus de cela.
Cardinal