INTRODUCTION

Top Abstract Introduction Appendix Appendix Appendix Appendix References

THIS PAPER DESCRIBES and summarizes a study of infant breathing by using data-analysis techniques derived from dynamic systems theory (DST). Such techniques have been useful for examining other complex physiological rhythms such as heart rate (7, 46, 47, 63, 67), electroencephalogram (2, 24, 26, 29, 36, 37), parathyroid hormone secretion (31), and opticokinetic nystagmus (45) and can distinguish variations that are random from those that are deterministic ("chaotic").

Any model that describes the control of respiration must account for two important behaviors, rhythm and pattern generation. There are two distinct approaches to the problem of rhythm: the first requires the existence of discrete "pacemaker" cells with intrinsic activity that drives other respiratory neurons. The output of various respiratory centers or pools of motoneurons is then organized by a pattern generator. A second approach implies that networks of cells with oscillatory behavior interact in a complex way to produce respiratory rhythms, which are either further organized by a pattern generator or might be self-organizing (11). The first model is, in effect, a linearly coupled model with pattern generation dependent on feedback mechanisms from the various components of the respiratory system, whereas the second is more likely to be complex. Advances in neurobiology have allowed recordings to be made from individual neurons and groups of neurons in the brain. Using these techniques, various studies have demonstrated that the concept of discrete respiratory centers made up of neurons with specific functions defined by the nature of a particular "center" is obsolete (11). Whereas there is organization of neurons into functional networks or pools, these are not necessarily anatomically discrete. Also, there are conflicting data in regard to the presence of a specific pattern generator. Given the complexity of the connections between the various groups of oscillating, respiratory-related neurons and the capacity for interactions between simple oscillating systems to produce complex behavior, we have argued that information about the organization of respiratory control can be determined by using DST. In essence, the argument that there is a simple pattern generator that coordinates the output from various respiratory centers is unnecessary if the output from interacting networks is dynamic and self-organizing. To examine such a hypothesis, respiration must first be adequately shown to be chaotic, then examined further to describe its dynamic structure.

Most studies of the dynamic behavior of biological systems have used fractal dimension estimation to establish that a system's behavior is chaotic or to classify distinct types of behavior by their complexity. Recent studies have suggested that respiration in humans is chaotic. If that is the case, then techniques derived from DST should allow the dynamic structure of respiratory behavior to be better described, thus improving our understanding of the control of breathing. However, these earlier studies have important limitations. Most studies have used the Grassberger and Procaccia algorithm (14, 15) for estimating fractal dimension, which is simple and easy to implement. Unfortunately, it is now recognized (19, 33) that this algorithm has some technical problems that can lead to misinterpretations of data. The most serious problems occur with small data sets or when the system incorporates a substantial noise component. The studies described below each suffer these limitations to a greater or lesser extent and none examined respiratory patterns in quiet sleep when breathing is most regular.

Donaldson (9) used DST to examine respiratory behavior in humans. He studied resting adults and, by estimating Lyapunov exponents, concluded that resting respiration is chaotic. However, this relatively simple approach could not distinguish a nonlinear dynamic system from linearly filtered noise.

Pilgram et al. (30) studied the dimensional complexity in the respiratory traces of infants during rapid-eye-movement (REM) sleep. They concluded that breathing during REM sleep is chaotic but did not study patterns of breathing during quiet sleep. In reaching their conclusion, Pilgram et al. also employed the method of surrogate data. Dimension-estimation algorithm used by Pilgram et al. failed to give reliable estimates for their noise-driven surrogates, and so the authors inferred that the control of respiration is possibly deterministic.

In another set of studies, Sammon and colleagues (38-42) examined dynamic respiratory behavior in rats. From their observations they concluded that, in anesthetized vagotomized rats, the respiratory system behaves as an oscillator with a single degree of freedom. With the vagus intact, however, respiratory behavior was more complex, exhibiting low-order chaos, which, the authors speculated, was due to feedback from various types of pulmonary afferent activity.

In this paper, we test for cyclic modulation in the amplitude of respiration. Cyclic fluctuation in the amplitude of respiration indicates the presence of additional complexity in the respiratory system and is consistent with the results of our correlation dimension calculations. Other authors (see, e.g., Ref. 64) have identified cyclic behavior in breath sizes. Waggener and colleagues (64) identify cyclic variation in the breath sizes of resting adults at simulated altitude.

The study reported here employs a new algorithm to determine fractal dimension that is technically more complex but is in practice more reliable, more robust under the restrictions of finite data, and less prone to misinterpretation. Estimates of fractal dimension were used as an aid to identify the dynamic system that produced the data we have measured.

From dimension estimations we conclude that the dynamics of breathing during quiet sleep are consistent with a large-scale, low-dimensional system with a substantial small-scale, high-dimensional component, i.e., a periodic orbit with a few (perhaps 2 or 3) degrees of freedom supplemented by smaller, more complex deterministic fluctuations.

The nature of the low-dimensional system was investigated further by constructing surrogate data, which enabled us to test whether the dynamics were consistent with linearly filtered noise or with a nonlinear dynamic system. Our results show clearly that, in almost all cases, the dynamics are best described as a low-dimensional, nonlinear dynamic system being driven by a high-dimensional noise source. In those cases when such a model is inconsistent with the data, the measured data have strong indications of nonstationarity; that is, the breathing patterns changed during the recording (for example, a sudden switch to periodic breathing occurred).

To describe a complex dynamic system as chaotic may not be helpful when the components that produce the structure cannot be described. Understanding these components is more important than simply classifying the system as chaotic. However, before the mechanism can be understood, it is important to determine basic properties of the system: Is it dynamic? Is it nonlinear? Is it chaotic?

We, therefore, investigated the structure of the dynamics we described in greater detail. We used a new form of linear modeling to expose a possible source of our multidimensional periodic orbit. By applying a reduced autoregressive modeling (RARM) algorithm to some data from sleeping infants we observe cyclic amplitude modulation (CAM) indicative of substantial structure in breath-to-breath variability. These calculations lead us to conclude that respiration in quiet sleep is most probably controlled by multidimensional oscillator supplemented by small-scale, complex (chaotic) behavior. Moreover, our results are inconsistent with a noise-driven, one-dimensional model of eupnea (or any other model that does not include substantial structure in breath-to-breath variability).

After a brief introduction to the new dimension estimation algorithm, we describe the experimental methodology, including a description of our surrogate data-generation methods and our RARM algorithm. Finally, we discuss the dimension calculations and the results of the hypothesis testing by using the surrogate data sets.

	METHODS

Using standard noninvasive inductive-plethysmography techniques, we obtained a measurement proportional to the cross-sectional area of the chest or abdomen, which is a gauge of the lung volume. The present study collected measurements of the cross-sectional area of the abdomen of infants during natural sleep. The study was approved by the Princess Margaret Hospital Ethics Committee.

Subjects

Data Collection

The only practical limitation on the length of time for which data could be collected is the period during which the infant remains asleep and still. The cross-sectional area of the lung varies with the position of the infant. However, in this study, we were interested only in the variation due to the breathing, and so we have been careful to avoid artifacts due to changes in position or band slippage. We have made observations of up to 2 h that are free from significant movement artifacts, although, typically, observations are in the range of 5-30 min. All 27 observations used to calculate dimension are between 240 and 360 s; those used to identify CAM are between 400 and 1,400 s.

In contrast to the study by Pilgram and colleagues (30) that examined breathing in REM sleep, we have studied infants in quiet sleep. From measurements of electroencephalogram, electromyogram, and electrooculogram, sleep stage was determined by using standard polysomnographic criteria (6). During quiet sleep, breathing often appears relatively regular. The possibly chaotic features of most interest are the small variations from this regular periodic behavior. Because we wished to observe such fine detail, we did not filter signals. The analog output of the respiratory plethysmograph (operating in its direct-current mode) has no built-in filtering. Filtering methods, such as linear filters and singular-value decomposition methods (30), can remove some features that we wished to observe. Furthermore, filtering (even to avoid aliasing) has been shown in some cases to lead to erroneous identification of chaos (25, 28, 34).

The abdominal signal is not necessarily proportional to lung volume. Takens' embedding theorem (56) (and, therefore, the methods of this paper; see APPENDIX A, Embedding) only requires a diffeomorphism (a smooth function) of a measurement of the system. Moreover, present methods are not capable of dealing well with multichannel data, and, therefore, use of both rib and abdominal signals to approximate actual lung volume is difficult. Of the available measurements, we found that the abdominal cross section was the easiest to measure experimentally.

The 27 observations used to calculate dimension were selected based on sleep state (quiet stage 3-4 sleep) and then on the basis of sufficient stationarity and a minimum of 4 min in length. From each of these observations, 240 s of stationary data (the 240 s that had the most stationary moving average) were used to calculate dimension. Data for RARM calculations were selected based on being at least 400 s in length and taken during uninterrupted quiet sleep.

	DATA ANALYSIS

In this study, we employed three main analysis methods: correlation dimension estimation, surrogate data analysis, and RARM. This section provides a description of these methods. The mathematical detail is deferred to APPENDIXES A-D. First, we discuss generalized dimension and correlation dimension estimation, and then we provide an overview of surrogate-data techniques and describe RARM procedures.

Glossary

C_N()

Correlation function

d_c, d_c(₀)

Correlation dimension

d_e, n

Embedding dimension

I(X )

An indicator function

MDL

Minimum description length

N

Length of a data set

P_j, P_j(·)

An orthogonal projection

Z⁺

Set of positive integers

Length scale

₀

Observation scale

Embedding lag

·

Euclidean normthe usual distance function

Dimension Estimation

Most applications of correlation dimension to physiological sciences have utilized the Grassberger and Procaccia algorithm (14, 15). However, this paper employs a new algorithm, which is technically more complex but is in practice more reliable and less prone to misinterpretation. Unlike previous estimation methods, this new algorithm recognizes that the dimension of an object (its structural complexity) may vary depending on how closely you examine it. Hence, the value of the estimate of correlation dimension may change with scale. It therefore offers a more informative and accurate estimate of dimension.

Computing correlation dimension d_c as a function of scale d_c(₀) can tell us much more about the structure of an object; for example, it can indicate the presence of large scale "periodic" motion and simultaneously detect smaller scale, higher dimensional chaotic motion. Quoting a single number as the correlation dimension of a data set ignores much of this information; in many respects, it produces an "average dimension." Plots of dimension as a function of scale are particularly important when complex physiological behavior is studied because they yield far more information than a single estimate at a fixed scale.

The estimation algorithm. The estimation algorithm used for the calculations in this paper is described in detail by Judd (19, 20), and an alternative treatment may be found in (for example) Ref. 17. The accuracy of this algorithm has been independently verified by Galka and colleagues (13). One important advantage of the new method is that it is possible to calculate error bars for dimension estimates. The confidence intervals on the dimension that the algorithm provides are dependent on the length of the time series.

Phenomenological Models

Surrogate-data techniques involve (implicitly or explicitly) building a model from data, as do RARM procedures. It is important to emphasize that the types of models we are describing are phenomenological rather than mechanistic. That is, these models are not described as a set of equations that have been motivated by physical laws or physiological insight. Rather, the models we describe here are fitted to experimental data to predict future behavior of the system or to simulate the system on a computer and produce an output indistinguishable from observations of the actual system.

Usually, a phenomenological model is selected from a broad abstract class of models by fitting data using some statistical criteria, that is, to find the most probable model within the class. Nonlinear models and reduced autoregressive models are particular examples of broad abstract classes. Mechanistic models include, for example, a model of the lungs as bellows or a physical model of blood-gas concentrations. In a mechanistic model, parameters and variables represent physical quantities, for example, the elasticity of the lung or blood oxygen concentration. In a phenomenological model, this is not the case. However, phenomenological models can provide clear and useful information about a system, such as the periodicities of the system.

Surrogate Data

To determine whether our results indicate the presence of anything more complicated than a noisy linear system, we employed the surrogate-data methods described by Theiler et al. (59). The principle of surrogate data is the following. One first assumes that the data come from some specific class of dynamic system, possibly fitting a parametric (phenomenological) model to the data. One then generates surrogate data from this hypothetical system and calculates various statistics of the surrogate data and original data. The surrogate data give the expected distribution of statistical values, and one can check that the original has a typical value. If the original data have atypical statistics, then we reject the hypothesis that the system that generated the original data is of the assumed class. One always begins with simple assumptions and progresses to more sophisticated models if the data are inconsistent with the surrogate data.

We expect our data to be most consistent with some type of nonlinear dynamic system. Before considering this type of surrogate, it is necessary to determine whether a simpler description of the data would be sufficient. To do this, we compared our data with surrogates generated by the traditional (linear) methods (see APPENDIX B, Traditional Surrogate Data). Many studies in biological sciences have employed these traditional surrogate methods (2, 31, 42, 63, 67). These methods determine whether experimental data are significantly different from specific (broad) categories of linear systems. In addition to these linear surrogate tests, we applied a new, more complicated, nonlinear surrogate test (50, 53). This method was used to determine whether the data are distinguishable from these generated by a broad class of nonlinear models (see APPENDIX B, Noise-Driven Nonlinear System Surrogates, and APPENDIX C).

Reduced Autoregressive Models

Our new method extends the traditional autoregressive (AR) model of order n [AR(n)], which predicts the next value in a time series as a weighted average of the last n values. We consider instead a RARM where any past values may be used to predict the upcoming value but only those that are important are used. To determine which past values are important, we employed Rissanen's minimum description length (MDL) criterion (35) using a modeling procedure described by Judd and Mees (21, 22); our implementation of these methods has been presented elsewhere (55) and is discussed elsewhere (48, 54). A computer software implementation of these techniques is available online (49) or from the first author. A brief description of the mathematical methods is presented in APPENDIX D. For now, let us assume that RARM can produce a model consisting only of those previous values that are useful in predicting future values. This is not necessarily a particularly good model in terms of predictionit is only an approximation to the breath-to-breath dynamics that we utilize to extract important information. Using this information, we deduce the period of quasi-periodic behavior in the time series from the temporal separation of the previous values. Hence, RARM can identify the period of CAM in much the same way as autocorrelation may, except our methods prove to be more sensitive and more discriminatory. Our method involves automated application of statistical procedures to time series and can be applied even when amplitude modulation is not apparent to the eye.

By reducing the original data to a breath-to-breath time series, we effectively rescale the time axes so that each breath is of equal length. However, CAM on a breath-to-breath basis does not (necessarily) suppress time-dependent dynamics. A hypothesized cyclic variation in breath duration may be related to a cyclic variation in breath amplitude (and, hence, evident in the breath amplitude time series). These could essentially be two separate observations of the same periodic behavior. The duration of a single breath is also far more difficult to measure accurately because of (relatively) long flat peaks and troughs.

	RESULTS

We first present our results from applying our dimension estimation algorithm. Next, we describe the results of our surrogate-data and RARM calculations.

Dimension Estimation

View larger version (27K): [in this window] [in a new window]   Fig. 1.   Correlation dimension estimates for 1 representative data set from each of the 10 subjects. Any data sets that produced dimension estimates dissimilar to those illustrated here are discussed in text (see Dimension Estimation). Plots are of scale (log 0) against correlation dimension, with confidence intervals shown as dotted lines (often indistinguishable from the estimate). Correlation dimension estimates were produced for embedding dimensions of 2, 3, 4, 5, 7, and 9 for all data sets, except for subjects 2, 4, and 7. Subjects 4 and 7 failed to produce an estimate for the 9 dimensional embedding. Subject 2 did not produce an estimate when embedded in 3 or 9 dimensions. All other dimension estimates are illustrated; higher embedding dimension produces larger correlation dimension.

The more or less linear increase in dimension with decreasing scale ₀, and the shift to higher dimensions as the embedding dimension is increased, are both indications that the system, or measurements, have a substantial component of high-dimensional dynamics, or noise, at small-to-moderate scales. The increase of dimension with decreasing scale is an obvious effect of high-dimensional dynamics or noise. The shifting to higher dimensions with increasing embedding dimension occurs because in a higher dimensional embedding the points "move away" from their neighbors and tend to become equidistant from each other, which in effect amplifies, or propagates, the small-scale, high-dimensional properties to large scales. (This effect is related to the counterintuitive fact that spheres in higher dimensions have most of their volume close to their surfaces rather than near their centers, as is the case in two and three dimensions.)

Some of the dimension estimates, particularly in two and three dimensions, produced curves that linearly increased for large-length scales but appeared to level off as length scale decreased. For most of the estimates we have computed, this is the case when the data are embedded in two dimensions. Furthermore, for these embeddings in two-dimensional space, the correlation dimension estimate seemed to approach two. This indicates that as we look closer at the data (that is, at a smaller length scale), the data appear to fill up all of our embedding space. For many of the dimension estimates (Fig. 1, subjects 7 and 9), the embedding in three dimensions also leveled at values slightly less than three. This behavior can be attributed to an attractor with correlation dimension of ~2.8-2.9. However, it is probably more likely that this too is simply due to the data "filling up" the three-dimensional space. This is consistent with the results of our false nearest neighbor calculations, which suggested that three- or four-dimensional space would be required to successfully embed the data.

There is one particular estimate that appeared to behave quite differently from all the others. Some of the curves of the estimates for subject 8 appeared to increase, decrease, and then increase again. This could indicate that, as we look more closely at the structure, there is some length scale for which the embedding structure seems to be relatively high in dimension, whereas by looking at an even smaller length scale the behavior has significantly lower dimension. These observations are supported by what we can observe directly from the data. This time series includes an episode of periodic breathingincreasing the complexity of the large-scale behavior (see Fig. 2). Similarly, some of the data sets for subject 2 include large sighs, causing the dimension estimate to increase at large scales (see Fig. 3).

View larger version (31K): [in this window] [in a new window]   Fig. 2.   Dimension estimate for subject 8, one of data sets used in our analysis. Periodic breathing caused the dimension estimates (dimension estimates used embedding dimensions of 2, 3, 4, 5, 7, and 9) at large scale to increase.

View larger version (26K): [in this window] [in a new window]   Fig. 3.   Dimension estimate for subject 2, one of our data sets along with the dimension estimates (shown are estimate with an embedding dimension of 2, 3, and 4). Note large sighs during recording and the corresponding increase in dimension estimate at moderate scale. Another data set from same infant exhibited similar behavior and produced a similar dimension estimate.

Finally, the remainder of the estimates (e.g., Fig. 1, subjects 1, 2, 4, 6, 7, 8, and 10) behaved in yet another manner. These estimates were approximately constant for a small range of large length scales and gradually increased over small length scales. The estimates at large length scales were generally about two to three, indicating that the large-scale behavior is slightly above two dimensional. The increase in dimension estimate for smaller length scales can again be attributed to either noise or high-dimensional dynamics. However, the scale of "small-scale structure" in the dimension estimates is at a larger scale than the instrumentation noise level. Typically, the smallest scale is log_e(₀)  2.5, a scale of ~5% of the attractor (e³  0.049787  0.05). The digitized signal will typically use at least 10 bits of the analog-to-digital converter (2¹⁰ = 1/1,024 < 0.001); other sources of instrumental error are certainly at levels <5%.

The approximately two-dimensional behavior is probably due to the inspiration-expiration cycle along with the breath-to-breath variation within that cycle. This is easily visualized as the orbit of a point around the surface of a torus. A dimension estimate of two could indicate that the attractor was any two-dimensional surface; the embedded data, however, have an approximately toroidal shape. In this motion, there are two characteristic cycles, first, the motion around the center of the torus and, second, a twisting motion around the surface. Our estimates slightly over two indicate that this behavior is complicated further by some other roughness over the surface of the torus. The shape of such an attractor would very closely resemble the textured surface of a doughnut.

Surrogate-Data Analysis

Linear surrogates. By comparing the value of dimension obtained from our data and surrogates consistent with each of these three null hypotheses, we were able to reject all three null hypotheses (see Fig. 4 for an example of such a calculation). These results are summarized in Table 1, which shows the number of SDs between the values of d_c(₀) for data and surrogate, for the value of log (₀), which gave the greatest difference. This is calculated over the range 2.5  log (₀)  0.5, and for d_e = 3, 4, 5. Data are from infants at 2 mo of age. The symbol n/a indicates that none of the surrogates produced convergent dimension estimate at any value of ₀. For each data set and each hypothesis test, there are three pairs of numbers. These three pairs of numbers are the results for d_e = 3, 4, 5, respectively. The first number is the number of SDs by which the mean value of dimension for the surrogates exceeded that for the data. The second number (in parentheses) is the value of log (₀) for which this occurred.

View larger version (30K): [in this window] [in a new window]   Fig. 4.   Linear surrogate calculations, an example of surrogate-data calculations for algorithms 0, 1, and 2. Here, we compared correlation dimension estimate for 1 of our data sets (solid line) and 30 surrogates (dotted lines). There is a clear difference between the correlation dimension of the data and that of the surrogates.

View larger version (61K): [in this window] [in a new window]   Fig. 5.   Surrogate data. Sections of 3 surrogates generated by traditional techniques: algorithms 0, 1, and 2, and a section of a surrogate data set generated from a radial-basis model. Also shown is a section of real data (top) used to generate these surrogates. There are obvious similarities between true data and the nonlinear surrogate, whereas other surrogates are obviously different.

Nonlinear surrogates. For each set of data, we have calculated its correlation dimension. Using a modeling algorithm described elsewhere (21, 22, 51), we constructed a radial-basis model of the data. From this model, we constructed 30 surrogate data sets. The surrogates were embedded in two, three, and four dimensions, by using the same embedding strategy as the true data set. We then calculated the correlation dimension curve for each set of surrogate data.

View larger version (38K): [in this window] [in a new window]   Fig. 6.   Nonlinear surrogate dimension estimates. Surrogate-data calculations for 2 of our data sets, embedded in 2, 3, and 4 dimensions. First set indicated a close agreement between data and surrogate. Second set of calculations indicated very clear distinction. Hence, the model of first data set is indistinguishable (according to correlation dimension) from a noise-driven periodic orbit, whereas model of the second fails to produce particularly strong similarities. Notice that, for almost any value of 0, comparison of the value of dc(0) of the data and the surrogates would also lead to these conclusions.

CAM Detected by Using RARM

Table 2 outlines the results of our calculations applied to 14 data sets from 14 infants. Data used for these calculations were recorded during quiet sleep. Subjects 1-10 are the same subjects as used for correlation dimension estimates. Data for subjects 1-6 were recorded during the same study as those used for dimension estimates. Data for subjects 7-12 were recorded at 4 mo of age, for subjects 13-14 at 6 mo. Respiratory rate is the average respiratory rate over the duration of the recording. Note that, although there was some variability in both the respiratory rate and the period expressed as a number of breaths, the period in seconds was relatively constant. In most cases, this period also falls within the range of periodic breathing. In subject 11, periodic breathing with cycle times from 13.5 to 15.5 s occurred during the same study.

	DISCUSSION

This study has confirmed that apparently regular breathing during quiet sleep is probably chaotic. This conclusion should be qualified. As Rapp (33) observed, to conclude that a phenomenon is chaotic is both difficult and often irrelevant. In real data sets, noise contamination will always increase the dimensional complexity of the data, and almost any experimental data will exhibit noninteger correlation dimension. Identification of apparently chaotic behavior is, however, a good first step in dynamic analysis. Our data are the most convincing evidence that respiratory variability in humans is deterministic and not random, due to noise. Furthermore, we have extended our observations and analyses to describe the dynamic structure of the system in greater detail. We have demonstrated the presence of a two-dimensional periodic system: one dimension corresponds to the regular breath-by-breath cycle of respiration, the other to a regular periodic modulation of breath amplitude (CAM), with cycle length approximating that of periodic breathing.

By applying a new dimension estimation algorithm, we have produced reliable, accurate, and informative dimension estimates. This method produces error bars which, for moderately long time series (such as ours), are reassuringly small. This, however, must not be taken to indicate that we were able to estimate the dimension of the dynamic system responsible for this behavior with similar confidence. As with any experiment, the accuracy of the results is limited by physical constraints such as observational error and noise interference. Even if the data can be collected in such a way that they accurately represented the volume of the subject's lung, there may be other phenomena contributing to that measurement, which are independent of the dynamic system we are attempting to observe. Hence, the tightness of the error bars we have calculated was due more to the quantity of the data than to lack of experimental error.

Our dimension estimate results indicate that, on a large scale, there is low-dimensional behavior, whereas the small-scale behavior is often dominated by very high-dimensional dynamics or noise (that is, extremely high-dimensional dynamics). Even though false nearest neighbor techniques suggest that we were embedding in high enough dimensions, there was still some small-scale behavior that filled the embedding space. The scale at which the embedding space is filled by the dynamics could indicate the level of experimental noise.

The most conclusive estimates from this study indicated that the structure of the attractor is likely to be similar to a torus with small-scale, very high-dimensional dynamics. Hence, at large length scales, the structure looked like the surface of a torus, whereas at smaller length scales dimension increased. This indicates that the attractor appears to be a torus with a very rough surface. The most important conclusion from these data is that this two-dimensional periodic system indicates two levels of periodicity. Hence, in addition to the periodic inspiration-expiration motion, it is likely that there was some cyclic breath-to-breath variation in breath depth.

By applying the method of surrogate data we demonstrated that the correlation dimension is related to the data from which we estimate it in a nontrivial way. The surrogates produced by algorithms 0, 1, and 2 are clearly inadequate. It is apparent that they should fail, and this was confirmed by our results. We have constructed our own surrogates using a nonlinear modeling process and compared surrogates and data to test the accuracy of the model. For the majority of our data sets, we found that the data and nonlinear surrogates were indistinguishable according to correlation dimension. For those data sets that were distinguishable from their surrogates, we found that there were several possible reasons for this. Usually, if the data were nonstationary, the model simply failed to produce surrogates that were close enough to the data. The model is stationary and periodic, whereas the data are not.

Occasionally, with nonstationary data, the model failed to produce even periodic surrogates. If this was the case, then the model had a stable fixed point. In these cases, the dimension estimates of data and surrogate were obviously different, and a better model is required. The fact that this modeling algorithm failed in cases where the data were not stationary is not particularly surprisingboth modeling and dimension estimation algorithms require stationarity. Perhaps, with improved modeling techniques, similar results could be obtained in these cases.

Even if both data and surrogate were stationary, the dimension estimates of the surrogates could still be different from those of the data. In all these cases, however, this has been found to be a problem with the level of dynamic noise introduced to the model to generate the surrogates. By changing the noise level, the dimension would also change, effectively moving the dimension estimate vertically. Because, in these cases, the shapes of the dimension estimate curves were approximately the same, by altering the noise level it was possible to produce surrogate estimates that were indistinguishable from the data. In all cases, however, the dynamic noise was substantially less than the model's root mean square prediction error. The root mean square prediction error is the noise level predicted by the modeling algorithm. However, this is the total noise and includes both dynamic and observational noise.

Dynamic noise and observational noise have a different effect on the correlation dimension estimates. Observational noise will increase the value of correlation dimension at length scales less than and equal to the noise level. It appears from our calculations that an increase in the level of dynamic noise increases the correlation dimension estimate equally across all length scales, effectively producing a vertical shift in the estimate. An increase in dynamic noise will certainly have a greater affect on d_c(₀) for larger ₀ than would a similar increase in observational noise. Assuming one has correctly identified the underlying deterministic dynamics, it may be possible to "tune" the level of dynamic noise so that surrogates and data have approximately the same value of correlation dimension estimates at moderate-to-large length scales and then alter the level of observational noise to tune the dimension estimates at small length scales. Hence, the level of noise required to be unable to reject the surrogate data is an indication of the relative proportion of dynamic and observational noise in the system. That is, we can distinguish between random behavior within the system (dynamic noise) and experimental error (observational noise).

Our calculations have demonstrated that the respiration of infants in quiet sleep can be modeled well with radial-basis methods. We conclude that the behavior is similar to a noise-driven periodic orbit. Furthermore, the radial-basis surrogates confirm that the periodic orbit associated with this behavior has at least two degrees of freedom, as suggested by our initial dimension estimates. Finally, our RARM calculations to detect CAM support these results.

The observation of CAM is intriguing. We serendipitously recorded periodic breathing from one infant. The cycle time of CAM (13.3-15.6 s) in the same infant corresponded almost exactly with that of the observed periodic breathing (13.5-15.5 s) as demonstrated in Fig. 2. The relationship to periodic breathing needs further investigation, but we believe that these two behaviors with identical cycle lengths (CAM and periodic breathing) are likely to be related and determined by similar factors, whatever they might be. These data support the hypothesis that oscillatory activity responsible for periodic breathing is ubiquitously present but masked during apparently regular breathing by the normal regular stimulation from respiratory motoneurons. Periodic breathing occurs when this normal regular drive is decreased [i.e., in infants when core body temperature is raised (18)]. The adoption of one particular physiological state, i.e., regular tonic respiration with CAM or periodic breathing is likely to be dependent on the environmental conditions and maturity of respiratory control as well as the presence of any pathological conditions.

In conclusion, this study has addressed the limitations of previous studies that have examined whether respiration is chaotic. We investigated children in quiet sleep when breathing appears most regular. Correlation dimension estimates are consistent with a chaotic system. Furthermore, unlike in most previous studies, we used surrogate data analyses to test whether the apparently chaotic behavior was due to linearly filtered noise. We found this unlikely and concluded that the simplest system consistent with our data is a noise-driven, nonlinear, radial-basis model. Our data are the most convincing evidence that respiratory variability in infants is deterministic and not random, due to noise. A recent study has demonstrated reduced variability of respiratory movements in infants, who subsequently died of sudden infant death syndrome (SIDS) [Schechtman et al. (43)]. This observation was retrospective but suggests that, because the variability that we have observed during quiet breathing is deterministic, then further study using DST could allow early identification of infants at risk of SIDS from simple measurements of respiratory patterns.

We speculate that since CAM is an important contributor to the complexity observed during quiet breathing, further studies might demonstrate distinct patterns of CAM in infants with respiratory control problems; for example, absence of CAM might explain the reduction in variability observed by Schechtman et al. (43) in infants who died of SIDS.

We believe that the techniques employed by us in this study could be used in the future for studying the complex dynamic behavior of respiration and its control in health and disease.