Surface EMG and muscle fatigue: multi-channel approaches to the study of myoelectric manifestations of muscle fatigue

In a broad view, fatigue is used to indicate a degree of weariness. On a muscular level, fatigue posits the reduced capacity of muscle fibres to produce force, even in the presence of motor neuron excitation via either spinal mechanisms or electric pulses applied externally. Prior to decreased force, when sustaining physically demanding tasks, alterations in the muscle electrical properties take place. These alterations, termed myoelectric manifestation of fatigue, can be assessed non-invasively with a pair of surface electrodes positioned appropriately on the target muscle; traditional approach. A relatively more recent approach consists of the use of multiple electrodes. This multi-channel approach provides access to a set of physiologically relevant variables on the global muscle level or on the level of single motor units, opening new fronts for the study of muscle fatigue; it allows for: (i) a more precise quantification of the propagation velocity, a physiological variable of marked interest to the study of fatigue; (ii) the assessment of regional, myoelectric manifestations of fatigue; (iii) the analysis of single motor units, with the possibility to obtain information about motor unit control and fibre membrane changes. This review provides a methodological account on the multi-channel approach for the study of myoelectric manifestation of fatigue and on the experimental conditions to which it applies, as well as examples of their current applications.


Introduction
Fatigue, in its broad sense, occurs in everyday life and can be described as a feeling of weakness, of muscle pain or a decrease in performance during physical or cognitive activities. When analyzed at the level of the neuromuscular system, fatigue is classically defined as an exercise-induced decline of the maximal force muscles can generate (Gandevia 2001, Enoka andDuchateau 2008).
During sustained submaximal contractions, both peripheral and central changes lead to alterations in the activity of muscles undergoing fatigue prior to task failure (Gandevia 2001, Barry andEnoka 2007, Enoka andDuchateau 2008). These changes, termed myoelectric manifestation of fatigue, can be assessed non-invasively from EMGs detected with a pair of surface electrodes positioned appropriately on the target muscle.
Many methods have been developed to characterize muscle fatigue from surface EMGs during isometric and dynamic contractions. Among them are the Fourier-based spectral estimators, parametric based spectral estimators, joint analysis of EMG spectrum and amplitude (JASA), time-frequency and time-scale methods, Wigner-Ville distribution, Choi-Williams distribution, time-varying autoregressive approaches, wavelet analysis, entropy, recurrence quantification analysis. Some recent works in literature give an exhaustive overview and critical analysis of the mentioned methods for fatigue evaluation with the bipolar electrodes (Cifrek et al 2009, Gonzàlez-Izal et al 2012, Merletti and Farina 2016; for this reason they will not be described in this review. This review specifically focuses on providing readers with a comprehensive view of the insights gained into the study of myoelectric manifestation of fatigue with multiple, surface electrodes; i.e. with the multi-channel approach.

Myoelectric manifestations of fatigue assessed with the multi-channel approach
The possibility of using multiple surface electrodes to sample EMGs from multiple regions of any given muscle has been attracting the interest of a progressively greater number of researchers and clinicians. Surface EMGs may be detected serially along the longitudinal or transverse axis of the muscle using linear arrays of surface electrodes (Wood et al 2001, Merletti et al 2003. Alternatively, the distribution of EMG amplitude and spectral descriptors may be assessed over the whole skin region covering the target muscle with bi-dimensional arrays of electrodes (Falla and Farina 2007, Gallina et al 2013b, Vieira et al 2015. Hereafter, multi-channel approach is used to indicate indistinctively the sampling of EMGs with multiple electrodes, regardless of whether arranged serially or bi-dimensionally. The simplest multi-channel approach consists in the manual positioning of individual electrodes over the muscle of interest (Sacco et al 2009, Staudenmann et al 2009. Even though it allows for the distribution of electrodes in different ways over any target muscle, the manual positioning of multiple surface electrodes may demand considerable time and is unlikely to provide a fixed inter-electrode distance. Moreover, the number of detection points is limited to the size of individual electrodes, which is typically relatively large; sampling from multiple regions of e.g. the small hand muscles with conventional pre-gelled electrodes is likely unviable. Embedding electrodes into flexible dielectric materials, on the other hand, has led to the development of an assortment of multi-channel detection systems, of varying size, dimension of electrodes and inter-electrode distances that adapt to almost any superficial muscles, from the external anal sphincter (Cescon et al 2011) to the facial muscles (Lapatki et al 2006). The technology considered for the design of such systems and the insights they provide into the neuromuscular system were the subject of reviews and text books published recently (Blok et al 2002, Drost et al 2006, Holobar and Farina 2014, Merletti and Farina 2016. In this section, we specifically indicate to readers how the multi-channel approach may contribute to advancing our knowledge on the neurophysiological mechanisms underpinning muscle fatigue, both at the global muscle level and at the level of individual motor units (MUs).

Global muscle level
In virtue of the relatively large inter-electrode distance, surface electrodes sample from a considerably larger muscle volume than intramuscular electrodes (Lynn et al 1978). With respect to EMGs collected from within the muscle, surface EMGs consequently provide a more representative view of the whole muscle electrical activity and are therefore more likely sensitive to myoelectric changes resulting from muscle fatigue (Rudroff et al 2008). The multi-channel approach further broadens the representation of muscle activity in the EMGs, opening new fronts for the study of muscle fatigue; it allows for: (i) the quantification of the propagation velocity, a physiological variable of marked interest to the study of fatigue (Merletti et al 1990); (ii) the assessment of regional, myoelectric manifestations of fatigue (Gallina et al 2011, Watanabe et al 2013. Both arguments are covered below.

Conduction velocity (CV).
Muscle fibre CV is defined as the propagation velocity of depolarization along the muscle fibre (Merletti et al 1990). Differently from EMG amplitude and spectral variables, CV is not only a mathematical descriptor but an important physiological variable directly related to fibre membrane properties, fibre diameter and fibre contractile characteristics (Andreassen and Arendt-Nielsen 1987). It is well established that changes in CV affect the shape of the measured motor unit action potential (MUAP) waveform and, finally, the amplitude and spectral variables of surface EMG (Lindstrom and Magnusson 1977, De Luca 1984, Arendt-Nielsen and Mills 1985. During sustained contractions, metabolic changes in the muscle affect the cellular mechanisms underlying ionic exchanges across the membrane and therefore the propagation of action potential (Allen et al 2008). These changes result in a progressive reduction of CV. This reduction is one of the main causes of the changes occurring in amplitude and spectral EMG variables during fatigue (Merletti et al 1990). Accessing this basic, physiological variable provides a relevant means to describe and quantify the fatigue processes taking place within the muscle. This section describes the methods for CV estimation, highlighting the important factors for its interpretation as an index of muscle fatigue.
2.1.1.1. Methods for CV estimation and implications for assessing fatigue. The estimation of muscle fibre CV from surface recordings is based on the measure of the delay between potentials recorded at fixed distances along the muscle fibres. In ideal conditions (i.e. when signals are identical, continuous and noise free) this delay can be determined with very simple methods. However, as described later, these conditions are hardly met in real recordings and the problem of delay estimation becomes more complex. Indeed, if signals are similar but not identical, the delay can be mathematically defined in several ways, depending on the criteria used to compare signal shapes (Merletti and Farina 2016).
There is an extensive literature concerning methods to measure the delay between two signals. The simplest approach is to identify latency between reference points (e.g. peak value or zero crossing) of the two signals (Andreassen and Arendt-Nielsen 1987). This method is highly affected by subtle variations in the shape of the potentials due to electrode-fibre misalignment or to noise that can markedly modify the time of occurrence of the reference points. This approach seems appropriate only in recordings with high signal-to-noise ratio, such as in electrically induced contractions and in experimental conditions where the high similarity between recorded potentials can be ensured. Another simple delay estimator, though more robust to noise and slight changes in signal shape, is the time shift that minimizes the least square error between the two detected signals. This delay definition corresponds to the maximum likelihood estimation in case of white additive noise . According to this delay definition, estimators based on the maximum of cross-correlation between pairs of signals have been proposed (Sollie et al 1985).
A common issue associated with the time domain methods applied to discrete signals is the temporal sampling. In such case, regardless of the method used, the delay estimation in time can be performed at best with a resolution of ±0.5 sampling interval (hundreds of µs for the sampling rates classically used for surface EMG, 1-2 kHz). Considering that the delay resulting from propagation is in the range of a few ms (for electrodes 1-2 cm apart), this limited temporal resolution may lead to relevant errors in the delay estimation , even if the signals are identical and sampled above the Nyquist frequency (see figure 1 in McGill and Dorfman (1984)). This problem can be addressed by increasing significantly (10-20 times) the sampling frequency, with the obvious consequences in terms of memory occupation, or by interpolating the cross correlation function around the peak (Wheeler and Smith 1988).
Alternative approaches consist in matching the signals in the frequency domain (i.e. to the Fourier transforms of the signals). In the frequency domain the delay is a continuous rather than a discrete variable . On this regard, spectral matching (McGill and Dorfman 1984) and phase delay (Hunter et al 1987) techniques have been proposed. The comparison between temporal and frequency methods showed that methods based on interpolated cross-correlation function have similar performance of spectral matching and that the phase delay method should be avoided in case of low signal-to-noise ratio .
The aforementioned methods have been originally developed to be applied to pairs of single differential or double differential signals. More recently, methods based on more than two channels arranged serially (linear arrays) or bi-dimensionally (grids of electrodes) have been proposed (Farina et al 2001, Farina and Merletti 2004, Farina et al 2004c, Grönlund et al 2005, Mesin 2015, Soares et al 2015. These methods have been shown to significantly improve CV estimates, in terms of variance of estimation, sensitivity of the estimate to electrode location and repeatability. Approaches based on image processing techniques applied to spatio-temporal pattern of signal propagation detected with grids of electrodes have been tested on simulated signals, providing promising results (Grönlund et al 2005, Mesin 2015. These methods may also be used to identify anatomical features of the muscle (and their change during the contraction) such as innervation location and muscle fibre orientation, critically affecting the estimation of EMG variables (Rainoldi et al 2000, Farina et al 2002b as described below.
Fatigue-related changes in CV can be assessed by applying the methods described above to the global EMGs (M-wave in case of electrically induced contractions or the interferential EMGs in case of voluntary contractions) or to single motor unit action potentials (MUAPs). When applied to the global signal, the estimated CV is an average value originating from the actual CV of all MUAPs detected in the considered time window; this average is weighted for the size of the detected MUAPs and is therefore biased towards the CV of MUs better represented in the surface potential, i.e. the largest and the most superficial ones. Indeed, considering that large MUs usually have higher CV (Andreassen and Arendt-Nielsen 1987), if MU territories are randomly scattered across the detection volume, it is more likely that larger MUs, with higher CV, generate larger potentials, thus biasing CV estimates and their changes with fatigue. Additional factors possibly influencing CV estimates are the presence of generation and extinction potentials, the inclination of fibres with respect to the detection system and the non-homogeneity of subcutaneous layers along the fibre length (Schneider et al 1991, Farina et al 2002b. Among these factors, misalignment between the fibres and the detection system seems to be particularly critical. It has been shown that when anisotropic detection systems (single differential or double differential) are used, a fibre-electrode misalignment frequently results in CV overestimation, depending on the position of the MU territory with respect to the electrodes and on the characteristics of the detection volume (Schneider et al 1991, Merletti et al 1999, Farina et al 2002b. Moreover, beside the cases where fibre-electrode alignment is overtly not viable, as for in-depth-pinnate muscles (Mesin et al 2011, Gallina et al 2013b (see also section 2.1.3), it is often difficult and time-consuming to identify the optimal alignment, which could change during the contraction. In addition, even in isometric and constant force contractions, relevant changes in fibre geometry may occur due to mechanical adaptations of muscle structures during fatigue (Kubo et al 2001, Shi et al 2007, Mitsukawa et al 2009. This is a possible additional confounding factor that, together with recruitment/de-recruitment of MU populations (Gazzoni et al 2001, Houtman et al 2003, Bawa and Murnaghan 2009, complicates the interpretation of CV changes even in well controlled conditions, such as isometric and constant force contractions. Interpretations may be further aggravated when CV is normalized with respect to the initial value. As described above, the CV estimation obtained from interferential EMGs is the mean value of the CV distribution of active MUs, weighted by the effect of the conduction volume. The possibility of accessing the profile of CV distribution with non-invasive techniques is certainly intriguing, as it provides the means of describing in more details the complex events occurring during fatigue likely hidden in the mean value. In 1987, Davies and Parker (1987) proposed the first method to extract features of the CV probability density function from surface EMGs. The method was based on the measurement of the cross-and auto-power spectra of surface EMGs detected by a pair of electrodes. More recently, different methods have been used to investigate the changes in CV distribution during fatigue , Lange et al 2002, Gazzoni et al 2004. All these approaches are based on the identification of MUAP occurrences in the interference signal and CV is then estimated for each identified MUAP. Methods substantially differ on the way MUAPs are extracted. The simplest approach was proposed by Lange and collaborators (Lange et al 2002); it is based on the identification of the occurrence of pairs of peaks in two bipolar EMGs. Epochs containing peaks separated by a delay compatible with a physiological CV are selected for CV estimation. This approach has been used to investigate the recruitment of MU populations during fatiguing contractions (Houtman et al 2003), fatigue in fibromyalgia (Klaver-Król et al 2012) and differences in CV distributions between sprinters and endurance athletes (Klaver-Król et al 2010). A more sophisticated method proposed by  applies a continuous wavelet transform to the double-differential surface EMG to determine the time location and to estimate the duration of waveforms that are candidate MUAPs. An algorithm based on beamforming is then applied to estimate the CV and to decide whether the detected waveform is a travelling potential (and therefore contributing to CV distribution) or not. An example of fatigue-induced changes in CV distribution in healthy biceps brachii estimated with this method is provided in Figure 1. (a) Global CV estimates and (b)-(d) CV distributions obtained from EMGs detected during a 30 s long constant force (60% MVC) isometric contraction of a healthy biceps-brachii. Global CV was estimated using the method described in McGill and Dorfman (1984) whereas CV distributions were calculated from 1s epochs using the method described in . It can be observed that global estimate and its changes are in agreement with those provided by the estimation of CV distribution. Adapted from  with permission. figure 1; the time course of mean muscle fibre CV estimated with spectral matching (McGill and Dorfman 1984) is also reported for comparison. This method was experimentally validated by comparing CV estimates with those obtained with spike trigger averaging technique (Farina et al 2002a) and was applied to study fatigue in elderly (Merletti et al 2002) and in subjects affected by carpal tunnel syndrome (Rainoldi et al 2008). , Gazzoni et al (2004 proposed a method for the extraction and classification of single MUAPs from surface EMG: MUAPs are first identified and extracted from the interference signal (segmentation step) and then clustered in MUs using a multichannel neural network. By estimating CV for all MUAPs identified in the segmentation step, Casale et al (2009) investigated the fatiguerelated changes in CV distribution in subject affected by fibromyalgia.
Methods for CV distribution seem to provide a powerful tool to obtain more insights into the changes of mean muscle CV. It is important to underline, however, that in all the methods described, CV is estimated from all the identified firings, thus providing the CV distribution of the identified MUAPs, not of the identified MUs (which would imply classification of identified MUAPs and a single value of CV estimation for each, identified MU). This means that MUs firing with higher rates contribute to the CV distribution with higher number of occurrences, thus biasing the CV distribution and its descriptors (e.g. standard deviation). For the same reason, MUs that are represented with higher energy on the surface potential tend to mask the contribution of low energy ones, thus biasing the CV distribution profile (Merletti et al 2002). Methods based on MU decomposition (Gazzoni et al 2004, De Luca et al 2006, Holobar and Zazula 2007, Kleine et al 2007, Holobar et al 2009, Negro et al 2016 have the potential of providing the actual MU CV distribution. However, at the present time, the number of MUs extracted by these methods (few tens) may not provide CV histograms sufficiently representative of the MUs within the electrodes' detection volume.

Changes in CV and spectral descriptors, do they provide the same information?
The effect of CV reduction on spectral and amplitude variables of surface EMG is well established (Lindstrom andMagnusson 1977, De Luca 1984). Lindström and Magnusson (Lindstrom and Magnusson 1977) demonstrated analytically how muscle fibre CV acts as a scaling factor for EMGs amplitude and its power spectrum. They showed that if a signal is scaled in time by a factor k, due to a change in CV (k < 1 for CV reduction), its power spectrum is compressed by the same factor. This compression in the frequency scale results in a change of spectral descriptors, such as mean frequency (MNF) and median frequency (MDF). The scaling affect amplitude variables as well; average rectified value (ARV) changes by a factor 1/k and root mean square (RMS) changes by a factor 1/sqrt(k). From these theoretical considerations, it follows that if the CV reductions were the only changes occurring during fatiguing contractions, scaling factors of CV and spectral parameters would be identical. In this case, the CV rate of decrease could be predicted by measuring the scaling factor (k) determined by the decrease of mean (or median) frequency of the power spectral density and a single channel EMG detection would provide a substantial description of fatigue processes taking place at the muscle level. Unfortunately, as explained in the following, this does not seem to be the case for most of the real scenarios. Indeed, as observed by several Authors, the rates of decrease of spectral variables are generally greater than those of CV, suggesting that CV reduction is not the only factor determining spectral changes in surface EMG (Bigland-Ritchie et al 1981, Merletti et al 1990, Hägg 1992, Krogh-Lund and Jørgensen 1993, Rainoldi et al 1999, Dimitrova and Dimitrov 2003, Farina et al 2006. The observed differences between the behavior of CV and spectral descriptors are generally more subtle, though often present, during electrically induced contractions and more pronounced during voluntary efforts (figure 2). A number of possible causes have been suggested to explain these experimental observations. The first argument in this regard is that, besides CV, local fatigue induces changes in other parameters affecting the shape of the potential, such as duration of intramuscular action potential (IAP,) amplitudes of IAP and negative and positive after-potential (Dimitrova and Dimitrov 2003). Changes in the statistical distribution of CV values of active MU have been also regarded as a factor influencing the differences in the time course of CV and MNF estimates during fatigue. Specifically, Merletti et al (1990) suggested that different profiles of CV decrease for different fibres within the same MU may lead to a widening of surface potential whose effect on MNF is more pronounced than that due to the sole CV reduction. This behavior, described for a single MU, becomes particularly relevant in electrically induced contractions as the recorded EMG (the M-wave) results from the synchronous contribution of several MUs with a wider range of CV and CV rates of change. When voluntary contractions are considered, MU recruitment and firing properties play an additional, relevant role. It has been demonstrated that the MUs recruited for compensating the fatigue-related reduction in muscle contractility may lead to an increase or decrease of frequency content of the interferential surface EMGs, depending on their position within the volume conductor, thus resulting in an unpredictable effect on spectral variables of surface EMGs (Farina et al 2006). Since CV estimations are theoretically insensitive to MU location within the muscle, MUs' recruitment may not affect equally CV estimates and spectral variables. The second consideration concerns the increased degree of MU synchronization associated with fatigue. Kleine et al (2001) showed that firing synchronization enhances the low frequency components of power spectrum thus reducing its median frequency. Modifications of spectral properties, also associated with an increase of surface EMG average rectified amplitude, do not seem however to affect significantly CV estimations (Farina et al 2002c). These considerations suggest a possible role of MU synchronization in the observed differences between relative changes of CV and spectral variables in voluntary, fatiguing contractions.
Overall, these observations indicate that fatigue-induced changes in spectral variables of global EMG cannot be considered a surrogate of CV. The effect of volume conductor properties and of fatigue-induced changes in the MU recruitment patterns (i.e. activation of additional MUs, de-recruitment of fatigued MUs and synchronization of active MUs) represents a relevant confounding factor for the interpretation of global variables, especially the spectral ones (Farina et al 2004a). This aspect is analyzed in more detail in section 2.2.
2.1.2. Regional variations of muscle activity. When collecting EMGs from a single muscle region, with the aim of e.g. studying muscle fatigue, one expects the recorded EMGs to be sensitive to neural changes within the whole muscle volume. Key results obtained from simulated and experimental EMGs detected with the multi-channel approach are though in contrast with the assumption of single, surface EMGs being representative of the whole muscle volume (Mesin et al 2011, Gallina et al 2013a, even for muscles as small as those in the forearm (Gallina and Botter 2013). Indeed, the amplitude of surface EMGs detected from a single muscle with electrodes positioned at different locations was observed to change by different extents in different circumstances, among which during fatiguing contractions (Zijdewind et al 1995, Falla and Farina 2007, Holtermann et al 2010, Gallina et al 2011. These pieces of evidence suggest the detection volume of surface electrodes, that is the muscle region within the pick-up volume of electrodes (Lynn et al 1978), is not as large as previously appreciated. In opposition to the detection of activity from muscles other than that of interest (i.e. crosstalk), surface EMGs detected locally may not reveal the changes occurring somewhere else within the muscle (Mesin et al 2011. This possibility, previously considered by others in applied studies (see page 20 in Nashner 1977 and page 489 in Joseph and Nightingale 1952), leads us into the question: do the myoelectric alterations signalling fatigue manifest equally on different skin regions covering the target muscle? Here we explore how the multi-channel approach has been used to address this question and highlight its added value to the study of muscle fatigue.
Before all, we would like to clarify what we mean by regional variations of muscle activity. The amplitude and the spectrum of surface EMGs collected from different muscle regions change in different ways. These changes have different meanings, depending on the direction along which they are assessed with respect to the muscle fibres. Considering for example a series of surface electrodes aligned parallel to the muscle fibres, the amplitude and bandwidth of EMGs differ between electrodes. Electrodes located nearby the innervation zone (i.e. the muscle end-plates) and the tendon regions provide single-differential EMGs with remarkably low amplitude and high frequency content (Rainoldi et al 2004). If detected in monopolar derivation, EMGs collected nearby the innervation zone and tendon regions show, conversely, high amplitude and high energy in the low frequencies (Kleine et al 2000). The amplitude and spectral descriptors of EMGs detected both in monopolar and single-differential derivations do not change across electrodes located between the innervation zone and tendon regions; the same action potential is detected by these electrodes. Explaining the reason for such dependence of EMG descriptors on anatomical factors and on the detection modality is beyond the scope of this review. The interested reader may refer to previous literature for a comprehensive view on this argument (Farina et al 2002a). Most importantly, our point here is that changes in surface EMGs collected at different regions along the longitudinal axis of muscle fibres should be not conceived as reflecting neural changes across the muscle. Even though aligning a series of consecutive electrodes parallel to the muscle fibres is a necessary condition for the estimation of CV, the EMGs detected by these electrodes may be not used to support inferences on regional variations in muscle activity; that is, changes in activity between different regions within individual muscles. Inferences on regional variations in muscle activity proceed only from surface EMGs collected by electrodes positioned at skin regions covering different muscle fibres and not different sections of the same group of fibres.
The distinction between regional variations in muscle activity and anatomically-induced changes in EMGs is clear when EMGs are collected with grids of electrodes. Figure 3 shows examples of single-differential action potentials of two MUs detected from the vastus medialis muscle of a single subject (Gallina and Vieira 2015). MUs were identified by decomposing monopolar EMGs detected by 128 (16 × 8) electrodes, with columns of electrodes aligned parallel to the muscle longitudinal axis (figure 3a). Given the vastus medialis fibres are oriented obliquely with respect to the muscle longitudinal axis (Smith et al 2009), columns of electrodes and fibres are not located in parallel directions. Consequently, action potentials of individual MUs are clearly represented in oblique directions (45°) within the grid (figure 3). Close inspection of traces shown in figure 3(b) reveals there is a phase opposition between action potentials detected from the first to the fourth column and those detected from the two right most columns of electrodes. Similarly, a phase opposition is observed between potentials detected at the first column and those detected from the second to the sixth column for the MU shown in figure 3(c). Moreover, consecutive potentials detected at either side from the region where phase opposition is observed clearly show a consistent delay between themselves. Such a consistent delay indicates the propagation and therefore the fibres' direction is oblique with respect to rows and columns of electrodes. The oblique representation of action potentials of these two MUs is better appreciated in RMS images created separately for each surface potential shown in figures 3(b) and (c). For both MUs, potentials with high RMS amplitude appear along the proximo-medial disto-lateral direction, coinciding with the oblique orientation of vastus medialis fibres (Gallina and Vieira 2015). Narrow regions with markedly low RMS amplitude interposed between regions with large RMS amplitude are clearly seen (i.e. the innervation zone location). Also, potentials with small RMS amplitude are seen at proximomedial extremity for MU 8 (figure 3(b)) and at the disto-lateral extremity for MU 10 (figure 3(c)), likely suggesting electrodes were located close to the fibres' ending. The changes in RMS amplitude observed in the proximo-medial disto-lateral direction result from the MU anatomy; the location of innervation zone and tendon regions in relation to electrodes. These are anatomically-induced changes in the surface EMGs. On the other hand, RMS images clearly show the action potentials of these two MUs are represented at somewhat distinct proximo-lateral disto-medial regions. Fibres of MU 8 are located ~6 cm proximo-laterally to the fibres of MU 10. Variations in surface EMGs collected at the proximo-lateral disto-medial direction, that is transverse to the muscle fibres' direction, are therefore indicative of regional variations in muscle activity.
Regional variations in muscle activity, studied with the multi-channel approach, provided two main insights into the myoelectric manifestations of muscle fatigue. The first insight concerns the manifestation of fatigue in the surface EMGs itself; namely the location within individual muscles where changes in surface EMGs are most expressive and thus where fibres are expected to most likely get fatigued. During a positional, endurance task, for instance, changes in the amplitude and spectrum of surface EMGs were observed to manifest locally within the trapezius muscle (Farina et al 2008c). Specifically, 11 subjects were asked to sustain their shoulder joint abducted at 90 degrees until failure while surface EMGs were collected from the trapezius muscle with a grid of 64 electrodes (13 × 5 arrangement, with columns aligned along the cranio-caudal direction; one missing electrode at the cranio-lateral edge of the grid). As shown in figure 4, at contraction start, EMGs with largest RMS amplitude were detected centrally. With the progression of contraction, the region where largest EMGs were detected changed systematically towards electrodes located cranially (see crossed circles in figure 4). At the contraction end (100% endurance time), EMGs with greatest amplitude were observed at the cranial third of the grid. Using a different EMG descriptor, an independent, contemporary study to that of Farina et al (2008a) further substantiated the occurrence of regional variations in muscle activity during fatiguing contractions. With grids of 130 electrodes (13 × 10 electrodes), Holtermann et al (2008) observed a local decrease in CV along the trapezius cranio-caudal direction while subjects sustained isometric contractions at different effort levels. Interestingly, changes in CV did not manifest consistently across subjects (figure 5). While CV decreased along the whole grid for some subjects (subject 1 in figure 5), for others it decreased mainly cranially or caudally (see subjects 2 and 3 in figure 5, respectively). These authors also reported local variations in CV for the biceps brachii muscle; the probability of identifying decreased CV during sustained, elbow flexion contraction differed along the muscle medio-lateral direction. Collectively, these results and results from other studies, assessing different muscles (Zijdewind et al 1995, Gallina et al 2011, Watanabe et al 2013, suggest fatigue is more likely to affect muscles locally rather than globally. The second insight provided by the multi-channel approach concerns the mechanisms causing or anticipating muscle fatigue. Variations in surface EMG descriptors shown in figures 4 and 5 refer to changes occurring in directions transverse to muscle fibres; i.e. they denote regional variations in muscle activity during fatigue. These regional variations in muscle activity have been interpreted as regional manifestations of muscle fatigue (Gallina et al 2011, Watanabe et al 2013. In the course of a fatiguing contraction, prior to a decline in muscle force (Holtermann et al 2008) or inability to sustain a given joint position (Farina et al 2008a), EMGs collected in different muscle regions change differently (figures 4 and 5). An interesting observation is the potential, functional relevance of such regional myoelectric manifestations of fatigue. Changes in the spatial distribution of EMG descriptors were observed to correlate positively with the duration along which subjects were able to sustain a constant force level (Falla andFarina 2007, Holtermann et al 2010). Interestingly, subjects exhibiting more expressive regional variations in muscle activity seem to have a greater potential for sustaining a given workload (Farina et al 2008a). Similar findings were observed for EMGs collected intramuscularly from the gastrocnemius muscle (McLean and Goudy 2004). The physiological mechanisms explaining the association between better endurance and greater variations in muscle activity are so far subject to speculation. One possibility is the optim isation of extracellular environment (i.e. expediting the removal of fatigue by-products) through the local modulation of blood flow. Kouzaki et al (2003), for example, reported a significant association between alternate activation of knee extensors and changes in blood volume during sustained, isometric contractions. As shown in figure 6, the amplitude of integrated EMGs collected from rectus femoris changes in an opposite fashion when compared to the amplitude of EMGs detected from the vastii. Most importantly, decreases/increases in rectus femoris activity are associated with increases/decreases in the muscle blood volume (see traces below integrated EMGs shown in figure 6). In an attempt to address the relevance of such alternate, synergistic activity, Kouzaki and Shinohara (2006) further assessed the association between alternate muscle activity and force production in fatiguing contractions; they indeed observed smaller decreases in knee extension torque for subjects showing more frequent periods of alternate muscle activity. Even though the association between modulation of blood flow and alternate activity within individual muscles seems yet untested in the literature, it is possible that regional modulation of blood flow contributes to explaining why greater variation in muscle activity correlates with better endurance (Falla andFarina 2007, Holtermann et al 2010).
Motor unit rotation or substitution is another possibility accounting for the notion that greater spatial variability in surface EMGs is associated with delayed fatigue. Both terminologies have different meanings (Westgaard and de Luca 1999). Substitution corresponds to the recruitment of MUs to replace active, tonically discharging MUs of slightly lower recruitment thresholds. Rotation, on the other hand, indicates different MUs are periodically de-and rerecruited in an alternate fashion throughout a sustained contraction. In spite of semantic differences, both substitution and rotation represent a mechanism for prolonging the ability of individual muscles to sustain a given force; the burden of fatigued MUs is taken up by fresh Figure 6. Typical example of knee torque and of EMGs integrated over 1 s epochs, after full-wave rectification, for rectus femoris and for vastus lateralis and medialis muscles, respectively from top to bottom. Traces denoting the blood volume profile for recuts femoris and vastus lateralis muscles are shown just below to the integrated EMGs. Data was collected during sustained, knee extension at 2.5% of maximal voluntary contraction (MVC) for 60 min. Scales for the surface EMG are expressed in relation to the maximal amplitude value obtained for integrated EMGs collected at maximal, isometric contractions. Modified from Kouzaki et al (2003) with permission.
MUs. Substitution and rotation have been observed for an assortment of muscles, small and large, spanning distal and proximal joints and in different situations, constant-force and-force varying contractions (Fallentin et al 1993, McLean and Goudy 2004, Bawa and Murnaghan 2009, Manning et al 2010. An example of MU rotation is shown in figure 7 for the extensor carpi radialis muscle. The first MU fires regularly at ~10 pps from the beginning of the figure. At about 10 min later, MU 2 starts to discharge tonically while unit 1 ceases firing. A few minutes later unit 1 resumes firing and the MU 2 falls silent and then starts firing again after ~8 min. The controversy with the concept of MU substitution and rotation is the apparent violation of the size principle (Westgaard and de Luca 1999). If MUs are recruited according to their size, which is strongly associated with their recruitment threshold, then, the recruitment of a higher threshold unit should be not accompanied by the derecruitment of an already active, tonically discharging unit. Alterations to the intrinsic properties of motor neurons (e.g. slowing inactivation of ionic channels, slowing of the Na + -K + pump and persistent inward currents) have been therefore suggested to explain the rotation/substitution of MUs (Bawa and Murnaghan 2009), as these changes would affect the motor neuron recruitment threshold; note e.g. that MUs 1 and 2 in figure 7 are recruited and derecrutied at different force levels throughout the contraction. Consistent with this view, Manning et al (2010) observed the excitability of a motor neuron increased progressively after it ceased firing and after substitution took place during a sustained contraction.
An alternative, or more likely complementary, view to rotation triggered by intrinsic properties is the spatially differentiated distribution of net excitatory input to motor neurons. EMGs collected intramuscularly from the small hand muscles (Zijdewind et al 1995) and from the large calf muscles (McLean and Goudy 2004) suggest indeed that motor neurons serving fibres in different muscle regions may participate in rotation during fatiguing contractions. Even though not at the level of single MUs, changes in the distribution of the amplitude of surface EMGs detected with the multi-channel approach with fatigue (Farina et al 2008c, Holtermann et al 2008 are consistent with the idea that rotation may extend to different regions within muscles; i.e. the regional rotation of MUs. This notion of regional rotation of MUs is not necessarily in contradiction with the size principle. Supposing all motor neurons serving an individual muscle receive the same synaptic input, the intrinsic changes in motor neuron properties triggering regional rotation would affect to different extents motor neurons serving different muscle regions, otherwise no regional variations in muscle activity would be expected. The possibility of regional rotation triggered exclusively by regional changes in motor neurons' intrinsic properties seems however unlikely as in such case no differential activation of distinct muscle regions would seem possible, unless MUs of different sizes were confined at specific muscle regions. For example, the selective recruitment of different trapezius regions, crucial for controlling different scapular movements (Holtermann et al 2009), would appear physiologically unviable under the assumption that all motor neurons serving the trapezius muscle receive the same input. Accepting, on the other hand, that the synaptic input may distribute differently across pools of motor neurons serving the same muscle (Desmedt and Godaux 1981), one may advance the spatial differentiation of net excitatory input could contribute to the regional rotation of MUs during fatigue. In this case, rotation and substitution occurring within pools of motor neurons receiving a common input would be triggered by intrinsic changes in the motor neuron excitability (Manning et al 2010), without violating the size principle (see page 1730 in Henneman 1968). Regional rotation, if verified, would most likely result from altered distribution of excitatory and inhibitory inputs to different pools of motor neurons. Although we understand the concept of regional rotation advanced here may appear speculative, we decided nevertheless to bring it to the readers' attention; the verification of regional rotation would advance our current knowledge on the determinants of muscle fatigue. Moreover, our current technology allows for devising experiments aimed at addressing this open issue. The multi-channel approach does seem to provide an additional view to myoelectric manifestations taking place within different sub-volumes of individual muscles prior and during fatigue, otherwise unobserved with the conventional, localized EMG sampling.

Considerations on muscle architecture.
With advances in imaging technology, changes in muscle structure in a number of situations, including fatigue, maybe studied in detail. This section is not aimed at describing the architectural alterations the muscle undergoes with fatigue. Here we are rather focused on providing readers with indications on how these structural changes may affect the surface EMGs. The relationship between changes in muscle architecture and changes in EMG descriptors is largely unexplored and, therefore, in the following paragraphs we briefly develop this argument based on our currently limited understanding.
A first note to point out is that different muscles have different architectures, each likely demanding a specific approach to the study of myoelectric manifestation of fatigue. Muscles may be classified according to the mechanical arrangement of their fibres with respect to the muscle tendon (Lieber and Fridén 2000). From a surface EMG perspective, a more relevant classification should however be based on the arrangement of muscle fibres with respect to the skin surface. Two categories of muscles emerge in this case: (i) muscles with fibres residing in planes parallel to the skin or skin-parallel fibreed muscles; (ii) muscles pinnate in the depth direction or in-depth pinnate muscles. In the former case, electromyographic features of crucial relevance for the study of muscle fatigue may be assessed (see section 2.1.1 CV). Care is however necessary to avoid potential misalignment between electrodes and fibres; in the case of misalignment, CV may be estimated after correcting for the actual fibres' direction with methods specifically designed for the estimation of muscle fibres' direction from the surface EMGs (Grönlund et al 2005, Gallina and Vieira 2015). More specifically, for skin-parallel fibreed muscles, surface EMGs may be detected from skin regions covering the same group of fibres (electrodes positioned parallel to muscle fibres) or different groups of fibres (electrodes positioned transverse to the muscle fibres). As commented earlier in the previous subsection, EMGs collected in parallel and transverse directions to the fibres of skin-parallel fibreed muscles have different meanings. Some examples of muscles that may fall within this category are trapezius, biceps brachii, brachioradialis, abductor pollicis brevis, abductor digit minimi and the vastus medialis and vastus lateralis. Even though the vastus muscles show some degree of pinnation in the depth direction (Blazevich et al 2006), action potential propagation can be clearly observed for each of these quadriceps' heads. From these skin-parallel fibreed muscles, regional variations in activity and in CV may be both assessed with the multi-channel approach in fatiguing contractions.
Most fibres within in-depth pinnate muscles, on the other hand, do not lay in planes parallel to the skin; according to the in-depth pinnate definition, fibres run from the deep to the superficial muscle region. The relative arrangement between fibres and electrodes for an in-depth pinnate muscle is shown in figure 8(a). This schematic representation illustrates a parasagittal view of electrodes and fibres for the gastrocnemius muscle, as observed with ultrasound imaging (Hodson-Tole et al 2013). Gastrocnemius fibres extend from the deep to the superficial aponeurosis (thin, black lines in figure 8b). Electrodes located at skin regions covering the muscle superficial aponeurosis are located over different fibres; from the sixth to the 15th pair of electrodes in figure 8(a). Considering the amplitude of single-differential (i.e. bipolar) action potentials decreases steeply with the distance between the action potentials themselves and the pair of surface electrodes (Lynn et al 1978), it follows that different pairs of electrodes sample from different gastrocnemius fibres. The contribution of action potentials travelling along proximal fibres is for example expected to be much greater to EMGs detected proximally than distally. Our group obtained indeed theoretical and experimental data suggesting that single-differential, surface EMGs collected from the gastrocnemius muscle are somewhat selective (Vieira et al 2010, Mesin et al 2011. It is our experience that action potentials of individual MUs appear clearly in surface EMGs detected by at most five consecutive electrodes (1 cm inter-electrode distance) positioned along the gastrocnemius proximo-distal axis (figure 8(b); figure 7 in Vieira et al (2010); figures 3 and 5 in Vieira et al (2011)). Given the gastrocnemius in-depth pinnate architecture, variations in the amplitude and spectrum of surface EMGs detected along the muscle longitudinal or transverse direction are thus unlikely affected by the anatomical factors (i.e. innervation zone and tendon regions) accounting for the variations in surface EMGs detected in directions parallel to the fibres of skin-parallel fibered muscles (see section 3.1.2). Conversely, differences between surface EMGs detected from different gastrocnemius region are most presumably due to regional variations in activity; that is, changes in activity between different, proximo-distal or medio-lateral, gastrocnemius regions. Although commented for the gastrocnemius muscle, these considerations potentially hold for other muscles with the in-depth pinnate architecture, among which tibialis anterior, biceps femoris and soleus.
Questions of relevance for the study of muscle fatigue therefore arise for in-depth pinnate muscles: is it appropriate to consider CV for the assessment of fatigue-induced changes in these muscles? (ii) if not, are other EMG descriptors sensitive to the myoelectric manifestations of fatigue? With respect to the first question, the most likely answer is no. CV is quantified as the ratio between inter-electrode distance and the delay estimated between surface EMGs detected by consecutive electrodes. Delay estimates are highly dependent on the alignment between electrodes and fibres; the more the fibres and electrodes are in parallel direction the closer the estimate reflects the actual delay resulting from action potential propagation (Farina et al 2002a; see also section 2.1.1.1). As discussed above, the parallel alignment between fibres and electrodes is not possible for most fibres within in-depth pinnate muscles. For the example shown in figure 8, only electrodes covering the very distal region of in-depth pinnate muscles maybe parallel to the muscle fibres (see the first five pairs of electrodes in figure 8(a)). As a consequence, propagation can be only appreciated for the most distal surface EMGs (figure 8b). Indeed, with a grid of 128 (16 × 8) electrodes, Gallina et al (2013a) obtained physiological estimates of CV only for electrodes positioned at the most distal gastrocnemius region. Physiological estimates of CV may be estimated from the distal region of the tibialis anterior muscle as well (Houtman et al 2003). Obtaining physiological CV estimates for the soleus muscle seems however unlikely (Rainoldi et al 2004). The point here is that obtaining physiological estimates of CV is not straightforward for these muscles and if physiological estimates are obtained, they reflect the CV of a small and perhaps unrepresentative muscle portion. From surface EMG, for example, ascertaining whether the CV estimates from the distal region of gastrocnemius and tibialis anterior correlates with estimates from the muscle proximal region is not viable. In spite of the potentially limited, yet unverified representativeness of CV estimates, myoelectric alterations with fatigue seem nevertheless to manifest in surface EMGs detected from in-depth pinnate muscles (Sadamoto et al 1983, Weist et al 2004, Mademli and Arampatzis 2005, Ishikawa et al 2006, Mitsukawa et al 2009. Recently, Mesin et al (2011) specifically assessed simulated and experimental signals to test for whether surface EMGs from the gastrocnemius muscle reflect the changes in amplitude and spectral features typically occurring during fatigue. Decreases in EMG median frequency were observed in both simulated and experimental conditions and for different degrees of in-depth pinnation (see figure 6 in Mesin et al (2011)). Interestingly, these authors observed greater variability for the amplitude and spectral descriptors of experimental than of simulated EMGs. It was suggested that factors not controlled for during simulations could have accounted for the greater variability observed in the experimental signals with fatigue (Mesin et al 2011), with the regional variations in activity within and between plantar flexors (Tamaki et al 1998, McLean andGoudy 2004) being the most likely candidate explaining these discrepancies. In agreement with this view, our group observed local changes in the median frequency of surface EMGs collected with a large grid of electrodes from the gastrocnemius muscle (figure 9) while subjects exerted, until exhaustion, bouts of isometric plantar flexion contractions (5 s contraction at 50% MVC, 5 s rest). Similarly, local changes in surface EMG with fatigue have been recently reported for the in-depth pinnate, rectus femoris muscle (Watanabe et al 2013). While the validity of CV estimates remains unverified, the amplitude and spectral descriptors of surface EMGs detected with the multi-channel approach seem highly sensitive to myoelectric manifestations of fatigue in in-depth pinnate muscles. With the progressive increase in popularity of multi-channel systems, we foresee the appearance of high-quality studies specifically designed to address relevant open issues relating changes in muscle architecture to changes in muscle activation with fatigue.
The multi-channel approach constitutes a promising tool for unveiling the causes and consequences of alterations in muscle structure with fatigue. Changes in muscle structure occurring during fatiguing contractions may affect EMGs directly and indirectly. By direct, we mean alterations in the EMG descriptors resulting directly from structural muscle changes. For example, Mesin et al (2011) have shown the distribution of EMG amplitude changes markedly with the degree of pinnation in the depth direction; the action potential of more pinnate fibres appears with greater amplitude and at a more localized skin region, centred nearby the fibres' superficial attachment (see figures 2-4 in Mesin et al (2011)). With ultrasound imaging, recent studies reported a significant alteration of the architecture of leg muscles with fatigue, even though with conflicting results. While some authors observed an increase in the pinnation angle of the leg muscle with fatigue (Kubo et al 2001, Brancaccio et al 2008, others did not (Ishikawa et al 2006, Mitsukawa et al 2009. In spite of methodological differences potentially accounting for disparities between studies (e.g. muscle studied, time when measurements were taken, etc), these results suggest the changes in muscle architecture may explain, at least in part, the myoelectric manifestations of muscle fatigue. How much architectural changes affect EMGs directly is currently an open issue. Structural alterations with fatigue may also affect EMGs indirectly; that is, changes in EMG descriptors may result from physiological responses triggered by architectural changes with fatigue. More specifically, for instance, Kubo et al (2001) observed an increased compliance of the tendon structures in the knee extensor muscles after a fatigue test, accompanied by a significant ~10 ms increase in the electromechanical delay (i.e. delay between the onset of moment development and the onset of myoelectric activity). Unfortunately, these authors did not report any results on variations in muscle activity with fatigue. It is however possible that changes in tendon compliance may lead to a redistribution of activity within the muscle, optimising force development and transmission. In line with this view, Avancini et al (2015) detected EMGs with greater ampl itude over a significantly larger proximo-distal gastrocnemius region when subjects performed isometric plantar flexion contractions with the knee at a flexed rather than extended position. Interestingly, the proximo-distal redistribution of myoelectric activity was not associated with proximo-distal differences in gastrocnemius architecture. Avancini et al (2015) suggested that Achilles tendon receptors could have mediated the more spatially diffused activation of gastrocnemius volume, which was expected to better compensate for the increased tendon compliance observed at more flexed knee positions (Herbert et al 2002, De Monte et al 2006. Notwithstanding the relevance of understanding whether changes in surface EMGs observed with fatigue may result from changes in muscle architecture rather than, or in addition to, from changes of physiological origin, so far no systematic reports on the relationship between muscle architecture and activation with fatigue were found. The recent development of electrodes transparent to ultrasound  does however provide a suitable means for addressing this issue. After all, it seems therefore the multi-channel approach has a great potential to provide a comprehensive view of electrophysiological muscle alterations with fatigue, regardless of anatomical differences between muscles.

Single motor unit level
As reported in the previous sections, surface EMGs has been often analysed as an interference signal and different mathematical descriptors have been used to track changes in surface EMGs during the development of muscle fatigue. More common techniques are power spectrum or time-frequency distribution (Lindstrom and Magnusson 1977, De Luca 1984, Merletti et al 1990, high-order power spectral moments , spike analysis (Gabriel et al 2007), fractal dimension (Ravier et al 2005), and several variations of these methods (Mesin et al 2009).
In the last years, the limitations of these approaches emerged. The main attention focused on the limitations that apply to inferences on EMG and movement control when the classical approach is used (Farina et al 2004a(Farina et al , 2014 but the same considerations are true in the study of muscle fatigue . The main reason for these limitations is that global signal features extracted from the interference EMGs are influenced by both central and peripheral properties of the neuromuscular system (Farina et al 2002a, Dimitrova andDimitrov 2003). Surface EMG is, in fact, the sum of action potentials generated in muscle fibres belonging to active MUs and consequently, it contains information about the synaptic inputs received by the motor neurons and about the muscle fibre electrical properties as shown in the following model (Farina et al 2014): where s(t) is the recorded EMG signal, M is the number of active MUs, ϕ i (t) is the shape of the action potential for the ith MU, and δ(t − t ij ) is the delta function representing a MU discharge at the jth time (t ij ) at which the ith MU discharges. From the model of surface EMG reported in equation (1) it is clear that it is not possible to extract, from the interference EMG signal, information only related either to the inputs from the central nervous system (represented by δ(t − t ij )) or to the peripheral MU properties (represented by ϕ i (t)). Amplitude variables, for example, depend on both the neural drive (increase when the discharge rates or the number of recruited MUs increase) and the peripheral properties (the energy of the action potentials increases with the decreasing of CV) without the possibility of distinguishing between the two effects. Similarly, frequency variables are influenced by MU recruitment, tissue filtering, and CV of muscle fibres (Dimitrov and Dimitrova 1998); CV estimated from surface EMGs is influenced by MU recruitment and discharge patterns, and by the fibre membrane properties of the active MUs (Gazzoni et al 2001, Farina et al 2004a, McGill and Lateva 2011. Figure 10a shows an example of EMGs detected from the biceps brachii during a 15% MVC isometric contraction held for 10 min: the recruitment of progressively larger MUs from the beginning to the end of the contraction can be recognized. Figure 10(b) shows the time course of global variables: ARV, MNF, and CV estimated from non-overlapping signal epochs of 1 s for three different contraction levels: 5%, 10% and 15% MVC. In the case of 15% MVC it is possible to clearly identify (approximately at second 300) an increase of CV after the initial decrease due to fatigue; the same trend even if less evident is shown by MNF. The increase of CV and MNF can be explained, on the basis of the visual inspection of the raw signals, with the recruitment of new larger MUs (with higher CV) in order to maintain the requested force level when force generation capability of initially recruited MUs decrease because of fatigue. This example highlights that global variables are influenced by both membrane properties and central MU control suggesting that caution must be used in the interpretation of their trend during fatiguing contractions. However, the only way to verify such an interpretation is to gain access to MU firing patterns.
In order to highlight the limitation of time-frequency methods in the study of fatigue, Farina et al (2014) studied, in simulated conditions, the effect of recruitment/derecruitment of MUs and changes of CV due to fatigue on surface EMG variables estimated with a timefrequency method. During simulated ramp contractions, MUs were recruited and de-recruited according to the size principle with higher threshold MUs having higher conduction velocities. Fatigue was simulated by the decreasing of CV of active MUs during the contraction and by the recruitment of new MUs. No association was found between the estimates of instantaneous MNF and the recruitment and derecruitment of MUs. Despite significant changes in the recruited MU pool, the instantaneous MNF of the surface EMGs remained relatively constant during each ramp. Moreover, the instantaneous MNF was lower after the fatiguing contraction despite the higher MU activity (figure 11). These results demonstrate that global spectral estimator of surface EMG is not sensitive to changes in motor unit population during fatigue.
On the basis of these pieces of evidence, study of myoelectric manifestations of fatigue is progressively moving from the analysis of global features extracted from the interference surface EMGs to the identification and analysis of the single MUs. The use of multi-channel detection systems together with advanced processing techniques allows separating the neural information (MU recruitment/de-recruitment and discharge characteristics) from the peripheral information (membrane properties and MU anatomy) in the surface EMGs as showed in figure 12 and described in the following.

Surface EMGs decomposition.
Because of the filtering effect of the volume conductor and the higher number of detected MUs, the identification of single MUs from the surface EMGs is more challenging than from intramuscular EMGs. In order to counteract these effects, early studies addressed the identification of single MUs in surface EMGs with two approaches: the spatial filtering and the spatial sampling. Spatial filtering was used to increase the spatial selectivity of the recording system and then reduce the number of MUs contributing to the signal (Reucher et al 1987a, 1987b, Disselhorst-Klug et al 1999. Spatial sampling (Masuda et al 1985, Blok et al 2002 was used to increase the possibility to discriminate the action potentials of different MUs (Blok et al 2002, Zwarts andStegeman 2003). In fact, the uniqueness of the skin representation of action potentials of single MUs is a necessary condition for decomposition. The larger is the number of recording locations, the lower is the Figure 11. Example of a simulation study of sensitivity of surface EMG power spectral frequencies to motor unit activity. A protocol was simulated including an up-down ramp contraction to 30% of MVC force, a 140 s contraction sustained at 30% MVC force, and a final ramp contraction similar to the 1st one (on the top). The model predicted a decrease in peak force for the 2nd ramp contraction after the fatiguing contraction and more motor unit activity during the 2nd ramp contraction. Instantaneous MNF was computed from the spectrogram of the surface EMG signals simulated during the 2 ramp contractions. Instantaneous MNF remains relatively constant during the ramps. Modified from Farina et al (2014) with permission. likelihood that the detected multi-channel action potentials are similar among MUs (Kleine et al 2007). In order to quantify this observation, Farina et al (2008b) examined the capacity of different electrode configurations (in terms of number, spatial configuration and spatial filtering) to discriminate the action potentials of single MUs in surface EMG recordings. Results showed that, in a population of simulated MUs, one bipolar recording allows the discrimination of less than 5% of the MUs (Farina et al 2008b). Similar results were obtained for experimentally detected MUs (Farina et al 2008b).
In some very specific experimental conditions (combining low contraction levels, spatial sampling, spatial filtering, and EMG feedback) single MUs can be extracted from surface EMGs recorded with multi-channel detection systems with simple processing techniques such as amplitude thresholding (Farina et al 2004a).
More advanced methods based on template matching, the same technique used for intramuscular EMG decomposition, have been proposed to identify few active MUs from the interference surface EMGs in less controlled conditions (Blok et al 2002, Gazzoni et al 2004, Kleine et al 2007. With this technique, MUAPs are detected in the signal using a threshold or a matching filter and then they are clustered on the basis of the shape of each identified action potential, taking advantage from the high number of recording locations to discriminate the action potentials of different MUs. The results allow to obtain information about the recruitment/derecruitment of MUs and to characterize muscle fibre membrane and anatomical properties of single MUs (Gazzoni et al 2004, Kleine et al 2007. With these techniques it is possible to analyse single MUs at relatively high contraction levels and the information provided can be used to evaluate changes in the membrane properties of single MUs during sustained contractions. The limitations of template-matching methods in the reconstruction of the complete firing pattern of active MUs motivated the development of surface EMG decomposition techniques based on blind source separation. Different methods have been developed and proposed in the last years in order to extract the complete firing pattern of MUs form interference EMGs (Holobar and Zazula 2003, De Luca et al 2006, Holobar et al 2009, Nawab et al 2010, Liu et al 2014, Negro et al 2016, Chen and Zhou 2016, Ning et al 2015. Once the discharge times of MUs are obtained from the complete decomposition, it is possible to extract the spatial-temporal representation of MUs by applying the spike-triggered averaging technique to the original surface EMGs. Figure 13 shows an example of complete decomposition of a signal recorded from the abductor pollicis brevis during repetitions of an isometric force ramp contraction. Different MUs, depicted in different colors, are active for different durations. Thus, they demonstrate different levels of fatigue. The average CV of MUs 1-7 decreased across different contractions whereas the CV of MUs 8 and 9 remained constant from the first to the 18th contraction when it start decreasing. Average MU discharge rates (DR) do not vary significantly over time (figure 13b). The action potentials of MUs 1-7 change significantly over the 27 contractions, while much smaller changes are observed for MUs 8 and 9. The same colour in figure 13 (a)-(c) denotes the same MU. In figure 13 (d), different colours denote different contractions.

2.2.2.
Analysis of fatigue at single motor unit level. Using the techniques described in the previous section, research activity is progressively focusing on the study of fatigue at the level of single MUs in order to clarify some open questions about peripheral and central adaptations to fatigue. It is in fact generally accepted that fatigue is related to (1) peripheral changes at the level of the muscle (such as changes at the neuromuscular junction and the accumulation of metabolites) which result in a reduced ability of muscles to generate force and that are referred to as 'peripheral fatigue' and (2) to central factors (such as impaired cortico-spinal transmission to the spinal cord, reduced motor neuron excitability and firing frequency) which result in a reduced drive of the central nervous system to the motor neurons and that are referred to as 'central fatigue' (Gandevia 2001, Enoka et al 2011. However, the sources of muscle fatigue and evidence of central fatigue are still debated (Contessa et al 2016).
By using simple threshold-based techniques, Farina et al (2004b) conducted a study to get a deeper insight into the effect of ischemia on muscle fatigue. Several studies showed a decrease of mean power spectral frequency and average CV due to ischemia and interpreted these results as due to the accumulation of metabolic byproducts in the muscle. However, as stated before, average CV is also affected by MU recruitment-derecruitment and by the distribution of discharge rates and it does not provide direct indication on the modifications of fibre membrane properties in single MUs. In Farina et al (2004b) the similar relative changes in CV and mean power spectral frequency observed at the muscle level and at the level of single MUs confirm that a change in CV is the main phenomenon responsible for the observed modifications in surface EMG variables during sustained contractions.
With template matching techniques Gazzoni et al analyzed the CV of single MUs before and after a low-force sustained contraction during which they were not activated (Gazzoni et al 2005). The study showed that sustained contraction determines a decrease in CV of MUs partly independently of their activation and that the decrease of CV of MUs not recruited during the sustained contraction was similar to the decrease in average CV of the recruited MU pool. The analysis of changes of single MU CV adds new knowledge for the interpretation of changes in global surface EMG variables during fatiguing contractions. For example, consider the trend of global variables during submaximal endurance contractions reported in figure 10, where the average muscle fibre velocity shows a decrease followed by an increase. This phenomenon has been demonstrated to be related to the progressive recruitment of new MUs in order to sustain the contraction level (Gazzoni et al 2001, Houtman et al 2003). However, considering that CV decreases also in non-active MUs (Gazzoni et al 2005), the observed increase in average CV is probably smaller than that expected under the assumption of recruitment of fresh MUs. Moreover, results reported in Gazzoni et al (2005) underline at least two additional factors in the interpretation of changes of surface EMG amplitude with fatigue. It has been observed that, during submaximal fatiguing contractions, EMG amplitude increases but remains significantly lower than maximum at the endurance time (Fuglevand Figure 13. Discharge patterns of nine MUs identified from surface EMG detected from the abductor pollicis brevis during 27 repetitions of ramp contractions (from 0% to 10% MVC) are shown in (a). Each dot indicates a single MU discharge at a given instant, whereas its relative vertical displacement codes the instantaneous MU discharge rate. Discharge rate, and CV estimated for each MU are reported in the panels (b) and (c) respectively. Panel (d) shows the progressive changes of MUAP shape due to fatigue for two representative MUs detected by one pair of electrodes (from Farina and Holobar (2016) with permission). et al 1993). Changed membrane properties of fibres in the muscle may result in modifications of the twitch force and then a decrease in discharge rates of newly recruited MUs for maintaining tetanic fusion with respect to the condition of fresh fibres. Moreover, decreased CV increases surface EMG amplitude cancellation (Keenan et al 2006, Farina et al 2008a. Even if methods for the analysis of single MU based on template matching techniques allow to get a deeper insight into peripheral changes of single MUs with fatigue, these methods are in general not able to separate action potentials of different MUs which overlap in time. Template matching methods therefore do not allow studying changes in central control due to fatigue because they provide an incomplete MU discharge patterns. Kallenberg and Hermens (2008) conducted a study in order to verify if the number of MUAPs per second obtained from an incomplete decomposition reflects the input of the central nervous system to a muscle. Results showed that such index reflects more selectively than frequency and amplitude global variables the central motor control at low contraction levels. However, the assessment of the complete firing patterns of MUs is needed to study in details the input of the central nervous system to a muscle.
Taking advantage from techniques recently developed and validated for the complete decomposition of EMGs (De Luca et al 2006, Farina et al 2010, Holobar and Farina 2014, Negro et al 2014 it is possible to investigate both central and peripheral changes during fatigue. McManus et al (2015) investigated global surface EMG variables and single motor unit properties (considered collectively for each subject) to gain insight into changes due to muscle fatigue immediately after a sustained fatiguing contraction and in the following recovery period. Changes in surface EMG reflected changes in MUAP waveforms and authors obtained a strong correlation between the MUAP amplitude and the surface EMG RMS and between the MUAP duration and the surface EMG MNF. From the analysis of MU recruitment and firing rates, McManus and colleagues observed reduced MU firing rates and increased central drive, concluding that MU recruitment is favored over rate coding to maintain force during fatigue.
On the basis of the observation that specific adjustments can occur for specific MUs within a single muscle (Carpentier et al 2001), Farina et al (2010) studied the relationship between the adjustments of both peripheral and central characteristics exhibited by low-threshold MUs and the relative proportion of time each MU was active during fatiguing ramp isometric contractions. The relationship between changes in muscle fibre CV and recruitment and derecruitment thresholds were investigated; this requires a precise identification of firing pattern of MUs (provided only by complete decomposition) and reliable estimation of membrane properties (obtained using the multi-channel approach for the estimation of MU CV). The results showed that adjustments differed among low-threshold MUs and changes depended on the relative duration that each motor unit was active during the task.
In a recent work, Contessa et al (2009) investigated the interaction between motor unit firing behavior and muscle force during fatigue. By applying surface EMG decomposition (Nawab et al 2010) the authors analyzed MU firing patterns during repeated voluntary submaximal isometric contractions of the vastus lateralis muscle. The observed increase of MU firing rates and decrease of recruitment thresholds of all MUs together with simulated findings suggest the excitation to the motoneuron pool adjusts motor unit firing behavior to compensate for the changing of muscle force twitch.

Conclusions and future perspectives
Much remains to be understood about the central and peripheral processes that contribute to muscle fatigue. This review analyzed the new knowledge and new issues open by more recent advances in surface EMG detection and processing, herein referred as multi-channel approach. Three main topics were considered and analysed: (i) the quantification of a physiological variable of marked interest to the study of fatigue, the propagation velocity, obtained with multi-channel detection systems and advanced algorithms; (ii) the assessment of regional, myoelectric manifestations of fatigue and the analysis of the relationship between changes in muscle architecture and changes in EMG descriptors; (iii) the analysis of surface EMGs at the level of single MUs with the possibility to obtain information about motor unit control and fibre membrane changes.
The multi-channel approach does seem to have great potential in providing a comprehensive view of electrophysiological muscle alterations taking place within different sub-volumes of individual muscles prior and during fatigue, that cannot be observed with the conventional EMG sampling. Moreover, it represents a promising tool for unveiling the causes and consequences of alterations in muscle structure with fatigue.
The possibility to extract information about individual MUs overcomes the difficulty of the interpretation of surface EMG global variables often prone to be erroneous. multi-channel approach can contribute to the understanding of the interactions between the nervous system and the muscle during fatiguing contractions.
However, many issues are still open and new ones emerged from the multi-channel approach. Among others the study of changes in fatigability following a treatment or training and the study of fatigue during dynamic contractions.
Addressing the first issue could open interesting perspectives in rehabilitation. For this purpose, the validation of currently available decomposition methods must be extended to verify whether they provide comparable results in different measurement sessions. A recent study (Martinez-Valdes et al 2015) assessed the intra and inter-session reliability of MU properties estimated from the decomposed surface EMGs. Results suggest that current techniques allow assessing central and peripheral motor unit properties within and across experimental sessions in a reliable way. This opens new possibilities in the evaluation of fatigue changes in longitudinal studies but additional research is required to confirm these results in a wider range of conditions.
With regards to the study of fatigue during dynamic contractions, in the last years, the scientific interest is more and more moving toward dynamic contractions also thanks to the availability of new wearable acquisition systems. The new insights obtained from the multichannel approach highlight the importance of accounting for anatomical and geometrical factors on the detected surface EMGs. The importance of considering these factors increases in dynamic contractions. Moreover, to study fatigue during dynamic contractions at the level of single MUs, new algorithms are required. To date, surface EMG decomposition algorithms have been validated and used in isometric contractions. Only few works are published in literature in this field (Holobar et al 2008, De Luca et al 2015. With the progressive increase in popularity of multi-channel systems, we foresee the appearance of high-quality studies specifically designed to address relevant open issues relating muscle fatigue.