Let x i n , 0 ≤ n < N t be the recorded signal from a target muscle during the i th trial; n is the sampling instant, where n = 0 is the start of a trial, N t is the number of data points from each trial. Let N r and N m be the number of samples in the rest- and move-phases of a trial, respectively, then N t = N r + N m . The time segments 0 ≤ n < N r and N r ≤ n < N t correspond to the rest-phase and move-phase of the trial, respectively.

Problem definition: To detect the presence of EMG in real-time in the move-phase of a trial using only the current and past EMG data x k , 0 ≤ k ≤ n from the start of the trial. Let y n represent the binary output of the sEMG detector at the current sampling instant n , y n = D x k k = 0 n p ∈ 0 1 (1)

where, D ∙ is the detector function that maps the sEMG signal x k k = 0 n , from the start of a trial to a binary output corresponding to the presence or absence of sEMG at the current time instant n ; p is the set of detector parameters that control the detector’s performance. The function D ∙ is often a complex mathematical operation consisting of a series of simpler operations performed on the sEMG data to produce the binary output. This binary output can be used as a simple on/off control of robotic assistance by severely affected patients to relearn movement initiation. ^{34 , 35}

The analysis of the different sEMG detectors was performed using simulated sEMG data. To generate this simulated sEMG, we assume that: •

the measurement noise has a fixed variance throughout the experiment,

These assumptions were made to evaluate the sEMG detectors under the conditions that: (a) the sEMG signal has a fixed signal-to-noise ratio (SNR) in the move-phase, and (b) all other intra- and inter-trial variabilities in the sEMG signal characteristics are minimized. The detectors that perform poorly under these ideal conditions will likely perform worse with real sEMG data from patients, since real sEMG data from severely affected stroke patients might occur in random bursts and is likely to have time-varying amplitude.

We simulated 100 trials of sEMG data with an individual trial duration of 13 seconds (8 s and 5 s for rest- and move-phases, respectively). The sampling frequency of the simulated signal was set to 1000 Hz. Three different sEMG signal generation models were employed in the current analysis – two phenomenological (Gaussian and Laplacian) models and one biophysical model.

The phenomenological models were based on the work of De Luca ³⁶ where a surface sEMG signal from a muscle activated at a fixed level can be treated as zero-mean white noise followed by a shaping filter (electrode properties); this model is widely accepted in the literature. ^{23 , 37 – 40} The exact probability density of the sEMG signal depends on the muscle activation level, with high levels of activation following a Gaussian distribution. ^{37 , 38} However, at low levels of muscle activation, sEMG signals have been reported to follow a distribution that lies between a Gaussian and a Laplacian distribution. ⁴¹ Therefore, to ensure that the detectors are tested with the appropriate signals, we generated data using both white Gaussian and white Laplacian signals, resulting in two phenomenological models.

The first step in this model is the generation of the zero-mean unit-variance white Gaussian and Laplacian noise e n . A step change in the signal variance, within a trial, at the transition between the rest- and move-phases of a trial was obtained by multiplying e n by σ n : e ̂ n = σ n e n , 0 ≤ n ≤ N t (2) σ 2 n = σ 0 2 0 ≤ n < N r σ 0 2 + σ 1 2 N r ≤ n < N t (3) where σ 0 2 and σ 1 2 are the noise and signal variances, respectively. The noise variance is always set to σ 0 2 = 1 in this analysis, and the signal variance is chosen based on the desired signal-to-noise ratio (SNR) in the move-phase. The signal e ̂ n is then zero-phase bandpass filtered (8 ^th order FIR bandpass filter with cut-off frequency 10 Hz and 450 Hz) to have a signal with spectral characteristics like a surface sEMG signal. x n = h sf n ∗ e ̂ n (4) where, h sf n is the impulse response of the bandpass or shaping filter, and x n is the generated sEMG signal that is used for the analysis. The signal-to-noise ratio (SNR) of this simulated sEMG signal in the move-phase is given by, SNR = 10 log σ 1 2 σ 0 2 dB . (5)

In addition to the phenomenological models, we also wanted to test the detectors on more realistic data based on the biophysics of the sEMG signal, accounting for the physiological origin of the electrical muscle activity and the recording electrode geometry. In this paper, the biophysical model proposed in Ref. 42 was employed to generate the simulated sEMG data. Assuming a linear, isotropic volume conduction model, a simple muscle geometry with parallel muscles fibres ignoring the effects due to the finite muscle fibre length, the sEMG recorded by a bipolar electrode configuration can be approximated using the following expression: sEMG t = ∑ q = 1 Q R q t ∗ D q t ∗ e q t ∗ p t D q t = ∑ m = 1 M q δ t − τ m (6) where Q is the number of motor units in the muscle, R q t is the impulse train signal arriving at the q th motor unit through its corresponding motor neuron, M q is the number of muscle fibers in the q th motor unit, e q t is the approximate electrode transfer function between the q th motor unit and the recording electrodes, and p t is the single fibre action potential, which is assumed to be the same for all fibres. The full details of the model can be found in Devasahayam. ⁴² with the associated parameters provided in Table S3 of the supplementary material.

The sEMG simulator developed by Devasahayam ⁴⁰ was employed in the current work to generate the simulated sEMG signals. ⁴² A bipolar surface electrode configuration with a 10 mm interelectrode distance was considered. The simulator takes in the muscle force level as its input and computes the corresponding firing pattern for the motor units. In the current study, the force levels from the muscle were set to 0N in the rest-phase (no muscle activation) and 10N in the move-phase (average firing rate of 16.4 Hz for the muscle). The simulator generated pure muscle activity x pure n recorded by the chosen electrode configuration. The force level for the muscle in the move-phase was chosen empirically to ensure that the temporal profile of the simulated sEMG signal x pure n 0 ≤ n ≤ N t visually resembled that of real surface sEMG signals. A zero-mean white Gaussian noise e n of fixed noise variance σ 0 2 was added to x pure n to introduce measurement noise. The noise variance was chosen based on the signal power ( σ 1 2 ) of x pure n in the move-phase to obtain a signal with the desired SNR: σ 0 2 = σ 1 2 SNR (7)

Following this, the noisy signal x pure n + e n is bandpass filtered (8 ^th order non-causal FIR filter) between 10 Hz and 450 Hz cut-off frequencies: x n = h sf n ∗ x pure n + e n (8) where h sf n is the impulse response of the bandpass shaping filter. The characteristics of sEMG generated from these three models are shown in Figure 1. We wanted to investigate the performance of the detectors under two conditions where the signal power was: (a) equal to noise power, and (b) less than the noise power. Thus, the current study employed two different SNRs of 0dB (signal power equals noise power) and -3dB (signal power is half the noise power).

The general structure for sEMG detectors proposed by Staude et al. ²³ is shown in Figure 2, which consists of three steps carried out sequentially to map the given real-time sEMG data into a binary output: 1.

Signal conditioning is the first step to improve sEMG signal quality for better detection, often involving high-pass filtering for movement artefact removal. Some detectors might employ additional filtering operations, such as adaptive whitening for stable sEMG amplitude estimation. ^{43 , 44} The conditioned signal is represented by x ~ n = S x k k = 0 n , 0 < n < N t (9) where S ∙ represents the mathematical operation performed by the signal conditioning step.

We note that each detector has a set of parameters associated with it. The current study compares the performance of 13 detector types reported in the literature which can be implemented in real-time, listed in Table 1. Each detector type fits into the general structure shown in Figure 2. The different parameters associated with these detector types are also provided in Table 2. A detailed description of the individual detector types and the algorithms for their implementation are provided in the extended data ( Table S1). All detector algorithms were obtained from the literature and implemented in MATLAB ⁵⁰ with appropriate modifications required for real-time detection.

The simulated sEMG data from the three different signal models and the two different SNRs were used to evaluate the performance of the different detector types. Each trial of sEMG signal (13 seconds long) was input to the different detectors to compute binary output indicating the presence or absence of sEMG signal. An optimal EMG detector designed for use in EMG-driven robot-assisted therapy should possess the capability to quickly identify the onset of sEMG, efficiently eliminate false positives, and consistently detect sEMG when it is present (lower false negative). Such a detector might be essential for maintaining the user’s motivation and sense of agency. These are computed from the output y n of each trial ( Figure 3), where the sEMG signal from each trial was analysed in the following three steps: 1.

the first three seconds (0–3s) of the rest-phase data is used for estimating the threshold h for detection: h = μ g + α σ g α ∈ 1 5 (12) where μ g and σ g are the mean and standard deviation of the test function in this period.

We defined a performance measure to compute a single number referred to as the cost of detection that considers the false positive rate r FP , the false negative rate r FN , and the detection latency ∆ t . Let c ≜ r FP r FN f ∆ t T be the cost vector associated with the output y n of the detector for a particular trial. We define the cost of detection C as the infinity norm of the cost vector c . C = max r FP , r FN , f ∆ t (13) providing the worst-case performance of the detector on the given trial.

The false positive rate r FP is defined as the ratio of the number of 1s in the detector output y n in the rest-phase of a trial, and the false negative rate r FN is defined as the ratio of the number of 0s in y n in the move-phase. 0 1 ∋ r FP ≜ 1 N r ∑ n = 0 N r − 1 y n 0 1 ∋ r FN ≜ 1 N t − N r ∑ n = N r N t − N r − 1 1 − y n (14)

From Figure 3, the detection latency is defined as the time delay from the start of the move-phase to when the detector output goes to 1: ∆ t = T s × min n − N r y n = 1 , N r ≤ n < N t (15) where T s is the sampling period of data in milliseconds, and ∆ t ∈ 0 5000 ms . The cost due to this latency is quantified by the function f ∆ t that maps ∆ t to a real number in the closed interval between 0 and 1: 0 1 ∋ f ∆ t = 0 ∆ t < 0 ms ∆ t 250 0 ms ≤ ∆ t < 250 ms 1 ∆ t ≥ 250 ms (16)

Latencies between 0 to 250 ms have linearly increasing costs while the ones above 250 ms are considered as bad as 250 ms. Based on the definitions of r FP , r FN , and f ∆ t , C ∈ 0 1 . A detector with a consistently lower cost of detection C would be considered a better detector.

A detector’s performance or cost is determined by the SNR of the input signal, the detector type, and its associated parameters. Thus, for a fixed SNR input signal, comparing two detector types must be done only after controlling for the influence of their corresponding detector parameters. In the current work, this was done by first choosing the optimal parameters for each detector type, before comparing different detector types. The optimal parameters for a detector type were selected by first splitting the 100 movement trials of simulated sEMG data of the three models (which was generated as explained above) into training and validation datasets with 50 trials each. This was done for both SNRs (0dB and -3dB) and for all three signal models (Gaussian, Laplacian, and biophysical). The training dataset was used to identify the optimal parameter values for the different detector types, i.e. the values of the parameter combination that consistently resulted in the least cost for the detector on the training dataset. The exact procedure is given in Algorithm 1 (end of the document), while details are provided in the extended data.

After identifying the optimal parameter combination for each detector type, the optimal parameter values were used to run the detector on the 50 trials of the validation dataset, which resulted in the validation cost set C val D = C i D i = 1 50 for the detector type D . The cost set from the different detector types were compared using two-way ANOVA with the detector type and signal SNR as the two factors for each of the three-signal model. The complete code for the analysis can be found here.

The 13 detector types were compared after choosing the optimal parameter set for each detector type using the training dataset of 50 trials for each of the six combinations of the three signal models and two SNRs. ⁵² This procedure is depicted in Figure 4 for the Modified Hodges detector for the 0 dB SNR Gaussian signal model, which shows the outcomes from the different steps in the optimization process described in Algorithm 1. Figure 4(a) shows the histograms of the cost C for the different parameter combinations, in light blue traces. These histograms are estimated from the cost values C i i = 1 50 obtained from the 50 trials in the training dataset for the different combinations of the detector parameters. The scatter plot of the median c med and inter-quartile range c iqr of these histograms are shown in Figure 4(b). The choice of the best parameter for the detector was determined to be the one with the least Euclidean norm c med 2 + c iqr 2 which is shown as the red circle in Figure 4(b); its corresponding histogram is shown in the thick red trace in Figure 4(a). Figure 4(c) shows the marginal histograms of the individual contributors r FP r FP f ∆ t to the cost C for the optimal parameter combination for the Modified Hodges detector. The values of the optimal parameters for the different detector types are listed in Table 2.

The performance of the different “optimal” detector types, i.e., detectors using the optimal parameter values, were compared using the 50 trials from the validation datasets. The boxplot of the performance of these different detector types for the three different signal models – Gaussian, Laplacian, and biophysical – are shown in Figure 5(a), (b), (c), respectively; each of these subplots displays the performance for the 0dB and -3dB SNRs in red and blue boxplots, respectively. Note that the order of the depiction of the different detectors is in terms of the increasing average cost across the three signal models and two SNRs; the detectors on the left are better than the ones on the right in an average sense. A two-way ANOVA on the effect of the detector type and SNR on performance revealed a significant difference between the detector types (p < 0.001) and SNRs (p < 0.001) for all three signal models. The test revealed a significant interaction between the factors for all three models (biophysical: p < 0.0001; Gaussian: p < 0.0001; Laplacian: p < 0.0001). These statistical results confirm the results shown in the boxplots in Figure 5, where the performance is different among the detector types, with consistently poorer performance for -3dB compared to 0dB. The costs for both 0dB and -3dB appear to be lower for the biophysical model compared to the Gaussian and Laplacian models.

Most detector types perform similarly except for the Sample entropy, Continuous wavelet transform (CWT), and Singular spectrum analysis (SSA) detector types which perform worse across the different signal models and SNRs. Among the other detector types – the Modified Hodges, the approximate generalise likelihood ratio test- Gaussian (AGLR-G), and the approximate generalise likelihood ratio test- Laplacian (AGLR-L) detectors – have almost similar costs for the different signal models and SNRs. The other detector types – Root mean square (RMS), Hodges, Bonato, Lidierth, Modified Lidierth, Teager Kaiser Energy Operator (TKEO), and Fuzzy entropy – have slightly higher costs for one or more specific signal models and SNRs. We note that the Fuzzy entropy detector performs very well on the biophysical signal model for both SNRs with an acceptable cost for more than 95% of the validation trials.

The choice of appropriate detector type(s) for use in robot-assisted therapy requires the specification of an acceptable cost of detection C accept . To this end, we specify the upper limits for the false positive rate, false negative rate, and the latency of detection as C accept = 0.2 , which corresponds to a detector with the following cost components: r FP ≤ 0.2 r FN ≤ 0.2 ∆ t ≤ 50 ms ⟺ c ∞ = C ≤ C accept = 0.2 (17)

We believe that these upper limits are a reasonable compromise among the three competing factors determining the cost. Any detector type with costs consistently lower than C accept would be deemed an appropriate detector for use in robot-assisted therapy. To determine the detector types with consistently lower costs than C accept , we computed the proportion ( r accept D ) of the 50 validation trials with acceptable C ≤ C accept for the detector D for the three signal models and two SNRs using r accept D = 1 50 ∑ i = 1 50 I C i D ≤ C accept (18) where C i D is the cost of detection for the i th validation trial for detector D . The value of r accept D for the different detector types is shown in Table 3, where the cells with r accept D ≥ 0.8 are highlighted. We can observe there that: 1.

All detectors perform poorly for the -3dB Laplacian signal model. The highest value of r accept D is 0.22 for this signal model, which interestingly is from the AGLR-L detector designed for the Laplacian signal. Many detectors perform a little better with higher r accept D values for the -3dB Gaussian and biophysical signal models.

Based on these observations, the Modified Hodges appears to be the most consistent detector for low SNR signal models, irrespective of the sEMG signal model. The two statistical detectors (AGLR-G, and AGLR-L) and the fuzzy entropy detectors provide similar but slightly lower performance than the Modified Hodges detector.