
Citation 
 Permanent Link:
 https://ufdc.ufl.edu/UFE0017765/00001
Material Information
 Title:
 New Approaches to Multiscale Signal Processing
 Creator:
 HU, JING ( Author, Primary )
 Copyright Date:
 2008
Subjects
 Subjects / Keywords:
 Brownian motion ( jstor )
Chaos ( jstor ) Chaos theory ( jstor ) Eggshells ( jstor ) Fractals ( jstor ) Radar ( jstor ) Seas ( jstor ) Signals ( jstor ) Time series ( jstor ) Trajectories ( jstor )
Record Information
 Source Institution:
 University of Florida
 Holding Location:
 University of Florida
 Rights Management:
 Copyright Jing Hu. Permission granted to University of Florida to digitize and display this item for nonprofit research and educational purposes. Any reuse of this item in excess of fair use or other copyright exemptions requires permission of the copyright holder.
 Embargo Date:
 7/12/2007
 Resource Identifier:
 659872455 ( OCLC )

Downloads 
This item has the following downloads:

Full Text 
NEW APPROACHES TO MULTISCALE SIGNAL PROCESSING
By
JING HU
A DISSERTATION PRESENTED TO THE GRADUATE SCHOOL
OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT
OF THE REQUIREMENTS FOR THE DEGREE OF
DOCTOR OF PHILOSOPHY
UNIVERSITY OF FLORIDA
2007
S2007 Jing Hu
To my family, teachers and friends
ACKENOWLED GMENTS
First of all, I would like to thank my advisor, Professor .Jianho Gao, for his invaluable
guidance and support. I was lucky to work with someone with so many original ideas
and such a sharp mind. Without his help, my Ph.D. degree and this dissertation would
have been impossible. During my years at the University of Florida, I was moved by his
magnificent personality. A famous Chinese adage ; .va, "One d .< as my advisor, entire life
as my father." Professor Gao highly qualifies.
I also want to thank my coninittee nienters (Professors SuShing C'I. in~ Yuguang
Fang, .Jose C. Principe, and K~eith D. White), for their interest in my work and valuable
aHulMy thanks of course also go to all the faculty, staff and my fellow students in
the Department of Electrical and Computer Engineering. I extend particular thanks to my
office mates .Jaentin Lee, Yi Z1n in_ and IUngsik K~in, for their helpful discussions.
I express my appreciation to all of my friends who offered encouragement. I specially
thank Qian Zhan, Xingyan Fan, Lingyu Gu and his wife (.Jing Zhang), Yannmin Xu, Yan
Liu, and Yana Liang. They gave me a lot of help during my early years in the University
of Florida. They taught me to cook and made my life more colorful.
Finally, I owe a great debt of thanks to my parents, grandparents, and my boyfriend,
.Jing Ai, for their pillarlike support and their unwavering relief and confidence in me
without this I do not think I would have made it this far!
TABLE OF CONTENTS
page
ACK(NOWLEDGMENTS ......... . . 4
LIST OF FIGURES ......... .. . 7
ABSTRACT ......... ..... . 9
CHAPTER
1 INTRODUCTION ......... . 10
1.1 Overview of Dynamical Systems . ..... .. 11
1.2 Examples of Multiscale Phenomena . .... .. 15
1.3 Brief Introduction to C'!s I... ......... .. 20
1.4 Brief Introduction to Fractal Geometry ... .. . . 21
1.5 Importance of Connecting ('!, I. .s and Random Fractal Theories .. .. 2:3
1.6 Importance of the Concept of Scale . ..... .. 24
1.7 Structure of this Dissertation ........ .. .. 24
2 GENERALIZATION OF CHAOS: POWERLAW SENSITIVITY TO INITIAL
CONDITIONS (PSIC) ......... ... 27
2.1 Dynamical Test for C'!s I... .... .. . .. .. 27
2.2 General Computational Framework for Powerlaw Sensitivity to Initial
Conditions (PSIC) .. . .. .. .. .. 3:3
2.3 C'll I l .terizing Edge of C'I!s I. by Powerlaw Sensitivity of Initial Conditions
(PSIC) .... ........ ............ :36
:3 CHARACTERIZING RANDOM FRACTALS BY POWERLAW SENSITIVITY
TO INITIAL CONDITIONS (PSIC) . ..... .. 41
:3.1 C'll I l .terizing 1/ fP Processes by Powerlaw Sensitivity to Initial Conditions
(PSIC) .... ........ .. .. ...... ..... 41
:3.1.1 C'll I l .terizing Fractional Brownian Motion (fBm) Processes by
Powerlaw Sensitivity to Initial Conditions (PSIC) .. .. .. .. 4:3
:3.1.2 C'll I l .terizing ON/OFF Processes with Pareto distributed ON and
OFF Periods by Powerlaw Sensitivity to Initial Conditions (PSIC) 50
:3.2 C'll I l .terizing Levy Processes by Powerlaw Sensitivity to Initial Conditions
(PSIC) .... ........ ............ 54
4 STITDY OF INTERNET DYNAMICS BY POWERLAW SENSITIVITY TO
INITIAL CONDITIONS (PSIC) ........ ... .. 60
4.1 Study of Transport Dynamics by Network Simulation .. .. .. .. 61
4.2 Complicated Dynamics of Internet Transport Protocols .. .. .. .. 69
5 STUDY OF SEA CLUTTER RADAR RETURNS BY POWERLAH
SENSITIVITY TO INITIAL CONDITIONS (PSIC) .. .. .. .. 76
5.1 Sea Clutter Data ........ ... .. 78
5.2 Non C'!s I..tic Behavior of Sea Clutter .... .. .. 79
5.3 Target Detection within Sea Clutter by Separating Scales .. .. .. .. 81
5.4 Target Detection within Sea Clutter by Powerlaw Sensitivity to Initial
Conditions (PSIC) ........ . .. 86
6 MITLTISCALE ANALYSIS BY SCALEDEPENDENT LYAPITNOV EXPONENT
(SDLE) ........... ......... .. 89
6.1 Basic Theory
6.2 Classification of Complex Motions.
6.2.1 C'! I ..~ Noisy (I I ..~ and Noiseinduced ChI ... .
6.2.2 Processes Defined by Powerlaw Sensitivity to Initial
(PSIC).
6.2.3 Complex Motions with Multiple Scaling Behaviors.
6.3 C'!I. .) ..terizing Hidden Frequencies
Conditions
7 CONCLUSIONS ......... ... .. 104
REFERENCES .. .......... ........... 105
BIOGRAPHICAL SK(ETCH ...... .. 115
LIST OF FIGURES
Figure page
11 Example of a trajectory in the phase space ..... .. .. 14
12 Giant ocean wave (tsunami). ......... .. 17
13 Example of sea clutter data. ......... .. 18
14 Example of Heart rate variability data for a normal subject. .. .. .. 19
21 Timedependent exponent A(k) vs. evolution time k curves for Lorenz system. .32
22 Timedependent exponent A(k) vs. evolution time k curves for the MackeyGlass
system. ......... ..... 33
23 Timedependent exponent A(t) vs. t curves for time series generated from logistic
map. ......... ............ .... 37
24 Timedependent exponent A(t) vs. t curves for time series generated from Henon
map. ............ ............... 40
31 White noise and Brownian motion . ...... .. 46
32 Several fBm processes with different H. . ... .. 48
33 Timedependent exponent A(k) vs. In k curves for fractional Brownian motion. 51
34 Example of ON/OFF processes. .... ... .. 51
35 Timedependent exponent A(k) vs. In k curves for ON/OFF processes with Pareto
distributed ON and OFF periods. ... ... .. 53
36 Examples for Brownian motions and Levy motions. .. .. .. 57
37 Timedependent exponent A(k) vs. In k curves for Levy processes. .. .. .. 58
41 Example of quasiperiodic congestion window size W(i) time series. .. .. .. 63
42 Examples of chaotic T(i) time series. . ..... .. .. 65
43 Complementary cumulative distribution functions (CCDFs) for the two chaotic
T(i) time series of Figs. 42(a,c). . ... ... .. 66
44 The time series T(i) extracted from a W(i) time series collected using net100
instruments over ORNLLSU connection. ..... .. .. 68
45 Time series for the congestion window size cwnd for ORNLGaTech connection. 72
46 Timedependent exponent A(k) vs. k curves for cwnd data corresponding to
Figf.45. ......... ..... 73
47 Timedependent exponent A(k) vs. In k curves for cwnd data corresponding to
Figf.45. ......... ..... 75
51 Collection of sea clutter data ......... ... 79
52 Typical sea clutter amplitude data. ........ .. .. 79
53 Two short segments of the amplitude sea clutter data severely affected by clipping. 80
54 Examples of the timedependent exponent A(k) vs. k curves for the sea clutter
data ........ ... . .... 81
55 Variations of the Lyapunov exponent A estimated by conventional methods vs.
the 14 range bins for the 10 HH measurements. .... .. .. 84
56 Variations of the scaledependent exponent corresponding to large scale vs. the
range bins for the 10 HH measurements. ..... .. .. 85
57 Examples of the A(k) vs. In k curves for the sea clutter data. .. .. .. .. 87
58 Variation of the P parameter with the 14 range bins. .. .. ... .. 87
59 Frequencies of the bins without targets and the bins with primary targets for
the HH datasets. .. ... .. .. 88
61 Scaledependent Lyapunov exponent A(e) curves for the Lorenz system and logistic
map ........ ... . .... 93
62 Scaledependent Lyapunov exponent A(e) curve for the MackeyGlass system. .. 94
63 The functions F(x) and G(x). ......... ... .. 97
64 Scaledependent Lyapunov exponent A(e) for the model described by Eq. 69. .97
65 Powerspectral density (PSD) for time series generated from Lorenz system. 100
66 Lorenz attractor. .. ... .. .. 101
67 Hidden frequency phenomenon of laser intensity data. ... .. .. .. 102
Abstract of Dissertation Presented to the Graduate School
of the University of Florida in Partial Fulfillment of the
Requirements for the Degree of Doctor of Philosophy
NEW APPROACHES TO MITLTISCALE SIGNAL PROCESSING
By
Jing Hu
May 2007
C'I I!r: Jianho Gao
Major: Electrical and Computer Engineering
Complex systems usually comprise multiple subsystems with highly nonlinear
deterministic and stochastic characteristics, and are regulated hierarchically. These
systems generate signals with complex characteristics such as sensitive dependence on
small disturbances, long nienory, extreme variations, and nonstationarity. Chaos theory
and random fractal theory are two of the most important theories developed for analyzing
complex signals. However, they have entirely different foundations, one being deterministic
and the other being random. To synergistically use these two theories, we developed a new
theoretical framework, powerlaw sensitivity to initial conditions (PSIC), to encompass
both chaos and random fractal theories as special cases. To show the power of this
framework, we applied it to study two challenging and important problems: Internet
dynamics and sea clutter radar returns. We showed that PSIC can readily characterize
the complicated Internet dynamics due to interplay of nonlinear AINID (additive increase
and niultiplicative decrease) operation of TCP and stochastic user behavior, and robustly
detect low observable targets from sea clutter radar returns with high accuracy.
We also developed a new niultiscale complexity measure, scaledependent Lyapunov
exponent (SDLE), that can he computed from short noisy data. The measure readily
classified various types of complex motions, and simultaneously characterized behaviors of
complex signals on a wide range of scales, including complex irregular behaviors on small
scales, and orderly behaviors, such as oscillatory motions, on large scales.
CHAPTER 1
INTRODUCTION
A dynamical system is one that evolves in time. Dynamical systems can be stochastic
(in which case they evolve according to some random process such as the toss of a coin)
or deterministic (in which case the future is uniquely determined by the past according
to some rule or mathematical formula). Whether the system in question settles down to
equilibrium, keeps repeating in cycles, or does something more complicated, dynamics are
what we use to analyze the behavior of systems.
Complex dynamical systems usually comprise multiple subsystems with highly
nonlinear deterministic and stochastic characteristics, and are regulated hierarchically.
These systems generate signals with complex characteristics such as sensitive dependence
on small disturbances, long memory, extreme variations, and nonstationarity [1]. A
stock market, for example, is strongly influenced by multi l1~ we decisions made by
market makers, as well as influenced by interactions of heterogeneous traders, including
intradwl~i traders, shortperiod traders, and longperiod traders, and thus gives rise to
highly irregular stock prices. The Internet, as another example, has been designed in a
fundamentally decentralized fashion and consists of a complex web of servers and routers
that cannot be effectively controlled or analyzed by traditional tools of queuing theory
or control theory, and give rise to highly bursty and multiscale traffic with extremely
high variance, as well as complex dynamics with both deterministic and stochastic
components. Similarly, biological systems, being heterogeneous, massively distributed, and
highly complicated, often generate nonstationary and multiscale signals. With the rapid
accumulation of complex data in life sciences, systems biology, nanosciences, information
systems, and physical sciences, it has become increasingly important to be able to analyze
multiscale and nonstationary data.
Multiscale signals behave differently depending upon which scale the data are looked
at. In order to fully characterize such complex signals, the concept of scale has to be
explicitly considered, and existing theories need to be used synergistically, instead of
individually, since on different scale ranges different theories may be most pertinent.
This dissertation aints to build an effective arsenal for niultiscale signal processing by
synergistically integrating approaches based on chaos theory and random fractal theory,
and going beyond, to complement conventional approaches such as spectral analysis and
machine learning techniques. To make such an integration possible, two important efforts
have been made:
* Powerlaw sensitivity to initial conditions (PSIC): We developed a new
theoretical framework for signal processing, to create a coninon foundation for chaos
theory and random fractal theory, so that they can he better integrated.
* Scaledependent Lyapunov exponent (SDLE): We developed an excellent
niultiscale measure, which is a variant of the finite size Lyapunov exponent (FSLE).
We proposed a highly efficient algorithm for calculating it, and showed that it can readily
classify different types of motions, including truly lowdintensional chaos, noisy chaos,
noiseinduced chaos, random 1/ fP and ccstable Levy processes, and complex motions with
chaotic behavior on small scales but diffusive behavior on large scales. The measure can
aptly characterize complex behaviors of real world niultiscale signals on a wide range of
scales.
In the rest of this chapter, we shall first present an overview of dynamics systems,
especially we shall point out why nonlinear systems are much harder to analyze than
linear ones. We then give a few examples of nmultiscale phenomena to show the difficulties
and excitement of niultiscale signal processing. After that, we shall briefly introduce some
background knowledge about chaos and fractal geometry. Then we discuss the importance
of connecting chaos theory and random fractal theory for characterizing the behaviors
of niultiscale signals on a wide range of scales. We also emphasize the importance of
explicitly incorporating the concept of scale in devising measures for characterizing
niultiscale signals. Finally, we outline the scope of the present study.
1.1 Overview of Dynamical Systems
There are two main types of dynamical systems: differential equations and iterated
maps (difference equations). Differential equations describe the evolution of systems
in continuous time, whereas iterated maps arise in problems where time is discrete. To
show how systems evolve, we take differential equations as examples. When confining our
attention to differential equations, the main distinction is between ordinary and partial
differential equations. For instance, the equation for a damped harmonic oscillator
m +1 b + kx = 0 (11)
is an ordinary differential equation, because it involves only ordinary derivatives dx/dt and
d2 ld2. That is, there is only one independent variable, the time t. In contrast, the heat
equation
Bu 82U
dt 8x2
is a partial differential equation: it has both time t and sapce x as independent variables.
Our concern in this proposal is with purely temporal behavior, and so we deal with
ordinary differential equations almost exclusively.
A very general framework for ordinary differential equations is provided by the system
(12)
Here the overdots denote differentiation with respect to t. Thus xi = dxi/dt. The
variables xl, x, might represent concentrations of chemicals in a reactor, populations
of different species in an ecosystem, or the positions and velocities of the planets in the
solar system. The functions fl, f, are determined by the problem at hand.
For example, the damped oscillator (Eq. 11) can be rewritten in the form of
(Eq. 12), thanks to the following trick: we introduce new variables xl = x and x2 x
Then xl = Z2. Hence the equivalent system (Eq. 12) is
xi = x2
b k
xa2=x 2 1i
m m
This system is said to be linear, because all the xi on the righthand side appear to
the first power only. Otherwise the system would be nonlinear. Typical nonlinear terms
are products, powers, and functions of the Xi, such aS xlx2, 1 3), Or COSZ2.
For example, the swinging of a pendulum is governed by the equation
x + sinx = 0,
where x is the angle of the pendulum from vertical, g is the acceleration due to gravity,
and L is the length of the pendulum. The equivalent system is nonlinear:
x z = x2
Z2 = L~smnxl.
Nonlinearity makes the pendulum equation very difficult to solve analytically. The
usual way around this is to fudge, by invoking the small angle approximation sinx e x for
x << 1. This converts the problem to a linear one, which can then be solved easily. But
by restricting to small x, we are throwing out some of the physics, like motions where the
pendulum whirls over the top. Is it really necessary to make such drastic approximation?
It turns out that the pendulum equation can be solved analytically, in terms of elliptic
functions. But there ought to be an easier way. After all, the motion of the pendulum is
simple: at low energy, it swings back and forth, and at high energy it whirls over the top.
There should be some way of extracting this information from the system directly using
geometric methods.
Here is the rough idea. Suppose we happen to know a solution to the pendulum
system, for a particular initial condition. This solution would be a pair of functions xl(t)
and x2(t), representing the position and velocity of the pendulum. If we construct an
I (t), X 2(t))
(x (0), X2(0))
Figure 11. Example of a trajectory in the phase space
abstract space with coordinates (xl,,Z2), then the solution (xl(t),2 xa1)) COTTOSponds to
a point moving along a curve in this space, as shown Fig. 11. This curve is called a
trajectory, and the space is called the phase space for the system. The phase space is
completely filled with trajectories, since each point can serve as an initial condition.
Our goal is to run this construction in reverse: given the system, we want to draw
the trajectories, and thereby extract information about the solutions. In many cases,
geometric reasoning will allow us to draw the trajectories without actually solving the
system!
The phase space for the general system (Eq. 12) is the space with coordinates
xl, x,. Because this space is ndimensional, we will refer to Eq. 12 as an ndimensional
system or an nthorder system. Thus a represents the dimension of the phase space.
As we have mentioned earlier, most nonlinear systems are impossible to solve
analytically. Why are nonlinear systems so much harder to analyze than linear ones?
The essential difference is that linear systems can be broken down into parts. Then each
part can be solved separately and finally recombined to get the answer. This idea allows a
fantastic simplification of complex problems, and underlies such methods as normal modes,
Laplace transforms, superposition arguments, and Fourier analysis. In this sense, a linear
system is precisely equal to the sum of its parts.
But many things in nature do not act this way. Whenever parts of a system interfere,
or cooperate, or compete, there are nonlinear interactions going on. Most of every w life
is nonlinear, and the principle of superposition fails spectacularly. If you listen to your two
favorite songs at the same time, you will not get double the pleasure! Within the realm of
physics, nonlinearity is vital to the operation of a laser, the formation of turbulence in a
fluid, and the superconductivity of Josephson junctions.
1.2 Examples of Multiscale Phenomena
Multiscale phenomena are ubiquitous in nature and engfineeringf. Some of them are,
unfortunately, catastrophic. One example is the Tsunami that occurred in south Asia
around C!~!sI it ~ 2004. Another example is the gigantic power outage occurred in North
America on August 14, 2003, which affected more than 4000 megawatts, was more than
300 times greater than mathematical models would have predicted, and cost between $4
billion and 71 billion, according to the U.S. Department of energy. Both events involved
scales that, on one hand, were so huge that our human being could not easily fathom,
and on the other hand, involved the very scales that were most relevant to individual life.
Below, we consider three more examples.
Multiscale phenomena in computer networks: Largescale communications
networks, especially the Internet, are among the most complicated systems that man
has ever made, with multiscales in terms of a number of aspects. Intuitively speaking,
these multiscales come from the hierarchical design of protocol stack, the hierarchical
topological architecture, and the multipurpose and heterogeneous nature of the Internet.
More precisely, there are multiscales in
*Time, manifested by the prevailing fract al, multifract al, and longrange dependent
properties in traffic,
* Space, essentially due to topology and geography and again manifested by scalefree
properties,
* State (e.g., queues, windows),
* Size (e.g., number of nodes, number of users).
Also, it has been observed that the failure of a single router may trigger routing
instability, which may be severe enough to instigate a route flap storm. Furthermore,
packets may be delivered out of order or even get dropped, and packet reordering is not
a pathological network behavior. As the next generation Internet applications such as
remote instrument control and computational steering are being developed, another facet
of complex niultiscale behavior is beginning to surface in terms of transport dynanxics.
The networking requirements for these next generation applications belong to (at least)
two broad classes involving vastly disparate time scales:
* High handwidths, typically multiples of 10Ghps, to support bulk data transfers,
* Stable handwidths, typically at much lower handwidths such as 10 to 100 1\bps, to
support interactive, steering and control operations.
Multiscale phenomena in sea clutter: Sea clutter is the radar backscatter
front a patch of ocean surface. The complexity of the signals contes front two sources,
the rough sea surface, sometimes oscillatory, sometimes turbulent, and the nmultipath
propagation of the radar backscatter. This can he appreciated by imagining how radar
pulses massively reflecting front the wavetip of Fig. 12. To be quantitative, in Fig. 13,
two 0.1 s duration sea clutter signals, sampled with a frequency of 1 K(Hz, are plotted
in Fig. 1:3(a,h), a 2 s duration signal is plotted in Fig. 1:3(c), and an even longer signal
(about 1:30 s) is plotted in Fig. 1:3(d). It is clear that the signal is not purely random,
since the waveform can he fairly smooth on short time scales (Fig. 1:3(a)). However, the
signal is highly nonstationary, since the frequency of the signal (Fig. 1:3(a,b)) as well as
the randoniness of the signal (Fig. 1:3(c,d)) change over time drastically. One thus can
perceive that naive Fourier analysis or deterministic chaotic analysis of sea clutter may
Figure 12. Giant ocean wave (tsunami). Suppose our field of observation includes the
wavetip of length scale of a few meters, it is then clear that the complexity of
sea clutter is mainly due to massive reflection of radar pulses from a wavy and
even turbulent ocean surface.
not e very,,,,,, useul From Fig 13(e,\ where Xt(m) is the nonoverlapping running mean of
X over block size m, and X is the sea clutter amplitude data, it can be further concluded
that neither autoregressive (AR) models (As an example, we give the definition for the
firstorder AR model, which is given by x, I = ax, + rl,, where the constant coefficient a
satisfies 0 / a < 1 and rlo is a white noise with zero mean.) nor textbook fractal models
can describe the data. This is because AR modeling requires exponentially decaying
autcorrelation (which,1 amounts, to Var(Xe(m)) ~ ml, or a Hurst parameter of 1/2. See
later chapters for the definition of the Hurst parameter), while fractal modeling requires
the variation between Vanr(Xe"m)) and m to follow a powerlaw~. However, neither behavior
is observed in Fig. 13(e). Indeed, albeit extensive work has been done on sea clutter,
the nature of sea clutter is still poorly understood. As a result, the important problem of
80 80.02 80.04 80.06
Time (sec)
80.08 8C
Time (sec)
8
log2 m
U 2U 4U bU0 tU 1UU 12U
Time (sec)
Figure 13.
Example of sea clutter data. (a,b) two 0.1 s duration signal; (c) a 2 s duration
sea clutter signal; (d) the entire sea clutter signal (of about 130 s); and (e)
log2 m2Var(X(m))] vs. log2 m, Where
Xim) {Xe(m) : X(m) (Xtmm+l + + Xtm)/m, t = 1, 2, } is the
nonoverlapping running mean of X = {Xt : t = 1, 2, } over block size m,
and X is the sea clutter amplitude data. To better see the variation of
Vawr(X~) with mVa(e) is multiplied by m2. When the autocorrelation
of the data decays exponentially fast (such as modeled by an autoregressive
(AR) process), V'ar(Xe(m")) ~ m1 Here Var(Xe(m)) decatys much faster. A
fractal process would have m2"Var(Xe"m)) ~ m ". However, th~is is not thle case.
Therefore, neither AR modeling nor ideal textbook fractal theory can be
readily applied here.
target detection within sea clutter remains a tremendous challenge. We shall return to sea
clutter later.
Multiscale and nonstationary phenomena in heart rate variability (HRV):
HRV is an important dynamical variable of the cardiovascular function. Its most salient
feature is the spontaneous fluctuation, even when the environmental parameters are
(b)
(d)
a,6
4
E
0 1 2 3 4 5 6 7
Index n4
x 104
100 110
(b) (c)
100
95
90
I80
2.66 2.67 2.68 2.69 2.7 2.71 3.4 .9 396 .7 398 .9
Index n x 104 Index n x 104
(d) (e)
104 105
S102 0
102 10 102 10
f f
Figure 14. Example of Heart rate variability data for a normal subject. (a) the entire
signal, (hsc) the segments of signals indicated as A and B in (a); (d,e) power
spectral density for the signals shown in (b,c).
maintained constant and no perturbing influences can he identified. It has been observed
that HRV is related to various cardiovascular disorders. Therefore, analysis of HRV is
very important in medicine. However, this task is very difficult, since HRV data are highly
complicated. An example is shown in Fig. 14, for a normal young subject. Evidently,
the signal is highly nonstationary and niultiscaled, appearing oscillatory for some period
of time (Figs. 14(h,d)), and then varying as a powerlaw for another period of time
(Figs. 14(c,e)). The latter is an example of the socalled 1/ f processes, which will be
discussed in depth in later chapters. While the niultiscale nature of such signals cannot
he fully characterized by existing methods, the nonstationarity of the data is even more
troublesome, since it requires the data to be properly segmented before further analysis
hy methods from spectral analysis, chaos theory, or random fractal theory. However,
automated segmentation of complex biological signals to remove undesired components is
itself a significant open problem, since it is closely related to, for example, the challenging
task of accurately detecting transitions from normal to abnormal states in physiological
data.
1.3 Brief Introduction to Chaos
Imagine we are observing an aperiodic, highly irregular time series. Can such a
signal arise from a deterministic system that can he characterized by only very few state
variables instead of a random system with infinite numbers of degrees of freedom? A
chaotic system is capable of just that! This discovery has such far reaching implications
in science and engineering that sometimes chaos theory is considered as one of the three
most revolutionary scientific theories of the twentiethcentury, along with relativity and
quantum mechanics.
At the center of chaos theory is the concept of sensitive dependence on initial
conditions: a very minor disturbance in initial conditions leads to entirely different
outcomes. An often used metaphor illustrating this point is that a sunny weather in New
York could be replaced by a rainy one some time in the future after a butterfly flaps
its wings in Boston. Such a feature contrasts sharply with the traditional view, largely
based on our experience with linear systems, that small disturbances (or causes) can only
generate proportional effects, and that in order for the degree of randomness to increase,
the number of degrees of freedom has to be infinite.
No definition of the term chaos is universally accepted yet, but almost everyone would
agree on the three ingredients used in the following working definition:
C'!s I. .s is periodic longterm behavior in a deterministic system that exhibits sensitive
dependence on initial conditions.
*"Aperiodic longterm b. !! lit.! I means that there are trajectories which do not settle
down to fixed points, periodic orbits, or quasiperiodic orbits as t oc. For practical
reasons, we should require that such trajectories are not too rare. For instance, we could
insist that there he an open set of initial conditions leading to periodic trajectories, or
perhaps that such trajectories should occur with nonzero probability, given a random
initial condition.
* "Deterministic" means that the system has no random or noisy inputs or parameters.
The irregular behavior arises from the system's nonlinearity, rather than from noisy
driving forces.
* "Sensitive dependence on initial conditions" means that nearby trajectories separate
exponentially fast (i.e., the system has a positive Lyapunov exponent).
C'!s I. .s is also commonly called a strange attractor. The term attractor is also difficult
to define in a rigorous way. We want a definition that is broad enough to include all
the natural candidates, but restrictive enough to exclude the imposters. There is still
disagreement about what the exact definition should he.
Loosely lo' I1:;0s an attractor is a set to which all neighboring trajectories converge.
Stable fixed points and stable limit cycles are examples. 1\ore precisely, we define an
attractor to be a closed set ,4 with the following properties:
1. A is an invariant set: any trajectory x(t) that starts in ,4 stays in ,4 for all time.
2. A attracts an open set of initial conditions: there is an open set U containing 24 such
that if x(0) E U, then the distance from x(t) to A tends to zero as t 00. This
means that ,4 attracts all trajectories that start sufficiently close to it. The largest
such U is called the basin of attraction of ,4.
3. A is minimal: there is no proper subset of ,4 that satisfies conditions 1 and 2.
Now we can define a strange attractor to be an attractor that exhibits sensitive dependence
on initial conditions. Strange attractors were originally called strange because they are
often fractal sets. ? ue ul .vis this geometric property is regarded as less important than
the dynamical property of sensitivity dependence on initial conditions. The terms chaotic
attractor and fractal attractor are used when one wishes to emphasize one or the other of
those aspects.
1.4 Brief Introduction to Fractal Geometry
Euclidean geometry is about lines, planes, triangles, squares, cones, spheres, etc. The
common feature of these different objects is regularity: none of them is irregular. Now let
us ask a question: Are clouds spheres, mountains cones, and islands circles? The answer is
obviously no. Pursuing answers to such questions, Alandelbrot has created a new branch of
science fractal geometry.
For now, we shall be satisfied with an intuitive definition of a fractal: a set that shows
irregular but selfsimilar features on many or all scales. Selfsimilarity means a part of
an object is similar to other parts or to the whole. That is, if we view an irregular object
with a microscope, whether we enlarge the object by 10 times, or by 100 times, or even
by 1000 times, we ah in find similar objects. To understand this better, let us imagine
we were observing a patch of white cloud drifting .l. li in the sky. Our eyes were rather
passive: we were staring more or less at the same direction. After a while, the part of the
cloud we saw drifted away, and we were viewing a different part of the cloud. Nevertheless,
our feeling remains more or less the same.
Mathematically, fractal is characterized by a powerlaw relation, which translates to
a linear relation in loglog scale. For now, let us again resort to imagination we were
walking down a wild and very j I__ d mountain trail or coastline. We would like to know
the distance covered by our route. Suppose our ruler has a length of e which could be
our stepsize, and different hikers may have different stepsizes a person riding a horse
has a huge stepsize, while a group of people with a little child must have a tiny stepsize.
The length of our route is
L = N(e) (13)
where NV(e) is the number of intervals needed to cover our route. It is most remarkable
that typically NV(e) scales with e in a powerlaw manner,
N~e) ED e 0(14)
with D being a noninteger, 1 < D < 2. Such a nonintegral D is often called the fractal
dimension to emphasize the fragmented and irregular characteristics of the object under
study.
Let us now try to understand the meaning of the nonintegral D. For this purpose, let
us consider how length, area, and volume are measured. A common method of measuring
a length, a surface area, or a volume is by covering it with intervals, squares, or cubes
whose length, area, or volume is taken as the unit of measurement. These unit intervals,
squares, and volumes are called unit boxes. Suppose, for instance, that we have a line
whose length is 1. We want to cover it by intervals (hoxes) whose length is e. It is clear
that we need NV(e) ~ e1 hoxes to completely cover the line. Similarly, if we want to cover
an area or volume by boxes with linear length e, we would need NV(e) ~ E2 to cover the
area, or NV(e) ~ e3 hoxes for the volume. Such D is called the topological dimension, and
takes on a value of one for a line, two for an area, and three for a volume. For isolated
points, D is zero. That is why a point, a line, an area, and a volume are called an 0 D,
1 D, 2 D, and 3 D objects, respectively.
Now let us examine the consequence of 1 < D < 2 for a j I__ d mountain trail. It is
clear that the length of our route increases as e becomes smaller, i.e., when e 0 L o.
T> he more concrete, let us visualize a race between the Hare and the Tortoise on a fractal
trail with D = 1.25. Assume that the length of the average step taken by the Hare is 16
times that taken by the Tortoise. Then we have
LHare 2 ,Lartoise
That is, the Tortoise has to run twice the distance of the Hare! Put differently, if you
were walking along a wild mountain trail or coastline and tired, slowing down your pace,
shrinking your steps, then you were in trouble! It certainly would be worse if you also got
lost!
1.5 Importance of Connecting Chaos and Random Fractal Theories
Multiscale signals behave differently depending upon which scale the data are
looked at. How can the behaviors of such complex signals on a wide range of scales
he simultaneously characterized? One strategy we envision is to use existing theories
synergistically, instead of individually. To make this possible, appropriate scale ranges
where each theory is most pertinent need to be identified. This is a difficult task, however,
since different theories may have entirely different foundations. For example, chaos theory
is mainly concerned about apparently irregular behaviors in a complex system that are
generated by nonlinear deterministic interactions of only a few numbers of degrees of
freedom, where noise or intrinsic randomness does not pIIl w any important role. Random
fractal theory, on the other hand, assumes that the dynamics of the system are inherently
random. Therefore, to fully characterize the behaviors of multiscale signals on a wide
range of scales, different theories, such as chaos and random fractal theories, need to be
integrated and even generalized.
1.6 Importance of the Concept of Scale
Complex systems often generate highly nonstationary and multiscale signals because
of nonlinear and stochastic interactions among their component systems and hierarchical
regulations imposed by the operating environments. Rapid accumulation of such complex
data in life sciences, systems biology, nanosciences, information systems, and physical
sciences, has made it increasingly important to develop complexity measures that
incorporate the concept of scale explicitly, so that different behaviors of signals on varying
scales can be simultaneously characterized by the same scaledependent measure. The
most ideal scenario is that a scaledependent measure can readily classify different types
of motions based on analysis of short noisy data. In this case, one can readily see that the
measure will not only be able to identify appropriate scale ranges where different theories,
including information theory, chaos theory, and random fractal theory, apply, but also be
able to automatically characterize the behaviors of the data on those scale ranges. The
importance of the concept of scale will be further illustrated in later chapters.
1.7 Structure of this Dissertation
We emphasized that to develop new approaches for multiscale signal processing, it is
most desirable to integrate and even generalize different theories. Recently, the defining
property for deterministic chaos, exponential sensitivity to initial conditions (ESIC) has
been generalized to powerlaw sensitivity to initial conditions (PSIC). In C'!s Ilter 2,
we first introduce the new concept of PSIC. Then we extend the concept of PSIC to
highdimensional case. We shall present a general computational framework for assessing
PSIC from real world data. We also apply the proposed computational procedure to
study noisefree and noisy logistic and Henon maps at the edge of chaos. we show that
when noise is absent, PSIC is hard to detect from a scalar time series. However, when
there is dynamic noise, motions around the edge of chaos, he it simply regular or truly
chaotic when there is no noise, all collapse onto the PSIC attractor. Hence, dynamic
noise makes PSIC observable. In (I Ilpter 3, we demonstrate that the new framework of
PSIC can readily characterize random fractal processes, including 1/ fP processes with
longrangecorrelations and Levy processes. From our study in C'!s Ilters 2 and 3, it is
clear that the concept of PSIC not just bridges standard chaos theory and random fractal
theory, but in fact provides a more general framework to encompass both theories as
special cases. To illustrate the power of the framework of PSIC, in OsI Ilpters 4 and 5, we
apply it to study two important but challenging problems: Internet dynamics and sea
clutter radar returns, respectively. Both are outstanding examples of complex dynamical
systems with both nonlinearity and randomness. The nonlinearity of Internet dynamics
comes from AllMD (additive increase and multiplicative decrease) operation of TCP
(Transmission Control protocol), while the stochastic component comes from the random
user behavior. Sea clutter is a nonlinear dynamical process with stochastic factors due
to interference of various wind and swell waves and to local atmospheric turbulence.
We shall show that the new theoretical framework of PSIC can effectively characterize
the complicated Internet dynamics as well as to robustly detect low observable targets
from sea clutter radar returns with high accuracy. In Chapter 6, we shall introduce
the basic theory of an excellent multiscale measure the scaledependent Lyapunov
exponent (SDLE), and develop an effective algorithm to compute the measure. To
understand the SDLE as well as appreciate its power, we apply it to classify various types
of complex motions, including truly lowdintensional chaos, noisy chaos, noiseinduced
chaos, processes defined by PSIC, and complex motions with chaotic behavior on small
scales but diffusive behavior on large scales. We also discuss how the SDLE can help
detect hidden frequencies in large scale orderly motions. Finally we suninarize our work in
CHAPTER 2
GENERALIZATION OF CHAOS: POWERLAW SENSITIVITY TO INITIAL
CONDITIONS (PSIC)
The fornialism of nonextensive statistical niechanics (NESM) [2] has found numerous
applications to the study of systems with longrangeinteractions [36] and niultifractal
behavior [7, 8], and fully developed turbulence [810], among many others. In order to
characterize a type of motion whose complexity is neither regular nor fully chaotic/randons
recently the concept of exponential sensitivity to initial conditions (ESIC) of deterministic
chaos has been generalized to powerlaw sensitivity to initial conditions (PSIC) [7, 11].
Mathematically, the formulation of PSIC closely parallels that of NESM. PSIC has been
applied to the study of deterministic 1D logisticlike maps and 2D Henon niap at the
edge of chaos [7, 1120], yielding considerable new insights. In this chapter, we first
briefly describe some background knowledge about chaos. In particular, we will discuss a
direct dynamical test for deterministic chaos developed by Gao and Zheng [21, 22]. This
method offers a more stringent criterion for detecting lowdintensional chaos, and can
simultaneously monitor motions in phase space at different scales. Then we introduce the
more generalized concept of PSIC, and extend the new concept to highdintensional case.
Specifically, we present a general computational framework for assessing PSIC in a time
series. We also apply the proposed computational procedure to study two model systems:
noisefree and noisy logistic and Henon maps at the edge of chaos.
2.1 Dynamical Test for Chaos
C'!s I. .s is also coninonly called a strange attractor. Being an I11 .. 1 ..1~", the
trajectories in the phase space are bounded. Being lI I.ng.", the nearby trajectories,
on the average, diverge exponentially fast. The latter property is characterized hv the
exponential sensitivity to initial conditions (ESIC), and can he niathentatically expressed
as follows. Let d(0) be the small separation between two arbitrary trajectories at time
0, and d(t) be the separation between them at time t. Then, for true lowdintensional
deterministic chaos, we have
d(t) ~ d(0)ex' (21)
where X1 is positive and called the largest Lyapunov exponent. Due to the boundedness
of the attractor and the exponential divergence between nearby trajectories, a strange
attractor typically is a fractal, characterized by a simple and elegant scaling law as defined
by Eq. 14. A nonintegral fractal dimension of an attractor contrasts sharply with the
integervalued topological dimension (which is 0 for finite number of isolated points, 1 for
a smooth curve, 2 for a smooth surface, and so on).
Conventionally, it has been assumed that a time series with an estimated positive
largest Lyapunov exponent and a nonintegral fractal dimension is chaotic. However, it has
been found that this assumption may not he a sufficient indication of deterministic chaos.
One counterexample is the socalled 1/ fP processes. Such processes have spectral density
S(f ) f (2H+1) (22)
where 0 < H < 1 is sometimes called the Hurst parameter. This type of processes
will be further discussed in much depth in ChI Ilpter 3. A trajectory formed by such a
process has a fractal dimension of 1/H. As we shall explain later, with most algorithms
of estimating the largest Lyapunov exponent, one obtains a positive number for the
"Lyapunov expei ill ~! Therefore, such processes could be interpreted as deterministic
chaos .
Due to the ubiquity of noise, experimental data are ahrl : corrupted by noise to some
degree. mathematically speaking, a noisy system, no matter how weak the noise is, has
infinite dimensions. Experimentally speaking, one would be more interested in a certain
range of finite scales. If the noise is very weak, then its influence on the dynamics may be
limited to very small scales, leaving the dynamics on finite scales deterministiclike.
The above discussions motivate us to employ a more sophisticated method, the
direct dynamical test for deterministic chaos [21, 22], to determine whether the data
under investigation is truly chaotic or not. The method offers a more stringent criterion
for lowdimensional chaos, and can simultaneously monitor motions in phase space at
different scales. The method has found numerous applications in the study of the effects
of noise on dynamical systems [2327], estimation of the strength of measurement noise
in experimental data [28, 29], pathological tremors [30], shearthickening surfactant
solutions [31], dilute sheared aqueous solutions [32], and serrated plastic flows [33].
Now let us briefly describe the direct dynamical test for deterministic chaos [21, 22].
Given a scalar time series, x(1), x(2),..., x(NV), (.III11.11. for convenience, that they have
been normalized to the unit interval [0,1]), one first constructs vectors of the following
form using the time delay embedding technique [3436]:
1K = [x(i), x(i + L), ..., x(i + (m 1)L)] (23)
with m being the embedding dimension and L the delay time. For example, when m = 3
and L = 4, we have VI = [x(1), x(5), x(9)], VW = [x(2), x(6), x(10)], and so on. For the
analysis of purely chaotic signals, m and L has to be chosen properly. This is the issue of
optimal embedding (see [21, 22, 37] and references therein). After the scalar time series is
properly embedded, one then computes the time dependent exponent A(k) curves
A(k) = n +kV (24)
with r < %~ Vy  < r + ar, where r and Ar are prescribed small distances. The angle
brackets denote ensemble averages of all possible pairs of (1M, Vj). The integer k, called
the evolution time, corresponds to time kbt, where 6t is the sampling time. Note that
geometrically (r, r + Ar) defines a shell, and a shell captures the notion of scale. The
computation is typically carried out for a sequence of shells. Comparing Eq. 24 with
Eq. 21, one notices that %~+k Vy,  phI i4 the role of d(t), while %~ Vy1  phi4 the role
of d(0). Figure 21(a) shows the A(k) vs. k curves for the chaotic Lorenz system driven by
stochastic forcing:
dx
=16(x y) + Drll(t)
= xz + 45.92x y + Drl2(t) (25)
= xy 4z + Dr/3(
with < rli(t) > 0, < i6i,(')> b(tt'), i, j = 1, 2, 3. Note that D2 is the variance of
the Gaussian noise terms, and D = 0 describes the clean Lorenz system. We observe from
Fig. 21(a) that the A(k) curves are composed of three parts. These curves are linearly
increasing for 0 < k < k,. They are still linearly increasing for k, < k < ky, but with a
slightly different slope. They are flat for k > k. Note that the slope of the second linearly
increasing part gives an estimate of the largest positive Lyapunov exponent [21, 22]. k,
is related to the time scale for a pair of nearby points (1M, Vyj) to evolve to the unstable
manifold of 1M or Vyj. It is on the order of the embedding window length, (m 1)L. kg
is the prediction time scale. It is longer for the A(k) curves that correspond to smaller
shells. The difference between the slopes of the first and second linearly increasing parts
is caused by the discrepancy between the direction defined by the pair of points (1K, Vyj)
and the unstable manifold of 1M or Vyj. This feature was first observed by Sato et al. [38]
and was used by them to improve the estimation of the Lyapunov exponent. The first
linearly increasing part can be made smaller or can even be eliminated by adjusting
the embedding parameters such as by using a larger value for m. Note that the second
linearly increasing parts of the A(k) curves collapse together to form an envelope. It is
this very feature that forms the direct dynamical test for deterministic chaos [21, 22].
This is because the A(k) curves for noisy data, such as white noise or 1/ f processes,
are only composed of two parts, an increasing (but not linear) part for k < (m 1)L
and a flat part [21, 22]. Furthermore, different A(k) curves for noisy data separate from
each other, hence an envelope is not defined. We also note, most other algorithms for
estimating the largest Lyapunov exponent is equivalent to compute A(k) for r < ro, where
ro is selected more or less arbitrarily, then obtain A(k)/kbt, for not too large k, as an
estimation of the largest Lyapunov exponent. With such algorithms, one can obtain a
'pm enLyapunov exponent for white noise and for 1/ f processes. However, the value of
the Lyapunov exponent estimated this way sensitively depends on the parameter ro chosen
in the computation, hence, typically is different for different researchers. We thus observe
a random element here! With these discussions, it should be clear why 1/ f processes may
be interpreted as chaotic.
One can expect that the behavior of the A(k) curves for a noisy chaotic system lies
in between that of the A(k) curves for a clean chaotic system and that of the A(k) curves
for white noise or for 1/ f processes. This is indeed so. An example is shown in Fig. 21(b)
for the noisy Lorenz system with D = 4. We observe from Fig. 21(b) that the A(k)
curves corresponding to different shells now separate. Therefore, an envelope is no longer
defined. This separation is larger between the A(k) curves corresponding to smaller shells,
indicating that the effect of noise on the sniallscale dynamics is stronger than that on the
largescale dynanxics. Also note that k, + k, is now on the order of the embedding window
length, and is almost the same for all the A(k) curves. With stronger noise (D > 4), the
A(k) curves will be more like those for white noise [21, 22].
Finally, we consider the MackeyGlass delay differential system [39]. When a=
0.2, b = 0.1, c = 10, F = 30, it has two positive Lyapunov exponents, with the largest
Lyapunov exponent close to 0.007 [21].
axr(t +F
dxr/dt = bxr(t) (26)
1 + r(t +F)"
Having two positive Lyapunov exponents while the value of the largest Lyapunov exponent
of the system is not much greater than 0, one might he concerned that it may be difficult
to deal with this system. This is not the case. In fact, this system can he analyzed as
straightforwardly as other dynamical systems including the Logistic nmap, Henon niap and
the Rossler system. An example of the A(k) vs. k curves for the chaotic AlackeyGlass
"tI
3I
O 50
200
200
Timedependent exponent A(k) vs. evolution time k curves for (a) clean and
(b) noisy Lorenz system. Six curves, from bottom to top, correspond to shells
(2(i+1)/2, 2i/2) With i 8, 9, 10, 11, 12, and 13. The sampling time for the
system is 0.03 sec, and embedding parameters are m = 4, L = 3. 5000 points
are used in the computation.
Figure 21.
100 150
0 50 100 150
Evolution timne k
3.5 
2.5
1.5
0.5
0 100 200 300 400 500 600
Evolution time k
Figure 22. Timedependent exponent A(k) vs. evolution time k curves for the
MackeyGlass system. Nine curves, from bottom to top, correspond to shells
(2(i+1)/2, 2i/2) With i 5, 6, and 13. The computation was done with
m = 5, L = 1, and 5000 points sampled with a time interval of 6.
delay differential system is shown in Fig. 22, where we have followed [21] and used
m = 5, L = 1, and 5000 points sampled with a time interval of 6. Clearly, we observe that
there exists a well defined common envelope to the A(k) curves. Actually, the slope of the
envelope estimates the largest Lyapunov exponent. This example illustrates that when
one works in the framework of the A(k) curves developed by Gao and Zheng [21, 22], one
does not need to be very concerned about nonuniform growth rate in highdimensional
systems .
2.2 General Computational Framework for Powerlaw Sensitivity to Initial
Conditions (PSIC)
To understand the essence of powerlaw sensitivity to initial conditions (PSIC),
let us first consider the 1D case and consider Eq. 27, where Ax(0) is the infinitesimal
discrepancy in the initial condition, Ax(t) is the discrepancy at time t > 0.
((t) =lim (27)
When the motion is chaotic, then
((t) = exlt (28)
((t) satisfies
Tsallis and coworkers [7, 11] have generalized Eq. 29 to
di = A,((t)V (210)
where q is called the entropic index, and X, is interpreted to be equal to K,, the
generalization of the K~olmogorovSinai entropy. Eq. 210 defines the PSIC in the 1D
case. Obviously, PSIC reduces to ESIC when q 1 The solution to Eq. 210 is
s'(t) =[1 +(1 q)Agt] 'il) (211)
When t is large and q / 1, ((t) increases with t as a powerlaw,
((t) Ct 7 ')(212)
where C = [(1 q) A,] 1 (V). For Eq. 212 to define an unstable motion with A, > 0, we
must have q < 1. Later we shall map different types of motions to different ranges of q.
To apply PSIC to the analysis of time series data, one can first construct a phase
space by constructing embedding vectors 1M as defined by Eq. 23. Eq. 211 can then be
generalized to highdimensional case [40],
((t) = Ilimo~l =lvo~l lv t 1 + (1 q)X()t (2 13
where  a V(0)  is the infinitesimal discrepancy between two orbits at time 0,  a V(t)  is
the distance between the two orbits at time t > 0, q is Ithe entropic index, and A i) is the
first qLyapunov exponent, corresponding to the powerlaw increase of the first principle
axis of an infinlitesimral ball in thle phase space. AI1 mayi not be equal to K,. TIhisi is
understood by recalling that for chaotic systems, the K~olmogorovSinai entropy is the sum
of all the positive Lyapunov exponents. We conjecture a similar relation may hold between
the qLyapunov exponents and K. When there are multiple unstable directions, then in
general A I) maty not be equal to Ki,. Whe~n t is large andc q / 1, Eq. 213 again gives a
powerlaw increase of ((t) with t. Note that under the above framework, the motion may
not be like fully developed chaos, thus not ergodic [41].
We now consider the general computational framework for PSIC. Given a finite
time series, the condition of av(0) 0 cannot be satisfied. In order to find the law
governing the divergence of nearby orbits, one thus has to examine how the neighboring
points, (1K, Vj), in the phase space, evolve with time, by forming suitable ensemble
averages. Notice that if (1MI, by) and (1%2, V/j2) arT tWO pairs of nearby points, when
 %1 Vjll  
I %2+t Vj2+tl CRllllOt be simply averaged to provide estimates for q anld Af ). Inl fact, it
would be most convenient to consider ensemble averages of pairs of points (1M, V4j) that
all fall within a very thin shell, rl < %~ g1 < r2, Where rl and T2 arT ClOSe. These
arguments so~ 1 that the timedependent exponent curves defined by Eq. 24 provide a
natural framework to assess PSIC from a time series. This is indeed so. In fact, we have
In I(t) a A(t) = In (+ +t) (214)
Now, by the discussions in Sec. 2.1, it is clear that PSIC is a generalization of ESIC:
as long as the A(k) curves from different shells form a linear envelope, then q = 1 and the
motion is chaotic. The next question is: Does PSIC also include random fractals as special
cases? The answer is yes. It will be given in the next chapter. There, we will also gain
a better undlerstanding of the meaning of A l. The re~st of this chapter shows examples
of characterizing the edge of chaos by the new framework of PSIC. Specifically, we study
time series generated from noisefree and noisy logistic and Henon maps.
2.3 Characterizing Edge of Chaos by Powerlaw Sensitivity of Initial
Conditions (PSIC)
Let us first examine the logistic map around the edge of chaos. Considering that
dynamic noise is often an important component of a system (e.g., spontaneous emission
noise in semiconductor lasers [26], and physiological noise in biological systems [30]), we
study both the deterministic and the noisy logistic map:
X,+1 = ax,(1 x,) + arl, (215)
where a is the bifurcation parameter and rlo is a white Gaussian noise with mean zero
and variance 1. The parameter a characterizes the strength of noise. For the clean system
(a = 0), the edge of chaos occurs at the accumulation point, a, = 3.569945672 We
shall study three parameter values, al = a, 0.001, a,, and a2 = a, + 0.001. When
noise is absent, al corresponds to a periodic motion with period 25, while a2 COTTOSponds
to a truly chaotic motion. We shall only study transientfree time series. In Figs. 23(ac),
we have plotted the A(t) vs. t curves for parameter values al, a,, and a2, TOSpectively.
We observe from Fig. 23(a) that the variation of A(t) with t is periodic (with period 16,
which is half of the period of the motion) when the motion is periodic. This is a generic
feature of the A(t) curves for discrete periodic attractors, when the radius of the shell is
larger than the smallest distance between two points on the attractor (when a periodic
attractor is continuous, A(t) can be arbitrarily close to 0). It has been found by Tsallis
and coworkers [11] that at the edge of chaos for the logistic map, ((t) is given by Eq. 211
with q m 0.2445. Surprisingly, we do not observe such a divergence in Fig. 23(b). In
fact, if one plots A(t) vs. In t, one only observes a curve that increases very slowly (similar
to that shown in Fig. 24(a)). The more interesting pattern is the one that is shown in
Fig. 23(c), where we observe a linearly increasingly A(t) vs. t curve. In fact, shown in
Fig. 23(c) are two such curves, corresponding to two different shells. Very interestingly,
the two curves collapse to form a common envelope in the linearly increasing part of the
6
(d)
1
05
6
(e)
O
4
1
O
0 1 23 456
In
Timedependent exponent A(t) vs. t curves for time series generated from the
noisefree logistic map with (a) al a, 0.001, where the motion is periodic
with period 2s, (b) a, 3.569945672 and (c) a2 a, + 0.001, where the
motion is chaotic. Plotted in (df) are A(t) vs. In t curves for the noisy logistic
map with o = 0.001. Very similar results were obtained when o = 0.0001.
Shown in (cf) are actually two curves, corresponding to two different shells.
104 pOintS Were used in the computation, with embedding parameters
m = 4, L = 1. However, so long as m > 1, the results are largely independent
of embedding. When m = 1, the A(t) curves are not smooth, and the
estimated 1/(1 q) value is much smaller than the theoretical value.
Figure 23.
curve. The slope of the envelope gives a good estimate of the largest positive Lyapunov
exponent. This is a generic feature of chaos [21, 22], as explained in Sec. 2.1. Since the
chaos studied here is close to the edge of chaos, the curves shown in Fig. 23(c) are less
smooth than those reported earlier [2123, 26].
Why cannot the theoretical prediction of PSIC at the edge of chaos for the logistic
map be observed from a clean time series? In a recently published very interesting and
2.0
1.5 
0.5
0 50 100 150 200 250
2.0 '
1.5
0.0 .
0 50 100 150 200 250
0 50 100 150 200 250
insightful paper, Beck [42] so__~ ; that dynamic noise may be of fundamental importance
to the Tsallis NESM. May dynamic noise pl ai similarly significant role for PSIC? As is
shown below, the answer is yes. In Figs. 23(df), we have shown the A(t) vs. In t curves
for the three parameters considered, with noise strength al = 0.001. In fact, shown
in each figure are two curves, corresponding to two different shells. They parallel with
each other. The slopes of those curves are about 1.20, close to the theoretical value of
1/(1 0.2445) a 1.32. While it is very satisfactory to observe PSIC at the edge of chaos,
it is more thrilling to observe the collapse of regular as well as chaotic motions onto the
PSIC attractor around the edge of chaos. This signifies the stableness of PSIC when there
is dynamic noise. It is important to emphasize that the results shown in Figs. 23(df)
are largely independent of the noise strength, so long as noise is neither too weak nor too
strong. For example, very similar results have been observed with a2 = 1/10 = 0.0001.
Before we move on to discuss the Henon map near the edge of chaos, let us comment
on the difference between the ESIC and the PSIC. When the ESIC is the case, the A(t)
vs. t curves are straight for T, < t < T,, where T, is a prediction time scale, and T, is a
time scale for the initial separation to evolve to the most unstable direction. When A(t)
vs. t curves corresponding to different shells are plotted together, they collapse together
for a considerable range of t. When the PSIC is the case, the A(t) vs. In t curves are
fairly straight, and the A(t) vs. In t curves corresponding to different shells separate and
parallel with each other. This feature is the fundamental reason that ensemble averages
are most conveniently formed by requiring neighboring pairs of points to all fall within
a thin shell. For example, if two distinct thin shells, described by rl < Xi Xj < T2 r
and T2 II Vy < T 3, arT joined together to form a single thicker shell, described by
rl < V~ Vy1 < T3, and one averages the separation between Xi and Xj before taking
the logarithm, then the resulting slope between A(t) and In t, assumed to be still linear,
will have to be smaller than that estimated from the two thin shells separately. At this
point it is also worth noting that the linearly increasing part of the A(t) vs. In t curves
may only contain a small interval of t. Actually, for the simple logistic nmap, A(t) saturates
around i s 20, when the radius of the shell is about 104 and the length of the time series
is 104. For experimnertal tilie series, the t for A(t) to saturate may be even smaller (we
suspect it would be about 10 samples irrespective of the sampling time, because of
the logarithmic timre scale). Tob get more accurate estimnates of the parametersi q and A ',
it may be better to fit ((t) using Eq. 213 when dealing with experimental time series.
In fact, we have found that if we do so, the estimated 1/(1 q) value from Figs. 1(df)
increases to about 1.27, much closer to the theoretical value of 1.32.
Next let us consider the Henon nmap, where ty., and ty, are white Gaussian noise, and
are uncorrelated with each other. The parameter a measures the strength of the noise.
r,z+l = 1 axr + Utz + wif, (77)
Un+l = br,z + aty,(iv)
Tirnakli [20] studied the deterministic version of this nmap at the edge of chaos, for
parameter values a,. = 1.40115518 b = O; a,. = 1.39966671 b = 0.001; and
a,. = 1.:386:37288 b = 0.01. We have studied both the noisefree and noisy nmap
for parameter values listed above, and found very similar results to those presented in
Fig. 2:3. In Figs. 24(a,h), we have plotted A(t) vs. In t curves for the noisefree and noisy
nmap for a,. = 1.:386:37288 b = 0.01. As can he observed clearly, Fig. 24(b) is very
similar to Figs. 2:3(df), while the slope for the curves in Fig. 24(a) is much smaller than
that in Fig. 24(b). Again, we have observed (but not shown by figures here, since they are
very similar to Figs. 2:3(df) and Fig. 24(b)) that dynamic noise makes the regular and
chaotic motions to collapse onto the PSIC attractor for parameters around those definingf
the edge of chaos, and that such transitions are largely independent of the noise strength.
To suninarize, we have described an easily intplenientable procedure for computationally
examining PSIC front a time series. By studying two model systems: noisefree and noisy
logistic and Henon maps near the edge of chaos, we have found that when there is no
2.0 (a)''""'""''"" 6 (b
1.5
S0.5 
0.0 1 2
0.5 1 1
1.0O ~......... ......... ......... ......... .......... ............ O .
0123456 0123456
Int Int
Figure 24. Timedependent exponent A(t) vs. t curves for time series generated from (a)
the noisefree Henon map with ac = 1.38637288 b = 0.01. Plotted in (b)
are two A(t) vs. In t curves (corresponding to two different shells) for the noisy
map with a 0.001. Similar results were obtained with a 0.0001. 104 pOintS
were used in the computation, with embedding parameters m = 4, L = 1.
noise, the PSIC attractor cannot be observed from a scalar time series. However, when
dynamic noise is present, motions around the edge of chaos, be it simply regular or
truly chaotic when there is no noise, all collapse onto the PSIC attractor. While these
examinations signify the ubiquity of PSIC, they also highlight the importance of dynamic
noise. The existence of the latter is perhaps the very reason that truly chaotic time series
can seldom be observed.
CHAPTER 3
CHARACTERIZING RANDOM FRACTALS BY POWERLAW SENSITIVITY TO
INITIAL CONDITIONS (PSIC)
In C'!s Ilter 2, we have discussed that PSIC is a generalization of ESIC: as long as
the A(t) curves from different shells form a linear envelope, then q = 1 and the motion
is chaotic. We have also shown that edge of chaos can be well characterized by PSIC.
In this chapter, we shall study two 1!! I r~~ types of random fractals: 1 fP processes with
longrangecorrelations and Levy processes. We shall show both analytically and through
numerical simulations that both types of processes can be readily characterized by PSIC.
3.1 Characterizing 1/ f Processes by Powerlaw Sensitivity to Initial
Conditions (PSIC)
A continuoustime stochastic process, X = {X(t), t > 0}, is said to be selfsimilar if
X(At) d AHX(t), t > 0 (31)
for A > 0, O < H < 1, where u~~denotes equality in distribution. H is called the
selfsimilarity parameter, or the Hurst parameter.
Before proceeding on, we note that, more rigorously speaking, processes defined by
Eq. 31 should be called selfaffine processes [43] instead of selfsimilar processes, since X
and t have to be scaled differently to make the function look similar. A simple physical
explanation for this is that the units for X and t are different. In this dissertation,
however, we will continue to call such processes selfsimilar, to follow the convention in
some of the engineering disciplines.
The following three properties can be easily derived from Eq. 31:
E[X(At)]
E [X(t)] Mean (32)
Var [X(At)]
Var[X(t)]= 2HVariance (33)
R, (t, s) = 2HAutocorrelation (34)
Note that Var[X(t)] = R,(t, t); hence, Eq. 33 can be simply derived from Eq. 34; we
have listed it as a separate equation for convenience of future reference. If we consider
only secondorder statistics, or if a process is Gaussian, Eq. 32 to Eq. 34 can be used
instead of Eq. 31 to define a selfsimilar process.
A very useful way of describing a selfsimilar process is by its spectral representation.
Strictly p.' I1:;14 the Fourier transform of X(t) is undefined, due to the nonstationary
nature of X(t). One can, however, consider X(t) in a finite interval, ;?i, O < t < T:
X (t, T) = t)0< t< T (35)
otherwise
and take the Fourier transform of X(t, T):
F(f, T) = lX(t)e2xjftd (36)
F(f, T)12df is the contribution to the total energy of X(t, T) from those components with
frequencies between f and f + df. The (average) power spectral density (PSD) of X(t, T)
is then
S(f T) = FfT
and the spectral density of X is obtained in the limit as T oo
S( f) = lim S( f ,T)
T>oo
Noting that the PSD for AHX(At) is
A2H limX(t2ft _2S/)
T>oo AT$ I 'o X.X)ZiLt 2 X"sfX
and that the PSD for AHX(At) and X(t) are the same [44], we have
S(f) = S(f/A)A2H1
The solution to the above equation is
S( f) ~ f P (37)
with
/7 = 2H + 1 (38)
Processes with power spectral densities as described by Eq. 37 are called 1/ f
processes. Typically, the powerlaw relationships that define these processes extend over
several decades of frequency. Such processes have been found in numerous areas of science
and engineering. Some of the older literatures on this subject can be found, for example,
in Press [45], Bak [46], and Wornell [47]. Some of the more recently discovered 1/ f
processes are in traffic engineering [4850], DNA sequence [5153], human cognition [54],
ambiguous visual perception [55, 56], coordination [57], posture [58], dynamic images [59,
60], and the distribution of prime numbers [61], among many others. It is further observed
that principle component analysis of such processes leads to powerlaw decaying eigenvalue
spectrum [62].
Two important prototypical models for 1/ f processes are the fractional Brownian
motion (fBm) processes and the ON/OFF intermittency with powerlaw distributed ON
and OFF periods. Below, we apply the concept of PSIC to both types of processes.
3.1.1 Characterizing Fractional Brownian Motion (fBm) Processes by Power
law Sensitivity to Initial Conditions (PSIC)
We first briefly introduce Brownian motion (Bm). Brownian motion is defined as a
(nonstationary) stochastic process B(t) that satisfies the following criteria:
1. All Brownian paths start at the origin: B(0) = 0.
2. For 0 < 1
independent.
3. For all (s, t) > 0, the variable B(t + s) B(t) is a Gaussian variable with mean 0 and
variance s.
4. B(t) is a continuous function of t.
We may infer from the above definition that the probability distribution function of B
is given by:
Pt { B(t + ) B(t) ]
This function also satisfies the scaling property:
P{[B(t + s) B(t)] < x} = P{B(At + As) B(At) < X1/2X (310)
In other words, B(t) and A1/2B(At) have the same distribution. Thus, we see from
Eq. 31 that Bm is a selfsimilar process with Hurst parameter 1/2.
Suppose we have measured B(tl), B(t2) 2~ 1 > 0. What can we ;?i about
B(s), 1 < s < 2? The answer is given by the Levy interpolation formula:
(S tl) (t2 s) 8 1l]~ 1/21
B(s) = B(tl) +[B(t2) 1 ~t
(t2 1l 2a 1
where W is a zeromean and unitvariance Gaussian random variable. The first two terms
are simply a linear interpolation. The third term gives the correct variance for B(s)B(t )
and B(t2) B 8), Which are s 1 and t2 8, TOSpectively.
A popular way of generating Bm is by the following random midpoint displacement
method, which is essentially an application of the Levy interpolation formula: suppose we
are given B(0) = 0 and B(1) as a sample of a Gaussian random variable with zero mean
and variance O.2. We can obtain B(1/2) from the following formula:
B(1/2) B(0) = n[B(1) B(0)] + D1
where D1 is a random variable with mean 0 and variance 22 2. D1 is simply the third
term of Eq. 311, and the coefficient 22 equals (t2 s) 8 1 2~t 1l), With tl = 0, t2 =
1, a = 1/2. Similarly, we have
B(1/4) B(0) = ,[B(1/2) B(0)] + D2
where D2 1S a random variable with zero mean and variance 23 2. We can apply the
same idea to obtain B(3/4). In general, the variance for D, is 2("+l)a2.
We note that a Brownian path is not a differentiable function of time. Heuristically,
this can be understood in this way: consider the variable, B(t + s) B(t), with variance
s. Its standard deviation, which is a measure of its order of magnitude, is ~ 2. Thus, the
derivative of B at t behaves like the limit Z/S = S1/2 aS s 0.
Although B(t) is almost surely not differentiable in t, symbolically one still often
writes
B(t) = w~(T)dv (312)
where w(t) is stationary white Gaussian noise, and extends the above equation to t < 0
through the convention,
JO~ Jt
It should be understood that integrals with respect to the differential element w(t)dt
should be interpreted more precisely as integrals with respect to the differential element
dB(t), in the RiemannStieltjes integral sense. This has profound consequence when
numerically integrating Eq. 312.
To see how, let us partition [0, t] into a equally spaced intervals, at = t/n. Since w(t)
are Gaussian random variables with zero mean and unit variance, in order for B(t) to have
variance t, we should have
i= 1
That is, the coefficient is (At)1/2 inStead of at, as one might have guessed. Typically,
at at. Thus, if one incorrectly uses at instead of (At)1/2, One is
severely underestimating the strength of the noise.
(a) White noise
(b) Brownian motion
Figure 31. White noise (a) and Brownian motion (b). Axes are arbitrary.
When the time is genuinely discrete, or the units of time can be arbitrary, one can
take at to be 1 unit, and Eq. 313 becomes
B(t) = w~7(i) (3 14)
i= 1
Eq. 314 provides perhaps the simplest method of generating a sample of Bm. An example
is shown in Fig. 31.
Eq. 314 is also known as a random walk process. A more sophisticated random
walk (or jump) process is given by summing up an infinite number of jump functions:
Je(t) = AiP(t ti), where P(t) is a unitstep function
1t>0
P(t) = O < (315)
and Ai, ti are random variables with Gaussian and Poisson distributions, respectively.
Next we discuss fractional Brownian motion (fBm) processes. A normalized fBm
process, Z(t), is defined as follows:
Z(t) WiH (t > 0) (316)
where W is a Gaussian random variable of zero mean and unit variance, and H is the
Hurst parameter. When H = 1/2, the process reduces to the standard Brownian motion
(Bm). It is trivial to verify that this process satisfies the defining Eq. 31 for a selfsimilar
stochastic process.
Fractional Brownian motion BH 1) is a Gaussian process with mean 0, stationary
increments, variance
E [(BH )2] t2H n7
and covariance:
E[B 8)H 2H 2PH _2HJ I
where H is the Hurst parameter. Due to its Gaussian nature, according to Eq. 32
to Eq. 34, the above three properties completely determine its selfsimilar character.
Fig. 32 shows several fBm processes with different H.
Roughly, the distribution for an fBm process is [63]:
BH~l TTii~ H/2 (319)
where F(t) is the gamma function:
0/OO"te~d 7
It is easy to verify that fBm satisfies the following scaling property:
P(BHt + ) BH~t I) = P 8H(X + ) BXH XHX} (320)
In other words, this process is invariant for transformations conserving the similarity
variable X/tH
(c) H=0.75
(d) H=0.90
Figure 32. Several fBm processes with different H.
The integral defined by Eq. 319 diverges. The more precise definition is
BH t) BH(0) =r( K1 7dB0 (321)
where the kernel K(t r) is given by [63]
K(t r) =l l (322)
The increment process of fBm, Xi = BH(( + )t BH i 1), i > i, Where at can
be considered a sampling time, is called fractional Gaussian noise (fGn) process. It is a
zero mean stationary Gaussian time series. Noting that
E(XXik) = E {[BH(Z + )t BH 8t~[ H ((i + 1 + k) at) BH ((i + k) at)]}
by Eq. 318, one obtains the autocovariance function y(k) for the fGn process:
q~~~~k)2` =. E(4ik/(X)= ( )2H 2k2H + k 12H, k >0 (3 23)
Notice that the expression is independent of at. Therefore, without loss of generality, one
could take at = 1. In particular, we have
y(1) = 222
Let us first note a few interesting properties of y(k):
(i) When H = 1/2, y(k) = 0 for k / 0. This is the wellknown property of white
Gaussian noise.
(ii) When 0 < H < 1/2, y(1) < 0.
(iii) When 1/2 < H < 1, y(1) > 0.
Properties (ii) and (iii) are often termed antipersistent and persistent correlations [43],
respectively.
Next, we consider the behavior of y(k) when k is large. Noting that when x
(1 + x)" a 1 + a~x + x
we have, when k > 1,
(k + 1)2H + k~ 12H = k~2H[(1 + 1/k)2H + (1 1/k~)2H] = k~2H[2 + 2H(2H 1)k2]
hence,
y(k) ~ H(2H 1)k2H2, Sk 00(324)
when H / 1/2. When H = 1/2, y(k) = 0 for k > 1, the Xi's are simply white noise.
Now we are ready to characterize 1/ f processes in the framework of PSIC. Recalling
that the defining property for a 1 fP process is that its variance increases with t as t2H
Irrespective of which embedding dimension is chosen, this property can be translated to
I(t) = tH (325)
Therefore ,
dt)= HtH1 (326)
Expressing t in terms of (, we have
= H( $ H327)
Comparing with the defining equation of PSIC (Eq. 210), we find that
1 2
q =1 1 (328)
H P1
X( ) = H= (329)
Noticing 0 < H < 1, from Eq. 328, we have oo < q < 0.
Computationally, the key is to demonstrate behaviors described by Eq. 325. This
can be readily done by calculating the time dependent exponent curves defined by
Eq. 214 (see also Eq. 24). In Fig. 33 we have shown a few examples for the fBm
processes with H = 0.25, 0.5, and 0.75. Clearly, the slopes of the straight lines correctly
estimate the H parameter used in the simulations.
3.1.2 Characterizing ON/OFF Processes with Pareto distributed ON and
OFF Periods by Powerlaw Sensitivity to Initial Conditions (PSIC)
ON/OFF intermittency is a ubiquitous and important phenomenon. For example,
a queueing system or a network can alternate between idle and busy periods, a fluid
flow can switch from a turbulent motion to a regular (called laminar) one. Fig. 34
shows an example for an ON/OFF traffic model with three independent traffic sources
Xi(t), i = 1, 3, where each Xi(t) is a 0/1 renewal process.
Let us denote an ON period by 1 and an OFF period by 0. We study ON/OFF trains
where ON and OFF periods are independent and both have heavytailed distribution.
Time
Figure 34. ON/OFF sources with NV = 3, X (t), X2(t), X3(t), and their summation S3 )
SON OFF ON OFF ON OFF
H = 0.25
H = 0.50
H = 0.75
___

6.5
2 2.5 3 3.5 4
Ink
Figure 33. Timedependent exponent A(k) vs. In k curves for fractional Brownian motion.
I I
X ,(t)
I I
I I
X 2t)
X3,t
S3(t
A heavytailed distribution is commonly expressed as
f (x) ~ x"1, x 00
Equivalently, one may write:
P X > x] ~ x". x 0 o 330)
This expression is often called the complementary cumulative distribution function
(CCDF) of a ]!. l.ir I.l:d distribution. In fact, it is more popular for expressing a
heavytailed distribution, since it emphasizes the tail of the distribution. It is easy to
prove that when p < 2, the variance and all higher than 2ndorder moments do not exist.
Furthermore, when p < 1, the mean also diverges.
When the powerlaw relation extends to the entire range of the allowable x, we have
the Pareto distribution:
P[X J > ~ I x]= >biU > 2> 2> (331)
where p and b are called the shape and the location parameters, respectively. It can be
proven [64] that when 1 < p < 2, the Hurst parameter of the ON/OFF processes with
Pareto distributed ON and OFF periods is given as
H = (3 p) /2 (332)
Our purpose here is to check whether the ON/OFF processes with Pareto distributed
ON and OFF periods can be characterized by PSIC. In Fig. 35 we have shown a few
examples for p = 1.2, 1.6, and 2.0. Noticing Eq. 332, we see that the slopes of the
straight lines again correctly estimate the H parameter. Also note that when 0 < p < 1,
1 < H < 1.5. Therefore, Eqs. (325) and (328) are still correct when 1 < H < 1.5. The
entire range of H = (3 p)/2, with 0 < p I 2, determines that 1 I q < 1/3 for Pareto
distributed ON/OFF processes.
4L = 2.0
4 = 1.6
4 = 1.2
2.8
3.4
.53 3.5 4 4.5
Ink
Figure 35. Timedependent exponent A(k) vs. In k curves for ON/OFF processes with
Pareto distributed ON and OFF periods.
Before ending this section, we would like to emphasize that the possibility of
misinterpreting 1/ f processes being deterministic chaos never occurs if one works
under the framework of PSIC. Under this framework, one monitors the evolution of two
nearby trajectories. If two nearby trajectories diverge exponentially fast, then the time
series is chaotic, if the divergence increases in a powerlaw manner, then the trajectories
belong to 1/ f processes. These ideas can be easily expressed precisely. Let  % Vy  
be the initial separation between two trajectories. This separation is assumed to be not
larger than some prescribed small distance r. After time t, the distance between the two
trajectories will be %~+t Vj+t. For true chaotic systems,
 K~+t Vj+t oc K 1, Vj ext
where A > 0 is called the largest Lyapunov exponent. On the other hand, for 1/ f
processes,
The latter equation thus provides another simple means of estimating the Hurst
parameter.
3.2 Characterizing Levy Processes by Powerlaw Sensitivity to Initial
Conditions (PSIC)
We now consider Levy processes, another important type of random fractal
models that have found numerous applications [6570]. There are two types of Levy
processes [71]. One is Levy flights, which are random processes consisting of many
independent steps, each step being characterized by a stable law, and consuming a
unit time regardless of its length. The other is Levy walks, where each step takes time
proportional to its length. A Levy walk can be viewed as sampled from a Levy flight
with a uniform speed. The increment process of a Levy walk, obtained by differencing
the Levy walk, is very similar to an ON/OFF train with powerlaw distributed ON
and OFF periods. Therefore, in the following, we shall not be further concerned about
it. We shall focus on Levy flights. Note that Levy flights, having independent steps,
are memoryless processes characterized by H = 1/2, irrespective of the value of the
exponent a~ characterizing the stable laws [71]. In other words, methods such as detrended
fluctuation analysis (DFA) [72] cannot be used to estimate the a~ parameter.
Now we define Levy flights more precisely. A (standard) symmetric acstable Levy
process {L,(t), t > 0} is a stochastic process with the following properties [73]:
1. L,(t) is almost surely 0 at the origin t = 0;
2. L,(t) has independent increments; and
3. L,(t) L,(s) follows an acstable distribution with characteristic function e(ts) lo~
where 0 < s < t < 00.
A random variable Y is called (strictly) stable if the distribution for CE is the same as
that for niggY
i= 1
where Yi, Y2, are independent random variables, each having the same distribution as
Y. From Eq. 333, one then finds that nVarY = n2/"VarY. For the distribution to be
valid, O < a~ < 2. When a~ = 2, the distribution is Gaussian, and hence, the corresponding
Levy flight is just the standard Brownian motion. When 0 < a~ < 2, the distribution is
heavytailed, P[X > x] ~ x", x 00o, and has infinite variance. Furthermore, when
0 < a~ < 1, the mean is also infinite.
At this point, it is worth spending a few paragraphs to explain how to simulate an
acstable Levy process. The key is to simulate its increment process, which follows an
acstable distribution. The first breakthrough for the simulation of stable distributions
was made by K~anter [74], who gave a direct method for simulating a stable distribution
with a~ < 1 and p = 1. The approach was later generalized to the general case by
C'I I1!.1wers et al. [75]. We first describe the algorithm for constructing a standard stable
random variable X ~ S(a~, p; 1).
Theorem Simulation of stable random variables: Let 8 and W be independent
with 8 uniformly distributed on (xr/2, xr/2), W exponentially distributed with mean 1.
For 0 < a~ < 2, 1 < p < 1, and 8o = ard i I[Sltan(;ac/2)]/a,
Z (cosina(00+8) cos[aeo+(a1)8] 1a/ 34
(cos000cos8)1/a
has a S(a~, p; 1) distribution.
For symmetric stable distribution with P = 0, the above expression can be greatly
simplified. To simulate an arbitrary stable distribution S(a~, P, y, 5), we can simply take
the following linear transform,
Y yZ + 6 a rJ 1 (335)
yZ+( S"qln n+6 a~= 1
where Z is given by Eq. 334.
We have simulated a number of symmetric acstable Levy motions with different
a~. Figure 36 shows two examples of Levy motions with a~ = 1.5 and 1. The difference
between the standard Brownian motions and Levy motions lies in that Brownian motions
appear everywhere similarly, while Levy motions are comprised of dense clusters connected
by occasional long jumps: the smaller the a~ is, the more frequently the long jumps appear.
Let us now pause for a while and think where Levy flightslike esoteric processes
could occur. One situation could be this: a mosquito head to a giant spider net and got
stuck; it struggled for a while, and luckily, with the help of a gust of wind, escaped. When
the mosquito was struggling, the steps it could take would be tiny; but the step leading
to its fortunate escape must be huge. As another example, let us consider (American)
football games. For most of the time during a game, the offense and defense would be
fairly balanced, and the offense team may only be able to proceed for a few yards. But
during an attack that leads to a touchdown, the hero getting the touchdown often "flies"
tens of yards somewhat he has escaped the defense. While these two simple situations
are only meant as analogies, remembering them could be helpful when reading research
papers searching for Levy statistics in various types of problems such as animal foraging.
At this point, we also wish to mention that Levyflightslike patterns have been used as a
type of screen saver for Linux operating systems.
The symmetric acstable Levy flight is 1/a~ selfsimilar. That is, for c > 0, the
processes {L,(ct), t > 0} and {cl/oL,(t), t > 0} have the same finitedimensional
distributions. By this argument as well as Eq. 333, it is clear that the length of the
(a) (b)
(c) (d)
(e) (f)
Figure 36. Examples for Brownian motions and Levy motions. 1D and 2D Brownian
motions (a,b) vs. Levy motions with a~ = 1.5 and 1 for (c,d) and (e,f),
respectively.
o a =1.0
4.5 0 a =1.5
5
5.5
6
06~.5 1 1.5 2 2.5 3
Ink
Figure 37. Timedependent exponent A(k) vs. In k curves for Levy processes.
motion in a time span of at, AL(At), is given by the following scaling:
AL(At) oc Atl/o (336)
This contrasts with the scaling for fractional Brownian motion:
ABHat oc OC (33H
We thus identify that 1/a ]I phi the role of H. Therefore,
q = 1 a~ (338)
A 0 = 1/al (339)
Noticing 0 < a~ < 2, from Eq. 338, we have 1 < q < 1.
We have applied the concept of PSIC to Levy processes with different a~. Examples
for Levy processes with a~ = 1 and 1.5 are shown in Fig. 37. The slopes of the straight
lines correctly estimate the values of a~ used in the simulations.
CHAPTER 4
STUDY OF INTERNET DYNA1\ICS BY POWERLAW SENSITIVITY TO INITIAL
CONDITIONS (PSIC)
Time series analysis is an important exercise in science and engineering. One of the
most important issues in time series analysis is to determine whether the data under
investigation is regular, deterministically chaotic, or simply random. Along this line, a
lot of work has been done in the past two decades. A hit surprisingly, quite often the
exponential sensitivity to initial conditions (ESIC) in the rigorous mathematical sense
cannot he observed in experimental data. While current consensus is to attribute this
fact to the noise in the data, we shall show that the more general concept of powerlaw
sensitivity to initial conditions (PSIC) provides an interesting framework for time series
analysis. We illustrate this by studying the complicated dynamics of Internet transport
protocols.
The Internet is one of the most complicated systems that man has ever made.
In recent years, two types of fascinating multiscale behaviors have been found in the
Internet. One is the temporaldomain fractal and multifractal properties of network traffic
(see [76, 77] and references therein), which typically represents .I_ _regation of individual
host streams. The other is the spatialdomain scalefree topology of the Internet (see [78]
and references therein). Also, it has been observed that the failure of a single router may
trigger routing instability, which may be severe enough to instigate a route flap storm [79].
Furthermore, packets may be delivered out of order or even get dropped, and packet
reordering is not a pathological network behavior [80]. As the next generation Internet
applications such as remote instrument control and computational steering are being
developed, it is becoming increasingly important to understand the transport dynamics
in order to sustain the needed control channels over widearea connections. Typically, the
transport dynamics over the Internet connections are a result of the nonlinear dynamics
of the widelydeploi II Transmission Control Protocol (TCP) interacting with the Internet
traffic, which is stochastic and often selfsimilar. This leads one to naturally expect the
transport dynamics to be highly complicated. To tackle the hard problem of studying
the transport dynamics, we first systematically carry out network simulation using a
coninonly used network simulator, the ns2 simulator (we have chosen the Tahoe version
of TCP, see http://www.isi. edu/nsnam/ns/), and then carefully study the dynamics of
the TCP congestion windowsize traces collected over the real Internet connections.
4.1 Study of Transport Dynamics by Network Simulation
By carrying out network simulation, we critically examine whether TCP dynamics
can he chaotic, and if yes, to identify a route to chaos. We shall also present Internet
measurements to complement the simulation results. Since the advent of the Internet,
it has been perceived that the dynamics of complicated coninunication networks may
be chaotic [81]. If transport dynamics are indeed chaotic, then steadystate analysis and
study of convergence to a steadystate will not he the sole focus in practice, particularly
for remote control and computational steering tasks. Rather, one should also focus on the
Il II! 1!1" nonconvergent transport dynanxics. Although recently a lot of research has
been carried out [8288] to better understand this issue, a lot of fundamental problems
remain to be answered. For example, the observation of chaos in [84] cannot he repeated,
which leads to the suspicion that the chaotic motions reported in [84] may correspond to
periodic motions with very long periods [89]. Fundamentally, we shall show, by developing
a 1D discrete nmap, that the network scenarios with only a few competing TCP flows
studied in [84] cannot generate chaotic dynanxics. However, with other parameter settings,
chaos is possible.
In order to examine the temporal evolution of a dynamical systent, one often
measures a scalar time series at a fixed point in the state space. Takens embedding
theorem [3436 ensures that the collective behavior of the dynamics described by the
dynamical variables coupled to the measured one can he conveniently studied by the
measured scalar time series. To illustrate the qualitative aspects of TCP dynamics by only
measuring a scalar time series, it is most effective to consider a single bottleneck link, as
such a link must be tightly coupled with the rest of the network. For simplicity, we shall
follow the setup of Veres and Boda [84] consisting of long lived TCP flows on a single link.
This setup has four parameters, the link speed C in Mbps, the link propagation delay d in
ms, the buffer size B in unit of packets of 1000 bytes, and the number of competing TCP
flows NV. We have found that NV acts as a critical bifurcation parameter. For example,
when C = 0.1, d = 10 ms, B = 10 and NV is small, we have observed only periodic
and quasiperiodic motions (Figs. 41(ad)). When NV is large, such as NV = 17, chaos is
observed (see Fig. 42(a,b)). We thus conclude that one route to chaos in TCP dynamics
is the wellknown quasiperiodic route. In order to show more diversified dynamical
behavior such as quasiperiodic motions with more than two incommensurate (i.e.,
independent) frequencies, below, however, we shall vary all the parameters. In particular,
we shall show that when NV is small, chaos cannot happen.
For each of the NV competing TCP flows, one can record its congestion window
size data. Let us only record one of them, and denote it by W(i). An example of
quasiperiodic W(i) is shown in Fig. 41(a), where C = 0.1 Mbps, delay d = 10
ms, B = 10, NV = 3. Its power spectral density (PSD) is shown in Fig. 41(b). We
observe many discrete sharp peaks. Note that when only around 105 points are used
to compute the PSD, one does not observe ]rn lw: peaks in the spectrum. Rather one
observes a broad spectrum, and hence, is tempted to interpret the W(i) time series
to be chaotic. Therefore, the W(i) time series is not ideally suited for the study of
(quasi)periodic motions with long periods. This motivates us to define a new time
series, T(i), which is the time interval between the onset of two successive exponential
backoff (i.e., multiplicative decrease) of TCP, as indicated in Fig. 41(a). The start of an
exponential backoff indicates triggering of a loss episode or heavy congestion. Hence, the
T(i) is related to the wellknown round trip time (RTT). In fact, it can be considered
to be more or less equivalent to the time interval between two successive loss bursts.
The T(i) time series for the W(i) of Fig. 41(a) is shown in Fig. 41(c), together with its
1015
1010
0 105
100
10
0 0.05 0.1
Frequency
10'" (d)
0 1000 2000 3000 4000 5000
index i
50
index i
Frequency
105
100
105
0 0.1 0.2 0.3 0.4 0.5
Frequency
Figure 41.
Example of quasiperiodic congestion window size IT(i) time series sampled by
an equal time interval of 10 ms (a) and its power spectral density (PSD) (b).
The parameters used are C = 0.1 Mbps, delay d = 10 ms, B = 10 packets,
where a packet is of size 1000 bytes, and the number of TCP streams NV = :3.
(c) The T(i) time series corresponding to (a) and (d) its PSD. (e) Another
T(i) time series corresponding to C = 0.5 Mbps, d = 10 ms, B = 10, NV = :3,
and (f) its PSD.
PSD in Fig. 41(d). We observe that this T(i) time series is a simple periodic sequence,
and hence, the IT(i) time series is quasiperiodic. However, T(i) time series can still be
quasiperiodic. An example is shown in Fig. 41(e), with C
0.5 Mbps, delay d
ms, B = 10 packets, NV = :3, together with its PSD is Fig. 41(f). Note that there are still
many discrete sharp peaks in Fig. 41(f). To appreciate how many independent frequencies
100
90
80
0 50 100 150
index i
Iruilu~
may be contained in the T(i) time series of Fig. 41(e), we have further extracted two new
time series, Til) (i), which is the interval time series between the successive local maxima
of the T(i) time series, and T(2) i), Which is the interval time series between the successive
local maxima of the Til) (i) time series. We have found that the T(2) i) iS Simply periodic.
Hence, we conclude that the W(i) time series is quasiperiodic with four independent
frequencies.
Before we proceed to discuss the chaotic TCP dynamics, we note that simple periodic
or quasiperiodic W(i) time series with less or more than four independent frequencies
can all be observed. In fact, we have found the "chaotic" congestion window size data
studied in [84] to be quasiperiodic with two independent frequencies. Furthermore, the
probability distributions for the constant, periodic, and quasiperiodic T(i) time series are
only composed of a few discrete peaks. This II__ ; that the observed quasiperiodic T(i)
time series constitutes a discrete torus.
Let us now turn to the discussion of chaotic TCP dynamics. Shown in Figs. 42(a,c)
are two irregular T(i) time series. The parameters C, d, and B used for Fig. 42(a) are the
same as those for Figs. 41(ad). Noting from Fig. 42 that T(i) often can be on the order
of 104 to 10s, hence, a not too long W(i) time series is even more illsuited for the study of
the underlying irregular dynamics.
To show that the T(i) time series in Figs. 42(a,c) is indeed chaotic, we employ
the direct dynamical test for deterministic chaos developed by Gao and Zheng [21, 22].
The A(k) curves for the T(i) time series of Figs. 42(a,c) are shown in Figs. 42(b,d),
where the four curves, from bottom to top, correspond to shells of sizes (2(i+1)/2, 2i/2)
i = 14, 15, 16, 17. we have simply chosen m = 6 and L = 1. Very similar curves have
been obtained for other choices of m and L. We observe a very well defined linear common
envelope at the lower left corner of Figs. 42(b,d). The existence of a common envelope
guarantees that a robust positive Lyapunov exponent will be obtained no matter which
shell is used in the computation. Hence, the two T(i) time series are indeed chaotic.
(c)
5
1
5
x14
8 5
2 15
xr 101
(d)
5
4
3
2
0 0
0 2000 4000 6000 0 10 20 30 40 50
Figure 42. Two examples of chaotic T(i) time series obtained with parameters (a)
C = 0.1, delay d = 10 ms, B = 10, and the number of TCP streams NV = 17
and (c) C = 0.1, d = 10 ms, B = 12, and NV = 19. The A(k) curves for (a) and
(c) are shown in (b) and (d), respectively.
We have remarked that the probability distributions for the T(i) time series of regular
motions are comprised of a few discrete peaks. What are the distributions for the chaotic
T(i) time series such as shown in Figs. 42(a,c)? They are powerlawlike, P(T < t) ~ t Y,
for almost two orders of magnitude in t, as shown in Figs. 43(a,b). The exponent y is not
a universal quantity, however.
We now develop a simple 1D map to describe the operation of TCP. TCP relies
on two mechanisms to set its transmission rate: flow control and congestion control.
10' 10'
10 a 10
103~ 103
141 14
102 103 104 105 102 103 104 105
t t
Figure 43. Complementary cumulative distribution functions (CCDFs) for the two chaotic
T(i) time series of Figs. 42(a,c).
Flow control ensures that the sender sends no more than the size of the receiver's last
advertised flowcontrol window based on its available buffer size while congestion control
ensures that the sender does not unfairly overrun the network's available bandwidth.
TCP implements these mechanisms via a flowcontrol window ( fwnd) advertised by the
receiver to the sender and a congestioncontrol window (cwnd) adapted based on inferring
the state of the network. Specifically, TCP calculates an effective window (ewnd), where
ewnd = min( fwnd, cwnd), and then sends data at a rate ewnd/RTT, where RTT is the
round trip time of the connection. Recently Rao and Chua [85] developed an analytic
model describing TCP dynamics. By assuming fwnd being fixed and the TCP's slow start
only contributing to transient behavior, Rao and Chua have developed the wup~date map
M~: [1, Wmax]i [1, Wmax],
I, + 1/l, if It, E1 RI, r tE R1
M~re)= < It < n /2"" if It, E R1, I tE R2
I /2"" if n, 6R2, I, 1 1 6
I, + 1/l, if no Ea R I, 1I R1 1
In region R1, there is no packet loss, as = 0, and the TCP is in additive increase mode. In
R2, there are ni > 0 packet losses, where ni may be a complicated function of time, and
the TCP is in multiplicative decrease model.
The above map (Eq. 41) can be simplified [90],
I / 2"' if n , > Wmax
Since fwnd is not assumed to be a constant, the modified map is more general than
the original Rao and Chua's map. To apply the above map to the study of a TCP
flow competing with NV 1 other TCP flows, one may lump the effect of the NV 1
other TCP flows as a background traffic. Hence, Wmax is determined by cwnd and the
background traffic, and typically is a complicated function of time, noticing that the
available bandwidth of a bottleneck link can vary with time considerably due to the
dynamic nature of background traffic. We have carried out simulation studies of Eq. 42
by assuming Wmax to be periodic and quasiperiodic, and have indeed observed chaos.
We are now ready to understand why chaos cannot occur when the number of
competing TCP flows is small, such as 2. To see this, let us assume that I, and no +
al, as well as I,. It and I, t + Al, t are all smaller than Wmax. It is then easy to
see that 1An' ,I t < An . That is, small disturbance decays. This means during the
additive increase phase of TCP, nearby trajectories contract instead of diverge. Hence,
the dynamics are stable. In order for the dynamics to be unstable so that the Lyapunov
exponent is positive, the transition from the additive increase phase to the multiplicative
decrease phase has to be fast. This can be true only when the number of competing TCP
streams is large. We thus conclude that the network scenarios considered in [84] cannot
generate chaotic TCP dynamics.
How relevant is our simulation result to the dynamics of the real Internet? Are the
irregular T(i) time series shown in Figs. 42(a,c) an artifact of Tahoe version of TCP
S10
400
200
0 CIVYUU yIYYUU"UU VRl V'IIIIYUII IIIVV 111 10
0 50 100 150 200 250 300 350 100 101 102 103
Figure 44. The time series T(i) extracted from a W(i) time series collected using net100
instruments over ORNLLSU connection. (a) the time series T(i), (b) the
complementary cumulative distribution function (CCDF) for the T(i) time
series.
we used or more generally of the ns2 simulator? To find an answer, we collected W(i)
measurements between ORNL and Louisiana State University at millisecond resolution
using the net100 instruments and computed T(i) as shown in Fig. 44(a), whose profile is
qualitatively quite similar to Fig. 42(a). In fact, it also follows a (truncated) powerlaw,
as shown in Fig. 44(b). The exponent for the powerlaw part is even smaller than the
time series of Fig. 42. The truncation in the powerlaw could be due to the shortness
of the T(i) time series. As we have pointed out, the T(i) time series is related to RTT.
Cottrell and Bullot have been experimenting with many advanced versions of TCP,
and also observed RTT time series very similar to these. We also have observed from
RTT and loss data measured on geographically dispersed paths on the Internet (with
a resolution of 1 sec) by researchers at the San Diego Supercomputer Center that the
probability distributions for the time interval between successive loss bursts typically
follow a powerlawlike distribution, with the exponent y also smaller than 1. Thus,
we have good reason to believe that the observed irregular T(i) time series is not an
artifact of the ns2 simulator, but reflects reality to some degree. In fact, the complicated
dynamics described here is much simpler than the actual dynamics of the Internet,
considering that there are all kinds of randomness in the Internet, especially that typical
TCP flows are not long lived [91], but vary considerably in durations determined by
specific applications. In the next section, we will analyze the actual dynamics of the
Internet measurements.
4.2 Complicated Dynamics of Internet Transport Protocols
As we have discussed earlier, it is very difficult to understand the dynamics of a
transport protocol over Internet connections. The 1! in r~ difficulties lie in two aspects.
First, it has been technologically hard to collect high quality measurements of network
transport variables on hli,' Internet connections. Secondly, it is very difficult to
analytically handle the interaction between the deterministic and stochastic parts of
transport dynamics. The deterministic dynamics are due to the highly nonlinear behavior
of TCP's Additive Increase and Multiplicative Decrease (AIMD) method emploix & for
congestion control. The stochastic component is due to the often selfsimilar traffic [92]
with which TCP competes for bandwidth and router buffers. Roughly speaking, the
interaction between the two is in terms of the additions to TCP congestion windowsize,
denoted by cwnd, in response to acknowledgments, and multiplicative decreases in
response to inferred losses (ignoring the initial slowstart phase). Any protocol, however
simple its dynamics are, is generally expected to exhibit apparently complicated
dynamics due to its interaction with the stochastic network traffic. TCP in particular
is shown (albeit in simulation) to exhibit chaotic or chaoslike behavior even when the
competing traffic is much simpler such as a single competing TCP [84] or User Datagram
Protocol (UDP) stream [85]. The difficulty is to understand the dynamics when both the
phenomena are in effect, and equally importantly to understand their impact on actual
Internet streams. Historically, the methods from standard chaos and stochastic theories
have been unable to offer much new perspective on the transport dynamics.
In this section, we utilize the recently developed net100 instruments (obtained front
http://www.net100. org) to collect high quality TCP cwnd traces for actual Internet
connections. We analyze the measurements using the concepts of time dependent exponent
curves [21, 22, 37]. Our major purpose is to elucidate how the deterministic and stochastic
components of transport dynamics interact with each other on the Internet. It is of
considerable interest to note that recently there have been several important works on
TCP dynamics with the purpose of improving congestion control [9396]. In fact, there
has been considerable effort in developing new versions and/or alternatives to TCP so
that the network dynamics can he more stable. Conventional v 0<~ Of analyzing the
network dynamics are unable to readily determine whether newer methods result in stable
transport dynamics or help in designing such methods. Our analysis shows two important
features in this direction:
(a) Randoniness is an integral part of the transport dynamics and must he explicitly
handled. In particular, the step sizes utilized for adjusting the congestion parameters
must he suitably varied (e.g. using Robbins1\onro conditions) to ensure the eventual
convergence [97]. This is not the case for TCP which utilizes fixed step sizes.
(b) The chaotic dynamics of the transport protocols do have a significant impact on
practical transfers, and the protocol design must he cognizant of its effects. In
particular, it might he worthwhile to investigate protocols that do not contain
dominant chaotic regimes, particularly for remote control and steering applications.
The traditional transport protocols are not designed to explicitly address the above two
issues, but justifiably so since their original purpose is data transport rather than finer
control of dynanxics.
In the rest of this section, we shall first briefly describe the ownd data studied here.
Then we shall briefly explain the analysis procedure and describe typical results of the
analysis.
We have collected a number of cwnd traces using single and two competing TCP
streams. This data was collected on two different connections front Oak Ridge National
Laboratory (ORNL) to Georgia Institute of Technology (GaTech) and to Louisiana State
University (LSU). The first connection has highbandwidth (OC192 at 10Gbps) with
relatively low backbone traffic and a roundtrip time of about 6 milliseconds. Four traces
were collected for this connection, two with a single TCP stream and the others with two
competing TCP streams. The second connection has much lower bandwidth (10 Mbps)
with higher levels of traffic and a round triptime of about 26 milliseconds. The sampling
time is approximately 1 millisecond, with an error of 10s of microseconds (due to a "busy
waitingt measurement loop). Eight traces were collected for the second connection. The
results based on these eight traces are qualitatively very similar to the ORNLGaTech
traces, and we only discuss the results for the latter. To appreciate better what these data
look like, we have plotted in Figs. 45(ad) segments of four datasets for ORNLGaTech
connection. Power spectral analysis of these data does not show any dominant peaks,
and hence, the dynamics are not simply oscillatory. Since our data was measured on the
Internet with Int ' background traffic, it is apparently more complicated and realistic
than the traces obtained by network simulation.
Next let us analyze cwnd data, x(i), i = 1, n, using the concepts of time
dependent exponent A(k) curves [21, 22, 37]. As has been discussed earlier, we first
need to employ the embedding theorem [3436, 98] to construct vectors of the form:
K = [x(i), x(i + L), ..., x(i + (m 1)L)], where m is the embedding dimension and L the
delay time. The embedding theorem [3436] states that when the embedding dimension
m is larger than twice the box counting dimension of the attractor, then the dynamics
of the original system can be studied from a single scalar time series. We note that the
embedding dimension used in [84] is only two, which has to be considered not large enough
(this may call for a closer examination of their conclusions).
Before we move on, we point out a few interesting features of the A(k) curves:
(i) For noise, only for k up to the embedding window size (m 1)L will A(k) increase.
Thus, whenever A(k) increases for a much larger range of k, it is an indication of
nontrivial deterministic structure in the data.
(ii) For periodic signals, A(k) is essentially zero for any k.
x 104
4
rl
qX104
(b) 1 TCP source

3
1
0
0
2000 4000 6000
2000 4000 6000
x 104
3~ i
X 104
41(d) 2 TCP sources
[Hdi n8 n, Ellrfli rfl
TCP sources
1
2000 4000
sample n
6000
2000 4000
sample n
6000
Figure 45. Time series for the congestion window size cwnd for ORNLGaTech
connection.
(iii) For quasiperiodic signals, A(k) is periodic with an amplitude typically smaller than
0.1, hence, for practical purposes, A(k) can be considered very close to 0.
(iv) When a chaotic signal is corrupted by noise, then the A(k) curves break themselves
away from the common envelope. The stronger the noise is, the more the A(k)
curves break away till the envelope is not defined at all. This feature has actually
been used to estimate the amount of both measurement and dynamic noise [99].
Now we are ready to compute and understand the A(k) curves for cwnd traces. The
set of A(k) curves, corresponding to Fig. 45 are plotted in Figs. 46. In the computations,
3 x 104 pOintS are used, and m
10, L = 1. The seven curves, from the bottom to
(a) 1 TCP source
3
1
0
0 50 100 150 200 250
4
(c 12 TC urces
3
2
1
0 50 100 150 200 250
50 100 150 200
Evolution time k
50 100 150 200
Evolution time k
250
Figure 46.
Timedependent exponent A(k) vs. k curves for cwnd data corresponding to
Fig. 45. In the computations, 3 x 104 pOintS are used, and m 10, L 1 .
Curves from the bottom to top correspond to shells with sizes (2(i+1)/2, 2i/2)
i = 2, 3, ,8.
top, correspond to shells of sizes (2(i+1)/2, 2i/2), i
2, 3, 8. We make the following
observations:
(i) The dynamics are complicated and cannot be described as either periodic or
quasiperiodic motions, since A(k) is much larger than 0.
(ii) The dynamics cannot be characterized as pure deterministic chaos, since in no case
can we observe a welldefined linear envelope. Thus the random component of the
dynamics due to competing network traffic is evident and can not simply be ignored.
(iii) The data is not simply noisy, since otherwise we should have observed that A(k) is
almost flat when k > (m 1)L. Thus, the deterministic component of dynamics
which is due to the transport protocol phI i an integral role and must be carefully
studied. The features (ii) and (iii) indicate that the Internet transport dynamics
contains both chaotic and stochastic components.
(iv) There are considerable differences between the data with only 1 TCP source and
with 2 competing TCP sources. In the latter case, the A(k) curves sharply rise when
k just exceeds the embedding window size, (m 1)L. On the other hand, A(k)
for Fig. 46(a) with only one TCP source increases much slower when k just exceeds
(m1)L. Also important is that the A(k) curves in Figs. 46(c,d) are much smoother
than those in Figs. 46(a,b). Hence, we can ;?i that the deterministic component of
the dynamics is more visible when there are more than 1 competing TCP sources
(along the lines of [84]).
Since the increasing part of the A(k) curves are not very linear, let us next examine if
the cwnd traces can be better characterized by the more generalized concept of powerlaw
sensitivity to initial conditions (PSIC). Figure 47 shows the A(k) vs. In k curves for
ORNLGaTech connection. Very interestingly, we have now indeed observed a bunch of
better defined linear lines, especially for small scales. In particular, the transport dynamics
with more than 1 competing TCP sources are better described by the concept of PSIC.
To summarize, by analyzing a number of high quality congestion windowsize data
measured on the Internet, we have found that the transport dynamics are best described
by the concept of PSIC, especially for small scales. It is interesting to further examine
how one might suppress the stochasticity of the network by executing more controls on
the network when making measurements, such as using a large number of competing TCP
sources together with a constant UDP flow. Our analysis motivates that both chaotic and
stochastic aspects be paid a close attention to in designing Internet protocols that are
required to provide the desired and tractable dynamics.
02 4
Timedependent exponent A(k) vs. In k curves for cwnd data corresponding to
Fig. 45. In the computations, 3 x 104 pOintS are used, and m 10, L 1 .
Curves from the bottom to top correspond to shells with sizes (2(i+1)/2, 2i/2)
i = 2, 3, ,8.
Figure 47.
CHAPTER 5
STITDY OF SEA CLITTTER RADAR RETITRNS BY POWERLAW SENSITIVITY TO
INITIAL CONDITIONS (PSIC)
Understanding the nature of sea clutter is crucial to the successful modeling of
sea clutter as well as to facilitate target detection within sea clutter. To this end, an
important question to ask is whether sea clutter is stochastic or deterministic. Since
the complicated sea clutter signals are functions of complex (sonietintes turbulent) wave
motions on the sea surface, while wave motions on the sea surface clearly have their
own dynamical features that are not readily described by simple statistical features,
it is very appealing to understand sea clutter hv considering some of their dynamical
features. In the past decade, Haykin et al. have carried out analysis of some sea clutter
data using chaos theory [100, 101], and concluded that sea clutter was generated hv an
underlying chaotic process. Recently, their conclusion has been questioned by a number of
researchers [102108]. In particular, Unsworth et al. [105, 106] have demonstrated that the
two main invariants used by Haykin et al. [100, 101], namely the in! I::!!Itsinl likelihood of
the correlation dimension alI lin II. and the "false nearest neighbors" are problematic in
the analysis of measured sea clutter data, since both invariants may interpret stochastic
processes as chaos. They have also tried an improved method, which is based on the
correlation integral of Grassherger and Procaccia [109] and has been found effective in
distinguishing stochastic processes front chaos. Still, no evidence of deternxinisni or chaos
has been found in sea clutter data.
To reconcile ever growing evidence of stochasticity in sea clutter with their chaos
hypothesis, recently, Haykin et al. [110] have II_tr 1.. that the nonchaotic feature of sea
clutter could be due to many types of noise sources in the data. To test this possibility,
1\kDonald and Dantini [111] have tried a series of lowpass filters to remove noise; but
again they have failed to find any chaotic features. Furthermore, they have found that
the coninonly used chaotic invariant measures of correlation dimension and Lyapunov
exponent, computed by conventional r .4, produce similar results for measured sea
clutter returns and simulated stochastic processes, while a nonlinear predictor shows little
intprovenient over linear prediction.
While these recent studies highly ell_ r that sea clutter is unlikely to be truly
chaotic, a number of fundamental questions are still unknown. For example, most of the
studies in [102106] are conducted by comparing measured sea clutter data with simulated
stochastic processes. We can ask: can the nonchaotic nature of sea clutter he directly
demonstrated without resorting to simulated stochastic processes'? Recognizing that
simple lowpass filtering does not correspond to any definite scales in phase space, can we
design a more effective method to separate scales in phase space and to test whether sea
clutter can he decomposed as signals plus noise'? Finally, will studies along this line he of
any help for target detection within sea clutter'?
In this chapter, we employ the direct dynamical test for deterministic chaos discussed
in Sec. 2.1 to analyze 280 sea clutter data measured under various sea and weather
conditions. The method offers a more stringent criterion for detecting lowdintensional
chaos, and can simultaneously monitor motions in phase space at different scales.
However, no chaotic feature is observed front any of these data at all the different scales
examined. But interestingly, we show that scaledependent exponent corresponding to
large scale appears to be useful for distinguishing sea clutter data with and without
targets. This II_ is that sea clutter may contain interesting dynamic features and that
the scaledependent exponent may be an important parameter for target detection within
sea clutter. Furthermore, we find that sea clutter can he conveniently characterized by the
new concept of powerlaw sensitivity to initial conditions (PSIC). We show that for the
purpose of detecting targets within sea clutter, PSIC offers a more effective framework.
Below, we shall first briefly describe the sea clutter data. Then we study the
sea clutter data by employing the direct dynamical test for lowdintensional chaos
developed by Gao and Zheng [21, 22]. And we show that the scaledependent exponent
corresponding to large scale can he used to effectively detect low observable targets within
sea clutter, while the conventional Lyapunov exponent fails in such a task. Finally, we
apply the new concept of PSIC to detect targets within sea clutter.
5.1 Sea Clutter Data
14 sea clutter datasets were obtained from a website maintained by Professor Simon
Haykin: http: //soma.ece.mcmaster. ca/ipix/dartmouth/datasets.html. The measurement
was made using the McMaster IPIX radar at Dartmouth, Nova Scotia, Canada. The
radar was mounted in a fixed position on land 2530 m above sea level, with the operating
(carrier) frequency 9.39 GHz (and hence a wavelength of about 3 cm). It was operated at
low grazing angles, with the antenna dwelling in a fixed direction, illuminating a patch
of ocean surface. The measurements were performed with the wave height in the ocean
varying from 0.8 m to 3.8 m (with peak height up to 5.5 m), and the wind conditions
varying from still to 60 km/hr (with gusts up to 90 km/hr). For each measurement, 14
areas, called antenna footprints or range bins, were scanned. Their centers were depicted
as B1, B2, B14 in Fig. 51. The distance between two .Il11 Il:ent range bins was 15 m.
One or a few range bins (;?,, Bi_l, Bi and Bi 1) hit a target, which was a spherical block
of styrofoam of diameter 1 m, wrapped with wire mesh. The locations of the three targets
were specified by their azimuthal angle and distance to the radar. They were (128 degree,
2660 m), (130 degree, 5525 m), and (170 degree, 2655 m), respectively. The range bin
where the target is strongest is labeled as the primary target bin. Due to drift of the
target, bins .Il11 Il:ent to the primary target bin may also hit the target. They are called
secondary target bins. For each range bin, there were 217 complex numbers, sampled with
a frequency of 1000 Hz. Amplitude data of two polarizations, HH (horizontal transmission,
horizontal reception) and VV (vertical transmission, vertical reception) were an~ llh. 1
Fig. 52 shows two examples of the typical sea clutter amplitude data without and
with target. However, careful examination of the amplitude data indicates that 4 datasets
are severely affected by clipping. This can be readily observed from Figs. 53(a, b). We
40 60 80 100 120 140
Time (sec)
11 11_1 11* I)
h : Antenna height
S: Grazing angle
R. : Range (distance from the radar)
B1 ~ B14 : Range bins
Target
BB..B. Bi B.+ ... B,
1 2 11 1 111
(Secondary)(Primary)(Secondary)
I~RRi
Figure 51. Collection of sea clutter data
25
20
o 15
E 10
5
(b)
E 1
0 20
rrl
Figure 52. Typical sea clutter amplitude data (a) without and (b) with target.
discard those 4 datasets, and >.1, lli. .. the remaining 10 measurements, which contain 280
sea clutter time series.
5.2 Non Chaotic Behavior of Sea Clutter
In this section, we employ the direct dynamical test for deterministic chaos developed
by Gao and Zheng [21, 22], to 2.1, lli. .. sea clutter data. The method offers a more
stringent criterion for lowdimensional chaos, and can simultaneously monitor motions
in phase space at different scales. The method has found numerous applications in the
0 20 40 60 80 100 120 140
Time (sec)
2.5 2.5
02 02
E 1.5 E 1.5
1 1
0.5 0.5
4.00 4.01 4.02 4.03 4.04 5.22 5.23 5.24 5.25 5.26
Time (sec) Time (sec)
Figure 5:3. Two short segments of the amplitude sea clutter data severely affected hv
clipping.
study of the effects of noise on dynamical systems [2:3, 24, 26, 27], estimation of the
strength of measurement noise in experimental data [28, 29], pathological tremors [:30],
shearthickening surfactant solutions [:31], dilute sheared aqueous solutions [:32], and
serrated plastic flows [:33]. In particular, this method was used by Gao et al. [107, 108] to
analyze one single set of sea clutter data. While chaos was not observed front that dataset,
no definite general conclusion could be reached, due to lack of large amount of data at
that time. Here we systematically study 280 sea clutter data measured under various sea
and weather conditions, and examine whether any chaotic features can he found front
these sea clutter data.
The explicit incorporation of scales in the Gao and Zheng's test [21, 22] enables us to
simultaneously monitor motions in phase space at different scales. We have systematically
analyzed 280 amplitude sea clutter time series measured under various sea and weather
conditions. However, no chaotic feature has been observed front any of these data at
all the scales examined. Typical examples of the A(k) vs. k curves for the sea clutter
amplitude data without targets are shown in Figs. 54(a, c, e) and the curves for the data
with targets shown in Figs. 54(h, d, f), respectively. We have simply chosen m = 6 and
4
3
4r
3
r2
1
0
4
3
1
(a)
5 10 15 20 25 30
(b)
5 10 15 20 25 30
(c)
5 10 15 20 25 30
(d)
5 10 15 20 25 30
Oi(e '" (f)c
0 5 10 15 20 25 30 0 5 10 15 20 25 30
k k
Figure 54. Examples of the timedependent exponent A(k) vs. k curves for the sea clutter
data (a, c, e) without and (b, d, f) with the target. Six curves, from bottom to
top, correspond to shells (2(i+1)/2, 2i/2) With i 13, 14, 15, 16, 17, and 18.
The sampling time for the sea clutter data is 1 msee, and embedding
parameters are m 6, L 1 217 data points are used in the computation.
L = 1. Very similar curves have been obtained for other choices of m and L. We have
not observed a common envelope at any scales. In fact, the results of Fig. 54 are generic
among all the 280 sea clutter data on~ li. .1 here. Hence, we have to conclude that none of
the sea clutter data is chaotic.
5.3 Target Detection within Sea Clutter by Separating Scales
Robust detection of low observable targets within sea clutter is a very important issue
in remote sensing and radar signal processing applications, for a number of reasons: (i)
identifying objects within sea clutter such as submarine periscopes, lowflying aircraft,
and missiles can greatly improve our coastal and national security; (ii) identifying small
marine vessels, navigation buoys, small pieces of ice, patches of spilled oil, etc. can
significantly reduce the threat to the safety of ship navigation; (iii) monitoring and
policing of illegal fishing is an important activity in the environmental monitoring. Due
to the rough sea surface, which is a consequence of energy transfer from small scale to
large scale of sea surface wave, and the multipath propagation of the radar backscatter,
sea clutter is often highly nonGaussian, even spiky, especially in heavy sea conditions.
Hence, sea clutter modelling is a very difficult problem, and a lot of effort has been made
to study sea clutter, both through the analysis of the distributions of sea clutter, including
Weibull [112], lognormal [113], K( [114, 115], and compoundGaussian distributions [116],
as well as using chaos theory [100104, 107, 108, 110, 111], wavelets [117], neural
networks [118, 119], waveletneural net combined approaches [120, 121], and the concept of
fractal dimension [122] and fractal error [123, 124]. However, no simple method has been
found to robustly detect low observable objects within sea clutter [125].
In this section, we examine whether the Lyapunov exponent estimated by conventional
methods and the scaledependent exponent may be helpful for target detection within sea
clutter.
We first check whether the Lyapunov exponent A estimated by conventional methods
can be used for distinguishing sea clutter data with and without targets. As pointed
out earlier, the conventional way of estimating the Lyapunov exponent is to compute
A(k)/kbt, where A(k) is defined as in Eq. 24, subject to the conditions that XiXj < r
and Xi+k Xj,,  < R, where r and R are prescribed small distance scales. The condition
Xi Xj  < r amounts to our smallest shell in computing the A(k) curves. The condition
Xi+k Xj+k I < R is presumably to set the time scale k smaller than the prediction
time scale k. Our analysis finds that the Lyapunov exponent estimated this way fails to
detect targets from the sea clutter data on~ lli. here. This motivates us to try a modified
approach, as described below.
The modified approach for estimating the Lyapunov exponent A uses the leastsquares
fit to the first few samples of the A(k) curve where the A(k) curve increases linearly or
quasilinearly. This is done for all the 280 sea clutter time series. To better appreciate the
variation of the Lyapunov exponent A among the 14 range hins of the sea clutter data,
we have first subtracted the A parameter of each hin by the minimum of A values for that
measurement, and then normalized the obtained A values by its maximum. The variations
of the A parameters with the 14 range hins for the 10 HH measurements are shown in
Figs. 55(a)(j), respectively, where open circles denote the range hins with the target,
and asteroids denote the hins without the target. The primary target hin is explicitly
indicated hv an arrow. Similar results have been obtained for the 10 VV measurements.
While Figs. 55(a) and (h) indicate that the primary target hin can he separated from bins
without the target, in general, we have to conclude that the Lyapunov exponent estimated
by methods equivalent or similar to conventional means cannot he used for distinguishing
sea clutter data with and without targets.
Next let us examine whether the scaledependent exponent corresponding to large
scale may be useful for detecting targets within sea clutter data. By scaledependent
exponent, we mean that if we use the leastsquares fit to the linearly or quasilinearly
increasing part of the A(k) curves of different shells, the slopes of those lines depend
on which shells are used for computing the A(k) curves. In other words, if we plot the
exponent estimated this way against the size of the shell, then the exponent varies
with the shell size or the scale. The scaledependent exponent has been used in the
study of noiseinduced chaos in an optically injected semiconductor laser model to
examine how noise affects different scales of dynamic systems [26]. We now focus on the
exponent corresponding to large shell or scale. Figures 56(a)(j) show the variations of
the scaledependent exponent corresponding to large scale with the 14 range hins for the
Primary
JTarget Bin
&z 05
0 5 10 15
(d)
05
0 5 10 15
1x
(f)~*
0 5 10 15
(c) J
0 5 10 15
(e) +
0 5 10 15
5 10 15
00 5 10 15 00 5v 10 `15
Bln Number Bin Number
Figure 55. Variations of the Lyapunov exponent A estimated by conventional methods vs.
the 14 range bins for the 10 HHI measurements. Open circles denote the range
bins with target, while denote the bins without target. The primary target
bin is explicitly indicated by an arrow.
same 10 HH measurements as studied in Fig. 55. We observe that the primary target bin
can be easily separated from the range bins without the target, since the scaledependent
exponent for the primary target bin is much larger than those for the bins without the
target.
It is thus clear that the scaledependent exponent corresponding to large scale is
very useful for distinguishing sea clutter data with and without targets. This II__ ;
that sea clutter may contain interesting dynamic features and that the scaledependent
exponent may be an important parameter for target detection within sea clutter. This
(a) Primar (b)
Target Bin
0505 /q
0g 5 10 15 0 5 10 15
(c) (d)
4z 05 0
0 5 10 15 0 5 10 15
(e) ? (f) I
A 05 /05
0 5 10 15 0O 5 10 15
(g) / (h)
Az 05 05
0 0
0 5 10 15 0 5 10 15
4z 05 / 05
O' 0
0 5 10 15 0 5 10 15
Bln Number Bin Number
Figure 56. Variations of the scaledependent exponent corresponding to large scale vs. the
range bins for the 10 HHI measurements. Open circles denote the range bins
with target, while denote the bins without target. The primary target bin is
explicitly indicated by an arrow.
example clearly signifies the importance of incorporating the concept of scale in a measure.
In fact, the concept of scale is only incorporated in the time dependent exponent A(k)
curves [21, 22] in a static manner. In ('! .pter 6, we shall see that when a measure
dynamically incorporates the concept of scale, it will become much more powerful.
Note that one difficulty of using the scaledependent exponent for target detection
within sea clutter is that we may need to choose a suitable scale for estimating the
exponent. This makes the method not easy to use, and it is hard to make the method
automatic.
5.4 Target Detection within Sea Clutter by Powerlaw Sensitivity to Initial
Conditions (PSIC)
In this section, we show how the concept of PSIC can he applied to detect target
within sea clutter data. Examples of the A(k) vs. In k curves for the range hins without
target and those with target of one measurement are shown in Figs. 57(a) and (b),
respectively, where the curves denoted by asteroids are for data without the target,
while the curves denoted by open circles are for data with the target. We observe that
the curves are fairly linear for the first a few samples. Also notice that the slopes of
the curves for the data with the target are much larger than those for the data without
the target. For convenience, we denote the slope of the curve by the parameter /9. To
better appreciate the variation of the /9 parameter with the range hims, we normalize /9
of each hin by the nmaxinial /9 value of the 14 range hins within the single measurement.
Fig. 58 shows the variation of the /9 parameter with the 14 range hims. It is clear that
the /9 parameter can he used to distinguish sea clutter data with and without target.
Interestingly, the feature shown in Fig. 58 is generically true for all the measurements. It
is worth pointing out that usually the A(k) vs. In k curves for different scales are almost
parallel, thus the estimated /9 parameter is relatively less dependent on the scale. This is a
very nice feature of this method.
Let us examine if a robust detector for detecting targets within sea clutter can he
developed based on the /9 parameter. We have systematically studied 280 time series
of the sea clutter data measured under various sea and weather conditions. To better
appreciate the detection performance, we have first only focused on hins with primary
targets, but omitted those with secondary targets, since sometimes it is hard to determine
whether a hin with secondary target really hits a target or not. After omitting the range
hin data with secondary targets, the frequencies for the /9 parameter under the two
hypotheses (the hins without targets and those with primary targets) for HH datasets
0.5'
0 0.5 1 1.5 2 2.5 3 3.5
In k
In k
Figure 57.
Examples of the A(k) vs. In k curves for range bins (a) without and (b) with
target. Open circles denote the range bins with target, while denote the bins
without target.
Primary
Target Bin
0.9
co. 0.85
0.75
x
0.7'
0 2 4
8
Bin Number
10 12 14
Figure 58.
Variation of the parameter with the 14 range bins. Open circles denote the
range bins with target, while denote the bins without target.
are shown in Fig. 59. We observe that the frequencies completely separate for the HH
datasets. This means the detection accuracy can be 1011' .
no targets
12 primary targets
10
a,
LL
n mm
t0.4 0.5 0.6 0.7 0.8 0.9 1
Figure 59. Frequencies of the bins without
the HH datasets.
targets and the bins with primary targets for
CHAPTER 6
1\ULTISCALE ANALYSIS BY SCALEDEPENDENT LYAPUNOV EXPONENT
(SDLE)
In C'!. Ilter 1, we have emphasized the importance of developing scaledependent
measures to simultaneously characterize behaviors of complex multiscaled signals on a
wide range of scales. In this chapter, we shall develop an effective algorithm to compute
an excellent multiscale measure, the scaledependent Lyapunov exponent (SDLE), and
show that the SDLE can readily classify various types of complex motions, including truly
lowdimensional chaos, noisy chaos, noiseinduced chaos, processes defined by powerlaw
sensitivity to initial conditions (PSIC), and complex motions with chaotic behavior on
small scales but diffusive behavior on large scales. Finally, we shall discuss how the SDLE
can help detect hidden frequencies in large scale orderly motions.
6.1 Basic Theory
The scaledependent Lyapunov exponent (SDLE) is a variant of FSLE [126]. The
algorithm for calculating the FSLE is very similar to the Wolf et al's algorithm [127].
It computes the average rfold time by monitoring the divergence between a reference
trajectory and a perturbed trajectory. To do so, it needs to define 1.. I.rest neighbors",
as well as needs to perform, from time to time, a renormalization when the distance
between the reference and the perturbed trajectory becomes too large. Such a procedure
requires very long time series, and therefore, is not practical. To facilitate derivation of
a fast algorithm that works on short data, as well as to ease discussion of continuous but
nondifferentiable stochastic processes, we cast the definition of the SDLE as follows.
Consider an ensemble of trajectories. Denote the initial separation between two
nearby trajectories by co, and their averagee separation at time t and t + at by et and et as,
respectively. Being defined in an average sense, et and et as can he readily computed from
any processes, even if they are nondifferentiable. Next we examine the relation between et
and et+nt, where at is small. When at 0 we have,
et+nt = etex(tt)nt (61)
where A(et) is the SDLE. It is given by
In et+nt In et
A(et) = (62)
Given a time series data, the smallest at possible is the sampling time 7r.
The definition of the SDLE so_~ ; Ma simple ensemble average based scheme to
compute it. A straightforward way would be to find all the pairs of vectors in the phase
space with their distance approximately e, and then calculate their average distance after a
time at. The first half of this description amounts to introducing a shell (indexed as k),
Ek  Vs Vy  < k Ek(63)
where 1%, Vj are reconstructed vectors, Ek, (the radius of the shell) and A~rk (the width
of the shell) are arbitrarily chosen small distances. Such a shell may be considered as
a differential element that would facilitate computation of conditional probability. To
expedite computation, it is advantageous to introduce a sequence of shells, k = 1, 2, 3, 
Note that this computational procedure is similar to that for computing the socalled
timedependent exponent (TDE) curves [21, 22, 37].
With all these shells, we can then monitor the evolution of all of the pairs of vectors
(1M, Vj) within a shell and take average. When each shell is very thin, by assuming that
the order of averaging and taking logarithm in Eq. 62 can be interchanged, we have
where t and at are integers in unit of the sampling time, and the angle brackets denote
average within a shell. Note that contributions to the SDLE at a specific scale from
different shells can be combined, with the weight for each shell being determined by the
number of the pairs of vectors (1M, Vyj) in that shell. In the following, to see better how
each shell characterizes the dynamics of the data on different scales, we shall plot the A(e)
curves for different shells separately.
In the above formulation, it is implicitly assumed that the initial separation, % 1 Vy l,
aligns with the most unstable direction instantly. For highdimensional systems, this is
not true, especially when the growth rate is nonuniform and/or the eigenvectors of the
Jacobian are nonnormal. Fortunately, the problem is not as serious as one might be
concerned, since our shells are not infinitesimal. When computing the TDE curves [21,
22, 37], we have found that when difficulties arise, it is often sufficient to introduce an
additional condition,
I i> m 1 L65
when finding pairs of vectors within each shell. Such a scheme also works well when
computing the SDLE. This means that, after taking a time comparable to the embedding
window (m 1)L, it would be safe to assume that the initial separation has evolved to the
most unstable direction of the motion.
Before proceeding on, we wish to emphasize the 1!! ri ~ difference between our
algorithm and the standard method for calculating the FSLE. As we have pointed out,
to compute the FSLE, two trajectories, one as reference, another as perturbed, have to
be defined. This requires huge amounts of data. Our algorithm avoids this by employing
two critical operations to fully utilize information about the time evolution of the data: (i)
The reference and the perturbed trajectories are replaced by time evolution of all pairs of
vectors satisfying the inequality (65) and falling within a shell, and (ii) introduction of
a sequence of shells ensures that the number of pairs of vectors within the shells is large
while the ensemble average within each shell is well defined. Let the number of points
needed to compute the FSLE by standard methods be NV. These two operations imply
that the method described here only needs about z/Vpoints to compute the SDLE. In the
following, we shall illustrate the effectiveness of our algorithm by examining various types
of complex motions.
6.2 Classification of Complex Motions
T> understand the SDLE as well as appreciate its power, we apply it to classify
various types of complex motions.
6.2.1 Chaos, Noisy Chaos, and Noiseinduced Chaos
Obviously, for truly lowdimensional chaos, A(e) equals the largest positive Lyapunov
exponent, and hence, must he independent of e over a wide range of scales. For noisy
chaos, we expect A(e) to depend on small e. T> illustrate both features, we consider the
chaotic Lorenz system with stochastic forcing described by Eq. 25. Figure 61(a) shows
five curves, for the cases of D = 0, 1, 2, 3, 4. The computations are done with 10000 points
and m = 4, L = 2. We observe a few interesting features:
For the clean chaotic signal, A(e) slightly fluctuates around a constant (which
numerically equals the largest positive Lyapunov exponent) when e is smaller than a
threshold value which is determined by the size of the chaotic attractor. The reason
for the small fluctuations in A(e) is that the divergence rate varies from one region of
the attractor to another.
When there is stochastic forcing, A(e) is no longer a constant when e is small, but
increases as y In e when the scale e is decreased. The coefficient ]* does not seem to
depend on the strength of the noise. This feature II__ is that entropy generation is
infinite when the scale e approaches to zero. Note that the relation of A(e) ~ y In e
has also been observed for the FSLE and the eentropy. In fact, such a relation can
he readily proven for the eentropy.
When the noise is increased, the part of the curve with A(e) ~ y 1n e shifts to the
right. In fact, little chaotic signature can he identified when D is increased beyond
:3. When noise is not too 1i~~e.: this feature can he readily used to quantify the
strength of noise.
Next we consider noiseinduced chaos. To illustrate the idea, we follow [2:3] and study
the noisy logistic map
r,w+1 = p~r, (1 r,w) + P,, O < r,z < 1 (66)
2.5 0.7
(a) LorenZ0. (b) Logistic
2 05 (1): slope =0.28
S0.4 h (2): slope = 0
D. = o0.
D = 3 0.1
~D=4
0.5 v 0
1.5 1 0.5 2 1.5 1 0.5 0
log ,E log ,E
Figure 61. Scaledependent Lyapunov exponent A(e) curves for (a) the clean and the
noisy Lorenz system, and (b) the noiseinduced chaos in the logistic map.
Curves from different shells are designated by different symbols.
where p is the bifurcation parameter, and P, is a Gaussian random variable with
zero mean and standard deviation o. In [23], we reported that at p = 3.74 and
o = 0.002, noiseinduced chaos occurs, and thought that it may be difficult to distinguish
noiseinduced chaos from clean chaos. In Fig. 61(b), we have plotted the A(et) for this
particular noiseinduced chaos. The computation was done with m = 4, L = 1 and 10000
points. We observe that Fig. 61(b) is very similar to the curves of noisy chaos plotted in
Fig. 61(a). Hence, noiseinduced chaos is similar to noisy chaos, but different from clean
chaos .
At this point, it is worth making two comments: (i) On very small scales, the effect
of measurement noise is similar to that of dynamic noise. (ii) The A(e) curves shown in
Fig. 61 are based on a fairly small shell. The curves computed based on larger shells
collapse on the right part of the curves shown in Fig. 61. Because of this, for chaotic
systems, one or a few small shells would be sufficient. If one wishes to know the behavior
of A on ever smaller scales, one has to use longer and longer time series.
1x 10
35
5
10 100
Figure 62. Scaledependent Lyapunov exponent A(e) curve for the MackeyGlass system.
The computation was done with ni = 5, L = 1, and 5000 points sampled with a
time interval of 6.
Finally, we consider the MackeyGlass delay differential system [39], defined by
Eq. 26. The system has two positive Lyapunov exponents, with the largest Lyapunov
exponent close to 0.007 [21]. Having two positive Lyapunov exponents while the value of
the largest Lyapunov exponent of the system is not much greater than 0, one might he
concerned that it may be difficult to compute the SDLE of the system. This is not the
case. In fact, this system can he analyzed as straightforwardly as other dynamical systems
including the Henon map and the Rossler system. An example of the A(e) curve is shown
in Fig. 62, where we have followed [21] and used ni = 5, L = 1, and 5000 points sampled
with a time interval of 6. Clearly, we observe a well defined plateau, with its value close to
0.007. This example illustrates that when computing the SDLE, one does not need to be
very concerned about nonuniform growth rate in highdimensional systems.
Up till now, we have focused on the positive portion of A(es). It turns out that when t
is large, A(e) becomes oscillatory, with mean about 0. Denote the corresponding scales by
ea, and call them the limiting or characteristic scales. They are the stationary portion of
et, and hence, they may still be a function of time. We have found that the limiting scales
capture the structured component of the data. This feature will be further illustrated
later, when we discuss identification of hidden frequencies. Therefore,the positive portion
of A(et) and the concept of limiting scale provide a comprehensive characterization of the
signals .
6.2.2 Processes Defined by Powerlaw Sensitivity to Initial Conditions (PSIC)
In C'!s Ilter 3, we have seen that powerlaw sensitivity to initial conditions (PSIC)
provides a common foundation for chaos theory and random fractal theory. Here, we
examine whether the SDLE is able to characterize processes defined by PSIC. Pleasingly,
this is indeed so. Below, we derive a simple equation relating the A, and q of PSIC to the
SDLE.
First, we recall that PSIC is defined by
( = lim = 1 ( ) )
Since et = av(0)s,, it is now more convenient to express the SDLE as a function of (c.
Using Eq. 62, we find that
(()=In $t+nt In (a(67
When at 0, 1+ (1 q)Agt > (1 q)XAft. Simplifying Eq. 67, we obtain
A((s) = Aq tq1 (68)
We now consider three cases:
*For chaotic motions, q = 1, therefore,
A((s) = Xq = COnStant
F'or 1/ fB noise, we, hae = 1_ (Eq. :328) and A," = H (Eq. :329). TIherefor~el
A((s) = Hg, H
For Le~vy fightsl weit have q = 1 a? (Eq. :338) and A,"l = 1/0~ (Eq. :339). Therefore,
6.2.3 Complex Motions with Multiple Scaling Behaviors
Some dynamical systems may exhibit multiple scaling behaviors, such as chaotic
behavior on small scales but diffusive behavior on large scales. To see how the SDLE can
characterize such systems, we follow Cencini et al. [128] and study the following map,
Xrn+l = [Xrn] + F(xrn [Xrnl) + aife (69)
where [r,z] denotes the integer part of r,t, r7< is a noise uniformly distributed in the interval
[1, 1], o is a parameter quantifying the strength of noise, and Fly) is given by
Fly) = (2 + A) y if y E [0, 1/2) (610)
(2 + A) y (1 + a) if y E (1/2, 1]
The map Fly) is shown in Fig. 63 as the dashed lines. It gives a chaotic dynamics with
a positive Lyapunov exponent In(2 + a). On the other hand, the term [:re] introduces a
random walk on integer grids.
It turns out this system is very easy to analyze. When a = 0.4, with only 5000 points
and m = 2, L = 1, we can resolve both the chaotic behavior on very small scales, and the
normal diffusive behavior (with slope 2) on large scales. See Fig. 64(a).
We now ask a question: Given a small dataset, which type of behavior, the chaotic
or the diffusive, is resolved first? To answer this, we have tried to compute the SDLE
with only 500 points. The result is shown in Fig. 64(b). It is interesting to observe that
the chaotic behavior can he well resolved by only a few hundred points. However, the
(30.6
S0.4
LL.
1
T
0.2'
0
Figure 63.
The function F(x) (Eq. 610) for a = 0.4 is shown as the dashed lines. The
function G(x) (Eq. 611) is an approximation of F(x), obtained using 40
intervals of slope 0. In the case of noiseinduced chaos discussed in the paper,
G(x) is obtained from F(x) using 104 HinerValS Of Slope 0.9.
100
102
100
102
clean data
noisy data
104 103 102 10 1
104 103 102
Scaledependent Lyapunov exponent A(e) for
(a) 5000 points were used, for the noisy case,
used.
Figure 64.
the model described by Eq. 69.
o = 0.001. (b) 500 points were
diffusive behavior needs much more data to resolve. Intuitively, this makes sense, since
the diffusive behavior amounts to a Brownian motion on the integer grids and is of much
higher dimension than the small scale chaotic behavior. Therefore, more data are needed
to resolve it.
We have also studied the noisy map. The resulting SDLE for a = 0.001 is shown in
Fig. 64(a), as squares. We have again used 5000 points. While the behavior of the SDLE
II _ r ;noisy dynamics, with 5000 points, we are not able to well resolve the relation
A(e) ~ In e. This indicates that for the noisy map, on very small scales, the dimension is
very high.
Map (69) can be modified to give rise to an interesting system with noiseinduced
chaos. This can be done by replacing the function Fly) in map (69) by G(y) to obtain
the following map [128],
xt+l = [Xt] + G(xt [xt]) + arlt (611)
where Orl is a noise uniformly distributed in the interval [1, 1], a is a parameter
quantifying the strength of noise, and G(y) is a piecewise linear function which approximates
Fly) of Eq. 610. An example of G(y) is shown in Fig. 63. In our numerical simulations,
we have followed Cencini et al. [128] and used 104 interValS Of Slope 0.9 to obtain G(y).
With such a choice of G(y), in the absence of noise, the time evolution described by
the map (611) is nonchaotic, since the largest Lyapunov exponent In(0.9) is negative.
With appropriate noise level (e.g., a = 104 or 103), the SDLE for the system becomes
indistinguishable to the noisy SDLE shown in Fig. 64 for the map (69). Having a
diffusive regime on large scales, this is a more complicated noiseinduced chaos than the
one we have found from the logistic map.
Before proceeding on, we make a comment on the computation of the SDLE from
deterministic systems with negative largest Lyapunov exponents, such as the map (611)
without noise. A transientfree time series from such systems is a constant time series.
Therefore, there is no need to compute the SDLE or other metrics. When the time series
contains transients, if the time series is sampled with high enough sampling frequency,
then the SDLE is negative. In the case of simple exponential decay to a fixed point, such
as expressible as ext, where A > 0, one can readily prove that the SDLE is A. Since such
systems are not complex, we shall not be further concerned about them.
6.3 Characterizing Hidden Frequencies
Defining and characterizing large scale orderly motions is a significant issue in many
disciplines of science. One of the most important types of large scale orderly motions is
the oscillatory motions. An interesting type of oscillatory motions is associated with the
socalled hidden frequency phenomenon. That is, when the dynamics of a complicated
system is monitored through the temporal evolution of a variable x, Fourier analysis of
x(t) may not II 1 any oscillatory motions. However, if the dynamics of the system is
monitored through the evolution of another variable, 11, z, then the Fourier transform of
z(t) may contain a welldefined spectral peak, indicating oscillatory motions. An example
is the chaotic Lorenz system described by Eq. 25. In Figs. 65 (ac), we have shown the
powerspectral density (PSD) of the x, y, z components of the system. We observe that the
PSD of x(t) and y(t) are simply broad. However, the PSD of z(t) shows up a very sharp
spectral peak. Recall that geometrically, the Lorenz attractor consists of two scrolls (see
Fig. 66). The sharp spectral peak in the PSD of z(t) of the Lorenz system is due to the
circular motions along either of the scrolls.
The above example illustrates that if the dynamics of a system contains a hidden
frequency that cannot be revealed by the Fourier transform of a measured variable
( 11, x(t)), then in order to reveal the hidden frequency, one has to embed x(t) to a
suitable phase space. This idea has led to the development of two interesting methods
for identifying hidden frequencies. One method is proposed by Ortega [129, 130], by
computing the temporal evolution of density measures in the reconstructed phase space.
Another is proposed by C'I. I i. et al. [131], by taking singular value decomposition of local
neighbors. The Ortega's method has been applied to an experimental time series recorded
10 "
104
102 104~
0 0.1 0.2 0.3 0.4 0.5 0 0.1 0.2 0.3 0.4 0.5
10 Frequency f 14Frequency f
1(b) (e)
106
,8100
4
104
1021 104
0 0.1 0.2 0.3 0.4 0.5 0 0.1 0.2 0.3 0.4 0.5
109Frequency f 14Frequency f
106 100
1031 1 104
0 0.1 0.2 0.3 0.4 0.5 0 0.1 0.2 0.3 0.4 0.5
Frequency f Frequency f
F'igur~e 65. P~owerspectral density (PSD>) for (a) x(t), (b) y(t), (c) z(t), (d) ei)(t), (e)
from a farinfrared laser in a chaotic state. The laser dataset can be downloaded from
the link http://wwwpsy ch. stanf ord .edu/~ andreas/Time Ser ies /Sant aFe .htm1. It
contains 10000 points, sampled with a time interval of 80 ns. Fig. 67(a) shows the laser
dataset. The PSD of the data is shown in Fig. 67(b), where one observes a sharp peak
around 1.7 MHz. Fig. 67(c) shows the PSD of the density time series, where one notes
an additional spectral peak around 37 kHz. This peak is due to the envelope of chaotic
pulsations, which is discernable from Fig. 67(a).

Full Text 
PAGE 1
1
PAGE 2
2
PAGE 3
3
PAGE 4
Firstofall,Iwouldliketothankmyadvisor,ProfessorJianboGao,forhisinvaluableguidanceandsupport.Iwasluckytoworkwithsomeonewithsomanyoriginalideasandsuchasharpmind.Withouthishelp,myPh.D.degreeandthisdissertationwouldhavebeenimpossible.DuringmyyearsattheUniversityofFlorida,Iwasmovedbyhismagnicentpersonality.AfamousChineseadagesays,\Onedayasmyadvisor,entirelifeasmyfather."ProfessorGaohighlyqualies.Ialsowanttothankmycommitteemembers(ProfessorsSuShingChen,YuguangFang,JoseC.Principe,andKeithD.White),fortheirinterestinmyworkandvaluablesuggestions.Mythanksofcoursealsogotoallthefaculty,staandmyfellowstudentsintheDepartmentofElectricalandComputerEngineering.IextendparticularthankstomyocematesJaeminLee,YiZheng,andUngsikKim,fortheirhelpfuldiscussions.Iexpressmyappreciationtoallofmyfriendswhooeredencouragement.IspeciallythankQianZhan,XingyanFan,LingyuGuandhiswife(JingZhang),YanminXu,YanLiu,andYanaLiang.TheygavemealotofhelpduringmyearlyyearsintheUniversityofFlorida.Theytaughtmetocookandmademylifemorecolorful.Finally,Ioweagreatdebtofthankstomyparents,grandparents,andmyboyfriend,JingAi,fortheirpillarlikesupportandtheirunwaveringbeliefandcondenceinmewithoutthisIdonotthinkIwouldhavemadeitthisfar! 4
PAGE 5
page ACKNOWLEDGMENTS ................................. 4 LISTOFFIGURES .................................... 7 ABSTRACT ........................................ 9 CHAPTER 1INTRODUCTION .................................. 10 1.1OverviewofDynamicalSystems ........................ 11 1.2ExamplesofMultiscalePhenomena ...................... 15 1.3BriefIntroductiontoChaos .......................... 20 1.4BriefIntroductiontoFractalGeometry .................... 21 1.5ImportanceofConnectingChaosandRandomFractalTheories ...... 23 1.6ImportanceoftheConceptofScale ...................... 24 1.7StructureofthisDissertation ......................... 24 2GENERALIZATIONOFCHAOS:POWERLAWSENSITIVITYTOINITIALCONDITIONS(PSIC) ................................ 27 2.1DynamicalTestforChaos ........................... 27 2.2GeneralComputationalFrameworkforPowerlawSensitivitytoInitialConditions(PSIC) ............................... 33 2.3CharacterizingEdgeofChaosbyPowerlawSensitivityofInitialConditions(PSIC) ...................................... 36 3CHARACTERIZINGRANDOMFRACTALSBYPOWERLAWSENSITIVITYTOINITIALCONDITIONS(PSIC) ........................ 41 3.1Characterizing1=fProcessesbyPowerlawSensitivitytoInitialConditions(PSIC) ...................................... 41 3.1.1CharacterizingFractionalBrownianMotion(fBm)ProcessesbyPowerlawSensitivitytoInitialConditions(PSIC) .......... 43 3.1.2CharacterizingON/OFFProcesseswithParetodistributedONandOFFPeriodsbyPowerlawSensitivitytoInitialConditions(PSIC) 50 3.2CharacterizingLevyProcessesbyPowerlawSensitivitytoInitialConditions(PSIC) ...................................... 54 4STUDYOFINTERNETDYNAMICSBYPOWERLAWSENSITIVITYTOINITIALCONDITIONS(PSIC) ........................... 60 4.1StudyofTransportDynamicsbyNetworkSimulation ............ 61 4.2ComplicatedDynamicsofInternetTransportProtocols ........... 69 5
PAGE 6
............... 76 5.1SeaClutterData ................................ 78 5.2NonChaoticBehaviorofSeaClutter ..................... 79 5.3TargetDetectionwithinSeaClutterbySeparatingScales .......... 81 5.4TargetDetectionwithinSeaClutterbyPowerlawSensitivitytoInitialConditions(PSIC) ............................... 86 6MULTISCALEANALYSISBYSCALEDEPENDENTLYAPUNOVEXPONENT(SDLE) ........................................ 89 6.1BasicTheory .................................. 89 6.2ClassicationofComplexMotions ....................... 92 6.2.1Chaos,NoisyChaos,andNoiseinducedChaos ............ 92 6.2.2ProcessesDenedbyPowerlawSensitivitytoInitialConditions(PSIC) .................................. 95 6.2.3ComplexMotionswithMultipleScalingBehaviors .......... 96 6.3CharacterizingHiddenFrequencies ...................... 99 7CONCLUSIONS ................................... 104 REFERENCES ....................................... 105 BIOGRAPHICALSKETCH ................................ 115 6
PAGE 7
Figure page 11Exampleofatrajectoryinthephasespace ..................... 14 12Giantoceanwave(tsunami). ............................. 17 13Exampleofseaclutterdata. ............................. 18 14ExampleofHeartratevariabilitydataforanormalsubject. ........... 19 21Timedependentexponent(k)vs.evolutiontimekcurvesforLorenzsystem. 32 22Timedependentexponent(k)vs.evolutiontimekcurvesfortheMackeyGlasssystem. ........................................ 33 23Timedependentexponent(t)vs.tcurvesfortimeseriesgeneratedfromlogisticmap. .......................................... 37 24Timedependentexponent(t)vs.tcurvesfortimeseriesgeneratedfromHenonmap. .......................................... 40 31WhitenoiseandBrownianmotion ......................... 46 32SeveralfBmprocesseswithdierentH. ....................... 48 33Timedependentexponent(k)vs.lnkcurvesforfractionalBrownianmotion. 51 34ExampleofON/OFFprocesses. ........................... 51 35Timedependentexponent(k)vs.lnkcurvesforON/OFFprocesseswithParetodistributedONandOFFperiods. .......................... 53 36ExamplesforBrownianmotionsandLevymotions. ................ 57 37Timedependentexponent(k)vs.lnkcurvesforLevyprocesses. ........ 58 41ExampleofquasiperiodiccongestionwindowsizeW(i)timeseries. ....... 63 42ExamplesofchaoticT(i)timeseries. ........................ 65 43Complementarycumulativedistributionfunctions(CCDFs)forthetwochaoticT(i)timeseriesofFigs. 42 (a,c). .......................... 66 44ThetimeseriesT(i)extractedfromaW(i)timeseriescollectedusingnet100instrumentsoverORNLLSUconnection. ...................... 68 45TimeseriesforthecongestionwindowsizecwndforORNLGaTechconnection. 72 46Timedependentexponent(k)vs.kcurvesforcwnddatacorrespondingtoFig. 45 ........................................ 73 7
PAGE 8
45 ........................................ 75 51Collectionofseaclutterdata ............................ 79 52Typicalseaclutteramplitudedata. ......................... 79 53Twoshortsegmentsoftheamplitudeseaclutterdataseverelyaectedbyclipping. 80 54Examplesofthetimedependentexponent(k)vs.kcurvesfortheseaclutterdata. .......................................... 81 55VariationsoftheLyapunovexponentestimatedbyconventionalmethodsvs.the14rangebinsforthe10HHmeasurements. .................. 84 56Variationsofthescaledependentexponentcorrespondingtolargescalevs.therangebinsforthe10HHmeasurements. ...................... 85 57Examplesofthe(k)vs.lnkcurvesfortheseaclutterdata. ........... 87 58Variationoftheparameterwiththe14rangebins. ............... 87 59FrequenciesofthebinswithouttargetsandthebinswithprimarytargetsfortheHHdatasets. ................................... 88 61ScaledependentLyapunovexponent()curvesfortheLorenzsystemandlogisticmap. .......................................... 93 62ScaledependentLyapunovexponent()curvefortheMackeyGlasssystem. .. 94 63ThefunctionsF(x)andG(x). ............................ 97 64ScaledependentLyapunovexponent()forthemodeldescribedbyEq. 6{9 . 97 65Powerspectraldensity(PSD)fortimeseriesgeneratedfromLorenzsystem. .. 100 66Lorenzattractor. ................................... 101 67Hiddenfrequencyphenomenonoflaserintensitydata. ............... 102 8
PAGE 9
Complexsystemsusuallycomprisemultiplesubsystemswithhighlynonlineardeterministicandstochasticcharacteristics,andareregulatedhierarchically.Thesesystemsgeneratesignalswithcomplexcharacteristicssuchassensitivedependenceonsmalldisturbances,longmemory,extremevariations,andnonstationarity.Chaostheoryandrandomfractaltheoryaretwoofthemostimportanttheoriesdevelopedforanalyzingcomplexsignals.However,theyhaveentirelydierentfoundations,onebeingdeterministicandtheotherbeingrandom.Tosynergisticallyusethesetwotheories,wedevelopedanewtheoreticalframework,powerlawsensitivitytoinitialconditions(PSIC),toencompassbothchaosandrandomfractaltheoriesasspecialcases.Toshowthepowerofthisframework,weappliedittostudytwochallengingandimportantproblems:Internetdynamicsandseaclutterradarreturns.WeshowedthatPSICcanreadilycharacterizethecomplicatedInternetdynamicsduetointerplayofnonlinearAIMD(additiveincreaseandmultiplicativedecrease)operationofTCPandstochasticuserbehavior,androbustlydetectlowobservabletargetsfromseaclutterradarreturnswithhighaccuracy. Wealsodevelopedanewmultiscalecomplexitymeasure,scaledependentLyapunovexponent(SDLE),thatcanbecomputedfromshortnoisydata.Themeasurereadilyclassiedvarioustypesofcomplexmotions,andsimultaneouslycharacterizedbehaviorsofcomplexsignalsonawiderangeofscales,includingcomplexirregularbehaviorsonsmallscales,andorderlybehaviors,suchasoscillatorymotions,onlargescales. 9
PAGE 10
Adynamicalsystemisonethatevolvesintime.Dynamicalsystemscanbestochastic(inwhichcasetheyevolveaccordingtosomerandomprocesssuchasthetossofacoin)ordeterministic(inwhichcasethefutureisuniquelydeterminedbythepastaccordingtosomeruleormathematicalformula).Whetherthesysteminquestionsettlesdowntoequilibrium,keepsrepeatingincycles,ordoessomethingmorecomplicated,dynamicsarewhatweusetoanalyzethebehaviorofsystems. Complexdynamicalsystemsusuallycomprisemultiplesubsystemswithhighlynonlineardeterministicandstochasticcharacteristics,andareregulatedhierarchically.Thesesystemsgeneratesignalswithcomplexcharacteristicssuchassensitivedependenceonsmalldisturbances,longmemory,extremevariations,andnonstationarity[ 1 ].Astockmarket,forexample,isstronglyinuencedbymultilayereddecisionsmadebymarketmakers,aswellasinuencedbyinteractionsofheterogeneoustraders,includingintradaytraders,shortperiodtraders,andlongperiodtraders,andthusgivesrisetohighlyirregularstockprices.TheInternet,asanotherexample,hasbeendesignedinafundamentallydecentralizedfashionandconsistsofacomplexwebofserversandroutersthatcannotbeeectivelycontrolledoranalyzedbytraditionaltoolsofqueuingtheoryorcontroltheory,andgiverisetohighlyburstyandmultiscaletracwithextremelyhighvariance,aswellascomplexdynamicswithbothdeterministicandstochasticcomponents.Similarly,biologicalsystems,beingheterogeneous,massivelydistributed,andhighlycomplicated,oftengeneratenonstationaryandmultiscalesignals.Withtherapidaccumulationofcomplexdatainlifesciences,systemsbiology,nanosciences,informationsystems,andphysicalsciences,ithasbecomeincreasinglyimportanttobeabletoanalyzemultiscaleandnonstationarydata. Multiscalesignalsbehavedierentlydependinguponwhichscalethedataarelookedat.Inordertofullycharacterizesuchcomplexsignals,theconceptofscalehastobe 10
PAGE 11
Thisdissertationaimstobuildaneectivearsenalformultiscalesignalprocessingbysynergisticallyintegratingapproachesbasedonchaostheoryandrandomfractaltheory,andgoingbeyond,tocomplementconventionalapproachessuchasspectralanalysisandmachinelearningtechniques.Tomakesuchanintegrationpossible,twoimportanteortshavebeenmade:Powerlawsensitivitytoinitialconditions(PSIC):Wedevelopedanewtheoreticalframeworkforsignalprocessing,tocreateacommonfoundationforchaostheoryandrandomfractaltheory,sothattheycanbebetterintegrated.ScaledependentLyapunovexponent(SDLE):Wedevelopedanexcellentmultiscalemeasure,whichisavariantofthenitesizeLyapunovexponent(FSLE).Weproposedahighlyecientalgorithmforcalculatingit,andshowedthatitcanreadilyclassifydierenttypesofmotions,includingtrulylowdimensionalchaos,noisychaos,noiseinducedchaos,random1=fandstableLevyprocesses,andcomplexmotionswithchaoticbehavioronsmallscalesbutdiusivebehavioronlargescales.Themeasurecanaptlycharacterizecomplexbehaviorsofrealworldmultiscalesignalsonawiderangeofscales. Intherestofthischapter,weshallrstpresentanoverviewofdynamicssystems,especiallyweshallpointoutwhynonlinearsystemsaremuchhardertoanalyzethanlinearones.Wethengiveafewexamplesofmultiscalephenomenatoshowthedicultiesandexcitementofmultiscalesignalprocessing.Afterthat,weshallbrieyintroducesomebackgroundknowledgeaboutchaosandfractalgeometry.Thenwediscusstheimportanceofconnectingchaostheoryandrandomfractaltheoryforcharacterizingthebehaviorsofmultiscalesignalsonawiderangeofscales.Wealsoemphasizetheimportanceofexplicitlyincorporatingtheconceptofscaleindevisingmeasuresforcharacterizingmultiscalesignals.Finally,weoutlinethescopeofthepresentstudy. 11
PAGE 12
dt2+bdx dt+kx=0(1{1) isanordinarydierentialequation,becauseitinvolvesonlyordinaryderivativesdx=dtandd2x=dt2.Thatis,thereisonlyoneindependentvariable,thetimet.Incontrast,theheatequation@u @t=@2u @x2 Averygeneralframeworkforordinarydierentialequationsisprovidedbythesystem _x1=f1(x1;;xn)..._xn=fn(x1;;xn):(1{2) Heretheoverdotsdenotedierentiationwithrespecttot.Thus_xi=dxi=dt.Thevariablesx1;;xnmightrepresentconcentrationsofchemicalsinareactor,populationsofdierentspeciesinanecosystem,orthepositionsandvelocitiesoftheplanetsinthesolarsystem.Thefunctionsf1;;fnaredeterminedbytheproblemathand. Forexample,thedampedoscillator(Eq. 1{1 )canberewrittenintheformof(Eq. 1{2 ),thankstothefollowingtrick:weintroducenewvariablesx1=xandx2=_x. 12
PAGE 13
1{2 )is_x1=x2_x2=b mx2k mx1: Forexample,theswingingofapendulumisgovernedbytheequationx+g Lsinx=0; Lsinx1: Itturnsoutthatthependulumequationcanbesolvedanalytically,intermsofellipticfunctions.Butthereoughttobeaneasierway.Afterall,themotionofthependulumissimple:atlowenergy,itswingsbackandforth,andathighenergyitwhirlsoverthetop.Thereshouldbesomewayofextractingthisinformationfromthesystemdirectlyusinggeometricmethods. Hereistheroughidea.Supposewehappentoknowasolutiontothependulumsystem,foraparticularinitialcondition.Thissolutionwouldbeapairoffunctionsx1(t)andx2(t),representingthepositionandvelocityofthependulum.Ifweconstructan 13
PAGE 14
Exampleofatrajectoryinthephasespace abstractspacewithcoordinates(x1;x2),thenthesolution(x1(t);x2(t))correspondstoapointmovingalongacurveinthisspace,asshownFig. 11 .Thiscurveiscalledatrajectory,andthespaceiscalledthephasespaceforthesystem.Thephasespaceiscompletelylledwithtrajectories,sinceeachpointcanserveasaninitialcondition. Ourgoalistorunthisconstructioninreverse:giventhesystem,wewanttodrawthetrajectories,andtherebyextractinformationaboutthesolutions.Inmanycases,geometricreasoningwillallowustodrawthetrajectorieswithoutactuallysolvingthesystem! Thephasespaceforthegeneralsystem(Eq. 1{2 )isthespacewithcoordinatesx1;;xn.Becausethisspaceisndimensional,wewillrefertoEq. 1{2 asanndimensionalsystemorannthordersystem.Thusnrepresentsthedimensionofthephasespace. Aswehavementionedearlier,mostnonlinearsystemsareimpossibletosolveanalytically.Whyarenonlinearsystemssomuchhardertoanalyzethanlinearones?Theessentialdierenceisthatlinearsystemscanbebrokendownintoparts.Theneachpartcanbesolvedseparatelyandnallyrecombinedtogettheanswer.Thisideaallowsa 14
PAGE 15
Butmanythingsinnaturedonotactthisway.Wheneverpartsofasysteminterfere,orcooperate,orcompete,therearenonlinearinteractionsgoingon.Mostofeverydaylifeisnonlinear,andtheprincipleofsuperpositionfailsspectacularly.Ifyoulistentoyourtwofavoritesongsatthesametime,youwillnotgetdoublethepleasure!Withintherealmofphysics,nonlinearityisvitaltotheoperationofalaser,theformationofturbulenceinauid,andthesuperconductivityofJosephsonjunctions. 15
PAGE 16
Also,ithasbeenobservedthatthefailureofasingleroutermaytriggerroutinginstability,whichmaybesevereenoughtoinstigatearouteapstorm.Furthermore,packetsmaybedeliveredoutoforderorevengetdropped,andpacketreorderingisnotapathologicalnetworkbehavior.AsthenextgenerationInternetapplicationssuchasremoteinstrumentcontrolandcomputationalsteeringarebeingdeveloped,anotherfacetofcomplexmultiscalebehaviorisbeginningtosurfaceintermsoftransportdynamics.Thenetworkingrequirementsforthesenextgenerationapplicationsbelongto(atleast)twobroadclassesinvolvingvastlydisparatetimescales:Highbandwidths,typicallymultiplesof10Gbps,tosupportbulkdatatransfers,Stablebandwidths,typicallyatmuchlowerbandwidthssuchas10to100Mbps,tosupportinteractive,steeringandcontroloperations. 12 .Tobequantitative,inFig. 13 ,two0.1sdurationseacluttersignals,sampledwithafrequencyof1KHz,areplottedinFig. 13 (a,b),a2sdurationsignalisplottedinFig. 13 (c),andanevenlongersignal(about130s)isplottedinFig. 13 (d).Itisclearthatthesignalisnotpurelyrandom,sincethewaveformcanbefairlysmoothonshorttimescales(Fig. 13 (a)).However,thesignalishighlynonstationary,sincethefrequencyofthesignal(Fig. 13 (a,b))aswellastherandomnessofthesignal(Fig. 13 (c,d))changeovertimedrastically.OnethuscanperceivethatnaiveFourieranalysisordeterministicchaoticanalysisofseacluttermay 16
PAGE 17
Giantoceanwave(tsunami).Supposeoureldofobservationincludesthewavetipoflengthscaleofafewmeters,itisthenclearthatthecomplexityofseaclutterismainlyduetomassivereectionofradarpulsesfromawavyandeventurbulentoceansurface. notbeveryuseful.FromFig. 13 (e),whereX(m)tisthenonoverlappingrunningmeanofXoverblocksizem,andXistheseaclutteramplitudedata,itcanbefurtherconcludedthatneitherautoregressive(AR)models(Asanexample,wegivethedenitionfortherstorderARmodel,whichisgivenbyxn+1=axn+n,wheretheconstantcoecientasatises06=jaj<1andnisawhitenoisewithzeromean.)nortextbookfractalmodelscandescribethedata.ThisisbecauseARmodelingrequiresexponentiallydecayingautocorrelation(whichamountstoVar(X(m)t)m1,oraHurstparameterof1/2.SeelaterchaptersforthedenitionoftheHurstparameter),whilefractalmodelingrequiresthevariationbetweenVar(X(m)t)andmtofollowapowerlaw.However,neitherbehaviorisobservedinFig. 13 (e).Indeed,albeitextensiveworkhasbeendoneonseaclutter,thenatureofseaclutterisstillpoorlyunderstood.Asaresult,theimportantproblemof 17
PAGE 18
Exampleofseaclutterdata.(a,b)two0.1sdurationsignal;(c)a2sdurationseacluttersignal;(d)theentireseacluttersignal(ofabout130s);and(e)log2[m2Var(X(m))]vs.log2m,whereX(m)=fX(m)t:X(m)t=(Xtmm+1++Xtm)=m;t=1;2;gisthenonoverlappingrunningmeanofX=fXt:t=1;2;goverblocksizem,andXistheseaclutteramplitudedata.TobetterseethevariationofVar(X(m)t)withm,Var(X(m)t)ismultipliedbym2.Whentheautocorrelationofthedatadecaysexponentiallyfast(suchasmodeledbyanautoregressive(AR)process),Var(X(m)t)m1.HereVar(X(m)t)decaysmuchfaster.Afractalprocesswouldhavem2Var(X(m)t)m.However,thisisnotthecase.Therefore,neitherARmodelingnoridealtextbookfractaltheorycanbereadilyappliedhere. targetdetectionwithinseaclutterremainsatremendouschallenge.Weshallreturntoseaclutterlater. 18
PAGE 19
ExampleofHeartratevariabilitydataforanormalsubject.(a)theentiresignal,(b,c)thesegmentsofsignalsindicatedasAandBin(a);(d,e)powerspectraldensityforthesignalsshownin(b,c). maintainedconstantandnoperturbinginuencescanbeidentied.IthasbeenobservedthatHRVisrelatedtovariouscardiovasculardisorders.Therefore,analysisofHRVisveryimportantinmedicine.However,thistaskisverydicult,sinceHRVdataarehighlycomplicated.AnexampleisshowninFig. 14 ,foranormalyoungsubject.Evidently,thesignalishighlynonstationaryandmultiscaled,appearingoscillatoryforsomeperiodoftime(Figs. 14 (b,d)),andthenvaryingasapowerlawforanotherperiodoftime(Figs. 14 (c,e)).Thelatterisanexampleofthesocalled1=fprocesses,whichwillbediscussedindepthinlaterchapters.Whilethemultiscalenatureofsuchsignalscannotbefullycharacterizedbyexistingmethods,thenonstationarityofthedataisevenmoretroublesome,sinceitrequiresthedatatobeproperlysegmentedbeforefurtheranalysis 19
PAGE 20
Atthecenterofchaostheoryistheconceptofsensitivedependenceoninitialconditions:averyminordisturbanceininitialconditionsleadstoentirelydierentoutcomes.AnoftenusedmetaphorillustratingthispointisthatasunnyweatherinNewYorkcouldbereplacedbyarainyonesometimeinthefutureafterabutteryapsitswingsinBoston.Suchafeaturecontrastssharplywiththetraditionalview,largelybasedonourexperiencewithlinearsystems,thatsmalldisturbances(orcauses)canonlygenerateproportionaleects,andthatinorderforthedegreeofrandomnesstoincrease,thenumberofdegreesoffreedomhastobeinnite. Nodenitionofthetermchaosisuniversallyacceptedyet,butalmosteveryonewouldagreeonthethreeingredientsusedinthefollowingworkingdenition: Chaosisaperiodiclongtermbehaviorinadeterministicsystemthatexhibitssensitivedependenceoninitialconditions.\Aperiodiclongtermbehavior"meansthattherearetrajectorieswhichdonotsettledowntoxedpoints,periodicorbits,orquasiperiodicorbitsast!1.Forpracticalreasons,weshouldrequirethatsuchtrajectoriesarenottoorare.Forinstance,wecould 20
PAGE 21
Chaosisalsocommonlycalledastrangeattractor.Thetermattractorisalsodiculttodeneinarigorousway.Wewantadenitionthatisbroadenoughtoincludeallthenaturalcandidates,butrestrictiveenoughtoexcludetheimposters.Thereisstilldisagreementaboutwhattheexactdenitionshouldbe. Looselyspeaking,anattractorisasettowhichallneighboringtrajectoriesconverge.Stablexedpointsandstablelimitcyclesareexamples.Moreprecisely,wedeneanattractortobeaclosedsetAwiththefollowingproperties: 1. 2. 3. Nowwecandeneastrangeattractortobeanattractorthatexhibitssensitivedependenceoninitialconditions.Strangeattractorswereoriginallycalledstrangebecausetheyareoftenfractalsets.Nowadaysthisgeometricpropertyisregardedaslessimportantthanthedynamicalpropertyofsensitivitydependenceoninitialconditions.Thetermschaoticattractorandfractalattractorareusedwhenonewishestoemphasizeoneortheotherofthoseaspects. 21
PAGE 22
Fornow,weshallbesatisedwithanintuitivedenitionofafractal:asetthatshowsirregularbutselfsimilarfeaturesonmanyorallscales.Selfsimilaritymeansapartofanobjectissimilartootherpartsortothewhole.Thatis,ifweviewanirregularobjectwithamicroscope,whetherweenlargetheobjectby10times,orby100times,orevenby1000times,wealwaysndsimilarobjects.Tounderstandthisbetter,letusimaginewewereobservingapatchofwhiteclouddriftingawayinthesky.Oureyeswereratherpassive:wewerestaringmoreorlessatthesamedirection.Afterawhile,thepartofthecloudwesawdriftedaway,andwewereviewingadierentpartofthecloud.Nevertheless,ourfeelingremainsmoreorlessthesame. Mathematically,fractalischaracterizedbyapowerlawrelation,whichtranslatestoalinearrelationinloglogscale.Fornow,letusagainresorttoimaginationwewerewalkingdownawildandveryjaggedmountaintrailorcoastline.Wewouldliketoknowthedistancecoveredbyourroute.Supposeourrulerhasalengthofwhichcouldbeourstepsize,anddierenthikersmayhavedierentstepsizesapersonridingahorsehasahugestepsize,whileagroupofpeoplewithalittlechildmusthaveatinystepsize.Thelengthofourrouteis whereN()isthenumberofintervalsneededtocoverourroute.ItismostremarkablethattypicallyN()scaleswithinapowerlawmanner, withDbeinganoninteger,1
PAGE 23
Nowletusexaminetheconsequenceof1
PAGE 24
24
PAGE 25
2 ,werstintroducethenewconceptofPSIC.ThenweextendtheconceptofPSICtohighdimensionalcase.WeshallpresentageneralcomputationalframeworkforassessingPSICfromrealworlddata.WealsoapplytheproposedcomputationalproceduretostudynoisefreeandnoisylogisticandHenonmapsattheedgeofchaos.weshowthatwhennoiseisabsent,PSICishardtodetectfromascalartimeseries.However,whenthereisdynamicnoise,motionsaroundtheedgeofchaos,beitsimplyregularortrulychaoticwhenthereisnonoise,allcollapseontothePSICattractor.Hence,dynamicnoisemakesPSICobservable.InChapter 3 ,wedemonstratethatthenewframeworkofPSICcanreadilycharacterizerandomfractalprocesses,including1=fprocesseswithlongrangecorrelationsandLevyprocesses.FromourstudyinChapters 2 and 3 ,itisclearthattheconceptofPSICnotjustbridgesstandardchaostheoryandrandomfractaltheory,butinfactprovidesamoregeneralframeworktoencompassboththeoriesasspecialcases.ToillustratethepoweroftheframeworkofPSIC,inChapters 4 and 5 ,weapplyittostudytwoimportantbutchallengingproblems:Internetdynamicsandseaclutterradarreturns,respectively.Bothareoutstandingexamplesofcomplexdynamicalsystemswithbothnonlinearityandrandomness.ThenonlinearityofInternetdynamicscomesfromAIMD(additiveincreaseandmultiplicativedecrease)operationofTCP(TransmissionControlprotocol),whilethestochasticcomponentcomesfromtherandomuserbehavior.Seaclutterisanonlineardynamicalprocesswithstochasticfactorsduetointerferenceofvariouswindandswellwavesandtolocalatmosphericturbulence.WeshallshowthatthenewtheoreticalframeworkofPSICcaneectivelycharacterizethecomplicatedInternetdynamicsaswellastorobustlydetectlowobservabletargetsfromseaclutterradarreturnswithhighaccuracy.InChapter 6 ,weshallintroducethebasictheoryofanexcellentmultiscalemeasurethescaledependentLyapunovexponent(SDLE),anddevelopaneectivealgorithmtocomputethemeasure.To 25
PAGE 26
7 26
PAGE 27
Theformalismofnonextensivestatisticalmechanics(NESM)[ 2 ]hasfoundnumerousapplicationstothestudyofsystemswithlongrangeinteractions[ 3 { 6 ]andmultifractalbehavior[ 7 8 ],andfullydevelopedturbulence[ 8 { 10 ],amongmanyothers.Inordertocharacterizeatypeofmotionwhosecomplexityisneitherregularnorfullychaotic/random,recentlytheconceptofexponentialsensitivitytoinitialconditions(ESIC)ofdeterministicchaoshasbeengeneralizedtopowerlawsensitivitytoinitialconditions(PSIC)[ 7 11 ].Mathematically,theformulationofPSICcloselyparallelsthatofNESM.PSIChasbeenappliedtothestudyofdeterministic1Dlogisticlikemapsand2DHenonmapattheedgeofchaos[ 7 11 { 20 ],yieldingconsiderablenewinsights.Inthischapter,werstbrieydescribesomebackgroundknowledgeaboutchaos.Inparticular,wewilldiscussadirectdynamicaltestfordeterministicchaosdevelopedbyGaoandZheng[ 21 22 ].Thismethodoersamorestringentcriterionfordetectinglowdimensionalchaos,andcansimultaneouslymonitormotionsinphasespaceatdierentscales.ThenweintroducethemoregeneralizedconceptofPSIC,andextendthenewconcepttohighdimensionalcase.Specically,wepresentageneralcomputationalframeworkforassessingPSICinatimeseries.Wealsoapplytheproposedcomputationalproceduretostudytwomodelsystems:noisefreeandnoisylogisticandHenonmapsattheedgeofchaos. 27
PAGE 28
where1ispositiveandcalledthelargestLyapunovexponent.Duetotheboundednessoftheattractorandtheexponentialdivergencebetweennearbytrajectories,astrangeattractortypicallyisafractal,characterizedbyasimpleandelegantscalinglawasdenedbyEq. 1{4 .Anonintegralfractaldimensionofanattractorcontrastssharplywiththeintegervaluedtopologicaldimension(whichis0fornitenumberofisolatedpoints,1forasmoothcurve,2forasmoothsurface,andsoon). Conventionally,ithasbeenassumedthatatimeserieswithanestimatedpositivelargestLyapunovexponentandanonintegralfractaldimensionischaotic.However,ithasbeenfoundthatthisassumptionmaynotbeasucientindicationofdeterministicchaos.Onecounterexampleisthesocalled1=fprocesses.Suchprocesseshavespectraldensity where0
PAGE 29
23 { 27 ],estimationofthestrengthofmeasurementnoiseinexperimentaldata[ 28 29 ],pathologicaltremors[ 30 ],shearthickeningsurfactantsolutions[ 31 ],diluteshearedaqueoussolutions[ 32 ],andserratedplasticows[ 33 ]. Nowletusbrieydescribethedirectdynamicaltestfordeterministicchaos[ 21 22 ].Givenascalartimeseries,x(1);x(2);:::;x(N),(assuming,forconvenience,thattheyhavebeennormalizedtotheunitinterval[0,1]),onerstconstructsvectorsofthefollowingformusingthetimedelayembeddingtechnique[ 34 { 36 ]: withmbeingtheembeddingdimensionandLthedelaytime.Forexample,whenm=3andL=4,wehaveV1=[x(1);x(5);x(9)],V2=[x(2);x(6);x(10)],andsoon.Fortheanalysisofpurelychaoticsignals,mandLhastobechosenproperly.Thisistheissueofoptimalembedding(see[ 21 22 37 ]andreferencestherein).Afterthescalartimeseriesisproperlyembedded,onethencomputesthetimedependentexponent(k)curves (k)=lnkVi+kVj+kk kViVjk(2{4) withrkViVjkr+r,whererandrareprescribedsmalldistances.Theanglebracketsdenoteensembleaveragesofallpossiblepairsof(Vi;Vj).Theintegerk,calledtheevolutiontime,correspondstotimekt,wheretisthesamplingtime.Notethatgeometrically(r;r+r)denesashell,andashellcapturesthenotionofscale.Thecomputationistypicallycarriedoutforasequenceofshells.ComparingEq. 2{4 withEq. 2{1 ,onenoticesthatkVi+kVj+kkplaystheroleofd(t),whilekViVjkplaystheroleofd(0).Figure 21 (a)showsthe(k)vs.kcurvesforthechaoticLorenzsystemdrivenby 29
PAGE 30
dt=16(xy)+D1(t)dy dt=xz+45:92xy+D2(t)dz dt=xy4z+D3(t)(2{5) with=0,=ij(tt0),i;j=1;2;3.NotethatD2isthevarianceoftheGaussiannoiseterms,andD=0describesthecleanLorenzsystem.WeobservefromFig. 21 (a)thatthe(k)curvesarecomposedofthreeparts.Thesecurvesarelinearlyincreasingfor0kka.Theyarestilllinearlyincreasingforkakkp,butwithaslightlydierentslope.Theyareatforkkp.NotethattheslopeofthesecondlinearlyincreasingpartgivesanestimateofthelargestpositiveLyapunovexponent[ 21 22 ].kaisrelatedtothetimescaleforapairofnearbypoints(Vi;Vj)toevolvetotheunstablemanifoldofViorVj.Itisontheorderoftheembeddingwindowlength,(m1)L.kpisthepredictiontimescale.Itislongerforthe(k)curvesthatcorrespondtosmallershells.Thedierencebetweentheslopesoftherstandsecondlinearlyincreasingpartsiscausedbythediscrepancybetweenthedirectiondenedbythepairofpoints(Vi;Vj)andtheunstablemanifoldofViorVj.ThisfeaturewasrstobservedbySatoetal.[ 38 ]andwasusedbythemtoimprovetheestimationoftheLyapunovexponent.Therstlinearlyincreasingpartcanbemadesmallerorcanevenbeeliminatedbyadjustingtheembeddingparameterssuchasbyusingalargervalueform.Notethatthesecondlinearlyincreasingpartsofthe(k)curvescollapsetogethertoformanenvelope.Itisthisveryfeaturethatformsthedirectdynamicaltestfordeterministicchaos[ 21 22 ].Thisisbecausethe(k)curvesfornoisydata,suchaswhitenoiseor1=fprocesses,areonlycomposedoftwoparts,anincreasing(butnotlinear)partfork(m1)Landaatpart[ 21 22 ].Furthermore,dierent(k)curvesfornoisydataseparatefromeachother;henceanenvelopeisnotdened.Wealsonote,mostotheralgorithmsforestimatingthelargestLyapunovexponentisequivalenttocompute(k)forr
PAGE 31
Onecanexpectthatthebehaviorofthe(k)curvesforanoisychaoticsystemliesinbetweenthatofthe(k)curvesforacleanchaoticsystemandthatofthe(k)curvesforwhitenoiseorfor1=fprocesses.Thisisindeedso.AnexampleisshowninFig. 21 (b)forthenoisyLorenzsystemwithD=4.WeobservefromFig. 21 (b)thatthe(k)curvescorrespondingtodierentshellsnowseparate.Therefore,anenvelopeisnolongerdened.Thisseparationislargerbetweenthe(k)curvescorrespondingtosmallershells,indicatingthattheeectofnoiseonthesmallscaledynamicsisstrongerthanthatonthelargescaledynamics.Alsonotethatka+kpisnowontheorderoftheembeddingwindowlength,andisalmostthesameforallthe(k)curves.Withstrongernoise(D>4),the(k)curveswillbemorelikethoseforwhitenoise[ 21 22 ]. Finally,weconsidertheMackeyGlassdelaydierentialsystem[ 39 ].Whena=0:2;b=0:1;c=10;=30,ithastwopositiveLyapunovexponents,withthelargestLyapunovexponentcloseto0.007[ 21 ]. 1+x(t+)cbx(t)(2{6) HavingtwopositiveLyapunovexponentswhilethevalueofthelargestLyapunovexponentofthesystemisnotmuchgreaterthan0,onemightbeconcernedthatitmaybediculttodealwiththissystem.Thisisnotthecase.Infact,thissystemcanbeanalyzedasstraightforwardlyasotherdynamicalsystemsincludingtheLogisticmap,HenonmapandtheRosslersystem.Anexampleofthe(k)vs.kcurvesforthechaoticMackeyGlass 31
PAGE 32
Timedependentexponent(k)vs.evolutiontimekcurvesfor(a)cleanand(b)noisyLorenzsystem.Sixcurves,frombottomtotop,correspondtoshells(2(i+1)=2;2i=2)withi=8,9,10,11,12,and13.Thesamplingtimeforthesystemis0.03sec,andembeddingparametersarem=4,L=3.5000pointsareusedinthecomputation. 32
PAGE 33
Timedependentexponent(k)vs.evolutiontimekcurvesfortheMackeyGlasssystem.Ninecurves,frombottomtotop,correspondtoshells(2(i+1)=2;2i=2)withi=5;6;,and13.Thecomputationwasdonewithm=5;L=1,and5000pointssampledwithatimeintervalof6. delaydierentialsystemisshowninFig. 22 ,wherewehavefollowed[ 21 ]andusedm=5;L=1,and5000pointssampledwithatimeintervalof6.Clearly,weobservethatthereexistsawelldenedcommonenvelopetothe(k)curves.Actually,theslopeoftheenvelopeestimatesthelargestLyapunovexponent.Thisexampleillustratesthatwhenoneworksintheframeworkofthe(k)curvesdevelopedbyGaoandZheng[ 21 22 ],onedoesnotneedtobeveryconcernedaboutnonuniformgrowthrateinhighdimensionalsystems. 2{7 ,wherex(0)istheinnitesimaldiscrepancyintheinitialcondition,x(t)isthediscrepancyattimet>0. x(0)(2{7) 33
PAGE 34
Tsallisandcoworkers[ 7 11 ]havegeneralizedEq. 2{9 to whereqiscalledtheentropicindex,andqisinterpretedtobeequaltoKq,thegeneralizationoftheKolmogorovSinaientropy.Eq. 2{10 denesthePSICinthe1Dcase.Obviously,PSICreducestoESICwhenq!1.ThesolutiontoEq. 2{10 is Whentislargeandq6=1,(t)increaseswithtasapowerlaw, whereC=[(1q)q]1=(1q).ForEq. 2{12 todeneanunstablemotionwithq>0,wemusthaveq1.Laterweshallmapdierenttypesofmotionstodierentrangesofq. ToapplyPSICtotheanalysisoftimeseriesdata,onecanrstconstructaphasespacebyconstructingembeddingvectorsViasdenedbyEq. 2{3 .Eq. 2{11 canthenbegeneralizedtohighdimensionalcase[ 40 ], jjV(0)jj=1+(1q)(1)qt1=(1q)(2{13) wherejjV(0)jjistheinnitesimaldiscrepancybetweentwoorbitsattime0,jjV(t)jjisthedistancebetweenthetwoorbitsattimet>0,qistheentropicindex,and(1)qistherstqLyapunovexponent,correspondingtothepowerlawincreaseoftherstprincipleaxisofaninnitesimalballinthephasespace.(1)qmaynotbeequaltoKq.Thisisunderstoodbyrecallingthatforchaoticsystems,theKolmogorovSinaientropyisthesum 34
PAGE 35
2{13 againgivesapowerlawincreaseof(t)witht.Notethatundertheaboveframework,themotionmaynotbelikefullydevelopedchaos,thusnotergodic[ 41 ]. WenowconsiderthegeneralcomputationalframeworkforPSIC.Givenanitetimeseries,theconditionofjjV(0)jj!0cannotbesatised.Inordertondthelawgoverningthedivergenceofnearbyorbits,onethushastoexaminehowtheneighboringpoints,(Vi;Vj),inthephasespace,evolvewithtime,byformingsuitableensembleaverages.Noticethatif(Vi1;Vj1)and(Vi2;Vj2)aretwopairsofnearbypoints,whenjjVi1Vj1jjjjVi2Vj2jj1,thentheseparationssuchasjjVi1+tVj1+tjjandjjVi2+tVj2+tjjcannotbesimplyaveragedtoprovideestimatesforqand(1)q.Infact,itwouldbemostconvenienttoconsiderensembleaveragesofpairsofpoints(Vi;Vj)thatallfallwithinaverythinshell,r1jjViVjjjr2,wherer1andr2areclose.TheseargumentssuggestthatthetimedependentexponentcurvesdenedbyEq. 2{4 provideanaturalframeworktoassessPSICfromatimeseries.Thisisindeedso.Infact,wehave ln(t)(t)=lnkVi+tVj+tk kViVjk(2{14) Now,bythediscussionsinSec. 2.1 ,itisclearthatPSICisageneralizationofESIC:aslongasthe(k)curvesfromdierentshellsformalinearenvelope,thenq=1andthemotionischaotic.Thenextquestionis:DoesPSICalsoincluderandomfractalsasspecialcases?Theanswerisyes.Itwillbegiveninthenextchapter.There,wewillalsogainabetterunderstandingofthemeaningof(1)q.TherestofthischaptershowsexamplesofcharacterizingtheedgeofchaosbythenewframeworkofPSIC.Specically,westudytimeseriesgeneratedfromnoisefreeandnoisylogisticandHenonmaps. 35
PAGE 36
26 ],andphysiologicalnoiseinbiologicalsystems[ 30 ]),westudyboththedeterministicandthenoisylogisticmap: whereaisthebifurcationparameterandnisawhiteGaussiannoisewithmeanzeroandvariance1.Theparametercharacterizesthestrengthofnoise.Forthecleansystem(=0),theedgeofchaosoccursattheaccumulationpoint,a1=3:569945672.Weshallstudythreeparametervalues,a1=a10:001,a1,anda2=a1+0:001.Whennoiseisabsent,a1correspondstoaperiodicmotionwithperiod25,whilea2correspondstoatrulychaoticmotion.Weshallonlystudytransientfreetimeseries.InFigs. 23 (ac),wehaveplottedthe(t)vs.tcurvesforparametervaluesa1,a1,anda2,respectively.WeobservefromFig. 23 (a)thatthevariationof(t)withtisperiodic(withperiod16,whichishalfoftheperiodofthemotion)whenthemotionisperiodic.Thisisagenericfeatureofthe(t)curvesfordiscreteperiodicattractors,whentheradiusoftheshellislargerthanthesmallestdistancebetweentwopointsontheattractor(whenaperiodicattractoriscontinuous,(t)canbearbitrarilycloseto0).IthasbeenfoundbyTsallisandcoworkers[ 11 ]thatattheedgeofchaosforthelogisticmap,(t)isgivenbyEq. 2{11 withq0:2445.Surprisingly,wedonotobservesuchadivergenceinFig. 23 (b).Infact,ifoneplots(t)vs.lnt,oneonlyobservesacurvethatincreasesveryslowly(similartothatshowninFig. 24 (a)).ThemoreinterestingpatternistheonethatisshowninFig. 23 (c),whereweobservealinearlyincreasingly(t)vs.tcurve.Infact,showninFig. 23 (c)aretwosuchcurves,correspondingtotwodierentshells.Veryinterestingly,thetwocurvescollapsetoformacommonenvelopeinthelinearlyincreasingpartofthe 36
PAGE 37
Timedependentexponent(t)vs.tcurvesfortimeseriesgeneratedfromthenoisefreelogisticmapwith(a)a1=a10:001,wherethemotionisperiodicwithperiod25,(b)a1=3:569945672,and(c)a2=a1+0:001,wherethemotionischaotic.Plottedin(df)are(t)vs.lntcurvesforthenoisylogisticmapwith=0:001.Verysimilarresultswereobtainedwhen=0:0001.Shownin(cf)areactuallytwocurves,correspondingtotwodierentshells.104pointswereusedinthecomputation,withembeddingparametersm=4;L=1.However,solongasm>1,theresultsarelargelyindependentofembedding.Whenm=1,the(t)curvesarenotsmooth,andtheestimated1=(1q)valueismuchsmallerthanthetheoreticalvalue. curve.TheslopeoftheenvelopegivesagoodestimateofthelargestpositiveLyapunovexponent.Thisisagenericfeatureofchaos[ 21 22 ],asexplainedinSec. 2.1 .Sincethechaosstudiedhereisclosetotheedgeofchaos,thecurvesshowninFig. 23 (c)arelesssmooththanthosereportedearlier[ 21 { 23 26 ]. WhycannotthetheoreticalpredictionofPSICattheedgeofchaosforthelogisticmapbeobservedfromacleantimeseries?Inarecentlypublishedveryinterestingand 37
PAGE 38
42 ]suggeststhatdynamicnoisemaybeoffundamentalimportancetotheTsallisNESM.MaydynamicnoiseplaysimilarlysignicantroleforPSIC?Asisshownbelow,theanswerisyes.InFigs. 23 (df),wehaveshownthe(t)vs.lntcurvesforthethreeparametersconsidered,withnoisestrength1=0:001.Infact,shownineachgurearetwocurves,correspondingtotwodierentshells.Theyparallelwitheachother.Theslopesofthosecurvesareabout1.20,closetothetheoreticalvalueof1=(10:2445)1:32.WhileitisverysatisfactorytoobservePSICattheedgeofchaos,itismorethrillingtoobservethecollapseofregularaswellaschaoticmotionsontothePSICattractoraroundtheedgeofchaos.ThissigniesthestablenessofPSICwhenthereisdynamicnoise.ItisimportanttoemphasizethattheresultsshowninFigs. 23 (df)arelargelyindependentofthenoisestrength,solongasnoiseisneithertooweaknortoostrong.Forexample,verysimilarresultshavebeenobservedwith2=1=10=0:0001. BeforewemoveontodiscusstheHenonmapneartheedgeofchaos,letuscommentonthedierencebetweentheESICandthePSIC.WhentheESICisthecase,the(t)vs.tcurvesarestraightforTatTp,whereTpisapredictiontimescale,andTaisatimescalefortheinitialseparationtoevolvetothemostunstabledirection.When(t)vs.tcurvescorrespondingtodierentshellsareplottedtogether,theycollapsetogetherforaconsiderablerangeoft.WhenthePSICisthecase,the(t)vs.lntcurvesarefairlystraight,andthe(t)vs.lntcurvescorrespondingtodierentshellsseparateandparallelwitheachother.Thisfeatureisthefundamentalreasonthatensembleaveragesaremostconvenientlyformedbyrequiringneighboringpairsofpointstoallfallwithinathinshell.Forexample,iftwodistinctthinshells,describedbyr1jjXiXjjjr2andr2jjViVjjjr3,arejoinedtogethertoformasinglethickershell,describedbyr1jjViVjjjr3,andoneaveragestheseparationbetweenXiandXjbeforetakingthelogarithm,thentheresultingslopebetween(t)andlnt,assumedtobestilllinear,willhavetobesmallerthanthatestimatedfromthetwothinshellsseparately.Atthispointitisalsoworthnotingthatthelinearlyincreasingpartofthe(t)vs.lntcurves 38
PAGE 39
2{13 whendealingwithexperimentaltimeseries.Infact,wehavefoundthatifwedoso,theestimated1=(1q)valuefromFigs.1(df)increasestoabout1.27,muchclosertothetheoreticalvalueof1.32. NextletusconsidertheHenonmap,wherexandyarewhiteGaussiannoise,andareuncorrelatedwitheachother.Theparametermeasuresthestrengthofthenoise. Tirnakli[ 20 ]studiedthedeterministicversionofthismapattheedgeofchaos,forparametervaluesac=1:40115518;b=0;ac=1:39966671;b=0:001;andac=1:38637288;b=0:01.Wehavestudiedboththenoisefreeandnoisymapforparametervalueslistedabove,andfoundverysimilarresultstothosepresentedinFig. 23 .InFigs. 24 (a,b),wehaveplotted(t)vs.lntcurvesforthenoisefreeandnoisymapforac=1:38637288;b=0:01.Ascanbeobservedclearly,Fig. 24 (b)isverysimilartoFigs. 23 (df),whiletheslopeforthecurvesinFig. 24 (a)ismuchsmallerthanthatinFig. 24 (b).Again,wehaveobserved(butnotshownbygureshere,sincetheyareverysimilartoFigs. 23 (df)andFig. 24 (b))thatdynamicnoisemakestheregularandchaoticmotionstocollapseontothePSICattractorforparametersaroundthosedeningtheedgeofchaos,andthatsuchtransitionsarelargelyindependentofthenoisestrength. Tosummarize,wehavedescribedaneasilyimplementableprocedureforcomputationallyexaminingPSICfromatimeseries.Bystudyingtwomodelsystems:noisefreeandnoisylogisticandHenonmapsneartheedgeofchaos,wehavefoundthatwhenthereisno 39
PAGE 40
Timedependentexponent(t)vs.tcurvesfortimeseriesgeneratedfrom(a)thenoisefreeHenonmapwithac=1:38637288;b=0:01.Plottedin(b)aretwo(t)vs.lntcurves(correspondingtotwodierentshells)forthenoisymapwith=0:001.Similarresultswereobtainedwith=0:0001.104pointswereusedinthecomputation,withembeddingparametersm=4;L=1. noise,thePSICattractorcannotbeobservedfromascalartimeseries.However,whendynamicnoiseispresent,motionsaroundtheedgeofchaos,beitsimplyregularortrulychaoticwhenthereisnonoise,allcollapseontothePSICattractor.WhiletheseexaminationssignifytheubiquityofPSIC,theyalsohighlighttheimportanceofdynamicnoise.Theexistenceofthelatterisperhapstheveryreasonthattrulychaotictimeseriescanseldombeobserved. 40
PAGE 41
InChapter 2 ,wehavediscussedthatPSICisageneralizationofESIC:aslongasthe(t)curvesfromdierentshellsformalinearenvelope,thenq=1andthemotionischaotic.WehavealsoshownthatedgeofchaoscanbewellcharacterizedbyPSIC.Inthischapter,weshallstudytwomajortypesofrandomfractals:1=fprocesseswithlongrangecorrelationsandLevyprocesses.WeshallshowbothanalyticallyandthroughnumericalsimulationsthatbothtypesofprocessescanbereadilycharacterizedbyPSIC. for>0,0
PAGE 42
3{3 canbesimplyderivedfromEq. 3{4 ;wehavelisteditasaseparateequationforconvenienceoffuturereference.Ifweconsideronlysecondorderstatistics,orifaprocessisGaussian,Eq. 3{2 toEq. 3{4 canbeusedinsteadofEq. 3{1 todeneaselfsimilarprocess. Averyusefulwayofdescribingaselfsimilarprocessisbyitsspectralrepresentation.Strictlyspeaking,theFouriertransformofX(t)isundened,duetothenonstationarynatureofX(t).Onecan,however,considerX(t)inaniteinterval,say,0
PAGE 43
with ProcesseswithpowerspectraldensitiesasdescribedbyEq. 3{7 arecalled1=fprocesses.Typically,thepowerlawrelationshipsthatdenetheseprocessesextendoverseveraldecadesoffrequency.Suchprocesseshavebeenfoundinnumerousareasofscienceandengineering.Someoftheolderliteraturesonthissubjectcanbefound,forexample,inPress[ 45 ],Bak[ 46 ],andWornell[ 47 ].Someofthemorerecentlydiscovered1=fprocessesareintracengineering[ 48 { 50 ],DNAsequence[ 51 { 53 ],humancognition[ 54 ],ambiguousvisualperception[ 55 56 ],coordination[ 57 ],posture[ 58 ],dynamicimages[ 59 60 ],andthedistributionofprimenumbers[ 61 ],amongmanyothers.Itisfurtherobservedthatprinciplecomponentanalysisofsuchprocessesleadstopowerlawdecayingeigenvaluespectrum[ 62 ]. Twoimportantprototypicalmodelsfor1=fprocessesarethefractionalBrownianmotion(fBm)processesandtheON/OFFintermittencywithpowerlawdistributedONandOFFperiods.Below,weapplytheconceptofPSICtobothtypesofprocesses. 1. AllBrownianpathsstartattheorigin:B(0)=0. 2. For0
PAGE 44
WemayinferfromtheabovedenitionthattheprobabilitydistributionfunctionofBisgivenby: 2sdu(3{9) Thisfunctionalsosatisesthescalingproperty: Inotherwords,B(t)and1=2B(t)havethesamedistribution.Thus,weseefromEq. 3{1 thatBmisaselfsimilarprocesswithHurstparameter1/2. SupposewehavemeasuredB(t1);B(t2);t2t1>0.WhatcanwesayaboutB(s);t1
PAGE 45
WenotethataBrownianpathisnotadierentiablefunctionoftime.Heuristically,thiscanbeunderstoodinthisway:considerthevariable,B(t+s)B(t),withvariances.Itsstandarddeviation,whichisameasureofitsorderofmagnitude,isp AlthoughB(t)isalmostsurelynotdierentiableint,symbolicallyonestilloftenwrites wherew(t)isstationarywhiteGaussiannoise,andextendstheaboveequationtot<0throughtheconvention,Zt0Z0t 3{12 Toseehow,letuspartition[0;t]intonequallyspacedintervals,t=t=n.Sincew(t)areGaussianrandomvariableswithzeromeanandunitvariance,inorderforB(t)tohavevariancet,weshouldhave Thatis,thecoecientis(t)1=2insteadoft,asonemighthaveguessed.Typically,t1;hence,(t)1=2t.Thus,ifoneincorrectlyusestinsteadof(t)1=2,oneisseverelyunderestimatingthestrengthofthenoise. 45
PAGE 46
Whitenoise(a)andBrownianmotion(b).Axesarearbitrary. Whenthetimeisgenuinelydiscrete,ortheunitsoftimecanbearbitrary,onecantakettobe1unit,andEq. 3{13 becomes Eq. 3{14 providesperhapsthesimplestmethodofgeneratingasampleofBm.AnexampleisshowninFig. 31 Eq. 3{14 isalsoknownasarandomwalkprocess.Amoresophisticatedrandomwalk(orjump)processisgivenbysummingupaninnitenumberofjumpfunctions:Ji(t)=Ai(tti),where(t)isaunitstepfunction andAi;tiarerandomvariableswithGaussianandPoissondistributions,respectively. 46
PAGE 47
whereWisaGaussianrandomvariableofzeromeanandunitvariance,andHistheHurstparameter.WhenH=1=2,theprocessreducestothestandardBrownianmotion(Bm).ItistrivialtoverifythatthisprocesssatisesthedeningEq. 3{1 foraselfsimilarstochasticprocess. FractionalBrownianmotionBH(t)isaGaussianprocesswithmean0,stationaryincrements,variance andcovariance: 2ns2H+t2Hjstj2Ho(3{18) whereHistheHurstparameter.DuetoitsGaussiannature,accordingtoEq. 3{2 toEq. 3{4 ,theabovethreepropertiescompletelydetermineitsselfsimilarcharacter.Fig. 32 showsseveralfBmprocesseswithdierentH. Roughly,thedistributionforanfBmprocessis[ 63 ]: (H+1=2)Zt(t)H1=2dB()(3{19) where(t)isthegammafunction:(t)=Z10yt1eydy Inotherwords,thisprocessisinvariantfortransformationsconservingthesimilarityvariablex=tH. 47
PAGE 48
SeveralfBmprocesseswithdierentH. TheintegraldenedbyEq. 3{19 diverges.Themoreprecisedenitionis (H+1=2)ZtK(t)dB()(3{21) wherethekernelK(t)isgivenby[ 63 ] TheincrementprocessoffBm,Xi=BH((i+1)t)BH(it);i1,wheretcanbeconsideredasamplingtime,iscalledfractionalGaussiannoise(fGn)process.ItisazeromeanstationaryGaussiantimeseries.NotingthatE(XiXi+k)=Ef[BH((i+1)t)BH(it)][BH((i+1+k)t)BH((i+k)t)]g
PAGE 49
3{18 ,oneobtainstheautocovariancefunction(k)forthefGnprocess: 2n(k+1)2H2k2H+jk1j2Ho;k0(3{23) Noticethattheexpressionisindependentoft.Therefore,withoutlossofgenerality,onecouldtaket=1.Inparticular,wehave(1)=1 222H2 (i) WhenH=1=2,(k)=0fork6=0.ThisisthewellknownpropertyofwhiteGaussiannoise. (ii) When00. Properties(ii)and(iii)areoftentermedantipersistentandpersistentcorrelations[ 43 ],respectively. Next,weconsiderthebehaviorof(k)whenkislarge.Notingthatwhenjxj1,(1+x)1+x+(1) 2x2 hence, whenH6=1=2.WhenH=1=2,(k)=0fork1,theXi'saresimplywhitenoise. Nowwearereadytocharacterize1=fprocessesintheframeworkofPSIC.Recallingthatthedeningpropertyfora1=fprocessisthatitsvarianceincreaseswithtast2H. 49
PAGE 50
Therefore, Expressingtintermsof,wehave ComparingwiththedeningequationofPSIC(Eq. 2{10 ),wendthat 2(3{29) Noticing0
PAGE 51
Timedependentexponent(k)vs.lnkcurvesforfractionalBrownianmotion. ON/OFFsourceswithN=3,X1(t),X2(t),X3(t),andtheirsummationS3(t). 51
PAGE 52
Thisexpressionisoftencalledthecomplementarycumulativedistributionfunction(CCDF)ofaheavytaileddistribution.Infact,itismorepopularforexpressingaheavytaileddistribution,sinceitemphasizesthetailofthedistribution.Itiseasytoprovethatwhen<2,thevarianceandallhigherthan2ndordermomentsdonotexist.Furthermore,when1,themeanalsodiverges. Whenthepowerlawrelationextendstotheentirerangeoftheallowablex,wehavetheParetodistribution: x;xb>0;2>0(3{31) whereandbarecalledtheshapeandthelocationparameters,respectively.Itcanbeproven[ 64 ]thatwhen12,theHurstparameteroftheON/OFFprocesseswithParetodistributedONandOFFperiodsisgivenas OurpurposehereistocheckwhethertheON/OFFprocesseswithParetodistributedONandOFFperiodscanbecharacterizedbyPSIC.InFig. 35 wehaveshownafewexamplesfor=1:2;1:6,and2.0.NoticingEq. 3{32 ,weseethattheslopesofthestraightlinesagaincorrectlyestimatetheHparameter.Alsonotethatwhen0<<1,1
PAGE 53
Timedependentexponent(k)vs.lnkcurvesforON/OFFprocesseswithParetodistributedONandOFFperiods. Beforeendingthissection,wewouldliketoemphasizethatthepossibilityofmisinterpreting1=fprocessesbeingdeterministicchaosneveroccursifoneworksundertheframeworkofPSIC.Underthisframework,onemonitorstheevolutionoftwonearbytrajectories.Iftwonearbytrajectoriesdivergeexponentiallyfast,thenthetimeseriesischaotic;ifthedivergenceincreasesinapowerlawmanner,thenthetrajectoriesbelongto1=fprocesses.Theseideascanbeeasilyexpressedprecisely.LetjjViVjjjbetheinitialseparationbetweentwotrajectories.Thisseparationisassumedtobenotlargerthansomeprescribedsmalldistancer.Aftertimet,thedistancebetweenthetwotrajectorieswillbejjVi+tVj+tjj.Fortruechaoticsystems,jjVi+tVj+tjj/jjViVjjjet
PAGE 54
65 { 70 ].TherearetwotypesofLevyprocesses[ 71 ].OneisLevyights,whicharerandomprocessesconsistingofmanyindependentsteps,eachstepbeingcharacterizedbyastablelaw,andconsumingaunittimeregardlessofitslength.TheotherisLevywalks,whereeachsteptakestimeproportionaltoitslength.ALevywalkcanbeviewedassampledfromaLevyightwithauniformspeed.TheincrementprocessofaLevywalk,obtainedbydierencingtheLevywalk,isverysimilartoanON/OFFtrainwithpowerlawdistributedONandOFFperiods.Therefore,inthefollowing,weshallnotbefurtherconcernedaboutit.WeshallfocusonLevyights.NotethatLevyights,havingindependentsteps,arememorylessprocessescharacterizedbyH=1=2,irrespectiveofthevalueoftheexponentcharacterizingthestablelaws[ 71 ].Inotherwords,methodssuchasdetrendeductuationanalysis(DFA)[ 72 ]cannotbeusedtoestimatetheparameter. NowwedeneLevyightsmoreprecisely.A(standard)symmetricstableLevyprocessfL(t);t0gisastochasticprocesswiththefollowingproperties[ 73 ]: 1. 2. 3. 54
PAGE 55
whereY1;Y2;areindependentrandomvariables,eachhavingthesamedistributionasY.FromEq. 3{33 ,onethenndsthatnVarY=n2=VarY.Forthedistributiontobevalid,0<2.When=2,thedistributionisGaussian,andhence,thecorrespondingLevyightisjustthestandardBrownianmotion.When0<<2,thedistributionisheavytailed,P[Xx]x;x!1,andhasinnitevariance.Furthermore,when0<1,themeanisalsoinnite. Atthispoint,itisworthspendingafewparagraphstoexplainhowtosimulateanstableLevyprocess.Thekeyistosimulateitsincrementprocess,whichfollowsanstabledistribution.TherstbreakthroughforthesimulationofstabledistributionswasmadebyKanter[ 74 ],whogaveadirectmethodforsimulatingastabledistributionwith<1and=1.TheapproachwaslatergeneralizedtothegeneralcasebyChambersetal.[ 75 ].WerstdescribethealgorithmforconstructingastandardstablerandomvariableXS(;;1). (cos0cos)1=ncos[0+(1)] hasaS(;;1)distribution. Forsymmetricstabledistributionwith=0,theaboveexpressioncanbegreatlysimplied.TosimulateanarbitrarystabledistributionS(;;;),wecansimplytake 55
PAGE 56
whereZisgivenbyEq. 3{34 WehavesimulatedanumberofsymmetricstableLevymotionswithdierent.Figure 36 showstwoexamplesofLevymotionswith=1:5and1.ThedierencebetweenthestandardBrownianmotionsandLevymotionsliesinthatBrownianmotionsappeareverywheresimilarly,whileLevymotionsarecomprisedofdenseclustersconnectedbyoccasionallongjumps:thesmallertheis,themorefrequentlythelongjumpsappear. LetusnowpauseforawhileandthinkwhereLevyightslikeesotericprocessescouldoccur.Onesituationcouldbethis:amosquitoheadtoagiantspidernetandgotstuck;itstruggledforawhile,andluckily,withthehelpofagustofwind,escaped.Whenthemosquitowasstruggling,thestepsitcouldtakewouldbetiny;butthestepleadingtoitsfortunateescapemustbehuge.Asanotherexample,letusconsider(American)footballgames.Formostofthetimeduringagame,theoenseanddefensewouldbefairlybalanced,andtheoenseteammayonlybeabletoproceedforafewyards.Butduringanattackthatleadstoatouchdown,theherogettingthetouchdownoften\ies"tensofyardssomewhathehasescapedthedefense.Whilethesetwosimplesituationsareonlymeantasanalogies,rememberingthemcouldbehelpfulwhenreadingresearchpaperssearchingforLevystatisticsinvarioustypesofproblemssuchasanimalforaging.Atthispoint,wealsowishtomentionthatLevyightslikepatternshavebeenusedasatypeofscreensaverforLinuxoperatingsystems. ThesymmetricstableLevyightis1=selfsimilar.Thatis,forc>0,theprocessesfL(ct);t0gandfc1=L(t);t0ghavethesamenitedimensionaldistributions.BythisargumentaswellasEq. 3{33 ,itisclearthatthelengthofthe 56
PAGE 57
ExamplesforBrownianmotionsandLevymotions.1Dand2DBrownianmotions(a,b)vs.Levymotionswith=1:5and1for(c,d)and(e,f),respectively. 57
PAGE 58
Timedependentexponent(k)vs.lnkcurvesforLevyprocesses. motioninatimespanoft,L(t),isgivenbythefollowingscaling: L(t)/t1=(3{36) ThiscontrastswiththescalingforfractionalBrownianmotion: BH(t)/tH(3{37) Wethusidentifythat1=playstheroleofH.Therefore, Noticing0<2,fromEq. 3{38 ,wehave1q<1. 58
PAGE 59
37 .Theslopesofthestraightlinescorrectlyestimatethevaluesofusedinthesimulations. 59
PAGE 60
Timeseriesanalysisisanimportantexerciseinscienceandengineering.Oneofthemostimportantissuesintimeseriesanalysisistodeterminewhetherthedataunderinvestigationisregular,deterministicallychaotic,orsimplyrandom.Alongthisline,alotofworkhasbeendoneinthepasttwodecades.Abitsurprisingly,quiteoftentheexponentialsensitivitytoinitialconditions(ESIC)intherigorousmathematicalsensecannotbeobservedinexperimentaldata.Whilecurrentconsensusistoattributethisfacttothenoiseinthedata,weshallshowthatthemoregeneralconceptofpowerlawsensitivitytoinitialconditions(PSIC)providesaninterestingframeworkfortimeseriesanalysis.WeillustratethisbystudyingthecomplicateddynamicsofInternettransportprotocols. TheInternetisoneofthemostcomplicatedsystemsthatmanhasevermade.Inrecentyears,twotypesoffascinatingmultiscalebehaviorshavebeenfoundintheInternet.Oneisthetemporaldomainfractalandmultifractalpropertiesofnetworktrac(see[ 76 77 ]andreferencestherein),whichtypicallyrepresentsaggregationofindividualhoststreams.TheotheristhespatialdomainscalefreetopologyoftheInternet(see[ 78 ]andreferencestherein).Also,ithasbeenobservedthatthefailureofasingleroutermaytriggerroutinginstability,whichmaybesevereenoughtoinstigatearouteapstorm[ 79 ].Furthermore,packetsmaybedeliveredoutoforderorevengetdropped,andpacketreorderingisnotapathologicalnetworkbehavior[ 80 ].AsthenextgenerationInternetapplicationssuchasremoteinstrumentcontrolandcomputationalsteeringarebeingdeveloped,itisbecomingincreasinglyimportanttounderstandthetransportdynamicsinordertosustaintheneededcontrolchannelsoverwideareaconnections.Typically,thetransportdynamicsovertheInternetconnectionsarearesultofthenonlineardynamicsofthewidelydeployedTransmissionControlProtocol(TCP)interactingwiththeInternettrac,whichisstochasticandoftenselfsimilar.Thisleadsonetonaturallyexpectthe 60
PAGE 61
81 ].Iftransportdynamicsareindeedchaotic,thensteadystateanalysisandstudyofconvergencetoasteadystatewillnotbethesolefocusinpractice,particularlyforremotecontrolandcomputationalsteeringtasks.Rather,oneshouldalsofocusonthe\transient"nonconvergenttransportdynamics.Althoughrecentlyalotofresearchhasbeencarriedout[ 82 { 88 ]tobetterunderstandthisissue,alotoffundamentalproblemsremaintobeanswered.Forexample,theobservationofchaosin[ 84 ]cannotberepeated,whichleadstothesuspicionthatthechaoticmotionsreportedin[ 84 ]maycorrespondtoperiodicmotionswithverylongperiods[ 89 ].Fundamentally,weshallshow,bydevelopinga1Ddiscretemap,thatthenetworkscenarioswithonlyafewcompetingTCPowsstudiedin[ 84 ]cannotgeneratechaoticdynamics.However,withotherparametersettings,chaosispossible. Inordertoexaminethetemporalevolutionofadynamicalsystem,oneoftenmeasuresascalartimeseriesataxedpointinthestatespace.Takensembeddingtheorem[ 34 { 36 ]ensuresthatthecollectivebehaviorofthedynamicsdescribedbythedynamicalvariablescoupledtothemeasuredonecanbeconvenientlystudiedbythemeasuredscalartimeseries.ToillustratethequalitativeaspectsofTCPdynamicsbyonlymeasuringascalartimeseries,itismosteectivetoconsiderasinglebottlenecklink,as 61
PAGE 62
84 ]consistingoflonglivedTCPowsonasinglelink.Thissetuphasfourparameters,thelinkspeedCinMbps,thelinkpropagationdelaydinms,thebuersizeBinunitofpacketsof1000bytes,andthenumberofcompetingTCPowsN.WehavefoundthatNactsasacriticalbifurcationparameter.Forexample,whenC=0:1,d=10ms,B=10andNissmall,wehaveobservedonlyperiodicandquasiperiodicmotions(Figs. 41 (ad)).WhenNislarge,suchasN=17,chaosisobserved(seeFig. 42 (a,b)).WethusconcludethatoneroutetochaosinTCPdynamicsisthewellknownquasiperiodicroute.Inordertoshowmorediversieddynamicalbehaviorsuchasquasiperiodicmotionswithmorethantwoincommensurate(i.e.,independent)frequencies,below,however,weshallvaryalltheparameters.Inparticular,weshallshowthatwhenNissmall,chaoscannothappen. ForeachoftheNcompetingTCPows,onecanrecorditscongestionwindowsizedata.Letusonlyrecordoneofthem,anddenoteitbyW(i).AnexampleofquasiperiodicW(i)isshowninFig. 41 (a),whereC=0:1Mbps,delayd=10ms,B=10,N=3.Itspowerspectraldensity(PSD)isshowninFig. 41 (b).Weobservemanydiscretesharppeaks.Notethatwhenonlyaround105pointsareusedtocomputethePSD,onedoesnotobservemanypeaksinthespectrum.Ratheroneobservesabroadspectrum,andhence,istemptedtointerprettheW(i)timeseriestobechaotic.Therefore,theW(i)timeseriesisnotideallysuitedforthestudyof(quasi)periodicmotionswithlongperiods.Thismotivatesustodeneanewtimeseries,T(i),whichisthetimeintervalbetweentheonsetoftwosuccessiveexponentialbacko(i.e.,multiplicativedecrease)ofTCP,asindicatedinFig. 41 (a).Thestartofanexponentialbackoindicatestriggeringofalossepisodeorheavycongestion.Hence,theT(i)isrelatedtothewellknownroundtriptime(RTT).Infact,itcanbeconsideredtobemoreorlessequivalenttothetimeintervalbetweentwosuccessivelossbursts.TheT(i)timeseriesfortheW(i)ofFig. 41 (a)isshowninFig. 41 (c),togetherwithits 62
PAGE 63
ExampleofquasiperiodiccongestionwindowsizeW(i)timeseriessampledbyanequaltimeintervalof10ms(a)anditspowerspectraldensity(PSD)(b).TheparametersusedareC=0:1Mbps,delayd=10ms,B=10packets,whereapacketisofsize1000bytes,andthenumberofTCPstreamsN=3.(c)TheT(i)timeseriescorrespondingto(a)and(d)itsPSD.(e)AnotherT(i)timeseriescorrespondingtoC=0:5Mbps,d=10ms,B=10,N=3,and(f)itsPSD. PSDinFig. 41 (d).WeobservethatthisT(i)timeseriesisasimpleperiodicsequence,andhence,theW(i)timeseriesisquasiperiodic.However,T(i)timeseriescanstillbequasiperiodic.AnexampleisshowninFig. 41 (e),withC=0:5Mbps,delayd=10ms,B=10packets,N=3,togetherwithitsPSDisFig. 41 (f).NotethattherearestillmanydiscretesharppeaksinFig. 41 (f).Toappreciatehowmanyindependentfrequencies 63
PAGE 64
41 (e),wehavefurtherextractedtwonewtimeseries,T(1)(i),whichistheintervaltimeseriesbetweenthesuccessivelocalmaximaoftheT(i)timeseries,andT(2)(i),whichistheintervaltimeseriesbetweenthesuccessivelocalmaximaoftheT(1)(i)timeseries.WehavefoundthattheT(2)(i)issimplyperiodic.Hence,weconcludethattheW(i)timeseriesisquasiperiodicwithfourindependentfrequencies. BeforeweproceedtodiscussthechaoticTCPdynamics,wenotethatsimpleperiodicorquasiperiodicW(i)timeserieswithlessormorethanfourindependentfrequenciescanallbeobserved.Infact,wehavefoundthe\chaotic"congestionwindowsizedatastudiedin[ 84 ]tobequasiperiodicwithtwoindependentfrequencies.Furthermore,theprobabilitydistributionsfortheconstant,periodic,andquasiperiodicT(i)timeseriesareonlycomposedofafewdiscretepeaks.ThissuggeststhattheobservedquasiperiodicT(i)timeseriesconstitutesadiscretetorus. LetusnowturntothediscussionofchaoticTCPdynamics.ShowninFigs. 42 (a,c)aretwoirregularT(i)timeseries.TheparametersC,d,andBusedforFig. 42 (a)arethesameasthoseforFigs. 41 (ad).NotingfromFig. 42 thatT(i)oftencanbeontheorderof104to105,hence,anottoolongW(i)timeseriesisevenmoreillsuitedforthestudyoftheunderlyingirregulardynamics. ToshowthattheT(i)timeseriesinFigs. 42 (a,c)isindeedchaotic,weemploythedirectdynamicaltestfordeterministicchaosdevelopedbyGaoandZheng[ 21 22 ].The(k)curvesfortheT(i)timeseriesofFigs. 42 (a,c)areshowninFigs. 42 (b,d),wherethefourcurves,frombottomtotop,correspondtoshellsofsizes(2(i+1)=2;2i=2),i=14;15;16;17.wehavesimplychosenm=6andL=1.VerysimilarcurveshavebeenobtainedforotherchoicesofmandL.WeobserveaverywelldenedlinearcommonenvelopeatthelowerleftcornerofFigs. 42 (b,d).TheexistenceofacommonenvelopeguaranteesthatarobustpositiveLyapunovexponentwillbeobtainednomatterwhichshellisusedinthecomputation.Hence,thetwoT(i)timeseriesareindeedchaotic. 64
PAGE 65
TwoexamplesofchaoticT(i)timeseriesobtainedwithparameters(a)C=0:1,delayd=10ms,B=10,andthenumberofTCPstreamsN=17and(c)C=0:1,d=10ms,B=12,andN=19.The(k)curvesfor(a)and(c)areshownin(b)and(d),respectively. WehaveremarkedthattheprobabilitydistributionsfortheT(i)timeseriesofregularmotionsarecomprisedofafewdiscretepeaks.WhatarethedistributionsforthechaoticT(i)timeseriessuchasshowninFigs. 42 (a,c)?Theyarepowerlawlike,P(Tt)t,foralmosttwoordersofmagnitudeint,asshowninFigs. 43 (a,b).Theexponentisnotauniversalquantity,however. Wenowdevelopasimple1DmaptodescribetheoperationofTCP.TCPreliesontwomechanismstosetitstransmissionrate:owcontrolandcongestioncontrol. 65
PAGE 66
Complementarycumulativedistributionfunctions(CCDFs)forthetwochaoticT(i)timeseriesofFigs. 42 (a,c). Flowcontrolensuresthatthesendersendsnomorethanthesizeofthereceiver'slastadvertisedowcontrolwindowbasedonitsavailablebuersizewhilecongestioncontrolensuresthatthesenderdoesnotunfairlyoverrunthenetwork'savailablebandwidth.TCPimplementsthesemechanismsviaaowcontrolwindow(fwnd)advertisedbythereceivertothesenderandacongestioncontrolwindow(cwnd)adaptedbasedoninferringthestateofthenetwork.Specically,TCPcalculatesaneectivewindow(ewnd),whereewnd=min(fwnd;cwnd),andthensendsdataatarateewnd=RTT,whereRTTistheroundtriptimeoftheconnection.RecentlyRaoandChua[ 85 ]developedananalyticmodeldescribingTCPdynamics.ByassumingfwndbeingxedandtheTCP'sslowstartonlycontributingtotransientbehavior,RaoandChuahavedevelopedthewupdatemapM:[1;Wmax]![1;Wmax], 66
PAGE 67
Theabovemap(Eq. 4{1 )canbesimplied[ 90 ], Sincefwndisnotassumedtobeaconstant,themodiedmapismoregeneralthantheoriginalRaoandChua'smap.ToapplytheabovemaptothestudyofaTCPowcompetingwithN1otherTCPows,onemaylumptheeectoftheN1otherTCPowsasabackgroundtrac.Hence,Wmaxisdeterminedbycwndandthebackgroundtrac,andtypicallyisacomplicatedfunctionoftime,noticingthattheavailablebandwidthofabottlenecklinkcanvarywithtimeconsiderablyduetothedynamicnatureofbackgroundtrac.WehavecarriedoutsimulationstudiesofEq. 4{2 byassumingWmaxtobeperiodicandquasiperiodic,andhaveindeedobservedchaos. WearenowreadytounderstandwhychaoscannotoccurwhenthenumberofcompetingTCPowsissmall,suchas2.Toseethis,letusassumethatwiandwi+wiaswellaswi+1andwi+1+wi+1areallsmallerthanWmax.Itistheneasytoseethatjwi+1jjwij.Thatis,smalldisturbancedecays.ThismeansduringtheadditiveincreasephaseofTCP,nearbytrajectoriescontractinsteadofdiverge.Hence,thedynamicsarestable.InorderforthedynamicstobeunstablesothattheLyapunovexponentispositive,thetransitionfromtheadditiveincreasephasetothemultiplicativedecreasephasehastobefast.ThiscanbetrueonlywhenthenumberofcompetingTCPstreamsislarge.Wethusconcludethatthenetworkscenariosconsideredin[ 84 ]cannotgeneratechaoticTCPdynamics. HowrelevantisoursimulationresulttothedynamicsoftherealInternet?AretheirregularT(i)timeseriesshowninFigs. 42 (a,c)anartifactofTahoeversionofTCP 67
PAGE 68
ThetimeseriesT(i)extractedfromaW(i)timeseriescollectedusingnet100instrumentsoverORNLLSUconnection.(a)thetimeseriesT(i),(b)thecomplementarycumulativedistributionfunction(CCDF)fortheT(i)timeseries. weusedormoregenerallyofthens2simulator?Tondananswer,wecollectedW(i)measurementsbetweenORNLandLouisianaStateUniversityatmillisecondresolutionusingthenet100instrumentsandcomputedT(i)asshowninFig. 44 (a),whoseproleisqualitativelyquitesimilartoFig. 42 (a).Infact,italsofollowsa(truncated)powerlaw,asshowninFig. 44 (b).TheexponentforthepowerlawpartisevensmallerthanthetimeseriesofFig. 42 .ThetruncationinthepowerlawcouldbeduetotheshortnessoftheT(i)timeseries.Aswehavepointedout,theT(i)timeseriesisrelatedtoRTT.CottrellandBullothavebeenexperimentingwithmanyadvancedversionsofTCP,andalsoobservedRTTtimeseriesverysimilartothese.WealsohaveobservedfromRTTandlossdatameasuredongeographicallydispersedpathsontheInternet(witharesolutionof1sec)byresearchersattheSanDiegoSupercomputerCenterthattheprobabilitydistributionsforthetimeintervalbetweensuccessivelossburststypicallyfollowapowerlawlikedistribution,withtheexponentalsosmallerthan1.Thus,wehavegoodreasontobelievethattheobservedirregularT(i)timeseriesisnotan 68
PAGE 69
91 ],butvaryconsiderablyindurationsdeterminedbyspecicapplications.Inthenextsection,wewillanalyzetheactualdynamicsoftheInternetmeasurements. 92 ]withwhichTCPcompetesforbandwidthandrouterbuers.Roughlyspeaking,theinteractionbetweenthetwoisintermsoftheadditionstoTCPcongestionwindowsize,denotedbycwnd,inresponsetoacknowledgments,andmultiplicativedecreasesinresponsetoinferredlosses(ignoringtheinitialslowstartphase).Anyprotocol,howeversimpleitsdynamicsare,isgenerallyexpectedtoexhibitapparentlycomplicateddynamicsduetoitsinteractionwiththestochasticnetworktrac.TCPinparticularisshown(albeitinsimulation)toexhibitchaoticorchaoslikebehaviorevenwhenthecompetingtracismuchsimplersuchasasinglecompetingTCP[ 84 ]orUserDatagramProtocol(UDP)stream[ 85 ].Thedicultyistounderstandthedynamicswhenboththephenomenaareineect,andequallyimportantlytounderstandtheirimpactonactualInternetstreams.Historically,themethodsfromstandardchaosandstochastictheorieshavebeenunabletooermuchnewperspectiveonthetransportdynamics. 69
PAGE 70
21 22 37 ].OurmajorpurposeistoelucidatehowthedeterministicandstochasticcomponentsoftransportdynamicsinteractwitheachotherontheInternet.ItisofconsiderableinteresttonotethatrecentlytherehavebeenseveralimportantworksonTCPdynamicswiththepurposeofimprovingcongestioncontrol[ 93 { 96 ].Infact,therehasbeenconsiderableeortindevelopingnewversionsand/oralternativestoTCPsothatthenetworkdynamicscanbemorestable.Conventionalwaysofanalyzingthenetworkdynamicsareunabletoreadilydeterminewhethernewermethodsresultinstabletransportdynamicsorhelpindesigningsuchmethods.Ouranalysisshowstwoimportantfeaturesinthisdirection: (a) Randomnessisanintegralpartofthetransportdynamicsandmustbeexplicitlyhandled.Inparticular,thestepsizesutilizedforadjustingthecongestionparametersmustbesuitablyvaried(e.g.usingRobbinsMonroconditions)toensuretheeventualconvergence[ 97 ].ThisisnotthecaseforTCPwhichutilizesxedstepsizes. (b) Thechaoticdynamicsofthetransportprotocolsdohaveasignicantimpactonpracticaltransfers,andtheprotocoldesignmustbecognizantofitseects.Inparticular,itmightbeworthwhiletoinvestigateprotocolsthatdonotcontaindominantchaoticregimes,particularlyforremotecontrolandsteeringapplications. Thetraditionaltransportprotocolsarenotdesignedtoexplicitlyaddresstheabovetwoissues,butjustiablysosincetheiroriginalpurposeisdatatransportratherthannercontrolofdynamics. Intherestofthissection,weshallrstbrieydescribethecwnddatastudiedhere.Thenweshallbrieyexplaintheanalysisprocedureanddescribetypicalresultsoftheanalysis. WehavecollectedanumberofcwndtracesusingsingleandtwocompetingTCPstreams.ThisdatawascollectedontwodierentconnectionsfromOakRidgeNationalLaboratory(ORNL)toGeorgiaInstituteofTechnology(GaTech)andtoLouisianaState 70
PAGE 71
45 (ad)segmentsoffourdatasetsforORNLGaTechconnection.Powerspectralanalysisofthesedatadoesnotshowanydominantpeaks,andhence,thedynamicsarenotsimplyoscillatory.SinceourdatawasmeasuredontheInternetwith\live"backgroundtrac,itisapparentlymorecomplicatedandrealisticthanthetracesobtainedbynetworksimulation. Nextletusanalyzecwnddata,x(i);i=1;;n,usingtheconceptsoftimedependentexponent(k)curves[ 21 22 37 ].Ashasbeendiscussedearlier,werstneedtoemploytheembeddingtheorem[ 34 { 36 98 ]toconstructvectorsoftheform:Vi=[x(i);x(i+L);:::;x(i+(m1)L)],wheremistheembeddingdimensionandLthedelaytime.Theembeddingtheorem[ 34 { 36 ]statesthatwhentheembeddingdimensionmislargerthantwicetheboxcountingdimensionoftheattractor,thenthedynamicsoftheoriginalsystemcanbestudiedfromasinglescalartimeseries.Wenotethattheembeddingdimensionusedin[ 84 ]isonlytwo,whichhastobeconsiderednotlargeenough(thismaycallforacloserexaminationoftheirconclusions). Beforewemoveon,wepointoutafewinterestingfeaturesofthe(k)curves: (i) Fornoise,onlyforkuptotheembeddingwindowsize(m1)Lwill(k)increase.Thus,whenever(k)increasesforamuchlargerrangeofk,itisanindicationofnontrivialdeterministicstructureinthedata. (ii) Forperiodicsignals,(k)isessentiallyzeroforanyk. 71
PAGE 72
TimeseriesforthecongestionwindowsizecwndforORNLGaTechconnection. (iii) Forquasiperiodicsignals,(k)isperiodicwithanamplitudetypicallysmallerthan0.1,hence,forpracticalpurposes,(k)canbeconsideredverycloseto0. (iv) Whenachaoticsignaliscorruptedbynoise,thenthe(k)curvesbreakthemselvesawayfromthecommonenvelope.Thestrongerthenoiseis,themorethe(k)curvesbreakawaytilltheenvelopeisnotdenedatall.Thisfeaturehasactuallybeenusedtoestimatetheamountofbothmeasurementanddynamicnoise[ 99 ]. Nowwearereadytocomputeandunderstandthe(k)curvesforcwndtraces.Thesetof(k)curves,correspondingtoFig. 45 areplottedinFigs. 46 .Inthecomputations,3104pointsareused,andm=10;L=1.Thesevencurves,fromthebottomto 72
PAGE 73
Timedependentexponent(k)vs.kcurvesforcwnddatacorrespondingtoFig. 45 .Inthecomputations,3104pointsareused,andm=10;L=1.Curvesfromthebottomtotopcorrespondtoshellswithsizes(2(i+1)=2;2i=2),i=2;3;;8. top,correspondtoshellsofsizes(2(i+1)=2;2i=2),i=2;3;;8.Wemakethefollowingobservations: (i) Thedynamicsarecomplicatedandcannotbedescribedaseitherperiodicorquasiperiodicmotions,since(k)ismuchlargerthan0. (ii) Thedynamicscannotbecharacterizedaspuredeterministicchaos,sinceinnocasecanweobserveawelldenedlinearenvelope.Thustherandomcomponentofthedynamicsduetocompetingnetworktracisevidentandcannotsimplybeignored. 73
PAGE 74
Thedataisnotsimplynoisy,sinceotherwiseweshouldhaveobservedthat(k)isalmostatwhenk>(m1)L.Thus,thedeterministiccomponentofdynamicswhichisduetothetransportprotocolplaysanintegralroleandmustbecarefullystudied.Thefeatures(ii)and(iii)indicatethattheInternettransportdynamicscontainsbothchaoticandstochasticcomponents. (iv) Thereareconsiderabledierencesbetweenthedatawithonly1TCPsourceandwith2competingTCPsources.Inthelattercase,the(k)curvessharplyrisewhenkjustexceedstheembeddingwindowsize,(m1)L.Ontheotherhand,(k)forFig. 46 (a)withonlyoneTCPsourceincreasesmuchslowerwhenkjustexceeds(m1)L.Alsoimportantisthatthe(k)curvesinFigs. 46 (c,d)aremuchsmootherthanthoseinFigs. 46 (a,b).Hence,wecansaythatthedeterministiccomponentofthedynamicsismorevisiblewhentherearemorethan1competingTCPsources(alongthelinesof[ 84 ]). Sincetheincreasingpartofthe(k)curvesarenotverylinear,letusnextexamineifthecwndtracescanbebettercharacterizedbythemoregeneralizedconceptofpowerlawsensitivitytoinitialconditions(PSIC).Figure 47 showsthe(k)vs.lnkcurvesforORNLGaTechconnection.Veryinterestingly,wehavenowindeedobservedabunchofbetterdenedlinearlines,especiallyforsmallscales.Inparticular,thetransportdynamicswithmorethan1competingTCPsourcesarebetterdescribedbytheconceptofPSIC. Tosummarize,byanalyzinganumberofhighqualitycongestionwindowsizedatameasuredontheInternet,wehavefoundthatthetransportdynamicsarebestdescribedbytheconceptofPSIC,especiallyforsmallscales.Itisinterestingtofurtherexaminehowonemightsuppressthestochasticityofthenetworkbyexecutingmorecontrolsonthenetworkwhenmakingmeasurements,suchasusingalargenumberofcompetingTCPsourcestogetherwithaconstantUDPow.OuranalysismotivatesthatbothchaoticandstochasticaspectsbepaidacloseattentiontoindesigningInternetprotocolsthatarerequiredtoprovidethedesiredandtractabledynamics. 74
PAGE 75
Timedependentexponent(k)vs.lnkcurvesforcwnddatacorrespondingtoFig. 45 .Inthecomputations,3104pointsareused,andm=10;L=1.Curvesfromthebottomtotopcorrespondtoshellswithsizes(2(i+1)=2;2i=2),i=2;3;;8. 75
PAGE 76
Understandingthenatureofseaclutteriscrucialtothesuccessfulmodelingofseaclutteraswellastofacilitatetargetdetectionwithinseaclutter.Tothisend,animportantquestiontoaskiswhetherseaclutterisstochasticordeterministic.Sincethecomplicatedseacluttersignalsarefunctionsofcomplex(sometimesturbulent)wavemotionsontheseasurface,whilewavemotionsontheseasurfaceclearlyhavetheirowndynamicalfeaturesthatarenotreadilydescribedbysimplestatisticalfeatures,itisveryappealingtounderstandseaclutterbyconsideringsomeoftheirdynamicalfeatures.Inthepastdecade,Haykinetal.havecarriedoutanalysisofsomeseaclutterdatausingchaostheory[ 100 101 ],andconcludedthatseaclutterwasgeneratedbyanunderlyingchaoticprocess.Recently,theirconclusionhasbeenquestionedbyanumberofresearchers[ 102 { 108 ].Inparticular,Unsworthetal.[ 105 106 ]havedemonstratedthatthetwomaininvariantsusedbyHaykinetal.[ 100 101 ],namelythe\maximumlikelihoodofthecorrelationdimensionestimate"andthe\falsenearestneighbors"areproblematicintheanalysisofmeasuredseaclutterdata,sincebothinvariantsmayinterpretstochasticprocessesaschaos.Theyhavealsotriedanimprovedmethod,whichisbasedonthecorrelationintegralofGrassbergerandProcaccia[ 109 ]andhasbeenfoundeectiveindistinguishingstochasticprocessesfromchaos.Still,noevidenceofdeterminismorchaoshasbeenfoundinseaclutterdata. Toreconcileevergrowingevidenceofstochasticityinseaclutterwiththeirchaoshypothesis,recently,Haykinetal.[ 110 ]havesuggestedthatthenonchaoticfeatureofseacluttercouldbeduetomanytypesofnoisesourcesinthedata.Totestthispossibility,McDonaldandDamini[ 111 ]havetriedaseriesoflowpasslterstoremovenoise;butagaintheyhavefailedtondanychaoticfeatures.Furthermore,theyhavefoundthatthecommonlyusedchaoticinvariantmeasuresofcorrelationdimensionandLyapunovexponent,computedbyconventionalways,producesimilarresultsformeasuredsea 76
PAGE 77
Whiletheserecentstudieshighlysuggestthatseaclutterisunlikelytobetrulychaotic,anumberoffundamentalquestionsarestillunknown.Forexample,mostofthestudiesin[ 102 { 106 ]areconductedbycomparingmeasuredseaclutterdatawithsimulatedstochasticprocesses.Wecanask:canthenonchaoticnatureofseaclutterbedirectlydemonstratedwithoutresortingtosimulatedstochasticprocesses?Recognizingthatsimplelowpasslteringdoesnotcorrespondtoanydenitescalesinphasespace,canwedesignamoreeectivemethodtoseparatescalesinphasespaceandtotestwhetherseacluttercanbedecomposedassignalsplusnoise?Finally,willstudiesalongthislinebeofanyhelpfortargetdetectionwithinseaclutter? Inthischapter,weemploythedirectdynamicaltestfordeterministicchaosdiscussedinSec. 2.1 toanalyze280seaclutterdatameasuredundervariousseaandweatherconditions.Themethodoersamorestringentcriterionfordetectinglowdimensionalchaos,andcansimultaneouslymonitormotionsinphasespaceatdierentscales.However,nochaoticfeatureisobservedfromanyofthesedataatallthedierentscalesexamined.Butinterestingly,weshowthatscaledependentexponentcorrespondingtolargescaleappearstobeusefulfordistinguishingseaclutterdatawithandwithouttargets.Thissuggeststhatseacluttermaycontaininterestingdynamicfeaturesandthatthescaledependentexponentmaybeanimportantparameterfortargetdetectionwithinseaclutter.Furthermore,wendthatseacluttercanbeconvenientlycharacterizedbythenewconceptofpowerlawsensitivitytoinitialconditions(PSIC).Weshowthatforthepurposeofdetectingtargetswithinseaclutter,PSICoersamoreeectiveframework. Below,weshallrstbrieydescribetheseaclutterdata.ThenwestudytheseaclutterdatabyemployingthedirectdynamicaltestforlowdimensionalchaosdevelopedbyGaoandZheng[ 21 22 ].Andweshowthatthescaledependentexponentcorrespondingtolargescalecanbeusedtoeectivelydetectlowobservabletargetswithin 77
PAGE 78
51 .Thedistancebetweentwoadjacentrangebinswas15m.Oneorafewrangebins(say,Bi1,BiandBi+1)hitatarget,whichwasasphericalblockofstyrofoamofdiameter1m,wrappedwithwiremesh.Thelocationsofthethreetargetswerespeciedbytheirazimuthalangleanddistancetotheradar.Theywere(128degree,2660m),(130degree,5525m),and(170degree,2655m),respectively.Therangebinwherethetargetisstrongestislabeledastheprimarytargetbin.Duetodriftofthetarget,binsadjacenttotheprimarytargetbinmayalsohitthetarget.Theyarecalledsecondarytargetbins.Foreachrangebin,therewere217complexnumbers,sampledwithafrequencyof1000Hz.Amplitudedataoftwopolarizations,HH(horizontaltransmission,horizontalreception)andVV(verticaltransmission,verticalreception)wereanalyzed. Fig. 52 showstwoexamplesofthetypicalseaclutteramplitudedatawithoutandwithtarget.However,carefulexaminationoftheamplitudedataindicatesthat4datasetsareseverelyaectedbyclipping.ThiscanbereadilyobservedfromFigs. 53 (a,b).We 78
PAGE 79
Collectionofseaclutterdata Figure52. Typicalseaclutteramplitudedata(a)withoutand(b)withtarget. discardthose4datasets,andanalyzetheremaining10measurements,whichcontain280seacluttertimeseries. 21 22 ],toanalyzeseaclutterdata.Themethodoersamorestringentcriterionforlowdimensionalchaos,andcansimultaneouslymonitormotionsinphasespaceatdierentscales.Themethodhasfoundnumerousapplicationsinthe 79
PAGE 80
Twoshortsegmentsoftheamplitudeseaclutterdataseverelyaectedbyclipping. studyoftheeectsofnoiseondynamicalsystems[ 23 24 26 27 ],estimationofthestrengthofmeasurementnoiseinexperimentaldata[ 28 29 ],pathologicaltremors[ 30 ],shearthickeningsurfactantsolutions[ 31 ],diluteshearedaqueoussolutions[ 32 ],andserratedplasticows[ 33 ].Inparticular,thismethodwasusedbyGaoetal.[ 107 108 ]toanalyzeonesinglesetofseaclutterdata.Whilechaoswasnotobservedfromthatdataset,nodenitegeneralconclusioncouldbereached,duetolackoflargeamountofdataatthattime.Herewesystematicallystudy280seaclutterdatameasuredundervariousseaandweatherconditions,andexaminewhetheranychaoticfeaturescanbefoundfromtheseseaclutterdata. TheexplicitincorporationofscalesintheGaoandZheng'stest[ 21 22 ]enablesustosimultaneouslymonitormotionsinphasespaceatdierentscales.Wehavesystematicallyanalyzed280amplitudeseacluttertimeseriesmeasuredundervariousseaandweatherconditions.However,nochaoticfeaturehasbeenobservedfromanyofthesedataatallthescalesexamined.Typicalexamplesofthe(k)vs.kcurvesfortheseaclutteramplitudedatawithouttargetsareshowninFigs. 54 (a,c,e)andthecurvesforthedatawithtargetsshowninFigs. 54 (b,d,f),respectively.Wehavesimplychosenm=6and 80
PAGE 81
Examplesofthetimedependentexponent(k)vs.kcurvesfortheseaclutterdata(a,c,e)withoutand(b,d,f)withthetarget.Sixcurves,frombottomtotop,correspondtoshells(2(i+1)=2;2i=2)withi=13,14,15,16,17,and18.Thesamplingtimefortheseaclutterdatais1msec,andembeddingparametersarem=6,L=1.217datapointsareusedinthecomputation. 54 aregenericamongallthe280seaclutterdataanalyzedhere.Hence,wehavetoconcludethatnoneoftheseaclutterdataischaotic. 81
PAGE 82
112 ],lognormal[ 113 ],K[ 114 115 ],andcompoundGaussiandistributions[ 116 ],aswellasusingchaostheory[ 100 { 104 107 108 110 111 ],wavelets[ 117 ],neuralnetworks[ 118 119 ],waveletneuralnetcombinedapproaches[ 120 121 ],andtheconceptoffractaldimension[ 122 ]andfractalerror[ 123 124 ].However,nosimplemethodhasbeenfoundtorobustlydetectlowobservableobjectswithinseaclutter[ 125 ]. Inthissection,weexaminewhethertheLyapunovexponentestimatedbyconventionalmethodsandthescaledependentexponentmaybehelpfulfortargetdetectionwithinseaclutter. WerstcheckwhethertheLyapunovexponentestimatedbyconventionalmethodscanbeusedfordistinguishingseaclutterdatawithandwithouttargets.Aspointedoutearlier,theconventionalwayofestimatingtheLyapunovexponentistocompute(k)=kt,where(k)isdenedasinEq. 2{4 ,subjecttotheconditionsthatkXiXjk
PAGE 83
ThemodiedapproachforestimatingtheLyapunovexponentusestheleastsquaresttotherstfewsamplesofthe(k)curvewherethe(k)curveincreaseslinearlyorquasilinearly.Thisisdoneforallthe280seacluttertimeseries.TobetterappreciatethevariationoftheLyapunovexponentamongthe14rangebinsoftheseaclutterdata,wehaverstsubtractedtheparameterofeachbinbytheminimumofvaluesforthatmeasurement,andthennormalizedtheobtainedvaluesbyitsmaximum.Thevariationsoftheparameterswiththe14rangebinsforthe10HHmeasurementsareshowninFigs. 55 (a)(j),respectively,whereopencirclesdenotetherangebinswiththetarget,andasteroidsdenotethebinswithoutthetarget.Theprimarytargetbinisexplicitlyindicatedbyanarrow.Similarresultshavebeenobtainedforthe10VVmeasurements.WhileFigs. 55 (a)and(h)indicatethattheprimarytargetbincanbeseparatedfrombinswithoutthetarget,ingeneral,wehavetoconcludethattheLyapunovexponentestimatedbymethodsequivalentorsimilartoconventionalmeanscannotbeusedfordistinguishingseaclutterdatawithandwithouttargets. Nextletusexaminewhetherthescaledependentexponentcorrespondingtolargescalemaybeusefulfordetectingtargetswithinseaclutterdata.Byscaledependentexponent,wemeanthatifweusetheleastsquaresttothelinearlyorquasilinearlyincreasingpartofthe(k)curvesofdierentshells,theslopesofthoselinesdependonwhichshellsareusedforcomputingthe(k)curves.Inotherwords,ifweplottheexponentestimatedthiswayagainstthesizeoftheshell,thentheexponentvarieswiththeshellsizeorthescale.Thescaledependentexponenthasbeenusedinthestudyofnoiseinducedchaosinanopticallyinjectedsemiconductorlasermodeltoexaminehownoiseaectsdierentscalesofdynamicsystems[ 26 ].Wenowfocusontheexponentcorrespondingtolargeshellorscale.Figures 56 (a)(j)showthevariationsofthescaledependentexponentcorrespondingtolargescalewiththe14rangebinsforthe 83
PAGE 84
VariationsoftheLyapunovexponentestimatedbyconventionalmethodsvs.the14rangebinsforthe10HHmeasurements.Opencirclesdenotetherangebinswithtarget,while*denotethebinswithouttarget.Theprimarytargetbinisexplicitlyindicatedbyanarrow. same10HHmeasurementsasstudiedinFig. 55 .Weobservethattheprimarytargetbincanbeeasilyseparatedfromtherangebinswithoutthetarget,sincethescaledependentexponentfortheprimarytargetbinismuchlargerthanthoseforthebinswithoutthetarget. Itisthusclearthatthescaledependentexponentcorrespondingtolargescaleisveryusefulfordistinguishingseaclutterdatawithandwithouttargets.Thissuggeststhatseacluttermaycontaininterestingdynamicfeaturesandthatthescaledependentexponentmaybeanimportantparameterfortargetdetectionwithinseaclutter.This 84
PAGE 85
Variationsofthescaledependentexponentcorrespondingtolargescalevs.therangebinsforthe10HHmeasurements.Opencirclesdenotetherangebinswithtarget,while*denotethebinswithouttarget.Theprimarytargetbinisexplicitlyindicatedbyanarrow. exampleclearlysigniestheimportanceofincorporatingtheconceptofscaleinameasure.Infact,theconceptofscaleisonlyincorporatedinthetimedependentexponent(k)curves[ 21 22 ]inastaticmanner.InChapter 6 ,weshallseethatwhenameasuredynamicallyincorporatestheconceptofscale,itwillbecomemuchmorepowerful. Notethatonedicultyofusingthescaledependentexponentfortargetdetectionwithinseaclutteristhatwemayneedtochooseasuitablescaleforestimatingtheexponent.Thismakesthemethodnoteasytouse,anditishardtomakethemethodautomatic. 85
PAGE 86
57 (a)and(b),respectively,wherethecurvesdenotedbyasteroidsarefordatawithoutthetarget,whilethecurvesdenotedbyopencirclesarefordatawiththetarget.Weobservethatthecurvesarefairlylinearfortherstafewsamples.Alsonoticethattheslopesofthecurvesforthedatawiththetargetaremuchlargerthanthoseforthedatawithoutthetarget.Forconvenience,wedenotetheslopeofthecurvebytheparameter.Tobetterappreciatethevariationoftheparameterwiththerangebins,wenormalizeofeachbinbythemaximalvalueofthe14rangebinswithinthesinglemeasurement.Fig. 58 showsthevariationoftheparameterwiththe14rangebins.Itisclearthattheparametercanbeusedtodistinguishseaclutterdatawithandwithouttarget.Interestingly,thefeatureshowninFig. 58 isgenericallytrueforallthemeasurements.Itisworthpointingoutthatusuallythe(k)vs.lnkcurvesfordierentscalesarealmostparallel,thustheestimatedparameterisrelativelylessdependentonthescale.Thisisaverynicefeatureofthismethod. Letusexamineifarobustdetectorfordetectingtargetswithinseacluttercanbedevelopedbasedontheparameter.Wehavesystematicallystudied280timeseriesoftheseaclutterdatameasuredundervariousseaandweatherconditions.Tobetterappreciatethedetectionperformance,wehaverstonlyfocusedonbinswithprimarytargets,butomittedthosewithsecondarytargets,sincesometimesitishardtodeterminewhetherabinwithsecondarytargetreallyhitsatargetornot.Afteromittingtherangebindatawithsecondarytargets,thefrequenciesfortheparameterunderthetwohypotheses(thebinswithouttargetsandthosewithprimarytargets)forHHdatasets 86
PAGE 87
Examplesofthe(k)vs.lnkcurvesforrangebins(a)withoutand(b)withtarget.Opencirclesdenotetherangebinswithtarget,while*denotethebinswithouttarget. Figure58. Variationoftheparameterwiththe14rangebins.Opencirclesdenotetherangebinswithtarget,while*denotethebinswithouttarget. areshowninFig. 59 .WeobservethatthefrequenciescompletelyseparatefortheHHdatasets.Thismeansthedetectionaccuracycanbe100%. 87
PAGE 88
FrequenciesofthebinswithouttargetsandthebinswithprimarytargetsfortheHHdatasets. 88
PAGE 89
InChapter 1 ,wehaveemphasizedtheimportanceofdevelopingscaledependentmeasurestosimultaneouslycharacterizebehaviorsofcomplexmultiscaledsignalsonawiderangeofscales.Inthischapter,weshalldevelopaneectivealgorithmtocomputeanexcellentmultiscalemeasure,thescaledependentLyapunovexponent(SDLE),andshowthattheSDLEcanreadilyclassifyvarioustypesofcomplexmotions,includingtrulylowdimensionalchaos,noisychaos,noiseinducedchaos,processesdenedbypowerlawsensitivitytoinitialconditions(PSIC),andcomplexmotionswithchaoticbehavioronsmallscalesbutdiusivebehavioronlargescales.Finally,weshalldiscusshowtheSDLEcanhelpdetecthiddenfrequenciesinlargescaleorderlymotions. 126 ].ThealgorithmforcalculatingtheFSLEisverysimilartotheWolfetal'salgorithm[ 127 ].Itcomputestheaveragerfoldtimebymonitoringthedivergencebetweenareferencetrajectoryandaperturbedtrajectory.Todoso,itneedstodene\nearestneighbors",aswellasneedstoperform,fromtimetotime,arenormalizationwhenthedistancebetweenthereferenceandtheperturbedtrajectorybecomestoolarge.Suchaprocedurerequiresverylongtimeseries,andtherefore,isnotpractical.Tofacilitatederivationofafastalgorithmthatworksonshortdata,aswellastoeasediscussionofcontinuousbutnondierentiablestochasticprocesses,wecastthedenitionoftheSDLEasfollows. Consideranensembleoftrajectories.Denotetheinitialseparationbetweentwonearbytrajectoriesby0,andtheiraverageseparationattimetandt+tbytandt+t,respectively.Beingdenedinanaveragesense,tandt+tcanbereadilycomputedfromanyprocesses,eveniftheyarenondierentiable.Nextweexaminetherelationbetweent
PAGE 90
where(t)istheSDLE.Itisgivenby Givenatimeseriesdata,thesmallesttpossibleisthesamplingtime. ThedenitionoftheSDLEsuggestsasimpleensembleaveragebasedschemetocomputeit.Astraightforwardwaywouldbetondallthepairsofvectorsinthephasespacewiththeirdistanceapproximately,andthencalculatetheiraveragedistanceafteratimet.Thersthalfofthisdescriptionamountstointroducingashell(indexedask), whereVi;Vjarereconstructedvectors,k(theradiusoftheshell)andk(thewidthoftheshell)arearbitrarilychosensmalldistances.Suchashellmaybeconsideredasadierentialelementthatwouldfacilitatecomputationofconditionalprobability.Toexpeditecomputation,itisadvantageoustointroduceasequenceofshells,k=1;2;3;.Notethatthiscomputationalprocedureissimilartothatforcomputingthesocalledtimedependentexponent(TDE)curves[ 21 22 37 ]. Withalltheseshells,wecanthenmonitortheevolutionofallofthepairsofvectors(Vi;Vj)withinashellandtakeaverage.Wheneachshellisverythin,byassumingthattheorderofaveragingandtakinglogarithminEq. 6{2 canbeinterchanged,wehave wheretandtareintegersinunitofthesamplingtime,andtheanglebracketsdenoteaveragewithinashell.NotethatcontributionstotheSDLEataspecicscalefromdierentshellscanbecombined,withtheweightforeachshellbeingdeterminedbythe 90
PAGE 91
Intheaboveformulation,itisimplicitlyassumedthattheinitialseparation,kViVjk,alignswiththemostunstabledirectioninstantly.Forhighdimensionalsystems,thisisnottrue,especiallywhenthegrowthrateisnonuniformand/ortheeigenvectorsoftheJacobianarenonnormal.Fortunately,theproblemisnotasseriousasonemightbeconcerned,sinceourshellsarenotinnitesimal.WhencomputingtheTDEcurves[ 21 22 37 ],wehavefoundthatwhendicultiesarise,itisoftensucienttointroduceanadditionalcondition, whenndingpairsofvectorswithineachshell.SuchaschemealsoworkswellwhencomputingtheSDLE.Thismeansthat,aftertakingatimecomparabletotheembeddingwindow(m1)L,itwouldbesafetoassumethattheinitialseparationhasevolvedtothemostunstabledirectionofthemotion. Beforeproceedingon,wewishtoemphasizethemajordierencebetweenouralgorithmandthestandardmethodforcalculatingtheFSLE.Aswehavepointedout,tocomputetheFSLE,twotrajectories,oneasreference,anotherasperturbed,havetobedened.Thisrequireshugeamountsofdata.Ouralgorithmavoidsthisbyemployingtwocriticaloperationstofullyutilizeinformationaboutthetimeevolutionofthedata:(i)Thereferenceandtheperturbedtrajectoriesarereplacedbytimeevolutionofallpairsofvectorssatisfyingtheinequality( 6{5 )andfallingwithinashell,and(ii)introductionofasequenceofshellsensuresthatthenumberofpairsofvectorswithintheshellsislargewhiletheensembleaveragewithineachshelliswelldened.LetthenumberofpointsneededtocomputetheFSLEbystandardmethodsbeN.Thesetwooperationsimplythatthemethoddescribedhereonlyneedsaboutp 91
PAGE 92
2{5 .Figure 61 (a)showsvecurves,forthecasesofD=0;1;2;3;4.Thecomputationsaredonewith10000pointsandm=4;L=2.Weobserveafewinterestingfeatures: Nextweconsidernoiseinducedchaos.Toillustratetheidea,wefollow[ 23 ]andstudythenoisylogisticmap 92
PAGE 93
ScaledependentLyapunovexponent()curvesfor(a)thecleanandthenoisyLorenzsystem,and(b)thenoiseinducedchaosinthelogisticmap.Curvesfromdierentshellsaredesignatedbydierentsymbols. whereisthebifurcationparameter,andPnisaGaussianrandomvariablewithzeromeanandstandarddeviation.In[ 23 ],wereportedthatat=3:74and=0:002,noiseinducedchaosoccurs,andthoughtthatitmaybediculttodistinguishnoiseinducedchaosfromcleanchaos.InFig. 61 (b),wehaveplottedthe(t)forthisparticularnoiseinducedchaos.Thecomputationwasdonewithm=4;L=1and10000points.WeobservethatFig. 61 (b)isverysimilartothecurvesofnoisychaosplottedinFig. 61 (a).Hence,noiseinducedchaosissimilartonoisychaos,butdierentfromcleanchaos. Atthispoint,itisworthmakingtwocomments:(i)Onverysmallscales,theeectofmeasurementnoiseissimilartothatofdynamicnoise.(ii)The()curvesshowninFig. 61 arebasedonafairlysmallshell.ThecurvescomputedbasedonlargershellscollapseontherightpartofthecurvesshowninFig. 61 .Becauseofthis,forchaoticsystems,oneorafewsmallshellswouldbesucient.Ifonewishestoknowthebehaviorofoneversmallerscales,onehastouselongerandlongertimeseries. 93
PAGE 94
ScaledependentLyapunovexponent()curvefortheMackeyGlasssystem.Thecomputationwasdonewithm=5;L=1,and5000pointssampledwithatimeintervalof6. Finally,weconsidertheMackeyGlassdelaydierentialsystem[ 39 ],denedbyEq. 2{6 .ThesystemhastwopositiveLyapunovexponents,withthelargestLyapunovexponentcloseto0.007[ 21 ].HavingtwopositiveLyapunovexponentswhilethevalueofthelargestLyapunovexponentofthesystemisnotmuchgreaterthan0,onemightbeconcernedthatitmaybediculttocomputetheSDLEofthesystem.Thisisnotthecase.Infact,thissystemcanbeanalyzedasstraightforwardlyasotherdynamicalsystemsincludingtheHenonmapandtheRosslersystem.Anexampleofthe()curveisshowninFig. 62 ,wherewehavefollowed[ 21 ]andusedm=5;L=1,and5000pointssampledwithatimeintervalof6.Clearly,weobserveawelldenedplateau,withitsvaluecloseto0.007.ThisexampleillustratesthatwhencomputingtheSDLE,onedoesnotneedtobeveryconcernedaboutnonuniformgrowthrateinhighdimensionalsystems. Uptillnow,wehavefocusedonthepositiveportionof(t).Itturnsoutthatwhentislarge,()becomesoscillatory,withmeanabout0.Denotethecorrespondingscalesby 94
PAGE 95
3 ,wehaveseenthatpowerlawsensitivitytoinitialconditions(PSIC)providesacommonfoundationforchaostheoryandrandomfractaltheory.Here,weexaminewhethertheSDLEisabletocharacterizeprocessesdenedbyPSIC.Pleasingly,thisisindeedso.Below,wederiveasimpleequationrelatingtheqandqofPSICtotheSDLE. First,werecallthatPSICisdenedbyt=limjjV(0)jj!0jjV(t)jj jjV(0)jj=1+(1q)(1)qt1=(1q) 6{2 ,wendthat Whent!0,1+(1q)qt(1q)qt.SimplifyingEq. 6{7 ,weobtain Wenowconsiderthreecases: 95
PAGE 96
3{28 )and(1)q=H(Eq. 3{29 ).Therefore,(t)=H1 3{38 )and(1)q=1=(Eq. 3{39 ).Therefore,(t)=1 128 ]andstudythefollowingmap, where[xn]denotestheintegerpartofxn,tisanoiseuniformlydistributedintheinterval[1;1],isaparameterquantifyingthestrengthofnoise,andF(y)isgivenby ThemapF(y)isshowninFig. 63 asthedashedlines.ItgivesachaoticdynamicswithapositiveLyapunovexponentln(2+).Ontheotherhand,theterm[xn]introducesarandomwalkonintegergrids. Itturnsoutthissystemisveryeasytoanalyze.When=0:4,withonly5000pointsandm=2;L=1,wecanresolveboththechaoticbehavioronverysmallscales,andthenormaldiusivebehavior(withslope2)onlargescales.SeeFig. 64 (a). Wenowaskaquestion:Givenasmalldataset,whichtypeofbehavior,thechaoticorthediusive,isresolvedrst?Toanswerthis,wehavetriedtocomputetheSDLEwithonly500points.TheresultisshowninFig. 64 (b).Itisinterestingtoobservethatthechaoticbehaviorcanbewellresolvedbyonlyafewhundredpoints.However,the 96
PAGE 97
ThefunctionF(x)(Eq. 6{10 )for=0:4isshownasthedashedlines.ThefunctionG(x)(Eq. 6{11 )isanapproximationofF(x),obtainedusing40intervalsofslope0.Inthecaseofnoiseinducedchaosdiscussedinthepaper,G(x)isobtainedfromF(x)using104intervalsofslope0.9. ScaledependentLyapunovexponent()forthemodeldescribedbyEq. 6{9 .(a)5000pointswereused;forthenoisycase,=0:001.(b)500pointswereused. 97
PAGE 98
Wehavealsostudiedthenoisymap.TheresultingSDLEfor=0:001isshowninFig. 64 (a),assquares.Wehaveagainused5000points.WhilethebehavioroftheSDLEsuggestsnoisydynamics,with5000points,wearenotabletowellresolvetherelation()ln.Thisindicatesthatforthenoisymap,onverysmallscales,thedimensionisveryhigh. Map( 6{9 )canbemodiedtogiverisetoaninterestingsystemwithnoiseinducedchaos.ThiscanbedonebyreplacingthefunctionF(y)inmap( 6{9 )byG(y)toobtainthefollowingmap[ 128 ], wheretisanoiseuniformlydistributedintheinterval[1;1],isaparameterquantifyingthestrengthofnoise,andG(y)isapiecewiselinearfunctionwhichapproximatesF(y)ofEq. 6{10 .AnexampleofG(y)isshowninFig. 63 .Inournumericalsimulations,wehavefollowedCencinietal.[ 128 ]andused104intervalsofslope0.9toobtainG(y).WithsuchachoiceofG(y),intheabsenceofnoise,thetimeevolutiondescribedbythemap( 6{11 )isnonchaotic,sincethelargestLyapunovexponentln(0:9)isnegative.Withappropriatenoiselevel(e.g.,=104or103),theSDLEforthesystembecomesindistinguishabletothenoisySDLEshowninFig. 64 forthemap( 6{9 ).Havingadiusiveregimeonlargescales,thisisamorecomplicatednoiseinducedchaosthantheonewehavefoundfromthelogisticmap. Beforeproceedingon,wemakeacommentonthecomputationoftheSDLEfromdeterministicsystemswithnegativelargestLyapunovexponents,suchasthemap( 6{11 )withoutnoise.Atransientfreetimeseriesfromsuchsystemsisaconstanttimeseries.Therefore,thereisnoneedtocomputetheSDLEorothermetrics.Whenthetimeseries 98
PAGE 99
2{5 .InFigs. 65 (ac),wehaveshownthepowerspectraldensity(PSD)ofthex;y;zcomponentsofthesystem.WeobservethatthePSDofx(t)andy(t)aresimplybroad.However,thePSDofz(t)showsupaverysharpspectralpeak.Recallthatgeometrically,theLorenzattractorconsistsoftwoscrolls(seeFig. 66 ).ThesharpspectralpeakinthePSDofz(t)oftheLorenzsystemisduetothecircularmotionsalongeitherofthescrolls. TheaboveexampleillustratesthatifthedynamicsofasystemcontainsahiddenfrequencythatcannotberevealedbytheFouriertransformofameasuredvariable(say,x(t)),theninordertorevealthehiddenfrequency,onehastoembedx(t)toasuitablephasespace.Thisideahasledtothedevelopmentoftwointerestingmethodsforidentifyinghiddenfrequencies.OnemethodisproposedbyOrtega[ 129 130 ],bycomputingthetemporalevolutionofdensitymeasuresinthereconstructedphasespace.AnotherisproposedbyChernetal.[ 131 ],bytakingsingularvaluedecompositionoflocalneighbors.TheOrtega'smethodhasbeenappliedtoanexperimentaltimeseriesrecorded 99
PAGE 100
Powerspectraldensity(PSD)for(a)x(t),(b)y(t),(c)z(t),(d)(x)1(t),(e)(y)1(t),and(f)(z)1(t). fromafarinfraredlaserinachaoticstate.Thelaserdatasetcanbedownloadedfromthelink 67 (a)showsthelaserdataset.ThePSDofthedataisshowninFig. 67 (b),whereoneobservesasharppeakaround1.7MHz.Fig. 67 (c)showsthePSDofthedensitytimeseries,whereonenotesanadditionalspectralpeakaround37kHz.Thispeakisduetotheenvelopeofchaoticpulsations,whichisdiscernablefromFig. 67 (a). 100
PAGE 101
Lorenzattractor. ToappreciatethestrengthoftheadditionalspectralpeakinFig. 67 (c),wehavetonotethattheunitsofthePSDinbothFigs. 67 (b,c)arearbitrarythelargestpeakistreatedas1unit,aswasdonebyOrtega.Infact,theactualenergyoftheoriginallaserintensitytimeseriesismuchgreaterthanthatofthedensitytimeseries.Therefore,theadditionalspectralpeakof37kHzin(c)isreallyquiteweak.Inotherwords,althoughtheOrtega'smethodisabletorevealhiddenfrequencies,itisnotveryeective.Infact,Chernetal.havepointedoutthattheeectivenessoftheirmethodaswellastheOrtega'smethodtoexperimentaldataanalysisremainsuncertain. Oneofthemajordicultiesoftheabovetwomethodsisconceptualbothmethodsarebasedonverysmallscaleneighbors,whilethephenomenonitselfislargescale.Recognizingthis,onecanreadilyunderstandthattheconceptoflimitingscale,developedinSec. 6.2.1 ,providesanexcellentsolutiontotheproblem.Inotherwords,onecansimplytaketheFouriertransformofthetemporalevolutionofthelimitingscale,andexpectadditionalwelldenedspectralpeaks,ifthedynamicscontainhiddenfrequencies.To 101
PAGE 102
Hiddenfrequencyphenomenonoflaserintensitydata.(a)laserintensitydata,(b)thepowerspectraldensity(PSD)fortheoriginallaserintensitytimeseries,(c)thePSDforthedensitytimeseriesofthereconstructedphasespacefromthelaserdata,and(d)thePSDforthelimitingscaleofthelaserdata.In(c)and(d),theembeddingparametersarem=4;L=1. 102
PAGE 103
65 (df)thePSDofthelimitingscalesofx(t);y(t)andz(t)oftheLorenzsystem.Weobserveverywelldenedspectralpeaksinallthreecases.Infact,nowthezcomponentnolongerplaysamorespecialrolethanthexandycomponents.Themethodissimilarlyamazinglyeectiveinidentifyingthehiddenfrequencyfromthelaserdata,asshowninFig. 67 (d):wenownotonlyobserveanadditionalspectralpeakaround37KHzin(d),butalsothatthispeakisevenmoredominantthanthepeakaround1.7MHz.Thereasonfortheexchangeofthisdominanceisduetofairlylargesamplingtimealthough80nsissmall,itisonlyabletosampleabout8pointsineachoscillation.Whentheembeddingwindow,(m1)L,becomesgreaterthan4,theNyquistsamplingtheoremisviolatedinthereconstructedphasespace,thenthispeakmayevendiminish. Itisimportanttonotethatthelimitingscaleisnotaectedmuchbyeithermeasurementnoiseordynamicnoise.Therefore,onecanexpectthatthehiddenfrequenciesrevealedbylimitingscaleswillnotbeaectedmuchbynoiseeither. 103
PAGE 104
Inthiswork,wehavedevelopedanewtheoreticalframework,powerlawsensitivitytoinitialconditions(PSIC),toprovidechaostheoryandrandomfractaltheoryacommonfoundation.Forpracticaluses,wehavedevelopedaneasilyimplementableprocedureforcomputationallyexaminingPSICfromascalartimeseries.Wehavefoundthatdierenttypeofmotions,suchaschaos,edgeofchaos,1=fprocesseswithlongrangecorrelationsandLevyprocesses,arereadilycharacterizedbytheframeworkofPSICindierentrangeofq.Specically,trulychaoticmotionsarecharacterizedbyq=1;edgeofchaosisusuallycharacterizedbyaqvaluelargerthan0butsmallerthan1;and1=fprocessesandLevyprocessesarecharacterizedbyqintherangeof
PAGE 105
[1] P.Fairley,\Theunrulypowergrid,"IEEESpectrum,vol.41,pp.2227,2004. [2] C.Tsallis,\PossiblegeneralizationofBoltzmannGibbsstatistics,"JStatPhys,vol.52,pp.479487,1988. [3] A.R.PlastinoandA.Plastino,\Stellarpolytropesandtsallisentropy,"Phys.Lett.A,vol.174,pp.384386,1993. [4] A.Lavagno,G.Kaniadakis,M.RegoMonteiro,P.Quarati,andC.Tsallis,\Nonextensivethermostatisticalapproachofthepeculiarvelocityfunctionofgalaxyclusters,"Astrophys.Lett.Commun.,vol.35,pp.449455,1998. [5] M.AntoniandS.Ruo,\ClusteringandrelaxationinHamiltonianlongrangedynamics,"Phys.Rev.E,vol.52,pp.23612374,1995. [6] V.Latora,A.Rapisarda,andC.Tsallis,\NonGaussianequilibriuminalongrangeHamiltoniansystem,"Phys.Rev.E,vol.64,pp.056134,2001. [7] M.L.LyraandC.Tsallis,\Nonextensivityandmultifractalityinlowdimensionaldissipativesystems,"Phys.Rev.Lett.,vol.80,pp.5356,1998. [8] T.ArimitsuandN.Arimitsu,\Tsallisstatisticsandfullydevelopedturbulence,"J.Phys.A,vol.33,pp.L235L241,2000. [9] C.Beck,\Applicationofgeneralizedthermostatisticstofullydevelopedturbulence,"PhysicaA,vol.277,pp.115123,2000. [10] C.Beck,G.S.Lewis,andH.L.Swinney,\MeasuringnonextensitivityparametersinaturbulentCouetteTaylorow,"Phys.Rev.E,vol.63,pp.035303,2001. [11] C.Tsallis,A.R.Plastino,andW.M.Zheng,\PowerlawsensitivitytoinitialconditionsNewentropicrepresentation,"Chaos,Solitons,&Fractals,vol.8,pp.885891,1997. [12] U.M.S.Costa,M.L.Lyra,A.R.Plastino,andC.Tsallis,\Powerlawsensitivitytoinitialconditionswithinalogisticlikefamilyofmaps:Fractalityandnonextensivity,"Phys.Rev.E,vol.56,245250,1997. [13] V.Latora,M.Baranger,A.Rapisarda,andC.Tsallis,\Therateofentropyincreaseattheedgeofchaos,"Phys.Lett.A,vol.273,pp.97103,2000. [14] E.P.Borges,C.Tsallis,G.F.J.Ananos,andP.M.C.deOliveira,Phys.Rev.Lett.,vol.89,pp.254103,2002. 105
PAGE 106
U.Tirnakli,C.Tsallis,andM.L.Lyra,\Circularlikemaps:sensitivitytotheinitialconditions,multifractalityandnonextensivity,"Eur.Phys.J.B,vol.11,pp.309315,1999. [16] U.Tirnakli,\Asymmetricunimodalmaps:Someresultsfromqgeneralizedbitcumulants,"Phys.Rev.E,vol.62,pp.78577860,2000. [17] U.Tirnakli,G.F.J.Ananos,andC.Tsallis,\GeneralizationoftheKolmogorovSinaientropy:logisticlikeandgeneralizedcosinemapsatthechaosthreshold,"Phys.Lett.A,vol.289,pp.5158,2001. [18] U.Tirnakli,\Dissipativemapsatthechaosthreshold:numericalresultsforthesinglesitemap,"PhysicaA,vol.305,pp.119123,2001. [19] U.Tirnakli,C.Tsallis,andM.L.Lyra,\Asymmetricunimodalmapsattheedgeofchaos,"Phys.Rev.E,vol.65,pp.036207,2002. [20] U.Tirnakli,\Twodimensionalmapsattheedgeofchaos:NumericalresultsfortheHenonmap,"Phys.Rev.E,vol.66,pp.066212,2002. [21] J.B.GaoandZ.M.Zheng,\Directdynamicaltestfordeterministicchaosandoptimalembeddingofachaotictimeseries,"Phys.Rev.E,vol.49,pp.38073814,1994. [22] J.B.GaoandZ.M.Zheng,\Directdynamicaltestfordeterministicchaos," [23] J.B.Gao,C.C.Chen,S.K.Hwang,andJ.M.Liu,\Noiseinducedchaos,"Int.J.Mod.Phys.B,vol.13,pp.32833305,1999. [24] J.B.Gao,S.K.Hwang,andJ.M.Liu,\Whencannoiseinducechaos?"Phys.Rev.Lett.,vol.82,pp.11321135,1999. [25] J.B.Gao,S.K.Hwang,andJ.M.Liu,\Eectsofintrinsicspontaneousemissionnoiseonthenonlineardynamicsofanopticallyinjectedsemiconductorlaser,"Phys.Rev.A,vol.59,pp.15821585,1999. [26] S.K.Hwang,J.B.Gao,andJ.M.Liu,\Noiseinducedchaosinanopticallyinjectedsemiconductorlaser,"Phys.Rev.E,vol.61,pp.51625170,2000. [27] J.B.Gao,W.W.Tung,andN.Rao,\NoiseInducedHopfBifurcationtypeSequenceandTransitiontoChaosintheLorenzEquations,"Phys.Rev.Lett.,vol.89,pp.254101,2002. [28] J.Hu,J.B.Gao,andK.D.White,\Estimatingmeasurementnoiseinatimeseriesbyexploitingnonstationarity,"Chaos,Solitons,&Fractals,vol.22,pp.807819,2004. 106
PAGE 107
C.J.Cellucci,A.M.Albano,P.E.Rapp,R.A.Pittenger,R.C.Josiassen,\Detectingnoiseinatimeseries,"Chaos,vol.7,pp.414422,1997. [30] J.B.GaoandW.W.Tung,\Pathologicaltremorsasdiusionalprocesses,"BiologicalCybernetics,vol.86,pp.263270,2002. [31] R.BandyopadhyayandA.K.Sood,\Chaoticdynamicsinshearthickeningsurfactantsolutions,"Europhys.Lett.,vol.56,pp.447453,2001. [32] R.Bandyopadhyay,G.Basappa,andA.K.Sood,\ObservationofchaoticdynamicsindiluteshearedaqueoussolutionsofCTAT,"Phys.Rev.Lett.,vol.84,pp.20222025,2000. [33] S.Venkadesan,M.C.Valsakumar,K.P.N.Murthy,andS.Rajasekar,\Evidenceforchaosinanexperimentaltimeseriesfromserratedplasticow,"Phys.Rev.E,vol.54,pp.611616,1996. [34] N.H.Packard,J.P.Crutcheld,J.D.Farmer,andR.S.Shaw,\Geometryfromatimeseries,"Phys.Rev.Lett.,vol.45,pp.712716,1980. [35] F.Takens,inDynamicalsystemsandturbulence,LectureNotesinMathematics,vol.898,pp.366,1981,editedbyD.A.RandandL.S.Young,SpringerVerlag,Berlin. [36] T.Sauer,J.A.Yorke,andM.Casdagli,\Embedology,"J.Stat.Phys.,vol.65,pp.579616,1991. [37] J.B.GaoandZ.M.Zheng,\Localexponentialdivergenceplotandoptimalembeddingofachaotictimeseries,"Phys.Lett.A,vol.181,pp.153158,1993. [38] S.Sato,M.Sano,andY.Sawada,\PracticalmethodsofmeasuringthegeneralizeddimensionandthelargestLyapunovexponentinhighdimensionalchaoticsystems,"Prog.Theor.Phys.,vol.77,pp.15,1987. [39] M.C.MackeyandL.Glass,\Oscillationandchaosinphysiologicalcontrolsystems,"Science,vol.197,pp.287288,1977. [40] J.B.Gao,W.W.Tung,Y.H.Cao,J.Hu,andY.Qi,\Powerlawsensitivitytoinitialconditionsinatimeserieswithapplicationstoepilepticseizuredetection,"PhysicaA,vol.353,pp.613624,2005. [41] J.P.EckmannandD.Ruelle,\Ergodictheoryofchaosandstrangeattractors,"Rev.Mod.Phys.,vol.57,pp.617656,1985. [42] C.Beck,\DynamicalFoundationsofNonextensiveStatisticalMechanics,"Phys.Rev.Lett.,vol.87,pp.180601,2001. [43] B.B.Mandelbrot,TheFractalGeometryofNature,SanFrancisco,Freeman,1982. 107
PAGE 108
D.Saupe,\Algorithmsforrandomfractals",in:H.Peitgen,D.Saupe(Eds.),TheScienceofFractalImages,SpringerVerlag,Berlin,1988,pp.71113. [45] W.H.Press,\Flickernoisesinastronomyandelsewhere",Commentson Astrophysics,vol.7.pp.103119,1978. [46] P.Bak,HowNatureWorks:theScienceofSelfOrganizedCriticality,Copernicus,1996. [47] G.M.Wornell,Signalprocessingwithfractals:awaveletbasedapproach,PrenticeHall,1996. [48] W.E.Leland,M.S.Taqqu,W.Willinger,andD.V.Wilson,\OntheselfsimilarnatureofEthernettrac(extendedversion)",IEEE/ACMTrans.onNetworking,vol.2,115,1994. [49] J.Beran,R.Sherman,M.S.Taqqu,andW.Willinger,\Longrangedependenceinvariablebitratevideotrac",IEEETrans.onCommun.,vol.43,pp.15661579,1995. [50] V.PaxsonandS.Floyd,\WideAreaTracThefailureofPoissonmodeling,"IEEE/ACMTrans.onNetworking,vol.3,pp.226244,1995. [51] W.LiandK.Kaneko,\LongRangeCorrelationandPartial1/fSpectruminaNonCodingDNASequence",Europhys.Lett.,vol.17,pp.655660,1992. [52] R.Voss,\Evolutionoflongrangefractalcorrelationsand1/fnoiseinDNAbasesequences",Phys.Rev.Lett.,vol.68,pp.38053808,1992. [53] C.K.Peng,S.V.Buldyrev,A.L.Goldberger,S.Havlin,F.Sciortino,M.SimonsandH.E.Stanley,\LongRangeCorrelationsinNucleotideSequences,"Nature,vol.356,pp.168171,1992. [54] D.L.Gilden,T.Thornton,andM.W.Mallon,\1/fnoiseinhumancognition",Science,vol.267,pp.18371839,1995. [55] Y.H.Zhou,J.B.Gao,K.D.White,I.Merk,andK.Yao,\PerceptualDominanceTimeDistributionsinMultistableVisualPerception,"BiologicalCybernetics,vol.90,pp.256263,2004. [56] J.B.Gao,V.A.Billock,I.Merk,W.W.Tung,K.D.White,J.G.Harris,andV.P.Roychowdhury,\Inertiaandmemoryinambiguousvisualperception,"CognitiveProcessing,vol.7,pp.105112,2006. [57] Y.Chen,M.Ding,andJ.A.S.Kelso,\LongMemoryProcesses(1/falphaType)inHumanCoordination",Phys.Rev.Lett.,vol.79,pp.45014504,1997. 108
PAGE 109
J.J.CollinsandC.J.DeLuca,\RandomWalkingduringQuietStanding,"Phys.Rev.Lett.,vol.73,pp.764767,1994. [59] V.A.Billock,\Neuralacclimationto1/fspatialfrequencyspectrainnaturalimagestransducedbythehumanvisualsystem,"PhysicaD,vol.137,pp.379391,2000. [60] V.A.Billock,G.C.deGuzman,andJ.A.S.Kelso,\Fractaltimeand1/fspectraindynamicimagesandhumanvision",PhysicaD,vol.148,pp.136146,2001. [61] M.Wolf,\1/fnoiseinthedistributionofprimes,"PhysicaA,vol.241,pp.493499,1997. [62] J.B.Gao,YinheCao,andJaeMinLee,\PrincipalComponentAnalysisof1=fNoise,"Phys.Lett.A,vol.314,pp.392400,2003. [63] B.B.MandelbrotandV.Ness,\FractionalBrownianmotions,fractionalnoisesandapplications,"SIAMRev.,vol.10,pp.422437,1968. [64] M.S.Taqqu,W.Willinger,andR.Sherman,Proofofafundamentalresultinselfsimilartracmodeling,ACMSIGCOMMComputerCommunicationReview,vol.27,pp.523,1997. [65] B.B.Mandelbrot,FractalsandScalinginFinance,NewYork:Springer,1997. [66] T.Solomon,E.Weeks,andH.Swinney,\ObservationofAnomalousDiusionandLevyFlightsinaTwoDimensionalRotatingFlow,"Phys.Rev.Lett.,vol.71,pp.39753979,1993. [67] T.Geisel,J.Nierwetberg,andA.Zacherl,\AcceleratedDiusioninJosephsonJunctionsandRelatedChaoticSystems,"Phys.Rev.Lett.,vol.54,pp.616620,1985. [68] S.Stapf,R.Kimmich,andR.Seitter,\ProtonandDeuteronFieldcyclingNMRrelaxometryofLiquidsinPorousGlasses:EvidenceofLevywalkStatistics,"Phys.Rev.Lett.,vol.75,pp.28552859,1993. [69] G.M.Viswanathan,V.Afanasyev,S.V.Buldyrev,E.J.Murphy,P.A.Prince,andH.E.Stanley,\Levyightsearchpatternsofwanderingalbatrosses,"Nature,vol.381,pp.413415,1996. [70] G.M.Viswanathan,S.V.Buldyrev,S.Havlin,M.G.E.daLuz,E.P.Raposo,andH.E.Stanley,\Optimizingthesuccessofrandomsearches,"Nature,vol.401,pp.911914,1999. [71] J.B.Gao,J.Hu,W.W.Tung,Y.H.Cao,N.Sarshar,andV.P.Roychowdhury,\Assessmentoflongrangecorrelationintimeseries:Howtoavoidpitfalls,"Phys.Rev.E,vol.73,pp.016117,2006. 109
PAGE 110
C.K.Peng,S.V.Buldyrev,S.Havlin,M.Simons,H.E.Stanley,andA.L.Goldberger,\Mosaicorganizationofdnanucleotides,"Phys.Rev.E,vol.49,pp.16851689,1994. [73] A.JanickiandA.Weron,\Canoneseestablevariablesandprocesses,"StatisticalScience,vol.9,pp.109126,1994. [74] M.Kanter,\Stabledensitiesunderchangeofscaleandtotalvariationinequalities,"AnnalsofProbability,vol.3,pp.697707,1975. [75] J.M.Chambers,C.L.Mallows,B.W.Stuck,\Methodforsimulatingstablerandomvariables,"JournaloftheAmericanStatisticaLAssociation,vol.71,pp.340344,1976. [76] S.ResnickandG.Samorodnitsky:Fluidqueues,on/oprocesses,andteletracmodelingwithhighlyvariableandcorrelatedinputs,inSelfSimilarTracandPerformanceEvaluation,eds,K.ParkandW.Willinger,JohnWiley&Sons(2000),pp.171192. [77] J.B.GaoandI.Rubin,\MultifractalmodelingofcountingprocessesofLongRangeDependentnetworkTrac,"ComputerCommunications,vol.24,pp.14001410,2001;\MultiplicativeMultifractalModelingofLongRangeDependent(LRD)TracinComputerCommunicationsNetworks,"J.NonlinearAnalysis,vol.47,pp.57655774,2001;\MultiplicativemultifractalmodelingofLongRangeDependentnetworktrac,"Int.J.Comm.Systems,vol.14,pp.783801,2001. [78] R.AlbertandA.L.Barabasi,\Statisticalmechanicsofcomplexnetworks,"Rev.ModernPhys.,vol.74,pp.4797,2002. [79] C.Labovitz,G.R.Malan,andF.Jahanian,\Internetroutinginstability,"IEEE/ACMTransactionsonNetworking,vol.6,pp.515528,1998. [80] J.C.R.Bennett,C.Partridge,andN.Shectman,\Packetreorderingisnotpathologicalnetworkbehavior,"IEEE/ACMTransactionsonNetworking,vol.7,pp.789798,1999. [81] A.ErramilliandL.J.Forys,\TracsynchronizationEectsinTeletracSystems,"Proc.ITC13,Copenhagen,1991. [82] R.Srikant,TheMathematicsofInternetCongestionControl,Birkhauser,2004. [83] C.Li,G.Chen,X.Liao,andJ.Yu,\Experimentalqueueinganalysiswithlongrangedependentpackettrac,"Chaos,Solitons,&Fractals,vol.19,pp.853868,2004. [84] A.VeresandM.Boda,\ThechaoticnatureofTCPcongestioncontrol,"Proc.ofInfocom,Piscataway,NJ,pp.17151723,2000. 110
PAGE 111
N.S.V.RaoandL.O.Chua,\Ondynamicsofnetworktransportprotocols,"Proc.WorkshoponSignalProcessing,Communications,ChaosandSystems,2002. [86] P.Ranjan,E.H.Abed,andR.J.La,\NonlinearinstabilitiesinTCPRED,"IEEE/ACMTransactionsonNetworking,vol.12,pp.10791092,2004. [87] N.S.V.Rao,J.GaoandL.O.Chua,\Ondynamicsofnetworktransportprotocols,"inComplexDynamicsinCommunicationNetworks,editedbyL.KocarevandG.Vattay,2004. [88] J.GaoandN.S.V.Rao,\ComplicatedDynamicsofInternetTransportProtocols,"IEEECommun.Lett.,vol.9,pp.46,2005. [89] D.R.Figueiredo,B.Liu,V.Misra,andD.Towsley,\OntheAutocorrelationStructureofTCPTrac,"ComputerNetworks,vol.40,pp.339361,2002;SpecialIssueon\AdvancesinModelingandEngineeringofLongRangeDependentTrac",ComputerSci.Tech.Report0055,Univ.ofMassachusetts,Nov.2002. [90] J.B.GaoandN.S.V.Rao,\SynchronizedOscillations,QuasiPeriodicity,andChaosinTCP,"Dept.ofElectricalandComputerEngineering,UniversityofFlorida,Technicalreport. [91] S.FloydandE.Kohler,inFirstWorkshoponHotTopicsinNetworks,2829October2002,Princeton,NewJersey,USA. [92] K.ParkandW.Willinger,SelfSimilarNetworkTracandPerformance Evaluation,Wiley,NewYork,2000. [93] V.FiroiuandM.Borden,\Astudyofactivequeuemanagementforcongestioncontrol,"Proc.ofInfocom,vol.3,pp.14351444,2000. [94] V.Misra,W.Gong,andD.Towsley,\Astudyofactivequeueforcongestioncontrol,"Proc.ofSigcomm,2000. [95] C.Hollot,V.Misra,D.Towsley,andW.Gong,\AcontroltheoreticanalysisofRED,"Proc.ofInfocom,vol.3,pp.15101519,2001. [96] P.Kuusela,P.Lassilaand,andJ.Virtamo,\StabilityofTCPREDcongestioncontrol,"inProc.ofITC17,pp.655666,Elsevier. [97] N.S.V.Rao,Q.WuandS.S.Iyengar,\Onthroughputstabilizatinofnetworktransport,"IEEECommun.Lett.,vol.8,pp.6668,2004. [98] T.Sauer,\ReconstructionofDynamicalSystemsfromInterspikeIntervals,"Phys.Rev.Lett.,vol.72,pp.38113814,1994. 111
PAGE 112
J.B.Gao,\Recognizingrandomnessinatimeseries,"PhysicaD,vol.106,pp.4956,1997. [100] S.HaykinandS.Puthusserypady,\Chaoticdynamicsofseaclutter,"Chaos,vol.7,pp.777802,1997. [101] S.Haykin,Chaoticdynamicsofseaclutter,JohnWiley,1999. [102] J.L.Noga,\BayesianStateSpaceModellingofSpatioTemporalNonGaussianRadarReturns,"Ph.Dthesis,CambridgeUniversity,1998. [103] M.Davies,\LookingforNonLinearitiesinSeaClutter,"IEERadarandSonarSignalProcessing,Peebles,Scotland,UK,July1998. [104] M.R.CowperandB.Mulgrew,\NonlinearProcessingofHighResolutionRadarSeaClutter,"Proc.IJCNN,vol.4,pp.26332638,1999. [105] C.P.Unsworth,M.R.Cowper,S.McLaughlin,andB.Mulgrew,\Falsedetectionofchaoticbehaviorinthestochasticcompoundkdistributionmodelofradarseaclutter,"Proc.10thIEEEWorkshoponSSAP,August2000,pp.296300. [106] C.P.Unsworth,M.R.Cowper,S.McLaughlin,andB.Mulgrew,\Reexaminingthenatureofseaclutter",IEEProc.RadarSonarNavig.,vol.149,pp.105114,2002. [107] J.B.GaoandK.Yao,\Multifractalfeaturesofseaclutter,"IEEERadarConference,LongBeach,CA,April,2002. [108] J.B.Gao,S.K.Hwang,H.F.Chen,Z.Kang,K.Yao,andJ.M.Liu,\Canseaclutterandindoorradiopropagationbemodeledasstrangeattractors?"The7thExperimentalChaosConference,SanDiego,USA,Augustpp.2529,2002. [109] P.GrassbergerandI.Procaccia,\Measuringthestrangenessofthestrangeattractor,"PhysicaD,vol.9D,pp.189208,1983. [110] S.Haykin,R.Bakker,andB.W.Currie,\Uncoveringnonlineardynamicsthecasestudyofseaclutter,"Proc.IEEE,vol.90,pp.860881,2002. [111] M.McDonaldandA.Damini,\Limitationsofnonlinearchaoticdynamicsinpredictingseaclutterreturns,"IEEProc.RadarSonarNavig.,vol.151,pp.105113,2004. [112] F.A.Fay,J.Clarke,andR.S.Peters,\Weibulldistributionappliedtoseaclutter,"Proc.IEEConf.Radar,London,U.K.,pp.101103,1977. [113] F.E.Nathanson,Radardesignprinciples,McGrawHill,pp.254256,1969. 112
PAGE 113
E.JakemanandP.N.Pusey,\AmodelfornonRayleighseaecho,"IEEETransAntennas&Propagation,vol.24,pp.806814,1976. [115] S.SayamaandM.Sekine,\Lognormal,logWeibullandKdistributedseaclutter,"IEICETrans.Commun.,vol.E85B,pp.13751381,2002. [116] F.Gini,A.Farina,andM.Montanari,\VectorsubspacedetectionincompoundGaussianclutter,PartII:performanceanalysis,"IEEETransactionsonAerospaceandElectronicSystems,vol.38,pp.13121323,2002. [117] G.DavidsonandH.D.Griths,\Waveletdetectionschemeforsmalltargetsinseaclutter,"ElectronicsLett.,vol.38,pp.11281130,2002. [118] T.BhattacharyaandS.Haykin,\NeuralNetworkbasedRadarDetectionforanOceanEnvironment,"IEEETrans.AerospaceandElectronicSystems,vol.33,pp.408420,1997. [119] N.Xie,H.Leung,andH.Chan,\Amultiplemodelpredictionapproachforseacluttermodeling,"IEEETran.GeosciRemote,vol.41,pp.14911502,2003. [120] C.P.Lin,M.Sano,S.Sayama,andM.Sekine,\Detectionofradartargetsembeddedinseaiceandseaclutterusingfractals,wavelets,andneuralnetworks,"IEICETrans.Commun.,vol.E83B,pp.19161929,2000. [121] C.P.Lin,M.Sano,S.Obi,S.Sayama,andM.Sekine,\DetectionofoilleakageinSARimagesusingwaveletfeatureextractorsandunsupervisedneuralclassiers,"IEICETrans.Commun.,vol.E83B,pp.19551962,2000. [122] T.Lo,H.Leung,J.LitvaandS.Haykin,\Fractalcharacterisationofseascatteredsignalsanddetectionofseasurfacetargets,"IEEProc.,vol.F140,pp.243250,1993. [123] C.P.Lin,M.Sano,andM.Sekine,\Detectionofradartargetsbymeansoffractalerror,"IEICETrans.Commun.,vol.E80B,pp.17411748,1997. [124] B.E.Cooper,D.L.Chenoweth,andJ.E.Selvage,\Fractalerrorfordetectingmanmadefeaturesinaerialimages,"ElectronicsLetters,vol.30,pp.554555,1994. [125] J.Hu,W.W.Tung,andJ.B.Gao,\Detectionoflowobservabletargetswithinseaclutterbystructurefunctionbasedmultifractalanalysis,"IEEETransactionsonAntennas&Propagation,vol.54,pp.135143,2006. [126] E.Aurell,G.Boetta,A.Crisanti,G.Paladin,andA.Vulpiani,\Predictabilityinthelarge:AnextensionoftheconceptofLyapunovexponent,"J.PhysicsA,vol.30,pp.126,1997. [127] A.Wolf,J.B.Swift,H.L.Swinney,andJ.A.Vastano,\DeterminingLyapunovexponentsfromatimeseries,"PhysicaD,vol.16,pp.285317,1985. 113
PAGE 114
M.Cencini,M.Falcioni,E.Olbrich,H.Kantz,andA.Vulpiani,\Chaosornoise:Dicultiesofadistinction,"Phys.Rev.E,vol.62,pp.427437,2000. [129] G.J.Ortega,\Anewmethodtodetecthiddenfrequenciesinchaotictimeseries,"Phys.Lett.A,vol.209,pp.351355,1995. [130] G.J.Ortega,\InvariantmeasuresasLagrangianvariables:Theirapplicationtotimeseriesanalysis,"Phys.Rev.Lett.,vol.77,pp.259262,1996. [131] J.L.Chern,J.Y.Ko,J.S.Lih,H.T.Su,andR.R.Hsu,\Recognizinghiddenfrequenciesinachaotictimeseries,"Phys.Lett.A,vol.238,pp.134140,1998. 114
PAGE 115
JingHureceivedtheB.S.andM.E.degreesfromtheDepartmentofElectronicsandInformationEngineeringfromHuazhongUniversityofScienceandTechnology(HUST),Wuhan,China,in2000and2002,respectively.InMay2007,ShewasawardedthePh.D.degreefromtheDepartmentofElectricalandComputerEngineeringattheUniversityofFlorida,Gainesville,Florida.Herresearchinterestsincludesignalandimageprocessing,informationscience,nonlineartimeseriesanalysis,biologicaldataanalysis,bioinformatics,andbiomedicalengineering. 115

