Multiplex Recurrence Networks from multi-lead ECG data
MMultiplex Recurrence Networks from multi-lead ECG data
Sneha Kachhara a) and G. Ambika b) Indian Institute of Science Education and Research (IISER) Tirupati, Tirupati - 517507,India (Dated: 28 August 2020)
We present an integrated approach to analyse the multi-lead ECG data using the framework of multiplex recurrence networks (MRNs). We explore how their intralayer andinterlayer topological features can capture the subtle variations in the recurrence patternsof the underlying spatio-temporal dynamics. We find MRNs from ECG data of healthycases are significantly more coherent with high mutual information and less divergencebetween respective degree distributions. In cases of diseases, significant differences inspecific measures of similarity between layers are seen. The coherence is affected most inthe cases of diseases associated with localized abnormality such as bundle branch block.We note that it is important to do a comprehensive analysis using all the measures to arriveat disease-specific patterns. Our approach is very general and as such can be applied in anyother domain where multivariate or multi-channel data are available from highly complexsystems.PACS numbers: 05.45.-a, 05.45.Ac, 87.18.-h, 05.90.+m, 05.45.Tp a) Electronic mail: [email protected] b) Electronic mail: [email protected] a r X i v : . [ phy s i c s . s o c - ph ] A ug he Electrocardiogram (ECG) is a record of the electrical activity of the heart in the formof a time series. The study of cardiac dynamics through ECG has gathered a lot of attentionin the nonlinear dynamics community, with attempts to justify the degree of chaos in thissystem as well as to identify anomalies in the case of a disease. Nonlinear time series meth-ods, sometimes coupled with Machine Learning approaches, have been employed to this endwith reasonable success. However, the interpretation of feature based classification studiesin terms of underlying dynamics is not explored, except for some particular ailments suchas Arrythmias and Chronic Heart Failure. Moreover, the ECG data comes as highly corre-lated multivariate data, from 3, 5 or 12 leads. Most of the study from the dynamics pointof view till now are on the average data over leads or on a single lead. In the present workwe aim to study multi-lead ECG within the framework of Multiplex Recurrence Networks,which highlights spatio-temporal features of the cardiac dynamics as reflected in ECG. Tothis end, we employ layer similarity/dissimilarity measures in addition to the standard com-plex network measures defined for multiplex complex networks. We include three levels ofstructural aspects of MRNs: the coarse structure in terms of links across layers, the inter-layer features in degree distributions, and the local micro-structures using local clusteringcoefficients. We show that the cardiac dynamics in the case of a disease manifests abnormal-ities in a multitude of ways and can be understood only by consolidated results from a set ofmeasures.I. INTRODUCTION Research related to chaotic behaviour in the cardiovascular system has seen consistent progressin the last two decades . While the initial attempts were aimed at establishing the presence ofdeterministic chaos in the cardiac system, techniques of nonlinear time series analysis and ma-chine learning are now being applied to understand cardiac dynamics . The focus has shifted tothe identification of signatures of altered dynamics (chaotic or not) in patients as compared tohealthy , leading to good progress in machine learning based diagnostics and early-warning tools.However, this approach relies heavily on the availability of huge databases of quality data to trainthe algorithms and provides little insight into the underlying dynamics that reflects the intricaciesrelated to cardiac malfunctions. 2he physiologically relevant signal in the context of the heart is the Electrocardiogram (ECG)which records the electrical activity of the heart . The data on Heart Rate Variability (HRV) onthe other hand, is the number of cardiac cycles per minute. ECG and HRV, both reflect related butdifferent aspects of the cardiovascular dynamics . While HRV data has the advantage that it iseasy to obtain, even with simple wearable devices, the full ECG waveform, on the other hand, isvery sensitive to noise, requires specific devices to record and is usually available only for shortduration. However, the ECG contains more information about the dynamics, and also at differenttimescales. As such, analyzing ECG from a dynamical systems’ perspective can be quite reward-ing and relevant . Studies in this direction indicate reduction of complexity in some diseases suchas CHF and Arrythmias , however the question of how exactly the dynamics is altered in thecase of specific diseases, and measures that reflect these abnormalities from a dynamical point ofview, remain relatively unexplored.We note that limitations such as short duration and non-stationarity may have restricted theapplication of tools of Nonlinear time series analysis to ECG data. In this context the methodof Recurrence Networks (RNs) proves to be useful in analysing short and non-stationary data. Ittransforms the recurrence pattern in the reconstructed dynamics from a given time series into acomplex network . In this context we have recently reported bimodality in the degree distribu-tion and scaling of link density with recurrence threshold as characteristic features of RNs fromECG .The multivariate data taken from a complex dynamical system has properties that can be un-covered only with a comprehensive approach of analysis. In the case of a multivariate time series,Eroglu et al. have proposed the framework of Multiplex Recurrence Networks in which layersof the network correspond to the different time series of the data. The patterns of connectionsinside each layer are governed by the dynamics as reflected in the corresponding time series ofthe original data. The framework is very recent but has been successful with coupled map latticesand regime changes in palaeobotanical data , oil-water spatial flow , and even for self-reportsof human experience (EMA and ESM data) .In this context, the multi-lead ECG provides multivariate data as it is recorded from the differ-3nt leads as discrete time series . There has been some progress in developing diagnostic toolsbased on 3-lead vectorcardiogram (VCG), which is constructed from the 12-lead ECG and pro-vides a succinct representation of phase relationships. However, in clinical practice 12-lead ECGis preferred by cardiologists to make diagnosis based on some or all of the leads, depending onthe nature of the disease.We start with the hypothesis that the lead-to-lead variations and the similarity underlying databetween pairs of leads are to be studied to understand the complexity of the spatio-temporal dy-namics of normal heart. Moreover, such a study can give information on any anomaly appearingin a specific set of leads which can be an evidence of a specific disorder. We show how this isfeasible by using the framework of multiplex recurrence networks (MRNs). We use the ECG datafrom the six precordial leads placed closest to the heart as the multivariate data. MRNs can beconstructed from multi-lead ECG such that recurrence pattern in the reconstructed dynamics withtime series from each lead form one layer of the multiplex, with correspondence between nodesacross layers that relate to the same instance in time in the reconstructed dynamics. The analysisof interlayer similarities in such MRNs can uncover patterns not apparent in a single layer alone,as well as provide deeper insight into localized anomalies.We begin by describing briefly the construction of RNs and hence MRNs from the 6 layersof RNs, and the network quantifiers defined on them. We study the MRNs at different levels ofcomplexity; from individual links in layers to the overall topological features. We base the analy-sis on three quantifiers: links, quantified by average edge overlap; degree distributions comparedbetween pairs of different layers; and distribution of simple cliques, quantified by average localclustering coefficients in different layers. We illustrate how specific diseases can cause differencesin measures of similarity between particular pairs of layers, reflecting the intricate variations inthe underlying spatio-temporal dynamics. II. DATA AND PROCESSING
We use the PTB dignostics database from Physionet . In total, data from 125 subjects areused in our analysis, 51 of which are healthy (HC) and the rest of them are from one of the four4iseases: Bundle Branch Block (BB, 15), Cardiomyopathy/Heart Failure (CM, 16), Dysrhythmia(DR, 14) and Myocardial Infarction with no secondary diagnosis (MI, 29). Each record originallyconsists of 12-lead ECG, mostly 60 seconds in duration, sampled at 1000 Hz (60,000 points). Wechoose data from the precordial leads v to v which are placed closest to the heart, for analysis.Each data is pre-processed ; filtered with 0.5-40 Hz, normalized as 0 to 1 and binned to 5000points. III. CONSTRUCTION OF MULTIPLEX RECURRENCE NETWORKS
The construction of Multiplex Recurrence Networks (MRNs) is done by first constructing therecurrence network from the time series of each lead. This procedure involves two main steps:embedding the time series in a phase space based on Taken’s embedding theorem , and identi-fying each point on the reconstructed phase space trajectory (or attractor) as node of the network,making connections between them based on their proximity in recurrences . The reconstructionof the attractor in phase space requires a delay time τ , and embedding dimension m that are tobe chosen appropriately for the data under analysis . For ECG data from each lead, we fixembedding dimension as 4 using the method of False Nearest neighbours (FNN) . The delaytime τ , i.e. the time at which auto correlation falls to 1/e, is chosen as the minimum value of such τ among the time series from the 6 leads . Then a recurrence threshold ( ε ) is chosen to constructthe links between nodes and to arrive at the adjacency matrix of the recurrence network . Thevariation of network measures with recurrence threshold itself is an indicator of the complexityof the underlying attractor, as reported recently . In the present work we set it to ε = 0.1 for thesake of uniformity, based on previous study.For the chosen ε , the recurrence matrix R is constructed such that if two points i and j on thereconstructed attractor lie within distance ε of each other, we set the corresponding matrix element R i j to be 1, and 0 otherwise. i.e. R i j = Θ (cid:0) ε − (cid:13)(cid:13) (cid:126) v i − (cid:126) v j (cid:13)(cid:13)(cid:1) (1)where (cid:126) v i and (cid:126) v j are the corresponding vectors of i and j in the phase space, and Θ is the Heavisidestep function . 5 a)(b) FIG. 1: Reconstructed attractors in 2-d (top panel) and the corresponding recurrence networks(bottom panel) from two leads 2 (left) and 4 (right) for two typical datasets, (a) healthy and (b)BB. The color of a node indicates its degree.6he adjacency matrix A of a Recurrence Network is obtained from R as: A i j = R i j − δ i j (2)where δ i j is the Kronecker delta function that is inserted to avoid self-links in the network. Thisconstruction results in an unweighted and undirected network of size N, where N is the number ofpoints on the reconstructed attractor.We display in figure 1 the reconstructed attractors and the corresponding recurrence networksfrom two different leads for two typical datasets, one of healthy and one of BB. The differencesbetween leads is obvious and more pronounced in the case of a disease. To capture such differ-ences, we construct the MRN by treating RN from each lead as a single layer of a multilayernetwork, by connecting the nodes in different layers that correspond to the same instant of time inthe reconstructed attractors. Thus constructed, the MRN can be represented by a supra-adjacencymatrix M , which consists of blocks of adjacency matrices from different layers as the diagonalelements and identity matrices as off-diagonal elements: M = A [ ] I N · · · I N I N A [ ] . . . ...... . . . . . . ... I N · · · I N A [ m ] (3)Thus it encapsulates properties of the multivariate data effectively, with display of its complexfeatures in the interlayer attributes. For a typical dataset, we show ECG data from each lead,the supra-adjacency matrix and the constructed MRN in figure 2. We analyze the MRNs withprimary focus on interlayer similarities and differences. The quantifiers used in the study arebriefly described in the next section and more details can be found elsewhere . IV. MEASURES FROM MULTIPLEX RECURRENCE NETWORKS
We expect a rich structure in multilayer networks because of both interlayer and intralayerconnections. For multiplex networks, there is node-to-node correspondence by construction thatgrants them an additional level of order. Keeping this ordered structure in mind, we select a fewmeasures that can capture their properties. For comparing RNs from two layers, we can considerfeatures at two different levels: local, like degree distribution and global, such as link density7IG. 2: Construction of Multiplex Recurrence Networks from ECG data. The upper panel showsECG time series and the corresponding supra-adjacency matrix, M . The lower panel shows MRNfrom the matrix. For clarity, the time series and M are shown only for the duration of 10 seconds,and the MRN for 100 nodes. The actual calculations are performed on entire time series, andcorresponding MRN has 5000 nodes in each layer.(LD) . There are other measures such as the average clustering coefficient (CC) which incor-porate both local and global structure . In addition to these, in the case of multiplex networks,we can effectively use a few other measures defined below which quantify the extent of similarityamong the layers. 8he average edge overlap, denoted by ω , is a consolidated quantifier of the overlap of links (oredges) considering all layers of the network. First proposed by Lacasa et al. , ω gives an estimateof the expected fraction of layers containing a link, defined as: ω = ∑ i ∑ j > i ∑ l A [ l ] i j m ∑ i ∑ j > i ( − δ , ∑ l A [ l ] i j ) (4)where δ i j stands for the Kronecker delta function, m is the number of layers, and A [ l ] i j correspondsto the layer l as defined in eq. (2). For a network with m layers, ω by definition, can take valuesin the range [1/m,1], ω being 1 only if links are identical in all layers. Thus for MRNs from ECGwith 6 layers, ω will be in [ , 1].The cosine similarity (CS) is based on the inner product of two vectors and in the presentcontext, in terms of vectors of the degrees of nodes D [ l ] and D [ l ] in two respective layers l and l , as: CS = D [ l ] . D [ l ] (cid:13)(cid:13) D [ l ] (cid:13)(cid:13) (cid:13)(cid:13) D [ l ] (cid:13)(cid:13) (5)CS is very useful in comparing similarities between layers, node-by-node .The index of dissimilarity (ID) calculates differences between two distributions . The distri-butions considered here are that of local clustering coefficients . For two layers l and l , ID isdefined as: ID = N N ∑ i = (cid:12)(cid:12)(cid:12) C [ l ] i − C [ l ] i (cid:12)(cid:12)(cid:12) (6)where C [ l ] i and C [ l ] i are the local clustering coefficients for the respective layers and N is the totalnumber of nodes.The Jensen-Shannon divergence is a measure to discern two distributions based on the entropyof mixing . In the context of multiplex networks, we can use the concept of Jensen-Shannondistance ( JSD ) to assess the similarity of the probability degree distributions in different layers.For two layers l and l , with probability degree distributions P ( k l ) and P ( k l ) respectively, wehave: JSD ( P ( k [ l ] ) || P ( k [ l ] )) = (cid:112) D ( P ( k [ l ] ) || M ) + D ( P ( k [ l ] ) || M ) M = P [ l ] + P [ l ] , is the point-wise mean and D represents Kullback-Leibler divergence .The interlayer mutual information is calculated on the respective degree distributions asfollows : I l l = ∑ k [ l ] ∑ k [ l ] P ( k [ l ] , k [ l ] ) ln P ( k [ l ] , k [ l ] ) P ( k [ l ] ) P ( k [ l ] ) (8)where P ( k [ l ] , k [ l ] ) is the joint distribution of probability of degrees k [ l ] in layer l and degree k [ l ] in layer l . I l l captures the topological similarity in the two layers.In addition, we use the Pearson Correlation Coefficient to compare distributions of local clus-tering coefficients in two layers. Each of the measures described above except ω can be averagedover all pairs of layers, which can then be used to compare two MRNs. Since the structure of MRNreflects the underlying spatio-temporal dynamics of the cardiac system, any specific variations inmeasures from two layers and statistically significant changes in a property across subjects of aclass will help to understand how diseases can alter heart dynamics and functions. V. AVERAGE MEASURES IN MRNS
The one-to-one correspondence of nodes between layers of multiplex networks enables us tocompare structures across individual layers on a very basic and concrete level, that of the indi-vidual links. In case of MRNs, links across layers can be associated with the recurrences in theunderlying dynamics that occur synchronously, making them even more relevant. We computethe average edge overlap, ω (eq. 4) for the different datasets and the average link density for the6 layers of MRN from every data set. The results for all the datasets used in the study are shownin figure 3. Each circle represents one subject, and the size of the circle is proportional to thevariance in link density (LD) across layers.We note from figure 3 that there is an apparent positive correlation between ω and average LD.This correlation is further explored using the correlation coefficient for all datasets of a category,as shown in table I. This positive correlation is not surprising since a network with a very high LDnaturally has high ω , merely because of increased number of links in all layers. However, their10IG. 3: Average edge overlap and average link density across layers of MRNs constructed fromECG data for different subjects. Each circle represents a subject and its size is proportional to thevariance in link density across layers. The colors represent different classes of patients, withhealthy in green.exact relationship depends on how the links are distributed inside the layers. We note that, thecorrelation is low for DR, and high for MI and BB.The results suggest that the MRNs from ECG in general have high degree of association amonglayers, leading to high values for ω . We note almost all the values for ω are higher than 0.4, muchhigher than the lower threshold of = . ω are high (as compared to a healthy person withthe same average LD). This leads to the conjecture that there are disproportionate changes in someof the layers in these cases. We will explore this possibility further with the interlayer similarity11 ategory edge overlap ( ω ) Average LD Correlation CoefficientHealthy 0.5845 ± ± ± ± ± ± ± ± ± ± TABLE I: Values of average edge overlap and average link density, for each category of ECGdata; along with correlation. The errors indicate standard deviation across subjects of a category.measures in section VI.
VI. INTER-LAYER SIMILARITIES OF DEGREES
In this section we discuss the measures to relate the degree distributions across layers. We firstdiscuss how degree of each node differs from layer to layer as captured by the Cosine Similarityor CS , and then, how the overall distribution of degrees differs from layer to layer through Jensen-Shannon Distance, JSD and Mutual Information, I .The CS by definition, is a local measure and captures node-to-node differences across layers.We compute the CS values taking leads pairwise as (1,2), (1,3) etc. for each dataset. We present theaverage values of CS computed pair wise for data of each category in figure 4, along with standarderror. We find that for layer pairs (1,2), (1,4), (2,4), (2,6) and (4,6), the values from healthy are thelowest and in general BB has the highest value in all layers.Similarly the values for JSD and I are calculated pairwise and presented in figure 5 (a) and (b)respectively. The values of I for pairs of layers (1,2), (1,4), (2,4), (2,6) and (4,6) are highest forhealthy.All the three measures are averaged over all pairs in each data set and then their means andstandard deviations for all datasets in each category are tabulated in table II. We see that JSD ishigher for BB and CM, while for DR and MI it stays close to the corresponding range for healthy.As for I , healthy have a high value for most of the pairs of layers.12IG. 4: Cosine Similarity of degrees among pairs of layers of MRNs constructed from ECG data.Each point represents the average value of CS for the corresponding pair of layers, across subjectsof a category. The colors indicate different categories as before with healthy in green. Theerrorbars indicate standard error. Category CS
JSD I
Healthy 0.1626 ± ± ± ± ± ± ± ± ± ± ± ± ± ± ± TABLE II: Interlayer similarity measures CS , JSD and I in degree distributions, averaged acrossevery pair of layers within each category of data sets. The errors indicate standard deviationacross subjects of the category.13 a)(b) FIG. 5: (a) Jensen-Shannon distance and (b) Mutual Information from degree distributions amongpairs of layers of MRNs constructed from ECG data. Each point represents the average value forthe corresponding pair of layers, across subjects of a category. The respective colors correspondto different categories as before, with healthy in green. The error bars indicate standard error.14IG. 6: Values of average Clustering Coefficient (CC) and average link density of differentsubjects. Each circle represents a subject and its size indicates variance in CC across layers forthat subject. The colors are representative of the category the subject belongs to, with healthy ingreen.
VII. INTERLAYER SIMILARITY AS REFLECTED IN LOCAL CC
In this section, we compute the average value of local Clustering Coefficients (CC) of eachlayer of the MRN in order to analyze the micro-structure. We take the average CC for all thelayers and plot them against average LD for every dataset, as shown in the figure 6. Each circlein the figure represents a subject and its size is proportional to the variance in average CC acrosslayers for that subject. We note that the magnitude of variance is low for healthy subjects. Also,the average CC and LD are correlated, as reflected in the correlation coefficient summarized intable III. Specifically, for DR, we observe that the values lie significantly away from the maindiagonal and hence differ from healthy and MI.Now to depict the dissimilarity in interlayer topology, we compute the Index of Dissimilarity15 a)(b)
FIG. 7: (a) Index of Dissimilarity ( ID ) and (b) Pearson Correlation Coefficient for local CCdistributions among pair of layers of MRNs constructed from ECG data. Each point in the plotrepresents average value of the corresponding measure for the pair of layers across differentsubjects of a category (color-coded by category).16 ategory ID Correlation CoefficientHealthy 0.0395 ± ± ± ± ± ± ± ± ± ± TABLE III: Interlayer similarity measures in local CC distributions, averaged across every pair oflayers of MRNs for each category. The errors indicate standard deviation across subjects of acategory.( ID ), which is a measure of the differences in the two distributions of local CC values in twolayers. The results are shown in figure 7 (a). We can also compare how correlated the distributionsare, as depicted using Pearson Correlation Coefficients in figure 7 (b). The average value for eachcategory over all datasets for all pairs of layers is presented in table III. We find that in general,ID is higher than healthy for all diseases, but most prominently so for CM. Correspondingly, thecorrelation is also low for CM. VIII. VARIATION IN MEASURES AMONG PAIRS OF LAYERS
The different measures computed for MRNs from multi-lead ECG data and presented in theabove sections are further analysed statistically and consolidated together in this section. For thiswe use the Welch’s t- test to compute a significance value for every measure for each categoryand pair of layers, so that we can understand how significant the differences in computed measuresare for each disease from healthy. The results are summarized in the form of significance matricesin figure 8. Each entry in a given matrix represents the significance value for that category ascompared to the corresponding measure for healthy, for the pair of layers indicated. A p-value< 0.05 is color coded (p > 0.05 is white) from blue to green in decreasing order such that greenindicates the least p-value, corresponding to the most significance.From the figure, we can conclude that different classes of diseases have different types of17 a) CS (degrees) (b) JSD (degree distribution)(c) I (degree distribution) (d) ID (CC) FIG. 8: Significance of computed measures for different pair of layers of the MRNs of eachdisease, as compared to healthy. Each off-diagonal matrix element corresponds to a particularpair of layers as indicated by its row and column number. The color code is as per the p-valuecomputed from Welch’s t-test. Statistically most significant differences are shown in green (lessso in blue), and most insignificant in white.variations in the cardiac dynamics and not all measures show differences for all diseases. For BBwe observe isolated elements in the matrix of significance for measures CS , JSD and ID . For CM,there seems to be a gradual change across adjacent layers, reflected strongly in JSD . For DR,there are few elements in green in all measures except in I , and for MI, significant difference isseen for multiple elements in both I and ID . No pair of layers show any significant differencein JSD for MI, hence the corresponding matrix is completely white. These results indicate that a18isease like BB could be localized, hence some layers are affected more than others, while DRand MI manifest in an integrated manner affecting all the leads. Since each measure encapsulatessimilarity of a different kind, we expect them to differ from disease to disease.
IX. CONCLUSIONS
The spatio-temporal features of the electrical excitation patterns of heart are complex and themultiplex recurrence networks is the ideal framework for understanding variations in their com-plexity. In this study, we analyse the multi-lead ECG data from 125 subjects by constructingmultiplex recurrence networks (MRNs) in order to discern patterns of variations in the underlyingcardiac dynamics. This includes 51 healthy subjects, and disease cases like Bundle Branch Block(BB), Cardiomyopathy (CM), Dysrhythmia (DR) and Myocardial Infarction (MI).The measures specific to multiplex networks, such as Edge overlap ( ω ) and Mutual Information( I ) are highly relevant for ECG data as they can highlight features of cardiac dynamics obscuredby inherent non-linearity and correlations in data. Measures such as Cosine Similarity ( CS ) andJensen-Shannon Distance ( JSD ) can provide additional insight into the topological differencesacross layers. Moreover, the differences in distribution of local clustering coefficients (CC), cap-tured by the Index of Dissimilarity ( ID ) and correlation coefficient, can provide more informationon the interlayer similarities in their micro-structures.The results on ω and average link density establish that MRNs from ECG have high similarityfrom layer-to-layer at the most basic level. ω for healthy is found between 0.4-0.8, and is gener-ally higher as compared to patients for the same average link density. Moreover, most cases ofpatients with high ω , show large variation in link density from layer to layer. Our results illustratethat healthy cardiac dynamics has less range of variations across layers, as compared to diseases.Most extreme cases are those of BB and CM, while DR and MI are mostly within the range forhealthy. The extreme values observed could be either due to some of the layers being very similar,or vice-a-versa. The value of CS , which measure variations in degrees of nodes across layers, hasthe highest value for BB, as compared to other diseases and healthy.19e extend the study to interlayer similarities as reflected in degree distributions and distribu-tions of local CCs. We find that the JSD is consistently high for BB and CM among all pairs oflayers, while for DR and MI the values are closer to healthy. On the other hand, I is highest forhealthy, on the average across all pairs of layers. In particular, we observe that the layer 4 differsmost from other layers in the case of BB and CM. For DR and MI, no particular layer or pair oflayers stand out in terms of differences in degree distributions. We note that the healthy do notdiffer much from layer to layer in terms of average local CCs. Thus, we can infer that there isan overall coherence in the healthy cardiac system which is hampered in case of a disease. Thereare some specific differences in the way different abnormalities manifest in cardiac dynamics. Alocalized anomaly, such as that of BB for example, will affect only some of the layers but in caseof MI, all the layers show significant differences.Our study is aimed at exploring the nature of variations in the underlying spatio-temporal dy-namics due to any type of disease as revealed through the measures computed from the multi-leadECG data. Such an understanding of the dynamics is basic to a proper knowledge of cardiacsystem and can provide insight for intelligent algorithms. Further studies in this direction cancombine deep learning approaches with dynamical systems theory. Moreover, the type of analysispresented here can be applied to multivariate data in other domains where traditional approachesof data analysis are insufficient. The framework of MRNs is not limited to time series data, as themeasures can be employed to analyze similarities in any real-world multiplex networks and it willbe interesting to see applications in more real data-based networks. X. DATA AVAILABILITY
The data that support the findings of this study are openly available in PTB diagnostics ECGdatabase at https://doi.org/10.13026/C28C71, ref. . ACKNOWLEDGMENTS
One of the authors, Sneha Kachhara, acknowledges financial support from the Council of Sci-entific and Industrial Research (CSIR), India through Senior Research Fellowship.20
EFERENCES L. Glass, “Introduction to controversial topics in nonlinear science: Is the normal heart ratechaotic?” Chaos (2009). E. M. Cherry, F. H. Fenton, T. Krogh-Madsen, S. Luther, and U. Parlitz, “Introduction to focusissue: Complex cardiac dynamics,” Chaos (2017). R. Acharya, S. Krishnan, J. A. Spaan, and J. Suri,
Advances in cardiac signal processing (Springer, 2007). T. B. Garcia, (Jones & Bartlett Publishers, 2013). U. R. Acharya, K. P. Joseph, N. Kannathal, C. M. Lim, and J. S. Suri, “Heart rate variability: areview,” Medical and Biological Engineering and Computing , 1031–1051 (2006). S. M. Shekatkar, Y. Kotriwar, K. Harikrishnan, and G. Ambika, “Detecting abnormality in heartdynamics from multifractal analysis of ECG signals,” Scientific reports , 1–11 (2017). S. Banerjee, S. K. Palit, S. Mukherjee, M. Ariffin, and L. Rondoni, “Complexity in congestiveheart failure: A time-frequency approach,” Chaos , 033105 (2016). Z. Qu, “Chaos in the genesis and maintenance of cardiac arrhythmias,” Progress in biophysicsand molecular biology , 247–257 (2011). Z. Qu, G. Hu, A. Garfinkel, and J. N. Weiss, “Nonlinear and stochastic dynamics in the heart,”Physics reports , 61–162 (2014). V. Tuzcu, S. Nas, T. Börklü, and A. Ugur, “Decrease in the heart rate complexity prior to theonset of atrial fibrillation,” Europace , 398–402 (2006). Y. Zou, R. V. Donner, N. Marwan, J. F. Donges, and J. Kurths, “Complex network approachesto nonlinear time series analysis,” Physics Reports , 1–97 (2019). R. V. Donner, Y. Zou, J. F. Donges, N. Marwan, and J. Kurths, “Recurrence networks—a novelparadigm for nonlinear time series analysis,” New Journal of Physics , 033025 (2010). S. Kachhara and G. Ambika, “Bimodality and scaling in recurrence networks from ECG data,”EPL (Europhysics Letters) , 60004 (2019). D. Eroglu, N. Marwan, M. Stebich, and J. Kurths, “Multiplex recurrence networks,” PhysicalReview E , 012312 (2018). Z.-K. Gao, W.-D. Dang, Y.-X. Yang, and Q. Cai, “Multiplex multivariate recurrence networkfrom multi-channel signals for revealing oil-water spatial flow behavior,” Chaos , 035809(2017). 21 F. Hasselman and A. M. Bosman, “Studying complex adaptive systems with internal states: Arecurrence network approach to the analysis of multivariate time-series data representing self-reports of human experience,” Frontiers in Applied Mathematics and Statistics , 9 (2020). S. Man, A. C. Maan, M. J. Schalij, and C. A. Swenne, “Vectorcardiographic diagnostic & prog-nostic information derived from the 12-lead electrocardiogram: Historical review and clinicalperspective,” Journal of electrocardiology , 463–475 (2015). R. Bousseljot, D. Kreiseler, and A. Schnabel, “Nutzung der ekg-signaldatenbank cardiodat derptb über das internet,” Biomedizinische Technik/Biomedical Engineering , 317–318 (1995). D. Kreiseler and R. Bousseliot, “Automatisierte ekg-auswertung mit hilfe der ekg-signaldatenbank cardiodat der ptb,” Biomedizinische Technik/Biomedical Engineering , 319–320 (1995). A. L. Goldberger, L. A. Amaral, L. Glass, J. M. Hausdorff, P. C. Ivanov, R. G. Mark, J. E.Mietus, G. B. Moody, C.-K. Peng, and H. E. Stanley, “Physiobank, physiotoolkit, and physionet:components of a new research resource for complex physiologic signals,” circulation , e215–e220 (2000). S. K. Berkaya, A. K. Uysal, E. S. Gunal, S. Ergin, S. Gunal, and M. B. Gulmezoglu, “A surveyon ECG analysis,” Biomedical Signal Processing and Control , 216–235 (2018). H. Kantz and T. Schreiber,
Nonlinear time series analysis , Vol. 7 (Cambridge university press,2004). G. Ambika and K. Harikrishnan, “Methods of nonlinear time series analysis and applications: Areview,” in
Dynamics and Control of Energy Systems (Springer, 2020) pp. 9–27. K. H. Kraemer, R. V. Donner, J. Heitzig, and N. Marwan, “Recurrence threshold selectionfor obtaining robust recurrence characteristics in different embedding dimensions,” Chaos ,085720 (2018). R. V. Donner, Y. Zou, J. F. Donges, N. Marwan, and J. Kurths, “Ambiguities in recurrence-basedcomplex network representations of time series,” Physical Review E , 015101 (2010). R. Hegger, H. Kantz, and T. Schreiber, “Practical implementation of nonlinear time series meth-ods: The tisean package,” Chaos , 413–435 (1999). S. Schinkel, O. Dimigen, and N. Marwan, “Selection of recurrence threshold for signal detec-tion,” The european physical journal special topics , 45–53 (2008). R. Jacob, K. Harikrishnan, R. Misra, and G. Ambika, “Uniform framework for the recurrence-network analysis of chaotic time series,” Physical Review E , 012202 (2016).22 N. Marwan, M. C. Romano, M. Thiel, and J. Kurths, “Recurrence plots for the analysis ofcomplex systems,” Physics reports , 237–329 (2007). P. Bródka, A. Chmiel, M. Magnani, and G. Ragozini, “Quantifying layer similarity in multiplexnetworks: a systematic study,” Royal Society open science , 171747 (2018). L. Lacasa, V. Nicosia, and V. Latora, “Network structure of multivariate time series,” Scientificreports , 15508 (2015). M. Newman,
Networks: an introduction (Oxford university press, 2010). L. Pardo,
Statistical inference based on divergence measures (CRC press, 2018). J. L. Myers, A. Well, and R. F. Lorch,