MEFFNet: Forecasting Myoelectric Indices of Muscle Fatigue in Healthy and Post-Stroke During Voluntary and FES-Induced Dynamic Contractions

Myoelectric indices forecasting is important for muscle fatigue monitoring in wearable technologies, adaptive control of assistive devices like exoskeletons and prostheses, functional electrical stimulation (FES)-based Neuroprostheses, and more. Non-stationary temporal development of these indices in dynamic contractions makes forecasting difficult. This study aims at incorporating transfer learning into a deep learning model, Myoelectric Fatigue Forecasting Network (MEFFNet), to forecast myoelectric indices of fatigue (both time and frequency domain) obtained during voluntary and FES-induced dynamic contractions in healthy and post-stroke subjects respectively. Different state-of-the-art deep learning models along with the novel MEFFNet architecture were tested on myoelectric indices of fatigue obtained during ${a}\text {)}$ voluntary elbow flexion and extension with four different weights (1 kg, 2 kg, 3 kg, and 4 kg) in sixteen healthy subjects, and ${b}\text {)}$ FES-induced elbow flexion in sixteen healthy and seventeen post-stroke subjects under three different stimulation patterns (customized rectangular, trapezoidal, and muscle synergy-based). A version of MEFFNet, named as pretrained MEFFNet, was trained on a dataset of sixty thousand synthetic time series to transfer its learning on real time series of myoelectric indices of fatigue. The pretrained MEFFNet could forecast up to 22.62 seconds, 60 timesteps, in future with a mean absolute percentage error of 15.99 ± 6.48% in voluntary and 11.93 ± 4.77% in FES-induced contractions, outperforming the MEFFNet and other models under consideration. The results suggest combining the proposed model with wearable technology, prosthetics, robotics, stimulation devices, etc. to improve performance. Transfer learning in time series forecasting has potential to improve wearable sensor predictions.


I. INTRODUCTION
D ETECTION of muscular fatigue is essential for sportspersons and athletes [1], [2], workers of firms [3] or construction sites [4], patients undergoing rehabilitation (such as stroke and spinal cord injury (SCI) patients) [5], etc. Overuse of fatigued muscles can lead to muscular injury in both healthy individuals and patients undergoing rehabilitation.It could also lead to absenteeism, poor turn-ups, and lack of enthusiasm for future training sessions.Therefore, there is a significant need to foresee the indicators of muscular fatigue for preventing its ill effects [6].
Muscle fatigue is often characterized as a reduction in maximum force or power output [7].In sports or in any physiotherapy treatment which involves repetitive movements either voluntary or FES-induced, it becomes difficult for both trainer/ coach/ physiotherapist and trainee/ patient to determine the optimum effort or optimum stimulation dose without causing fatigue or muscular injury.Athletes can maintain high-end performance by preventing injuries such as of anterior cruciate ligament with proper fatigue management [8].Previous studies [9], had discussed how training on fatigued muscles lead to muscle strain injuries and may lead to further distressing effects in healthy individuals such as changes in control and coordination, postural stability, and kinematics (particularly decreased range of motion and lifting velocity) [10].Patients undergoing low back pain rehabilitation may display healthy equivalent kinematic performance even with residual impairment in strength [11].However, the unevaluated fatigue and continued activity may lead to pain and alteration in neuromuscular control and coordination system [10], [12].In stroke rehabilitation [13], the post-stroke adversaries involve central and peripheral mechanisms which limit the activity of daily living.To overcome the post-stroke challenges, different rehabilitation modalities are implemented.One such modality is the use of FES to stimulate paralysed muscles.However, the benefits of FES in post-stroke individuals are hindered by the detrimental effects of stimulation-induced muscle fatigue [14], [15].
Non-assessment of muscular fatigue and status of muscular health leads to over-exercise or overstimulation of fatigued muscles causing an adverse effect rather than benefiting the concerned person.Since, the transition from non-fatigued to fatigued state of muscles is gradual and varies from person to person [16], muscular fatigue assessment in advance can therefore help in preventing fatigue related injuries.
Muscular fatigue assessment, classification, and prediction are predominantly obtained using surface electromyography (sEMG) signals.A multitude of information about muscle health can be obtained from the sEMG signals [17], [18].It being a reliable source of fatigue assessment for decades, showed advantages over ultrasound-based, accelerometerbased, and torque-based assessments in terms of portability, accuracy, and feasibility respectively [6].Usually, the sEMG based assessment of fatigue is done using myoelectric indices [18].With recent advancements and the existence of wearable sensors the fatigue assessments through these indices have become simpler.
Although the data from wearables are accessible throughout the entire duration of a training/therapy/exercise session, achieving an instantaneous estimation/assessment of fatigue is challenging since the transition to fatigue is not instantaneous rather progressive [19].Therefore, continuous monitoring along with futuristic estimations are important in determining the failure point in muscle fatigue assessments [20].
Previous works that had done forecasting of sEMG features in healthy subjects were limited to isometric contractions-induced muscle fatigue [21], [22].The shallow network using Normalized Least Mean Squares (NLMS) algorithm and the stacked dilated convolution neural network (CNN) had been proved effective in forecasting the features of trunk muscle fatigue up to a forecasting window of 25 seconds in healthy subjects [21], [22].Whereas in rehabilitation domain, such as for post-stroke and spinal cord injury population, muscle fatigue prediction and estimation are associated with FES-induced contractions, where stimulation-induced muscle fatigue is a big performance deterrent [23].In this regard, a past study utilized a nonlinear ARX -recurrent neural network (NARX-RNN) based model to estimate torque from evoked EMG feature-mean absolute value (MAV) for muscle fatigue tracking [23] with an application to customize the patient-specific stimulation parameters under muscle fatigue.However, till now the studies have covered only the hindsight of muscle fatigue, and a noticeable gap exists in terms of foreseeing muscle fatigue in FES-induced contractions.Since, evoked EMG can manifest both electrical and mechanical behaviors of a muscle [24], it can be used to extract its time and frequency domain features to indicate muscle fatigue [18].These features are advantageous as they reflect various information about the muscle status directly [18] unlike the joint angle and torques which measure only the biomechanical outputs produced by muscles under fatigue.
The features of sEMG which could be forecasted, such as the time domain indices of fatigue [18]-root mean square (RMS), average rectified value (ARV), and frequency domain indices of fatigue [18]-mean frequency (MNF), median frequency (MDF), usually follow rising (time domain indices) or falling trend (frequency domain indices) along with seasonality (periodically varying component) and irregularities during fatigue progression in repetitive dynamic contractions.
This typical characteristic (trend and seasonality) could be seen in both voluntary [25] and FES-induced contractions [26].The presence of dependencies of future values on past values and non-linear trend in these signals during dynamic contractions make accurate predictions challenging for long time horizons [27].Additionally, the periodicity of seasonal components in these indices depends upon the time duration within which one cycle/repetition is completed.Therefore, the seasonal component of sEMG signals evolve dynamically with repetitions and fatigue in both voluntary [25] and FES-induced contractions [26].Hence, seasonal decomposition techniques, which assumes constant seasonality in a time series, cannot be applied to reduce time series complexity [28], [29].In these scenarios, where the time series is complex, forecasting models such as Autoregressive Integrated Moving Average (ARIMA) and Exponential Smoothing (ETS) turns out to be less efficient than machine learning and deep learning models that can handle complex non-stationary time series [30].
However, the various machine learning and deep learning architectures are problem specific, such as, stacked dilated causal CNN for forecasting of isometric trunk-muscle-fatigue features for lower back pain management [21], forecasting of exercise-induced fatigue indices from inertial measurement unit (IMU) signals through generative adversarial network (GAN)-based networks [31] etc.Their viability in forecasting of myoelectric indices of fatigue pertaining to dynamic voluntary and FES-induced contractions is still unknown.Therefore, these available deep learning architectures were studied in this work from the perspective of cross-domain implementation of deep learning models along with the novel MEFFNet which was specifically designed for forecasting of myoelectric indices of muscle fatigue (both time and frequency domain) during voluntary and FES-induced dynamic contractions for healthy and post-stroke subjects respectively.
The proposed MEFFNet model takes into consideration the evolving seasonal pattern along with the changing trend observed in the myoelectric indices of fatigue.It consists of a convolutional subnetwork and a bidirectional long short-term memory (BiLSTM) subnetwork combined with an attention mechanism.The presence of convolutional subnetwork extracts the important features from the time series while the attention mechanism gives preference to the most important features.The subsequent BiLSTM subnetwork learns the long term dependencies among the features to provide the future predictions in the subsequent layers.
Additionally, the MEFFNet model was trained on a large dataset of synthetic time series and was named pretrained MEFFNet.The pretrained MEFFNet uses a transfer learning strategy by being trained on a synthetic dataset.It was subsequently used to predict myoelectric fatigue indices to improve the forecasting performance compared to MEFFNet.A break of 20 minutes was given before proceeding further with another weight.The maximum flexion angle and pace at which the activity was performed were subject-specific.
The major contributions of the present study are: • A deep learning forecasting model MEFFNet is developed which is aimed to forecast myoelectric indices of fatigue that could inherit short or long range periodic component making it useful in both voluntary contractions and FES-induced contractions.
• The performance of MEFFNet is compared with the existing forecasting models in the muscle fatigue domain, such as LSTM, BiLSTM, CNN, stacked dilated causal CNN [21], and GAN-based architectures [31].
• Furthermore, the pretrained MEFFNet has been introduced with an objective to implement transfer learning in deep networks to forecast time series with small sample size, thereby capturing the complex patterns and helping in reducing overfitting.
• An ablation study on the input and forecasting window lengths has been presented which helped in deciding the optimum number of timesteps that could be forecasted.

A. Participants and Experimental Protocols
Task: Elbow flexion and extension movements were chosen for this study since it being the most common upper limb exercise for healthy individuals and post-stroke patients [32], [33].Ethical clearance for these experiments was approved by the All India Institute of Medical Sciences ethics committee (Ref.No. IEC-299/07.05.2021).
1) Voluntary Dynamic Contractions: Sixteen healthy, righthanded male subjects (mean ± SD, age: 30.69 ± 3.14 years) with no neurological disorders and muscular injuries, participated in the experiment.They performed non-stop biceps curls (elbow flexion and extension) under four weight conditions (1 kg, 2 kg, 3 kg, and 4 kg) while seated on a chair (Fig. 1 (a)).Each weight condition involved a single 7-minute trial with a 20-minute break before proceeding with other weight for removing the effect of fatigue from previous weight lift [34].Participants were motivated to target at least 20 cycles per minute to induce fatigue; however, this criterion was not made stringent [34], [35].One cycle was defined as one elbow flexion and extension performed while carrying weight from the reference position to the maximum flexion angle (subjectspecific) and back to the reference position (Fig. 1 (b)).A cycle was not delimited within time through any visual or verbal cue, rather the subjects were instructed to maintain their initial angle of maximum flexion.Participants who could not perform activity till 7 minutes with the highest weight (4 kg) were removed during initial screening.During the activity electromyography data were recorded from ten muscles using sEMG sensors.These muscles were: brachioradialis, biceps brachii, and brachialis (the main flexors of the elbow joint); triceps short head and triceps long head (the main extensors of the elbow joint [36], [37]); and deltoid medial, deltoid anterior, deltoid posterior, trapezius, and anconeus (the synergist muscles).The elbow joint angle was measured using a goniometer.
2) FES-Induced Dynamic Contractions: Seventeen poststroke and sixteen healthy subjects (other than those who participated in voluntary dynamic contractions) were recruited to perform FES-induced elbow flexion in vertical plane while sitting on a comfortable chair (as shown in Fig. 2 (a)) under three different stimulation patterns/templates as shown in Fig. 3.These patterns were customized rectangular, trapezoidal, and muscle-synergy based which were adopted from our previous work [26], where these patterns were applied to the subjects to explore the differences in the stimulation-induced muscle fatigue.Detailed information on these stimulation patterns, Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

TABLE I SUBJECT DESCRIPTION
Fig. 3.This figure shows all the stimulation patterns used for stimulating biceps brachii and brachioradialis.The first stimulation pattern is a customized rectangular stimulation pattern with an initial ramp-up stage followed by a 2 s hold duration; the second is a trapezoidal stimulation pattern with equal ramp-up and ramp-down time, it contained a 2 s peak stimulation hold time too, while the third one is a muscle-synergy-based stimulation pattern which also contained a 2 s peak stimulation hold time followed by a two-tier decreasing pattern.The peak current to achieve maximum elbow flexion was subject specific.
their construction, and pattern-dependent differences on myoelectric indices of muscle fatigue can be found here [26].
The seventeen post-stroke patients were selected from sub-acute phase (significant functional recovery is observed during this phase [38]).These subjects complied with the initial screening.The subject description for post-stroke patients and healthy participants can be found in Table I.These experiments were approved by the All India Institute of Medical Sciences ethics committee (Ref.No. IEC-299/07.05.2021).Participants' biceps brachii and brachioradialis muscles were stimulated simultaneously to flex their elbow vertically [26].Doublet stimulation pulses, which are two closely spaced electrical pulses (8 ms apart) delivered in rapid succession for additional muscle twitches [39], were delivered with a pulse width of 200 µs [40], [41] at 40 Hz frequency [39], [41] to generate forceful and tetanic contractions which are desirable in weak muscles of the paralysed limb [39].Each stimulation pulse was symmetric and biphasic to eliminate the undesirable charge stored in the tissues [39].
The maximum stimulation current required to achieve full elbow flexion was determined through ramp signal.The participants were exposed to a gradually increasing ramp current on both biceps brachii and brachioradialis muscles, and the current level at which full flexion was achieved was determined as the maximum stimulation current.However, the participants could bear higher current levels but the initial level that caused full flexion was chosen as the maximum stimulation current.This maximum current varied across the participants as shown in Table I.For each participant the stimulation current was the same for biceps brachii and brachioradialis muscles.By further increasing the stimulation current, the level barely causing pain was recorded as the pain threshold.The current injection was limited to 35 mA to ensure safety [26], [40].The maximum stimulation current and pain threshold were assessed for each patient/participant prior to the FES sessions.
Participants were administered with conventional stimulation patterns/templates-customized rectangular and trapezoidal as well as muscular synergy-based stimulation template tailored according to their respective maximum stimulation current to achieve full elbow flexion.The various stimulation templates are shown in Fig. 3.The sequence in which different templates were administered was randomly generated via program.The goal of each stimulation pattern was to induce maximum elbow flexion.35 cycles of elbow flexion were completed during each stimulation template as per the time flow shown in Fig. 2 (b).A 10 s inter-cycle break was asserted between consecutive cycles of elbow flexion.After completion of 35 cycles under a particular stimulation template, a 1-hour 30-minute break was given to reduce muscular fatigue before proceeding towards the next stimulation template.FES stimulation parameters were regulated by MATLAB 2019 b.The participants did not hold any weight in this experiment.

B. Instrumentation and Data Acquisition
The sEMG data during voluntary and FES-induced dynamic contractions were acquired using DELSYS Trigno Wireless EMG sensors at 2000 Hz and the elbow angle was recorded using goniometer and data logger (Biometrics Ltd., UK) at 2000 Hz.For FES-induced contractions, two pairs of oval shaped RehaTrode FES electrodes (4 cm × 6.4 cm) were used to stimulate the flexor muscles.Each target muscle was connected to one pair of FES electrodes placed 20 mm apart.The EMG system, goniometer, and FES device were time synchronized with a common start and stop trigger using DELSYS Trigger Module and MATLAB 2019 b.Skin preparation before placement of sEMG sensors, goniometer, and FES electrodes was done by hair removal and cleaning with ethanol.The sEMG sensors and the FES electrodes placement were done using standard protocols [42], [43].Stimulations were given by Hasomed Rehastim2 under expert's guidance.The evoked EMG and angular displacements were continuously and simultaneously acquired while stimulations were delivered.

C. Data Processing and Feature Extraction
The data processing was done on MATLAB 2019b.For voluntary dynamic contractions, the raw sEMG data was filtered using 4 th order zero phase lag bandpass Butterworth filter in the frequency range [20,450] Hz to remove low frequency motion artifacts and high frequency noise.The goniometer data was filtered using 2 nd order low pass IIR filter at 10 Hz to remove high frequency noise.The pace at which the activity was performed was subjective and hence goniometer data was used to segment each cycle in the sEMG data.For FES-induced dynamic contractions, the evoked EMG (eEMG) was contaminated by the stimulation artifacts of magnitude several times higher than the M-waves [44].In our previous work [26], empirical mode decomposition (EMD) was found to be more efficient in preserving signal content than blanking window, comb filters, and band-pass filters, which were used for filtering evoked EMG.The EMD based filtering involved decomposing of the eEMG signal into a series of intrinsic mode functions (IMFs).Through selective removal of unwanted IMFs, EMD was used for stimulation artifact removal.First and foremost, the raw eEMG was filtered between [20,150] Hz using 4 th order Butterworth bandpass filter.Subsequently, IMFs were obtained and were inspected for the specific IMFs that contained artifacts and noise.The IMFs related to stimulation artifacts contained harmonics of stimulation frequency (40 Hz), while IMFs related to noise contained high frequency components.The undesired IMFs were removed, and the remaining ones were used to reconstruct the eEMG signal.The percentage reconstruction error of EMD for raw electrically evoked EMG was negligibly small.The EMD filtering methodology [44] has been explained mathematically in detail in the supplementary material.
The time domain features [18] -ARV and RMS, and frequency domain features [18] -MNF and MDF were extracted from the filtered sEMG data on 0.5 s window which moved with an overlapping percentage of 25 (explained in supplementary).Extracted features from voluntary contractions and FES induced contractions are shown in Fig. 4 and Fig. 5 respectively.

D. Techniques for Deep Learning 1) Windows Slicing Techniques and Variables for Window
Size Effect: The most crucial aspect while making a time series forecasting model is the choice of input and forecasting window lengths, which needs optimization to deliver efficient forecasting results.Previous work [31] had shown that an input to forecasting horizon ratio of 1500: 80 timesteps generated optimum forecasting and fatigue classification results with IMU in dynamic contractions.This ratio was taken as a cue in this analysis.However, the input window length was not increased above 600 timesteps, i.e., 225.125 s or say 3.75 minutes to meet with the demands of dynamic contractions, which Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.2) Different Time Series Forecasting Models: To identify the most suitable forecasting model, the following models were tested using real-time series data of myoelectric indices derived from spontaneous and FES-induced dynamic contractions.To begin, a Vanilla LSTM model was utilized, comprising a single LSTM layer with 650 hidden units, ReLU activation, time distributed layer, and a subsequent dense layer with neurons equal to the length of the forecasting horizon.Further, LSTM and BiLSTM models each with two layers and 648 hidden units, ReLU activation, time distributed layer, and subsequent dense layer was implemented separately.Next, a CNN model was constructed, comprising two Conv1D layers with filter sizes of 256 and 648, kernel size of 2, stride of 1, and ReLU activation followed with global pooling, time distributed layer, and dense layer.Next a GAN-based forecasting framework was also explored where the generator was based on a transformer, and the critic consisted of a network composed of three convolution layers followed by a dense layer, as detailed in the prior work [31].Additionally, a stacked dilated causal CNN [21] was used.The model consisted of stacked dilated causal convolution layers where the number of convolution layers varied with the input training length.For a training window of n timesteps the number of convolution layers were found as log 2 n.The dilation factor progressed in a manner of 2 n−1 for n convolution layers.The kernel size and filter size were chosen as 4 and 2 respectively [21].
a) MEFFNet: MEFFNet is an attention based deep CNN-BiLSTM model consisting of a convolutional subnetwork and a Bi-LSTM subnetwork.The convolutional subnetwork comprised of a 1-D convolutional layer of 256 filters followed by batch normalization, max pooling, and attention layer.This subnetwork extracted the important features from the time series which can handle the evolving seasonal component in a time series.Its output was fed into the Bi-LSTM subnetwork having two hidden layers of 750 units, batch normalization, dropout, and attention layer which was then passed on to a fully connected network of dense layer which yielded a sequenced output vector of forecasting window length.The network architecture can be found in Fig. 7. Detailed description about the layers can be found in the supplementary material.
b) Synthetic dataset: The synthetic data was constructed using mathematically calculated trends, seasonality, and white noise [45].These time series had variable increasing, decreasing, and mixed global trends.For seasonal component the waveform dataset available in public domain on the Math-Works website was also made part of the time series generation.Additionally, synthetic seasonal pattern mimicking an original time series was generated by averaging the seasonal component and adding white noise to it.The periodicity of seasonal component in these generated time series varied from 5 to 200 timesteps (with intervals of 5).This way the synthetic time series mimicked the real time series of myoelectric indices whose average periodicity was found to vary between [6 (2.375 s) and 18 (6.875s) timesteps] in voluntary contractions and [33 (12.5 s) and 162 (60.875 s) timesteps] in FES-induced contractions.Furthermore, the trends of the original time series were extracted using moving average and were multiplied with trends extracted using other methods such as exponential smoothing, holt-winter's method, linear fit, and polynomial fit.The signal to noise ratio varied from 5 dB (poor) to 20 dB (good).Thus, different trends of increasing, decreasing, and mixed global trends were created.These time series were also compared amongst each other and any two waveforms having coherence more than 0.65 and Spearman's rank correlation more than 0.85 was discarded.This way a dataset of 60,000 synthetic time series was formed.
c) Pretraining MEFFNet on synthetic dataset: Unlike traditional time series forecasting approaches that rely on training a network from scratch on each time series, pretrained MEFFNet was constructed by training MEFFNet on a large dataset of synthetic time series with 1000 epochs (under early stopping mechanism).This pretraining was done to improve the prediction metric and training convergence via transfer learning when applied on the real time series data.Transfer learning is useful when the dataset under consideration is limited, and the likelihood of overfitting increases as the number of epochs increases [46].

E. Forecasting Performance Metrics
The model's performance was evaluated with mean absolute percentage error between original and forecasted signal values according to the equation (1).
Mean absolute percentage error formula is given by: where, A t is the actual value of the time series at timestamp t, F t is the forecasted value of the time series at timestamp t, n is the forecasted window length.

F. Statistical Analysis on the MAPE Values of Various Models
For voluntary dynamic contractions, performance of the models (through MAPE values pooled across features) has been statistically compared using one way repeated measures ANOVA and post-hoc tests (Tuckey Kramer test).Similarly, for FES-induced contractions, the performances of the models were compared using one way repeated measures ANOVA and post-hoc tests (Tuckey Kramer test) for statistically significant differences.

G. Statistical Analysis on the MAPE Values of Pretrained MEFFNet
In the following analyses the subjects were a random factor and were pooled together [47].The statistical tests were performed on the MAPE values generated by pretrained MEFFNet for the time series of myoelectric indices under voluntary dynamic contractions and FES-induced dynamic contractions.The normality of the distribution of MAPE values was checked through the Kolmogorov-Smirnov test.
• Since the seasonal components of myoelectric indices under voluntary contractions had smaller periodicity as compared to their FES-induced counterparts, hence, all MAPE values obtained under voluntary dynamic contractions were compared with those obtained under FES-induced dynamic contractions through the Kruskal Wallis test.

1) Voluntary Dynamic Contractions:
• Since frequency domain indices were more prone to noises and fluctuations, one way ANOVA was performed across the myoelectric indices-ARV, RMS, MNF, and MDF where all subjects and muscles were pooled.Tuckey Kramer test was used for post-hoc analysis to locate the source of difference among the groups.
• In another analysis, MAPE values of time domain indices (ARV and RMS) were pooled together in case of no Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.significant differences [47] and were compared among elbow flexors, elbow extensors, and synergist muscles through the Kruskal Wallis test.A post-hoc analysis by Dunn test was performed to investigate the source of difference.Similarly, the MAPE values of frequency domain indices MNF and MDF were pooled together in absence of significant differences and compared across elbow flexors, elbow extensors, and synergist muscles through the Kruskal Wallis test and Dunn test.
• The forecasting performance of pretrained MEFFNet was studied (through MAPE values pooled across features) for different weight levels (1 kg, 2 kg, 3 kg, 4 kg) using the one way repeated measures ANOVA.
2) FES-Induced Dynamic Contractions: • The MAPE values obtained from a stimulation pattern were compared across myoelectric indices-ARV, RMS, MNF, and MDF through one way ANOVA and post-hoc Tuckey Kramer test.
• MAPE values of myoelectric indices under a stimulation pattern were pooled when found not significant [47] and were compared across customized rectangular, trapezoidal, and muscle synergy-based stimulation templates through one way ANOVA and post-hoc Tuckey Kramer test.
• The forecasting performance of pretrained MEFFNet was studied for healthy and post-stroke group through the Kruskal Wallis test.Note: The statistical level of significance was α = 0.05.All pairwise comparisons were done with the Bonferroni correction, where p-value = α/n, and n= number of comparisons among the groups.

III. RESULTS
The results were obtained with TensorFlow 2.13.0 and python 3.10.12versions respectively.

A. Window Size Effects
1) Voluntary Dynamic Contractions: The seasonal (periodic) component in the time series had an average periodicity of 8.66 samples which means that on an average the participant had taken 3.375 seconds to complete one flexion and extension while carrying weight (feature extracted from sliding window of 0.5 seconds with 25% overlapping on sEMG).Since in each time series the seasonal component and trend varied, hence, an ablation study was required for studying the size effects of input and forecasting window for each of the myoelectric indices.The input window length varied in the range of (60, 120, 240, 480, 600), while the forecasting window length varied along (12,48,60,100,200).The ablation study on window effect size for MEFFNet can be found in Fig. 8 for voluntary contractions.For most models, the input window length of 600 samples and the forecasting horizon of 60 samples produced the optimal results (the results of other models are shown in supplementary material Fig. A).In this study, longer forecasting horizons were desired in models, thus horizons of 60, 100, 200 samples were prioritized more than 12 and 48 samples even though the errors were less.
2) FES-Induced Dynamic Contractions: The seasonal (periodic) component of the myoelectric indices of fatigue in the chosen time series had an average periodicity of 66.65 samples across participants which means that an FES-induced elbow flexion was completed in 25.118 seconds.Therefore, the input window length varied along (100, 140, 280, 480, 600), while the forecasting window varied along (12,48,60,100,200).The ablation study on window effect size for MEFFNet can be found in Fig. 9 for FES-induced contractions.The input window of 600 samples and forecasting horizon of 60 samples yielded the optimal results for MEFFNet and other models (the results of other models are shown in supplementary material Fig. B).

B. Comparison of Different Time Series Forecasting Models
Seven distinct deep learning architectures (including the novel MEFFNet) and a pretrained version of MEFFNet were studied in this work.The pretrained version of MEFFNet outperformed all other models for both voluntary and FES-induced contractions as seen in Fig. 10 and Fig. 11 respectively.For an input window of 600 timesteps and a  forecasting horizon of 60 timesteps, pretrained MEFFNet outperformed the other models with MAPE of 15.99 ± 6.48% and 11.93 ± 4.77% under voluntary (Fig. 10) and FES-induced (Fig. 11) contractions respectively.

C. Forecasting Performance
The forecasted samples (in blue color) given by the pretrained MEFFNet for the first testing window of the ARV time series taken from an exemplar trial of FES-induced dynamic contractions is shown in Fig. 12 (a).The ARV time series was obtained from the eEMG of the biceps brachii muscle of a participant under muscle synergy-based stimulation pattern, as shown in Fig. 12  The mean absolute percentage error in the training and testing samples was found to be 3.21% and 3.60% respectively.The forecasting performance of pretrained MEFFNet for voluntary and FES-induced dynamic contractions is further shown in Fig. 13 under myoelectric manifestation of muscle fatigue.The time domain indicator-RMS was observed as increasing whereas the frequency domain indicator -MNF was found to be decreasing in subfigures (a) and (c) of Fig. 13, under voluntary contractions.For FES-induced contractions (since the stimulation parameters were not modified according to physiological requirements in real-time) these indicators were found to be decreasing with time indicating muscle fatigue, as can be seen in subfigures (b) and (d) of Fig. 13.

D. Statistical Analysis on the MAPE Values of Various Models
For voluntary dynamic contractions, statistically significant difference (p-value < 0.05, F-value: 3581, one way repeated Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.measures ANOVA) was observed among the MAPE values (averaged across features) of different models.Tuckey Kramer test for post-hoc analysis revealed that each model's performance differed significantly (α/28 = 0.00178, p<0.00178,Bonferroni corrected) from others, also shown in Fig. 10 (b).Similarly, for FES-induced dynamic contractions, statistically significant difference (p-value < 0.05, F-value: 3601.62,one way repeated measures ANOVA) was observed among the MAPE values (averaged across features) of different models.Tuckey Kramer test for post-hoc analysis revealed that each model's performance differed significantly (α/28= 0.00178, p<0.00178,Bonferroni corrected) from others, except GAN and MEFFNet which did not show statistically significant difference (p=0.067) in their performance, also shown in Fig. 11 (b).

E. Statistical Analysis on the MAPE Values of Pretrained MEFFNet
• The MAPE values obtained in FES-induced contractions (11.• There were no significant differences observed in pretrained MEFFNet's performance among different weight levels (p=0.0568,one way repeated measures ANOVA).
• Further, the difference between the MAPE values of myoelectric indices (pooled together) were insignificant (p-value > 0.05, one way ANOVA) across the three different stimulation patterns (customized rectangular, trapezoidal, and muscle synergy-based).
• The performance of pretrained MEFFNet across the subjects (healthy and post-stroke) in FES-induced dynamic contractions has been shown in Fig. 14 (a).The performance of pretrained MEFFNet was compared between the healthy and post-stroke group as shown in Fig. 14
IV. DISCUSSION In this study a deep learning network MEFFNet is presented which was pretrained on a large dataset of synthetic time series to forecast 60 samples of myoelectric indices of fatigue produced during voluntary and FES-induced contractions.The 60 samples of forecasting window estimated to 22.62 s (myoelectric indices were extracted from the sEMG using a sliding window of 500ms with 25% overlapping).

A. Window Size Effects
Large input window length generally accommodates many periods of any feature signal which can help in forecasting non-stationary signals.It might also help the model to learn the trend and seasonal repetitions more efficiently.On the other hand, the forecasting window length is dependent both on the input window length as well as the parameters that a model learns [48].The long-term dependencies in a non-stationary signal are difficult to locate, hence, long duration forecasting without increasing the errors is difficult to achieve unless the number of parameters of the model is increased [31], [48].
From the ablation study on window effect size (Fig. 8 and Fig. 9), low MAPE values were found for ratios of 600: 60, 480: 60 and 480: 48 for input to forecasting window lengths.In case of FES-induced contractions each cycle of elbow flexion activity varied anywhere between 12.5 s to 60.875 s (minimum and maximum stimulation duration per cycle as observed in a healthy and post-stroke participant respectively)whereas in voluntary contractions the same activity was achieved in 2.375 s to 6.875 s z(as observed under different weights-1 kg, 2 kg, 3 kg, and 4 kg.Therefore, choosing 600 samples in input improved the overall model's performance as it accommodated the need of the neural network to learn long term dependencies.Since a larger forecasting horizon is mostly a desirable aspect, hence forecasting horizon was chosen as 60 timesteps for further training and analysis.In this study, desirable performance of models was obtained with the minimum ratio of 10:1 for input window length to output window length.The pretrained MEFFNet performed best for 600:60 in both voluntary and FES-induced contractions among all other models.

B. Comparison of Different Time Series Forecasting Models
The performance of all seven models under consideration along with the pretrained version of MEFFNet is shown in Fig. 10 and Fig. 11.From Fig. 10 (b), regarding voluntary dynamic contractions, it could be inferred that each model's performance differed significantly (p<0.00178)from others.However, in Fig. 11 (b), related to FES-induced dynamic contractions, all models differed in their performances significantly (p<0.00178)except for MEFFNet and GAN (p>0.05).
The models which could generate a lower MAPE (less than 25%) in both voluntary and FES-induced contractions were considered as the models of choice.In this criterion, MEFFNet and GAN performed better than the other models with MAPE values (across all myoelectric indices) of 23.04 ± 6.39% and 24.58 ± 7.79% in voluntary contractions and 24.04 ± 5.80% and 24.87 ± 6.66% in FES-induced contractions respectively.Stacked dilated causal CNN couldn't perform well according Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.Hence, stacked dilated causal CNN was not chosen as the good performer unlike MEFFNet and GAN.Other models such as CNN, BiLSTM, LSTM, and Vanilla LSTM produced MAPEs much higher than 30%.To improve the forecasting performance of the chosen models prior training on a large synthetic dataset was done.For this, MEFFNet was chosen over GAN based on the MAPE to prediction-time ratio of the test samples, which was found to be 18.97% per 15ms for MEFFNet and 19.87% per 28ms for GAN.The performance of MEFFNet improved significantly (p-value < 0.05, t-test) in its pretrained version by a drop of 7.05% and 12.12% in MAPE values for voluntary and FES-induced contractions respectively.For a pretrained MEFFNet on synthetic data, it became advantageous to predict future timesteps on real datasets.This could be seen from the MAPE distribution of pretrained MEFFNet and other models for voluntary (Fig. 10) and FES-elicited contractions (Fig. 11).

C. Pretrained MEFFNet
For the pretrained MEFFNet, MAPE values in FES-induced contractions were 4.06% lower than those in voluntary contractions (p-value < 0.05, Kruskal Wallis Test).The high fluctuations in frequency domain myoelectric indices (quantified using autocorrelation as shown in Fig. C of supplementary material) may have resulted in these higher MAPE values [17], [18] in voluntary dynamic contractions.Whereas in FES-induced contractions the frequency domain indices were relatively less affected by fluctuations (quantified using autocorrelation as shown in Fig. D of supplementary material).This attribute could be related to the properties of both compound motor unit action potential (CMUAP) (in case of voluntary contractions) and M-waves (in case of FESinduced contractions).The CMUAPs are more non-stationary than the M-waves of FES-induced contractions [17], [18].Additionally, vEMG signals are occasionally contaminated by physiological, extrinsic, or inadvertent background spikes [49].Therefore, despite having large periodicity in the seasonal component, eEMG myoelectric fatigue indices displayed relatively lower fluctuations (and high autocorrelation, supplementary Fig. D) than vEMG ones which can also be observed in Fig. 4 and Fig. 5.
The pretrained MEFFNet with an MAPE of 15.99 ± 6.48% was better than the other models (Fig. 10) for voluntary contractions across sixteen healthy subjects, ten muscles, and all myoelectric indices.However, this overall performance was affected by the MAPE values of the frequency domain indices (18.17 ± 6.60%) which were significantly higher (p-value < 0.008, Tuckey Kramer post-hoc test) than that of the time domain indices (13.81 ± 5.57%).The higher MAPE values of frequency domain indices could be traced from the high fluctuations present in the signal.This could be due to the muscular vibrations and noises which gets more pronounced in the frequency domain resulting in high fluctuations (and low autocorrelation, supplementary Fig. C) in these time series [50], [51].
Further, the MAPE values of the time domain myoelectric indices (pooled together) were compared among the three muscle groups-flexors (14.78 ± 5.69%), extensors (14.08 ± 5.45%), and synergist muscles (12.05 ± 5.42%) using one-way ANOVA.The MAPE values of the synergist muscles were significantly lower (p-value < 0.016, Bonferroni corrected, Tuckey Kramer post-hoc test) than the flexor muscles and extensor muscles.The consistent seasonal components (as shown through autocorrelation in Fig. E of the supplementary material) in the time series of synergist muscles could have caused lower MAPE than the flexors and extensors.Additionally, the synergist muscles are activated for a shorter period than the prime movers [52], hence, they may have been less affected by the muscle fatigue resulting in more consistent seasonal pattern (high autocorrelation, supplementary Fig. E).
It is interesting to note that for voluntary dynamic contractions, no significant differences were observed in the performance of pretrained MEFFNet among the weightlifting conditions (1 kg, 2 kg, 3 kg, 4 kg).
For FES-induced contractions across all subjects (poststroke and healthy), muscles (biceps brachii and brachioradialis), and myoelectric indices, the MAPE was found to be 11.93 ± 4.77% which was better than the other models (Fig. 11).The performance of pretrained MEFFNet was similar (p-value > 0.5, one way ANOVA) across all myoelectric indices under any stimulation pattern.Since the electrically evoked EMG are quasi-deterministic and quasi-periodic, they are less prone to background physiological noises and random fluctuations unlike the volitional EMG [18].This leads to better spectral estimates (MNF and MDF) in eEMG than vEMG.This might have resulted in the homogenous performance of pretrained MEFFNet across all myoelectric indices obtained from eEMG which are governed by the periodicity of consistent and regular seasonal pattern (measured through autocorrelation in Fig. D supplementary material) and trend in the signal.
Since, evoked EMG and its features depend on stimulation pattern, hence, the performance of pretrained MEFFNet was tested using the MAPE values of myoelectric indices (pooled together) across the three stimulation patterns-customized rectangular, trapezoidal, and muscle synergy-based.No statistically significant differences were observed across the three different stimulation patterns (p-value > 0.05, one-way ANOVA) for time domain and frequency domain myoelectric indices.This suggests that the pretrained MEFFNet performed evenly for the different time series of myoelectric indices evoked under different stimulation patterns which is advantageous for the network.Thus, for FES-induced contractions, the performance of pretrained MEFFNet was non-discriminatory towards stimulation templates and myoelectric indices.
Additionally, for FES-induced contractions, the performance of pretrained MEFFNet was not affected by population type (healthy and post-stroke).All the mentioned observations show that in the present study, the performance of pretrained MEFFNet was affected only by the type of contractions (voluntary vs FES-induced) and by the type of feature (time domain vs frequency domain), as in case of voluntary contractions, frequency domain features could have captured noise prominently and reflected its characteristics in poor forecasting performance of the model.The various weight levels, stimulation patterns, and population type did not affect pretrained MEFFNet's performance as observed in the present study.
Limitations: The study proposes a method for forecasting fatigue indices in voluntary and FES-induced dynamic contractions, however, any claims on generalisation of results can only be made after a proper investigation through experiments on population from all demographics with proper representation from both the genders.A further study in this regard will be envisioned.
An extended discusiion on challenges and way out, applications and future scope can be found in the supplementary material.

V. CONCLUSION
The present study investigated the prospects of implementing pretrained deep learning model MEFFNet for forecasting the myoelectric indices of muscle fatigue in both voluntary as well as FES-induced dynamic contractions targeting healthy and post-stroke subjects at the same time.The pretrained MEFFNet was found to perform forecasting of 60 timesteps which corresponded to 22.62 s with an MAPE value of 15.99 ± 6.48% and 11.93 ± 4.77% in voluntary and FES elicited dynamic contractions respectively.The quantized version of this network could be combined with wearable devices to evaluate the performance in real time.The presented study has potential application in other clinical population such as spinal cord injury, multiple sclerosis, and cerebral palsy patients, where muscle fatigue is a prevalent and detrimental phenomenon.In subsequent work, we will investigate the possibility of extending the scope of the current study to encompass a clinical population that is both large and diverse.

Fig. 1 .
Fig. 1.(a) This figure shows the experiment setup for the biceps curl activity (that is a full elbow flexion-extension) with 1 kg, 2 kg, 3 kg, and 4 kg weights to induce fatigue.sEMG sensors were synchronized with goniometer using DELSYS Trigger Module.(b) This figure shows one cycle of the biceps curl, here, cycle means one complete elbow flexion and extension.Participants were requested to repeatedly perform biceps curls for 7 minutes under each weight (1 kg, 2 kg, 3 kg, 4 kg).A break of 20 minutes was given before proceeding further with another weight.The maximum flexion angle and pace at which the activity was performed were subject-specific.

Fig. 2 .
Fig. 2. (a) This figure shows the experimental setup for elbow flexion in vertical plane.(b) This figure shows the protocol followed in every session of electrically induced elbow flexion.There was a total of 35 cycles recorded.Between the two cycles there was a 10 s off period.Three different stimulation patterns (customized rectangular, trapezoidal, muscle synergy-based) were tested with an inter-session gap of 1 hour 30 minutes on both healthy and post-stroke participants.

Fig. 4 .
Fig. 4. In this figure, (a) shows the volitional EMG, whereas (b) and (c) shows the time domain and the frequency domain myoelectric indices of fatigue, RMS and MNF respectively.The signals were obtained during voluntary dynamic contractions done by a participant bearing 1 kg weight.The indices shown in (b) and (c) as time series, were obtained on EMG shown in (a).

Fig. 5 .
Fig. 5.In this figure, (a) shows the evoked EMG, whereas (b) and (c) shows the time domain and the frequency domain myoelectric indices of fatigue, RMS and MNF respectively.The signals were obtained during FES induced dynamic contractions from a participant under customized rectangular stimulation pattern.The indices shown in (b) and (c) as time series, were obtained on EMG shown in (a).
are high intensity yet short durational when compared with the isometric contractions, without much compromise in accuracy of the forecasted results.Since predicting multiple steps into the future is desirable for wearables, an ablation study was performed with different combinations of input and forecasting window lengths.Fig.6(a) shows the deep learning training procedure schematically whereas Fig. 6 (b) shows the window slicing scheme for the input and the target window length for the deep learning model.Each time series was divided into training to testing-ratio of 70:30 and was fed into the model in the form of successive input and forecasting windows.During training and testing, the model was given real data and not the forecasted data from the previous instance.

Fig. 6 .
Fig. 6.(a) This figure shows the forecasting framework where previous data stored as data buffer is split into input and target windows which are used to train deep learning model.Once trained the deep learning model can start forecasting on real time data incoming to the model from timestamp t m .(b) This figure shows window slicing scheme or say how the time series data has been split into input and forecasting windows successively for training a deep learning model.It also shows the training and testing methodology.For better understanding the ratio of input to forecasting window has been kept as 10:1 which reflects the findings of the present work.A time series has been distributed into training to testing ratio of x: y.The x timesteps or samples were fed into the model, each model input was an input window which consisted of 20 timesteps or samples.Whereas the model output also known as the forecasting window consisted of 2 timesteps or samples.Once the model got trained it was used to forecast future timesteps, however the forecasted timesteps were also compared with the real values in that period which resulted towards the forecasting error.

Fig. 7 .
Fig.7.Network architecture of MEFFNet: the input vectors shaped according to input window length were passed through a 1 D convolution layer for feature extraction followed by a max pooling layer, the resultant was passed through an attention layer for selection of the most prominent features by giving it more attention.The resultant of the attention block was fed into the bidirectional LSTM network of two hidden layers followed by an attention, flatten and fully connected dense layer.The output of the last dense layer corresponds to the forecasting window length.

Fig. 8 .
Fig. 8.This figure shows the mean absolute percentage error (MAPE), averaged across subjects, conditions, and muscles, as a function of the input window and forecasting window for MEFFNet while forecasting myoelectric indices of fatigue (ARV, RMS, MNF, and MDF) obtained during voluntary dynamic contractions.The positive standard deviations are represented in black star ( * ), while the negative standard deviations are shown in pink filled circles ( o ).

Fig. 9 .
Fig. 9.This figure shows the mean absolute percentage error (MAPE), averaged across subjects, conditions, and muscles, as a function of the input window and forecasting window for MEFFNet while forecasting myoelectric indices of fatigue (ARV, RMS, MNF, and MDF) obtained during FES-induced dynamic contractions.The positive standard deviations are represented in black star ( * ), while the negative standard deviations are shown in pink filled circles ( o ).

Fig. 10 .
Fig. 10.(a) The boxplots show the distribution of mean absolute percentage error (MAPE) for each model under the four myoelectric indices of fatigue (ARV, RMS, MNF, and MDF).Each plot contains data from all ten muscles involved in biceps curl (voluntary dynamic contractions) activity of sixteen healthy subjects under all weights: 1 kg, 2 kg, 3 kg, 4 kg.The pretrained MEFFNet outperforms other models under consideration with MAPE of 15.99 ± 6.48%.(b) The statistical comparison across the models reveals significant differences between each pair.The presence of significant difference can be observed by the blue star placed on top of each model which is compared with the model having red square on top.

Fig. 11 .
Fig. 11.(a) The boxplots show the distribution of mean absolute percentage error (MAPE) for each model under the four myoelectric indices of fatigue (ARV, RMS, MNF, and MDF).Each plot contains data from two stimulated muscles (biceps brachii and brachioradialis) under three different stimulation patterns (customized rectangular, trapezoidal, and muscle synergy-based) involved in FES induced elbow flexion activity in seventeen post-stroke and sixteen healthy participants.The pretrained MEFFNet outperforms other models under consideration with MAPE of 11.93 ± 4.77%.(b) The statistical comparison across the models reveals significant differences between each pair.The presence of significant difference can be observed by the blue star placed on top of each model which is compared with the model having red square on top.The absence of blue star over MEFFNet when compared with GAN (red square) represents absence of any significant difference, p = 0.067.
(a) and (b).The corresponding actual samples (in red color) are shown in Fig. 12 (a).The performance of pretrained MEFFNet during training and testing is shown in Fig. 12 (b).In Fig. 12 (c), the percentage error generated in each forecasted sample is shown in light blue and orange colors during training and testing respectively.

Fig. 12 .
Fig. 12.(a) This figure represents the first testing window for ARV time series taken from an exemplar trial of FES-induced dynamic contractions under muscle synergy-based stimulation pattern (biceps brachii muscle), where 600 samples of ARV shown in red color were used to forecast next 60 samples.The forecasted samples shown in blue color are compared with the real values also shown with red color (b) This figure shows the actual values of ARV time series (in red color) and the forecasted values using pretrained MEFFNet (in purple color during training and blue color during testing).(c) This figure shows the percentage error on the forecasted values demonstrated in (b) for both training and testing.The mean absolute percentage value of the errors obtained in (c) was found to be 3.21% for training samples and 3.60% for testing samples.

Fig. 13 .
Fig. 13.This figure shows the time domain indicator of muscle fatigue RMS and frequency domain indicator of muscle fatigue MNF for voluntary (under 1 kg weight) in (a) and (c) and for FES-induced contractions (upon administration of muscle synergy-based stimulation pattern) in (b) and (d) respectively.The frequency domain feature shows declining values with consecutive cycles for both voluntary and FES-induced contractions.The time domain feature had shown increasing trend for voluntary dynamic contractions and decreasing trend for FES-induced contractions.All the myoelectric indices have been closely predicted by pretrained MEFFNet.The actual values of these time series of features (in red color) and the forecasted values using pretrained MEFFNet (in purple color during training and blue color during testing) show the predicting ability of pretrained MEFFNet.

Fig. 14 .
Fig. 14.(a) This figure shows the forecasting performance of pretrained MEFFNet across the healthy (subjects-1 to 16) and post-stroke (subjects-17 to 33) group for all myoelectric indices obtained under three different stimulation patterns for biceps brachii and brachioradialis muscles.(b) This figure represents the statistical comparison between the healthy and post-stroke group tested using the Kruskal Wallis test (p > 0.05, χ 2 = 3.7) for all conditions and muscles pooled.
Statistically significant difference (p-value < 0.05, Fvalue: 108.65, one way ANOVA) was observed among the MAPE values of ARV (13.79 ± 5.38%), RMS (13.83 ± 5.75%), MNF (18.16 ± 6.67 %) and MDF (18.19 ± 6.54 %) for subjects and muscles pooled.Tuckey The MAPE values of time domain indices (ARV and RMS pooled) obtained from pretrained MEFFNet were compared among the flexor, extensor, and synergist muscles through the Kruskal Wallis test.The significant difference (p-value < 0.05, χ 2 :23.67) was observed among these three groups.The post-hoc analysis by Dunn test revealed that the differences (p-value < 0.016, Bonferroni corrected) was due to the group of synergist muscles having MAPE values (12.05 ± 5.42%) which was less than those of extensor (14.08 ± 5.45%) and flexor (14.78 ± 5.69%) muscles.No significant difference was observed between the MAPE values of extensor and flexor muscles.Therefore, the time domain indices obtained from the synergist muscles were less prone to forecasting error than the flexor and extensor muscles.However, no statistically significant difference was observed in the MAPE values of frequency domain myoelectric indices (MNF and MDF pooled) when compared among these muscle groups.