# Graph-Based Multivariate Multiscale Dispersion Entropy: Efficient
Implementation and Applications to Real-World Network Data

John Stewart Fabila-Carrasco, Chao Tan, and Javier Escudero J.S.F.C is with
the School of Informatics, University of Edinburgh, Edinburgh, EH8 9LE, UK
(e-mail<EMAIL_ADDRESS>C.T. is with the School of Electrical and
Information Engineering Tianjin University Tianjin, 300072, China (e-mail:
tanchao@tju.edu.cn). J.E is with the School of Engineering, Institute for
Imaging, Data and Communications, University of Edinburgh, Edinburgh, EH9 3FB,
UK (e-mail J.E. Javier.Escudero@ed.ac.uk). The work of J.S.F.C. and J.E. was
supported by Leverhulme Trust through a Research Project under Grant
RPG-2020-158. For the purpose of open access, the author has applied a
Creative Commons Attribution (CC BY) licence to any Author Accepted Manuscript
version arising from this submission.

###### Abstract

We introduce Multivariate Multiscale Graph-based Dispersion Entropy
($\operatorname{mvDE_{G}}$), a novel, computationally efficient method for
analysing multivariate time series data in graph and complex network
frameworks, and demonstrate its application in real-world data.
$\operatorname{mvDE_{G}}$ effectively combines temporal dynamics with
topological relationships, offering enhanced analysis compared to traditional
nonlinear entropy methods. Its efficacy is established through testing on
synthetic signals, such as uncorrelated and correlated noise, showcasing its
adeptness in discerning various levels of dependency and complexity. The
robustness of $\operatorname{mvDE_{G}}$ is further validated with real-world
datasets, effectively differentiating various two-phase flow regimes and
capturing distinct dynamics in weather data analysis. An important advancement
of $\operatorname{mvDE_{G}}$ is its computational efficiency. Our optimized
algorithm displays a computational time that grows linearly with the number of
vertices or nodes, in contrast to the exponential growth observed in classical
methods. This efficiency is achieved through refined matrix power calculations
that exploit matrix and Kronecker product properties, making our method faster
than the state of the art. The significant acceleration in computational time
positions $\operatorname{mvDE_{G}}$ as a transformative tool for extensive and
real-time applications, setting a new benchmark in the analysis of time series
recorded at distributed locations and opening avenues for innovative
applications.

###### Index Terms:

Multivariate Entropy, Dispersion Entropy, Time Series Analysis, Graph signals,
Networks, Real-world data.

##  I Introduction

Real-world data, particularly from systems like industrial and weather
phenomena, is often characterized by dynamic and complex behaviours that
exhibit transient, chaotic, and nonlinear dynamics over time, accompanied by a
mixture of randomness and uniformity in spatial distribution [22, 30]. These
systems are fundamental to various phenomena across environmental,
technological, and scientific domains, influencing everything from water
distribution [32] to advanced manufacturing processes [8]. Gaining a deep
understanding of the underlying mechanisms of these systems is crucial.
Despite substantial advancements in the field, a notable challenge persists in
the development of sophisticated mathematical tools. These tools are essential
not only for an accurate analysis of specific applications like two-phase flow
and weather dynamics but also for understanding and analysing general network
behaviours across a wide array of domains, including any form of networked or
graph-based data, addressing both academic research needs and practical
applications [34].

Multivariate entropy techniques are essential tools in the analysis of data
consisting of multiple time series, or channels [2, 31]. These techniques
offer a perspective of complexity that extends beyond the scope of univariate
analysis [1, 27]. Among these, classical Multivariate Dispersion Entropy
($\operatorname{mvDE}$) is particularly noteworthy for its robust approach in
quantifying the complexity and dynamics of multivariate data [5]. Its wide-
ranging applications, from biomedicine to fault diagnosis, demonstrate its
significance in various research fields [4]. In comparison to other entropy
methods such as the Multivariate Sample Entropy [1], $\operatorname{mvDE}$ has
better performance, more stability, and is particularly effective for shorter
time series [5].

The landscape of data analysis is rapidly evolving, increasingly pivoting
towards the utilization of graph-based structures and complex networks [15,
24]. This shift is driven by technological advancements that have vastly
increased our capacity to collect and analyse data across diverse fields,
ranging from environmental monitoring to public health systems [6, 25]. Graph-
based methods are gaining prominence as they offer a powerful framework to
uncover complex relationships and dependencies within this data, which
traditional tabular data formats and analysis techniques might not fully
capture. These methods leverage the unique topological properties of graph-
structured data, enabling a more nuanced understanding of the intricate
interplay within these varied systems [24]. The development of Permutation
Entropy and Dispersion Entropy tailored for graph data marks a significant
advancement, facilitating analyses that incorporate these crucial topological
dimensions [14, 13, 15].

Moreover, the approach we have developed for Permutation and Dispersion
Entropy on graph structures has influenced subsequent methodologies, as
evidenced by recent advancements like the extension of Bubble Entropy to graph
signals [11], directly utilising our foundational work. This not only
demonstrates the robustness of our methods but also indicates that our
approach is pioneering paths for further innovations in the analysis of
network data. However, while these methods provide valuable insights into the
topological dimension of the data, they often do not concurrently consider the
temporal dimension. This limitation is crucial, especially in applications
where understanding the interplay between topological (or spatial) structure
and temporal dynamics is essential for a comprehensive analysis, such as
climate science [26] or industrial applications [10].

For example, in two-phase flow analysis, various entropy metrics have been key
to understanding the system dynamics. Sample Entropy has been used to measure
regularity in fluid dynamics time-series [4]. Permutation Entropy helps in
detecting changes in flow patterns due to its sensitivity to temporal
variations [29]. Despite their utility, these metrics often do not fully
capture the combined spatial and temporal complexities of two-phase flow
systems. This limitation highlights the need for more advanced tools like
$\operatorname{mvDE_{G}}$, which are capable of integrating both topological
and temporal data for a more comprehensive analysis.

###  I-A Contributions

We introduce the Multiscale Multivariate Graph-based Dispersion Entropy
($\operatorname{mvDE_{G}}$), an advanced method for analysing data on complex
networks. Differing from classical Dispersion Entropy [27], which focuses on
temporal information, and Dispersion Entropy for Graphs [14], which emphasizes
graph structure, $\operatorname{mvDE_{G}}$ integrates both temporal and
topological dimensions for a more comprehensive analysis.

Furthermore, the proposed algorithmic implementation notably increases
computational efficiency, making $\operatorname{mvDE_{G}}$ suitable for
processing both short and long time series. Our optimized algorithm displays a
computational time that grows linearly with the number of vertices or nodes,
in contrast to the exponential growth observed in classical methods.
$\operatorname{mvDE_{G}}$ demonstrates improved performance over classical
$\operatorname{mvDE}$ in various test scenarios, including synthetic and real-
world datasets. The efficiency of $\operatorname{mvDE_{G}}$ can be extended to
other graph-based entropy metrics, notably, the principles underlying
$\operatorname{mvDE_{G}}$ could be seamlessly adapted to improve the
efficiency of the Multivariate Multiscale Permutation Entropy [12].
Specifically, algorithms that use Cartesian Graph products combined with our
adjacency matrix approach, as detailed in previous works on Permutation
Entropy [13] and extended to Bubble Entropy [11], will find our method
straightforwardly applicable. This adaptability ensures that our efficient
approach can significantly benefit a wide range of graph-based entropy
calculations, fostering advancements across various fields.

In applying $\operatorname{mvDE_{G}}$ to real-world data including two-phase
flow systems and weather phenomena, we have successfully delineated distinct
entropy profiles for various complex patterns, providing acute insights into
their dynamics. The precision of the method is especially notable in
distinguishing different flow regimes and weather patterns, revealing clear
entropy demarcations at lower scales. These results not only demonstrate the
robustness of $\operatorname{mvDE_{G}}$ in handling intricate data from
industrial processes and environmental monitoring but also underscore its
potential as an indispensable tool in advanced time series analysis.

###  I-B Graph Theory Notation

A _graph_ $G$ is defined as $G=(\mathcal{V},\mathcal{E})$, where
$\mathcal{V}=\\{1,2,3,\dots,n\\}$ is the set of vertices, and
$\mathcal{E}\subset\\{(i,j):i,j\in\mathcal{V}\\}$ is a set of edges [7]. The
adjacency matrix $\mathbf{A}$, a $N\times N$ matrix, represents the
connections between vertices, with $\mathbf{A}_{ij}=1$ if
$(i,j)\in\mathcal{E}$ and $0$ otherwise. For weighted graphs,
$\mathbf{A}_{ij}$ represents the weight of the connection between $i$ and $j$.
Our approach is applicable to any weighted graph, whether directed or
undirected, by appropriately adjusting the adjacency matrix $\mathbf{A}$.

The _Cartesian product_ of two graphs ${G}=(\mathcal{V},\mathcal{E})$ and
${G}^{\prime}=(\mathcal{V}^{\prime},\mathcal{E}^{\prime})$, denoted by
${G}\square{G}^{\prime}$, is a graph defined as:

  1. 1.

$\mathcal{V}({G}\square{G}^{\prime})=\mathcal{V}\times\mathcal{V}^{\prime}=\\{(v,v^{\prime})\mid
v\in\mathcal{V}\text{ and }v^{\prime}\in\mathcal{V}^{\prime}\\}$.

  2. 2.

$(v,v^{\prime})$ and $(u,u^{\prime})$ are adjacent in ${G}\square{G}^{\prime}$
if and only if either: $v=u$ and $v^{\prime}$ is adjacent to $u^{\prime}$ in
${G}^{\prime}$, or $v^{\prime}=u^{\prime}$ and $v$ is adjacent to $u$ in
${G}$.

##  II Multivariate Multiscale Graph-based Dispersion Entropy
($\operatorname{mvDE_{G}}$)

The algorithm $\operatorname{mvDE_{G}}$ enhances dispersion entropy by
integrating topological and temporal data dimensions. It builds on the
topological approach in [14] and the temporal focus of [5], offering a more
thorough analysis. However, the initial computation of large matrix powers in
$\operatorname{mvDE_{G}}$ presents a computational challenge, which is
addressed in Sec. III through an efficient implementation strategy. The
algorithm consists of two main steps: 1) a coarse-graining process, and 2) the
calculation of $\operatorname{mvDE_{G}}$ at each scale $\tau$.

Consider a multivariate signal
$\mathbf{X}=\left\\{x_{k,i}\right\\}_{k=1,2,\ldots,p}^{i=1,\ldots,N}$, where
$p$ is the number of channels and $N$ is the number of observations per
channel.

###  II-A Coarse-Graining Process for Multivariate Signals

For each channel, the original signal is divided into non-overlapping segments
of length $\tau$, named scale factor. Next, for each channel, the average of
each segment is calculated to derive the coarse-grained signals as follows:

$z_{k,i}^{(\tau)}=\frac{1}{\tau}\sum_{b=(i-1)\tau+1}^{i\tau}x_{k,b},1\leq
i\leq\left\lfloor\frac{N}{\tau}\right\rfloor=L,1\leq k\leq p\;;$

where $L$ denotes the length of the coarse-grained signal. While this study
uses a straightforward coarse-graining approach, various alternative coarse-
graining methods are also explored in the literature [4]. In the second step,
$\operatorname{mvDE_{G}}$ is calculated for each coarse-grained signal.

###  II-B Graph-Based Multivariate Dispersion Entropy

To analyse the inter-channel interactions within $\mathbf{X}$, we utilize an
adjacency matrix $I_{p}\in\mathbb{R}^{p\times p}$. This matrix, $I_{p}$,
encapsulates the connectivity or interaction patterns between the channels. It
can be predefined based on system-specific knowledge, or fully connected graph
for simplicity. Alternatively, $I_{p}$ can be inferred from the data itself,
employing various graph learning techniques such as in [33]. This flexibility
in defining $I_{p}$ allows for tailored analysis of the multivariate signal,
taking into account the unique interaction dynamics present in different
datasets or systems. The $\operatorname{mvDE_{G}}$ is defined through the
following steps:

  1. 1.

Embedding Matrix Construction: Given an embedding dimension $2\leq
m\in\mathbb{N}$ and a class number $c\in\mathbb{N}$, the embedding matrix
$\textbf{Y}\in\mathbb{R}^{N\times m}$ is formed as
$\textbf{Y}=[\textbf{y}_{0},\textbf{y}_{1},\cdots,\textbf{y}_{m-1}]$. Each
column $\textbf{y}_{k}=D\mathbf{A}^{k}\textbf{v}$, where $\mathbf{A}$ is the
adjacency matrix of the graph ${\overrightarrow{P}_{N}}\square I_{p}$, $D$ is
a diagonal normalization matrix, and v is a vectorized form of the
multivariate signal $\mathbf{X}$.

  2. 2.

Mapping to Classes: A mapping function
$M\colon\mathbb{R}\rightarrow\mathbb{N}_{c}$, where
$\mathbb{N}_{c}=\\{1,2,\dots,c\\}$, is applied element-wise to $\mathbf{Y}$ to
classify each entry into one of $c$ classes. This results in a matrix
$M(\mathbf{Y})\in\mathbb{N}_{c}^{N\times m}$.

  3. 3.

Dispersion Patterns and Relative Frequencies: Each row of $M(\mathbf{Y})$,
termed an embedding vector, corresponds to a unique dispersion pattern. The
relative frequency of each pattern $\pi\in\Pi$ is calculated across the
dataset.

  4. 4.

Entropy Computation: The $\operatorname{mvDE_{G}}$ is the normalized Shannon
entropy of the dispersion patterns, computed as:

$\operatorname{mvDE_{G}}(\mathbf{X},m,L,c)=-\dfrac{1}{\log(c^{m})}\sum_{\pi\in\Pi}p(\pi)\ln
p(\pi).$

This approach effectively captures the complex interactions within
multivariate data, leveraging both the temporal sequence and the underlying
graph structure. For embedding matrix construction details, refer to [14]. The
primary computational challenge involves computing powers of the large
$Np\times Np$ matrix $\mathbf{A}_{{\overrightarrow{P}_{N}}\square I_{p}}$,
representing the graph’s adjacency matrix. We address this challenge and its
solutions in Sec. III.

##  III Efficient Algorithm Implementation

Utilizing matrix properties and the Kronecker product [18], this section
introduces an efficient $\operatorname{mvDE_{G}}$ implementation to overcome
computational challenges of prior graph signal entropy metrics [13].

Graph | Model | Notation | $\lvert\mathcal{V}\rvert$ | Adjacency
---|---|---|---|---
Directed Path | Time | ${\overrightarrow{P}_{N}}$ | $N$ | $\mathbf{A}_{\overrightarrow{P}_{N}}$
Interaction Matrix | Topology | ${I_{p}}$ | $p$ | $\mathbf{A}_{I_{p}}$
TABLE I: Overview of the graphs and matrices.

The adjacency matrix of ${\overrightarrow{P}_{N}}\square I_{p}$ is given by:

$\mathbf{A}_{{\overrightarrow{P}_{N}}\square
I_{p}}=\mathbf{A}_{\overrightarrow{P}_{N}}\otimes\mathbb{I}_{p}+\mathbb{I}_{N}\otimes\mathbf{A}_{I_{p}}$
(1)

where $\otimes$ denotes the Kronecker product and $\mathbb{I}_{N}$ is the
$N\times N$ identity matrix. Let
$X=\mathbf{A}_{\overrightarrow{P}_{N}}\otimes\mathbb{I}_{p}$ and
$Y=\mathbb{I}_{N}\otimes\mathbf{A}_{I_{p}}$. We demonstrate that $B$ and $C$
commute:

$\displaystyle BC$
$\displaystyle=(\mathbf{A}_{\overrightarrow{P}_{N}}\otimes\mathbb{I}_{p})(\mathbb{I}_{N}\otimes\mathbf{A}_{I_{p}})=\mathbf{A}_{\overrightarrow{P}_{N}}\mathbb{I}_{N}\otimes\mathbb{I}_{p}\mathbf{A}_{I_{p}}$
$\displaystyle=\mathbf{A}_{\overrightarrow{P}_{N}}\otimes\mathbf{A}_{I_{p}}=\mathbb{I}_{N}\mathbf{A}_{\overrightarrow{P}_{N}}\otimes\mathbf{A}_{I_{p}}\mathbb{I}_{p}$
$\displaystyle=(\mathbb{I}_{N}\otimes\mathbf{A}_{I_{p}})(\mathbf{A}_{\overrightarrow{P}_{N}}\otimes\mathbb{I}_{p})=CB\;.$

Consequently, the $m$-th power matrix in Eq. 1 is:

$\displaystyle(\mathbf{A}_{{\overrightarrow{P}_{N}}\square I_{p}})^{m}$
$\displaystyle=\sum_{k=0}^{m}\binom{m}{k}(\mathbf{A}_{\overrightarrow{P}_{N}}\otimes\mathbb{I}_{p})^{k}(\mathbb{I}_{N}\otimes\mathbf{A}_{I_{p}})^{m-k}$
$\displaystyle=\sum_{k=0}^{m}\binom{m}{k}((\mathbf{A}_{\overrightarrow{P}_{N}})^{k}\otimes\mathbb{I}_{p})(\mathbb{I}_{N}\otimes(\mathbf{A}_{I_{p}})^{m-k})$
$\displaystyle=\sum_{k=0}^{m}\binom{m}{k}((\mathbf{A}_{\overrightarrow{P}_{N}})^{k}\mathbb{I}_{n})\otimes(\mathbb{I}_{p}(\mathbf{A}_{I_{p}})^{m-k})$
$\displaystyle=\sum_{k=0}^{m}\binom{m}{k}(\mathbf{A}_{\overrightarrow{P}_{N}})^{k}\otimes(\mathbf{A}_{I_{p}})^{m-k}\;.$

Finally, the adjacency matrix of $\overrightarrow{P}_{N}$,
$\mathbf{A}_{\overrightarrow{P}_{N}}$, and its $k$-th power for $1\leq k\leq
N$ are computed as:

$\mathbf{A}_{\overrightarrow{P}_{N}}^{k}(i,j)=\begin{cases}1&\text{if
}j=i+k\text{ and }j\leq N,\\\ 0&\text{otherwise},\end{cases}$

where
$\mathbf{A}_{\overrightarrow{P}_{N}}(i,j)=\mathbf{A}_{\overrightarrow{P}_{N}}^{1}(i,j)$
for $k=1$.

This efficient implementation reduces the computation of the $m$-th power of a
$Np\times Np$ matrix to a sum of Kronecker products involving smaller $p\times
p$ matrices. This approach effectively addresses the high computational cost
associated with graph-based entropy metrics, a major challenge noted in
previous studies [13, 14].

##  IV Synthetic Signals

Synthetic signals like multivariate $1/f$ noise and White Gaussian Noise (WGN)
are essential for validating multivariate multiscale entropy techniques.
Multivariate $1/f$ noise, common in natural phenomena, represents complex,
correlated systems and contrasts with the high irregularity of WGN [5, 9, 20].
This contrast enables a detailed evaluation of entropy methods, as these
signals provide benchmarks to assess how entropy measures handle different
complexities and predictabilities in signal patterns, vital for applications
[34].

###  IV-A Uncorrelated Noise Analysis

We generate a trivariate time series $\mathbf{F}(p)$, where
$p\in\\{0,1,2,3\\}$ represents the number of channels with $1/f$ noise, and
$3-p$ channels with WGN, i.e.,
$\mathbf{F}(p)=\\{\underbrace{\text{WGN},\ldots,\text{WGN}}_{p\text{
times}},\underbrace{1/f,\ldots,1/f}_{3-p\text{ times}}\\}\;.$

Here, $p$ allows a systematic transition from all channels representing
independent $1/f$ noise to all channels being independent WGN. For algorithms
on generating these multivariate signals, see [5, 20, 9].

We applied the $\operatorname{mvDE}$ and $\operatorname{mvDE_{G}}$ to $40$
independent realizations of uncorrelated trivariate WGN and $1/f$ noise, each
with $15,000$ sample points. The mean and standard deviations (SD) of the
results for $\operatorname{mvDE}$ and $\operatorname{mvDE_{G}}$ are shown for
$m=4$ and $c=6$ in both cases in Figs. 1 and 2, respectively. These $m$ and
$c$ values are consistent with literature recommendations [4, 27].

Figure 1: Mean value and standard deviation (SD) for the classical
multivariate Dispersion Entropy $\operatorname{mvDE}$.

(a)

(b)

Figure 2: Entropy value comparison across scale factors for
$\operatorname{mvDE_{G}}$, the mean and SD with (a) identity adjacency graph,
indicating no prior inter-channel correlation, and (b) correlation matrix
graph, reflecting inherent inter-channel dependencies.

In both Figs. 1 and 2, the x-axis represents the Scale Factor, while the
y-axis delineates the normalised Entropy Value. The results from both methods
are strikingly similar, with the entropy values of trivariate WGN signals
being higher than those of other trivariate time series at low scale factors.
However, the entropy values for coarse-grained trivariate $1/f$ noise signals
are almost constant along the temporal scale factor, while those for coarse-
grained WGN signal monotonically decrease with increasing scale factors. This
indicates that multivariate WGN time series contain information primarily at
small temporal scale factors, whereas trivariate $1/f$ noise signals exhibit
higher complexity across multiple time series. These observations align with
literature findings [9].

In Fig. 2(a), we consider the case of a zero adjacency graph ($I_{p}$ is a
$3\times 3$ matrix with all zero entries) since the noise is uncorrelated and
theoretically, there is no interaction between time series. The results are
similar to $\operatorname{mvDE}$, reflecting the primarily temporal nature of
the information. However, the simulated time series exhibit slight
correlations in practice. When these correlations are calculated and used as
$I_{p}$ in Fig. 2(b), the entropy values of $\operatorname{mvDE_{G}}$ are
slightly lower, indicating a decrease in complexity due to the estimated
correlations. These results are consistent with the notion that $1/f$ noise is
structurally more complex than multivariate WGN.

In single-channel scenarios, multivariate entropy metrics reduce to their
univariate counterparts. Thus, classical Dispersion Entropy [27] and
Dispersion Entropy for Graph Signals [14] yield identical results. In
multivariate cases, while both $\operatorname{mvDE}$ and
$\operatorname{mvDE_{G}}$ produce different entropy values, they similarly
characterize dynamics. However, $\operatorname{mvDE_{G}}$ distinguishes itself
with its computational efficiency. Unlike $\operatorname{mvDE}$, where the
number of patterns increases exponentially with parameters $m$ and $p$ without
considering the graph structure, $\operatorname{mvDE_{G}}$ incorporates graph
structure linearly, significantly reducing the pattern count and computational
load. This efficiency, coupled with the ability to integrate topological
information, makes $\operatorname{mvDE_{G}}$ a superior choice for complex
data analysis, especially in multichannel scenarios, as in Sec. V.

###  IV-B Correlated Noise Analysis

In this section, we explore the behaviour of multivariate methods on
correlated noise. The experiment involves generating four sets of time series,
each with distinct correlation structures: (1) Four Uncorrelated WGN Series,
(2) Two Correlated and Two Uncorrelated Series (a mix of correlated and
uncorrelated series), (3) Two Pairs of Correlated Series (Two distinct pairs,
each internally correlated), (4) Three Correlated and One Uncorrelated Series
and (5) Four Correlated Series.

For each set of time series, we computed both $\operatorname{mvDE}$ and
$\operatorname{mvDE_{G}}$ using $N=500$, $m=4$, and $c=6$. The mean entropy
values and their SD were determined over $40$ simulations. In these
computations, the underlying graph ${G}$ was based on the theoretical
correlation matrix. When using the estimated correlation matrix $\hat{{G}}$,
derived from the data, the entropy values obtained for each simulation were
remarkably consistent, indicating the robustness of $\operatorname{mvDE_{G}}$
in capturing the dynamics of the time series even with small variations in the
graph structure.

Additionally, our focus on shorter time series ($N=500$) demonstrates
$\operatorname{mvDE_{G}}$’s robustness with limited data, a frequent challenge
where traditional methods struggle. This underscores the method’s adaptability
to scenarios with constrained data, as exemplified in Sec. VI. Comparable
results are observed with larger sample sizes.

Figure 3: Correlated noise entropy across various structures, ranging from
uncorrelated to fully connected across scale factors.

The results, as shown in Fig. 3, reveal interesting insights into the dynamics
of correlated noise. In both $\operatorname{mvDE}$ and
$\operatorname{mvDE_{G}}$ methods, the entropy values for WGN decrease as the
scale factor increases, indicating a reduction in complexity. However, a key
observation in mvDE is the overlapping of entropy values across all sets of
time series and scales, suggesting that mvDE struggles to differentiate
between various degrees of correlation and uncorrelated dynamics.

In contrast, our proposed method, $\operatorname{mvDE_{G}}$, not only presents
distinct mean entropy values for each set of time series but also avoids the
overlapping of standard deviations. This distinct separation in the entropy
values and their variances demonstrates that $\operatorname{mvDE_{G}}$ is more
adept at characterizing the dynamics of time series with varying degrees of
correlation. This capability makes $\operatorname{mvDE_{G}}$ a more robust
tool for analysing complex multivariate time series data, especially in
scenarios where understanding the interplay between correlation and signal
complexity is crucial.

Figure 4: Entropy for varying degrees of correlated noise.

Further, we explore the dynamics of correlated noise with varying degrees of
correlation. For this purpose, we generate three time series for each set,
with correlation values set at 0.95, 0.75, 0.55, 0.35, and 0.15. The
parameters used are $m=4$, $c=6$, and the analysis is based on 40
realizations. The results, depicted in Fig. 4, include both the mean and SD of
the entropy values.

The analysis reveals that while both methods show a decrease in entropy values
with increasing scale factor, the classical $\operatorname{mvDE}$ fails to
distinguish the dynamics influenced by varying degrees of correlation, as
evidenced by the overlapping SDs. In contrast, $\operatorname{mvDE_{G}}$
successfully differentiates between these dynamics, with minimal overlap in
the SDs across different correlation levels. This distinction is not only
pivotal for understanding complex signal interactions but also highlights the
superior computational efficiency of $\operatorname{mvDE_{G}}$, especially for
shorter time series. The following section will further elaborate on the
computational advantages of our proposed method.

##  V Computational Time Comparison

This section compares the computational efficiency of classical
$\operatorname{mvDE}$ and $\operatorname{mvDE_{G}}$. To assess the
computational performance, we utilized uncorrelated multivariate WGN time
series of varying lengths, ranging from $100$ to $10,000$ sample points,
across $10$ channels. The experiments were conducted on a PC running MATLAB
R2023b, equipped with 32.0 GB RAM and an Intel(R) Core(TM) i7-10610U CPU @
1.80GHz, 2.30 GHz processor.

Figure 5: Computational times comparison for $\operatorname{mvDE_{G}}$ against
classical $\operatorname{mvDE}$, for different data lengths.

The classical $\operatorname{mvDE}$ method exhibited significant limitations
in processing large channels. For a dataset with $\sim 2,000$ sample points,
$m=5$ and $p=8$ channels, MATLAB encountered an overflow RAM error. This
limitation is problematic in scenarios requiring the analysis of high-
dimensional data or multiple channels. In contrast, $\operatorname{mvDE_{G}}$
showcased computational efficiency under analogous conditions. It successfully
processed data without encountering the memory constraints that impeded the
classical method. As depicted in Fig. 5, $\operatorname{mvDE_{G}}$
consistently outperformed $\operatorname{mvDE}$ in terms of computational
time, especially in datasets with a larger number of sample points and
channels.

The classical $\operatorname{mvDE}$ method computes $(N-m+1)\binom{mp}{m}$
dispersion patterns. This combinatorial approach, while making
$\operatorname{mvDE}$ robust to very short data sequences by generating a
multitude of patterns even from a limited time series, also introduces a
computational complexity. The quantity of such dispersion patterns highlights
the exponential growth in the number with increasing $m$ and $p$, leading to
increases in both computational time and resource usage.

The superior performance of $\operatorname{mvDE_{G}}$ is primarily attributed
to its algorithmic design. The computational requirement of
$\operatorname{mvDE_{G}}$ is significantly lower, as it is limited to
processing at most $(N-m)p$ dispersion patterns and also the algorithm
optimizes the computation of $Np\times Np$ matrix powers. This difference in
computational requirements underscores the advanced algorithmic design of
$\operatorname{mvDE_{G}}$, making it a preferable choice for applications
involving large-scale multivariate time series analysis even if the results
are similar. For instance, in a dataset with $N=2000$, $p=8$, and $m=5$, the
classical $\operatorname{mvDE}$ method needs to compute approximately
$1.3\times 10^{9}$ dispersion patterns. In contrast, $\operatorname{mvDE_{G}}$
requires only about $1.6\times 10^{3}$ patterns, demonstrating a substantial
reduction in computational complexity.

This efficiency, combined with its robustness for short time series, positions
$\operatorname{mvDE_{G}}$ as a superior tool in time series analysis.

##  VI Real-World Applications

This section demonstrates the application of $\operatorname{mvDE_{G}}$ to
meteorological data and industrial two-phase flow systems. The complex nature
of weather data [16] and the dynamic characteristics of industrial flows [29]
provide ideal test to showcase the efficacy and adaptability of
$\operatorname{mvDE_{G}}$.

###  VI-A Weather Data on Ground Station Graphs

####  VI-A1 Data Description

We analysed temperature data from Brittany ground stations (January 2014),
provided by the French national meteorological service [17]. A graph was
constructed where vertices represent stations, and edges’ weights are based on
the Gaussian kernel of the Euclidean distances between stations. This setup
allows for a detailed examination of spatial interactions and temperature
variability across the region. The weighting formula between two vertices $i$
and $j$ is as follows:

$\mathbf{A}_{ij}=\begin{cases}\exp\left(\frac{-d(i,j)^{2}}{2\sigma_{1}^{2}}\right)&\text{if
}d(i,j)\leq\sigma_{2}\\\ 0&\text{otherwise.}\end{cases}$ (2)

where $d(i,j)$ denotes the Euclidean distance between stations $i$ and $j$,
with parameters $\sigma_{1}^{2}=5.1^{8}$ and $\sigma_{2}=10^{5}$, as per [17].

Figure 6: Geographical distribution of Brittany ground station, with selected
cities highlighted for detailed analysis. Figure 7: Representative time series
for wind velocity, temperature, and rainfall recorded at the Brittany stations
over a 30-day period, as identified in Fig. 6.

Fig. 6 depicts Brittany’s map with sensor locations as black nodes. The edges
are computed using Eq. 2. Three cities (vertices) are marked in different
colours (red, green, and blue), and their corresponding time series for wind
velocity, temperature, and rainfall are shown in Fig. 7. Each time series
consists of $744$ samples, representing hourly data for $31$ days. The entropy
values for $\operatorname{mvDE_{G}}$ were computed using the distance graph
for $m=4$, $c=6$, and the results are presented in Fig. 8.

Figure 8: Entropy of meteorological parameters: wind speed, maximum wind
speed, temperature, and rainfall.

####  VI-A2 Results and Discussion

The entropy analysis distinctly identifies patterns in meteorological data.
Rainfall data (purple) consistently show low entropy values across scales,
akin to the characteristics of $1/f$ noise, indicating a consistent, but low,
level of complexity. In comparison, temperature data (blue) display slightly
higher, yet stable entropy values, pointing to a different pattern of
complexity. Wind data (mean and max wind speed) demonstrate the greatest
complexity, with a decrease at lower scales.

$\operatorname{mvDE_{G}}$ effectively captures the varying characteristics of
meteorological elements, especially the complex and fluctuating patterns of
wind [23]. Influenced by factors like local terrain and atmospheric
conditions, wind patterns pose significant challenges for accurate prediction
and modelling. In contrast, temperature changes tend to be more predictable,
following daily and seasonal cycles. Rainfall, subject to variability,
generally follows recognizable weather patterns, showing less sudden changes
compared to the more unpredictable nature of wind.

The entropy analysis distinctly differentiates meteorological elements, with
higher entropy in wind data emphasizing the challenges in its prediction and
modelling – a critical factor for industries reliant on accurate wind
forecasts. In contrast, the more stable and lower entropy patterns in rainfall
data underscore $\operatorname{mvDE_{G}}$’s capacity to effectively delineate
the varied complexities of meteorological phenomena.

####  VI-A3 Comparison with other state-of-the-art methods

Despite the relatively small size of our dataset, consisting of data from $37$
cities each with $744$ samples, the computational demand of the classical
$\operatorname{mvDE}$ is substantial. Specifically, it requires the processing
of approximately $3.4\times 10^{11}$ patterns. This high computational load,
largely due to each city being treated as a separate channel, renders the
classical $\operatorname{mvDE}$ impractical for execution on standard personal
computers because it leads to memory errors. . In contrast, our
$\operatorname{mvDE_{G}}$ method drastically reduces this requirement to only
about $2.7\times 10^{4}$ patterns. This significant decrease in computational
complexity not only makes $\operatorname{mvDE_{G}}$ feasible for processing
large datasets on typical computing setups but also extends its applicability
to a broader range of real-world scenarios.

Due to the impracticality of applying classical $\operatorname{mvDE}$ to our
dataset, we resorted to using the univariate version of the Multiscale
Dispersion Entropy and Sample Entropy for analysing the weather data. We use
Sample Entropy due to its widespread application in data analysis [21]. This
approach limits our analysis to the temporal information of each sensor,
disregarding the spatial relationships between them. The mean entropy results
from both methods are presented in Fig. 9. While these methods are effective
in detecting the dynamics and yield lower entropy values for rainfall data,
they fail to distinguish the complex dynamics between temperature and wind
data. Additionally, Sample Entropy becomes undefined [4] for scale factors
$>6$.

(a)

(b)

Figure 9: Comparison of entropy analysis of weather data with state-of-the-art
methods: (a) Dispersion Entropy [27], (b) Sample Entropy [4].

###  VI-B Two-phase Flow Data

####  VI-B1 Data Description

We analyse data from two-phase flow experiments conducted with Electrical
Resistance Tomography (ERT), a technique that leverages the varying
conductivities of mediums to reconstruct the conductivity distribution within
a pipe. This is achieved by measuring boundary voltages between electrodes,
with a constant electrical current of 50 kHz serving as the excitation signal
and a data acquisition rate of 120 frames per second [30, 29]. The
$16-$electrode ERT setup captures $16\times 13=208$ voltage data points. To
manage the redundancy in these measurements, feature vectors $V_{Ri}$ are
extracted from each electrode, reducing the dimensionality of the data. This
process involves averaging the voltage measurements relative to a baseline
condition (when the pipe is full of water), calculated as
$V_{Ri}=\frac{1}{13}\sum_{j=1}^{13}(V_{ij}-V_{i{j_{0}}})/V_{i{j_{0}}}$, where
$V_{ij}$ is the measured voltage value and $V_{i{j_{0}}}$ is the baseline
voltage for electrode $i$ [30]. In the experiments, gas and water were mixed
at velocities ranging from $0.4\,\text{m/s}$ to $2.9\,\text{m/s}$ for water
and $0.06\,\text{m/s}$ to $5.64\,\text{m/s}$ for gas. This produced diverse
flow patterns, which includes $105$ distinct experiments conducted under
various gas and water flow rates, each generated around $1400$ time samples.

Two-phase flow experiments usually require data dimensionality reduction for
computational analysis, often by averaging feature vectors $V_{Ri}$ from a
16-electrode ERT system into four aggregated time series, as per [30]. In
contrast, $\operatorname{mvDE_{G}}$ processes the entire dataset directly,
eliminating the need for initial data reduction. This capability allows
$\operatorname{mvDE_{G}}$ to efficiently handle the full dataset, including
raw, unfiltered data, demonstrating its robustness against noise.

####  VI-B2 Results and Discussion

For our analysis, we applied the $\operatorname{mvDE_{G}}$ technique,
selecting $m=4$ and $c=6$. This choice of the graph was based on the fact that
each $V_{Ri}$ represents an unweighted average of nodes, suggesting the use of
a complete graph as the underlying structure. We found consistent results
across various underlying graphs derived from the data, as well as with
different $m$ and $c$ values, and even when analysing cleaned data. The
entropy values obtained for the six flow patterns were plotted against the
scale factor in Fig. 10.

A key element of our $\operatorname{mvDE_{G}}$ analysis is assessing entropy
variations across scales, offering vital insights into the dynamics of each
flow pattern. For example, the consistent entropy profile in the bubbly flow
regime suggests uniform complexity and predictability, while the entropy
shifts in slug and annular flows reflect changes in system complexity and
dynamics at various scales. These fluctuations are linked to alterations in
flow structure and phase interactions with changing scale factors.

Figure 10: Entropy values for the two-phase flow data.

#### Bubbly Flow

The bubbly flow regime is marked by the lowest entropy values, indicative of a
structured and predictable system. This pattern results from small bubbles’
limited impact on selected electrodes, intensified by gravity’s tendency to
concentrate bubbles at the top of the electrode array. Therefore, the sensors
exhibit regular readings, leading to lower entropy values [30]. As the scale
factor increases, the influence of these bubbles becomes more noticeable,
causing slight entropy variations. Our $\operatorname{mvDE_{G}}$ method
effectively captures this nuanced impact, especially in electrodes more
directly exposed to the gas phase. The overall flattened entropy profile
across scales in this regime reflects its stable and predictable nature.

#### Stratified Flow

Stratified flow is characterized by a low and stable entropy profile,
indicative of a system with reduced complexity. This flow regime, marked by a
clear separation between the gas and liquid phases, exhibits minimal
fluctuations across all $16$ electrodes. The $\operatorname{mvDE_{G}}$
method’s ability to delineate between closely related flow patterns,
particularly stratified and bubbly flows, is evident in lower scales where
traditional methods might falter [30]. This capability to separate stratified
flow from other flow patterns, even at lower scales, highlights the robustness
and precision of the $\operatorname{mvDE_{G}}$ method in analysing complex
two-phase flow dynamics.

#### Slug Flow

Slug flow displays distinctive entropy values, with higher values at lower
scales, indicative of the pronounced influence of large gas bubbles. This
complexity is most evident in the initial scales, as the large bubbles’
dynamics generate unique patterns. As shown in Fig. 10, a decreasing trend in
entropy values emerges beyond scale factor 10, reflecting a change in flow
complexity across scales. $\operatorname{mvDE_{G}}$ adeptly differentiates
slug flow from other flow regimes, particularly at lower scales, underscoring
its capacity to capture the intricate dynamics unique to this type of flow.
The entropy profile, with its initial rise and subsequent fall, mirrors the
varying influence of large gas bubbles within slug flow, revealing their
significant impact on the flow’s overall complexity.

#### Plug Flow

In contrast to slug flow, plug flow exhibits lower entropy values,
particularly at smaller scales. This is attributed to the complex dynamics of
smaller bubbles, which produce more variation and an increase in entropy
values as the scale factor increases. However, a notable transition occurs
around scale factor 10, where the entropy values for plug flow diverge
significantly from those of slug flow. This divergence is indicative of the
transition from plug to slug flow, characterized by the coalescence of smaller
bubbles into larger ones. The transition from plug flow to slug flow is
distinct from bubbly flow across all scales considered. With increasing gas
velocity, a complexity transition is observed at larger scales, and the
transition scale tends to shift towards lower scales. This shift, as depicted
in Fig. 10, suggests that scale 10 is particularly sensitive to the transition
from plug flow to slug flow, as also supported by similar findings in the
literature (see [30]).

#### Churn Flow

Churn flow, occupying an intermediate position between slug and annular flows,
is indicative of the system’s increasing complexity. This flow regime is
characterized by entropy values that surpass those of slug flow but remain
below the peaks of annular flow, embodying a transitional state with elements
of both adjacent patterns. [28] considered churn flow as a developing slug
flow, while [19] described it as intermediate between slug and annular flows.
This characterization aligns with our entropy analysis, where churn flow’s
entropy values are consistently higher than slug flow but lower than annular
flow across all scale factors. The churn regime, known for its chaotic nature
and complexity [3], is accurately reflected in the high entropy values
obtained from our analysis. This phase serves as a critical indicator of the
system’s transition towards the more chaotic nature characteristic of annular
flow.

#### Annular Flow

Annular flow is characterized by the highest entropy values, indicative of its
complex and turbulent nature. This regime involves a rapidly moving gas core
surrounded by a liquid film with embedded gas bubbles and droplets,
contributing to its high complexity [29]. In horizontal annular flow, gravity
causes an uneven liquid film distribution, creating surface waves. These
waves, disrupted by gas shear forces, produce droplets in the gas core, adding
to the regime’s intricacy [30]. $\operatorname{mvDE_{G}}$ effectively captures
the dynamic interplay between these elements, highlighting the chaotic
characteristics of annular flow. The entropy values initially rise and then
stabilize at higher scales, reflecting a transition from variable to more
consistent complexity within this regime.

Our $\operatorname{mvDE_{G}}$ analysis captures the unique entropy profiles
for each flow pattern in two-phase flow, providing insights into the
underlying dynamics. The method’s sensitivity to nuanced dynamics is evident
in the clear separation of entropy values for each regime, particularly at
lower scales. Additionally, $\operatorname{mvDE_{G}}$ allows us to utilize the
total information from all 16 electrodes without the need to reduce them to a
smaller set of vertices, thus preserving the richness of the data. This
capability contrasts with the $\operatorname{mvDE}$ approach, where attempting
to use all $16$ electrodes would lead to significant computational
difficulties. Indeed, even with only 8 channels and a similar size of the
sample of two-phase flow data, the state-of-the-art methods, as detailed in
Sec. V, have been observed to produce out-of-memory errors due to their
computational demands. In this way, $\operatorname{mvDE_{G}}$ not only offers
more detailed insights but also ensures computational feasibility when
handling extensive sensor data.

##  VII Conclusion

This paper presented the Multivariate Multiscale Graph-based Dispersion
Entropy ($\operatorname{mvDE_{G}}$), a new method for analysing multivariate
time series on graphs and networks. $\operatorname{mvDE_{G}}$ uniquely
combines temporal dynamics with topological structure, enhancing analysis
beyond traditional entropy methods.

$\operatorname{mvDE_{G}}$ has proven effective in differentiating complex
patterns in multivariate time series, excelling particularly with short series
where conventional methods falter. Its application to real-world datasets,
including weather and two-phase flow data, highlights its versatility and
robustness in practical scenarios.

A key innovation of $\operatorname{mvDE_{G}}$ is its computational efficiency,
a significant leap forward from existing methods. By utilizing matrix
properties and the Kronecker product, $\operatorname{mvDE_{G}}$ transitions
computational time growth from the exponential rates typical of classical
methods to a linear growth with respect to the number of vertices or nodes,
making it apt for extensive and real-time analyses.

In summary, $\operatorname{mvDE_{G}}$ stands as a powerful tool in
multivariate time series analysis, adept at unraveling complex spatial-
temporal interplays in various fields.

## References

  * [1] M. U. Ahmed and D. P. Mandic.  Multivariate multiscale entropy.  Phys. Rev. E - Stat. Nonlinear, Soft Matter Phys., 84(6):061918, Dec 2011.
  * [2] M. U. Ahmed and D. P. Mandic.  Multivariate multiscale entropy analysis.  IEEE Signal Processing Letters, 19(2):91–94, 2011.
  * [3] A. Ansari, N. Sylvester, C. Sarica, O. Shoham, and J. Brill.  A comprehensive mechanistic model for upward two-phase flow in wellbores.  SPE Production & Facilities, 9(02):143–151, 1994.
  * [4] H. Azami and J. Escudero.  Coarse-graining approaches in univariate multiscale sample and dispersion entropy.  Entropy, 20(2):138, 2018.
  * [5] H. Azami, A. Fernández, and J. Escudero.  Multivariate multiscale dispersion entropy of biomedical time series.  Entropy, 21(9):1–21, 2019.
  * [6] S. Boccaletti and et al.  Complex networks: Structure and dynamics.  Physics Reports, 424(4-5):175–308, 2006.
  * [7] J. A. Bondy.  Theory with Applications.  1982\.
  * [8] Z. Chai and C. Zhao.  Enhanced random forest with concurrent analysis of static and dynamic nodes for industrial fault classification.  IEEE Trans Industr Inform, 16(1):54–66, 2019.
  * [9] M. Costa, A. L. Goldberger, and C.-K. Peng.  Multiscale entropy analysis of complex physiologic time series.  Physical Review Letters, 89(6):068102, 2002.
  * [10] W. Dang and et al.  A novel deep learning framework for industrial multiphase flow characterization.  IEEE Trans Industr Inform, 15(11):5954–5962, 2019.
  * [11] K. Dong and D. Li.  An information theoretic approach to analyze irregularity of graph signals and network topological changes based on bubble entropy.  Fluctuation and Noise Letters, 2024.
  * [12] J. S. Fabila-Carrasco, C. Tan, and J. Escudero.  Multivariate permutation entropy, a cartesian graph product approach.  In 2022 30th European Signal Processing Conference (EUSIPCO), pages 2081–2085. IEEE, 2022.
  * [13] J. S. Fabila-Carrasco, C. Tan, and J. Escudero.  Permutation entropy for graph signals.  IEEE Trans. Signal Inf. Process. over Networks, 8:288–300, 2022\.
  * [14] J. S. Fabila-Carrasco, C. Tan, and J. Escudero.  Dispersion entropy for graph signals.  Chaos, Solitons & Fractals, 175:113977, Oct 2023.
  * [15] Q. Fan and et al.  Graph multi-scale permutation entropy for bearing fault diagnosis.  Sensors, 24(1):56, 2023.
  * [16] K. Fraedrich.  Estimating the dimensions of weather and climate attractors.  Journal of Atmospheric Sciences, 43(5):419–432, 1986.
  * [17] B. Girault.  Stationary graph signals using an isometric graph translation.  In 2015 EUSIPCO, pages 1516–1520. IEEE, 2015.
  * [18] A. Graham.  Kronecker Products and Matrix Calculus with Applications.  Courier Dover Publications, 2018.
  * [19] G. Hewitt.  Annular two-phase flow.  Elsevier, 2013.
  * [20] A. Humeau-Heurtier.  Multivariate generalized multiscale entropy analysis.  Entropy, 18(11):411, 2016.
  * [21] Z. Huo, M. Martínez-García, Y. Zhang, R. Yan, and L. Shu.  Entropy measures in machine fault diagnosis: Insights and applications.  IEEE Transactions on Instrumentation and Measurement, 69(6):2607–2620, 2020.
  * [22] R. H. Jhaveri, S. V. Ramani, G. Srivastava, T. R. Gadekallu, and V. Aggarwal.  Fault-resilience for bandwidth management in industrial software-defined networks.  IEEE Transactions on Network Science and Engineering, 8(4):3129–3139, 2021.
  * [23] M. Laib, J. Golay, L. Telesca, and M. Kanevski.  Multifractal analysis of the time series of daily means of wind speed in complex regions.  Chaos, Solitons & Fractals, 109:118–127, 2018.
  * [24] A. Ortega, P. Frossard, J. Kovačević, J. M. Moura, and P. Vandergheynst.  Graph signal processing: Overview, challenges, and applications.  Proceedings of the IEEE, 106(5):808–828, 2018.
  * [25] B. Prasse and P. Van Mieghem.  Network reconstruction and prediction of epidemic outbreaks for general group-based compartmental epidemic models.  IEEE Transactions on Network Science and Engineering, 7(4):2755–2764, 2020.
  * [26] K. Qiu, X. Mao, X. Shen, X. Wang, T. Li, and Y. Gu.  Time-varying graph signal reconstruction.  IEEE Journal of Selected Topics in Signal Processing, 11(6):870–883, 2017.
  * [27] M. Rostaghi and H. Azami.  Dispersion entropy: A measure for time-series analysis.  IEEE Signal Process. Lett., 23(5):610–614, May 2016.
  * [28] Y. Taitel, D. Barnea, and A. Dukler.  Modelling flow pattern transitions for steady upward gas-liquid flow in vertical tubes.  AIChE Journal, 26(3):345–354, 1980.
  * [29] C. Tan, H. Jia, G. Liang, X. Wang, W. Niu, and F. Dong.  Combinational multi-modality tomography system for industrial multiphase flow imaging.  IEEE Trans Instrum Meas, 2023.
  * [30] C. Tan, J. Zhao, and F. Dong.  Gas-water two-phase flow characterization with electrical resistance tomography and multivariate entropy analysis.  ISA Transactions, 55, 2015.
  * [31] X. Wang and et al.  Variational embedding multiscale diversity entropy for fault diagnosis of large-scale machinery.  IEEE Transactions on Industrial Electronics, 69(3):3109–3119, 2021\.
  * [32] Z. Wei, A. Pagani, G. Fu, I. Guymer, W. Chen, J. McCann, and W. Guo.  Optimal sampling of water distribution network dynamics using graph fourier transform.  IEEE Transactions on Network Science and Engineering, 7(3):1570–1582, 2019.
  * [33] F. Xia and et al.  Graph learning: A survey.  IEEE Transactions on Artificial Intelligence, 2(2):109–127, 2021\.
  * [34] R. Zhou and et al.  Edm-fuzzy: An euclidean distance-based multiscale fuzzy entropy technology for diagnosing faults of industrial systems.  IEEE Trans Industr Inform, 17(6), 2020.