Repositorio Institucional de Documentos

Proyecto Fin de Carrera Ingenier´ıa en Inform´atica CLASIFICACI´ ON DE LOS ESTADOS DE REPOSO, PREPARACI´ ON Y EJECUCI´ ON DEL MOVIMIENTO DEL BRAZO A TRAV´ ES DEL ELECTROENCEFALOGRAMA Autor: Miguel Rodrigo G´omez Directores: Luis Montesano del Campo Javier M´ınguez Zafra Departamento de Inform´atica e Ingenier´ıa de Sistemas Zaragoza, Septiembre 2011 Agradecimientos En primer lugar, me gustar´ıa agradecer a mi director del proyecto, Luis Montesano, y mi codirector, Javier M´ınguez, por todo su esfuerzo, dedicaci´on, paciencia y por lo much´ısimo que he aprendido a lo largo de todo el tiempo dedicado al proyecto, adem´as por haberme dado la oportunidad de participar en la elaboraci´on de un art´ıculo para la comunidad cient´ıfica. En segundo lugar, dar las gracias a todos los trabajadores de la empresa BitBrain: Mar´ıa, Soraya, M´onica, Isabel, Escolano, Dani, Marco, Edu y Sergio, por todo el feedback proporcionado, as´ı como por esos grandes momentos en su compa˜n´ıa, que dieron lugar a grandes ideas para este proyecto. Y al Grupo de Rob´otica de la universidad de Zaragoza por los datos suministrados, en especial a Mauricio, por el tiempo dedicado a ense˜narme el trabajo realizado. A todos los viejos amigos de C´aceres, en especial David Rivas, por haberme hecho la persona que soy ahora, ya que sin su apoyo y su amistad no habr´ıa llegado al final del camino. A los nuevos amigos de Zaragoza, los cu´ales han hecho que mi estancia aqu´ı me permitiera seguir adelante y dar lo mejor de m´ı mismo. A todas las personas que de una manera u otra, han estado a mi lado y que me haya podido olvidar de ellas: lo siento, mi memoria ya no es lo que era. A mi familia, por no haber perdido la fe y haber apostado por m´ı, contra viento y marea. Y, por ´ultimo, a Jara Gonz´alez ´ I˜nigo, cambi´o mi vida: le dio un significado mayor del que ya ten´ıa y me hizo coger el impulso necesario para afrontar cualquier reto del camino, sabiendo que lo afrontaremos juntos. 3 Resumen Un inter´es creciente ha aparecido, en los ´ultimos a˜nos, en torno al uso de ”interfaces cerebro computador” aplicadas a terapias de rehabilitaci´on, para ayudar a personas con discapacidad motora a mover sus miembros. Una de las aplicaciones m´as importantes dentro de este contexto, es la de detectar las intenciones del paciente y, gracias a ellas, ajustar mejor las terapias de rehabilitaci´on, as´ı como, un control m´as preciso de pr´otesis rob´oticas. El objetivo del proyecto ser´a detectar los estados de reposo, preparaci´on o anticipaci´on y ejecuci´on del movimiento del brazo, utilizando la se˜nal de electroencefalograma obtenida durante el proceso. Para ello, a partir de la se˜nal original, estudiaremos c´omo distintos filtros espaciales lineales, filtros frecuenciales, rectificado y suavizado de la se˜nal, act´uan en la mejora de la separabilidad de las tres clases. Con las mejores caracter´ısticas seleccionadas, trataremos de detectar, por un lado, los estados por parejas: reposo y preparaci´on, reposo y movimiento y preparaci´on y movimiento, y, por otro lado, las tres clases de manera simult´anea, para lo cu´al entrenaremos dos clasificadores lineales: ”lineal discriminant analysis”,LDA, y ”shrinkage”. Los resultados de la clasificaci´on, en el caso de detectar los estados en parejas, fueron cercanos al 65 % de datos bien clasificados, alcanzando en algunos sujetos una cota cercana al 75 %. En el caso de la clasificaci´on de los tres estados simult´aneamente conseguimos una tasa de aciertos en torno al 55 % de promedio, obteniendo en algunos casos valores por encima del 60 %. Como conclusi´on, hay que destacar que el procesado de la se˜nal para seleccionar las caracter´ısticas que mejoren la separabilidad de las clases, da unos resultados aceptables, puesto que ´este es un campo de investigaci´on reciente y en el cu´al se est´a investigando mucho en la actualidad. La mejora de las t´ecnicas de selecci´on de caracter´ısticas permitir´a en el futuro, desarrollar unas terapias de rehabilitaci´on que se ajusten de manera personal a cada paciente, adem´as lograr´a que la integraci´on de pr´otesis rob´oticas a la persona sea m´as ´ıntima y eficaz de lo que es en la actualidad. 6 ´ Indice general 1. Introducci´on 13 1.1. Contexto..................................... 13 1.2. EstadodelArte................................. 14 1.3. Objetivos .................................... 14 1.4. Organizaci´on de la memoria . . . . . . . . . . . . . . . . . . . . . . . . . . 16 1.5. Herramienta................................... 17 2. Protocolo de Experimentaci´on 19 2.1. Introducci´on................................... 19 2.2. Descripci´on del experimento . . . . . . . . . . . . . . . . . . . . . . . . . . 19 2.3. Obtenci´on de los datos de EEG . . . . . . . . . . . . . . . . . . . . . . . . 21 3. Procesado de la se˜nal 25 3.1. Introducci´on................................... 25 3.2. Filtros espaciales lineales . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26 3.2.1. Descripci´on ............................... 26 3.2.2. CAR................................... 27 3.2.3. Filtros Bipolares . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 3.2.4. Filtros Laplacianos . . . . . . . . . . . . . . . . . . . . . . . . . . . 27 3.2.5. Selecci´on del filtro . . . . . . . . . . . . . . . . . . . . . . . . . . . 28 3.3. Selecci´on de Filtros y Canales . . . . . . . . . . . . . . . . . . . . . . . . . 31 3.4. Rectificado y Suavizado . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33 3.5. Se˜nalseleccionada................................ 34 4. Clasificaci´on 37 4.1. Introducci´on................................... 37 4.2. Metodolog´ıa................................... 37 4.2.1. Validaci´on cruzada de los datos . . . . . . . . . . . . . . . . . . . . 37 4.2.2. An´alisis de las Componentes Principales (PCA) . . . . . . . . . . . 38 7 4.2.3. Normalizaci´on.............................. 39 4.2.4. Clasificadores .............................. 39 5. Resultados 41 5.1. Parejas de clases independientes . . . . . . . . . . . . . . . . . . . . . . . . 41 5.1.1. Shrinkage ................................ 41 5.1.2. LDA................................... 42 5.1.3. Interpretaci´on de los resultados . . . . . . . . . . . . . . . . . . . . 42 5.2. Resultados de las tres clases simult´anemente . . . . . . . . . . . . . . . . . 43 5.2.1. Shrinkage ................................ 43 5.2.2. LDA................................... 43 5.2.3. Interpretaci´on de los resultados . . . . . . . . . . . . . . . . . . . . 44 6. Fases del proyecto 45 6.1. Descripci´on de las fases . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45 6.2. DiagramadeGantt............................... 47 7. Cierre del proyecto 49 7.1. Conclusiones................................... 49 7.2. Trabajofuturo ................................. 50 8 ´ Indice de figuras 1.1. Desarrollo del proyecto . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.1. Descripci´on del experimento . . . . . . . . . . . . . . . . . . . . . . . . . . 20 2.2. Trayectorias ejecutadas por los usuarios . . . . . . . . . . . . . . . . . . . . 20 2.3. Disposici´on de los electrodos en la cabeza . . . . . . . . . . . . . . . . . . . 21 2.4. Segmentaci´on de los canales . . . . . . . . . . . . . . . . . . . . . . . . . . 23 3.1. estimador r2para el Sujeto 1 con los datos raw . . . . . . . . . . . . . . . 29 3.2. estimador r2para el Sujeto 1 con filtro bipolar FC1 para Reposo - Preparaci´on y Reposo - Movimiento y filtro bipolar P3 para Preparaci´on - Movimiento ................................... 30 3.3. Estimador r2para todos los canales del sujeto 4 . . . . . . . . . . . . . . . 30 3.4. Tipos de se˜nales en el dominio del tiempo para el sujeto 1, de arriba abajo: filtrado espacial, frecuencial, rectificado, suavizado . . . . . . . . . . . . . . 34 6.1. Diagrama de Gantt con las fases del proyecto . . . . . . . . . . . . . . . . 47 9 1.4. ORGANIZACI ´ ON DE LA MEMORIA 1. Introducci´on Figura 1.1: Desarrollo del proyecto 1.4. Organizaci´on de la memoria La presente memoria est´a estructurada en una serie de cap´ıtulos de manera que abarquen todo el proceso seguido durante la realizaci´on del proyecto. Debido a la imposibilidad de plasmar todos los aspectos t´ecnicos en los distintos cap´ıtulos, se har´an diferentes referencias a la bibliograf´ıa utilizada, as´ı como una serie de anexos que ampl´ıen la informaci´on de los cap´ıtulos. El primer cap´ıtulo expone el contexto en el que se realiz´o el proyecto, adem´as de una descripci´on del estado del arte y la herramienta utilizada para codificar los programas. El segundo cap´ıtulo describe el protocolo del experimento para obtener los datos del EEG. El tercero describe el procesamiento de la se˜nal para la mejora de sus caracter´ısticas: uso de filtros espaciales, frecuenciales, rectificado y suavizado de la se˜nal, y cu´al fue la selecci´on de caracter´ısticas para cada sujeto. El cuarto cap´ıtulo trata sobre el proceso de clasificaci´on: reducci´on de la dimensionalidad de la se˜nal, normalizaci´on y entrenamiento de los clasificadores. El quinto expone los resultados que se obtuvieron, tanto clasificando en parejas de estados, como utilizando la se˜nal con las tres clases a la vez. El sexto enumera las fases seguidas para la realizaci´on del trabajo, mostr´andolo adem´as de manera gr´afica con el diagrama de Gantt utilizado. Por ´ultimo, el s´eptimo cap´ıtulo muestra las conclusiones a las que se han llegado, as´ı como las l´ıneas de trabajo futuro en el campo de la investigaci´on sobre el reconocimiento de los tres estados. 16 1. Introducci´on 1.5. HERRAMIENTA Debido a la densidad del trabajo fue necesario incluir una serie de anexos: el primer anexo describir´a el estimador r2ocoeficiente de determinaci´on m´as al detalle. En el segundo anexo analizaremos en profundidad los aspectos te´oricos de los clasificadores lineales utlizados. En el siguiente anexo, expondremos todos los resultados obtenidos en el caso de la clasificaci´on simult´anea de los tres estados. Mostraremos en el cuarto anexo cu´ales han sido los programas desarrollados y, por ´ultimo, en el quinto adjuntamos el art´ıculo cient´ıfico que se envi´o y se aprob´o para la trig´esimo tercera conferencia anual de ingenier´ıa en medicina y biolog´ıa (IEEE EMBC) con toda la investigaci´on realizada en este proyecto, as´ı como los resultados y conclusiones obtenidas. 1.5. Herramienta Seleccionamos el lenguaje de programaci´on Matlab para codificar todo el procesado de la se˜nal y los clasificadores, porque es el lenguaje que se est´a utilizando en la empresa BitBrain Technologies para desarrollar el proyecto CORBYS. La versi´on elegida fue, en un primer momento, la revisi´on r2010a, cuando apareci´o la nueva actualizaci´on se sustituy´o por ´esta. 17 1.5. HERRAMIENTA 1. Introducci´on 18 Cap´ıtulo 2 Protocolo de Experimentaci´on 2.1. Introducci´on Este cap´ıtulo describe el experimento realizado por el Grupo de Rob´otica de la Universidad de Zaragoza, para la captura de la se˜nal de EEG producida por el movimiento del brazo. Con los datos de la se˜nal, se pretende detectar los tres estados asociados al movimiento: el sujeto se encuentra en reposo, el sujeto se prepara para desplazar el brazo y el sujeto inicia la trayectoria. 2.2. Descripci´on del experimento El objetivo del experimento es grabar el EEG, en el proceso de mover el brazo derecho, y, con estos datos, decodificar los estados de reposo, preparaci´on y ejecuci´on asociados al movimiento. Para ello, seis sujetos fueron sentados en una c´omoda silla frente del instrumental de investigaci´on, ver Figura 2.1, y se les pidi´o que movieran el brazo derecho desde la posici´on de reposo hasta una de las ocho localizaciones posibles. Una vez que hubieran tocado la posici´on deseada deber´ıan volver a mover el brazo hasta la localizaci´on origen. Este proceso lo iniciaba la persona voluntariamente cada 7’5 segundos aproximadamente de media(habi´endose registrado un m´ınimo de 2’8 y un m´aximo de 9,7 segundos). Mientras que el sujeto estaba realizando el desplazamiento al lugar elegido, se le instaba a que no parpadee o que lo haga lo menos posible, sugiri´endole que mirase fijamente a un punto. Entre repetici´on y repetici´on, pod´ıan parpadean todo lo que quisieran. Todo el proceso se repiti´o a lo largo de cinco bloques de cinco minutos cada uno, durante los cuales cada uno de los usuarios hizo un n´umero diferente de repeticiones o ”trials”, como se puede ver en la Tabla 2.1. Durante el tiempo que duraba cada repeticion, se iba registrando una se˜nal de sincronizaci´on que nos cu´ando es el instante en el que se comenz´o a ejecutar el movimiento. Adem´as de la obtenci´on de la se˜nal de EEG, se grababan las diversas trayectorias 19 2.2. DESCRIPCI ´ ON DEL EXPERIMENTO 2. Protocolo de Experimentaci´on que realizaba el sujeto para ser utilizados en la decodificaci´on de las trayectorias a partir del EEG para otros proyectos, ver Figura 2.2 Figura 2.1: Descripci´on del experimento SUJETO REPETICIONES 1 196 2 98 3 157 4 75 5 135 6 130 Tabla 2.1: N´umero de repeticiones o ”trials” por cada sujeto Figura 2.2: Trayectorias ejecutadas por los usuarios 20 2. Protocolo de Experimentaci´on 2.3. OBTENCI ´ ON DE LOS DATOS DE EEG 2.3. Obtenci´on de los datos de EEG La actividad de EEG fue grabada con un sistema gTec (2 gUSBamp amplificadores con 28 electrodos situados de acuerdo al sistema 10/10 internacional como se puede ver en la Figura 1, con la toma a tierra en FPz y tomando como referencia el l´obulo de la oreja. Para registrar el movimiento de los musculos que rodean a los ojos se utiliz´o un Electrooculograma o EOG, ya que dicho movimiento interfiere en la se˜nal de EEG. Los movimientos tanto verticales como horizontales fueron grabados para su posterior eliminaci´on de la se˜nal del electroencefalograma. Los dos tipos de se˜nales, tanto el EEG como el EOG, son digitalizadas con una frecuencia de sampleo de 256 Hz, a las cuales se les aplica un filtro de primera clase y un filtro de paso bajo a 60 Hz. Los 28 electrodos o canales reciben los siguientes nombres: Fp1, Fp2, F7, F8, F3, Fz, F4, T7, T8, C3, Cz, C4, P7, P8, P3, Pz, P4, O1, O2, FC1, FC2, FC5, FC6, CPz, CP1, CP2, CP3, CP4 y su disposici´on, a lo largo del cr´aneo de los sujetos, era: Figura 2.3: Disposici´on de los electrodos en la cabeza 21 2.3. OBTENCI ´ ON DE LOS DATOS DE EEG 2. Protocolo de Experimentaci´on Adem´as de las medidas de EEG, con la se˜nal de sincronizaci´on que indica en qu´e momento se inicia el movimiento, el sistema tambi´en grab´o marcadores de tiempo correspondientes al inicio del movimiento del brazo, quedando la se˜nal en todos los casos organizada de la siguiente manera: tenemos un marcador negativo hasta que llegamos al marcador 0 que nos indica que es entonces cuando el sujeto ha comenzado a mover el brazo. Esta informaci´on se obtuvo teniendo tanto en el brazo, como en el punto destino, botones que registraban esta actividad. A trav´es de los marcadores del tiempo, la se˜nal se dividi´o en tres segmentos de 300 ms cada uno, uno por cada clase, del siguiente modo: los 300 ms iniciales de la se˜nal para el reposo, los 300 ms desde el instante marcado como 0 corresponder´an al estado de ejecuci´on del movimiento y los 300 ms previos al instante 0 ser´an considerados como anticipaci´on, como se puede ver en la Figura 2.4. Como resultado tendremos la se˜nal en crudo, o ”raw data”, la cu´al est´a organizada de la siguiente manera: X∈RCXT , siendo Tel n´umero de muestras temporales, una vez que nos hemos quedado con los 300 ms por clase (en nuestro caso en particular, como estaba realizado el experimento y la frecuencia de sampleo de las muestras, tenemos que son 77 muestras por clase), y Cel n´umero de canales. 22 2. Protocolo de Experimentaci´on 2.3. OBTENCI ´ ON DE LOS DATOS DE EEG Figura 2.4: Segmentaci´on de los canales 23 2.3. OBTENCI ´ ON DE LOS DATOS DE EEG 2. Protocolo de Experimentaci´on 24 Cap´ıtulo 3 Procesado de la se˜nal 3.1. Introducci´on Como la se˜nal de EEG presenta ruido, as´ı como otros procesos que no son de inter´es para este proyecto, es imprescindible estudiar las caracter´ısticas o ”features” para poder clasificar adecuadamente los estados de reposo, preparaci´on y movimiento. Con el estudio de las caracter´ısticas, podemos diferenciar cu´ales son las m´as significativas y seleccionarlas para la clasificaci´on. Con las m´as importantes seleccionadas, ser´a el momento de comenzar un proceso de realce de dichas propiedades, lo cu´al redundar´a en una mejor clasificaci´on de las clases. De acuerdo con estudios anteriores, la informaci´on m´as relevante es la conocida como ”desincronizaciones o sincronizaciones relacionadas con el evento” o ”event-related potentials” (ERD/ERS) de las ´areas motoras [2], las cuales se encuentran en las frecuencias comprendidas entre los 12 y los 15 Hz, as´ı como una importante actividad en las bajas frecuencias correspondientes a los ”potenciales corticales lentos” o ”slow cortical potentials”, los cuales se relacionan con la preparaci´on motora [3]. A priori, esta informaci´on deber´ıa ser suficiente por s´ı sola para poder distinguir claramente los tres estados, puesto que hasta la fecha los ERD/ERS han sido utilizados para clasificar reposo y movimiento y los potenciales corticales lentos lo fueron para separar preparaci´on o anticipaci´on de la ejecuci´on del movimiento. Los resultados obtenidos con estas ”features” no dieron buenos resultados, por lo que tuvimos que recurrir a un refinado de la se˜nal, que nos mejorara la separabilidad de las clases, para lo que se aplicaron distintos filtros: filtros espaciales, como pueden ser laplacianos, bipolares o CAR, y filtros frecuenciales, en el rango de frecuencias en los cuales se encuentre la mayor diferencia estad´ıstica entre clases utilizando el estimador r2, conoci- 25 3.3. SELECCI ´ ON DE FILTROS Y CANALES 3. Procesado de la se˜nal tad´ıstica en las frecuencias bajas. Este ´ultimo responde, probablemente, a los potenciales bajos relacionados con el movimiento como se describe en [14]. Como vemos en la Figura 3.3 para el peor sujeto, las frecuencias significativas son totalmente distintas: entre reposo y preparaci´on surge en la parte media de la banda β, mientras que entre reposo y movimiento aparece en la banda α, igual que entre preparaci´on y movimiento. En la Tabla 3.2 se puede ver qu´e filtros y qu´e rango de frecuencias frhan sido seleccionados para cada uno de los seis sujetos. SUJETO REPPREP REP - MOV PREPMOV 1 paso banda 20 - 30 hz paso banda 20 - 30 hz paso bajo 10 hz 2 paso banda 25 - 30 hz paso bajo 10 - 20 hz paso banda 10 - 20 hz 3 paso banda 20 - 30 hz paso bajo 5 hz paso banda 20 - 30 hz 4 paso banda 20 - 25 hz paso banda 10 - 15 hz paso banda 10 - 15 hz 5 paso banda 10 - 15 hz paso banda 15 - 25 hz paso banda 17 - 20 hz 6 paso banda 5 - 10 hz paso bajo 10 hz paso bajo 10 hz Tabla 3.2: Tipo de filtro butterworth para cada sujeto El proceso de selecci´on de canales es an´alogo al de selecci´on del filtro espacial y de las frecuencias: a partir del coeficiente de determinaci´on, podemos ver en qu´e canales, adem´as del C3, son los que presentan un valor num´erico m´as alto y por lo tanto, van a ser los m´as significativos a la hora de poder separar las clases de la mejor manera posible. En las Figuras 3.2 y 3.3 podemos ver de una manera visual cu´ales son los canales seleccionados para los sujetos 1 y 4. En la Tabla 3.3 podemos ver todos los canales seleccionados para cada sujeto, en todas ellas debe aparecer C3 por su importancia para la investigaci´on. 32 3. Procesado de la se˜nal 3.4. RECTIFICADO Y SUAVIZADO SUJETO REPPREP REP - MOV PREPMOV 1 T8, C3, C4, P7, P8, P3, Pz, P4, 01, 02, FC2, FC6, CPZ, CP1, CP2, CP3, CP4 T8, C3, C4, P7, P8, P3, Pz, P4, 01, 02, FC2, FC6, CPZ, CP1, CP2, CP3, CP4 C3, Cz, FC1, FC5, CP1, CP3 2 FP1, F7, F3, T7, C3, P7, FC5 F3, C3, P7, CP1, CP3 F3, C3, P7, CP1, CP3 3 F3, C3, P7, FC1, CP1, CP3 F3, C3, P7, FC1, CP1, CP3 F3, C3, P7, FC1, CP1, CP3 4 F7, F3, T7, T8, C3, C4, P7, O1, FC5, CP3 F7, F3, T7, T8, C3, P7, FC1, FC5, FC6 F7, F3, T7, T8, C3, P7, P3, FC1, FC5 5 FP1, FP2, F7, F8, F3, T7, T8, C3, P7, FC5, FC6, CP3 T8, C3, C4, P7, P8, P3, CP1, CP2, CP3, CP4 F3, Fz, C3, P3, FC1 6 F7, F3, C3, Cz, FC1, FC5, CP1, CP3 C3, C4, P3, P4, FC2, FC6, CP2, CP3 C3, Cz, C4, P3, Pz, P4, O1, O2, FC2, CPz, CP1, CP2, CP3 Tabla 3.3: Canales seleccionados para cada sujeto 3.4. Rectificado y Suavizado El ´ultimo paso antes de entrenar el clasificador consisti´o en procesar los canales seleccionados en el dominio del tiempo nuevamente. Primero, aplicamos un filtro ”Butterworth”, ver Ecuaci´on 3.8 en el espectro de la frecuencia seleccionado fr. Diversos estudios demuestran que una buena manera de mejorar la estabilidad de la se˜nal e incrementar la diferencia estad´ıstica entre clases consiste en rectificar y suavizar la se˜nal [15]. El rectificado consiste en hacer positiva toda la se˜nal, es decir, el valor absoluto de toda la se˜nal. Por otro lado, la t´ecnica de suavizado utiliz´o una media m´ovil simple (SMA) utilizando los n = 3 datos anteriores. 33 3.5. SE ˜ NAL SELECCIONADA 3. Procesado de la se˜nal 3.5. Se˜nal seleccionada Por ´ultimo, igual que en el caso de la elecci´on los distintos tipos de caracter´ısticas, se utiliz´o un test de r2para seleccionar cu´al de los cuatro tipos de se˜nales es la mejor para la clasificaci´on. Los cuatro tipos de se˜nales son: se˜nal filtrada espacialmente, se˜nal filtrada espacialmente y frecuencialmente, se˜nal filtrada espacialmente, frecuencialmente y rectificada y, por ´ultimo, la se˜nal con todo el procesado indicado: filtrada espacialmente, frecuencialmente, rectificada y suavizada. En la Figura 3.4 se puede ver las cuatro se˜nales para el sujeto 1 y cu´ales se seleccionaron a tenor de lo mostrado por el coeficiente de determinaci´on. Para ver de manera gr´afica cu´ales fueron los canales seleccionados, las bandas de frecuencias y el tipo de se˜nal seleccionada para cada uno de los sujetos consultar el Anexo B. Figura 3.4: Tipos de se˜nales en el dominio del tiempo para el sujeto 1, de arriba abajo: filtrado espacial, frecuencial, rectificado, suavizado En la Tabla 3.4 se puede ver cu´al de los cuatro tipos de se˜nales se seleccionaron para la clasificaci´on. 34 3. Procesado de la se˜nal 3.5. SE ˜ NAL SELECCIONADA SUJETO REPPREP REP - MOV PREPMOV 1 Suavizada Suavizada Filtrado frecuencial 2 Filtrado espacial Filtrado espacial Filtrado espacial 3 Filtrado espacial1Filtrado frecuencial Filtrado espacial1 4 Filtrado espacial Suavizado Filtrado espacial 5 Suavizado Filtrado frecuencial Suavizado 6 Filtrado espacial Filtrado espacial Filtrado espacial Tabla 3.4: Se˜nal seleccionada para cada sujeto 1Emp´ıricamente se demostr´o que para clasificar las tres clases por separado es mejor usar la se˜nal suavizada 35 3.5. SE ˜ NAL SELECCIONADA 3. Procesado de la se˜nal 36 Cap´ıtulo 4 Clasificaci´on 4.1. Introducci´on Con las caracter´ısticas finales de la se˜nal seleccionadas, el ´ultimo paso del proceso consisti´o en el entrenamiento de dos clasificadores lineales multiclase: LDA y ”shrinkage” [16], [17], [18], [19], [20], [21]. Como el n´umero de caracter´ısticas de la se˜nal es parecido, en el mejor de los casos, o mucho mayor, en el peor, se utiliz´o ”shrinkage” como estimador de la matriz de covarianza de las clases [16], [17]. Como se disponen de pocas muestras de la se˜nal, y m´as si los comparamos con el n´umero de caracter´ısticas, se necesita una validaci´on de estos datos. Para ello, se utilizar´a la t´ecnica de validaci´on cruzada de los datos, que nos dividir´a en dos grupos de datos: ”entrenamiento” y”test” un n´umero determinado de veces. Con los datos separados en entrenamiento y validaci´on, les aplicaremos un an´alisis de las componentes principales o PCA para reducir la dimensionalidad de los datos, despu´es normalizaremos los datos y, por ´ultimo, entrenaremos los clasificadores con los datos de ”entrenamiento” y validaremos con los datos de ”test” 4.2. Metodolog´ıa 4.2.1. Validaci´on cruzada de los datos Necesitamos dividir los datos en datos para entrenar los clasificadores y datos para probarlos, pero disponemos de un conjunto peque˜no de muestras de los datos, se hace necesario un proceso de validaci´on de los datos, en nuestro caso usamos el m´etodo de validaci´on cruzada de datos. La validaci´on cruzada o ”cross-validation” es una pr´actica estad´ıstica que consiste en un grupo de datos, se subdividen en dos subconjuntos: uno 37 4.2. METODOLOG´ IA 4. Clasificaci´on de entrenamiento y otro de validaci´on. ´ Este proceso se repite un n´umero ”n” de veces obteniendo as´ı la validaci´on del modelo y de los datos. Para este proyecto se utiliza la variaci´on conocida como ”ten-fold cross validation” que consiste en separar los datos en dos conjuntos: el primer conjunto de entrenamiento con el 90 % de los datos y el conjunto de validaci´on con el 10 % restante. ´ Este proceso se repite un total de 10 veces, con lo que todos los datos han sido usados como entrenamiento y como validaci´on. 4.2.2. An´alisis de las Componentes Principales (PCA) El an´alisis de las componentes principales, en ingl´es ”principal component analysis” (PCA), es una t´ecnica estad´ıstica que sirve para reducir la dimensionalidad de los datos. Esta t´ecnica es utilizada para estudiar la variabilidad de las componentes de los datos y ordenarlas de mayor a menor importancia. PCA elabora una serie de transformaciones lineales que seleccionan unos nuevos sistemas de coordenadas en los cuales las componentes ordenadas de mayor a menor variabilidad ocupan los nuevos ejes de manera ordenada. Tenemos una matriz Xcon mfilas y ncolumnas: XT∈Rnxm (4.1) Calculamos su descomposici´on en valores singulares: X=W∗Σ∗VT(4.2) Con Wlos autovectores de X∗XT, Σ ∈Rmxn una matriz diagonal regular y Vlos autovectores de XT∗Xy siendo WT L, los Lprimeros autovectores de W, la matriz con la dimensionalidad reducida ser´ıa: Y=WT L∗X(4.3) En el caso del proyecto, con unas matrices de datos muy grandes fue necesario reducir la dimensionalidad de los datos para que los clasificadores dieran resultados v´alidos. Entonces, con el 90 % de los datos para validaci´on o ”entrenamiento” aplicamos PCA para quedarnos con el 99 % de la variabilidad de la se˜nal. Una vez que ten´ıamos la matriz de transformaci´on para quedarnos con el 99 % de los datos la utiliz´abamos para reducir las 38 4. Clasificaci´on 4.2. METODOLOG´ IA dimensiones de los datos de ”entrenamiento” y de ”test” por separado, para que no ubiera interferencia entre los datos. 4.2.3. Normalizaci´on Con los datos de validaci´on obten´ıamos el m´aximo y el m´ınimo, siendo estos valores 1 y 0. Luego aplic´abamos esta normalizaci´on con esos valores a los datos de validaci´on, para que ´estos estuvieran en la misma escala, pudiendo haber valores de m´as de 1 o menores que 0 puesto que puede que el m´aximo global de la se˜nal se encontrara en los datos de validaci´on y no en los de entrenamiento. Si tenemos que xmax =max(entrenamiento) , xmin =min(entrenamiento), ymax = 1, ymin = 0, entonces tenemos que: valor normalizado =(ymax −ymin)∗(xi−xmin) (xmax −xmin) + ymin ,∀xi∈Xentrenamiento ∨Xtest (4.4) 4.2.4. Clasificadores El paso final de la consisti´o en elegir qu´e m´etodos de clasificaci´on se iban a utilizar y c´omo se iban a clasificar los datos. Los m´etodos que seleccionamos fueron dos clasificadores lineales, ampliamente utilizados en investigaci´on sobre BCI: el an´alisis lineal discriminante o LDA y el m´etodo de estimaci´on de las matrices de covarianza derivado del LDA, conocido como ”shinkage”. LDA (ver Anexo A.2.2) asume que los datos est´an distribuidos siguiendo una normal con media µky covarianza Σk,∀k∈1..N, siendo N el n´umero de clases distintas. En nuestro caso k puede tomar valores 1 (reposo), 2 (preparaci´on) y 3 (anticipaci´on). La particularidad de este m´etodo es que asume que todas las clases tienen la misma covarianza: Σk= Σ. La matriz de covarianza estimada es la matriz de covarianza emp´ırica. Cada muestra pertenecer´a a la clase kque minimice la funci´on: δk(X) = (X−ˆµk)0∗ˆ Σ−1∗(X−ˆµk) + ln|ˆ Σ| − 2∗lnπk(4.5) El estimador est´andar de la matriz de covarianza es la covarianza emp´ırica. Con datos que tienen muchas dimensiones pero pocas muestras en comparaci´on, como es nuestro caso, este estimador se vuelve impreciso porque el n´umero de par´ametros desconocidos que hay que aproximar es cuadr´atico con respecto al n´umero de dimensiones, por ello se 39 4.2. METODOLOG´ IA 4. Clasificaci´on hizo necesario usar otro m´etodo: ”shrinkage”. El m´etodo ”shrinkage”, ver Anexo A.2.3, parte del mismo principio, pero surge para contrarrestar el ”bias” o sesgo sistem´atico que surge en la estimaci´on de la matriz de covarianza. Para solucionar esto, la matriz de covarianza Σ, es reemplazada por la estimaci´on, e Σ: e Σ(γ) = (1 −γ)∗Σ + γ∗ν∗I(4.6) Siendo γ∈[0,1] el par´ametro de afinaci´on, νel promedio de la traza de los autovalores dividido entre la dimensi´on de los datos e ”I” la matriz identidad. Con γ= 1 asumimos matrices de covarianza esf´ericas y con γ= 0 estamos en el caso de LDA. Los dos tipos de clasificaciones hechas fueron: clasificando entre dos clases (reposo - preparaci´on, reposo - movimiento y preparaci´on - movimiento), utilizando el vector de caracter´ısticas formado por el par de se˜nales en cuesti´on, y clasificando las tres clases a la vez simult´aneamente usando la concatenaci´on de las caracter´ısticas seleccionadas para cada par de clases. 40 Cap´ıtulo 5 Resultados 5.1. Parejas de clases independientes 5.1.1. Shrinkage SUJETO REPPREP REP - MOV PREPMOV 1 69’59 67’86 80’78 2 67’30 81’01 75’01 3 72’01 72’57 68’47 4 58’63 55’39 63’18 5 57’83 58’85 54’09 6 66’92 60’85 65’38 PROMEDIO 65’38 66’08 67’81 Tabla 5.1: Porcentaje de aciertos por parejas usando ”shrinkage” En la T6 mostramos que el promedio de aciertos de todos los sujetos se sit´ua en torno al 66 %, teniendo picos en momentos puntuales que llegan a alcanzar el 80 %. Podemos comprobar que los individuos que mejor han respondido han alcanzado un porcentaje de clasificaci´on por encima del 71 %, teniendo un caso intermedio en el sujeto 6 que llega al 65 %. En el otro extremo, los sujetos que peor funcionan no alcanzan el 60 % de los datos bien clasificados. 41 6.2. DIAGRAMA DE GANTT 6. Fases del proyecto 48 Cap´ıtulo 7 Cierre del proyecto 7.1. Conclusiones Este proyecto trata sobre la investigaci´on en la clasificaci´on de los estados de reposo, preparaci´on y ejecuci´on del movimiento del brazo derecho a trav´es del electroencefalograma. Siguiendo un procesado de la se˜nal, ampliamente conocido en el campo de la neuro-ciencia y el BCI , se ha aplicado a seis sujetos de manera dispar: tres sujetos responden muy bien al procesado, que son conocidos como ”responders”, y otros tres que no responden al mismo sistema. Esto se puede apreciar ya en el test de r2, pues hay muchas dificultades para encontrar canales con una diferencia estad´ıstica significativa en el canal C3. El hecho de que haya ”non-responders” no responde a una causa clara y com´un, a veces es debido a fallos en la interpretaci´on del experimento por parte del sujeto, la aparici´on de muchos artefactos por la tensi´on inherente al propio experimento o porque su estructura mental no se ajusta a lo que se considera ”normalidad”, es decir, que el canal C3 no es el m´as relevante en el movimiento del brazo derecho, es decir, que el punto de mayor diferencia no se encuentra justo donde lo tiene la mayor´ıa de las personas. Por otro lado, cuando clasificamos los estados por parejas entrenando tres clasificadores para cada pareja, obtenemos en los sujetos que s´ı que responden al tratamiento de la se˜nal unos porcentajes de acierto superior al 72 %, llegando a tener picos del 80 %, lo que nos quiere decir que el procesado es correcto en los casos de los sujetos ”responders” y que con investigaciones m´as exhaustivas y el estudio de las se˜nales y filtros m´as complejos se conseguir´an en el futuro unos porcentajes a´un mejores, lo cu´al sentar´a las bases de una integraci´on definitiva de las pr´otesis rob´oticas, de manera m´as ´ıntima, en las personas con miembros amputados, as´ı como unas terapias de rehabilitaci´on personalizadas. 49 7.2. TRABAJO FUTURO 7. Cierre del proyecto Por otro lado, el caso de clasificar las tres clases a la vez, este es un proyecto pionero en el campo de investigaci´on, pues no hab´ıa nada de importancia publicado hasta la fecha. Como vemos en los resultados el porcentaje de aciertos para los ”responders” sufre un baj´on de un 10 % pero a´un as´ı est´a por encima del 60 % llegando incluso a obtener un 70 % en momentos puntuales. Este proyecto viene a demostrar que seleccionando las caracter´ısticas de la se˜nal de EEG, es razonable pensar que se pueden llegar a alcanzar cotas de un 70 % de ejemplos bien clasificados, lo que hoy por hoy es un buen resultado en este campo de investigaci´on, porque es una ciencia relativamente nueva y en la que continuamente se est´an descubriendo nuevos datos. El seguir estudiando en este campo para la mejora de las caracter´ısticas servir´a para en un futuro ser capaces de personalizar la pr´otesis, ya que ´esta adaptar´a los par´ametros de la misma al amputado, as´ı como tener unos mejores porcentajes en la clasificaci´on har´a que las propias pr´otesis respondan mejor a las intenciones de sus usuarios, incrementando as´ı sus prestaciones. Tanto la investigaci´on llevada a cabo en la selecci´on de caracter´ısticas, como los resultados y conclusiones obtenidos, fueron plasmados en un art´ıculo que fue enviado y aprobado para su presentaci´on en la conferencia IEEE EMBC (Engineering in Medicine and Biology Conference) en Agosto de 2011. 7.2. Trabajo futuro El trabajo a corto y medio plazo sobre este tema, tanto por el Grupo de Rob´otica de la Universidad de Zaragoza, como BitBrain Technologies, tiene tres vertientes fundamentales: la selecci´on autom´atica de caracter´ısticas, eliminaci´on de las marcas de la se˜nal y la clasificaci´on de la se˜nal a lo largo del tiempo. El primer aspecto de futuro consistir´a en reducir la intervenci´on del terapeuta a la hora de seleccionar las caracter´ısticas de manera visual o num´erica, como realizamos en la actualidad, para ello se buscar´a la clasificaci´on y selecci´on de las caracter´ısticas m´as importantes de la manera m´as autom´aticamente posible, as´ı que de la propia obtenci´on de la se˜nal ya habr´an descartado las propiedades menos relevantes de la misma, que minan la capacidad de clasificaci´on, al ser comunes para los tres estados de la se˜nal y s´olo aportan redundancia en los datos. El segundo, la eliminaci´on de las marcas de la se˜nal, consistir´a en eliminar la marca 50 7. Cierre del proyecto 7.2. TRABAJO FUTURO que presenta la se˜nal en la que indicamos que empieza el movimiento y que est´a marcado como instante cero, a partir del cu´al las marcas temporales son positivas mientras que antes lo han sido negativas. En esto hay que investigar m´as para poder ser capaces, a partir del estudio de la propia forma fisiol´ogica de la se˜nal, ´esta nos indique cu´al es el momento exacto en el que se empieza a ejecutar el movimiento, y no recurrir a sensores colocados en los brazos. El ´ultimo punto de investigaci´on futura trata sobre la clasificaci´on de la se˜nal a lo largo del tiempo. En este proyecto, la se˜nal se separaba en los tres estados cogiendo 300 ms de tiempo para cada uno de ellos, pero en un futuro, se investigar´a en que a partir de la se˜nal continua ser capaces de decir en cu´al de los tres estados se encuentra la se˜nal, para lo cu´al se empez´o a trabajar en el campo de los ”Modelos Ocultos de Markov” para este mismo proyecto, pero por falta de tiempo no se consiguieron resultados concluyentes y por eso no se ha reflejado en esta memoria. A t´ıtulo personal, considero que en los experimentos, el hecho de que al mover el brazo se tuviera que volver luego a la posici´on de origen puede provocar que el estado de ejecuci´on de movimiento de ida se solape con la intenci´on del movimiento de vuelta del brazo haciendo que la clasificaci´on del estado de ejecuci´on del movimiento no sea tan precisa como debiera ser. En este tema, ser´ıa conveniente realizar nuevos experimentos que obligaran a mantener la mano en el destino una vez que se haya efectuado el movimiento. 51 7.2. TRABAJO FUTURO 7. Cierre del proyecto 52 Bibliograf´ıa [1] K. P. Tee, C. Guan, K. Keng Ang, K. Soon Phua, C. Wang, H. Zhang. ”Augmenting cognitive processes in robot-assisted motor rehabilitation”,Proceedings of the 2nd Biennial IEEE/RAS-EMBS International Conference on Biomedical Robotics and Biomechatronics, Scottsdale, USA, 2008. [2] G. Pfurtscheller, F. Lopez da Silva. ”Event-related synchronization”,Handbook of Electroencephalography and Clinical Neurophysiology, 1999 [3] M. Kutas, E. Donchin ”Preparation to respond as manifested by movement-related brain potentials”,Brain Research, 1980, n´umero 95 [4] C. Nauper, G. Pfurtscheller. ”Event-related dynamics of cortical rhythms: frequencyspecific features and functional correlates”,International Journal of Psychophysiology, 2001, p´ags. 41-58 [5] Valerie Morash, Ou Bai, Stephen Furlani, Peter Lin, Mark Hallet. ”Classifying eeg signals preceding right hand, left hand, tongue and right foot movements and motor imageries”,Clinical Neurophysiology, vol. 119, no. 11, p´ags. 2570 - 2578, 2008 [6] B. O. Peters, G. Pfurtscheller, G. Edlinger. ”Automatic differentiation of multichannel eeg signals”,IEEE Transations on Biomedical Engineering, vol. 48, no. 1, 2001 [7] B. Blankertz, G. Dornhege, C. Schafer, R. Krepki, J. Kohlmorgen, K. R. Muller, F. Losch, V. Kunzmann, G. Curio. ”Boosting bit rates and error detection for the classification of fast-paced motor commands based on single-trial eeg analysis”,IEEE Transactions Neural Systems Rehabilitation Engineering, vol. 11, no. 127, 2003 [8] M. Krauledat, G. Dornhege, B. Blankertz, F. Losch, G. Curio, K. R. Muller. ”Improving speed and accuray of brain-computer interfaces using readiness potential features”,Engineering in Medicine and Biology Society, p´ags. 4511 - 4515, 2004 [9] Everitt, B. S. ”Cambridge Dictionary of Statistics”,Cambridge University Press, 2a ed., 2002 53 BIBLIOGRAF´ IA BIBLIOGRAF´ IA [10] Benjamin Blankertz, Ryota Tomioka, Steven Lemm, Motoaki Kawanabe, Klaus- Robert M¨uller. ”Optimizing Spatial Filters for Robust EEG Single-Trial Analysis”, IEEE signal processing magazine, Vol. XX, 2008 [11] Paul L. Nunez, Ramesh Srinivasan, Andrew F. Westdorp, Ranjith S. Wijesinghe, Don M. Tucker, Richard B. Silberstein, Peter J. Cadusch ”EEG coherency I: statistics, reference electrode, volume conduction, Laplacians, cortical imaging, and interpretation at multiple scales”,Electroencephalography and clinical Neurophysiology, no. 103, p´ags. 499 - 515 [12] Herbert Ramoser, Johannes M¨uller-Gerking, and Gert Pfurtscheller. ”Optimal Spatial Filtering of Single Trial EEG During Imagined Hand Movement”,IEEE Transactions on rehabilitation engineering, vol. 8, no. 4, 2000 [13] G. Pfurtscheller, C. Neuper ”Event-related synchronization of mu rhythm in the EEG over the cortical hand area in man”,Neuroscience Letters, vol. 147, no. 1, p´ags. 93 - 96, 1994 [14] S. Waldert, H. Preissl, E. Demandt, C. Braun, N. Birbaumer, A. Aertsen, C. Mehring. ”Hand movement direction decoded from meg and eeg”,J. Neuroscience, vol. 28, no. 4, p´ags. 1000 - 1008, 2008 [15] Carmen Vidaurre. ”Clasificaci´on de componentes frecuenciales”, Machine Learning Laboratory, Berlin Institute of Technology. [16] Benjamin Blankertz, Steven Lemm, Matthias Treder, Stefan Haufe, Klaus-Robert M¨uller. ”Single-trial analysis and classification of ERP components - A tutorial”, NeuroImage. 2010 [17] J. B. Copas ”Regression , prediction and shrinkage”,Journal of the Royal Statistical Society. Series B (Methodological, vol. 45, no. 3, p´ags. 311 - 354, 1983 [18] Jerome H. Friedman. ”Regularized Discriminant Analysis”,Journal of the American Statistical Association, vol. 84, no. 405, p´ags. 165 - 175, 1989. [19] Jia Li. ”Linear Discriminant Analysis”, Department of Statistics, The Pennsylvania State University. [20] Jia Li. ”Regularized Discriminant Analysis and Reduced-Rank LDA”, Department of Statistics, The Pennsylvania State University. [21] Christopher M. Bishop. ”Machine Learning and Pattern Recognition”, Ed. Springer, 2aed., p´ags 101-135 54 BIBLIOGRAF´ IA BIBLIOGRAF´ IA [22] Julianne Sch¨afer, Korbinian Strimmer. ”A Shrinkage Approach to Large-Scale Covariance Matrix Estimation and Implications for Functional Genomics”,Statistical Applications in Genetics and Molecular Biology, vol. 4, art. 32, 2005. [23] Eduardo L´opez-Larraz, I˜naki Iturrate, Luis Montesano, Javier M´ınguez ”Real-Time Recognition of Feedback Error-Related Potentials during a Time-Estimation Task”, International Conference of the IEEE Engineering in Medicine and Biology Society, 2010. 55