- Email: [email protected]
The analysis of pulmonary and cardiovascular parameters with a mathematical model of external respiration
MATHEMATICAL
BIOSCIENCES
13, 369-389
(1972)
369
The Analysis of Pulmonary and Cardiovascular Parameters with a Mathematical Model of External Respiration * J. A. KELMAN
AND
L. D. HOMER?
Division of Biological and Medical Sciences Brown University, Providence, Rhode Island 02912
Communicated
by F. Grodins
ABSTRACT The behavior of the external respiratory system was described by a set of differential equations. The system was solved by digital computer, to express the relationship between stationary cycle blood gas tensions and fixed cardiac and respiratory parameters. Dynamic respiratory maneuvers were simulated for different values of these parameters to investigate the feasibility of estimating cardiac output, respiratory dead space, and ventilation-perfusion distribution, from flow and composition of expired air, using non-linear regression techniques.
INTRODUCTION
Christiansen, Douglas, and Haldane (1914) made the original attempt to estimate blood gas parameters with a bloodless technique utilizing respiratory data. Kim, Rahn, and Farhi (1966) have developed a more sophisticated method of estimating venous and arterial PCOZ from respiratory maneuvers, and recently Bickel, Diener, and Brammell (1970) have devised a computerized data gathering system to estimate cardiac output using this technique and the Fick Equation. The physiological basis for these and other studies in this area, is the rapid response to changes in cardiac output, of the partial pressure of gases at the mouth. The qualitative events are often obvious, and the major problems in the investigations involve the quantitative predictions. The accuracy of these predictions may be expected to depend on the exactness of the assumptions on which they are based, as well as the amount and accuracy of the data analyzed. Any set of assumptions contains the concept of a model, and their validity is limited to the extent that the model can accurately simulate the real system. It is our intention to apply the techniques developed for * This work supported in part by ONR Contract NOOO14-67-A-0191-0013. t Reprint requests should be addressed to: L. D. Homer, Naval Medical Research Institute National Naval Medical Center, Bethesda, Md. 20014 U.S.A. Copyright 0 1972 by American Elsevier Publishing Company, Inc.
370
J. A. KELMAN
AND L. D. HOMER
pulmonary control system modeling, to the analysis and estimation of pulmonary and cardiovascular parameters. The first mathematical model for mammalian respiration was proposed by Gray (1945). It consisted of a system of equations which could be solved algebraically for steady state values. An advanced model for the determination of steady state respiratory information has been devised by West (1969), including complex non-differential equations describing blood gas dissociation and ventilation-perfusion phenomena, which are solved by numerical means. Grodins, Gray, Schroeder, Norins, and Jones (1954) published the first mathematical description of a respiratory model that allowed the solution of transient values, using an analog computer to predict respiratory parameters arising in response to changes in carbon dioxide inhalation. Grodins’ model has since served as the basis for many other analyses of the respiratory system, involving analog and digital simulation, which are discussed in a recent review by Yamamoto and Raub (1967). Many of the ideas in our model are based on Grodins’ work, although applied towards a different goal. In the bulk of the respiratory models published to date, the major emphasis has been placed on the “closed loop” system, with the implementation of different “ventilation equations” in an effort to simulate respiratory feedback control. We are interested here, in an open-loop respiratory model with no internal control. All the parameters can be set independently, so that the quantitative effect of changes in a single parameter, on respiratory function, can be isolated and evaluated independently of changes in the other parameters of the system. By this kind of simulation, we can specify those non-pulmonary parameters that it is possible to estimate by the measurement of expired air, and suggest numerical procedures for doing so. This publication includes the system of equations comprising our model, the methods of digital solution of these equations, and some preliminary results indicating the quantitative changes in blood gas tensions due to changes in cardiac output, ventilation, and ventilation-perfusion distribution.
Tl IEORY
In our analysis of gas exchange, the numerical model constructed contains ten homogenous compartments (Fig. 1). The gas flow between a single dead space compartment and two alveolar compartments is controlled by forced variation in alveolar volume. Blood circulation is a closed loop which involves two pulmonary capillary compartments feeding into a single arterial compartment. The arterial blood enters the venous blood
MODEL OF EXTERNAL
RESPIRATION
371
compartment, where blood can equilibrate with the tissue space, and is then returned to the pulmonary capillary beds. The sizes of the three tissue compartments are equal to the volumes of distribution of oxygen, nitrogen, and carbon dioxide in turn, and the blood returning to the lungs is assumed to have the same content as each component has in its tissue space. The values used for the volumes of distribution may not be exact, but have negligible effects on the data we are interested in with the time periods under consideration. As this system is perfusion-limited at the alveolus as well as in the peripheral circulation, blood leaving each set of pulmonary capillaries is assumed to have the same partial pressure as its associated alveolar compartment.
FIG. 1. Flow diagram representing compartmental analysis.
the pulmonary-cardiovascular
model in terms of
The nitrogen tension is related to gas content by Henry’s Law, as it involves only dissolved gas; while oxygen and carbon dioxide tensions are related to content by more complex equations. The oxygen content includes dissolved gas as well as oxygen bound to hemoglobin, while carbon dioxide is present in blood as dissolved gas, as serum bicarbonate, and as carbamino compounds bound to blood protein. The distribution of CO2 is also dependent on hemoglobin content and saturation. The pH of the
372
J. A. KELMAN
AND L. D. HOMER
blood, in each circulatory compartment, follows the Henderson-Hasselbalch equation depending on PCOz and blood bicarbonate (Comroe, 1965). By including the effect of pH on the oxyhemoglobin dissociation curve (the Bohr effect) oxygen content, carbon dioxide content, and pH must be solved for simultaneously in an implicit routine. The set of differential equations that mathematically define the above relationships, is presented here as an ordered set, following the logical sequence suggested by the flow diagram. The equations are written using Fortran conventions for arithmetic. The symbols i, j, 1,2,02, C02, and N2 at the end of a variable name are to be read as subscripts. The subscript i refers to alveolus 1 or 2, and the subscript j refers to oxygen, carbon dioxide, or nitrogen. All the symbols used in the equations are defined in Appendix I. We define gas flow to alveolar compartment i from the dead space, as VDOTi, positive on inspiration. This value is determined by the respiration rate (RATE) and the tidal volume (amp = tidal volume/2.), as well as the distribution of flow term FVi (the fraction of gas flow directed to alveolus i). VDOTi is also adjusted for gas exchange in the alveolus, with the diffusion terms Dij. The quantity VDOTi must be adjusted for gas exchange across the alveolar membrane, and is therefore not identical to the rate of change of alveolar volume i [d(Vai)/dt in Eq. IV]. In effect, Eq. 1 adjusts gas flow at the mouth and at the main bifurcation of the pulmonary airways, so that the alveolar volumes each fluctuates around a fixed value (the “resting volume”). The degree of fluctuation is determined by the fraction of flow to each alveolus. Normal respiratory gas flow, at time t, is approximated with a fixed sine wave forcing function: VDOTi = FVi*2*pi*RATE*amp*sin(2*pi*t*RATE) + [(273. + 37.)/273.]*(760./PTO)*(Di02 + DiC02
+ DiN2)
(I)
Two equations are used to calculate the dead space composition, one for inhalation and one for exhalation. The fraction of gas j in the dead space is Fdj, Foj refers to gas j fraction in the inspired air, and Faij represents gas j fraction in alveolar volume i. VD is the dead space volume. Then for inhalation : d(Fdj)/dt
= (VDOTl
+ VDOT2)*(Foj
- Fdj)/VD
(Ira)
and for exhalation: d(Fdj)/dt
= VDOTl*(Fdj
- FaIj)/VD + VDOT2*(Fdj
- Fa2j)/VD
(IIb)
Separate inhalation and exhalation equations are used also in determining alveolar composition. Gij refers to the quantity of gas j in alveolus i at BTP, rather than gas fraction. PTO is the ambient pressure, and Dij
MODEL OF EXTERNAL
RESPIRATION
373
the alveolar capillary transfer of gas j in alveolus converted to BTP in (IIIa) and (IIIb): d(Gij)/dt d(Gij)/dt
= VDOTi*Fdj = VDOTi*Faij
- [(273. + 37.)/273]*(760./PTO)*Dij
(IIIa)
- [(273. + 37.)/273.]*(760./PTO)*Dij
(IIIb)
The change in total alveolar Equation III over j gases: d(Vai)/dt
= VDOTi
i, in ml. at STP which is
volume
(Vai) is determined
by summing
- ((273. + 37.)/273.)*(760./PTO)r(Di02 + DiC02 + DiN2)
(IV)
or what is the same, d(Vai)/dt
= FVi*2*pi*RATE*amp*sin(2*pi*RATE*t)
Diffusion (D), in ml/min at STP, is calculated from cardiac output (CO), fraction of flow to alveolar compartment i (FBQi), end pulmonary capillary concentration (Cc), and venous return concentration (Cv) : Dij = CO*FQBi*(Ccij
- Cvj)
(V)
Arterial content of gas j (Carj) is altered by influx from two capillary compartments, and outflow into one venous-tissue compartment. Vart refers to arterial compartment volume : d(Carj)/dt
= (l/Vart)*(CO*FQBl*Cclj + CO*FQB2*Cc2j
- CO*Carj)
(VI)
In the same manner, venous-tissue gas contents are calculated from arterial inflow, outflow to the pulmonary capillaries, and metabolic flow (METAB), which represents gas production and consumption in tissue metabolism. METABC02 divided by METAB defines the overall respiratory quotient. Each gas j, has its own volume of distribution (VDISj) : d(Cvj)/dt
= (l/VDISj)*(CO*Carj
- CO*Cvj
The relationship between the partial pressure in the blood presents the most complex problems of this model. Blood gas tension (Pcj) is set in compartment equal to the partial pressure in the Pcij = Paij
- METABj)
(VII)
and the content of gases in the numerical solution the pulmonary capillary alveolus (Paij): (VIII)
The calculation of gas content in the blood from the partial pressure follows from three different routines for nitrogen (IX), oxygen (X), and carbon dioxide (XI). However, as indicated previously, Eqs (X) and (XI) must be considered together. The blood N2 content is directly proportional to nitrogen tension, as described by Henry’s Law: CciN2 = kN2*PciN2
(IX)
314
J. A. KELMAN
AND L. D. HOMER
The blood oxygen content is related to oxygen pressure through the oxy-hemoglobin dissociation curve, which in our calculations is approximated by the Hill Equation: P50 = PPRIM + 30.*(7.4 - pH) SAT = (PciO2/P50)**2.6/((PciO2/P50)**2.6 CciO2 = SAT*HBMX + k02*Pci02 The carbon dioxide content is computed and depends on bicarbonate concentration tration (the Haldane effect), as well as the and Bart, 1967). The method essentially is repeatedly until the calculated value on the value for CciCO2 used in the computation
+ 1.)
(Xa) (Xb)
(Xc) by solving an implicit function, and oxy-hemoglobin concenCO2 tension (Grodins, Buell, to compute the left side of XI left is within 0.1 percent of the on the right.
CciCO2 = BHC03 + 0.375*(HBMX - SAT*HBMX) - (0.16 + 2.3*HBMX)*(log[(CciC02 - kC02*PciC02)/(.01rPciCO2)] - 0.14) + kC02*PciC02
(XI)
Because of the Bohr effect in Eq. (Xa), pH becomes a prerequisite for the solution of Eqs. (X) and (XI). The pH of plasma may be calculated with the Henderson-Hasselbalch Equation Hions = K’*((kC02*PciC02)/(CciC02 pH = - log(Hions)
- kC02tPciC02))
(Xlla)
(Xllb) (Hions represents the concentration of hydrogen ions). As pH appears in Eq. (X), the Eqs. (X)-(X11) must all be solved as an implicit routine, which includes an independent implicit equation (XI). Once again, the method of solution is merely to solve Eqs. (X)-(X11) repeatedly. To complete the mathematical model, alveolar partial pressures are calculated from gas contents and total alveolar volume using the definition of partial pressure Paij = (Gij/Vai)*(PTO - 47.) (XIII) The set of differential equations is solved numerically by a RungeKutta method (Davis, 1960), which generates an updated value for each variable as time is incremented by about 1/15th of a respiratory cycle. Gas partial pressure in expired air, and gas tension in the arterial and venous-tissue compartments, are values that are not required for the solution of the set of differential equations. However, these variables are the familiar ways of expressing the composition of expired air, and the calculation is described below. The dead space gas fractions are converted to pressure by, Pmj = Fdj+(PTO In order to calculate
- 47.)
the blood gas tensions,
in the arterial
(XIV) and venous
MODEL OF EXTERNAL
375
RESPIRATION
compartments, Eqs. (X)-(X11) must again be used in an implicit loop. However, in this context the method of solution must be inverted to generate partial pressure from blood content of the gases, and in this use, both (X) and (XI) become implicit subroutines. The model, as presented here, was solved using a Fortran program with an IBM 360/67 computer. RESULTS
The results generated by this simulation of the pulmonary system, were obtained under stationary cycle conditions, and under transient conditions following some perturbation in the system at stationary cycle. The “normal” values used in the model are defined in Appendix II, except as noted in the text. The usual procedure is to set the input parameters (inspired gas fractions, ventilation rate, cardiac output, etc.), then run the simulation until the system approaches a stationary cycle (approximately 70 respiratory cycles, using artificially small volumes of distribution), yielding the values used in Fig. 3, 4, 5, and 6. The data presented in Figs. 2, 7, 8, 9, and 10 are obtained by instituting a timed respiratory maneuver after the stationary conditions are achieved, and after readjusting the volumes of distribution of the gases to normal. EXPIRED FIND RRTERIRL CFIRBBN DIOXIDE OURING NORMRL RESPIFlFiTICIN
TENSIONS
g s
FIG. 2. Expired carbon dioxide partial pressure and arterial carbon dioxide tension, during normal sinusoidal respiration.
376
J. A. KELMAN
AND L. D. HOMER
Figure 2 demonstrates the fluctuations in arterial carbon dioxide tension, and expired carbon dioxide partial pressure, that are due to the respiratory cycle. As shown, the smaller regular variations in arterial PCO, have a phase lag in comparison to the larger regular variations in expired PCOI. We have found that the difference between mean arterial PCO, and end tidal volume expired PCO, is highly sensitive to the magnitude of the tidal volume, at a constant total ventilation. STEADY
STRTE
\IS CRRDIFlC
6LpJl30
GRS
TENSIONS
OUTPUT
8 F
(rJ PRO2 Q PRCW
FIG. 3. Arterial and venous blood gas tensions values of cardiac output.
at stationary
h
PVO2
x
PVC02
cycle, for different
Figure 3, a plot of blood gas tension against cardiac output, illustrates the phenomena described by the Fick Equation (Eq. V). Assuming perfusion limited alveolar gas exchange (Eq. VIII), arterial PO2 and PC02 are unaffected by variation in cardiac output, and change occurs only in the venous blood gases. The gas content difference between arterial and venous blood (a-v difference) is proportional to cardiac output, while the gas tension difference is non-linear, as depicted, due to the non-linearity of gas-blood dissociation equations for oxygen and carbon dioxide. The feature of our system that is most different from other dynamic pulmonary models, is the inclusion of two alveolar compartments to allow ventilation/perfusion (V/Q) distribution studies. Imbalance in V/Q distribution is a major cause of clinical hypoxemia, resulting in less efficient gas exchange in pathological situations (West, 1965). This is illustrated
MODEL OF EXTERNAL
377
RESPIRATION
in Fig. 4, where blood gas tensions are plotted against increasing V/Q imbalance. The overall ventilation-perfusion ratio remains constant, as the total ventilation and the cardiac output are fixed. The scale on the abscissa refers to an increase in the percentage of total ventilation distributed to compartment 1, and to a reciprocal increase in the percentage of total perfusion (cardiac output) distributed to compartment 2. A value of 10 percent therefore, defines a situation in which 60 percent of the ventilation and 40 percent of the perfusion are directed to alveolus 1, while 40 percent of the ventilation and 60 percent of the perfusion are directed to alveolus 2. STEFlDY STATE BLCJ00 GRS TENSIONS ‘45 SHIFT IN V/Q OISTRIEWTION 8 frJ PA02 c) PflCO2 & PVi32 x
3o.m
PERCENT
%-IlFT
IN
V/Q
ART?km
FIG. 4. Arterial and venous blood gas tensions at stationary percent shifts in ventilation-perfusion distribution.
l-S.00
PVC02
wJ.cil
cycle, for different
The major effect of ventilation perfusion distribution imbalance appears in the arterial oxygen tension which falls 60 mm. of mercury over the range studied. Venous P02, which is a measure of oxygen delivered to the tissues, falls less than 14 mm. of mercury in this range, and arterial and venous carbon dioxide tensions show an increase of approximately 18 mm. of mercury from normal to the pathological state described as a 30 percent shift. For the values in Fig. 4, the mean volume of each alveolar compartment remains constant. As the cardiac output in Fig. 4 remains fixed at 4 liters per minute, the arterial-venous oxygen content difference also remains constant, from the 2.5
378
J. A. KELMAN
AND L. D. HOMER
Fick Equation. The non-linearity of the oxygen tension a-v difference, is due to the sigmoid shaped oxy-hemoglobin dissociation curve [Eq. (X)]. In Fig. 5, blood gas tensions are plotted against the volume of the anatomic dead space compartment, with fixed tidal volume and respiratory rate. As expected, arterial PO2 falls, while arterial and venous PC02 rise, as the dead space volume increases (decreasing alveolar ventilation). The venous P02, however, increases slightly, with increasing dead space. This phenomenon is due to the Bohr Effect, where the increasing venous PC02 lowers the pH and shifts the oxyhemoglobin curve to the right [Eq. (Xa)]. As a result, although the venous oxygen content falls as the dead space volume increases, oxygen becomes less tightly bound to hemoglobin, resulting in an elevated oxygen tension. STERDY STRTE BLCIOO GAS TENSIONS VS ANFlTClMIC DEAO SPfKE B 2
(TJPRO2
i:
FIG. 5. Arterial and venous blood gas tensions at stationary cycle, for different values of anatomic dead space.
The discovery that ventilation perfusion imbalance has a much greater effect on arterial oxygen tension than on arterial carbon dioxide tension; and that a decrease in alveolar ventilation has an approximately equal effect on arterial PO2 and PC02, is not unexpected. However, a mathematical model of the lung can be very valuable in assessing the relative contribution of changes in these two parameters to observed clinical pathology. Figure 6 is an example of the type of nomogram which can be
MODEL OF EXTERNAL
379
RESPIRATION
constructed with the lung model presented here, from stationary cycle data as described in Figs. 3-5. Following the general form of Olszowka and Farhi (1969), we have plotted arterial oxygen tension against arterial carbon dioxide tension for different preset values of alveolar ventilation (defined as percent of a standard Va : Vastd = 4725. ml/min), and different preset shifts in ventilation-perfusion distribution. The dashed lines are lines of constant V/Q distribution, and the solid lines represent constant alveolar ventilation. By entering patient arterial blood gases in this nomogram, values can be obtained differentiating between the two kinds of pathology. Figure 6, of course, is an oversimplification of a complicated process. The “normal” values can all be set for the individual patient; and in addition to the data collected for breathing room air, nomograms generated with other inspired gas mixtures can be used to improve the technique. FIRTERIRL LINES 8 Y
PW
VS
RRTERIFIL
PC02
OF CONSTF(NT RLVEOLFW VENT~LRTION
RN0 CONSTRNT VENTILRTION
PERFUSION
SHIFT
g d
: ui n
gB zj 2 % 5
\
\
-. -.
.
z ri Y
‘. .
20s Shl’+
z + %m
.%.w
30.00
JS.po PRCOP
w.m tnn no
-45.00
s0.m
11.00
s0.m
FIG. 6. Arterial carbon dioxide tension plotted against arterial oxygen tension, at stationary cycle. A nomogram is constructed with lines of equal alveolar ventilation (solid), and lines of constant shift in ventilation-perfusion distribution (dashed).
With on-line usage, this system could control tidal volume, frequency, and gas content of clinical respirators in order to maximize the efficiency of ventilation therapy. Figures 7 and 8 represent the simulation of respiratory maneuvers initiated after the steady cycle values described by the preceding figures
380
J. A. KELMAN
AND L. D. HOMER
have been achieved. In Fig. 7, expired PO2 and PC02 are plotted under conditions where rebreathing expired air is started at time 0, at the end of a normal inspiration. Figure 8 represents expired PO2 and PC02 during a timed prolonged expiration, begun after a larger than normal inspiration. The data in these figures were generated with the standard set of parameters described in Appendix II, i.e. cardiac output of 4 liters per minute, anatomic dead space of 175 ml., and equal ventilation and perfusion distribution to the alveolar compartments. To evaluate the effect of parameter changes on these expiratory curves, the respiratory maneuver simulation can be repeated for different sets of parameter values. PFIRTIOL
PRESSURES
DURING
REBRERTHING
,1/5 0.“7
FIG.
0.15
0.92
T I KE
OF
0.30
IMI N1JlES)
EXPIRED
0.37
GRS
0.45
rJ.52
0.60
MODEL OF EXTERNAL
381
RESPIRATION
PRRTIRL
PRESSURES
GURING
PROLONGED
OF
EXPIRED
GRS
EXPIRflTION
x
FIG. 8. Expired carbon dioxide and oxygen partial pressures during a prolonged expiration maneuver. The arrows in this figure, as in Fig. 7, represent the points at which the expired carbon dioxide pressure equals the arterial and the venous carbon dioxide tensions. EXPIRED DURING 8
STRNDRRO
OXYGEN PROL[1NGEO FlNO
PRRTIAL
PRESSURES
EXPIRATION
INCREMENTEO
CRROICI-PULNONARY
PRRRMETERS
z [3 STRNORRO a
CRROIRC
A OERO x
V/Q
VRLUES OUTPUT=
SPRCE= StilFT=
IOL/flIN
SOML 20
PERC.
FIG. 9. Expired oxygen partial pressure curves during prolonged expiration. The curves are constructed from the standard case by single parameter changes in cardiac output, anatomic dead space, and ventilation-perfusion distribution.
382
J. A. KELMAN
AND L. D. HOMER
from the expired gases, occur at different times for oxygen and carbon dioxide. The differences, in general, are detectable with present instrumentation. EXPIRED OURlNG D 9 D
STANORRD
CRRECid PROL0NGED AN0
DIC?XiOE
PflRTIRL
PRESSURES
EXPIFWTI0N
INCREMENTED
CRRDICI-PULMONRRY
PWiRiETERS
IRC
CUT?UT=
SHIFT=
20
IOL /,YilJ PERC.
FIG. 10. Expired carbon dioxide partial pressure curves during prolonged expiration, with single parameter changes as in Fig. 9. DISCUSSION
In the set of mathematical equations representing our analysis of the pulmonary system, there are several controversial areas. The question of perfusion-limited gas exchange at the alveolus as opposed to diffusion limited gas exchange, has never been conclusively answered. The published models to date have most commonly assumed the perfusion limited case, in which pulmonary venous blood gas tensions are identical to alveolar partial pressures (Yamamoto and Raub, 1967). We have used this assumption in Eq. (V), applying it to blood leaving the two alveolar compartments. An alternate approach to gas exchange across the alveolar-capillary membrane, would be a gradient equation where diffusion of gas j across alveolar capillary interface i, depends on the pressure gradient of the gas across the interface, and its permeability. This approach can be used to simulate pathological situations such as interstitial fibrosis and alveolar proteinosis, where a diffusion defect may contribute to the alveolararterial (A-a) differences in gas tension. In the model presented here,
MODEL OF EXTERNAL
RESPIRATION
383
including Eq. (V), A-a differences occur only because of V/Q imbalance. We are investigating the similarity between the effects of pathology in V/Q distribution and pathology in alveolar permeability. Pulmonary dead space includes three components, anatomic, stratified, and alveolar. The anatomic dead space consists of the respiratory airways which do not undergo gas exchange, alveolar dead space is a convention used to describe the effects of alveolar ventilation-perfusion imbalance, and the stratified dead space represents the area between the terminal airways and the alveoli which is incompletely mixed. The gas concentrations in the individual alveoli, and the anatomic dead space volume, are homogeneous from diffusion mixing, while the stratified component is characterized by gas concentration gradients (Strieder, Barnes, Levine, and Kazemi, 1970). Our lung model allows only for alveolar dead space as expressed by V/O shift, and anatomic dead space as a well-mixed compartment [Eq. (II)]. This defines the extreme case in which the stratified contribution is taken to be negligible. We have found this to be a satisfactory approximation, although a more advanced model could include stratified flow (Sikand, Cerretelli, and Farhi, 1966), as well as still gas diffusion and convective flow (Wagner, Latham, Brinkman, and Filley, 1969). The treatment of V/Q distribution can also have several interpretations. While the distribution of perfusion is straightforward, the distribution of ventilation is associated with the choice of resting lung volume. We have kept the mean volume in each alveolar compartment constant and equal, with fluctuations around this mean based on FVi, the fraction of ventilation directed towards each alveolus i. This method is equivalent to decreasing the functional residual capacity (FRC) of alveolus i, as its ventilation increases. The phenomenon of sequential emptying of unequally ventilated alveoli is not a part of this model, as it has been shown that the relative contribution of different zones of the lung to expiratory flow remains constant (Farhi, 1969). The results presented in this paper have two functions, the stationary cycle data represents a simulation of pulmonary-cardiovascular pathology, while the transient state data represents a first attempt at parameter estimation on the basis of pulmonary function studies. Further research in the area of parameter estimation is possible with these results. For the measurement of cardiac output, an established bloodless (i.e. pulmonary) technique is that of Kim, Rahn, and Farhi (1966) as applied by Bickel et al. (1970). This method involves the estimation of arterial and venous PC02 from expired air gases on prolonged expiration. Using our model, we have simulated the respiratory maneuver described by Kim, Rahn, and Farhi (inhalation of a tidal volume larger than normal followed by steady exhalation for 20 to 30 sec.), and have used their methods of numerical
384
J. A. KELMAN
AND L. D. HOMER
analysis to calculate cardiac output, for different cases of cardio-pulmonary pathology. The estimated cardiac outputs come within 20 percent of the value preset in the programmed simulation, and even greater accuracy is possible (10 percent of preset values), with other choices of mathematical operations to convert PC02 to CO2 content. The successful application of this technique is indicative of internal consistency in the model, and also suggests a generalized approach to the problem of parameter estimation. Accuracy in Kim’s estimation procedure is limited to the degree that changes in cardiac output (or any other parameter being studied) will produce changes in recorded variables. There is the further limitation of ambiguity, as several combinations of cardiac output, dead space volume, and ventilation perfusion distribution may be capable of producing a single expired PO2 or PC02 value during prolonged expiration. Kim, Rahn, and Farhi use expired oxygen and carbon dioxide simultaneously to decrease this ambiguity, because of the different behavior of the two gases during expiration, as described in connection with Fig. 9 and 10. Data from the oxygen curve is used to calculate instantaneous respiratory quotient, the value of which determines the choice of the two points on the carbon dioxide curve which best correspond to the arterial and venous CO2 tensions (indicated in Fig. 8). Calculations based on two points from an expiratory flow curve, however, cannot be as accurate, or as free of ambiguity, as the use of the entire curve. Therefore, a logical extension of the techniques previously mentioned, is parameter estimation based on curve fitting of expired gas data with values produced by computer simulation. For purposes of parameter estimation, the lung model can be redefined in terms of a parameter set Bj, j = 1, k, to be estimated; a set of observed independent variables Xg, g = 1, p; and a set of dependent variables Yh, h = 1, m, to be observed. Tn this notation, B would include cardiac output, dead space volume, and V/Q distribution; X would include tidal volume, respiratory rate, gas flow, etc., * and Y would consist of observed Fd02, FdC02, and FdN2 at the mouth. If Yhbar is the value of expired partial fraction of gas h calculated by the model, then: Yhbar = fh(X1, X2,. . . Xp, Bl, B2,. . ., Bk) By generating a series of Yhbari, i = 1, n (n increments over time for some simulated respiratory maneuver), these values can be compared to the set of observed Y for each gas. Using programmed algorithms for leastsquares estimation of non-linear parameters (Horwitz and Homer, 1970; Marquardt, 1963) the set of Bj can be altered until a best fit of the model generated data (Ybar) to the observed data (Y) is achieved. The computer generated data can be simultaneously fitted to several observed variables, i.e. P02, PC02, and PN2, to increase the precision of the estimation and
MODEL
OF EXTERNAL
RESPIRATION
385
reduce the interdependence between the parameters. We can define a matrix (X), in which each element element Xij is the derivative of Ybar (for a single gas) with respect to parameter j, at time increment i: j = 1,k. i = 1,n; (X)nxk = (dYbari/dBj) If these initial calculations are generalized for three gases, the cumulative (X’X) inverse matrix can be constructed and normalized. These matrices, for a rebreathing maneuver, are presented in Table I. TABLE 1 PARAMETER TOMIC DEAD
21s ANA3 IS VENTILATION-
1 IS CARDIAC OUTPUT,PARAMETER
SPACE, AND
PARAMETER
PERFUSION DISTRIBUTION SHIFT(PERCENT)
Cumulative (X’X) 0.514E 09 - 0.394E 07 - 0.207E 05
inverse matrix: 0.394E 07 0.136E 07 0.104E 04
- 0.207E 05 - 0.104E 04 0.309E 01
0.149E 00 O.lOOE 01 - 0.505E 00
- 0.520E 00 - 0.505E 00 O.lOOE 01
-
Normalized: O.lOOE 01 - 0.149E 00 - 0.520E 00
-
The (X’X) inverse matrix can provide information on two aspects of parameter estimation, variance, and co-variance. The diagonal elements of this matrix are indications of the error expected in the parameter estimation. If we can simulate observed data with model generated data to a standard deviation of less than two percent, which is a realistic goal for this kind of curve fitting, the cardiac output can be estimated to within 12 percent, the dead space can be estimated to within 14 percent, and the ventilation perfusion shift can be estimated to within 8 percent. However, these values are based on the assumption that the estimation of each parameter is completely independent. This is obviously incorrect, and a measure of the interdependence of the estimations is given by the off-diagonal elements in the normalized (X’X) inverse matrix. For this maneuver (rebreathing), the co-variance between cardiac output and V/Q distribution is approximately 50 percent, while the co-variance between cardiac output and dead space is less than 15 percent. Therefore, 50 percent of the error in estimating cardiac output is due to error in the estimation of V/Q distribution, and 15 percent due to the dead estimation error. If these calculations are repeated for a prolonged expiratory maneuver, the standard errors are increased
386
J. A. KELMAN
AND L. D. HOMER
for the three parameters, while the total co-variance (sum of individual co-variances) decreases. Therefore, a combination of respiratory maneuvers may eventually provide the best estimation procedure. This system could be implemented to give on-line results, if analog to digital conversion of patient data were used. We hope to report on experimental confirmation of our model in the near future. APPENDIX I
Character
Subscript
BHC03 co c
i
j (gas)
ar
Cl, c2 V D
i
FV
i
FQB
i
F
i
j
192 132
192 j (gas)
al, a2 cl
G
I
J
(gas> I,2 HBMX
Definition
Units
standard bicarbonate content cardiac output concentration of gas j
ml CO,/ ml blood 37” ml/min ml gas/ml blood
in the arterial compartment in pulmonary capillary compartment 1 and 2 in the mixed venous compartment gas diffusion of gas j from alveolus 1 and 2 ventilation distribution fraction to alveolus 1 and 2 perfusion distribution fraction to alveolus 1 and 2 gas fraction of gas j in alveolus 1 and 2 in the dead space compartment in the inspired air gas content of gas j in alveolus 1 and 2 maximum binding potential of 0, to Hb in blood
ml gas/min
dimensionless
dimensionless
dimensionless
ml
ml 02/ml
blood
MODEL OF EXTERNAL
Character
RESPIRATION
Subscript
Hions K
_i
(gas) K’ PERM
j (gas)
P 50
P Prim PTO P
i
j (gas)
al, a2 CI, c2 m RATE SAT t VDOT
i 192
i aI, a2 art d VDIS
j
(gas)
387 Definition
hydrogen ion concentration solubility coefficient for gas j in blood dissociation content for carbonic acid permeability coefficient of gas j the 0, tension at 50% saturation after correction for pH the P,, at pH = 7.4 the ambient pressure tension or partial pressure of gas j in alveolus 1 and 2 in pulmonary capillary compartment 1 and 2 at the month respiratory rate hemoglobin saturation time gas flow to alveolus 1 and 2 (uncorrected for diffusion across the alveolus capillary interface) volume of alveolus compartment 1 and 2 of arterial blood compartment of dead space compartment volume of distribution of gas j
Units moles/liter ml gas (STP)/ ml blood moles/liter ml gas
mm-mm Hg mm Hg
mm Hg mm Hg mm Hg
cycles/min percent min ml/min
ml
ml
388 APPENDIX
Normal
J. A. KELMAN
AND L. D. HOMER
II
parameter
values :
BHC03 = 0.547 co = 4000. Fo02 = 0.21 FoC02 = 0. FoN2 = 0.79 HBMX = 0.201 k02 = 0.0000282 kC02 = 0.0006732 kN2 = 0.0000172 K’ = 7.95*10**(7) METAB = 250. METABC02 = 200. PPRIM = 26.5 (National Research Council, 1944) PTO = 760. RATE = 14. tidal volume = 400. ml. avg. Val = avg. Va2 = 1300. (therefore the FRC = 2400.) Vart = 50. VD = 175. VDIS02 = 5000. VDISC02 = 21000. VDTSN2 = 42000.
REFERENCES 1 R. G. Bickel, C. F. Diener, and H. L. Brammell, Aerospace Med. 41, 203 (1970). C. G. Douglas, and J. S. Haldane, J. Physiol. (London) 48,244 (1914). 3 J. H. Comroe, Physiology of respiration, Yearbook Medical Publishers, Chicago (1965), pp. 157-181. 4 H. T. Davis, Introductio,l to non-linear differential and integral equations, United States Atomic Energy Commission (1960), pp. 482-486. 5 L. E. Farhi, in Circulatory and respiratory mass transport (Wolstenholme and Knight, eds.), Little, Brown, and Co., Boston (1969), pp. 277-297. 6 J. S. Gray, The multiple factor theory of respiratory regulation. Army Air Force
2 J. Christiansen,
School of Aviation Medicine Proj. Rept. No. 386 (I, 2, 3; 1945). 7 F. S. Grodins, J. S. Gray, K. R. Schroeder, A. L. Norins and R. W. Jones, J. Appl. Physiol. 7, 283 (1954). 8 F. S. Grodins, J. Buell, and A. J. Bart, J. Appl. Physiol. 22, 260 (1967). 9 D. L. Horwitz and L. D. Homer, Analysis of Biomedical Data by Time-sharing computers, 1, Non-linear regression analysis, Naval Medicine Research Institute
Research Rept. No. 25 (1970).
MODEL OF EXTERNAL
RESPIRATION
389
10 T. S. Kim, H. Rahn, and L. E. Farhi, J. Appl. PhysioL 21, 1338 (1966). 11 D. B. Marquardt, .Z. Sot. In&f. Math. 11, 431 (1963). 12 National Research Council, Oxygen dissociation curves for human blood at 37 degrees Centigrade, in The handbook of respiration data in aviation, Washington D.C. (1944). 13 A. J. Olszowka and L. E. Farhi, J. Appl. Physiol. 26, 141 (1969). 14 R. Sikand, P. Cerretelli, and L. E. Farhi, .Z. Appl. Physiol. 21, 1331 (1966). 15 D. J. Strieder, B. A. Barnes, B. W. Levine, and H. Kazemi, .Z. Appt. Physiol. 29,486 (1970). 16 W. W. Wagner, Jr., L. P. Latham,
P. D. Brinkman, and G. F. Filley, Science 163, 1210 (1969). 17 J. B. West, Ventilution/bloodflow andgus exchange, Blackwell Scientific Publications, Oxford (1965), pp. 14-15. 18 W. S. Yamamoto and W. F. Raub, Computers and Biom. Res. 1, 65 (1967).
BIOSCIENCES
13, 369-389
(1972)
369
The Analysis of Pulmonary and Cardiovascular Parameters with a Mathematical Model of External Respiration * J. A. KELMAN
AND
L. D. HOMER?
Division of Biological and Medical Sciences Brown University, Providence, Rhode Island 02912
Communicated
by F. Grodins
ABSTRACT The behavior of the external respiratory system was described by a set of differential equations. The system was solved by digital computer, to express the relationship between stationary cycle blood gas tensions and fixed cardiac and respiratory parameters. Dynamic respiratory maneuvers were simulated for different values of these parameters to investigate the feasibility of estimating cardiac output, respiratory dead space, and ventilation-perfusion distribution, from flow and composition of expired air, using non-linear regression techniques.
INTRODUCTION
Christiansen, Douglas, and Haldane (1914) made the original attempt to estimate blood gas parameters with a bloodless technique utilizing respiratory data. Kim, Rahn, and Farhi (1966) have developed a more sophisticated method of estimating venous and arterial PCOZ from respiratory maneuvers, and recently Bickel, Diener, and Brammell (1970) have devised a computerized data gathering system to estimate cardiac output using this technique and the Fick Equation. The physiological basis for these and other studies in this area, is the rapid response to changes in cardiac output, of the partial pressure of gases at the mouth. The qualitative events are often obvious, and the major problems in the investigations involve the quantitative predictions. The accuracy of these predictions may be expected to depend on the exactness of the assumptions on which they are based, as well as the amount and accuracy of the data analyzed. Any set of assumptions contains the concept of a model, and their validity is limited to the extent that the model can accurately simulate the real system. It is our intention to apply the techniques developed for * This work supported in part by ONR Contract NOOO14-67-A-0191-0013. t Reprint requests should be addressed to: L. D. Homer, Naval Medical Research Institute National Naval Medical Center, Bethesda, Md. 20014 U.S.A. Copyright 0 1972 by American Elsevier Publishing Company, Inc.
370
J. A. KELMAN
AND L. D. HOMER
pulmonary control system modeling, to the analysis and estimation of pulmonary and cardiovascular parameters. The first mathematical model for mammalian respiration was proposed by Gray (1945). It consisted of a system of equations which could be solved algebraically for steady state values. An advanced model for the determination of steady state respiratory information has been devised by West (1969), including complex non-differential equations describing blood gas dissociation and ventilation-perfusion phenomena, which are solved by numerical means. Grodins, Gray, Schroeder, Norins, and Jones (1954) published the first mathematical description of a respiratory model that allowed the solution of transient values, using an analog computer to predict respiratory parameters arising in response to changes in carbon dioxide inhalation. Grodins’ model has since served as the basis for many other analyses of the respiratory system, involving analog and digital simulation, which are discussed in a recent review by Yamamoto and Raub (1967). Many of the ideas in our model are based on Grodins’ work, although applied towards a different goal. In the bulk of the respiratory models published to date, the major emphasis has been placed on the “closed loop” system, with the implementation of different “ventilation equations” in an effort to simulate respiratory feedback control. We are interested here, in an open-loop respiratory model with no internal control. All the parameters can be set independently, so that the quantitative effect of changes in a single parameter, on respiratory function, can be isolated and evaluated independently of changes in the other parameters of the system. By this kind of simulation, we can specify those non-pulmonary parameters that it is possible to estimate by the measurement of expired air, and suggest numerical procedures for doing so. This publication includes the system of equations comprising our model, the methods of digital solution of these equations, and some preliminary results indicating the quantitative changes in blood gas tensions due to changes in cardiac output, ventilation, and ventilation-perfusion distribution.
Tl IEORY
In our analysis of gas exchange, the numerical model constructed contains ten homogenous compartments (Fig. 1). The gas flow between a single dead space compartment and two alveolar compartments is controlled by forced variation in alveolar volume. Blood circulation is a closed loop which involves two pulmonary capillary compartments feeding into a single arterial compartment. The arterial blood enters the venous blood
MODEL OF EXTERNAL
RESPIRATION
371
compartment, where blood can equilibrate with the tissue space, and is then returned to the pulmonary capillary beds. The sizes of the three tissue compartments are equal to the volumes of distribution of oxygen, nitrogen, and carbon dioxide in turn, and the blood returning to the lungs is assumed to have the same content as each component has in its tissue space. The values used for the volumes of distribution may not be exact, but have negligible effects on the data we are interested in with the time periods under consideration. As this system is perfusion-limited at the alveolus as well as in the peripheral circulation, blood leaving each set of pulmonary capillaries is assumed to have the same partial pressure as its associated alveolar compartment.
FIG. 1. Flow diagram representing compartmental analysis.
the pulmonary-cardiovascular
model in terms of
The nitrogen tension is related to gas content by Henry’s Law, as it involves only dissolved gas; while oxygen and carbon dioxide tensions are related to content by more complex equations. The oxygen content includes dissolved gas as well as oxygen bound to hemoglobin, while carbon dioxide is present in blood as dissolved gas, as serum bicarbonate, and as carbamino compounds bound to blood protein. The distribution of CO2 is also dependent on hemoglobin content and saturation. The pH of the
372
J. A. KELMAN
AND L. D. HOMER
blood, in each circulatory compartment, follows the Henderson-Hasselbalch equation depending on PCOz and blood bicarbonate (Comroe, 1965). By including the effect of pH on the oxyhemoglobin dissociation curve (the Bohr effect) oxygen content, carbon dioxide content, and pH must be solved for simultaneously in an implicit routine. The set of differential equations that mathematically define the above relationships, is presented here as an ordered set, following the logical sequence suggested by the flow diagram. The equations are written using Fortran conventions for arithmetic. The symbols i, j, 1,2,02, C02, and N2 at the end of a variable name are to be read as subscripts. The subscript i refers to alveolus 1 or 2, and the subscript j refers to oxygen, carbon dioxide, or nitrogen. All the symbols used in the equations are defined in Appendix I. We define gas flow to alveolar compartment i from the dead space, as VDOTi, positive on inspiration. This value is determined by the respiration rate (RATE) and the tidal volume (amp = tidal volume/2.), as well as the distribution of flow term FVi (the fraction of gas flow directed to alveolus i). VDOTi is also adjusted for gas exchange in the alveolus, with the diffusion terms Dij. The quantity VDOTi must be adjusted for gas exchange across the alveolar membrane, and is therefore not identical to the rate of change of alveolar volume i [d(Vai)/dt in Eq. IV]. In effect, Eq. 1 adjusts gas flow at the mouth and at the main bifurcation of the pulmonary airways, so that the alveolar volumes each fluctuates around a fixed value (the “resting volume”). The degree of fluctuation is determined by the fraction of flow to each alveolus. Normal respiratory gas flow, at time t, is approximated with a fixed sine wave forcing function: VDOTi = FVi*2*pi*RATE*amp*sin(2*pi*t*RATE) + [(273. + 37.)/273.]*(760./PTO)*(Di02 + DiC02
+ DiN2)
(I)
Two equations are used to calculate the dead space composition, one for inhalation and one for exhalation. The fraction of gas j in the dead space is Fdj, Foj refers to gas j fraction in the inspired air, and Faij represents gas j fraction in alveolar volume i. VD is the dead space volume. Then for inhalation : d(Fdj)/dt
= (VDOTl
+ VDOT2)*(Foj
- Fdj)/VD
(Ira)
and for exhalation: d(Fdj)/dt
= VDOTl*(Fdj
- FaIj)/VD + VDOT2*(Fdj
- Fa2j)/VD
(IIb)
Separate inhalation and exhalation equations are used also in determining alveolar composition. Gij refers to the quantity of gas j in alveolus i at BTP, rather than gas fraction. PTO is the ambient pressure, and Dij
MODEL OF EXTERNAL
RESPIRATION
373
the alveolar capillary transfer of gas j in alveolus converted to BTP in (IIIa) and (IIIb): d(Gij)/dt d(Gij)/dt
= VDOTi*Fdj = VDOTi*Faij
- [(273. + 37.)/273]*(760./PTO)*Dij
(IIIa)
- [(273. + 37.)/273.]*(760./PTO)*Dij
(IIIb)
The change in total alveolar Equation III over j gases: d(Vai)/dt
= VDOTi
i, in ml. at STP which is
volume
(Vai) is determined
by summing
- ((273. + 37.)/273.)*(760./PTO)r(Di02 + DiC02 + DiN2)
(IV)
or what is the same, d(Vai)/dt
= FVi*2*pi*RATE*amp*sin(2*pi*RATE*t)
Diffusion (D), in ml/min at STP, is calculated from cardiac output (CO), fraction of flow to alveolar compartment i (FBQi), end pulmonary capillary concentration (Cc), and venous return concentration (Cv) : Dij = CO*FQBi*(Ccij
- Cvj)
(V)
Arterial content of gas j (Carj) is altered by influx from two capillary compartments, and outflow into one venous-tissue compartment. Vart refers to arterial compartment volume : d(Carj)/dt
= (l/Vart)*(CO*FQBl*Cclj + CO*FQB2*Cc2j
- CO*Carj)
(VI)
In the same manner, venous-tissue gas contents are calculated from arterial inflow, outflow to the pulmonary capillaries, and metabolic flow (METAB), which represents gas production and consumption in tissue metabolism. METABC02 divided by METAB defines the overall respiratory quotient. Each gas j, has its own volume of distribution (VDISj) : d(Cvj)/dt
= (l/VDISj)*(CO*Carj
- CO*Cvj
The relationship between the partial pressure in the blood presents the most complex problems of this model. Blood gas tension (Pcj) is set in compartment equal to the partial pressure in the Pcij = Paij
- METABj)
(VII)
and the content of gases in the numerical solution the pulmonary capillary alveolus (Paij): (VIII)
The calculation of gas content in the blood from the partial pressure follows from three different routines for nitrogen (IX), oxygen (X), and carbon dioxide (XI). However, as indicated previously, Eqs (X) and (XI) must be considered together. The blood N2 content is directly proportional to nitrogen tension, as described by Henry’s Law: CciN2 = kN2*PciN2
(IX)
314
J. A. KELMAN
AND L. D. HOMER
The blood oxygen content is related to oxygen pressure through the oxy-hemoglobin dissociation curve, which in our calculations is approximated by the Hill Equation: P50 = PPRIM + 30.*(7.4 - pH) SAT = (PciO2/P50)**2.6/((PciO2/P50)**2.6 CciO2 = SAT*HBMX + k02*Pci02 The carbon dioxide content is computed and depends on bicarbonate concentration tration (the Haldane effect), as well as the and Bart, 1967). The method essentially is repeatedly until the calculated value on the value for CciCO2 used in the computation
+ 1.)
(Xa) (Xb)
(Xc) by solving an implicit function, and oxy-hemoglobin concenCO2 tension (Grodins, Buell, to compute the left side of XI left is within 0.1 percent of the on the right.
CciCO2 = BHC03 + 0.375*(HBMX - SAT*HBMX) - (0.16 + 2.3*HBMX)*(log[(CciC02 - kC02*PciC02)/(.01rPciCO2)] - 0.14) + kC02*PciC02
(XI)
Because of the Bohr effect in Eq. (Xa), pH becomes a prerequisite for the solution of Eqs. (X) and (XI). The pH of plasma may be calculated with the Henderson-Hasselbalch Equation Hions = K’*((kC02*PciC02)/(CciC02 pH = - log(Hions)
- kC02tPciC02))
(Xlla)
(Xllb) (Hions represents the concentration of hydrogen ions). As pH appears in Eq. (X), the Eqs. (X)-(X11) must all be solved as an implicit routine, which includes an independent implicit equation (XI). Once again, the method of solution is merely to solve Eqs. (X)-(X11) repeatedly. To complete the mathematical model, alveolar partial pressures are calculated from gas contents and total alveolar volume using the definition of partial pressure Paij = (Gij/Vai)*(PTO - 47.) (XIII) The set of differential equations is solved numerically by a RungeKutta method (Davis, 1960), which generates an updated value for each variable as time is incremented by about 1/15th of a respiratory cycle. Gas partial pressure in expired air, and gas tension in the arterial and venous-tissue compartments, are values that are not required for the solution of the set of differential equations. However, these variables are the familiar ways of expressing the composition of expired air, and the calculation is described below. The dead space gas fractions are converted to pressure by, Pmj = Fdj+(PTO In order to calculate
- 47.)
the blood gas tensions,
in the arterial
(XIV) and venous
MODEL OF EXTERNAL
375
RESPIRATION
compartments, Eqs. (X)-(X11) must again be used in an implicit loop. However, in this context the method of solution must be inverted to generate partial pressure from blood content of the gases, and in this use, both (X) and (XI) become implicit subroutines. The model, as presented here, was solved using a Fortran program with an IBM 360/67 computer. RESULTS
The results generated by this simulation of the pulmonary system, were obtained under stationary cycle conditions, and under transient conditions following some perturbation in the system at stationary cycle. The “normal” values used in the model are defined in Appendix II, except as noted in the text. The usual procedure is to set the input parameters (inspired gas fractions, ventilation rate, cardiac output, etc.), then run the simulation until the system approaches a stationary cycle (approximately 70 respiratory cycles, using artificially small volumes of distribution), yielding the values used in Fig. 3, 4, 5, and 6. The data presented in Figs. 2, 7, 8, 9, and 10 are obtained by instituting a timed respiratory maneuver after the stationary conditions are achieved, and after readjusting the volumes of distribution of the gases to normal. EXPIRED FIND RRTERIRL CFIRBBN DIOXIDE OURING NORMRL RESPIFlFiTICIN
TENSIONS
g s
FIG. 2. Expired carbon dioxide partial pressure and arterial carbon dioxide tension, during normal sinusoidal respiration.
376
J. A. KELMAN
AND L. D. HOMER
Figure 2 demonstrates the fluctuations in arterial carbon dioxide tension, and expired carbon dioxide partial pressure, that are due to the respiratory cycle. As shown, the smaller regular variations in arterial PCO, have a phase lag in comparison to the larger regular variations in expired PCOI. We have found that the difference between mean arterial PCO, and end tidal volume expired PCO, is highly sensitive to the magnitude of the tidal volume, at a constant total ventilation. STEADY
STRTE
\IS CRRDIFlC
6LpJl30
GRS
TENSIONS
OUTPUT
8 F
(rJ PRO2 Q PRCW
FIG. 3. Arterial and venous blood gas tensions values of cardiac output.
at stationary
h
PVO2
x
PVC02
cycle, for different
Figure 3, a plot of blood gas tension against cardiac output, illustrates the phenomena described by the Fick Equation (Eq. V). Assuming perfusion limited alveolar gas exchange (Eq. VIII), arterial PO2 and PC02 are unaffected by variation in cardiac output, and change occurs only in the venous blood gases. The gas content difference between arterial and venous blood (a-v difference) is proportional to cardiac output, while the gas tension difference is non-linear, as depicted, due to the non-linearity of gas-blood dissociation equations for oxygen and carbon dioxide. The feature of our system that is most different from other dynamic pulmonary models, is the inclusion of two alveolar compartments to allow ventilation/perfusion (V/Q) distribution studies. Imbalance in V/Q distribution is a major cause of clinical hypoxemia, resulting in less efficient gas exchange in pathological situations (West, 1965). This is illustrated
MODEL OF EXTERNAL
377
RESPIRATION
in Fig. 4, where blood gas tensions are plotted against increasing V/Q imbalance. The overall ventilation-perfusion ratio remains constant, as the total ventilation and the cardiac output are fixed. The scale on the abscissa refers to an increase in the percentage of total ventilation distributed to compartment 1, and to a reciprocal increase in the percentage of total perfusion (cardiac output) distributed to compartment 2. A value of 10 percent therefore, defines a situation in which 60 percent of the ventilation and 40 percent of the perfusion are directed to alveolus 1, while 40 percent of the ventilation and 60 percent of the perfusion are directed to alveolus 2. STEFlDY STATE BLCJ00 GRS TENSIONS ‘45 SHIFT IN V/Q OISTRIEWTION 8 frJ PA02 c) PflCO2 & PVi32 x
3o.m
PERCENT
%-IlFT
IN
V/Q
ART?km
FIG. 4. Arterial and venous blood gas tensions at stationary percent shifts in ventilation-perfusion distribution.
l-S.00
PVC02
wJ.cil
cycle, for different
The major effect of ventilation perfusion distribution imbalance appears in the arterial oxygen tension which falls 60 mm. of mercury over the range studied. Venous P02, which is a measure of oxygen delivered to the tissues, falls less than 14 mm. of mercury in this range, and arterial and venous carbon dioxide tensions show an increase of approximately 18 mm. of mercury from normal to the pathological state described as a 30 percent shift. For the values in Fig. 4, the mean volume of each alveolar compartment remains constant. As the cardiac output in Fig. 4 remains fixed at 4 liters per minute, the arterial-venous oxygen content difference also remains constant, from the 2.5
378
J. A. KELMAN
AND L. D. HOMER
Fick Equation. The non-linearity of the oxygen tension a-v difference, is due to the sigmoid shaped oxy-hemoglobin dissociation curve [Eq. (X)]. In Fig. 5, blood gas tensions are plotted against the volume of the anatomic dead space compartment, with fixed tidal volume and respiratory rate. As expected, arterial PO2 falls, while arterial and venous PC02 rise, as the dead space volume increases (decreasing alveolar ventilation). The venous P02, however, increases slightly, with increasing dead space. This phenomenon is due to the Bohr Effect, where the increasing venous PC02 lowers the pH and shifts the oxyhemoglobin curve to the right [Eq. (Xa)]. As a result, although the venous oxygen content falls as the dead space volume increases, oxygen becomes less tightly bound to hemoglobin, resulting in an elevated oxygen tension. STERDY STRTE BLCIOO GAS TENSIONS VS ANFlTClMIC DEAO SPfKE B 2
(TJPRO2
i:
FIG. 5. Arterial and venous blood gas tensions at stationary cycle, for different values of anatomic dead space.
The discovery that ventilation perfusion imbalance has a much greater effect on arterial oxygen tension than on arterial carbon dioxide tension; and that a decrease in alveolar ventilation has an approximately equal effect on arterial PO2 and PC02, is not unexpected. However, a mathematical model of the lung can be very valuable in assessing the relative contribution of changes in these two parameters to observed clinical pathology. Figure 6 is an example of the type of nomogram which can be
MODEL OF EXTERNAL
379
RESPIRATION
constructed with the lung model presented here, from stationary cycle data as described in Figs. 3-5. Following the general form of Olszowka and Farhi (1969), we have plotted arterial oxygen tension against arterial carbon dioxide tension for different preset values of alveolar ventilation (defined as percent of a standard Va : Vastd = 4725. ml/min), and different preset shifts in ventilation-perfusion distribution. The dashed lines are lines of constant V/Q distribution, and the solid lines represent constant alveolar ventilation. By entering patient arterial blood gases in this nomogram, values can be obtained differentiating between the two kinds of pathology. Figure 6, of course, is an oversimplification of a complicated process. The “normal” values can all be set for the individual patient; and in addition to the data collected for breathing room air, nomograms generated with other inspired gas mixtures can be used to improve the technique. FIRTERIRL LINES 8 Y
PW
VS
RRTERIFIL
PC02
OF CONSTF(NT RLVEOLFW VENT~LRTION
RN0 CONSTRNT VENTILRTION
PERFUSION
SHIFT
g d
: ui n
gB zj 2 % 5
\
\
-. -.
.
z ri Y
‘. .
20s Shl’+
z + %m
.%.w
30.00
JS.po PRCOP
w.m tnn no
-45.00
s0.m
11.00
s0.m
FIG. 6. Arterial carbon dioxide tension plotted against arterial oxygen tension, at stationary cycle. A nomogram is constructed with lines of equal alveolar ventilation (solid), and lines of constant shift in ventilation-perfusion distribution (dashed).
With on-line usage, this system could control tidal volume, frequency, and gas content of clinical respirators in order to maximize the efficiency of ventilation therapy. Figures 7 and 8 represent the simulation of respiratory maneuvers initiated after the steady cycle values described by the preceding figures
380
J. A. KELMAN
AND L. D. HOMER
have been achieved. In Fig. 7, expired PO2 and PC02 are plotted under conditions where rebreathing expired air is started at time 0, at the end of a normal inspiration. Figure 8 represents expired PO2 and PC02 during a timed prolonged expiration, begun after a larger than normal inspiration. The data in these figures were generated with the standard set of parameters described in Appendix II, i.e. cardiac output of 4 liters per minute, anatomic dead space of 175 ml., and equal ventilation and perfusion distribution to the alveolar compartments. To evaluate the effect of parameter changes on these expiratory curves, the respiratory maneuver simulation can be repeated for different sets of parameter values. PFIRTIOL
PRESSURES
DURING
REBRERTHING
,1/5 0.“7
FIG.
0.15
0.92
T I KE
OF
0.30
IMI N1JlES)
EXPIRED
0.37
GRS
0.45
rJ.52
0.60
MODEL OF EXTERNAL
381
RESPIRATION
PRRTIRL
PRESSURES
GURING
PROLONGED
OF
EXPIRED
GRS
EXPIRflTION
x
FIG. 8. Expired carbon dioxide and oxygen partial pressures during a prolonged expiration maneuver. The arrows in this figure, as in Fig. 7, represent the points at which the expired carbon dioxide pressure equals the arterial and the venous carbon dioxide tensions. EXPIRED DURING 8
STRNDRRO
OXYGEN PROL[1NGEO FlNO
PRRTIAL
PRESSURES
EXPIRATION
INCREMENTEO
CRROICI-PULNONARY
PRRRMETERS
z [3 STRNORRO a
CRROIRC
A OERO x
V/Q
VRLUES OUTPUT=
SPRCE= StilFT=
IOL/flIN
SOML 20
PERC.
FIG. 9. Expired oxygen partial pressure curves during prolonged expiration. The curves are constructed from the standard case by single parameter changes in cardiac output, anatomic dead space, and ventilation-perfusion distribution.
382
J. A. KELMAN
AND L. D. HOMER
from the expired gases, occur at different times for oxygen and carbon dioxide. The differences, in general, are detectable with present instrumentation. EXPIRED OURlNG D 9 D
STANORRD
CRRECid PROL0NGED AN0
DIC?XiOE
PflRTIRL
PRESSURES
EXPIFWTI0N
INCREMENTED
CRRDICI-PULMONRRY
PWiRiETERS
IRC
CUT?UT=
SHIFT=
20
IOL /,YilJ PERC.
FIG. 10. Expired carbon dioxide partial pressure curves during prolonged expiration, with single parameter changes as in Fig. 9. DISCUSSION
In the set of mathematical equations representing our analysis of the pulmonary system, there are several controversial areas. The question of perfusion-limited gas exchange at the alveolus as opposed to diffusion limited gas exchange, has never been conclusively answered. The published models to date have most commonly assumed the perfusion limited case, in which pulmonary venous blood gas tensions are identical to alveolar partial pressures (Yamamoto and Raub, 1967). We have used this assumption in Eq. (V), applying it to blood leaving the two alveolar compartments. An alternate approach to gas exchange across the alveolar-capillary membrane, would be a gradient equation where diffusion of gas j across alveolar capillary interface i, depends on the pressure gradient of the gas across the interface, and its permeability. This approach can be used to simulate pathological situations such as interstitial fibrosis and alveolar proteinosis, where a diffusion defect may contribute to the alveolararterial (A-a) differences in gas tension. In the model presented here,
MODEL OF EXTERNAL
RESPIRATION
383
including Eq. (V), A-a differences occur only because of V/Q imbalance. We are investigating the similarity between the effects of pathology in V/Q distribution and pathology in alveolar permeability. Pulmonary dead space includes three components, anatomic, stratified, and alveolar. The anatomic dead space consists of the respiratory airways which do not undergo gas exchange, alveolar dead space is a convention used to describe the effects of alveolar ventilation-perfusion imbalance, and the stratified dead space represents the area between the terminal airways and the alveoli which is incompletely mixed. The gas concentrations in the individual alveoli, and the anatomic dead space volume, are homogeneous from diffusion mixing, while the stratified component is characterized by gas concentration gradients (Strieder, Barnes, Levine, and Kazemi, 1970). Our lung model allows only for alveolar dead space as expressed by V/O shift, and anatomic dead space as a well-mixed compartment [Eq. (II)]. This defines the extreme case in which the stratified contribution is taken to be negligible. We have found this to be a satisfactory approximation, although a more advanced model could include stratified flow (Sikand, Cerretelli, and Farhi, 1966), as well as still gas diffusion and convective flow (Wagner, Latham, Brinkman, and Filley, 1969). The treatment of V/Q distribution can also have several interpretations. While the distribution of perfusion is straightforward, the distribution of ventilation is associated with the choice of resting lung volume. We have kept the mean volume in each alveolar compartment constant and equal, with fluctuations around this mean based on FVi, the fraction of ventilation directed towards each alveolus i. This method is equivalent to decreasing the functional residual capacity (FRC) of alveolus i, as its ventilation increases. The phenomenon of sequential emptying of unequally ventilated alveoli is not a part of this model, as it has been shown that the relative contribution of different zones of the lung to expiratory flow remains constant (Farhi, 1969). The results presented in this paper have two functions, the stationary cycle data represents a simulation of pulmonary-cardiovascular pathology, while the transient state data represents a first attempt at parameter estimation on the basis of pulmonary function studies. Further research in the area of parameter estimation is possible with these results. For the measurement of cardiac output, an established bloodless (i.e. pulmonary) technique is that of Kim, Rahn, and Farhi (1966) as applied by Bickel et al. (1970). This method involves the estimation of arterial and venous PC02 from expired air gases on prolonged expiration. Using our model, we have simulated the respiratory maneuver described by Kim, Rahn, and Farhi (inhalation of a tidal volume larger than normal followed by steady exhalation for 20 to 30 sec.), and have used their methods of numerical
384
J. A. KELMAN
AND L. D. HOMER
analysis to calculate cardiac output, for different cases of cardio-pulmonary pathology. The estimated cardiac outputs come within 20 percent of the value preset in the programmed simulation, and even greater accuracy is possible (10 percent of preset values), with other choices of mathematical operations to convert PC02 to CO2 content. The successful application of this technique is indicative of internal consistency in the model, and also suggests a generalized approach to the problem of parameter estimation. Accuracy in Kim’s estimation procedure is limited to the degree that changes in cardiac output (or any other parameter being studied) will produce changes in recorded variables. There is the further limitation of ambiguity, as several combinations of cardiac output, dead space volume, and ventilation perfusion distribution may be capable of producing a single expired PO2 or PC02 value during prolonged expiration. Kim, Rahn, and Farhi use expired oxygen and carbon dioxide simultaneously to decrease this ambiguity, because of the different behavior of the two gases during expiration, as described in connection with Fig. 9 and 10. Data from the oxygen curve is used to calculate instantaneous respiratory quotient, the value of which determines the choice of the two points on the carbon dioxide curve which best correspond to the arterial and venous CO2 tensions (indicated in Fig. 8). Calculations based on two points from an expiratory flow curve, however, cannot be as accurate, or as free of ambiguity, as the use of the entire curve. Therefore, a logical extension of the techniques previously mentioned, is parameter estimation based on curve fitting of expired gas data with values produced by computer simulation. For purposes of parameter estimation, the lung model can be redefined in terms of a parameter set Bj, j = 1, k, to be estimated; a set of observed independent variables Xg, g = 1, p; and a set of dependent variables Yh, h = 1, m, to be observed. Tn this notation, B would include cardiac output, dead space volume, and V/Q distribution; X would include tidal volume, respiratory rate, gas flow, etc., * and Y would consist of observed Fd02, FdC02, and FdN2 at the mouth. If Yhbar is the value of expired partial fraction of gas h calculated by the model, then: Yhbar = fh(X1, X2,. . . Xp, Bl, B2,. . ., Bk) By generating a series of Yhbari, i = 1, n (n increments over time for some simulated respiratory maneuver), these values can be compared to the set of observed Y for each gas. Using programmed algorithms for leastsquares estimation of non-linear parameters (Horwitz and Homer, 1970; Marquardt, 1963) the set of Bj can be altered until a best fit of the model generated data (Ybar) to the observed data (Y) is achieved. The computer generated data can be simultaneously fitted to several observed variables, i.e. P02, PC02, and PN2, to increase the precision of the estimation and
MODEL
OF EXTERNAL
RESPIRATION
385
reduce the interdependence between the parameters. We can define a matrix (X), in which each element element Xij is the derivative of Ybar (for a single gas) with respect to parameter j, at time increment i: j = 1,k. i = 1,n; (X)nxk = (dYbari/dBj) If these initial calculations are generalized for three gases, the cumulative (X’X) inverse matrix can be constructed and normalized. These matrices, for a rebreathing maneuver, are presented in Table I. TABLE 1 PARAMETER TOMIC DEAD
21s ANA3 IS VENTILATION-
1 IS CARDIAC OUTPUT,PARAMETER
SPACE, AND
PARAMETER
PERFUSION DISTRIBUTION SHIFT(PERCENT)
Cumulative (X’X) 0.514E 09 - 0.394E 07 - 0.207E 05
inverse matrix: 0.394E 07 0.136E 07 0.104E 04
- 0.207E 05 - 0.104E 04 0.309E 01
0.149E 00 O.lOOE 01 - 0.505E 00
- 0.520E 00 - 0.505E 00 O.lOOE 01
-
Normalized: O.lOOE 01 - 0.149E 00 - 0.520E 00
-
The (X’X) inverse matrix can provide information on two aspects of parameter estimation, variance, and co-variance. The diagonal elements of this matrix are indications of the error expected in the parameter estimation. If we can simulate observed data with model generated data to a standard deviation of less than two percent, which is a realistic goal for this kind of curve fitting, the cardiac output can be estimated to within 12 percent, the dead space can be estimated to within 14 percent, and the ventilation perfusion shift can be estimated to within 8 percent. However, these values are based on the assumption that the estimation of each parameter is completely independent. This is obviously incorrect, and a measure of the interdependence of the estimations is given by the off-diagonal elements in the normalized (X’X) inverse matrix. For this maneuver (rebreathing), the co-variance between cardiac output and V/Q distribution is approximately 50 percent, while the co-variance between cardiac output and dead space is less than 15 percent. Therefore, 50 percent of the error in estimating cardiac output is due to error in the estimation of V/Q distribution, and 15 percent due to the dead estimation error. If these calculations are repeated for a prolonged expiratory maneuver, the standard errors are increased
386
J. A. KELMAN
AND L. D. HOMER
for the three parameters, while the total co-variance (sum of individual co-variances) decreases. Therefore, a combination of respiratory maneuvers may eventually provide the best estimation procedure. This system could be implemented to give on-line results, if analog to digital conversion of patient data were used. We hope to report on experimental confirmation of our model in the near future. APPENDIX I
Character
Subscript
BHC03 co c
i
j (gas)
ar
Cl, c2 V D
i
FV
i
FQB
i
F
i
j
192 132
192 j (gas)
al, a2 cl
G
I
J
(gas> I,2 HBMX
Definition
Units
standard bicarbonate content cardiac output concentration of gas j
ml CO,/ ml blood 37” ml/min ml gas/ml blood
in the arterial compartment in pulmonary capillary compartment 1 and 2 in the mixed venous compartment gas diffusion of gas j from alveolus 1 and 2 ventilation distribution fraction to alveolus 1 and 2 perfusion distribution fraction to alveolus 1 and 2 gas fraction of gas j in alveolus 1 and 2 in the dead space compartment in the inspired air gas content of gas j in alveolus 1 and 2 maximum binding potential of 0, to Hb in blood
ml gas/min
dimensionless
dimensionless
dimensionless
ml
ml 02/ml
blood
MODEL OF EXTERNAL
Character
RESPIRATION
Subscript
Hions K
_i
(gas) K’ PERM
j (gas)
P 50
P Prim PTO P
i
j (gas)
al, a2 CI, c2 m RATE SAT t VDOT
i 192
i aI, a2 art d VDIS
j
(gas)
387 Definition
hydrogen ion concentration solubility coefficient for gas j in blood dissociation content for carbonic acid permeability coefficient of gas j the 0, tension at 50% saturation after correction for pH the P,, at pH = 7.4 the ambient pressure tension or partial pressure of gas j in alveolus 1 and 2 in pulmonary capillary compartment 1 and 2 at the month respiratory rate hemoglobin saturation time gas flow to alveolus 1 and 2 (uncorrected for diffusion across the alveolus capillary interface) volume of alveolus compartment 1 and 2 of arterial blood compartment of dead space compartment volume of distribution of gas j
Units moles/liter ml gas (STP)/ ml blood moles/liter ml gas
mm-mm Hg mm Hg
mm Hg mm Hg mm Hg
cycles/min percent min ml/min
ml
ml
388 APPENDIX
Normal
J. A. KELMAN
AND L. D. HOMER
II
parameter
values :
BHC03 = 0.547 co = 4000. Fo02 = 0.21 FoC02 = 0. FoN2 = 0.79 HBMX = 0.201 k02 = 0.0000282 kC02 = 0.0006732 kN2 = 0.0000172 K’ = 7.95*10**(7) METAB = 250. METABC02 = 200. PPRIM = 26.5 (National Research Council, 1944) PTO = 760. RATE = 14. tidal volume = 400. ml. avg. Val = avg. Va2 = 1300. (therefore the FRC = 2400.) Vart = 50. VD = 175. VDIS02 = 5000. VDISC02 = 21000. VDTSN2 = 42000.
REFERENCES 1 R. G. Bickel, C. F. Diener, and H. L. Brammell, Aerospace Med. 41, 203 (1970). C. G. Douglas, and J. S. Haldane, J. Physiol. (London) 48,244 (1914). 3 J. H. Comroe, Physiology of respiration, Yearbook Medical Publishers, Chicago (1965), pp. 157-181. 4 H. T. Davis, Introductio,l to non-linear differential and integral equations, United States Atomic Energy Commission (1960), pp. 482-486. 5 L. E. Farhi, in Circulatory and respiratory mass transport (Wolstenholme and Knight, eds.), Little, Brown, and Co., Boston (1969), pp. 277-297. 6 J. S. Gray, The multiple factor theory of respiratory regulation. Army Air Force
2 J. Christiansen,
School of Aviation Medicine Proj. Rept. No. 386 (I, 2, 3; 1945). 7 F. S. Grodins, J. S. Gray, K. R. Schroeder, A. L. Norins and R. W. Jones, J. Appl. Physiol. 7, 283 (1954). 8 F. S. Grodins, J. Buell, and A. J. Bart, J. Appl. Physiol. 22, 260 (1967). 9 D. L. Horwitz and L. D. Homer, Analysis of Biomedical Data by Time-sharing computers, 1, Non-linear regression analysis, Naval Medicine Research Institute
Research Rept. No. 25 (1970).
MODEL OF EXTERNAL
RESPIRATION
389
10 T. S. Kim, H. Rahn, and L. E. Farhi, J. Appl. PhysioL 21, 1338 (1966). 11 D. B. Marquardt, .Z. Sot. In&f. Math. 11, 431 (1963). 12 National Research Council, Oxygen dissociation curves for human blood at 37 degrees Centigrade, in The handbook of respiration data in aviation, Washington D.C. (1944). 13 A. J. Olszowka and L. E. Farhi, J. Appl. Physiol. 26, 141 (1969). 14 R. Sikand, P. Cerretelli, and L. E. Farhi, .Z. Appl. Physiol. 21, 1331 (1966). 15 D. J. Strieder, B. A. Barnes, B. W. Levine, and H. Kazemi, .Z. Appt. Physiol. 29,486 (1970). 16 W. W. Wagner, Jr., L. P. Latham,
P. D. Brinkman, and G. F. Filley, Science 163, 1210 (1969). 17 J. B. West, Ventilution/bloodflow andgus exchange, Blackwell Scientific Publications, Oxford (1965), pp. 14-15. 18 W. S. Yamamoto and W. F. Raub, Computers and Biom. Res. 1, 65 (1967).