Identifying controllable cortical neural markers with machine learning for adaptive deep brain stimulation in Parkinson's disease.

The identification of oscillatory neural markers of Parkinson's disease (PD) can contribute not only to the understanding of functional mechanisms of the disorder, but may also serve in adaptive deep brain stimulation (DBS) systems. These systems seek online adaptation of stimulation parameters in closed-loop as a function of neural markers, aiming at improving treatment's efficacy and reducing side effects. Typically, the identification of PD neural markers is based on group-level studies. Due to the heterogeneity of symptoms across patients, however, such group-level neural markers, like the beta band power of the subthalamic nucleus, are not present in every patient or not informative about every patient's motor state. Instead, individual neural markers may be preferable for providing a personalized solution for the adaptation of stimulation parameters. Fortunately, data-driven bottom-up approaches based on machine learning may be utilized. These approaches have been developed and applied successfully in the field of brain-computer interfaces with the goal of providing individuals with means of communication and control. In our contribution, we present results obtained with a novel supervised data-driven identification of neural markers of hand motor performance based on a supervised machine learning model. Data of 16 experimental sessions obtained from seven PD patients undergoing DBS therapy show that the supervised patient-specific neural markers provide improved decoding accuracy of hand motor performance, compared to group-level neural markers reported in the literature. We observed that the individual markers are sensitive to DBS therapy and thus, may represent controllable variables in an adaptive DBS system.

Introduction

Deep brain stimulation (DBS) has developed into a standard therapy for treating refractory stages of Parkinson’s disease (PD). The large number of DBS systems routinely implanted nowadays are relatively simple from a technical perspective: they uninterruptedly deliver high-frequency stimulation pulse trains 24 h a day. Stimulation is applied to the target area—e.g., the subthalamic nucleus (STN), or the internal globus pallidus—without taking into account the current motor state of the patient, sleep-/wake-cycles, or other contextual information. Changes to the stimulation parameters—such as pulse width, amplitude, frequency, or the choice of electrode contacts—can be made only by a trained expert. This limits the number of adjustment sessions to a few per year. While this fitting strategy may be sufficient for adapting the system to long-term variations of PD-induced changes in the patient’s state—which take place over months and years—it is not sufficient to react upon varying daily conditions or changes in even smaller temporal scales (Carron et al., 2013, Little et al., 2013). Furthermore, it has been shown that patients undergoing such continuous DBS therapy are prone to motor and neuropsychiatric side effects (Witt et al., 2012, Little et al., 2016, Castrioto et al., 2014).

Adaptive DBS systems provide alternatives to continuous DBS therapy by allowing online adaptation of stimulation parameters as a function of motor state surrogates. An adaptive system seeks to identify pathological states in real time and to adjust the stimulation strategy accordingly, in order to improve treatment and ameliorate DBS-induced side effects (Khobragade et al., 2015, Little et al., 2013, Little et al., 2016, Tinkhauser et al., 2017, Priori et al., 2013). The motor state is accessible via external sensors which measure, for example, muscle tone or tremor (Herron et al., 2016, Graupe et al., 2010). Alternatively, surrogates of the motor state may be provided by neural markers (NMs) (Carron et al., 2013, Priori et al., 2013, Beudel and Brown, 2016). It has been proposed to extract NMs from local field potentials (LFPs) of subthalamic target structures like the STN (Little et al., 2013, Priori et al., 2013, Whitmer et al., 2012). The majority of these NMs are bandpower features derived from oscillatory processes in the beta-band ( Hz), showing a significant correlation with common PD motor symptoms, such as bradykinesia and rigidity (Kühn et al., 2009, Kühn et al., 2008, Blumenfeld and Brontë-Stewart, Dec. 2015, Whitmer et al., 2012, Neumann et al., 2017).

Due to the strong functional coupling of cortical structures with subcortical ganglia, the former may also provide NMs (Deniau et al., 2010). In addition to providing information about broader brain areas, cortical signals are further away from the stimulation contacts and thus, are less affected by DBS-induced artifacts (Neumann et al., 2019). In this regard, several electrocorticography studies have found that alpha-, beta-, and gamma-band oscillations may serve as NMs (Swann et al., 2018, Kondylis et al., 2016). Likewise, non-invasive recordings may provide similar information, as reported using magnetoencephalographic (Airaksinen et al., 2012, Cao et al., 2017) and electroencephalographic (EEG) signals (Soikkeli et al., 1991, Weiss et al., 2015). In this regard, studies have found a negative correlation between beta-band and alpha-band power with respect to clinical scores (Melgari et al., 2014, He et al., 2017), such as the unified PD rating scale (UPDRS) (Ramaker et al., Sept. 2002, Goetz et al., 2008), whereas other studies have reported a positive correlation between theta-band power and Hoehn & Yahr stage (Goetz et al., 2004, He et al., 2017). Furthermore, effects of levodopa and DBS upon theta-, alpha-, beta-, and mu-band activity have also been observed (Brown and Marsden, 1999, Melgari et al., 2014, Jech et al., 2006, Whitmer et al., 2012). The cortical and subcortical NMs reported in the above studies stem from group-level analyses, which largely make use of physiologically motivated hypotheses of the underlying PD mechanisms. Such group-level NMs may have limited use in adaptive DBS systems because of the intrinsic heterogeneity of PD and its sub-classes (van Rooden et al., 2010, van Rooden et al., 2011), which indicates that one-NM-fits-all may not exist.

In contrast, we suggest to utilize a bottom-up approach pursuing personalized NMs based on brain signals of an individual. A similar approach to achieve this has been recently proposed by Yao and colleagues (Yao et al., 2020), where a machine learning model was used to decode Parkinsonian tremor from STN-LFP signals. For high dimensional signals, as those typically measured from the cortex by means of EEG, the machine learning methods can be borrowed from the field of brain-computer interfaces. Such algorithms are capable of decoding intentions from brain signals for the purpose of controlling assistive devices or communication (Blankertz et al., 2011, Tangermann et al., 2012) and monitoring ongoing brain states, e.g, for forecasting hand motor performance (Meinel et al., 2016). As suggested in the outlook papers by Neumann and colleagues (Neumann et al., 2019 and Hell and colleagues Hell et al., 2019), an adaptive system based on NMs extracted via machine learning methods, offer a highly personalized treatment that could significantly improve the therapy. While such adaptive DBS systems can not be expected to generalize well across patients, an adaptive DBS treatment should not entail any disadvantages for the individual patient, and is compatible with standard clinical practices that pursue individual adaptation of stimulation parameters according to the individual semiology of the patients (Wagle Shukla et al., 2017, 2017.).

In this regard, we present a novel approach for a supervised, machine learning based identification of cortical NMs for hand motor performance in PD patients, and offer general guidelines for the utilization of this approach in a prospective adaptive DBS system. Even though our contribution is focused on the identification of NMs, it is framed around the characterization of the interactions between the three main phenomena relevant to an adaptive system, i.e., the motor state, its NMs, and the DBS-induced effects upon such NMs. Each of these interactions is individually addressed in three separate tasks:

Task T1 – Extraction of DBS-induced changes of hand motor performance: We determine the multidimensional behavioral characteristics of hand movement, as captured during the execution of a fine-motor control task termed the copy-draw test (Castaño-Candamil et al., 2019). Additionally, task T1 provides a motor score that condensates such information in a single uni-dimensional variable.

Task T2 – Decoding of motor performance from brain signals: The identification of NMs in itself is based on machine learning methods designed to extract neural sources contained in EEG recordings, which maximally correlate with the copy-draw test performance identified in task T1.

Task T3 – Decoding of DBS condition from brain signals: This auxiliary task makes sure that NMs identified in T2 not only correlate with motor performance but are also controllable by DBS—a prerequisite for using those NMs to control an adaptive DBS system. This is achieved by identifying those NMs from task T2 which contain discriminative information about the DBS condition.

Methods

The copy-draw test – a fast hand motor assessment tool

In a recent study Castaño-Candamil et al. (2019), we have presented the copy-draw test, inspired by Prichard et al. (2014). It allows fast extraction of features related to complex and continuous hand and arm movements. The experimental setup of our study is sketched in Fig. 1. At the beginning of a copy-draw trial, a trace composed of pseudo-letters is displayed on screen. The patient’s task is to copy-draw the trace trajectory “as fast and as accurately as possible” within a predefined time limit.

Fig. 1

Left: Illustrative EEG layout used in sessions . The EEG electrodes avoided incision sites used for DBS electrode implantation. In sessions , a full EEG layout was used. Right: The experimental setup and design of the copy-draw test, as composed by (a) single trial trace consisting of three (b) trace atoms, (c) get-ready box, and (d) starting point of the trace.

We extracted two categories of features from every trace drawn: directed and undirected features. For each category, kinematic and precision features were computed: The kinematic features describe the average velocity, acceleration, and jerk per trial. The precision feature is obtained as a point-wise Euclidean distance between the trace drawn and the template trace. The difference between directed and undirected features is that for the former, each point of the trace is put into one of eight 45° bins according to the instantaneous direction in the drawing plane. In total, four undirected features and 32 directed features characterize each trial.

A scalar motor performance score for the e-th trial, , should now be obtained from the concatenation of all features described above. We assume that the highest variability of performance in the copy-draw test is induced by switching between DBS-on and DBS-off. Then, a dimensionality reducing projection can be determined by linear discriminant analysis (LDA), as , where is a vector containing the movement features for the e-th trial. Using information about the DBS condition as 0/1 class labels, projection vector is obtained by the shrinkage-regularized linear discriminant analysis as described in Blankertz et al. (2011). The task T1, i.e., the extraction of motor performance under DBS-on/-off, is realized by this projection.

Data-driven identification of controllable neural markers

EEG can contain oscillatory NMs hidden by noise and overlapping confounds. In the following, we describe a machine learning framework typically used in brain-computer interfaces, which allows to recover such oscillatory NMs on a single trial basis.

Let be the EEG data segmented into a window of interest e, corresponding to the -th trial of the copy-draw test, recorded at time samples with sensors. The power of an underlying narrow-band neural source, acting as motor performance surrogate, can be recovered as , where is a so-called spatial filter. To obtain , consider the optimization problemwhere matrices define some desired properties of the underlying subspace spanned by . It can be shown that a solution to this optimization problem can be obtained by solving the generalized eigenvalue problem:where is the corresponding eigenvalue.

In order to solve tasks T2 and T3, we select matrices and , such that we obtain the properties aimed at by the source power co-modulation (SPoC) and the common spatial patterns (CSP) algorithms, as detailed in the following.

Spatial filters for task T2

The algorithm SPoC (Dähne et al., 2014) implements a supervised regression approach aimed at estimating a spatial filter that maximizes the correlation between —in a narrow frequency band—and the one-dimensional motor performance score , obtained using task T1. This is achieved by defining , i.e., the -weighted sample covariance matrix, and , i.e., sample covariance matrix. Using data of all the copy-draw test trials available in a single session, both matrices are determined as and , respectively. The decoding performance obtained in task T2 is evaluated with the Pearson correlation coefficient .

Spatial filters for T3

The CSP (Ramoser et al., 2000) algorithm is an established supervised method for classification. It seeks a one or multiple spatial filters that extract sources which maximize the contrast in between DBS-on and DBS-off trials. This can be achieved by defining and , where and are the sample covariance matrices of the respective EEG signals. Decoding performance for T3 is evaluated using classification accuracy.

Interpretation of a spatial filter for tasks T2 and T3

One of the main advantages provided by CSP and SPoC is their spatial interpretability. However, as pointed out by Haufe et al. (2014), a spatial filter should not be interpreted directly. Instead, the corresponding spatial pattern should be used. The pattern corresponding to can be computed as .

Experimental setup

Patients

The study was approved by the local ethics committee at the University Medical Center Freiburg and was conducted according to the Declaration of Helsinki. All patients provided written informed consent prior to participation.

Seven PD patients (, …,) were bilaterally implanted in the STN with directional DBS electrodes (Vercise Cartesia™, Boston Scientific). Implantation trajectory and target location were planned using information from pre-operative computer tomography and magnetic resonance imaging and based on intraoperative microelectrode recordings. The dorsolateral motor area of the STN was targeted and the leads were placed such that the dead space between the segmented electrodes of the implanted lead was located at the optimal intraoperative stimulation site. Typically, this resulted in a further contact (ring) above the STN in the field of Forel. The lead locations were intraoperatively checked with lateral x-rays and post-operatively using helical computer tomography and fusion back to the planning data (using magnetic resonance imaging).

The DBS leads were kept externalized during three days following the implantation surgery. On day four, leads were internalized and the pulse generator was implanted. In addition to a short familiarization session with the copy-draw test, each patient participated in five EEG experimental sessions: one pre-session executed a day before lead implantation, three sub-acute sessions () on the three days following lead implantation, and one subchronic session () during the first inpatient control visit approximately four weeks after implantation surgery. Throughout the text, we use the notation to denote session b executed with patient a. For example, indicates the session executed with patient 5, three days after lead implantation. A description of the all available sessions can be found in Table 1.

Table 1

Description of patients and copy-draw sessions analyzed. Unless indicated otherwise, DBS was performed with a pulse width of 60 μs and frequency of 130 Hz. The STN contralateral to the dominant hand was stimulated during DBS-on, unless bilateral stimulation was applied, which is indicated. Medication intake times are provided relative to the experiment start. Scores of unified PD rating scale part III were obtained—unless indicated otherwise—in medication off state and several days prior to the implantation surgery. Tremor (T), bradykinesia (B), and rigor (R) scores obtained under DBS-on are given relative to this DBS-off baseline, ranging from 0 (no improvement), to +3 (complete symptom suppression), according to the average of following hand motor items of the unified PD rating scale part III: 3 (rigor); 4, 5 and 6 (bradykinesia); 15 and 16 (tremor). n/a marks unavailable information. Please refer to the main text for an explanation of the session codes.

Patient	Session	Age	Sex	Medication	DBS parameters	Trial duration	N trials	Preimplantation UPDRS-Part III	T	B	R
S_1	+2	55	M	normal acting l-Dopa at +26 min and fast acting l-Dopa at +43 min	3 mA	7.5 s	144	25 (med-on)	0	+2	+2
	c			fast acting l-Dopa at +15 min	2.4 mA, bi		216		n/a
S_2	+2	54	F	medication off	6 mA	6 s	144	18	+2	+2	+2
	c						168		+3	0	+2
S_3	c	56	M	medication off	6 mA	7.5 s	120	41	n/a
S_4	+2			retarded l-Dopa at −150 min	6.5 mA	6 s	168	23	0	0	+1
	+3	44	M	retarded l-Dopa at −330 min	5 mA		168		0	0	+1
	c			normal acting l-Dopa combined withCOMT-Inhibitor (Entacapone) at −210 min	3.5 mA, 30 μs		156		0	+1	+2
S_5	+2			retarded l-Dopa at −300 min	6.5 mA	6 s	144	33	0	+1	+1
	+3	69	M	retarded l-Dopa at −210 min	6.5 mA		144		+1	+1	+2
	c			normal acting l-Dopa −660 min	6.0 mA		192		0	+1	+1
S_6	+2	63	M	retarded l-Dopa at −330 min	3 mA	9 s	144	35	+1	+2	0
	c			medication off	3.5 mA		192		0	+2	+1
S_7	+2			n/a	4.5 mA	9 s	168	34	0	0	+2
	+3	55	M	retarded l-Dopa at −330 min	5 mA	6 s	168		n/a
	c			normal acting l-Dopa −660 min	3.5 mA	7.5 s	192		0	+2	+2

DBS parameters

The stimulation amplitude was adjusted for each patient individually to a level of 0.5 mA below side-effects onset but maximally at 6.5 mA. The electric field possessed an omnidirectional bipolar pattern, with the top unsegmented and adjacent segmented electrodes used as anode and cathode, respectively. During the sessions in the acute phase, DBS was performed unilaterally due to hardware limitations. During the subchronic phase, bilateral stimulation was possible but we opted for unilateral stimulation for the sake of consistency with the acute sessions. Patients were not explicitly informed about the DBS condition, however, a blinding probably was not obtained. The reason for this is that transient paresthesias were reported by most patients every time the DBS condition changed. As a limitation of the study, it is important to point out that also the experimenters were not blinded to the DBS condition.

Time structure of the copy-draw test

A copy-draw session is constituted by a series of blocks, with each block comprising 12 copy-draw trials (see Section 2.1). Self-timed pauses between trials and a longer break after each block were offered. Within a block, DBS was either on or off and alternated between blocks. In each session, at least ten copy-draw blocks were executed. Depending on the fatigue level reported by the patient, additional blocks could be performed. The trial duration was adjusted for each patient in the range between 6 s and 9 s, see Table 1. For further details about preprocessing of the behavioral data, please refer to Castaño-Candamil et al. (2019).

Electroencephalographic signals

Acquisition

EEG signals were recorded from passive Ag/AgCl electrodes (EasyCap GmbH, Germany) placed according to the extended 10–20 system. Impedance was kept below 20 k. All channels were referenced against the nose. The EEG signals were registered by BrainAmp DC amplifiers (Brain Products GmbH, Germany) at a sampling rate of 2 kHz, with an analog lowpass filter of 250 Hz applied before digitization.

For sessions , patient-specific fronto-central-parietal areas of the scalp were unreachable for EEG recording due to the surgical wounds. Consequently, EEG signals were recorded from a number of electrodes that varied from patient to patient. An example of a resulting EEG channel layout is depicted in Fig. 1. Fig. 2 gives a general overview of the following data processing pipeline.

Fig. 2

Processing pipeline of the EEG data. The lower row represent that main pipeline, the upper row is performed for identifying biological artifactual epochs and ICA components (by visual inspection), which are removed from the data in the main pipeline.

Preprocessing

Raw EEG signals were downsampled at 200 Hz and then segmented into epochs relative to the trial start ( s, determined by the patient touching the starting point of the trace, see Fig. 1). The epoch length is defined from  s until the end of each trial. Obtained trials were further segmented into non-overlapping s sub-epochs, where the ones exceeding a peak-to-peak amplitude of 80  in the  Hz band were marked as artifacts for later rejection. These thresholds are capable of detecting biological artifacts, such as eye movements and eye blinks, typically found in EEG recordings. Finally, independent component analysis (ICA) was carried out in the  Hz band on the non-artifactual sub-epochs. Then, the EEG data was projected onto the ICA space after components representing eye and muscle artifacts had been removed manually. All data analysis shown was performed with the MNE-Python package (Gramfort et al., 2013).

The sole purpose of the sub-epoching is to facilitate the identification and rejection of short time segments that might be contaminated by biological artifacts. To derive NMs, however, data of several non-contaminated sub-epochs acquired in the same trial were aggregated. Using this aggregated data, the covariance matrices necessary for SPoC and CSP were estimated.

Frequency filtering for decoding of motor performance (task T2) and DBS condition (task T3)

As SPoC and CSP require narrow-band signals, we performed a random search over the hyperparameters central frequency and frequency width , for each session individually. For each frequency band, preprocessed data were filtered using a 5-th order butterworth filter. During the random search, a total of 600 configurations were tested per session, with varying in the range  Hz and in  Hz. The frequency filtered data is then used to compute the spatial filters , according to the methods shown in the previous section.

For each hyperparameter configuration, the decoding accuracy was determined using cross-validation. For each cross-validation fold, two blocks of the copy-draw test were as validation set and the remaining blocks formed the training set. The reported score corresponds to the average across all possible block permutations. The NMs discussed in the next session correspond to configurations performing significantly above chance level. The latter was obtained by cross-validation, using a bootstrapping procedure with 200 label shuffles and a significance level of .

Results

In Sections 4.1, 4.2 we present the results for the extraction of motor performance (task T1), the decoding of motor performance (task T2), and for the decoding of the DBS condition (task T3). In Section 4.3, we relate the obtained motor decoding accuracy to the contrast of motor performance induced by the stimulation. Later in Section 4.4, we provide introspection about the obtained NMs by analyzing their corresponding spatial, temporal, and spectral features. Finally, Section 4.5 compares our data-driven NMs against group-level NMs.

Task 1: Extraction of hand motor performance labels

Fig. 3 provides the area under the ROC curve (AUC) to describe the performances obtained in task T1. The figure also shows the distribution of motor scores z separately per session. We observed a separability between DBS conditions above chance level for all of the sessions. This demonstrates that DBS on/off has induced a modulation of the hand motor performance, that the copy-draw test is sensitive to capture this modulation, and that a linear classification model is capable to project this modulation to a meaningful scalar value.

Fig. 3

Hand motor performance z extracted in task T1 and the area under the ROC curve (AUC) describing the separability between DBS-on/-off. For all sessions, AUC scores are significantly above chance level, tested using a bootstrapping procedure with 200 label shuffles at a significance level of (uncorrected).

Task 2 and Task 3: Decoding of motor performance and of the DBS condition from EEG signals

For all 16 sessions, Fig. 4 shows the accuracy achieved for the decoding of motor performance (task T2) and the DBS condition (task T3) across different hyperparameters and . The chance levels are given by solid lines in corresponding colors. By different circle sizes in Fig. 4 we indicate the robustness of a central frequency to provide statistically significant NMs under varying band widths . Specifically, this robustness is provided by the percentage of components achieving an statistically significant decoding performance.

Fig. 4

Per patient and session, accuracies (y-axes) achieved for the decoding of motor performance (green), and of DBS-condition (pink) are provided relative to varying central frequencies (x-axes). Accuracies for different hyperparameter configurations of and have been grouped into 2 Hz bins and averaged. Circle size represents the percentage of evaluations for a given 2 Hz bin which achieve an accuracy above chance level. Chance level is indicated by solid lines.

Defining bins robust as where at least 40 % of the components showed significant decoding accuracy, we can also observe that almost all sessions showing robust frequency bands for the decoding of motor performance also allowed to decode the DBS condition above chance level using the same frequency bands.

Finally, to determine if the NMs identified are able to detect changes of motor state while DBS is either on or off only (a minimum prerequisite in any aDBS pipeline), we analyze the decoding performance for task T2 according to the homogeneity of the stimulation condition in the test sets during the crossvalidation procedure. Specifically, we categorize each test set according to the stimulation conditions contained in it: only DBS-on, only DBS-off, or mixed. Fig. 5 shows the distribution of performances according to this categorization criterion.

Fig. 5

Distribution of decoding performance in the regression task (task T2), according to the stimulation condition of the trials in the test folds. In some specific cases, there is an apparent difference in performance according to the homogeneity of the stimulation condition in the test set, however, a global trend cannot be identified.

It can be seen that for two sessions, namely and , the decoding performance in mixed test folds is superior than that of folds with homogeneous stimulation condition. This observation calls for caution on the interpretation of these two sessions, since it might indicate that the high Pearson correlation achieved could be given by the discrete step-wise increase of motor performance due to changes in the stimulation. In some other cases, for example and , DBS-off folds delivered a low decoding performance. This suggests that the corresponding NMs would only be suitable for decoding hand motor performance during DBS-on. A reversed case is observed for sessions and , where NMs provided an improved decoding accuracy under DBS-off, compared to DBS-on folds.

Low performance variation explains low decoding accuracy

As observed in the previous section, not only the informative frequency bands varied largely between sessions, but also the maximum decoding accuracy obtained. To examine this closer, Fig. 6 shows the relationship between the motor performance separability induced by DBS-on/-off—quantified as the AUC obtained during the extraction of motor performance in task T1—and the maximum accuracy obtained for the decoding of motor performance via SPoC in task T2. The linear relation in this scatter plot suggests that it is necessary to have DBS-induced contrast of motor performance in the first place, in other to obtain NMs for task T2.

Fig. 6

DBS-induced changes in motor performance (AUC during extraction of motor performance) vs. the highest accuracy achieved in the decoding of motor performance from EEG signals.

Spatio-temporal signatures of informative neural sources

Unfortunately the incomplete EEG layouts used in sessions impede an interpretation of spatial patterns obtained from these sessions, regardless of the subspace decoding method (SPoC or CSP). Components obtained for sessions , however, are based on the full EEG layout. We performed a qualitative analysis of those NMs which were obtained for the random search in tasks T2 and T3 and which were statistically significant. Their spatial patterns are shown in Fig. 7. We do not consider them to reflect stimulation artefacts of the DBS, as they lack the specific spatial (and frequency) characteristics of stimulation artifacts explained in the Appendix A. (See Fig. 8.).

Fig. 7

Informative components of neural origin obtained for the decoding of motor performance and of DBS condition. Each panel contains: frequency spectrum of the spatially filtered data (top-left); time-locked dynamics of the spatially and frequency-filtered data wrt. to the beginning of the trial (bottom-left); spatial pattern of the component (right). Vertical black lines in the bottom-left plot mark the beginning and end of the copy-draw trial. Unless indicated, components shown were found informative for both motor performance and DBS condition. For , the spatial patterns of components indicated biological artifactual origins and thus, are not shown. Note that the polarity of the spatial patterns shown is arbitrary.

Fig. 8

Continued: informative components of neural origin obtained for the decoding of motor performance and DBS condition.

For each component we show (1) the spatial pattern , (2) the frequency spectrum of the projected data and (3) the band-power dynamics (”envelope”) of the component’s source, time-locked to the trial start. The latter is computed by the magnitude of the Hilbert transform of the spatially and temporally filtered signals.

We observed that the topological origin of the components varies across patients. Specifically, we found centro-parietal components for , and , frontal components for and , and an occipital component for . The topological features of the patterns match the expected spectral characteristics thereof, i.e., occipital components stemming from alpha-band frequencies ( and ), frontal components from the theta-band ( and ), and central components from the beta-band () which may origin from sensory-motor regions.

Trial-locked dynamics could also be identified: an event-related synchronization of the alpha-band (enhanced by DBS-on) can be observed for . For the same session, also a DBS-enhanced beta-band event-related desynchronization was found. Similarly, event-related desynchronization effects were observed for all the patterns of . End-of-trial artifacts are present in and , very likely unrelated to any neural activity of interest and probably stemming from muscle activity or eye movements. However, it is important to mention that only the interval between and the end of the trial was considered for solving tasks T2 and T3 and therefore, these biological artifacts did not play any role in the decoding performance reported.

Decoding accuracy of individual NMs compared to group-level NMs

We also investigated if NMs obtained with our supervised approach are capable of providing more information about hand motor performance than group-level NMs can provide using a fixed frequency band, extracted under an unsupervised framework. The latter were estimated using the spatio-spectral decomposition algorithm (Haufe et al., 2014, Nikulin et al., 2011). This algorithm delivers a spatial filter that enhances the signal to noise ratio of the pre-defined frequency band compared to neighboring bands, disregarding motor labels z.

Fig. 9 provides three scatter plots. They compare the decoding accuracy between the best spatio-spectral decomposition and best individually determined NMs (our approach) that represent cortical sources in the theta ( Hz), alpha ( Hz), and beta ( Hz) band.

Fig. 9

Comparison between motor performance decoding using spatio-spectral decomposition of the fixed theta-, alpha-, and beta-bands (x-axis) vs. the overall best performing data-driven NMs identified with our approach (y-axis), not limited to any frequency band. Sessions above the dashed line indicate better performance of the data-driven NMs. Statistical significance is indicated by filled markers and was determined with Bonferroni corrected Wilcoxon signed-rank test at p-value of 0.05.

In at least 12 out of 16 sessions per frequency band analyzed, we obtained a statistically significant superior decoding accuracy of motor performance by our data-driven, patient-specific NMs. The remaining three to four sessions did not show any significant performance difference between patient-specifc NMs and the group-level NMs.

Discussion

In the following, we will first discuss the observed results in the context of tasks T1, T2, and T3 as formulated in the Introduction section. Then we will discuss the potential impact of the proposed machine learning approach upon future adaptive DBS systems.

Task 1 – Extraction of motor performance

Our approach is built around the copy-draw test, which is designed to capture hand motor performance in a fast and repeatable way. We have shown that a scalar performance score extracted from the copy-draw test is modulated by DBS. We also found that DBS-induced contrast in this score correlates with how well we can decode it from brain signals alone. In sessions where DBS induce only a weak behavioral contrast, then we can assume that most of the variance observed in labels z for these sessions is either caused by external noise sources independent of the DBS, or by neural processes not captured by the EEG.

A caveat of our approach is that we disregard the possible presence of a washout effect upon the motor features, despite the relatively short intervals (approximately every three to five minutes) after which DBS was toggled between on and off. From the perspective of the machine learning model trained with labels derived from task T1, the washout effect can be considered as some form of label noise, as some trials at the beginning of DBS-off blocks were executed with a residual effect of the previous DBS-on condition. This label noise increases the complexity of the learning tasks in T2 and T3 even further, however, it does not erode the interpretability and validity of the results obtained. Finally, it is important to mention that to this point there is no proof that an enhanced performance in the copy-draw test translates into a clinically improved movement pattern of the patient, although we think this might well be possible. A comparison of our functionally specialized hand task with UPDRS-III motor scores has not been performed here but shall be addressed in future research.

Task 2 – Decoding of motor performance from brain signals

In 13 out of the 16 sessions analyzed, we were able to identify NMs of hand motor performance. In eleven of them, the supervised session-specific data-driven NMs achieved a significantly better decoding accuracy compared to other NMs extracted with an unsupervised learning method from predefined frequency bands, i.e., enhanced beta-, alpha-, and theta-band power. This suggest that our data-driven approach should be the preferred strategy.

The NMs identified by SPoC reveal a prominent variability across different sessions, even within the same patient. On the one hand, variability within a subject might have its cause in the strong non-stationary dynamics elicited by the stun effect, which can last up to several months upon motor performance after lead implantation (Chen et al., 2006, Eusebio and Brown, 2009, Piña-Fuentes et al., 2017). However, since research regarding the impact of stun effects upon cortical rhythms is missing, its influence upon the NMs cannot be conclusively claimed. On the other hand, variability between subjects might be attributed to different symptomatological groups across subjects, a characteristic intrinsic to PD.

We also have shown that our approach allows for the decoding of motor performance on rather short time scales, e.g., when switching between DBS conditions every few tens of seconds. For adaptive DBS this is very valuable, since PD-symptoms are characterized by complex multi-timescale non-stationary dynamics. The latter are symptom-depended, but can govern fluctuations of the motor state in the order of seconds—e.g., for tremor and rigor—up to several hours—e.g., for bradykinesia (Cooper et al., 2013).

Task 3 – Difference in components relevant for decoding of motor performance and DBS condition

Most frequency bands that were found informative to regress the motor performance were also informative for the classification of DBS condition. A few exceptions only showed components informative solely for DBS condition decoding. They may be interpreted as underlying neural processes that are affected by DBS, but that do not interact with hand motor performance. On the contrary, components informative solely for the decoding of hand motor performance may stem from neural processes not affected by DBS.

Having in mind the design of an adaptive system, it is not only important to identify neural features correlating with motor performance, but these neural features must be required also to be susceptible to DBS modulation. Otherwise, a NM describing the motor state that cannot not be modulated by DBS, which would prevent to build a controllable adaptive DBS system.

Interpretability of identified NMs

For the subchronic sessions , we were able to characterize the NMs in terms of their spatial patterns. Here, we observed NMs that focussed not only over central electrodes—and thus, possibly from processes in the sensori-motor cortex—, but also components with a focus on occipital, parietal, and temporal areas. The latter suggest, that the decoding delivered information about neural processes not directly involved in motor execution, but which may modulate motor control, e.g., visual perception processes relevant for the hand motor task.

In two of the sessions, we determined frontal theta-band components. These NMs may reflect the influence of DBS upon cognitive processes via a functional connection between STN and the frontal cortex as reported in the literature (Zavala et al., 2014, Frank et al., 2007, Aiello et al., 2017). This hypothesis is supported by the lateralization of these components wrt. to the stimulated hemisphere. The modulation of oscillatory components in the theta-band has been observed upon processing negative feedback, realizing errors, or following an event which indicates need for increased control (Ishii et al., 1999, Cavanagh et al., 2012, Cavanagh and Frank, 2014). Certainly the copy-draw test could cause such events.

During the analysis of decoded NMs, we came across many components dominated by higher frequent oscillations ( Hz), which probably represent muscle activity. This non-neural contributions should not automatically be condemned, as they may potentially contribute valuable information for designing adaptive systems. Specifically, the electromyographic signals of cervical muscles may provide short term information about the muscle tone, which could be a highly valuable non-neural biomarker for PD. Finally, in the present analysis, the scope has been limited to the bandpower of cortical features. However, other phenomena such as connectivity indices shall also be investigated, as these less-explored features have been deemed promising (Bočková and Rektor, 2019, Geraedts et al., 2018).

Significance and implications for adaptive DBS systems

A closed-loop adaptive system, implemented as in Fig. 10, can only perform well if the controller is provided with information about the motor state of an individual patient. To avoid continuous behavioral testing, an estimate of the motor state can be provided by NMs. Our contribution tackles the problem of NMs individualization. Our results suggest that the NMs decoded by machine learning methods could contribute more meaningful information to the controller of an adaptive DBS system about the ongoing motor state than a fixed power feature can provide. This proposition is in line with recent findings of Tan and colleagues (Tan et al., 2019) in the context of adaptive DBS for essential tremor patients. Our NMs obtained with task T2 can be exploited to realize adaptive stimulation strategies in the same way as fixed power features are used in the existing literature. Thus, they can be integrated into control strategies that have already been validated, e.g., threshold-based or proportional controllers (Little et al., 2013, Rosa et al., 2015). Our NMs could also be deployed in combination with more complex control strategies such as coordinated reset (Adamchic et al., 2014), phase-specific stimulation (Cagnan et al., 2016), or with controllers realized by reinforcement learning models (Kumar et al., 2016).

Fig. 10

In an adaptive DBS system, motor performance is estimated by a decoder for neural markers (red) from brain signals. Its output informs a controller, which decides on the ongoing stimulation intensity. Using the neural marker information, the controller does not require any active polling of motor performance—this is done during calibration only (in yellow).

However, several caveats and limitations require further attention on the way to a future clinical deployment. First, in our present contribution, we focused solely on hand motor performance, and only in a hand writing context. This clearly covers a very small range of motor and non-motor processes which can be affected by PD. An expanded motor and non-motor test battery (possibly adapted to the typical activities of a patient) could provide a remedy.

Second, EEG is an impracticable recording modality for a prospective chronic adaptive system. We can only speculate, if invasive recordings techniques such as electrocorticography signals could provide comparable or better information to drive a controller.

Third, we observed a high NM variability even within subjects. This points at the necessity of a framework that allows for permanent updates of NMs to cope with other causes of non-stationarity.

Conclusion

The use of supervised machine learning methods in combination with a fast hand motor test reveals patient-specific bandpower NMs hidden in noisy EEG recordings. The components obtained are sensitive to DBS, or informative about the level of hand motor control, or combine both characteristics. Compared to state-of-the-art NMs derived from group-level studies, these data-driven NMs are more informative in the sense that motor performance can be decoded with higher accuracy. Thus, they provide a preferable surrogate of the patient’s (motor) state—even on short time scales—and should be considered for the design of individualized adaptive DBS systems. For chronically implanted systems, however, it needs to be explored if this information can also be accessed using, for example, electrocorticography signals, as implanted electrodes are preferred for reasons of usability, compared to EEG signals. With our contribution, we ultimately seek to contribute to a roadmap towards an individualized data-driven adaptive system based on cortical information that goes beyond the exclusive use of LFP signals from basal ganglia.

CRediT authorship contribution statement

Sebastián Castaño-Candamil: Conceptualization, Methodology, Software, Formal analysis, Visualization, Writing - original draft. Tobias Piroth: Methodology, Resources, Writing - review & editing. Peter Reinacher: Methodology, Resources, Writing - review & editing. Bastian Sajonz: Methodology, Resources, Writing - review & editing. Volker A. Coenen: Conceptualization, Methodology, Resources, Writing - review & editing, Supervision. Michael Tangermann: Conceptualization, Methodology, Resources, Writing - original draft, Writing - review & editing, Supervision.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.