Heartbeat dynamics in adrenergic blocker treated conscious beagle dogs

Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

Contents lists available at ScienceDirect

Journal of Pharmacological and Toxicological Methods j o u r n a l h o m e p a g e : w w w. e l s ev i e r. c o m / l o c a t e / j p h a r m t ox

Original article

Heartbeat dynamics in adrenergic blocker treated conscious beagle dogs Dingzhou Li, Alan Y. Chiang ⁎, Christine A. Clawson, Bradley W. Main, Derek J. Leishman Global Statistical Sciences and Toxicology, Lilly Research Laboratories, Eli Lilly and Company, P.O. Box 708, Greenfield, IN 46140, USA

A R T I C L E

I N F O

Article history: Received 28 February 2008 Accepted 2 May 2008 Keywords: Detrended fluctuation analysis Multifractality ECG Heart rate variability Spectral analysis

A B S T R A C T Introduction: Adrenergic blockade as a treatment for chronic heart failure (CHF) has proved effective, but its pharmacological mechanism on CHF remains unclear. In the past two decades, studies on heart rate variability (HRV) have reported that CHF patients generally have a reduced temporal complexity in heart rate variability. On the other hand, adrenergic blockers have been shown to restore such complexity. Fractal analysis is a novel and efficient tool to explore the adrenergic blockade effect on HRV. This paper applies the detrended fluctuation analysis (DFA) and multifractal DFA (MF-DFA) methods in an attempt to understand the effect of adrenergic blockade on cardiac dynamics in conscious beagle dogs. Methods: DFA and MF-DFA analysis are conducted on RR interval data generated from telemetry instrumented dogs receiving a combination of 15 mg/kg nadolol and 5 mg/kg phenoxybenzamine orally administered at the 22nd and 34th hour in a parallel design (n = 12). All dogs had approximately 48 h of beat-to-beat heart rate measurements recorded in the left ventricle. Complexity measures for heartbeat series are compared between the blocker and vehicle group. We also compute traditional statistics for HRV and spectral parameters and examine their correlation with fractal analysis. Results: When compared to the vehicle group, the adrenergic blocker group had: 1) longer RR intervals (p = 0.02) and lower beat-to-beat variability (p = 0.04); 2) decreased low frequency (LF) and high frequency (HF) power (p = 0.03), and higher LF-to-HF ratio; 3) larger middle-range scaling exponents (p b 0.01); 4) broader multifractal spectra (p = 0.03) with higher dominant singularity indices (p = 0.02). Discussion: Our results show that 1) adrenergic blockade alters the sympathovagal balance; 2) adrenergic blockers enhance the complexity of the cardiac dynamics; 3) the adrenergic blockade effect on cardiac dynamics is primarily the attenuation of small fluctuations in RR intervals. Fractal analysis also has the potential to be applied to early QT diagnosis. © 2008 Elsevier Inc. All rights reserved.

1. Introduction Adrenergic blockade as a treatment for chronic heart failure (CHF) has proved effective in both improving the left ventricular function and reversing pathological states in the heart (e.g., Eichhorn, 1998; Waagstein, Hjalmarson, Varnauskas, & Wallentin, 1975). However, little is known to date about the exact pharmacological mechanism of adrenergic blockade on CHF. In the past two decades, studies on heart rate variability (HRV) have reported that CHF patients in general have a reduced temporal complexity in their heartbeats (e.g., Goldberger et al., 2002; Malik, Farrell, Cripps, & Camm, 1989; Peng, Havlin, Stanley, & Goldberger, 1995), while administration of adrenergic blockers is able to restore such complexity (e.g., Chiu, Chan, Chu, & Lin, 2007). Measurements of HRV have also been used to assess autonomic tone in a variety of clinical and nonclinical conditions. Many drugs act directly or indirectly on the autonomic nervous system, and heart rate variability analysis has been shown to be a useful tool to explore the influence of various agents on sympathetic and parasympathetic activity (Task Force, 1996). A useful measure of HRV is SDNN, the standard deviation of all ⁎ Corresponding author. Tel.: +1 317 433 3049; fax: +1 317 651 9964. E-mail address: [email protected] (A.Y. Chiang). 1056-8719/$ – see front matter © 2008 Elsevier Inc. All rights reserved. doi:10.1016/j.vascn.2008.05.005

normal-to-normal (NN) interbeat intervals over 24 h. Other common measures of HRV under time domain include SDSD, the standard deviation of the successive differences of the RR intervals, and SDANN, the standard deviation of the RR interval calculated over short periods (usually 5 min). Though variations in heart rate may be evaluated by a number of methods, the use of spectral techniques to quantify short-term heart rate fluctuations on the order of seconds to minutes has helped define the autonomic contributions to beat-to-beat control of heart rate. When HRV is analyzed in the frequency domain, the variance of the signal is decomposed into its sinusoidal components at various underlying frequencies. This can be visualized in a plot of the HRV power spectrum in which the area under the curve at each frequency is equal to the variance explained at that frequency. Power spectral analysis of HRV has been applied to study the diurnal variation of autonomic nervous activity in dogs (Matsunaga et al., 2001) and other species (Kuwahara, Hiraga, & Nishimura, 1998). However, the HRV analysis remains a challenge because the signals obtained are under constantly varying conditions. Traditional methods mentioned above assume that the statistical properties of the signals are stationary functions of time. This is not true for the RR intervals — e.g., the statistical properties of the heart rate change when a subject

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

Fig. 1. Illustration of DFA and MF-DFA procedures. The slope of the bottom panel would be the scaling exponent, h(q).

Fig. 2. Typical time series and their DFA plots. From top to bottom: white noise, Brownian motion, and 1/f noise.

119

120

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

Fig. 4. Response to isoproterenol challenge in previous studies with and without betaadrenergic blockade presented as change in index of contractility (dP/dtmax). Fig. 3. Response to phenylephrine challenge in previous studies with and without alpha-adrenergic blockade presented as change in systolic blood pressure.

rises to a standing position or goes to sleep. This means that the statistical properties of those RR intervals are not constant in time; they are not served up neatly as independent, identically distributed (i.i.d) random variables. It can be easily shown that two heart rate time series with identical values of their means and standard deviations can have dramatic differences in their dynamics. Therefore, traditional methods based upon statistical properties such as mean and standard deviation may give spurious results and hence lack power of differentiating healthy and pathological states, not to mention detecting the drug effect. Even the spectral analysis can lead to misleading results in the presence of nonstationarity as temporal information is masked in the calculation of Fourier coefficients. Thus more reliable and robust methods are needed for the analysis of nonstationary and nonlinear time series. Peng et al. (1995) proposed a novel analysis of HRV, referred to as the detrended fluctuation analysis (DFA), that addresses the issue of nonstationarity. The scaling properties of HRV are self-similar on a scale of minutes to hours. In other words, the heartbeat time series have similar statistical properties on a wide range of different time scales, which is termed fractality in the nonlinear dynamics and chaos theory. The DFA algorithm quantifies fractal-like correlation properties by calculating the scaling property of the root-mean-square fluctuation of the integrated and detrended time series (Fig. 1). The results of DFA can be used to quantify the dynamics of the signals. For example, The 1/f noise has an exponent value of h = 1.0, the white noise h = 0.5, and the Brownian noise h = 1.5 (Fig. 2). The white noise suggests an uncorrelated random time series, and the Brownian noise is the integration of the white noise. The 1/f noise can be interpreted as a compromise between the completely unpredictable white noise and the much smoother landscape of Brownian noise. There has been ample evidence (e.g., Amaral et al., 2001; Peng et al., 1995) that heartbeat series possess the behaviors of typical 1/f noise. It has also been shown that increased randomness of heart rate patterns, i.e., a high degree of sinus arrhythmia of non-respiratory origin, is associated with a worse prognosis (Stein, Domitrovich, Hui, Rautaharju, & Gottdiener, 2005). The theoretical relationship between DFA and weighted power spectra has been demonstrated by Willson, Francis, Wensel, Coats, and Parker (2002). Kantelhardt et al. (2002) shows that fluctuations in physiological signals at different times usually demonstrate different scaling properties. It is very rare for the DFA plot to be a nearly straight line. Hence, one

or a few fractal exponents are no longer sufficient to characterize the long heartbeat series (usually lasting at least 24 h). A more useful method in this case is the multifractal DFA (MF-DFA) that examines a spectrum of scaling properties for the span of the time series. The advantage of MF-DFA is that it no longer assumes the existence of a single scaling properties (“monofractal”) in the time series of interest. Rather, it automatically searches for and calculates the properties and richness of a variety of structures in a time series. Therefore, a fuller comprehension of temporal variability in the heart rate is allowed by this approach, which is a more suited model for cardiovascular physiology as endogenous and exogenous factors make it constantly change on both small and large time scales. In addition, the problem of selecting an appropriate time window for HRV analysis is avoided because MFDFA explores a range of temporal resolutions in the heartbeat in a selfiterating fashion. This is especially useful when we do not have a priori information on the pharmacokinetic–pharmacodynamic aspects of compounds that have been targeted to the cardiovascular system. In general, a broad spectrum of multifractality implies a higher degree of complexity in the heartbeat dynamics. On the other hand, a narrow

Fig. 5. Systolic blood pressure in 30-minute windows.

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

121

2. Methods 2.1. Experimental design Twelve beagle dogs were fasted and an antibiotic regimen was initiated during surgery. General anesthesia was induced with a short acting barbiturate, intravenously to effect. Anesthesia was maintained with isoflurane inhalational anesthetic in oxygen delivered through a cuffed endotracheal tube. A Data Sciences International® (DSI) telemetry device (model number: TL11M3-D70-PCTP) was surgically

Fig. 6. Conventional HR quantities in 30-minute windows: (a) mean RR interval; (b) SDANN. The vertical dashed lines mark the times of two oral administrations more readable.

spectrum in the multifractality has been found in subjects with heart failure (Ivanov et al., 1999). In recent years multifractal analysis has been extensively employed to various types of time series such as seismic signals (Telesca, Lapenna, & Macchiato, 2005), financial data (Jiang & Zhou, 2007), and physiological signals (Amaral et al., 2001; West, Latka, Glaubic-Latka, & Latka, 2003). In this paper, we apply the DFA and MF-DFA, in addition to traditional methods, in an attempt to understand the dynamics of complex physiologic fluctuations arising from a preclinical cardiovascular study on adrenergic blockers. From these results, we evaluated the influence of adrenergic blockade on the heartbeat dynamics. Table 1 Group comparison of conventional HR quantities 1st administration RR-I (ms)

2nd administration

SDANN (ms) Skewness RR-I (ms)

Blocker 661.3 ± 97.5 156.6 ± 58.5 Vehicle 569.7 ± 71.3 195.4 ± 57.0 p-value 0.02⁎ 0.15 ⁎ p-values less than 0.05.

SDANN (ms) Skewness

1.09 ± 1.10 613.4 ± 110.3 110.0 ± 70.1 0.21 ± 0.33 538.6 ± 67.8 180.7 ± 58.6 0.002⁎ 0.12 0.04⁎

1.04 ± 1.26 0.25 ± 0.35 0.002⁎

Fig. 7. Average histograms of RR intervals: (a) before the treatment; (b) after the first administration; (c) after the second administration.

122

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

placed in each dog's abdominal cavity with the implant cannulae inserted into the left ventricular chamber of the heart and an abdominal artery. The left ventricular catheter was inserted by an abdominal approach through the diaphragm for the LVP placement in the left ventricle of the heart and no thorocotomy was performed. ECG leads were placed subcutaneously in the thorax. The transmitter was secured to the wall of the abdomen. The abdomen was closed in several layers. The 12 dogs were randomized into two groups. Adrenergic antagonists Nadolol (beta blocker) and Phenoxybenzamine (alpha blocker) were administered orally twice at doses of 12 mg/kg and

5 mg/kg, respectively, to 6 animals (Blocker Group). Previous studies were conducted to determine appropriate doses and frequency for adequate alpha and beta-adrenergic blockade. Phenylephrine was administered to challenge alpha-adrenergic blockade and isoproterenol was administered to challenge beta-adrenergic blockade. We determined that sufficient alpha and beta-adrenergic blockade was achieved with oral doses of 5 mg/kg phenoxybenzamine and 12 mg/kg nadolol, respectively, every 12 h. The contractility and systolic blood pressure results in previous studies are shown in Figs. 3 and 4. Fig. 3 shows that 4 h after dosing with Phenoxybenzamine (1.5 mg/kg), the phenylephrine challenge produced

Fig. 8. Spectral parameters in 30-minute windows. (a) raw LF power; (b) raw HF power; (c) normalized LF power; (d) normalized HF power; (e) LF:HF ratio. The vertical dashed lines mark the times of two oral administrations.

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

an attenuated increase in systolic blood pressure. This attenuation did not persist through 10 h. However, after two oral doses of Phenoxybenzamine (5 mg/kg 12 h apart), the response to phenylephrine challenge was reduced at both 2 and 10 h demonstrating adequate alpha-adrenergic blockade. Fig. 4 shows that after two oral doses of Nadolol (12 mg/kg 12 h apart), there was essentially no response to the isoproterenol challenge at 2 or 10 h after dosing demonstrating adequate beta-adrenergic blockade. Therefore, at least two oral doses of each compound were necessary to achieve consistent block as measured by lack of response to challenge agents. The systolic blood pressure results in this study that will be presented later also justify the need of two administrations. The vehicle treatment (100% ethyl alcohol, USP for Nadolol and 10% ethyl alcohol, USP with propylene glycol for Phenoxybenzamine) was administered at the same dose volume to the rest of 6 animals (Vehicle Group). Continuous data collection occurred from approximately 08:30 AM through approximately 06:00 AM the following day. Left ventricular pressure and systemic arterial blood pressure were recorded using the Data Sciences International and Ponemah Physiology Platform data acquisition and analysis systems. The peak value of the first derivative of left ventricular pressure (dP/dtmax) was used as an index of left ventricular inotropic state. Systolic, diastolic, and mean arterial pressures and arterial pulse pressure were derived by the data acquisition system from the arterial pressure signal. Heart rate was obtained from the left ventricular pressure signal. The ventricular and aortic pressure signals were digitized at a sampling rate of 500 and 250 Hz, respectively. The raw pressure signals from all dogs were digitized utilizing two Ponemah Physiology Platform 16-channel data acquisition and analysis systems. The systems reside within microcomputers in the cardiovascular laboratory at Eli Lilly and Company, Greenfield, IN. The data acquisition software was calibrated electronically at the beginning of each data collection period using values provided by the manufacturer of the DSI telemetry devices. Binary files containing the raw data were initially stored on the computer's hard disk and were then transferred to a qualified server following the completion of each data collection period. 2.2. Spectral analysis For spectral analysis, the RR interval series is first divided into 5minute non-overlapping windows. Summary statistics, including mean and standard deviation (SDANN), and power spectral parameters are calculated for each window. For the purpose of this illustration, heart rates higher than 400 beats per minute (bpm) or lower than 30 bpm are excluded from this analysis. We note that about 95% of the heartbeats in this experiment were actually in the range of 50–130 bpm. These exclusion criteria allow for dynamic heart rate increase in a transient state immediately after dosing. In addition, we plan to conduct new studies in the future with compounds that could stimulate the sympathetic tone, where the heart rate would usually be higher. A relatively higher threshold would allow for a consistent comparison. On examining the qualitative results, these cutoffs appeared to be not sensitive after we tried several sets of inclusion criteria, including 50– 200 bpm and 50–150 bpm. The high frequency (HF) and low frequency (LF) regions are defined as: 0.04–0.1 Hz for LF power and 0.1–0.6 Hz for HF power, which are consistent with results reported in the literature (Frazier, Moser, & Stone, 2001; Harada, Abe, Shiotani, Hamada, & Horii, 2005; Mastunaga et al., 2001). To eliminate artifacts at very low frequency (VLF; 0.003–0.04 Hz), we also calculate the normalized HF and LF by dividing the raw HF and LF power by the total power minus VLF power. It has been suggested that the ratio of LF-to-HF power mirrors the balance between sympathetic and vagal modulation of the heart rate (Task Force, 1996). Therefore, we also evaluate the LF-to-HF ratio in each window. We then average all the parameters over every six consecutive windows to obtain a 30-minute moving average. All the subsequent statistical analyses are applied to the derived moving average.

123

2.3. Detrended fluctuation analysis and multifractal detrended fluctuation analysis Due to the nonstationarity properties of a time series (e.g., circadian rhythms over 12 h and local drift in the order of several minutes), unless it is a complete random sequence, there are usually discernable trends on various time scales. By fitting a regression line on the trend, we can estimate the generic “residual” fluctuations in the time series. This is the basic procedure of the DFA method proposed by Peng et al. (1995). It is intuitively easy to expect that the average fluctuation in a time series is an increasing function of the time scale over which fluctuations are averaged (Fig. 1). However, it is observed that the average fluctuations in the RR intervals obey a power-law like function of the time scale, which is ranging from several beats to tens of thousand beats (Peng et al.,1995). Mathematically, such relationship is denoted by F ðnÞ~nh

ð1Þ

where F(n) is the average fluctuation on a given time scale (in units of heartbeats), n. The presence of a power-law is a signature of selfsimilarity in that fluctuations on all time scales are dictated by a single exponent h, meaning that there does not exist a characteristic time scale that stands out from others. Namely, fluctuations on any two time scales can be made very alike by a trivial scale transformation. In general, smoother time series are characterized by larger exponents h. The above is usually a simplified picture in assuming that all fluctuations, regardless of their magnitude, belong to a single family characterized by the exponent h. In reality, the behavior of the fluctuations depends not only on its time scale but on its magnitude. We should refrain from inferring that small and large fluctuations act equally in response to changes in conditions. In the framework of MFDFA (Kantelhardt et al., 2002), by the introduction of an order parameter q, different weights were added to fluctuations so as to reflect whether it is small or large fluctuations that dominate at a given time point. When taking into consideration the magnitude-based weighting parameters, we retrieve the power-law like dependence on time scale (i.e., selfsimilarity). The relationship is similar to Eq. (1): Fq ðnÞ~nhðqÞ

ð2Þ

where Fq(n) is the average qth order fluctuation on time scale n, and q can be any real number. When q is more positive (negative), larger (smaller) fluctuations contribute more in the calculation of h(q). Thus by plotting h(q) vs. q, we can examine the properties of fluctuations of different magnitudes. It is also worth noting that in practice Eq. (2) only holds in a finite range of q, which in this study is from −10 to 10. It is easy to observe that the DFA procedure proposed by Peng et al. (1995) is equivalent to the special case of MF-DFA at q = 2, and the scaling exponent, h, in their work is equal to h(2) in the MF-DFA framework. However, MF-DFA is not limited to q = 2 but rather examines h(q) over a range of q's, which provides more comprehensive understanding of time series for the aforementioned reasons. When a time series is characterized by a single scaling exponent h, it is said to be of “monofractality”; when it possesses a spectrum of Table 2 Group comparison of spectral parameters LF (ms2/Hz)

HF (ms2/Hz)

LF norm(%)

HF norm(%)

LF:HF

1st administration: Blocker 2951.0 ± 5081.0 Vehicle 1632.7 ± 888.0 p-value 0.40

22558 ± 22031 33120 ± 21097 0.17

13.3 ± 7.4 7.5 ± 4.2 0.03⁎

85.0 ± 7.4 91.6 ± 4.2 0.03⁎

0.17 ± 0.12 0.09 ± 0.06 0.03⁎

2nd administration: Blocker 1793.1 ± 3415.1 Vehicle 1620.0 ± 1634.8 p-value 0.15

13481 ± 18247 29128 ± 20060 0.03⁎

15.3 ± 7.4 9.0 ± 6.4 0.03⁎

84.8 ± 7.4 91.0 ± 6.4 0.03⁎

0.20 ± 0.02 0.11 ± 0.10 0.05⁎

⁎ p-values less than 0.05.

124

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

dimension associated with the index, f(α), indicates the abundance of such index in the span of the time series. A pronounced dependence of f(α) on α, as well as a wide range of α, is an indication of multifractality in the time series. The relationship between f(α) and α is usually a concave-downward spectrum with a peak at a certain value αmod and a width αmax − αmin, where αmax and αmin are maximal and minimal singularity indices in the spectrum, respectively. The modal (i.e. most probable) property of fluctuations is thus characterized by αmod, and the complexity in the whole spectrum by αmax − αmin. 2.4. Statistical analysis Mixed models with repeated measures (MMRM) is the primary statistical analysis. For mean, SDANN, and spectral parameters, values measured during the last 30 min before the first administration serve as the baseline variable in the model. Other factors include time and group (i.e., vehicle or blockers). For DFA and multifractal metrics, baseline values are calculated from about 40,000 heartbeats before the treatment. Responses are based upon 40,000 heartbeats after each administration. The significance level is set to p ≤ 0.05. 3. Results 3.1. Systolic blood pressure Fig. 5 shows time courses of mean systolic blood pressure in 30minute segments (average of six 5-minute windows). Before the treatment, the vehicle group and adrenergic blocker group had comparable systolic (134.5 ± 9.2 mm Hg vs. 138.1 ± 5.7 mm Hg for the vehicle and adrenergic blocker group, respectively; p = 0.48). After the first oral administration systolic blood pressure were not significantly affected, and the difference between the vehicle group and adrenergic blocker group remained non-significant (132.7 ± 19.3 mm Hg vs. 123.1 ± 4.8 mm Hg; p = 0.17). After the second administration systolic blood pressure was decreased in the adrenergic blocker group, which resulted in a significant difference between the two groups (129.7 ± 12.4 mm Hg vs. 109.5 ± 9.0 mm Hg; p = 0.02). 3.2. Mean RR interval and SDANN

Fig. 9. DFA plots of the vehicle and blocker groups: (a) before the treatment; (b) after the first administration; (c) after the second administration. The vertical dashed lines delineate the boundaries of different ranges.

h(q), it is said to be of “multifractality”. Alternatively, we can calculate the singularity index spectrum to quantify multifractality. Details of correspondence between the scaling exponent spectrum (i.e., h(q)) and the singularity index spectrum can be found in the original work (Kantelhardt et al., 2002). Briefly, a singularity index α measures the leading behavior of the time series at a given time point, and a noninteger value of α is evidence of fractality in the time series. The

Fig. 6 shows time courses of mean and SDANN of RR intervals in 30minute segments (average of six 5-minute windows). Table 1 summarizes statistics for both groups. The adrenergic blocker group had larger mean RR intervals and lower SDANN after each oral administration. It is generally agreed that adrenergic blockade can reduce the heart rate by suppressing the sympathetic tone (Eichhorn,1998; Nickerson,1949), which is consistent with our findings. However, there is a marked difference between the behavior of the mean RR interval and its beat-to-beat variability measured by SDANN. The increase in the RR interval (i.e., decrease in the heart rate) was significant after the first administration of the blockers, and the treatment effect diminished after the second administration. In contrast, the decrease in SDANN was not significant until after the second administration. Note that the error bars on the mean RR interval plot measure animal-to-animal variation, which is completely different from beat-to-beat variability. To further illustrate the group difference in mean and variability, histograms of RR intervals are shown in Fig. 7. Table 3 Group comparison of DFA slopes in short, middle, and long range 1st administration

Blocker Vehicle p-value

2nd administration

Short

Middle

Long

Short

Middle

Long

0.73 ± 0.11 0.73 ± 0.20 0.99

0.58 ± 0.08 0.50 ± 0.10 0.16

1.07 ± 0.04 1.03 ± 0.03 0.18

0.80 ± 0.17 0.89 ± 0.14 0.34

0.66 ± 0.11 0.40 ± 0.10 b0.01⁎

1.03 ± 0.05 0.99 ± 0.02 0.26

⁎ p-values less than 0.05.

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

Fig. 6(b) shows that SDANN in the vehicle group continued to increase after the second administration. This is because the second administration was only 1 h through the dark cycle (i.e., 18:00) when respiratory sinus arrhythmia became pronounced as a result of the dominant vagal tone. Large fluctuations in the heart rate are common during this stage, and hence SDANN is usually larger. On the other hand, the increase in SDANN was much lower in the blocker group, implying a perturbation on the neuroautnomic control by adrenergic blockade.

125

3.3. Spectral analysis Fig. 8 shows the time course of various spectral parameters in 30minute segments for approximately 48 h. Table 2 shows statistical comparison of the parameters between groups. After the first administration, normalized LF power was higher in the blocker group (13.3% vs. 7.5%), whereas normalized HF power was lower in the blocker group (85.0% vs. 91.6%). As a result, the LF-to-HF ratio was

Fig. 10. MF-DFA plots of the vehicle and blocker groups: (a), (b): before treatment; (c), (d): after the first administration; (e), (f): after the second administration.

126

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

higher in the blocker group (0.17 vs. 0.09). Similar results are also found after the second administration. In addition, the raw HF power was also lower in the blocker group after the second dosing. The dominance of HF over LF in both groups reflects strong vagal tones in dogs. However, there is a strong group difference in the relative distribution of LF and HF power, which indicates that adrenergic blockade altered sympapthovagal balance. 3.4. Detrended fluctuation analysis Fig. 9 shows (on a log–log scale) the average fluctuation, F2(n), versus the box size measured by number of beats, n. We calculate the slopes of the log F2(n) versus log(n) plot in three different ranges: short range (n = 4–20 beats), middle range (n = 20–50 beats), and long range (n = 50– 10,000 beats) for the following reasons. First, it can be observed from the DFA plot that the relationship is rarely a straight line across all the different box sizes. Rather we should fit the plot in a piece-wise fashion. Second, some authors (Willson et al., 2002) have shown that different parts on the DFA plot may correspond to spectral properties in different frequency ranges, and the cut-off value n and cut-off frequency f have a roughly reciprocal relationship: n ≈r− /f where r− is the mean heart rate in beats per second (bps). Based upon the definitions of LF and HF regimes and an estimate of r− ≈ 2 bps, we derived n = 20 and 50 as the cut-off values for the three different ranges. Group comparison of the DFA slope in different ranges are summarized in Table 3. As can be seen, the only difference was the higher middle-range slope in the blocker group (0.66 vs. 0.40). By recalling the relationship between DFA and spectral analysis, we know that middlerange correlation of heart rate fluctuations roughly corresponds to the LF power, which is believed to be a characteristic of sympathetic control. The higher slope value means that the heartbeat series on this scale (i.e. from 20 to 50 beats) is smoother in the blocker group. Not only is the group difference interesting, the exact values of slopes are also informative in disclosing the behavior of heartbeat series in each group. A slope value of 0.5 is the boundary between so-called antipersistent and persistent time series. In persistent time series (h N 0.5), larger (smaller) values are more likely to be followed by larger (smaller) values, while occurrences of large and small values tend to alternate in antipersistent time series (h b 0.5). The vehicle group had an average middle-range slope lower than 0.5, which is a characteristic of antipersistence. This is clearly attributed to respiratory sinus arrhythmia in dogs. In contrast, the blocker group had a slope larger than 0.5, meaning the heartbeats in this group is more like persistent time series. Hence large (small) RR intervals are more likely followed by large (small) intervals in the blocker group. Therefore, adrenergic blockade resulted in a qualitatively different structure in the heartbeats. 3.5. Multifractal detrended fluctuation analysis Fig. 10 shows the h(q) vs. q and f (α) vs. α calculated from 40,000 heart beats before and after each oral administration. Recall that these quantities are not independent but can be derived from each other. However, each plot is particularly useful for a certain aspect of multifractality. The scaling exponents, h(q)s, are shown in Fig. 10a, c, and e. A dependence of h(q) on q is indication of multifractality in the time series and hence high level of complexity in the dynamics. As can be seen from the figures, h(q) is dependent of q in both groups regardless of the treatment. It is evidence of multifractality in the cardiac dynamics of dogs. The h(q) dependence became stronger in the blocker group only after the second administration. The group difference also varies with q. When q is positive (negative), large (small) fluctuations dominate in the calculation of h(q). In general, small (large) h(q) characterizes larger (smaller) fluctuation in the scaling properties. Therefore, small fluctuations in the heat rate became even smoother after the second administration of adrenergic blockers, whereas large fluctuations were not significantly affected by the treatment.

Table 4 Group comparison of multifractal spectra 1st administration

Blocker Vehicle p-value

2nd administration

Peak

Width

Peak

Width

0.92 ± 0.05 0.85 ± 0.03 0.19

0.67 ± 0.21 0.64 ± 0.09 0.82

0.95 ± 0.11 0.84 ± 0.02 0.02⁎

0.92 ± 0.38 0.62 ± 0.12 0.03⁎

⁎ p-values less than 0.05.

The singularity index spectra (Fig. 10b, d, and f) are readily characterized by two parameters: the index at which f(α) peaks, αmod, and the width of spectrum defined by the difference between the maximum and minimum singularity indices, αmax and αmin. Statistical comparison of the peak index and width are shown in Table 4. The singularity index spectrum was wider in the blocker group after the second dosing (0.92 vs. 0.62), which is evidence of increased multifractality in the blocker group. The blocker group also has a higher value of modal singularity index in the multifractal spectrum (0.95 vs. 0.84), suggesting heartbeat series in this group appeared slightly smoother for the most of time. Thus, we conclude that adrenergic blockade had both enhanced the richness in the structure and changed the dominant behavior of the heartbeat series. 4. Discussion Under physiological conditions, complex dynamical behavior of the heart rate series has shown erratic fluctuations and patchiness. These irregularities have been historically ignored in conventional preclinical studies that focus on averaged quantities. In fact, these fluctuations are often inadvertently labeled as “noise” to distinguish them from the true “signal” of interest. By extending methods developed in nonlinear dynamics such as DFA and MF-DFA, hidden scaling structures can be disclosed in the physiologic fluctuations. The emergence of these structures is a result of the balance between two or more antagonistic regulations on the physiology. Such balance positions the cardiac system to a dynamical equilibrium around a “critical point” in that heartbeat fluctuations that are temporally remote are strongly correlated, contrary to what we would expect from a steady and stable system. Therefore, cardiovascular drug safety assessment based upon small segments of ECG or heart rate would be inadequate. Alteration of such balance can be due to either endogenous adaptation, exogenous pharmacological intervention, or both. Our findings raise the possibility that understanding the origin of such temporal structures and their alternations may elucidate certain basic features of heart rate control mechanisms, and have practical value in clinical and preclinical monitoring. Our results show that adrenergic blockade has much more profound impact on the cardiac dynamics than merely decreasing the heart rate in beagle dogs. First, the beat-to-beat variability was markedly reduced by the blocker treatment. In previous work, we demonstrated a reduction in beat-to-beat variability in the presence of increased heart rates from sympathetic stimulation (Chiang, Li, & Leishman, 2007). What we have found in this study is the opposite in terms of the relationship between heart rate and beat-to-beat variability. It appears that the modulation of heartbeats by sympathetic regulation can be more complicated than expected and may involve the antagonistic vagal tone. It would also be imprudent to conclude that there exists a rigid relationship between heart rate and its variability. Furthermore, we found by means of DFA and MF-DFA that the adrenergic blockers primarily changed the behavior of small fluctuations in RR intervals. The blocker group demonstrates a broader multifractal spectrum in the heartbeats implying an increased complexity in the cardiac dynamics, as well as an overall smoothing effect on the appearance of heartbeat series. The spectral analysis results also support the smoothing effect on heartbeat series by adrenergic blockade. A systematic study demonstrated that the 1/f-type noise and the monofractal structure in the heart rate series is likely to manifest the

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

balance between parasympathetic and sympathetic nervous systems (Struzik, Hayano, Sakata, Kwak, & Yamamoto, 2004). Hence, changes in the fractal structure can be indications of altered autonomic nervous regulation. There has been evidence that carvedilol, a type of beta blockers, can restore multifractality in patients with congestive heart failure (Chiu et al., 2007). Amaral et al. (2001) finds that healthy human subjects treated with metoprolol, a β1-selective antagonist, had decreased multifractality compared to placebo. Our findings in dogs tend to support increased multifractality caused by adrenergic blockers. Nevertheless, we realize that there are several major differences in the subjects as well as methods in the work of Amaral et al. (2001) from ours. First, the subjects in their study were humans, and it is common knowledge that dogs have much a stronger vagal tone than humans. Second, the drug effect on multifractality in dogs became significant at the second administration in our study, whereas the human subjects in their study only had a single dosing. The pharmacokinetics and pharmacodynamics can be very different for adrenergic blockade in humans from those in dogs. Third, the blocker in their study was primarily a beta blocker while our treatment was composed of a mix of alpha and beta blockers. A recent study in rats has shown a distinction between alpha- and beta-blockade in their abilities to affect the nonlinearity and complexity in the cardiac dynamics (Beckers, Verheyden, Ramaekers, Swynghedauw, & Aubert, 2006). Our study features two administrations about 12 h apart. We observed that different measures were affected at different times. Conventional measurements such as mean and spectral parameters changed shortly after the first administration of the blockers. In contrast, changes in fractality (DFA and MF-DFA) were only significant after the second dosing. While such difference can be partly attributed to the timing of administration in the circadian cycle, it is possible that different quantities, as they measure different aspects of the cardiac dynamics, respond to extraneous interventions differently. In our future studies, we aim to distinguish the circadian effect from repeated dosing effect by allowing for additional days of experiments. In doing so, we hope to gain insights into which quantities are to use when we assess cardiovascular safety in the context of repeated dosing. To our knowledge, our work is the first effort to analyze the cardiac dynamics of beagle dogs treated with adrenergic blockers by several different methods including spectral analysis and multifractality analysis. By applying methods that explore the multifaceted dynamics in the heartbeats, we hope to gain more comprehensive and precise understanding of how adrenergic blockers modulate the heartbeat dynamics on different time scales. Our findings demonstrate that DFA and MF-DFA can circumvent the issues from nonstationarity and nonlinearity and hence provide more reliable conclusions regarding cardiac dynamics as well as neuroautnomic regulations. As one of the next steps, we would examine the scenario where the adrenergic blocker treatment is followed by sympathetic tone stimulation, either by medication or physical routines, and compare the cardiac dynamics between this scenario and an untreated condition. By controlling the heart rate in both the vehicle and blocker group, we aim to explore the adrenergic blockade effect on more profound aspects of cardiac dynamics when the overall heart rate does not change. It is worth noting that fractal analysis also has the potential for applications on other types of cardiovascular investigations such as QT prolongation. Sanbe et al. (2005) found a higher scaling exponent of the QT intervals in transgenic cTnIR146G rabbits when compared to the non-transgenic control group, which demonstrates the ability of DFA to detect pathological change in the repolarization phase before the onset of more obvious ECG abnormalities. A recent study by Kazachenko, Astashev, and Grinevich (2007) shows that activities of both KCa channels in kidney Vero cells and of KV channels in mollusc neurons can be characterized as multifractal processes. It has been conjectured that increased complexity in the QT series, partly characterized by multifractality in ion-channel mechanisms, could lead to increased propensity of Torsade de Pointes (TdP) due to

127

increased heterogeneity of repolarization. Persistent increase in QT intervals as a precursor to TdP can be appropriately captured by monitoring the leading singularity index, αmod even when the QT increase itself is sometimes not obviously shown on the cardiogram. Singularity indices larger than 0.5 signify transition to the regime of persistent QT prolongation. Thus fractal analysis would provide additional diagnostic values when QT prolongation cannot be easily detected using conventional methods. Alteration of dominant scaling exponents along with multifractal characteristics in repolarization would also enable us to detect the effect of interventions such as hERG blockers by examining multifractal spectra of QT intervals. Acknowledgements The authors would like to thank Stephanie Adams for providing detailed information on experimental procedures, Dr. Hui-Rong Qian for reviewing the manuscript and providing suggestions on presenting the statistical results, and Dr. Ary Goldberger for informative discussions on multifractality in HRV. References Amaral, L. A. N., Ivanov, P. Ch., Aoyagi, N., Hidaka, I., Tomono, S., Goldberger, A. L., et al. (2001). Behavior-independent features of complex heartbeat dynamics. Physical Review Letters, 86, 6026−6029. Beckers, F., Verheyden, B., Ramaekers, D., Swynghedauw, B., & Aubert, A. (2006). Effects of autonomic blockade on non-linear cardiovascular variability indices in rats. Clinical and Experimental Pharmacology and Physiology, 33(5–6), 431−439. Chiang, A. Y., Li, D., & Leishman, D. J. (2007). Quantifying long-range heart rate variability. Journal of Pharmacological and Toxicological Methods, 56, e63. Chiu, K., Chan, H., Chu, S., & Lin, T. (2007). Carvedilol can restore the multifractal properties of heart beat dynamics in patients with advanced congestive heart failure. Autonomic Neuroscience, 132, 76−80. Eichhorn, E. (1998). Restoring function in failing hearts: The effects of beta blockers. American Journal of Medicine,, 104, 163−169. Frazier, S. K., Moser, D. K., & Stone, K. S. (2001). Heart rate variability and hemodynamic alterations in canines with normal cardiac function during exposure to pressure support, continuous positive airway pressure, and a combination of pressure support and continuous positive airway pressure. Biological Research for Nursing, Vol. 2, 167−174. Goldberger, A. L., Amaral, L. A. N., Hausdorff, J. M., Ivanov, P. Ch., Peng, C. -K., & Stanley, H. E. (2002). Fractal dynamics in physiology: Alterations with disease and aging. Proceedings of the National Academy of Sciences of the United States of America, 99(Suppl 1), 2466−2472. Harada, T., Abe, J., Shiotani, M., Hamada, Y., & Horii, I. (2005). Effect of autonomic nervous function on QT interval in dogs. Journal of Toxicological Sciences, 30, 229−237. Ivanov, P. Ch., Amaral, L. A. N., Goldberger, A. L., Havlin, S., Rosenblum, M. B., Struzik, Z., et al. (1999). Multifractality in human heartbeat dynamics. Nature (London), 399, 461−465. Jiang, Z., & Zhou, W. (2007). Scale invariant distribution and multifractality of volatility multipliers in stock markets. Physica A, 381, 343−350. Kantelhardt, J. W., Zschiegner, S. A., Koscielny-Bunde, E., Havlin, S., Bunde, A., & Stanley, H. E. (2002). Multifractal detrended fluctuation analysis of nonstationary time series. Physica A, 316, 87−114. Kazachenko, V. K., Astashev, M. E., & Grinevich, A. A. (2007). Multifractal analysis of K+ channel activity. Biochemistry (Moscow) Supplementary Series AMembrane and Cell Biology, 1, 169−175. Kuwahara, M., Hiraga, A., & Nishimura, T. (1998). Power spectral analysis of heart rate variability in a horse with atrial fibrillation. Journal of Veterinary Medical Science, 60, 111−114. Malik, M., Farrell, T., Cripps, T., & Camm, A. J. (1989). Heart rate variability in relation to prognosis after myocardial infarction: Selection of optimal processing techniques. European Heart Journal, 10, 1060−1074. Matsunaga, T., Harada, T., Mitsui, T., Inokuma, M., Hashimoto, M., Miyauchi, M., et al. (2001). Spectral analysis of circadian rhythms in heart rate variability of dogs. Americal Journal of Veterinary Research, 62, 37−42. Nickerson, M. (1949). The pharmacology of adrenergic blockade. Pharmacological Reviews, 1, 27−101. Peng, C. K., Havlin, S., Stanley, H. E., & Goldberger, A. L. (1995). Quantification of scaling exponents and crossover phenomena in nonstationary heartbeat time series. Chaos, 5, 82−87. Sanbe, A., James, J., Tuzu, V., Nas, S., Martin, L., Gulick, J., et al. (2005). Transgenic rabbit model for human troponin I-based hypertrophic cardiomyopathy. Circulation, 111, 2330−2338. Stein, P. K., Domitrovich, P. P., Hui, N., Rautaharju, P., & Gottdiener, J. (2005). Sometimes higher heart rate variability is not better heart rate variability: Results of graphical and nonlinear analyses. Journal of Cardiovascular Electrophysiology, 16, 954−959. Struzik, Z. R., Hayano, J., Sakata, S., Kwak, S., & Yamamoto, Y. (2004). 1/f scaling in heart rate requires antagonistic autonomic control. Physical Review E, 70(5), 050901.

128

D. Li et al. / Journal of Pharmacological and Toxicological Methods 58 (2008) 118–128

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. Circulation, 93, 1043−1065. Telesca, L., Lapenna, V., & Macchiato, M. (2005). Multifractal fluctuations in seismic interspike series. Physica A, 354, 629−640. Waagstein, F., Hjalmarson, A., Varnauskas, E., & Wallentin, I. (1975). Effect of chronic beta-adrenergic receptor blockade in congestive cardiomyopathy. British Heart Journal, 37, 1022−1036.

West, B., Latka, M., Glaubic-Latka, M., & Latka, D. (2003). Multifractality of cerebral blood flow. Physica A,, 318, 453−460. Willson, K., Francis, D. P., Wensel, R., Coats, A. J. S., & Parker, K. H. (2002). Relationship between detrended fluctuation analysis and spectral analysis of heart rate variability. Physiological Measurement, 23, 385−401.