El análisis de datos caóticos
Asociación
EURATOM-CIEMAT
Para Fusión
• Esta sección trata del análisis de datos con características de caos
• Antes de comenzar, hay que enfatizar que los experimentos de investigación de
caos son experimentos muy controlados, donde se ha puesto mucho empeño en
eliminar el ruido y las perturbaciones externas, y generalmente de sistemas
relativamente sencillos (imanes, haces de láser, etc.).
• Las mismas técnicas de análisis se aplican también a la turbulencia, pero hay que
tener en cuenta que la turbulencia es mucho más compleja que un simple sistema
caótico.
• Por ejemplo, un péndulo con 3 grados de libertad puede ser un sistema caótico.
• En un fluido, sin embargo, hay N partículas en un volúmen dado, y todas tienen
varios grados de libertad (6: para su posición y su velocidad, y posiblemente más).
N es inmensamente grande (número de Avogadro), y las interacciones entre
partículas aumentan la complejidad.
• Por tanto, un sistema turbulento (dimensión casi ∞)es algo muy diferente a los
sistemas caóticos sencillos típicos que se estudian en la teoría del caos (dimensión
O(1)).
• Sin embargo, se pueden aplicar técnicas similares de análisis (hasta cierto punto).
1
El análisis de datos caóticos
Asociación
EURATOM-CIEMAT
Para Fusión
• El análisis de datos caóticos requiere nuevos métodos de análisis. Los métodos
usuales para sistemas lineales (especialmente Fourier) sirven de muy poco para
sistemas no-lineales.
• Los datos que uno obtiene de un sistema caótico/turbulento suelen ser
incompletos:
• Una serie temporal en un punto proveniente de un sistema con una
extensión en 3 dimensiones
• Una variable de las N que constituyen el sistema.
• El objetivo del análisis es:
• Clasificar el sistema
• Predecir (aunque sea de manera estadística) su comportamiento futuro
• Construir un modelo del sistema y estimar sus parámetros a partir de los
datos
• Si se logra construir un modelo, podemos decir que hemos “entendido” el
sistema.
2
Reconstrucción del
espacio de fase
Asociación
EURATOM-CIEMAT
Para Fusión
• El espacio de fase tiene tantas dimensiones como sean necesarias para definir el
estado de un sistema completamente (cada estado corresponde a 1 sólo punto).
• Considera un fluido en un volúmen V. Su movimiento se describe con la ecuación
de Navier-Stokes. El espacio de fase tiene esencialmente infinitas dimensiones (6N,
donde N es del orden del número de Avogadro).
• ¿Cómo es posible que se pueda obtener alguna información valiosa de medidas
escalares en uno o varios puntos? ¿No serían necesarías infinitas medidas?
• Debido a que el sistema es disipativo (generalmente hablando), el movimiento no
cubre todo el espacio de fases, sino que se limita esencialmente a mover a lo largo de
un atractor de dimensión mucho inferior. En algunos casos, esto puede reducir la
dimensión efectiva del sistema a valores manejables. En fórmula:
F (u ) 
u ( x , t )
Sistema
Determinante
D  det 
 F (u ( x , t ))

dinámico
del Jacobiano
t
 u 
• D = 1: sistema Hamiltoniano (conserva volúmen, por ejemplo líneas de campo
con B = 0; caso muy particular e interesante)
• D < 1: sistema disipativo (la situación usual)
• D > 1: sistema no físico, imposible porque gana energía; o sistema no cerrado
3
Programa
Asociación
EURATOM-CIEMAT
Para Fusión
En lo que sigue, vamos a tratar la manera de extraer información sobre el sistema
caótico subyacente, suponiendo que sólo tenemos unas medidas escalares (en uno
o varios puntos) a nuestra disposición. Los pasos a seguir son:
• Hallar la dimensión de embedding (del espacio que contiene el atractor; la cual
tiene una dimensión superior a la dimensión fractal del atractor) para poder
reconstruir el atractor.
• Determinar los invariantes de la dinámica del sistema (dimensiones fractales,
exponentes de Lyapunov).
• Construcción de modelos para el sistema
• Separación de la señal del ruido
4
Uso de coordenadas retrasadas
Asociación
EURATOM-CIEMAT
Para Fusión
• Partimos de una señal escalar: s(n) = s(t0 + nts), donde ts es el intervalo de
muestreo.
• La reconstrucción de la órbita en el espacio de fases requiere pintar s(n) frente a
ds(n)/dt, d2s(n)/dt2, etc.
• El cálculo de derivadas discretas involucra diferencias finitas:
ds(n)/dt ≈ s(n+1)–s(n); d2s(n)/dt2 ≈ s(n+1)–2s(n)+s(n-1)
• Cada diferencia de un orden superior involucra valores del índice más alejados.
• Pero también cada diferencia finita introduce un error, y amplifica el ruido (la
diferencia finita es como un filtro paso-alto).
• Para evitar eso, se ha introducido el concepto de “coordenadas retrasadas”
(delay-time coordinates):
• En lugar de usar s(n), ds(n)/dt, d2s(n)/dt2, etc. para contruir el espacio de fases, se
usa s(n), s(n+T), s(n+2T), etc.
• Ventajas: muy fácil y no introduce errores.
• Problema: hay que determinar el valor de T que da el mejor resultado (no
necesariamente es igual al intervalo de muestreo).
5
Uso de coordenadas retrasadas
Asociación
EURATOM-CIEMAT
Para Fusión
Usando coordenadas retrasadas, la “órbita” en el “espacio de fase” está dada por:
y (n)   s(n), s(n  T ), s(n  2 T ), ..., s(n  (d  1)T )
donde d es la dimensión del espacio embedding (por determinar).
• Ejemplo: el atractor de Lorenz:
6
Embedding dimension
Asociación
EURATOM-CIEMAT
Para Fusión
• La elección de la dimensión de embedding d, que identificaremos con dE, no es
muy crítica.
• La literatura sugiere que cualquier valor para dE que es mayor que 2 dA servirá
para “desenvolver” el atractor, donde dA es la dimensión (fractal) del atractor.
• Otra cuestión que uno puede plantearse es: admitiendo este procedimiento como
válido, ¿cómo sabemos que la órbita y(n) tiene relación alguna con el sistema
físico que estudiamos?
• Lo sabemos con seguridad porque la relación entre las variables físicas u(n) y
las coordenadas retrasadas y(n) consiste en un cambio de variables lineal y suave;
es decir podemos estudiar y(n) en lugar de u(n) para obtener los mismos
invariantes del movimiento (dimensiones fractales, exponentes de Lyapunov,
etc.).
• Este es un resultado fundamental para todo el estudio de fractales (Mañé,
Takens, Sauer).
7
Elección del tiempo de retraso
Asociación
EURATOM-CIEMAT
Para Fusión
• La elección del tiempo de retraso para las coordenadas retrasadas es importante:
• Demasiado corto: no hay información independiente (las dos coordenadas
no forman una base)
• Demasiado largo: ausencia total de correlación (la órbita es aleatoria)
• A primera vista, este argumento sugiere el uso de la correlación lineal (elegir el
tiempo de decorrelación – tiempo en que la correlación ha caído con un factor e).
• Sin embargo, la correlación lineal no toma en cuenta una gran parte de la
información que podemos obtener del sistema (en particular, la parte no-lineal).
• Por tanto, es mejor usar una generalización de la correlación lineal a sistemas
no-lineales: la información mutua.
8
La información mutua
Asociación
EURATOM-CIEMAT
Para Fusión
• Considera 2 sistemas, A y B.
• Se efectúan medidas y se obtienen 2 series de medidas ak y bk.
• Ahora se pueden calcular las distribuciones de probabilidad PA(a) yPB (b), así
como la distribución de probabilidad conjunta o bidimensional PAB (a,b) – la
probabilidad de obtener un valor a simultáneamente con un valor b.
• La cantidad de información binaria (bits) que podemos saber/predecir sobre el
sistema A midiendo sólo el sistema B, y conociendo las distribuciones de
probabilidad, es:
 PAB (a i , b k ) 
I AB (a i , b k )  log 2 

PA ( a i )PB ( b k ) 
• La cantidad promedia de información mutua es el promedio sobre todas las
combinaciones ai y bk posibles:
I AB (T ) 
 PAB ( ai , b k )I AB (a i , b k )
ai ,b k
9
La información mutua
Asociación
EURATOM-CIEMAT
Para Fusión
• Para medidas de un sólo sistema, ai es s(i) y bk es s(k) = s(i+T).
• En otras palabras, la información mutua promedia es:
N
I (T ) 

n 1
 P ( s(n), s(n  T )) 
P (s( n), s(n  T )) log 2 

P (s(n )) P (s( n  T )) 
• Máximo en retraso 0
• Cae a cero para retrasos grandes
• Retraso óptimo: es un punto
intermedio; un punto donde se
alcanza una cierta independencia
estadística para poder construir una
base
“Reglas de oro”:
• Posición del primer mínimo ó
• Posición donde el valor ha caído
hasta 1/5 del máximo.
10
La dimensión del espacio de
embedding
Asociación
EURATOM-CIEMAT
Para Fusión
• La elección de un espacio suficientemente grande para representar el atractor es
fundamental por la siguiente razón.
• Si proyectamos un objeto (el atractor) que tiene una cierta dimensión en un
espacio demasiado pequeño, entonces algunos puntos del objeto que realmente
están muy alejados unos de otros aparecerán como cercanos (“vecinos
equivocados”):
• órbita resuelta correctamente en
un espacio suficientemente
grande
• órbita proyectada en un espacio
demasiado pequeño: el objeto se
intersecta a sí mismo donde no
debe y parece más complicado
de lo que es en realidad
11
Determinar el espacio de
embedding
Asociación
EURATOM-CIEMAT
Para Fusión
• Además de la indicación proporcionada por la Información Mutua, disponemos de
métodos más sofisticados para establecer cuál es el espacio de embedding correcto.
• Estos métodos se fundamentan en la idea que alguna propiedad del atractor
cambia al pasar de una dimensión inferior a una superior, pero cuando se alcanza la
dimensión de embedding ya no varía más. Las dimensiones “sobrantes” se quedan
esencialmente vacíos.
• Los siguientes son dos ejemplos de métodos de este tipo:
• El primero consiste en ejecutar el algoritmo de Grassberger y Procaccia para la
determinación de la dimensión fractal (u otras características del atractor) para
varias dimensiones de embedding. Cuando la dimensión obtenida se estabiliza, se
ha llegado a la dimension embedding.
• El segundo consiste en determinar si los puntos, cercanos en una dimensión, son
lejanos en la siguiente. Cuando ya no hay más puntos de este tipo, se ha llegado a la
dimension embedding. Ahora describiremos esta técnica.
12
Asociación
EURATOM-CIEMAT
Para Fusión
Vecinos cercanos equivocados
• La técnica consiste en aumentar la dimensión de embedding y observar los
puntos que, cercanos en dimensión d, ya no lo están en dimensión d+1.
• En dimensión d, cada punto:
y (k )   s(k ), s( k  T ), ..., s(k  (d  1)T )
tiene un vecino (el que está más cerca), digamos y(j). La distancia entre ellos es:
R d ( k )  y (k )  y ( j )
Al pasar a la siguiente dimensión, la distancia cambia a:
R d 1 ( k )  Rd (k )  s(k  dT )  s( j  dT )
2
2
2
Si este cambio de distancia (normalizado) es mayor que un cierto umbral, entonces
decimos que en dimensión d los puntos eran “vecinos falsos”, ya que se han
alejado en dimensión d+1. En otras palabras, si se cumple:
s( k  dT )  s( j  dT )
Rd ( k )
 RT
los “vecinos” son “falsos”.
13
Vecinos cercanos equivocados
Asociación
EURATOM-CIEMAT
Para Fusión
• El valor del umbral RT resulta no ser muy crítico, un valor entre 10 y 50 funciona
generalmente.
• Ejemplo con el atractor de Lorenz (la escala vertical da el porcentaje de “vecinos
falsos” entre todos los puntos):
Obsérvese que el valor verdadero
de la dimensión del atractor de
Lorenz es 2.06, así que este
algoritmo detecta muy bien que a
partir de la dimensión 3 el atractor
“cabe” en el espacio de
embedding.
14
Efecto de ruido en el algoritmo
de vecinos falsos
Asociación
EURATOM-CIEMAT
Para Fusión
• Lo que llamamos “ruido” es una señal con un espectro de Fourier plano (o muy
ancho) y un tiempo de autocorrelación muy corto.
• Sin embargo, lo mismo podemos decir de muchas señales caóticas, y por tanto eso
no nos permite distinguir entre las dos.
• Otra característica del ruido es que es una señal con una dimensionalidad muy
alta.
• En el algoritmo de vecinos cercanos, esto tiene como consecuencia que en cada
nueva dimensión los puntos se siguen alejando unos de otros.
• Esto sí nos proporciona una manera de distinguir el ruido de un atractor con una
dimensionalidad baja. Introduciendo el tamaño del atractor RA, decidimos que un
punto pertenece al ruido (y es un vecino falso) si
R d 1  2 R A
(el “2” es un poco arbitrario). Como estimación del tamaño del atractor RA
podemos tomar simplemente el valor RMS de la señal.
15
Asociación
EURATOM-CIEMAT
Para Fusión
Invariantes
• Los invariantes son propiedades estadísticas del movimiento que no dependen de
condiciones iniciales, la órbita exacta, o perturbaciones pequeñas (ruido).
• Ejemplo: la dimensión fractal.
• Considera la regla dinámica
y (k  1)  F ( y (k ))
que nos indica cómo evoluciona un punto hacia el siguiente en la órbita.
• La “densidad de puntos” r en un punto x del espacio es, por definición:
r( x ) 
1
N
 ( x  y ( k ))

N
k 1
• El uso de la función delta () es conveniente para la matemática, pero en la
práctica se reemplazará por una probabilidad de encontrar un punto en una
caja/división del espacio.
• Cualquier función g(x) integrado sobre la densidad de puntos es invariante.
Demostración:
16
Asociación
EURATOM-CIEMAT
Para Fusión
Invariantes
 g( x ) r ( x )dx

g( F ( x )) r ( x )dx 
1
g 

1
N
N
 g( x )  g( F ( x ))
N
 g( y ( k ))
k 1
N

k 1
g( y (k  1)) 
1
N
 g( y ( N
 1))  g( y (1))   g
para N grande
• Por tanto, los promedios de g(x) son invariantes (para órbitas suficientemente
largas).
• Expresando r en términos de la función de probabilidad:
N
g 

i 1
g( x i ) p( x i )
p( x i ) 
n( x i )
 n( x i )
i
donde xi es una división (elemento) del espacio, y n(xi) el número de puntos en ella.
17
Asociación
EURATOM-CIEMAT
Para Fusión
Dimensiones generalizadas
• De nuevo vamos a calcular la dimensión, ahora desde la perspectiva de invariantes.
• Tomemos como función g un momento de la densidad, es decir:
g( x )  r ( x )
p
• Por tanto, una medida del volumen ocupado es: g
• La definición (diferencial) de la dimensión geométrica es:
1/ p
D  lim
r 0
log( Volúmen )
log( Tamaño )
 lim
r 0
log g
1/ p
log r
 lim
r0
1
q1
log
r
log r
q
 Dq
con q = p+1.
• Así llegamos a la misma dimensión generalizada que antes.
• Ya sabemos: D0 es la dimensión del algorimo de conteo de cajas (mide el volúmen
geométrico) y D1 es la dimensión que toma en cuenta la densidad de puntos, o la
probabilidad de encontrar puntos, en cada caja (mide la “masa” de las cajas, y se
llama dimensión de información).
• En cambio, D2 mide la probabilidad de encontrar parejas de puntos en una caja, y
por eso se llama dimensión de correlación.
18
Asociación
EURATOM-CIEMAT
Para Fusión
Estimación numérica de la
dimensión
• La dimensión D2 da lugar a un algoritmo cómodo para estimar la dimensión de
datos experimentales (el algoritmo de Grassberger y Procaccia). El algoritmo cuenta
cuántos puntos están separados por menos de r:
C 2 (r ) 
 dx r ( x )n(r , x ) 
1
N ( N  1)
N
  r 
y ( j )  y (i )

i j
donde  es la función de Heaviside ((x) = 0 para x < 0, y (x) = 1 para x ≥ 0).
• Para hacer un estimado de la dimensión de un objeto con una dimensión d, se
necesitan muchos puntos (del orden de ad, donde a es un número entre 2 y 3). Por
tanto, en la práctica es imposible medir dimensiones por encima de 20 o 30.
• La dimensión se estima mediante un dibujo de log C(r) frente a log r.
• La definición matemática dice que el valor correcto corresponde al límite r  0,
pero cerca de r = 0 hay poca estadística y por tanto mucho ruido.
• Por tanto, se extrae la dimensión de una parte recta en la gráfica logarítmica doble.
• Usualmente, D0, D1 y D2 son todos muy similares en la práctica.
19
Información dinámica
Asociación
EURATOM-CIEMAT
Para Fusión
• La dimensión del sistema es muy importante, pero no proporciona información
sobre la dinámica del sistema.
• Considera un volúmen elemental (un hipercubo ed).
• Deja evolucionar el sistema. El hipercubo se transforma en otro:
• En unas direcciones crece, en otras se contrae.
• Las razones entre los tamaños de la caja en cada dirección finales y originales se
llaman Números de Lyapunov: L1, L2,..., Ld.
• Ya que el volúmen de la caja debe disminuir para movimientos caóticos, sabemos:
d
L
j
1
j 1
20
Asociación
EURATOM-CIEMAT
Para Fusión
Los exponentes de Lyapunov
• Se ordenan los números Lj tal que L1 > L2 >...> Ld.
• La definición de los exponentes de Lyapunov es: lj = log Lj.
• Los exponentes de Lyapunov también dan lugar a una estimación de la dimensión:
K
DL  K  l
1
K 1
l
i 1
K 1
K
i
donde K es tal que:
l
i 1
i
 0 y
l
i
 0
i 1
• El primer exponente de Lyapunov mide la velocidad de separación de 2 puntos
(proporcional a el1t), y por tanto nos da una información muy importante del
“grado de caoticidad”.
21
Cálculo de los exponentes de
Lyapunov
Asociación
EURATOM-CIEMAT
Para Fusión
• El cálculo de los exponentes de Lyapunov para datos experimentales procede
en varios pasos.
• Primero se elige un punto de una órbita, y unos puntos cercanos.
• Se siguen todos estos puntos en su desarrollo a lo largo de sus respectivas
órbitas durante un tiempo suficientemente largo.
• Se construye un mapa (una transformación matricial) que, en una
aproximación polinómica de bajo orden, lleva los puntos iniciales a los puntos
finales.
• De este mapa se determina cuál es el Jacobiano. Los valores propios del
Jacobiano nos dan los exponentes de Lyapunov.
• Se repite esto para muchos puntos iniciales. Así se obtienen muchos estimados
de los exponentes de Lyapunov, que se promedian para obtener el valor final.
• El algoritmo es algo complicado y requiere mucha atención en el caso de
presencia de ruido. Sin embargo, en algunos casos se han podido extraer unos
exponentes con mucha precisión.
22
Predecir
Asociación
EURATOM-CIEMAT
Para Fusión
• En el cálculo de los exponentes de Lyapunov, que describen el desarrollo del
sistema dinámico, se contruyeron mapas locales del movimiento.
• Estos mapas se pueden aplicar a la predicción del movimiento, ya que relacionan
nuevos puntos y con anteriores: y = F(yprevio).
• Se puedo optimizar la predicción minimizando el error de predicción sobre el
conjunto de datos del que se dispone.
• Nótese que el mapa puede ser global (un sólo mapa para todo el espacio) o local
(un mapa que depende de la posición).
• Se puede diseñar un mapa específicamente para el intervalo T de predicción
requerido, o para un paso temporal D para luego calcular una predicción para el
intervalo T ejecutando la predicción para D T/D veces.
23
Asociación
EURATOM-CIEMAT
Para Fusión
Predecir
• En lugar de polynomios sencillos, uno puede utilizar polinomios ortonormales
especiales para construir el mapa. Si la ortogonalidad de los polinomios f(x) se
determina respecto a la densidad de puntos sobre el atractor, es decir:
 dx r ( x )f
a
( x )f b ( x )   ab
entonces la determinación de los parámetros de ajuste se simplifica mucho,
eliminando la necesidad de hacer ajustes (el procedimiento es muy similar a una
transformada).
• Otra forma de obtener un ajuste es usando redes neuronales. Teóricamente la
calidad de la predicción con esta técnica es muy buena, pero el proceso de ajuste
de parametros es muy no-lineal y costoso.
• Ningún tipo de predicción escapa al problema fundamental de predictibilidad
limitada que supone la dinámica caótica.
24
Ruido
Asociación
EURATOM-CIEMAT
Para Fusión
• La cuestión de los datos “caóticos” contaminados con ruido es un problema
complicado.
• No se pueden aplicar las técnicas de filtración habituales del análisis lineal, ya que
tanto la señal (caótica) como el ruido tienen un espectro ancho.
• Un método bastante eficaz consiste en, otra vez, construir mapas polinómicos de
las órbitas en un entorno con radio r de un punto dado:
y (n  1)  g ( y ( n))  A  By (n )  C y ( n) y (n)  ...
donde A, B y C son matrices, y comparar la predicción polinómica con el valor
siguiente medido.
• Se buscan pequeñas correcciones de la órbita, y, que minimicen:
e  y (n  1)   y (n  1)  g ( y (n )   y (n ))
2
• De esta forma es posible “limpiar” una señal contaminada hasta cierto punto.
• Un punto a considerar es que el filtraje habitual de señales experimentales puede
cambiar por completo las características caóticas de la señal medida.
25
Descargar

Análisis de datos caóticos, dimensión fractal