Journal of Biosciences and Medicines, 2014, 2, 19-23
Published Online June 2014 in SciRes.
How to cite this paper: Trifonov, M. and Rozhkov, V. (2014) Age-Related Changes in Probability Density Function of Pair-
wise Euclidean Distances between Multichannel Human EEG Signals. Journal of Biosciences and Medicines, 2, 19-23.
Age-Related Changes in Probability Density
Function of Pairwise Euclidean Distances
between Multichannel Human EEG Signals
Mikhail Trifonov*, Vladimir Rozhkov
Sechenov Institute of Evolutionary Physiology and Biochemistry of the Russian Academy of Sciences,
Saint-Petersburg, Russia
Email: *mtri fon ov@ma ,
Received April 2014
The probability density functions (pdfs) and the first order structure functions (SF’s) of the pair-
wise Euclidean distances between scaled multichannel human EEG signals at different time lags
under hypoxia and in resting state at different ages are estimated. It is found that the hyper gam-
ma distribution is a good fit for the empirically derived pdf in all cases. It means that only two pa-
rameters (sample mean of EEG Euclidean distances at a given time lag and relevant coefficient of
variation) may be used in the approximate classification of empirical pdf’s. Both these parameters
tend to increase in the first twenty years of life and tend to decrease as healthy adults getting old-
er. Our findings indicate that such age-related dependence of these parameters looks like as age-
related dependence of the total brain white matter volume. It is shown that 15 min hypoxia (8%
oxygen in nitrogen) causes a significant (about 50%) decrease of the mean relative displacement
EEG value that is typical for the rest state. In some sense the impact of the oxygen deficit looks like
the subject getting older during short-term period.
EEG Developme nt, Hypoxia, Probability Density Functions, Hyper Gamma Distribution
1. Introduction
The progress in information theory, nonlinear dynamics, deterministic chaos theory, and random fractal theory
caused a wave of researches where the analysis of complexity EEG signals is done on the base on the using of
various complexity measures derived from them. Since the foundations of these theories are fundamentally dif-
ferent one can get a variety of complexity measures concerning the same EEG process. Detailed examination of
a number of such measures given in [1] shows that their variations with time are either similar or reciprocal, but
behaviors some of them are counter-intuitive and puzzling. The attempt to understand such behaviors is done in
[1] through a new multiscale complexity measure of EEG. Despite of very promising findings in these studies
Corresponding author.
M. Trifonov, V. Rozhkov
some of complexity measures adopted from nonlinear dynamics and chaos theory domain should be used with
caution for EEG classification, since the brain isn’t completely deterministic and the stochasticity may influence
its function in some cases [2]. This means that it is reasonable just now to avoid any speculations about what
types of deterministic and/or stochastic processes govern the EEG signals and to use some characteristics of
these signals that are insensitive to any process underlying them. The one-dimensional probability density func-
tions (pdf’s) of the human EEG relative displacements may be used as one of such characteristics. It is not suffi-
cient to infer the EEG dynamics but it is enough to capture some of its features.
2. Methods
2.1. EEG Data Collection
The eye closed resting state EEG data were recorded from 46 healthy subjects (12 adults and 34 children of
school-age and preschool-age) with 16 Ag/AgCl scalp electrodes placed according to the international 10-20
system over both hemispheres at a sampling frequencies fs of 185 Hz and 250 Hz. Three adults were also en-
gaged in night time EEG recording and one was involved in 15 minutes EEG recording under hypoxia (8%
oxygen in nitrogen). All school-age children took part in a longitudinal study of the EEG activity started at age
8.8 - 11.5 in 2005 and ended at age of 16.3 - 17.4 in 2011. The artifact free epochs selected for our analysis vary
in length from a several seconds to one minute.
2.2. EEG Data Analysis
Let all EEG data samples be represented as m-dimensional vectors (m-dimensional time series)
( )( )( )( )
{ }
,, ,T
XtX tXtXt= …
whe r e m = 16 is a number of channels (electrodes), Xj(t) is the signal amplitude on the channel j at the discrete
integer valued time moment t = 1, 2,…, N, N is the number of samples received, and the superscript T denotes the
matrix transpose operation. Since signals reveal significant spread of amplitude values from subject to subject
and for different sleep stages within a single sleep recording the original EEG data was centered by subtracting
their mean in every channel first and then scaled by the [det(R)](1/2m), where det denotes determinant, R =
E[δXδXT] is the sampling covariance matrix, δX = X - E[X], and E denotes statistical expectation. As the result
any new vector Y(t) = δX(t)/[det(R)](1/2m) has the same generalized variance independently on the subject since
the determinant of the covariance matrix Σ = EYδYT] is equal here to 1. Geometrically the quantity [det(Σ)]1/2
determines the volume of the confidence ellipsoid for any particular confidence level and the scaling proposed
here makes the distributions of any vectors Y to be equivalent in the sense that they occupy the same hypervolume
in the m-space. It means that the ellipsoids with different orientations and different semi-axes but having the same
generalized variance will be considered here as equivalent.
The vector sequence (Y(1), Y(2 ) , …, Y(N)) specifies a personal EEG trajectory in initial m-di me nsio na l space.
Our aim is to estimate the probability density functions (pdf’s) of the EEG relative displacements ΔYτ (the Euc-
lidean distances between scaled multichannel human EEG signals pairs at given time lag τ) along this trajectory,
defined by
2 1/2
( )()(( )())[]
YYtYtY tY t
∆= +−=+−
The idea to use such pdf’s for describing and discriminating EEG patterns was inspired by the approach [3]
proposing to reduce the shape matching problem to the comparison of probability distributions of the distance
between two random points on a surface provide a robust method for discriminating between some classes of 3D
objects. Since the relative displacement ΔYτ value, being m-dimensional Euclidean distance, doesnt depend on
how the axes of the initial space are chosen and how many of principal components exhaust a given amount of
the total variance there is no need to care about the proper state space reconstruction for pdf’s deriving. But un-
like artificial objects having rather stable shapes, the EEG trajectory forms a “living” shape that may evolve in
time, i.e. may change its geometrical characteristics in dependence on epoch selection during a given EEG re-
cording. Our aim is not only to analyze the empirically derived pdf’s but also to find some types of theoretical
distributions that can fit them. Here we will follow the general assumption borrowed from [4] “that if the me-
M. Trifonov, V. Rozhkov
chanism (experiment) to generate the samples is the same, then the distribution type that describes the datasets
will also be the same”. In our case it means that we need to identify a single type of theoretical distribution that
can fit the different datasets of {ΔYτ} by altering its parameters. Most likely the best candidate for such distribu-
tion may be a hyper gamma distribution proposed in [5] that is written here as
( )exp()fY YY
ττ τ
∆=⋅∆− ∆
( )
γαβν α
= Γ
, Γ(∙) is the gamma function, α and ν are shape parameters, and β is the scale parame-
ter. It follows from the fact that the “random EEG” relative displacement ΔYτ (obtained from m-dimensional
vector X(t) whose components are independent normal random variables having zero mean and variance σ2) is
distributed according to the chi distribution that is a special case of (2) when ν = m, α = 2, and β = 0.25/σ2.
The pdf’s f(ΔYτ) are obtained in the standard way by fitting an empirical probability functions f0(ΔYτ)
represented by normalized histograms. The histogram bin size Wτ is selected here as the half-sum of the
23.73 ()WN
= −
calculated according to [6] and [7] suggestions respect-
tively, where IQRτ is the interquartile range defined as the 75th percentile minus the 25th percentile, and στ is the
standard deviation of the distribution. The number of bins hτ is defined as entire(max(ΔYτ)/Wτ) + 1.
The fitting procedure is based here on the estimation of parameters ν, α, and β by the moment method. Since
any covariance matrix Σ for our EEG data under analysis has full column rank, we set ν = m. Then the estimates of
α is defined by the squared coefficient of variation (CVτ) according to the following equation
([])((3) )((1) )[((2) )]1EY mmm
σαα α
∆ =Γ+Γ+Γ+−
Getting the estimate of α using the Equation (3) one can obtain the estimate of β according to the following
[((2) )][[]((1) )]mEY m
βα α
= Γ+∆⋅Γ+
The fitting procedure is limited here by the fitting range (FRτ) endpoints that are defined as 1st percentile and
the 99th percentile respectively. The criterion of fitting quality evaluation is based on the relative error of the
first three moments and the percentage root mean square difference (PRD) suggested in [8]
( )
0, ,0,
100%( )( )(( ))
kk k
τττ τ
= =
=⋅∆ −∆∆
3. Results
Examples of the empirical pdf f0(ΔY1) and its hyper gamma fit f(ΔY1) for one of the subjects are shown in Figure
1(a). In this case the value of N is equal to 8060 that corresponds to about 43.6 second EEG recording with the
sampling frequency fs = 185 Hz. Time lag τ = 1 corresponds to time interval of 0.0054 sec and W1 = 0.2. For
comparison the pdf’s for the “random” EEG with σn = 1 and for the original EEG with randomly permuted data
samples within each channel are shown on the same Figure 1(a). Visually the single type of theoretical f(ΔY1)
(a) (b )
Figure 1. (a) The examples of the empirical pdf’s f0(ΔY1) and hyper gamma fits f(ΔY1) for the real (blue), real, but randomly
permuted (green), and “random EEG” data (red) and (b) relevant structure functions SFτ.
M. Trifonov, V. Rozhkov
provides a reasonably good fit for all three empirical f0(ΔY1). The maximum value of PRD1 here refers to the
case of the original EEG and is less than 8.8% while the maximum relative error in the first three moments does
not exceed 0.5%.
In general there is a noticeable inter-individual variation in PRDτ value even for the same sample size N but in
average PRDτ decreases proportional to N0.17 as N increases. As usual the relatively high PRDτ values occur
when the empirical pdf’s f0(ΔYτ) have a rather long tales. It might mean that they can represent a mixture of hy-
per gamma distributions with different time scales. Since the mean PRDτ value here is about 9.1% and visually
any empirical pdf is fitted well enough with the single hyper gamma pdf, we restrict our attention to this case.
According to Equations (3) and (4) it is enough to analyze here only two parameters-E[ΔYτ] and στ/E[ΔYτ].
The first one defines the individual first order structure function (SFτ) that rapid increases for small τ as τδ (δ <
1) leading to shifting the f0(ΔYτ) from f0(ΔY1) to the right. The CVτ also reveals the individual power law depen-
dence on τ for small τ. But starting with some τ* value the SFτ begins to oscillate around its sill level defined by
the variance of stochastic component of ΔYτ. It causes the oscillatory behavior of f0(ΔYτ) that depends on the
subject’s individuality and is most conspicuous in the resting state. It was found that for all our subjects in this
state the oscillation frequency lies in the alpha range (see Figure 1(b)). For comparison, both SFτ and CVτ cor-
responding to the “random” EEG do not depend on τ and for σn = 1 are equal to 5.57 and 0.99 respectively. The
SFτ and CVτ for the original EEG with randomly permuted data samples within each channel do not depend on τ
as well.
The age dependence of E[ΔY1] in the eye-closed resting state is shown in Figure 2(a). The original E[ΔY1]
values getting at fs = 185 Hz were multiplied by 0.74 to approximately adopt them to fs = 250 Hz. The blue dots
show such adopted data representing in fact the upper bounds for E[ΔY1] that actually could be in this case. The
presented data reflect high variability of E[ΔY1] both within some subjects and between different subjects. Nev-
ertheless there is some tendency of E[ΔY1] increasing in the first twenty years of life and of E[ΔY1] decreasing as
adults getting older. This tendency can be approximately described by log-normal curves in Figure 2(a) de-
picted by colored lines. We are going to define more exactly the type of curve fitting at a later time.
Having conducted multiple longitudinal studies, we find that the CV1 reveals noticeable inter-individual dif-
ferences being rather stable within an individual subject over long period of time. The analysis of 10 min resting
state, following 15 min hypoxia and 15 min after hypoxia study did not show significant variation in CV 1 as
well. The relative standard deviation of CV1 for this case is about 4.4%.
The time-specific dependence of E[ΔY1] in 15 min hypoxia (red dots) was followed by 10 min resting state
(blue dots) is shown in Figure 2(b). The hypoxia causes a dramatic decreasing of E[ΔY1] relative to the resting
state. During the last minute of the hypoxia the initial E[ΔY1] value is reduced by 50%. In some sense the impact
of the oxygen deficit looks like as the subject get older during short-term period (see age-dependence of E[ΔY1]
on Figure 2(a)). During the 15 min period of time after the hypoxia E[ΔY1] tends to reach the value inherent in
the resting state. It was found that changes of E[ΔY1] have a strong positive correlation (r = 0.89) with changes
in oxygen saturation.
(a) (b )
Figure 2. (a) Age dependence of E[ΔY1] in the eye-closed resting state and (b) time-specific dependence of E[ΔY1] in 15
min hypoxia (red dots) was followed by 10 min resting state (blue dots) and 15 min post-hypoxic recovery (green dots).
The different color lines on the panel (a) corresponding to log-normal curves and the green line and curve on the panel (b)
are used for a fitting of the empirical data.
M. Trifonov, V. Rozhkov
4. Conclusions and Future Work
The results of this study allow us to make a preliminary conclusion that one-dimensional pdf’s of EEG relative
displace ments can be used for understanding of the real EEG dynamics in various functional states and different
subjects. To a first approximation, in each case the empirically derived pdf are fitted quite well by the single
hyper gamma distribution. It means that only two parameters (sample mean of EEG relative displacements and
coefficient of variation) may be taken into account. Both these parameters exhibit subject’s individuality. The
first one reveals age and state dependence while the second one stays rather stable for a given subject over long
period of time except sleep stages. It is interesting to note that age-related dependence of E[ΔY1] looks like as
age-related dependence of the total brain white matter volume given in [9]. In addition the non-linear age effect
on E[ΔY1] is consistent with the suggestion that during late childhood period there is a shift of topological or-
ganization of brain white matter toward a more randomized configuration [10].
In our future research we are going to analyze much more EEG records to investigate in details the age and
time-specific dependence of the parameters mentioned above. The next interesting aspect of our future work is
the improvement of the fitting quality of the empirically derived pdf of EEG relative displacements by using
distribution that consists of a mixture of hyper gamma distribution.
The authors would like to express their sincere thanks to Prof. M. N. Tsitseroshin, Dr. E. A. Panasevich, and Dr.
S. S. Bekshaev for providing EEG data.
[1] Gao, J., Hu, J. and Tung, W. (2011) Complexity Measures of Brain Wave Dynamics. Cognitive Neurodynamics, 5,
[2] Masquelier, T. (2013) Neural Variability, or Lack There of. Frontiers in Computational Neuroscience, 7, 1-7.
[3] Osada, R., Funkhouser, T., Chazelle, B. and Dobkin, D. (2002) Shape Distributions. ACM Transactions on Graphics,
21, 807-832.
[4] Nikolova, N.D., Toneva-Zheynova, D., Kolev, K. and Tenekedjiev, K. (2013) Monte Carlo Statistical Tests for Identity
of Theoretical and Empirical Distributions of Experimental Data. In: Chan, V., Ed., Theory and Applications of Monte
Carlo Simulations, InTech, 1-26.
[5] Suzuki, E. (1964) Hyper Gamma Distribution and Its Fitting to Rainfall Data. Papers in Meteorology and Geophysics,
15, 31-51.
[6] Izenman, A.J. (1991) Recent Developments in Nonparametric Density Estimation. Journal of the American Statistical
Association, 86, 205-224.
[7] Scott, D.W. (2004) Multivariate Density Estimation and Visualization. P apers / Humboldt-Universität Berlin, Center for
Applied Statistics and Economics (CASE), No. 2004, 16.
[8] Alfaouri, M., Daqrouq, K., Abu-Isbeih, I.N., Khalaf, E.F., Al-Qawasmi, A-R. and Al-Sawalmeh, W. (2009) Quality
Evaluation of Reconstructed Biological Signals. American Journal of Applied Science, 5, 187-193.
[9] Sowell, E.R., Peterson, B.S., Thompson, P.M., Welcome, S.E., Henkenius, A.L. and Toga, A.W. (2003) Mapping Cor -
tical Change across the Human Life Span. Nature Neuroscience, 6, 309-315.
[10] Chen, Z., Liu, M., Gross, D.W. and Beaulieu, C. (2013) Graph Theoretical Analysis of Developmental Patterns of the
White Matter Netwo rk. Frontiers in Human Neuroscience, 7, 1-13.