


ORIGINAL ARTICLE 

Year : 2015  Volume
: 5
 Issue : 3  Page : 177180 


Estimation of a single motor unit's threshold and activation range, a study on patients with muscular disorders
Nammam Ali Azadi, Daem Roshani
Kurdistan Research Center for Social Determinants of Health; Department of Epidemiology and Biostatistics, Medical School, Kurdistan University of Medical Sciences, Sanandaj, Islamic Republic of Iran
Date of Submission  03Dec2014 
Date of Acceptance  23May2015 
Date of Web Publication  15Sep2015 
Correspondence Address: Daem Roshani Department of Epidemiology and Biostatistics, Medical School, Kurdistan University of Medical Sciences, Pasdaran Street, 66177.13446 Sanandaj Islamic Republic of Iran
Source of Support: This work was supported financially by Kurdistan University of Medical Sciences, Sanandaj, Iran., Conflict of Interest: None declared.  Check 
DOI: 10.4103/2229516X.165380
Abstract   
Background: In clinical neurophysiology, threshold tracking studies are used to evaluate the functionality of a muscle through studying the functionality of its motor units (MUs) that govern the muscle. The functionality of an MU can be quantified by estimation of its excitability properties via MU's stimulusresponse curve. In this study, we aim to develop a modelbased approach to estimate MU's threshold mean and its activation range as indications of MU's excitability. This is a different approach from routine strategies in neurophysiology, which are mostly subjective. Methods: To assess the excitability of a single MU, needle electromyography examination was used to obtain the axonal activity of that MU. To improve estimation, the examination was repeated several times on individuals. Replication of experiment introduces serial correlation between observations. We account for this correlation by using a mixedeffects model. We investigate the appropriateness of classical logistic mixedeffects model and its Bayesian formulation for estimation purpose. Results: Both classical and Bayesian models can obtain a reliable estimation of MU's threshold. However, we found Bayesian approach to provide a better estimate of MU's activation range. Moreover, if data contain outliers both classical and Bayesian methods are vulnerable to some extent. Conclusions: Compared to the classical approach, Bayesian method is more flexible in dealing with overdispersion and provides more robust estimation of MU's parameters.
Keywords: Electromyography, logistic mixedeffects model, motor unit, muscular disorder, needle electromyography, neurophysiology
How to cite this article: Azadi NA, Roshani D. Estimation of a single motor unit's threshold and activation range, a study on patients with muscular disorders. Int J App Basic Med Res 2015;5:17780 
How to cite this URL: Azadi NA, Roshani D. Estimation of a single motor unit's threshold and activation range, a study on patients with muscular disorders. Int J App Basic Med Res [serial online] 2015 [cited 2021 Dec 7];5:17780. Available from: https://www.ijabmr.org/text.asp?2015/5/3/177/165380 
Introduction   
Needle electromyography (NEMG) examination is a routine clinical technique to examine the functionality of muscles. The contraction of muscles releases a microvolt level electrical signals that can be captured through a needle electrode inserted into the target muscle of the subject's tissue. By quantification of these bioelectrical signals, NEMG technique enables clinicians to study the nerve and muscular disorders.
A key concept in clinical EMG is the socalled motor units (MUs). A MU is formed by an alpha motor neuron in the spinal cord, its long protruding axon, and the tens to thousands of muscle fibers that this neuron innervates.^{[1],[2]} In patients with some extend of muscular disability, the muscle activity is limited due to nonfunctionality of a group of MUs. In testing muscle activities, the aim is to gather information from MUs' activity using conventional methods to understand the extent of loss (for future therapies and evaluation of disease progression). Of particular interest is to obtain a full description of MU's functionality in this evaluation.^{[3]}
In the laboratory, in the event of muscular disorder, the interest is to isolate a single MU to study its excitability. The excitability of a MU is reflected in the intensity of the electrical stimulus that is necessary to evoke an action potential (electrical signal) in MU's axon. Therefore, in the socalled threshold tracking studies, a full description of a MU is obtained by observing MU's responses to a range of stimulus intensities over which it displays a stochastic behavior (MU activation range).^{[4]} The activation range is spanned around the unknown MU's threshold so that it starts with a low stimulus intensity that is not large enough in the magnitude to innervate the nerve that supplies the MU, and ends with a stimulus intensity at which the nerve is always activated. MU's threshold is usually defined as the stimulus intensity that elicits a response to 50% of the stimuli.^{[5]} Generally, a MU fires (response) if the applied stimulus exceeds its threshold, otherwise it remains off (no response).^{[6]} In this sense, the MU responses are binary observations.
As an indication of the MU's excitability behavior, a threshold tracking study leads to MU's stimulusresponse curve (a sigmoid curve). [Figure 1] displays an example of stimulusresponse curve obtained from one NEMG experimentation. The horizontal scale represents the stimulus intensity values, varied from 59.9 mA to 69 mA. For each stimulus intensity, the vertical scale represents the number of successful MU's responses (firing) out of 50 attempts.  Figure 1: A typical stimulusresponse curve showing the general excitability behavior of a single motor unit. The motor unit (MU) starts with no response at lower stimulus intensities and ends with always active at high stimulus intensities. The mean threshold is located somewhere close to the middle of xaxis where MU response to an applied stimulus is totally unpredictable. Activation range is the range of stimulus intensities on xaxis
Click here to view 
The sshape curve in [Figure 1] is known as MU's stimulusresponse curve. In practice, MU's threshold (mean of the curve) and its activation range (slope of the curve) are calculated from the stimulusresponse curve and interpreted as MU's excitability properties. Thus, one important question is how one can obtain efficient estimates of these parameters. For estimation of MU's threshold, various strategies are available.^{[5]} However, the judgments are mainly subjective, and statistical methods are seldom used.^{[6]}
In this paper, our goal is to investigate how statistical methods can be merely employed for estimation purpose of both MU's threshold and activation range parameters. We use a real NEMG data to apply our proposed models. NEMG data are noisy thus a great care is needed to model the heterogeneity of data properly. This issue is discussed in Section 3 and our strategy to tackle this problem is introduced. Our models include classical approaches, which are readily available from standard statistical packages, and a more sophisticated approach that it's implementing requires experience and deeper statistical skills (Bayesian models).
The paper is organized as follows: NEMG experiment is introduced in Section 2. At Section 3, our strategy for estimating MU's threshold and activation range parameters is introduced. Results are given in Section 4. The paper ends with a discussion in Section 5.
A brief introduction to needle electromyography experiment
The NEMG is a routine technique in diagnostics of muscular disorders. In this technique, experimenter inserts a needle electrode into the target muscle of the subject's tissue to acquire microwavelet electrical signals. In the NEMG experiment we used in this study, the needle was manipulated to isolate a single MU. A potential was considered acceptable when it demonstrated marked "allornone" behavior and could be clearly separated out from other potentials throughout the recording session. When visual inspection indicated that such a stable potential was obtained, stimulus intensity was increased until no probabilistic inactivity occurred during 50 stimuli. From that level, which marks the end of the MU's activation curve, stimulus value was decreased in the smallest steps allowed by the equipment until the muscle fiber responded to none of 50 stimuli. At each step, the number of MU's responses to the 50 stimuli was recorded. The experiment was repeated several times at couple of stimulus voltages.
[Table 1] shows an example dataset of this experiment with six replications on a subject suffering from amyotrophic lateral sclerosis. At each replication, the MU's responses have been observed at a few distinct stimulus intensities (s). At each stimulus intensity, the observation (y), is the number of MU's firings out of 50 stimuli.  Table 1: The results of parameter estimations under Bayesian and nonBayesian (classical) logit mixedeffects model using data from [Table 2]
Click here to view 
In a typical NEMG experiment such as [Table 1], MU's responses can be treated as binomially distributed observations.
Procedure
Motor unit's threshold and activation range parameters
In this section, we develop our statistical model to estimate MU's parameters. Since MU's observations were binary identities, NEMG observations modeled under a logistic regression function. However, there was one major issue with experimentation to tackle; NEMG experiment was carried out over a long period of time (several minutes) where during this period the needle inserted to evoke MU moves slightly altering the stimulus voltage and the firing rate. This introduces a new source of variation into the data, the socalled overdispersion or extrabinomial variation. Because of the problem of overdispersion, to improve estimates experiment was repeated on an individual. In statistics overdispersion on binomially distributed data can be handled using a logistic mixedeffects model.^{[7]}
To introduce the notation, let's denote y_{it} as the number of MU's responses for a given stimulus intensity s_{it} at replication i = 1, 2, …r and time t = 1, 2, …n_{i} so that y_{it}❘s_{it} ~ Bernoulli(p_{it}). Consider the following model;
[Inline: 1]
Where p_{it} is the probability the MU to fire at time t and replication i, m_{i} and b_{i} are threshold and slope (activation range) parameters of the MU in replication i respectively. Model 1 assumes that both threshold and slope parameters are subject to change when we move between replicated NEMG experiments. That is, b_{i} and m_{i} are random effects distributed normally with mean 0 and variance [Inside: 1] and [Inside: 2] respectively. The error terms ɛ_{it} ˜ N(0, σ^{2}) are added to the model to account for overdispersion.
The main criticism of the classical model 1, however, is that it underestimates the uncertainty around the estimates. Bayesian approach solves this inadequacy by the formulation of the beliefs available prior to experimentation. This is the main motivation to consider Bayesian methods as an alternative choice with the hope to improve parameter estimations. For the sake of simplicity, we skip introducing mathematical Bayesian formulation for model 1. Readers can refer to ^{[8]} for details of maths and other justifications of Bayesian approach.
Results   
Weaklyinformative priors were used for Bayesian models. [Table 1] displays the results of parameter estimations of two models. Note that [Inside: 3] and [Inside: 4] are indications of variability between MU's threshold and slope parameters at six replicates and σ^{2} can be used as an indication of overdispersion. The R statistical package ^{[9]} was used to fit model 1 under classical model and WinBUGS software ^{[10]} for Bayesian paradigm. The result of MU's threshold and its activation range parameters for both models are given by [Table 1].
Note that classical and Bayesian estimates are similar but have different interpretations. For instance, under model 1, the value 67.54 is the maximum likelihood estimation of the threshold and 0.38 is its standard deviation whereas under model 2, [Inside: 5] is the posterior mean of the threshold parameter distribution and 1.08 is the posterior standard errors.
Out of six NEMG replications presented in [Table 2], for unknown reasons replication 6 was not completed. To address the effect of this unexpected disturbance on parameter estimations by two statistical methods, this replication was removed from the analysis. The results are given by [Table 3].  Table 2: A MU activity recorded by needle EMG experiment for an individual at six replications
Click here to view 
 Table 3: The results of parameter estimations under classical and Bayesian logit mixedeffects models after excluding replication 6 from data in [Table 2]
Click here to view 
Discussion and Conclusion   
From the results presented in [Table 1], the immediate conclusion is the similarity between classical logit mixed model and its Bayesian estimations. However, a mild difference can be observed in estimation of the activation range and variance components between two models. The larger posterior estimation (or maximum likelihood estimator) of variance between thresholds ([Inside: 6]) can be due to the last replication with threshold mean estimation of 71.52, markedly different from others. Activation range of this replication (and hence its threshold) is significantly different from others [Table 2]. Note that NEMG experiment for this replication is not complete and has been terminated for unknown reason/s. In general, under Bayesian approach, the estimation of variance components ([Inside: 7] and [Inside: 8]) is larger than that of classical method but estimation of [Inside: 9], an indication of overdispersion, is smaller. This means Bayesian model can handle overdispersion in the NEMG data more efficiently. However, both approaches suffer to some extend from overdispersion ([Inside: 10] means no overdispersion).
As expected, after discarding replication 6 [Table 3], the variation among thresholds was decreased significantly for both models. Overdispersion was improved under Bayesian model but no improvement achieved using classical model. Furthermore, note that for classical logit model, [Inside: 11]. As Punheiro and Bates ^{[11]} noted, this variance should not be interpreted as no variation between activation range in different replications, it might be an indication of small variation between them that is not large enough to be captured by classical model. In other words, the variations have to be large in order to be estimated by the classical model. Otherwise it fails to estimate this variation appropriately.
Generally speaking, using a real NEMG experiment, the main conclusion is that both approaches can achieve a robust estimation of threshold parameter, but they dealt with estimation of slope parameters (activation range) in different manners. Bayesian approach appeared to provide better estimate of slope parameter, mainly due to better handling of overdispersion. However, when outlier presented in data (such as replication 6), both classical and Bayesian method were vulnerable though latter model showed better performance. It should be noted to avoid exaggerating our finding in better performance of Bayesian approach for overdispersed data, a comprehensive simulation study is recommended to demonstrate the reliability of both approaches in a broad spectrum of data. This can be a topic for future research.
Acknowledgments   
Authors would like to thank Joleen Blok at Erasmus MC, Netherlands, for making needle EMG data available to us.
References   
1.  Ritchie J. The Axon: Structure, Function and Pathophysiology. Oxford: Oxford University Press; 1995. 
2.  Ignacio RC, Luis GU, Armando MT. Potential duration: Measurement and significance. In: Ajeena IM, editor. Advances in Clinical Neurophysiology. Rijeka, Croatia: InTech Online Publishers; 2012. p. 13360. 
3.  Abdelmaseeh M, Poupart P, Smith B, Stashuk DW. Muscle categorization using quantitative needle electromyography: A 2stage Gaussian mixture model based approach. In ICMLA, '12 Proceedings of the 2012 11th International Conference on Machine Learning and Applications', Vol 1, IEEE. Computer Society Washington, DC, USA; 2012. p. 54853. 
4.  Bostock H, Cikurel K, Burke D. Threshold tracking techniques in the study of human peripheral nerve. Muscle Nerve 1998;21:13758. 
5.  Stein RB, Yang JF. Methods for estimating the number of motor units in human muscles. Ann Neurol 1990;28:48795. 
6.  Ridall PG, Pettitt AN, Henderson RD, McCombe PA. Motor unit number estimation – A Bayesian approach. Biometrics 2006;62:123550. 
7.  Zhang H, Lu N, Feng C, Thurston SW, Xia Y, Zhu L, et al. On fitting generalized linear mixedeffects models for binary responses using different statistical packages. Stat Med 2011;30:256272. 
8.  Azadi NA. Modeling the Excitability Behaviour of a Single Motor Unit and Experimental Design for the Surface EMG Test. Ph.D. Thesis. Lancaster, UK: Mathematics and Statistics Department, Lancaster University; 2011. 
9.  R Development Core Team. R: A Language and Environment for Statistical Computing, Reference Index version 2.15.1. R Foundation for Statistical Computing, Vienna, Austria; 2012. Available from: http://www.Rproject.org. [Last accessed on 2014 Nov 25]. 
10.  Lunn DJ, Thomas A, Best N, Spiegelhalter D. Winbugs – A bayesian modelling framework: Concepts, structure, and extensibility. Stat Comput 2000;10:32537. 
11.  Punheiro JC, Bates DM. Mixedeffects Models in S and SPlus. New York, USA: Springer; 2000. 
[Figure 1]
[Table 1], [Table 2], [Table 3]
