diff --git a/text/thesis/02MaterialsAndMethods.tex b/text/thesis/02MaterialsAndMethods.tex index 57b8032..fce3a3b 100644 --- a/text/thesis/02MaterialsAndMethods.tex +++ b/text/thesis/02MaterialsAndMethods.tex @@ -27,26 +27,25 @@ \label{fig:10-20} \end{figure} \subsubsection{Frequency bands} - EEG can be divided into several bands. According to \cite{Gerrard07} we can propose following bands: + EEG can be divided into several bands according to the frequency of the waves. Usually this also corresponds to the form of waves and allows insights in the mental state from rest to highly aroused. According to \cite{Gerrard07} we can propose following bands: \begin{itemize} \item Delta: 1-4 Hz \item Theta: 4-7 Hz \item Alpha: 7.5-12.5 Hz (depending also on age) \item Beta: 13-20Hz \end{itemize} - They were named almost alphabetically concerning their discovery, why they seem assorted sometimes. %TODO: naja + There are different definitions of the limits of the bands, as we only use them for rough estimation we stick on these. For more exact results an analysis of wave patterns would be necessary. \subsection{Power estimation} \subsubsection{EEG} To use data from EEG one way is to analyze the occurring frequencies and their respective power.\\ To gain these from the continuous signal there are different methods. The intuitive approach would be to use Fourier transformation however the Fourier transform does not need to exists for a continuous signal. So we used power spectral density (PSD) estimation. \subsubsection{Power spectral density estimation} - The PSD is the power per frequency. Power here refers to the square of the amplitude. %TODO: formulation,fft + The PSD is the power per frequency. Power here refers to the square of the amplitude.\\ If the Fourier transform is existing, PSD can be calculated from it e.g. as periodogram. If not it has to be estimated. One way to do so is parametrized with an Autoregressive model(AR). Here one assumes that the there is a correlation of the spectral density between $p$ consecutive samples and the following one. This leads to an equation with only $p$ parameters which can be estimated in different ways. We used Burg's method (\texttt{pburg} from \matlab{} library).\\ In Figure~\ref{fig:psd} we see the difference between autoregressive \texttt{pburg} and periodogram \texttt{pwelch} PSD estimation. \begin{figure} \includegraphics[width=\textwidth]{psd.png} - \caption{PSD with Autoregressive model and FFT respectively\protect\footnotemark} - + \caption{PSD with FFT or an Autoregressive model respectively\protect\footnotemark} \label{fig:psd} \end{figure} \footnotetext{The signal was unfiltered EEG data from channel \textit{Cz} second run of second session with subject AO} @@ -64,7 +63,6 @@ In the 2000s there began a movement using new techniques to record ultrafast and infraslow brainwaves (above 50Hz and below 1Hz). These were found to have some importance (cf. \cite{Vanhatalo04}).\\ Also in predicting movements there was found some significance in low frequency as was done by Liu et al. (\cite{Liu11}) and Antelis et al. (\cite{Antelis13}) for example. Antelis et al. found correlations between hand movement and low frequency signal of $(0.29,0.15,0.37)$ in the dimensions respectively.\\ Lew et al. (\cite{Lew14}) state low frequencies are mainly involved in spontaneous self-induced movement and can be found before the movement starts. By this they may be a great possibility to lower reaction time of neuroprostheses for example. - %TODO: more details (idea, possible explanantion) \subsection{Filtering} Filtering of the recorded EEG signal is necessary for different reasons. First there are current relics from 50Hz current. These can be filtered out with bandstop filters.\\ Secondly we need to concentrate on the interesting frequencies (for classical EEG 1-50Hz). This is done by applying lowpass or highpass filters respectively. This is necessary because the PSD of lower frequency is a lot higher than that of higher frequencies. The relation $$PSD(f)=\frac{c}{f^\gamma}$$ holds for constants $c$ and $\gamma$ (\cite{Demanuele07}).\\ @@ -73,7 +71,7 @@ When using muscles they are contracted after an signal via an efferent nerve activates them. Contraction of muscles also releases measurable energy which is used for Electromyography (EMG). There are intramuscular applications of EMG but we only used surface EMG.\\ From surface EMG activity of muscles can be estimated however not very precisely without repetition. Since the muscles used for arm movements are quite large in our setting EMG allows relatively precise estimations of underlying muscle activity. - EMG is mainly developed for diagnostic tasks. However it is also applicable in science to track muscle activity. %TODO? + EMG is mainly developed for diagnostic tasks. However it is also applicable in science to track muscle activity. \subsection{Synergies} \label{back:synergies} Movement of the arm (and other parts of the body) are under-determined meaning with given trajectory there are different muscle contractions possible. One idea how this problem could be solved by our nervous system are synergies. Proposed by Bernstein in 1967 (\cite{Bernstein67}) they describe the goal of the movement (e.g. the trajectory) instead of controlling single muscles. This would mean however that predicting the activity of single muscles from EEG is harder than predicting a synergy which in turn determines the contraction of muscles.\\ @@ -91,15 +89,14 @@ \label{fig:pca} \end{figure} In Figure~\ref{fig:pca} we see the eigenvectors of the data. The longer vector is the principal component the shorter one is orthogonal to it and explains the remaining variance. The second component here also is the component which explains least variance, since most variance is orthogonal to it. - %TODO? formula, ... \subsection{NMF} \label{mat:nmf} In some applications Non-negative Matrix Factorization (NMF) is preferred over PCA (cf. \cite{Lee99}). This is because it does not learn eigenvectors but decomposes the input into parts which are all possibly used in the input. When seen as matrix factorization PCA yields matrices of arbitrary sign where one represents the eigenvectors the other the specific mixture of them. Because an entry may be negative cancellation is possible. This leads to unintuitive representation in the first matrix.\\ NMF in contrast only allows positive entries. This leads to \qq{what is in, is in} meaning no cancellation which in turn yields more intuitive matrices. The first contains possible parts of the data, the second how strongly they are represented in the current input.\\ The formula for NMF is $$Input\approx \mathbf{WH}$$ - where Input is $n\times m$, $W$ is $n\times r$ and $H$ is $r\times m$ with $r<<\min\{m,n\}$. So $\mathbf{WH}$ is only an approximation of the input however with significant lower dimension - the number of Synergies used.\\ %TODO: formulation - The factorisation is learnt with an update rule that my be chosen. We used the \matlab{} default, an alternating least squares (ALS) algorithm. It can be described as in algorithm \ref{alg:als} (cf. \cite{Berry07}).\\ + where Input is $n\times m$, $W$ is $n\times r$ and $H$ is $r\times m$ with $r<<\min\{m,n\}$. So $\mathbf{WH}$ is only an approximation of the input however with significant lower dimension - the number of Synergies used.\\ + The factorization is learned with an update rule that my be chosen. We used the \matlab{} default, an alternating least squares (ALS) algorithm. It can be described as in algorithm \ref{alg:als} (cf. \cite{Berry07}).\\ \begin{algorithm} \begin{algorithmic} \State initialize $\mathbf{W}$ randomly @@ -241,9 +238,9 @@ \label{mm:design} The data used for this work were mainly recorded by Farid Shiman, Nerea Irastorza-Landa, and Andrea Sarasola-Sanz for their work (\cite{Shiman15},\cite{Sarasola15}). We were allowed to use them for further analysis.\\ There were 9 right-handed subjects with an average age of 25 (variance 6.67, minimum 20, maximum 28). Three female and 6 male subjects were tested. All the tasks were performed with the dominant right hand.\\ - To perform was a center-out reaching task to one of four targets (see \ref{fig:experimentalDesign}) while 32 channel EEG, at least% + To perform was a center-out reaching task to one of four targets (see \ref{fig:experimentalDesign}) while 32 channel EEG, 6 channel% \footnote{\texttt{'AbdPolLo', 'Biceps', 'Triceps', 'FrontDelt', 'MidDelt'} and \texttt{'BackDelt'} were recorded for every subject, others only in some. Only the 6 channels tracked in every session were used} % - 6 channel surface EMG and 7 DOF kinematics were tracked. + surface EMG and 7 DOF kinematics were tracked. \begin{figure} \centering \includegraphics{experimentalDesign.jpg} @@ -385,8 +382,7 @@ The kinematic record was started after the EEG recording. In synchronization channel\footnote{cf. Table~\ref{tab:channelNames}} there is a peak when kinematic recording is started. This was used to align movement with EEG and EMG data. In addition we adjusted the kinematic data to the EMG window and shift to be able to use corresponding data for the same time step. This was done by summing all differences (for movement) or by calculating the mean position in the time window.\\ Size of this data is same as EMG and Synergies in length but has only three features per time step since we used only 3D positioning ($x,y$ and $\theta$) of the hand and no information about the fingers. \section{Data Analysis} - %\subsection{Overview} - In Figure~\ref{fig:overview} is shown what is predicted from where. + Figure~\ref{fig:overview} shows the steps of our work. EEG, EMG and positions were recorded, Synergies and velocities were calculated from them. To check the performance of our methods the relations between them were predicted. \begin{figure} \centering \input{pictures/overview.tikz} diff --git a/text/thesis/04Discussion.tex b/text/thesis/04Discussion.tex index 75f38a1..4d71a2e 100644 --- a/text/thesis/04Discussion.tex +++ b/text/thesis/04Discussion.tex @@ -56,8 +56,8 @@ %TODO \subsection{Prediction via Synergies} Of course the prediction via Synergies is a bit worse than direct prediction, since the machine learning techniques could do the same dimensionality reduction and also much more.\\ - Nevertheless we also get good correlations when predicting from synergies meaning the model may match the reality. %TODO formulierung? - In addition synergies are predicted significantly ($p<0.05$) better from EEG than EMG. So the representation as synergies probably matches the representation in the brain better. This could mean that the controlling of a prostheses should be done via synergies - representing the representation in the brain and being easier to implement than a prosthesis listening to 32 EEG channels. + This decrease however is not large which suggests that synergies are a valid step in between.\\ + In addition the prediction of synergies from EEG are significantly ($p<0.05$) better than the prediction of EMG. So the representation as synergies probably matches the representation in the brain better. This could mean that the controlling of a prostheses should be done via synergies - representing the representation in the brain and being easier to implement than a prosthesis listening to 32 EEG channels. \subsection{Comparison with EMG} The results show that the dimensionality reduction from 6 dimensional EMG to 3 dimensional Synergies (here via autoencoder) does not cost much information when predicting velocities and positions.\\ For velocities there is no significant difference and even for positions the mean only differs about $0.03$ (EMG: $0.23$, Autoencoder: $0.20$). diff --git a/text/thesis/Acd.tex b/text/thesis/Acd.tex index a2801bb..2e63d85 100644 --- a/text/thesis/Acd.tex +++ b/text/thesis/Acd.tex @@ -1,2 +1,45 @@ \chapter{Contents on CD} \label{app:cd} +\begin{forest} + for tree={ + font=\ttfamily, + grow'=0, + child anchor=west, + parent anchor=south, + anchor=west, + calign=first, + edge path={ + \noexpand\path [draw, \forestoption{edge}] + (!u.south west) +(7.5pt,0) |- node[fill,inner sep=1.25pt] {} (.child anchor)\forestoption{edge label}; + }, + before typesetting nodes={ + if n=1 + {insert before={[,phantom]}} + {} + }, + fit=band, + before computing xy={l=15pt}, + } +[CD + [matlabCode + [topoplot + [plotOneSubjectOneDay.m] + [topoplot\_martin.m] + [topoplot\_Wrapper.m] + ] + [T + [O + [DO] + ] + ] + ] + [bashscripts + [noSynergies.bash] + [runOnEciton.bash] + [runOnEcitonPos.bash] + ] + [thesis + [thesis.pdf] + ] +] +\end{forest} diff --git a/text/thesis/Bfunctions.tex b/text/thesis/Bfunctions.tex index 21fbf23..47b5560 100644 --- a/text/thesis/Bfunctions.tex +++ b/text/thesis/Bfunctions.tex @@ -22,6 +22,16 @@ For each data point it is decided whether it belongs to movement or not according to the given threshold in EMG activity. If not there is no movement so the class is set to 0 (rest).\\ If there is movement it is decided whether there should be according to the given classification. If not the old class is applied also for this point. If yes it is checked whether the movement just started (up to now class was 0). If the movement just started, 1 second before is taken out (pause \true) or half second before is classified same, 0.5s to 1s is dropped (pause \false). + \subsection{\texttt{generateTrainingData.m}} + \label{code:generate} + In this function the transformation of EEG and EMG signals is done as PSD with Burg's method for EEG and waveform length for EMG. Additionally velocities are computed from kinematic data. + \subsection{\texttt{generateTrainingDataPos.m}} + Same as \ref{code:generate} but kinematic data is used as is as positions. + \subsection{\texttt{myDownsample.m}} + \label{code:myDown} + \texttt{myDownsample} takes the given number of equidistant samples to represent the input. + \subsection{\texttt{namesAndNumbers.m}} + \texttt{namesAndNubers} returns names and numbers of subjects and runs according to the file given (created by \ref{code:run.bash}) \subsection{\texttt{readAll.m}} \label{code:readAll} This is the central function for the acquisition of data. @@ -35,15 +45,47 @@ Finally the kinematics and the synergies are generated matching the size of EMG data. All is then saved under the given path as \texttt{.mat} file. \subsection{\texttt{readAllPos.m}} Same as \ref{code:readAll} but using position instead of movement from kinematics data. - \subsection{\texttt{generateTrainingData.m}} - \label{code:generate} - %TODO - \subsection{\texttt{generateTrainingDataPos.m}} - Same as \ref{code:generate} but using position instead of movement from kinematics data. + \subsection{\texttt{readEMGsig.m}} + Reads only the EMG signal (used for \ref{code:noSyn}). + \subsection{\texttt{shiftingDownsample}} + Shifts \texttt{myDownsample} \ref{code:myDown} over a signal. + \subsection{\texttt{shiftingKin.m}} + Calculates velocity in a window + \subsection{\texttt{shiftingMean.m}} + Downsamples by using the mean for a window. + \subsection{\texttt{shiftingPburg.m}} + Applies Pburg to a window. + \subsection{\texttt{shiftingPos.m}} + Calculates position for a window. If there is no data the next existent data will be used. + \subsection{\texttt{shiftingSum.m}} + Downsamples by using the sum in a window. + \subsection{\texttt{waveformLength.m}} + Calculates the waveform length for the given signal. + \subsection{\texttt{waveformLengthAll.m}} + Calculates the waveform length for the whole EMG signal calling \texttt{waveformLength.m} on defined windows \section{Data Analysis} \subsection{\texttt{correlation2.m}} Calculates the squared correlation between each corresponding columns of two matrices.\\ This is used as a measure for fit when comparing the number of synergies. + \subsection{\texttt{kFoldCV.m}} + \label{code:kfoldCV} + \texttt{kfoldCV} does a k-fold Cross-validation of parameter c for SVM. + \subsection{\texttt{kFoldRidge.m}} + \label{code:kFoldRidge} + \texttt{kFoldRidge} does a k-fold Cross-validation of parameter lambda for Ridge. + \subsection{\texttt{noSynergies.m}} + \label{code:noSyn} + Plots performance (squared elementwise correlation) for 1 to 6 synergies for Autoencoder, PCA and NNMF (is called by \ref{code:noSyn.bash}).\\ + Whether the figure is shown or the data is saved can be adjusted by altering the corresponding flag. + \subsection{\texttt{ridgeCorrelation.m}} + Calculates the correlation with original data, according coefficients and the testing data for the prediction. + \subsection{\texttt{ridgeCV.m}} + \label{code:ridgeCV} + Runs a ridge regression to predict 'goal' from 'data'. Therefore calls kFoldRidge.m \ref{code:kFoldRidge}. + \subsection{\texttt{ridgeCVvia.m}} + Same as \ref{code:ridgeCV} but has an intermediate step (usually synergies). + \subsection{\texttt{svmEciton.m}} + Runs a SVM to classify the data. Cross validation is done by calling \texttt{kfoldCV} \ref{code:kfoldCV}. \section{Evaluation} \subsection{\texttt{evaluation.m}} Script for composing the file containing the results (\texttt{evaluation.mat}). This file is used for the evaluations after some manual post-processing (mainly renaming). @@ -51,3 +93,16 @@ Evaluations run on results from SVMs. \subsection{\texttt{evaluationCorrelations.m}} Evaluations run on results from ridge regressions. + \subsection{\texttt{plotCM.m}} + Plots a confusion matrix. Axes are labeled according to the classes of EEG and EMG. +\section{Plots} + \subsection{\texttt{mySaveFigure.m}} + \texttt{mySaveFigure.m} save a given figure with default size at given filename. +\section{Bash Scripts} +\label{code:noSyn.bash} +\label{code:run.bash} +\section{miscellaneous} + \subsection{\texttt{psdPlot.m}} + Plots PSD with pburg and pwelch. + \subsection{\texttt{pickFromStruct.m}} + Allows to pick $i$th entry of a struct (without knowing the name). diff --git a/text/thesis/thesis.tex b/text/thesis/thesis.tex index f3726cd..e0b0500 100644 --- a/text/thesis/thesis.tex +++ b/text/thesis/thesis.tex @@ -8,6 +8,7 @@ \documentclass[twoside,12pt,a4paper]{book} %\usepackage{reportpage} +\usepackage{forest} \usepackage{natbib} \usepackage{epsf} \usepackage{graphics, graphicx} @@ -29,8 +30,6 @@ \usetikzlibrary{calc,positioning,shadows,arrows,automata} -\newcommand{\qq}[1]{``#1''} -\renewcommand{\exp}[1]{\ensuremath{\cdot10^{#1}}} \textwidth 14cm \textheight 22cm @@ -56,6 +55,8 @@ \newcommand{\matlab}{\textsc{Matlab}} \newcommand{\true}{\texttt{true}} \newcommand{\false}{\texttt{false}} +\newcommand{\qq}[1]{``#1''} +\renewcommand{\exp}[1]{\ensuremath{\cdot10^{#1}}} \begin{document} @@ -116,7 +117,7 @@ \newpage -\section*{Acknowledgments} +\section*{Acknowledgments} %TODO : morgen nochmal lesen \addcontentsline{toc}{section}{Acknowledgments} \noindent First of all I would like to thank my supervisor Martin Spüler who always had an open ear for problems and gave lots of useful input. %TODO diff --git a/usedMcode/myDownsample.m b/usedMcode/myDownsample.m index 8113abd..14141e7 100644 --- a/usedMcode/myDownsample.m +++ b/usedMcode/myDownsample.m @@ -1,5 +1,5 @@ function [out]=myDownsample(in,noSamples) -% myDownsample samples in as noSamples equidistant values +% myDownsample samples as noSamples equidistant values stepSize=max(size(in))/noSamples; out=in(fix((0:noSamples-1)*stepSize)+1); end \ No newline at end of file diff --git a/usedMcode/shiftingPos.m b/usedMcode/shiftingPos.m index 901c555..c8701a1 100644 --- a/usedMcode/shiftingPos.m +++ b/usedMcode/shiftingPos.m @@ -4,7 +4,7 @@ for j=1:size(kinPerSec,1) tmp=mean(kin(kin(:,1)>(j-1)*shiftEMG*1000 & kin(:,1)<=(j-1)*shiftEMG*1000+windowEMG*1000,2:4)); i=j; - while isnan(tmp) %interval is empty + while isnan(tmp) %interval is empty, next existen data will be used tmp=mean(kin(kin(:,1)>i*shiftEMG*1000 & kin(:,1)<=i*shiftEMG*1000+windowEMG*1000,2:4)); i=i+1; end