Characterization of T Wave Amplitude, Duration and Morphology Changes During Hemodialysis: Relationship With Serum Electrolyte Levels and Heart Rate

— Objective: Chronic kidney disease affects 7 more than 10% of the world population. Changes in serum 8 ion concentrations increase the risk for ventricular ar- 9 rhythmias and sudden cardiac death, particularly in end- 10 stage renal disease (ESRD) patients. We characterized how 11 T wave amplitude, duration and morphology descriptors 12 change with variations in serum levels of potassium and 13 calcium and in heart rate, both in ESRD patients and in 14 simulated ventricular ﬁbers. Methods: Electrocardiogram 15 (ECG) recordings from twenty ESRD patients undergoing 16 hemodialysis (HD) and pseudo-ECGs (pECGs) calculated 17 from twenty-two simulated ventricular ﬁbers at varying 18 transmural heterogeneity levels were processed to quantify 19 T wave width ( T w ), T wave slope-to-amplitude ratio ( T S/A ) 20 and four indices of T wave morphological variability based 21 on time warping ( d w , d NL w , d a and d NL a ). Serum potassium 22 and calcium levels and heart rate were measured along 23 HD. Results: d NL a was the marker most strongly correlated high inter-patient variability was observed in the pattern of such relation-29 ships. This variability, accentuated during the ﬁrst time 30 points, was reproduced in the simulations and shown to be 31 inﬂuenced by differences in transmural heterogeneity. Con-32 clusion: Changes in serum potassium and calcium levels 33 and in heart rate strongly affect T wave descriptors, par-34 ticularly those quantifying morphological variability. Sig-35 niﬁcance: ECG markers have the potential to be used for 36 monitoring serum ion concentrations in ESRD patients.

ships. This variability, accentuated during the first HD time 30 points, was reproduced in the simulations and shown to be 31 influenced by differences in transmural heterogeneity. Con- 32 clusion: Changes in serum potassium and calcium levels 33 and in heart rate strongly affect T wave descriptors, par- 34 ticularly those quantifying morphological variability. Sig- 35 nificance: ECG markers have the potential to be used for 36 monitoring serum ion concentrations in ESRD patients. 37 Index Terms-Calcium, heart rate, hemodialysis, potas- 38 sium, time warping, transmural heterogeneity, T wave mor- 39 phology, ECG, in silico modeling. 40 I. INTRODUCTION 41 C HRONIC kidney disease represents a global health burden, 42 with an estimated 10% of the population being affected. 43 All stages of this disease, but particularly the late ones, are asso- 44 ciated with increased mortality and decreased quality of life [1].  ] during the HD session: the first one at the HD onset 123 and the next three samples every hour during the HD session 124 (Fig. 1, h 0 to h 3 in red). The 5th blood sample was collected 125 at the end of the HD (h 4 , at minute 215 or 245, depending 126 on the patient). A 6th blood sample, h 5 , was taken after 48 127 hours, immediately before the next HD session, but was not 128 analyzed as part of this work. The study protocol was approved 129 by the ethical committee (CEICA ref. PI18/003) and all patients 130 signed the informed consent form. Table I shows the population 131 characteristics.

B. ECG Pre-Processing 133
Pre-processing of ECG signals from ESRD patients included 134 band-pass filtering (0.5-40 Hz) to remove baseline wander as 135 well as muscular and powerline noise. A wavelet-based single-136 lead delineation method was used for QRS detection and wave 137 delineation of each of the twelve leads [20]. 138 Principal component (PC) analysis was spatially applied to 139 the T waves of the eight independent leads [21] to enhance the 140 T wave energy. The coefficients defining the PC transformation 141 were obtained from the eigenvectors of the 8 × 8 inter-lead auto-142 correlation matrix estimated by including all segmented T waves 143 within a 10-minute window at the end of the HD session, as this 144 is the time when the patient was discharged from hospital with 145 restored serum ion levels, thus being an acceptable reference for 146 ambulatory monitoring. The first PC computed by projecting the 147 ECG recording was used for subsequent ECG analysis, as it is 148 the transformed lead where the T waves have maximal energy, 149 thus allowing better morphological characterization. 150 The T waves in the first PC were delineated using the single-151 lead delineation algorithm described in [20]. The

179
For the patients' ECGs, reference T waves were calculated from 180 the MWTW at the end of the HD session.

181
The T wave for a given HD time point was expressed 182 as f s (t s ) = [f s (t s (1)), . . ., f s (t s (N s ))] and the reference  Fig. 2 (a) shows f r and f s , with their respective 188 time domains, t r and t s . Let γ(t r ) be the warping function that 189 relates t r and t s , such that f s (γ(t r )) denotes the time-domain 190 warping of f s (t s ) using γ(t r ). The square-root slope function 191 (SRSF) transformation was used to find the optimal warping 192 function by warping the SRSFs of the original T waves [22]. 193 This transformation is defined as: The optimal warping function was determined as the one mini-195 mizing the SRSF amplitude difference: shows reference (blue) and investigated (red) T waves obtained from an ECG segment during HD. Panel (b) shows the warped T waves, which have the same duration while keeping the original amplitude. Panel (c) depicts the warped T waves after normalization by their L2-norms. The area (cyan region) between both T waves in panel (d) represents d w , which quantifies the total amount of warping. The black solid line is the linear regression function γ * l (t r ) best fitted to γ * (t r ). The marker d NL w quantifies the non-linear warping by computing the area of the dashed grey region between γ * (t r ) and γ * l (t r ).
A dynamic programming algorithm was used to obtain the 197 function γ * (t r ) that optimally warps f r (t r ) into f s (t s ). This 198 function is shown in Fig. 2 (d). The warped T wave, f s (γ * (t r )), 199 is shown in Fig. 2 (b), together with the reference T wave, f r (t r ). 200 The descriptor d w , shown in Fig. 2 (d), was used to quantify 201 the level of warping required to optimally align the T waves 202 f s (t s ) and f r (t r ): where s d = N u r n=1 (γ * (t r (n)) − t r (n)) + N r n=N u r +1 (t r (n) − 204 γ * (t r (n))) is used to account for the sign, with N u r denoting 205 the number of samples in the T wave upslope.

206
The amplitude descriptor d a was computed from the area 207 contained between f r (t r ) and f s (γ * (t r )) normalized by the 208 L2-norm of f r (t r ), thus quantifying amplitude differences after 209 time warping the two T waves: where s a = N r n=1 (f s (γ * (t r (n))) − f r (t r (n))) is used to ac-211 count for the sign.

212
The warping parameter d w has a positive sign if the analyzed 213 T wave is globally widened during the warping procedure to 214 fit the reference T wave, and a negative sign if the T wave 215 is compressed. In the amplitude domain, d a is positive if the 216 warped T wave has larger amplitude than the reference T wave, 217 and negative if the T wave has smaller amplitude.

218
The marker d w incorporates information from the linear and 219 non-linear warping required to fit the two T waves in the time 220 domain. The non-linear component of d w can be quantified as: where γ * l (t r ) (black line in Fig. 2 (d)) was derived by linearly fitting γ * (t r ) through the least absolute residual method.

223
The marker d NL a was defined by computing the L 2 norm of 224 the difference between L 2 -normalized versions of f r (t r ) and 225 f s (γ * (t r )): The set of all morphology-based T wave markers analyzed in 227 this study included: whereX andȲ are the sample means. , [26]. The correlation coefficient after removing 248 the effects of Z in both X and Y was calculated as: The correlation coefficient between X and Y after removing the 250 effects of two variables Z 0 and Z 1 was calculated as: where Z 0 , To test for significant differences in Different proportions of endocardial, midmyocardial and epi-273 cardial cells were simulated in a total of 22 combinations with 274 10% variations in the proportion of each cell type: endocardial 275 layer ranging from 10% to 50%, midmyocardial layer from 10% 276 to 50% and epicardial layer from 20% to 80%. We used the 277 notation Cuvw, where C stands for the word "case" and u, v 278 and w denote the first digit of the proportions of endocardial, 279 midmyocardial and epicardial cells, respectively (e.g. C334 280 represents the case with 30%, 30% and 40% of endocardial, 281 midmyocardial and epicardial cells, respectively).

282
A train of 10 stimuli was applied to the first cell of each 283 fiber with a basic cycle length of 1000 ms and amplitude equal 284 to 1.5 times the diastolic threshold. The initial state for each 285 simulation was pre-calculated from a single cell simulation, 286 where the values of the model state variables after 1000 paced 287 beats were considered as representative of the cell at steady 288 state. To compute electrical propagation, a finite element-based 289 software [32] was used with a time step of 0.01 ms and space 290 discretization of 0.01 cm.

291
Unipolar pECGs were computed as described in previous 292 studies [29] using the expression: where is a constant proportional to the ratio of intracellular and 294 extracellular conductivities, V (x, y, z) is the transmembrane 295 potential and r(x, y, z) is the distance between each source 296 point (x, y, z) in the 1D fiber and the virtual electrode (x, y, z ) 297 located, in this study, 2 cm away from the epicardium in the fiber 298 direction: The last pECG beat of each simulated condition was delin-317 eated using the same delineation method mentioned above [20].

318
The time-, amplitude-and morphology-based T wave de-319 scriptors of section II-C were measured over those pECGs.
where Y c;a is the average value of the T wave marker Y from  (Fig. 3). A decreasing trend 387 of d w and d NL a from the maximum to the minimum level of 388 [K + ] was observed in all the simulated cases (panels j and n). 389 Monotonic trends of T w , T S/A , d a and d NL w were observed in most 390 of the simulated cases (panels i, l, m and n).

391
The bottom panels in Fig. 4           shown in Table III, [10], [42]. In this study, 536 we could observe such behavior in some of the ESRD patients' 537 recordings, as illustrated in the bottom panels of Fig. 3. However, 538 a decrease in T wave amplitude could not be consistently mea-539 sured for all patients, but large inter-individual variability was 540 noted in the relationship between [K + ] and T wave amplitude. 541 Other studies have analyzed the effects of [K + ] changes on 542 the width, slope and amplitude-to-slope ratio of the T wave 543 as well as the ratio of the T wave amplitude to the R wave 544 amplitude [14], [15], [17], [43], [44]. The main limitation of 545 these descriptors is that, even if some of them may show a high 546 degree of correlation with the level of [K + ], their changes cannot 547 be exclusively attributed to [K + ] variations, as confirmed in our 548 study by including in the analysis additional confounders like 549 variations in [Ca 2+ ] or HR.

550
Regarding the analysis of the T wave shape, a morphology 551 combination score (MCS) based on T wave asymmetry, flatness 552 and notching [45]- [47] has been used to analyze its relationship 553 with [K + ] in a primary care population [48]. A clear association 554 between MCS and [K + ] could only be found among individuals 555 with [K + ] in the range 2-4.1 mM, but not among those with 556 [K + ] in the range 4.2-6 mM. In ESRD patients, we found 557 that morphological variability, specifically quantified by our 558 analyzed T wave marker d NL a was closely related to serum [K + ] 559 in a wide range of values, covering both hyper-and hypokalemic 560 values.

561
As for the effects of [Ca 2+ ] variations on the ECG, a recent 562 large-scale study has found that low [Ca 2+ ] values are asso-563 ciated with clinically relevant QT prolongation in the general 564 population [7]. In chronic patients undergoing HD, changes 565 in [Ca 2+ ] have been found to be negatively correlated with 566 changes in the last part of the ECG repolarization measured by 567 the T-peak to T-end interval [49]. In this study, we showed that 568 the full repolarization duration measured by T w indeed presents 569 an inverse relationship with [Ca 2+ ] after removing the effects 570 of other confounders. Nevertheless, such a relationship between 571 T w and [Ca 2+ ] was not as strong as that of other markers like 572 d w reflecting temporal variations in T wave morphology.  ECG repolarization and, in particular, T wave morphology [22],

631
[54]- [59]. Our study confirms these observations on the impact 632 of transmural heterogeneities on T wave width, amplitude and 633 shape characteristics, not only at physiological electrolyte con-634 centrations but also at high and low [K + ] and [Ca 2+ ] levels and at 635 different heart rates. Even if transmural heterogeneities can con-636 tribute to inter-individual differences in the T wave response to 637 electrolyte and HR variations, other ventricular heterogeneities, 638 like interventricular, apicobasal or anteroposterior, may play a 639 relevant role, which should be assessed in further studies.

640
Our results on the sensitivity of T wave morphological mark-641 ers with respect to variations in transmural heterogeneities, and 642 more specifically to the proportion of epicardial cells within 643 the ventricular wall, are aligned with computational findings 644 presented by Janusek et al. [54], which demonstrated the influ-645 ence of epicardial cells on the development of T wave alternans, 646 a form of repolarization variability [54]. The contribution of 647 variations in the midmyocardial layer to T wave morphology 648 has been shown in a recent study too [56].

650
This study investigated 20 ECG recordings of ESRD patients 651 during an HD session, with 5 blood samples available along HD. 652 Future studies should investigate application of the proposed 653 methods to larger numbers of patients and, if possible, with more 654 available blood samples during the full 48-hour ECG recording. 655 This would allow more robust assessment of the relationship 656 between changes in T wave markers and specific variations 657 in [K + ], [Ca 2+ ] or HR, potentially using nonlinear regression 658 statistical techniques [60], [61].  [64]. 663 In the present study, [Mg 2+ ] was not investigated due to the 664 unavailability of serum [Mg 2+ ] levels.

665
Our electrophysiological simulations considered human 666 transmural ventricular fibers. Future research is aimed at ex-667 tending the investigations of the present study to include sim-668 ulations in bi-ventricular models embedded in patient-specific 669 torso models, from which more realistic ECGs can be computed. 670 This research will additionally allow exploring the role of other 671 types of ventricular heterogeneities, on top of transmural ones, 672 on the T wave response to electrolyte and HR variations.