Discrete Takagi-Sugeno Fuzzy Models Applied to Control the Knee Joint Movement of Paraplegic Patients

In this manuscript, a method for designing Takagi-Sugeno (T-S) fuzzy discrete-time regulators based on linear matrix inequalities (LMIs) is proposed to control the variation of the knee joint angle movement of paraplegic patients through electrical stimulation. A simple method for discretizing nonlinear systems described by T-S fuzzy models is used. The control strategy is applied for a paraplegic volunteer and a healthy one. The results and analysis show that the controlled system attended the design specifications for small values of the sample time considered for the discretization.


I. INTRODUCTION
Electrical stimulation has been used in conjunction with control systems to control the movement of the legs of paraplegic patients. For example, in [1], a Mamdani fuzzy controller was employed.
In this manuscript, a Takagi-Sugeno (T-S) fuzzy controller [2] is proposed to control the leg position of a paraplegic patient. The goal is to control the knee joint angle variation through electrical stimulation at the quadriceps muscle. The mathematical modeling of the shin-foot complex is proposed in [3]. This method relates the applied pulse width to the torque generated at the knee joint. In the control action, the leg moves out of its static angular position 0 • to the angular position of 30 • . When the electrical stimulation is disconnected, the leg returns to the initial angular position of 0 • . In the control action, the leg is brought to the desired angular position, maintained at a balanced state, and then returned to the initial position.
The use of the model developed in [3] is adequate for the design of the controller, considering that the relationship between the torque and the width of the electric pulse applied to the muscle is nonlinear and complex, and its mathematical model is difficult to obtain [1], [3].
The associate editor coordinating the review of this manuscript and approving it for publication was Xiaoli Luan.
In [4], [5], and [6], function electrical stimulation and continuous-time T-S fuzzy models were for the first time used to study and simulate the control of leg position for a paraplegic patient [2]. In the problem studied in [2], for a given operation point, the maximum, and the minimum values of the nonlinear function were considered in a desired interval.
This method allows the exact representation of a class of nonlinear systems through T-S fuzzy models. Other operation points can be considered in the design procedure presented in [7].
The concept of parallel distributed compensation (PDC) [8] is used in the design of fuzzy regulators to stabilize nonlinear systems described by fuzzy models. The idea is to design a compensator for each rule of the fuzzy model, where a controller is associated with each rule. The resulting global fuzzy regulator, which is nonlinear, in general, is a fuzzy combination of the individual linear regulators. The PDC offers a procedure to design a regulator for each T-S fuzzy local model, where each control rule is designed from the correspondent T-S local model of the plant.
The controller design can be realized using linear matrix inequalities (LMIs) [8]. The software Matlab, with the LMI control toolbox [9], can solve the LMI-described problems when they are feasible. Relaxed conditions can be used [10], but in this manuscript, LMIs presented in [8], which are less relaxed, were used, giving good results. In [4], [5] and [6], LMI-based design control systems are used in conjunction with electrical stimulation to vary the knee joint angle of patients.
Most of the dynamic systems are described by continuoustime models. However, instead of implementing analogic controllers, implementation using a digital dispositive is preferable as it features greater flexibility and lower cost. In nonlinear controller's design, considering the plant described by T-S fuzzy models, several results have been published, as in [11]. The discretization of nonlinear systems is not always easy. In [12], [13] and [14], linear digital controllers designs are described.
A technique to discretize linear systems with emulation is suggested in [12]. This technique can also be used to discretize T-S fuzzy nonlinear systems in a way less complicated than the method presented in [11]. In Theorem 1, the validity of the discretization methodology with T-S fuzzy models to design discrete-time controllers for the rehabilitation of paraplegic patients using electrical stimulation is shown in a way more simplified than in [15].
This method and its analysis with simulations show the validity. To the knowledge of the authors, scientific studies on the application of discrete time T-S controllers to control the leg position of paraplegic patients have not yet been performed.
Recently, few articles have been published in this area, such as [16] and [17]. In [18] and subsequent articles [19]- [22], and [23] published by our research group, some theories are employed to control knee joint movement using neuromuscular electrical stimulation. In [24], an article was published describing control application in paraplegic patients.
In the present study, to improve the control techniques of previous works with T-S fuzzy models [17], [25], an embedded fuzzy logic was implemented in a digital signal processor (DSP), with closed loop.
Other articles on the human muscles command issue have been recently published, such as [16], [26]- [32], and [33]. The model of [16] uses a linear and a nonlinear block of the state space model to apply variable structure control. The presented model considers factors such as fatigue, spasticity, and retraining effects.
Models [4] through [17] work with nonlinear state space modeling. Although our model does not consider all the effects considered in [16], the T-S fuzzy models utilized in the design follows the recommendations of [3] to reduce the fatigue effect by the variation of pulse width and the stimulus application time instructions with their appropriate reapplication intervals. Therefore, the implemented fuzzy logic has a good response to the idealized design [4].

II. DYNAMIC MODEL USED IN THE CONTROL OF THE LEG POSITION OF A PATIENT
The mathematical model of the leg employed in this manuscript is proposed in [3]. This model relates the applied pulse width with the torque generated on the knee joint. In the modeling [3], the leg was considered as an open kinematic system composed by two stiff segments: the thigh and the shin-foot complex, as shown in Fig. 1. From [3], the equilibrium equation, around the knee joint is given as where • J is the inertial moment of the shin-foot complex; • θ is the knee angle between the shin and the thigh on sagittal plane; •θ =θ v is the knee angular velocity; •θ =θ v is the knee angular acceleration; • m is the mass of shin-foot complex; • g is the gravitational acceleration; • l is the distance between the knee and the shin-foot complex mass center; • B is the viscous friction coefficient; • M s is the torque due to the stiffness component; • M a is the knee active torque generated by electrical stimulation.
The stiffness moment is defined as: where λ and E are the coefficients of the exponential terms, and ω is the elastic rest angle of the knee. In [3], it was observed that the torque applied to the muscle (M a ) and the electrical stimulation pulse width (P) can be adequately related by the transfer function described in (3), where G and τ are positive constants obtained from parametric identification when a stimulus is applied to quadriceps and angular variation is read, as mentioned in [3]: VOLUME 8, 2020 FIGURE 2. Graphics of the exact functionf 21 (x 1 ) and its approximation using an 11 th -order Taylor series, as done in [17].
Considering (1)-(3), the resulting state space representation is as follows [4], [5]: where by definition, θ v0 is the desired final position and (θ v ,θ v , M a ) = (θ v0 , 0, M a0 ) is an equilibrium point of this system. The goal of the tracking control is to assure that this equilibrium point of the controlled system is asymptotically stable and thus The functionf 21 (x 1 ) is a nonlinear parameter of the system and is written as Expanding (5) in Taylor series, as done in [17], the term x 1 at the denominator can be removed, avoiding the indetermination problem at x 1 = 0. From Fig. 2, one can observe that the 11 th -order Taylor series gives a good approximation of the exact curve of the nonlinear function at the interval If a more exact representation is necessary, the Taylor series order can be increased.

III. GENERAL FUZZY T-S REPRESENTATION
Certain classes of nonlinear systems with T-S fuzzy models can be exactly represented using the method described in [2]. According to this construction method, the local models are obtained in the function of the operation region.
For i = 1, 2, . . . , r, the obtained i-th rule of the continuous-time T-S fuzzy model is expressed below: The resulting fuzzy model is the ponderated mean of the local models, given bẏ where: ,

IV. DISCRETIZATION PROBLEM DESCRIPTION
In digital controller design, the conventional methodology can be used for discretization of linear time-invariant plants, as described below:ẋ If an input u(t) is generated by a digital computer, followed by a digital-analog converter, where T is the sample time, one can define: For equation (7), considering (8), the discretization of the model, according to [12], is given by Therefore, from (9), one has: In (10), if the sample time is small enough, one can despise the high-order terms of the Taylor series, and, consequently, = I. An intuitive idea to discretize the continuous T-S fuzzy model is to directly apply the discretization method for linear time-invariant systems at the consequent parts of the fuzzy rules [11].
For i = 1, 2, . . . , r, the i-th rule of the discrete-time T-S fuzzy model correspondent to the continuous-time local T-S model (6) is shown below: As shown in [11]: Considering the fuzzy model of (11), one has, for definition, M i j , j=1,2,. . . ,p as the fuzzy set j of the rule i; the premise variables x 1 (kT ), . . . , x p (kT ), and the membership function µ i j (x j (kT )) of the fuzzy set M i j is given as Finally, one has: where: From these definitions, the discretized form for the T-S model is found: The analysis of (17), as well as (12) can be found in [11], where the nonlinear behavior, and the fuzzy rules are considered. Here, the proposed discretization method is presented. Consider the continuous-time T-S fuzzy model, given bẏ Theorem 1: Consider the continuous-time T-S fuzzy model, described in (18). Thus, for a small enough sample time T, (18) can be represented by the discrete-time T-S fuzzy model, presented in (17), where for i = 1, 2, . . . , r. Proof: Suppose that the sample time T is small enough, the following approximation holds: So, from (18): Considering (16), (20), and (21) note that, The proof is concluded.

V. NONLINEAR DYNAMIC MODEL OF THE LEG OF THE PARAPLEGIC PATIENT USING TAKAGI-SUGENO FUZZY MODELING
From the continuous model presented in (4), one can find a discrete state space model. The design of the continuous case can be seen in [4] and [5]; the design considers the operation points at θ v0 = 30 • and is based on the T-S modeling proposed in [2]. Consequently, (4) is decomposed into local models, and its membership functions, considering the maximum and minimum values of the nonlinear functioñ , are a 211 = 5.5506 and a 212 = −35.9035. Thus, the membership functions are formalized as [8] where the local models (4) are: From Theorem 1, the discretized local fuzzy models with sample time T are given by: VOLUME 8, 2020 Finally, the discrete T-S fuzzy model is obtained, using local models (24) and the membership functions (23): α i (x(kT ))(L i x(kT )+V i u(kT )), whose coefficients α 1 (x(kT )) and α 2 (x(kT )) are constant into the intervals [kT , kT + T ), but change their values at each sampling instant.
Remark 2: According the T-S fuzzy theory in [2], the number of local models depends on the number of nonlinear parameters. In paraplegic patient case described in (4), this equationing contains one nonlinearity. Thus, as in [2], one defines two local models. In (24), the terms a 211 and a 212 are the maximum and the minimum values off 21 (x(kT )), respectively.

VI. REGULATORS WITH TAKAGI-SUGENO FUZZY MODELS
Similarly to the continuous case, the PDC [8] in discrete fuzzy regulators design can stabilize nonlinear systems described by fuzzy models. The idea is to design a compensator for each fuzzy rule. For each rule, there exists an associated controller. The resulting global fuzzy regulator, which is nonlinear in general, is a fuzzy combination of each individual linear regulator. The PDC offers a procedure to design a regulator for each T-S fuzzy model, where each control rule is designed from the corresponding plant T-S model rule. The designed fuzzy regulator shares the local controllers, each one given by The fuzzy global regulator is given by: where α = [α 1 , . . . , α r ] T . So, the closed-loop system is composed in the form below: The goal of the fuzzy regulator design is to determine the local feedback gains F i at the consequent parts for each fuzzy rule i. Therefore, from (27), the matrices Q ij are defined as Then, from (16), (23), and (25), note that Therefore, where, for definition: For the controllers' design, the concepts about LMI in the discrete-time case were used.

3)
Q ij +Q ji 2 T P Q ij +Q ji 2 − P ≤ 0, for i < j. Replacing (28) at items 1, 2, and 3, one gets the LMIs to be implemented in a computational program to determine the controller gains using relaxed stability conditions. Therefore, the optimization problem to determine the controller law is to find matrices X = X T > 0, Y = Y T ≥ 0, V i , and M i , for i = 1, . . . , r, satisfying [8]: where s is the number of rules that can be active [8]. When the LMIs are feasible, the positive definite matrix P and the controller gains F i are given by:

B. RELAXED STABILITY CONDITIONS USING CONTROLLER DESIGN WITH DECAY RATE
Similarly to items 1, 2, and 3 of Subsection VI-A, [8] presents LMIs with relaxed conditions that increase the feasibility area of the system. Defining X = P −1 , M i = F i X, and Y = XWX, the fuzzy control design, considering decay rate (that is related with the setting time of the system), for the discrete-time case is to find P = P T > 0 and W = W T ≥ 0 such that 1) P > 0, P = P T ; 2) Q T ii PQ ii − σ 2 P + (s − 1)W < 0, for i = 1, 2, . . . , r; 3) where β = σ 2 is the decay rate, 0 ≤ β < 1, and s is the number of rules that can be active. For instance, the design to solve the movement control of paraplegic patients contains two active rules. In the following, the LMIs that solve items 1, 2, and 3 above using Schur complement [34] are presented. From [8], these LMIs are  If (33) and (34) are feasible for i, j = 1, 2, . . . , r and X = X T > 0, the controller gains can be found by the equation F i = M i X −1 , i = 1, 2, . . . , r.

C. INPUT CONSTRAINTS
Consider a known initial condition x(0). From [8], the constraint ||u(kT )|| 2 ≤ µ is guaranteed for all time t = kT ≥ 0 if the LMIs To guarantee the input constraint, LMIs (35) and (36) need to be added to the LMIs that guarantee stability, given in (31) and (32), or decay rate, given in (33) and (34).
For the equilibrium point study case, conservative LMIs were used to solve this problem, with a constant matrix P of the Lyapunov function V (x) = x T Px for this modeling. When a time-variant matrix P was used, more relaxed LMIs were obtained, as given in [35].
Therefore, for this modeling, the LMIs in articles [36]- [39], and [40] can adequately solve the conditions with a time-variant matrix P.

A. THEORIC CASE SIMULATION BASED ON THE MODEL GIVEN IN [3]
The system described in [3] presents two local models, as shown in (24), with the parameters given in Table 2: At T-S modeling, with θ v0 = 60 • , LMIs (33)-(36) are feasible for X = X T > 0, for T = 0.001s, β = 0.999, and µ = 0.005. The conservative and relaxed forms of the LMIs for resolution presented in [8] sufficiently solve this kind of problem, but the methods presented in [10], with more relaxed conditions, can be used. The gains F 1 and F 2 are In Fig. 3 the simulation results are shown and one can observe that the approximation (22) is valid. In Table 1, one can see the expression of matrix in (10), for the sample time values T = 0.1 s, T = 0.01 s, and T = 0.001 s. As the value decreases, the matrix gets closer to the identity matrix. Therefore, the discretization error is negligible if the sample time is very small. For greater periods, imprecisions occur due to the approximation of the equation model (21); that is, system (22) corresponds to a discrete emulation of nonlinear systems. As T increases, the discretization error may affect the system stability. For the system considered in this subsection, the maximum period such that LMIs (33)- (36) are feasible is T = 0.5214 s. Observe that this maximum value of T is much larger than the pulse width applied to the patient. Thus, due to VOLUME 8, 2020   [21].

TABLE 4.
Identification of the parameters [21], according to the model in [3]. the nature of the system, the sampling time is much smaller than the maximum value for feasibility of LMIs (33)- (36).
The simulation result shows the following: a steady-state angular position of 30 • , null angular velocity, and active torque of 4.6068 Nm. Therefore, the design result followed the stability and performance specifications required for the discrete system.
In the text with β, knowing that its value is between 0 and 1, as shown in [8], the authors observed in the simulations that the variation of β from 0.9 to 0.9999 produced different con-  troller types; therefore, the most adequate value was chosen as β = 0.999.
In [3], the authors suggest methods of experimentally obtaining the parameters of interest. In this manuscript, the same parameters in [3] were adopted, as given in Table 2, where θ v ,θ v , M a = (θ v0 , 0, M a0 ) is an equilibrium point.
In Table 3, the characteristics of the volunteers are presented.
The parametric identification methodology presented in [3] is applied for the volunteers. The anthropometric parameters of the volunteers are given in Table 4. Considering the parameters given in Table 4, each system presents two local models, as shown in (24), with the parameters given in Table 4.    For the paraplegic patient, the frequency used was 50 Hz, because it is an intermediate frequency that recruits both fast fibers as well as slow fibers.
Then, the discrete-time T-S fuzzy controller is computed such that items 1, 2, and 3 of Subsection VI-B are feasible; that is, LMIs (33)-(36) is feasible for X = X T > 0.
Considering β = 0.99 and T = 0.001 s, the gains for the healthy volunteer (H1) for a frequency of 77 Hz are  Fig. 4 shows a photo of the neuromuscular electrical stimulator used. It has eight channels and can provide current of up to 140 mA [20], [21].
The patient remained sitting on the instrumented chair, implemented at the Advanced Control Robotics and Biomedical Engineering Laboratory, and the chair allowed the application of the stimulus on quadriceps muscle, which is shown in Fig. 5. A bar with an electrogoniometer, an accelerometer, and a gyroscope gives the measures of the variable signals θ, θ v , and M a [19]- [21], [24].
Considering the state feedback T-S design, the controller was developed in Simulink, as shown in Fig. 6. This figure contains the controller generated from LMIs (33)- (36) using the Matlab Control Toolbox of Mathworks R12.
The controller design implemented at Simulink was compiled in embedded codes for utilization on the digital signal processor DSP TMS320F28335 kit, from Texas Instruments [19]- [21]. The procedure for parameter estimation, according to [3], as well as the T-S fuzzy control design, is shown in [21].
Two cases were analyzed to validate the proposed model, for one healthy, and one paraplegic subjects.   Although the transient behavior for the simulated and the implemented systems seem to be different, the controllers used were the same. This is possibly because the used model does not consider the following cases: • instrumented shaft, since that the paraplegic patients show stiffness, because of involuntary spasms, depending on the instant and the physiologic conditions of the patient; • the fiber recruitment delay. In future works, it is advisable to consider the spasm and delay in the recruitment of fibers, since these factors make the initial movement a bit slower, as mentioned by Lynch and Popovic in [16]. The shaft dynamics also must be considered.
Although some initial differences in the transient behavior occurred, after a certain time, the implemented model had the same behavior as the simulated system. VOLUME 8, 2020 For volunteer H1, repeatability tests were made for the fuzzy controller, considering the muscle relaxation time between one application and another. Fig. 9 shows the muscular behavior in the repeatability test using a frequency of 77 Hz.
As seen in the Figures 7, 8, and 9, one can observe that for the paraplegic patient, the frequency used was 50 Hz, and for the healthy person, the frequency was 77 Hz. For the paraplegic patient, 50 Hz was used because it is an intermediate frequency that enables the recruitment of both fast and slow fibers. This can be explained by the fact that before the accident, the paraplegic patient was triathlete and the healthy person had a large circumference on the femoral thigh and a sedentary lifestyle. The paraplegic patient maintains the physical activities, and he has a great muscle tone.
The frequencies of 50 and 77 Hz were experimentally optimized. For the example, with the healthy person, 70 Hz was used for the first test, 75 Hz for the second test, and 80 Hz for the third test. The best results occurred between 75 and 80 Hz; thus, 77 Hz, which is an intermediary point between these values, was adopted.
From the simulation results, one can observe that for values of T below 0.01s, the idea of emulation was evident, where the discretized nonlinear system worked as a linear discretized system, with = I , despising the superior-order terms and minimizing the approximation error. Therefore, T = 0.001 s was adopted. A more complex model was studied by [11]; this theory greatly simplifies the problem.

VIII. CONCLUSION
A simple discretization method of nonlinear systems, described through continuous-time T-S fuzzy models, is proposed, and the method is demonstrated to be valid. The LMIs, with relaxed conditions proposed in [8], which consider stability, decay rate, and input constraints, were used for the design of discrete-time T-S regulators. Other LMIs, with more relaxed conditions can also be used [10]. Digital simulations show the validity of the proposed method at discrete-time nonlinear control of the leg position of a paraplegic patient through electrical stimulation. The T-S fuzzy controller presented excellent results when it was applied to the volunteers, with a final angle near 30 • , and a transient behavior similar that obtained in the simulation results.