## Abstract

We examine the potential to treat unstable ventilatory control (seen in periodic breathing, Cheyne-Stokes respiration, and central sleep apnea) with carefully controlled dynamic administration of supplementary CO_{2}, aiming to reduce ventilatory oscillations with minimum increment in mean CO_{2}. We used a standard mathematical model to explore the consequences of phasic CO_{2} administration, with different timing and dosing algorithms. We found an optimal time window within the ventilation cycle (covering ∼1/6 of the cycle) during which CO_{2} delivery reduces ventilatory fluctuations by >95%. Outside that time, therapy is dramatically less effective: indeed, for more than two-thirds of the cycle, therapy increases ventilatory fluctuations >30%. Efficiency of stabilizing ventilation improved when the algorithm gave a graded increase in CO_{2} dose (by controlling its duration or concentration) for more severe periodic breathing. Combining gradations of duration and concentration further increased efficiency of therapy by 22%. The (undesirable) increment in mean end-tidal CO_{2} caused was 300 times smaller with dynamic therapy than with static therapy, to achieve the same degree of ventilatory stabilization (0.0005 vs. 0.1710 kPa). The increase in average ventilation was also much smaller with dynamic than static therapy (0.005 vs. 2.015 l/min). We conclude that, if administered dynamically, dramatically smaller quantities of CO_{2} could be used to reduce periodic breathing, with minimal adverse effects. Algorithms adjusting both duration and concentration in real time would achieve this most efficiently. If developed clinically as a therapy for periodic breathing, this would minimize excess acidosis, hyperventilation, and sympathetic overactivation, compared with static treatment.

- heart failure
- Cheyne-Stokes respiration
- inhaled carbon dioxide

the presence of oscillations of respiratory and cardiac parameters with an approximate period of 1 min in patients with congestive heart failure (CHF) indicates adverse prognosis (6, 8, 16, 23, 25, 30). This periodic breathing (PB) is primarily associated with oscillations of arterial CO_{2} (7, 10), where the troughs in CO_{2} engender apneas and/or hypopneas (40). These respiratory oscillations are accompanied by cardiovascular fluctuations (9, 11, 12), which may lead to adverse effects in patients with CHF (15, 17, 18).

Although PB in CHF patients is associated with increased mortality (8, 16), there is no clinical evidence suggesting that abolition or reduction of PB would improve this outcome in CHF patients. However, the potential improvement of quality of life alone, from improved sleep, may make the development and investigation of novel therapies worthwhile. Modeling studies can facilitate such developments and inform subsequent clinical investigations to test their relevance.

The use of inhaled CO_{2} to stabilize PB has been investigated since the early 1980s (5). Several other papers since then have reported that the administration of constant concentrations of inhaled CO_{2} elevates the troughs in CO_{2} oscillations (2, 24, 32, 39), preventing apneas, provided that sufficient CO_{2} is administered. However, the systemic consequences of the resulting hyperventilation and sympathetic overactivation have made this a less attractive treatment option (1, 24). While this form of static CO_{2} delivery reduces the oscillations in ventilatory parameters, the associated arousals are not reduced, with no overall improvement in sleep quality (33). This may be due to the excess CO_{2} directly stimulating the cortex or reducing the threshold for cortical arousal.

A fixed increase in inspired CO_{2} [inspired CO_{2} fraction (Fi_{CO2})] mandates an increase in ventilation and, therefore, increases the work of breathing. This is often commented on in previous studies (1, 2, 5, 24, 32, 33, 35, 38, 39), but not always quantified. However, it can be easily calculated because, for constant net exhaled volume of CO_{2}, the narrowed gap between expired and inspired CO_{2} must be compensated for by a reciprocal increase in alveoli ventilation (*Eq. 1*). (1) where ΔV̇ is the change in alveolar ventilation; V̇co_{2} is a constant volume of body production of CO_{2}; and Δ(Fet_{CO2}-Fi_{CO2}) is the change in difference of end-tidal CO_{2} fraction (Fet_{CO2}) to Fi_{CO2}. Using this calculation, which provides only a lower limit on the increment in ventilation, we see that ventilation is raised by up to 96% (24–96%). We have summarized calculated values of increase in ventilation due to inspired CO_{2} in previous studies (Supplementary Table S2; the online version of this article contains supplementary data). In fact, alveolar ventilation would have to have increased by even more than this, because the increase in ventilatory work is likely to increase V̇co_{2}. Lorenzi et al. concluded: “Although CO_{2} inhalation abolishes CSR [Cheyne-Stokes respiration]-CSA [central sleep apnea] and reduces the frequency of arousals, it also increases ventilation and would therefore augment the energy and blood flow demands of the respiratory muscles in the face of low cardiac output. Therefore, it is unlikely that CO_{2} inhalation would be a useful long-term therapy for CSR-CSA in patients with CHF” (24).

In this study, we use mathematical modeling to investigate treating PB using optimally graded doses of inspired CO_{2} through dynamic delivery.

We aim to deliver targeted CO_{2} to counteract the troughs in end-tidal CO_{2} oscillations, thereby limiting the total quantity of CO_{2} administered. This might minimize the oscillations in the respiratory gases without producing the observed large increases in end-tidal CO_{2} and ventilation. We hypothesized that a small dose of inspired CO_{2} would reduce PB, if administered dynamically at an optimal phase and amplitude of ventilatory oscillation. This modeling approach allows a study of the consequences of different algorithms on the resultant control patterns of reproducibly unstable ventilatory control system.

## METHODS

### Regulation of Ventilation

Ventilation can be considered to be regulated through a feedback control loop of two physiological mechanisms: the chemoreflex and pulmonary gas exchange, which produce the controller and plant gain, respectively. The controller gain is the change in ventilation due to a change in end-tidal CO_{2}, and plant gain is a change in end-tidal CO_{2} as a result of a change in ventilation (20).

### Potential Steady State of the Feedback Control System

The dependence of Fet_{CO2} on steady-state alveolar ventilation (plant gain) is an inversely proportional relationship, because the product of their means of ventilation and Fet_{CO2} must equal the mean metabolic production of CO_{2} (provided that the total body production of CO_{2} remains constant). This relationship (*Eq. 2*) produces a hyperbolic curve (the isometabolic curve) when CO_{2} is plotted against ventilation. (2)

The dependence of ventilation on end-tidal CO_{2} (controller gain) is determined by the chemoreflex response curve. This is shown as a linear relationship (*Eq. 3*) in Fig. 1, implying a constant value of chemoreflex gain, but can also be nonlinear (26). (3)

The intercept of these two response curves is where the ventilation control system is potentially in equilibrium (Fig. 1).

The stability of this feedback control system depends on how it responds to small perturbations from the equilibrium state. If, after a small perturbation away from the equilibrium point, the system responses cause it to diverge further away from the equilibrium point, the system is unstable. In contrast, if the system returns to equilibrium, it is stable (26).

### The Dynamic Treatment Model

We modified a previously published iterative model (26). In its original form, the model maps system behavior following a small, transient perturbation in CO_{2} from a steady-state condition and hence predicts whether the cardiorespiratory control system will be stable or unstable for a wide range of both linear and nonlinear chemoreflex response curves.

The development of the model described here allows the introduction of inspired CO_{2} at any desired time and dose during the simulation and allows the effects of this on the ventilation pattern to be quantified. Hence we were able to identify those treatment strategies that are likely to be of clinical benefit.

### The Modified Model

The iteration of the model begins with the instantaneous CO_{2} value at the steady-state position at time *t*_{o} followed by an introduction of a small, transient, positive perturbation on the CO_{2} (Fig. 2). The instantaneous ventilation value V̇ at time *t* (i.e., *t*_{o} + Δ*t*, where Δ*t* is 1 s in the model) is then calculated from a previous value of CO_{2} (C_{t-δ}) using a linear function of the chemoreflex gain. The new value of CO_{2} for time *t* is calculated using a difference equation below (*Eq. 4*), derived from a previously published analytic equation (13). (4) where ΔC is the incremental change in end-tidal CO_{2} concentration during the small time increment Δ*t*, Vl is a constant lung volume, V̇ and C are instantaneous ventilation and end-tidal CO_{2} variables, respectively, at time *t*, C̄ is the mean alveolar CO_{2} fraction, Q̇ is a constant cardiac output, and β is a constant of solubility of CO_{2} in blood.

The new value of V̇ is then tested according to a configured treatment regime for whether the model should administer CO_{2} therapy or not. If ventilation is unstable and the value of V̇ (both in amplitude and phase of oscillation) is in the configured treatment zone (as described in *The Treatment Regime* section below), then a corresponding programmed CO_{2} treatment dose is added to the new value of the end-tidal CO_{2}.

In the treatment model, the instantaneous CO_{2} at time *t* (C_{t}) is given by, (5) Here C_{(t-Δt)} is the instantaneous CO_{2} at the previous time point (*t*-Δ*t*). ΔC is the component of change in CO_{2} that is caused by normal gas exchange processes, but assuming zero inspired CO_{2}, and is obtained from *Eq. 4*. CO_{2} therapy is the component of change in CO_{2} that is caused by inspired CO_{2} and is calculated from the configured values of “peak concentration” and treatment phase window “duration”. This value is dependent on the amplitude and phase of the ventilation cycle at time *t*, by an algorithm discussed in more detail below, in *The Treatment Regime* section.

The new end-tidal CO_{2} value is then used to calculate the new value of ventilation.

In the model, the following cardiorespiratory parameters were held constant to standard values for all of our simulations: cardiac output, 3.5 l/min; metabolic production of CO_{2}, 0.2 l/min; lung volume, 5 liters; linear chemoreflex response curve with slope of 1,200 l·min^{−1}·atm^{−1}; chemoreflex delay, 20 s; solubility of CO_{2} in blood, 5 l·l^{−1}·atm^{−1}; start point of CO_{2} equal to steady-state CO_{2}, 4%; start point of ventilation equal to steady-state ventilation, 5 l/min; and positive CO_{2} perturbation, 0.001% 10 s after start. It has been shown that these parameters can be assumed constant near the equilibrium point, since they are not critical in the genesis of ventilatory instability (13). The constant values are based on previously published figures that would be typical for a patient with heart failure and PB (4, 14, 22, 26, 28, 36, 37).

The CO_{2} therapy is configured to automatically adjust to the amplitude of oscillation and the phase of ventilation during a treatment episode. This allowed us to dynamically configure the therapy (in concentration, duration, and delivery phase) to test the most effective treatment profile.

### Severity of Ventilatory Oscillations

In this study, the severity of ventilatory oscillations was measured by calculating the relative standard deviation of ventilation over a 15-min simulation period following a 5-min initial transient interval. The smaller the relative standard deviation of ventilation, the more stable the breathing.

### The Treatment Regime

Visually, a cycle of ventilation can be represented as a clock face, with peak ventilation (hyperpnea) at 12 o'clock and trough ventilation (hypopnea) at 6 o'clock (Fig. 3). Increasing amplitude of oscillation is represented by an increase in distance away from the center of the clock. The treatment envelope is represented as a sector of the clock, where a longer duration of therapy is a proportionally larger sector, and where the concentration of therapy can be represented by the degree of shading: the darker the treatment area, the greater the concentration of CO_{2} delivered (Fig. 4).

The angular distribution of the concentration follows a profile of [(1 − cosine)/2] curve, where it rises from zero at the beginning of the treatment episode, peaking in the midangle of the total duration, and then falls back to zero at the end of the treatment episode (Fig. 4). We reasoned that this profile would be best suited to counteract the sinusoidal-like shapes of the troughs of the CO_{2} oscillations.

During the dynamic therapy, the treatment is potentially available during the whole 20 min of the simulation. However, a dose of inspired CO_{2} is delivered during a cycle of ventilation, only if the amplitude of the oscillatory ventilation is >0.6 l/min during the configured treatment phase window. (This allows situations of only trivial oscillation to be left untreated, automatically). In a single treatment episode of a ventilation cycle, the actual dose of the delivered inspired CO_{2} is dependent on the amplitude of ventilatory oscillation during that cycle and the length of the treatment duration in the cycle.

For oscillations between 0.6 and 6.0 l/min, the dose of treatment increases from zero quadratically up to the programmed peak concentration and duration. The quadratic increase in dose is achieved by simultaneously increasing both the concentration and the treatment duration (per cycle) linearly, in response to increasing amplitude of ventilation. PB with oscillation >6 l/min is treated with the programmed maximum peak concentration and duration, as shown in Fig. 4. A three-dimensional representation of the treatment envelope is shown in Fig. 5. The configuration of the model allows the amount of CO_{2} treatment to vary from cycle to cycle, according to the pattern of ventilatory oscillations that develops.

In general, when oscillations are more prominent, the algorithm delivers CO_{2}, either at a higher peak concentration, or for a longer proportion of the PB cycle, or both. The relationship between oscillation amplitude and configured peak concentration and duration can be programmed at the start of each simulation and can be altered at any stage of the treatment simulation, if dose adjustment is required.

### Use of Real-Time Fourier Analysis to Detect Phase Within the PB Cycle

PB sometimes has a sinusoidal envelope of ventilation, but more often has a nonsinusoidal shape, with the lower part flat, representing apnea.

In our model, we monitor the ratio of instantaneous to mean ventilation. The model uses a Fourier transform method to provide a sinusoidal fit to a wide range of oscillatory ventilation waveforms. The software determines the phase in the cycle of PB. This is achieved by performing a sliding Fourier transform on the preceding ventilation data, acquired during an approximately 1-min period, that can be adjusted to match the observed PB cycle length. The sinusoidal fit to the ventilation is then derived by constructing a sine wave with the amplitude and phase of the lowest frequency oscillatory component in the Fourier transform (i.e., that which is closest in cycle time to the period of ventilation signal transformed). This fitting technique has the advantage that it is not adversely affected by oscillatory cycles containing periods of apnea. It will also be unaffected by breath-to-breath fluctuations of the respiratory signal at frequencies much higher than that of PB.

Figure 6 shows an example of the sinusoidal fitting process satisfactorily applied to simulated noisy PB with periods of apnea. The phase of PB is detected validly.

## RESULTS

This model allowed the effect of dynamically timed therapy to be simulated, and the effect of changes in the therapy protocol to be studied. With no treatment programmed, when the model was run, it developed a pattern of PB, as shown in Fig. 7. The period of oscillations was ∼60 s.

We then implemented one of a range of treatment regimes and repeated the simulation, from the same starting point as the no-treatment simulation. The resulting ventilatory pattern was documented and compared with the pattern resulting from no-treatment simulation.

The treatment regimes all contained a zone of therapy (duration), which occurred at a certain phase in the PB cycle, and had a certain peak concentration of administered CO_{2}. The regimes differed in their *1*) timing of center point (midangle) of the treatment episode within the PB cycle (expressed as a “treatment phase”, with 0° representing the time of peak ventilation); *2*) duration of treatment episode within the PB cycle (expressed as “duration” an angle within a circle, with, for example, 90° representing treatment for one-fourth of the cycle); and *3*) peak concentration of administered CO_{2} [expressed as “peak concentration” a fraction between 0 and 1 (100% CO_{2})].

During a treatment episode (i.e., in the therapy zone), the concentration rises from zero to the peak concentration and back to zero again, following a smooth sigmoid distribution of a cosine curve, with peak at the treatment phase, which is the center point of the duration.

### Study 1: Existence of an Optimum Phase Within PB Cycle for CO_{2} Therapy

We used the model to vary the phase of ventilation at which the CO_{2} was delivered, while keeping peak concentration and treatment duration constant. We set the duration of therapy to be one-half the PB cycle (180° out of the 360° cycle), and peak CO_{2} concentration to be 1%. We repeatedly ran the model from the same baseline starting conditions (Fig. 7) with treatment phases varying in steps of 10° from 180° before peak ventilation to 180° after peak ventilation. To assess efficacy, we measured the ratio of the size of ventilatory oscillations (standard deviation of ventilation) with and without treatment. Values < 1.0 indicate reduced amplitude of PB.

We found that the CO_{2} therapy is most effective when delivered at a specific range of treatment phases, in this case between −30° and +40° around peak ventilation (0°), as shown in Fig. 8. The severity of PB is sharply reduced (by >95%) when the treatment phase approaches the optimum range from either end. Outside this phase range, the CO_{2} therapy actually worsens the PB (increased by >30%).

This particular range of optimum treatment phase is dependent on other parameters, and, if these other parameters are changed, the range of optimum treatment phase also changes.

### Study 2: Effect of CO_{2} Dose on Ventilatory Stability

We investigated the effect, on the ventilatory pattern, of administering peak concentrations of CO_{2}, ranging from 0 (no treatment) to 10%. The duration of CO_{2} therapy was kept the same for each simulation at 180°, as in *study 1*. The treatment phase was also kept at peak ventilation for each simulation; this was shown in *study 1* to be within the optimum range of phases for reducing PB for the given baseline conditions.

Initially, we titrated the peak concentration value of administered CO_{2} from zero in small increments to assess when sufficient CO_{2} had been delivered to affect the ventilatory pattern. When the concentration of administered CO_{2} was <0.01%, there was little or no effect on the PB. Small increments of the peak concentration from 0.01% started to affect the ventilatory pattern. Thereafter, peak concentrations between 0.1 and 0.2% reduced the PB sharply. Beyond this point, further increases in concentration resulted in only minimal reductions in PB (Fig. 9).

For any concentration of 0.2% or above, there was a time window during which treatment could be administered and yield a >90% reduction in PB. The most effective abolition of PB was observed with the highest peak CO_{2} dose (oscillations reduced by >98%). However, this high concentration has disadvantages: high concentration treatments reduce the margin of error for the optimum treatment phase (−20° to +30° for 10% CO_{2} vs. −30° to +40° for 2% CO_{2}) and would be associated with greater systemic side effects from the CO_{2}. When higher peak concentration treatment (10% CO_{2}) is delivered outside the optimum phase window, PB worsens by >300%, whereas at lower peak concentration treatment (2% CO_{2}), PB worsens by <33%. We, therefore, used 2% CO_{2} for subsequent simulations, as it provided the best balance of attenuation of ventilatory oscillations with the widest margin of error. Figure 9 demonstrates how the magnitude of residual PB varies with increasing concentration of administered CO_{2} therapy per ventilation cycle (treatment phase = 0°, duration = 180°).

### Study 3: Effect of Duration of CO_{2} Therapy on Ventilatory Pattern

We studied the effect of duration of treatment on the ventilatory pattern. The simulation was run for various treatment durations between 0° and 360° (where supplemental CO_{2} treatment continues throughout the PB cycle, following the cosine curve). For each simulation, the peak concentration and treatment phase of administered CO_{2} were kept the same (2% and 0°, respectively). These were shown to be within the optimal range for reducing PB.

We titrated the duration of CO_{2} therapy from 0° in small increments to assess when sufficient CO_{2} had been delivered to affect the ventilatory pattern. When the duration of CO_{2} delivery was <0.5°, there was little or no effect on the PB. Small increments of duration between 0.5° and 12° reduced the PB sharply. Further increase of the duration beyond this point had little further effect (Fig. 10).

For any duration of treatment of 12° or above, there was a time window during which treatment could be administered and yield >90% reduction in PB. The most effective abolition of PB was observed with the longest treatment duration of 360° (oscillations reduced by >98%). However, this long duration has disadvantages: it allows a smaller margin of error for the optimum treatment phase (20° less for duration of 360° compared with 180°). Furthermore, treatment with 360° duration worsens PB by ∼50%, if delivered outside the optimum treatment phase, compared with <33% increase in degree of PB when applying at 180° duration. In view of these factors, we used a treatment duration of 180° for subsequent simulations, which we felt to be optimal.

### Study 4: Interchangeability of CO_{2} Therapy Concentration and Duration

We found that a decrease in the concentration of CO_{2} administered could be compensated for by an increase in the duration of therapy. In *studies 2* and *3*, we showed that the effect of peak concentration on PB has a pattern similar to that of the effect of treatment duration. Therefore, it should be possible to use several combinations of peak concentration and duration to deliver a similar dose of CO_{2} therapy, to produce an equivalent effect.

We have presented a three-dimensional representation of efficacy of treatment in the online supplement (Supplementary Fig. S2). The relationship between the peak concentration and duration can be expressed as: (6)

Linear increase or decrease of peak concentration and duration simultaneously will change the resultant dose quadratically. This has its benefits in treating PB with moderate-to-severe amplitude oscillations, as we show below.

### Study 5: The Effect of Grading Peak Concentration and Duration Simultaneously When Treating Moderate to Severe (0.6–6.0 l/min) Amplitude Oscillation

Therapy is not administered for low-amplitude PB (<0.6 l/min), while, for high-amplitude PB (severe) (>6.0 l/min), it is more effective to deliver uniform duration and concentration to reduce the amplitude of PB. However, for moderate-to-severe amplitude (between 0.6 and 6.0 l/min) PB, the optimal dose of CO_{2} is no longer uniform but is quadratically related to the amplitude of the PB. The quadratic dose relation is achieved by simultaneously grading the duration and concentration delivery linearly (Fig. 11). This increased the efficacy of therapy by 22%, compared with linearly graded dosing. Grading either parameter linearly in isolation produced a narrower spectrum of optimum treatment range at these moderate-to-severe amplitudes, as shown in Fig. 11 diagrams.

### Study 6: Comparison of the Effectiveness of Static and Dynamic CO_{2} Delivery on PB

We programmed the model to deliver a static CO_{2} therapy at all times to the baseline configuration (Fig. 7) and hence established the constant CO_{2} concentration that optimally reduced PB (Fig. 12).

We performed the simulation with small increments of CO_{2} concentration of 0.05% from 0 to 1% and in increments of 0.5% from 1 to 10%. We found that the optimum concentration for the given configuration with static therapy is ∼2%. This therapeutic concentration reduced the PB by >99%. However, this static treatment increased the average end-tidal CO_{2} from 4 to 4.3%, and the average ventilation from 5 to 8.4 l/min, over the 20-min treatment period. Therefore, abolishing ventilatory oscillations using static CO_{2} delivery is associated with significant increase in average ventilation (Fig. 12).

We studied the effect of using a dynamic treatment at a peak concentration of 2% (the same value as for this static treatment) on the average end-tidal CO_{2} and average ventilation. We also found the static concentration that gave the same reduction in PB (97%), as this dynamic treatment, and calculated the effect on the average end-tidal CO_{2} and average ventilation. The increases in mean end-tidal CO_{2} of the dynamic and static treatments were 0.0005 and 0.1710 kPa, respectively. The increases in mean ventilation of the dynamic and static treatments were 0.005 and 2.015 l/min, respectively (for the same 97% reduction in PB). The total CO_{2} delivered during the static treatments were >80% more than the amount delivered during dynamic treatment (Table 1).

#### Slow steady increase in inspired CO_{2}.

We modeled the effect of a slow steady increase in Fi_{CO2} on the ventilatory pattern, for comparison with our dynamic CO_{2} delivery method, by gradually increasing inspired CO_{2} from 0 to 5% over 2 h of simulated time.

We started with the baseline configuration, but, after the initial period, during which the PB had become fully developed and persistent, we introduced a series of 0.25% increments in inspired CO_{2} at 5-min intervals, from 0.00 to 5.00%.

We found that inspired concentrations >0.75% were large enough to affect the ventilatory pattern. By the time inspired CO_{2} had risen to 1.75%, the residual PB oscillations were <1/100 of the baseline size. However, this entailed a substantial increase in ventilation. At an optimum inspired CO_{2} value of 1.75%, average ventilation had already increased by >50%.

### Study 7: The Effect of Other Parameters on the Optimum Treatment Phase of Ventilation

We studied the effect of other parameters on the optimal treatment window: cardiac output, chemoreflex gain, lung volume, metabolic production of CO_{2}, and chemoreflex delay. Instead of the fixed values used in *studies 1*–*6*, in this study, we took values for these parameters lower and higher (cardiac output 2.5–4.5 l/min, chemoreflex gain 1,100–1,800 l·min^{−1}·atm^{−1}, lung volume 3–6 liters, metabolic production of CO_{2} 100–400 ml/min and chemoreflex delay 17–40 s) (13).

For reasons of brevity, we have not shown here the results of every combination of these values, or of other values, but the findings can be summarized by observing that whenever the system configuration is unstable, there turns out to be an optimum window of therapy. The size and precise timing of the optimal window varies with system configuration, but is generally close to peak ventilation. For the extremes of the cardiac output range quoted above, the optimal delivery phase varied by 10° (0° at 2.5 l/min to 10° at 4.5 l/min). A similar variation was obtained for the extremes of chemoreflex gain, metabolic production of CO_{2}, and lung volume.

Chemoreflex delay has a more dramatic effect, as might be expected. As chemoreflex delay ranges from 17 to 40 s, the period of oscillation increases proportionately (by a factor of ∼2.5), and the optimal treatment phase moves from −10 to +80° referenced to peak ventilation. This change in optimal treatment phase corresponds to a change in time of a magnitude similar to the increase in chemoreflex delay. Therefore, in a clinical device, it may be possible to predict the optimum treatment phase region by measuring the period of the ventilatory oscillation and inferring from this the approximate value of the chemoreflex delay.

Further sensitivity analysis of the model is described in the online supplement.

### Study 8: Peripheral and Central Chemoreceptors

In the previous sections, we have concentrated on a single chemoreceptor response for clarity in introducing the principles of the dynamic administration of CO_{2}. To investigate the response of the algorithm on multiple chemoreceptors, we reconfigured the model with two notional groups of chemoreceptors with individual group chemoreflex gain and delay. The first set of chemoreceptors was programmed with a gain of 180 l·min^{−1}·atm^{−1} and had one of three short time delays (15, 20, or 25 s). The second set of chemoreceptors was programmed with a gain of 1,020 l·min^{−1}·atm^{−1} and had one of three long time delays (80, 100, or 120 s). This corresponds to a split into a ratio of 15 to 85% of the gains used in our simulations using a single chemoreflex gain. This ratio is typical of published values (34). We ran the simulation with various combinations of the individual set chemoreflex delay times. We incorporated a low-pass filter in line with the second set of receptors to account for the sluggish response of central chemoreceptors.

In all cases, we found that, whenever the system is unstable, there turns out to be a window of optimum CO_{2} administration to effect a reduction in PB. Similar findings were obtained with more equal distribution and when three receptor groups were simulated.

### Study 9: Effect on System Stability of Large-amplitude, Long-period Initial Driving Oscillations

To establish whether the results above are peculiar to the particular transient, small-amplitude, initial perturbation that we tried, we repeated the simulations with large-amplitude, long-period initial driving oscillations. The driving oscillations lasted for 5 min, had cycle time of 1 min, and were of a size sufficient to generate apnea at the nadir of ventilation. The treatment was then applied as above, and, after a further 10 min, the resulting size of oscillation was calculated over a 15-min period. An example is shown in Fig. 13*A*.

The results for optimum treatment phase are shown in Fig. 13*B*. In general, they are the same as those found with single, transient, small-amplitude, initial perturbations. Figure 13*B* also shows that the treatment delivered during the apneic phase (close to 180°) did not affect the PB.

## DISCUSSION

We have devised and studied a novel mathematical model to deliver optimally timed dynamic CO_{2} therapy to ameliorate PB. Using our model, we were able to control the timing, concentration, and duration of CO_{2} delivered, as a function of the phase and amplitude of ventilatory oscillations. We have consistently shown that there is an optimum phase window in which a small dose of CO_{2} therapy can improve PB. If therapy is delivered outside the optimum range of treatment phase, PB worsens. Importantly, the model allows the prediction of this margin of error.

We considered several factors to determine the optimum treatment values, namely, the maximum PB reduction, with minimum dose delivered and maximum margin of error to remain in the optimum treatment window. For example, in Figs. 9 and 10, a treatment with concentration of 0.8% and duration of 40° may appear as effective as 2% concentration and 180° duration. However, these graphs show the response of the system only when treatment is delivered at the ideal time. In practice, it is useful to have an algorithm that works, even if treatment is not delivered at the ideal time, giving a wider window of opportunity to treat effectively. The window of optimum treatment phase is only one-half the size (−10 to +20°), with the former configuration than with the latter (−30 to +30°). The latter configuration is, therefore, preferable. Furthermore, we have identified that the ideal pattern for this algorithm to work would be with both the concentration and duration (length of treatment episode per ventilation cycle) increasing linearly (quadratic increase of dose), with increasing degree of oscillation.

Most relevant to potential clinical therapies, this study showed that the dynamic treatment increased neither the average end-tidal CO_{2}, nor the average ventilation, from the baseline value. This is in contrast to the traditional static treatment, where our model shows these parameters to be increased by physiologically significant amounts. Our investigation, in *study 6*, demonstrates that the dynamic treatment requires a much smaller dose of CO_{2}, compared with static treatment that reduces PB by the same level.

Previous clinical studies have used static CO_{2} to treat PB (1, 5, 24, 32, 33, 35, 39), based on the rationale that patients with PB have arterial CO_{2} levels below their apneic threshold. With static supplementary CO_{2} administration, the apneas are eliminated, as corroborated by our *study 6*. However, the mean ventilation increases, as demonstrated in Supplemental Table S2 of the online supplement. In contrast, our dynamic therapy aims to target only the troughs in end-tidal CO_{2} that drive hypoventilation. The targeted inspired CO_{2} will potentially counteract a transient reduction in arterial Pco_{2}, preventing it from falling below the apneic threshold, without grossly increasing its mean value. Our proposed system aims not only to prevent apneas, but also to attenuate hypopneas. By minimizing end-tidal CO_{2} oscillations, ventilatory oscillations may potentially be reduced or abolished without significantly increasing mean ventilation.

There are some patients in whom mean ventilation is low: in these, our concern to minimize the increment in mean ventilation is not justified, and so the proposed therapeutic model need not be considered particularly advantageous.

We believe that the reason why PB worsens when treatment timing is outside the optimum region is that the inspired CO_{2} contributes to augmenting hypercapnia instead of counteracting the transient reduction of CO_{2}.

One potential advantage of algorithms such as the one we present here is that it has a relatively short “memory”. Sleep is known to be characterized by discontinuous changes in system characteristics (20, 21, 22), which could be seriously problematic for a therapy strategy that relied on static therapy, or therapy aimed at a particular CO_{2} or ventilation target. Our dynamic cycle-detecting approach, however, automatically adjusts to new system characteristics as they arise. Moreover, if system characteristics become stable, therapy automatically ceases.

Previous static CO_{2} administration regimes have improved breathing in CSR, but not sleep quality, as reported by Steens et al. (32), Andreas et al. (1) and Szollosi et al. (33). The latter two groups suggested that the main cause for this lack of improvement in sleep quality is elevated sympathetic activity, due to hypercapnia. If this is so, the dynamic treatment suggested by our model may potentially enable reduction of PB without such arousals.

Thomas et al. (35) showed that CO_{2} concentrations as low as 0.5% were an effective adjunct to continuous positive airway pressure in the treatment of severe complex sleep apneas, improving sleep quality as well as respiratory disturbances. Although our model corroborates the effectiveness of low-concentration CO_{2} therapy, it is not possible, from their work, to determine which proportion of the synergistic effects of CO_{2} and continuous positive airway pressure resulted in these beneficial effects. Our baseline simulation differs from this study in that we used levels of cardiac output compatible with those found in CHF patients (as in the clinical studies cited above), whereas they treated patients without heart failure.

Although our simulations showed that the approach of gradual increase in inspired CO_{2} abolished PB at low inspired CO_{2} level, it also showed that this increased the average ventilation hugely. This demonstrates that we still face similar problems of elevated ventilation with gradual CO_{2} increase as we would with static CO_{2} administration.

Added dead space has been shown to stabilize PB (19). However, this is achieved at the expense of increased ventilation similar to static CO_{2} treatment.

Two previous groups have produced ingenious and simple devices, which arrange to deliver added inspired CO_{2} during the patient's hyperventilation phase to increase end-tidal CO_{2} (3, 31). With these systems, therapy always centers around a phase near peak ventilation with a fixed time delay due to (circuit and physiological) dead space. As our model suggests, the period of PB can vary such that significant changes in the phase of optimal CO_{2} delivery away from peak ventilation are required for varying chemoreflex delays. Such variation cannot be achieved with these systems, but is obtainable with a model-mediated approach.

In both cited methods, once the fresh gas flow is set in its application to control breathing, the circuit could interpret increases in baseline ventilation as hyperventilation and would deliver CO_{2}. Our proposed method would treat only periodic fluctuations in ventilations and not changes in the baseline value.

### Study Limitations

This study is a mathematical simulation, and not a clinical study. Although this can be considered a limitation, it has the advantage of our being able to guarantee a constant underlying control system instability. This ensures that differences between one run of the simulation and another can be attributed to differences in the therapy algorithm and are not due to random environmental fluctuations. Clinical studies paralleling this one in terms of variety of algorithms tested would be difficult, because biological variability would necessitate large numbers of replicate experiments to give reasonable certainty that observed differences were due to differences between algorithms.

Another limitation of this study is that, to limit the number of simulations to show the aims of the study, we concentrated on the simple linear response of the chemoreflex gain, although the model can be configured to different curve responses, as shown in Manisty et al. (26). In a linear response chemoreflex gain, its slope near the steady-state CO_{2} is the predominant factor to determine system stability (27).

In this study, only a single example of an unstable respiratory control system was studied in great depth, i.e., a single set of constant values describing chemoreflex response characteristics, chemoreflex delay, lung volume, cardiac output, and metabolic production of CO_{2}. However, in *study 7*, we varied each parameter between its extreme ranges, while maintaining our standard set of values for the others.

We have presented in an online supplement a further sensitivity analysis of the model, in which we have shown how it reacts to combinations of extreme changes to the various physiological parameters that are likely to be found in practice. These sensitivity analyses enable us to appreciate the relative importance of the various physiological parameters as sources of variations to the model outputs. Although these analyses demonstrate that the extreme range of unstable system configurations gives qualitatively similar responses to simulated treatments as those using our standard physiological parameter set, they do not give the detailed behavior of the system for the full spectrum of the input states.

One other limitation of this study is that, for simplicity and brevity, we have not included the possible effect of upper airway resistance in the model. However, in reality, a fall in CO_{2} may well result in the rise of upper airway resistance, which will reduce ventilation in a synergistic manner with the same direction of effect as the chemoreceptor-mediated fall in ventilatory effort. Therefore, we generally expect the overall behavior of the system will be similar, although the instability of a system with upper airway resistance sensitive to CO_{2} will be more enhanced (29).

In this study, we have used only a limited number of chemoreceptors for simplicity and brevity. Future development of the model should include a large number of chemoreceptors with individual gain and delay times.

This is a “basic science” study, applying respiratory and control system modeling. Previous mathematical models have increased our understanding of the pathophysiology of PB. Similarly, mathematical modeling provides insights into the construction of therapeutic algorithms. This study demonstrates that it is feasible to develop an algorithm that could work in real time (since our CO_{2}-admininstering software module uses only real-time input of modeled ventilation and outputs only a desired inspired CO_{2} concentration). We recognize that the next stage of this research is to apply this modeling system to some real cases (with various unpredictable apnea-hyperpnea patterns) to further tune the timing and dosage of CO_{2} delivery and design an experimental study. However, this study gives a good background on the desirable characteristics of the affecting parameters (optimum dose and treatment timing) that an experimental setup would need.

### Conclusions

Although static concentrations of CO_{2} administration can ameliorate PB, very much smaller quantities of CO_{2} can achieve the same degree of amelioration if it can be administered dynamically with careful timing and dosing. Algorithms adjusting both duration and concentration in real time appear to have the greatest potential to do this with minimal unwanted elevation in average CO_{2} levels.

With appropriate development of real-time delivery technology, it might be possible to develop clinical therapies for PB diseases, which use such a dynamic approach, with only very minor increments in systemic CO_{2} levels and, therefore, essentially avoiding the undesirable physiological effects of CO_{2} administration, such as hyperventilation or sympathetic overactivation.

## GRANTS

Y. Mebrate was supported by the Royal Brompton and Harefield National Health Service Trust and St. Mary's Coronary Flow Trust (P11209). C. H. Manisty was supported by the Wellcome Trust (077049/Z/05/Z). K. Willson was supported by Foundation for Circulatory Health and the Coronary Flow trust. D. P. Francis was supported by the British Heart Foundation (FS/04/079).

- Copyright © 2009 the American Physiological Society