Non‐invasive estimation of muscle fibre size from high‐density electromyography

Because of the biophysical relation between muscle fibre diameter and the propagation velocity of action potentials along the muscle fibres, motor unit conduction velocity could be a non‐invasive index of muscle fibre size in humans. However, the relation between motor unit conduction velocity and fibre size has been only assessed indirectly in animal models and in human patients with invasive intramuscular EMG recordings, or it has been mathematically derived from computer simulations. By combining advanced non‐invasive techniques to record motor unit activity in vivo, i.e. high‐density surface EMG, with the gold standard technique for muscle tissue sampling, i.e. muscle biopsy, here we investigated the relation between the conduction velocity of populations of motor units identified from the biceps brachii muscle, and muscle fibre diameter. We demonstrate the possibility of predicting muscle fibre diameter (R2 = 0.66) and cross‐sectional area (R2 = 0.65) from conduction velocity estimates with low systematic bias (∼2% and ∼4% respectively) and a relatively low margin of individual error (∼8% and ∼16%, respectively). The proposed neuromuscular interface opens new perspectives in the use of high‐density EMG as a non‐invasive tool to estimate muscle fibre size without the need of surgical biopsy sampling. The non‐invasive nature of high‐density surface EMG for the assessment of muscle fibre size may be useful in studies monitoring child development, ageing, space and exercise physiology, although the applicability and validity of the proposed methodology need to be more directly assessed in these specific populations by future studies.


Introduction
The motor unit consists of a motor neuron and the group of muscle fibres it innervates (Heckman & Enoka, 2004;Liddell & Sherrington, 1925). Skeletal muscles are composed of muscle fibres, whose heterogeneity results from the broad range of properties of both the innervating motor neuron and the structural and biochemical composition of the fibres (Duchateau & Enoka, 2022;Heckman & Enoka, 2004).
Since the pioneering work of Adrian & Bronk (1928), different techniques have been developed to record motor unit activity in humans either invasively or non-invasively (Duchateau & Enoka, 2011). At present, intramuscular EMG recordings remain the classic means for investigating the activity and properties of individual motor units. Nevertheless, advances in electrode fabrication, biosignal processing, and modelling over the last 20 years have contributed to expanding the opportunities to study motor unit activity in humans non-invasively . For instance, the recent introduction of high-density surface EMG has greatly helped to overcome some of the limitations of intramuscular EMG recordings Merletti et al., 2008;Negro et al., 2016), making it possible to decode and monitor the concurrent activity of large and representative populations of motor units non-invasively Enoka, 2019).
On the other hand, the muscle biopsy technique provides an invaluable opportunity to study several characteristics of skeletal muscle tissue and muscle fibre properties, including fibre type composition and cross-sectional area (Bergström, 1962;Henriksson, 1979). Since the seminal work by Charriere & Duchenne (1865), methodological advances in biopsy sampling and processing of the specimen have contributed to spreading the use of the muscle biopsy procedure beyond its original context of clinical diagnosis of neuromuscular disorders, to become the gold standard technique for the evaluation of fibre morphology within general physiology and exercise science (Ekblom, 2017;Folland & Williams, 2007). For instance, the use of muscle biopsy has been fundamental for our understanding of muscle function, dysfunction and adaptability (Jones et al., 1989;McDonagh & Davies, 1984;Pette & Staron, 1997;Schiaffino & Reggiani, 2011).
However, although safe and with a low rate of complications, the biopsy method requires the insertion of a needle through the skin and muscle fascia, which exposes the individual to some risks both during and after the procedure (Ekblom, 2017;Tarnopolsky et al., 2011). The poor reproducibility of this technique for fibre area measurements (Blomstrand et al., 1984;Bottinelli, 2001;Folland & Williams, 2007;Lexell & Taylor, 1989) and the fact that it provides morphological information dissociated from the neural control of the sampled muscle fibres  have also been well-documented. Thus, it would be enticing to identify alternative physiological parameters that could inform us about the size of muscle fibres belonging to individual motor units non-invasively. J Physiol 601.10 In this respect, there have been previous attempts at associating EMG-derived parameters to muscle fibre size (Blijham et al., 2006;Methenitis et al., 2016;Pope et al., 2016). One potential candidate to indirectly estimate fibre size is muscle fibre conduction velocity, which can be measured from surface EMG. The inverse of the average latency between surface EMG signals recorded from adjacent electrodes aligned in the direction of the muscle fibres, provides an estimate of the average propagation velocity of the action potentials of several concurrently active motor units (Farina & Merletti, 2004). Alternatively, by decomposing the surface EMG signal and separating action potentials of individual motor units Holobar & Zazula, 2007), it is possible to estimate the motor unit conduction velocity (Farina et al., 2001(Farina et al., , 2002. Although both parameters can be estimated non-invasively from surface EMG, a major advantage of motor unit conduction velocity over muscle fibre conduction velocity estimates is that the former provide information about individual motor unit properties and not the average properties of all active motor units. Motor unit conduction velocity is considered a 'size principle parameter' , due to its linear association with motor unit recruitment threshold (Andreassen & Arendt-Nielsen, 1987;Del Vecchio et al., 2017). This linear association is due to the direct relation between muscle fibre diameter and conduction velocity, observed with invasive techniques in animal models (Hakansson, 1956) and humans (Blijham et al., 2006;Methenitis et al., 2016). Therefore, the possibility of estimating the conduction velocity for a large population of motor units non-invasively would theoretically allow an estimate of the distribution of muscle fibre sizes alternative to muscle biopsies.
By combining high-density surface EMG  with the gold standard technique for muscle tissue sampling (i.e. muscle biopsy), here we systematically investigated the relation between the conduction velocity of large populations of motor units identified non-invasively during voluntary contractions, and muscle fibre size, i.e. fibre diameter and cross-sectional area, directly derived from muscle biopsies. We studied this relation in the biceps brachii muscle of a heterogeneous population consisting of healthy untrained and chronically resistance-trained subjects, in order to have a sample of widely varying muscle fibre size. From this analysis, we then propose a new methodology to estimate the average fibre diameter and cross-sectional area from surface EMG processing without the need for muscle biopsies.
The main purpose of this study was to determine how motor unit conduction velocity and muscle fibre size (fibre diameter and cross-sectional area) were related in a heterogeneous population of healthy adults, by considering inter-participant variability of both parameters. We hypothesized that the conduction velocity, assessed in the present study from large populations of concurrently active motor units, would have a strong association with fibre diameter, thus becoming a novel non-invasive marker of muscle fibre size. We reveal that it is possible to accurately transform motor unit conduction velocity values into a distribution of fibre diameters, therefore opening new avenues for research in ageing, training and neuromuscular disease.

Overview of the study
Participants visited the laboratory on three occasions for measurements of the elbow flexors, and specifically biceps brachii muscle, of the non-dominant arm (defined as the non-writing hand). The familiarization session (visit 1) involved habituation with the experimental procedures and the performance of a series of maximum and submaximal isometric contractions with the elbow flexors. Participants' physical activity was assessed with the International Physical Activity Questionnaire (IPAQ, short format; Craig et al., 2003), supplemented with a resistance training experience questionnaire and follow-up individual discussion. These measurements ensured that only participants that fit the criteria for the two groups (described below) were retained in the study. The main neuromuscular assessment (visit 2), 7-10 days after visit 1, involved the simultaneous recording of elbow flexor force with isometric dynamometry, and high-density surface EMG for motor unit conduction velocity estimation from the biceps brachii muscle during both maximum voluntary contractions and submaximal force-matching ramp contractions. Prior to high-density surface EMG, participants underwent B-mode ultrasonography to quantify the muscle-electrode distance, i.e. an indicator of subcutaneous adipose tissue thickness, and its influence on conduction velocity estimation. Approximately, 7-8 days after visit 2, participants reported to the laboratory for the muscle biopsy (visit 3), which was taken in the region of EMG recordings, for the assessment of fibre diameter and cross-sectional area. Participants were instructed not to participate in demanding physical activity (48 h) and to avoid caffeine (24 h) consumption before each session. An overview of the study is shown in Fig. 1.

Participants
The participants involved in the present study were the same as in our recent publication , which compared motor unit behaviour during voluntary contractions in chronically resistance-trained Figure 1. Overview of the study and experimental setup A, muscle biopsy was taken from the biceps brachii muscle from the same region where electromyographic recordings were made. B, transverse serial cross sections were obtained using a cryotome before incubation with a primary antibody for myosin heavy chain I. C, images were captured using a fluorescence microscope. Fibre cross-sectional area was assessed for ∼200 different fibres per participant. Representative examples of type I and II fibre cross-sectional area measured from one resistance-trained (RT) and one untrained (UT) participant are displayed. D, experimental configuration for the main neuromuscular assessment. Two bi-dimensional electrode grids were placed over the biceps brachii muscle. E, representative example of a linearly increasing isometric ramp contraction to 70% MVF (force signal in black), with concomitant monopolar recordings from 13 electrodes from one column of the grid (displayed in different colours). F, colour-coded raster plots of the individual motor units (n = 15) identified from biceps brachii muscle for one representative participant. G, motor unit action potential waveforms from the individual pulse trains of the identified motor units were obtained via the spike triggered averaging technique. Clear propagation of a single motor unit action potential from the innervation zone (IZ, in red) both to the proximal and to the distal tendon regions along one electrode column can be observed. H, motor unit conduction velocity (MUCV) was estimated with a validated multichannel maximum likelihood algorithm. A J Physiol 601.10 representative example of a low threshold motor unit (LT MU, in red) and of a high threshold motor unit (HT MU, in blue) is shown. The action potentials of the HT MU (5.04 m s −1 ) propagate with a greater velocity than that of the LT MU (4.27 m s −1 ). [Colour figure can be viewed at wileyonlinelibrary.com] (RT) athletes versus a cohort of healthy untrained (UT) individuals. Here we report data from motor unit conduction velocity estimated from high-density surface EMG and muscle fibre size parameters derived directly from muscle biopsies.
Thirty healthy young males volunteered to take part to this cross-sectional study and provided written informed consent prior to their participation. As general requisites for inclusion, volunteers were required to be aged 18−40 years, without any underlying health issues, with no history of traumatic upper body injury and/or surgery and with no self-reported history of androgenic-anabolic steroid use. Specific inclusion criteria for the RT cohort were extensive history of upper limb resistance training (≥2 sessions week −1 for ≥10 months year −1 , for ≥3 years) and elbow flexion isometric maximum voluntary torque >90 N m. Specific inclusion criteria for the UT group were no history of regular upper limb resistance training and no involvement in systematic physical training in the 18 months prior to their participation in the study. Considering the aim of this current analysis, participants were treated as one single heterogeneous group, presenting a wide variability in terms of elbow flexor maximal voluntary isometric force, muscle fibre size and anatomy (Table 1). One participant was excluded from the analyses because the biopsy sample did not yield sufficient useful muscle tissue. Thus, here we present the results for twenty-nine participants ( Table 1).
The experiments were carried out at Loughborough University (UK), the protocols and procedures were approved by the Loughborough University Ethics Review Sub-Committee (R17 -P174) and conformed to the requirements of the Declaration of Helsinki, except for registration in a database.
Experimental procedures Neuromuscular assessment. In visit 2, after a standardized warm-up , participants performed three to four maximal voluntary contractions with the elbow flexors of their non-dominant arm, and were instructed to 'pull as hard as possible' for 3−5 s. Rest between contractions was ≥30 s, and strong verbal encouragement to overcome the peak force of the previous contraction was displayed with a horizontal cursor. The highest instantaneous force recorded during any maximal voluntary contraction was defined as elbow flexor maximum voluntary isometric force (MVF), and was used as reference to define the submaximal force levels.
Five minutes after the maximal voluntary contractions, participants performed eight submaximal ramp contractions: two contractions to each force level of 15%, 35%, 50%, and 70%MVF in a randomized order. The contractions involved a linear ascending ramp phase increasing force from rest to the specified target force level, a plateau phase of 10 s (for 15−35%MVF trials) or 5 s (for 50−70%MVF trials) of constant force production at the target force level. The rate of force increase was kept constant at 10%MVF s −1 for all contractions, and trials were separated by 3−5 min of rest. During the aforementioned tasks, participants were instructed to match as precisely as possible a visual template displayed on a monitor positioned at a 1 m distance from the participants' eyes, for 5 s before and throughout each contraction. The real-time force output was superimposed on the template for visual feedback.
Force recording. Participants' neuromuscular force production was recorded whilst they were seated in a custom-built isometric elbow flexion dynamometer, consisting of a rigid strength-testing chair, which was adjusted to the participants' height and upper arm length. Participants were seated in an upright position (hip joint angle at ∼90°) with their non-dominant shoulder at ∼90°of flexion (i.e. perpendicular to the trunk) and slightly horizontally abducted at ∼10°, with the elbow at ∼70°of flexion (0°= full extension), and with forearm in half-supination ∼45°(∼0°= anatomical position) . The wrist was tightly strapped to an adjustable wrist brace in series with a calibrated S-beam strain gauge (Force Logic, Swallowfield, UK), positioned perpendicular to forearm rotation. Participants were firmly strapped to the seat back across the pelvis, chest and shoulder to limit extraneous movements.
The analog force signal was amplified (×200), digitized at 2048 Hz using the same acquisition system as that used for high-density surface EMG (EMG-Quattrocento, OT Bioelettronica, Turin, Italy), and recorded with the software OT BioLab Ver. 2.0.6352.0 (OT Bioelettronica). Real-time force visual feedback as well as the ramp templates were provided with the computer software Spike2 (CED, Cambridge, UK).
Ultrasound recording. Prior to electrode placement, ultrasonographic (US) images of biceps brachii muscle were captured at rest using B-mode ultrasonography (EUB-8500; Hitachi Medical systems UK Ltd, Wellingborough, UK) and a 92 mm, 5−10 MHz linear-array transducer (EUP-L53L). Resting US images were captured with the participants standing in an upright position with their non-dominant shoulder at ∼90°of abduction (i.e. humerus perpendicular to the trunk), with elbow fully extended and forearm supinated. Participants were instructed to stand still and relax their arm, elbow and shoulder, which were supported by adjustable padding. A surgical pen was used to delineate the profiles of the two biceps brachii muscle bellies (i.e. short and long head) firstly identified through palpation, and a central line running along the length of each belly at 50% of the medio-lateral width was drawn. Narrow echo-absorbent markers were placed on the dermal surface (perpendicular to the length of the humerus) from the muscle-tendon junction of biceps brachii (i.e. 0 cm) at 5 cm intervals (i.e. 5, 10, 15 cm) along the length of the muscle bellies. Water soluble transmission gel was then applied above each muscle belly to optimize US image detection along the entire length of the muscle. The transducer, also coated with water soluble transmission gel, was gradually moved with minimal pressure applied on the dermal surface along the centre line of each belly (i.e. short head and long head) from the distal biceps brachii muscle-tendon junction to the proximal end of the muscle. A single sweep of US images was performed and recorded (see below) for each muscle belly.
US measurement was used to quantify the muscle-electrode distance, which was adopted as an indicator of subcutaneous adipose tissue thickness interposed between the muscle belly and the recording electrodes. The echo-absorbent markers were adopted to align three separate images (i.e. 0−5 cm, 5−10 cm and 10 cm+). This allowed measurement of muscle-electrode distance, i.e. distance from the surface of the skin to the muscle fascia, at 1 cm intervals along the length of biceps brachii muscle. Specifically, muscle-electrode distance was measured starting at −2 cm distal relative to the placement of the distal edge of the electrode grid, and then every 1 cm until 13 cm above the starting location (electrode grid length of 10.9 cm). This approach was used to at least partly account for distal movement of the underlying muscle beneath the electrode grid, between the US measurement (i.e. elbow fully extended) and elbow flexion contractions (i.e. elbow at ∼70°of flexion). Furthermore, in order to assess the influence of muscle-electrode distance on motor unit conduction velocity estimation at the individual subject level, only the average muscle-electrode distance values measured from the short head of biceps brachii were correlated with average conduction velocity values (MUCV AVERAGE ) recorded from the muscle region underlying the electrode grid.
The video output from the US was transferred to a computer and US images were recorded through Ezcap video capture software (EZCAP. Forward Video Technology Co. Ltd. Shenzen, PRC). US image analyses were performed with the Tracker image analysis software (Ver. 5.0.6 http://www.physlets.org/tracker/). High-density surface EMG recording. Myoelectrical activity was recorded from the biceps brachii muscle during the submaximal isometric contractions using two high-density surface adhesive grids with 64 electrodes each (13 rows × 5 columns, 1 mm electrode diameter, 8 mm inter-electrode distance; GR08MM1305, OT Bioelettronica; Fig. 1B). This electrode configuration formed an array of 128 electrodes covering the short and long heads of the biceps brachii. To ensure an optimal placement and orientation of the electrodes, an experienced investigator identified the two bellies of the biceps brachii long and short heads through US images, and marked their profiles with a surgical marker. The portion of the skin identified using the criteria defined above was shaved, gently abraded and cleansed with alcohol wipes. The electrode grids were positioned with the participant's non-dominant arm in the testing configuration (∼70°of elbow flexion). The grids were centred between the proximal and distal tendons of the biceps brachii muscle (average distance of the distal side of the array from the antecubital fossa, 5.9 ± 1.2 cm), with the electrode columns aligned with the assumed direction of muscle fascicles Casolo, Nuccio et al., 2020;Del Vecchio, Negro, Falla et al., 2018). Disposable bi-adhesive foam layers, whose cavities were filled with conductive paste (SpesMedica, Battipaglia, J Physiol 601.10 Italy) to optimize the skin-to-electrode contact, were used to attach the electrode grids to the skin overlying the biceps brachii. The large number of electrodes covering the biceps brachii muscle bellies allowed the identification of the innervation zone in the proximal portion of the muscle, and the visual selection of recording electrodes with propagating motor unit action potentials along the entire length of the grids. Reference electrodes were located on the ulna styloid process (i.e. main ground electrode) and the radial styloid process (i.e. high-density grid references).
High-density surface EMG signals were acquired in monopolar derivation, amplified (×150), band-pass filtered (10-500 Hz) and digitized at a sampling rate of 2048 Hz using a multichannel 16-bit analog-to-digital acquisition system (EMG-Quattrocento, OT Bioelettronica). The signals were recorded with the software OT BioLab Ver. 2.0.6352.0 (OT Bioelettronica) and synchronized with the force signal at source.
Muscle biopsy. During visit 3 (8.4 ± 5.4 days from visit 2), muscle biopsies were taken from the biceps brachii under local anaesthesia (1% lidocaine) using the Weil-Blakesley conchotome technique (Ekblom, 2017). Muscle biopsies were taken from the distal anterior biceps brachii (average distance of the biopsy from the antecubital fossa, 6.6 ± 0.9 cm) in the same region the electromyographic recordings were made. Muscle samples were dissected of any visible connective tissue and fat, and then immediately embedded in a mounting medium (Tissue-Tek O.C.T. Compound, Sakura Finetek Europe, Alphen aan den Rijn, The Netherlands) and frozen in liquid nitrogen-cooled isopentane. Samples were then stored at −80°C for immunohistochemistry analysis.

Data analysis
Force analysis. The acquired force signal was converted to newtons (N) and low-pass filtered (4th order, zero-lag Butterworth filter, 15 Hz cut-off frequency). The offset was removed by correcting for the effect of gravity. Only one of the two isometric ramp contractions at each force target (15, 35, 50, 70%MVF) was analysed for each participant. The selection criteria was the lowest deviation of the force trajectory from the given template Nuccio et al., 2020).
High-density surface EMG decomposition. Monopolar electromyographic signals were initially band-pass filtered (2nd order, zero-lag Butterworth filter, 20−500 Hz cut-off frequencies) and visually inspected for noise and artifacts. Channels with noise or artifacts were excluded from further analysis. The filtered signals were decomposed into individual motor unit pulse trains with the Convolution Kernel Compensation decomposition algorithm implemented in the DEMUSE tool software (Ver. 5.01; The University of Maribor, Slovenia) (Holobar & Zazula, 2007;. This decomposition approach accurately identifies motor unit discharge times over a broad range of volitional forces and has been extensively validated using experimental (Holobar et al., 2010) or simulated  signals. After the automatic identification of individual pulse trains, the decomposition accuracy was first assessed by calculating the pulse-to-noise ratio for each motor unit. The pulse-to-noise ratio metric is correlated with both the sensitivity and false alarm rate in the identification of motor unit discharges . The identified motor unit pulse trains from both grids of electrodes were visually inspected and manually edited by an experienced investigator to check for false positives and false negatives (Hug et al., 2021) according to available guidelines . Only motor units exhibiting a reliable discharge pattern with a pulse-to-noise ratio >30 dB (sensitivity >90%, false alarm rates <2%) and/or interspike intervals <2s were retained. Furthermore, only motor units extracted from the grid of electrodes, which showed the highest number of reliably identified motor units, were included in the following analyses. In all cases (29/29), this was the grid of electrodes located on the short head of the biceps brachii .
For each retained motor unit, basic properties were extracted. Specifically, the normalized recruitment threshold and the discharge rate in the ascending phase of the ramp contractions were computed for all identified motor units. Normalized motor unit recruitment threshold was defined as the percentage of force (%MVF) generated by the elbow flexors at which the first motor unit action potential was discharged. Motor unit recruitment discharge rate was calculated as the average of the first 20 discharges (i.e. action potentials) of each motor unit during the ascending phase of the ramp. This number of discharge timings was chosen to minimize the effects of interspike interval variations both on the computation of motor unit discharge rate in the ascending phase of the ramp contractions and on the estimation of motor unit conduction velocity Farina et al., 2002). For clarity, data across the four different target forces (15, 35, 50, 70%MVF) for each subject were collapsed to produce overall motor unit properties (i.e. motor unit recruitment threshold, discharge rate and conduction velocity), irrespective of contraction level .
Motor unit conduction velocity estimation. Since the decomposition algorithm identifies only individual pulse trains but not the waveforms of the corresponding multi-channel motor unit action potentials, these were obtained via the spike triggered averaging technique Farina et al., 2002). The multichannel motor unit action potential waveforms were estimated by averaging the recorded monopolar (raw) EMG signals over 15 ms intervals (action potential duration), using the first 20 discharge times of each motor unit as triggers Farina et al., 2002).
Double differential derivations were computed by differentiating the averaged monopolar motor unit action potentials along the columns of the electrode grids, and were used for motor unit conduction velocity estimation (Farina et al., 2002). A minimum of three up to a maximum of eight double differential channels from the same electrode column were visually inspected and selected for conduction velocity estimation for each identified motor unit. The visual selection of the highest number of double differential EMG channels yielding the clearest propagation of motor unit action potentials along the electrode columns with minimal change in action potential waveform (from innervation zone to the distal tendon region) and with the highest cross-correlation between consecutive channels (≥0.70) is the most accurate approach for EMG channel selection and for motor unit conduction velocity estimation Del Vecchio et al., 2017;Farina et al., 2002). Once the channels were selected, motor unit conduction velocity was estimated through a validated multichannel maximum likelihood algorithm that ensures the accurate calculation of propagation velocity in single motor units with a low standard deviation (<0.1 m s −1 ) (Farina et al., 2001) and estimation errors (2-3%) (Farina et al., 2002). Firstly, motor unit conduction velocity was estimated for each of the 701 motor units identified from biceps brachii muscles of all the participants. Secondly, the MUCV AVERAGE of the identified and active motor units in the ascending phase of the ramp contractions was computed for each participant.
Muscle fibre analysis. Transverse serial cross-sections (8 μm thick) of muscle tissue were obtained using a cryotome and placed onto poly-l-lysine-coated glass slides. Sections were fixed for 10 min in 3.7% paraformaldehyde at room temperature and blocked with Tris-buffered saline containing 5% goat serum, 2% bovine serum albumin and 0.2% Triton X-100 for 1 h at room temperature. Serial muscle sections were then incubated with a primary antibody for myosin heavy chain I (A4.951, Developmental Studies Hybridoma Bank, Iowa City, IA, USA) and diluted 1:200 in the blocking solution for 1 h at room temperature. Sections were then incubated for 2 h at room temperature with an appropriate secondary antibody consisting of goat anti-mouse Alexa Fluor 488 (A11029, Fisher Scientific, Pittsburgh, PA, USA) diluted 1:500 and wheat germ agglutinin Alexa Fluor 350 Conjugate (W11263, Fisher Scientific) diluted 1:20 in the blocking solution. Following incubation, coverslips were mounted with Fluoromount aqueous mounting medium (F-4680, Sigma-Aldrich, St Louis, MO, USA).
Images were captured at ×20 magnification using a fluorescence microscope (Leica DM2500; Leica Microsystems, Wetzlar, Germany). All fibre image analyses were performed using Fiji (ImageJ) software (Schindelin et al., 2012), and the investigator was blinded to the participant group allocation. Any fibres that were clearly oblique or not transverse to the longitudinal axis of the fibre were excluded from the analyses (Mizuno, 1991). Fibre cross-sectional area was assessed by manually drawing around the perimeter of each fibre for ∼200 different fibres per participant. In a small number of participants (n = 4), only 145−165 fibres were analysed for area due to limited number of clear images/fibre perimeters.
Estimation of muscle fibre size from motor unit conduction velocity. First, by assuming that muscle fibres are of uniform cylindrical shape (i.e. circular in cross-section), we converted the single muscle fibre cross-sectional area values (n = 5619) obtained from muscle biopsies into measured muscle fibre diameters, as follows: Measured fibre diameter (μm)

=
Fibre cross sectional area/π × 2 (1) Second, we converted single conduction velocity values estimated from the overall population of identified motor units across all subjects into estimated muscle fibre diameters with the equation of Nandedkar & Stålberg (1983), as follows: This equation was based on preliminary evidence that reported a biophysical proportionality between the size of a biological conductor, such as a muscle fibre (or axon), and the propagation velocity of action potentials in both animal (Hakansson, 1956) and human skeletal muscle (Stålberg, 1966). The two parameters K 1 and K 2 have been estimated in previous studies, but in conditions different from the current study (e.g. with intramuscular EMG signals and with electrical stimulation of the fibres). Therefore, they needed to be identified (calibrated) based on the current data. For calibration, we fitted the estimated fibre diameter data derived from Eqn (2) with the experimentally derived measured fibre diameter data (Eqn (1)) and identified the optimal values of the parameters K 1 and K 2 that minimized the mean square error between the estimated and measured distributions of fibre diameters. This led to the following equation: Estimated fibre diameter (μm) = Motor unit conduction velocity m s −1 − 2.9 /0.021 Equation (3) was derived for the first time from inter-subject observations of both measured fibre diameters from muscle biopsies and diameters estimated from motor unit conduction velocity. To validate the proposed procedure in terms of predictability, we then adopted a leave-one-out calibration. The leave-one-out approach consisted in estimating the parameters in Eqn (2) using all dataset except for one subject and repeating the procedure for all subjects. Equation (2) with the parameter values optimized on all subjects except for one was used to estimate the fibre diameter distribution from the conduction velocity distribution for the left-out subject. In this way, the estimates of fibre diameter for each subject were obtained without the use of the biopsy data of the subject, but rather using a database of motor unit conduction velocity and biopsy data from a separate group of subjects. This analysis provided a direct indication of the quality of estimates of fibre diameter directly from motor unit conduction velocity data. From the leave-one-out estimates of fibre diameter of all subjects, correlation and regression analyses were performed to assess the relation between measured and estimated fibre diameter at the individual subject level.
Lastly, we also explored the relation between muscle fibre cross-sectional area estimated from conduction velocity and muscle fibre cross-sectional area directly measured from muscle biopsies, with the same leave-one-out procedure as described above for fibre diameter. Specifically, the estimated fibre diameter obtained from Eqn (3) was converted into estimated fibre cross-sectional area as follows: Estimated fibre cross sectional area μm 2 = π ×(Estimated f ibre diameter/2) 2 (4) Further correlation and regression analyses were performed to assess the relation between measured and estimated mean fibre cross-sectional area.
Statistical analysis. The Shapiro-Wilk test was used to check the normality of the distribution of data for all the variables analysed. In this respect, all of the variables analysed exhibited a normal distribution.
The relative frequency histograms of both individual fibre cross-sectional area values ( Fig. 2A) and individual motor unit conduction velocity estimates (Fig. 2B) were composed of 41 bins with a ∼500 μm 2 width and 44 bins with ∼0.07 m s −1 width.
An independent Student's t test was used to assess differences between average fibre type I and type II cross-sectional area (Fig. 2C) and between average low-threshold and high-threshold motor unit conduction velocity (Fig. 2D).
The relation between conduction velocity and normalized recruitment threshold (%MVF) of all the motor units identified for each participant was assessed with Pearson's product-moment correlation coefficient (r). The coefficient of determination (R 2 ) of the linear regressions was computed for each participant and averaged across participants as an index of prediction power ( Fig. 2E and F).
The relation between MUCV AVERAGE and muscle-electrode distance was assessed with Pearson's product-moment correlation coefficient (r), at the whole population level (Fig. 2G). Similarly, the relation between MUCV AVERAGE and muscle fibre type composition (i.e. % type II fibres by area) was assessed with Pearson's product-moment correlation coefficient (r).
Individual motor unit conduction velocity estimates of all identified motor units were converted into estimated fibre diameters by Eqn (2). Individual fibre cross-sectional area values were converted into measured fibre diameters by Eqn (1). The relative frequency histograms of both estimated and measured fibre diameters were composed of 170 bins with a ∼2 μm diameter increment (Fig. 3A). Similarly, the relative frequency histograms of both estimated and measured fibre cross-sectional area were composed of 81 bins with a ∼250 μm 2 area increment (Fig. 3D). The Kolmogorov-Smirnov test was adopted to compare the estimated and the measured fibre diameter/cross-sectional area distributions of data.
The relations between estimated mean fibre diameter/cross-sectional area and measured mean fibre diameter/cross-sectional area were assessed with Pearson's product-moment correlation coefficient (r) at the individual subject level. The coefficient of determination (R 2 ) of the linear regressions was computed as an index of prediction power ( Fig. 3B and E). The agreement between fibre diameter/cross-sectional area derived from conduction velocity estimates (estimated mean muscle fibre diameter/cross-sectional area) and fibre diameter/cross-sectional area derived from muscle biopsies (measured mean muscle fibre diameter/fibre cross-sectional area) was assessed with Bland-Altman analysis. The average absolute bias across subjects, i.e. the average of the absolute difference between measured . C, bar plots displaying the average type I and type II fibre CSA measured from muscle biopsies for each participant (n = 29). For each participant, average type II fibre CSA was significantly greater than that of type I fibres (mean ± SD; P = 0.015; independent Student's t test). D, bar plot displaying the mean MUCV of low-threshold (LT MU) and high-threshold (HT MU) motor units for each participant (n = 29). The mean MUCV of HT MU was significantly higher than that of LT MU (mean ± SD; P < 0.001; independent Student's t test). E and F, MUCV regression lines plotted as a function of normalized recruitment threshold (MU RT) for the 29 participants (in E, resistance-trained, n = 16; in F, untrained individuals, n = 14). Each dot represents a single MUCV estimate and each colour represents a different participant. MUCV linearly increased as a function of normalized MU RT in all participants (P < 0.001 in all cases; linear regression analysis). The average coefficient of determination (R 2 ) across participants is reported in the upper left corner. G, correlation analysis between average muscle-electrode distance and MUCV AVERAGE . Each green dot represents a different participant (n = 29). The coefficient of correlation (r) is displayed in the upper left corner (P = 0.498; Pearson's product-moment correlation analysis).
[Colour figure can be viewed at wileyonlinelibrary.com] J Physiol 601.10

Figure 3. Relations between estimated mean and measured mean muscle fibre diameter and cross-sectional area
A and D, relative frequency histograms of muscle fibre diameters (MFD) (A) and muscle fibre cross-sectional area (fibre CSA) (D), when considering all individual data from the 29 participants. In both graphs, the distributions of measured MFD (A) and measured fibre CSA (D), directly derived from muscle biopsies, are shown in red, whereas the distributions of estimated MFD (A) and estimated fibre CSA (D), derived from motor unit conduction velocity, are shown in blue. Measured and estimated MFD distributions A were composed of 170 bins with a ∼2 μm width; measured and estimated fibre CSA distributions (D) were composed of 81 bins with a ∼250 μm 2 width. Differences between estimated and measured MFD/fibre CSA distributions were assessed with the Kolmogorov-Smirnov Test. B and E, inter-subject relations between estimated mean MFD and measured mean MFD (n = 29) (B), and estimated mean fibre CSA and measured mean fibre CSA (n = 29) (E). Regression analyses confirmed that the estimated mean MFD can predict the majority of the variance of the measured mean MFD (R 2 = 0.65; P < 0.001; linear regression analysis), at the individual subject level (B). Similarly, the estimated mean fibre CSA can predict the majority of the variance of the measured mean fibre CSA (R 2 = 0.66; P < 0.001; linear regression analysis), at the individual subject level (E). Each green dot represents a different participant. C and F, Bland-Altman analyses showing absolute differences between estimated mean MFD and measured mean MFD (n = 29) (C), and between estimated mean fibre CSA and measured mean fibre CSA in absolute values (n = 29) (F). A mean absolute bias in estimated mean MFD compared with measured mean MFD of 1.9 ± 9.2 μm (C) and in estimated mean fibre CSA compared with measured mean fibre CSA of 256 ± 1307 μm 2 (F) was detected (mean ± SD; Bland-Altman analysis). LLOA, lower limit of agreement; ULOA, upper limit of agreement.
[Colour figure can be viewed at wileyonlinelibrary.com] and estimated fibre diameter/cross-sectional area for each participant, was determined between estimated mean fibre diameter/cross-sectional area and measured mean fibre diameter/cross-sectional area. Similarly, the average relative bias across subjects, i.e. the average of the relative difference between measured and estimated fibre diameter/cross-sectional area for each participant, was computed.
The absolute and relative errors for each individual were determined by calculating the difference between estimated mean fibre diameter/cross-sectional area and measured mean fibre diameter/cross-sectional area (normalized to measured mean cross-sectional area for relative error), before averaging across participants.
All statistical analyses were performed with the software GraphPad Prism, Ver. 9.2.0 (GraphPad Software, San Diego, CA, USA) and MATLAB, Ver. 2021 (The MathWorks, Natick, MA, USA). Statistical significance was set at α = 0.05 for all tests. Results are expressed as means ± SD.

Results
We report a systematic assessment of the relation between action potentials propagation velocity of individual motor units, i.e. motor unit conduction velocity, here estimated non-invasively from large and representative populations of motor units during isometric voluntary contractions, and the size of muscle fibres, i.e. fibre diameter and cross-sectional area, directly derived from muscle biopsies, in an heterogeneous population of healthy individuals.

High-density EMG decomposition and motor unit properties
Myoelectrical signals recorded during voluntary ramp contractions were decomposed into discharge instants of individual motor units   (Fig. 1D-F ;Holobar & Zazula, 2007). After visual inspection and manual editing (see 'High-density surface EMG decomposition' in Methods), when considering all the 29 participants and the four contraction levels performed by each participant, we identified a total of 701 unique motor units within the 29 biceps brachii short-head muscles. The average number of identified motor units per participant was 24.2 ± 8.6 (range 11-47).
The average normalized motor unit recruitment threshold, i.e. the force value relative to maximum (%MVF) at which the first motor unit action potential was discharged, was 31.5 ± 4.0%MVF (range 24-42.7). Based on their normalized recruitment thresholds, the identified motor units were clustered into low-threshold (range of recruitment thresholds: 0−30%MVF) and high-threshold motor units (range of recruitment thresholds: 50−70%MVF) . A total of 298 and 76 motor units were classified as low-threshold and high-threshold, respectively. The average number of identified low-threshold and high-threshold motor units per participant was 10.3 ± 4.3 (range 2-22) and 2.6 ± 2.2 (range 1-9), respectively.
For each participant, the average recruitment motor unit discharge rate, i.e. the mean of the instantaneous discharge rates of the first 20 action potentials in the ascending phase of the ramp contractions, was 17.4 ± 2.7 pulses s −1 (range 14.3-23.2).

Motor unit conduction velocity
The conduction velocity of the identified motor units was estimated with a validated multi-channel maximum likelihood algorithm (see 'Motor unit conduction velocity estimation' in Methods). Figure 1G shows a representative example of one motor unit with action potentials propagating from the innervation zone (in red) in both directions along an electrode column of the grid. Figure 2B shows the relative frequency histogram for conduction velocity values of all the identified motor units (n = 701; range 3.77-6.02 m s −1 ). In agreement with previous observations (Del Hogrel, 2003), the distribution was unimodal and did not highlight distinct clusters of conduction velocity values. The average conduction velocity (MUCV AVERAGE ) of each participant from all their active motor units in the ascending phase of the ramp contraction was 4.76 ± 0.33 m s −1 (range 4.34-5.48 m s −1 ).
High-threshold motor units displayed significantly higher conduction velocity compared to that of low-threshold motor units (5.20 ± 0.40 m s −1 vs. 4.53 ± 0.35 m s −1 , P < 0.001, Fig. 2D). As expected, the faster propagation velocity of action potentials in high-threshold compared to low-threshold motor units was further confirmed by the significant linear relation between conduction velocity and recruitment threshold, which was observed in each of the 29 participants. Individual participant R 2 values for the relationship between conduction velocity and recruitment threshold ranged from 0.59 to 0.93 (P < 0.001, in all participants) with an average of 0.78 ± 0.09. The relations between conduction velocity and recruitment threshold for the 29 participants are shown in Fig. 2E (resistance-trained, n = 16) and Fig. 2F (untrained individuals, n = 14).
In order to assess the influence of subcutaneous adipose tissue thickness (i.e. volume conductor) on conduction velocity estimates, we also assessed the inter-individual association (i.e. n = 29) between the mean muscle-electrode distance, measured with B-mode J Physiol 601.10 ultrasonography from the biceps brachii short head (see 'Ultrasound recording' in Methods), and MUCV AVERAGE . The absence of a significant correlation between the two parameters (P = 0.498, Fig. 2G) indicated that motor unit conduction velocity values estimated from high-density surface EMG decomposition were not influenced by the subcutaneous adipose tissue interposed between the muscle belly and the recording site. This rules out the possibility that adipose tissue would be a confounding factor for the results.

A non-invasive interface to estimate muscle fibre size
Single motor unit conduction velocity values were converted into estimated fibre diameters and subsequently related to the measured fibre diameters directly obtained from muscle biopsies. By relating the estimated fibre diameters with actual measured fibre diameters, we developed a non-invasive interface which allowed the prediction of the mean measured fibre diameter and of fibre cross-sectional area at the individual subject level from an EMG-derived neuromuscular parameter (see 'Estimation of muscle fibre diameters from motor unit conduction velocity' in Methods). Figure 3A shows the relative frequency histograms of measured and estimated fibre diameters, when pooling all analysed muscle fibres (n = 5619, in red) from all participants. Measured fibre diameters ranged from 37.9 to 162.8 μm, with an average value of 90.8 ± 19.6 μm. Similarly, when considering all identified motor units (n = 701, in blue), estimated fibre diameters ranged from 41.6 to 148.7 μm, with an average value of 90.8 ± 19.7 μm, following a calibration (see 'Estimation of muscle fibre size from motor unit conduction velocity' in Methods). The estimated fibre diameter distribution was very similar to that of the measured fibre diameter, which showed that the calibration procedure allowed almost perfect mapping of motor unit conduction velocity into fibre diameter values. The similarity between the measured and estimated fibre diameter distributions was further confirmed by the Kolmogorov-Smirnov test (D = 0.050, P = 0.086). It should be noted that, while the agreement of mean and standard deviation of the two distributions was expected because of the calibration, the non-significant difference between the entire distributions of motor unit conduction velocity and fibre diameter indicated similar shapes of the two histograms that could not be due only to the calibration (which was a linear transformation).
The distributions of motor unit conduction velocity and fibre diameter obtained from the entire dataset provide an indication of the relation between propagation velocity and fibre diameter in our dataset. Equation (3) in Methods provides the linear transformation to obtain conduction velocity from fibre diameter and vice versa, according to the estimates of conduction velocity obtained by the proposed technique. The next step was to test whether the relation between the two variables could only be obtained for the full population or whether it was robust enough to allow for the prediction of the distribution of fibre diameters given the measured distribution of motor unit conduction velocity for individual subjects. Basically, we wanted to study the predictability of fibre diameter distribution for each subject when the motor unit conduction velocity was measured. For this analysis, we obtained the linear relation between motor unit conduction velocity and fibre diameter for each subject by a leave-one-out procedure. According to this procedure, the optimal fitting parameters for Eqn (2) (see Methods) were obtained using the entire dataset but excluding the subject for whom the estimates were derived. Once the parameters were determined using data from all subjects but one, the fibre diameter distribution for the left-out subject was estimated using Eqn (2) with the optimal parameters obtained from the remaining dataset. This was repeated by leaving out all subjects, one at a time, so that estimated fibre diameter distributions were obtained for each subject without using the corresponding subject-specific biopsy data. This procedure corresponded to estimating fibre diameter distributions from a general association with motor unit conduction velocity obtained from a training/calibration set of data. After this operation, we compared the measured fibre diameter distributions from each subject with those estimated from the motor unit conduction velocity distributions to which the subject-specific linear transformation derived from all other subjects was applied. In this way, the calibration was obtained from data not used for the analysis of the association between estimated and measured fibre diameter distributions at the individual subject level. This approach provided information on the predictability of fibre diameter at the individual subject level, following a separate calibration on a representative sample of individuals. With this procedure, Pearson's correlation analysis revealed a strong and significant correlation between measured and estimated mean fibre diameter across individuals (r = 0.81, 95% CI = 0.62 0.90, P < 0.001). Regression analysis further confirmed that estimated mean fibre diameter was able to predict the majority of the variance of the measured mean fibre diameter (R 2 = 0.65, P < 0.001, Fig. 3B). This indicated that the fibre diameter could be estimated using only the motor unit conduction velocity data with a calibration on a separate dataset. For new subjects, it would be therefore possible to estimate fibre diameter distributions without needing the biopsy data. According to the high level of predictability obtained on an individual subject basis, the parameters estimated for each subject with the leave-one-out procedure varied in a very small range when changing the left-out subject.
Bland-Altman plots are shown in Fig. 3C for absolute values and revealed a mean absolute bias in estimated mean fibre diameter compared with measured mean fibre diameter of 1.9 ± 9.2 μm across all subjects. The mean relative bias across subjects was 2.2 ± 10.4%. Moreover, the mean absolute and relative errors in estimated mean fibre diameter compared with measured mean fibre diameter across subjects were 7.3 ± 5.8 μm and 8.1 ± 6.8%, respectively.
Lastly, we also explored the relation between muscle fibre cross-sectional area estimated from motor unit conduction velocity and the muscle fibre cross-sectional area directly measured from muscle biopsies. Firstly, estimated fibre diameter data from all participants were converted into estimated fibre cross-sectional area (see Eqn (4) in Methods). Figure 3D shows the relative frequency histograms of measured and estimated fibre cross-sectional area, when pooling individual data from all participants (n = 701). Similarly to what was observed for fibre diameter, the estimated fibre cross-sectional area distribution (6782 ± 2829 μm 2 ) was very similar to that of the measured fibre cross-sectional area (6772 ± 2912 μm 2 ), as shown by the almost identical means and standard deviations. The similarity between the measured and estimated cross-sectional area distributions was further confirmed by the Kolmogorov-Smirnov test (D = 0.049, P = 0.090). As for fibre diameter, this was expected because of the global calibration. Moreover, at the individual subject level, the regression analysis revealed a strong and significant association between measured and estimated mean fibre cross-sectional area (R 2 = 0.66, P < 0.001, Fig. 3E) (this was also obtained with the leave-one-out procedure, as discussed above).
In contrast to the findings for fibre size (fibre diameter and cross-sectional area) and motor unit conduction velocity, there was no association between conduction velocity and muscle fibre type composition (i.e. percentage of type II fibres by area; r = 0.20, P = 0.302), as well as between fibre cross-sectional area and muscle fibre type composition (i.e. percentage of type II fibres by area; r = 0.31, P = 0.105).
Bland-Altman plots are shown in Fig. 3F for absolute values and indicate a mean absolute bias in estimated mean fibre cross-sectional area compared with measured mean fibre cross-sectional area of 256 ± 1307 μm 2 across subjects. The mean relative bias across subjects was 4.0 ± 20.3%. Moreover, the mean absolute and relative errors in estimated mean fibre cross-sectional area compared with measured mean fibre cross-sectional area across subjects were 1037 ± 814 μm 2 and 15.9 ± 14.4%, respectively.
These results suggest a good degree of accuracy in predicting important indicators of muscle fibre size, i.e. the mean fibre diameter and mean fibre cross-sectional area, with our non-invasive neuromuscular interface.

Discussion
We demonstrated that a non-invasive neuromuscular interface allows the prediction of the mean muscle fibre size for a heterogeneous population of healthy individuals with a relatively low margin of individual error (≤16%) and very low systematic bias (≤4%). The strong and positive relations observed between the mean fibre diameter/cross-sectional area estimated from motor unit conduction velocity and the mean fibre diameter/cross-sectional area directly derived from muscle biopsies (R 2 = 0.65 and R 2 = 0.66, respectively), show experimentally that motor unit conduction velocity may be adopted as a novel strong index of muscle fibre size in healthy individuals. Our results provide substantial novelty about the potential adoption of an EMG-derived, and fully non-invasive, physiological parameter to estimate muscle fibre size with a relatively low margin of error.
The main finding of the current analysis is that we revealed the possibility of transforming motor unit conduction velocity data into estimated measures of muscle fibre size (i.e. estimated fibre diameters/cross-sectional area), which in turn showed a good degree of association with the actual muscle fibre size measured directly by muscle biopsy. Regression analysis demonstrated the possibility to predict the mean measured fibre diameter (Fig. 3B) from an EMG-derived parameter across subjects, with a relatively low bias and error. Since fibre cross-sectional area is generally considered as a more precise indicator of fibre size, because of the typical irregular shape of skeletal muscle fibres (Mizuno, 1991), we also assessed whether the proposed equation allowed an accurate prediction of measured fibre cross-sectional area. Similarly, our results confirmed the possibility of predicting the mean measured fibre cross-sectional area (Fig. 3E) from an EMG-derived parameter across subjects, with a relatively low bias and error.
From a methodological standpoint, the proposed equation which describes the linear relation between motor unit conduction velocity and fibre diameter (Eqn (3)) presents normalization values (constants) that are somewhat different from those observed by previous studies on both healthy (Nandedkar & Stålberg, 1983) and pathological (Blijham et al., 2006) individuals. These discrepancies may be due to several methodological differences. First, in these previous studies the conduction velocity has been estimated with intramuscular EMG (e.g. concentric needle electrode), during electrically elicited stimulations in resting conditions (Blijham et al., 2006), J Physiol 601.10 or it has been mathematically derived from simulations of intramuscular recordings (Nandedkar & Stålberg, 1983). In the first case, the high selectivity associated with intramuscular EMG recordings combined with electrical stimulations may have hindered the sampling of motor units in a representative way, a condition which may have affected the apparent relations between the conduction velocity and fibre diameter (Blijham et al., 2006). Moreover, the estimates of conduction velocity from intramuscular EMG are generally lower, in the order of ∼1.0 m s −1 , and more variable than the estimates obtained from the surface EMG (Zwarts, 1989). In the second case, the adoption of simulated intramuscular EMG recordings for conduction velocity estimation did not allow a systematic assessment of the relation with fibre diameter based on intrasubject observations of both measurements from the same muscle (Nandedkar & Stålberg, 1983), as was done in the current analysis. It is worth considering also that conduction velocity estimates of earlier studies (Blijham et al., 2006;Methenitis et al., 2016;Nandedkar & Stålberg, 1983) were likely less accurate than those obtained in the present study, not only due to the different recording technique adopted, but also to the different algorithms applied (Farina et al., 2001(Farina et al., , 2002. Individual motor unit conduction velocity values lay within the physiological ranges in all cases (Del Vecchio, Negro, Falla et al., 2018;Farina et al., 2000;Hogrel, 2003) and displayed a unimodal distribution (Fig. 2B). Our results are in agreement with previous studies on biceps brachii (Hogrel, 2003;Hogrel et al., 2008) and other muscles , which did not observe distinct classes of fibres when estimating their motor unit conduction velocity. Similarly, the distribution of fibre cross-sectional area derived from muscle biopsies was continuously distributed, despite fibre typing ( Fig. 2A). As suggested by , the continuous distribution of motor unit properties (i.e. conduction velocity, fibre cross-sectional area) may be a prerequisite for a smooth generation of force in the full recruitment range and for energy control.
In agreement with previous reports on biceps brachii and other muscles (Andreassen & Arendt-Nielsen, 1987;Casolo, Farina et al., 2020;Del Vecchio, Negro, Falla et al., 2018;Nuccio et al., 2020;Sadoyama et al., 1988), we found a strong association between motor unit conduction velocity and recruitment threshold (Fig. 2E and 2F). This indicated the progressive recruitment of motor units with increasing propagation velocities, and further confirmed that motor unit conduction velocity can be considered as a robust size principle parameter (Andreassen & Arendt-Nielsen, 1987;Del Vecchio et al., 2017;. Another important factor to consider that relates to the methodology adopted to estimate motor unit conduction velocity and to derive fibre diameters in the current study is the observed independence of conduction velocity from the subcutaneous adipose tissue interposed between the muscle belly and the recording site (Fig. 2G). Our results contrast with those of previous studies (Cescon et al., 2008;Nordander et al., 2003), which reported that muscle fibre conduction velocity estimated from the interference EMG was influenced by the subcutaneous fat layer. Indeed, in the current study we found that the conduction velocity of individual motor units was not influenced by the subcutaneous fat layer. Thus, the proposed methodology has the potential to reliably estimate fibre diameter in populations with very different distributions of both muscle fibre size and subcutaneous adipose tissue.
In this study we found motor unit conduction velocity to be associated with muscle fibre size, whilst being unrelated to fibre type composition. These findings indicate that fibre size, independent from fibre type composition, determines motor unit conduction velocity. A possible explanation that could justify the existence of a biophysical relation between conduction velocity and muscle fibre size is the lower cytoplasmic resistance of larger muscle fibres, which results in faster propagation of a motor unit action potential along their sarcolemma compared to smaller muscle fibres. The lower cytoplasmic resistance, along with the specific electrophysiological features of larger muscle fibres (Clausen, 1986(Clausen, , 2003Zhang et al., 2006), which directly influence the excitability of the sarcolemma, seems to justify the faster propagation of action potentials along the sarcolemma of bigger muscle fibres (i.e. fibres with greater diameter or cross-sectional area), generally innervated by larger α-motor neurons and hence generally belonging to higher-threshold motor units, compared to that of smaller muscle fibres, generally innervated by smaller α-motor neurons and hence generally belonging to lower-threshold motor units.
There are some limitations within the current study that should be recognized. Firstly, we established a relation between fibre diameter and motor unit conduction velocity, an EMG decomposition-derived parameter which reflects the average propagation velocity of action potentials along muscle fibres innervated by individual motor neurons. Consequently, this relation could not be assessed on individual muscle fibres because it is technically impossible to estimate conduction velocity and diameter from one and the same fibre in vivo. Furthermore, fibre diameters were indirectly derived from fibre cross-sectional area measurements by assuming that the muscle fibre is straight and circular. However, this appears to be the case in fetal and infant skeletal muscle, whereas as the fibres become larger with increasing age, the fibre areas are known to become less rounded and somewhat irregular in cross-section (Saltin & Gollnick, 1983). Moreover, although the adoption of high-density grids of electrodes ensured a representative spatial sampling of motor unit behaviour and properties that allowed the overcoming of the limitations of intramuscular EMG recordings , we have to consider that the biopsy-derived fibre diameters were likely affected by sampling errors, since only a small proportion of muscle fibres were assessed. Therefore, both fibre diameter and conduction velocity values might have been affected by a sampling error, whose magnitude is difficult to quantify. It has also to be noted that in order to assess the validity and generalizability of the proposed interface to estimate muscle fibre size, future studies assessing the relations between the estimated and measured fibre size indicators on an independent sample are needed. Closely related to this point, is also the fact that the proposed methodology was developed on data collected in healthy adult males (either trained or untrained) presenting a wide variability in the distribution of data. As a consequence, despite our belief in the potential widespread applicability of this novel neuromuscular interface for the non-invasive estimation of muscle fibre size in vivo, future investigations on a wider sample of specific populations (e.g. children, elderly, neuromuscular patients), including both males and females but displaying more homogeneous characteristics (more skewed data), are needed in order to test its applicability and validity. Similarly, the applicability of the proposed methodology to detect changes in muscle fibre size following specific longitudinal interventions should be directly tested in future investigations.
In addition, in the current study motor unit conduction velocity was estimated from a fusiform (non-pennate) muscle, i.e. biceps brachii, whose specific fibre orientation tends to be parallel with the longitudinal axis of the muscle. Therefore, in this muscle, the alignment of the electrode grid with muscle fibre orientation, a prerequisite for reliable conduction velocity estimation, is relatively straightforward. Nevertheless, more caution is needed when motor unit conduction velocity has to be estimated in pennate muscles (i.e. vastus medialis and lateralis or tibialis anterior) Del Vecchio et al., 2017;Nuccio et al., 2020).
Lastly, it seems legitimate to analyse the strengths and weaknesses of the proposed electrophysiological method to estimate muscle fibre size with respect to existing muscle imaging techniques. The major difference between available non-invasive muscle imaging techniques (e.g. conventional magnetic resonance imaging, computed tomography scanning, ultrasound imaging) and the proposed electrophysiological method is that these facilitate the measurement of whole muscle size at the macroscopic level (e.g. volume, area and thickness, respectively) . Conversely, the proposed neuromuscular interface offers the possibility of estimating muscle fibre (and motor unit) size at the microscopic level, which at present is only possible through immunohistochemistry from muscle biopsy samples. A major advantage of the proposed methodology is the fact that it allows measurement of the characteristics of the voluntarily activated motor units, and their concurrent neural (e.g. recruitment threshold, discharge rate) and morphological (e.g. muscle fibre size) properties in a unique manner. In this respect, novel advances in magnetic resonance imaging systems and approaches (e.g. diffusion-weighted imaging, diffusion tensor imaging, q-space imaging), which allow the obtaining of high-resolution images, have shown the potential to access and quantify accurately the muscle tissue structure at the microscopic level (i.e. at the fibre level; Hata et al., 2019;Kalia et al., 2017). However, their application is still confined to a few research laboratories and clinical settings, since it is limited by operation costs and availability, and requires highly specialized technicians to acquire and analyse images, in addition to the fact that there is no standardized protocol for image analysis and tissue segmentation. Whilst in the future advanced non-invasive imaging techniques may become common and replace the traditional more invasive means of quantifying muscle fibre size, we believe that at present the proposed neuromuscular interface offers a valid cost-and time-efficient alternative to traditional available methods for muscle fibre size assessment in vivo.
Overall, our results demonstrate the potential of using an EMG decomposition-derived physiological parameter, i.e. motor unit conduction velocity, to reliably estimate muscle fibre size in healthy adults in a fully non-invasive way. Here we confirmed that motor unit conduction velocity can be considered not only as a reliable and valid indicator of the progressive recruitment of motor units, i.e. a size-principle parameter, but also as a novel strong biomarker of muscle fibre size. We believe that the proposed neuromuscular interface may therefore open new perspectives in the use of a fully non-invasive means for muscle fibre size estimation not only in healthy populations, but also in other fields where the non-invasive and painless determination of muscle fibre diameter or cross-sectional area becomes a priority.