Real-time estimation of FES-induced joint torque with evoked EMG

Background Functional electrical stimulation (FES) is a neuroprosthetic technique for restoring lost motor function of spinal cord injured (SCI) patients and motor-impaired subjects by delivering short electrical pulses to their paralyzed muscles or motor nerves. FES induces action potentials respectively on muscles or nerves so that muscle activity can be characterized by the synchronous recruitment of motor units with its compound electromyography (EMG) signal is called M-wave. The recorded evoked EMG (eEMG) can be employed to predict the resultant joint torque, and modeling of FES-induced joint torque based on eEMG is an essential step to provide necessary prediction of the expected muscle response before achieving accurate joint torque control by FES. Methods Previous works on FES-induced torque tracking issues were mainly based on offline analysis. However, toward personalized clinical rehabilitation applications, real-time FES systems are essentially required considering the subject-specific muscle responses against electrical stimulation. This paper proposes a wireless portable stimulator used for estimating/predicting joint torque based on real time processing of eEMG. Kalman filter and recurrent neural network (RNN) are embedded into the real-time FES system for identification and estimation. Results Prediction results on 3 able-bodied subjects and 3 SCI patients demonstrate promising performances. As estimators, both Kalman filter and RNN approaches show clinically feasible results on estimation/prediction of joint torque with eEMG signals only, moreover RNN requires less computational requirement. Conclusion The proposed real-time FES system establishes a platform for estimating and assessing the mechanical output, the electromyographic recordings and associated models. It will contribute to open a new modality for personalized portable neuroprosthetic control toward consolidated personal healthcare for motor-impaired patients.


Background
Around 90 million people currently suffer from spinal cord injury (SCI) worldwide, and 85 thousand people each year survive the traumatic SCI and prepare to spend an average of 40 years or more in a wheelchair [1]. Various rehabilitation techniques and intelligent healthcare systems emerged in order to assist disabled people, releasing the heavy works from caregivers and reducing

France
Full list of author information is available at the end of the article economical/social cost [2,3]. Functional electrical stimulation (FES) is one of the neuroprosthetic techniques that can actively restore lost motor functions for motorimpaired subjects by delivering electrical current pulses to their paralyzed or weakened muscles. It may reduce patients' dependence on assistive devices and caregivers [4,5]. Under FES, motor units are excited synchronously [6] and action potentials are generated at motor neurons' axons level or muscles' fibers level, resultantly causing muscle contractions. The contractions transfer the generated force through the tendon and actuate related joints to achieve movement. FES has been applied to SCI patients, aiming to functionally assist their daily movements or training their muscles undergoing atrophy. For lower limb movements, FES is mainly dedicated to assist locomotion such as walking [7], standing [8], and correcting drop foot [9]. For upper limb movements, FES is used to help the patients to reach and grasp [10], hold [11], or lift target objects [12].
When utilizing FES, i.e. synchronous activation, the resultant joint torque was found to be highly correlated to evoked electromyography (eEMG), which brings possibility of modeling relationship between eEMG and joint torque [13][14][15][16]. Achieving accurate joint torque control using FES, is dependent on the modeling and estimation quality [17,18]. Additionally, explicit involvement of eEMG can distinguish the actual moment elicited by stimulus and by external forces. As a result, the issue of eEMG-based torque modeling and estimation play a critical role for further implementation of joint moment control by FES.
Volitional EMG based joint movement estimation approaches are widely used on healthy subjects and good movement prediction results are reported [19,20]. Different from able-bodied subjects [21], SCI patients are not able to generate volitional EMG due to damage of the descending efferent pathways at the spinal level. FES provides a way to synchronously recruit the motor units in place of the central nervous system. The recorded eEMG is then the summation of synchronous action potentials called M-wave. This M-wave can give pieces of information about muscle fatigue so that its tracking becomes possible [14,22]. Some representative works were proposed to predict muscle torque under metabolic fatigue during long and interrupted stimulation with multiple fatigue-recovery cycles [23][24][25], aiming at modeling longterm fatigue effect of stimulated muscles when metabolic fatigue is prominent. Different from the aforementioned works on prediction of muscle torque with long recovery from metabolic consumption, our previous works were focusing on FES-induced torque offline predictions with eEMG on SCI patients both with surface and implanted FES using Kalman filter as torque predictor [22,26], in order to predict FES-induced joint torque under potential early-middle stage fatigue without configuring recovery phase. Furthermore, we also compared different estimation methods to track muscle fatigue and improve its performance in a offline manner [27].
However, when transferring laboratory-environmentvalidated approaches into realistic clinical FES applications, real-time FES systems are required [7,28,29] and technical points should be considered. For instance, in real-time context, data acquisition, signal processing and estimation algorithms should be embedded and executed within time intervals defined by the close loop control, otherwise the effective real-time FES can not be achieved and the lag will degrade the control severely and even may induce instabilities. All the estimation approaches validated in offline are ultimately needed to be implemented in real-time FES systems to be assessed regarding their usability in a daily environment.
In this study, we develop a real-time system for estimating FES-induced torque based on eEMG for clinical scenarios. Efficient implementation of recursive Kalman filter and recurrent neural network (RNN) are used as the online estimators for identifying eEMGto-torque relationship. Prediction results on three ablebodied subjects and three SCI patients are presented to show the performance of real-time eEMG based joint torque estimation system with Kalman filter and RNN. This framework would also open a new modality for personalized portable neuroprosthetic usage as it can capture the subject-specific contraction dynamics by the tracking of muscle activation to joint moment relationship.

Methods
The system tests were conducted on three able-bodied subjects and three SCI patients upon their consent with a protocol (approved by Nîmes Ethics Committee, France, 2013). In this experiment, we abide by the norms of Helsinki Declaration of 1975, revised in 2000. The SCI patients' characteristics are shown in Table 1, and all the SCI patients have never experienced FES training until performing the experiment. Preliminary estimation results on one healthy subject through pre-testing the wireless FES system were presented in [30]. The system setup for the torque estimation experiment is shown in Fig. 1. The system consists of the wireless stimulator [31] (emanating from a transfer of our technology to Vivaltis Inc., Montpellier, France) as shown in Figs. 2 and 3, eEMG (MP100, Biopac Systems Inc., Santa Barbara, CA, USA), torque acquisition devices (Biodex 3, Shirley Corp.,  The dynamometer is measuring the ankle torque along with EMG measurement and the stimulation is applied through wireless stimulator to TA muscle NY, USA), and a laptop computer with the MATLAB (version 2012a) interface for remote control of wireless stimulator and real-time data acquisition/processing. The subject/patient was seated on the chair with the ankle at 90 degrees, while the joint center was aligned with the axis of a calibrated dynamometer. The shank was adjusted with the knee joint at 40 degrees. The foot was strapped to the pedal to transmit ankle torque to the dynamometer and to allow the optimal recording of isometric ankle torque.
The bipolar AgCl EMG electrodes were positioned over the muscle belly in the direction of muscle fiber with 20mm interelectrode spacing. The reference electrode was placed on the patella of contra lateral leg [27]. Electrical current pulses were delivered to the Tibialis Anterior (TA) or Medial Gastrocnemius (MG) muscle group with surface electrodes placed. Raw eEMG of TA or MG muscle group and ankle joint torque were recorded, amplified (gain 1000) and sampled at a frequency f samp = 4096 Hz Fig. 2 The wireless stimulator is a regulated-current dual-channel stimulator with the features as follows. Maximal electrical current is 100 mA, stimulation frequency ranges between 1 to 1000 Hz, and electrical polarity can be configured. All parameters are dynamically and remotely adjustable The stimulation frequency f stim was set between 30 and 40 Hz leading to the loop execution time between 25 ms and 33 ms. Truncation effect would appear at high stimulation frequency (up to 80 Hz) making the end of M-waves overlaps with the next forthcoming stimulation pulse. Such truncation effect becomes more obvious as fatigue level increases. The maximum pulse width (PW) of the stimulator was limited to 350 μs. During all FES experiments on both healthy subjects and SCI patients, the trapezoidal stimulation train of PW was used in all sessions for all individuals for torque estimation/prediciton. In different sessions with different individuals, the plateau amplitudes of trapezoidal PW trains were randomly and specifically generated with range limited to 50 to 100 % of maximum PW. The suitable stimulation intensities were found to be from 20 to 50 mA generally and specific to each patient. The stimulation frequency for healthy subjects was chosen as 40 Hz and for SCI patients was chosen as 30 Hz, since SCI patients may be more easy to be fatigue under a higher stimulation frequency. The stimulation intensity was varying from different one to different one, and it was also subject-specific for producing stable torque. For healthy subjects, the stimulation intensity was ranging from 20 to 40 mA with their TA muscles stimulated. For patients, the stimulation intensity is ranging from 35 to 50 mA, and we have used the protocol to choose the muscle which is giving stronger response, to avoid the torque measurement problem with too small torque scale which may happen for the patients. Thus, the stimulation intensity is fixed and PW modulation is performed to have different levels of muscle contractions to induce regular joint movement. The test session included 2 phases: identification and prediction phases. Each sequence contained trapezoidal trains consisting of 2 s stimulation (0.5 s ramp-up, 1 s plateau and 0.5 s ramp-down) and 2 s rest, with 18±2 rest cycles. During identification phase, the plateau stimulation PW of each trapezoidal train was increased gradually with step size of from 40 to 100 % of the maximum PW. After identification, the plateau stimulation PW was randomly determined within 50 to 100 % of maximum PW in the prediction phase. The stimulation artifacts in the raw eEMG were removed by blanking with a window size at 10 ms. Further processing are: 1) Calculation of differences between adjacent raw eEMG samples recorded during one loop (the number of samplings is around f samp /f stim ); 2) Comparison of these differences with a threshold, if the differences are larger than the threshold value, set the corresponding adjacent raw eEMG samples in the template window as zero. During every loop the mean absolute value (MAV) of eEMG and mean value of torque were computed based on the f samp /f stim raw eEMG and raw torque sampling respectively, and the MAVs of eEMG were averaged to be smoother (with a about 0.8s time length sliding window). The envelope of M-wave would be specifically different from session/subject to session/subject according to stimulation intensity, pulse width, and even stimulation electrode location. Therefore we computed the mean absolute value of M-wave sampling values recorded in a fixed time window in a batch and macroscopic manner in order to make them available for input of model identification.
Two online estimators are presented for the proposed torque estimation system. The identical models with the same orders are explored for each individual's eEMG based torque estimation. In our previous work [27], we used these two estimators to track FES-induced torque with eEMG in offline scenario and found both of them have promising performances. In this work, we online implement them toward real-time clinical trials on SCI patients to validate and evaluate their further performances. Both estimators are computed online and have been designed to be embedded, i.e., to be later implemented in the real-time architecture of the control unit. Based on previous off-line analysis and validation [27], both of the estimators adopt NARX-type architecture to describe the muscle dynamics as shown in Fig. 4.
Both estimators can guarantee their execution time within 1/f stim (between 25 to 33 ms) in every stimulation loop under MATLAB2012a interface environment to enable their real-time utilization. About more details on the estimation approaches, please refer to Appendix.   Table 2, the prediction performances are more stable for healthy subjects as compared to performances of SCI patients. Although EMG signals from SCI patient are more reliable and easy to identify as no volitional contraction that may change MAV appears, muscles of SCI patients become weak to make deleted motor units inducing less progressive and stable recruitment. Moveover, the extension levels of the lesion can be drastically different from patient to patient, estimators are needed to offer more patient-specific adjustment on SCI patients than on healthy ones. It can be seen as one of the main reasons why clinical trials with SCI patients are mandatory to really validate the FES system. Again observed from Fig. 8, we could see that the NARX-RNN can deal with SCI patients with different muscle strength and clinical symptom. Patient P1 had much stronger muscle strength and put up spasticity effect in certain posture during the experiment, while patient P2 had weaker muscle strength and produced ankle joint torque of lower level.

Results with discussion
In the current work, we have performed torque prediction tasks in isometric condition which is at neutral position of ankle angle is at 90 degree. This study is specifically focusing on the muscle response changes under FES, which can directly reflect time-variant muscle force changes. In dynamic motion situation, the torque can be altered as joint angle changes, so in non-isometric case muscle force-length and force-velocity relation should be involved. However, the potential muscle fatigue is originated from the muscle force only, thus torque prediction based on eEMG is necessary in isometric case before considering non-isometric case which needs measuring moment arm and estimating muscle force-length-velocity relation.
The concept of this work is to obtain personalized models for FES control purpose. When strong metabolic muscle fatigue arrives, FES control would not be feasible in any case in any FES system, as the muscle can not respond correctly to the stimulation. Thus, this study is targeted for early-middle fatigue stage. On FES trials among different days, even for the same subject, real muscle activation dynamics can differ a lot due to the change of stimulation contact situations. The stimulation contact situation could be drastically different due to stimulation electrode re-location. If the reality of stimulation is altered, it is necessary to make model re-identification. The proposed method can handle the different days case as well by applying this quick identification. Even with the good modeling, if the reality of stimulation is altered, we would need to make model re-identification in any case. Therefore, the quick on-line robust estimator such as RNN, could be very practical in the real world application in FES.
As compared with previous works on validation of eEMG based methods for FES-induced torque prediction from offline acquired data [22,27], this work focuses on applying the offline-validated torque estimators into a clinical application by considering realistic implementation factors, so as to condense suitable FES approaches in specific scenarios. It will contribute to transferring offline-validated FES techniques developed in laboratory environment to clinical applications which are mostly in real-time progress. It is worth mentioning here that, concerning adaptation to wider stimulation paradigms such as stimulation intensity and stimulation frequency as new inputs of model will lead to larger model complexity and longer identification process, and eventually may degrade online performance. It is not realistic in a clinical context and not necessary to provide excessive movements. However few extensions can be considered with dual different frequencies to produce more torque in specific conditions with limited models' complexity increase.

Conclusion
In this paper, a real-time FES-induced torque estimation system based on a wireless portable stimulator, is proposed and developed. Kalman filter and RNN are used online to identify and then predict torque output based on eEMG recordings only. Experiment validation on 3 able-bodied subjects and 3 SCI patients show the promising performance of the proposed real-time FES system. The real-time FES-induced torque estimation system can contribute to have personalized quantitative evaluation of muscle response under FES, which is first useful for clinical diagnostics to estimate the exact mechanical response by using eEMG signal. Further, it can be useful for the implementation of real-time closed-loop FES torque control, which is previously implemented in wired stimulators and EMG system for validation [18]. Attempting to get rid of torque sensor all the time will improve portability and mobility of ambulatory autonomous closed loop FES system. The next step on which we work is to perform wireless eEMG and embedded processing on a portable controller toward personalized neuroprosthesis, and more statistic analysis is to be addressed on the real-time FES system across a diverse population of SCI patients and/or under different stimulation regimes.

Kalman filter estimator
In our previous work, a polynomial Hammerstein model (PHM) was used as the estimator model which could represent muscle contraction. We defined contraction dynamics as the relationship between muscle activation and torque production [14,18]: where y(k) denotes model response output and u(k) denotes model input, a i , b i , and c j denote the model parameters to be identified, w(k) is the zero mean Gaussian white noise. In case of identifying relationship between stimulus and eEMG, u(k) denotes normalized stimulus PW and y(k) denotes MAV of eEMG; In case of identifying relationship between eEMG and torque, u(k) denotes MAV of eEMG and y(k) can denote joint torque information obtained through analog voltage signal (with the baseline-value effect eliminated).
The state space form of such PHM is described as follows: transforms the previous state vector x(k − 1) into the current state vector x(k). Coefficient B(k) denotes coupled linear and nonlinear parameter matrix associated with input u(k − 1).
Recursive Kalman filter with a forgetting factor λ ∈ [ 0.9 1] [22] is adopted as the online estimator in the realtime environment. One is estimator between stimulus and eEMG, and the other is between eEMG and torque [18]. Here, the recursive Kalman filter contains 2 parts, first one is a priori estimate phase and the second one is a posteriori update phase.

P(k) ← (I − K(k)H(k))P(k)
Such recursive Kalman filter can be implemented with the real-time measurement and acquisition system. It is worth noting that the actual obtained torque τ (k) (Nm) is scaled according to its linear relationship to the measured analog voltage y(k) (V) exported from the analog port of the Biodex chair, where the linear scaling relationship between torque and measured analog signal is τ (k) = 8.64y(k).

Online RNN-based estimator
The recursive NARX-RNN model which online describes muscle contraction dynamics is as follows: where k denotes the number of stimulation loops, and torque τ (k) can be as the acquired analog torque signal without scaling conversion or actual torque amplitude.
The parameters updating at every stimulation loop is governed by: . . .
where operator • denotes Hadamard product or Schur product between vectors. The predicted torque from k ≥ k id is represented by: From k > k id the torque model is no longer updated with the last estimated parameters at k = k id .
Actually, for RNN model (2) it might not be necessary to update its parameters at each stimulation loop k, although it can be implemented in an online recursive way like Kalman filter. In the identification phase, which is ended at k = k id , the latest parameters are estimated for the prediction phase. In this sense, the identification could only be done once at the stimulation loop k = k id . From this point of view, rather than updating its states every loop like Kalman filter, RNN could reduce a lot of computational load. Meanwhile, RNN may need less tuning the key parameters empirically in a trial and error manner like Kalman filter, which makes it more convenient for personalized FES systems. Moreover, Kalman filter can fail to work in non-normalization scale situation, but RNN may not suffer from such problem.