Evaluating Inhibition of Motoneuron Firing From Electromyogram Data to Assess Vestibular Output Using Vestibular Evoked Myogenic Potentials S. R. Prakash,1,2,3 Barbara S. Herrmann,1,2,4 Rupprecht Milojcic,2,3 Steven D. Rauch,2,3 and John J. Guinan, Jr.1,2,5
Objectives: Vestibular evoked myogenic potentials (VEMPs) are due to vestibular responses producing brief inhibitions of muscle contractions that are detectable in electromyographic (EMG) responses. VEMP amplitudes are traditionally measured by the peak to peak amplitude of the averaged EMG response (VEMPpp) or by a normalized VEMPpp (nVEMPpp). However, a brief EMG inhibition does not satisfy the statistical assumptions for the average to be the optimal processing strategy. Here, it is postulated that the inhibition depth of motoneuron firing is the desired metric for showing the influence of the vestibular system on the muscle system. The authors present a metric called “VEMPid” that estimates this inhibition depth from the EMG data obtained in a usual VEMP data acquisition. The goal of this article was to compare how well VEMPid, VEMPpp, and nVEMPpp track inhibition depth.
(EMG) potential of the corresponding muscle. The inhibitory strength of such a volley has been quantified by the change in the average EMG, taking into account that the EMG response is proportional to the level of muscle contraction (Widmer & Lund 1989). Although averaging is the optimal signal processing technique when measuring a fixed signal in an independent random noise, these conditions are not satisfied for measurements of averages of briefly inhibited EMGs. The brief inhibition of an ongoing EMG corresponds better to a signal that modulates a noise than to a signal that is added to an independent noise. Thus, a signal processing method that is different than simple averaging may produce better results for quantifying such an inhibition. Here, we present a signal processing method for quantifying a brief inhibition in the context of a clinically used vestibular test, but it can be used to quantify a brief inhibition or a subthreshold excitation in ongoing responses of any motoneuron pool. Activation of the vestibular system by brief sounds or by vibrational stimuli produces a brief inhibition of various muscles. The inhibitory signal from vestibular system cannot be measured directly, but its effect on the EMG of a target muscle can be measured, and this effect is used as a clinical index of vestibular sensitivity. In particular, the vestibular evoked myogenic potential (VEMP), a surface EMG potential, is used to assess the functioning of the saccule (Welgampola & Colebatch 2005; Curthoys 2010; Rosengren et al. 2010). Saccule excitation activates vestibular reflex pathways that inhibit motoneurons that innervate the ipsilateral sternocleidomastoid (SCM) muscle (Colebatch et al. 1994; Colebatch 2001; Murofushi and Kaga 2009). Although VEMP responses can be measured from many head and neck muscles, in this study we concentrate on the response to sound measured from the ipsilateral SCM, which is sometimes called the cervical VEMP (cVEMP) response. Previous studies on VEMP potentials focused on the average stimulus-locked VEMP waveform, particularly its peak to peak (PP) amplitude, which we refer to as “VEMPpp.” Obtaining a VEMP metric by averaging the rectified EMG signal has also been tried, but it is not as good as simple averaging (Rosengren et al. 2010). VEMPpp amplitudes scale approximately linearly with SCM contraction strength although saturation of VEMPpp amplitude has been reported for some types of SCM contractions (Bogle et al. 2013). To minimize the effects of different contraction strengths, VEMPs are commonly “normalized” by dividing the VEMPpp by a measure of the ongoing EMG response, for example, by the average amplitude of the rectified EMG. This yields a normalized VEMPpp or an “nVEMPpp.” VEMP responses arise because the saccule responds to sound and through neural pathways produces inhibition of SCM motoneurons (Colebatch & Rothwell 2004). Measures such as
Design: To find a robust method to compare VEMPid, VEMPpp, and nVEMPpp, realistic physiological models for the inhibition of VEMP EMG signals were made using VEMP data from four measurement sessions on each of the five normal subjects. Each of the resulting 20 EMGproduction models was adjusted to match the EMG autocorrelation of an individual subject and session. Simulated VEMP traces produced by these models were used to compare how well VEMPid, VEMPpp, and nVEMPpp tracked model inhibition depth. Results: Applied to simulated and real VEMP data, VEMPid showed good test-retest consistency and greater sensitivity at low stimulus levels than VEMPpp or nVEMPpp. For large-amplitude responses, nVEMPpp and VEMPid were equivalent in their consistency across subjects and sessions, but for low-amplitude responses, VEMPid was superior. Unnormalized VEMPpp was always worse than nVEMPpp or VEMPid. Conclusions: VEMPid provides a more reliable measurement of vestibular function at low sound levels than the traditional nVEMPpp, without requiring a change in how VEMP tests are performed. The calculation method for VEMPid should be applicable whenever an ongoing muscle contraction is briefly inhibited by an external stimulus. Key words: cVEMP, Vestibulocollic reflex.
Electromyogram,
EMG,
oVEMP,
VEMP,
(Ear & Hearing 2015;36;591–604)
INTRODUCTION A brief inhibitory volley delivered to a motoneuron pool can be detected by the change it produces in the electromyographic 1 Harvard-MIT Division of Health Science and Technology, Speech and Hearing Bioscience and Technology Program, Massachusetts Institute of Technology, Cambridge, Massachusetts, USA; 2Department of Otology and Laryngology, Harvard Medical School, Boston, Massachusetts, USA; 3 Department of Otolaryngology, Massachusetts Eye and Ear Infirmary, Boston, Massachusetts, USA; 4Department of Audiology, Massachusetts Eye and Ear Infirmary, Boston, Massachusetts, USA; and 5Eaton-Peabody Laboratories, Massachusetts Eye and Ear Infirmary, Boston, Massachusetts, USA.
0196/0202/2015/365-0591/0 • Ear & Hearing • Copyright © 2015 Wolters Kluwer Health, Inc. All rights reserved • Printed in the U.S.A. 591
592
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
the VEMPpp provide only an indirect measure of this inhibition. A metric that would be a step closer to the saccule’s response to the stimulus is the decrement in the motor-unit spike rate. With this in mind, we propose the “Inhibition Depth,” defined as the percentage reduction in the composite spike rate of SCM motoneurons produced by the activation of the saccule, as an index of motoneuron inhibition (Prakash 2009). The modulation of motoneuron firing rate was also identified by Lütkenhöner and Basel (2012) as being the key variable to quantify the saccular effect on the motoneuron pool. For clinical use, we developed a single-number index for VEMP inhibition that is termed the VEMP inhibition depth or “VEMPid” (Prakash 2009). The goal for this index was that it would represent the inhibition depth and that for a given stimulus, subjects with similar saccular sensitivity would have similar VEMPids. This means that the VEMPid metric should be constant even though the VEMP average waveforms may be different due to different muscle contractions, different placements of the VEMP electrodes, or other nonvestibular factors. If this robustness is achieved, VEMPid would be a clinically significant measure of saccule function and a useful indicator of pathological changes in the saccule. The VEMPid is intended as an alternative to VEMPpp and nVEMPpp. VEMPid can be calculated from the EMG data acquired using the already well-established experimental protocols for recording VEMP responses (Curthoys 2010; Rosengren et al. 2010), that is, VEMPid is derived from the set of individual peristimulus EMG waveforms (called VEMP “traces”). In this article, we test VEMPid versus VEMPpp and nVEMPpp to see how well these metrics track the motoneuron inhibition depth. Inhibition depth cannot be directly measured, so to make this comparison we develop a model of VEMP production that includes the inhibition depth as an explicit variable, and we test how well VEMPid, VEMPpp, and nVEMPpp track the model inhibition depth across various conditions. Tracking inhibition depth does not guarantee clinical usefulness, which ultimately must be determined from tests on a large body of clinical data. Nonetheless, determining how well these metrics track inhibition depth provides a first comparison of these metrics and indicates relative performance of the metrics under well-understood conditions. Our test results show that for high response levels, VEMPid and nVEMPpp are approximately equivalent but for low response levels, VEMPid does better than nVEMPpp. However, for the present implementation of VEMPid to work, a recognizable VEMP waveform needs to be recorded (typically from a response at a high sound level), which limits VEMPid’s current value to making sensitive lowlevel measurements only when a good higher-level measurement is available.
METHODS Methods Overview We begin by explaining in detail what we mean by inhibition depth and the approximations needed to get a single-number metric of inhibition depth. We then consider how to obtain a sensitive, linear metric for the VEMP component in each VEMP trace. With these in hand, we describe how VEMPid is computed from the VEMP traces that are obtained during a standard run used to obtain VEMP traces. For comparison, VEMPpp was computed as the PP value of the waveform average of the VEMP
traces, and nVEMPpp was computed by normalizing VEMPpp by the average of the rectified EMG from the VEMP traces. Although the comparison of VEMPid with VEMPpp and nVEMPpp is to be done using a model, we want the comparison to be robust, that is, to apply to data from a variety of subjects and conditions. To do a robust overall comparison, we created a family of models and compared the VEMP metrics using results from each model. Each model had the same overall structure but was individualized to different muscle-specific properties derived from actual EMG data. To get these data, we used five subjects, and on each measured VEMPs in four separate sessions. We constructed the family of VEMP models using anatomical and physiological data, as well as modeling procedures that have been developed in many previous articles. Although some of the techniques we used to estimate muscle properties are different than those used in previous models, our resulting models are basically similar to previous VEMP models (Lütkenhöner et al. 2010, 2011; Lütkenhöner & Basel 2011, 2012). Armed with a family of models that have inhibition depth as an explicit variable, we ran the various models with different values of inhibition depth, thereby creating model VEMP traces that we term “synthetic VEMP data.” From each set of synthetic VEMP data, we computed VEMPid, VEMPpp, and nVEMPpp. These were then compared for how well they tracked the inhibition depth set in the models.
Defining VEMPid and Obtaining VEMPid From VEMP Traces Data on the time course of VEMP motoneuron inhibition come from several sources. Click-induced inhibition was measured in individual SCM fibers in humans by Colebatch and Rothwell (2004). They found a mean duration of inhibition (inhibition was defined as a reduction greater than 50% in two consecutive histogram bins) of 3.6 msec with longer inhibitions for louder clicks. Lütkenhöner and Basel (2012) used a model analysis of VEMP data to calculate a quantity they called the motor-unit “Rate Modulation.” The Rate Modulation had two components, a brief inhibition which roughly corresponds to VEMP waves p13-n23 and a longer, lower-amplitude excitation which roughly corresponds to VEMP waves n34-p44. They concluded that the first component corresponds to the vestibular-mediated motor-unit inhibition that produces the VEMP response of interest here, whereas the second component may have additional input from cochlear responses and is ignored. Lütkenhöner and Basel (2012) noted that although their calculated inhibitory wave was like a Gaussian pulse with a halfwidth of 1.56 msec, their methods could not distinguish this shape from an equivalent rectangular pulse of inhibition. Because we are primarily interested in deriving a single statistic that reflects the saccule’s inhibition, we have defined the inhibition depth as a rectangular pulse of inhibition that is equivalent to the early inhibitory VEMP wave of Lütkenhöner and Basel (2012). We chose 6 msec as the duration of the inhibition depth pulse. Motor unit firing versus time for a 6-msec pulse of inhibition is illustrated in Figure 1. Inhibition with a fixed 6-msec duration is a simplification based on studies of inhibitory postsynaptic potentials evoked by the vestibulocollic reflex (Kushiro et al. 1999), the poststimulus-time histograms of single VEMP motor units (Colebatch & Rothwell 2004) and
593
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
Overall Firing Rate
-20
Saccule-Induced Inhibition Inhibition Depth Motor Unit Spike Sequences
0
20
40
60
TIME AFTER SOUND ONSET (ms)
Fig. 1. Motor-unit firing versus time for a 6-msec pulse of inhibition shown in two ways. These simulated data were produced by assuming that the saccule was excited at time zero and inhibited sternocleidomastoid (SCM) motoneurons in the 6-msec period between the dashed lines. Top: Model results showing the summed firing rate of all motor units that contribute to the SCM electromyogram. The saccule-induced inhibition can be seen by the 6-msec period in which a fraction (the inhibition depth) of spikes is inhibited, reducing the overall firing rate. Bottom: Spike trains from groups of SCM motor units with spikes missing during the interval of inhibition.
the duration of the initial component of the motor unit Rate Modulation of Lütkenhöner and Basel (2012). Our choice of a 6-msec window is longer than the halfamplitude values of Colebatch and Rothwell (2004), and Lütkenhöner and Basel (2012), but the exact duration chosen is not important as long as it is constant and short compared with the VEMP waveform. With the inhibition depth duration set to a fixed 6 msec, the inhibition depth statistic should be viewed as providing a measure of the integral over time of the motor-unit Rate Modulation. The integral is the inhibition depth times inhibition period, so that a 5% inhibition for 6 msec is equivalent to a 10% inhibition for 3 msec, a 20% inhibition for 1.5 msec, and so on.
of VEMP traces acquired at the highest stimulus intensity in the session. We then selected a time window that included all, or almost all, of the VEMP response, for example, 12 to 37 msec for a 500-Hz tone burst. The VEMP template for a session is the resulting waveform segment scaled to have unit energy. The template correlation value (TCV) for a trace is defined as the correlation (or the inner product) of the template with the corresponding part of the VEMP trace (see Fig. 2), which results in a single number for each response trace. One powerful advantage of using the template correlation method is that any ensemble of N traces can be represented as a set of N TCVs, and the size and variability of the VEMP response can be examined by computing the statistics of the TCVs. To understand why the TCV provides a preferred statistic, one must understand the limitations of PP VEMP measurements. VEMPpp is obtained by simple waveform averaging. Waveform averaging is appropriate when each trace consists of the same “true” response embedded in a random noise that is independent of the signal. Many auditory-evoked potentials fit this model, for instance, the auditory brainstem response. However, in the case of a VEMP measurement, the underlying physiology reveals that the “measurement noise” primarily originates in the random nature of the motor-unit activity, and the saccule’s response is encoded by the modulation of this random process. VEMP noise and signal are therefore neither additive nor independent. While averaging may capture some aspects of the signal involved, it ignores other information available in the individual traces. The TCV does a better job of capturing this information. Although the TCV provides a sensitive measure of the amount of VEMP signal present in each VEMP trace, it was not found to be proportional to inhibition depth (Prakash 2009). A metric that is approximately proportional to the inhibition depth TCV= ((a1
Template Window
Quantifying Individual VEMP Responses
a1 a2 a3 Template
Averaged VEMP
VOLTAGE
To obtain a VEMP metric, we chose to quantify how much of the VEMP response signal is present in each trace. The conventional PP amplitude metric for VEMP has two important limitations: (1) it is not linear, therefore the PP amplitude of an average of VEMP traces cannot be derived from the PP amplitudes of the individual VEMP traces, and (2) subject-specific differences in the VEMP waveform morphology strongly influence the PP amplitudes. The use of a linear metric overcomes these limitations and allows computationally powerful linear signal processing algorithms to be applied. Suitable linear metrics could be defined in a number of ways; here, we used one that we call the “template correlation” method. This method resembles a “matched filter,” which is the optimum linear filter for the detection of a known signal (the “template”) in additive noise. Although the VEMP is not a signal in additive noise, the template correlation method, like a matched filter, provides a way of quantifying that part of the signal that resembles the template. Also, because it measures from a wide duration of the waveform (instead of just the two points used by PP), this method provides a greater signal to noise ratio and therefore a greater sensitivity than a PP measurement. To apply this method, we need a template. We formed a VEMP template for each subject and session by averaging a set
b1) + ( a2 b2) + ... +( an bn))
b1 b2 b3 VEMP single trace
0
10
20 30 40 TIME (ms)
50
60
Fig. 2. Cartoon showing a calculation of a template correlation value (TCV). The TCV provides a metric of the amount of template-like signal that is contained in a vestibular evoked myogenic potential (VEMP) single trace. The top trace is the average of a set of traces from the highest-level stimulus and illustrates the VEMP response used for getting the “Template.” The “Template Window” shows the time period selected as containing the VEMP response and the resulting Template is shown as the thick-line part of the trace. To get the TCV, the voltage at each time point in the template is multiplied by the voltage at the corresponding time point in the VEMP single trace, and all these products are summed as shown by the equation at top. To illustrate the correspondences for the first three points in each waveform, the circled regions in the waveforms are shown expanded at right. At the end of the TCV calculation, the points in the Template are normalized by (i.e., divided by) the overall r.m.s. value of the template. For illustration purposes, the VEMP single trace has a very low noise compared with a typical real VEMP single trace.
594
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
is the mean TCV divided by the standard deviation of the TCV. We call this metric the “muscle modulation amplitude” (MMA) because it reflects the modulation of the muscle EMG by the saccule’s response to sound. Because the MMA is approximately proportional to the inhibition depth, an estimate of the inhibition depth can be obtained by multiplying the MMA by a constant, “KMMA.” The proportionality factor between inhibition depth and MMA was found to be approximately 0.2 (with both expressed as percent), that is, KMMA = 0.2 (a fraction, not a percentage).
The Procedure for Calculating VEMPid MMA, and from that VEMPid, can be computed from an ensemble of measured VEMP traces by the following steps: (1) Get the VEMP template. To do this, compute the conventional VEMP average waveform from the traces recorded at the highest stimulus level in a measurement session, and select the time window that contains the VEMP response. This time window will depend on the stimulus used and on the subject’s response, for example, for 500-Hz tone bursts, a window from 12 to 37 msec may be used (see Discussion). Normalize the windowselected waveform segment by dividing each point by the overall root mean square (rms) value of the waveform segment. This normalized waveform segment is the template for the current subject and session. (2) Get the TCV for each trace. That is, compute the correlation of each trace with the template waveform. The correlation is obtained by multiplying each sample of the template by the corresponding sample of the trace and adding together the resulting products (Fig. 2). The result of the template correlation operation on an ensemble of N traces is a set of N numbers (the TCVs); (3) compute the MMA as the mean TCV divided by the standard deviation of the N TCVs. (4) The VEMPid is the MMA multiplied by 0.2. The method described earlier requires that every VEMP trace be recorded and stored; however, there is an alternate method that needs less storage. If, during the recording, running sums are maintained of the TCVs and the squared TCVs, then the MMA and VEMPid can be computed from these sums at the end of the recording. Note that with both methods, the TCV calculation requires using a VEMP template waveform for the session. If all the traces are stored, then the selection of the template waveform and the computation of VEMPid can be done off-line. It is important to note that the experimental setup and procedures for obtaining a VEMPid are the same as those used in obtaining a conventional VEMPpp or nVEMPpp. The placement of VEMP electrodes, the tensing of the muscle, and the application of a sound stimulus are all exactly the same. All that differs is the processing of the resulting traces to calculate TCVs and from them VEMPid.
VEMP Data Acquisition Five healthy adult subjects (3 females, 2 males, ages 20 to 60 years) were recruited for VEMP measurements. Subjects were screened for normal hearing (pure-tone air and bone conduction thresholds at or below 20 dB HL for the octave frequencies 250 to 8000 Hz with no air-bone gaps greater than 10 dB at 250, 500, or 1000 Hz). In addition, subjects were screened for good general health and to exclude a history of vestibular, auditory, neurological, or musculoskeletal disease. An initial VEMP
test was conducted after the Massachusetts Eye and Ear Infirmary (MEEI) standard clinical protocol to ensure that the subjects had a detectable VEMP response at the highest stimulus intensity used (500-Hz tone bursts with a peak sound pressure equivalent to that of a 123-dB SPL tone, i.e., 123-dB peSPL; see Rauch et al. 2004). To measure peSPL, the sound level meter was a Bruel and Kjaer 2203 using a B&K 4144 1” microphone and a NBS-9A Coupler. Informed consent was obtained from each subject in accordance with MEEI and Massachusetts Institute of Technology Institutional Review Boards. The system and procedures routinely used at MEEI for VEMP testing were used to deliver stimuli and to display the running VEMP average and the rms EMG level to the tester. A stimulus train consisted of at least 500 tone bursts (500-Hz frequency, 8-msec duration, two cycle rising and declining phases, both shaped by a Blackman window) of alternating polarity, presented 13 per second for a total duration of 38 to 47 sec. A repetition rate of 13/sec was used to enable us to gather a large amount of data in each session. During tests, each analog stimulus and EMG response was acquired, digitized, and stored using a separate computer to allow for off-line processing of individual responses (traces). Each trace, which consisted of the approximately 77 msec of EMG after the onset of each tone burst, was anti-aliased (i.e., low-pass filtered at 2000 Hz), sampled at 25 kHz, and stored. Each subject was tested in four identical 1.5-hr sessions, spaced approximately 1 week apart. At the start of a session, the SCM was palpated to identify the approximate midpoint (belly) and the sternal attachment (tendon). The EMG was recorded from the right SCM using conductive-gel ECG electrodes in a single differential configuration with the active electrode at the belly, the reference at the tendon, and the ground on the forehead. The stimuli were delivered by a supra-aural headphone (Telephonics TDH-49) to the right ear. The electrodes and headphones remained in place for the duration of the session. Subjects were tested in a sitting position and contracted the SCM muscle by turning the head so that the chin was over the left shoulder. The tester monitored the ongoing rms EMG level and provided verbal feedback to the subject whenever the level fell below 200 μV. During each session, three sets of recordings were made with stimuli of 103, 108, 113, 118, and 123 dB peSPL (the value momentarily reached at the peak of the sound). The subject was instructed to maintain maximal SCM contractions for the first and third recording sets and to maintain half-maximal SCM contractions for the middle recording set. Between recordings, the subject rested (typically for a few minutes, but more if the subject showed signs of fatigue) and was encouraged to relax the neck and turn the head forward but without detaching the electrodes or headphones. We focus mainly on the data from the maximum-strength recordings: for each subject, from the four sessions with two recordings in each session, we obtained eight sets of over 500 traces each, at each of the five stimulus levels. The average VEMP waveforms from the second stimulation at 123-dB peSPL from each session of each subject are shown in Figure 3. Several of these subjects did not yield VEMP waveforms with the classical VEMP shape even though they were measured between electrodes on the belly and tendon of the SCM muscle. We note that others have found VEMP waveforms that were reversed relative to the normal VEMP morphology and that individual SCM muscles sometimes have more than one actionpotential initiation point which will also affect the VEMP morphology (Falla et al. 2002; Wei et al. 2013). We deemed these
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
A
B
C
D
E
Fig. 3. Vestibular evoked myogenic potential (VEMP) waveforms across subjects and sessions. Each column shows data from one subject, and each row shows data from one recording session. Each waveform was from a 500-trace average from a recording at maximum voluntary contraction with 123 dB peSPL, 500-Hz tone bursts. The shaded regions show the time windows that typically contained the VEMP responses. Positive at the muscle belly is plotted up.
waveforms as suitable for our purposes because clinical measurements at MEEI sometimes had VEMP waveforms like these (although not nearly as often as shown in Fig. 3), and an objective of our study was to do a robust test of the VEMP metrics, which means to determine how well the VEMP metrics track inhibition depth even when the VEMP waveform is not the classic shape.
Computational Model Framework To produce ensembles of synthetic VEMP traces, we used a computational model structure made up of two components
595
(Fig. 4). The first component, the vestibulocollic reflex model, summarizes the influence of the vestibular system on the motor system as the resulting inhibition depth of motoneuron firing. A key assumption is that the vestibular system influences the motor system in the same way across subjects and that differences in the stimulus to the vestibular system and in vestibular sensitivity are translated to differences in the resulting motoneuron inhibition depth. Thus, this part of the model has the same structure for all subjects but varies in inhibition depth. The second component, the EMG-production model, simulates the production of the surface EMG signal arising from muscle activity due to voluntary contraction effort. The EMG-production model takes into account the muscle organization of the subject, as well as the electrode placement and recording session characteristics. It is therefore assumed to vary across subjects and sessions. This dual model structure allowed making multiple models that were tailored to a variety of individual subjects and sessions and provided a range of model variation that enabled a robust comparison of how well each VEMP metric tracked inhibition depth. Most features of our model were adapted from previous modeling work (see References), but the use of a fixed, simple vestibulocollic reflex model followed by an EMG-production model tailored to the specific subject and session is unique to this article.
The EMG-Production Model The EMG model simulates the production of the SCM surface EMG, as illustrated in Figure 4. This model is similar to the ones used in prior studies that examined either surface EMG produced by a muscle and/or motor-unit statistics related to muscle force (Clamann 1969; Brody et al. 1974; De Luca 1979; Kukulka & Clamann 1981; Duchêne & Hogrel 2000; Taylor et al. 2003; McGill 2004; Zhou & Rymer 2004). In our EMG model, the muscle is considered to be a collection
Fig. 4. The physiologic model for the production of synthetic VEMP traces. The Vestibular Model (top left) includes the saccular response to sound and the transmission + conversion of this response to motor-neuron inhibition which is summarized by the inhibition depth. The Motor Drive (bottom left) represents the subject’s descending control of motor-neuron spike generation and ultimately muscle contraction strength. The spike inhibition box represents the deletion of a percentage (the inhibition depth) of motor-neuron spikes in the brief time interval T. The volume conductor box shows the production of each motor unit’s contribution to the electromyographic (EMG) as the convolution of each spike with a motor-unit action-potential (MUAP) surface response. This part of the model was adjusted for the motor system of each subject and session by choosing a MUAP that mimicked the autocorrelation of the EMG from that subject and session (a MUAPss). One execution of the model produced a single synthetic VEMP trace. An ensemble of synthetic traces (analogous to the set of real VEMP traces obtained in a VEMP average run from a subject) was produced by executing an identical model N times, each with a different randomization of the spike generation and spike inhibition on each motor unit. See the text for further details.
596
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
of motor units, with each unit contributing a fraction of the overall surface EMG (Burke, 1981). The size of a motor unit’s EMG contribution depends on the number of muscle fibers in the unit and the distance of the unit from the EMG electrodes. Based on anatomical and physiological studies (Galea et al. 1991; Lille & Bauer 1994; Kamibayashi & Richmond 1998; McComas 1998; Routal & Pal 2000; Hamilton et al. 2004), the number of motor units and the total number of muscle fibers were set to 250 and 50,000, respectively, for the SCM muscle. The EMG-production characteristics of all motor units were assumed to be similar, a simplification that was assumed in previous VEMP models (e.g., Wit & Kingma 2006; Lütkenhöner & Basel 2011). The fiber counts in the motor units were exponentially distributed (Fuglevand et al. 1993). The threshold of spiking (i.e., the motor-drive value at which a motor unit begins to fire) and the spike rate (i.e., average number of action potentials per unit time) of each unit were determined by a recruitment rule based on Henneman’s “size principle,” which recruits from the smallest to the largest motor unit (Henneman 1957; Buchthal & Schmalbruch 1980). Each motor unit was assumed to generate a random train of motor-unit action potentials (MUAPs) at a rate determined by the corticospinal motor drive. The interspike intervals of each motor unit were independent and Gaussian distributed (Clamann 1969) with a mean value equal to the inverse of the spike rate assigned by the recruitment rule and a standard deviation of 20% of the mean (negative and zero interspike intervals were excluded). The spike trains from different motor units were statistically independent. Because the EMG-production model simulates the randomness inherent in motor-neuron firing, the resulting spike times vary across individual simulation runs, even when all parameters in the EMG model are constant. This randomness in spike timing is a main contributor to the randomness of the synthetic VEMP traces, as well as to the variability of the VEMP averages seen clinically. Because the recordings from our VEMP subjects showed that EMG properties can vary across subjects and sessions (Fig. 3), we adapted the EMG-production model to the characteristics of each subject and session. Although each motor unit of a subject produces a unique pattern at the EMG electrodes, we sought to mimic only the overall behavior of an individual muscle. Therefore, for simplicity, we assumed that for a given subject and session, all motor units in the muscle produced the same surface response waveform shape (Wit & Kingma 2006; Lütkenhöner & Basel 2011) but different amplitudes depending on the motor-unit size and position relative to the electrodes. We refer to the surface response waveform as the MUAP of the subject and session or “MUAPss.” We obtained the relevant EMG characteristics for each subject and session from their EMG autocorrelation functions. The superposition of uncorrelated spike sequences from many motor units has an approximately flat spectrum, so that the shape of the autocorrelation of each EMG waveform must arise largely from the temporal structure of the MUAPss waveform. The EMG model was customized to each subject and session by creating individual MUAPss waveforms that produced synthetic-EMG autocorrelation functions similar to the experimental EMG autocorrelation from that subject and session. For the experimental EMG autocorrelation, we used the EMG obtained at the lowest (103 dB peSPL) stimulus intensity (i.e., the waveform from the whole time period of all the traces from the
VEMP recording at 103 dB peSPL). Although some subjects appeared to show a small VEMP response at 103 dB peSPL, this would have little influence on the EMG autocorrelation because the VEMP response duration is short compared with the overall trace interval, and, as will be seen later, even at the highest sound level, our analysis indicates that the motor units are inhibited only by a small percentage during a VEMP response. For each subject and session, we created a MUAPss waveform that adequately matched the experimentally obtained EMG autocorrelation function. The computation of a waveform from its autocorrelation is an ill-posed inverse problem (detailed for EMGs in Brody et al. 1974). The problem was made tractable by constructing each MUAPss function as a weighted sum of three Gaussian pulses, each with three parameters: width, amplitude, and relative time-shift. With one pulse set to a shift of zero and an amplitude of one, each MUAPss was described by seven free parameters. The pulses were further constrained to have the integral of a MUAPss be zero so the EMG would have zero mean. For each MUAPss, the pulse parameters were determined by an iterative steepest-descent algorithm to minimize the total squared error between the autocorrelations of the MUAPss waveform and the experimental EMG. Using the methods outlined earlier, MUAPss waveforms were computed for each of our 4 × 5 subject/session combinations (Fig. 5). The fitting error (area in the difference between model and experimental autocorrelations) was less than 1.1% of the area under the autocorrelation function. Note that the goal of computing the MUAPss is to match statistics of the simulated EMG to the statistics of the experimental (i.e., actual) EMG. An infinite number of different waveforms can yield the same autocorrelation function (Papoulis & Pillai 2002). Our method computes one choice that matches the autocorrelation function for each MUAPss. It is not expected that this waveform matches the experimental MUAP. The point is to provide a variety of models (20 models) that reflect the experimental EMG variation and to use the models to compare how well the three VEMP metrics track VEMP inhibition depth.
The Motor Drive Each EMG model receives a motor drive that is expressed as a percentage of maximum voluntary contraction. At a motor drive of 100%, all motor units are active at their maximum firing rates. For other motor drives, motor-unit rates were assigned by the recruitment rule. In each EMG-production model, there is a scaling parameter that specifies the amplitude of the MUAPss waveform. The rms value of the synthetic EMG depends on both the scaling parameter and the motor drive. To unambiguously determine the values of these two parameters, the scale parameter was set such that the synthetic EMG at a motor drive of 95% had the same rms value as the largest experimental EMG value recorded during the session. Ninety-five percentage represents an estimate of how close the subject came on their strongest SCM activation to a full motor activation. During the simulations used to compare the three VEMP metrics, the motor drive was varied to mimic the variations in motor drive in the corresponding subject and session. The variation in the motor drive for each run from a subject was obtained from the power spectrum of the overall (40 sec or more) EMG time record. Each run showed an increase in energy at a very-low,
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
Autocorrelation
NORMALIZED AMPLITUDE (sessions offset vertically)
Subject A
Motor Unit Action Potential
Subject B
597
vestibular model generates a signal that modulates the spike rates of the motor units in the EMG-production model (Fig. 4, bottom). The vestibular model produces an inhibition signal that deletes spikes that fall within a time interval T relative to the acoustic stimulus, with a probability of deletion determined by the inhibition depth. The thinned spike sequences are then convolved with the MUAPss to generate each motor unit’s contribution to the surface EMG trace.
Generating Synthetic VEMP Traces Each synthetic VEMP trace was generated by a single execution of the VEMP model consisting of (1) the vestibular and EMG models with specific values of motor-drive and Inhibition Depth parameters, and (2) the EMG-production model with a MUAPss corresponding to one subject and session, all implemented in MATLAB. As in a run in a real VEMP average, a typical ensemble of simulated VEMP traces consisted of a large number, N, of traces (e.g., 200 to 500) obtained by executing the model N times. In each simulation of N traces, the inhibition depth was constant and the motor drive was varied to match the variation of that subject and session. Note that even when the model parameters remained the same, successive executions of the model differed in the probabilistic sequence of motor-unit spikes and in which spikes were inhibited. Because deleted spikes were evenly distributed across interval T, and each synthetic trace was the convolution of a spike sequence and the spike’s surface EMG response in that subject and session (the MUAPss), an approximation of the mean of the traces is the convolution of the MUAPss and a 6-msec rectangular Inhibition Depth change in the motor-unit firing rate. This approximation yields the average VEMP but does not capture the variation across trials, which is an important component of a VEMP measurement, one that influences the accuracy of individual runs.
Subject C
Subject D
Subject E
RESULTS -10
0
10
-10
0
TIME re PEAK (ms)
10
Fig. 5. The electromyographic (EMG) autocorrelations (left) and the corresponding motor-unit action potentials (MUAPss) (right) for each subject and session. At left, each experimental autocorrelation is shown by a solid line, and the synthetic autocorrelation reconstructed from the MUAPss at right is shown as a dashed line offset lower for visibility. The MUAPss reconstruction does not yield a latency, and so the MUAPs at right are centered at zero. The autocorrelations at left, by their nature, are centered at zero. Because there is no useful information in the amplitudes of the autocorrelations or MUAPsss, they are shown on a normalized (dimensionless) scale with the peak values all scaled to the same amplitude.
subject-specific frequency (
Objectives: Vestibular evoked myogenic potentials (VEMPs) are due to vestibular responses producing brief inhibitions of muscle contractions that are detectable in electromyographic (EMG) responses. VEMP amplitudes are traditionally measured by the peak to peak amplitude of the averaged EMG response (VEMPpp) or by a normalized VEMPpp (nVEMPpp). However, a brief EMG inhibition does not satisfy the statistical assumptions for the average to be the optimal processing strategy. Here, it is postulated that the inhibition depth of motoneuron firing is the desired metric for showing the influence of the vestibular system on the muscle system. The authors present a metric called “VEMPid” that estimates this inhibition depth from the EMG data obtained in a usual VEMP data acquisition. The goal of this article was to compare how well VEMPid, VEMPpp, and nVEMPpp track inhibition depth.
(EMG) potential of the corresponding muscle. The inhibitory strength of such a volley has been quantified by the change in the average EMG, taking into account that the EMG response is proportional to the level of muscle contraction (Widmer & Lund 1989). Although averaging is the optimal signal processing technique when measuring a fixed signal in an independent random noise, these conditions are not satisfied for measurements of averages of briefly inhibited EMGs. The brief inhibition of an ongoing EMG corresponds better to a signal that modulates a noise than to a signal that is added to an independent noise. Thus, a signal processing method that is different than simple averaging may produce better results for quantifying such an inhibition. Here, we present a signal processing method for quantifying a brief inhibition in the context of a clinically used vestibular test, but it can be used to quantify a brief inhibition or a subthreshold excitation in ongoing responses of any motoneuron pool. Activation of the vestibular system by brief sounds or by vibrational stimuli produces a brief inhibition of various muscles. The inhibitory signal from vestibular system cannot be measured directly, but its effect on the EMG of a target muscle can be measured, and this effect is used as a clinical index of vestibular sensitivity. In particular, the vestibular evoked myogenic potential (VEMP), a surface EMG potential, is used to assess the functioning of the saccule (Welgampola & Colebatch 2005; Curthoys 2010; Rosengren et al. 2010). Saccule excitation activates vestibular reflex pathways that inhibit motoneurons that innervate the ipsilateral sternocleidomastoid (SCM) muscle (Colebatch et al. 1994; Colebatch 2001; Murofushi and Kaga 2009). Although VEMP responses can be measured from many head and neck muscles, in this study we concentrate on the response to sound measured from the ipsilateral SCM, which is sometimes called the cervical VEMP (cVEMP) response. Previous studies on VEMP potentials focused on the average stimulus-locked VEMP waveform, particularly its peak to peak (PP) amplitude, which we refer to as “VEMPpp.” Obtaining a VEMP metric by averaging the rectified EMG signal has also been tried, but it is not as good as simple averaging (Rosengren et al. 2010). VEMPpp amplitudes scale approximately linearly with SCM contraction strength although saturation of VEMPpp amplitude has been reported for some types of SCM contractions (Bogle et al. 2013). To minimize the effects of different contraction strengths, VEMPs are commonly “normalized” by dividing the VEMPpp by a measure of the ongoing EMG response, for example, by the average amplitude of the rectified EMG. This yields a normalized VEMPpp or an “nVEMPpp.” VEMP responses arise because the saccule responds to sound and through neural pathways produces inhibition of SCM motoneurons (Colebatch & Rothwell 2004). Measures such as
Design: To find a robust method to compare VEMPid, VEMPpp, and nVEMPpp, realistic physiological models for the inhibition of VEMP EMG signals were made using VEMP data from four measurement sessions on each of the five normal subjects. Each of the resulting 20 EMGproduction models was adjusted to match the EMG autocorrelation of an individual subject and session. Simulated VEMP traces produced by these models were used to compare how well VEMPid, VEMPpp, and nVEMPpp tracked model inhibition depth. Results: Applied to simulated and real VEMP data, VEMPid showed good test-retest consistency and greater sensitivity at low stimulus levels than VEMPpp or nVEMPpp. For large-amplitude responses, nVEMPpp and VEMPid were equivalent in their consistency across subjects and sessions, but for low-amplitude responses, VEMPid was superior. Unnormalized VEMPpp was always worse than nVEMPpp or VEMPid. Conclusions: VEMPid provides a more reliable measurement of vestibular function at low sound levels than the traditional nVEMPpp, without requiring a change in how VEMP tests are performed. The calculation method for VEMPid should be applicable whenever an ongoing muscle contraction is briefly inhibited by an external stimulus. Key words: cVEMP, Vestibulocollic reflex.
Electromyogram,
EMG,
oVEMP,
VEMP,
(Ear & Hearing 2015;36;591–604)
INTRODUCTION A brief inhibitory volley delivered to a motoneuron pool can be detected by the change it produces in the electromyographic 1 Harvard-MIT Division of Health Science and Technology, Speech and Hearing Bioscience and Technology Program, Massachusetts Institute of Technology, Cambridge, Massachusetts, USA; 2Department of Otology and Laryngology, Harvard Medical School, Boston, Massachusetts, USA; 3 Department of Otolaryngology, Massachusetts Eye and Ear Infirmary, Boston, Massachusetts, USA; 4Department of Audiology, Massachusetts Eye and Ear Infirmary, Boston, Massachusetts, USA; and 5Eaton-Peabody Laboratories, Massachusetts Eye and Ear Infirmary, Boston, Massachusetts, USA.
0196/0202/2015/365-0591/0 • Ear & Hearing • Copyright © 2015 Wolters Kluwer Health, Inc. All rights reserved • Printed in the U.S.A. 591
592
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
the VEMPpp provide only an indirect measure of this inhibition. A metric that would be a step closer to the saccule’s response to the stimulus is the decrement in the motor-unit spike rate. With this in mind, we propose the “Inhibition Depth,” defined as the percentage reduction in the composite spike rate of SCM motoneurons produced by the activation of the saccule, as an index of motoneuron inhibition (Prakash 2009). The modulation of motoneuron firing rate was also identified by Lütkenhöner and Basel (2012) as being the key variable to quantify the saccular effect on the motoneuron pool. For clinical use, we developed a single-number index for VEMP inhibition that is termed the VEMP inhibition depth or “VEMPid” (Prakash 2009). The goal for this index was that it would represent the inhibition depth and that for a given stimulus, subjects with similar saccular sensitivity would have similar VEMPids. This means that the VEMPid metric should be constant even though the VEMP average waveforms may be different due to different muscle contractions, different placements of the VEMP electrodes, or other nonvestibular factors. If this robustness is achieved, VEMPid would be a clinically significant measure of saccule function and a useful indicator of pathological changes in the saccule. The VEMPid is intended as an alternative to VEMPpp and nVEMPpp. VEMPid can be calculated from the EMG data acquired using the already well-established experimental protocols for recording VEMP responses (Curthoys 2010; Rosengren et al. 2010), that is, VEMPid is derived from the set of individual peristimulus EMG waveforms (called VEMP “traces”). In this article, we test VEMPid versus VEMPpp and nVEMPpp to see how well these metrics track the motoneuron inhibition depth. Inhibition depth cannot be directly measured, so to make this comparison we develop a model of VEMP production that includes the inhibition depth as an explicit variable, and we test how well VEMPid, VEMPpp, and nVEMPpp track the model inhibition depth across various conditions. Tracking inhibition depth does not guarantee clinical usefulness, which ultimately must be determined from tests on a large body of clinical data. Nonetheless, determining how well these metrics track inhibition depth provides a first comparison of these metrics and indicates relative performance of the metrics under well-understood conditions. Our test results show that for high response levels, VEMPid and nVEMPpp are approximately equivalent but for low response levels, VEMPid does better than nVEMPpp. However, for the present implementation of VEMPid to work, a recognizable VEMP waveform needs to be recorded (typically from a response at a high sound level), which limits VEMPid’s current value to making sensitive lowlevel measurements only when a good higher-level measurement is available.
METHODS Methods Overview We begin by explaining in detail what we mean by inhibition depth and the approximations needed to get a single-number metric of inhibition depth. We then consider how to obtain a sensitive, linear metric for the VEMP component in each VEMP trace. With these in hand, we describe how VEMPid is computed from the VEMP traces that are obtained during a standard run used to obtain VEMP traces. For comparison, VEMPpp was computed as the PP value of the waveform average of the VEMP
traces, and nVEMPpp was computed by normalizing VEMPpp by the average of the rectified EMG from the VEMP traces. Although the comparison of VEMPid with VEMPpp and nVEMPpp is to be done using a model, we want the comparison to be robust, that is, to apply to data from a variety of subjects and conditions. To do a robust overall comparison, we created a family of models and compared the VEMP metrics using results from each model. Each model had the same overall structure but was individualized to different muscle-specific properties derived from actual EMG data. To get these data, we used five subjects, and on each measured VEMPs in four separate sessions. We constructed the family of VEMP models using anatomical and physiological data, as well as modeling procedures that have been developed in many previous articles. Although some of the techniques we used to estimate muscle properties are different than those used in previous models, our resulting models are basically similar to previous VEMP models (Lütkenhöner et al. 2010, 2011; Lütkenhöner & Basel 2011, 2012). Armed with a family of models that have inhibition depth as an explicit variable, we ran the various models with different values of inhibition depth, thereby creating model VEMP traces that we term “synthetic VEMP data.” From each set of synthetic VEMP data, we computed VEMPid, VEMPpp, and nVEMPpp. These were then compared for how well they tracked the inhibition depth set in the models.
Defining VEMPid and Obtaining VEMPid From VEMP Traces Data on the time course of VEMP motoneuron inhibition come from several sources. Click-induced inhibition was measured in individual SCM fibers in humans by Colebatch and Rothwell (2004). They found a mean duration of inhibition (inhibition was defined as a reduction greater than 50% in two consecutive histogram bins) of 3.6 msec with longer inhibitions for louder clicks. Lütkenhöner and Basel (2012) used a model analysis of VEMP data to calculate a quantity they called the motor-unit “Rate Modulation.” The Rate Modulation had two components, a brief inhibition which roughly corresponds to VEMP waves p13-n23 and a longer, lower-amplitude excitation which roughly corresponds to VEMP waves n34-p44. They concluded that the first component corresponds to the vestibular-mediated motor-unit inhibition that produces the VEMP response of interest here, whereas the second component may have additional input from cochlear responses and is ignored. Lütkenhöner and Basel (2012) noted that although their calculated inhibitory wave was like a Gaussian pulse with a halfwidth of 1.56 msec, their methods could not distinguish this shape from an equivalent rectangular pulse of inhibition. Because we are primarily interested in deriving a single statistic that reflects the saccule’s inhibition, we have defined the inhibition depth as a rectangular pulse of inhibition that is equivalent to the early inhibitory VEMP wave of Lütkenhöner and Basel (2012). We chose 6 msec as the duration of the inhibition depth pulse. Motor unit firing versus time for a 6-msec pulse of inhibition is illustrated in Figure 1. Inhibition with a fixed 6-msec duration is a simplification based on studies of inhibitory postsynaptic potentials evoked by the vestibulocollic reflex (Kushiro et al. 1999), the poststimulus-time histograms of single VEMP motor units (Colebatch & Rothwell 2004) and
593
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
Overall Firing Rate
-20
Saccule-Induced Inhibition Inhibition Depth Motor Unit Spike Sequences
0
20
40
60
TIME AFTER SOUND ONSET (ms)
Fig. 1. Motor-unit firing versus time for a 6-msec pulse of inhibition shown in two ways. These simulated data were produced by assuming that the saccule was excited at time zero and inhibited sternocleidomastoid (SCM) motoneurons in the 6-msec period between the dashed lines. Top: Model results showing the summed firing rate of all motor units that contribute to the SCM electromyogram. The saccule-induced inhibition can be seen by the 6-msec period in which a fraction (the inhibition depth) of spikes is inhibited, reducing the overall firing rate. Bottom: Spike trains from groups of SCM motor units with spikes missing during the interval of inhibition.
the duration of the initial component of the motor unit Rate Modulation of Lütkenhöner and Basel (2012). Our choice of a 6-msec window is longer than the halfamplitude values of Colebatch and Rothwell (2004), and Lütkenhöner and Basel (2012), but the exact duration chosen is not important as long as it is constant and short compared with the VEMP waveform. With the inhibition depth duration set to a fixed 6 msec, the inhibition depth statistic should be viewed as providing a measure of the integral over time of the motor-unit Rate Modulation. The integral is the inhibition depth times inhibition period, so that a 5% inhibition for 6 msec is equivalent to a 10% inhibition for 3 msec, a 20% inhibition for 1.5 msec, and so on.
of VEMP traces acquired at the highest stimulus intensity in the session. We then selected a time window that included all, or almost all, of the VEMP response, for example, 12 to 37 msec for a 500-Hz tone burst. The VEMP template for a session is the resulting waveform segment scaled to have unit energy. The template correlation value (TCV) for a trace is defined as the correlation (or the inner product) of the template with the corresponding part of the VEMP trace (see Fig. 2), which results in a single number for each response trace. One powerful advantage of using the template correlation method is that any ensemble of N traces can be represented as a set of N TCVs, and the size and variability of the VEMP response can be examined by computing the statistics of the TCVs. To understand why the TCV provides a preferred statistic, one must understand the limitations of PP VEMP measurements. VEMPpp is obtained by simple waveform averaging. Waveform averaging is appropriate when each trace consists of the same “true” response embedded in a random noise that is independent of the signal. Many auditory-evoked potentials fit this model, for instance, the auditory brainstem response. However, in the case of a VEMP measurement, the underlying physiology reveals that the “measurement noise” primarily originates in the random nature of the motor-unit activity, and the saccule’s response is encoded by the modulation of this random process. VEMP noise and signal are therefore neither additive nor independent. While averaging may capture some aspects of the signal involved, it ignores other information available in the individual traces. The TCV does a better job of capturing this information. Although the TCV provides a sensitive measure of the amount of VEMP signal present in each VEMP trace, it was not found to be proportional to inhibition depth (Prakash 2009). A metric that is approximately proportional to the inhibition depth TCV= ((a1
Template Window
Quantifying Individual VEMP Responses
a1 a2 a3 Template
Averaged VEMP
VOLTAGE
To obtain a VEMP metric, we chose to quantify how much of the VEMP response signal is present in each trace. The conventional PP amplitude metric for VEMP has two important limitations: (1) it is not linear, therefore the PP amplitude of an average of VEMP traces cannot be derived from the PP amplitudes of the individual VEMP traces, and (2) subject-specific differences in the VEMP waveform morphology strongly influence the PP amplitudes. The use of a linear metric overcomes these limitations and allows computationally powerful linear signal processing algorithms to be applied. Suitable linear metrics could be defined in a number of ways; here, we used one that we call the “template correlation” method. This method resembles a “matched filter,” which is the optimum linear filter for the detection of a known signal (the “template”) in additive noise. Although the VEMP is not a signal in additive noise, the template correlation method, like a matched filter, provides a way of quantifying that part of the signal that resembles the template. Also, because it measures from a wide duration of the waveform (instead of just the two points used by PP), this method provides a greater signal to noise ratio and therefore a greater sensitivity than a PP measurement. To apply this method, we need a template. We formed a VEMP template for each subject and session by averaging a set
b1) + ( a2 b2) + ... +( an bn))
b1 b2 b3 VEMP single trace
0
10
20 30 40 TIME (ms)
50
60
Fig. 2. Cartoon showing a calculation of a template correlation value (TCV). The TCV provides a metric of the amount of template-like signal that is contained in a vestibular evoked myogenic potential (VEMP) single trace. The top trace is the average of a set of traces from the highest-level stimulus and illustrates the VEMP response used for getting the “Template.” The “Template Window” shows the time period selected as containing the VEMP response and the resulting Template is shown as the thick-line part of the trace. To get the TCV, the voltage at each time point in the template is multiplied by the voltage at the corresponding time point in the VEMP single trace, and all these products are summed as shown by the equation at top. To illustrate the correspondences for the first three points in each waveform, the circled regions in the waveforms are shown expanded at right. At the end of the TCV calculation, the points in the Template are normalized by (i.e., divided by) the overall r.m.s. value of the template. For illustration purposes, the VEMP single trace has a very low noise compared with a typical real VEMP single trace.
594
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
is the mean TCV divided by the standard deviation of the TCV. We call this metric the “muscle modulation amplitude” (MMA) because it reflects the modulation of the muscle EMG by the saccule’s response to sound. Because the MMA is approximately proportional to the inhibition depth, an estimate of the inhibition depth can be obtained by multiplying the MMA by a constant, “KMMA.” The proportionality factor between inhibition depth and MMA was found to be approximately 0.2 (with both expressed as percent), that is, KMMA = 0.2 (a fraction, not a percentage).
The Procedure for Calculating VEMPid MMA, and from that VEMPid, can be computed from an ensemble of measured VEMP traces by the following steps: (1) Get the VEMP template. To do this, compute the conventional VEMP average waveform from the traces recorded at the highest stimulus level in a measurement session, and select the time window that contains the VEMP response. This time window will depend on the stimulus used and on the subject’s response, for example, for 500-Hz tone bursts, a window from 12 to 37 msec may be used (see Discussion). Normalize the windowselected waveform segment by dividing each point by the overall root mean square (rms) value of the waveform segment. This normalized waveform segment is the template for the current subject and session. (2) Get the TCV for each trace. That is, compute the correlation of each trace with the template waveform. The correlation is obtained by multiplying each sample of the template by the corresponding sample of the trace and adding together the resulting products (Fig. 2). The result of the template correlation operation on an ensemble of N traces is a set of N numbers (the TCVs); (3) compute the MMA as the mean TCV divided by the standard deviation of the N TCVs. (4) The VEMPid is the MMA multiplied by 0.2. The method described earlier requires that every VEMP trace be recorded and stored; however, there is an alternate method that needs less storage. If, during the recording, running sums are maintained of the TCVs and the squared TCVs, then the MMA and VEMPid can be computed from these sums at the end of the recording. Note that with both methods, the TCV calculation requires using a VEMP template waveform for the session. If all the traces are stored, then the selection of the template waveform and the computation of VEMPid can be done off-line. It is important to note that the experimental setup and procedures for obtaining a VEMPid are the same as those used in obtaining a conventional VEMPpp or nVEMPpp. The placement of VEMP electrodes, the tensing of the muscle, and the application of a sound stimulus are all exactly the same. All that differs is the processing of the resulting traces to calculate TCVs and from them VEMPid.
VEMP Data Acquisition Five healthy adult subjects (3 females, 2 males, ages 20 to 60 years) were recruited for VEMP measurements. Subjects were screened for normal hearing (pure-tone air and bone conduction thresholds at or below 20 dB HL for the octave frequencies 250 to 8000 Hz with no air-bone gaps greater than 10 dB at 250, 500, or 1000 Hz). In addition, subjects were screened for good general health and to exclude a history of vestibular, auditory, neurological, or musculoskeletal disease. An initial VEMP
test was conducted after the Massachusetts Eye and Ear Infirmary (MEEI) standard clinical protocol to ensure that the subjects had a detectable VEMP response at the highest stimulus intensity used (500-Hz tone bursts with a peak sound pressure equivalent to that of a 123-dB SPL tone, i.e., 123-dB peSPL; see Rauch et al. 2004). To measure peSPL, the sound level meter was a Bruel and Kjaer 2203 using a B&K 4144 1” microphone and a NBS-9A Coupler. Informed consent was obtained from each subject in accordance with MEEI and Massachusetts Institute of Technology Institutional Review Boards. The system and procedures routinely used at MEEI for VEMP testing were used to deliver stimuli and to display the running VEMP average and the rms EMG level to the tester. A stimulus train consisted of at least 500 tone bursts (500-Hz frequency, 8-msec duration, two cycle rising and declining phases, both shaped by a Blackman window) of alternating polarity, presented 13 per second for a total duration of 38 to 47 sec. A repetition rate of 13/sec was used to enable us to gather a large amount of data in each session. During tests, each analog stimulus and EMG response was acquired, digitized, and stored using a separate computer to allow for off-line processing of individual responses (traces). Each trace, which consisted of the approximately 77 msec of EMG after the onset of each tone burst, was anti-aliased (i.e., low-pass filtered at 2000 Hz), sampled at 25 kHz, and stored. Each subject was tested in four identical 1.5-hr sessions, spaced approximately 1 week apart. At the start of a session, the SCM was palpated to identify the approximate midpoint (belly) and the sternal attachment (tendon). The EMG was recorded from the right SCM using conductive-gel ECG electrodes in a single differential configuration with the active electrode at the belly, the reference at the tendon, and the ground on the forehead. The stimuli were delivered by a supra-aural headphone (Telephonics TDH-49) to the right ear. The electrodes and headphones remained in place for the duration of the session. Subjects were tested in a sitting position and contracted the SCM muscle by turning the head so that the chin was over the left shoulder. The tester monitored the ongoing rms EMG level and provided verbal feedback to the subject whenever the level fell below 200 μV. During each session, three sets of recordings were made with stimuli of 103, 108, 113, 118, and 123 dB peSPL (the value momentarily reached at the peak of the sound). The subject was instructed to maintain maximal SCM contractions for the first and third recording sets and to maintain half-maximal SCM contractions for the middle recording set. Between recordings, the subject rested (typically for a few minutes, but more if the subject showed signs of fatigue) and was encouraged to relax the neck and turn the head forward but without detaching the electrodes or headphones. We focus mainly on the data from the maximum-strength recordings: for each subject, from the four sessions with two recordings in each session, we obtained eight sets of over 500 traces each, at each of the five stimulus levels. The average VEMP waveforms from the second stimulation at 123-dB peSPL from each session of each subject are shown in Figure 3. Several of these subjects did not yield VEMP waveforms with the classical VEMP shape even though they were measured between electrodes on the belly and tendon of the SCM muscle. We note that others have found VEMP waveforms that were reversed relative to the normal VEMP morphology and that individual SCM muscles sometimes have more than one actionpotential initiation point which will also affect the VEMP morphology (Falla et al. 2002; Wei et al. 2013). We deemed these
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
A
B
C
D
E
Fig. 3. Vestibular evoked myogenic potential (VEMP) waveforms across subjects and sessions. Each column shows data from one subject, and each row shows data from one recording session. Each waveform was from a 500-trace average from a recording at maximum voluntary contraction with 123 dB peSPL, 500-Hz tone bursts. The shaded regions show the time windows that typically contained the VEMP responses. Positive at the muscle belly is plotted up.
waveforms as suitable for our purposes because clinical measurements at MEEI sometimes had VEMP waveforms like these (although not nearly as often as shown in Fig. 3), and an objective of our study was to do a robust test of the VEMP metrics, which means to determine how well the VEMP metrics track inhibition depth even when the VEMP waveform is not the classic shape.
Computational Model Framework To produce ensembles of synthetic VEMP traces, we used a computational model structure made up of two components
595
(Fig. 4). The first component, the vestibulocollic reflex model, summarizes the influence of the vestibular system on the motor system as the resulting inhibition depth of motoneuron firing. A key assumption is that the vestibular system influences the motor system in the same way across subjects and that differences in the stimulus to the vestibular system and in vestibular sensitivity are translated to differences in the resulting motoneuron inhibition depth. Thus, this part of the model has the same structure for all subjects but varies in inhibition depth. The second component, the EMG-production model, simulates the production of the surface EMG signal arising from muscle activity due to voluntary contraction effort. The EMG-production model takes into account the muscle organization of the subject, as well as the electrode placement and recording session characteristics. It is therefore assumed to vary across subjects and sessions. This dual model structure allowed making multiple models that were tailored to a variety of individual subjects and sessions and provided a range of model variation that enabled a robust comparison of how well each VEMP metric tracked inhibition depth. Most features of our model were adapted from previous modeling work (see References), but the use of a fixed, simple vestibulocollic reflex model followed by an EMG-production model tailored to the specific subject and session is unique to this article.
The EMG-Production Model The EMG model simulates the production of the SCM surface EMG, as illustrated in Figure 4. This model is similar to the ones used in prior studies that examined either surface EMG produced by a muscle and/or motor-unit statistics related to muscle force (Clamann 1969; Brody et al. 1974; De Luca 1979; Kukulka & Clamann 1981; Duchêne & Hogrel 2000; Taylor et al. 2003; McGill 2004; Zhou & Rymer 2004). In our EMG model, the muscle is considered to be a collection
Fig. 4. The physiologic model for the production of synthetic VEMP traces. The Vestibular Model (top left) includes the saccular response to sound and the transmission + conversion of this response to motor-neuron inhibition which is summarized by the inhibition depth. The Motor Drive (bottom left) represents the subject’s descending control of motor-neuron spike generation and ultimately muscle contraction strength. The spike inhibition box represents the deletion of a percentage (the inhibition depth) of motor-neuron spikes in the brief time interval T. The volume conductor box shows the production of each motor unit’s contribution to the electromyographic (EMG) as the convolution of each spike with a motor-unit action-potential (MUAP) surface response. This part of the model was adjusted for the motor system of each subject and session by choosing a MUAP that mimicked the autocorrelation of the EMG from that subject and session (a MUAPss). One execution of the model produced a single synthetic VEMP trace. An ensemble of synthetic traces (analogous to the set of real VEMP traces obtained in a VEMP average run from a subject) was produced by executing an identical model N times, each with a different randomization of the spike generation and spike inhibition on each motor unit. See the text for further details.
596
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
of motor units, with each unit contributing a fraction of the overall surface EMG (Burke, 1981). The size of a motor unit’s EMG contribution depends on the number of muscle fibers in the unit and the distance of the unit from the EMG electrodes. Based on anatomical and physiological studies (Galea et al. 1991; Lille & Bauer 1994; Kamibayashi & Richmond 1998; McComas 1998; Routal & Pal 2000; Hamilton et al. 2004), the number of motor units and the total number of muscle fibers were set to 250 and 50,000, respectively, for the SCM muscle. The EMG-production characteristics of all motor units were assumed to be similar, a simplification that was assumed in previous VEMP models (e.g., Wit & Kingma 2006; Lütkenhöner & Basel 2011). The fiber counts in the motor units were exponentially distributed (Fuglevand et al. 1993). The threshold of spiking (i.e., the motor-drive value at which a motor unit begins to fire) and the spike rate (i.e., average number of action potentials per unit time) of each unit were determined by a recruitment rule based on Henneman’s “size principle,” which recruits from the smallest to the largest motor unit (Henneman 1957; Buchthal & Schmalbruch 1980). Each motor unit was assumed to generate a random train of motor-unit action potentials (MUAPs) at a rate determined by the corticospinal motor drive. The interspike intervals of each motor unit were independent and Gaussian distributed (Clamann 1969) with a mean value equal to the inverse of the spike rate assigned by the recruitment rule and a standard deviation of 20% of the mean (negative and zero interspike intervals were excluded). The spike trains from different motor units were statistically independent. Because the EMG-production model simulates the randomness inherent in motor-neuron firing, the resulting spike times vary across individual simulation runs, even when all parameters in the EMG model are constant. This randomness in spike timing is a main contributor to the randomness of the synthetic VEMP traces, as well as to the variability of the VEMP averages seen clinically. Because the recordings from our VEMP subjects showed that EMG properties can vary across subjects and sessions (Fig. 3), we adapted the EMG-production model to the characteristics of each subject and session. Although each motor unit of a subject produces a unique pattern at the EMG electrodes, we sought to mimic only the overall behavior of an individual muscle. Therefore, for simplicity, we assumed that for a given subject and session, all motor units in the muscle produced the same surface response waveform shape (Wit & Kingma 2006; Lütkenhöner & Basel 2011) but different amplitudes depending on the motor-unit size and position relative to the electrodes. We refer to the surface response waveform as the MUAP of the subject and session or “MUAPss.” We obtained the relevant EMG characteristics for each subject and session from their EMG autocorrelation functions. The superposition of uncorrelated spike sequences from many motor units has an approximately flat spectrum, so that the shape of the autocorrelation of each EMG waveform must arise largely from the temporal structure of the MUAPss waveform. The EMG model was customized to each subject and session by creating individual MUAPss waveforms that produced synthetic-EMG autocorrelation functions similar to the experimental EMG autocorrelation from that subject and session. For the experimental EMG autocorrelation, we used the EMG obtained at the lowest (103 dB peSPL) stimulus intensity (i.e., the waveform from the whole time period of all the traces from the
VEMP recording at 103 dB peSPL). Although some subjects appeared to show a small VEMP response at 103 dB peSPL, this would have little influence on the EMG autocorrelation because the VEMP response duration is short compared with the overall trace interval, and, as will be seen later, even at the highest sound level, our analysis indicates that the motor units are inhibited only by a small percentage during a VEMP response. For each subject and session, we created a MUAPss waveform that adequately matched the experimentally obtained EMG autocorrelation function. The computation of a waveform from its autocorrelation is an ill-posed inverse problem (detailed for EMGs in Brody et al. 1974). The problem was made tractable by constructing each MUAPss function as a weighted sum of three Gaussian pulses, each with three parameters: width, amplitude, and relative time-shift. With one pulse set to a shift of zero and an amplitude of one, each MUAPss was described by seven free parameters. The pulses were further constrained to have the integral of a MUAPss be zero so the EMG would have zero mean. For each MUAPss, the pulse parameters were determined by an iterative steepest-descent algorithm to minimize the total squared error between the autocorrelations of the MUAPss waveform and the experimental EMG. Using the methods outlined earlier, MUAPss waveforms were computed for each of our 4 × 5 subject/session combinations (Fig. 5). The fitting error (area in the difference between model and experimental autocorrelations) was less than 1.1% of the area under the autocorrelation function. Note that the goal of computing the MUAPss is to match statistics of the simulated EMG to the statistics of the experimental (i.e., actual) EMG. An infinite number of different waveforms can yield the same autocorrelation function (Papoulis & Pillai 2002). Our method computes one choice that matches the autocorrelation function for each MUAPss. It is not expected that this waveform matches the experimental MUAP. The point is to provide a variety of models (20 models) that reflect the experimental EMG variation and to use the models to compare how well the three VEMP metrics track VEMP inhibition depth.
The Motor Drive Each EMG model receives a motor drive that is expressed as a percentage of maximum voluntary contraction. At a motor drive of 100%, all motor units are active at their maximum firing rates. For other motor drives, motor-unit rates were assigned by the recruitment rule. In each EMG-production model, there is a scaling parameter that specifies the amplitude of the MUAPss waveform. The rms value of the synthetic EMG depends on both the scaling parameter and the motor drive. To unambiguously determine the values of these two parameters, the scale parameter was set such that the synthetic EMG at a motor drive of 95% had the same rms value as the largest experimental EMG value recorded during the session. Ninety-five percentage represents an estimate of how close the subject came on their strongest SCM activation to a full motor activation. During the simulations used to compare the three VEMP metrics, the motor drive was varied to mimic the variations in motor drive in the corresponding subject and session. The variation in the motor drive for each run from a subject was obtained from the power spectrum of the overall (40 sec or more) EMG time record. Each run showed an increase in energy at a very-low,
PRAKASH Et al. / EAR & HEARING, VOL. 36, NO. 5, 591–604
Autocorrelation
NORMALIZED AMPLITUDE (sessions offset vertically)
Subject A
Motor Unit Action Potential
Subject B
597
vestibular model generates a signal that modulates the spike rates of the motor units in the EMG-production model (Fig. 4, bottom). The vestibular model produces an inhibition signal that deletes spikes that fall within a time interval T relative to the acoustic stimulus, with a probability of deletion determined by the inhibition depth. The thinned spike sequences are then convolved with the MUAPss to generate each motor unit’s contribution to the surface EMG trace.
Generating Synthetic VEMP Traces Each synthetic VEMP trace was generated by a single execution of the VEMP model consisting of (1) the vestibular and EMG models with specific values of motor-drive and Inhibition Depth parameters, and (2) the EMG-production model with a MUAPss corresponding to one subject and session, all implemented in MATLAB. As in a run in a real VEMP average, a typical ensemble of simulated VEMP traces consisted of a large number, N, of traces (e.g., 200 to 500) obtained by executing the model N times. In each simulation of N traces, the inhibition depth was constant and the motor drive was varied to match the variation of that subject and session. Note that even when the model parameters remained the same, successive executions of the model differed in the probabilistic sequence of motor-unit spikes and in which spikes were inhibited. Because deleted spikes were evenly distributed across interval T, and each synthetic trace was the convolution of a spike sequence and the spike’s surface EMG response in that subject and session (the MUAPss), an approximation of the mean of the traces is the convolution of the MUAPss and a 6-msec rectangular Inhibition Depth change in the motor-unit firing rate. This approximation yields the average VEMP but does not capture the variation across trials, which is an important component of a VEMP measurement, one that influences the accuracy of individual runs.
Subject C
Subject D
Subject E
RESULTS -10
0
10
-10
0
TIME re PEAK (ms)
10
Fig. 5. The electromyographic (EMG) autocorrelations (left) and the corresponding motor-unit action potentials (MUAPss) (right) for each subject and session. At left, each experimental autocorrelation is shown by a solid line, and the synthetic autocorrelation reconstructed from the MUAPss at right is shown as a dashed line offset lower for visibility. The MUAPss reconstruction does not yield a latency, and so the MUAPs at right are centered at zero. The autocorrelations at left, by their nature, are centered at zero. Because there is no useful information in the amplitudes of the autocorrelations or MUAPsss, they are shown on a normalized (dimensionless) scale with the peak values all scaled to the same amplitude.
subject-specific frequency (
