|
|
||||||||
1 Cardiovascular Division, Department of Internal Medicine, and Department of Molecular Physiology and Biological Physics, and Cardiovascular Research Center, University of Virginia Health Sciences Center, Charlottesville, Virginia 22908; and 2 Medical Automation Systems, Charlottesville, Virginia 22903
| |
ABSTRACT |
|---|
|
|
|---|
Entropy, as it relates to dynamical systems, is the rate of information production. Methods for estimation of the entropy of a system represented by a time series are not, however, well suited to analysis of the short and noisy data sets encountered in cardiovascular and other biological studies. Pincus introduced approximate entropy (ApEn), a set of measures of system complexity closely related to entropy, which is easily applied to clinical cardiovascular and other time series. ApEn statistics, however, lead to inconsistent results. We have developed a new and related complexity measure, sample entropy (SampEn), and have compared ApEn and SampEn by using them to analyze sets of random numbers with known probabilistic character. We have also evaluated cross-ApEn and cross-SampEn, which use cardiovascular data sets to measure the similarity of two distinct time series. SampEn agreed with theory much more closely than ApEn over a broad range of conditions. The improved accuracy of SampEn statistics should make them useful in the study of experimental clinical cardiovascular and other biological time series.
probability; nonlinear dynamics
| |
INTRODUCTION |
|---|
|
|
|---|
NONLINEAR DYNAMICAL ANALYSIS is a powerful approach to understanding biological systems. The calculations, however, usually require very long data sets that can be difficult or impossible to obtain. Pincus (21, 22) devised the theory and method for a measure of regularity closely related to the Kolmogorov entropy, the rate of generation of new information, that can be applied to the typically short and noisy time series of clinical data. This family of statistics, named approximate entropy (ApEn), is rooted in the work of Grassberger and Procaccia (10) and Eckmann and Ruelle (5) and has been widely applied in clinical cardiovascular studies (3, 6-8, 12-19, 24, 26, 29, 32-34).
The method examines time series for similar epochs: more frequent and more similar epochs lead to lower values of ApEn. Informally, given N points, the family of statistics ApEn(m, r, N ) is approximately equal to the negative average natural logarithm of the conditional probability that two sequences that are similar for m points remain similar, that is, within a tolerance r, at the next point. Thus a low value of ApEn reflects a high degree of regularity. Importantly, the ApEn algorithm counts each sequence as matching itself, a practice carried over from the work of Eckmann and Ruelle (5) to avoid the occurrence of ln (0) in the calculations. This step has led to discussion of the bias of ApEn (22, 23, 27). In practice, we find that this bias causes ApEn to lack two important expected properties. First, ApEn is heavily dependent on the record length and is uniformly lower than expected for short records. Second, it lacks relative consistency. That is, if ApEn of one data set is higher than that of another, it should, but does not, remain higher for all conditions tested (22). This shortcoming is particularly important, because ApEn has been repeatedly recommended as a relative measure for comparing data sets (22-24).
To reduce this bias, we have developed and characterized a new family of statistics, sample entropy (SampEn), that does not count self-matches. SampEn is derived from approaches developed by Grassberger and co-workers (2, 9-11). SampEn(m, r, N ) is precisely the negative natural logarithm of the conditional probability that two sequences similar for m points remain similar at the next point, where self-matches are not included in calculating the probability. Thus a lower value of SampEn also indicates more self-similarity in the time series. In addition to eliminating self-matches, the SampEn algorithm is simpler than the ApEn algorithm, requiring approximately one-half as much time to calculate. SampEn is largely independent of record length and displays relative consistency under circumstances where ApEn does not.
Cross-ApEn is a recently introduced technique for analyzing two related time series to measure the degree of their asynchrony (20, 28). Cross-ApEn is very similar to ApEn in design and intent, differing only in that it compares sequences from one series with those of the second. Because it does not compare a series with itself, bias from self-matches does not arise. A potential problem, however, remains in the necessity for each template to generate a defined, nonzero probability. Thus each template must find at least one match for m + 1 points, or a probability must be assigned to it according to a "correction" strategy. We tested the effect of two extremes of correction strategies on cross-ApEn analysis. We find that cross-ApEn analysis lacks relative consistency, and conclusions about relative synchrony of pairs of time series depend on the unguided selection of analysis schemes. Cross-SampEn, on the other hand, is defined as long as one template finds a match, and we find that cross-SampEn remains relatively consistent for conditions where cross-ApEn does not.
| |
THEORY |
|---|
|
|
|---|
ApEn reports on similarity in time series. We employ the terminology and notation of Grassberger and Procaccia (10), Eckmann and Ruelle (5), and Pincus (21) in describing techniques for estimating the Kolmogorov entropy of a process represented by a time series and the related statistics ApEn and SampEn. The parameters N, m, and r must be fixed for each calculation. N is the length of the time series, m is the length of sequences to be compared, and r is the tolerance for accepting matches. It is convenient to set the tolerance as r × SD, the standard deviation of the data set, allowing measurements on data sets with different amplitudes to be compared. Throughout this work, all time series have been normalized to have SD = 1.
We proceed as follows: For a time series of N points, {u( j): 1
j
N } forms the N
m + 1 vectors xm(i) for
{i|1
i
N
m
+ 1}, where xm(i) = {u(i + k):
0
k
m
1} is the vector of m
data points from u(i) to
u(i + m
1). The distance between two
such vectors is defined to be d[x(i),
x( j)] = max {|u(i + k)
u( j + k)|:0
k
m
1},
the maximum difference of their corresponding scalar components. Let
Bi be the number of vectors
xm( j) within r of
xm(i) and let
Ai be the number of vectors
xm + 1( j ) within
r of xm + 1(i).
Define the function Cmi(r) = (Bi)/(N
m + 1).
In calculating
Cmi(r), the vector xm(i) is called the
template, and an instance where a vector
xm( j) is within r
of it is called a template match.
Cmi(r) is the probability that any vector
xm( j) is within r
of xm(i). The function
m(r) = (N
m + 1)
1
N
m + 1i = 1 lnt[Cmi(r)]
is the average of the natural logarithms of the functions
Cmi(r). Eckmann and Ruelle (5) suggest approximating the entropy of the
underlying process as limr
0
limm
limN
[
m(r)
m + 1(r)]. Because of these
limits, this definition is not suited to the analysis of the finite and
noisy time series derived from experiments.
Pincus (21) saw that the calculation of
m(r)
m + 1(r) for fixed parameters
m, r, and N had intrinsic interest as a measure of regularity and complexity. He defines the related parameter ApEn(m, r) = limN
[
m(r)
m + 1(r)], which for finite data
sets is estimated by the statistic ApEn(m, r, N ) =
m(r)
m + 1(r). Algebraic manipulation
reveals that ApEn(m, r,N) = (N
m + 1)
1
N
m + 1i = 1 ln
[Cmi(r)]
(N
m)
1
N
mi = 1 ln[Cm + 1i(r)]. When N is large,
ApEn(m, r, N ) is approximately equal to
(N
m)
1
N
mi = 1 [
ln (Ai/Bi)],
the average over i of the negative natural logarithm of the
conditional probability that
d[xm + 1(i), xm + 1( j)]
r given that
d[xm(i), xm( j)]
r. The
difference between ApEn(m, r, N ) and this
average logarithmic probability is less than
(N
m + 1)
1
ln (N
m + 1), which is <0.05 for N
m + 1
90 and <0.02 for N
m + 1
283. They differ because there are N
m + 1 templates of length m, but only N
m
templates of length m + 1 (27). Thus, although the quantity
CmN
m + 1 (r) is defined,
Cm + 1N
m + 1 (r) is not, simply because the vector
um + 1(N
m + 1)
does not exist. For practical purposes, ApEn can be thought of as the
negative natural logarithm of the probability that sequences that are
close for m points remain close for an additional point. Because conditional probabilities lie between 0 and 1, the parameter ApEn(m, r) is a positive number of infinite range. For
finite N, however, the largest possible value of
ApEn(m, r, N ) is

m + 1(r)
ln
(N
m)
1
, so
ApEn(m, r, N )
ln (N
m).
ApEn(m, r, N) is biased and suggests more similarity than is
present.
It is important to note that ApEn takes a template-wise approach to
calculating this average logarithmic probability, first calculating a
probability for each template. The ApEn algorithm thus requires that
each template contribute a defined, nonzero probability. This
constraint is overcome by allowing each template to match itself.
Formally, because
d[xm(i),
xm(i)] = 0
r, the
ApEn algorithm counts each template as matching itself, a practice we
will refer to as self-matching. This ensures that the functions
Cmi(r)
are > 0 for all i, thereby avoiding the occurrence of ln (0)
in the calculation. Thus ApEn(m, r, N ) is
certain to be defined under all circumstances. Pincus and Goldberger
(23, 27) note that, as a consequence of this practice,
ApEn(m, r, N ) is a biased statistic.
Formally, a statistic is biased if its expected value is not equal to
the parameter it estimates. ApEn statistics are biased, because the
expected value of ApEn(m, r, N ) is less
than the parameter ApEn(m, r) (27).
i, such that
d[xm + 1(i), xm + 1( j)]
r by Ai. The ApEn algorithm thus
assigns to the template xi(m) a
biased conditional probability of
(Ai + 1)/(Bi + 1),
which is always greater than the unbiased
Ai/Bi. In the
limit as N approaches infinity, Ai
and Bi will generally be large, making the
biased and unbiased probabilities asymptotically equivalent. Therefore,
this bias is evident only for the analysis of finite data sets and is a
characteristic of the statistic
ApEn(m, r, N ), rather than the
parameter ApEn(m, r). For a finite N, however,
the result is that ApEn(m, r, N ) is
biased toward lower values of ApEn and returns values below those
predicted by theory.
The largest deviation occurs when a large proportion of templates have
Bi = Ai = 0, since these templates are assigned a conditional probability of 1, corresponding to perfect order. Furthermore, the difference between the
biased and unbiased conditional probabilities assigned to individual templates makes the calculation sensitive to record length in a way
that depends on the conditional probability. Suppose that the unbiased
conditional probability is known and denote it by CP. For a given
um(i) let
Bi denote the number of template matches
without counting self-matches. The original algorithm estimates CP as
(1 + Ai)/(1 + Bi) = (1 + Bi × CP)/(1 + Bi).
The fractional error of this relative to CP is Err = {[(1 + Bi × CP)/
(1 + Bi)]
CP} /CP). To find
the value of Bi necessary to keep the
fractional error below a threshold Errmax, we isolate Bi from the inequality Err
Errmax, yielding Bi
[1
CP(Errmax + 1)]/(CP × Errmax).
For independent, identically distributed (iid) random numbers,
obtaining Bi matches of length
m requires a data set containing, on average,
Bi/(CP)m templates.
For iid random numbers and m = 2, estimating CP = 0.368 [ApEn(m, r, N ) = 1] within
Errmax = 0.05 requires Bi
33 and a data set of >240 points. Estimating CP = 0.135 [ApEn(m, r, N ) = 2] with similar
resolution requires Bi
127 and
>6,900 points.
The most straightforward way to eliminate the bias would be to remove
self-matching from the ApEn algorithm, leaving it otherwise unaltered.
However, without the inclusion of self-matches, the ApEn algorithm is
not defined unless
Cm + 1i(r) > 0 for every i. Removing self-matches would make ApEn
statistics highly sensitive to outliers; if there were a single
template that matched no other vector, ApEn could not be calculated
because of the occurrence of ln (0). For the set of uniform random
numbers shown in Fig. 1A, without
self-matches, ApEn(1, r, N ) would not have been
defined for r
0.63 and N < 4,000;
ApEn(2, r, N ) would have required r
1 for N < 4,000. Thus, for many practical applications, self-matches cannot simply be excluded when ApEn statistics are calculated.
|
Can the bias be corrected?
It is suggested that this bias can be reduced with a family of
estimators
by defining
Cmi,
(r) = (N
m + 1)
1 {
+ number
of j
i such that
d[xm(i),
xm( j)]
r}
(27). Then
m
(r) = (N
m + 1)
1
N
m + 1i = 1 Cmi,
(r)
and ApEn
(m, r, N ) =
m
(r)
m + 1
(r). For large N, this is very close to estimating the conditional probabilities by
(
+ Ai)/(
+ Bi).
When
< 1, the error due to counting self-matches is reduced, but
this algorithm would still report a conditional probability of 1 in the
event that only self-matches are encountered. Another strategy would be
to define
ApEn
2,
1(m, r,
N ) =
m
2(r)
m + 1
1(r)
where
1 and
2 are small and
1
2. This is similar to estimating the
conditional probabilities by
(
1 + Ai)/
(
2 + Bi). Thus a
template reporting only self-matches would be assigned a probability of
1/
2. This remains a biased estimate of
ApEn, with bias toward
ln (
1/
2)
when a large proportion of templates report only self-matches. If
1/
2 = (N
m + 1)
1 and
Ai = Bi = 0, the bias is toward the lowest observable probability. This is the
opposite bias of the original formulation of ApEn, where the bias is
toward the highest possible probability. Still another approach would
be to use the estimators
1 and
2 but
incorporate them into the calculation only when Ai or Bi = 0, respectively. This would estimate the probabilities as
Ai/Bi and would
simply redefine Ai =
1 when
Ai = 0 and Bi =
2 when Bi = 0.
1/
2 = CP, but CP is not known beforehand.
Thus no family of estimators
that minimizes bias can be chosen a priori.
SampEn statistics have reduced bias. We developed SampEn statistics to be free of the bias caused by self-matching. The name refers to the applicability to time series data sampled from a continuous process. In addition, the algorithm suggests ways to employ sample statistics to evaluate the results, as explained below.
There are two major differences between SampEn and ApEn statistics. First, SampEn does not count self-matches. We justified discounting self-matches on the grounds that entropy is conceived as a measure of the rate of information production (5), and in this context comparing data with themselves is meaningless. Furthermore, self-matches are explicitly dismissed in the later work of Grassberger and co-workers (2, 9, 11). Second, SampEn does not use a template-wise approach when estimating conditional probabilities. To be defined, SampEn requires only that one template find a match of length m + 1. We began from the work of Grassberger and Procaccia (10), who defined Cm(r) = (N
m + 1)
1
N
m + 1i = 1 Cmi(r),
the average of the Cmi(r)
defined above. This differs from
m(r)
only in that
m(r) is the average
of the natural logarithms of the
Cmi(r). They suggest approximating the Kolmogorov entropy of a process represented by a time series by limr
0
limn
limN
ln [Cm + 1(r)/
Cm(r)]; self-matches are
counted, and
Cm + 1(r)/Cm(r) = (N
m + 1)
N
mi = 1 Ai/(N
m)
N
m + 1i = 1 Bi.
In this form, however, the limits render it unsuitable for the analysis
of finite time series with noise. We therefore made two alterations to
adapt it to this purpose. First, we followed their later practice in
calculating correlation integrals (2, 9, 11) and did not consider
self-matches when computing
Cm(r). Second, we considered
only the first N
m vectors of length m,
ensuring that, for 1
i
N
m,
xm(i) and xm + 1(i) were defined.
We defined
Bmi(r)
as (N
m
1)
1 times the
number of vectors xm( j)
within r of xm(i),
where j ranges from 1 to N
m, and j
i to exclude self-matches. We then defined
Bm(r) = (N
m)
1
N
mi = 1 Bmi(r).
Similarly, we defined Ami(r)
as (N
m
1)
1 times the
number of vectors
xm + 1( j) within
r of xm + 1(i),
where j ranges from 1 to N
m
( j
i), and set
Am(r) = (N
m)
1
N
mi = 1 Ami(r).
Bm(r) is then the probability
that two sequences will match for m points, whereas
Am(r) is the probability that
two sequences will match for m + 1 points. We then defined the
parameter SampEn(m, r) = limN
{
ln [Am(r)/
Bm(r)]}, which is estimated
by the statistic SampEn(m, r, N ) =
ln [Am(r)/Bm(r)].
Where there is no confusion about the parameter r and the length m of the template vector, we set B = {[(N
m
1)(N
m)]/2}
Bm(r) and A = {[(N
m
1)(N
m)]/2}
Am(r), so that B is the
total number of template matches of length m and A is
the total number of forward matches of length m + 1. We note that A/B = [Am(r)]/Bm(r)],
so SampEn(m, r, N ) can be expressed as
ln (A/B).
The quantity A/B is precisely the conditional
probability that two sequences within a tolerance r for
m points remain within r of each other at the next
point. In contrast to ApEn(m, r, N ),
which calculates probabilities in a template-wise fashion, SampEn(m, r, N ) calculates the negative
logarithm of a probability associated with the time series as a whole.
SampEn(m, r, N ) will be defined except
when B = 0, in which case no regularity has been detected, or
when A = 0, which corresponds to a conditional probability of 0 and an infinite value of SampEn(m, r, N ).
The lowest nonzero conditional probability that this algorithm can report is
2[(N
m
1)(N
m)]
1.
Thus, the statistic SampEn(m, r, N )
has ln (N
m) + ln (N
m
1)
ln (2) as an upper bound,
nearly doubling ln (N
m), the dynamic range of
ApEn(m, r, N ).
Confidence intervals inform the implementation of SampEn statistics.
SampEn is not defined unless template and forward matches occur and is
not necessarily reliable for small numbers of matches. We have reviewed
SampEn(m, r, N ) calculation as a process
of sampling information about regularity in the time series and used sample statistics to inform us of the reliability of the calculated result. For example, say that we find B template matches,
allowing for no more than B forward matches, A of which
actually occur. We assign a value of 1 to the A forward matches
and a value of 0 to the (B
A) potential forward
matches that do not occur and compute the conditional probability
measured by SampEn as the average of this sample of 0s and 1s. For
operational purposes, we will assume that the sample averages follow a
Student's td distribution, where d
is the number of degrees of freedom. We can then say with 95%
confidence that the "true" average conditional probability of the
process is within
SDt(B
1,0.975)/
of the sample average, where SD is the sample standard deviation and
tB
1,0.975 is the upper 2.5th
percentile of a t distribution with B
1 degrees of
freedom (31). The size of the confidence intervals depends on the
number B and the number of forward matches. Informally, large
confidence intervals around
SampEn(m, r, N ) indicate that there are
insufficient data to estimate the conditional probability with
confidence for that choice of m and r. In addition,
confidence intervals allow standard statistical tests of the
significance of differences between data sets.
Cross-ApEn and cross-SampEn measure asynchrony.
Cross-ApEn is a recently introduced technique for comparing two
different time series to assess their degree of asynchrony or
dissimilarity (20, 28). The definition of cross-ApEn is very similar to
ApEn. Given two time series of N points
{u( j): 1
j
N } and
{v( j):
1
j
N }, form the vectors
xm(i) = {u(i + k):
0
k
m
1} and
ym(i) = {v(i + k):
0
k
m
1}. The distance between two
such vectors is defined as
d[xm(i),
ym( j)] = max {|u(i + k)
v( j + k)|:
0
k
m
1}. Define
Cmi(r)(v
u) as the number of ym( j)
within r of xm(i) divided by (N
m + 1), then define
m(r)(v
u) = (N
m + 1)
1
N
m + 1i = 1 ln [Cmi(r)
(v
u)], and
cross-ApEn(m, r, N )(v
u) =
m(r)(v
u)
m + 1(r)(v
u).
This is identical to the definition of the statistic ApEn, except the
templates are chosen from the series u and compared with
vectors from v. There is thus a directionality to this
analysis, and we will call the series that contributes the templates
the template series and the series with which they are compared the target series. We can then refer to
cross-ApEn(m, r, N )
(target
template).
u) can equal 0, and there is no assurance that
cross-ApEn(m, r, N )(v
u)
will be defined. Second, there is a "direction dependence" of
cross-ApEn analysis. We defined cross-SampEn to avoid both potential problems.
Cross-ApEn is not always defined. As noted above, cross-ApEn does not include self-matching and, thus, does not inherently suffer from the same bias as ApEn. A potential problem, however, remains in the necessity for each template to generate a defined, nonzero conditional probability. Thus each template must find at least one match for m + 1 points, or a probability must be assigned to it. No guidelines have been suggested for handling this potential difficulty. Cross-SampEn, on the other hand, requires only that one pair of vectors in the two series match for m + 1 points.
The family of MIX(P) stochastic processes (21) provided a testing ground for cross-ApEn. Informally, the MIX(P) time series of N points, where P is between 0 and 1, is a sine wave, where N × P randomly chosen points have been replaced with random noise. We calculated cross-ApEn(1, r, 250) for the pair [MIX(Q)
MIX(P)] and its direction
conjugate [MIX(P)
MIX(Q)] for 16 realizations of each of the 6 combinations of P = 0.1, 0.2, 0.3 and Q = 0.5, 0.7 over a range of values of r from 0.01 to 1.0. Cross-ApEn(1, r, 250)
[MIX(Q)
MIX(P)] was not defined
for any of the 96 pairs for r
0.16 and was defined for
all of them only for r
0.50. Cross-ApEn (1, r, 250)
[MIX(P)
MIX(Q)] was not defined
for any values of r
0.32 and was defined for all pairs only
for r = 1.0.
To broaden the conditions for which cross-ApEn was defined, we
introduced a correction factor into its algorithm. To avoid ln (0)
whenever
Cmi(r)(v
u) = 0 or
Cm + 1i(r)(v
u) = 0, we redefined them to be positive and nonzero. Rather than include
these factors into the calculation of each
Cmi(r)(v
u) and
Cm + 1i(r)(v
u), we minimized their impact on the overall calculation by including them
only when needed to ensure that cross-ApEn is defined. The cost of this
adjustment, however, is the introduction of bias. For this reason, we
tested cross-ApEn with two different correction strategies. The first
strategy, which we called bias 0, was very similar to self-matching; we
simply set any
Cmi(r)(v
u) or
Cm + 1i(r)(v
u) = 1 if it would otherwise have been 0. Thus, if a template matched no
other at all, it would be assigned a conditional probability of 1, as
with the original description of ApEn. If, however,
Cmi(r)(v
u)
0 but
Cm + 1i(r)(v
u) = 0, we redefined
Cm + 1i(r)(v
u) = (N
m)
1, so that the
probability assigned would be the lowest observable, nonzero
probability given the nonzero value of
Cmi(r)(v
u).
The second approach, which we called bias max, also only modified the
functions
Cmi(r)(v
u) and
Cm + 1i(r) (v
u) that would otherwise have been 0. Here,
Cmi(r) (v
u) was redefined to be 1 when it
would otherwise have been 0 and
Cm + 1i(r)(v
u) was redefined to be
(N
m + 1)
1, as for bias 0.
The only difference between the two strategies is that bias max assigns
to a template yielding no matches at all a probability of
(N
m)
1, the lowest nonzero
probability allowed by the length of the time series. Thus bias 0 sets
the bias toward a cross-ApEn value of 0 in the absence of any matches,
whereas bias max sets the bias toward the highest observable value of
cross-ApEn.
Cross-ApEn is direction dependent; cross-SampEn is not.
Because of the logarithms inside the summation,
m(r)(v
u)
will not generally be equal to
m(r)(u
v).
Thus
cross-ApEn(m, r, N )(v
u)
and its direction conjugate cross-ApEn(m, r, N )(u
v)
are unequal in most cases.
u)
as (N
m)
1 times the number of
vectors ym( j) within
r of xm(i), where
j ranges from 1 to N
m. We then defined
Bm(r)(v
u) = (N
m)
1
N
mi = 1 Bmi(r)(v
u).
Similarly, we set
Ami(r)(v
u)
as (N
m)
1 times the number of
vectors ym + 1( j)
within r of
xm + 1(i), where j ranges from 1 to N
m. We then defined
Am(r)(v
u) = (N
m)
1
N
mi = 1 Ami(r)(v
u).
Finally, we set
cross-SampEn(m, r, N )(v
u) =
ln {[Am(r)(v
u)]/
[Bm(r)(v
u)]}.
Examining this definition for direction dependence, we see that
(N
m)
Bmi(r)(v
u) is the number of vectors from v within r of the
ith template of the series u. Summing over the
templates, we see that
N
mi = 1 (N
m)
Bmi(r)(v
u)
simply counts the number of pairs of vectors from the two series that match within r. The number of pairs that match is clearly
independent of which series is the template and which is the target.
Because the last summation is equal to
(N
m)2
Bm(r)(v
u),
it follows that
Bm(r)(v
u)
is also direction independent, implying that
cross-SampEn(m, r, N )(v
u) = cross-SampEn(m, r, N )(u
v).
We further noted that cross-SampEn will be defined provided that
Am(r)(v
u)
0. Cross-SampEn, on the other hand, requires only that one pair of
vectors in the two series match for m + 1 points.
We calculated cross-SampEn(1, r, 250) for the same
realizations of the
[MIX(Q)
MIX(P)] pairs used above
to test cross-ApEn and over the same range of r. In contrast to
cross-ApEn(1, r, 250)
[MIX(P)
MIX(Q)] and
cross-ApEn(1, r, 250)
[MIX(Q)
MIX(P)], cross-SampEn(1, r, 250)
[MIX(P)
MIX(Q)], which is
identical to cross-SampEn(1, r, 250)
[MIX(Q)
MIX(P)], was defined for
all 96 pairs over the entire range of r considered.
ApEn and SampEn can be calculated analytically for series of random numbers. ApEn and SampEn derive from formulas suggested to estimate the Kolmogorov entropy of a process represented by a time series. At their root, each is a measurement of the conditional probability that two vectors that are close to each other for m points will remain close at the next point. There are several models, including sets of iid random numbers, for which the theoretical values of the parameters ApEn(m, r)(21) and SampEn(m, r) can be calculated. We show here the case of uniform random numbers, for which the theoretical values of ApEn and SampEn are nearly identical.
The expected value of the key probability can be calculated analytically for series of iid numbers based only on their probabilistic distribution. The numbers' independence implies that the probability that two randomly selected sequences within rSD of each other for the first m points will remain within r at their next points is simply the probability that any two points will be within a distance rSD of each other. For random numbers with density function p(x) and standard deviation SD, the expression for this probability is
+

p(t)
[
t + rSDt
rSD p(x) dx] dt. Because the parameter
ApEn(m, r) measures the average of the negative
logarithm of this probability, it is expected to return a value of
+

p(t)
ln [
t + rSDt
rSD p(x) dx] dt (21).
SampEn(m, r), on the other hand, measures the
logarithm of the conditional probability and is thus expected to return
a value of
ln {
+

p(t)
[
t + rSDt
rSD p(x) dx] dt}.
Because these expressions for the expected values of the parameters
ApEn and SampEn depend solely on r and the probabilistic character of the data, uniform random numbers provide a benchmark for
testing the estimation of ApEn over a range of parameters. In
particular, ApEn and SampEn are expected to give identical results for
uniformly distributed random numbers. Figure 1A shows a
histogram of the random numbers used and an excerpt from the sequence
(inset). Taking the derivative of

+

p(t)
ln [
t + rSDt
rSD p(x) dx] dt with respect to
ln (r) reveals that, for small r (in this case
r
1), ApEn(m, r) decreases in proportion
to ln (r) for small r. This result generalizes to
random numbers with other distributions provided that their density
functions are smoothly continuous and bounded. The expected values are
shown as the straight line in Fig. 1B. For random numbers with
other distributions, the theoretical values of ApEn and SampEn will
differ slightly. This is explained by noting that because
ln is a
convex function, Jensen's inequality (4) applied to each of the
integral expressions implies that the ApEn is expected to return a
value greater than or equal to SampEn.
| |
RESULTS AND DISCUSSION |
|---|
|
|
|---|
SampEn agrees with theory more closely than ApEn. For most processes, ApEn statistics are expected to have two properties. First, the conditional probability that sequences within r of each other remain within r should decrease as r decreases and the criterion for matching becomes more stringent. In other words, ApEn(m, r, N ) should increase as r decreases (22, 27). This expected property is demonstrated in Fig. 1B by the straight line, which plots the theoretically predicted values of ApEn and SampEn. Second, ApEn should be independent of record length. The plot of theoretical values of ApEn and SampEn in Fig. 1B illustrates this expectation. Because of their similarity, we expect SampEn statistics to exhibit similar properties. It has been suggested that record lengths of 10m-20m should be sufficient to estimate Ap- En(m, r) (27).
We tested ApEn and SampEn statistics on uniform, iid random numbers, because the results could be compared with the analytically calculated expected values. Figure 1, B and C, shows the performance of ApEn(2, r, N ) and SampEn(2, r, N ) on uniform iid random numbers. SampEn(2, r, N ) very closely matches the expected results for r
0.03 and
N
100 (Fig. 1, B and C), whereas
ApEn(2, r, N ) differs markedly from expectations
for N < 1,000 and r < 0.2. Figure
2 shows SampEn and ApEn as functions of
r (m = 2) for three sets of uniform random numbers
consisting of 100, 5,000, and 20,000 points. SampEn statistics for
r = 0.2 are in agreement with theory for much shorter
data sets. We investigated the general applicability of SampEn and
ApEn statistics to random numbers with other distributions. The
analysis of numbers with Gaussian, exponential, and
-distributions with parameter
= 1,2, ... , 10
gave results essentially identical to those shown in
Fig. 1.
|
SampEn shows relative consistency where ApEn does not.
A critically important expected feature of ApEn is relative consistency
(22, 27). That is, it is expected that, for most processes, if
ApEn(m1, r1)(S)
ApEn(m1, r1)(T ),
then
ApEn(m2, r2)(S),
ApEn(m2, r2)(T ).
That is, if record S exhibits more regularity than record
T for one pair of parameters m and r, it is
expected to do so for all other pairs. Graphically, plots of ApEn as a
function of r for different data sets should not cross over one
another. The determination that one set of data exhibits greater
regularity than another can be made only when this condition is met. We
tested this expectation using 1,000-point realizations of the
MIX(P) process, where the degree of order could be specified. Figure 3A shows that the ApEn
statistics of the less-ordered MIX(0.9), which has, on average, few
matches of length m for a given template (small
Bi) for small r, rises as a
function of r and crosses over a plot of ApEn statistics of the
more-ordered MIX(0.1). For r < 0.05, one would conclude
incorrectly that MIX(0.9) was more ordered than MIX(0.1). Thus relative
consistency does not hold for ApEn statistics.
|
Are SampEn statistics relatively consistent? We have shown that SampEn statistics appear to be relatively consistent over the family of MIX(P) processes, whereas ApEn statistics are not. Although we believe that relative consistency should be preserved for processes for which probabilistic character is understood, we see no general reason why ApEn or SampEn statistics should remain relatively consistent for all time series and all choices of parameters.
We propose a general, but by no means exhaustive, explanation for this phenomenon. SampEn is, in essence, an event-counting statistic, where the events are instances of vectors being similar to one another. When these events are sparse, the statistics are expected to be unstable, which might lead to a lack of relative consistency. Recall that SampEn(m, r, N ) is less than or equal to ln (B), the natural logarithm of the number of template matches. Suppose SampEn(m, r, N )(S) < SampEn(m, r, N ) (T ) and that the number of T's template matches, BT, is less than the number of S's template matches, BS, which would be consistent with T displaying less order than S. Provided that AT and AS, the number of forward matches, are relatively large, both SampEn statistics will be considerably lower than their upper bounds. As r decreases, BT and AT are expected to decrease more rapidly than BS and AS. Thus, as BT becomes very small, SampEn(m, r, N )(T ) will begin to decrease, approaching the value ln (BT), and could cross over a graph of SampEn(m, r, N )(S), where or while BS is still relatively large. Furthermore, as the number of template matches decreases, small changes in the number of forward matches can have a large effect on the observed conditional probability. Thus the discrete nature of the SampEn probability estimation could lead to small degrees of crossover and intermittent failure of relative consistency, and we cannot say that SampEn will always be relatively consistent. We have shown, however, that SampEn is relatively consistent for conditions where ApEn is not, and we have not observed any circumstance where ApEn maintains relative consistency and SampEn does not.One source of the residual bias in SampEn is correlation of
templates.
Although more consonant with theory than ApEn, we found that SampEn
statistics deviated from predictions for very short data sets. For
105 sets of Gaussian random numbers with m = 2, and
r = 0.2, we found that the deviation was <3% for record
lengths >100 but as high as 35% for sets of 15 points. Figure
3C shows the biased results of SampEn(2, 0.2, N )
for the range of 4
N
100. We suspected that the integral
expressions for the parameters ApEn(m, r) and SampEn(m, r) could not be used as expected values of
the statistics ApEn(m, r, N ) and
SampEn(m, r, N ) under all conditions,
because the expressions relied on the assumption that the templates
were independent of one another. As N decreases and m
increases, however, a larger proportion of templates are comprised of
overlapping segments of the record and are thus not independent.
Because of this correlation, results might deviate from these
predictions for short data sets. We thus tested the hypothesis that the
majority of the bias results from nonindependence of the templates.



p(w)
(
w + rSDw
rSD p(x)
{
x + rSDx
rSD p( y)
[
y + rSDy
rSD p(z) dz] dy} dx) dw]/
(


p(w) {
w + rSDw
rSD p(x)
[
x + rSDx
rSD p( y) dy] dx} dw),
where p(z) is the Gaussian probability density function. For r = 0.2, this was evaluated by numerical
integration and found to be 0.137. For 106 sets of 4 iid
Gaussian numbers, we calculated SampEn(2, 0.2, 4) to be 0.140. Thus
the observed and expected results nearly agree, indicating that bias
due to nonindependence of templates accounts for most of the observed
deviation from the conditional probability of 0.112 expected for
independent templates in this case.
A second way to test the hypothesis is to calculate SampEn statistics
for one set of data under two conditions: with and without overlapping
templates. The hypothesis predicts that the results should not match.
For this test, we chose the case of m = 2 and N = 6, the shortest record containing two nonoverlapping templates of length
m + 1 = 3. We can represent each set as
{a, b, c, x, y, z}.
For 106 sets of six Gaussian random numbers, we calculated
the conditional probability that (a, b, c)
was within r of (x, y, z) given that their first two points were close, thus calculating the probability for
pairs of disjoint templates. The result was 0.111, very close to the
expected 0.112 for independent templates. For the same number sets, we
then calculated the average value of SampEn(2, 0.2, 6), and the
result was 0.094. Thus the two results do not match, in support of the
hypothesis. We conclude from this analysis that the statistics
SampEn(m, r, N ) are not completely
unbiased under all conditions and that the bias of SampEn for very
small data sets is largely due to nonindependence of templates.
One method for removing this bias would be to partition the time series
{uj|1
j
N } into the m + 1 sets of neighboring, disjoint
vectors of length m + 1 Xi = {[ui + k(m + 1), u
i + k(m + 1) + 1, ... , u i + k
(m + 1) + m]: 0
k
[N
(m + i)]/(m + 1)},
where i, the initial point of the first template, ranges from 1 to m + 1. The conditional probability that vectors close for
m points remain close at the next point would be calculated for
each of the sets of vectors Xi and then
averaged. Because this calculation compares only disjoint templates, it
will not suffer from the bias introduced by nonindependent templates.
This truly unbiased approach has the potentially severe limitation of
reducing the number of possible template matches and enlarging the
confidence intervals about the SampEn estimate. Because this bias
appears to be present only for very small N, the disjoint
template approach does not appear necessary in usual practice.
Cross-SampEn shows relative consistency where cross-ApEn does not. As noted above, an essential feature of the measures of order is their relative consistency. That is, if one series is more ordered than another, it should have lower values of ApEn and SampEn for all conditions. We can extend this idea to cross-ApEn and cross-SampEn; if a pair of series is more synchronous than another pair, it should have lower values of cross-ApEn and cross-SampEn statistics for all conditions tested.
We tested the ability of cross-ApEn and cross-SampEn to distinguish between MIX(0.1) and the less-ordered MIX(0.6) processes. The strategy was to compare each with the intermediate MIX(0.3) process. The expected result was that the [MIX(0.3), MIX(0.1)] pair should appear more ordered than the [MIX(0.3), MIX(0.6)] pair, because MIX(0.1) is significantly more ordered than MIX(0.6). That is, cross-ApEn(2, r, 250) [MIX(0.3)
MIX(0.1)]
should be less than cross-ApEn(2, r, 250)
[MIX(0.3)
MIX(0.6)], and
cross-SampEn(2, r, 250) [MIX(0.3)
MIX(0.1)]
should be less than cross-SampEn(2, r, 250)
[MIX(0.3)
MIX(0.6)].
We tested this prediction bidirectionally, that is, MIX(0.3) served as
the template series for one analysis and as the target series for the
next, and over a range of tolerances r by using the bias-max
and bias-0 strategies for ensuring that cross-ApEn was defined. The
results are shown in Fig. 4. MIX(0.3) was
the template series in Fig. 4, C and E, and the target
series in Fig. 4, D and F.
|