## Abstract

This study explored the effects of gender and aging on the complexity of cardiac pacemaker activity. Electrocardiogram signals were studied in normal women (*n* = 240) and men (*n* = 240) ranging in age from 40 to 79 yr. Nonlinear analysis of short-term resting R-R intervals was performed using the correlation dimension (CD), approximate entropy (ApEn), and largest Lyapunov exponent (LLE). Evidence of nonlinear structure was obtained by the surrogate data test. CD, ApEn, and LLE were negatively correlated with age. Despite similar means and SDs of the R-R intervals, women had a significantly higher CD, ApEn, and LLE compared with men in the age strata of 40–44 and 45–49 yr. CD and ApEn were strongly (*r* > 0.71) correlated with low- and high-frequency components. We conclude that the resting cardiac pacemaker activity of women is more complex than that of men in middle age, and the gender-related difference diminishes after the age of 50 yr. The higher complexity implies a more comprehensive neural modulation.

- chaos
- discharge pattern
- gender difference
- autonomic nervous system

neuronal activity is generated by a composite deterministic system with some stochastic input, and this system can be either linear or nonlinear (6, 24, 25, 29). On the basis of assumptions of linear models, power spectral analysis has been applied the to biomedical sciences and is particularly suitable for periodic and stationary signals (15, 35, 36). On the other hand, nonlinear models have gradually gained in popularity as a number of experimental signals are difficult to delineate by traditional linear methods (6, 16,19). Of the developed algorithms for nonlinear analysis, the correlation dimension (CD) (8), approximate entropy (ApEn) (21), and largest Lyapunov exponent (LLE) (33) are three frequently used algorithms for the dynamic analysis of biological signals including neuronal discharges (24) and electrocardiograms (ECGs) (17). They are often related to measurements of complexity or “chaos.” Nonlinear analysis of heart rate (HR) variability (HRV) is especially popular in clinical research for providing better detection sensitivity and specificity for some pathological states (22). In addition to clinical applications, nonlinear analysis of HRV may also provide a window on the coding of neural activity in humans by analyzing the ECG signals originating from cardiac pacemaker cells.

Gender and age are two basic factors continuously affecting body functions. These two factors are also correlated with disease categories and even life expectancy (23). Sexual dimorphism in many aspects of neural structures and functions has been noted (3, 9). With regard to cardiac pacemaker activity, there is still debate as to whether HR dynamics differ between women and men. The mean and SD of the R-R interval appear to be similar in both sexes (14). Women have been reported to have lower (28), similar (2), and higher (14) high-frequency power (HF) and similar low-frequency power (LF) of HRV compared with men. Despite differences in experimental setups, we wondered whether the confusing results could be due to the nonlinear characteristics of pacemaker activity (24), which may cause large variations when analyzed by traditional linear methods. Such variations might lead to severe interference if the sexual dimorphism in sinus node discharge is relatively small.

The overall objective of the present study was to determine whether the complexity or chaos of cardiac pacemaker activity differs between women and men. Because aging may be a major determinate of HR dynamics, we systemically studied the effect of aging on nonlinear properties and on gender-related differences. Experiments were performed in healthy humans in the supine rest position. We also compared our results with standard frequency-domain methods to measure parasympathetic and sympathetic regulation of the HR.

## MATERIALS AND METHODS

#### Study sample and experimental setup.

The investigation conforms with the principles outlined in the Declaration of Helsinki. In total, 480 normal volunteers (240 women and 240 men), aged 40–79 yr, were randomly enrolled in this study. The age distribution is similar to those of the individual female (59.5 ± 0.7 yr) and male (59.6 ± 0.7 yr) populations. They were evenly distributed into 8 strata based on 5-yr age intervals, and each age stratum contained 30 women and 30 men. We excluded subjects with the following conditions that affect cardiovascular fluctuations (26): hypotension (systolic pressure <110 mmHg or diastolic pressure <60 mmHg), hypertension (systolic pressure >140 mmHg or diastolic pressure >90 mmHg), diabetic neuropathy, an implanted cardiac pacemaker, frequent occurrence of atrial fibrillation, premature atrial or ventricular contractions, or other forms of arrhythmia. There was no significant difference in mean arterial pressure between the male and female populations. Furthermore, no patients were receiving medication reported to influence cardiovascular fluctuations, such as hypnotics or autonomic blockers. Informed consent was obtained from each participant.

#### Processing of ECG signals.

The analytic algorithms for ECG signals have been detailed in previous investigations (14, 34). Because 5 min has been suggested by the Task Force of the European Society of Cardiology and the North American Society of Pacing and Electrophysiology as the data length for a short-term HRV analysis (26), a precordial ECG was taken for 5 min in the day time (8:00–16:00) while each subject lay quietly and breathed normally in a quiet environment. The raw ECG signals were recorded using an analog-to-digital converter with a sampling rate of 256 Hz. The digitized ECG signals were stored on a hard disk for off-line analysis. Signal acquisition, storage, and processing were performed on IBM personal computer-compatible computers. The computer program was written in Pascal (Borland Pascal 7.0).

The HRV analysis algorithm was designed according to the recommended procedures (26). In the QRS identification procedure, the computer first detected all peaks of the digitized ECG signals using a spike detection algorithm (12) similar to general QRS detection algorithms. Parameters such as the amplitude and duration of all spikes were measured so that their means and SDs could be calculated as standard QRS templates. Each QRS complex was then identified, and each ventricular premature complex or noise was rejected according to its likelihood in standard QRS templates. The R point of each valid QRS complex was defined as the time point of each heartbeat, and the interval between two R points (the R-R interval) was estimated as the interval between current and latter R points. In the R-R interval rejection procedure, a temporary mean and SD of all R-R intervals were first calculated for standard reference. Each R-R interval was then validated as follows: if the standard score of an R-R value exceeded 4, it was considered erroneous or nonstationary and was rejected. The average percentile of R-R rejection according to this procedure was 0.9%.

#### Time-domain and frequency-domain analyses.

The mean (RR) and SD of the validated R-R values were estimated. R-R values were subsequently resampled at the rate of 7.11 Hz to accomplish the continuity in time domain. It should be noted that the resampling procedure was only applied to this section but not applied to the following nonlinear analyses. Frequency-domain analysis was performed using the nonparametric method of fast Fourier transform (FFT). The direct current component was deleted, and a Hamming window was used to attenuate the leakage effect (13). For each time segment (288 s; 2,048 data points), our algorithm estimated the power spectral density based on the FFT. The resulting power spectrum was corrected for attenuation resulting from sampling and the Hamming window (14). The power spectrum was subsequently quantified into various frequency-domain measurements as defined previously (14, 26), including LF (0.04–0.15 Hz), HF (0.15–0.4 Hz), and the LF-to-HF (LF/HF) ratio.

#### Correlation dimension.

We estimated the CD according to the method developed by Grassberger and Procaccia (8, 32). Each reconstructed phase space point serves as the center of a step-by-step growing hypersphere. A sequence of radii (*r*) is defined, and the number of points within the resulting hypersphere [the correlation integral C(*r*)] is determined. C(*r*) versus *r* is plotted as a log-log diagram. Within a certain range, there is a linear relation between log C(*r*) and log *r*. The slope of this relation equals the CD. To apply this technique to HRV research, we made a technical improvement of the CD estimation algorithm to accomplish an automatic analysis. Traditionally, log C(*r*) and/or the slope are plotted versus log *r*. Our computer algorithm additionally plots the linearity as a function of log*r* (Fig. 1). Within the range of interest (between 4% and 8% of the mean value of the time series in this study), the algorithm locates the *r* at which the linearity exhibits the largest value. It then selects the corresponding slope as the CD. This procedure ensured the highest reliability for quantifying the CD. Because this procedure is fully automatic and no manual operation is needed, the experimental data can be analyzed efficiently and objectively. According to Mansier and colleagues (17), the data number (*N*) should be >10^{CD/2}, whereas the embedding dimension (*m*) has been suggested to be ≥3. In the present study, we computed CD with*N* = 300 and *m* = 4.

#### Approximate entropy.

ApEn is a nonlinear measure of variability, expressing the logarithmic likelihood that data points of a certain deviation, represented by*r*, over a defined number of observations, represented by*m*, retain the same distance between incremental comparisons. More regularity is associated with small values of ApEn and greater irregularity with large values. ApEn was calculated using the algorithm of Pincus and co-workers (21, 22). The ApEn calculation has three parameters (*N*, *m*, *r*). In the ApEn calculation, *N* and *m* were fixed (*N* = 300 and *m* = 2). The *r*value for ApEn was calculated as 4% of the mean value of the time series.

#### Lyapunov exponent.

The LLE was estimated via the fixed-time evolution program developed by Wolf and colleagues (32, 33). A *d*-dimensional system has *d* Lyapunov exponents. If it is a nonlinear chaotic system, at least one exponent is positive. The magnitude of the Lyapunov exponent is directly related to the degree of dependence on initial conditions. Wolf's algorithm determines the LLE. Several parameters are involved in this algorithm, on which a successful determination depends. According to Mansier and colleagues (17), *N* should be >10^{CD} and*m* should be >2 × CD + 1. Thus *N*= 300 and *m* = 9 were used to calculate LLE in this study. A reference trajectory and a near-test trajectory in reconstructed space are observed over time, and their mutual distance is determined. If the distance becomes greater than a fixed value (SCALMX), the test trajectory is exchanged for a new trajectory, which lies closer to the reference trajectory. Points that lie closer than a constant (SCALMN) are not suitable as starting points for test trajectories, because the distance between the two trajectories may not exceed the noise level. In the present study, *N* was equal to 300, and the parameters SCALMX and SCALMN were set to 2% and 0.5%, respectively, of the mean value of the time series.

#### Surrogate data test.

To obtain statistical evidence as to whether the HRV presented a nonlinear structure, the estimation of CD and ApEn were repeated with surrogate data (27). We obtained the surrogate data based on randomizing the phases of the Fourier transform of the time series (27, 31). Ten surrogate data sets of each time series were calculated, which all share the same mean, variance, and power spectrum of the original data. The value ρ = (m_{surr} −*x*
_{exp})/SD_{surr} was then estimated, where *x*
_{exp} is the result (the calculated CD or ApEn) of the applied algorithm on the experimental data and m_{surr} and SD_{surr} are the mean and SD of the results on the 10 surrogate data series, respectively. A value of ρ > 2 allows one to reject the null hypothesis (that the result is due to a stochastic linear process) and is indicative of nonlinear dynamics in the data (11, 27). This computation was performed for CD and ApEn but, because of the limitations described in the original paper of surrogate data method (27), not for LLE.

#### Statistical methods.

LF, HF, and LF/HF of the linear analysis were logarithmically transformed to correct for the skewness of the distribution (14). Correlations between parameters were assessed using Pearson's correlation coefficient. We considered a good or strong correlation between two variables to be indicated by*r*
^{2} > 0.50 (or *r* > 0.71), because over 50% of the change in one variable can be explained by a change in the other. Differential effects of the two sexes and the eight age strata on HRV parameters were compared using two-way ANOVA. When indicated by a significant *F*-statistic, regional differences were isolated using post hoc comparisons with Fisher's least-significant-difference test. Comparisons between the two sets of data were performed with an unpaired Student's *t*-test. Statistical significance was assumed for *P* < 0.05. Values are expressed as means ± SE.

## RESULTS

Surrogate data test revealed that the mean ρ values of CD and ApEn were 7.04 and 5.02, respectively. Because the mean ρ values are >2, the null hypothesis of linearly correlated noise can be rejected. Women in their 40s (40–49 yr old) had a higher CD of R-R intervals in the reconstructed phase space (Fig. 1
*A*) than did the men in their 40s (Fig. 1
*B*). As age increased, the CD of both sexes decreased toward a similar level (Fig. 1, *C*and *D*).

Figure 2 provides a typical example demonstrating how gender and aging influence ApEn. It is obvious that R-R values of women in their 40s (Fig. 2
*A*) are more scattered than those of the men in their 40s (Fig. 2
*B*) in the third dimension (RR_{n}
_{+2}) of the reconstructed phase space under the condition that the two-dimensional scattergrams (RR_{n}, RR_{n}
_{+1}) were similar in both sexes. Sexual dimorphism eventually diminished in the elderly (Fig. 2,*C* and *D*), whose R-R values significantly converged in the third dimension.

Although *m* used for the analysis of LLE was eventually chosen to be 9, we specially offered a reconstructed picture with an*m* of 2 for ease in demonstrating the effects of gender and aging on LLE (Fig. 3). The trajectories of women in their 40s diverged more widely than did those of men at the same ages (Fig. 3, *A* and *B*), indicating that women had a larger LLE. However, such a sex-related difference was hard to identify in older age groups (Fig. 3, *C* and*D*), where the trajectories diverged less in both sexes.

#### Differential effects of gender and aging on nonlinear properties of HRV.

Before we describe the quantitative results of the nonlinear HRV analyses, it should be mentioned that there were no significant differences in means or SDs of R-R intervals between male and female populations. We found that our female population had a higher CD and ApEn (Table 1). Linear regression analysis showed that the three nonlinear measurements of HRV were statistically (*P* < 0.001) but weakly (*r* ≈ −0.4) correlated with age (Table2). When analyzed with the two genders and eight age strata, two-way ANOVA detected significant effects of age on SD, CD, ApEn, and LLE and significant effects of gender on CD, ApEn, and LLE (Fig. 4). The interaction between age and gender was significant for ApEn and LLE but not for CD. When we used SD, LF, HF, and LF/HF as covariates to depurate CD, ApEn, and LLE from the linear effects of SD and these spectral components, the effects of both age and gender on these nonlinear indexes were still significant. For women, significant changes from the age stratum of 40 yr were not detected until the age strata of 55 yr for SD and 50 yr for CD, ApEn, and LLE. For men, changes were not detected until the age strata of 60 yr for SD, CD, ApEn, and LLE. When the effect of gender was analyzed at each age stratum, it was noted that women exhibited a greater CD, ApEn, and LLE in the age strata of 40–44 and 45–49 yr. All differences disappeared in age strata ≥ 50 yr.

#### Relationship between nonlinear and linear measurements of HRV.

Table 2 compares the nonlinear and linear measurements of HRV. The three nonlinear measurements were only weakly correlated with RR (*r* ≤ 0.22). CD and ApEn but not LLE were strongly (*r* ≥ 0.72) correlated with SD, LF, and HF. Of all the linear indexes considered, HF was best correlated with the nonlinear indexes. None of the nonlinear indexes, however, was strongly correlated with LF/HF (*r* ≤ 0.20). ApEn was well correlated with individual CD and LLE (*r* ≥ 0.82; Table3).

## DISCUSSION

Our short-term HRV analyses reveal that woman's HR is characterized by a higher CD, ApEn, and LLE compared with that of man in the middle age, indicating a more complex signal broadcasting from woman's cardiac pacemaker. In the aspect of neurophysiology, a more complex neural signal is always accompanied with a more involved neural network. On the contrary, a very simple firing pattern can only be observed in an isolated neuron (6, 36). Thus the higher complexity of women's HR dynamics implies that the female heart is modulated more comprehensively by the autonomic nervous system (ANS), especially the rapid vagal influence (1), although such modulation is not strong enough to produce an evident change in the mean and SD of R-R intervals. With regard to the effect of age, previous studies have revealed that aging is accompanied with a decrease of complexity in either cardiac pacemaker activity (7) or midbrain neural activity (6). Our older subjects, either women or men, had a lower CD, ApEn, and LLE, indicating a lower degree of neural modulation to the cardiac pacemaker. The lower complexity may be due to a decrease in the ANS potency resulting from the normal aging process. The higher complexity of HR dynamics observed in women before the age of menopause is especially noteworthy because it may be related to the lower cardiac mortality and longer life expectancy of women. Such changes in physiological functions can hardly be detected and studied by classical physiological techniques such as calculations of the mean and SD of R-R intervals.

With regard to neuronal discharge, it has been noted that the spike pattern may not be satisfactorily delineated by linear analysis in some circumstances. The importance of nonlinear theories and related techniques were therefore highlighted (6, 16). On the other hand, the ECG signal is one of the strongest electrical signals on the surface of the human body. It ultimately reflects the electric activity of pacemaker cells in the sinus node. The study of ECGs is therefore not only of clinical interest but also important for exploring neural coding problems, especially in the ANS. As with other aspects of the nervous system, cardiac pacemaker activity has been proposed to consist of a combination of linear and nonlinear properties. The former properties have been broadly studied (14,26), but the latter ones were hard to appreciate until the recent development of these analytic techniques.

CD, ApEn, and LLE are three nonlinear indexes that have been applied to the analysis of cardiac pacemaker activity. Grassberger and Procaccia (8) developed an algorithm for calculating CD, which is similar or even equal to the fractal dimension. The concept of entropy, as it applies to signals like R-R intervals, is to quantify the repetition of patterns in that signal (18). Larger values of entropy correspond to greater apparent randomness or irregularity, whereas smaller values correspond to more instances of recognizable patterns in the data. Pincus and co-workers (21) published an algorithm of ApEn for evaluating entropy that has been applied to ECG signals. This index has been shown to be sensitive to some pathological conditions such as sudden infant death syndrome (20) and autonomic disturbances (18). Lyapunov exponents provide information in regard to the sensitivity of a system to initial conditions. The LLE can be estimated via the fixed-time evolution program developed by Wolf and colleagues (33). The magnitude of the Lyapunov exponent is directly related to the degree of dependence on initial conditions.

When we analyze the relationship between each of the three nonlinear indexes, a good correlation was noted between CD and ApEn (*r* = 0.82). Such an observation may be due to the close mathematical origins of CD and ApEn (8, 21). However, ApEn is simpler to calculate mathematically and requires less calculation time. LLE originates from a different model and is computed by a special computer algorithm (33). Because of the complexity of the algorithm as well as the high *m* required (17), calculating LLE eventually took the most time in this study. The relationship between LLE and ApEn was also good (*r* = 0.83). All three of these nonlinear indexes can clearly detect the effect of gender and aging on cardiac pacemaker activity (Fig. 4).

Previous studies regarding the nonlinear analysis of HR dynamics chiefly focused on its capability to discriminate pathological states (20, 22) or senile changes (7) from normal conditions. Because illness and aging may lead to a prominent change in ANS function, which in turn leads to a prominent change in HR dynamics, such alterations can usually be detected by traditional linear methods. Thus the theoretical advantage of nonlinear over linear indexes cannot be highlighted. With regard to the effect of gender, studies using linear analysis techniques have led to controversial results. For example, women have been reported to have a lower (28), similar (2), and higher (10, 14) HF to that of men. In our observations, neither the time-domain indexes nor some of the frequency-domain indexes (10, 14) could detect the effect of gender. In contrast, all of these nonlinear indexes could detect that effect. Such an observation indicates that these nonlinear techniques may be proper tools to study the ambiguous gender-related differences in HR dynamics.

Linear regression analysis revealed strong correlations between HF and CD and between HF and ApEn. Vagal control provides the fastest beat-to-beat influence on R-R intervals from the central nervous system (1). At the present time, HF has been well accepted as a noninvasive index for cardiac vagal function (26). Thus it is very likely that the values of CD and ApEn are closely related with the vagal function (4). Such an observation is also compatible with that of a previous animal study (4).

#### Perspectives.

Women have a longer life expectancy than men in almost every country of the world. It has been proposed that gender-related differences in cardiovascular function are responsible. The exact mechanism underlying these effects is unclear at the present time. According to our data, however, it can be deduced that changes in the nonlinear properties may reflect effects of gender and aging on cardiac parasympathetic regulation. The gender-related difference in autonomic regulation function is especially noteworthy. Women before menopause have a lower risk of heart diseases than do men. The “cardioprotective” effect of estrogen has been proposed for this. But the linkage between the sex hormone and the heart is still ambiguous. Because a potentiated vagal function may protect the heart from tachyarrhythmias after ischemic heart disease and decrease the mortality rate (5, 30), the dominance of vagal function in middle-aged women may produce a protective effect against lethal tachyarrhythmias. In our study, 71% women were postmenopausal; the mean age of menopause was 47.7 yr, and only 6% of postmenopausal women received hormone replacement therapy. Thus the loss of gender-related differences in the nonlinear indexes after 50 yr of age supports the hypothesis that estrogen has a facilitating effect on cardiac parasympathetic regulation as revealed by the nonlinear content of HR dynamics.

Nonlinear analysis techniques may reveal new aspects of HR dynamics. Mechanisms underlying the effects of gender and aging on nonlinear properties of cardiac pacemaker activity and related clinical applications warrant further investigations. Because the algorithms may be directly incorporated into current HRV analysis protocols, CD, ApEn, and LLE have a high potential for studying gender-related differences in cardiac pacemaker activity and related autonomic regulation.

## Acknowledgments

This study was supported by Veteran General Hospital-Taipei Grant VGH90-427-2 and National Science Council (Taiwan) Grant NSC-90-2320-B-320-010.

## Footnotes

Address for reprint requests and other correspondence: T. B. J. Kuo, Institute of Neuroscience, Tzu Chi Univ., 701, Sect. 3, Chung Yan Rd., Hualien 970, Taiwan (E-mail:tbjkuo{at}mail.tcu.edu.tw).

The costs of publication of this article were defrayed in part by the payment of page charges. The article must therefore be hereby marked “

*advertisement*” in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.June 21, 2002;10.1152/ajpheart.00169.2002

- Copyright © 2002 the American Physiological Society