^{1}Institute for Bioengineering of Catalonia, The Barcelona Institute of Science and Technology, Barcelona, Spain^{2}Automatic Control Department (ESAII), Barcelona East School of Engineering (EEBE), Universitat Politècnica de Catalunya, Barcelona, Spain^{3}Institute of Innovative Health Technologies, Ernst-Abbe-Hochschule Jena, Jena, Germany^{4}Centro de Investigación Biomédica en Red de Bioengenieria, Biomateriales y Nanomedicina, Madrid, Spain

Cardiovascular diseases are one of the most common causes of death; however, the early detection of patients at high risk of sudden cardiac death (SCD) remains an issue. The aim of this study was to analyze the cardio-vascular couplings based on heart rate variability (HRV) and blood pressure variability (BPV) analyses in order to introduce new indices for noninvasive risk stratification in idiopathic dilated cardiomyopathy patients (IDC). High-resolution electrocardiogram (ECG) and continuous noninvasive blood pressure (BP) signals were recorded in 91 IDC patients and 49 healthy subjects (CON). The patients were stratified by their SCD risk as high risk (IDC_{HR}) when after two years the subject either died or suffered life-threatening complications, and as low risk (IDC_{LR}) when the subject remained stable during this period. Values were extracted from ECG and BP signals, the beat-to-beat interval, and systolic and diastolic blood pressure, and analyzed using the segmented Poincaré plot analysis (SPPA), the high-resolution joint symbolic dynamics (HRJSD) and the normalized short time partial directed coherence methods. Support vector machine (SVM) models were built to classify these patients according to SCD risk. IDC_{HR} patients presented lowered HRV and increased BPV compared to both IDC_{LR} patients and the control subjects, suggesting a decrease in their vagal activity and a compensation of sympathetic activity. Both, the cardio -systolic and -diastolic coupling strength was stronger in high-risk patients when comparing with low-risk patients. The cardio-systolic coupling analysis revealed that the systolic influence on heart rate gets weaker as the risk increases. The SVM IDC_{LR} vs. IDC_{HR} model achieved 98.9% accuracy with an area under the curve (AUC) of 0.96. The IDC and the CON groups obtained 93.6% and 0.94 accuracy and AUC, respectively. To simulate a circumstance in which the original status of the subject is unknown, a cascade model was built fusing the aforementioned models, and achieved 94.4% accuracy. In conclusion, this study introduced a novel method for SCD risk stratification for IDC patients based on new indices from coupling analysis and non-linear HRV and BPV. We have uncovered some of the complex interactions within the autonomic regulation in this type of patient.

## Introduction

According to the 2015 update of the heart disease and stroke statistics of the American Heart Association (Mozaffarian et al., 2015), ~325,000 cases of sudden cardiac death (SCD) occurred in the United States in that year, and it is the cause of 15–20% of mortality worldwide (Saour et al., 2017). The implantable cardioverter defibrillator (ICD) is commonly recommended in patients that are at high risk of suffering SCD, and the risk of SCD is halved when one is implanted, although the presence or absence of an ICD implant has no significant influence over the rate of death itself (Kober et al., 2016).

The implantation of an ICD is suggested in patients with an ejection fraction (EF) lower or equal to 35%. Currently, there is no effective means to stratify SCD risk in patients with EF above the risk threshold, who constitutes at least 70% of the patients who will suffer SCD (Chugh, 2017). Additionally, the effectiveness of ICD therapy is time dependent, making a prediction of the duration of the treatment desirable for the purpose of optimizing costs. Therefore, there is still a need for additional predictors to identify patients with idiopathic dilated cardiomyopathy (IDC) who have an increased risk of SCD and who benefit from ICD implantation (Bristow et al., 2004; Duray et al., 2006).

Previous studies have hypothesized various strategies to assess SCD risk based on clinical tests and data, imaging techniques, and signal processing methods, among others. The corrected QT interval was tested in elderly subjects and was associated with SCD risk (Panikkath et al., 2011); others analyzed T-wave inversions, wide QRS-T angle, and the left bundle branch block and found it prolonged the ability of QRS to predict all-cause mortality, including SCD (Aro et al., 2011, 2012a,b), but despite being useful for prediction, it is unable to predict individual risk. Sympathetic dominance of the autonomic nervous system in conjunction with pro-arrhythmic processes increases the probability of SCD in patients with ventricular fibrillation problems (Schwartz et al., 1992).

Studies related to cardiac imaging have analyzed the myocardial scar, including the peri-infarction border zone, to stratify SCD risk. The quantification of the total myocardial scar was explored, proving to be superior to the left ventricular ejection fraction (LVEF) in predicting an appropriate ICD therapy in cardiomyopathy patients (Bertini et al., 2012). Other studies demonstrated that the image-based analysis of myocardial scars can contribute to the decision to implant ICDs in SCD patients (Schmidt et al., 2007; Roes et al., 2009; Wu et al., 2012).

Heart rate variability (HRV) has been studied as a measurement of autonomic tone. Higher vagal tone activity is related to increased spontaneous variations in heart rate, and multiple non-linear techniques have been applied to study it (Wolf et al., 1978; Task-Force-of-the-European-Society-of-Cardiology-the-North-American-Society-of-Pacing-Electrophysiology, 1996; Ebrahimzadeh et al., 2014; Wu et al., 2014; Fujita et al., 2016). Lower levels of indices related to this variability have been associated with patients at SCD risk, regardless of their LVEF (La Rovere et al., 1998, 2003). Another index explored is heart rate turbulence, a consistent phenomenon in low-risk ischemic heart disease patients, as a measure of autonomic function. It is capable of predicting SCD-related mortality by assessing the absence of this behavior (Schmidt et al., 1999). The dynamics of the cardiovascular system behave in a highly complex way through the interplays of different linear and non-linear subsystems (Voss et al., 2009). Changes in blood pressure are reflected in changes in heart-rate regulation, and vice versa (Cohen and Taylor, 2002).

Several linear and non-linear time series analysis approaches have been developed for the quantitative analysis of the cardiovascular system in bivariate ways. However, linear approaches might be insufficient to quantify non-linear structures and the complexity of physiological systems. Therefore, approaches from non-linear time-series analysis seem to be more suited to capture complex interactions between time series, and are able to quantify direct interrelationships such as the nonlinear influence of blood pressure on heart rate. These coupling approaches are used to quantify direct and indirect relationships, as well as causal and non-causal relationships between time series, providing deeper insights into alterations of the cardiovascular system and leading to improved knowledge of the interacting regulatory mechanisms under different physiological and pathophysiological conditions. These approaches represent promising tools for detecting multivariate information flows (Schulz et al., 2013a). Some studies are based on the analysis of these interactions through the application of bivariate coupling methods.

For example, the directional cardiovascular interactions on young healthy subjects were assessed by bivariate and multivariate coupling measures, finding that bivariate measures better quantify the information transferred between indices, while trivariate better reflects the existence and delay of directed interactions (Javorka et al., 2017). Information decomposition measurements, in terms of variance or entropy, were explored to assess information dynamics in cardiovascular networks (Faes et al., 2017), analyzing the heart period, the systolic blood pressure, and the respiratory activity. The authors concluded that these measures of information transfer and information modification are better assessed through entropy-based and variance-based methods. Other work explored the polysomnographic recordings and the finger blood pressure measurements in healthy subjects in order to investigate the differences between the wake-sleep states in the heart period and the systolic blood pressure coupling (Silvani et al., 2008). They found that at low frequencies there are differences between these states in human subjects. Additionally, the complexity and causality of the interactions of cardiovascular variability series was assessed through linear model-based and non-linear model free techniques (Porta et al., 2014), and deducing that model free methods provide additional insights compared to the simpler linear model based approaches.

The behavior of cardiovascular coupling differs based on physiological conditions. Consequently, we hypothesize that the relationships between the cardiac and vascular systems will be different between the IDC patients at a high and low risk of SCD. Therefore, the aim of this study is to analyze the suitability of cardiovascular couplings for risk stratification in these patients. We propose the characterization of these interactions through features extracted from ECG and blood pressure signals, to better describe the complex dynamics of the cardiovascular interaction and identify new indices of cardiovascular risk.

## Database

The signals analyzed in this study are part of German Autonomic Regulation Trial (ART) study, oriented to evaluate risk predictors of SCD and to improve the risk stratification in IDC patients. Signals from 220 IDC patients were recorded in two hospitals: the Friedrich-Schiller University Hospital in Jena (44 patients) and the Franz-Volhard Clinic in Berlin (176 patients). All participants provided their written informed consent to a protocol approved by the local ethics committee of the two hospitals. This study complies with the Declaration of Helsinki. The study was approved by the local Ethics Committee (No. 0986-11/02).

In the acquisition protocol, the ECG signals (32 bit resolution, 1,600 Hz sampling frequency) and blood pressure (32 bit, 500 Hz) were synchronously recorded for 30 min. The ECG recordings were performed with a Porti system (TMSi BV, Netherlands), and the blood pressure recordings with the Portapress NIBP monitor (TNO Biomedical Instrumentation, Amsterdam, the Netherlands). All patients were recorded under resting conditions in supine position.

The inclusion criteria for the IDC patients were LVEF <45% and/or fractional shortening <25%, and the left ventricle end-diastolic diameter (LVEDD) >117%. Additionally, the New York heart association index (NYHA) of each patient was included. The exclusion criteria were systemic hypertension at rest, coronary artery disease, congenital heart disease, pericardial diseases, valvular heart diseases, systemic disease known to cause dilated cardiomyopathy, chronic alcoholism, sustained ventricular tachycardia, atrial fibrillation, diabetes mellitus, renal failure, and active permanent cardiac pacemaker as well as patients without sinus rhythm (Voss et al., 2010, 2012b).

Afterwards, an additional exclusion criteria was applied to discard patients with comorbidities and confounding factors influencing the autonomic regulation system. 129 of the 220 patients were excluded due to the following reasons: 36 with paced rhythm, 34 patients suffered coronary artery disease, 32 presented a high (>5%) percentage of ectopic beats or artifacts, 19 with atrial fibrillation, four because of technical problems, two suffered from hypertrophic non-obstructive cardiomyopathy, one patient with arrhythmogenic right ventricular cardiomyopathy, and another patient who was clinically unstable due to acute decomposition (Voss et al., 2012b). Finally, a total number of 91 IDC patients (21 female, 70 male) were investigated in this work.

After a median follow-up period of 28 months (range: 17–38 months), the patients were classified into two groups according to their SCD risk. The group of patients that remained in stable physical condition were considered as low risk for cardiovascular sudden death (*IDC*_{LR}). The remaining patients, who either died of SCD or needed resuscitation because of a life-threatening tachyarrhythmia, were categorized as high risk for cardiac sudden death (*IDC*_{HR}). None of these patients died from a non-cardiac disease. Additionally, 49 healthy subjects (30 male, 19 female; aged 46 ± 14 years) were used as control group (CON). Table 1 presents the baseline clinical information of the IDC patients.

**Table 1**. Baseline clinical information of ART database IDC patients (median and interquartile range).

In order to characterize autonomous regulation of the cardiovascular systems, the following time series from ECG (RR interval) and BP signals were extracted using in-house software (programming environment Delphi 3 and MATLAB^{®} R2011b):

- Time series of heart rate consisting of successive beat-to-beat intervals [BBI, tachogram, (ms)] from the ECG signal.

- Time series consisting of the maximum successive end-systolic blood pressure amplitude values over time in relation to the previous R-peak [SBP, systogram, (mmHg)] from the BP signal.

- Time series consisting of the minimum successive end-diastolic blood pressure amplitude values over time in relation to the previous R-peak [DBP, diastogram, (mmHg)] from the BP signal.

All extracted time series (BBI, SBP, DBP) were filtered by applying an adaptive variance estimation algorithm (Wessel et al., 2000) to remove and interpolate seldom occurring ventricular premature beats and artifacts (e.g., movement, electrode noise, and extraordinary peaks) to obtain normal-to-normal beat time series (NN). To obtain synchronized time series, BBI, SBP, and DBP were resampled using a linear interpolation method (2 Hz).

## Methods

The BBI, SBP, and DBP time series were extracted using algorithms based on zero crossings and different thresholds. These data were evaluated through several non-linear characterization techniques such as high resolution joint symbolic dynamics (HRSJD), segmented Poincaré plot analysis (SPPA) and normalized short-time partial directed coherence (NSTPDC).

### High Resolution Joint Symbolic Dynamics

The joint symbolic dynamics (JSD) method (Baumert et al., 2002) is based on the analysis of dynamic processes by the means of symbols. Considering BBI, SBP, and DBP time series, ** X **is defined as a bivariate sample vector that contains two out of the three time series, for all their possible combinations (BBI-SBP, BBI-DBP, and SBP-DBP), as expressed in Equation (1),

being *N* the total number of samples.

In JSD, the increments between two successive values of the temporal series are coded as “1” and the decrements and equilibriums are coded as “0,” these increments and decrements are considered in relation to a threshold *l*. The vector ** X **can be transformed into the symbolic vector

**using the rules set out in Equation (2), were the threshold**

*S**l*is set equal to 0 (Baumert et al., 2002; Giraldo et al., 2015).

Sequences of symbols are considered words of length *k*. The words are then arranged in a vector matrix ** W**. Particularly, for

*k*= 3 (

*S*

_{n},

*S*

_{n+1},

*S*

_{n+2}) an 8 ×8 vector matrix can be derived, taking values from (000, 000)

^{T}to (111, 111)

^{T}.

In order to obtain JSD indices that are more robust against noise, fluctuations and artifacts, the comparison threshold should be something other than 0. There are several advantages for choosing a non-zero threshold. For instance, a state will be generated that will help to distinguish between small and large changes in the system's variability. It is also possible to differentiate between decrements and equilibrium because both states are no longer coded with the same symbols. Lastly, the number of word types including “0” will not be the most prevalent within the ** W **matrix (Schulz et al., 2013b, 2015b). The high-resolution joint symbolic dynamics (HRJSD) method implements a three symbol JSD after setting a threshold

*l*. The increment states are coded as “2,” the decrements are coded as “1” and the equilibrium states are coded as “0.” With this technique, the transformation from

**to**

*X***varies as is shown in Equation (3).**

*S*The new space is composed of combinations of 27 different possible types of words (from 000 to 222) and a total of 729 indices. All the word types were grouped into eight pattern families, transforming the vector matrix ** W **into a vector matrix family (

*W*_{f}

**)**. These indices were analyzed by their occurrence probabilities. Afterwards, these indices were grouped according to their family description (Table 2). These pattern families represent different interactions between the branches of the autonomic regulation system, leading to indices with statistically sufficient occurrence probabilities (Schulz et al., 2015b). Additionally, the Shannon entropy (Rundle et al., 2019) was calculated for all the proposed families to assess the complexity of the coupling.

The thresholds applied were 5 ms and 1 mmHg on the BBI and the blood pressure time series, respectively. These threshold values were successfully applied in earlier work (Schulz et al., 2013b). The threshold level using spontaneous baroreflex sensitivity, in contrast to other thresholds, is the most suitable for highlighting different specific cardiovascular coupling patterns.

### Segmented Poincaré Plot Analysis

The Poincaré plot analysis (PPA) is used to quantify self-similarity in processes by plotting the data into a higher dimensional state space. Considering a time series *X*(*n*) = *x*_{1}, *x*_{2}, *x*_{3}, …, *x*_{n}, the Poincaré plot is obtained by plotting *X*(*n*) vs. *X*(*n*+1). The typical representation is an elongated scatter of plots through the line of identity with all the points whose values are near the mean placed toward the center. Long- and short-term variability indices can be by fitting an ellipse to the shape of the plot and measuring the dispersion along the minor (SD1) and major (SD2) axis of the ellipse, corresponding to their standard deviation (Rodriguez et al., 2017).

Although PPA is a non-linear characterization method, indices SD1 and SD2 can be correlated with the linear behavior of the system, hence making them a suboptimal way to explore information about the non-linear part of the process (Seeck et al., 2011; Voss et al., 2012a).

The Segmented Poincaré plot analysis (SPPA) is an enhanced pseudo-phase space quantification method that yields indices that also represent the non-linear information of the system. In SPPA the SD1 and SD2 indices are calculated similarly to PPA. Then, the scatter points are rotated α degrees around the main focus of the plot as defined in Equation (4),

The *x* and *y* axis correspond to *X*(*t*) and *X*(*t* + 1) values, $\overline{{X}_{n}}\text{}$and $\text{}\overline{{X}_{n+1}}$ are the mean values of the original and shifted ** X **time series respectively, and

*z*is the axis of the rotation (Voss et al., 2010, 2012a; Seeck et al., 2011). In this case we define α = 45 degrees in order to simplify the estimating procedure of the SD1/SD2 adapted probability. Afterwards, a 12 ×12 rectangular grid is drawn for the plot. The size of the rectangles (height, width) is adapted based on the SD1 (row) and SD2 (column) values.

Finally, for each rectangle in position (*i, j*) the single probability (ρ_{ij}) is calculated considering the number of points contained by the total number of points in the series. Afterwards, the probabilities of each row (ρ_{ri}) and column (ρ_{cj}) are calculated as the sum of their single probabilities, as shown in Equation (5):

### Normalized Short-Time Partial Directed Coherence

The directed coherence method (DC) describes how and whether two complex physiological signals are functionally connected (Faes et al., 2012). The DC method studies the relative structural relationships between the systems by breaking down their interactions into *feedback* and *feedforward* aspects.

The partial directed coherence method (PDC) determines the either direct or indirect causality between the systems analyzed. The PDC is limited to work on stationary signals and is unable to yield information about the partial correlative short-time interaction properties (Baccala and Sameshima, 2001).

The NSTPDC is able to manage non-stationary signals by evaluating their dynamic coupling changes and detecting their level and directions in multivariate and complex dynamic systems (Adochiei et al., 2013; Schulz et al., 2015a).

Normalized short-time partial directed coherence based on an *m*-dimensional multichannel auto-regressive model (MAR) process with model order *p* to determine Granger causality in the frequency domain. For the selection of the optimal model order p_{opt} of the AR(*p*) model and for the estimation of its coefficients, the stepwise least squares algorithm (Neumaier and Schneider, 2001) and the Schwarz's Bayesian Criterion (SBC) were applied (Schneider and Neumaier, 2001). NSTPDC is based on the time-variant partial directed coherence approach (tvPDC, π_{xy}(*f, n*)) providing information about the partial correlative short-time interaction properties of non-stationary signals, with *f* as the frequency and *n* the number of windows (Milde et al., 2011).

To quantify the coupling direction between two time series, ** X **and

**(e.g., BBI and SBP: with**

*Y**x*

_{BBI}and

*y*

_{SBP}) with the covariate

*z*(e.g., DBP with

*z*

_{DBP}), a coupling factor (CF) was introduced. CF was obtained by dividing the mean value π

_{xy}(

*f, n*) by the mean value of π

_{yx}(

*f, n*), defined as

These results were normalized to become a specific set of values leading to the (normalized) factor NF representing the coupling direction, given by

A normalization procedure was applied to CF leading to the normalized factor (NF). NF determinates the direction of the causal connections between the investigated time series (*x*_{BBI} and *y*_{SBP}) as a function of frequency *f*. NF takes the following values: NF = {−2, −1, 0, 1, 2}. Strong unidirectional coupling is indicated if NF is −2 or 2 (where −2 denotes *y*_{SBP} as the driver), bidirectional coupling if NF = −1 or 1 (−1 denotes *y*_{SBP} as the driver), and an equal influence in both directions and/or no coupling if NF = 0 in respect to coupling strengths (if both area indices reveal equal values that are larger than zero an equal influence in both directions is present, on the other hand, if both area indices reveal equal values but are zero there is no coupling).

For determining the coupling strength between two time series, e.g., *x*_{BBI} and *y*_{SBP} with covariate (*z*_{DBP}), the areas (A_{BBI → SBP(DBP)}, A_{SBP → BBI(DBP)}, [a.u.]) generated in space by CF were estimated in each window within the frequency band *f* = 0–2 Hz, and afterwards, averaged. A_{BBI → SBP(DBP)} and A_{SBP → BBI(DBP)} ranges between 0 to 1 [0, 1]. Hereby, 1 indicates that all causal influence originating from time series ** X **are directed toward (arrows: → ) time series

**.**

*Y*In order to take advantage of the aspect of stationarity and scale-invariance for NSTPDC analyses, a normalization (zero mean and unit variance) of the time series BBI, SBP, and DBP was performed (Schulz et al., 2015a). Therefore, each sample *i* of the BBI-, SBP-, or DBP- time series **X** = {*x*_{i}, *i* = 1, …*N*} and **Y** = {*y*_{i}, *i* = 1, …*N*} with *N* as the maximal number of samples *i* (temporal index) was first normalized by subtracting the mean of $\overline{x}$ and, then divided by the standard deviation (std) of ** X **or

**, respectively. The normalized time series x**

*Y*_{norm}and y

_{norm}(zero mean and unit variance) were thus obtained as seen in (8)

### Dual Sequence Method

A widely spread method to investigate the spontaneous baroreflex sensitivity (BRS) is the sequence method (Bertinieri et al., 1988). The BRS was obtained scanning the SBP and BBI time series for sequences of three or more successive heart beats in which a progressive increase (or decrease) in SBP was followed (with a one-beat delay). The dual sequence method (DSM) was developed to improve the analysis of the baroreflex sensitivity (Malberg et al., 1999, 2002). A fluctuation is defined when there are variations greater than 1 mmHg (increasing or decreasing) in SBP and greater than 5 ms in BBI values. The highest slope of every sequence was taken for the linear regression. The slopes of the regression lines between the SBP and BBI sequences were taken as an index for local BRS [ms/mmHg] and calculated in each recording. The DSM is based on standard sequence methods, the enhancement lies in the analysis of two different kinds of BBI response: bradycardic (an increase in SBP causes an increase in BBI) and tachycardic fluctuations (a decrease in SBP causes a decrease in BBI), whereas only the bradycardic fluctuations represented the classical spontaneous baroreflex sensitivity (*bslope*). The analysis of the tachycardic fluctuations (*tslope*) provides additional information about autonomous cardiovascular regulation (Parati et al., 2000).

### Heart Rate and Blood Pressure Variability Standard Indices

The HRV and the blood pressure variability (BPV) were quantified using the standard indices from time and frequency domain. The following indices were calculated from the BBI, SBP, and DBP time series, with NN being the normal-to-normal intervals:

- The mean of BBI, SBP, and DBP time series (BBI_meanNN, SBP_meanNN and DBP_meanNN, respectively).

- The standard deviation of BBI, SBP, and DBP time series (BBI_sdNN, SBP_sdNN, DBP_sdNN).

- The proportion derived by dividing NN50 (number of pairs of adjacent NN intervals differing by more than 50 ms in the entire recording) by the total number of NN intervals (BBI_PNN50, SBP_PNN50, DBP_PNN50).

- The square root of the mean squared differences of BBI, SBP, and DBP time series (BBI_rmssd, SBP_rmssd, DBP_rmssd).

- The power of the low frequency components (0.04–0.15 Hz) of BBI, SBP, and DBP time series (BBI_LF, SBP_LF, DBP_LF).

- The power of the high frequency components (0.15–0.4 Hz) of BBI, SBP, and DBP time series (BBI_HF, SBP_HF, DBP_HF).

- The ratio between the low and high frequency power components of BBI, SBP, and DBP time series (BBI_LF/HF, SBP_LF/HF, DBP_LF/HF).

### Feature Extraction

A total number of 621 indices were extracted to analyze the interaction between the coupling across all signals (ECG and blood pressure signals). The “cd” label represents the cardiac and diastolic coupling, “cs” label indicates the cardiac and systolic coupling, and “ds” label the diastolic and systolic coupling.

The distribution of these indices is as follows: 264 indices from the HRJSD, 216 from the SPPA, 97 from the JSD, 12 from the PPA, 9 from the NSTPDC, 21 from the standard HRV and BPV indices, and 2 from the DSMs. Summary descriptions of the indices are shown in Table 3.

**Table 3**. Coupling indices extracted from high-resolution joint symbolic dynamics, segmented Poincaré plot analysis, and normalized short-time partial directed coherence.

## Statistical Analysis

In order to reduce dataset dimensionality, Mann Whitney non-parametric statistical test was used to determine the statistical significance of the indices obtained in the characterization process. The results were analyzed for different levels of significance, including the Bonferroni criterion, considering:

Additionally, a correlation analysis was performed on the statistically significant indices. The ones with high correlation (ρ ≥ 0.7) and relative lower significance were discarded.

The leave-one-out cross-validation procedure was used to validate the results. The classification results are presented in terms of accuracy (Acc), sensitivity (Sn), specificity (Sp), and area under the curve (AUC).

## Classification Techniques

The objective of the support vector machines (SVM) method is to find a higher dimensional space where a classification problem can be more easily solved than in the original space. The vectors that defines the hyperplane are called support vectors. This technique allows separating groups of data that are not originally possible using linear classifiers (Giraldo et al., 2015).

Being ** X** = {

*x*

_{1}, …,

*x*

_{L}},

*x*∈ ℝ a given set of data vectors and

**= {**

*Y**y*

_{1}, …,

*y*

_{L}} their corresponding labels, the SVM function, defined as a linear discriminant function, known as hyperplane, is given can be defined by Equation (9),

where *w* is the normal vector to the hyperplane. The function *K*(*x*_{i}*y*_{i}) is the Kernel function that will shape the hyperplane and α_{i} and *b* define the efficiency of the classifier on the optimal values. In this study we evaluated the Gaussian, Laplace, and ANOVA kernels.

The Gaussian kernel is useful to model radially distributed data, and is defined in Equation (10), where σ is a penalization term,

Laplace kernel is a less σ influenced version of the Gaussian kernel, given by Equation (11)

The ANOVA kernel works well on multidimensional support vector regression models (Stitson et al., 1999) and is defined by Equation (12), where σ and *d* are the optimization indices of its function,

The classification problem is solved by maximizing the margin while minimizing the training error. Using the Lagrange multipliers method, a dual formulation can be obtained (Cortes and Vapnik, 1995) (Equation 13),

where *C* is a penalty parameter. Besides the scale of *C* having no direct meaning, as its value increases, the penalty assigned to errors is stronger, narrowing the decision boundary (Ben-Hur et al., 2008).

Each feature was scaled and normalized (zero mean and unit variance) in order to avoid scaling biases. For each iteration of features, the model was built by optimizing the value of *C* for each of the kernels considered, by iterating different values of σ and *d*. The indices that showed statistical differences and low correlation were used in pairs to build several SVM models. The accuracy of each model was then calculated and the one with the higher value was chosen as optimal for each type of kernel.

## Results

The calculated indices were used to analyze the cardiovascular coupling in 91 IDC patients and 49 healthy subjects. Four different comparisons were performed:

• The high risk IDC patients (*IDC*_{HR}) vs. the low risk IDC patients(*IDC*_{LR}).

• The *IDC* patients vs. the CON subjects.

• The high risk IDC patients (*IDC*_{HR}) vs. the CON subjects.

• The low risk IDC patients (*IDC*_{LR}) vs. the CON subjects.

**IDC**_{HR} Patients vs. **IDC**_{LR} Patients Compared With CON Subjects

We obtained statistically significant differences in the symbolic dynamic analysis, for both the cardio-diastolic (BBI–DBP) and diastolic-systolic (DBP-SBP) couplings. Figures 1, 2 present an example of the three-dimensional plots of the word distribution density matrix of the couplings, using the HRJSD method from (Figures 1A,C, 2A,C) IDC_{LR} and (Figures 1B,D, 2B,D) IDC_{HR} patients, respectively.

**Figure 1**. Three-dimensional plots of the word distribution density matrix using the HRJSD method (single word probabilities, word families), from a patient of **(A,C)** *IDC*_{LR} and **(B,D)** *IDC*_{HR}, respectively, for cardio-diastolic coupling.

**Figure 2**. Three-dimensional plots of the word distribution density matrix using the HRJSD method (single word probabilities, word families), from a patient of **(A,C)** *IDC*_{LR} and **(B,D)** *IDC*_{HR}, respectively, for diastolic-systolic coupling.

Figure 3 shows an example of the Poincaré plot method applied to a patient for each analyzed group, considering the systogram from (Figure 3A) a CON subject, (Figure 3B) a IDC_{LR} patient, and (Figure 3C) a IDC_{HR} patient.

**Figure 3**. Systolic blood pressure Poincaré plot analysis results from **(A)** a CON subject, **(B)** a *IDC*_{LR} patient, and **(C)** a *IDC*_{HR} patient.

Figure 4 represents the averaged NSTPDC applied on the BBI, SBP, and DBP time series couplings, for (Figure 4A) the CON, (Figure 4B) IDC_{LR}, and (Figure 4C) IDC_{HR} groups.

**Figure 4**. Averaged NSTPDC plots for cardiovascular coupling analyses for **(A)** the CON, **(B)** the *IDC*_{LR}, and **(C)** the *IDC*_{HR} group. Arrows indicating the causal coupling direction from one time series to another time series, e.g., SYS←BBI, indicating the causal link from BBI to SYS. Coupling strength ranges from blue (0, no coupling) to red (1, maximum coupling) where BBI are beat-to-beat intervals, and SYS are successive end-systolic blood pressure amplitude values over time.

Figure 5 presents the relationships between all three analyzed groups, where the arrows represent the coupling direction and the arrow thickness indicating the coupling strength. In addition, the level of statistical significance has been represented (*p* ≤ 0.01).

**Figure 5**. Graphical representation of the entire cardiovascular coupling structure (coupling strengths and directions) among the cardiac (BBI), systolic blood pressure (SYS), and diastolic blood pressure (DIA) systems, when comparing **(A)** CON vs. IDC_{LR}, **(B)** CON vs. IDC_{HR}, and **(C)** IDC_{HR} vs. IDC_{LR}. The arrows direction indicates the causal coupling direction and the thickness the coupling strength. *Indicates that the NSTPDC index associated to the coupling represented by the arrow is statistically significant when **p* ≤ 0.01; ***p* ≤ 0.001; ****p* ≤ 0.0000167.

When comparing the *IDC*_{HR} and *IDC*_{LR} groups, 96 indices presented statistically significant differences, corresponding to: 76 indices from the HRJSD, 12 from the SPPA, and 7 from the NSTPDC. After the correlation analysis, a total of 36 statistically significant indices were chosen for the classification process. Some of the most relevant indices, expressed in a mean value and 95% confidence interval, are shown in Table 4.

**Table 4**. Mean value and 95% confidence interval of the most significant indices when comparing the IDC_{HR} and IDC_{LR} patients.

### IDC Patients vs. CON Subjects

When comparing IDC patients vs. CON subjects, 261 indices presented statistically significant differences between the two groups: 169 from the HRJSD, 76 from the SPPA, and 13 from NSTPDC. A total of 82 indices remained after discarding the highly correlated indices. Table 5 shows the most relevant of them.

**Table 5**. Mean value and 95% confidence interval of the most significant indices comparing IDC patients and CON subjects.

**IDC**_{HR} and **IDC**_{LR} Patients vs. CON Subjects

When comparing *IDC*_{HR} patients vs. CON subjects and *IDC*_{LR} patients vs. CON subjects, the differences were found in 182 and 247 indices, respectively. After the correlation analysis, 61 and 85 remaining indices were chosen for the classification step. A summary of the most relevant indices is shown in Tables 6, 7.

**Table 6**. Mean value and 95% confidence interval of the most significant indices comparing the IDC_{HR} patients and CON subjects.

**Table 7**. Mean value and 95% confidence interval of the most significant indices comparing the IDC_{LR} patients and CON subjects.

### HRV and BPV Standard Indices and DSM Results

When the HRV and BPV standard indices and the DSM were evaluated in the IDC_{HR} vs. IDC_{LR} comparison, no statistically significant indices were found. In the remaining comparisons, seven indices presented statistical significances: two from the sequence analysis and five from the HRV and BPV standard indices. A summary of these results is shown in Table 8.

**Table 8**. Mean value, 95% confidence interval and *p*-value of the HRV and BPV standard indices and the dual sequence method across all comparisons.

### Classification Results

After the SVM classification step, the PPAs_SD1/SD2 and HRJSDds_LU1-P indices were deemed the optimal choices for the Laplace kernel SVM model, achieving an accuracy of 98.9% and an AUC of 0.96 for the *IDC*_{HR} vs. *IDC*_{LR} comparison. The HRJSDds-E0d and NSTPDCcs_Area_c-s allowed us to classify *IDC* patients from the CON group with an accuracy of 93.6 % and an AUC of 0.94 using the Laplace kernel. Meanwhile, the HRJSDds_LD1-V and the NSTPDCcs_Area_c-s were able to discriminate between *IDC*_{HR} patients and the CON group with an accuracy of 96.9% and AUC of 0.95 applying the Gaussian kernel. Finally, the HRJSDds-E1d and the NSTPDCcs_Area_c-s were found to be the best indices to classify *IDC*_{LR} patients from the CON group, obtaining 89.6 % accuracy and a 0.85 AUC with the Laplace kernel. The classification plots and the results are shown in Figure 6 and Table 9, respectively.

**Figure 6**. SVM classification plots: **(A)** IDC vs. CON (Laplace kernel), **(B)** IDC_{HR} vs. IDC_{LR} (Laplace kernel), **(C)** IDC_{HR} vs. CON (Gaussian kernel), and **(D)** IDC_{LR} vs. CON (Laplace kernel).

**Table 9**. Accuracy (Acc), sensitivity (Sn), specificity (Sp), and area under the curve (AUC) obtained with the best SVM model for each comparison.

### Cascaded Risk Stratification

In order to consider the typical clinical case in which the original condition of the subject is unknown, a cascade model was developed to classify this new subject using the label SDC risk, comparing the IDC vs. CON and IDC_{HR} vs. IDC_{LR} models (Figure 7). The general idea of this model is to classify the subject in either CON, LR, or HR without any prior labeling. The first step is to decide if a subject is an IDC patient or not using the IDC vs. CON model. Afterwards, those classified as IDC patients are analyzed based on the LR vs. HR model to predict their risk level.

**Figure 7**. Cascade model structure. The subject is evaluated in the IDC vs. CON model, and resulting IDC patients' level of risk is evaluated by means of the HR vs. LR model.

The cascade model achieved 94.4% accuracy. From a One vs. All approach, the results for the CON group were 86.1% sensitivity and 92.8% specificity, from the perspective of the HR group, the model yielded a sensitivity of 93.2 and specificity of 94.5%, and for the LR group, the model achieved 98.8% sensitivity and 89.0% specificity.

## Discussion

The aim of this study was to find indices capable of stratifying SCD risk in IDC patients. In order to achieve this, we derived indices extracted from BBI (ECG) and systolic and diastolic blood pressure (BP) temporal series using linear and non-linear univariate and bivariate (coupling analysis) techniques. The indices extracted from these techniques revealed patterns that behave differently in patients at high risk of SCD. These indices (mainly from coupling analyses) were used to train several SVM models in order to classify the subjects based on different levels of SCD risk.

Our findings revealed that the cardio-diastolic coupling (NSTPDCcd_NF) is bidirectional with the diastolic activity as a driver in IDC_{LR} patients. Additionally, the coupling strength from cardiac activity over both the systolic and diastolic blood pressure (NSTPDCcs_Area_c->s, NSTPDCcs_Area_c->d) gets stronger when the patients are at high risk, suggesting that the cardiac activity is significantly dominating the blood pressure. This type of relationship has been observed before in congestive heart failure patients (Marinazzo et al., 2006).

Our results also revealed that there is a significant increase in systolic pressure activity as a response to the alternant activity of diastolic pressure in patients at risk of SCD (HRJSDds_P-LU1). Similar patterns were also observed in the cardio-diastolic coupling (HRJSDcd_LD1-P), suggesting that the deterioration of autonomic regulation is more severe in patients at high risk of SCD. Earlier studies indicated that symmetric patterns in the HRJSD could be related to baroreflex-like response patterns. This suggests that this kind of behavior is also more pronounced in patients at high risk (Baumert et al., 2005; Schulz et al., 2016).

The Poincaré plot analysis revealed that the patients from the HR group have higher short-term systolic blood pressure deviation than the patients from the LR group and the CON subjects. This pathological behavior is also reflected in the PPAs_SD1/SD2 index, which is less balanced as the illness progresses. The diastolic blood pressure behaved in a similar way: the short-term diastolic blood pressure deviation was significantly lower in HR patients and their PPAd_SD1/SD2 was less balanced as well, indicating higher BPV in patients with critical conditions. In earlier studies (Tatasciore et al., 2013; Ribeiro et al., 2017), higher short-term BPV was associated with several cardiac maladies such as left ventricular systolic dysfunction and atherosclerosis, amongst others. This BPV behavior was also present in sinoaortic denervated cats (Di Rienzo et al., 1991). Additionally, baroreflex effectiveness has been studied in paraplegic subjects (Castiglioni et al., 2007), and it their BPV, compared to control subjects, was found to be higher.

The aforementioned patterns were also present when the IDC patients were compared to the CON group. In general, the indices found in the PPA suggest that short-term BPV is higher in patients with pathological conditions (Mancia et al., 1983; Mehlum et al., 2018). The PPA indices regarding the BBI revealed that short-term deviation and the short- and long-term deviation ratio of the heart rate is higher in the CON group, suggesting that the HRV in the CON group is higher than in IDC patients. Additionally, the BBI_rmssd index suggests that the HRV is lower in IDC patients compared with the CON subjects (Shaffer and Ginsberg, 2017). It is known that reduced HRV is a predictor of an adverse prognosis in patients with cardiac disease (Task-Force-of-the-European-Society-of-Cardiology-the-North-American-Society-of-Pacing-Electrophysiology, 1996; Gang and Malik, 2003). Several studies have related low HRV with heart failure (Dekker et al., 2000; Musialik-Lydka et al., 2003; Sandercock and Brodie, 2006; Patel et al., 2017). Additionally, an increase complexity in BBI randomness and a lower fractal-like behavior have been associated with SDC (Huikuri et al., 2000).

The DSM revealed that both the tslope and bslope are significantly lower in IDC patients compared to CON subjects, and lower values of these sequences are associated with baroreflex dysfunction (Di Rienzo et al., 2001; La Rovere et al., 2008). In addition, a trend emerged in these indices when the patients were compared by their level of risk: IDC_{HR} patients showed lower tslopes and bslopes.

The segmented approach of the Poincaré Plot analysis revealed that some patterns in the cardiovascular coupling are more common in HR patients. There was a significantly higher concentration of these patterns in column 5 and 8 in both cardio-diastolic and cardio-systolic couplings in patients at higher risk of SDC, indicating a lower variability in their baroreflex activity compared with patients at low risk. SPPA patterns (SPPAcd_Column_5, SPPAcd_Column_8) occurred more frequently and were more concentrated in low risk patients, suggesting that these patients present a higher HRV compared to patients at higher risk of SCD. The viability of HRV as a reliable predictor of SCD in IDC patients was questioned in an earlier study (Grimm et al., 2003), which did not support the hypothesis that HRV is a reliable predictor of SCD in IDC patients. However, their results were based on a time domain analysis of HRV only, whereas the significant indices analyzed in our work come primarily from non-linear methods. These characterization methods are more suitable for describing the non-linear behavior of HRV in pathological conditions.

The HRJSD results suggest that patients at high risk adapt less frequently to changes in blood pressure, reflected in the lower presence of decreasing patterns of BBI (E0) in response to decreasing patterns in diastolic blood pressure (LD1). This may be a result of the vagal response, causing less frequent parasympathetic activity, leading to a less effective control of the blood pressure, and consequently the heart rate, in patients at higher risk conditions.

These results were consistent when CON subjects were compared with IDC patients: decreasing patterns (E0) in diastolic BPV were reflected in decreasing heart rates at higher frequencies in CON subjects. The DIA_LF/HF was higher in the patients when compared to the CON group. Higher levels of this index reflects efferent sympathetic activity (Mccraty and Shaffer, 2015). Additionally, a higher prevalence of unchanging patterns (E1) was found in the IDC patients compared with the CON group, indicating that changes in blood pressure are frequently not reflected in changes in heart rate in patients with pathological conditions.

In addition, the HRJSDcd-E2d and HRJSDcd-LD1d indices showed that steady (E2) and low decreasing (LD1) diastolic blood pressure patterns, independent from all BBI patterns, are more frequent in CON subjects. These indices suggest a worsened circulatory homeostasis in IDC patients and support the idea of the influence of baroreflex activity in pathological conditions (Kishi, 2018).

Earlier studies (Gavish et al., 2008; Schillaci and Pucci, 2010) have stated that the relationship between systolic and diastolic blood pressure should be coherent: if one of them increases, the other is expected to increase as well. The results of the systolic-diastolic coupling revealed that patterns that are opposing in nature (sLU1-dP, sLD1-dV) are more frequent in the HR group. This suggests that the relationship between systolic and diastolic blood pressure loses linearity as the pathological condition worsens.

The coupling strength of the cardiovascular and diastolic-systolic couplings are stronger in pathological conditions. In addition, the symmetric patterns of the diastolic-systolic coupling activity were less frequent in patients. This may be caused by the effect of the autonomic regulation mechanisms in pathological conditions (Floras, 2009; Kishi, 2012).

To summarize, our results suggest that there is a gradual loss of HRV as SCD risk increases and, at the same time, BPV increases alongside SCD risk. There is great controversy about the prognostic value of linear time and frequency domain HRV indices for risk stratification among this type of patient (Grimm et al., 2003; Hohnloser et al., 2003; Minamihaba et al., 2003; Valencia et al., 2009; Voss et al., 2010, 2012a). The results of this work support the idea that the commonly used techniques for analysis of the time and frequency domain of HRV is not suitable for risk stratification. However, the combination of non-linear HRV analysis and linear as well as non-linear coupling analysis seems to be a promising tool for risk assessment in IDC patients (Voss et al., 2012b; Valencia et al., 2013; Fischer and Voss, 2014). The processes involved in the circulatory homeostasis are by nature non-linear. Therefore, the differences between the stages of this process can be more adequately revealed by the quantification of the signal properties rather than by the assessment of their magnitude.

We hypothesize that a dysfunction of the vagal activity, and in general of the baroreflex mechanism as a whole, prevents the body from maintaining circulatory homeostasis correctly. This reduction in vagal activity and increase in the sympathetic influence exposes the cardiovascular system to frequent states of stress that contribute to the worsening of the condition over time. The gradual worsening of the heart rate and blood pressure variability's in the different SDC risk stages considered here supports this assumption. This kind of impairment has also been associated with other cardiac pathological conditions like ventricular fibrillation during myocardial ischemia (La Rovere et al., 2001). Abnormal sympathetic neural firing has been associated with SCD and the genesis of ventricular tachyarrhythmias (Esler et al., 1995). A similar behavior has been observed in elder mice (Freeling and Li, 2015), who presented a reduced baroreflex bradycardic response compared to younger mice.

This study has some limitations that are important to consider. The average age in the CON group was lower than in the IDC group, the influence of age in the study of HRV has been widely studied in earlier work (Voss et al., 2015). However, this limitation does not affect our results for HR vs. LR comparisons. A higher number of indices were analyzed than the number of patients in the database. Consequently, in order to minimize problems due to possible overfitting, we stablished different levels of statistical significance to reduce the dataset dimensionality, including the Bonferroni-Holm correction criterion.

The characterization of linear and non-linear coupling could be also analyzed based on linear and non-linear Granger causality in time- and frequency domains (Schulz et al., 2013a). The objective of this study was to evaluate the general behavior of the underlying coupling, through the average of all features (windows) over time was performed. Nevertheless, it is uncertain if another time-invariant time domain measure based on Granger causality would be more appropriate and would have more discriminative power (Porta and Faes, 2016). It would be interesting for an ongoing study.

Another limitation of this work is related to comorbidities and confounding factors influencing the autonomic regulation system. Therefore, these exclusion criterions make risk stratification not yet applicable to every patient.

## Conclusion

With this study we suggest that indices from coupling analysis and non-linear HRV and BPV can contribute to the development of risk stratification in IDC patients. We have introduced a novel cascade model that successfully classified subjects into different levels of SCD risk (CON, IDC_{LR}, and IDC_{HR}). Further, this study allowed us to uncover, for the first time, some of the complex interactions that take place within autonomic regulation, leading to a more accurate modeling, and interpretation of these processes in pathological conditions. Our findings revealed that there is a gradual decrease of HRV and an increase of BPV as the SCD risk increases. We conclude that patients at high risk of SCD can no longer maintain circulatory homeostasis consistently, leading to states of stress that worsens the condition over time.

However, these results should be validated with a greater number of patients, especially in the high-risk group. Therefore, the results presented in this work are more of a hypothesis-generating nature than confirmatory.

## Ethics Statement

This study are part of German Autonomic Regulation Trial (ART) study, oriented to evaluate risk predictors of SCD and to improve the risk stratification in IDC patients. Signals from 220 IDC patients were recorded in two hospitals: the Friedrich-Schiller University Hospital in Jena (44 patients) and the Franz-Volhard Clinic in Berlin (176 patients). Each study was performed following the recommendations of the Declaration of Helsinki and patient consent was granted. The study was approved by the local Ethics Committee (No. 0986-11/02).

## Author Contributions

BG and AV: conceptualization. JR, SS, BG, and AV: methodology, validation, writing—original draft preparation, writing—review and editing, and final version approval. BG, JR, and SS: software and data analysis. AV: data acquisition.

## Funding

This work was partly supported by grants from the German Federal Ministry for Economic Affairs and Energy (KF 2447309KJ4, ZF4485201SB7), by the Secretariat of Universities and Research of the Department of Economy and Knowledge of the Government of Catalonia (Consolidated research group GRC 2017 SGR 1770), by the CERCA Program/Generalitat de Catalunya, and by the Spanish Ministry of Economy and Competitiveness (DPI2015-68820-R MINECO/FEDER), and the Spanish Ministry of Science, Innovation and Universities (RTI2018-098472-B-I00).

## Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## References

Adochiei, F., Schulz, S., Edu, I., Costin, H., and Voss, A. (2013). A new normalised short time PDC for dynamic coupling analyses. *Biomed. Tech.* 58 (Suppl. 1). doi: 10.1515/bmt-2013-4167

Aro, A. L., Anttonen, O., Tikkanen, J. T., Junttila, M. J., Kerola, T., Rissanen, H. A., et al. (2011). Intraventricular conduction delay in a standard 12-lead electrocardiogram as a predictor of mortality in the general population. *Circ. Arrhythm. Electrophysiol.* 4, 704–710. doi: 10.1161/CIRCEP.111.963561

Aro, A. L., Anttonen, O., Tikkanen, J. T., Junttila, M. J., Kerola, T., Rissanen, H. A., et al. (2012a). Prevalence and prognostic significance of T-wave inversions in right precordial leads of a 12-lead electrocardiogram in the middle-aged subjects. *Circulation* 125, 2572–2577. doi: 10.1161/CIRCULATIONAHA.112.098681

Aro, A. L., Huikuri, H. V., Tikkanen, J. T., Junttila, M. J., Rissanen, H. A., Reunanen, A., et al. (2012b). QRS-T angle as a predictor of sudden cardiac death in a middle-aged general population. *Europace* 14, 872–876. doi: 10.1093/europace/eur393

Baccala, L. A., and Sameshima, K. (2001). Partial directed coherence: a new concept in neural structure determination. *Biol. Cybern.* 84, 463–474. doi: 10.1007/PL00007990

Baumert, M., Baier, V., Truebner, S., Schirdewan, A., and Voss, A. (2005). Short- and long-term joint symbolic dynamics of heart rate and blood pressure in dilated cardiomyopathy. *IEEE Trans. Biomed. Eng.* 52, 2112–2115. doi: 10.1109/TBME.2005.857636

Baumert, M., Walther, T., Hopfe, J., Stepan, H., Faber, R., and Voss, A. (2002). Joint symbolic dynamic analysis of beat-to-beat interactions of heart rate and systolic blood pressure in normal pregnancy. *Med. Biol. Eng. Comput.* 40, 241–245. doi: 10.1007/BF02348131

Ben-Hur, A., Ong, C. S., Sonnenburg, S., Scholkopf, B., and Ratsch, G. (2008). Support vector machines and kernels for computational biology. *PLoS Comput Biol.* 4:e1000173. doi: 10.1371/journal.pcbi.1000173

Bertini, M., Schalij, M. J., Bax, J. J., and Delgado, V. (2012). Emerging role of multimodality imaging to evaluate patients at risk for sudden cardiac death. *Circ. Cardiovasc. Imaging* 5, 525–535. doi: 10.1161/CIRCIMAGING.110.961532

Bertinieri, G., Di Rienzo, M., Cavallazzi, A., Ferrari, A. U., Pedotti, A., and Mancia, G. (1988). Evaluation of baroreceptor reflex by blood pressure monitoring in unanesthetized cats. *Am. J. Physiol.* 254, H377–H383. doi: 10.1152/ajpheart.1988.254.2.H377

Bristow, M. R., Saxon, L. A., Boehmer, J., Krueger, S., Kass, D. A., De Marco, T., et al. (2004). Cardiac-resynchronization therapy with or without an implantable defibrillator in advanced chronic heart failure. *N. Engl. J. Med.* 350, 2140–2150. doi: 10.1056/NEJMoa032423

Castiglioni, P., Di Rienzo, M., Veicsteinas, A., Parati, G., and Merati, G. (2007). Mechanisms of blood pressure and heart rate variability: an insight from low-level paraplegia. *Am. J. Physiol. Regul. Integr. Comp. Physiol.* 292, R1502–1509. doi: 10.1152/ajpregu.00273.2006

Chugh, S. S. (2017). Sudden cardiac death in 2017: spotlight on prediction and prevention. *Int. J. Cardiol.* 237, 2–5. doi: 10.1016/j.ijcard.2017.03.086

Cohen, M. A., and Taylor, J. A. (2002). Short-term cardiovascular oscillations in man: measuring and modelling the physiologies. *J. Physiol.* 542, 669–683. doi: 10.1113/jphysiol.2002.017483

Cortes, C., and Vapnik, V. (1995). Support vector networks. *Mach. Learn.* 20, 273–297. doi: 10.1007/BF00994018

Dekker, J. M., Crow, R. S., Folsom, A. R., Hannan, P. J., Liao, D., Swenne, C. A., et al. (2000). Low heart rate variability in a 2-minute rhythm strip predicts risk of coronary heart disease and mortality from several causes: the ARIC Study. Atherosclerosis risk in communities. *Circulation* 102, 1239–1244. doi: 10.1161/01.CIR.102.11.1239

Di Rienzo, M., Parati, G., Castiglioni, P., Omboni, S., Ferrari, A. U., Ramirez, A. J., et al. (1991). Role of sinoaortic afferents in modulating BP and pulse-interval spectral characteristics in unanesthetized cats. *Am. J. Physiol.* 261, H1811–H1818. doi: 10.1152/ajpheart.1991.261.6.H1811

Di Rienzo, M., Parati, G., Castiglioni, P., Tordi, R., Mancia, G., and Pedotti, A. (2001). Baroreflex effectiveness index: an additional measure of baroreflex control of heart rate in daily life. *Am. J. Physiol. Regul. Integr. Comp. Physiol.* 280, R744–R751. doi: 10.1152/ajpregu.2001.280.3.R744

Duray, G., Israel, C. W., and Hohnloser, S. H. (2006). Recent primary prevention implantable cardioverter defibrillator trials. *Curr. Opin. Cardiol.* 21, 15–19. doi: 10.1097/01.hco.0000198978.55637.b5

Ebrahimzadeh, E., Pooyan, M., and Bijar, A. (2014). A novel approach to predict sudden cardiac death (SCD) using nonlinear and time-frequency analyses from HRV signals. *PLoS ONE* 9:e81896. doi: 10.1371/journal.pone.0081896

Esler, M. D., Thompson, J. M., Kaye, D. M., Turner, A. G., Jennings, G. L., Cox, H. S., et al. (1995). Effects of aging on the responsiveness of the human cardiac sympathetic nerves to stressors. *Circulation* 91, 351–358. doi: 10.1161/01.CIR.91.2.351

Faes, L., Erla, S., and Nollo, G. (2012). Measuring connectivity in linear multivariate processes: definitions, interpretation, and practical analysis. *Comput. Math. Methods Med.* 2012:140513. doi: 10.1155/2012/140513

Faes, L., Porta, A., Nollo, G., and Javorka, M. (2017). Information decomposition in multivariate systems: definitions, implementation and application to cardiovascular *Networks* 19, 1–28. doi: 10.3390/e19010005

Fischer, C., and Voss, A. (2014). Three-dimensional segmented poincaré plot analyses SPPA3 investigates cardiovascular and cardiorespiratory couplings in hypertensive pregnancy disorders. *Front. Bioeng. Biotechnol.* 2:51. doi: 10.3389/fbioe.2014.00051

Floras, J. S. (2009). Sympathetic nervous system activation in human heart failure: clinical implications of an updated model. *J. Am. Coll. Cardiol.* 54, 375–385. doi: 10.1016/j.jacc.2009.03.061

Freeling, J. L., and Li, Y. (2015). Age-related attenuation of parasympathetic control of the heart in mice. *Int. J. Physiol. Pathophysiol. Pharmacol.* 7, 126–135.

Fujita, H., Acharya, R., Sudarshan, V., Ghista, D., Sree, V., Lim Wei Jie, L., et al. (2016). Sudden cardiac death (SDC) prediction based on nonlinear heart rate variability features and SCD index. *Appl. Soft Comput.* 43, 210–519. doi: 10.1016/j.asoc.2016.02.049

Gang, Y., and Malik, M. (2003). Heart rate variability analysis in general medicine. *Indian Pacing Electrophysiol. J.* 3, 34–40.

Gavish, B., Ben-Dov, I. Z., and Bursztyn, M. (2008). Linear relationship between systolic and diastolic blood pressure monitored over 24 h: assessment and correlates. *J. Hypertens.* 26, 199–209. doi: 10.1097/HJH.0b013e3282f25b5a

Giraldo, B. F., Rodriguez, J., Caminal, P., Bayes-Genis, A., and Voss, A. (2015). Cardiorespiratory and cardiovascular interactions in cardiomyopathy patients using joint symbolic dynamic analysis. *Conf. Proc. IEEE Eng. Med. Biol. Soc.* 2015, 306–309. doi: 10.1109/EMBC.2015.7318361

Grimm, W., Herzum, I., Muller, H. H., and Christ, M. (2003). Value of heart rate variability to predict ventricular arrhythmias in recipients of prophylactic defibrillators with idiopathic dilated cardiomyopathy. *Pacing Clin. Electrophysiol.* 26, 411–415. doi: 10.1046/j.1460-9592.2003.00060.x

Hohnloser, S. H., Klingenheben, T., Bloomfield, D., Dabbous, O., and Cohen, R. J. (2003). Usefulness of microvolt T-wave alternans for prediction of ventricular tachyarrhythmic events in patients with dilated cardiomyopathy: results from a prospective observational study. *J. Am. Coll. Cardiol.* 41, 2220–2224. doi: 10.1016/S0735-1097(03)00467-4

Huikuri, H. V., Makikallio, T. H., Peng, C. K., Goldberger, A. L., Hintze, U., and Moller, M. (2000). Fractal correlation properties of R-R interval dynamics and mortality in patients with depressed left ventricular function after an acute myocardial infarction. *Circulation* 101, 47–53. doi: 10.1161/01.CIR.101.1.47

Javorka, M., Krohova, J., Czippelova, B., Turianikova, Z., Lazarova, Z., Javorka, K., et al. (2017). Basic cardiovascular variability signals: mutual directed interactions explored in the information domain. *Physiol. Meas.* 38, 877–894. doi: 10.1088/1361-6579/aa5b77

Kishi, T. (2012). Heart failure as an autonomic nervous system dysfunction. *J. Cardiol.* 59, 117–122. doi: 10.1016/j.jjcc.2011.12.006

Kishi, T. (2018). Disruption of central antioxidant property of nuclear factor erythroid 2-related factor 2 worsens circulatory homeostasis with baroreflex dysfunction in heart failure. *Int. J. Mol. Sci.* 19:646. doi: 10.3390/ijms19030646

Kober, L., Thune, J. J., Nielsen, J. C., Haarbo, J., Videbaek, L., Korup, E., et al. (2016). Defibrillator implantation in patients with nonischemic systolic heart failure. *N. Engl. J. Med.* 375, 1221–1230. doi: 10.1056/NEJMoa1608029

La Rovere, M. T., Bigger, J. T. Jr., Marcus, F. I., Mortara, A., and Schwartz, P. J. (1998). Baroreflex sensitivity and heart-rate variability in prediction of total cardiac mortality after myocardial infarction. ATRAMI (Autonomic Tone and Reflexes After Myocardial Infarction) Investigators. *Lancet* 351, 478–484. doi: 10.1016/S0140-6736(97)11144-8

La Rovere, M. T., Pinna, G. D., Hohnloser, S. H., Marcus, F. I., Mortara, A., Nohara, R., et al. (2001). Baroreflex sensitivity and heart rate variability in the identification of patients at risk for life-threatening arrhythmias: implications for clinical trials. *Circulation* 103, 2072–2077. doi: 10.1161/01.CIR.103.16.2072

La Rovere, M. T., Pinna, G. D., Maestri, R., Mortara, A., Capomolla, S., Febo, O., et al. (2003). Short-term heart rate variability strongly predicts sudden cardiac death in chronic heart failure patients. *Circulation* 107, 565–570. doi: 10.1161/01.CIR.0000047275.25795.17

La Rovere, M. T., Pinna, G. D., and Raczak, G. (2008). Baroreflex sensitivity: measurement and clinical implications. *Ann. Noninvasive Electrocardiol.* 13, 191–207. doi: 10.1111/j.1542-474X.2008.00219.x

Malberg, H., Wessel, N., Hasart, A., Osterziel, K. J., and Voss, A. (2002). Advanced analysis of spontaneous baroreflex sensitivity, blood pressure and heart rate variability in patients with dilated cardiomyopathy. *Clin. Sci.* 102, 465–473. doi: 10.1042/cs1020465

Malberg, H., Wessel, N., Schirdewan, A., Osterziel, K. J., and Voss, A. (1999). [Dual sequence method for analysis of spontaneous baroreceptor reflex sensitivity in patients with dilated cardiomyopathy]. *Z. Kardiol.* 88, 331–337. doi: 10.1007/s003920050294

Mancia, G., Ferrari, A., Gregorini, L., Parati, G., Pomidossi, G., Bertinieri, G., et al. (1983). Blood pressure and heart rate variabilities in normotensive and hypertensive human beings. *Circ. Res.* 53, 96–104. doi: 10.1161/01.RES.53.1.96

Marinazzo, D., Pellicoro, M., and Stramaglia, S. (2006). Nonlinear parametric model for Granger causality of time series. *Phys. Rev. E Stat. Nonlin. Soft Matter Phys.* 73:066216. doi: 10.1103/PhysRevE.73.066216

Mccraty, R., and Shaffer, F. (2015). Heart rate variability: new perspectives on physiological mechanisms, assessment of self-regulatory capacity, and health risk. *Glob. Adv. Health Med.* 4, 46–61. doi: 10.7453/gahmj.2014.073

Mehlum, M. H., Liestol, K., Kjeldsen, S. E., Julius, S., Hua, T. A., Rothwell, P. M., et al. (2018). Blood pressure variability and risk of cardiovascular events and death in patients with hypertension and different baseline risks. *Eur. Heart J.* 39, 2243–2251. doi: 10.1093/eurheartj/ehx760

Milde, T., Schwab, K., Walther, M., Eiselt, M., Schelenz, C., Voss, A., et al. (2011). Time-variant partial directed coherence in analysis of the cardiovascular system. A methodological study. *Physiol. Meas.* 32, 1787–1805. doi: 10.1088/0967-3334/32/11/S06

Minamihaba, O., Yamaki, M., Tomoike, H., and Kubota, I. (2003). Severity in myocardial dysfunction contributed to long-term fluctuation of heart rate, rather than short-term fluctuations. *Ann. Noninvasive Electrocardiol.* 8, 132–138. doi: 10.1046/j.1542-474X.2003.08207.x

Mozaffarian, D., Benjamin, E. J., Go, A. S., Arnett, D. K., Blaha, M. J., Cushman, M., et al. (2015). Heart disease and stroke statistics−2015 update: a report from the American Heart Association. *Circulation* 131, e29–322. doi: 10.1161/CIR.0000000000000152

Musialik-Lydka, A., Sredniawa, B., and Pasyk, S. (2003). Heart rate variability in heart failure. *Kardiol. Pol.* 58, 10–16.

Neumaier, A., and Schneider, T. (2001). Estimation of parameters and eigenmodes of multivariate autoregressive models. *ACM Trans. Math. Softw.* 27, 27–57. doi: 10.1145/382043.382304

Panikkath, R., Reinier, K., Uy-Evanado, A., Teodorescu, C., Hattenhauer, J., Mariani, R., et al. (2011). Prolonged Tpeak-to-tend interval on the resting ECG is associated with increased risk of sudden cardiac death. *Circ. Arrhythm. Electrophysiol.* 4, 441–447. doi: 10.1161/CIRCEP.110.960658

Parati, G., Di Rienzo, M., and Mancia, G. (2000). How to measure baroreflex sensitivity: from the cardiovascular laboratory to daily life. *J. Hypertens.* 18, 7–19. doi: 10.1097/00004872-200018010-00003

Patel, V. N., Pierce, B. R., Bodapati, R. K., Brown, D. L., Ives, D. G., and Stein, P. K. (2017). Association of holter-derived heart rate variability parameters with the development of congestive heart failure in the cardiovascular health study. *JACC Heart Fail.* 5, 423–431. doi: 10.1016/j.jchf.2016.12.015

Porta, A., and Faes, L. (2016). Wiener-granger causality in network physiology with applications to cardiovascular control and neuroscience. *Proc. IEEE* 104, 282–309. doi: 10.1109/JPROC.2015.2476824

Porta, A., Faes, L., Bari, V., Marchi, A., Bassani, T., Nollo, G., et al. (2014). Effect of age on complexity and causality of the cardiovascular control: comparison between model-based and model-free approaches. *PLoS ONE.* 9:e89463. doi: 10.1371/journal.pone.0089463

Ribeiro, A. H., Lotufo, P. A., Fujita, A., Goulart, A. C., Chor, D., Mill, J. G., et al. (2017). Association between short-term systolic blood pressure variability and carotid intima-media thickness in ELSA-brasil baseline. *Am. J. Hypertens.* 30, 954–960. doi: 10.1093/ajh/hpx076

Rodriguez, J., Voss, A., Caminal, P., Bayes-Genis, A., and Giraldo, B. F. (2017). Characterization and classification of patients with different levels of cardiac death risk by using Poincaré plot analysis. *Conf. Proc. IEEE Eng. Med. Biol. Soc.* 2017, 1332–1335. doi: 10.1109/EMBC.2017.8037078

Roes, S. D., Borleffs, C. J., Van Der Geest, R. J., Westenberg, J. J., Marsan, N. A., Kaandorp, T. A., et al. (2009). Infarct tissue heterogeneity assessed with contrast-enhanced MRI predicts spontaneous ventricular arrhythmia in patients with ischemic cardiomyopathy and implantable cardioverter-defibrillator. *Circ. Cardiovasc. Imaging* 2, 183–190. doi: 10.1161/CIRCIMAGING.108.826529

Rundle, J. B., Giguere, A., Turcotte, D. L., Crutchfield, J. P., and Donnellan, A. (2019). Global seismic nowcasting with shannon information entropy. *Earth Space Sci.* 6, 191–197. doi: 10.1029/2018EA000464

Sandercock, G. R., and Brodie, D. A. (2006). The role of heart rate variability in prognosis for different modes of death in chronic heart failure. *Pacing Clin. Electrophysiol.* 29, 892–904. doi: 10.1111/j.1540-8159.2006.00457.x

Saour, B. J. Y, Philbin, E, Tan, H., Nguyen, D., O'brien, J., Sidhu, S., et al. (2017). Following theraoeutic hypothermia in patients with sudden cardiac death due to ventricular tachycardia or fibrillation. *Int. J. Clin. Med.* 8, 293–305. doi: 10.4236/ijcm.2017.85028

Schillaci, G., and Pucci, G. (2010). The dynamic relationship between systolic and diastolic blood pressure: yet another marker of vascular aging? *Hypertens. Res.* 33, 659–661. doi: 10.1038/hr.2010.95

Schmidt, A., Azevedo, C. F., Cheng, A., Gupta, S. N., Bluemke, D. A., Foo, T. K., et al. (2007). Infarct tissue heterogeneity by magnetic resonance imaging identifies enhanced cardiac arrhythmia susceptibility in patients with left ventricular dysfunction. *Circulation* 115, 2006–2014. doi: 10.1161/CIRCULATIONAHA.106.653568

Schmidt, G., Malik, M., Barthel, P., Schneider, R., Ulm, K., Rolnitzky, L., et al. (1999). Heart-rate turbulence after ventricular premature beats as a predictor of mortality after acute myocardial infarction. *Lancet* 353, 1390–1396. doi: 10.1016/S0140-6736(98)08428-1

Schneider, T., and Neumaier, A. (2001). Algorithm 808: ARfit—a matlab package for the estimation of parameters and eigenmodes of multivariate autoregressive models. *ACM Trans. Math. Softw.* 27, 58–65. doi: 10.1145/382043.382316

Schulz, S., Adochiei, F. C., Edu, I. R., Schroeder, R., Costin, H., Bar, K. J., et al. (2013a). Cardiovascular and cardiorespiratory coupling analyses: a review. *Philos. Trans. A Math. Phys. Eng. Sci.* 371:20120191. doi: 10.1098/rsta.2012.0191

Schulz, S., Bär, K. J., and Voss, A. (2015a). Analyses of heart rate, respiration and cardiorespiratory coupling in patients with schizophrenia. *Entropy* 17, 483–501. doi: 10.3390/e17020483

Schulz, S., Bolz, M., Bar, K. J., and Voss, A. (2016). Central- and autonomic nervous system coupling in schizophrenia. *Philos. Trans. Math. Phys. Eng. Sci.* 374:20150178. doi: 10.1098/rsta.2015.0178

Schulz, S., Haueisen, J., Bar, K. J., and Andreas, V. (2015b). High-resolution joint symbolic analysis to enhance classification of the cardiorespiratory system in patients with schizophrenia and their relatives. *Philos. Trans. Math. Phys. Eng. Sci.* 373:20140098. doi: 10.1098/rsta.2014.0098

Schulz, S., Tupaika, N., Berger, S., Haueisen, J., Bär, K. J., and Voss, A. (2013b). Cardiovascular coupling analysis with high-resolution joint symbolic dynamics in patients suffering from acute schizophrenia. *Physiol. Meas.* 34, 883–901. doi: 10.1088/0967-3334/34/8/883

Schwartz, P. J., La Rovere, M. T., and Vanoli, E. (1992). Autonomic nervous system and sudden cardiac death. Experimental basis and clinical observations for post-myocardial infarction risk stratification. *Circulation* 85, I77–91.

Seeck, A., Baumert, M., Fischer, C., Khandoker, A., Faber, R., and Voss, A. (2011). Advanced Poincaré plot analysis differentiates between hypertensive pregnancy disorders. *Physiol. Meas.* 32, 1611–1622. doi: 10.1088/0967-3334/32/10/009

Shaffer, F., and Ginsberg, J. P. (2017). An overview of heart rate variability metrics and norms. *Front. Public Health* 5:258. doi: 10.3389/fpubh.2017.00258

Silvani, A., Grimaldi, D., Vandi, S., Barletta, G., Vetrugno, R., Provini, F., et al. (2008). Sleep-dependent changes in the coupling between heart period and blood pressure in human subjects. *Am. J. Physiol. Regul. Integr. Comp. Physiol.* 294, R1686–R1692. doi: 10.1152/ajpregu.00756.2007

Stitson, M., Gammerman, A., Vapnik, V., Vovk, V., Watkins, C., and Weston, J. (1999). *Support Vector Regression With ANOVA Decomposition Kernels.* Cambridge, MA: MIT Press.

Task-Force-of-the-European-Society-of-Cardiology-the-North-American-Society-of-Pacing-Electrophysiology (1996). Heart rate variability: standards of measurement, physiological interpretation and clinical use. Task force of the european society of cardiology and the North American society of pacing and electrophysiology. *Circulation* 93, 1043–1065.

Tatasciore, A., Zimarino, M., Tommasi, R., Renda, G., Schillaci, G., Parati, G., et al. (2013). Increased short-term blood pressure variability is associated with early left ventricular systolic dysfunction in newly diagnosed untreated hypertensive patients. *J. Hypertens.* 31, 1653–1661. doi: 10.1097/HJH.0b013e328361e4a6

Valencia, J. F., Vallverdu, M., Porta, A., Voss, A., Schroeder, R., Vazquez, R., et al. (2013). Ischemic risk stratification by means of multivariate analysis of the heart rate variability. *Physiol. Meas.* 34, 325–338. doi: 10.1088/0967-3334/34/3/325

Valencia, J. F., Vallverdu, M., Schroeder, R., Voss, A., Vazquez, R., Bayes De Luna, A., et al. (2009). Complexity of the short-term heart-rate variability. *IEEE Eng. Med. Biol. Mag.* 28, 72–78. doi: 10.1109/MEMB.2009.934621

Voss, A., Fischer, C., Schroeder, R., Figulla, H. R., and Goernig, M. (2010). Segmented Poincaré plot analysis for risk stratification in patients with dilated cardiomyopathy. *Methods Inf. Med.* 49, 511–515. doi: 10.3414/ME09-02-0050

Voss, A., Fischer, C., Schroeder, R., Figulla, H. R., and Goernig, M. (2012a). Lagged segmented Poincaré plot analysis for risk stratification in patients with dilated cardiomyopathy. *Med. Biol. Eng. Comput.* 50, 727–736. doi: 10.1007/s11517-012-0925-5

Voss, A., Goernig, M., Schroeder, R., Truebner, S., Schirdewan, A., and Figulla, H. R. (2012b). Blood pressure variability as sign of autonomic imbalance in patients with idiopathic dilated cardiomyopathy. *Pacing Clin. Electrophysiol.* 35, 471–479. doi: 10.1111/j.1540-8159.2011.03312.x

Voss, A., Schroeder, R., Heitmann, A., Peters, A., and Perz, S. (2015). Short-term heart rate variability–influence of gender and age in healthy subjects. *PLoS ONE* 10:e0118308. doi: 10.1371/journal.pone.0118308

Voss, A., Schulz, S., Schroeder, R., Baumert, M., and Caminal, P. (2009). Methods derived from nonlinear dynamics for analysing heart rate variability. *Philos. Trans. Math. Phys. Eng. Sci.* 367, 277–296. doi: 10.1098/rsta.2008.0232

Wessel, N., Voss, A., Malberg, H., Ziehmann, C., Voss, H., Schirdewan, A., et al. (2000). Nonlinear analysis of complex phenomena in cardiological data. *Z. Herzschr. Elektrophys* 11, 159–173. doi: 10.1007/s003990070035

Wolf, M. M., Varigos, G. A., Hunt, D., and Sloman, J. G. (1978). Sinus arrhythmia in acute myocardial infarction. *Med. J. Aust.* 2, 52–53.

Wu, K. C., Gerstenblith, G., Guallar, E., Marine, J. E., Dalal, D., Cheng, A., et al. (2012). Combined cardiac magnetic resonance imaging and C-reactive protein levels identify a cohort at low risk for defibrillator firings and death. *Circ. Cardiovasc. Imaging* 5, 178–186. doi: 10.1161/CIRCIMAGING.111.968024

Keywords: idiopathic dilated cardiomyopathy, heart rate variability, blood pressure variability, coupling analysis, sudden cardiac death, risk stratification

Citation: Rodriguez J, Schulz S, Giraldo BF and Voss A (2019) Risk Stratification in Idiopathic Dilated Cardiomyopathy Patients Using Cardiovascular Coupling Analysis. *Front. Physiol.* 10:841. doi: 10.3389/fphys.2019.00841

Received: 07 February 2019; Accepted: 19 June 2019;

Published: 09 July 2019.

Edited by:

Zbigniew R. Struzik, The University of Tokyo, JapanReviewed by:

Luca Faes, University of Palermo, ItalyPhilip Thomas Clemson, Lancaster University, United Kingdom

Copyright © 2019 Rodriguez, Schulz, Giraldo and Voss. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Steffen Schulz, Steffen.Schulz@eah-jena.de