Acta Polytechnica doi:10.14311/AP.2013.53.0895 Acta Polytechnica 53(6):895–900, 2013 © Czech Technical University in Prague, 2013 available online at http://ojs.cvut.cz/ojs/index.php/ap ROBUST CONTROL OF END-TIDAL CO2 USING THE H∞ LOOP-SHAPING APPROACH Anake Pomprapaa,∗, Berno Misgelda, Verónica Sorgatoa, André Stollenwerkb, Marian Waltera, Steffen Leonhardta a Philips Chair for Medical Information Technology, Helmholtz Institute for Biomedical Engineering, RWTH Aachen University, Pauwelsstrasse 20, D-52074 Aachen, Germany b Embedded Software Laboratory, Chair of Computer Science 11, RWTH Aachen University, Ahornstrasse 55, D-52074 Aachen, Germany ∗ corresponding author: pomprapa@hia.rwth-aachen.de Abstract. Mechanically ventilated patients require appropriate settings of respiratory control variables to maintain acceptable gas exchange. To control the carbon dioxide (CO2) level effectively and automatically, system identification based on a human subject was performed using a linear affine model and a nonlinear Hammerstein structure. Subsequently, a robust controller was designed using the H∞ loop-shaping approach, which synthesizes the optimal controller based on a specific objective by achieving stability with guaranteed performance. For demonstration purposes, the closed-loop control ventilation system was successfully tested in a human volunteer. The experimental results indicate that the blood CO2 level may indeed be controlled noninvasively by measuring end-tidal CO2 from expired air. Keeping the limited amount of experimental data in mind, we conclude that H∞ loop-shaping may be a promising technique for control of mechanical ventilation in patients with respiratory insufficiency. Keywords: Closed-loop ventilation, embedded system, system identification, H∞ loop-shaping control design, biomedical application, control of etCO2. 1. Introduction In an intensive care unit, patients with respiratory insufficiency due to lung diseases, injury or undergo- ing a surgical procedure, may require support from a mechanical ventilator to maintain appropriate gas exchange: oxygenation and carbon dioxide (CO2) elim- ination from blood circulation [1]. Concerning CO2 exchange, if an abnormal value of CO2 pressure in arterial blood (PaCO2) persists for a long period of time, this can cause an imbalance of the pH value, which may be life threatening. Therefore, it is es- sential to regulate PaCO2 during ventilation therapy in order to avoid both hypercapnia and hypocapnia. Assuming good diffusion conditions (as in a healthy lung), PaCO2 can be approximated by end-tidal CO2 (etCO2) or CO2 partial pressure at end of expiration, which can be noninvasively measured from the exhaled air. Currently, no accurate mathematical model of the cardiopulmonary system is available that allows to estimate etCO2. Therefore, in this work, we pro- pose to define a model structure and to parametrize this model using system identification [2]. To this end, we have simplified the problem and assumed a single-input single-output (SISO) system. Minute ventilation (MV ; in L/min) was used to control the value of etCO2 [3–7]. To identify the parameters of a Hammerstein model, our results from a grey box identification are presented. For this, we assumed a nonlinear steady-state (or static) characteristic of the controlled plant (the patient) with, in addition, some linear time dynamics. Note that the nonlinear Ham- merstein model is composed of a serial interconnection between a linear time-invariant system and a static nonlinearity; this is also classified as a block-oriented structure [8, 9]. The advantage of the block-oriented structure is that it is able to represent input and output multiplicities; this makes it well suited for application on a cardiopulmonary system, due to its similar behavior of input and output multiplicities [10]. In addition, a linear affine model was identified in or- der to compare this with the nonlinear Hammerstein modeling results. After model validation, a robust controller was designed by an H∞ loop-shaping approach [11, 12]. Note that this method guarantees closed-loop stability whilst offering performance and robustness trade-offs [13]. Our goal for the H∞ loop-shaping approach was to tune the singular value of the open-loop gain or the open-loop transfer function gain, to increase the bandwidth of the system and to eliminate steady-state errors. The advantage of this method is to add perfor- mance possibilities while obtaining an exact solution in the H∞ optimal sense. For this design method, we present simulation results when evaluating the control performance based on various conditions of model uncertainty. Finally, a patient-in-the-loop ventilation system was connected to a human volunteer (first author) to test the etCO2 control algorithm in vivo. The remainder of this article is organized as follows. 895 http://dx.doi.org/10.14311/AP.2013.53.0895 http://ojs.cvut.cz/ojs/index.php/ap Anake Pomprapa, Berno Misgeld et al. Acta Polytechnica metabolism Patient Mechanical Ventilator Sensor (CO SMO+)2 • FiO2 • PIP • PEEP • RR • I:E ratio pH• etCO2 Pulse Oximeter Gas analyser dSPACE Control Unit C A N B U S R S 2 3 2 • SpO2 Medical Panel PC • FiO2 • MV µC µC Figure 1. Configuration of the system for closed-loop ventilation. Section 2 describes a unified approach for system identification including a physiological description; here, all the models are parametrized and validated. Section 3 presents the robust control design based on the H∞ loop-shaping technique and the simulation results of the controller under various conditions of model uncertainty. Section 4 presents a discussion of this work and the conclusions are presented in Section 5. 2. System modeling 2.1. System configuration The proposed closed-loop system is composed of a medical panel PC for process monitoring, user inter- face and data storage, a mechanical ventilator (VEN- TIlogic LS, Weinmann Geräte für Medizin GmbH, Germany), a capnograph (CO2SMO+, Philips GmbH, Germany) for etCO2 monitoring, a MicroAutoBox II dSpace control unit, and two ARM-based micro- controllers. The system configuration is presented in Figure 1. The data communication protocol was de- signed based on the CAN (Controller Area Network) protocol. CAN-Bus is a serial fieldbus, which allows additional devices to be connected to the system ar- chitecture using this topological arrangement. A data transfer rate of 1 Mbit/s can be obtained and collision avoidance between the messages can be achieved for all connected devices, based on priority assignment. Therefore, the proposed closed-loop system is suit- able for real-time automatic control of mechanical ventilation. 2.2. Static nonlinearity To serve as an example, system identification was conducted based on one male volunteer with healthy homogeneous lungs and a body mass index within normal range (i.e. 23.3 kg/m2). This person was con- nected to a mechanical ventilator, which was operating in pressure-controlled mode. All ventilation settings were manually adjusted to extract cardiopulmonary information from the subject. Based on the static characteristics of the patient, etCO2 is a nonlinear function of MV [4]. Its static 8 10 12 14 16 18 20 22 24 20 25 30 35 MV [L/min] e tC O 2 [m m H g ] Measured data from a human subject Estimated curve of the relationship Figure 2. Static nonlinearity between etCO2 and MV . nonlinearity is the so-called “metabolic hyperbola” (Figure 2). By extracting CO2 information, MV input was increased stepwise from 10 L/min to 25 L/min, and etCO2 was measured at steady state. We could indeed confirm that etCO2 is a decreasing function in terms of MV input at steady state. Hence, the more MV applied to the lungs, the less etCO2 can be measured from the subject. This relationship is important and can be employed for controlling etCO2 with MV input. A mathematical description of the nonlinearity may approximate the nonlinearity as a parabolic equation as provided in (1). etCO2 = N[MV ] = a · MV 2 + b · MV + c (1) where N[MV ] denotes a nonlinear function of MV . In this particular example, a = 0.05, b = −2.55 and c = 52.80 are the best parameters according to a least-squares fit. 2.3. Linear affine model As an initial estimation, a simple dynamic model can be applied to complex input-output relationship of the system in order to evaluate to what extent such a simple model can represent the real system (Figure 3). Such an affine model is formulated by a linear combi- nation form of primitive variables, provided in (2): y(k) = yo + p∑ i=1 a(i) ·y(k−i) + q∑ j=1 b(j) ·u(k−j) (2) where y(k) represents etCO2 at the sampling point k, and u(k) denotes the input ∆P at the sampling k. The sampling time for this study was 4.28 s. Also, p and q are finite order parameters with p ≥ q and yo representing a constant so called offset. The unknown parameters (yo,a(i) and b(j)) can be estimated by using a least-squares algorithm based on the input and output measurements. The following algorithm is used to estimate the unknown parameters. According to the experimental data and the model structure from (2) for the linear affine model, it can 896 vol. 53 no. 6/2013 Robust Control of End-Tidal CO2 15 20 25 30 35 40 e tC O 2 [m m H g ] measured etCO 2 1 st order model 2 nd order model 2 nd order model with zero 1 st order Hammerstein model 0 100 200 300 400 500 600 700 800 900 1000 0 5 10 Time [s] D P [h P a ] Figure 3. Input-output relationship used for system identification with the results of model identification. be expressed in terms of a vector and matrices, as in (3): Y = X ·β + ε, (3) where Y = [ y(ko + 1) y(ko + 2) · · · y(M) ]T , XT =   1 1 · · · 1 y(ko) y(ko + 1) · · · y(M − 1) ... ... ... y(ko −p) y(ko + 1 −p) · · · y(M − 1 −p) u(ko) u(ko + 1) · · · u(M − 1) ... ... ... u(ko −q) u(ko + 1 −q) · · · u(M − 1 −q)   , β = [ yo a1 · · · ap b1 · · · bp ]T for ko ≥ p and ko ≥ q and ε is an unknown disturbance vector. The vector β or the unknown parameters can then be estimated using the ordinary least-squares algorithm that minimizes the sum of squared errors provided in (4) [2]: β = (XT · X)−1 · XT · Y. (4) In order to relax the muscles involved in respiratory breathing, a fixed respiratory rate (RR) is introduced at 14 bpm to allow the subject to minimize work of breathing [5]. Since the mechanical ventilator was set to pressure-controlled mode, MV was obtained by a multiplication of tidal volume (VT ) and respira- tory rate (RR). Because of the fixed RR of 14 bpm, the input for this particular system is transformed from MV to the difference in the driving pressure (∆P = PIP −PEEP) between inspiration and expira- tion pressure, and has a direct influence on VT . 2.4. Nonlinear Hammerstein model Providing a representation of signal flow using a block- oriented structure, the Hammerstein model comprises a static nonlinearity N[·] at the input u(k) cascaded with a linear dynamic model H(z): N[•] H(z) u(k) v(k) y(k) The static nonlinearity maps the input u(k) into the intermediate variable v(k) (see (5)) and the linear dynamic model maps the intermediate variable v(k) into the output y(k) (see (6)). The model can be represented by v(k) = N[u(k)], (5) y(k) = p∑ i=1 a(i) ·y(k − i) + q∑ j=1 b(j) ·v(k − j). (6) Rearranging (5) and (6), the Hammerstein model can be described as shown in (7). y(k) = p∑ i=1 a(i) ·y(k− i) + q∑ j=1 b(j) ·N[u(k−j)] (7) It can be seen from (7) that the Hammerstein model is very similar to the linear dynamic model. Because the qualitative behavior of the transient response is entirely determined by the discrete transfer function of the linear subsystem H(z), it can be used as an alternative to the linear model. This model can ex- hibit input multiplicities if the static nonlinearity is in the form of input multiplicities. According to the experimental data and the model structure, the un- known parameters (a(i) and b(j)) can be estimated by a constrained least-squares algorithm, as provided in (3) and (4). 2.5. Evaluation of model structure Both the linear and nonlinear model structures are used to describe this system. The performance index used in our evaluation was obtained by root mean square error (RMSE). Table 1 presents the results of the comparisons, divided into an estimated dataset and a validated dataset for the different model struc- tures. Based on the validation dataset, the first-order Hammerstein model provides the best result of all the listed models. Nevertheless, the first-order lin- ear model also offers the best result of the RMSE evaluation of all the linear models. In addition, qualitative comparison of the selected models is provided in Figure 3. In the following sec- tion, the design of the H∞ controller and simulation re- sults are conducted using the first-order linear model. It should be emphasized that capnography has an accuracy of ±2 mmHg for values in the 0–40 mmHg range, ±5 % for values in the 41–70 mmHg range, and ±8 % for values in the 71–150 mmHg range [14]. 897 Anake Pomprapa, Berno Misgeld et al. Acta Polytechnica Estimated Validated RMSE RMSE 1st Linear Affine 2.2475 2.2880 2nd Linear Affine 2.2116 2.2988 2nd Affine w. zero 2.1597 2.4093 1st Hammerstein 2.1988 1.6709 2nd Hammerstein 2.1680 1.7804 2nd Ham. w. zero 2.1351 1.8085 Table 1. Root mean squared error (RMSE) evalua- tion of the different model structures. r y K Gp - - N M -1 u DN DM Figure 4. The H∞ control structure with left coprime factorization of a plant G. 3. Design of robust controller 3.1. Loop-shaping design using H∞ This method combines the principle of Bode’s sensitiv- ity integral [15] and the H∞ optimization technique by minimizing the H∞ norm in the presence of uncer- tainty. In designing the H∞ controller, both stability and performance are taken into account, with bounded differences between the nominal model and the real nonlinear plant. Given a nominal discrete-time model of a plant G, it can be represented using a normalized left coprime factorization (LCF) G = M−1 ·N, (8) where M and N are coprime matrices in RH∞ (Fig- ure 4): A perturbed model associated with the LCF rep- resentation of the plant G is given by (9), where perturbations are assumed to be bounded, is given by Gp = (M + ∆M)−1 · (N + ∆N),∥∥∆M ∆N∥∥∞ < � (9) The objective is to find a robust controller K that stabilizes Gp and minimizes (10). γmin = ∥∥∥∥K(I −GK)−1M−1(I −GK)−1M−1 ∥∥∥∥ ∞ (10) The solution of this H∞ norm problem can be solved using the algorithm proposed by McFarlane and Glover [11]. Define [A,B,C, 0] to be a state space 10 -3 10 -2 10 -1 10 0 -20 -10 0 10 20 30 40 50 60 70 80 Singular Values Frequency [rad/s] S in g u la r V a lu e s [ d B ] Plant Singular Value Desired Singular Value Designed Loop Gain Singular Value s( )G s( )GK Figure 5. Shaping the singular value of the single- input single-output (SISO) system. minimal realization of plant G. Then, the subopti- mal H∞ controller K can be computed by discrete algebraic Riccati equations [16]. The loop-shaping objective is to design a robust controller K, so that σ(GK) � 1 or |GK(jω)| � 1 (for a case of SISO) at low frequencies (minimizing the effect of output disturbances) and σ(GK) � 1 or |GK(jω)| � 1 at high frequencies (minimizing the effect of sensor noise and providing robustness for additional uncertainty). The singular value of the open loop gain GK is shaped based on this design criterion. Based on the singular value of the plant in Figure 5, the integral action has been chosen to ensure zero steady-state error. In addition, the cut-off frequency is designed to be 0.28 rad/s in comparison with a very low cut-off frequency at 0.03 rad/s of the plant. In this way, the bandwidth is increased by a factor of approximately 10 times. Referring to [11], γ = 1.86 or � = 0.5376 indicates the allowable proportional uncertainty in N and M of approximately 50 % in the crossover frequency range of the shaped plant. 3.2. Step response simulation A step response, shown in Figure 6, represents the con- trol performance of the H∞ loop-shaping controller for both the linear affine and the nonlinear Hammerstein models. The etCO2 is controlled to have a target of 34 mmHg before t = 0 s and is forced to 35 mmHg after t = 0 s. In both models, rise time is 14 s and settling time is within 60 s. The response has no steady-state error. The robust H∞ loop-shaping controller gives a good transient and steady-state performance for the linear affine and the nonlinear Hammerstein models. For all model parameters, a satisfied step response can also be achieved with a parameter uncertainty of about 12 %. Therefore, the controller is robust in coping with model and parameter uncertainty, and with disturbance. 898 vol. 53 no. 6/2013 Robust Control of End-Tidal CO2 33 34 35 36 e tC O 2 [m m H g ] H ∞ loop-shaping controller for linear affine model 0 10 20 30 40 50 60 33 34 35 36 Time [s] e tC O 2 [ m m H g ] H loop-shaping controller for Hammerstein model ∞ -10-20 Figure 6. Step response of the H∞ loop-shaping controller for the linear affine and the Hammerstein models. 0 20 40 60 80 100 120 140 160 180 200 220 32 34 36 38 40 42 Time [s] e tC O 2 [m m H g ] Figure 7. Details of controller performance obtained from a human volunteer. 3.3. Evaluation of closed-loop control ventilation After simulating the controller performance, closed- loop control ventilation is implemented and tested on a human volunteer. The result of the control performance is presented in Figure 7. Robustness with a step response can be achieved for the required etCO2 at 35 mmHg under real dis- turbance and measurement uncertainty. Within ap- proximately 60 s, the target etCO2 is satisfied. The response lies in an acceptable range for clinical ap- plication. It should be emphasized that a relatively good performance can be obtained by the robust H∞ loop-shaping controller. 4. Discussion For an abnormal lung condition, like acute respiratory distress syndrome (ARDS), etCO2 does not corre- spond to PaCO2 and therefore invasive measurement of PaCO2 is required. Therefore, our focus of iden- tification and control design for etCO2 is only valid for patients after treatment with the Open Lung re- cruitment maneuver [17], or for patients whose etCO2 reading appropriately reflects the true PaCO2 (as, for example, in chronic obstructive pulmonary disease) [18]. In such cases, the mean value of CO2 pressure in arterial blood (PaCO2) is approximately 9 mmHg higher than the value of etCO2 and the required value of etCO2 should be adapted based on the calibration with PaCO2, which can be obtained from blood gas analysis. The second-order model with one zero is correlated to the pharmacological two-compartment model as proposed in [7]. The etCO2 represents the output from the lung compartment and is one of the state variables in the model. Zero position relies on gas transport from tissues to the lung. Our findings of poles and zero positions based on animal experiments correlate well with the findings derived from 18 patients [7]. One pole is located near the origin of the unit circle and another pole is near the point (1, 0) in the unit circle. The identification parameters are subject to measur- ing errors due to the limitations of capnography: its accuracy for etCO2 monitoring is ±2 mmHg for the 0– 40 mmHg range, ±5 % for the 41–70 mmHg range, and ±8 % for the 71–150 mmHg range, and the resolution is 1 mmHg [14]. Therefore, the identified parameters will not perfectly reflect the underlying parameters of the plant. In other words, parameter uncertainty also exists because of the limitations related to the accuracy and resolution of the measuring device itself. Considering the input applied to the system, MV shows a better result compared to inverse minute ventilation (IMV ) for both the linear affine and the Hammerstein model. The design problem is simplified to be a SISO system by regarding MV as the input and etCO2 as the output. In a pressure-controlled ventilation mode, tidal volume cannot be directly ad- justed. Thus, a pressure difference should be changed in order to meet the required tidal volume. The Hammerstein model provided better numerical results compared with the affine model, especially for the validation dataset (Table 1). Note that the Hammerstein model has been successfully applied in several other biomedical applications [19], including the stretch reflex EMG [20] and heart rate regulation [21]. In our clinical application, the complex nonlinear cardiopulmonary system can be better modeled by a Hammerstein model than by an affine model. It should be noted that the block-oriented NARMAX models [8], which offer a modeling for output multiplicities, did not give acceptable results when they were tested; therefore, those results are not presented or discussed here. Regarding patient safety, MV that is too low leads to low oxygenation and a risk of mortality. On the other hand, extremely high MV carries a high risk of trauma and a possibility of lung damage. Thus, ac- tuator saturation should be introduced in our system design with the aid of an anti-windup technique and should be considered for future research work. 5. Conclusion In a clinical application, etCO2 is required to be feedback-controlled to a certain value to minimize the risk of hypercapnia or hypocapnia. To realize this task, we propose a model-based approach by 899 Anake Pomprapa, Berno Misgeld et al. Acta Polytechnica identifying the model parameters of a complex non- linear cardiopulmonary system using a block-oriented structure with the linear affine and the Hammerstein models. A robust control design was implemented using the H∞ loop-shaping approach based on the derived affine model. The simulation results indicate a good control performance of the H∞ loop-shaping controller for both the linear affine and the nonlinear Hammerstein models, including possible parameter variations up to 12 %. Finally, for demonstration pur- poses, the controller was tested in a task to control etCO2 in a healthy volunteer and a positive result was achieved with the robust H∞ loop-shaping controller. Acknowledgements The authors acknowledge financial support from the Ger- man Federal Ministry of Science and Education (Bun- desministeriums für Bildung und Forschung – BMBF) within the Oxivent project under grant number 16SV5605. References [1] A. Pomprapa, D. Schwaiberger, B. Lachmann, S. Leonhardt. Predictive 3D model of CO2 elimination for optimal pressure-controlled ventilation. Am J Respir Crit Care Med 185(3):A1717, 2012. [2] L. Ljung. System Identification Theory for the User. Prentice Hall, New Jersey, 2nd edn., 1999. [3] F. W. Chapman, J. C. Newell, R. J. Roy. A feedback controller for ventilatory therapy. Ann Biomed Eng 13:359–372, 1985. [4] R. Rudowski, C. Spanne, G. Matell. Computer simulation of a patient end tidal CO2 controller system. Comput Meth Prog Bio 28:243–248, 1989. [5] A. Pomprapa, M. Walter, C. Goebel, et al. L1 adaptive control of end-tidal co2 by optimizing the muscular power for mechanically ventilated patients. 9th IFAC Symposium on Nonlinear Control Systems 2013. [6] S. Mersmann. Smartcare: automatizing clinical guidelines. Biomed Tech 54:283–288, 2009. [7] J. Hahn, G. Dumont, M. Anersmino. System identification and closed-loop control of end-tidal CO2 in mechanically ventilated patients. IEEE Trans Inf Technol Biomed pp. 1–9, 2012. [8] M. Pottmann, R. K. Pearson. Block-oriented NARMAX models with output multiplicities. AIChE Journal 44(1):131 – 140, 1998. [9] R. K. Pearson, M. Pottmann. Gray-box identification of block-oriented nonlinear models. J Process Contr 10(4):301–315, 2000. [10] A. Pomprapa, R. Pikkemaat, H. Luepschen, et al. Self-tuning adaptive control of the arterial oxygen saturation for automated open lung maneuvers. 3 Dresdner Medizintechnik-Symposium 2010. [11] D. C. McFarlane, K. Glover. A loop-shaping design procedure using H∞-bounded uncertainty. IEEE Trans Autom Control 37(6):759–769, 1992. [12] R. A. Hyde. H∞ aerospace control design: A VSTOL flight application. Springer Verlag, 1995. [13] D. C. McFarlane, K. Glover. Robust controller design using normalized coprime factor plant descriptions. Springer Verlag, 1989. [14] Novametrix Medical Systems Inc. CO2SMO+ respiratory profile monitor service manual model 8000. Catalog No 9758-90-02 (Rev 2) p. 67, 2001. [15] H. W. Bode. Network analysis and feedback amplifier design. D. Van Nostrand Inc., 1945. [16] D. W. Gu, P. H. Petkov, M. M. Konstantinov. Formulae for discrete H∞ loop-shaping design procedure controllers. 15th IFAC World Congress 15(1):351, 2002. [17] J. J. Haitsma, R. A. Lachmann, B. Lachmann. Open lung in ARDS. Acta Phamacol Sin 12:1304–1307, 2003. [18] M. Kartal, E. Goksu, O. Eray, et al. The value of etCO2 measurement for COPD patients in the emergency. Eur J Emerg Med 18:9 – 12, 2011. [19] I. W. Hunter, M. J. Korenberg. The identification of nonlinear biological systems: Weiner and Hammerstein cascade models. Biol Cybern 55:135–144, 1986. [20] D. T. Westwick, R. E. Kearney. Identification of a Hammerstein model of the stretch reflex EMG using separable least squares. IEEE EMBS 3:1901 – 1904, 2000. [21] S. W. Su, S. Huang, L. Wang, et al. Nonparametric Hammerstein model based model predictive control for heart rate regulation. IEEE EMBS pp. 2984 – 2987, 2007. 900 Acta Polytechnica 53(6):895–900, 2013 1 Introduction 2 System modeling 2.1 System configuration 2.2 Static nonlinearity 2.3 Linear affine model 2.4 Nonlinear Hammerstein model 2.5 Evaluation of model structure 3 Design of robust controller 3.1 Loop-shaping design using H infty 3.2 Step response simulation 3.3 Evaluation of closed-loop control ventilation 4 Discussion 5 Conclusion Acknowledgements References