Multiple model adaptive control of neuromuscular blockade: Design guidelines and clinical cases

Multiple Model AdaptiveControl of neuromuscular blockade: Design guidelines and clinical cases Teresa Mendonça, Hugo Magalhães, João M. Lemos, Paula Rocha and Simão Esteves Abstract—The high level of uncertainty of the dynamic response of patients subject to anaesthesia motivates the use of adaptivecontrol methods. This paper proposes an approach based on Switched Multiple Model AdaptiveControl (SMMAC) to tackle this problem in what concerns the control of the neuromuscular blockade level. It is shown howto design the different elements of the SMMACcontroller,enhancing the importance of the observer polynomial, that is shown to be instrumental to stabilize the loop. Clinical results using atracurium as blocking agent are reported, thereby illustrating the application of the proposed approach in actual clinical practice. Index Terms—Anaesthesia, Neuromuscular Blockade, Computer Control, AdaptiveControl, Biomedicine, control applications. I.INTRODUCTION Feedback control for drug dosing in clinical pharmacology is receiving increasing attention [1], [2]. Besides circumventing the tedious, imprecise and sometimes low performance procedures associated to open-loop control, feedback control can also significantly advance the medical understanding of the effects of pharmacology agents and anesthetics and provide progress in drug delivery systems. Significant examples [1] are provided by the closed-loop control of the cardiovascular function and automated anaesthesia. During surgical procedures patients are usually under general anaesthesia, defined as the lack of response and recall to noxious stimuli, reflected in loss of conscience, pain insensitivity and muscle paralysis. Several approaches to the neuromuscular blockade control problem as well as an introduction to this issue can be seen in [3], [6], [4] [7] and [5]. Muscle relaxant drugs are frequently given during surgical operations. The non-depolarizing types of muscle relaxant act by blocking the neuromuscular transmission (NMT), thereby producing muscle paralysis. The level of muscle relaxation is measured from an evoked EMG at the hand by electrical stimulation of the adductor policies muscle to supra-maximal train-of-four stimulation of the ulnar nerve. In aclinical environment the measurement of the neuromuscular blockade level corresponds to the first single response (T1%)calibrated by areference twitch, This paper was written in relation with the NISIS EC network. T.Mendonça and H. Magalhães are with the DMA, FCUP,Rua do Campo Alegre, 687, 4169-007 Porto, Portugal. {tmendo,hfmagalh }@fc.up.pt J.M. Lemos is with the INESC-ID/IST,Rua Alves Redol, 9, 1000-029 Lisboa, Portugal. [email protected] P.Rocha is with the DM, UA, Campo de Santiago, 3810-193 Aveiro, Portugal. [email protected] Simão Esteves is with the Anestesiology Department, Hospital Geral de Santo António, Largo Prof. Abel Salazar,4050 Porto, Portugal. obtained by defining asupra-maximal stimulation current. This measuring process is prone to the occurrence of outliers, aproblem dealt with in [8]. The control of the neuromuscular blockade provides agood illustration of the main features and inherent constraints associated with the control of physiological variables. It is characterized by avery high degree of uncertainty in the dynamics of the system due both to interpatient variability as well as time variations. This suggests that multiple-models based control techniques would provide asuitable solution [9]–[11]. Such ascheme (incorporating several modifications in order to accommodate the specific characteristics of the problem) has been developed in [12], using multiple controllers constructed from abank of models which replicate the observed variability of the dynamic responses to the muscle relaxant [13]. In [14] restriction techniques (controller localization) relying on arobustness condition havebeen proposed. These techniques yield good results but suffer from the major drawback of assuming a strong assumption on the plant dynamics. The main contribution of this paper consists in providing guidelines for designing the Switched Multiple Model AdaptiveControl algorithm for neuromuscular blockade using atracurium,and to validate them with actual clinical cases. Aspecial emphasis is placed in the selection of the observer dynamics, for which guidelines for the selection of the characteristic polynomial are provided. The paper is organized as follows: after the introduction which motivates and formulates in general terms the problem to be solved, the model for neuromuscular blockade is described in section 2. The switching multiple model control algorithm used is described in section 3. Section 4describes results, both simulations and clinical trials. Finally,section 5draws conclusions. II.NEUROMUSCULAR BLOCKADEMODEL In order to design the bank of local controllers upon which the switched multiple model algorithm relies. Furthermore, this model reveals the structure of the system to control, thereby providing insight on design options. The dynamic response of the neuromuscular blockade for atracurium may be modelled as shown in Figure 1[7], [13]. Alinear pharmacokinetic model (block 1of Figure 1), described by the following linear system of state equations        ˙x1(t)=−λ1x1(t)+a1u(t) ˙x2(t)=−λ2x2(t)+a2u(t) cp(t)=X2 i=1 xi(t), (1) Proceedings of the European Control Conference 2007 Kos, Greece, July 2-5, 2007 WeA01.3 ISBN: 978-960-89028-5-5 2541 1 λ λ + s 2 τ τ / 1 /1 + s 3 Hill equation 4 u(t)cp(t)c(t)ce(t)r(t) 2 2 1 1 λλ + + +s a s a Fig. 1. Block diagram of the neuromuscular blockade model. relates the drug infusion rate u(t)[µg kg−1min−1]with the plasma concentration cp(t)[µg ml−1], where xi(i= 1,2) are state variables (implicitly defined by (1)) and ai [kgml−1], λi[min−1](i=1,2) are patient dependent parameters. The physiological basis of the model described by equation (1) consists of assuming two plasma compartments (central and peripheric) both communicating with each other. Alinear second order model (blocks 2and 3of Figure 1), described by the cascade of two first order systems, written as ˙c(t)=−λc(t)+λcp(t)(2) and ˙ce(t)=−1/τce(t)+1/τc(t),(3) is assumed to relate cp(t)with the concentration in the effect compartment, ce(t)[µg ml−1]. Here, c(t)is an intermediate variable and λ[min−1], τ[min]are patient dependent parameters. It is remarked that standard models developed for atracurium [?], [15] do not consider the block 3. As shown in [13], the inclusion of the extra delay associated to τallows abetter replication of the observed experimental responses. Finally the pharmacodynamic effect, that relates ce(t)to the induced pharmacodynamic response, r(t)[%], may be modelled by the Hill equation [15](block 4of Figure 1), r(t)=100Cγ 50/(Cγ 50 +cγ e(t)),(4) where the parameters C50 [µg ml−1]and γ(adimensional) are also patient-dependent. The variable r(t),normalized between 0and 100, measures the level of the neuromuscular blockade, 0corresponding to full paralysis and 100 to full muscular activity. III.SWITCHING CONTROL In recent years adaptivecontrol approaches based on supervisory switching control havebeen proposed to deal with systems presenting high level of uncertainty.This section describes the basic structure of switching control and presents two different solutions to overcome robustness issues related to the implementation of such acontrol scheme. It is assumed that adiscretization procedure has taken place, and hence from nowon the variable trepresents discrete time instants. A. Basic structure The basic structure of asupervisor based switched multiple model controller is seen in Figure 2, as described in [9], [10], [14]. In Figure 2, ydenotes the sensor measure, ref denotes the reference, ecthe control error and Pthe plant to be controlled. Abank of controllers Cj,j=1,..., Nis designed to match the plant models Mj.This set of models is assumed to “cover” all the possibilities of the actual plant P. In order to select at each time which controller best matches P,the principle according to which the best model performance implies the best controller performance is applied. One possibility for evaluating the model performance is to compare the output yjof each model Mjwith the process output y.Another possibility is to construct estimators Ej based on each of the models Mj,see subsection III-C. In either case an error ejis produced, which is measured through aperformance index(PI) πj,j=1,..., N,computed by lowpass filtering according to πj(t+1) =λππj(t)+(1 −λπ)ej(t),(5) where tis discrete time and λπis aparameter to adjust. The switching logic block SL selects the indexσof the controller to apply to the plant. This selection is given by the value of jcorresponding to the least value of πj,but ensuring that, once acontroller is applied to the plant, it remains so for a minimum amount of time τD.This is the so called “dwell time” condition, which prevents high frequencycommutation among controllers and prevents instabilities that could occur due to too fast switching [9]. An integrator common to all blocks ensures bumpless transfer between different controllers [9] ¡discrete control transfer function Cj(z)= ¯ Cj(z)z∆t z−1,∆tsampling time¢.Outlier removal is performed with aBayesian filter,according to the techniques described in [8]. P 1 − ∆ z tz C1 CN SL M1/E1 MN/EN PI PI σ ref r 1 π N π r r + + - - u + -v ec 1 e N e Fig. 2. Switched multiple model control strategy. WeA01.3 2542 B. Local controller design Each controller Cjin the bank of local controllers has a PID structure with the integrator separated, parameterized as v(t)=kpÃ1 Ti +z−1 z∆t+Tdµz−1 z∆t¶2!ec(t),(6) u(t)=z∆t z−1v(t).(7) The gains kp,Tiand Tdare designed according to a dominant-pole placement rule [16] by applying the expressions: n=b1w1−b2w2+ζa2w2, a=a2−a1+ζb2 αn, b=b2w1−b1w2−ζa2w1 αw1w2n, Ti=−a 2ra2 4−b, (8) c=−ζa2+b2−b1, d=(a2+ζb2)µw2αTi−1 w2Ti¶, e=−a1µw1αTi−1 w1Ti¶, kp=ζ c+d+e,(9) Td=αTi,(10) where G(jwi)=ai+jbi,i=1,2are two points of the plant Nyquist curve, ζis the damping ratio and αis a priori fixed. Typical points of the Nyquist curveare those which correspond to the plant phase frequencies −180◦and avalue close, but inferior to it, and acommon value for α is α=4.The abovemethod for designing PIDs proves to be an adequate one, even for non-minimum phase plants and provides astep response with little overshoot. C. Observer Dynamics As done in [9], the estimate ˆyj(t)of y(t)is made by including an observer polynomial. Let each model Mjbe represented by the ARX model Aj(q−1)yj(t)=Bj(q−1)u(t)+¯ej(t),(11) in which Aj(q−1)=1 +Xna i=1 aj,iq−i,(12) Bj(q−1)=Xnb i=1 bj,iq−i(13) are polynomials in the unit delay operator q−1,with Aj monic and Bjof fixed degrees naand nbrespectively,for all j=1,..., N.In order to associate to (11) astate space model, define (nb>1) s(t)= [yj(t). . . yj(t−na+1) u(t−1) ... u(t−nb+1)]T.(14) For nb=1,s(t)=[yj(t). . . yj(t−na+1)].In the case of an ARX model, i.e.when the disturbance acting on (11) is white noise, s(t)is a(nonminimal) state associated to (11), whose evolution is described by the state-space model s(t+1) =Φjs(t)+Γju(t)+G¯e(t+1) yj(t)=Hs(t)(15) in which Φj= "−aj,1... −aj,nabj,2. . . bj,nb Ina−10na−1×10na−1×nb−20na−1×1 0... 0 0 . . . 0 0 0nb−2×1Inb−20nb−2×1# (16) Γj=[bj,101×na−11 01×nb−2]T,(17) with the “1” in Γjin the na+1th position (nb>1) and G= [1 01×na+nb−2]T,H=[1 01×na+nb−2].The pair (Φj,Γj) is reachable and (H,Φj)is reconstructible (the state can be obtained from past data) but not observable (the state cannot be recovered from future data). The advantage of using this type of state consists in the fact that it readily provides ashared state realization, i.e.,astate-space realization in which the state is common to all models Mj/Ejin Figure 2. Consider nowthe problem of designing astate observer to (15) or equivalently,aprediction to (11). For that sake add Aoyj(t)−Ajyj(t)to both sides of (11) to conclude that the model may be represented by yj(t)=(Ao−Aj)1 Ao yj(t)+Bj Ao u(t)+1 Ao ¯ej(t),(18) where Ao(q−1)=1+Pna i=1 ao,iq−iis astable monic polynomial with the same degree as the Aj’s, which will hereafter be referred as the “observer polynomial”. Define an estimator Ejfor the plant by ˆyj(t)=(Ao−Aj)1 Ao y(t)+Bj Ao u(t).(19) The idea behind this is that if the process model would coincide with Mj,ywould coincide with yj.Therefore y(t)−ˆyj(t)=yj(t)−ˆyj(t)=1 Ao ¯ej(t)(20) and the estimate ˆyjwould be asymptotically accurate provided that Aois chosen to be Hurwitz. Indeed, the tracking error is alowpass filtering of the disturbances. If adeadbeat ("fast") observer is selected (all the roots of Aoat the origin) no filtering occurs. As confirmed by the experiments below,amuch better choice consists in making the observer dynamics much slower,by selecting the roots of Aoin the real segment between 0and 1,but away from 0. D. State space interpretation The estimator (19) can be interpreted in terms of astate observer as follows. Define η(t)=1 Ao(q−1)y(t),(21) ν(t)=1 Ao(q−1)u(t).(22) WeA01.3 2543 With the state defined by xη(t)=[η(t)η(t−1) . . . η(t−na+1)]T,(23) the state-space realization of (21) is seen to be xη(t)=¯ Φoxη(t−1) +¯ Γoy(t),(24) where ¯ Φo=·−ao,1−ao,2. . . −ao,na Ina−10na−1×1¸,(25) ¯ Γo=[1 01×na−1]T.(26) For (22), asimilar reasoning holds. The state defined by xν(t)=[ν(t)ν(t−1) . . . ν(t−na+1)]T(27) yields xν(t)=¯ Φoxν(t−1) +¯ Γou(t).(28) In turn, defining the state of the reconstructor as sf(t)=£xη(t)Txν(t)T¤T,(29) equation (24) together with (28) yield sf(t)=Φosf(t−1) +Γoyy(t)+Γouu(t),(30) in which Φo=·¯ Φo0na×na 0na×na¯ Φo¸,(31) Γoy =£¯ ΓT o01×na¤T,(32) Γou =£01×na¯ ΓT o¤T.(33) From (19) ˆyj(t)=Hj,osf(t−1) (34) with Hj,o given by Hj,o =[ao,1−aj,1...ao,na −aj,nabj,1...bj,nb01×na−nb] (35) The vector sf(t)corresponds to afiltering of s(t)by the observer dynamics. This dynamics may be chosen in order to provide the controller with robustness properties with respect to un-modelled plant dynamics. Apossibility consists in choosing adead-beat observer,in which case ao,i =0for i=1,...,nand s(t)is exactly recovered after nsteps. However,as already mentioned, the dead-beat observer is an inconvenient choice in the case of neuromuscular blockade where the possible dynamics are such that this may lead to awrong selection of the best local controller. IV.RESULTS Toachieve a high level of neuromuscular blockade in a short time, abolus of atracurium is always administered in the beginning of asurgery.After the administration of the bolus,the level of the muscular blockade increases very quickly (the variable r,that measures muscular activity decreases), and full muscle paralysis is induced in afew minutes. Following that initial period, the control objective is to followaspecific reference profile with afinal target value ref ≡ref0.The value of the reference profile is initially fixed at alowlevel (typically 2.5%) during the first 30 minutes. It is then gradually increased to the final value (typically 10%). Abank of N=100 non-linear dynamic models Mj,j=1, . . . , Nwith the same structure as described in section II was generated using the probabilistic model for atracurium [13]. Furthermore, for each Mj,acontroller Cj was tuned using the dominant-pole placement rule described in subsection III-B. This provides acontroller bank which covers awide range of behaviors. Anumber of simulations are shown hereafter which illustrate the results obtained with switching control. A. Base line examples In this set of simulations, abase line is established and the observer polynomial is selected with all the roots at the origin (no filtering). Figure 3illustrates the results for P=M69,where the pair (M69,C69)has been removed from the model-controller bank. As seen in Figure 3the resulting behavior is poor. 0 50 100 150 200 250 0 5 10 15 20 r(t) % Model 69 r ref 0 50 100 150 200 250 0 5 10 15 20 25 u(t) µg kg−1 min−1 0 50 100 150 200 250 −8 −6 −4 −2 0 2 4 6 Time (minutes) ec=ref−r 0 50 100 150 200 250 0 20 40 60 80 100 Time (minutes) Switching Signal 75 Fig. 3. Results obtained for P=M69 with all the roots of the observer polynomial selected at the origin (no filtering). B. Observer dynamics effect If the roots of the observer polynomial are chosen to be stable and real but not at the origin, the performance clearly improves as expected. Twoexamples may be seen in Figures 4and 5 with roots of Aoat {0.85,0.85,0.85,0.85}and {0.92,0.92,0.92,0.92}respectively. WeA01.3 2544 0 50 100 150 200 250 0 5 10 15 20 r(t) % Model 69 r ref 0 50 100 150 200 250 0 5 10 15 20 25 u(t) µg kg−1 min−1 0 50 100 150 200 250 −8 −6 −4 −2 0 2 Time (minutes) ec=ref−r 0 50 100 150 200 250 0 20 40 60 80 100 Time (minutes) Switching Signal 44 Fig. 4. Results obtained for P=M69 with observer dynamics where the roots of Aoare {0.85,0.85,0.85,0.85}. 0 50 100 150 200 250 0 5 10 15 20 r(t) % Model 69 r ref 0 50 100 150 200 250 0 5 10 15 20 25 u(t) µg kg−1 min−1 0 50 100 150 200 250 −8 −6 −4 −2 0 2 Time (minutes) ec=ref−r 0 50 100 150 200 250 0 20 40 60 80 100 Time (minutes) Switching Signal 39 Fig. 5. Results obtained for P=M69 with observer dynamics where the roots of Aoare {0.92,0.92,0.92,0.92}. An extensivesimulation study using each model of the bank for mimicking the real plant Pwas conducted. This study,based on the mean square tracking error (MSE) and initial overshoot value, proved the superior performance of the controller designed when Aoroots ro i∈[0.65,0.95] ,i= 1,...,na=4(order of the neuromuscular blockade model). All roots of Aowere made equal. Figure 6shows the MSE as afunction of the roots ro iof the observer polynomial Ao for the worst case observed, i.e.,for P=M69. C. Clinical results The clinical test of the abovealgorithm was performed as part of astudy to evaluate and compare the performance of different control strategies, including SMMAC, approved by the Ethics Committee of Hospital Geral de Santo António (Porto, Portugal). In this framework, patients with health features levels Ito IV according to the American Society 00.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 3 6 9 12 15 Roots of Ao Mean Square Error 0.6 0.7 0.8 0.9 0 0.05 0.1 0.15 Fig. 6. Mean square error (MSE) for P=M69 as afunction of the roots ro iof the observer polynomial Ao.The upper right corner shows azoom of the MSE for 0.65 ≤ro i≤0.95. of Anesthesiology (ASA) haveundergone electivesurgery with automatic control of neuromuscular blockade. Anesthesia was induced with intravenous fentanyl and propofol.After calibration of the NMT module [17] of the DatexAS/3 Anaesthesia Monitor,a500 µg kg−1min−1 bolus of atracurium was administered. Anesthesia was maintained with AIR or N2O/O2,propofol infusion or sevoflurane and fentanyl as needed. Atracurium (1mg/ml)was delivered by acomputer controlled syringe driver (B|Braun Perfusor compact S) [18]. The controller is implemented in aportable battery operated computer that receives the muscle relaxation level signal via the RS 232C port, from the DatexAS/3 and updates the diffusion rate delivered by the pump with asampling rate of 20 s.The control algorithms were programmed in C-code, automatically generated from MATLAB code. Figure 7provides aviewof apatient in the operating room with the set-up for automatic neuromuscular blockade control (NMB). The NMB sensor is seen on the left hand of the patient. The syringe is seen mounted on the computer controller drive(left of the picture, of green color) as well as the portable computer (on the right). Figures 8and 9showclinical results with observer dynamics where the roots of Aoare {0.85,0.85,0.85,0.85} and {0.92,0.92,0.92,0.92}respectively.The last situation (with the observer roots at 0,92)has clearly asuperior performance, both in tracking the reference and in reducing the fluctuations of the manipulated variable. Again this illustrates the advantage of aslower observer. V.CONCLUSIONS The control of neuromuscular blockade during anaesthesia, through the continuous infusion of amuscle relaxant, is characterized by avery high degree of uncertainty in the dynamics of the system, due both to inter-patient variability as well as time variations. In order to deal with these features adaptation methods using multiple model based switching control was used. The selection of the observer dynamics was WeA01.3 2545 Fig. 7. Aviewof apatient in the operating room with the set-up for automatic neuromuscular blockade control. 0 50 100 150 200 250 0 20 40 60 80 100 r(t) % 0 50 100 150 200 250 0 5 10 15 20 25 30 u(t) µg kg−1 min−1 0 50 100 150 200 250 0 20 40 60 80 100 Time (minutes) Filtered r(t) % Filtered r ref 0 50 100 150 200 250 0 20 40 60 80 100 Time (minutes) Switching Signal 65 end of infusion Fig. 8. Clinical results obtained with observer dynamics where the roots of Aoare {0.85,0.85,0.85,0.85}. shown to be acrucial issue for achieving good performance and even stabilizing the loop. As expected from theory,fast observers yield oscillating behaviors and may even lead to wrong aselection of the local controller.Another important aspect is the way local controllers are designed which reflects in the overall performance of SMMAC. In this paper,the local controllers are PIDs tuned according to adominant pole placement rule. The clinical results confirmed the simulations and present agood performance. The algorithm was able to stabilize the patient dynamics without imposing anyapriori restrictions on the bank of models/controllers. REFERENCES [1] J. Baileyand W.Haddad, Drug dosing in clinical pharmacology: paradigms, benefits, and challenges, IEEE Control Syst. Mag.,vol. 25(2), 2005, pp 35-51. [2] Araki, M. and E. Furutani (2005). Computer control of physiological states of patients under and after surgical operation. Annual Reviews in Control,29:229-236. 0 50 100 150 200 0 20 40 60 80 100 r(t) % 0 50 100 150 200 0 5 10 15 20 25 u(t) µg kg−1 min−1 0 50 100 150 200 0 20 40 60 80 100 Time (minutes) Filtered r(t) % Filtered r ref 0 50 100 150 200 0 20 40 60 80 100 Switching Signal Time (minutes) 1 end of infusion Fig. 9. Clinical results obtained with observer dynamics where the roots of Aoare {0.92,0.92,0.92,0.92}. [3] Jaklitsch, R. and D. Westenskov (1987). Amodel based self adjusting two-phase controller for vecuronium-induced muscle relaxation during anesthesia. IEEE Trans. Biom. Eng.,BME-34:583-594. [4] Kansanaho, M. and K. Olkkola (1996). Performance assessment of an adaptivemodel-based feedback controller: Comparison between atracurium,mivacurium,rocuronium and vecvuronium.J.Clin. Monitor., Comput.,13(4):217-224. [5] Mason, D.; J. Ross, N. Edwards, D. Linkens and C. Reilly (1999). Self-learning fuzzy control with temporal knowledge for atracuriuminduced neuromuscular block during surghery.Comput Biomed. Res., 32(3):187-197. [6] D. Linkens, Intelligent Control in Biomedicine,Taylor and Francis, 1994. [7] T.Mendonça and P.Lago, PID control strategies for the automatic control of neuromuscular blockade, Control Eng.Pract.,vol. 6(10), 1998, pp 1225-1231. [8] J. M. Lemos and H. Magalhães and T.Mendonça and R. Dionísio, Control of neuromuscular blockade in the presence of sensor faults, IEEE Trans. Biomed. Eng.,vol. 52(11), 2005, pp 1902-1911. [9] S. Morse, Supervisory control of families of linear set-point controllers –Part 1: Exact matching, IEEE Trans. Automat. Contr.,vol. 41(10), 1996, pp 1413-1431. [10] K. Narendra and J. Balakrishnan, Adaptivecontrol using multiple models, IEEE Trans. Automat. Contr.,vol. 42(2), 1997, pp 171-187. [11] Mosca, E. and T.Agnoloni (2001). Inference of candidate loop performance and data filtering for switching supoervisory control. Automatica,37:527-534. [12] A. Neves and T.Mendonça and P.Rocha, “Tracking by switching control: acase study”, in Proc. 4th Portuguese Conference on Automatic Control (CONTROLO’00),Guimarães, Portugal, 2000, pp. 61-65. [13] P.Lago and T.Mendonça and L. Gonçalves, “On-line autocalibration of aPID controller of neuromuscular blockade”, in Proc. of the 1998 IEEE International Conference on Control Applications,Trieste, Italy, 1998, pp. 363-367. [14] T.Mendonça and P.Lago and H. Magalhães and A. Neves and P. Rocha, “On-line multiple model switching control implementation: A case study”, in Proc. 10th Mediterranean Conference on Control and Automation (MED’02),Lisbon, Portugal, 2002. [15] B. Weatherleyand S. Williams and E. Neill, Pharmacokinetics, Pharmacodynamics and Dose-Response Relationships of Atracurium Administered i. v., Br.J.Anaesth.,vol. 55, 1983, pp 39s-45s. [16] K. Åström and T.Hägglund, Automatic tuning of PID controllers, Instrument Society of America, 1988. [17] The Datex-Ohmeda Website: http://www.datex-ohmeda.com/clinical/cw_prev_01_article5.htm [18] The B|Braun Website: http://www.bbraunusa.com/products/infusion/perfusor_c1.html WeA01.3 2546