Design of Type-1 and Interval Type-2 Fuzzy PID Control for Anesthesia Using Genetic Algorithms

This paper presents the automatic drug administration for the regulation of bispectral (BIS) index in the anesthesia process during the clinical surgery by controlling the concentration target of two drugs, namely, propofol and remifentanil. To realize the automatic drug administration, real clinical data are collected for 42 patients for the construction of patients’ models consisting of pharmacokinetic and pharmacodynamic models describing the dynamics reacting to the input drugs. A nominal anesthesia model is obtained by taking the average of 42 patients’ models for the design of control scheme. Three PID controllers are employed, namely linear PID controller, type-1 (T1) fuzzy PID controller and interval type-2 (IT2) fuzzy PID controller, to regulate the BIS index using the nominal patient’s model. The PID gains and membership functions are obtained using genetic algorithm (GA) by minimizing a cost function measuring the control performance. The best trained PID controllers are tested under different scenarios and compared in terms of control performance. Simulation results show that the IT2 fuzzy PID controller offers the best control strategy regulating the BIS index while the T1 fuzzy PID controller comes the second.


Introduction
Anesthesia is a reversible state of people who are temporarily lack of consciousness in the purpose of undergoing a surgery without pain.General process of anesthesia consists of induction, maintenance, emergence, and re-covery [1].The depth of anesthesia (DOA) represents the level of consciousness [2], which is to be controlled through the process of anesthesia for varieties of control targets such as steady error, settling time, and overshoot.In general, control strategies can be categorized into two classes: open-loop control and closed-loop control.In open-loop control, based on the knowledge and experience, anesthetists manually adjust the drug dosage to maintain the DOA assisted by some clinical indices of patients.In closed-loop control, the drug dosage is automatically adjusted according to some indices of DOA, which makes control input continuous and responsive [3].Closed-loop control is also expected to avoid over-and under-dosage and suppress the adverse effect of interindividual differences [4].Despite of its advantages, the stability of closed-loop control needs to be ensured due to the automated process without supervision [5].
To achieve the stabilization based on closed-loop control, the following components are required [6]: 1) a patient model, with an output as the index of DOA; 2) a controller for stabilization, such as proportional-integraldifferential (PID) controllers.Therefore, to begin with, an estimated patient model is required to represent real patients supporting the control design and performance evaluation.A widely employed mathematical model consists of a linear pharmacokinetic (PK) model and a nonlinear pharmacodynamic (PD) model [7].The PK model uses the drug dosage as the input and the drug concentration as the output.In what follows, the PD model exploits the drug concentration as the input and exports the index of DOA.
As stated above, the control input of the overall patient model is the drug dosage.According to the types of drugs, the anesthesia can be classified to inhalational and intravenous anesthesia.Inhalational anesthesia has been used since the mid-19th century with for instance nitrous oxide, isoflurane, and halothane as commonly used drugs [8].In contrast, intravenous anesthesia with propofol and remifentanil as commonly used drugs [9] has the following advantages [10]: separate provision of anesthesia from ventilation, reduced atmospheric pollution, rapid and clear-headed recovery, and so on.While inhalational anesthesia is still frequently used for children, intravenous anesthesia becomes more and more popular due to the rapid and safe transition [11].
Main difficulties involved in modeling are the determination of PD model parameters and the selection of the output index.In general, the parameters of PK and PD models need to be determined beforehand.For the PK model, parameters can be estimated depending on sex, age, and weight of patients.Nonetheless, for the PD model, it is not possible to estimate the parameters for certain patient.Accordingly, it requires the controller to be robust in a domain of PD model parameters [12].Although PK model parameters are varied with patients and PD model parameters are even changing for one patient, a general index can be designed to evaluate the DOA for all patients [13].The bispectral (BIS) index [14] is an extensively accepted index to measure the DOA.It is a univariate dimensionless parameter from 0 to 100, which is devised from a large set of electroencephalogram (EEG) data [15].
With regard to the controllers for stabilization, several control problems are faced for the closed-loop control of anesthesia: stability, which is the basic control objective; robustness, which overcomes the uncertainty of PD model parameters, measurement noise, surgical stimulation, and so forth; adaptiveness, which makes the controller adaptive to different patients rather than only one patient.Commonly used controllers include PID controllers, model-based controllers, and knowledge-based controllers [6].
To close the loop of control systems, classical PID controller is a doubtless option due to its successful applications in other areas.For stabilization, it is necessary to tune the parameters of PID controllers for a certain patient model.While the tuning process is an empirical manner, some tuning schemes were even developed to guarantee the robustness [16] [17].Combined with patient model identification from the induction phase of anesthesia [18], the adaptiveness can be further achieved.In addition, combined with genetic algorithm (GA), the parameters of PID controllers can also be online optimized [19].
Despite of various causes of the changes in PK and PD model parameters, the model-based controller relies on the current model which reflects the patient's current pharmacological behavior.For this reason, the patient model needs to be updated through the overall process of anesthesia.This online adaptation can be achieved by several approaches such as Kalman filter algorithm [20], Bayesian-based adaptive control [21], and adaptive genetic fuzzy clustering algorithm [19].Alternatively, the variability of PK and PD models between individuals as well as surgical stimulation and anesthetic-analgesic interaction can be explicitly considered offline such that the stability is mathematically guaranteed [22].Additionally, a Lyapunov-based adaptive controller was developed to attain partial asymptotic regulation [23].
Unlike model-based controllers, knowledge-based controllers do not require a known mathematical model.In fuzzy logic controllers, for example, decisions are made based on fuzzy rules predefined by expert knowledge and experience.Compared with linear PID controllers and model-based controller, knowledge-based controllers are easier to implement without tuning PID parameters or mathematical derivation.A hybrid control scheme, namely fuzzy-PID control, was developed to combine the merits of both control strategies [24].Moreover, a multivariable neural-fuzzy controller was proposed to simultaneously administrate both propofol and remifentanil [25].However, one problem of knowledge-based controllers is that the interaction of each piece of knowledge makes the controller less transparent and makes it difficult to achieve adaptive control [6].To deal with this problem, a direct adaptive interval type-2 (IT2) fuzzy logic controller was proposed for multivariable anesthesia systems [26] and a genetic fuzzy logic controller was developed to adjust fuzzy rules using GA [19].
Due to the advantages of closed-loop control, in this paper, we aim to implement the closed-loop control for anesthesia model.The performance of controllers is ensured by simulation and optimization.For the anesthesia model used in this paper, while existing anesthesia model is utilized, the parameters of PD model are specially obtained from clinical data collected from 42 patients.Following the same design procedure, these tailor-made controllers can be redesigned for other patients as long as clinical data have been collected.In the anesthesia model, the co-administration of both propofol and remifentanil is investigated.To realize the drug administration, we propose two groups of strategies: two fuzzy PID controllers and one fuzzy PID controller with scaling factors.For each group of strategies, linear PID controller, type-1 (T1) fuzzy PID controllers and IT2 fuzzy PID controllers are designed.To draw a distinction from existing IT2 fuzzy PID controllers, we combine IT2 fuzzy PID controllers with GA.All PID gains, scaling factors and parameters of membership functions are optimized by GA in an offline manner subject to a performance index (cost function) which quantifies the performance of the controllers.A BIS training profile which sets different local targets for regulation considering the real anesthesia situation is employed for training purposes.The trained PID control strategies will be tested by a testing profile to verify their performance towards unseen working conditions.Comparisons are made among all controllers to demonstrate the characteristics of each control strategy.
Following sections are organized in this sequence: in Section 2, a multivariable anesthesia model used in this paper is described; in Section 3, the control background of linear PID controller, T1 fuzzy PID controllers and IT2 fuzzy PID controllers is introduced; in Section 4, we design all the control strategies and overall procedure and necessary information of the simulation; in Section 5, simulation results are provided and comparison and analysis are carried out; finally in Section 6, a conclusion is drawn.

Multivariable Anesthesia Modeling
In this section, the anesthesia model with two drugs is introduced.As shown in Figure 1, the anesthesia model consists of the target controlled infusion (TCI) system, the PK model, and the PD model, where are effect-site concentrations, and the subscript prop and remi are for propofol and remifentanil, respectively.Details of each module are described in the following subsections.The clinical protocol for obtaining the parameters of PD model is also described.

Pharmacokinetic Modeling
PK model is applied to reckon the relation from the drug infusion rate and the amount of drugs to plasmatic concentration ( ) Cp and effect-site concentration ( ) Ce .In this paper, we utilize the PK model from Marsh [27] and Minto [28] for propofol and remifentanil, respectively.A three-compartmental structure of these PK models is presented in Figure 2 which can be formulated by the following dynamic system: is the drug infusion rate; ( ) is the rate constant for distribution and elimination of the anesthetic drug; ( ) ( ) is a scalar describing the delay on time of the plasmatic-effect-site drug concentration equilibration.

Pharmacodynamic Modeling
After Cp and Ce are obtained from PK model, PD model maps Ce to the BIS value which reflects the anesthetic condition of patients.In general, PD model is established by a sigmoid max E model named as Hill curve for single drug.Based on the Hill curve, a response surface model was proposed by Minto et al. [30] to characterize drugs' interaction for multi-drugs.The PD model for representing propofol and remifentanil can be formulated as follows [31]: where 0 E is the clinical effect when no drug is infused; max E is the maximum drug effect; γ is the steep- ness of the concentration-response curve;θ is the ratio of the two drugs defined in Equation ( 3

Ce
, β and γ are pa- tient-dependent parameters which are estimated by clinical data explained in the following subsection.

Target Controlled Infusion System
TCI system transfers the plasmatic concentration targets ( ) Cpt to corresponding drug infusion rates.Practical- ly, the infusion can be realized by Alaris PK syringe pumps which incorporates the TCI system.These pumps are able to infuse adequate amount of drugs according to predefined Cpt .Furthermore, they can predict the realistic concentration based on PK model.
The TCI system consists of a bolus-elimination-transfer (BET) infusion scheme designed from the PK model so the Cpt is achieved.The initial bolus, loading dose (LD), is calculated according to Equation ( 6) where ( ) V l is the volume of central compartment represented on Figure 2. The LD bolus will ensure the Cpt is achieved however it needs to be followed by a continuous infusion, ( ) r t , defined by Equation ( 7) to ensure Cpt is maintained due to elimination and compartmental transfer of the drug [32].10 12 13 e e , where ( ) t s is the continuous time; ( ) is the rate constant for distribution and elimination of the anesthetic drug defined in Section 2.1.
TCI systems are also able to deal with Cpt changes.When Cpt increases an initial additional loading dose (ADDLD), Equation (6) and Equation (7) will become the following equations to maintain the new ( ) where τ is the time of Cpt change; ( ) The new continuous infusion rate formula takes into account the amount of drug already present in the peripheral compartments and which will be redistributed.
In the case when the Cpt decreases, no additional bolus is necessary, however the infusion of anesthetic drug is halted until NCpt is achieved on central compartment due to elimination of drug presented on the PK model.Once this NCpt is achieved the continuous infusion restarts according to Equation (9) where τ is the time NCpt is achieved on central compartment through the elimination process.

Clinical Protocol
The parameters of PD model in this paper are determined by clinical data obtained from 42 elderly patients who had major vascular surgery of the lower limb.The information of these patients is summarized as follows (in the format of "mean ± standard deviation"): 35 males and 7 females, age 72 ± 8 years old, height 169 ± 10 cm and weight 73 ± 14 kg.
The process of collecting data from these patients is briefly described here.Before induction, all sensors were settled and a radial arterial line was inserted to receive baseline readings on all clinical monitors.Two Alaris PK syringe pumps were employed to implement the total intravenous anesthesia (TIVA) with propofol [27] and remifentanil [28]  were adjusted by anesthetists according to readings of different clinical apparatuses to maintain the BIS value within [40,60] for an appropriate DOA.Among these readings, cardiac output (CO) was monitored and stabilized at or near pre-induction level by anesthetists through all periods of the surgery.

PID Controllers
In this Section, 3 types of PID controllers, namely linear PID controller, T1 fuzzy PID controller and IT2 fuzzy PID controller, are introduced.These 3 types of PID controllers are employed to regulate the output BIS index of the anesthesia model.

Linear PID Controller
A linear PID controller is shown in Figure 3, which consists of 3 elements, namely proportional, integral and derivative blocks.The output of the linear PID controller in discrete time is given as where , 0,1, 2, , proportional, integral and derivative gains, respectively, to be determined.

Type-1 Fuzzy PID Controller
In view of the linear PID controller [33], as the proportional, integral and derivative gains are constant, it is not able to handle well a highly nonlinear system.It motivates the use of T1 fuzzy PID controller [34] [35] of which the gains are changing according to the operating domains.By applying different sets of gains in different operating domains, more appropriate PID controller is employed to deal with the nonlinear system resulting in an improvement of control performance.
A T1 fuzzy PID control system is shown in Figure 4, which consists of a T1 fuzzy PID controller and a patient's model (detailed in Figure 1) connected in a closed loop.Unlike the linear PID controller having a set of constant gains, the T1 fuzzy PID controller has a fuzzy inference system providing a set of feedback gains through a reasoning process according to the operating condition.
The behavior of the fuzzy inference system is governed by a set of fuzzy rules of the following format: Rule : IF is AND AND is THEN , where ( ), 1, 2, , , is the fuzzy term corresponding to the linguistic variable j x in the th i rule;  is a positive integer; ( ) ( ) is the output of the fuzzy inference  system; and i y is the singleton membership function corresponding to the th i rule.The inferred output is gi- ven as where , where 0 p > denotes the number of rules; is the normalized grade of membership; is the grade of membership corresponding to the fuzzy term The output of the fuzzy inference system Equation ( 12) is employed to replace the gains of the linear PID controller turning it to become a T1 fuzzy PID controller.More precisely, 3 fuzzy inference systems are required to implement a T1 fuzzy PID controller.The outputs of the 3 fuzzy inference systems will be employed as P K , I K and D K .Consequently, the PID gains are no longer constant but dependent on the operating condition characterized by the membership functions.

Interval Type-2 Fuzzy PID Controller
Type-2 (T2) fuzzy sets [36]- [39] demonstrate a superior characteristic handing uncertainties compared with the T1 fuzzy sets.Uncertainties are captured by the lower and upper membership functions which form the footprint of uncertainty (FOU).A T2 fuzzy inference system can be considered as a set of infinite number of T1 fuzzy inference systems.Consequently, T2 fuzzy inference system is able to outperform the T1 fuzzy inference system in terms of reasoning and generalization capability.In general, the defuzzification process for general T2 fuzzy sets is computational demanding.By using IT2 fuzzy sets [40], the computational demand can be significantly reduced.
By employing interval fuzzy sets for the T1 fuzzy inference system Equation ( 12), it becomes an IT2 fuzzy inference system.The behavior of an IT2 fuzzy inference system is described by a set of rules of the following format: Rule : IF is AND AND is THEN , , , and denote the lower grade of membership, upper grade of membership, lower membership function and upper membership function, respectively.Using the center-ofset type reducer, the inferred output of the IT2 fuzzy inference system is given as , where and l r y y can be obtained using the Karnik-Mendel (KM) algorithms [41].The final defuzzifized output is given by ( ) ( ) The IT2 fuzzy PID control system [42] can be represented by Figure 4 as well of which the PID gains are given by the outputs of 3 IT2 fuzzy inference systems.

Simulation Design
In this section, the simulation environment is presented.First, we design the training profiles and testing profiles which all controllers need to deal with.The purpose of selected target profiles is first described.Then different control strategies are provided for comparison.For fuzzy control strategies, the approach to determine fuzzy rules is offered.Finally, the procedure of using GA for optimization is described, and a performance index is designed as the cost function for optimization.

Target Profiles
A training profile as shown in Figure 5, which is a series of BIS values required for an appropriate DOA, is employed for the training of the PID control strategies using GA.It defines the target BIS values in different periods that a PID control strategy has to achieve.After the training, a testing profile as shown in Figure 6 is employed to verify if the trained PID control strategies are able to control the BIS value to reach some targets subject to different operating conditions.
Clinically, the induction process starts with ( ) Since the induction period is short compared with the maintenance period, there is not much difference during induction in terms of control performance between various control strategies.Thus, we use two linear PID controllers for two drugs to drive the BIS from 98 to 50.By trial and error, the PID gains are predefined, and it is guaranteed that ( ) BIS 1000 is  As for the recovery process, since it is not allowed take the drug out of humans, the cease of controllers is the only and fastest option.These will not exist overshoot either because the target BIS is the maximum value.Therefore, all controllers will perform the same at this stage.It is not necessary to design controllers for the recovery process which is thus not included in the target profiles.It is noted that although the induction process is not under comparison either, we cannot ignore it.The reason is that it is difficult to find all model parameters describing the state of ( ) On the other hand, the parameters for ( ) are known and it is easy to start from this initial condition.

Control Strategies
In this paper, we aim to compare different control strategies to achieve a better control performance.All six control strategies are listed in Table 1.These six cases can be separated into two groups: two-controller Cases 1 to 3 and one-controller Cases 4 to 6.In two-controller cases, we control each drug by an independent PID controller offering control outputs α .Apparently, parameters of these two divided controllers are constrained by the ratio of these two factors.Theoretically, this constraint leads to conservativeness and make the performance worse than the ones of two-controller cases.Nevertheless, one-controller cases have less number of parameters to be determined resulting in lower computational demand, faster convergence in the training process and lower implementation cost for the PID control strategies.In each of these two groups, we have the following control strategies: linear PID controllers, T1 fuzzy PID controllers, and IT2 fuzzy PID controllers.For fuzzy PID controllers, the block diagram of BIS index regulation using two controllers and one controller with scaling factors are shown in Figure 7 and Figure 8, respectively.

Fuzzy Rules
T1 and IT2 fuzzy PID controllers are discussed in Subsections 3.2 and 3.3, respectively.While the number of fuzzy rules and the shape of membership function are predefined, the parameters of input and output membership functions are optimized by GA.In this paper, the shape of membership functions is defined as triangular shape for simplicity.The triangular IT2 membership functions are shown in Figure 9.For each input IT2 membership functions, the lower and upper membership functions are characterized by seven points 1 7 to p p which are to be optimized by GA.As for T1 fuzzy PID controller, each input membership function is characterized by three points 1 3 to p p which are to be optimized by GA.For both T1 and IT2 fuzzy PID controllers, their output membership functions are T1 and IT2 singleton membership functions whose values are to be determined by GA.
In fact, from the determination of membership functions and fuzzy rules, we can find the relation between these control strategies: solutions from linear PID controllers can be implemented by T1 fuzzy PID controllers, and solutions from T1 fuzzy PID controllers can be implemented by IT2 fuzzy PID controllers.In other words, linear PID controller is a subset of T1 fuzzy PID controller, and T1 fuzzy PID controller is a subset of IT2 fuzzy PID controller.By adding some constraints on the membership functions and consequents in fuzzy rules, IT2 fuzzy PID controllers can be reduced to T1 fuzzy PID controllers, and T1 fuzzy PID controllers can be reduced to linear PID controllers.
Three rules are employed for each fuzzy inference system.
. More rules can be used to partition the universe of discourse for a better result.However, it will lead to slower convergence of training and more computational burden.Additionally, the rate of change of ( ) can also be treated as another linguistic variable together with ( ) . Likewise, it will cause more fuzzy rules and parameters which increase the computational burden and defers the convergence.From the above discussion, triangular membership functions and three rules are employed for both T1 and IT2 fuzzy PID controllers.In the GA optimization, the triangular shape and sequence of membership functions of three rules should be guaranteed.For the sequence, specifically, the same point in membership functions corresponding to linguistic terms N, Z and P should be in ascending order (except 7 p ).For example,    Another condition needed to be ensured is that at least one rule is fired for = +∞ .It can be found that solutions from the first approach can be implemented by the second approach.In other words, the second approach is more general than the first approach.
Hence, we adopt the second method.For that reason, rules N and P are two shoulder-shape membership functions which can be treated as a special case of triangular membership functions.In this approach, it is guaranteed that at least one rule is fired for are temporarily employed in the optimization.

Parameters Optimization
The PID gains, scaling factors and membership functions are optimized by GA subject to a cost function reflecting the control performance.Due to the disparity of each training, GA is run for 10 times for each case of PID control strategies as shown in Table 1.The best set of parameters for each control strategy for each run of GA is recorded for further comparison, analysis and practical application.Statistical information including the worst, mean and best costs and the standard deviation for the 10 runs are collected.Among the 10 runs for each PID control strategy, the best set of parameters given by the best cost is used to implement the corresponding PID controller.During the optimization, the lower and upper bounds (LB and UB) of PID gains and scaling factors are listed in Table 2, which are determined by trial and error for good control performance.Due to a large number of variables and the highly nonlinear cost function (defined in the following subsection 4.5), especially for fuzzy PID controllers, GA may not be able to reach the global optimal solution.To speed up the training process, we define the initial population for fuzzy PID controllers based on our knowledge on the PID control strategy.It is known that linear PID controller is a subset of T1 fuzzy PID controller, and T1 fuzzy PID controller is a subset of IT2 fuzzy PID controller.It is thus reasonable that the best set of PID gains for linear PID controllers obtained by GA is employed as the initial PID gains for T1 fuzzy PID controllers for all rules.As a result, initially, the T1 fuzzy controller is equivalent to a linear PID controller.Similarly, the best set of PID gains and membership function parameters obtained for T1 fuzzy PID controller is employed as the initial set of PID gains and membership function parameters for IT2 fuzzy PID controller.That is to say, we utilize the "knowledge" on the PID control strategy such that GA starts with a considerable well initial condition.In this way, although there is still the same number of variables to be trained, GA only needs to utilize the additional parameters and new structures (fuzzy rules and membership functions) provided by fuzzy PID controllers in order to obtain a better cost value.

Performance Index
The performance index is used to judge whether the control objective is achieved and show the merits of various control strategies.Different performance indices can be selected such as settling time, overshoot, steady error, mean absolute error (MAE), and mean square error (MSE).For the GA optimization, a well-defined performance index is required as the cost function.All parameters of controllers are optimized according to the cost value provided by the cost function.In this paper, we present the following index J mainly based on MAE:  26), the first term is MAE which aims at minimizing the difference between the current BIS value and its target.Unlike settling time, overshoot and steady error, MAE and MSE record the er-  ror information through all simulation periods which reflect more comprehensive properties of performance.
The reason we choose MAE instead of MSE is that MSE amplifies the effect of large errors and then the optimization tries to minimize large errors which results in oscillation.On the other hand, MAE keeps the original weights on large and small errors leading to a mild and smooth response.
The three terms after the first term of ( 26) are the rate of change of ( ) ( ) remi k Cpt t , respectively, which are designed to reduce the oscillation for a gentle progress.The last term is to reduce the difference of the concentration of two drugs because they are equally important for different purposes during the anesthesia.Serious bias to either of these two drugs will not provide an adequate anesthesia.

Simulation Results
In this section, we implement the simulation on the control of anesthesia using the anesthesia model in Figure 1.Different control strategies in Table 1 are applied to regulate the BIS using the training profiles shown in Figure 5 for training and their control performance are verified by testing profiles shown in Figure 6.PID gains, parameters of scaling factors and membership functions are optimized by GA according to the cost function Equation (26).Comparisons of performance are made between the six cases.It should be noted that the simulation was carried out in discrete time.The PK model in Equation ( 1) was discretized by zero-order-hold (ZOH) assumption.The TCI system implemented has also been adjusted to perform under discrete PK model.
A real-coded GA available from Matlab global optimization toolbox is employed for the training.The control parameters of GA are listed in Table 3.
We define the parameters in Equation ( 26) as follows: 1000, . It is worth to mention that although the weight 1 λ is not the largest, the MAE is still the main contribution for the total cost.The bounds of parameters in Table 2 are adopted for the training process of GA.By running GA 10 times for each case in Table 1, we obtain the statistical information of the cost J as shown in Table 4.
Referring to Table 4, comparing with Cases 1 to 3, Case 3 offers the best cost and Case 2 comes to the second which coincides with the theory.The same rank can be found in Cases 4 to 6. Comparing Cases 1 and 4, twocontroller case is better than one-controller case, which can be also proved by comparing between Cases 2 and 5, and Cases 3 and 6.The reason is that one-controller case has the constraint on the ratio of control signal between two drugs while two-controller case does not have such constraint.One-controller case is a subset of two-controller case.Although one-controller case performs worse, it has less number of parameters to be determined and thus reduces the computational demand and implementation cost of PID control strategy.
In spite of the "Best" cost, the "Std" is descending from linear PID controllers to IT2 fuzzy PID controller except Case 4. It provides the information that the convergence gets easier for fuzzy PID controllers.The reason is that we use the best set of parameters of previous cases as the initial populations of following training.It also indicates that the improvement becomes more and more difficult, and thus 10 times running offer similar results.The "Std" for two-controller cases is larger than the corresponding one-controller cases, which implies that the convergence for two-controller cases is harder due to more variables to be trained.
The best sets of parameters and membership functions are shown in Tables 5-14.Corresponding membership functions are exhibited in Figures 10-13.From the obtained membership functions, it can be summarized that points and the maximum change of ( ) as an estimated domain where the PID controllers work suggesting that the fuzzy blending should start from around 20  ± .For the IT2 membership functions of IT2 fuzzy PID controllers in Figure 11 and Figure 13 after training, it was found that the lower and upper membership functions are very close to each other such as membership functions of N and Z in Figure 11.The values of points ( ) to p p associated with the IT2 membership functions are shown in Table 9 and Table 14.
By applying the best sets of parameters to the training profile of BIS, the time responses of BIS and corresponding drug concentration information are shown in Figures 14-19 for Cases 1 to 6, respectively.By applying the best sets of parameters to the testing profile of BIS, the ones for the testing profile are exhibited in Figures 20-25   +∞ lower and upper green lines, the performance will be well acceptable.
It can be seen from the figures that all cases of PID control strategies succeed in tracking the BIS targets and maintaining BIS values.In terms of control performance, while the linear PID controllers have less steady error, fuzzy PID controllers have less overshoot and settling time.Since there is no rigorous requirement on small steady error, fuzzy PID controllers are preferable that they are able to keep the BIS values within the regions bounded by the green lines with less overshoot and settling time.Both training and testing profiles demonstrate  the same characteristics as discussed for the corresponding control strategies.Recalling a common controller with two linear PID controllers and predefined parameters is used to achieve the induction process for all cases, the proposed PID control strategies kick in at the start of the maintenance process.The switching point of induction and maintenance process occurs at 1000 seconds which results in a small vibration.
The lower panels of Figures 20-25 show the behavior of drug concentration.It can be seen that both Cp and Ce of propofol and remifentanil are adjusted by controllers according to target profiles.More importantly, they are adjusted in advance to decrease the overshoot.This is mainly because all controllers are based on PID control strategy, whose abilities are inherited.In addition, the total amount of remifentanil is larger than propofol,                 which illustrates that the remifentanil is more favorable to the stabilization of BIS than the propofol.This might be due to the faster action time of remifentanil comparatively to propofol.Although the total amount of remifentanil is larger, the cost function takes the difference of the amount of two drugs into consideration such that serious bias to remifentanil is avoided.Even though the slightly bias towards remifentanil can be noticed overall the propofol and remifentanil Ce reached are within acceptable clinical boundaries even when the BIS target was set at 20, a value inferior to the safe interval of 40 to 60. Overall from a clinical point of view the results and performance obtained are acceptable, considering the lack of significant overshoot, the stability of the controller response and the ability of achieving the desired targets on a short time with clinically acceptable Ce .The controller also provided the expected response when the value of BIS target increases, and it stops infusion of anesthetic drugs but only restarts its infusion in other to avoid an overshoot.
When the PID control strategies with the best set of parameters (offering the best cost of J ) are employed to regulate the BIS value using the testing profile in Figure 6, their cost values are listed in Table 15.The rank of which is identical to the sequence of cost values from the training profile.It indicates that the advantage of fuzzy PID controller over linear PID controller exists not only in the training profile, but also in various target profiles.

Conclusion
In this paper, drug administration for anesthesia has been realized by various PID control strategies.We have first constructed a multivariable anesthesia model with propofol and remifentanil based on 42 patients' clinical data.Simulation has been conducted to regulate the output BIS value governed by this model using six control strategies including linear PID controllers, T1 fuzzy PID controllers and IT2 fuzzy PID controllers.These six strategies are separated into two groups, namely two controllers and one controller with scaling factors, to handle the co-administration of two drugs.Parameters of these controllers have all been optimized by GA subject to a performance index quantitatively measuring the control performance.In order to make the control sophisticated, target profiles have been designed for BIS regulation.Different target profiles have been utilized to test and verify the performance of controllers.It has been demonstrated that the IT2 fuzzy PID controllers offer the best performance, and T1 fuzzy PID controllers come to the second.In the future, further investigation can be carried out to control a general and non-parametric multivariable anesthesia model rather than the model with parameters obtained from specific group of patients.

Figure 1 .
Figure 1.A block diagram of multivariable anesthesia model.

Figure 2 .
Figure 2. Structure of a three-compartmental PK model with effect-site compartment [29].

(
the amount of anesthetic drug in the second compartment at the time of target change; of anesthetic drug in the third compartment at the time of target change.

Figure 3 .
Figure 3.A block diagram of linear PID controller.

Figure 4 .
Figure 4.A block diagram of fuzzy PID control system.

Figure 6 .
Figure 6.Testing profile.around50 and the maintenance process starts from In one-controller cases, however, the control output of the PID controller is divided into two control outputs by two scaling factors r α and p

P are points 2 p
for membership functions corresponding to lin- guistic terms N, Z and P, respectively.

Figure 7 .
Figure 7. BIS index regulation using two fuzzy PID controllers.

Figure 8 .
Figure 8. BIS index regulation using one fuzzy PID controller with scaling factors.

Figure 9 .
Figure 9.An example of IT2 membership functions.Dashed line: lower membership function.Dotted line: Upper membership function.Gray area: footprint of uncertainty.
functions, as long as they intersect each other, there exists one rule to be fired.On the left-hand side of membership function N and the right-hand side of membership function P, there are two approaches to ensure this condition.One is defining

5 .
be around 20 related to the training profiles we define in Figure Considering that the only linguistic variable is ( )

Figure 10 .
Figure 10.Membership functions for two T1 fuzzy PID controllers.Dashed line: membership function N. Dotted line: membership function Z. Solid line: membership function P.

Figure 11 .
Figure 11.Membership functions for two IT2 fuzzy PID controllers.Dashed line: lower membership functions.Dotted line: upper membership functions.Solid line: the shoulder of membership functions.

Figure 12 .
Figure 12.Membership functions for one T1 fuzzy PID controller with T1 fuzzy scaling factors.Dashed line: membership function N. Dotted line: membership function Z. Solid line: membership function P.

Figure 13 .
Figure 13.Membership functions for one IT2 fuzzy PID controllers with IT2 fuzzy scaling factors.Dashed line: lower membership functions.Dotted line: upper membership functions.Solid line: the shoulder of membership functions.

Figure 14 .
Figure 14.BIS and drug concentration for training profile by two PID controllers.

Figure 15 .
Figure 15.BIS and drug concentration for training profile by two T1 fuzzy PID controllers.

Figure 16 .
Figure 16.BIS and drug concentration for training profile by two IT2 fuzzy PID controllers.

Figure 17 .
Figure 17.BIS and drug concentration for training profile by one PID controller with scaling factors.

Figure 18 .
Figure 18.BIS and drug concentration for training profile by one T1 fuzzy PID controller with T1 fuzzy scaling factors.

Figure 19 .
Figure 19.BIS and drug concentration for training profile by one IT2 fuzzy PID controller with IT2 fuzzy scaling factors.

Figure 20 .
Figure 20.BIS and drug concentration for testing profile by two PID controllers.

Figure 21 .
Figure 21.BIS and drug concentration for testing profile by two T1 fuzzy PID controllers.

Figure 22 .
Figure 22.BIS and drug concentration for testing profile by two IT2 fuzzy PID controllers.

Figure 23 .
Figure 23.BIS and drug concentration for testing profile by one PID controller with scaling factors.

Figure 24 .
Figure 24.BIS and drug concentration for testing profile by one T1 fuzzy PID controller with T1 fuzzy scaling factors.

Figure 25 .
Figure 25.BIS and drug concentration for testing profile by one IT2 fuzzy PID controller with IT2 fuzzy scaling factors.

Table 1 .
Six cases of PID control strategies.
t denotes the desired BIS target value, ( ) The premise of each rule takes only one linguistic variable, i.e., only [  ]

Table 2 .
Lower and upper bounds of parameters. .

Table 3 .
Control parameters of GA.

Table 4 .
The cost J from running GA 10 times.

Table 5 .
Best set of parameters for two PID controller.

Table 6 .
Best set of parameters for two T1 fuzzy PID controllers.

Table 7 .
Best set of parameters for two T1 fuzzy PID controllers.

Table 8 .
Best set of parameters for two IT2 fuzzy PID controllers.

Table 9 .
Best set of membership functions for two IT2 fuzzy PID controllers.

Table 10 .
Best set of parameters for one PID controller.

Table 11 .
Best set of parameters for one T1 fuzzy PID controller with T1 fuzzy scaling factors.

Table 12 .
Best set of membership functions for one T1 fuzzy PID controller with T1 fuzzy scaling factors.

Table 13 .
Best set of parameters for one IT2 fuzzy PID controller with IT2 fuzzy scaling factors.

Table 14 .
Best set of membership functions for one IT2 fuzzy PID controller with IT2 fuzzy scaling factors.

Table 15 .
The cost J for the testing profile.