Computer Simulation Tests of Feedback Error Learning Controller with IDM and ISM for Functional Electrical Stimulation in Wrist Joint Control

Feedforward controller would be useful for hybrid Functional Electrical Stimulation (FES) system using powered orthotic devices. In this paper, Feedback Error Learning (FEL) controller for FES (FEL-FES controller) was examined using an inverse statics model (ISM) with an inverse dynamics model (IDM) to realize a feedforward FES controller. For FES application, the ISM was tested in learning off line using training data obtained by PID control of very slow movements. Computer simulation tests in controlling wrist joint movements showed that the ISM performed properly in positioning task and that IDM learning was improved by using the ISM showing increase of output power ratio of the feedforward controller. The simple ISM learning method and the FEL-FES controller using the ISM would be useful in controlling the musculoskeletal system that has nonlinear characteristics to electrical stimulation and therefore is expected to be useful in applying to hybrid FES system using powered orthotic device.


Introduction
Functional electrical stimulation (FES), which applies electric current or voltage pulses to peripheral nerves and muscles, is a method of restoring or assisting motor functions lost by the spinal cord injury or the cerebrovascular disease.FES has been found to be effective clinically, especially in controlling paralyzed upper limbs [1][2][3].For restoring lower limb functions, the hybrid FES system, which uses an orthosis with FES, has been accepted as one of practical methods [4,5].
In the recent years, powered orthotic devices or robotic exoskeletons have been focused on an assist or rehabilitation of lower limb functions [6,7].Therefore, the hybrid FES system is also expected to be realized with powered orthotic devices.In such system, cooperative control between FES and powered orthosis will be necessary.Feedforward control scheme would be useful for controlling fast movements of lower limbs in tracking to movements developed by the powered orthosis because control performance of a feedback controller is limited by large time delay and time constant in responses of electrically stimulated muscles.However, complex, time-consuming adjustment of many parameters of the feedforward controller such as creating stimulation data for a lot of muscles and time-varying properties of the musculoskeletal system make it difficult to use practically the feedforward FES controller in clinical application.
The Feedback Error Learning (FEL) proposed by Kawato et al. [8,9] can realize a feedforward controller by learning inverse dynamics of controlled object.The FEL will be useful in FES control because it can learn nonlinear characteristics of the musculoskeletal system to electrical stimulation and can remove the problem of manual adjustment of controller parameters by medical staffs in applying to various subjects that have different characteristics of the musculoskeletal system.
In order to apply the FEL, a feedback controller is required.The multichannel feedback FES controller has to solve the ill-posed problem in regulating stimulation intensities because the number of stimulated muscles is larger than that of controlled joint angles.The feedback FES controller based on the Proportional-Integral-Derivative (PID) control algorithm that we developed could provide a way of solving the ill-posed problem [10,11].In our previous work, the FEL controller for FES (FEL-FES controller) using the PID controller was found to be feasible in controlling 1-Degree-Of-Freedom (1-DOF) of wrist joint movement (dorsi-and palmar flexions) stimulating 2 muscles [12].
The FEL-FES controller makes it possible to use both the feedforward and feedback controllers, which is an advantage for the cooperative control between FES and powered orthosis in the hybrid FES system.Therefore, we performed preliminary test to expand the previous FEL-FES controller into controlling 2-DOF movements stimulating 4 muscles through computer simulation.However, the previous FEL-FES controller had a problem in learning the inverse dynamics.That is, learning the inverse dynamics model (IDM) in the previous FEL-FES controller sometimes failed.
Since a major problem in applying the FEL to FES is inappropriate learning of the IDM in FES control, a modification of the FEL-FES controller was discussed through computer simulation before testing with human subjects and applying the controller to hybrid FES system in this paper.In the previous FEL-FES controller, the IDM was only used for the feedforward controller since learning an inverse statics model (ISM) was not easy in clinical applications of FES because of difficulty in acquiring training data, while the FEL controller by Kawato was composed of the ISM and the IDM.
In this paper, in order to include the ISM into the feedforward controller, a simple measurement method of training data for the ISM was introduced considering FES applications.The ISM learning and the modified FEL-FES controller including the ISM were examined in wrist joint movement control by computer simulations in order to be compared to our previous work.

Feedback Error Learning Controller for FES
2.1.Outline.A block diagram of the feedback error learning controller for FES examined in this study is shown in Figure 1.The sum of output stimulation intensities from feedforward controllers (ISM and IDM) and a feedback controller is applied to each muscle after adding offset (threshold value of electrical stimulation intensity) and clipping out with the limiter to prevent excessive stimulation.
The PID controller outputs positive and negative values of stimulation intensity for each muscle to cancel out the difference between the desired joint angle (θ d ) and the actual angle (θ) during movement control.The outputs were also used in IDM learning on line.
Two three-layered artificial neural networks (ANNS) were used for ISM and IDM.The IDM and the ISM output positive values of stimulation intensity to each muscle calculated from the desired joint angle (θ d ), while the IDM uses the first and second derivatives of the desired angle.The  ISM is trained off line before IDM learning, and then the IDM is done on line using outputs of the feedback controller.

Feedforward Controller.
The structure of ANN for the IDM is shown in Figure 2. The input data of the desired joint angle and its first and second derivatives at continuous 6 times, from t to t + 5, (50 ms interval) in the directions of dorsi/palmar flexion (θ d1 , θd1 , and θd1 ) and radial/ulnar flexion (θ d2 , θd2 , and θd2 ) were given simultaneously.Outputs were stimulation intensities to 4 muscles.Therefore, the numbers of neurons in the IDM were 36 for the input layer and 4 for the output layer.That for the hidden layer was 18, which was determined based on our previous results [12].
The output of each neuron in the hidden and the output layers was defined as where x i represents outputs of the neurons in the previous layer, w i is the connection weight from neurons in the previous layer, c is the bias term, and i is the index of the neuron in the previous layer.The output function f (x) of the neuron is the sigmoid function ( The IDM was trained on line by the error backpropagation algorithm [13,14] using outputs of the PID controller.ANN connection weights are changed to reduce total error, E, as follows: where I desired and I IDM are desired stimulation intensity and stimulation intensity of the IDM, respectively.ε is the learning speed coefficient that has effect on convergence speed of learning.I desired − I IDM is approximated by stimulation intensity of the PID controller, I PID .
The ISM was trained off line before the IDM learning by using the error backpropagation algorithm.The threelayered ANN that had 2 neurons for the input layer, 18 and 4 for the hidden and the output layers, was used for the ISM.The ISM and the PID controller output stimulation during control for IDM learning, although outputs of the ISM were not used for IDM learning.

Feedback Controller.
The following PID control algorithm was used in the FEL-FES controller as the feedback controller: where the error vector e(n) is defined as difference between desired and measured joint angle vectors at time n.The PID parameter matrices K P , K I , and K D were determined by modifying the Chien, Hrones, and Reswick (CHR) method, and their elements were expressed as follows [10]: where L i and T i are the latency and the time constant of the step response of muscle i, when the response is approximated to the first order delay with latency.Δt is the sampling period.In case that a muscle has two or more functions (j shows index of the function), the delay time and the time constant obtained for every components in a movement were averaged, respectively.The coefficient m − i j corresponds to a reciprocal of the steady state gain of the system, which is calculated as an element of a generalized inverse matrix of a transformation matrix M. The matrix M transforms change of stimulation intensity vector into change of joint angle vector.Calculation method of the coefficient m − i j is shown in Appendix A.

Computer Simulation Tests
The FEL-FES controller including the ISM was tested in controlling 2-DOF movements of the wrist joint.The muscles to be stimulated were the extensor carpi radialis longus/brevis (ECRL/ECRB), the extensor carpi ulnaris (ECU), the flexor crpi radialis (FCR) and the flexor crpi ulnaris (FCU).The ECRL and the ECRB were assumed to be one muscle group (ECR) because of difficulty in selective stimulation to them in experiments using surface electrodes that we performed [10].
For computer simulation tests of learning the ISM and the IDM and of control performance, a musculoskeletal model of the upper limb was developed.In brief, muscle force F CE produced by electrical stimulation was described by the Hill type muscle model with nonlinear length-force relationship k(l) and nonlinear velocity-force relationship h(v), which included muscle activation level a m (s) determined by nonlinear recruitment characteristics with dynamics to applied electrical stimulation (refer to Appendix B for details).That is, (7) where s, l, and v were normalized stimulation intensity, muscle length and contraction velocity, respectively.F max showed a constant of maximum muscle force.Active torque τ CE produced by electrical stimulation was calculated by muscle force F CE and moment arm r f (θ).That is, Moment arm r f (θ) was represented by an approximated polynomial equation as a nonlinear function of joint angle θ for each movement developed by each muscle [15].Six different subject models were prepared, in which the difference between 6 subjects was represented by adjusting mainly parameters of recruitment characteristics based on step responses and input-output (stimulus intensity-joint angle) relationships of the muscles measured on 6 neurologically intact subjects.In this study, ISM learning was carried out off line using training data that consisted of stimulation intensities to 4 muscles and 2 joint angles.A set of training data was obtained by the tracking control of very slow movements using the PID controller.Figure 3 shows target trajectories of the tracking controls to obtain the training data set.The cycle period was 30 s for all trajectories.In Figure 3 6 s, were used for all trajectories.Six cycles were included in one control trial for IDM learning.Three sets of initial values of ANN connection weights were prepared, which were random small values that did not have effect on movements at the 1st control trial (before IDM learning).Therefore, a total of 45 learning tasks were tested on 6 subject models with all controllers (without ISM, using ISM trained with 2 trajectories, and using ISM trained with 4 trajectories).
Iteration number of IDM learning was fixed at 50.

Results
The ISM was evaluated by feedforward control of positioning.Target position for the control was set by a pair of dorsi/palmar flexion and radial/ulnar flexion angles at every In the case of using 2 target trajectories for obtaining training data (ISM-2), the error did not reduce around the center of the target trajectory and at positions between training data.As for the 4 trajectories for training data (ISM-4), the errors were small inside the largest target trajectory.Larger target joint angles outside the largest trajectory could not be controlled appropriately with both ISM-2 and ISM-4.Figure 6 shows average errors in open loop control of the positioning for ISM-2 and ISM-4.There was no large difference in the error between 6 subject models.Positioning errors shown in Figure 6(a) are for evaluation including targets outside the largest trajectory, and those in Figure 6(b) show those excluding targets outside the largest trajectory.Average positioning errors inside the largest trajectory (Figure 6(b)) were smaller than those in Figure 6(a).Figure 6(b) suggests that positioning in the radial/ulnar flexion was not trained sufficiently with the ISM-2.Figure 7 indicates an example of control result of the FEL-FES controller using the ISM with the IDM.The IDM was trained during the tracking control.The first cycle period of 5 s, which was set for moving to the start position of tracking control, was not used in the IDM learning.Before IDM learning (the 1st control trial), the ISM and the PID controller performed tracking control without the IDM.After IDM learning (the 50th control trial), the FEL-FES controller could perform good tracking with very small outputs of the PID controller.
In order to evaluate performance of the FEL-FES controller, mean error (ME) and power ratio (PR) shown in the following equations were calculated in each learning task:  where, e(n) represents the error between target joint angle and the resulted one at time n.N is the number of sampled data.P FF (n) and P FB (n) represent the output power of the feedforward and the feedback controllers, respectively.The ME was calculated for each movement direction, and the PR was done for each muscle.Average values of ME are shown in Figure 8.The controllers using the ISM decreased the error at the 1st control trial (before IDM learning).Especially, the ME was very small for slow movement control.All 3 controllers performed good tracking control after the IDM learning (the 50th control trial).There was no difference in ME after the IDM learning between ISM-2 and ISM-4 and also between with and without the ISM.
The power ratio, PR, gives us information of IDM learning.Figure 9 shows average value, the minimum and the maximum values of the PR.The FEL-FES controller using the IDM and the ISM achieved larger average value and larger minimum value of the PR than those of the previous controller before and after IDM  learning.After IDM learning, the minimum value of PR was greatly improved by using the ISM.There was no difference in those improvements between ISM-2 and ISM-4.

Discussion
The off line ISM learning was effectively achieved with the small number of measurements of training data.For practical clinical application, small number of measurements and short period of control time for acquiring the training data are required to avoid muscle fatigue and burden to patients.Therefore, training data acquired from feedback FES control of very slow continuous movements can be useful in ISM learning for FES.
Increasing the number of target trajectories to obtain training data may be required for learning the ISM of the musculoskeletal system that has nonlinear characteristics.However, if the ISM is mainly used to improve learning performance of the IDM, it is possible to decrease the  number of measurements of training data because there was no large difference between ISM-2 and ISM-4.On the other hand, target positions that had larger joint angles outside the largest trajectory could not be controlled appropriately as seen in Figure 5.This was a natural result because those targets were outside the training data.Since the control   than 84% of the number of muscle outputs showed the increase of PR for movements with the cycle period of 2 s.For movements with the cycle period of 3 s and 6 s, it was more than 65% and more than 40%, respectively.These results show that IDM learning was improved in most of learning tasks.For evaluating the improvement of IDM learning, the large PR rate that was defined as the percentage of the number of muscle outputs that had PR larger than 80% was calculated (Figure 10).The large PR rate was also improved by using ISM, especially for fast movement control.These results suggest that the FEL-FES controller using the ISM can be effective to realize a feedforward controller by learning nonlinear characteristics of the musculoskeletal system to electrical stimulation.For practical applications of the FEL to FES, an effective method of IDM learning will be needed, because the musculoskeletal system has nonlinear characteristics and also has hysteresis characteristics.
The FEL-FES controller using the ISM made better control with small values of ME at the first control trial for IDM learning as expected (Figure 8(a)).Since the difference in ME between with and without the ISM was not so large, the feedback controller was considered to perform well.However, control performance of the feedback FES controller sometimes deteriorated in tracking control because of nonlinear characteristics of the musculoskeletal system to electrical stimulation [16] although the feedback controller has been shown to perform properly [10,11].Therefore, the ISM is expected to become useful in controlling before IDM learning.
After IDM learning, all controllers showed small values of ME with no significant difference between the controllers (Figure 8(b)).However, the controllers using ISM resulted in larger average and minimum values of PR than those of the controller without the ISM (Figure 9(b)).This suggests that the PID controller had effect on decreasing errors for the controller without the ISM even after IDM learning while the feedforward controller worked mainly in the controllers

S S i max
Stimulus intensity to muscle i Joint angle of function j (deg) using ISM.Therefore, there is a possibility that the controller without ISM has a problem in movement control of the musculoskeletal system that has nonlinear characteristics.

Conclusions
Feedback error learning (FEL) controller using the ISM with the IDM was applied to FES control.The FEL-FES controller was examined in controlling 2-DOF movements of the wrist joint through computer simulation.In order to train the ISM in FES application, training data were acquired by controlling very slow movements with the PID controller.The ISM trained off line using the training data obtained by the simple measurement method was found to perform properly in the positioning task.The output power ratio of the feedforward controller in the FEL-FES controller was increased by using the ISM showing improvement of IDM learning.The FEL-FES controller using ISM would be useful in realizing feedforward controller for controlling musculoskeletal system that has nonlinear characteristics to electrical stimulation and therefore expected to be useful in applying to hybrid FES system.
the input-output relationship of the musculoskeletal system was represented approximately by using experimentally determined constant matrix M: In case of controlling 2-DOF movements stimulating 4 muscles, the following equation is obtained: where Δθ 1 and Δθ 2 show change of joint angles of dorsi/palmar flexion and radial/ulnar flexion, respectively.ΔS i means change of stimulation intensity to muscle i.The matrix M is not the square matrix in general because the number of muscles stimulated is larger than that of degree-of-freedom of movement controlled.Therefore, the generalized inverse matrix of the matrix M, M − , was calculated.That is, Since there are many generalized inverse matrices for M, the generalized inverse matrix M − has to be determined uniquely.
Here, after changing negative sign of m i j into positive one, the calculation of the generalized inverse matrix can be solved as the quadratic programming problem using (A.5) as the objective function under the constraints shown by (A.6) and (A.
This type of the quadratic programming problem can be converted to the linear programming problem by the Wolfe's algorithm [18].The unique solution of such linear programming problem can be obtained after the finite number of iterative calculations by the simplex method [18].That is, a set of positive values of m − i j minimizing the value L can be calculated under the condition of MM − = I after changing negative sign of m i j into positive one.Finally, the sign of m − i j was changed to negative sign based on the sign of m i j .

B. Musculoskeletal Model for FES Control
In this study, the 2-DOF wrist joint movements (dorsi/ palmar flexions and radial/ulnar flexions) were controlled stimulating the flexor carpi radialis (FCR), the flexor carpi ulnaris (FCU), the extensor carpi radialis longus/brevis (ECRL/B), and the extensor carpi ulnaris (ECU).Since the four stimulated muscles also relate to forearm or elbow movements, the skeletal model structure of the upper extremity was constructed in order to represent elbow flexion/extension, forearm pronation/supination, and wrist dorsi/palmar flexions and radial/ulnar flexions as shown in Figure 12.The shoulder joint was designed to be fixed at arbitrary angles of flextion/extention and rotation.The 15 muscles relating these movements as the agonist were included as listed in Table 1.Some muscles were also modeled as the synergistic muscles for other movements.
The musculoskeletal model to predict responses of electrically stimulated muscles is outlined in Figure 13.Muscle force F CE produced by electrical stimulation was described by the Hill type muscle model including muscle activation level determined by electrical stimulation a m (s), lengthforce relationship k(l), velocity-force relationship h(v), and maximum muscle force F max .That is, where s, l, and v were normalized stimulation intensity, muscle length, and contraction velocity, respectively.Active torque τ CE produced by electrical stimulation was calculated by muscle force F CE and moment arm r f (θ).That is, for each movement developed by each muscle [15].For example, the moment arm for the wrist dosri/palmar flexion and elbow flexion/extension was described by the following equation: where s c , s h , x c , and y c were constants.Electrical stimulation was expressed in normalized stimulation intensity s.The muscle activation a m was described by the following dynamics using the recruitment property with different two time constants, t r and t f [20]: The length-force relationship k(l) was described by the following equation.l o means optimum muscle length [21]: The velocity-force relationship h(v) during shortening and lengthening of muscle was modeled.v max shows maximum contraction velocity [21,22]: The maximum muscle force produced by electrical stimulation F max was determined by PCSA (physiological cross-sectional area) as follows [15]: (B.8) The passive viscoelastic element developed passive torque τ P calculated by the following equation for each joint movement [23].The range of motion was also represented by this property: where θ and ω were joint angle and angular velocity, respectively.Constants k 0 , b 0 , k 1 , and k 2 were determined for each joint movement.

Figure 1 :
Figure 1: Feedback error learning controller tested in this study.The inverse statics model (ISM) and inverse dynamics model (IDM) were used as the feedforward controller.

Figure 2 :
Figure 2: Structure of ANN for IDM used in the FEL-FES controller.

Figure 3 :Figure 4 :
Figure 3: Target joint angle trajectories to obtain training data for ISM learning.Cycle period was 30 s for all trajectories.

2
deg in the range of 20 deg in dorsi-and palmar flexions and in the range of 16 deg in radial and ulnar flexions.An example of the evaluation result of the ISM is shown in Figure 5.

Figure 5 :
Figure 5: An example of evaluation results of ISM in positioning control (model A).

Figure 6 :
Figure 6: Evaluation results of ISM learning in positioning control.

Figure 7 :
Figure 7: An example of control result by the FEL-FES controller using the ISM and the IDM.(model C, center at palmar position, cycle period of 3 s).

Figure 8 :
Figure 8: Average values of the mean error (ME) in tracking control by FEL-FES controllers.Error bar shows the standard deviation.

Figure 9 :
Figure 9: Average values of the power ratio (PR) in tracking control by FEL-FES controllers.Error bar shows the minimum and the maximum values of the PR.

Figure 10 :
Figure 10: Large PR rate after IDM learning for each FEL-FES controller.

Figure 11 :
Figure 11: Outline of determination of gain of the musculoskeletal system.The gray solid line shows the approximated linear line of the input-output relationship of the muscle.

Figure 12 :
Figure 12: Skeletal model structure of the upper limb.
r f (θ) = a 0 + a 1 θ + a 2 θ 2 + a 3 θ 3 + a 4 θ 4 + a 5 θ 5 (B.3)where a 0 ∼ a 5 were parameters for each movement of each muscle.Each element of the F CE is described in the following.The nonlinear recruitment property of electrically stimulated muscle u(s) was modeled by the following[19]:u(s) = s c tanh{s h (s − x c )} + y c (B.4)

Table 1 :
Muscles included into the model.
2)Moment arm r f (θ) was represented by an approximated polynomial equation as a nonlinear function of joint angle θ