Complexity

Volume 2018, Article ID 7269494, 15 pages

https://doi.org/10.1155/2018/7269494

## A Wavelet-Based Correlation Analysis Framework to Study Cerebromuscular Activity in Essential Tremor

^{1}Through-life Engineering Services Centre, Cranfield University, Bedfordshire MK43 0AL, UK^{2}Computational and Software Techniques in Engineering, Cranfield University, Bedfordshire MK43 0AL, UK^{3}Cixi Institute of Biomedical Engineering, Ningbo Institute of Industrial Technology, Chinese Academy of Sciences, Ningbo, China^{4}Information Engineering, Zhejiang University of Technology, Hangzhou 310023, China^{5}Department of Neurosciences, Sheffield Teaching Hospitals, NHS Foundation Trust, Royal Hallamshire Hospital, Sheffield S10 2JF, UK

Correspondence should be addressed to Yifan Zhao; ku.ca.dleifnarc@oahz.nafiy and Yitian Zhao; nc.ude.tib@oahz.naitiy

Received 12 February 2018; Revised 23 April 2018; Accepted 22 May 2018; Published 3 July 2018

Academic Editor: Nadia Mammone

Copyright © 2018 Yifan Zhao et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### Abstract

*Objective*. Deep brain stimulation (DBS) provides dramatic tremor relief in patients with severe essential tremor (ET). Typically, the VIM nucleus is the most effective brain area to target for high-frequency electrical stimulation in these patients. Correlation analysis between electrical local field potential (LFP) recordings from the thalamic DBS leads and electrical muscle activity from the contralateral tremulous limb has become an attractive practical tool to interpret the LFPs and their association with the tremulous clinical manifestations. Although functional connectivity analysis between brain electrical recordings and electromyographic (EMG) signals from the tremor has been of interest to an increasing number of engineering researchers, there is no well-accepted tailored framework to consistently characterise the association between thalamic electrical recordings and the tremorogenic EMG activity. *Methods*. This paper proposes a novel framework to address this challenge, including an estimation of the interaction strength using wavelet cross-spectrum and phase lag index while demonstrating the statistical significance of the findings. *Results*. Consistent results were estimated for single and multiple trials of consecutive or partially overlapping epochs of data. The latter approach reveals a substantial increase on the range of statistically significant dynamic low-frequency interrelationships while decreasing the dynamic range of high-frequency interactions. *Conclusion*. Results from both simulation and real data demonstrate the feasibility and robustness of the proposed framework. *Significance*. This study offers the proof of principle required to implement this methodology to uncover VIM thalamic LFP-EMG interactions for (i) better understanding of the pathophysiology of tremor; (ii) objective selection of the DBS electrode contacts with the highest strength of association with the tremorogenic EMG, a particularly useful feature for the implementation of novel multicontact directional leads in clinical practice; and (iii) future research on DBS closed-loop devices.

#### 1. Introduction

Essential tremor (ET) is defined as a neurological disorder that causes involuntary abnormal repetitive shaking. This shaking can appear in different parts of the body such as the hands, forearms, or head [1]. ET is a common movement disorder affecting around four out of 100 adults over 40 years of age. It is considered to be a centrally driven tremor, and the constituents of the network comprise anatomical areas of the physiological motor system [2]. Treatment of ET is mainly based on pharmacotherapy and surgery for medically refractory cases [3]. The main surgical approach currently in use consists of continuous deep brain stimulation (DBS) through implantation of a depth electrode in the area of the ventral intermediate (VIM) nucleus of the thalamus, a key area in the neuronal loop generating the tremor [3, 4]. The depth electrodes after DBS surgery emit continuous electrical signals coming from the neurostimulator or pacemaker [5]. To confirm electrode position within the VIM thalamic nucleus, in addition to neuronavigation techniques to target this specific anatomical area, some centres perform microelectrode recordings during surgery to detect tremor-related electric neuronal bursts and find the ideal position for electrode implantation [6]. Through the routinely used macroelectrodes in DBS surgery, neuronal network electrical activity can be recorded from the VIM thalamus in the form of local field potentials (LFPs) a few days after surgery.

#### 2. Material

In this study, we use recordings from two patients (aged 64 and 53, both female) who underwent DBS surgery for medically refractory ET. Ethics approval for use of patients’ EEGs and LFPs for the development of new quantitative EEG (qEEG) methods was obtained both from the University of Sheffield and the NHS ethics committees (SMBRER207 and 11/YH/0414). The multichannel Natus Quantum Amplifier (Optima Medical Ltd.) at a sampling rate of 16,384 Hz was used for all EEG/LFP/EMG polygraphy recordings (analogue bandwidth 0.01–4000 Hz).

In our institution (Royall Hallamshire Hospital, Sheffield Teaching Hospitals, NHS Foundation Trust), every patient undergoing DBS surgery for tremor is offered, five days after implantation of the depth electrodes, a comprehensive electrophysiological analysis to include LFPs from the VIM thalamus, electroencephalography (EEG), and electromyographic (EMG) polygraphy recordings, measuring the electrical activity from the muscles, to look at the correlation between the tremor and the thalamic network oscillations. The cortical scalp EEG recordings are not used in this work. There are four electrode contacts on each macroelectrode placed in the thalamus in close proximity to each other (contacts 0, 1, 2, and 3) through which bipolar recordings of LFPs can be obtained. Many centres use an empirical approach in selecting the ideal electrode contact to stimulate after DBS surgery to obtain the ideal tremor suppression. This approach can be time-consuming, and it will become more of a problem as new multicontact directional leads make their way into clinical practice [7]. Additionally, it seems that the LFPs offer very significant information which in the future could be used to optimise clinical outcomes in closed-loop systems [8]. However, before this becomes possible, an objective measure of the strength of association between the VIM thalamic network and the tremor recorded on EMG, or the equivalent mechanical oscillations through accelerometers, is required. One practical solution is based on the correlation analysis between data coming from the DBS leads (LFP) and electrical muscular activity (EMG) of the tremulous limb. Corticomuscular functional connectivity is defined as the interaction between the electrical activity of the cortex in the brain or the thalamus in this instance and the electrical activity recorded from various muscles. During the last few years, a few methods have been developed to understand this type of functional connectivity [9–11].

From the engineering perspective, the approaches to study connectivity between two signals can be classified as time- and frequency-domain-based methods. One of the classic linear measures to estimate similarity in time domain is cross-correlation that has been used to study ET [12]. It measures the degree of connectivity when a time series is shifted from other reference series. This method is also useful to measure the time lag between two signals. However, this method usually assumes that the system is linear and stationary. Another branch to quantify the correlation or causality between signals in time domain is by mutually predicting selected observable measurements based on multivariate modelling, where the best-established methods are based on the Granger causality test [13] that has been used to study the correlation between EEG signals [14]. They are based on parametric modelling. A full and unbiased model is therefore required, which can be a challenge to achieve due to the limited knowledge on the human brain. Other nonparametric methods include mutual information [15] and transfer entropy [16, 17], which are model-free but usually require larger datasets or averaging over many realisations to mitigate the effects of noise. In the scope of frequency-based methods, cross-spectrum allows determining the connection between two stationary signals in terms of frequency [18]. When computing complex cross-spectrum, results can be divided in cospectrum (in-phase connectivity) and quad-spectrum (out-of-phase connectivity). Coherence [19, 20] uses a normalisation of cross-spectrum values that takes unit value when total linear phase relationship is detected or zero value for independent signals. Both methods are easy-to-use and computationally inexpensive, but the investigated correlation must be stationary and linear and there is no-coupling information among various frequencies. Cross-bispectrum [21] is used to detect the quadratic phase coupling (QPC) between frequency components of two target signals. One pitfall of this algorithm is that bispectrum is affected not only by nonlinearity but also by non-Gaussian data. Thus, it requires the data to present Gaussian distribution to make sure to detect nonlinearity. A suitable approach to overcome these non-phase-coupling peaks is through cross-bicoherence [22]. However, this method only favours the strongly phase-coupled signals, the ones that show QPC interactions.

In addition to the frequency- and time-based approaches, efforts have been made with the wavelet domain to study connectivity of brain networks. This approach aims to address limitations of the above methods on tackling dynamic systems by providing time-resolving value with accurate locality. Meanwhile, it is a model-free (nonparametric) measure, which reduces the requirement of a priori knowledge of the underlying model. Wavelet coherence has attracted increased interest on studying brain-related disorders. Jeong et al. [23] proposed to use wavelet energy and wavelet coherence as EEG biomarkers to distinguish Parkinson’s disease and Alzheimer’s disease. A wavelet coherence-based clustering of EEG signals has been developed to estimate the brain connectivity in absence epileptic patients [24]. It has also been applied to studies on autism [25], traumatic brain injury [26], schizophrenia [27], and poor sleep quality [28]. However, there is very limited research focusing on its application on ET. Furthermore, for the approaches based on wavelet coherence to understand brain connectivity, there are no well-accepted complete frameworks to understand cerebromuscular connectivity systematically.

Addressing these challenges, this paper develops a new wavelet-based correlation analysis framework combined by the estimation of connectivity strength, significance test, and phase-delay characterisation. It aims to better understand cerebromuscular interactions in a structured manner. It should be noted that this framework has the prospect to be applied on other applications of connectivity analysis, such as EEG-EEG and EEG-EMG.

#### 3. Methods

Model-free (nonparametric) measures are chosen in this paper to study the correlation between LFPs and EMG because this biological system is so complex that it would require a substantial number of parameters and computation time to build a satisfactory parametric model. Additionally, there is no well-accepted analytical model to start with. A common deficiency when applying biomedical signal-processing tools is the assumption of stationarity. A typical practice for the estimation of spectrum distribution is segmenting long records of data and averaging calculations over segments. The main disadvantage of this practice is the inability of having time-resolved values. There is therefore no time resolution, and the dynamic behaviour of neuronal interactions cannot be revealed. One solution to overcome this issue is the use of wavelets.

Wavelet transformation makes a decomposition of a time series into a frequency-time domain. It uses convolution of a mother wavelet and its scaled and shifted versions. Among all the possible mother wavelets, the *Morlet* wavelet provides a good balance between frequency and time resolution and has been widely used in EEG and EMG research [29]. This study is restricted to this wavelet. Equation (1) corresponds with a normalised *Morlet* wavelet where the dimensionless frequency is set as 6 and dimensionless time is denoted by . Considering an observed time-serial , in continuous wavelet transform (CWT), the mother wavelet is modified by varying the scale , as shown in (2), so that with a discretised time domain of time step .

##### 3.1. Wavelet Coherence

Coherence is one of the most widely used methods for measuring linear interactions. It is based on the Pearson correlation coefficient used in statistics but in frequency and time domain. It measures the mean resultant vector length (or consistency) of the cross-spectral density between two signals. Its squared value varies from 0 to 1, meaning low and high linear frequential correlation. During this study, coherence is used as a reference standard for comparison to other methods. The wavelet formulation of coherence between two signals, and , and in the frequency and time domain, can be formulated as where is the wavelet cross-spectrum between and and and are the corresponding autospectrums. Working with two single signals (single realisation) usually requires using a smoothing operator (see operator in (4)), and ergodicity properties should be assumed [30].

If multiple trials of both signals are available, square coherence can be estimated using (5), which is used in this study. where

The number of trails is denoted by .

To assert significant values of wavelet coherence, the statistical methodology stablished by Gallego et al. [31] is employed. It estimates a threshold based on the null hypothesis () of independency by analytically calculating the statistical distribution of coherence. Specifically, assumes that both signals are independent Gaussian variables. Under , with a specified probability , where is calculated as

In this paper, the parameter is set as a fixed value of , equivalent to a 95% of confidence interval.

##### 3.2. Wavelet Cross-Spectrum

Wavelet cross-spectrum (WCS) also provides information about linear synchronisation, but its values are not normalised as in wavelet coherence. Its calculation is written in (6) for multiple trials. Its values seem difficult to be interpreted by statisticians and complementary plots, such as coherence or coquadrature, that can help understand the frequential relationships between signals [9]. That is why few neurological connectivity studies have used this method. However, this paper proposes to use an appropriate significant test that allows WCS results to be judged with more clarity by including the variance of each signal.

It has been proven by Bigot et al*.* [32] that including the variance spectrum data (autospectrum) and benefits the interpretation of WCS in contrast to wavelet coherence. Its performance is based, analogously to coherence, on a threshold calculation under a null hypothesis. However, this makes the statistical procedure using a combination of parametric and nonparametric estimations. It has parametric characteristics since it stablishes the signals as independent Gaussian vectors. But its distribution is considered to show a zero-mean value and a general covariance matrix and estimated by the sampled data. Hence, it does not make any parametric assumption on the covariance. The statistical test stablishes, under and , that . The threshold , for each time and frequency point, is defined as
where and are the largest eigenvalues of the empirical covariance matrix of the time series *x* and *y*, respectively. The symbol means the energy of the wavelet, but since wavelet-normalised values are considered, this variable is omitted. This paper considers a significant probability of 95% with .

##### 3.3. Phase Lag Characterisation

Apart from the quantification of the linear interaction strength and associated significance, it is also important to measure the time or phase lag between signals at a certain frequency. It is particularly important for studying corticomuscular interactions as they carry significant time delays.

The most straightforward approach would be calculating the time lag from the phase information in the complex data of the WCS results with a simple fraction where is the phase difference between both signals at a specific time and frequency and is the sample rate. However, volume conduction can cause the coherence and the phase-locking value to spuriously increase, and (9) is not effective to deal with this kind of common noise sources. To overcome this problem, a measure called the weighted phase lag index (WPLI) algorithm [33] that describes the consistency of the phase difference or time lag is used. This value will inform if the phase lag existent between both signals is consistent. This method not only highlights frequency bands where the phase difference is constant, indicating a strong linear relationship, but also penalises those synchronisations where time lag is close to zero, thus avoiding possible artefacts coming from common noise sources. The weight can be illustrated by Figure 1, which clearly shows the WPLI weights cross-spectra according to the magnitude of the imaginary component of the cross-spectrum. Cross-spectra around the real axis contribute to a less extent than cross-spectra around the imaginary axis. Their values are calculated as where is the imaginary part of the argument and is the cross-spectrum value of the th trial of the total number of trials (). The value is normalised from −1 to 1. As indicated by Figure 1, a value of 1 suggests that the phase lag is 90° while a value of 0 corresponds to 0° or 180°. It should be noted that this process should be applied before the significance test.