A cross species thermoelectric and spatiotemporal analysis of alternans in live explanted hearts using dual voltage-calcium fluorescence optical mapping

Objective. Temperature plays a crucial role in influencing the spatiotemporal dynamics of the heart. Electrical instabilities due to specific thermal conditions typically lead to early period-doubling bifurcations and beat-to-beat alternans. These pro-arrhythmic phenomena manifest in voltage and calcium traces, resulting in compromised contractile behaviors. In such intricate scenario, dual optical mapping technique was used to uncover unexplored multi-scale and nonlinear couplings, essential for early detection and understanding of cardiac arrhythmia. Approach. We propose a methodological analysis of synchronized voltage-calcium signals for detecting alternans, restitution curves, and spatiotemporal alternans patterns under different thermal conditions, based on integral features calculation. To validate our approach, we conducted a cross-species investigation involving rabbit and guinea pig epicardial ventricular surfaces and human endocardial tissue under pacing-down protocols. Main results. We show that the proposed integral feature, as the area under the curve, could be an easily applicable indicator that may enhance the predictability of the onset and progression of cardiac alternans. Insights into spatiotemporal correlation analysis of characteristic spatial lengths across different heart species were further provided. Significance. Exploring cross-species thermoelectric features contributes to understanding temperature-dependent proarrhythmic regimes and their implications on coupled spatiotemporal voltage-calcium dynamics. The findings provide preliminary insights and potential strategies for enhancing arrhythmia detection and treatment.


Introduction
Temperature plays a pivotal role in regulating the behavior of excitable biological media, exerting precise control over both intra-and extracellular ionic dynamics (Fedorov et al 2008, Crampin et al 2020).In the context of mammalian hearts, alterations in temperature have been demonstrated to trigger irregular cardiac rhythms, such as tachycardia and fibrillation, as evidenced by the elevated incidence of cardiac arrest during the winter months (Katsouyanni et al 2008, Törö et al 2010).Despite these findings, moderate hypothermia is generally recommended as a therapeutic approach to mitigate injuries associated with cardiac arrest and anoxic brain damage (Bigelow et al 1950, Reuler 1978, Mattu et al 2002, Moore et al 2011).However, it is important to note that hypothermia also brings about challenges such as hypotension, edema, and heightened sympathetic activity during the cooling and rewarming phases, making the determination of an optimal temperature range for clinical interventions a subject of ongoing debate.
Temperature-induced electrical instabilities are typically quantified using AP-based parameters, including duration (APD), amplitude (APA), and conduction velocity (CV) (Crozier 1926, Bjørnstad et al 1993, Gizzi et al 2017).These metrics are significantly influenced by thermal variations, particularly in the context of period-doubling bifurcations and complex spatiotemporal alternans patterns (Fenton et al 2013).Alternans, which refer to the beat-to-beat oscillation of these features, have been extensively studied in the literature, both in terms of time and spatial characteristics (Pastore et al 1999, Watanabe et al 2001, Gizzi et al 2013).Despite numerous contributing factors, the bidirectional interaction between transmembrane voltage (V) and Ca ions is a critical element in the genesis of alternans (Shiferaw andKarma 2006, Edwards andBlatter 2014).This interaction occurs through V-regulating Ca transients via L-type channels, which, in turn, affect APD through Ca-sensitive membrane currents, as recently uncovered in the case of Ca-activated K + channels (Kanaporis and Blatter 2023).Furthermore, cell-based alternans mechanisms can be V-driven (Qu et al 2010, Huang et al 2020) or Ca-driven (Edwards andBlatter 2014, Kanaporis andLothar 2015), and they have been demonstrated to serve as precursors to life-threatening phenomena, such as spiral and scroll wave reentries leading to cardiac fibrillation.Alternans patterns, whether APD-driven or Ca-driven, can exhibit spatial in-phase behavior (concordant alternans-CA), where the entire tissue alternates between long/high and short/low phases, or out-of-phase behavior (discordant alternans-DA), where neighboring subregions alternate between long and short phases or high and low phases with transitions occurring via non-alternating areas (nodal lines) (Hörning et al 2017).Complex behaviors and higher-order rhythms have been confirmed through various spectral analyses (Loppini et al 2022), and advances in imaging techniques (Laurita and Singal 2001, Marina-Breysse et al 2021, Uzelac et al 2022), that allow simultaneous recordings of voltage and calcium dynamics, have enabled novel multi-field and multi-scale characterization of spatiotemporal alternans phenomena (figure 1).
The degree of calcium alternans is often defined in terms of transient Ca amplitude (CaA), as Ca accumulates within specific organelles (Qian et al 2001).Several alternative indicators have been proposed to quantify this phenomenon, including (i) the alternans ratio (AR), which quantifies large and small CaA for various threshold values (Wu and Clusin 1997, Choi and Salama 2004, Visweswaran et al 2013, Kanaporis and Lothar 2015, Martinez-Hernandez et al 2022); and (ii) the normalized calcium transient amplitude (∆CaT), calculated as the maximum amplitude value between two consecutive beats (Choi and Salama 2004, Herron et al 2012, Uzelac et al 2017, Hernández-Romero et al 2019).It is worth noting that optical mapping techniques often result in signals with a suboptimal signal-to-noise ratio (SNR), as extensively discussed in the literature (Mironov et al 2006).As a result, the calculation of these indices is susceptible to significant errors.This study will analyze signals with a minimum SNR value of approximately 2 dB, particularly for Ca traces that accumulate during pacing.Such a low SNR may hinder the accurate identification of the beat onset and its relative amplitude (figure 1(a)).
In this context, the present study introduces a novel method for detecting alternans based on integral features, such as the area under the curve, which have been employed in limited prior studies, specifically focusing on V traces (Kanaporis and Blatter 2023).In particular, we extend the scope of this method to encompass concurrent recording of epicardial fluorescence signals for both voltage and calcium acquired during pacing-down restitution protocols at varying thermal states (from normal body temperature, i.e. 37 • C, to mild We conducted a cross-species thermoelectric methodological analysis involving rabbit and guinea pig hearts within the temperature range of [37 ÷ 25] • C, with 4 • C steps.Our findings exhibit that integral indices serve two key purposes: (i) simplifying signal processing complexity while corroborating established temperature-dependent alternans patterns (Gizzi et al 2017), and (ii) enabling the construction of comprehensive alternans maps for the entire ventricle, even in case of noisy signals.Additionally, we conduct spatiotemporal correlation analyses, extending the notion of characteristic spatial length to Ca signals, a previously unexplored aspect in the literature.The proposed analysis involved comparing the statistical properties related to CV measurements for the two coupled signals.Our results show that an effective combination of diverse integral information improves predictive capabilities concerning alternans onset and period-doubling bifurcation, shedding light on the critical influence of varying thermal states in the initiation of arrhythmias (Fenton et al 2013, Loppini et al 2019, 2021).
The manuscript is organized as follows.In section 2, we describe experimental protocols, processing of recorded fluorescence data, and the methods adopted to quantify alternans, CV, and statistical spatiotemporal correlation.In section 3, we present the results of the analysis on three cases of study: a rabbit, a guinea pig, and a pathological human tissue.Finally, section 4 is devoted to discussions, comments, and future perspectives.

Experiment description
Fluorescent dual optical mapping technique was performed on ex-vivo epicardical surface of New Zealand white rabbit (female, 2-3 kg, n = 1) and guinea-pig (n = 1) whole-heart and for human right ventricle endocardical surface (n = 1) under approved IACUC from Georgia Tech.Simultaneous calcium handling and electrical wave propagation with high spatiotemporal resolution was obtained (figure 1(a)-(raw data)).The optical system was previously described (Fenton et al 2009, Ji et al 2017, Uzelac et al 2017) and composed by six high-power LEDs for excitation, three of them used for voltage imaging and the others for calcium one, both coupled with OD4 650/20 nm and OD4 550/20 nm excitation filters, respectively.A custom-designed apparatus was used to regulate LEDs operations and intensity.Fluorescent images were captured by a single camera (Photometric Evolve 128 EMCCD) with a time-multiplexing method that switch the different excitation LEDs to record voltage and calcium dynamics at spatial resolution of 128×128 pixels at 500 fps digitalized with a 16-bit dynamic range A/D.Tissue preparation and stimulation protocol are referred to (Uzelac et al 2017): animals were anesthetized with ketamine/xylazine/ace-promazine and then injected with heparin.After 5 min, euthanasia was induced with pentobarbital.Hearts were then quickly removed via a left thoracotomy and perfused retrogradely via the aorta with cardioplegic solution.Then the hearts were immersed in a chamber kept at 37.0 ± 0.3 • C and perfused with Tyrode's solution at a pressure of about 60 mmHg maintained by a peristaltic pump.Contraction was avoided by means of Blebbistatin (Fedorov et al 2007) at concentration of 1.8 uM and 2.5 uM for rabbit and guinea pig hearts, respectively.For the optical mapping system, the heart was stained with the voltage-sensitive dye JPW-6003 and Ca-sensitive dye Rhod-2.
External bipolar stimuli were delivered to the heart's apex with a full pacing-down restitution protocol.This involved delivering current pulses while progressively decreasing the PCL, as illustrated in figure 1(b).The PCL reductions started from physiological values (approximately 1 second) and extended down to faster arrhythmic intervals (around 200 milliseconds), a range dependent on the specific bath temperature conditions and species.For each PCL, around 100-150 stimuli were delivered to allow the system to reach steady state.The programming sequence was synchronized with the internal camera trigger clock using an Arduino (Uno R3) so that each pacing stimulus was delivered at a known time.
Experiments were conducted at four thermal states, i.e. from 37 • C to 25 • C with a step of 4 • C, in which the heart was placed in an oval, water-jacketed, and heated bath filled with Tyrode's solution at a specific temperature.A temperature-controlled reservoir, equipped with a heat pump, was directly connected with the bubble trap (Radnoti) and the thermal bath, creating a closed-loop system that ensures thermal isolation.Precisely, the temperature of the Tyrode was measured inside the bubble trap using a temperature probe (ThermoFisher) with an accuracy of 0.05 • C and well-calibrated.To minimize heat loss, the bubble trap was connected to the heart by a tube with a small diameter of 2 mm and a total length of 30 cm.This configuration enabled a precise temperature control of the whole apparatus within ±1 • C, considering a room temperature of 24 • C, so that temperature gradients were negligible.The temperature was further adjusted to match the steady-state temperature measured inside the bubble trap.The bath temperature was considered as a reference since the flow rate was relatively low.
The analyzed human heart was obtained from the recipients of orthotopic heart transplantation at Emory University Hospital.For human tissue preparation and protocol we refer to Uzelac et al (2022), Iravanian et al (2023).

Data processing
Optical mapping signals were analyzed by a customized and semi-automatized MATLAB (MathWorks, 2022) algorithm.Raw data were separated in V and Ca signals and then processed to eliminate the pixel-by-pixel drift, by subtracting a linear interpolation of minimum values and the intrinsic noise, using a low-pass Gaussian spatial filter (window = 7, variance = 2) combined with a moving mean average filter (5 samples for both types of traces).Since Ca recordings are noisier than V ones (figure 1(a)), signals were further filtered via stacking (Uzelac and Fenton 2015), where each image was summed up at specific equidistant time interval, i.e. each PCL.Finally, a pixel-by-pixel normalization was performed.
Time series with artifacts and low SNR were eliminated using a custom binarized mask that identifies the regions of interest (Ω in figure 1) for subsequent analysis.To detect Ω, the algorithm designates background pixels as zero (excluded pixels) through an edge detection technique employing Sobel filters (Kanopoulos et al 1988).Additionally, pixels with peaks in the low-frequency spectrum (0.1÷1.15 Hz), indicative of photobleaching or motion, are identified due to their low-frequency content.For both guinea pigs and human tissues, an additional criterion was applied based on the mean amplitude of the envelope of the Hilbert transform (thr = 0.18 for guinea pigs, 0.25 for human data).These thresholds were chosen for each experimental dataset to provide an optimal trade-off between including pixels potentially affected by motion artifacts while minimizing errors in alternans detection.It is important to note that thresholds are dataset-specific, given the inherent differences between experiments, including variations in heart conditions and dye concentrations.This setting was essential to increase the number of informative pixels within the mask.The resulting masks varied by PCL and signal type for the human heart due to increased noise.Elsewhere, a uniform mask was used across all recordings and thermal conditions.

Alternans measures
Spatiotemporal analysis of V-Ca alternans has been traditionally quantified by subtracting two consecutive APDs-CaDs (figure 1(b)), thus building two iterative maps: and identifying alternans according to the thresholding mechanism described in the previous section.In our analysis, we used thr APD ≡ thr CaD = 5 ms including recording time resolution and analyses artifatcs.This threshold value ensures that rapid changes in the signals are adequately detected, considering the raw signals Whereas the mentioned techniques are widely recognized and used for V alternans detection, the same is not directly applicable to Ca signals due to (1) intrinsic noise, (2) accumulation features, and (3) high space-time variability.As highlighted in figure 2, identifying a unique threshold is unreachable or inconsistent among the traces, becoming even more pronounced at lower bath temperatures or in diseased cases.
As a novel contribution, we introduce the integral area under the curve as a parameter of interest for the detection of V-Ca alternans, defining the two ratio indices: where, AV n and ACa n are the V and Ca areas for the nth beat, while AV n−1 and ACa n−1 represent the corresponding areas for the preceding beat (n − 1).These integrals are calculated for each pixel (x, y) ∈ Ω of the mapped field (figure 1(b)).
To enable cross-species comparison, we normalized each dataset using the minimum and maximum area values observed across all temperatures for each species.Subsequently, for each pixel, the normalized area sets are sorted, ensuring that the odd-numbered beat (n − 1) exhibits a larger area than the even-numbered beat (n).Such a list is then separated into two distinct subgroups: odd (AV n−1 , ACa n−1 ) and even (AV n , ACa n ).Finally, to generate global restitution curves, time and space average was applied to each subgroup.The associated spatial distribution is assessed by calculating the standard deviation of the averages, as described in the following section.
A dedicated sensitivity analysis was performed to find optimal threshold values thr ARV and thr ARCa for alternans onset detection: we defined the bifurcation point (PCL BP ) in the restitution curve as the PCL that presents an AR between mean values for odd and even beats Accordingly, we obtained a threshold value of 12% for alternans detection, aligning with documented trends in the literature (Gizzi et al 2017).In particular, a propensity for alternans onset is obtained at elevated values of PCL with decreasing temperature.For a detailed description of the methodology, we refer to section 1 supplementary materials.It is worth noticing that AV-DI and ACa-DI restitution curves are not shown since it is not straightforward to compute the DI integral value for Ca signals.
To cope with well-established APD-based analyses (Gizzi et al 2013), we consider a red-white-green color scale to visualize alternans maps on whole heart geometries highlighting complex temperature-dependent patterns.Multiple regions alternate in-phase (spatially concordant alternans) or out-of-phase (spatially discordant alternans).Figure 1(c) shows the 3D-plot of spatially concordant Ca-alternans (region a), and the transition to an out-of-phase region (b) passing through a nodal line (solid black trace).Inequalities in equation ( 2) were also used to calculate the percentage of the alternating pixels inside the mask.Accordingly, alternans onset-both concordant and discordant-was identified at that PCL for which at least 10% of the pixels display alternating area values.

Conduction Velocity
Conduction velocity in cardiac tissue is known to be temperature-dependent, creating complex spatially DA patterns (Fenton et al 2013, Gizzi et al 2017).CV restitution, as the underlying mechanism of this pathological mechanism, was rigorously proven in both experimental and theoretical studies (Cao et al 1999, Pastore et al 1999, Loppini et al 2019), and further generalized varying temperature (Loppini et al 2021).
In the present contribution, we computed CV for both voltage and calcium excitation, developing a custom MATLAB algorithm that: (i) select four time windows (i.e.four beats); (ii) identify wavefront isochrones at 70% of depolarization for consecutive time intervals (figure 1(b)) calculated as V t (P j ) = 0.7 for a generic point (P j ) and time t; (iii) iteratively select pairs of pixels on two consecutive isochrones whose distance vector is normal to the direction of propagation; (iv) iteratively evaluate the local CV for each pair of pixels and consecutive isochrones over each time window.Pixels inside a previously detected square box (51 × 51 pixels) were chosen for the calculation as representative of the whole tissue behavior.Median values were used to build CV-restitution curves.For further description, we refer to figure S4, supplementary materials.
In addition, we conducted a second comparative benchmark analysis in which CV was calculated via the circle method (Siles-Paredes et al 2022): local CV was calculated as a function of differences among local activation time t (with a threshold of 0.2) across the chords of a circle (radius r = 8 pixels) centered in each (x, y) position within the selected box, with an orientation angle θ, so as CV(θ) = 2r/∆t.The circle method evaluates the resulting CV by performing a geometrically weighted average of CV(θ) over an angular range Θ ∈ [θ − θ/2, θ + θ/2].We refer to Siles-Paredes et al (2022) for a detailed discussion of equations.Comparable trends are obtained between median CV-restitution curves.

Spatiotemporal correlation analysis
Correlation features in both V and Ca data were investigated to advance our predictive power for period-doubling bifurcation and generalize what previously found in Loppini et al (2019Loppini et al ( , 2021)).We calculated the characteristic spatial length (decay length L 0 ) for four consecutive beat cycles, i.e. a depolarization front and one repolarization back within the 2D mapped field (figure 1(b)).Specifically, Spearman's rank correlation coefficient (Corder and Foreman 2014) was evaluated among all pixels separated by a distance⃗ r within a square box (spatial resolution of 0.04 cm/pixel for rabbit and 0.03 cm/pixel for guinea pig tissue) in the mapped domain: where, Y is a generic notation for ranks of the two types of signal value at time t, positions ⃗ x (Y A ) and ⃗ x +⃗ r (Y B ); ⟨•⟩ t stands for the time average computed over one beat cycle, and σ is the standard deviation.Computing L 0 for four consecutive beats, as in previous studies (Loppini et al 2019) allows highlighting higher-order rhythms.Further, we defined the normalized characteristic length ( L0 ) for the maximum L 0 value (L + 0 ) using the average value (L 0 ) over the four beats: Accordingly, we identified alternans onset when a first net deviation from zero in L0 was present.Characteristic lengths during sustained fibrillation were further computed, verifying known values for V (Loppini et al 2019) and extending new observations for Ca at different thermal states.

Statistical analysis
To quantify the cross-species V-Ca alternans degree, we performed Wilcoxon signed rank tests (Hollander et al 2013), in the context of alternating phases and nodal lines.The test aimed to examine whether significant differences existed in the cross-species alternans patterns between paired observations.In addition, the AR metric in equation ( 2) was employed as a parameter of interest for the comparative analysis.To obtain the AR-distributions, we (i) selected a common PCL across species that corresponds to the 70% of the maximal PCL in each restitution curve-the minimum value to capture alternans regime-, (ii) we separated pixels that belong to nodal lines from alternating ones, and (iii) we computed the averaged AR values within each zone in the spatial domain for each analyzed couple of beats.The resulting AR distributions were visually presented using box plots, providing a clear overview of alternans patterns for each experimental condition.

Data representation
Restitution curves of APDs and area values are reported as mean ± standard deviation.Tables with observation numbers are provided in appendix (table 1).

Results
This section comprehensively analyzes the spatiotemporal features associated with the onset and progression of alternans patterns under PCL restitution.We compare various electrophysiological indices, including restitution curves and alternans maps.Additionally, we enhance the characterization of calcium dynamics across different thermal conditions.Our exploratory analysis will compare the following indices for each species at different bath temperatures: (1) Ca-V area restitution curves in relation to APD; (2) Ca-V AR-distribution spatial alternans maps; (3) Ca-V characteristic length calculation; (4) CV restitution curves; (5) cross-species statistical comparison.
To emphasize noteworthy temperature-dependent phenomena, we employ a color-coding scheme: transmembrane voltage is represented in a blue-magenta scale, Ca levels are depicted in a red-black scale, and alternans maps are visualized using a color spectrum spanning from green to red.

Restitution curves and alternans maps
In figure 3, we present surface-averaged restitution curves for APD (a) and AV, ACa values (b) considering the rabbit tissue dataset.The results show known trends in the literature with restitution curves decreasing and bifurcations occurring at lower PCLs.Besides, significant differences can be noted between AV and ACa values, in particular when temperature variations are considered.Notably, at the baseline temperature of  2)): white pixels correspond to nodal lines, red to positive phase, and green to negative one (black identifies out-of-mask or noisy pixels).Number of observations shown in appendix.37 • C, neither APD nor alternation in AV is observed; in other words, no bifurcation can be detected.However, the baseline ACa-based restitution curve shows a clear bifurcation around PCL = 200 ms.This feature becomes more pronounced as temperature decreases, highlighting the role of calcium in sustaining cardiac arrhythmias, mainly when the substrate is highly heterogeneous, i.e. at lower temperatures.It is worth noticing that ACa restitution curves present less dispersed distributions, emphasizing the underlying nonlinear coupled dynamics among V and Ca measured via fluorescence optical mapping.Such an effect may be linked to the alteration of ionic kinetics and time-constant changes that affect the loading time history of Ca within subcellular reservoirs.
According to the expected temperature-dependent APD trends, area restitution curves show increasing values as temperature decreases (Bjørnstad et al 1993, Filippi et al 2014, Loppini et al 2019).In addition, we note that the amount of dispersion decreases consistently, i.e. the standard deviation in our calculation.This reduction is attributed to the use of integral indicators rather than relying solely on traditional duration-based metrics.Notably, at lower temperatures, specifically at 29 • C and 25 • C, the APD restitutions display a wide dispersion, while the area amplitudes tend to cluster closely around their respective mean values.
Figure 3(c) presents spatial alternans maps for specific PCLs within rabbit hearts at four different temperatures, characterized in terms of duration (APD and CaD) and AR (ARV and ARCa), as described by equation ( 2).These ARV and ARCa patterns serve to reinforce the observations elucidated in panels 3(a) and (b).Interesting, the formation of alternans initially occurs in localized SCA clusters, subsequently diffusing throughout the tissue and ultimately evolving into SDA.Compared to the more coherent Ca-phase behavior, the V-signals display significant heterogeneity, often manifesting one or more nodal lines at lower PCLs.In addition, alternans occurred earlier under non-physiological temperatures, consistent with APD/CaD alternans maps (Fenton et al 2013, Loppini et al 2021).
Comparing alternans maps generated by both methods, we observe that the integral index expands the number of pixels analyzed (as further detailed below), yielding similar spatial patterns.Such enhanced performance is particularly evident for Ca, with a more consistent visualization of complex patterns under non-physiological thermal states.
Figures 4(a) and (b) shows guinea pig restitution curves considering the same structure as discussed for the rabbit dataset.Also, we observe similar features as temperature decreases, supporting the adoption of integral area indices to strongly reduce intrinsic noise and highlight the onset and evolution of bifurcations, SCA, and SDA patterns.For the two species, area restitution curves show rising values and larger dispersion as temperature decreases, indicating higher beat-to-beat spatiotemporal variability at low temperatures.Besides, figure 4(c) corroborates and expands previous observations in figure 3(c), where the number of usable pixels between the two methodologies results much higher-almost the whole ventricle-for the ACa.
Figure 5 compares the cross-species percentage of ARV and ARCa alternating pixels, offering a clear and comprehensive vision of the studied temperature-dependent spatiotemporal dynamics.Notably, the number of alternating pixels at low temperatures strongly increases at low PCLs.Furthermore, the Ca signals exhibit a higher percentage compared to the V ones, which aligns with the findings presented in figures 3(c) and 4(c).

Conduction Velocity restitution curves
Conduction velocity restitution curves are computed for V traces, CV V (PCL), at four temperatures as shown in figure 6 for rabbit and guinea pig, respectively.In general, we measured CV V ∈ 10 ÷ 70 cm s −1 and found decreasing values of CV at low PCLs and temperatures consistently with previous studies (Mironov et al 2008, Loppini et al 2021).Further, CV V presents steeper averaged slopes and reduced dispersion at physiological temperature.Such a trend highlights a larger dispersion of repolarization and a consequent propensity to develop alternans at lower temperatures.
From a cross-species perspective, we notice that rabbit tissue shows comparable CV V at 29 • C and 25 • C, which is not the case for guinea pig.Such an observation suggests the presence of intricate spatiotemporal coupling mechanisms between voltage and calcium dynamics (CV Ca are shown in figure S5, supplementary materials).These findings emphasize the need for further investigation into subcellular mechanisms.
The robustness of the area method is supported by observing CV V oscillations at low PCLs in both tissue samples just before conduction block or fibrillation occurs.CV oscillations are associated, in general, with complex spatial alternans patterns and alternating characteristic spatial lengths, as discussed in the following section, calling for advanced analytical techniques.

Spatiotemporal correlation analysis
Statistical correlation length analysis unveils spatiotemporal synchronization patterns and merges information derived from CV and APD restitution curves (or AV/ACa, as investigated in this study) into a unique composite metric.For instance, the correlation length L 0 decreases at low temperatures and PCLs. Figure 7 shows L 0 values for rabbit ((a), (b)) and guinea pig ((d), (e)) across four temperatures, further providing the normalized characteristic length L0 (equation ( 4)) ((c), (f)).Our cross-species analyses indicate that L 0 exhibits consistent average values across the species tested.V-data show a range of L 0 within the interval of 35 to 3 cm (figures 7(a) and (d)); similarly, Ca-data range from 28 to 2 cm (figures 7(b) and (e)).Notably, the behavior of rabbit V signals deviate from these general trends.For example, at 25 • C, L 0 partially overlaps with the values observed at 29 • C.An in-depth analysis was also conducted by comparing APD and CV values to investigate such a condition, revealing a consistent increase in APD for a constant CV.
Concerning the normalized characteristic length (figures 7(c) and (f)), we observed a net increase in L0 oscillations over four beats, indicative of a rise in beat-to-beat variability as temperature and PCL decrease.In particular, we identified alternans onset when the deviation from the baseline is L0 > 20% (see section 3.4 and appendix), thereby supporting and extending the findings of previous studies (Filippi et al 2014, Loppini et al 2021): the increasing oscillations of L0 at low PCLs result from the nonlinear spatiotemporal behavior associated with developing discordant alternans during paced-down restitution protocols.
We highlight low L 0 values for Ca signals to support our findings.These results follow and expand what was already known in the literature for V traces, demonstrating that a higher likelihood of conduction block and reentry arises as characteristic length decreases.The calculation of L 0 during sustained fibrillation resulted in L 0 ∼ 1 cm for voltage and ∼0.9 cm for calcium (see figure S6, supplementary materials).

Detection of alternans onset
Identification of alternans onset, quantified in terms of PCL as temperature decreases, is crucial to understand spatiotemporal evolution of the aforementioned pro-arrhythmic mechanisms.Besides classical analysis on the slope of restitution curves (Nolasco and Dahlen 1968), in this study we discuss three different alternans onset detection methods: • (1 • ) PCL fulfilling the condition AR ⩾ 12% for averaged odd-even area values; • (2 • ) PCL where a net deviation of L0 from baseline ( L0 > 20%) was observed; • (3 • ) percentage of alternating pixels with cut-off = 10% (see figure 5).Consistently with previous trends (table 2, appendix), Ca alternans onset often occurs first, over the two species, and appears at higher PCLs as temperature decreases.At lower temperatures, critical nonlinear effects may lead to inverse trends, e.g. at 29 • C alternans onset (2 • method) is higher than at 25 • C.

Cross-species statistical comparison
Figure 8 shows representative box plots of AR values for each species, temperature, and signal type.Statistically significant differences (p < 0.05) appear, resulting in distinct AR distributions during alternating phases between calcium and voltage for each temperature and species.Notably, differences are observed solely in rabbit nodal lines.Additionally, increased sample variability is evident in guinea pigs compared to rabbits.Performing the same statistical tests and considering signal type as a differentiating factor, i.e. all Vor Ca-distribution across all temperatures, statistically significant differences are obtained only for Ca-distributions at 37 • C and 33 • C. In contrast, alternating phases consistently exhibit significant differences across all conditions.

Human ventricular representative case
Human endocardium was analyzed at two selected temperatures, i.e. 37 • C and 29 • C, to test the integral-based method in diseased scenarios toward a comprehensive cross-species analysis.Figures 9(a) and (b) provides restitution curves showing APD, ARV, and ARCa.We observe similar trends as in previous sections and recent human studies (Koller et al 2005, Selvaraj et al 2007), i.e. steeper slopes of restitution curves and higher amplitudes at non-physiological temperatures.Evidence is provided using integral indices emphasizing a low V standard deviation, a comprehensive analysis of Ca traces, and highlighting significant differences between odd and even beats.
The diseased case appears evident by V-Ca spatial alternans maps (figure 9(c)).Alternating pixels do not show clustering effects as in healthy rabbit and guinea pig, further confirming the complex multi-scale role of Ca dynamics (Qian et al 2001, Clusin 2003, Díaz et al 2004, Pruvot et al 2004, Aistrup et al 2009).In this sense, the lack of clustering suggests differences or failure in sub-cellular and inter-cellular calcium handling, reflected in tissue-scale asynchronous alternans.In particular, interspecies variability in sarcoplasmic calcium content (Díaz et al 2004) or sarcolemmal conductances over species could be implied in loss of clustering.Further, a decrease in gap-junction permeability in the pathological human heart may lead to a lack of calcium oscillation synchronization in nearby areas (Qian et al 2001).It is worth mentioning that loss of clustering in restricted tissue areas was also due to low SNR in human traces.However, by proper filtering (see section 2.2) and supervised inspection of our post-analysis, we verified that loss of clustering was not an artifact due to noise.Indeed, nearby pixels consistently show different alternans patterns (see figure S7, supplementary materials).

Discussion
In this study, we investigated the spatiotemporal activation response of whole ventricles from different species (including rabbit, guinea pig, and diseased human heart tissue), applying a full pacing-down restitution protocol.A dual optical mapping technique was adopted to record synchronous voltage and calcium fluorescence signals at four temperatures.A thorough analysis across species compared duration-based methods with more robust integral spatiotemporal approaches.
Our main findings can be recast as follows: (i) calcium alternans onset anticipates voltage in a wide range of temperatures; (ii) the area method provides a comprehensive index for alternans detection as it takes into account both changes in amplitude and duration of each type of signals, giving results consistent with predictions based on generalized spatiotemporal correlation analyses.

Temperature effects on voltage-calcium bidirectional coupled dynamics
Cardiac alternans refer to a pathologic alternation into V-Ca signals, and it has been observed in various arrhythmogenic conditions, including hypothermia (Egorov et al 2012).Alternans onset and development in response to temperature changes still require a thorough understanding to gain predictability and controllability of arrhythmic phenomena.voltage-calcium bidirectional coupled dynamics play a critical role in alternans genesis via V-driven and/or Ca-driven mechanisms (Hirayama et al 1993, Pruvot et al 2004, Shiferaw and Karma 2006).The second hypothesis, in which alternans originates from a primary disruption of calcium cycling then followed by APD alternans (Visweswaran et al 2013, Uzelac et al 2017), is here investigated as a function of thermal tissue state.Data suggest that Ca-alternans arise at higher PCLs (section 3.4), and it sustains even in the absence of V-alternans, e.g. at baseline physiological temperature for rabbit dataset (figure 3(b)).Temperature affects calcium handling mechanisms, as extensively described in terms of voltage, changing signal morphology, gradually increasing duration, area values, and CV (figures 3(b), 4(b) and 9(b)), and (figure 6), slowing down as PCL and temperature reduce, triggering increasingly complex spatial alternans patterns (figures 3(c) and 4(c)).Furthermore, differences in CV values support characteristic length analyses calling for advanced investigations of the underlying electronic couplings.
It is worth noticing that features and trends described above are comparable across species, raising the question of common underlying processes marking alternans onset and evolution as early predictors of cardiac arrhythmias.

Integral indices as a general robust method
Optical mapping techniques have prompted the quantitative study of local activation and repolarization timings, as well as calcium handling, with high spatiotemporal resolution in different Langendorff-perfused heart species (Herron et al 2012, Crocini et al 2014, Hernández-Romero et al 2019, Uzelac et al 2022, Marchal et al 2023).Output signals, however, can be highly noisy, especially concerning calcium traces.In the case of Ca-alternans, the small-transient peak could merge with the large-transient descending part, making it challenging to detect Ca-transient duration or amplitude.
There are good reasons why APD, Ca amplitudes, or durations are measured during alternans: Ca amplitudes can influence L-type Ca inactivation or NCX, and APD is a reasonable estimation of refractoriness related to conduction blocks and wavebreaks.In this scenario, the area under the curve could merge amplitude and duration features in one synthetic index, returning generalized information on spatiotemporal alternans dynamics.

Limitations and future perspectives
The present thermoelectric cross-species spatiotemporal analysis is limited to one specimen per species: a healthy rabbit and guinea pig whole heart, a diseased human ventricle.Given augmented statistic significance, we expect to broaden the study by considering a larger number of tissue samples and categorizing genetic factors (e.g.age and sex) (Haq et al 2023) as well as different pathological states.Of note, the effect of anti-arrhythmogenic drugs should be accounted for, e.g. late sodium current blockers, gap junction openers, and beta-blockers (Larson et al 2022, Song et al 2023), in combination with thermal bath changes, usually not considered in optical studies.
Clinically, arrhythmias precursors generally result in pulsus (mechanical) and T-wave (electrical) alternans, pathological conditions that can lead to ventricular fibrillation and sudden cardiac death (Chudin et al 1999, Pruvot et al 2004).In such a scenario, simultaneous V-Ca optical mapping features combined with ECG ones-with a cause-effect relation-at different thermal states may help in understanding appropriate drug therapy.Besides, the methodology can be applied to other biological systems and organs, such as the gastrointestine tract, where the voltage-calcium couplings and electro-mechanical feedback are crucial for proper organ functioning and disease detection (Calder et al 2022).
Finally, we foresee the integration and assimilation of such novel spatiotemporal features (i.e.ARCa spatial heterogeneities and calcium spectral information) into thermo-electric cardiac models (Pullan et al 2005, Loppini et al 2018, Ruiz-Baier et al 2019, Barone et al 2020a, 2020b), also considering recent developments in uncertainty quantification approaches (Peirlinck et al 2019, Pathmanathan et al 2020, Brandstaeter et al 2021, Nitzler et al 2022) that account for statistical indices associated with restitution curves and spatiotemporal alternans maps.However, the differences in tissue dynamics could also be due to scar distributions, which are highly variable from one heart to the next.Accordingly, an accurate study of the optical ultrastructure (Loppini et al 2022) could be further associated with the proposed integral analyses toward a more comprehensive investigation of the underlying mechanisms extracted from fluorescence optical signals.Preliminary results in this direction are provided in figure S8, supplementary materials.

Conclusions
This work supports the temperature-dependent nature of pro-arrhythmic regimes in cardiac electrophysiology, which affect spatiotemporal voltage and calcium dynamics.Through concurrent optical fluorescence mapping of whole-heart voltage and calcium signals and employing a variety of methodological approaches for the post-processing analysis, this study finds that an integral-based index is a robust and reliable strategy for quantifying intricate spatial accumulation patterns and morphology alterations resulting from macroscopic temperature changes and repeated electrical stimulation.This approach could combine traditional methodologies for alternans detection while facilitating signal pre-processing complexities, e.g.data filtering, beat selection, signal oversampling, and identification of the correct threshold value for beat duration calculation.In such a perspective, extrapolating experimental characteristics across multiple species, varying thermal conditions, and correlating them with length-scale transitions determined through statistical correlation matrices offers a promising avenue to elucidate novel multiscale and nonlinear couplings-not limited to the sole cardiac context.Such insights finally aim to improve our predictive capabilities of early detection of cardiac arrhythmias.

Figure 1 .
Figure 1.Schematic diagram of the whole-heart signal analysis, the rabbit case is adopted as representative example.(a) Simultaneous optical mapping recordings for voltage (V) and intracellular calcium (Ca).(b) Filtered recordings for a generic pixel in the 2D mapped field; definition of signal duration (APD and CaD) and area (AV and ACa) values at each PCL.A thresholding technique (thrAPD,CaD = 5 ms, thrARV,ARCa = 0.12) is shown to quantify alternans according to (1b) and (2).Zoomed view of a selected box (dotted lines) with isochrones (white lines) at 70% of depolarization with representative pairs of pixels (grey points) along the contour normal direction to wavefront propagation (white arrows).(c) Spatial pattern of Ca area alternans.Color code refers to alternans phases in terms of ARCa (2).On the right, a spatiotemporal representation from the (a)-(b) line time series highlights transition phases: spatially concordant Ca-alternans (region a-red), out-of-phase region (b-green) passing through a nodal line (solid black trace).

Figure 2 .
Figure 2. Time traces (2 beats shown) of synchronous normalized voltage (left) and calcium (right) signals from one common pixel between traces at four different thermal states (organized columnwise) for rabbit dataset (n = 1).Common PCL values (ms) for V and Ca signals at each temperature are shown on each time series.

Figure 3 .
Figure 3. Restitution curves of rabbit epicardium (n = 1) for (a) APD, (b) AV, and ACa at four temperatures (37 ÷ 25 • C) organized column-wise (n. of observations: 98060).Data are reported as mean values (points) ± standard deviation (shaded) for two consecutive beats: long-short APD and odd-even area values.Uniform axes are adopted, comparing amplitudes and PCL.(c) Spatial area alternans maps: ∆APD and ∆CaD (left), ARV and ARCa (right), organized row-wise for selected PCLs.AR alternans are detected at 12% thresholding (equation (2)): white pixels correspond to nodal lines, red to positive phase, and green to negative one (black identifies out-of-mask or noisy pixels).Number of observations shown in appendix.

Figure 4 .
Figure 4. Restitution curves of guinea pig epicardical surface (n = 1) for (a) APD, (b) V (n. of obs: 112840) and Ca (n. of obs: 63750) areas at four temperatures organized column-wise.Data are reported as mean (points) ± standard deviation (shaded) for two consecutive beats: long-short APD, odd-even area values.Uniform axes are adopted, comparing amplitudes and PCL.(c) Spatial area alternans map: ∆APD and ∆CaD (left), ARV and ARCa (right), organized row-wise for selected PCLs.AR alternans are detected at 12% thresholding (equation (2)): white pixels correspond to nodal lines, red to positive phase, and green to negative one (black identifies out-of-mask or noisy pixels.Number of observations are provided in appendix.

Figure 5 .
Figure 5. Bar plot of the percentage of alternating pixels for V (blue) and Ca (red) at each PCL and temperature for rabbit (a) and guinea pig (b) dataset.The percentage is computed considering the ratio between the number of pixels whose AR ⩾ 12% (1 − Aeven/A odd ) and the total amount of in-mask pixels.Dotted lines refer to 10% as the cut-off value for alternans onset detection (see table 2, appendix).

Figure 7 .
Figure 7. Characteristic length (L0) for rabbit (a), (b) and guinea pig (d), (e) cross-comparing V-Ca signals at each PCL.Solid lines (a), (b), (d), (e) represent L0 mean values (n. of obs: 4) over four beats; minimum and maximum values are shaded.Normalized L0 values (4) are provided for V (solid blue) and Ca (solid red) traces: rabbit (c) and guinea pig (f) datasets are shown column-wise at four temperatures.Uniform axes are provided for cross-comparison.

Figure 8 .
Figure 8. Temperature-dependent cross-species (rabbit-guinea pig) box plots comparison for voltage (left) and calcium (right) AR nodal lines (white) and alternating phases (grey for spatially concordant alternans, red-green for spatially discordant alternans) mean values.Solid red lines are the median values of each distribution, while asterisks represent p < 0.05.The dotted grey lines correspond to the assumed threshold value of 12% (equation (2)).

Figure 9 .
Figure 9. Restitution curves of human endocardium: (a) APD, (b) V-Ca area at 37 • C and 29 • C organized column-wise.Curves show mean (points) ± standard deviation values (shaded) for long and short beats, APD, and odd and even beats, AV and ACa.Uniform axes are provided for comparison.(c) Spatial alternans map for V (left) and Ca (right) at selected PCLs.The alternating phase was identified as positive or negative (|AR(x, y)| ⩾ 12%) with a red-green scale otherwise white for non-alternating pixels, and black for out-of-mask ones.Number of observations is provided in appendix.

Table 1 .
Summary tables of the number of observations for each dataset.