AJP - Heart Calcium Transients and Cell-Sarcomere
HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS
 QUICK SEARCH:   [advanced]


     


Am J Physiol Heart Circ Physiol 281: H2714-H2730, 2001;
0363-6135/01 $5.00
This Article
Right arrow Abstract Freely available
Right arrow Full Text (PDF)
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Right arrow Citation Map
Services
Right arrow Email this article to a friend
Right arrow Similar articles in this journal
Right arrow Similar articles in ISI Web of Science
Right arrow Similar articles in PubMed
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Citing Articles
Right arrow Citing Articles via HighWire
Right arrow Citing Articles via ISI Web of Science (15)
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Mukkamala, R.
Right arrow Articles by Cohen, R. J.
Right arrow Search for Related Content
PubMed
Right arrow PubMed Citation
Right arrow Articles by Mukkamala, R.
Right arrow Articles by Cohen, R. J.
Vol. 281, Issue 6, H2714-H2730, December 2001

A forward model-based validation of cardiovascular system identification

Ramakrishna Mukkamala and Richard J. Cohen

Harvard-Massachusetts Institute of Technology Division of Health Sciences and Technology, Massachusetts Institute of Technology, Cambridge, Massachusetts 02139


    ABSTRACT
TOP
ABSTRACT
INTRODUCTION
CARDIOVASCULAR SYSTEM...
FORWARD MODEL DESCRIPTION
FORWARD MODEL-BASED EVALUATION...
FORWARD MODEL VALIDATION
FORWARD MODEL-BASED EVALUATION
DISCUSSION
APPENDIX A
APPENDIX B
APPENDIX C
REFERENCES

We present a theoretical evaluation of a cardiovascular system identification method that we previously developed for the analysis of beat-to-beat fluctuations in noninvasively measured heart rate, arterial blood pressure, and instantaneous lung volume. The method provides a dynamical characterization of the important autonomic and mechanical mechanisms responsible for coupling the fluctuations (inverse modeling). To carry out the evaluation, we developed a computational model of the cardiovascular system capable of generating realistic beat-to-beat variability (forward modeling). We applied the method to data generated from the forward model and compared the resulting estimated dynamics with the actual dynamics of the forward model, which were either precisely known or easily determined. We found that the estimated dynamics corresponded to the actual dynamics and that this correspondence was robust to forward model uncertainty. We also demonstrated the sensitivity of the method in detecting small changes in parameters characterizing autonomic function in the forward model. These results provide confidence in the performance of the cardiovascular system identification method when applied to experimental data.

beat-to-beat model; mathematical modeling; cardiovascular modeling; autonomic nervous system; hemodynamics


    INTRODUCTION
TOP
ABSTRACT
INTRODUCTION
CARDIOVASCULAR SYSTEM...
FORWARD MODEL DESCRIPTION
FORWARD MODEL-BASED EVALUATION...
FORWARD MODEL VALIDATION
FORWARD MODEL-BASED EVALUATION
DISCUSSION
APPENDIX A
APPENDIX B
APPENDIX C
REFERENCES

WHEN ONE CONSIDERS MODELING the cardiovascular system, one usually envisions constructing a model based on physical principles that is capable of generating realistic data. This type of modeling approach, which we refer to as forward modeling, is useful for enhancing one's understanding of cardiovascular physiology. One may also consider an inverse modeling approach, in which models are built from measured data. This type of modeling approach is referred to as system identification when system dynamics or memory is being considered. System identification may potentially provide a powerful means for the intelligent patient monitoring of cardiovascular function. Rather than simply recording hemodynamic signals, the signals are mathematically analyzed so as to provide a dynamical characterization of the physiological mechanisms responsible for coupling them.

To this end, we (26, 27, 29) have previously developed a cardiovascular system identification method for the analysis of short-term, beat-to-beat fluctuations in heart rate signals derived from the surface electrocardiogram (ECG) and noninvasively measured arterial blood pressure (ABP) and instantaneous lung volume (ILV) signals. The method provides a dynamical characterization of the important autonomic and mechanical mechanisms that couple the measured fluctuations specifically in terms of impulse responses. The method also provides power spectra of perturbing noise sources, which represent the residual variability in each of the signals not attributed to the variability generated through the impulse responses.

We have evaluated the cardiovascular system identification method with experimental data collected during conditions of pharmacological autonomic blockade, postural changes, and diabetic autonomic neuropathy (27, 29). We found that these three conditions altered the impulse response and power spectra of perturbing noise source estimates in a manner consistent with known physiological mechanisms. This suggests, but does not directly demonstrate, the validity of the estimated system dynamics.

To evaluate directly the system dynamics estimated by the method, it is necessary to be able to establish the impulse responses and power spectra of perturbing noise sources in a manner independent of system identification. These impulse responses and power spectra of perturbing noise sources may then be regarded as the gold standards against which the corresponding estimates from system identification may be compared. Ideally, one would validate the method based on experimental data. However, the establishment of gold standard impulse responses and power spectra of perturbing noise sources would require extreme experimental conditions that would be virtually impossible to implement in practice. Consider, for example, the establishment of an impulse response, which would require, according to mathematical definition, applying an arbitrarily narrow, unit-area input to the appropriate point of the cardiovascular system and measuring the output response of interest while all other perturbations to the output are held constant. Moreover, even if this type of experiment could be implemented, it is unreasonable to expect that the cardiovascular state and system operating point would be precisely maintained, which renders the meaning of such a comparison to be questionable.

By contrast, the method could be readily evaluated on a theoretical basis with a forward model of the cardiovascular system, because the gold standard impulse responses and power spectra of perturbing noise sources would either be precisely known or easily determined for the relevant cardiovascular state and system operating point. That is, the gold standards could be readily established according to their mathematical definitions. A forward model would also provide a powerful means to analyze the sensitivity of the cardiovascular system identification method. That is, we would be able to determine precisely how much the dynamical properties of the forward model would have to be altered before we would see a corresponding change in the estimates. The major limitation of this type of evaluation is that the results would be only as meaningful as the extent to which the forward model coincides with the actual cardiovascular system. This limitation can be attenuated by determining the robustness of the estimates over a set of models that reflect the uncertainty in the relevant properties of the forward model. We note that the general concept of analyzing inverse modeling algorithms based on forward models of the cardiovascular system is not novel. For example, investigators (9, 16) have utilized complex forward models of the systemic circulation to evaluate the estimation of lumped parameters representing systemic arterial resistance [Ra(t)] and systemic arterial compliance.

The principal goal of this study is to present a theoretical validation of our previously developed cardiovascular system identification method. To realize this aim, this study includes the following five components: 1) a brief review of the cardiovascular system identification method, which has been previously described in detail (26, 27, 29); 2) a description of a forward model of the human cardiovascular system; 3) the procedure for evaluating the system identification method against the forward model, which includes the establishment of gold standards based on their mathematical definitions; 4) the validation of the forward model in terms of power spectra of beat-to-beat variability; and 5) the evaluation of the impulse response and power spectra of perturbing noise source estimates derived by applying system identification to beat-to-beat variability generated from the forward model against the corresponding gold standards in terms of accuracy, robustness, and sensitivity. Note that the final component amounts to assessing the equivalence of the impulse responses and power spectra of perturbing noise sources determined from executing the forward model under two conditions: 1) resting conditions by applying system identification to analyze beat-to-beat variability, and 2) extreme conditions established by the literal meaning of these quantities. Because the extreme conditions cannot be implemented in practice, this study seeks to determine whether the cardiovascular system identification method is able to estimate reliably the physiologically important impulse responses and power spectra of perturbing noise sources during experimental conditions that could be readily implementable.


    CARDIOVASCULAR SYSTEM IDENTIFICATION REVIEW
TOP
ABSTRACT
INTRODUCTION
CARDIOVASCULAR SYSTEM...
FORWARD MODEL DESCRIPTION
FORWARD MODEL-BASED EVALUATION...
FORWARD MODEL VALIDATION
FORWARD MODEL-BASED EVALUATION
DISCUSSION
APPENDIX A
APPENDIX B
APPENDIX C
REFERENCES

Figure 1 illustrates the model upon which the cardiovascular system identification method is based. The model includes five physiological coupling mechanisms relating ECG-derived heart rate signals, ABP, and ILV: circulatory mechanics, heart rate (HR) baroreflex, ILVright-arrowHR, ILVright-arrowABP, and sinoatrial (SA) node.


View larger version (21K):
[in this window]
[in a new window]
 
Fig. 1.   Mean (solid traces) and standard deviation (dotted traces) of cardiovascular system identification estimates obtained from a standing human subject breathing randomly. Figure adapted from raw data described in Ref. 34.

Circulatory mechanics represents the mechanical feedforward effects of a pulsatile heart rate (PHR) signal [defined to be a train of unit-area impulses occurring at the times of ventricular contraction (R wave)] on the pulsatile ABP waveform. HR baroreflex represents the autonomically mediated feedback effects of fluctuations in ABP on fluctuations in an HR tachogram, which is derived from the R-R intervals (3). ILVright-arrowHR represents the autonomically mediated coupling between respiration and HR, which is responsible for governing respiratory sinus arrhythmia. ILVright-arrowABP represents the mechanical effects of respiration on ABP due to the alterations in venous return and the filling of intrathoracic vessels and heart chambers associated with the changes in intrathoracic pressure. SA node maps HR (the net autonomic input to the sinoatrial node) to PHR (the onset times of each ventricular contraction) through an "integrate and fire" device referred to as the integral pulse frequency modulation (IPFM) model (3). Unlike the other four physiological coupling mechanisms, SA node is predefined and not identified from the measured signals.

The model also incorporates two perturbing noise sources, NHR and NABP, which are determined from analysis of the measured signals. NHR represents the fluctuations in HR not caused by fluctuations in ABP or ILV. Such fluctuations may result, for example, from autonomically mediated perturbations driven by cerebral activity. NABP represents fluctuations in ABP not caused by PHR or fluctuations in ILV. Such ABP fluctuations may result, for example, from fluctuations in Ra(t) as tissue beds adjust local vascular resistance to match local blood flow to demand.

To obtain a complete characterization of each of the coupling mechanisms over their physiological range of frequencies, we employ a broadband excitation protocol, in which subjects breathe according to a random sequence of auditory tones (4). We have found that when the subjects are breathing according to this protocol and are otherwise at rest that the fluctuations are sufficiently small and stationary such that the coupling mechanisms may be characterized by linear time-invariant (LTI) transfer functions around a given system operating point (11). We specifically represent each of the coupling mechanisms with autoregressive moving average difference equations, a highly flexible subclass of LTI models, whose parameters may be conveniently identified with the analytic methods of linear least-squares estimation. Further details of the method are provided in (26, 29).

Figure 1 includes an example of the cardiovascular system identification results computed from a standing, healthy subject breathing randomly, in which ABP and ILV are noninvasively measured with a continuous blood pressure monitor (model 2300 Finapres; Ohmeda; Englewood, CO) and a Respitrace System two-belt chest-abdomen inductance plethysmograph (Ambulatory Monitoring; Ardsley, NY) (34). The identified transfer functions are presented in their time domain form of impulse responses. For example, the HR baroreflex impulse response estimate indicates that HR would immediately decrease and then return to baseline if a unit-area impulse of ABP were applied at time 0 while the other inputs to HR (ILV and NHR) were set to 0. The estimated perturbing noise sources in the figure are depicted in their frequency domain form of power spectra. For example, the power spectrum of NHR indicates that ABP and ILV do not fully account for HR variability within ~0.05 Hz.


    FORWARD MODEL DESCRIPTION
TOP
ABSTRACT
INTRODUCTION
CARDIOVASCULAR SYSTEM...
FORWARD MODEL DESCRIPTION
FORWARD MODEL-BASED EVALUATION...
FORWARD MODEL VALIDATION
FORWARD MODEL-BASED EVALUATION
DISCUSSION
APPENDIX A
APPENDIX B
APPENDIX C
REFERENCES

The forward model of the human cardiovascular system to which we apply the cardiovascular system identification method includes three major components: a pulsatile heart and circulation, a short-term regulatory system, and resting physiological perturbations. The pulsatile heart and circulation is a nonlinear, lumped parameter model (R. Mukkamala, D. A. Sherman, R. G. Mark, and R. J. Cohen, unpublished observations). The short-term regulatory system consists of three subcomponents: an arterial baroreflex, a cardiopulmonary baroreflex, and a direct neural coupling mechanism between respiration and HR. The resting physiological perturbations also include three subcomponents: respiratory activity, the autoregulation of local vascular beds, and 1/f HR fluctuations. We describe all of these components and subcomponents and include parameter values below. See APPENDIX A for a description of forward model implementation.

Pulsatile Heart and Circulation

The nonlinear, lumped parameter model of the pulsatile heart and circulation is illustrated in Fig. 2 in terms of its electrical circuit analog, in which charge is analogous to blood volume (Q, ml), current, to blood flow rate (q, ml/s), and voltage, to pressure (P, mmHg). The model consists of six compartments, which represent the left and right ventricles (l and r), systemic arteries and veins (a and v), and pulmonary arteries and veins (pa and pv). Each compartment consists of a conduit for viscous blood flow, which is characterized by either a linear or nonlinear resistance (R) and a volume storage element, which is characterized by either a linear or nonlinear compliance (C) with an associated unstressed volume (Q0). The reference (ref) pressure is atmospheric pressure (or ground) for the systemic compartments and intrathoracic (th) pressure for the ventricular and pulmonary compartments. The compliances of the ventricles vary periodically over time (t) according to the model of Suga and colleagues (38, 39) and are responsible for driving the flow of blood. The four ideal diodes represent the ventricular inflow and outflow valves and ensure unidirectional blood flow. We have demonstrated (R. Mukkamala, D. A. Sherman, R. G. Mark, and R. J. Cohen, unpublished observations) that this model behaves reasonably in terms of pulsatile waveforms and limiting static terminal pressure-flow rate properties of the systemic circulation and heart-lung unit. However, by itself, the model cannot generate the short-term, low-frequency hemodynamic variability that the cardiovascular system elicits during resting conditions. To account for this variability, it is necessary to incorporate a short-term regulatory system and resting physiological perturbations with the model here.


View larger version (12K):
[in this window]
[in a new window]
 
Fig. 2.   Electrical circuit analog of the human pulsatile heart and circulation model. Here, current is analogous to blood flow rate (q, ml/s), whereas voltage is analogous to pressure (P, mmHg). The model includes 6 compartments representing the left and right ventricles (l and r), systemic arteries and veins (a and v), and pulmonary arteries and veins (pa and pv). Each compartment consists of a linear or nonlinear resistance (R) and a linear or nonlinear and time (t)-varying compliance (C) with an associated unstressed volume. Note that each nonlinear element is denoted with a box. The reference pressure is intrathoracic (th) pressure for the ventricular and pulmonary compartments. The 4 ideal diodes represent the ventricular inflow and outflow valves.

Short-Term Regulatory System

Arterial baroreflex. The short-term regulatory system includes a baroreflex system, which is based on a previously developed model (13). This model, which is fully described here, conceptualizes the arterial (A) baroreflex arc as the feedback system illustrated in the block diagram of Fig. 3. The feedback system is aimed at tracking a setpoint (sp) pressure through the following sequence of events. The baroreceptors sense Pa(t) and then relay this pressure via autonomic afferent fibers to the autonomic centers in the brain. Here, the autonomic nervous system (ANS) compares the deviation between the sensed pressure and P<UP><SUB>a</SUB><SUP>sp</SUP></UP> with zero and then sends commands via autonomic efferent fibers to adjust four parameters of the pulsatile heart and circulation such that the ensuing Pa(t) is near P<UP><SUB>a</SUB><SUP>sp</SUP></UP>. The four adjustable parameters are instantaneous HR [F(t)], ventricular contractility [in terms of C<UP><SUB>l,r</SUB><SUP>es</SUP></UP>(t), the differential compliances that the ventricles assume at their unstressed volumes during end systole] (R. Mukkamala, unpublished communications), systemic venous unstressed volume [Q<UP><SUB>v</SUB><SUP>0</SUP></UP>(t)], and Ra(t). The ANS specifically adjusts these parameters based on the history of Pa(t- P<UP><SUB>a</SUB><SUP>sp</SUP></UP> according to the following nonlinear, dynamical map
X(t)=<LIM><OP>∫</OP></LIM> [G<SUB>p</SUB><IT>p</IT>(<IT>&tgr;</IT>)<IT>+G</IT><SUB>s</SUB><IT>s</IT>(<IT>&tgr;</IT>)]<IT>S</IT> (1)

× atan <FENCE><FR><NU><IT>P<SUB>a</SUB></IT>(<IT>t−&tgr;</IT>)<IT>−P</IT><SUP>sp</SUP><SUB><IT>a</IT></SUB></NU><DE><IT>S</IT></DE></FR></FENCE><IT> d&tgr;+X</IT><SUP>sp</SUP>
where X may represent any of the four controllable parameters, the atan function (parameterized by the constant S) imposes arterial baroreflex saturation, p(t) and s(t) are unit-area effector mechanisms, which, respectively, represent the fast, parasympathetic limb of the ANS and the slower, sympathetic limb (see Fig. 4), and Gp and Gs reflect the respective weighting values of the effector mechanisms.


View larger version (7K):
[in this window]
[in a new window]
 
Fig. 3.   Block diagram of the feedback system depicting the arterial baroreflex arc. The feedback system is aimed at tracking a setpoint (sp) pressure through the following sequence of events. The baroreceptors sense Pa(t). The autonomic nervous system (ANS) compares the deviation between the sensed pressure and P<UP><SUB>a</SUB><SUP>sp</SUP></UP> with 0 and then responds by adjusting 4 parameters of the pulsatile heart and circulation to keep the resulting Pa(t) near P<UP><SUB>a</SUB><SUP>sp</SUP></UP>. The four adjustable parameters are instantaneous heart rate [F(t)], ventricular contractility [C<UP><SUB>l,r</SUB><SUP>es</SUP></UP>(t)], systemic venous unstressed volume [Q<UP><SUB>v</SUB><SUP>0</SUP></UP>(t)], and systemic arterial resistance [Ra(t)]. The ANS controls these parameters through a static saturation mapping, followed by convolution with the impulse responses in Fig. 4 (see Eq. 1).



View larger version (11K):
[in this window]
[in a new window]
 
Fig. 4.   Unit-area effector mechanisms representing the fast, parasympathetic limb p(t) (A) and slower, sympathetic limb s(t) (B). These effector mechanisms characterize the dynamical properties of the block labeled ANS in Fig. 3. The signals p(t) and s(t) are derived from experimental data in Ref. 5 (figure adapted from Ref. 13).

To map F(t) to the times of onset of ventricular contraction, an IPFM model of the sinoatrial node is also incorporated. This model integrates F(t) (in units of beats/s) over time until the integral is equal to one. At this point, systole is initiated through the ventricular variable compliance model (with the systolic duration determined from the preceding cardiac cycle length), the integral is set to zero, and the integration is repeated.

Cardiopulmonary baroreflex. The baroreflex system also includes a cardiopulmonary (CP) baroreflex arc, which is conceptualized with a feedback diagram analogous to Fig. 3. However, the sensed pressure here is defined to be the effective right atrial transmural pressure [P<UP><SUB>“ra”</SUB><SUP>tr</SUP></UP>(t) = P"ra"(t)-Pth(t)] of the pulsatile heart and circulation model (see Fig. 2), where P<UP><SUB>ra</SUB><SUP>tr</SUP></UP> is the effective right artrial transmural pressure, Pra is right atrial pressure, and Pth is thoractic pressure.

Direct neural coupling mechanism between respiration and HR. The short-term regulatory system model also includes a direct neural coupling mechanism between respiration and HR. This coupling mechanism is implemented in terms of an LTI impulse response, which maps fluctuations in ILV to fluctuations in F(t). In accordance with previous studies (29, 40) with experimental data, the impulse response is defined here by a linear combination of s(t) and p(t), each of which are advanced in time by 1.5 s. The model ILV signal [Qlu(t)] that is necessary for the implementation of this mechanism is described below.

Parameter values. The numerical values of the parameters characterizing the baroreflex system are taken from (13 and T. Heldt, E. B. Shim, R. D. Kamm, and R. G. Mark, unpublished observations) (see below for exceptions) and are provided in Table 1. Note that the cardiopulmonary baroreflex control of F(t) and c<UP><SUB>l,r</SUB><SUP>es</SUP></UP>(t) are both neglected in the model due to either controversial data or a lack of available data. Although the alpha - and beta -sympathetic sublimbs are both characterized by s(t), the table specifies the particular sublimb corresponding to each s(t). This permits the modeling of, for example, the effects of propranolol, which selectively blocks the beta -sympathetic sublimb.

                              
View this table:
[in this window]
[in a new window]
 
Table 1.   Summary of the parameter values characterizing the arterial and cardiopulmonary baroreflexes

The weighting values for the direct neural coupling mechanism are taken from (29) and are as follows (in ml): -0.0002 beats/s for s(t) and 0.00012 beats/s for p(t). These values indicate that upon inspiration, there is a withdrawal of parasympathetic activity, followed by a withdrawal of beta -sympathetic activity.

Resting Physiological Perturbations

Respiratory activity. The resting physiological perturbations model includes both metronome (met) and random-interval respiratory activity represented in terms of Qlu(t). Because ILV strongly resembles a sinusoid during metronomic breathing, Qlu(t) is mathematically represented as follows
Q<SUB>lu</SUB>(<IT>t</IT>)<IT>=</IT>Q<SUB>fr</SUB><IT>+</IT><FR><NU>Q<SUB>t</SUB></NU><DE>2</DE></FR> <FENCE>1<IT>−</IT>cos <FENCE><FR><NU>2<IT>&pgr;t</IT></NU><DE><IT>T</IT><SUB>r</SUB></DE></FR></FENCE></FENCE> (2)
where Qfr is the functional residual volume of the lungs, Qt is the tidal volume, and Tr is the respiratory period. The model here may also be considered to represent ILV during normal, spontaneous breathing, which, at least to first approximation, appears to be sinusoidal as well.

In the random-interval breathing protocol, the period of each respiratory cycle is determined from the outcome of an independent probability experiment governed by a modified exponential probability density function, in which the respiratory period ranges from 1-15 s with a mean of 5 s (4). The subjects are allowed to control their tidal volume so that normal ventilation is not disrupted. Each respiratory cycle during random-interval breathing is therefore modeled to be one period of a sinusoid with the period determined from the experimental outcome of the probability experiment and the tidal volume chosen such that the alveolar (alv) ventilation rate is identical to that from metronome breathing. This amounts to the following mathematical representation for random-interval Qlu(t)
Q<SUB>lu</SUB>(<IT>t</IT>)<IT>=</IT>Q<SUB>fr</SUB><IT>+</IT><FR><NU><A><AC>q</AC><AC>˙</AC></A><SUP>met</SUP><SUB>alv</SUB><IT>T</IT><SUB>r</SUB>(<IT>j</IT>)<IT>+</IT>Q<SUB>ds</SUB></NU><DE>2</DE></FR> <FENCE>1<IT>−</IT>cos <FENCE><FR><NU>2<IT>&pgr;</IT>(<IT>t−t<SUB>j</SUB></IT>)</NU><DE><IT>T</IT><SUB>r</SUB>(<IT>j</IT>)</DE></FR></FENCE></FENCE><IT>,</IT> (3)

<IT>t<SUB>j</SUB>≤t<t<SUB>j</SUB>+T</IT><SUB>r</SUB>(<IT>j</IT>)
where Qds is the dead space in the airways (air), the argument j denotes the jth respiratory cycle, whereas the subscript j denotes the commencement of the jth respiratory cycle, and
T<SUB>r</SUB>(j)=1−<FR><NU>1</NU><DE>0.2083</DE></FR> ln [1<IT>−x</IT>(<IT>j</IT>)(1<IT>−e</IT><SUP>−2.92</SUP>)] (4)
with x(j) being the outcome of the jth independent probability experiment defined by a uniform distribution ranging from 0-1.

To account for the mechanical effects of respiratory activity on intrathoracic pressure, we incorporate the simple model of ventilation illustrated in Fig. 5 in terms of its electrical circuit analog. The electrical components may be interpreted similarly to those in Fig. 2 by considering air here rather than blood. Hence, the resistor may be thought of as a conduit for airflow between the atmosphere and the lungs whereas the capacitor may be interpreted as an air volume container representing the lung compartment, which is parameterized by an unstressed volume in addition to a linear compliance. The governing differential equations of this model provide a precise mapping from Qlu(t) to Pth(t) as well as to Palv(t), which influences qpa(t) (R. Mukkamala, D. A. Sherman, R. G. Mark and R. J. Cohen, unpublished observations).


View larger version (7K):
[in this window]
[in a new window]
 
Fig. 5.   Model of ventilatory mechanics in terms of its electrical circuit analog, in which charge is analogous to air volume (Q) and voltage is analogous to air pressure (P). The resistor (R) represents the airway (air) between the atmosphere (atm) and the lungs (lu), whereas the capacitor (C) represents the lung compartment. Pth(t) and Palv(t), intrathoracic and alveolar pressures, respectively.

Autoregulation of local vascular beds. Although respiratory-induced hemodynamic fluctuations are fairly well understood, the mechanisms responsible for hemodynamic variability at frequencies below ~0.1 Hz have not been adequately elucidated. It has been previously demonstrated that ABP fluctuations are not determined by HR fluctuations at these lower frequencies, and it has therefore been postulated that fluctuations in Ra(t) caused by the autoregulation of local vascular beds are responsible for these ABP fluctuations, which, in turn, induces HR fluctuations through the baroreflex (1, 24). However, little is known about the system dynamics characterizing autoregulatory processes except that these processes are relatively slow with a characteristic time constant from ~5 to 20 s (6). We therefore include into the model of resting physiological perturbations a zeroth-order model of autoregulatory processes. This model represents these processes as a stochastic, exogenous disturbance to Ra(t) [nRa(t)] defined by a Gaussian white noise process with zero mean and sigma 2 variance that is bandlimited to 0.1 Hz.

1/f HR fluctuations. Figure 1 demonstrates that ABP and ILV are not the only factors responsible for perturbing HR at frequencies within ~0.05 Hz. However, the other perturbing factors, which may include cerebral activity impinging on the ANS, are essentially unknown. It is well known that HR fluctuations demonstrate fractal behavior over at least four decades of frequency, in which the highest frequency decade is within the frequency band of interest here (7, 37). We therefore include in the resting physiological perturbations model an exogenous, stochastic 1/f disturbance to F(t) [nF(t)] The disturbance is created by passing a zero-mean, Gaussian white noise process with variance lambda 2 through a cascade combination of two LTI filters. One filter is of unit-area and is designed to have a 1/f magnitude squared frequency response over four decades of frequency starting at 10-4 Hz (12, 20, 36). The other filter is defined by a linear combination of s(t) (beta -sympathetic) and p(t) (with arbitrary weighting values of 12 beats/s each), because HR fluctuations are almost exclusively mediated by the ANS (2). Importantly, nF(t) as well as nRa(t) are considered to be unmeasurable quantities in the forward model.

Parameter values. The parameter values characterizing the models of Qlu(t) and ventilation are provided in Table 2. These values are taken from (21, 41).

                              
View this table:
[in this window]
[in a new window]
 
Table 2.   Summary of the parameter values characterizing the models of Qlu(t) and ventilation

Both sigma 2 and lambda 2, which characterize nRa(t) and nF(t), respectively, are considered to be free parameters whose values are chosen such that the model low-frequency hemodynamic variability matched that determined from a set of experimental data previously obtained from 12 healthy standing humans breathing at a fixed rate of 0.25 Hz (10). The data set consists of simultaneous, noninvasive measurements of ECG and ABP (Finapres), as well as uncalibrated ILV derived from the ECG for each subject. The specific procedure for choosing the parameter values based on this data set is as follows. We computed the mean ± SD of the power in three nonoverlapping frequency bands for HR and ABP over the group of 12 subjects (see Table 3) by using autoregressive spectral analysis (22). The free parameters sigma 2 and lambda 2 were then tuned during fixed-rate breathing excitation at 0.25 Hz, such that the power in these frequency bands for F(t) and Pa(t) were near their respective human values. To satisfy sufficiently this matching procedure, we also had to increase the weighting values of the parasympathetic and beta -sympathetic effector mechanisms by 33% from their original values given in (13). However, this increase did not move these values outside of their physiological range, as established by experimental data (14). On the basis of this procedure, we set sigma  = 0.16 mmHg · s · ml-1 and lambda  = 0.0035 beats/s.

                              
View this table:
[in this window]
[in a new window]
 
Table 3.   Power in three low-frequency bands of HR and ABP as determined from the forward model and experimental human data during fixed-rate breathing at 0.25 Hz

On the basis of this procedure, the hemodynamic variability generated by the forward model may be thought of as representing an "average" subject. We did not tune the model to represent an individual subject because there is tremendous intersubject variability in experimental hemodynamic spectra. If we had chosen to represent an individual subject, the results of this study would be biased toward the chosen subject.


    FORWARD MODEL-BASED EVALUATION PROCEDURE
TOP
ABSTRACT
INTRODUCTION
CARDIOVASCULAR SYSTEM...
FORWARD MODEL DESCRIPTION
FORWARD MODEL-BASED EVALUATION...
FORWARD MODEL VALIDATION
FORWARD MODEL-BASED EVALUATION
DISCUSSION
APPENDIX A
APPENDIX B
APPENDIX C
REFERENCES

Gold Standard Cardiovascular System Identification Results

To evaluate the performance of the cardiovascular system identification method applied to data generated from the forward model, it is necessary to establish, in a manner independent of system identification, the impulse responses and power spectra of perturbing noise sources in Fig. 1, which characterize the forward model. These impulse responses and power spectra of perturbing noise sources may then be regarded as the gold standard cardiovascular system identification results against which the estimates may be compared.

We establish the gold standards according to the literal mathematical meaning of the impulse responses and the power spectra of the perturbing noise sources. In particular, the impulse response is defined to represent the output response of a system to an arbitrarily narrow, unit-area input. The impulse response also completely characterizes the input-output relationship of the physiological coupling mechanisms here because of our LTI assumption. Hence, in establishing each of the gold standard impulse responses, we apply an impulse input to the desired point in the forward model and measure the output response of interest while holding all other perturbations to the output constant. For example, to establish the gold standard ILVright-arrowABP impulse response, we measure the ABP response to an impulse input of ILV while holding HR constant and setting the exogenous disturbance to Ra(t) to zero.

By definition, the perturbing noise source represents the residual variability in a measured output not due to the fluctuations in the other measured inputs. So, in establishing each of the gold standard power spectra of perturbing noise sources, we first measure the output of interest in the forward model while holding the two measured input fluctuations constant and then compute the power spectrum of the output. For example, to establish the gold standard power spectrum of NABP, we first measure ABP in the forward model while holding HR and ILV constant and then compute the power spectrum of the resulting ABP.

The specific procedures for establishing each of the gold standard impulse responses and power spectra of perturbing noise sources in Fig. 1 are described in APPENDIX B. This description includes the signal processing necessary for the estimated and gold standard impulse responses to be sampled identically.

Root-Normalized-Mean-Squared Error

To provide a compact means for evaluating the similarity of each of the estimates with respect to its corresponding gold standard, we utilize a statistic, which we refer to as the root-normalized-mean-squared error (RNMSE). Let the estimate be vector x with mean vector x and covariance matrix Lambda x (which reflects the uncertainty in the estimate) and its corresponding gold standard be vector x0. The RNMSE is then defined as follows
RNMSE (<B>x</B><IT>, </IT><B>x</B><SUB>0</SUB>)<IT>=</IT><RAD><RCD><FR><NU><IT>E</IT>[(<B>x</B><IT>−</IT><B>x</B><SUB>0</SUB>)<SUP><IT>T</IT></SUP>(<B>x</B><IT>−</IT><B>x</B><SUB>0</SUB>)]</NU><DE><IT>E</IT>(<B>x</B><SUP><IT>T</IT></SUP><SUB>0</SUB><B>x</B><SUB>0</SUB>)</DE></FR></RCD></RAD><IT>·</IT>100<IT>%</IT> (5)

<IT>=</IT><RAD><RCD><FR><NU>trace (<B>&Lgr;</B><SUB><IT>x</IT></SUB>)</NU><DE><B>x</B><SUP><IT>T</IT></SUP><SUB>0</SUB><B>x</B><SUB>0</SUB></DE></FR><IT>+</IT><FR><NU>(<B><A><AC>x</AC><AC>ˆ</AC></A></B><IT>−</IT><B>x</B><SUB>0</SUB>)<SUP><IT>T</IT></SUP>(<B><A><AC>x</AC><AC>ˆ</AC></A></B><IT>−</IT><B>x</B><SUB>0</SUB>)</NU><DE><B>x</B><SUP><IT>T</IT></SUP><SUB>0</SUB><B>x</B><SUB>0</SUB></DE></FR></RCD></RAD><IT>·</IT>100<IT>%</IT>
where E( · ) is the expectation operator and trace( · ) is the operator, which sums the diagonal elements of its matrix argument. This statistic may be thought of as representing the percentage of error of the estimate with respect to its gold standard. In the case of the power spectral estimates, in which no measure of uncertainty is available, the matrix Lambda x is assumed to be a matrix of zeros.

Monte Carlo Simulations

Although the impulse response estimates include a measure of uncertainty in terms of a covariance matrix, this uncertainty measure is only an estimate itself (31). To account more accurately for estimation error variance, we report the mean ± SD for the impulse response estimates as well as for the power spectral estimates and RNMSE statistics as determined from 20 different realizations of forward model generated data.


    FORWARD MODEL VALIDATION
TOP
ABSTRACT
INTRODUCTION
CARDIOVASCULAR SYSTEM...
FORWARD MODEL DESCRIPTION
FORWARD MODEL-BASED EVALUATION...
FORWARD MODEL VALIDATION
FORWARD MODEL-BASED EVALUATION
DISCUSSION
APPENDIX A
APPENDIX B
APPENDIX C
REFERENCES

The forward model-based evaluation of the cardiovascular system identification method is only meaningful provided that the low-frequency power spectral content of the signals analyzed by the method is realistic. Because the forward model is tuned to represent what one may think of as an "average" subject, we initially considered demonstrating that the model spectra resemble the group average spectra of the 12 healthy humans breathing at a fixed rate of 0.25 Hz (see Parameter Values in Resting Physiological Perturbations). However, averaging tended to smear out the spectral peaks, which were usually not centered at the same frequencies for each subject. We therefore demonstrate that the model spectra resemble what one may find experimentally through a comparison with the spectra from one of the healthy subjects. This does not seem unreasonable if we keep in mind that the power in three nonoverlapping frequency bands of the model spectra quantitatively match that computed from the group average of the 12 subjects (see Table 3).

Figure 6 illustrates that the power spectra for ILV, HR, and ABP at frequencies below the mean HR obtained from the forward model during fixed-rate breathing excitation at 0.25 Hz indeed resemble spectra obtained from one of the 12 subjects. Differences at the respiratory frequency in Fig. 6 and in Table 3 may be attributed to the fact that the subject here did not precisely follow the fixed-rate breathing protocol as well as discrepancies in tidal volume, which was not measured. Note that although the forward model was tuned to represent the "average" subject, it could have been tuned to match more precisely the subject here.


View larger version (19K):
[in this window]
[in a new window]
 
Fig. 6.   Power spectra of instantaneous lung volume (ILV; normalized to unity power), HR, and arterial blood pressure (ABP) at frequencies below the mean HR from the forward model and a single, standing human during fixed-rate breathing at 0.25 Hz. The model spectra represent the average over 10 different realizations of forward model-generated data, in which qalv = 70 ml/s.

Figure 6 also demonstrates that the model spectra exhibit a spectral peak at ~0.07 Hz, which is near the center frequency of the "posture" peak in humans (32). Because the exogenous disturbances within 0.1 Hz are broadband (see Autoregulation of Local Vascular Beds and 1/f HR Fluctuations), this peak implies a system resonance in the forward model. This result supports the simple computer simulation of de Boer (14), which demonstrated that the "posture" peak could be due to a system resonance. de Boer specifically implicated the system resonance to the arterial Ra(t) baroreflex arc. However, based on a few simple experiments, in which the weighting values of the open-loop effector mechanisms were varied, we have found that the resonance peak is substantially diminished in the absence of the arterial baroreflex feedback pathways responsible for controlling Q<UP><SUB>v</SUB><SUP>0</SUP></UP>(t) as well as Ra(t). That is, the arterial baroreflex controlling alpha -sympathetic activity in general is involved in eliciting the resonance peak.

It is also necessary to demonstrate that the cross spectra between each of these signals are realistic. This essentially amounts to demonstrating that the transfer functions, which characterize the couplings between the signals are consistent with those determined from experimental data. However, this has been virtually taken care of by the very construction of the model, which was based largely on physiological findings published in the literature.

Although it is only necessary for the purposes of this study that the low-frequency spectral content of ILV, HR, and ABP be reasonable, we assume that the low-frequency content of the remaining model signals are accounted for as well by virtue of reasonably representing these three signals. Figure 7 gives some credence to this assumption by illustrating that the model spectrum of left ventricular flow rate (LVFR) at frequencies below the mean HR resembles the corresponding spectrum measured with a Doppler ultrasound technique (17) from an individual subject breathing at a fixed rate of 0.25 Hz and tilted upright (30° with respect to the supine posture).


View larger version (16K):
[in this window]
[in a new window]
 
Fig. 7.   Power spectra of left ventricular flow rate (LVFR) at frequencies below the mean HR from the forward model and a single human tilted upright 30° with respect to the supine posture during fixed-rate breathing at 0.25 Hz. The model spectrum represents the average over 10 different realizations of forward model-generated data, in which qalv = 70 ml/s.


    FORWARD MODEL-BASED EVALUATION
TOP
ABSTRACT
INTRODUCTION
CARDIOVASCULAR SYSTEM...
FORWARD MODEL DESCRIPTION
FORWARD MODEL-BASED EVALUATION...
FORWARD MODEL VALIDATION
FORWARD MODEL-BASED EVALUATION
DISCUSSION
APPENDIX A
APPENDIX B
APPENDIX C
REFERENCES

Accuracy Analysis

We applied the cardiovascular system identification method to the beat-to-beat variability generated from the forward model during random-interval breathing excitation. Figure 8 illustrates that there is a close agreement between the resulting cardiovascular system identification estimates and the corresponding gold standards. These results indicate the validity of the cardiovascular system identification method to the extent that the forward model coincides with the actual cardiovascular system.


View larger version (21K):
[in this window]
[in a new window]
 
Fig. 8.   Mean (solid traces) and standard deviation (dotted traces) of cardiovascular system identification estimates along with their respective gold standards (dashed traces), which characterize the forward model.

Robustness Analysis

Perhaps the major source of uncertainty in the relevant properties of the forward model involves the system dynamics responsible for controlling HR. We therefore consider here the robustness of the HR baroreflex, ILVright-arrowHR, and NHR estimates over a set of forward models, which reflect some of our uncertainty in these system dynamics.

Cardiopulmonary HR baroreflex. We first consider how these estimates would be altered if a cardiopulmonary HR baroreflex were included in the forward model. Inclusion of this reflex may be achieved simply by adjusting its weighting values to nonzero values (see Table 1). To proceed with this analysis, the gold standard NHR power spectrum, which is no longer valid in the presence of a cardiopulmonary HR baroreflex, must be redefined (see APPENDIX C for the specific procedure).

Figure 9, A-C, illustrates the RNMSE results of the HR baroreflex, ILVright-arrowHR, and NHR estimates as a function of the ratio of the static gains of the cardiopulmonary to arterial HR baroreflexes. The range of this ratio is determined from published experimental data (15, 35). Note that a negative ratio indicates a Bainbridge type of cardiopulmonary baroreflex, whereas a positive ratio indicates a type of cardiopulmonary baroreflex, which contributes to ABP regulation. The results indicate that the ILVright-arrowHR impulse response estimate is most significantly and adversely influenced by the presence of the cardiopulmonary HR baroreflex. This result implies that low-frequency effective right atrial transmural pressure fluctuations are determined largely by ILV fluctuations.


View larger version (15K):
[in this window]
[in a new window]
 
Fig. 9.   Root-normalized-mean-squared error (RNMSE) results (mean ± SD) for HR baroreflex (A), ILVright-arrowHR (B), NHR, (C) and ILVright-arrowHR (D) with respect to the redefined gold standard (see APPENDIX C) as a function of the ratio of the static gains of the cardiopulmonary to arterial HR baroreflexes.

From Fig. 9A, we see that the HR baroreflex impulse response estimate is not affected by a positive ratio of the static gains of the cardiopulmonary to arterial HR baroreflexes. However, for a negative ratio, the estimate is adversely affected but not as significantly as the ILVright-arrowHR impulse response estimate (consider RNMSE standard deviation). We hypothesize that this discrepancy between positive and negative ratios is due to the size of low-frequency ABP fluctuations. That is, in contrast to the positive ratio, which results in tighter ABP regulation, the negative ratio leads to sufficiently large low-frequency ABP fluctuations such that they contribute somewhat to the generation of low-frequency effective right atrial transmural pressure fluctuations in the model.

Figure 9B illustrates that the ILVright-arrowHR estimate becomes progressively worse with respect to its corresponding gold standard with increasing absolute static gain of the cardiopulmonary HR baroreflex. However, this does not necessarily imply that the estimate is meaningless for large absolute static gain values. In fact, the estimate is quite meaningful regardless of the static gain value of the cardiopulmonary HR baroreflex provided that we redefine the gold standard such that it encompasses the dynamics of both the direct neural coupling mechanism and the cardiopulmonary HR baroreflex (see APPENDIX C for the specific procedure). Figure 9D shows that the RNMSE results for the ILVright-arrowHR estimate with respect to the new gold standard is essentially independent of the ratio of the static gain values of the cardiopulmonary to arterial HR baroreflexes. The results here suggest that when considering the ILVright-arrowHR impulse response identified from experimental data, it is probably more accurate to interpret the estimate analogously to the gold standard as redefined here.

Arterial baroreflex saturation. The S parameter in Eq. 1 characterizes the degree of arterial baroreflex saturation in the forward model. This parameter provides an upper bound on the deviation between sensed ABP and its setpoint pressure so that the manipulated variables cannot be controlled to arbitrarily high or low values. The maximum deviation permitted by this parameter in the forward model is ~28 mmHg. This pressure is approximately equal to the range of ABP, over which the atan mapping in Eq. 1 differs from true linearity by ~10% (that is, the linear ABP range). However, experimental data from one study (25) indicate that the linear ABP range may actually be only ~10 mmHg. We therefore consider the robustness of the HR baroreflex, ILVright-arrowHR, and NHR estimates against more narrow linear ABP ranges.

Figure 10 shows the RNMSE results of the HR baroreflex, ILVright-arrowHR, and NHR estimates as a function of the linear ABP range of the arterial baroreflex. These results indicate that only the HR baroreflex estimate is significantly and adversely affected by increasing degrees of saturation. Furthermore, the deviation of the RNMSE from its nominal value only becomes significant when the extent of saturation is no longer substantiated by the experimental data in (25). Hence, provided that arterial baroreflex saturation is the only significant nonlinearity, then the assumption that the fluctuations in cardiovascular system identification data are small enough such that the couplings between the fluctuations are related linearly seems quite tenable.


View larger version (13K):
[in this window]
[in a new window]
 
Fig. 10.   RNMSE results (means ± SD) for HR baroreflex (A), ILVright-arrowHR (B), and NHR (C) as a function of the linear ABP range (as defined in the text) of the arterial baroreflex.

Sensitivity Analysis

The ultimate potential of the cardiovascular system identification method is to provide a clinician with an intelligent, sensitive tool for monitoring a patient's cardiovascular state over time so as to guide therapy. An analysis of the resolving power of the method is germane to the realization of this potential. We therefore consider here an analysis of the sensitivity of the method particularly in terms of detecting changes in autonomic activity, as reflected by the weighting values for p(t) and beta -sympathetic s(t) in the forward model, through the HR baroreflex, ILVright-arrowHR, and NHR estimates.

Figure 11, A and B, shows the sensitivity results to changes in autonomic function in terms of the estimated versus actual percentage of change in the static gains (sum of the weighting values) of the autonomically mediated impulse responses, where the percentage of change here is with respect to the nominal static gain values. These results indicate that a change of 25% in the static gain of the gold standard HR baroreflex impulse response is required to detect reliably a change in the static gain of the corresponding estimate. On the other hand, a 25% change in the static gain of the gold standard ILVright-arrowHR impulse response is not sufficient to detect reliably a change in the static gain of the corresponding estimate.


View larger version (20K):
[in this window]
[in a new window]
 
Fig. 11.   Sensitivity results (solid trace; means ± SD) in terms of estimated versus actual percentage of change in the static gain of ILVright-arrowHR (A), static gain of HR baroreflex (B), absolute peak amplitude of ILVright-arrowHR (C), and absolute peak amplitude of HR baroreflex (D). The percentage of change is with respect to the nominal values and the dashed trace is the identity line.

On the basis of our experience with experimental human data (27, 29), we have found that the absolute peak amplitude of the impulse response estimates seems to be quite reliable in detecting changes in autonomic function. Figure 11, C and D, shows the sensitivity results in terms of the estimated versus actual percentage of change in the absolute peak amplitude of the autonomically mediated impulse responses. Because we varied the static gains of the impulse responses simply via scaling, the actual percentage of change in the static gain of the impulse responses is equal to the actual percentage of change in the absolute peak amplitude. The results indicate substantial improvement in sensitivity with respect to the ILVright-arrowHR estimate and some improvement with respect to the HR baroreflex estimate. That is, as small as a 10% change in the absolute peak amplitude of the gold standard ILVright-arrowHR impulse response is sufficient to detect reliably a change in the absolute peak amplitude of the corresponding estimate.

The substantial improvement in the sensitivity of the ILVright-arrowHR impulse response may be explained by considering the ratio of the spectra of ILV to the spectra of the unmeasured HR disturbance (signal-to-noise ratio as a function of frequency) in the model. This ratio is smallest at frequencies near 0 Hz because of the 1/f character of the unmeasured disturbance as well as the reduced ILV spectral content near 0 Hz due to the probability density defining the random-interval respiratory pattern (see Respiratory Activity). The static gain of the ILVright-arrowHR impulse response estimate, which reflects the 0 Hz estimate, is therefore not reliable. However, the ratio is larger at higher frequencies, and the absolute peak amplitude, which reflects the wider bandwidth parasympathetic filter, is therefore quite reliable. On the other hand, the spectral content of ABP is sufficiently significant with respect to that of the unmeasured 1/f HR fluctuations at frequencies near 0 Hz such that the static gain of the HR baroreflex impulse response is fairly reliable.

Figure 12 illustrates the sensitivity results in terms of the estimated and actual percentage of change in total, low-frequency (0-0.15 Hz), and high-frequency (0.15-0.4 Hz) power of the NHR spectra. These results emphasize the superior reliability of the high-frequency components of the estimates due to the 1/f character of the unmeasured HR disturbance. Note that the high-frequency power of NHR is just as sensitive a measure of parasympathetic function as the absolute peak amplitude of ILVright-arrowHR. However, we believe that this absolute peak amplitude may be a better measure of parasympathetic function, because the perturbing noise source NHR, unlike ILVright-arrowHR, is not normalized for inputs, e.g., cerebral activity.


View larger version (19K):
[in this window]
[in a new window]
 
Fig. 12.   Sensitivity results (solid trace; mean ± SD) in terms of estimated versus actual percentage of change in the total power (A), low-frequency (LF) power (0.0-0.15 Hz) (B), and high-frequency (HF) power (C) (0.15-0.4 Hz) of the NHR power spectrum. The percentage of change is again with respect to the nominal values, and the dashed trace again denotes the identity line.


    DISCUSSION
TOP
ABSTRACT
INTRODUCTION
CARDIOVASCULAR SYSTEM...
FORWARD MODEL DESCRIPTION
FORWARD MODEL-BASED EVALUATION...
FORWARD MODEL VALIDATION
FORWARD MODEL-BASED EVALUATION
DISCUSSION
APPENDIX A
APPENDIX B
APPENDIX C
REFERENCES

Previous Forward Models

A few forward models of the cardiovascular system have been previously developed for the purposes of eliciting and analyzing low-frequency hemodynamic variability. de Boer (14) presented a beat-to-beat model, which incorporated windkessel and Starling properties with arterial baroreflex control of HR and systemic arterial resistance. The model was perturbed with sinusoidal respiration and stochastic disturbances to HR and pulse pressure to generate resting hemodynamic fluctuations. Madwed et al. (23) developed a model to study the generation of Mayer waves. Their model included windkessel properties, constant stroke volume, arterial baroreflex control of HR and systemic arterial resistance, and sinusoidal respiration. Cavalcanti (8) presented a model that incorporated windkessel properties with nonlinear arterial baroreflex control of HR and assumed stroke volume and systemic arterial resistance to be constant. Despite the absence of exogenous perturbations, the model was demonstrated to elicit HR fluctuations due to the nonlinear delayed structure of closed-loop control.

None of these simple models are nearly as comprehensive as our forward model in terms of accounting for pulsatility, the short-term regulatory system, and resting physiological perturbations. Our forward model thus provides a more complete means for studying the mechanisms responsible for eliciting beat-to-beat variability compared with these oversimplified models. Consider, for example, our more thorough analysis of the physiological mechanisms responsible for the system resonance phenomenon with respect to the study of de Boer (14) (see FORWARD MODEL VALIDATION).

Forward Model Limitations

The mechanisms responsible for eliciting low-frequency hemodynamic variability, particularly at frequencies below ~0.1 Hz, have not been fully elucidated. Nevertheless, we constructed a forward model of the human cardiovascular system that is capable of generating reasonable low-frequency hemodynamic variability and is largely consistent with previous experimental findings. However, the forward model presented here is not the unique model for simultaneously realizing this behavior and being consistent with experimental findings. That is, we could alter the properties of the forward model and still generate realistic low-frequency hemodynamic variability. For example, the resonance peak at ~0.07 Hz could also be elicited with a stochastic, exogenous disturbance to Q<UP><SUB>v</SUB><SUP>0</SUP></UP>(t), which may reflect, for example, fast-acting hormonal loops (19). In the forward model, this disturbance would induce ABP fluctuations, and consequently HR fluctuations, through stroke volume fluctuations. This is consistent with experimental findings that have hypothesized systemic arterial resistance fluctuations, but have not precluded stroke volume fluctuations, in the genesis of low-frequency ABP and HR fluctuations (1, 24). The forward model as a general beat-to-beat variability model should therefore be considered with caution. Importantly, however, the forward model as a means for evaluating our cardiovascular system identification method is more reliable by virtue of the realistic autospectra of ILV, HR, and ABP at frequencies below the mean HR (Fig. 6) and the realistic cross spectra between these signals (Fig. 8).

Robustness to Other Forward Model Uncertainties

Although we evaluated how the presence of a cardiopulmonary HR baroreflex and varying degrees of arterial baroreflex saturation would affect the autonomically mediated estimates, we did not address the effects of nonstationarities or other types of nonlinearities on these estimates. However, the stationarity property of the forward model seems quite tenable given that the data are considered over short time periods during stable