Scaling and correlation properties of RR and QT intervals at the cellular level

We study complex scaling properties of RR and QT intervals of electrocardiograms (ECGs) with their equivalences at the cellular level, that is, inter-beat intervals (IBI) and field potential durations (FPD) of spontaneously beating human-induced pluripotent stem cell-derived cardiomyocyte (hiPSC-CM) aggregates. Our detrended fluctuation analysis and Poincaré plots reveal remarkable similarities between the ECG and hiPSC-CM data. In particular, no statistically significant difference was found in the short- and long-term scaling exponents α1 and α2 of RR and QT intervals and their cellular equivalences. Previously unknown scaling properties of FPDs of hiPSC-CM aggregates reveal that the increasing scaling exponent of QT intervals as a function of the time scale, is an intrinsic feature at the cellular level.


Validation of DFA implementation
We use publicly available MIT NSR database [1] to validate our implementation of detrended fluctuation analysis (DFA) [2] algorithm against the PhysioNet's DFA software package [1]. Using PhysioNet's algorithm and official annotations, RR and QT end intervals were calculated for each recordings. Each time series was filtered with the same preprocessing procedure described in the Method section. The preprocessed RR and QT intervals were then used as inputs to our own DFA implementation and PhysioNet DFA program. In both implementations, following input parameters were used: • linear detrending.
• maximum window size = N /4, where N is the length of the time series.
Both implementations produce fluctuations F (s) as a function of the scale s. Scaling exponents, the slope in the log-log plot, are calculated by a linear regression function, scipy.stats.linregress, in the scale range of 4-30 beats.
For 18 RR interval time series, mean ± standard deviation of the scaling exponents α computed with our implementation was 0.998±0.142, and that of scaling exponents computed using PhysioNet DFA was 0.999 ± 0.148. Root mean square error (RMSE) and normalised RMSE (NRMSE), i.e., RMSE normalised by the mean of the PhysioNet results, were computed to make comparison between results of two implementations: RMSE = 0.01880 and NRMSE = 0.01881. Similarly, for 18 QT interval time series, our implementation yielded α = 0.693 ± 0.126 and PhysioNet DFA α = 0.692 ± 0.124, with RMSE = 0.00433 and NRMSE = 0.00626.
The minor numerical differences arise in selecting a set of window sizes or scales s, at which fluctuations F (s) are computed. PhysioNet DFA selects a set of window sizes such that the ratio between successive window sizes is 2 1/8 ranging from predefined minimum window size to maximum window size. On the other hand, our implementation chooses number of window sizes as an input and generates a set of window sizes that are equally spaced in log scale and round them. Therefore, as fluctuations are computed at different set of scales, numerical deviations arise when computing the slope.
An important difference between the two implementations is in handling the window sizes that do not divide the time series evenly. While PhysioNet DFA discards the end of the data that is not covered by non-overlapping windows (there is also an option of using sliding windows), our implementation follows the algorithm illustrated by Kantelhardt et al. [3], in which the time series is divided into non-overlapping windows also starting from the end so that no data is neglected. Therefore, the different procedures introduce the deviations in the fluctuation.
As supported by low RMSE and NRMSE between the results of two implementations, the numerical deviations in the results do not affect our analysis, thus showing that our implementation produces valid DFA results.

Distributions of the measures
We report the distributions of the measures presented in Table 1 and 2 of the manuscript. Table  1 contains Poincaré plot related measures, such as mean, ratio of SD1 to SD2, and Pearson's correlation coefficient r for RR, QT intervals and IBIs and FPDs, of which the histograms are presented in Fig. S1 -S3. Table 2 contains scaling exponents from DFA. Histograms of scaling exponents α 1 and α 2 of RR, QT, IBI, and FPD are presented in Fig. S4.
Distributions are checked in order to make a valid comparison of the presented measures between ECG and hiPSC-CM data. The variables from IBI and FPD have higher variances in general compared to RR and QT. That is, equal variance is not satisfied between compared variables; therefore, Welch's t-test has been employed. Normality of the variables were also checked with Shapiro-Wilk test. When the normality requirement is not met, instead of t-test, non-parametric Wilcoxon ranksum test was used.

Standard HRV measures of IBI and FPD
Here we report the standard time-and frequency-domain measures that are often used to assess HRV. We chose the measures defined in Table S1. The results are summarised with min-max, median, and interquartile range for each measure in Table S2 and S3. The results were obtained with the HRV toolkit software package from PhysioNet [1]. IBI and FPD time series were preprocessed according to the Methods section in the manuscript and used as input to the HRV toolkit program.