PMCCPMCCPMCC

Search tips
Search criteria 

Advanced

 
Logo of nihpaAbout Author manuscriptsSubmit a manuscriptNIH Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
 
J Clin Neurophysiol. Author manuscript; available in PMC Nov 11, 2010.
Published in final edited form as:
PMCID: PMC2978662
NIHMSID: NIHMS85821
A Model Predicting Optimal Parameters for Deep Brain Stimulation in Essential Tremor
SE Cooper,1 AM Kuncel,2 BR Wolgamuth,1 A Rezai,3 and WM Grill2
1Department of Neurology, Center for Neurological Restoration Cleveland Clinic Foundation, Cleveland, OH
2Department of Biomedical Engineering, Duke University, Durham, NC
3Department of Neurosurgery, Center for Neurological Restoration Cleveland Clinic Foundation, Cleveland, OH
Corresponding Author: Scott E. Cooper PhD, MD, Cleveland Clinic Foundation / S31, 9500 Euclid Ave., Cleveland, OH 44195, Phone: 216.445.4745, Fax: 216.444.9401, Email: coopers2/at/ccf.org
In 9 patients with essential tremor (14 thalami), we varied frequency, voltage, and pulsewidth of thalamic deep brain stimulation (DBS), and quantified postural tremor. Low frequency stimulation aggravated tremor; the effect increased with increasing voltage. High frequency stimulation had a U-shaped relation to voltage, with minimum tremor at an optimal voltage characteristic of the individual thalamus and increases in voltage beyond the optimum reduced tremor suppression.
Based on the hypothesis that tremor response to DBS resulted from two competing processes, we successfully modeled the relationship of tremor to voltage and frequency of stimulation using a mathematical model. The optimum voltage predicted by the model agreed with the empirically measured value. Moreover, the model made accurate predictions at high stimulation frequency based on measurements made at low stimulation frequency.
Our results indicate there is an optimal voltage for tremor suppression by thalamic DBS in most patients with ET. The optimum varies across patients, and this is related to electrode position. A mathematical model based on ‘competing processes’ successfully predicts optimum voltage in individual patients. This supports a competing processes model of DBS effects.
Keywords: deep brain stimulation, movement disorders, essential tremor, surgery, thalamus, thalamic stimulation
Despite the effectiveness of deep brain stimulation (DBS), the mechanisms by which it acts remain subject to debate (Dostrovsky and Lozano, 2002; McIntyre et al., 2004; Benabid et al., 2005). The present study seeks to elucidate DBS mechanisms by quantifying the relationship between stimulation parameters and stimulation effects. We studied DBS of the Ventral Intermediate (Vim) nucleus of the thalamus for treatment of essential tremor (ET) as this afforded a readily quantified symptom promptly responsive to changes in stimulation.
Previously (Kuncel et al., 2006), using general linear regression models and model selection methods, we determined that the stimulation parameters most affecting tremor were voltage and frequency, while pulsewidth was least important. However, the models could not predict the clinical response, because the models were subject-specific: i.e. a different model for each patient. Moreover, the models in our previous analysis were parametric, meaning that the model structure did not attempt to incorporate information about physical processes underlying DBS.
The present paper analyzes the changes in tremor with changes in DBS voltage and frequency using a nonlinear, structural model of Vim DBS. Structural models incorporate information about underlying physical processes. They are more parsimonious than parametric ones and predict more accurately, with fewer model parameters, outside the range of data to which they were originally fit. Structural models offer an opportunity to test hypotheses: If a given hypothesis, incorporated into a structural model, yields more accurate and more parsimonious predictions, then the hypothesis is supported. In this paper we present a nonlinear model based on a “competing processes” hypothesis of the mechanisms of Vim DBS. The model supports the competing processes hypothesis and can predict the clinical response to stimulation.
Details of the experiments have been previously published (Kuncel et al., 2006) and are reviewed briefly below.
1. Subjects
We recruited subjects spanning a wide range of clinical outcomes from patients with essential tremor refractory to medical management who underwent DBS surgery at the Cleveland Clinic. Nine subjects (Table 1) participated. Four were tested on left and right sides, and one before and after his left electrode was replaced in a different location, for a total of 14 thalami. All experiments were approved and monitored by the Cleveland Clinic Institutional Review Board.
Table 1
Table 1
Subject Characteristics
Characteristics of subjects and of individual thalami studied. Subjects are identified by numbers, thalami by letters. Note that, in four subjects two thalami (left and right) were studied. In subject 1, the left thalamus was restudied after the stimulating electrode was revised to a more medial position: the two experiments are designated as thalami B & C. Age the age in years at time of testing. Postop: months, at time of testing, since implantation of the electrode being tested. Contacts: the electrode contacts clinically determined to be most effective, and used in these experiments. The electrode has 4 evenly spaced contacts numbered, from ventral to dorsal 0 to 3, any one or more of which can be either positive or negative during stimulation. In addition, the positive contact can be at a distant low impedance site (the metal case of the implantable pulse generator (IPG), implanted subcutaneously in the chest), effectively delivering monopolar stimulation to the brain. For example, 0−,3+ means contact 0 was negative, contact 3 positive; 1−, C + means contact 1 was negative, the case positive (monopolar stimulation with contact 1). %Improved: % decrease in Fahn-Tolosa-Marin (Fahn et al., 1988) exam for the contralateral arm, preop vs. postop-stim-ON. Global the subject’s global impression of the results of surgery on a scale of poor-fair-good-very good-excellent. This was influenced by other factors besides %Improved (see "Notes" column).
2. Experimental Design
We varied frequency, pulsewidth, and voltage of stimulation and used the stimulating contacts found most effective in the course of clinical adjustments, outside the research protocol*. Varying these three parameters over their full range resulted in about 60 different parameter setting combinations to be tested in each experiment.
In phase I of each experiment, for random points in the frequency-pulsewidth plane, we increased voltage until side effects (usually paresthesias) approached an uncomfortable level. In this way, a restricted set of tolerable test points was generated. Data for analysis were obtained from phase II when the test points established in phase I were presented in random order.
Tremor was quantified using a 3-axis linear accelerometer (Models CLX04LP3 and CLX01LF Crossbow Technology, San Jose CA) on the back of the hand contralateral to the stimulator, while the subject held an empty plastic cup " almost but not quite touching the lips " as if drinking.
3. Data Analysis
Accelerometer signals were digitised (PCI-6025E, National Instruments, Austin, TX) and analyzed in MATLAB (MathWorks, Natick MA). The 3 Cartesian axes were combined as a vector sum; the amplitude of the acceleration vector was then analyzed as a scalar time series in the frequency domain (power spectral density (psd) function, Welch's method, 1-second Hanning window)
We used the total power within a window around the tremor frequency as a measure of tremor severity (e.g. see Fig 1B), in units of decibels, i.e. 10 times the common logarithm of the ratio of tremor power for the test trial to the mean tremor power on control (no stimulation) trials for that thalamus.
Fig 1
Fig 1
Tremor as a function of stimulation voltage of DBS of the Vim thalamus with different combinations of stimulation frequency and pulsewidth. Data from one thalamus.
4. Constructing & Fitting the Model
Details of the model development are given in the Appendix, incluing rationale for the assumptions incorporated into to the model. Briefly, our data (see Results) suggested that total tremor reflected the net effect of a tremor-suppressing process plus a tremor-aggravating process, and that the effect of pulsewidth was small, compared to voltage and frequency, in keeping with Kuncel et al., 2006. We therefore represented tremor as the sum of two terms, representing the two processes. This resulted in the following equation:
equation M1
Where T is tremor, v is stimulation voltage and f is stimulation frequency. The remaining quantities are model parameters fit to the data.
Model fitting was done in the R statistical programming language (R Development Core Team, 2003). Control trials (trials without stimulation) were considered to be zero volts stimulation amplitude and zero Hz stimulation frequency. We linearly scaled the tremor from each thalamus to lie in the range ± 1 before fitting the model, in order to normalize the dynamic range of tremor, preserve the shape of the response surface, as well as any inter-patient or inter-thalamus variation in that shape, while removing inter-patient and inter-thalamus variation in absolute tremor severity in the absence of stimulation. The tremor data for each thalamus were scaled according to the formula Ti’=Ti/max(abs(Ti)), where Ti is a vector of all tremor measurements in thalamus i in dB, Ti’ is the same vector in normalized log units, abs(v) is a vector of the absolute values of vector v, max(v) is the largest element in vector v, and “/” is elementwise division of a vector by a scalar. This transformation preserves zero as the value of tremor in the absence of stimulation.
Model parameters were fit to the data using a two-stage algorithm. Global parameters (Table 2) were those constrained to have the same value for all thalami and were fit to the entire dataset using a Gauss-Newton multidimensional iterative gradient descent algorithm (R’s nlm() function). The calculation of mean squared error for this algorithm incorporated an internal Gauss-Newton fit of the single "local" parameter (Table 3) which was allowed to vary across thalami (R's nls() function).
Table 2
Table 2
Global model parameters
Table 3
Table 3
Local model parameter
1. Frequency-Voltage Response Surfaces
Tremor from a typical thalamus is plotted vs. voltage in Fig 1A with each frequency-pulsewidth combination represented by a different letter. The points fell into two clusters, dependent on stimulation frequency, but not pulsewidth. For low stimulation frequencies (lower case letters), tremor increased with increasing voltage. In contrast, at high stimulation frequencies, tremor exhibited an asymmetrical U-shaped relation to voltage, decreasing up to an optimal voltage (here, about 2 volts) at which tremor was minimal, and then increasing again at higher voltages. The U-shaped relation to voltage was seen in all thalami, with statistically significant separation between low- and high-frequency tremor-voltage curves, as we have previously reported (Kuncel et al., 2006). The lack of apparent effect of pulsewidth is also consistent with our earlier report. Our initial analysis focused on voltage and frequency alone, and the effect of pulsewidth is analyzed further below (Results, Part 3).
Data from two more individual thalami are plotted in Fig. 2A,B. These response surfaces can be viewed as topographic maps; latitude and longitude are frequency and voltage, while the height of the terrain is the amplitude of tremor. There was a "valley" (minimum) at upper left, and a "peak" (maximum) at lower right (i.e. the response surface was a double rising ridge system (Box and Draper, 1987).) The zero contour ("sea level") corresponds to baseline tremor, in the absence of stimulation. The minimum at upper left represented a therapeutic effect (reduced tremor at high stimulation frequency), the maximum, at lower right represented worsening of tremor at low stimulation frequency. The bottom of the "valley" (minimum) was an optimum, in the sense that tremor was minimal (therapeutic effect maximal) there.
Fig 2
Fig 2
Contour plots of tremor as a function of stimulation voltage and stimulation frequency of DBS of the Vim thalamus in two subjects with essential tremor. Part A,B: experimental data; Part C,D: model predictions with values of β1 taken from Table (more ...)
All thalami had response surfaces which were qualitatively similar, with a long, narrow valley at high frequency-moderate voltage, and a peak at low frequency-high voltage. However, there were quantitative differences among thalami. Some had a tremor minimum which was deep, with a minimum tremor at higher voltage (the optimal voltage); in those thalami, the tremor maximum was smaller. Other thalami had a shallow tremor minimum at lower optimal voltage. In those thalami, the tremor maximum was larger and tremor at high-frequency high-voltage stimulation could exceed baseline (unstimulated) levels. All 14 thalami appeared to fall at different points along that spectrum. Thus the amount of tremor suppression (depth of the tremor minimum) was inversely related to the amount of tremor aggravation at low stimulation frequency and both quantities were related to the optimal voltage.
2. Model of Frequency-Voltage Response Surfaces
The nonlinear competing processes model reproduced the variation of tremor responses across thalami by varying one of the parameters, β1 (Fig 2 C,D). This was not true of any of the other parameters. Low values of β1 corresponded to "good-responders" like thalamus D in Fig 2 with a deep optimum at a high stimulation voltage; high values of β1 corresponded to "poor-responders" like thalamus E in Fig 2, with a shallower minimum at lower stimulation voltage.
We fit the model with β1 allowed to vary between thalami (‘local’ parameter), while the other parameters were constrained to be the same for the entire dataset (‘global’ parameters). With this model, we were able to account for 46% of the variance in our data. We compared our nonlinear model to a linear (multiple polynomial regression) model similar to the ones in our earlier paper (Kuncel et al., 2006). The comparison model included frequency, frequency2, pulsewidth, pulsewidth2, voltage, voltage2, all main effects, and all interactions. The resulting model has a very large number of parameters (64 in all) and would rarely be used in practice. It is probably possible to construct a reduced linear model with fewer terms and similar goodness of fit, but the full model represents a highly conservative basis for comparison: with so many parameters it should to fit the data better than our nonlinear model, which has fewer parameters unless the structure of the latter corresponded better with the underlying physical processes.
We found that the linear parametric comparison model fit our data less well than our nonlinear structural model (34% as compared to 46% variance accounted for), despite having more free parameters. This supports our hypothesis that the nonlinear model captures something true about the physical processes underlying Vim DBS. Additional model validation appears in Results: section 5.
"Global" model parameters (those which are the same for all thalami), with 95% confidence intervals estimated by bootstrap resampling (948 samples drawn with replacement from 948; 10 replicates). For units, see Appendix. Values of the "local" model parameter β1 are given in Table 3
Values of model parameter β1, which varied between individual thalami. β1 has units of V−1 (see Appendix). For a comparison of tremor with different values of β1, see Fig 2.
3. Effect of Pulsewidth
Our earlier study (Kuncel et al., 2006) and the graphical analysis (see Results: Part 1, above) suggested that pulsewidth exerted relatively little influence over tremor. We re-asessed this, using the nonlinear model, by subtracting predicted from observed tremor to obtain residual variance not accounted for by frequency and voltage.
In some thalami, it appeared (Fig 3) that the model underestimated tremor aggravation (negative residual) at low frequency/short pulsewidth and underestimated tremor suppression (positive residual) at high frequency/short pulsewidth. This suggested that short pulsewidths were less effective than long pulsewidths in both aggravating and suppressing tremor. To assess this, we constructed a linear mixed model with the residual variance from the main model as the dependent variable, frequency & pulsewidth as regressors, and thalamus as the random factor. As expected, main effects were nonsignificant. However, the pulsewidth × frequency interaction term was borderline significant (p=0.043), supporting a weak dependence of stimulation effect on pulsewidth.
Fig 3
Fig 3
Contour plots of residual tremor not accounted for by the model as a function of stimulation frequency and stimulation pulsewidth in two thalami.
4. Effect of Electrode Location
We found that much of the variation in response surfaces between thalami was captured by the single parameter β1. In one subject whose electrode was surgically revised, changing the location of the electrode changed the value of β1. This suggested that β1 reflected electrode location.
To examine this, we computed the location of the active (negative) electrode contact used in each experiment from stereotaxic coordinates recorded at the time of surgery. (When there were 2 adjacent negative contacts, the location was taken as their midpoint.) We computed Z (dorsoventral) and Y (anteroposterior) coordinates relative to posterior commissure (PC), and the X (mediolateral) coordinate relative to the lateral ventricle wall. β1 was weakly corellated with X and Y, with lower (better clinical response) values located at more medial and more anterior locations; however, only the X (mediolateral) effect achieved statistical significance (multiple linear regression, p=0.03). This was true regardless of whether the absolute coordinates were used or the coordinates scaled by AC-PC distance.
5. Prediction of Optimal Stimulation Voltage
We tested our model by predicting a clinically meaningful quantity, namely optimum voltage at which tremor was minimized (Vopt). For each thalamus, we chose, for determination of Vopt, a frequency greater than 100 that was tested with a wide range of voltages. In 11 thalami, we were able to use 160 or 170 Hz. In one thalamus, we used 130 Hz. In 2 thalami the available range of voltages was too restricted by side effects for an unambiguous Vopt to be determined. Fig 4A shows the positive correlation between observed Vopt and the model prediction (linear regression, R2=0.59 slope significantly different from zero at p=0.0036). This supports the validity of the model.
Fig 4
Fig 4
Correlation between the optimal stimulation voltage predicated by the competing processes model and the optimal stimulation voltage determined from the experimental data for 12 of 14 thalami. Two thalami were excluded because an optimum voltage could (more ...)
As a more rigorous test, we split the dataset into two parts, fit the model to one part, and validated it against the other. Specifically, we fit β1 for each thalamus using only test points with frequency ≤ 100 Hz (Fig 4B). Using those values, predicted Vopt was again well correlated with the actual values measured from the data at frequencies > 100 Hz (linear regression, R2=0.65, slope significantly different from zero at p=0.0016).
In this study we modeled the quantitative relationship between Vim DBS stimulation parameters and the tremor response in essential tremor. Unlike our earlier patient-specific models (Kuncel et al., 2006), this model can predict the tremor response to varying stimulation parameters: It incorporated a “competing processes” hypothesis, and made more accurate predictions, with fewer free model parameters, than a conventional model lacking that hypothesis. Our model accounted for 46% of the variance in the data, which is far from perfect: however our goal was not perfect prediction but, rather, to test the “competing processes” hypothesis. Unlike other models (Ushe et al., 2004), our model accounts for both frequency and voltage-dependence, as well as inter-thalamus variation related to electrode position via the thalamus-specific parameter β1.
1. Competing Processes Model of DBS Effects
We found that, for therapeutically effective stimulation frequencies, the relation of tremor to stimulation voltage was U-shaped with an optimum voltage that differed between subjects and thalami. As voltage increased above the optimum, tremor suppression was decreased, and at higher voltages, stimulation aggravated tremor.
This is easily overlooked in clinical observations limited to patients with well-placed electrodes and moderate stimulation voltages. We detected it because we applied a quantitative analysis to voltage varied across a wide range, in patients with a wide range of electrode placements. We did not detect any sharp division into categories of well- and poorly-placed electrodes. Rather, the U-shaped tremor-voltage relationship was seen in all patients, with only quantitative differences according to electrode location.
We also found, with others (Ohye and Narabayashi, 1979; Bejjani et al., 2000) that tremor aggravation occurred not only with high voltage, but also with low frequency stimulation.
We propose a ‘competing processes’ model of Vim DBS to account for tremor suppression and both types of tremor-aggravation. Under this model, DBS exerts a tremor-suppressing effect due to activation of a region near the electrode (presumably Vim nucleus, or, possibly, Vop or prelemniscal radiations), and a tremor-aggravating effect due to current spread to a more distant site (perhaps internal capsule). At voltages less than the optimal voltage Vopt, increasing voltage activates a larger volume of the tremor-suppressing region. At voltages greater than Vopt, increasing voltage results in current spread to the more distant tremor-aggravating region. The two processes have qualitatively similar dependence on frequency, but quantitative differences such that, at low frequency, tremor-aggravation dominates.
The competing processes model arose as a way to account for our observations while being economical in the number of assumptions; in fact, the model makes only one novel assumption (see Appendix), namely, the notion that Vim DBS has simultaneous competing effects on tremor. We tested that hypothesis by incorporating it into our mathematical model of Vim DBS effects on tremor, then tested the model in two ways. First, we demonstrated that our model predicted Vim DBS effects on tremor more accurately, with fewer free parameters, than a multiple polynomial regression model. The point is not that our model was more accurate than the multiple regression model (although it was). The point is that building certain assumptions into our model increased its accuracy, without sacrificing parsimony, implying that the assumptions were correct. Second, we validated the model by fitting it to a subset of the data (low-frequency stimulation) and demonstrating that it could then accurately predict a separate subset of the data (high-frequency stimulation). The success of the model by both tests supports the competing processes hypothesis on which it was built.
The physical mechanisms of tremor-aggravation and -suppression is not directly addressed by our data. We have elsewhere proposed a mechanism for frequency-dependence of tremor-suppression (Grill et al., 2004). One possible physical mechanism by which tremor-aggravation might occur is by producing uncomfortable side effects, such as paresthesias. Discomfort would then aggravate tremor because it was stressful or anxiogenic. In its simplest form, this hypothesis is incompatible with our data, because tremor aggravation was prominent at low stimulation frequencies, which did not produce uncomfortable paresthesias. This hypothesis could still apply, however, to high-frequency / high voltage stimulation, if we postulate a second, separate mechanism of tremor aggravation for low-frequency stimulation. When tremor increased with low-frequency stimulation, in agreement with others (Ohye and Narabayashi, 1979; Bejjani et al., 2000) we observed a jerky quality apparently different from “natural” tremor,. However, we saw the same with tremor-aggravation with high-frequency high-voltage stimulation; moreover we found that with stimulation-aggravated tremor, the power spectral peak did not broaden or shift (Fig 1B) indicating that the movement was still a rhythmical oscillation at the same frequency. A separate mechanism of tremor aggravation by low-stimulation frequency could be incorporated into our model by splitting the tremor-aggravation term into two terms, one operative at low stimulation frequency, and the other at high frequency. At present, however, this seems a needless complication because the increased model complexity is not required to fit the data.
Another possible mechanism of tremor aggravation would be ‘driving’ the essential tremor oscillator. The amplitude of an oscillation can be increased by supplying energy to the oscillator at a susceptible frequency. For a linear oscillator, the forcing function must oscillate at the natural frequency of the oscillator, or a harmonic, or subharmonic thereof: for example, 2 Hz stimulation might aggravate a 4 Hz tremor, and 5 Hz stimulation a 5 Hz tremor, but not vice versa. In its simplest form, this hypothesis is incompatible with our data and that of Ushe et al., 2004 because the relationship of tremor to stimulation frequency is smooth, without sharp peaks and valleys corresponding to resonant frequencies. However, a nonlinear oscillator may exhibit more complex behavior: it is possible that a nonlinear essential tremor oscillator could exhibit a susceptibility to driving corresponding to the smooth frequency-dependence in our model (Fig 5 A).
Fig 5
Fig 5
Tremor aggravation and suppression components of the model are dependent on stimulation frequency and voltage. Parameters used to generate these plots were taken from Table 2 and Table 3 (thalamus D).
We did not investigate the effects of stimulation on tremor frequency. Inspection of individual power spectra did not reveal a frequency difference between spontaneous and stimulation-aggravated tremor (e.g. Fig 1B). However, our experiments were not designed to detect such an effect. Experiments concentrating on tremor-frequency effects might shed light on the mechanism by which activation of the tremor-aggravating region aggravates tremor.
2. Clinical Aspects
The present study was primarily intended to elucidate DBS mechanisms, rather than producing data of direct utility to clinicians. Our model could, in principle, be used to aid clinical stimulator adjustments: clinicians could test the effects of low-frequency stimulation, and use the results to predict the optimal voltage for high frequency stimulation, thereby avoiding extensive testing at high frequency. This might be advantageous in some patients, since low frequency stimulation is less apt to cause uncomfortable paresthesias. In practice, however, clinical Vim DBS adjustment for ET is rarely difficult enough to justify the weight of mathematical apparatus developed here.
We suggest that in practice the most important clinical implications derive from the qualitative shape of the response surface. We found that the frequency-voltage response surface was a ridge system (Box and Draper, 1987) oriented approximately parallel to the frequency axis. This means that, in the part of the frequency-voltage regime where net tremor-suppression takes place, variations in voltage produced a much larger change in tremor suppression than variations in frequency. This suggests that adjustments should begin at a fixed frequency, with more attention to fine-tuning voltage than frequency. If tremor suppression at the optimal voltage is inadequate, further voltage increase is likely to be counterproductive. Increasing frequency may gain additional benefit, although we (Kuncel et al., 2006) and others (Ushe et al., 2004)(O'Suilleabhain et al., 2003) have not found this to yield large improvements in tremor. Increased pulsewidth might also gain additional benefit, but we found the effect of pulsewidth to be quite small, in agreement with O'Suilleabhain et al., 2003 and our own earlier work (Kuncel et al., 2006).
This study examined a population of tremor patients that was heterogeneous with respect to electrode placement and clinical outcome. The findings of tremor aggravation by low-frequency stimulation and a U-shaped relationship between tremor and voltage with high-frequency stimulation were consistent across all patients. Additionally, we detected a correlation of a patient-specific characteristic (electrode location) with a model parameter (beta1). The ability to detect correlation is greater when there is more variance in the independent variable, so the clinical heterogeneity helped us detect the correlation, and restricting the population to patients with optimal electrode placement would reduce the strength of the correlation. For the clinician faced with a more heterogeneous patient population, the correlation may be useful for clinical decision-making, e.g. knowing how much improvement in symptoms can be achieved by optimizing stimulation parameters helps to decide whether to reposition the electrode.
3. Limitations of the Model
Our model captures much of the variation in our data in a mathematically tractable way. However, it has limitations. First it is still not a complete model of physical processes. The values of parameters have been fit to achieve a match to the experimental data, not to directly measured physical quantities like distance, electrical impedance, etc. Nonetheless, the present model is a step towards the goal of a complete structural model. For example, we believe a physical model will have to incorporate two different sites where stimulation acts, at different distances from the electrode.
Second, our model does not address effects of varying the active contact on the electrode used for stimulation. We have conducted a separate series of experiments in which we varied this parameter (Kuncel et al., MS in preparation), however, in the experiements reported here, in order to keep the number of stimulator setting combinations to be tested within manageable limits, we conducted all experiements using the electrode contacts previously established, by the treating physician, to be the most effective.
Third, in order to predict Vopt for an individual thalamus, the model needs an estimate of β1 for that thalamus. This requires, at a minimum, tremor measurements made at low stimulation frequency (though in future, it might be possible to estimate β1 from information about the electrode location). This can be regarded as either a defect or a virtue of the model: On the one hand, for maximum accuracy, one (though only one) of the model parameters needs to be fit individually. On the other hand, inter-thalamus variation in the response to stimulation is a real (and clinically important) phenomenon, and any model that fails to take it into account is, in some sense, incomplete.
Finally, our model accounts for less than 100% of the variance in our data. Some of the remaining variance may be truly random, but some may reflect as yet unidentified factors such as, clinical & demographic variables. We believe the present model will aid in uncovering such factors. By modeling, predicting, and then subtracting out dominant effects, subtler effects of other factors may be uncovered.
ACKNOWLEGEMENTS
Supported by NIH R01 NS-40894. We are indebted to Dr. Erwin Montgomery for advice and encouragement.
APPENDIX: MATHEMATICAL MODEL OF DBS EFFECTS
1 ASSUMPTIONS OF THE MODEL
1.1 Novel Assumptions
The model makes only one truly novel assumption, namely that thalamic DBS produces not one, but two competing effects, one suppressing, and the other aggravating tremor. In a sense, this is not a novel assumption at all, since both tremor suppression and tremor aggravation (Ohye and Narabayashi, 1979; Bejjani et al., 2000) have been described with thalamic DBS. However, previously, it was assumed that stimulation caused either tremor suppression or tremor aggravation, according to whether stimulation frequency was high or low.
Because we observed tremor aggravation also with high frequency stimulation as well (Kuncel et al., 2006) such that the same stimulation frequency caused either tremor suppression or tremor aggravation depending on stimulation voltage, we postulated that stimulation always causes both effects, regardless of frequency, and that either an increase or decrease of tremor was seen according to which effect predominated. We supposed that the reason tremor suppresion gave way to to tremor aggravation as voltage increased was that current spread extended beyond a local tremor-suppressing region, and began to affect a more distant tremor-aggravating region.
1.2 Commonplace Assumptions
In order to translate this hypothesis into a mathematical model, we needed to make several other more straightforward, less novel assumptions.
1.2.1
First, we assumed that both the tremor-aggravating and tremor-suppressing regions responded to stimulation amplitude and stimulation frequency in the same way. We did not postulate a special, adhoc property of responding more to low- than to high-frequency stimulation for one of the two regions based on hypothetical resonance effects of neuronal oscillators. Rather, we assumed that both regions behaved like most nervous tissue, i.e. the effect of stimulation increased with increasing stimulation voltage and increased with increasing stimulation frequency.
1.2.2
Next we assumed that the effect of stimulation was zero at zero volts stimulation amplitude. This should require no justification.
1.2.3
Next, we assumed that stimulation effect "saturated," or plateaued at high voltage. It should go without saying that stimulation effect plateaus at some voltage. Stimulation effect at 6V is certainly different from 1V, and effect at 16V may be different from 11V, but it is hardly likely that 56V differs from 51V by as much as 6V from 1V. The only real question is: at what voltage does the effect plateau – and the model makes no assumptions about this (it is determined when the model is fit to the data).
1.2.4
Next we assumed that stimulation effect "saturated," or plateaued at high frequency. The argument here is the same as for saturation at high voltage.
1.2.5
Finally, we assumed that stimulation effect also saturated at low frequency. This is justified by the fact that frequency is the reciprocal of period. The effect of 100 Hz stimulation (0.01 sec between pulses) certainly differs from that of 10 Hz stimulation (0.1 second between pulses), which may well differ from 1 Hz stimulation (1.0 second between pulses). It is hardly likely, though, that the effect of one pulse per second differs as much from one pulse every 10 seconds (0.1 Hz), still less that one pulse per 10 seconds differs so much from one pulse per hour, and so on. If, as we suppose, stimulation effect decreases monotonicly with stimulation frequency (1.2.1) and reaches zero at 0 Hz (1.2.3) then it is only reasonable to suppose that it approaches zero asymptoticly as stimulation frequency decreases.
2 MATHEMATICAL FORM OF THE MODEL
Assumption 1.1 gives:
equation M2
Equation 1
Where
θ = tremor aggravation (dimensionless),
ϕ = tremor suppression (dimensionless).
T = net tremor (dB)
a,b = parameters reflecting relative strength of tremor aggravation and tremor suppression (dB).
Assumptions 1.2.1 and 1.2.2 give:
equation M3
Equation 2
(See Fig 5A,B)
Where
v = stimulation voltage (V).
ψ = a quantity determining the steepness of the curve (V−1).
Substituting x for θ and ϕ in equation 1 gives:
equation M4
Equation 3
(See Fig 5C)
Dependence on stimulation frequency, satisfying assumptions 1.2.1, 1.2.3, 1.2.5, & 1.2.6 is incorporated by:
equation M5
Equation 4
(See Fig 5D,E)
Where
f is stimulation frequency (Hz)
α is the asymptotic maximum of ψ at high f
β is the asymptotic minimum of ψ at low f
w is the frequency at which ψ is midway between α and β (Hz)
s determines the "sharpness" of the sigmoidal curve (dimensionless). For low values, ψ is "spread out" along the frequency axis; for high values, ψ approximates a step function.
Substituting for ψ in eq 3 gives:
equation M6
Equation 5
Based on assumption 1.2.1, we kept the model as parsimonious as possible by assigning the same value to equivalent global parameters in the tremor-aggravating and tremor-suppressing terms, where possible. Thus, s1=s2 & w1=w2. Allowing these parameters to vary independently in the two terms did not improve the fit to the data. Obviously, not all the parameters can be the same in both terms (otherwise the predicted T would always be zero). The structure of the model required the remaining parameters (a,b,α1 & α2, β1 & β2) to vary independently in order to reproduce the experimental observations.
The complete model is thus:
equation M7
Equation 6
Footnotes
*The algorithm in use at that time for thalamic stimulation at our institution was to start with bipolar 0− 3+. For inadequate therepeutic effect, despite adjustment of voltage, frequency, and pulsewitdth a second cathode was added, and then, if necessary, the geometry was changed to monopolar (with one cathode), and then, if necessary, a second monopolar cathode was added. To reduce side effects (usually paresthesias), the distance between the positive and negative contacts was reduced and the most distal negative contact was moved up (e.g., from 0− 3+ to 1−3+).
  • Bejjani BP, Arnulf I, Vidailhet M, Pidoux B, Damier P, Papadopoulous S, Bonnet A-M, Cornu P, Dormont D, Agid Y. Irregular jerky tremor, myoclonus, and thalamus: A study using low-frequency stimulation. Movement Disorders. 2000;15(5):919–924. [PubMed]
  • Benabid AL, Wallace B, Mitrofanis J, Xia R, Piallat B, Chabardes S, Berger F. A putative generalized model of the effects and mechanism of action of high frequency electrical stimulation of the central nervous system. Acta Neurol Belg. 2005;105(3):149–157. [PubMed]
  • Box GEP, Draper NR. Empirical Model-Building and Response Surfaces. New York: John Wiley and Sons; 1987.
  • Cleveland WS. Robust locally weighted regression and smoothing scatterplots. Journal of the American Statistical Association. 1979;74:829–836.
  • Dostrovsky JO, Lozano AM. Mechanisms of deep brain stimulation. Mov Disord. 2002;17(Suppl 3):S63–S68. [PubMed]
  • Fahn S, Tolosa E, Marín C. Clinical rating scale for tremor. In: Jankovic J, Tolosa E, editors. Parkinson's Disease and Movement Disorders. Baltimore: Urban & Schwarzenberg; 1988. pp. 225–234.
  • Grill WM, Snyder AN, Miocinovic S. Deep brain stimulation creates an informational lesion of the stimulated nucleus. Neuroreport. 2004;15(7):1137–1140. [PubMed]
  • Kuncel AM, Cooper SE, Wolgamuth BR, Clyde MA, Snyder SA, Montgomery EBJ, Rezai AR, Grill WM. Clinical response to varying the stimulus parameters in deep brain stimulation for essential tremor. Mov Disord. 2006;21(11):1920–1928. [PubMed]
  • McIntyre CC, Grill WM, Sherman DL, Thakor NV. Cellular effects of deep brain stimulation: model-based analysis of activation and inhibition. Journal of Neurophysiology. 2004;91(4):1457–1469. [PubMed]
  • O'Suilleabhain PE, Frawley W, Giller C, RBD Tremor response to polarity, voltage, pulsewidth and frequency of thalamic stimulation. Neurology. 2003;60:786–790. [PubMed]
  • Ohye C, Narabayashi H. Physiological study of presumed ventralis intermedius neurons in the human thalamus. Journal of Neurosurgery. 1979;50(3):290–297. [PubMed]
  • R Development Core Team, The. R: A language and environment for statistical computing. Vienna: R Foundation for Statistical Computing; 2003.
  • Ushe M, Mink JW, Revilla FJ, Wernle A, Schneider Gibson P, McGee-Minnich L, Hong M, Rich KM, Lyons KE, Pahwa R, Perlmutter JS. Effect of stimulation frequency on tremor suppression in essential tremor. Mov Disord. 2004;19(10):1163–1168. [PubMed]