Impaired phase synchronization of motor-evoked potentials reflects the degree of motor dysfunction in the lesioned human brain.

The functional corticospinal integrity (CSI) can be indexed by motor-evoked potentials (MEP) following transcranial magnetic stimulation of the motor cortex. Glial brain tumors in motor-eloquent areas are frequently disturbing CSI resulting in different degrees of motor dysfunction. However, this is unreliably mirrored by MEP characteristics. In 59 consecutive patients with diffuse glial tumors and 21 healthy controls (CTRL), we investigated the conventional MEP features, that is, resting motor threshold (RMT), amplitudes and latencies. In addition, frequency-domain MEP features were analyzed to estimate the event-related spectral perturbation (ERSP), and the induced phase synchronization by intertrial coherence (ITC). The clinical motor status was captured including the Medical Research Council Scale (MRCS), the Grooved Pegboard Test (GPT), and the intake of antiepileptic drugs (AED). Motor function was classified according to MRCS and GPT as no motor deficit (NMD), fine motor deficits (FMD) and gross motor deficits (GMD). CSI was assessed by diffusion-tensor imaging (DTI). Motor competent subjects (CTRL and NMD) had similar ERSP and ITC values. The presence of a motor deficit (FMD and GMD) was associated with an impairment of high-frequency ITC (150-300 Hz). GMD and damage to the CSI demonstrated an additional reduction of high-frequency ERSP (150-300 Hz). GABAergic AED increased ERSP but not ITC. Notably, groups were indistinguishable based on conventional MEP features. Estimating MEP phase synchronization provides information about the corticospinal transmission after transcranial magnetic stimulation and reflects the degree of motor impairment that is not captured by conventional measures.

INTRODUCTION

Transcranial magnetic stimulation (TMS) to the primary motor cortex (M1) produces complex, multiphasic corticospinal volleys (CSV). These volleys are composed of an early direct (D‐) wave generated by the activation of the corticospinal neurons and followed by a series of indirect (I‐) waves generated by GABAergic transsynaptic input from distinct intracortical circuits (Chen et al., 2008; Groppa et al., 2012; Rossini et al., 2015). Changing the direction of the electric field activates different circuits and I‐waves within 20–40 ms after TMS (Di Lazzaro et al., 2001; Merton & Morton, 1980). CSV are integrated on the level of a single spinal alpha motoneuron which subsequently activates associated muscle fibers, the so‐called single motor unit (SMU) (Rossini et al., 2015; Rothwell, Thompson, Day, Boyd, & Marsden, 1991). In invasive electromyographic (EMG) recordings, CSV are mirrored by an increase of firing probabilities of SMU (Day et al., 1989; Hanajima et al., 1998; Merletti, Holobar, & Farina, 2008; Sakai et al., 1997). Contrary, in surface EMG, the net activity of all stimulated SMU is depicted as a motor‐evoked potential (MEP) blurring the insight to the CSV.

Temporal dispersion of the CSV is supposed to result in phase cancelation of individual SMUs (Bestmann & Krakauer, 2015; Groppa et al., 2012; Rossini et al., 2015). Thus, temporal synchrony of the CSV determines the latency, shape and amplitude of the MEPs and contributes to the high trial‐to‐trial variance of the MEP (Bestmann & Krakauer, 2015; Groppa et al., 2012; Rossini et al., 2015). While corticospinal conduction varies physiologically due to fluctuations of cortical (Khademi, Royter, & Gharabaghi, 2018; Naros, Lehnertz, Leão, Ziemann, & Gharabaghi, 2020) and spinal excitability (van Elswijk et al., 2010; Naros et al., 2020), lesions within the corticospinal system (e.g., brain tumors) can cause additional impairment of CSV transmission (Cirillo, Calabro, & Perez, 2016; Hallett, 2000; Kobayashi & Pascual‐Leone, 2003; Stinear et al., 2007). However, conventional MEP characteristics (i.e., latency and amplitude) are often ambiguous, resulting in an inconsistency between clinical (e.g., motor impairment) and electrophysiological findings (Machetanz et al., 2021; Mirchandani et al., 2020; Picht et al., 2012; Sollmann et al., 2017). Conversely, there is only limited data relating the temporal precision of CSV to the actual motor performance (Machetanz, Gallotti, et al., 2021; Machetanz et al., 2021).

Time‐frequency analysis enables the transformation of the surface EMG signal into a sum of sine waves estimating their magnitude and phase within the complete signal. The information depicted by the spectral representation of the surface EMG is often considered a global measure of SMU activity (Farina, Merletti, & Enoka, 2004). Thus, the decipherment of the MEP in the time‐frequency domain could extract information about the underlying CSV. In fact, first studies have shown that the MEP is characterized by a synchronization of EMG oscillations up to 500 Hz mirroring physiological SMU discharge patterns (Machetanz, Gallotti, et al., 2021; Machetanz, Wiesinger, et al., 2021).

Against this background, we have hypothesized that glial brain tumors might affect corticospinal transmission by temporal dispersion (i.e., phase asynchrony) of CSV (Machetanz, Gallotti, et al., 2021; Machetanz, Wiesinger, et al., 2021). The objective of the present study was to elaborate whether the phase synchronization of MEPs relates to the clinical status in patients with glial brain tumors.

METHODS

Participants

This prospective study enrolled 80 participants (i.e., 21 healthy subjects and 59 consecutive patients with glial brain tumors in motor eloquent areas) undergoing a navigated transcranial magnetic stimulation (nTMS) examination. In patients, nTMS mapping was performed for clinical purposes prior to brain surgery. Healthy subjects were included as a control group (CTRL). Patients were classified into three categories by experienced neurosurgeons based on their clinical status. The Medical Research Council Scale (MRCS) and the Grooved Pegboard Test (GPT) were used to determine their motor status and dexterity of the contralateral upper limb, respectively. The patients were considered to have a fine motor deficit (FMD) with MRCS = 5 and a GPT score below the normative cohort, and to have a gross motor deficit (GMD) with MRCS < 5. All other patients were considered to have no motor deficit (NMD). Details of the clinical and demographic participants characteristics are depicted in Table 1. The study was approved by the local ethics committee of the Medical Faculty of the Eberhard Karls University Tuebingen. All participants gave written informed consent.

TABLE 1

Patients' clinical and imaging characteristics

	CTRL	NMD	FMD	GMD
	Healthy subjects	No motor deficit	Fine motor deficits	Gross motor deficits
	n = 21	n = 28	n = 15	n = 16
Age	31.6 ± 11.7	45.0 ± 15.8	50.7 ± 17.3	56.5 ± 16.8	H = 20.75, p < .001
					Kruskal–Wallis
Gender (f:m)	14:7	11:17	5:10	9:7	X 2 = 5.49, p = .139
					X 2‐test
Height (cm)	173 ± 9	173 ± 13	172 ± 8	170 ± 9	H = 0.94, p = .816
					Kruskal–Wallis
Weight (kg)	65.1 ± 15.0	78.9 ± 18.7	81.0 ± 14.0	77.1 ± 14.5	H = 14.24, p = .003
					Kruskal–Wallis
Diagnosis
HGG	–	20	10	15	X 2 = 3.83, p = .147
LGG	–	8	5	1	X 2‐test
AED
Yes	–	18	6	7	X 2 = 0.42, p = .812
No	–	10	9	9	X 2‐test
Tumor size (cm3)	–	49.3 ± 51.7	38.7 ± 40.3	37.7 ± 36.9	H = 0.24, p = .887
					Kruskal–Wallis
MAParea (cm2)	11.9 ± 4.9	10.1 ± 6.0	12.3 ± 4.9	10.9 ± 5.4	H = 2.10, p = .553
					Kruskal–Wallis
DTI
Mean FA	–	0.50 ± 0.05	0.49 ± 0.07	0.44 ± 0.05	H = 5.65, p = .059
Mean ADC (10−4 mm2/s)	–	7.64 ± 0.76	8.04 ± 0.53	8.44 ± 0.47	H = 6.23, p = .044
		n = 12	n = 10	n = 10	Kruskal–Wallis

Note: Values in bold are significant p‐ values.

Abbreviations: ADC, apparent diffusion coefficient; AED, antiepileptic drugs; CTRL, control group; DTI, diffusion tensor imaging; FA, fractional anisotropy; FMD, patients with fine motor deficits; GMD, patients with gross motor deficits; HGG, high grade glioma; LGG, low grade glioma; NMD, patient with no motor deficits.

Magnetic resonance imaging

All patients and healthy subjects received anatomical MR imaging using a 1.5 T magnetic resonance imaging (MRI) unit (Skyra/Prisma‐fit/Aera, Siemens Healthineers, Erlangen, Germany) with an 8‐channel head coil (contrast‐enhanced T1‐weighted mprage sequence, isovoxel 1, TR/TE: 2300/2.29). In 32/59 patients, additional diffusion‐tensor imaging (DTI) was performed with a single‐shot spin echo at a b‐value of 1,000 s/mm2 along 12–64 geometric directions. The anatomical MRI data set was imported to the nTMS system (Nexstim Eximia, version 3.2.2, Helsinki, Finland) for cortical mapping. An experienced neurosurgeon (G.N.) delineated manually the tumor on the bases of individual MR images using the Mricro software (https://www.nitrc.org/projects/mricron). The individual tumor volume was noted and the tumor mask was saved for further analysis.

Navigated transcranial magnetic stimulation

The cortical mapping procedure has been described previously in detail (Kraus et al., 2016; Kraus & Gharabaghi, 2015; Leão, Naros, & Gharabaghi, 2020; Mathew, Kübler, Bauer, & Gharabaghi, 2016): Prior to the mapping, patients' and subjects' anatomical T1‐weighted MR images were co‐registered to the participant's head with a registration error of <2 mm. nTMS mapping was performed with a biphasic figure‐8 coil (eXimia®, Nexstim, Helsinki, Finland). After determining the “hotspot” yielding the largest motor‐evoked potential (MEP) from the contralateral abductor pollicis brevis muscle (APB) and the first dorsal interosseous muscle (FDI), the resting motor threshold (RMT) of the FDI was obtained, defined as the minimum stimulus intensity that resulted in a MEP > 50 μV in at least 5/10 trials. In case of insufficient quality of the FDI recording, the APB muscle was selected for RMT detection. The orientation of the induced current in the brain was posterior–anterior (PA) for the first phase and anterior–posterior (AP) for the second phase of the stimulus. The orientation of the electric field, calculated on the basis of the individual MRI of each subject by the eXimia software, was kept perpendicular to the mapped sulcus. Subsequently, the cortex was mapped with 110% RMT starting at the primary motor cortex and then extending around this spot to cover the primary motor cortex, somatosensory cortex, and premotor cortex. As glial brain tumors are often blurring the gyral and sulcal anatomy in the MRI, mapping was extended until no MEP responses were obtainable (Figure 1a). Thus, an average of 133.6 ± 69.3 (51–336) stimuli were applied per patient and map. Stimulation sites were visualized on the surface at a depth of 25–30 mm. Coordinates of the stimulation sites were automatically saved by the eXimia software for later analysis. The Matlab function “convhull.m” was used to calculate the convex hull enveloping the x and y coordinates of the stimulation sites. The area (in cm2) of the convex hull represents the extend of the mapping area (MAParea).

FIGURE 1

nTMS results. (a) exemplary data of a characteristic nTMS map. Colored dots represent nTMS coordinates eliciting a MEP (≥100 μV). In contrast, gray dots indicate spots with no MEPs. (b) Patient's nTMS results were used as a seed for deterministic DTI fiber tracking. (c) There was a high intrasubject variability of MEP time‐series (first row). However, most participant showed a biphasic EMG response within 20–40 ms after TMS. Time‐frequency decomposition revealed a group‐dependent power increase (ERSP, second row) as well as a phase synchronization (ITC, third row) between 37 and 500 Hz with a local maximum around 100–150 Hz. Notably, EMG phase synchronization was detectable prior to the MEP onset latency. Black triangles represent the mean group values for Lat0 (i.e., MEP onset latency), Lat1 (i.e., latency of the maximum positive deflection of the MEP) and Lat2 (latency of the minimum negative deflection of the MEP). For further statistical analysis, time‐frequency representation of the MEP was divided in six quadrants (ERSP1‐6) and the mean ERSP and ITC values were calculated for each quadrant. CTRL, control group; FMD, patients with fine motor deficits; GMD, patients with gross motor deficits; MEP, motor‐evoked potentials; NMD, patients with no motor deficit

EMG recordings and data analysis

During nTMS mapping, EMG of the contralesional APB and FDI were recorded with the eXimia EMG amplifier system (3 kHz sampling rate, band‐pass filter of 10–500 Hz) using Ag/AgCl wet gel surface electrodes (AmbuNeuroline 720, Ambu GmbH, Germany). Data analysis was performed with custom‐written scripts based on MATLAB (Mathworks Ltd., R2017a) and its open source toolboxes EEGlab (Delorme & Makeig, 2004) and Fieldtrip (Oostenveld, Fries, Maris, & Schoffelen, 2011). EMG data was segmented into epochs from −100 to +100 ms relative to the TMS pulse. There was no filtering or data preprocessing except of linear detrending of the epochs. Generally, the FDI muscle was selected for analysis. The selection of nTMS trials (i.e., MEP+ trials) was based on a three‐step approach: (1) Trials with peak‐to‐peak EMG amplitudes of ≥100 μV during 20–40 ms after TMS were selected; (2) Subsequently, trials with MEP onset latencies ≤15 ms or ≥30 ms and trials with pre‐stimulus baseline EMG activation exceeding 3 SDs of mean EMG baseline (i.e., −100 to 0 ms in relation to TMS) were removed as artifacts. (3) All trials and patients were controlled visually. The remaining MEP+ trials (31.8 ± 23.3 [9-135]) were included for further analysis, all other trials were classified as MEP‐. In case of a bad signal‐to‐noise ratio or a number of artifacts higher than the average, the APB muscle was chosen for further analysis. A Matlab‐based custom‐written script was used to automatically detect several time series characteristics of the MEP: Amp (i.e., peak‐to‐peak amplitude), Lat0 (i.e., MEP onset latency), Lat1 (i.e., latency of the maximum positive deflection of the MEP) and Lat2 (latency of the minimum negative deflection of the MEP). The time‐frequency analysis of the MEP was performed on the basis of a Morlet wavelet approach with a fixed wavelet length of 40 ms (as implemented by the newtimef function of the EEGlab toolbox) (Delorme & Makeig, 2004). This approach resulted in a spectral resolution of 1 Hz (37–500 Hz) and a temporal resolution of 0.333 ms (−79.333 to 79.333 ms relative to the TMS pulse). Event‐related spectral perturbations (ERSP) were calculated (in dB) and trial‐wise normalized to the baseline spectrum (−79.3 to −10 ms relative to the TMS pulse) (Grandchamp & Delorme, 2011). Transforming a signal to the time‐frequency domain (by wavelet or Fourier transformation) splits the time‐series signal into a sum of oscillations. It results in coefficients Z (t,f) for each combination of time point t and frequency f for each trial i. These coefficients are complex numbers with real A(t,f) and imaginary B(t,f) components that represent the magnitude and the phase of an oscillation with the frequency f. The intertrial coherence (ITC) is a measure of how consistent the oscillatory phase is across an ensemble of trials. The ITC is bounded between 0 and 1, with 1 being perfect intertrial coherence (i.e., exactly the same phase on every trial), and 0 being no intertrial coherence. ITC is calculated by

By dividing the complex Fourier coefficients by their absolute values, they are normalized on a unit circle. These normalized values are averaged over trials i. The ITC represents the absolute value of the complex means (Delorme & Makeig, 2004). To visualize the divergence of intertrial phase synchronization, we realigned trials according to the individual Lat1 (correcting for height‐related conductance differences) and calculated the amplitude‐weighted phase (AWP) for each trial by with Amax(f) being the maximum amplitude of the oscillation within the trial.

DTI tractography

After nTMS mapping, the coordinates of MEP+ were exported as DICOM from the Nexstim software and imported into the BrainLab iPlan 3.0 software. These coordinates were used to reconstruct a cortical region of interest (ROI) (Krieg et al., 2012; Machetanz et al., 2019). The ROI was fused to the anatomical T1‐weighted MRI and DTI dataset. In addition to the cortical ROI, a subcortical ROI was placed in the caudal pons based on the color‐coded fractional anisotropy (FA) map (Machetanz et al., 2019; Rosenstock et al., 2017). The CST was detected using a fiber length of 110 mm and an FA value corresponding to 75% of the individual FA threshold impeding any fiber detection (Frey et al., 2012; Krieg et al., 2012; Machetanz et al., 2019). Mean FA and apparent diffusion coefficient (ADC) values of the resulting CST were noted as imaging surrogates of its integrity.

Statistics

Statistical evaluation was performed using SPSS (IBM SPSS Statistics for Windows, Version 25.0, Armonk, NY: IBM Corp.) and custom‐written Matlab scripts including the FieldTrip toolbox and Matlab statistics toolbox. Group effects on clinical (e.g., age, gender, height, weight, diagnosis, AED intake), imaging (i.e., FA and ADC values, tumor size, MAParea) as well as electrophysiological characteristics (RMT, number of TMS trials, MEP amplitudes and latencies) were evaluated by nonparametric Kruskal–Wallis and X 2‐tests when applicable. Correlation analyses were based on Pearson's and Spearman's correlation coefficients. We used a cluster‐based permutation test implemented in the Fieldtrip toolbox (ft_freqstatistics.m) as a first explorative data analysis of the time‐frequency MEP representation when comparing two different cohorts (i.e., CTRL vs. NMD or AED+ vs. AED−). In general, time‐frequency data is characterized by high spectro‐temporal dimensionality. After transforming into the time‐frequency domain, the MEP signal is presented with a large number of (time, frequency)‐samples. A pairwise statistical comparison would result in a major multiple comparisons problem (MCP). The cluster‐based test statistics can help to reduce the dimensionality of the data and to solve the MCP by determining significant spectro‐temporal clusters. In the present study, each (time, frequency)‐sample was compared across the different conditions by an unpaired two‐tailed t‐test (ft_statfun_indepsamplesT.m). t‐Values were thresholded at the 2.5‐th and the 97.5‐th quantiles for a two‐sided test. Selected samples were clustered in connected sets on the basis of temporal and spectral adjacency. Cluster‐level statistics were then calculated by taking the sum of the t values within every cluster, and the resultant maximum summed t‐values were used to compute the statistical comparisons. The significance probability was calculated using a Monte‐Carlo method (with 1,000 permutations). By randomizing the data, the reference distribution of the maximum of summed cluster t values was acquired. Clusters from the original data were considered to be significant if their summed cluster t values were below or above the 2.5‐th and the 97.5‐th quantiles of the reference distribution (representing a two‐sided test with an alpha level 5%). As this approach is disputable from a statistical point of view (Maris & Oostenveld, 2007), we sought to confirm the results of the cluster‐based analysis. Based on the clusters depicted by the first analysis and on previous findings (Machetanz, Gallotti, et al., 2021; Machetanz, Wiesinger, et al., 2021), we divided the time‐frequency MEP representation depending on two time bins (15 ≤ t1 ≤ 30 ms and 30 ≤ t2 ≤ 45 ms) and three frequency bins (37 ≤ f1 ≤ 150 Hz, 150 ≤ f2 ≤ 300 Hz and 300 ≤ f3 ≤ 500 Hz) in six quadrants (Q): Q1(t1,f1), Q2(t1,f2), Q3(t1,f3), Q4(t2,f1), Q5(t2,f2), Q6(t2,f3) for additional statistical analysis. Differences in ERSP and ITC values in these quadrants were evaluated by a Kruskal–Wallis test. Every Kruskal–Wallis analysis of variance was followed up by Dunn's pairwise post hoc test with Bonferroni‐correction for MCP. Finally, as the effect of AED intake and the motor status (i.e., NMD, FMD, GMD) seemed to overlay, we confirmed our results by an additional multivariate analysis of variance (MANOVA) including the ERSP and ITC values of each spectro‐temporal quadrants as dependent variables and the AED intake (yes vs. no) and the patient's motor status (NMD, FMD, GMD) as independent variables. The MANOVA was followed up by univariate ANOVAs and a post hoc test with Bonferroni‐correction for MCP. Results are shown as mean ± SD.

RESULTS

Patient characteristics

The present study includes 59 consecutive patients with glial brain tumors and variable motor status (28 NMD, 15 FMD, and 16 GMD) and a control group of 21 healthy subjects (CTRL). Patients were significantly older and of larger weight than the healthy control group. However, there were no significant group differences in body height. Within the patient group, there were no differences in diagnosis and tumor size between NMD, FMD and GMD. The extent of the TMS‐based cortical mapping (MAParea) was comparable in all cohorts. 31/59 (53%) patients were on antiepileptic drugs (AED) at the time of the nTMS examination. AED intake was distributed equally between groups. All patients were treated with the GABAergic drug Levetiracetam. Patients with GMD demonstrated a significant impairment of the CST as measured by lower FA and higher ADC values in DTI in comparison to the other patient groups. Patients' characteristics are summarized in Table 1.

Time‐series analysis

nTMS cortical mapping was performed in all healthy subjects and patients in the same way. There were no significant group differences in the number of applied TMS pulses or the number of MEP+ trials (Table 2). Kruskal–Wallis analysis of variance showed a significant group effect on the resting motor threshold (RMT, Table 2). Dunn's post hoc test proved a significant higher RMT in patients with FMD in comparison to CTRL (p = .015, Bonferroni‐corrected).

TABLE 2

Patients' electrophysiological characteristics

	CTRL	NMD	FMD	GMD	Kruskal–Wallis
RMT (%)	34.4 ± 6.9	39.5 ± 9.1	44.5 ± 10.4	40.8 ± 10.2	H = 9.12, p = .028
No. of trials	135 ± 67	127 ± 66	151 ± 77	128 ± 73	H = 0.98, p = .827
MEP+	37 ± 17	34 ± 21	42 ± 21	41 ± 35	H = 1.86, p = .602
MEP−	92 ± 54	90 ± 56	105 ± 67	84 ± 46	H = 0.39, p = .943
Amp (μV)	649 ± 473	494 ± 317	365 ± 119	330 ± 160	H = 5.38, p = .146
Lat0 (ms)	23.8 ± 1.5	24.9 ± 2.2	25.8 ± 2.4	24.9 ± 3.0	H = 5.54, p = .136
Lat1 (ms)	28.0 ± 1.5	29.0 ± 2.1	29.8 ± 2.1	29.0 ± 3.0	H = 5.83, p = .120
Lat2 (ms)	32.0 ± 1.8	33.1 ± 2.4	34.2 ± 2.4	33.4 ± 3.3	H = 5.33, p = .149
Lat1‐0 (ms)	4.2 ± 0.5	4.1 ± 0.8	4.0 ± 0.6	4.1 ± 0.7	H = 1.95, p = .582
Lat2‐0 (ms)	8.2 ± 1.2	8.3 ± 1.3	8.7 ± 1.2	8.5 ± 1.4	H = 1.31, p = .728
Lat2‐1 (ms)	3.9 ± 1.2	4.1 ± 1.0	4.7 ± 1.0	4.4 ± 1.1	H = 2.59, p = .460
ERSP a	27.1 ± 6.1	25.5 ± 4.9	23.1 ± 4.7	19.0 ± 4.9	H = 19.53, p < .001
ERSP1 b	34.6 ± 7.2	29.7 ± 10.6	26.9 ± 9.6	22.3 ± 7.1	H = 15.91, p = .001
ERSP2 c	33.6 ± 6.1	26.2 ± 10.4	22.9 ± 7.9	20.30 ± 6.6	H = 21.97, p < .001
ERSP3 d	23.9 ± 5.8	18.2 ± 7.9	16.2 ± 6.0	16.1 ± 5.2	H = 15.06, p = .002
ERSP4 e	29.5 ± 8.3	29.6 ± 8.7	31.3 ± 7.1	23.7 ± 8.5	H = 5.93, p = .115
ERSP5 f	28.1 ± 8.3	26.0 ± 8.0	27.7 ± 4.3	19.2 ± 5.9	H = 14.03, p = .003
ERSP6 g	20.0 ± 7.7	17.9 ± 6.2	20.7 ± 4.3	15.8 ± 5.9	H = 4.63, p = .201
ITC a	0.54 ± 0.07	0.52 ± 0.09	0.42 ± 0.09	0.44 ± 0.06	H = 14.38, p = .003
ITC1 b	0.87 ± 0.07	0.84 ± 0.09	0.75 ± 0.19	0.76 ± 0.14	H = 11.33, p = .010
ITC2 c	0.65 ± 0.13	0.61 ± 0.16	0.46 ± 0.16	0.49 ± 0.17	H = 15.47, p = .001
ITC3 d	0.37 ± 0.12	0.34 ± 0.14	0.26 ± 0.07	0.31 ± 0.11	H = 7.82, p = .050
ITC4 e	0.79 ± 0.09	0.80 ± 0.10	0.72 ± 0.10	0.73 ± 0.14	H = 8.03, p = .045
ITC5 f	0.52 ± 0.12	0.52 ± 0.14	0.37 ± 0.14	0.40 ± 0.05	H = 19.60, p < .001
ITC6 g	0.32 ± 0.08	0.31 ± 0.09	0.23 ± 0.06	0.28 ± 0.08	H = 11.59, p = .009

Note: Values in bold are significant p‐values.

15 ≤ t ≤ 45 ms and 30 ≤ f ≤ 500 Hz.

15 ≤ t ≤ 30 ms and 30 ≤ f ≤ 150 Hz.

15 ≤ t ≤ 30 ms and 150 ≤ f ≤ 300 Hz.

15 ≤ t ≤ 30 ms and 300 ≤ f ≤ 500 Hz.

30 ≤ t ≤ 45 ms and 30 ≤ f ≤ 150 Hz.

30 ≤ t ≤ 45 ms and 150 ≤ f ≤ 300 Hz.

30 ≤ t ≤ 45 ms and 300 ≤ f ≤ 500 Hz.

In the time‐series domain, MEPs were characterized by several amplitude and latency measures and showed a large intrasubject and intersubject variability in morphology (Figure 1). In summary, there were no significant group differences of the conventional MEP amplitudes and latencies (Table 2). RMT was related to tumor volume (r = 0.27, p = .015; Pearson's), however, there was no correlation between tumor volume and MEP time‐series characteristics (i.e., Amp, Lat0, Lat1, Lat2, Lat1‐0, Lat2‐0. Lat2‐1; p > .05; Pearson's). AED intake did not affect time‐series parameter in patients (i.e., RMT, Amp, Lat0, Lat1, Lat2, Lat1‐0, Lat2‐0. Lat2‐1; p > .05, Kruskal–Wallis).

Time‐frequency analysis

Time‐frequency decomposition of the MEP transfers the EMG signal into a sum of sine waves and estimates their contribution to the complete signal. The transformation of the MEP revealed an ERPS projecting to a frequency band between 37 and 500 Hz with a local maximum around 100–150 Hz. At the same time, there was a high ITC. Notably, phase synchronization of the EMG was detectable up to 10 ms prior to the MEP onset (Figure 1c). Cluster‐based significance analysis showed a significant ERSP reduction for all patients in comparison to healthy subjects. Notably, there was no significant differences between NMD and FMD. The occurrence of a paresis (GMD), however, led to an additional ERSP reduction for frequencies between 150 and 300 Hz in comparison to other patients (NMD and FMD) (Figure 2a). For the ITC, there was no significant difference when comparing CTRL to patients without motor deficits (i.e., NMD). In contrast. There was a significant decrease of phase synchronization for patients with motor deficits (FMD and GMD) in comparison to CTRL and NMD (Figure 2b).

FIGURE 2

Cluster‐based analyses of group difference in ERSP power and phase synchronization. Group differences in ERSP (a) and ITC (b) values were analyzed using a cluster‐based permutation analysis (correcting for multiple comparisons). There was a significant reduction of high frequency (>150 Hz) ERSP for all patient groups in comparison to the control group. Within the patients groups an additional ERSP reduction was observed for patients with GMD. There was no significant ERSP difference between NMD and FMD. In contrast, decreased ITC phase synchronization of high frequencies (>150 Hz) accompanied the occurrence of motor deficits (FMD and GMD). Notably, NMD and FMD patients were distinguishable based on the impaired phase synchronization for frequencies >150 Hz. Significant time‐frequency bins are outlined in red (p < .05, cluster corrected). CTRL, control group; FMD, patients with fine motor deficits; GMD, patients with gross motor deficits; ITC, intertrial coherence; NMD, patients with no motor deficit

This explorative analysis was confirmed when comparing the overall (i.e., for 15 ≤ t ≤ 45 ms and 37 ≤ f ≤ 500 Hz) ERSP and ITC values (H = 19.53, p < .001 and H = 14.38, p = .003; Kruskal–Wallis). Dunn's post hoc pairwise test (Bonferroni‐corrected) determined lower ERSP values in patients with gross motor deficits (GMD) in comparison to CTRL and NMD (Figure 3a). In contrast, the occurrence of a motor impairment (FMD and GMD) led to an overall reduction of ITC (Figure 3b).

FIGURE 3

Group effects on time‐ and frequency‐domain parameters. (a) Bar plots show the overall ERSP (for 15 ≤ t ≤ 45 ms and 37 ≤ f ≤ 500 Hz) for the different cohorts. There was a significant reduction of overall ERSP for patients with GMD in comparison to the other cohorts. Horizontal lines represent significant group differences (p < .05, Dunns' post hoc test, Bonferroni‐corrected). (b) In contrast, a significant reduction of overall ITC (for 15 ≤ t ≤ 45 ms and 37 ≤ f ≤ 500 Hz) was observed for patients with motor deficits (FMD and GMD). Notably, NMD and FMD were distinguishable by the impaired phase synchronization. Horizontal lines represent significant group differences (p < .05, Dunns' post hoc test, Bonferroni‐corrected). (c) For further elaboration of different frequency bands, we divided time‐frequency representation of the MEP in six quadrants (ERSP1‐6 and ITC1‐6, see Figure 1). Radar plots show group mean values ± SD of each quadrant as dots and circles, while lines indicate significant pairwise group differences (p < .05, Dunns' post hoc test, Bonferroni‐corrected). Groups with no motor deficit (CTRL and NMD) are depicted in on the left side, groups with motor deficits (FMD and GMD) are illustrated on the right side of each radar plot. There was a reduction of high frequency (>150 Hz) ERSP of the first period of the MEP (i.e., 15–30 ms after TMS) in all brain tumor patients when compared to healthy subjects. (d) An ERPS reduction for 150–300 Hz oscillations during the second phase of the MEP (i.e., 30–45 ms after TMS) was associated with GMD. In contrast, there was a decrease of high frequency (>150 Hz) phase synchrony (ITC) in the early (e) and late (f) MEP period associated with the occurrence of motor deficits (FMD and GMD). Notably, there was no difference between CTRL and NMD. CTRL, control group; FMD, patients with fine motor deficits; ERSP, event‐related spectral perturbation; GMD, patients with gross motor deficits; ITC, intertrial coherence; MEP, motor‐evoked potentials; NMD, patients with no motor deficit

To evaluate the effect of different frequency bands, we divided the time‐frequency MEP representation in six quadrants depending on two time‐bins (15–30 ms, 30–45 ms) and three frequency‐bins (<150 Hz,150–300 Hz, 300–500 Hz) (Figure 1c). There was a significant group effect on the ERSP of all frequency bands during the early period of the MEP (ERSP1‐3, Table 2). Dunn's post hoc pairwise test (Bonferroni‐corrected) indicated a reduction of ERSP1‐3 of all patients groups in comparison to healthy subjects (Figure 3c). In contrast, the significant group effect on ERSP5 was attributed to a power reduction of 150–300 Hz frequencies during the later period of the MEP (i.e., for 30 ≤ t ≤ 45 ms) in patients with GMD (Figure 3d). Notably, NMD and FMD patients were indistinguishable based on their ERSP values.

For the ITC, there was a significant group effect in all frequency bands and throughout the MEP cycle (Table 2). While there was no impairment of the intertrial phase synchronization in motor competent subjects and patients (CTRL and NMD), ITC values dropped significantly in patients with motor deficits (GMD and FMD). Notably, patients with no motor deficit (NMD) were distinguishable from patients with fine motor deficits (FMD) by lower ITC values in frequencies >150 Hz (Figure 3e,f).

Considering the patients' motor status as ordinal variable (i.e., NMD, FMD, GMD), there was a significant (Spearman's) correlation between the motor status and ERSP1 (r = −0.30, p = .021), ERSP2 (r = −0.26, p = .046), ERSP4 (r = −0.27, p = .039), ERSP5 (r = −0.40, p = .002) and ITC5 (r = −0.26, p = .048). In contrast, there was no significant correlation between the motor status and conventional MEP features (i.e, RMT, Amp, Lat0, Lat1, Lat2, Lat1‐0, Lat2‐0. Lat2‐1; p > .05, Spearmen's).

Intertrial phase synchronization

To visualize the divergence of intertrial phase synchronization, we realigned the amplitude‐weighted phase (AWP) of each trial according to the individual Lat1 (correcting for height‐related conductance differences) (Figure 4). This visualization indicates an excellent phase alignment of the 100 Hz oscillation for all groups (Figure 4a,b). In contrast, 300 Hz oscillations showed a good phase alignment only in CTRL and NMD. Phase synchronization of 300 Hz oscillations appeared delayed in NMD in comparison to CTRL. Furthermore, phase alignment decreased in patients with motor deficits (FMD and GMD) (Figure 4c,d).

FIGURE 4

Phase synchronization in relation to the MEP cycle. Amplitude‐weighted phase (AWP) values of exemplary 100 Hz (a, b) and 300 Hz (c, d) oscillations were realigned to the individual Lat1 (first positive MEP peak) (t = 0 ms). (a, c) Each pixel row corresponds to the AWP values of one trial, smoothed with a 20‐trial boxcar window. (b, d) The average of the trials is shown in a and c, respectively. For comparison, mean Lat0 (MEP onset latency) and Lat2 (first negative MEP peak) are indicated as vertical lines. This visualization indicates an excellent phase alignment of the 100 Hz oscillation for all groups (a, b). In contrast, there was a delayed (NMD, FMD, and GMD) and decreased (FMD and GMD) phase synchronization of 300 Hz oscillations in glioma patients (FMD and GMD) (c, d). FMD, patients with fine motor deficits; GMD, patients with gross motor deficits; MEP, motor‐evoked potentials; NMD, patients with no motor deficit

Impact of antiepileptic medication

Patients were divided into two groups according to the presence of AED medication to evaluate its effect on the time‐frequency representation. Remarkably, patients with AEDs had a significant higher ERSP for a time period of 30–45 ms after TMS (Figure 5a). However, there was no AED effect on the phase synchronization of MEPs (Figure 5b).

FIGURE 5

Effect of AED on time‐frequency representation of MEPs. Cluster‐based permutation analysis proved that GABAergic AED intake significantly increased ERSP power (a) for high frequencies (>100 Hz) in the late MEP period (30–45 ms after TMS) without affecting phase synchronization (b) in comparison to patients without AED. Significant time‐frequency bins are outlined in red (p < .05, cluster corrected). AED, antiepileptic drug intake; ERSP, event‐related spectral perturbation; MEP, motor‐evoked potentials

Subsequently, we applied a multivariate analysis of variance (MANOVA) on the ERSP (ERSP1‐6) and ITC (ITC1‐6) values to evaluate the effect of AED intake (yes vs no) while controlling for the motor status MS (i.e., NMD, FMD, GMD). There was a statistically significant difference in time‐frequency parameter of MEP based on both the AED intake (F [12,42] = 2.41, p = .018; Wilk's Λ = 0.593, partial η 2 = 0.41) and the MS (F [24,84] = 2.28, p = .003; Wilk's Λ = 0.367, partial η 2 = 0.39). Notably, there was no significant multivariate AED x MS interaction (F [24,84] = 1.12, p = .341; Wilk's Λ = 0.575, partial η 2 = .24). Follow‐up ANOVAs confirmed an increase in ERSP4‐6 after AED intake (Table 3). Notably, there was a significant main effect of the motor status on the same ERSP values (ERSP4‐6) but no AED x MS interaction (Table 3). Pairwise post hoc tests (Bonferroni‐corrected) on MS subgroups (i.e., NMD, FMD and GMD) confirmed a power reduction for GMD in comparison to FMD (ERSP4, p = .020; ERSP5, p < .001; ERSP6: p = .033) and/or NMD (ERSP4, p = .109; ERSP5, p = .007; ERSP6: p = 1.000). Additionally, follow‐up ANOVAs showed a significant main effect of MS on ITC2 and ITC5 (Table 3). Pairwise post hoc tests (Bonferroni‐corrected) on MS subgroups (i.e., NMD, FMD and GMD) confirmed an impairment of phase synchronization for both FMD (p = .013) and GMD (p = .013) in comparison to NMD. There was no significant difference in ITC5 when comparing FMD and GMD (p = 1.000).

TABLE 3

ANOVA follow‐up analyses

Source		Type III sum of squares	df	Mean square	F	Significance	Partial eta squared
AED	ERSP1	184.62	1	184.62	2.31	.135	.04
	ERSP2	166.45	1	166.45	2.38	.129	.04
	ERSP3	25.37	1	25.37	0.56	.457	.01
	ERSP4	577.91	1	577.91	10.36	.002	.16
	ERSP5	445.53	1	445.53	12.71	.001	.19
	ERPS6	258.44	1	258.44	10.41	.002	.16
	ITC1	0.02	1	0.02	1.09	.302	.02
	ITC2	0.02	1	0.02	0.78	.381	.01
	ITC3	0.00	1	0.00	0.03	.857	.00
	ITC4	0.01	1	0.01	0.98	.326	.02
	ITC5	0.00	1	0.00	0.05	.819	.00
	ITC6	0.00	1	0.00	0.67	.416	.01
MS	ERSP1	459.83	2	229.91	2.88	.065	.10
	ERSP2	286.70	2	143.35	2.05	.139	.07
	ERSP3	22.60	2	11.30	0.25	.779	.01
	ERSP4	477.15	2	238.58	4.28	.019	.14
	ERSP5	581.78	2	290.89	8.30	.001	.24
	ERPS6	179.98	2	89.99	3.62	.033	.12
	ITC1	0.07	2	0.03	2.36	.104	.08
	ITC2	0.20	2	0.10	4.08	.023	.13
	ITC3	0.05	2	0.02	1.68	.196	.06
	ITC4	0.07	2	0.03	2.82	.069	.10
	ITC5	0.19	2	0.10	6.57	.003	.20
	ITC6	0.04	2	0.02	3.08	.054	.10
AED × MS	ERSP1	477.14	2	238.57	2.98	.059	.10
	ERSP2	404.88	2	202.44	2.89	.064	.10
	ERSP3	145.86	2	72.93	1.62	.208	.06
	ERSP4	146.63	2	73.32	1.31	.277	.05
	ERSP5	75.84	2	37.92	1.08	.346	.04
	ERPS6	52.27	2	26.14	1.05	.356	.04
	ITC1	0.20	2	0.10	6.80	.002	.20
	ITC2	0.17	2	0.08	3.47	.038	.12
	ITC3	0.01	2	0.01	0.33	.720	.01
	ITC4	0.06	2	0.03	2.59	.085	.09
	ITC5	0.04	2	0.02	1.27	.290	.05
	ITC6	0.02	2	0.01	1.90	.160	.07

Notes: A multivariate analysis of variance (MANOVA) was applied on the event‐related spectral perturbation (ERSP1‐6) and intertrial coherence (ITC1‐6) to evaluate the effect of antiepileptic drug (AED) intake (yes vs no) while controlling for the motor status (MS). There was a multivariate significance for both the AED intake (F [12,42] = 2.41, p = .018; Wilk's Λ = 0.593, partial η 2 = 0.41) and the MS (F [24,84] = 2.28, p = .003; Wilk's Λ = 0.367, partial η 2 = 0.39) without an AED × MS interaction (F [24,84] = 1.12, p = .341; Wilk's Λ = 0.575, partial η 2 = 0.24). The table shows the results of the follow‐up ANOVA. Significant values (p < .05) are highlighted in bold.

DISCUSSION

Corticospinal volleys (CSV) after transcranial magnetic stimulation (TMS) induce a cyclic motor‐evoked potential (MEP). Precise corticospinal transmission of CSV predicts a high trial‐to‐trial phase coherence (i.e., ITC) and a high magnitude (i.e., ERSP) of the induced EMG response. In turn, temporal dispersion will reduce spectral and phase synchronization (Ding & Simon, 2013; Makeig, Debener, Onton, & Delorme, 2004). Time‐frequency transformation of the MEP indicates a high ITC of a broad frequency band after TMS. In glial brain tumor patients, the occurrence of motor deficits was associated with a disruption of this high frequency phase alignment, in particular for frequencies between 150 and 300 Hz. In contrast, an ERSP reduction of this frequency band indicated a gross motor dysfunction (i.e., paresis). Notably, there were no significant differences in the conventional time‐series MEP parameter (i.e., MEP peak latencies and amplitudes).

TMS mediates its effect by transsynaptic activation of I‐waves seen as rhythmic CSV of up to 500 Hz (Opie & Semmler, 2021; Ziemann, 2020). The direction of the electric field determines the activation of different intracortical networks. A posterior‐lateral to anterior‐medial (PA) current, induced by a monophasic coil, preferentially activates early I‐waves (i.e., I1‐3) (Di Lazzaro et al., 2003). When the current is reversed in an anterior‐medial to posterior‐lateral (AP) direction, there is preferential activation of late I‐waves (i.e., I3‐5) (Hamada, Murase, Hasan, Balaratnam, & Rothwell, 2013; Di Lazzaro et al., 2001). In the present study a biphasic coil was used for cortical mapping. Biphasic coils have been shown to be beneficial in nTMS mapping as they produce more reliable cortical maps (Groppa et al., 2012; Rossini et al., 2015). The orientation of the induced current was PA for the first phase and AP for the second phase of the stimulus. In this configuration, early I‐waves of the PA pulse will initiate the MEP, whereas the later I‐waves of the AP pulse will arrive during the later period of the MEP cycle (Opie & Semmler, 2021; Ziemann, 2020).

Early and late I‐wave differ significantly in their physiology (Opie & Semmler, 2021). I‐wave recruitment has been attributed to several mechanisms, that is, intrinsic membrane properties of corticospinal neurons (Rusu, Murakami, Ziemann, & Triesch, 2014), or from backpropagating calcium action potentials from the corticospinal apical dendrites (Ugawa, Rothwell, & Paulus, 2020). The contributions of GABAergic inhibitory interneurons, however, have been specifically attributed to late I‐wave generation (Di Lazzaro, Rothwell, & Capogna, 2018; Di Lazzaro & Rothwell, 2014; Ziemann, 2020). Notably, GABAergic and AP‐TMS sensitive late I‐wave circuits have been shown to contribute to the fine‐tuning of muscle activation (Federico & Perez, 2017; Jo, Di Lazzaro, & Perez, 2018; Jo & Perez, 2019; Zoghi, Pearce, & Nordstrom, 2003). We hypothesize, that glial brain tumors might provoke the observed effects on MEP time‐frequency representation by impairing transsynaptic recruitment of cortico‐motoneuronal cells in the late I‐wave generation (see phase asynchrony of >150 Hz). In line with previous studies (Federico & Perez, 2017; Jo et al., 2018; Jo & Perez, 2019; Zoghi et al., 2003), this impairment was associated with the occurrence of FMD. In contrast, structural damage to horizontal or vertical white matter tracts (i.e., the CST as measured by the DTI tractography) might result in a general reduction of I‐waves represented by an ERSP reduction.

In cortical physiology, phase synchronization, that is, the alignment of the phase of the ongoing activity to an external input, is considered as a fundamental processing mechanism in event‐related potential (ERP) generation after an input to a neuronal network (Canavier, 2015; Ding & Simon, 2013; Makeig et al., 2004; Palva, Palva, & Kaila, 2005; Shah et al., 2004). In general, ERPs are suggested to evolve from neuronal background activity by (a) adding (i.e., evoking) stimulus‐locked activity recruited by the stimulus and/or (b) aligning the phase of ongoing oscillations to the stimulus (i.e., phase resetting) (Ding & Simon, 2013; Makeig et al., 2004). Phase alignment has been shown to provide a timing frame for gating sensory input in the context of visual and auditory processing, attention, spatial navigation, or memory (Voloh & Womelsdorf, 2016). So far, there is little theoretical workup on the relation of MEP to these models. We conjecture that after biphasic TMS to the motor cortex early I‐waves initiate (i.e., evoke) the MEP while late I‐waves might align the ongoing EMG phase. Temporal dispersion of CSV will reduce spectral (i.e., ERSP) and phase (i.e., ITC) synchronization of both components (Ding & Simon, 2013; Makeig et al., 2004). However, ITC is much more sensitive to stimulus‐synchronized activity than power (Ding & Simon, 2013; Makeig et al., 2004). In line, it is the ITC discriminating different degrees of motor dysfunction in the present study.

Notably, AED intake (i.e., levetiracetam) was associated with an increase of high frequencies ERSP without improving ITC. While the exact mechanism of levetiracetam is unclear, the GABAergic effect is suggested to restore the lost balance between the excitatory and inhibitory systems (Luz Adriana et al., 2018). Our data indicates that the GABAergic effect expands the late I‐wave network (Di Lazzaro et al., 2018; Di Lazzaro & Rothwell, 2014; Ziemann, 2020) without improving its temporal precision.

Conventional MEP time‐series parameters (i.e., onset latency, amplitude and RMT) did not depict any differences between healthy subjects and brain tumor patients or between the patient groups. It is well known that MEP time‐domain characteristics are ambiguous in brain tumor patients. MEP after TMS to the lesioned hemisphere are often indistinguishable in intrasubject (i.e., non‐lesioned hemisphere) or intersubject (i.e., healthy subjects) comparison studies (Picht et al., 2012; Sollmann et al., 2017). In contrast, time‐frequency MEP parameter have been shown to exert a higher intertrial reliability than conventional MEP characteristics (Machetanz, Gallotti, et al., 2021; Machetanz, Wiesinger, et al., 2021). Thus, time‐frequency analysis might increase the sensitivity of TMS in detecting an affection of corticospinal transmission (Machetanz, Gallotti, et al., 2021; Machetanz, Wiesinger, et al., 2021). This is an important aspect in the neurosurgical treatment of brain tumors as preoperative affection of the CST is a major predictor of the postoperative motor outcome (Rosenstock et al., 2017; Rosenstock et al., 2021). Additionally, time‐frequency‐based motor mapping could assist in delineating the “cortical origin” of the CST. There is good evidence that information derived from TMS‐based tractography helps to prevent postoperative deficits after neurosurgical interventions, such as brain (Rosenstock, Grittner, et al., 2017; Rosenstock et al., 2021; Sollmann et al., 2018) or brainstem surgeries (Zdunczyk, Roth, Picht, & Vajkoczy, 2021). Furthermore, postoperative information about the CST status could improve the prediction of functional outcome after brain surgery (Ille et al., 2021; Rosenberg et al., 2010). In line, time‐frequency characteristics of MEP have been shown to improve intraoperative neuromonitoring for prevention and predict functional outcome of postoperative motor deficits after neurosurgical interventions (Hu, Luk, Lu, Holmes, & Leong, 2001; Hu, Luk, Lu, & Leong, 2003).

Limitations of the study

There are several limitations of the study that should be explicitly addressed. Brain tumor patients are a very inhomogeneous cohort concerning molecular tumor characteristics, tumor location, tumor size and the resulting affection of cortico‐cortical and corticospinal connections. Additionally, the mass effect of brain tumors might change the anatomical microstructure without directly interfering with but shifting cortico‐cortical and corticospinal fibers. Thus, TMS might result in different MEP features in glioma patients due to slight variations of the magnetic field in relation to the underlying neuronal elements stimulated by TMS (Bashir, Perez, Horvath, & Pascual‐Leone, 2013). Controlling TMS location and orientation in space will not dissolve this problem completely. Additional studies are needed to evaluate these factors and to clarify whether time‐frequency analysis can help to address and distinguish these effects. Maybe, the inclusion of information about patient's white matter microstructure (e.g., DTI) might help to explore this problem in the future.

CONCLUSION

Brain lesions such as tumors affect the temporal precision of CSV in patients with impaired motor performance. Phase synchronization analysis provides a more robust technique for detecting temporal distortion of CSV than standard MEP analyses, thus, providing more precise information about the status of the motor system.

CONFLICT OF INTEREST

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

AUTHOR CONTRIBUTIONS

Georgios Naros was responsible for the conception and design, data acquisition, analysis and interpretation as well as writing and reviewing of the manuscript. Kathrin Machetanz contributed to the acquisition, analysis and interpretation of data and writing of the first draft. Maria Teresa Leao, Sophie Wang, and Marcos Tatagiba contributed to the data acquisition, interpretation of data and the review and critique of the final manuscript. Alireza Gharabaghi contributed to the interpretation of data and reviewing of the manuscript.

ETHICS STATEMENT

The study was approved by the local ethics committee of the Medical Faculty of the Eberhard Karls University Tuebingen.

PATIENT CONSENT

All participants gave written informed consent.

PERMISSIONS

There is no material from other sources.