Repositorio Institucional de Documentos

Proyecto Final de Carrera Ingenier´ıa de Telecomunicaciones Curso 2010-2011 Dise˜no de una Base Normativa de EEG Basal en Adultos M´onica Aguilar Herrero Marzo de 2011 Director: Javier M´ınguez Zafra Departamento de Inform´atica e Ingenier´ıa de Sistemas Centro Polit´ecnico Superior Universidad de Zaragoza Si el cerebro fuese tan sencillo para que pudi´eramos comprenderlo nosotros ser´ıamos tan simples que no lograr´ıamos entenderlo. Albert Einstein iv Dise˜no de una Base Normativa de EEG Basal en Adultos RESUMEN El electroencefalograma constituye una de las t´ecnicas de adquisici´on no invasivas m´as ampliamente utilizadas para el registro de actividad cerebral. El uso de Bases Normativas de EEG resulta de gran inter´es para determinar un patr´on de referencia en cuanto a se˜nal EEG se refiere. En el presente proyecto se dise˜na, desarrolla y valida una Base Normativa de EEG basal en adultos, con el objetivo de automatizar el diagn´ostico, tratamiento y evaluaci´on de las enfermedades neurol´ogicas. Actualmente la problem´atica existente en el desarrollo de una Base Normativa de EEG viene determinada por la alta variabilidad interpersonal de la se˜nal EEG. Este aspecto ha sido estudiado mediante la relaci´on de esta variabilidad con factores fisiol´ogicos llevando a cabo una estratificaci´on de sujetos clasificados por dichos factores. Esta t´ecnica ha resultado ´util en el caso de Bases Normativas de adultos, pero no en el caso de Bases Normativas infantiles. Este proyecto plantea una alternativa mediante la estratificaci´on de la Base Normativa de forma expl´ıcita sobre los cambios significativos de la propia se˜nal EEG. Adem´as, se han elaborado unas herramientas matem´aticas que permiten aplicar la Base Normativa en el diagn´ostico de una patolog´ıa neurol´ogica y en el an´alisis de la evoluci´on de un sujeto bajo una terapia establecida. La ventaja del conjunto de t´ecnicas desarrolladas frente a las existentes sobre Bases Normativas radica en la automatizaci´on tanto del proceso de inclusi´on de un sujeto a la base, como en la aplicaci´on de las herramientas que determinan la normalidad en la se˜nal EEG, evitando cualquier intervenci´on humana en el dise˜no de la Base Normativa. El conjunto de procesos desarrollados se centran en favorecer la capacidad de discriminaci´on entre sujetos sanos y patol´ogicos, en funci´on de unos patrones extra´ıdos de la se˜nal EEG, y analizar estad´ısticamente estos patrones para determinar su robustez frente a la variabilidad interpersonal. La b´usqueda de los patrones de referencia se orient´o hacia descriptores de la energ´ıa en las bandas del espectro de potencias de la se˜nal EEG los cuales se ha demostrado su relaci´on con diferentes patolog´ıas neurol´ogicas. Finalmente, se lleva a cabo la validaci´on global de las t´ecnicas dise˜nadas tomando la Base Normativa desarrollada como una prueba diagn´ostica bajo un contexto cl´ınico, mediante el uso de datos EEG cuasi-sint´eticos. v ´ Indice 1. Introducci´on 2 1.1. Contexto..................................... 2 1.2. Motivaci´on y trabajo relacionado . . . . . . . . . . . . . . . . . . . . . . . 2 1.3. Objetivo y Alcance del proyecto . . . . . . . . . . . . . . . . . . . . . . . . 4 1.4. Organizaci´on .................................. 6 2. Registro y procesado de se˜nal EEG 7 2.1. Introducci´on................................... 7 2.2. Electroencefalograma.............................. 7 2.3. Descripci´on de los datos y Registro . . . . . . . . . . . . . . . . . . . . . . 8 2.4. Procesadodese˜nal ............................... 9 2.4.1. Filtrado de artefactos . . . . . . . . . . . . . . . . . . . . . . . . . . 9 2.4.2. Extracci´on de patrones . . . . . . . . . . . . . . . . . . . . . . . . . 10 3. Fiabilidad 14 3.1. Introducci´on................................... 14 3.2. M´etodo1:Test-Retest ............................. 15 3.2.1. Resultados................................ 16 3.3. M´etodo2:Split&Half.............................. 17 3.3.1. Resultados................................ 18 3.4. Ampliaci´on de se˜nal EEG . . . . . . . . . . . . . . . . . . . . . . . . . . . 22 vii ´ INDICE ´ INDICE 3.4.1. Resultados................................ 23 4. An´alisis Estad´ıstico 27 4.1. An´alisis estad´ıstico del Periodograma . . . . . . . . . . . . . . . . . . . . . 27 4.2. Transformaci´on ................................. 29 4.3. AlgoritmoGen´etico............................... 33 4.3.1. Dise˜no del Algoritmo Gen´etico . . . . . . . . . . . . . . . . . . . . . 34 4.3.1.1. Codificaci´on.......................... 34 4.3.1.2. Par´ametros .......................... 35 4.3.1.3. Operadores .......................... 36 4.3.2. Funci´onobjetivo ............................ 37 4.3.3. Resultados................................ 39 5. Clusterizaci´on 40 5.1. Introducci´on................................... 40 5.2. Clusterizaci´on y M´etodos . . . . . . . . . . . . . . . . . . . . . . . . . . . . 41 5.2.1. Distancia de similitud . . . . . . . . . . . . . . . . . . . . . . . . . 42 5.3. Resultados.................................... 43 6. Validaci´on 46 6.1. Estudio de la Fiabilidad . . . . . . . . . . . . . . . . . . . . . . . . . . . . 46 6.1.1. Resultados................................ 47 6.2. Estudio del An´alisis del Diagn´ostico . . . . . . . . . . . . . . . . . . . . . 48 7. Conclusiones y l´ıneas futuras 50 A. Desarrollo 53 A.1.Hitosdelproyecto................................ 53 A.2.DiagramadeGantt............................... 54 B. Sistema BCI 55 viii ´ INDICE ´ INDICE B.1. Introducci´on. Sistema 10-20 . . . . . . . . . . . . . . . . . . . . . . . . . . 55 B.2.Artefactos.................................... 57 B.3.Infraestructura ................................. 60 B.3.1.Electrodosactivos............................ 60 B.3.2.Amplificador .............................. 60 B.3.3.Gorro .................................. 61 B.3.4.Geles................................... 61 B.3.5.Otros................................... 62 B.3.6. Software: BCI2000 . . . . . . . . . . . . . . . . . . . . . . . . . . . 62 C. An´alisis estad´ıstico del PSD 63 D. M´etodos 64 D.1. Mezcla de gaussianas (GMM) . . . . . . . . . . . . . . . . . . . . . . . . . 64 D.2. Distancias en los m´etodos de clustering . . . . . . . . . . . . . . . . . . . . 65 D.3. Z-test...................................... 65 D.4. Contraste para la diferencia de medias . . . . . . . . . . . . . . . . . . . . 65 D.5. Distancia Bhattacharyya . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66 ix 1. Introducci´on 1.3 Objetivo y Alcance del proyecto derivaciones y diferentes estados funcionales (estado de reposo, tarea activa...) pero cambia entre individuos. Adem´as se atribuye una posible dependencia del FEG con aspectos fisiol´ogicos del individuo (como la edad), y por tanto con propiedades de maduraci´on del volumen conductor del cerebro, como las caracter´ısticas del hueso, la geometr´ıa del cr´aneo y las conductancias de la piel y el cr´aneo. En consecuencia, el desarrollo de una Base Normativa requiere o bien un estudio de la variabilidad del EEG causado por motivos fisiol´ogicos (los cuales no han sido completamente demostrados actualmente, como la edad, la geometr´ıa del cr´aneo, etc), o bien una estratificaci´on de sujetos de la Base Normativa que tenga en cuenta expl´ıcitamente los cambios significativos en su EEG. Hasta la fecha todas las Bases Normativas realizan una estratificaci´on por edades [5] (aunque se ha visto que es ´util en sujetos maduros dado que el EEG es relativamente invariable este criterio no se satisface en el caso de Bases Normativas de ni˜nos) y otros autores han ajustado estos patrones dependiendo adem´as de la regi´on del cerebro analizada [6]. A estos inconvenientes se suma que la contribuci´on del FEG sobre la varianza total del EEG apenas alcanza el 42 % de los datos corregidos por la edad. Por tanto, actualmente en muchos de los estudios de EEG se trata de disminuir esta variabilidad que existe entre individuos de la misma edad [7]. M´as all´a de la restricci´on mencionada, los problemas principales asociados a la estratificaci´on del EEG en base a la edad se resume en: No linealidad de las relaciones entre el EEG y la edad del individuo. Heteroscedasticidad3de la predicci´on del error en la no linealidad mencionada. Vulnerabilidad en errores de calibraci´on. Todos estos aspectos apuntan a que la variabilidad interpersonal es una limitaci´on en la construcci´on de las Bases Normativas actuales. Este proyecto propone una alternativa a la estratificaci´on de se˜nal EEG por edades, la cual est´a basada en la aplicaci´on de una t´ecnica de clusterizaci´on de individuos en funci´on de un conjunto de caracter´ısticas de se˜nal extra´ıdas sobre el EEG. Como consecuencia, la poblaci´on normal queda dividida en clusters caracterizados por propiedades de la se˜nal. La gran ventaja de la t´ecnica es que permite la automatizaci´on completa del proceso de generaci´on de la Base Normativa. A su vez, realiza una estratificaci´on expl´ıcita basada en los cambios significativos del EEG. 1.3. Objetivo y Alcance del proyecto El objetivo de este proyecto fin de carrera es el dise˜no, desarrollo y validaci´on de una Base Normativa de EEG basal de 84 adultos sanos (comprendidos entre 18 y 30 a˜nos) en condiciones de adquisici´on de ojos cerrados. En esta direcci´on se ha dise˜nado un proceso innovador de creaci´on de una Base Normativa mediante el cual los sujetos que componen 3Heteroscedasticidad: Se produce cuando la varianza de las perturbaciones no es constante a lo largo de las observaciones. 4 1. Introducci´on 1.3 Objetivo y Alcance del proyecto Figura 1.2: Clusterizaci´on de una Base Normativa la Base Normativa se dividen en clases aplicando un algoritmo param´etrico construido sobre caracter´ısticas del EEG. A su vez se han creado las herramientas autom´aticas para comparar un individuo frente a la Base Normativa, y determinar si se encuentra dentro de una normalidad definida. Si el sujeto no concuerda con la Base Normativa, es posible determinar de que manera se desv´ıa de la normalidad y la magnitud de esta desviaci´on. La ventaja del proceso creado es que es autom´atico con lo que la inclusi´on de un individuo es directa y evita la intervenci´on humana en el proceso. La figura 1.3 muestra el proceso desarrollado y las etapas que se han llevado a cabo en este proyecto. Las cuales se describen a continuaci´on: el desarrollo de la Base Normativa comienza con la etapa de adquisici´on de EEG, la cual no ha formado parte del presente proyecto. En concreto, se ha trabajado con datos de se˜nal EEG recogidos por la empresa Nova Tech EEG (NTE). A esta etapa le sigue un tratamiento de se˜nal espec´ıfico para la extracci´on de las caracter´ısticas en las que sustenta la Base Normativa. Estos par´ametros son analizados mediante la t´ecnica de Split&Half para estudiar su fiabilidad. Posteriormente, se realiza un an´alisis estad´ıstico de dichas caracter´ısticas con el objetivo de obtener una distribuci´on estad´ıstica com´un a todas ellas, determin´andose la distribuci´on gaussiana como funci´on objetivo. La siguiente etapa define un m´etodo de clusterizaci´on de los sujetos en clases siguiendo un criterio de similitud en cuanto a caracter´ısticas de se˜nal EEG. Por ´ultimo, se concluye con una etapa de validaci´on, en la cual se analiza la sensibilidad de la Base Normativa y las herramientas de medida dise˜nadas como una prueba de diagn´ostico. 5 1. Introducci´on 1.4 Organizaci´on Figura 1.3: Procedimientos para el desarrollo y validaci´on de una base normativa de EEG. 1.4. Organizaci´on El presente documento est´a estructurado en siete cap´ıtulos, siendo este primero la introducci´on; el cap´ıtulo dos presenta los conceptos b´asicos sobre el electroencefalograma adem´as detalla el registro y procesado de se˜nal EEG llevado a cabo para la extracci´on de los patrones que caracterizan la Base Normativa; el cap´ıtulo tres detalla los estudios realizados sobre la fiabilidad de los datos adquiridos; el cap´ıtulo cuatro describe el an´alisis estad´ıstico realizado sobre la distribuci´on de las muestras utilizadas; el cap´ıtulo cinco expone el proceso de segmentaci´on realizado sobre los sujetos de la Base Normativa seg´un criterios de similitud en par´ametros de la se˜nal EEG; el cap´ıtulo seis describe la validaci´on de las herramientas de medida dise˜nadas junto con la Base Normativa y por ´ultimo en el cap´ıtulo siete se recopila la valoraci´on sobre el trabajo realizado, y se exponen las posibles extensiones del proyecto. 6 2. Registro y procesado de se˜nal EEG 2.1. Introducci´on El presente cap´ıtulo describe la se˜nal empleada en el desarrollo de este proyecto, su adquisici´on y procesado. La secci´on 2.2 presenta una introducci´on del electroencefalograma, describiendo el origen fisiol´ogico de la se˜nal EEG y de sus caracter´ısticas temporales y espectrales. A continuaci´on, la secci´on 2.3 describe los datos EEG utilizados en este proyec- to y el registro llevado a cabo. Y por ´ultimo, en la secci´on 2.4 se detalla el tratamiento de se˜nal realizado, una primera etapa de filtrado de artefactos y un procesado m´as concreto para la extracci´on de los par´ametros que recogen la Base Normativa de este proyecto. 2.2. Electroencefalograma El electroencefalograma (EEG) es el registro no invasivo de la actividad el´ectrica cerebral mediante electrodos colocados en la superficie del cuero cabelludo. Al conjunto de la actividad el´ectrica se le conoce como ritmo debido a su comportamiento oscilatorio. Concretamente la actividad de una sola neurona cortical no puede ser medida debido al grosor de las capas de tejido que componen la corteza cerebral. Sin embargo la actividad conjunta de millones de neuronas corticales, localizadas a una profundidad de varios mil´ımetros, producen un campo el´ectrico suficientemente fuerte para ser medido en la superficie cerebral, siendo esta profundidad dependiente de la potencia de la fuente neuronal. El campo el´ectrico es principalmente generado por corrientes que fluyen durante excitaciones sin´apticas1de las dendritas, la excitaci´on post-sin´aptica. La diversidad de ritmos EEG es inmensa y depende, entre otras cosas, del estado mental del sujeto (grado de atenci´on, de relajaci´on, etc). Los ritmos est´an convencionalmente caracterizados por: La amplitud de EEG: est´a relacionada con el grado de sincronizaci´on con el cual las neuronas corticales interaccionan. Una excitaci´on sincronizada se produce 1Sinapsis: uni´on intercelular producida por una descarga qu´ımica que origina una corriente el´ectrica en la membrana de la c´elula sin´aptica (emisora) una vez que este impulso alcanza el extremo del ax´on la propia neurona segrega neurotransmisores. 7 2. Registro y procesado de se˜nal EEG2.3 Descripci´on de los datos y Registro cuando las se˜nales originadas por cada neurona se suman de forma coherente en fase produciendo como resultado una se˜nal de amplitud elevada. Por el contrario, una excitaci´on as´ıncrona de las neuronas provoca un EEG irregular de baja amplitud. La amplitud es variable, de 20-60 µV (50 µV por t´ermino medio). La frecuencia de EEG: viene definida por la actividad proveniente del T´alamo. Esta parte del cerebro est´a formada por neuronas que poseen propiedades marcapasos, es decir, tienen la capacidad intr´ınseca para generar un patr´on r´ıtmico. En la actualidad, el registro de EEG se lleva a cabo mediante la colocaci´on de un gorro sobre el cuero cabelludo, en el cual se colocan un conjunto de electrodos. La disposici´on de los electrodos sobre el cuero cabelludo sigue el est´andar descrito por el Sistema Internacional 10-20, ver en anexo B.1. Cada par de electrodos registra las variaciones de potencial el´ectrico definiendo un canal de EEG. Cada una de estas se˜nales son amplificadas mediante el uso de una amplificador diferencial previo a la etapa de digitalizaci´on. La figura 2.2 representa una muestra de encefalograma compuesta por 16 canales. Tal y como se observa, la se˜nal es segmentada en fragmentos temporales denominados epochs (divididos por l´ıneas verticales en la figura). Figura 2.1: (Izquierda) Ejemplo de se˜nal EEG adquirida sobre 16 posiciones del cuello cabelludo. Cada fila representa un canal. (Derecha) Representaci´on visual en forma de mapa de color de la zona sombreada en la imagen izquierda 2.3. Descripci´on de los datos y Registro Los datos empleados en este proyecto fueron adquiridos por la empresa NTE.´ Estos se componen de la se˜nal EEG registrada en 84 adultos sanos, con edades comprendidas entre los 18 y 30 a˜nos. La selecci´on de los sujetos sigui´o unos criterios de exclusi´on basados en el historial psiqui´atrico, tanto de los sujetos como de sus familiares, adem´as de un an´alisis sobre posibles abusos de drogas, discapacidades f´ısicas o lesiones cerebrales. El registro 8 2. Registro y procesado de se˜nal EEG 2.4 Procesado de se˜nal se realiz´o durante una sesi´on (o trial) continua de unos 3 a 5 minutos de duraci´on, en la cual el sujeto se encontraba en condiciones de reposo, con los ojos cerrados y sentado c´omodamente en una habitaci´on silenciosa, con luz tenue. El protocolo de montaje junto con el registro de se˜nal EEG fue llevado a cabo por NTE (en el anexo B.3 se detalla de forma m´as extensa la infraestructura de montaje utilizado). Para la adquisici´on de se˜nal EEG se utilizaron 19 electrodos cuya colocaci´on sigui´o el est´andar prescrito por el Sistema Internacional 10-20.(FP1, FP2, F7, F3, Fz, F4, F8, T3, C3, Cz, C4, T4, T5, P3, Pz, P4, T6, O1, O2). Se realiz´o un montaje bipolar, en el cual, la se˜nal de referencia es adquirida en los l´obulos de las orejas, denominado ‘linked-ears reference’. Se comprob´o que la impedancia de los electrodos estuviera por debajo de 5 Kohmios asegurando una relaci´on se˜nal a ruido apropiada. El amplificador diferencial utilizado es el sistema de adquisici´on- 24 NeuroSearch empleando un conversor de 12-bit Anal´ogico-Digital (Lexicor Medical techology, Inc., bouder, CO). El cual permite un filtrado anal´ogico para eliminar la interferencia debido al ruido ambiental que proviene de las l´amparas fluorescentes y otros dispositivos que emiten ruido a trav´es de ondas de 60Hz, utilizando un filtro elimina-banda tipo Notch. En el sistema de digitalizaci´on (previo a la conversi´on Anal´ogico-Digital) se realiz´o un prefiltrado mediante un filtro antialiasing: filtro paso bajo con frecuencia de corte de 50Hz. Se seleccion´o una fMde sampleo de 128 Hz para garantizar que la fMes superior a dos veces la m´axima frecuencia presente en la se˜nal continua y as´ı satisfacer el Teorema de Muestreo de Nyquist: fM>2fmax (2.1) 2.4. Procesado de se˜nal 2.4.1. Filtrado de artefactos La se˜nal EEG se caracteriza por su reducida magnitud (del orden de las decenas de uV), por lo que es f´acilmente contaminada por ruido e interferencias. De tal manera, un aspecto especialmente importante en el tratamiento de se˜nal EEG consiste en la eliminaci´on de ruido y artefactos para mejorar la relaci´on se˜nal a ruido (en ingl´es Signal to Noise Ratio SNR o S/N). Un correcto filtrado de artefactos permite, adem´as de eliminar se˜nales no deseadas, minimizar la variabilidad interpersonal del EEG. Este es un aspecto cr´ıtico en el desarrollo de una Base Normativa, dado que en definitiva una BN agrupa se˜nales EEG provenientes de diferentes sujetos con el objetivo de definir un patr´on de referencia. Un artefacto es un registro de actividad no cerebral que aparece dentro del EEG alterando en gran medida la se˜nal adquirida debido a que su amplitud suele ser mucho mayor que la propia actividad cerebral. Los artefactos se presentan inevitablemente en todo registro de EEG pudiendo tener un origen biol´ogico, instrumental o medioambiental. Para m´as informaci´on ver el anexo ??. El filtrado de artefactos requiere de un importante conocimiento y experiencia. Este proceso consiste en un exhaustivo reconocimiento de artefactos mediante la ayuda de un software de visualizaci´on, procediendo manualmente 9 2. Registro y procesado de se˜nal EEG 2.4 Procesado de se˜nal a eliminar aquellos epochs afectados por un posible artefacto. Como consecuencia, la longitud de se˜nal EEG original se ve sensiblemente reducida. El filtrado de artefactos realizado por Nova Tech EEG se compone de dos fases. La primera de ellas, consiste en una inspecci´on visual del EEG digitalizado en una pantalla de alta-resoluci´on para identificar los epochs que contengan artefactos visibles. A continuaci´on, se procede a un filtrado m´as exhaustivo mediante la aplicaci´on Eureka 2, ayuda a la identificaci´on de artefactos m´as sutiles la cual necesita siempre la intervenci´on humana para la eliminaci´on. 2.4.2. Extracci´on de patrones Una Base Normativa recoge un conjunto de patrones para caracterizar el EEG de un sujeto sano. La extracci´on de patrones relevantes es particularmente crucial cuando el prop´osito es dise˜nar un sistema que permita discriminar entre un EEG ’normal’ y uno patol´ogico, como es una de las aplicaciones de este proyecto. Una vez la se˜nal EEG est´a libre de artefactos se prosigue con un tratamiento m´as espec´ıfico para la extracci´on de caracter´ısticas de la se˜nal. Los patrones m´as utilizados en las investigaciones mediante EEG se clasifican en patrones espectrales y patrones de conectividad. Los primeros se definen a partir del an´alisis del espectro de frecuencias en una ventana temporal y describen la intensidad del campo electromagn´etico de esa localizaci´on. Los patrones de conectividad son medidas m´as complejas que consideran caracter´ısticas espacio-temporales, que implican varias localizaciones simult´aneamente y describen la intensidad de la conexi´on entre esas regiones cerebrales. Los patrones de se˜nal EEG m´as comunes desde un punto de vista cl´ınico son los patrones espectrales de frecuencia. El tratamiento frecuencial de EEG comienza con el c´alculo de la densidad espectral de potencias (PSD). La PSD de la se˜nal EEG se ha llevado a cabo mediante el estimador espectral de potencia, denominado periodograma. El c´alculo del periodograma asume que la se˜nal es estacionaria3y peri´odica, es decir, que se puede identificar un patr´on que se va repitiendo con cierta periodicidad (a estos segmentos se les denomina epochs). Se seleccion´o una longitud de epoch igual a un segundo, coincidiendo con la duraci´on elegida en la limpieza de artefactos. Adem´as, esta longitud es usada habitualmente en anteriores estudios que analizan la estacionariedad del EEG. Debido a la propiedad no estacionaria del EEG, se asume en su an´alisis que el EEG est´a compuesto por segmentos consecutivos con propiedades estacionarias. El procesado se realiz´o sobre las 19 se˜nales EEG filtradas de artefactos por el proceso descrito en la subsecci´on anterior. La estimaci´on del PSD se llev´o a cabo mediante el periodograma modificado, en el cual cada una de las se˜nales es enventanada en segmentos de un segundo. Se utiliz´o la ventana de Tukey, o coseno alzado truncado, cuya expresi´on se define para una longitud L: 2Aplicaci´on para el an´alisis y filtrado de artefactos en la se˜nal EEG, de Nova Tech EEG. 3Se˜nal estacionaria: se˜nal cuyos par´ametros estad´ısticos permanecen constantes sobre el tiempo. 10 2. Registro y procesado de se˜nal EEG 2.4 Procesado de se˜nal Wtukey(i)) =            1 2(1 + cos(2π α[x−α/2])) 0 ≤x < αL 2 1αL 2≤x<L(1 −α 2) 1 2(1 + cos(2π αL[x−1 + α/2])) (1 −α 2)≤x≤L Wtukey(i) = 1 2(1+cos π Tt) (2.2) El valor de longitud de ventana L se tom´o igual a 128 puntos, los cuales representan un segundo de la se˜nal EEG en tiempo muestreada a una frecuencia de 128 Hz. El c´alculo de la Transformada Discreta de Fourier (Discrete Fourier Transform, DFT) se llev´o a cabo mediante la aplicaci´on de zero-padding4, lo que permiti´o emplear un total de 2048 puntos en el c´alculo de la DFT. Obteni´endose una resoluci´on frecuencial de: Mf=Fs N=128Hz 2048muestras = 0,0625Hz/muestra (2.3) La densidad espectral de potencias del EEG se caracteriza por la forma representada en la figura 2.2. Aproximadamente el 98 % de la potencia en la se˜nal se distribuye dentro del intervalo 0.5 y 30Hz, por lo que ´este intervalo es el analizado en este proyecto. Figura 2.2: Periodograma t´ıpico de una se˜nal EEG, el cual representa la distribuci´on de la potencia de la se˜nal sobre el dominio frecuencial. A trav´es del c´odigo de colores se representan las bandas de frecuencia t´ıpicas del espectro de EEG: Delta, Theta, Alpha y Beta. En general, la densidad espectral de potencia de EEG es fuertemente variable de forma intrapersonal, es decir, la distribuci´on de potencia en el dominio frecuencial var´ıa dependiendo de cada persona. Por este motivo, no resulta ´util en t´erminos de comparaci´on tomar un patr´on o (par´ametro) de se˜nal EEG extra´ıdo de una ´unica frecuencia o bin del espectograma ,sino que usualmente se agrupan en bandas frecuenciales. Tradicionalmente 11 2. Registro y procesado de se˜nal EEG 2.4 Procesado de se˜nal Delta Theta Alfa Beta Intervalo [0.5-3.5 Hz] [3.5- 7 Hz] [7-13 Hz] [13-30 Hz] Tabla 2.1: Definici´on bandas frecuenciales del espectro en la se˜nal EEG. se han definido cuatro bandas frecuencia, mostradas en la tabla 2.1. A pesar de que estas bandas de frecuencia han sido estudiadas durante d´ecadas y la gran mayor´ıa de trabajos con el EEG est´an descritos en base a ella, no existe consenso absoluto respecto a los l´ımites exactos de estos ritmos frecuenciales. Como resultado, existe pluralidad en la definici´on tanto del n´umero de bandas frecuenciales como de las frecuencias lim´ıtrofes. En este proyecto, se emplean las bandas definidas anteriormente, las cuales han sido definidas en anteriores Bases Normativas [1]. Para determinar los par´ametros de se˜nal que recogen la Base Normativa de este proyec- to, se estudiaron los diferentes patrones patol´ogicos conocidos en el EEG [8]. Finalmente, se determin´o por su relevancia en la caracterizaci´on mayoritaria de las patolog´ıas estudiadas, la potencia absoluta. La potencia absoluta (PA) se define como la cantidad media de potencia (en uV 2) en una banda frecuencial y es calculada como el ´area representada bajo la curva de PSD en el intervalo definido por dicha banda. En la figura 2.3 se representa el ´area correspondiente a la potencia absoluta de cada banda frecuencial. En consecuencia, Figura 2.3: Potencias Absolutas frecuenciales sobre el periodograma de EEG cada epoch de un canal de EEG es analizado en el dominio frecuencial para la extracci´on de la potencia absoluta en cada una de las bandas: delta, theta, alpha y beta. Como resultado se obtienen 4 conjuntos de valores por cada c´alculo de potencia en los canales EEG, es decir, 19*4=76 variables de potencia por sujeto. En las pr´oximas secciones se tomar´a la siguiente nomenclatura: 4Zero-Padding: consiste en rellenar con ceros la se˜nal en tiempo a partir de su ´ultimo valor hasta alcanzar el tama˜no deseado para mejorar la resoluci´on espectral en su Transformada de Fourier. 12 2. Registro y procesado de se˜nal EEG 2.4 Procesado de se˜nal Pi δ: Variable de potencia en la banda δ, del canal i. Pi θ: Variable de potencia en la banda θ, del canal i. Pi α: Variable de potencia en la banda α, del canal i. Pi β: Variable de potencia en la banda β, del canal i. {PEEG }={P1 δ, P1 θ, P1 α, P1 β, ..., Pn β}: Representa el conjunto de las variables de potencia absoluta, con n=19 canales. 13 3. Fiabilidad 3.3 M´etodo 2: Split&Half cada una de las 76 distribuciones de potencia que caracterizan a un participante, ser´ıa la siguiente: denotando Er=rc xx(i)−rc xx(i−1) la diferencia del coeficiente de correlaci´on ¯rxx como resultado de computar Split&Half con iei−1 iteraciones sobre una variable banda-canal denotada por x: E(i) = p19[E2 r+E2 r+E2 r+E2 r] (3.6) E(i) = p76 ∗E2 r(3.7) Siendo Emax(i)=0.005: Er<r0,0052 76 = 0,0081 (3.8) Es decir, este valor de error impone una cota m´axima de la diferencia del coeficiente de correlaci´on ¯rxx como resultado de computar Split&Half con iei−1 iteraciones sobre una variable banda-canal denotada por xpor ebajo de las mil´esimas. A continuaci´on se muestran los resultados del c´alculo de la fiabilidad. La figura 3.5 representa el valor del par´ametro E(i) en funci´on del n´umero de iteraciones del algoritmo Split &Half superponiendo de la soluci´on de cada participante. Las l´ıneas rojas representan el valor de error m´aximo admisible Emax = 0,005. ´ Este valor se cumple entre las 800 y 1000 iteraciones para el conjunto de participantes. Por tanto, se toma como par´ametro M del algoritmo de la ecuaci´on 3.4 el valor de 1000 iteraciones. En la figura 3.6 se muestra los resultados de fiabilidad de cada participante. En la cual se represnta la distribuci´on de los valores de rxx resultante en las 76 variables de potencia absoluta. Se concluye que los resultados de fiabilidad se encuentran en torno a rxx = 0,7. A diferencia de los resultados anteriores, se consigue aumentar significativamente el valor de fiabilidad, obteni´endose una fiabilidad ajustada con la literatura. Sin embargo, los resultados obtenidos no alcanzan un nivel de fiabilidad admisible en aplicaciones cl´ınicas. Por tanto, se requiere de una t´ecnica que favorezca esta fiabilidad. Te´oricamente la fiabilidad viene limitada gravemente por la duraci´on de la se˜nal analizada. En particular, el presente an´alisis trabaja con datos de longitud muy reducida (entre 60 y 120 muestras, epochs de un segundo). En consecuencia, la longitud de la se˜nal EEG es un valor cr´ıtico y determinante en lo que a fiabilidad se refiere. 20 3. Fiabilidad 3.3 M´etodo 2: Split&Half Figura 3.5: Representa el valor del par´ametro E(i) en funci´on el n´umero de iteraciones i, el cual define el error cometido en el computo del algoritmo Split &Half al tomar i−1 iteraciones en vez de i. Las l´ıneas verticales rojas representan el valor de error m´aximo admisible Emax = 0,0005. En ella se ha superpuesto la funci´on E(i) junto con la cota Emax para cada uno de los 84 participantes. Figura 3.6: Representa los resultados de fiabilidad al aplicar Split &Half con un valor del par´ametro Migual a 1000 iteraciones. Los resultados se representan mediante el diagrama estad´ıstico denominado Box-Plot. El cual es un gr´afico representativo la distribuci´on de la variable rho calculada sobre cada participante. La caja central indica el rango en el que se concentra el 50 % central de los datos, sus extremos son el primer y tercer cuartil de la distribuci´on y la l´ınea central en la caja es la mediana. Los extremos delimitan el 95 % central de los datos, y las muescas en rojo representan los datos considerados como espurios en la distribuci´on. 21 3. Fiabilidad 3.4 Ampliaci´on de se˜nal EEG 3.4. Ampliaci´on de se˜nal EEG Debido a una restricci´on en la longitud de la se˜nal EEG adquirida y por tanto del n´umero de epochs, los resultados de fiabilidad obtenidos no fueron los deseados. Por tanto, fue necesario encontrar una estrategia para aumentar la duraci´on de la se˜nal EEG de cada uno de los participantes a˜nadiendo epochs ‘virtuales’con el objetivo de mejorar los valores de fiabilidad intra-personal. Un primer dise˜no contemplaba el uso de ventanas temporales deslizantes de un segundo de duraci´on aplicando un porcentaje de solapamiento entre ventanas consecutivas. De este manera, se obtiene un valor de potencia por cada ventana aplicada. Al fijar la longitud de ventana junto con el porcentaje de solapamiento entre ellas, el n´umero de ventanas resultantes queda limitado por la duraci´on de la se˜nal temporal en segundos. Este primer dise˜no fue abandonado porque el n´umero de epochs ’virtuales’ resultantes quedaba limitado a un n´umero insuficiente para su aplicaci´on. Se dise˜n´o una segunda t´ecnica donde los l´ımites descritos por el m´etodo anterior quedaban resueltos. Para ello, se a˜nadi´o una reordenaci´on aleatoria de la se˜nal original segmentada en epochs de un segundo de duraci´on. Un epoch ’virtual’ resulta de unir dos mitades de epochs consecutivos de la se˜nal original reordenada, en consecuencia se obtienen n−1 epochs ’virtuales’ sobre una se˜nal de duraci´on n. Este procedimiento se recomputa hasta conseguir una se˜nal con una fiabilidad deseada calculada mediante la t´ecnica de Split & Half. En la Figura 3.7 se representa el algoritmo dise˜nado, el cual est´a descrito por el siguiente algoritmo iterativo: 1. Reordenar la secuencia original: permutar los epochs iniciales aleatoriamente. 2. Agrupar pares de epochs consecutivos de dos en dos. Ejemplo: [EPOCHi−EPOCHi+1]. 3. Dividir cada pareja en cuatro fragmentos cortando cada epoch por la mitad. Ejemplo: [EPOCHi(a)−EPOCHi(b)−EPOCHi+1(a)−EPOCHi+1(b)]. 4. Crear el nuevo epoch ‘virtual’seleccionando aleatoriamente entre las siguientes casos: uniendo la segunda mitad del primer ıepoch con la primera mitad del segundo, la primera mitad del primer epoch con la segunda mitad del segundo. Ejemplo: [EPOCHi(b)−EPOCHi+1(a]) o bien [EPOCHi(a)−EPOCHi+1(b)]. 5. Volver al punto 1 hasta alcanzar un valor de fiabilidad estipulado, en cada una de las caracter´ısticas extra´ıdas por separado. Este procedimiento no altera la distribuci´on de la se˜nal en el tiempo, debido a que la adquisici´on se produce en un estado continuo de reposo. Adem´as, como resultado del 22 3. Fiabilidad 3.4 Ampliaci´on de se˜nal EEG filtrado de artefactos se introducen discontinuidades de tiempo, debido a la eliminaci´on de epochs alterados con artefactos, obteni´endose una secuencia de epochs con discontinuidades de tiempo entre ellos. 3.4.1. Resultados En el dise˜no del algoritmo se determin´o como criterio de parada el m´ınimo valor de fiabilidad obtenido en las potencias de las cuatro bandas frecuenciales de los 19 canales que componen los par´ametros de EEG. Se fij´o un valor del coeficiente de correlaci´on m´ınimo de rxx = 0,8 (umbral m´ınimo para estudios cl´ınicos). Es decir, en cada una de las 76 se˜nales de potencias que caracterizan el EEG se obtuvo una fiabilidad de al menos el 80 %. En la figura 3.8 se muestra el resultado de fiabilidad obtenido por participante en cada una de las 76 distribuciones en particular. Como resultado se ha obtenido una fiabilidad media de rxx = 0,894. A continuaci´on, se describen estudios relacionados con la aplicaci´on de este m´etodo. El algoritmo se aplic´o en cada participante por separado, como resultado se obtuvo unas se˜nales ampliadas diferentes para da participante. La figura 3.10 muestra la ampliaci´on resultante para el conjunto de participantes. El m´etodo utilizado fue dise˜nado con el objetivo de mejorar fiabilidad. ´ Este debe cumplir el prop´osito para el cual fue dise˜nado, sin modificar la distribuci´on de la se˜nal original. En un primer an´alisis, se comprob´o gr´aficamente la variaci´on que sufre la media para diferentes longitudes de se˜nal. En la figura 3.11 se muestra esta variaci´on para longitudes N=56, 606 y 1101 epochs. Adem´as, se realiz´o un an´alisis estad´ıstico sobre la distribuci´on de las caracter´ısticas antes y despu´es del proceso de ampliaci´on de epochs, mediante el an´alisis estad´ıstico K-S. La gr´afica 3.9 muestra el n´umero de distribuciones de potencia por participante que mantienen su distribuci´on. Se observa unos resultados positivos, de las 76 distribuciones de potencia que componen una participante en la mayor´ıa de participantes se mantiene su 90 % de sus distribuciones. En resumen, en este cap´ıtulo se ha alcanzado unos resultados de fiabilidad sobre los datos de potencias en las bandas frecuenciales : Delta, Theta, Alpha y Beta, para los 19 canales adquiridos por encima del 80 %. Estos resultados demuestran una fiabilidad intrapersonal aceptables en un contexto cl´ınico. 23 3. Fiabilidad 3.4 Ampliaci´on de se˜nal EEG Figura 3.7: Algoritmo para a˜nadir epochs virtuales al conjunto inicial. 24 3. Fiabilidad 3.4 Ampliaci´on de se˜nal EEG Figura 3.8: Resultados de Split &Half para todo el conjunto de participantes. Cada fila representa la potencia banda-canal, de arriba a bajo: P1 D,P1 T,P1 A,P1 B,P2 D...P1 B9 Figura 3.9: Representa el n´umero de variables de las 76 que componen: {PEEG }que mantienen su distribuci´on estad´ıstica tras la aplicaci´on del m´etodo de ampliaci´on mediante epochs ’virtuales’ en cada uno de los participantes. La l´ınea horizontal en verde muestra la cota m´axima de esta gr´afica situada en el valor 76 25 3. Fiabilidad 3.4 Ampliaci´on de se˜nal EEG Figura 3.10: Muestra la longitud de la se˜nal EEG en n´umero de epochs adquirida en dada uno de los sujetos, junto con el n´umero de epochs resultantes tras la aplicaci´on del m´etodo de ampliaci´on mediante epochs ’virtuales’, empleado para conseguir un valor de fiabilidad de Split &Half superior a 0.8. Figura 3.11: Representa la la media en las cuatro variables de potencia en cada uno de los 19 canales, de las siguientes tres se˜nales EEG: el EEG original adquirido en el participante n´umero uno de duraci´on 56 epochs; la se˜nal anterior ampliada con una longitud de 606 epochs y la anterior se˜nal ampliada hasta 1101 epochs. De este modo, se puede estudiar si la media de las distribuciones var´ıan mediante la t´ecnica de ampliaci´on mediante epochs ’virtuales’. 26 4. An´alisis Estad´ıstico Una vez alcanzados unos niveles aceptables sobre la fiabilidad en las caracter´ısticas extra´ıdas en la se˜nal EEG, el siguiente proceso consiste en analizar estad´ısticamente los patrones espectrales extra´ıdos. El objetivo es conocer la funci´on de distribuci´on de las mismas y facilitar en la medida que sea posible los siguientes procesos. En secci´on 4.1 se describe la distribuci´on de probabilidad de los datos de partida. Las potencias absolutas en las diferentes bandas frecuenciales. Para ello se realiza un an´alisis estad´ıstico del periodograma, y en particular para cada una de las potencias absolutas, concluyendo que una distribuci´on que no son modeladas a trav´es de ninguna distribuci´on conocida. Con el objetivo de facilitar los c´alculos posteriores y alinearse con los actuales est´andares en Bases Normativas [14] es ajustar las distribuciones de potencias a la distribuci´on gaussiana. El proceso de aproximaci´on de cada una de las distribuciones hacia una gaussiana se lleva a cabo en dos etapas: la primera consiste en una transformaci´on de los datos de partida (secci´on 4.2), y la segunda un filtrado de los datos que acerque la distribuci´on a la gaussianidad mediante un algoritmo gen´etico dise˜nado ad hoc (secci´on 4.3). Como resultado, se obtiene para cada participante, un conjunto de 76 distribuciones normales unidimensionales, (definidas por los 19 canales, y 4 bandas frecuenciales). 4.1. An´alisis estad´ıstico del Periodograma Cuando se hace uso de una se˜nal estoc´astica como el EEG, una de las principales cuestiones es identificar la funci´on de densidad de distribuci´on (PDF, Probability Density Function) que permita una caracterizaci´on adecuada de la se˜nal. En particular, en este ep´ıgrafe se analiza la PDF de las caracter´ısticas extra´ıdas a partir de la densidad espectral de potencia del EEG. Estas caracter´ısticas est´an formadas por las potencias absolutas y relativas de las cuatro bandas de frecuencias δ,θ,αyβ. Tal y como se describe en la secci´on 2.4.2, una potencia absoluta es el ´area representada bajo la curva de PSD en una banda frecuencial. El c´alculo de la PA se ha llevado a cabo integrando de forma trapezoidal la funci´on PSD en un intervalo de frecuencias. El espectro de potencia (PSD) es no negativo por definici´on y seg´un estimadores espectrales de se˜nales estoc´asticas no est´a normalmente distribuido. Ha sido demostrado 27 4. An´alisis Estad´ıstico 4.1 An´alisis estad´ıstico del Periodograma te´oricamente que la parte real e imaginaria de la DFT de una se˜nal y(t) estacionaria pueden ser considerados independientes y modelados como una distribuci´on gaussiana de media cero. Por tanto, el m´odulo de la DFT (periodograma) puede ser modelado como la suma de estas variables al cuadrado. Bajo esta asunci´on cada bin del periodograma |Y(λ, k)|2est´a distribuido por una funci´on Chi-Cuadrado con dos grados de libertad, χ2(2) (Para una informaci´on m´as detallada consultar el anexo C). Esta distribuci´on puede ser interpretada como una distribuci´on exponencial: f|Y(λ,k)|2(x) = U(x) σ2 Y(λ, k)exp(−x/σ2 Y(λ, k)) (4.1) Debido a la ampliaci´on de la se˜nal EEG realizada mediante epochs virtuales, se estudi´o de forma experimental si los bins de los periodogramas resultantes cumpl´ıan la distribuci´on caracterizada por una funci´on χ2(2). Para ello, se hizo uso de la prueba de bondad de ajuste de Kolmog´orov-Smirnov (K-S) mediante la cual se comprob´o la distribuci´on χ2(2) en cada bin de frecuencia de 0.5 a 30 Hz. Una vez caracterizados los bins de la PSD, se procede a estudiar la distribuci´on estad´ıstica de los patrones espectrales en particular, es decir, de cada una de las potencias en las cuatro bandas frecuenciales calculadas en cada canal EEG. Tal y como se ha indicado anteriormente, cada una de las potencias calculadas, resultan de integrar de forma trapezoidal la curva PSD en una banda frecuencial determinada. La distribuci´on de las potencias no se conoce de forma directa. En consecuencia, el objetivo siguiente se centra en modelar las distribuciones de potencia mediante una distribuci´on conocida. Con el objetivo de facilitar los c´alculos posteriores y alinearse con los actuales est´andares en Bases Normativas se ha optado por ajustar las distribuciones de potencias a la distribuci´on gaussiana. En la pr´actica es muy frecuente buscar la gaussianidad en datos resultantes de una medici´on. Debido a la premisa b´asica de que la mayor´ıa de los fen´omenos tienen una distribuci´on normal. El teorema Central de L´ımite es una evidencia te´orica sobre este supuesto cuando se suman un conjunto de variables independientes. Aunque no es aplicable en el caso de potencias porque son definidas positivas. Dadao que no se conoce la distribuci´on real de las potencias, se eligi´o la gaussianidad como distribuci´on objetivo a alcanzar, por su simplicidad y porque facilita la aplicaci´on de numerosas herramientas matem´aticas (notar que el segundo objetivo en la construcci´on de una Base Normativa es proporcionar las herramientas matem´aticas para poder comparar sujetos frente a la Base Normativa). Adem´as, la importancia de aproximar par´ametros de se˜nal EEG hacia la normalidad ha sido enfatizada anteriormente por varios autores [15],[16]. Para alcanzar la distribuci´on gaussiana se estudi´o transformar los datos originales, como se ver´a en la siguiente secci´on. 28 4. An´alisis Estad´ıstico 4.2 Transformaci´on 4.2. Transformaci´on El asunto de alcanzar distribuciones normales en un conjunto de observaciones de distribuci´on desconocida, ha sido usualmente resuelto mediante la transformaci´on de los mismos. Esta opci´on tiene la propiedad de que las conclusiones obtenidas con los datos transformados tambi´en se aplican a los datos originales siempre y cuando la transformaci´on sea invertible. Este proceso requiere el conocimiento a priori sobre la distribuci´on de partida para determinar el tipo de transformaci´on m´as adecuada. En primer lugar, se analiz´o gr´aficamente el comportamiento estad´ıstico de las distribuciones de potencias en cada banda frecuencial. A modo de ejemplo, la Figura 4.1 representa el histograma de la potencia absoluta en la banda Delta en el canal F3 de un participante, de la se˜nal original (izquierda), y de la se˜nal ampliada (derecha). La distribuci´on representada se caracteriza por una fuerte asimetr´ıa positiva (las frecuencias m´as altas se encuentran en el lado izquierdo de la media, mientras que en el derecho hay frecuencias m´as peque˜nas) Figura 4.1: El histograma de la izquierda representa la distribuci´on de la potencia absoluta en la banda δ calculada sobre el canal F3 de la se˜nal EEG adquirida en el sujeto n´umero 4. El histograma de la derecha representa la distribuci´on de la potencia en la banda δcalculada sobre la se˜nal ampliada mediante epochs ’virtuales’ sobre el canal F3 del sujeto n´umero 4. En conclusi´on a los resultados obtenidos gr´aficamente, la funci´on de transformaci´on aplicada sobre estas distribuciones de datos, tiene que corregir la fuerte asimetr´ıa positiva con el objetivo de aproximar la distribuci´on a una sim´etrica como la normal. En la literatura existen diversas funciones de transformaci´on usadas en par´ametros de la se˜nal EEG [17]. Para el caso de distribuciones de datos con asimetr´ıa positiva se utilizan las 29 4. An´alisis Estad´ıstico 4.3 Algoritmo Gen´etico 4.3.1.3. Operadores Los operadores de un Algoritmo Gen´etico son empleados para obtener la poblaci´on de la siguiente generaci´on a partir de la poblaci´on actual. ´ Estos son: la selecci´on, el cruce y la mutaci´on. A continuaci´on se describe el dise˜no llevado a cabo para cada uno de ellos. Selecci´on: La selecci´on consiste en la elecci´on de los individuos de la poblaci´on para efectuar el operador de cruce y posteriormente la mutaci´on. En este dise˜no, la selecci´on viene determinada por el concepto de elitismo, que consiste en la idea de que los (Nelite) mejores individuos de la actual generaci´on pasen a la siguiente generaci´on sin modificarse por los operadores de cruce o mutaci´on. Por el contrario, los (Nmalos) peores individuos de la poblaci´on no pasar´an a la siguiente generaci´on. Sobre el resto de la poblaci´on se seleccionan de forma aleatoria cada pareja de padres para llevar a cabo el operador de cruce. Cruce: Crea una generaci´on de individuos nuevos a partir de la informaci´on de sus ancestros: padre y madre. El individuo resultante se forma aplicando la siguiente t´ecnica: de forma alternada se copian un ´ındice del padre y un ´ındice de la madre, comenzando por el extremo izquierdo de la secuencia del padre y por el extremo derecho de la secuencia de la madre. De tal forma se va completando la secuencia del hijo resultante, quedando la primera mitad de la secuencia del padre junto con la segunda mitad de la secuencia de la madre (Figura 4.7). Notar que la secuencia resultante no admite la repetici´on del mismo ´ındice, en el caso de intentar copiar un ´ındice ya incluido en la secuencia del hijo, se procede a copiar el ´ındice consecutivo. Figura 4.7: Operador de cruce que determina a partir de dos individuos de la poblaci´on seleccionada un individuo nuevo. Para ello se toma la informaci´on de los dos individuos y se copia para el nuevo hijo, tomando la primera mitad de la secuencia del padre y la primera mitad de la secuencia de la madre. Mutaci´on: Se aplica a cada hijo de manera individual y consiste en el intercambio de posiciones de dos ´ındices seleccionados aleatoriamente (Figura 4.8). El n´umero de mutaciones que se aplica a cada individuo se selecciona de forma aleatoria en cada iteraci´on (Pindices). Todos los hijos resultantes del cruce tienen la misma probabilidad de mutaci´on, la cual viene dada por el par´ametro global Nmutaci´on. Notar que los operadores de cruce y mutaci´on no alteran la longitud de la secuencia que codifica a un individuo, ´unicamente reordenan los ´ındices de la secuencia original 36 4. An´alisis Estad´ıstico 4.3 Algoritmo Gen´etico Figura 4.8: Operador de mutaci´on que se realiza sobre los hijos resultantes del operador de cruce. El individuo mutado resulta del intercambio de una pareja de ´ındices seleccionados aleatoriamente. El n´umero de mutaciones que se realiza sobre el mismo sujeto viene determinado aleatoriamente en cada iteraci´on del Algoritmo por Pindices. de epochs. El operador de cruce explota las buenas propiedades de los individuos, y sus efectos decrecen con la convergencia del Algoritmo Gen´etico. Por el contrario, la mutaci´on permite escapar de ´optimos locales explorando el espacio de b´usqueda. 4.3.2. Funci´on objetivo La funci´on objetivo representa la capacidad de adaptaci´on para cada individuo devolviendo un n´umero real proporcional a su nivel de adaptaci´on. El Algoritmo Gen´etico dise˜nado tiene como objetivo obtener gaussianidad en las distribuciones de potencias que caracterizan a un sujeto de la Base Normativa, mediante el filtrado de un conjunto de epochs de la se˜nal EEG. En primer lugar, es necesario definir un filtro de ’espurios’ (desde el punto de vista de la gaussianidad) y a continuaci´on dise˜nar el c´alculo de la funci´on objetivo. Dado que un individuo es codificado por una secuencia de ´ındices que determina la ordenaci´on de los epochs de la se˜nal EEG, el filtrado de epochs es aplicado directamente sobre cada individuo de la poblaci´on. El filtro es representado por una m´ascara centrada en la secuencia de ´ındices (Figura 4.9). La longitud de la m´ascara es constante en cada generaci´on de individuos, cuyo valor es igual a la longitud de la secuencia que codifica a cada individuo N(igual al n´umero de epochs de la se˜nal EEG) menos el valor de %Esp*N, el cual representa el n´umero de epochs que se desean filtrar. Por tanto, a cada lado de la m´ascara se encuentra la mitad de los ´ındices sin seleccionar N∗%Espurios/2. El filtrado de ´ındices conlleva un filtrado de epochs sobre la se˜nal EEG. A partir de cada epoch se ha determinado 76 variables que representan sobre cada canal la potencia absoluta en las cuatro bandas frecuenciales. Por lo tanto, un epoch determina una variable multidimensional de 76 dimensiones. En esta sentido, la b´usqueda de gaussianidad queda reflejado por la b´usqueda de una distribuci´on multigaussiana de 76 dimensiones. Debido al elevado n´umero de dimensiones y por tanto a la complejidad que conlleva alcanzar esta distribuci´on, se decidi´o simplificar el problema en 76 variables unidimensionales, sin violar la dependencia que existe entre ellos. Es decir, un epoch eliminado sobre una distribuci´on 37 4. An´alisis Estad´ıstico 4.3 Algoritmo Gen´etico Figura 4.9: Filtro sobre la secuencia de ´ındices. Est´a determinado por una m´ascara central sobre la secuencia. El n´umero de ´ındices que no son seleccionados sobre la secuencia viene determinado por N∗%Espurios, quedando en los extremos de la m´ascara la mitad de ´estos. unidimensional de potencias conlleva la eliminaci´on de ´este sobre las 75 distribuciones marginales restantes. De este modo se consigue alcanzar gaussianidad simult´aneamente en cada una de las distribuciones marginales. Para cuantificar la gaussianidad sobre la secuencia de epochs filtrados se hace uso del test de bondad Anderson-Darling, el cual devuelve un p-valor que indica la aproximaci´on a una distribuci´on gaussiana. Un p-valor igual o superior a 0.05 indica gaussianidad a un nivel de significancia del 95 %. El test Anderson-Darling es aplicado en cada una de las 76 distribuciones unidimensionales que determina la secuencia de epochs filtrados que resultan al aplicar la m´ascara sobre una secuencia de ´ındices. La funci´on objetivo se dise˜n´o con el prop´osito de maximizar la gaussianidad de las distribuciones unidimensionales de potencias, que se traduce a maximizar los 76 p-valores resultantes del test Anderson-Darling aplicado sobre cada una de las variables de potencia. En un primer dise˜no, se determin´o como funci´on objetivo la suma de los 76 p-valores para maximizar la gaussianidad de cada una simult´aneamente. Esta funci´on maximiza el resultado global ponderando cada p-valor de forma equitativa. Sin embargo, el objetivo de este Algoritmo no es aumentar la suma de los p-valores, sino que cada uno de ellos alcance al menos un valor igual a 0.05 (umbral que determina la gaussianidad al 95 % de significancia). Por tanto, se dise˜n´o otro funci´on que contemplara este valor umbral para cada uno de los p-valores. La funci´on objetivo a maximizar consiste en la suma de los p-valores por separado potenciando aquellas distribuciones con un p-valor por debajo del valor umbral de 0.05. Para ello se dise˜n´o la funci´on objetivo representada en la Figura 4.10, en la cual aquellas variables unidimensionales gaussianas (con un p-valor≥0,05) se establece su p-valor igual a un valor constante P. Por el contrario si la distribuci´on no es gaussiana (su p-valor <0,05) se mantiene su p-valor inicial. El valor de Pse fij´o por el valor m´aximo que puede tomar la funci´on objetivo si todas las distribuciones menos una han alcanzado la gaussianidad, es decir 0.05*3.75, en la pr´actica se tom´o P= 4. La convergencia del Algoritmo Gen´etico se determina cuando se obtiene una soluci´on ´optima. Por tanto se define un criterio de parada definido por la gaussianidad de las 76 variables de potencias unidimensionales, es decir, cuando el resultado de la funci´on objetivo es ≥4∗76. Como resultado se obtiene una secuencia de epochs ’virtuales’ que definen 76 distribuciones de potencias modeladas por una distribuci´on gaussiana unidimensional. 38 4. An´alisis Estad´ıstico 4.3 Algoritmo Gen´etico Figura 4.10: Funci´on Objetivo obtenida por el an´alisis de gaussianidad en cada una de las 76 distribuciones de potencias. La funci´on viene determinada por la suma de los p-valores resultantes al aplicar el test de gaussianidad Anderson-Darling. Si un p-valor es ≥0,05 se fija su valor a 4, si por el contrario es <0,05 su valor permanece constante. 4.3.3. Resultados El dise˜no del Algoritmo Gen´etico sufri´o diferentes cambios como resultado del estudio de sus par´ametros. En particular se estudi´o la influencia del par´ametro %Esp, por su relevancia en el factores como: la velocidad de evoluci´on del algoritmo o tiempo de procesado de cada iteraci´on. En las pruebas iniciales se estudiaron valores de %Esp iguales al 15 %, 25 %, 35 % y 45 %. Adem´as, se analiz´o la estabilidad del Algoritmo Gen´etico. Para ello se recomput´o el algoritmo sobre el mismo sujeto, y se compararon las distribuciones obtenidas mediante el an´alisis estad´ıstico K-S, para comprobar la consistencia del mismo. Como resultado se obtuvo las mismas distribuciones gaussianas en las 76 variables de potencias, por tanto se concluy´o que el Algoritmo Gen´etico dise˜nado ad hoc resultaba consistente. Como resultado final se obtuvo una convergencia en los Algoritmos Gen´eticos en 61 sujetos de los 84 que forman la Base Normativa. Los Algoritmos Gen´etico que no han conseguido converger durante la ejecuci´on del algoritmo durante tres meses. Para estos sujetos se podr´ıa aplicar una modificaci´on del algoritmo, como podr´ıa ser la modificaci´on del par´ametro Nmutaci´on que representa la probabilidad de mutaci´on, debido a que aumentando este porcentaje se evita el estado de stagnation, problema que ocurre cuando la soluci´on se estanca en un m´ınimo local. De este modo, se asegura una suficiente exploraci´on del espacio de b´usqueda favoreciendo la convergencia del Algoritmo Gen´etico. Esta modificaci´on queda determinada para un trabajo futuro de este proyecto. 39 5. Clusterizaci´on 5.1. Introducci´on El presente cap´ıtulo propone el desarrollo de t´ecnicas de clusterizaci´on para resolver el problema actual que se enfrentan las Bases Normativas sobre la variabilidad interpersonal del EEG. Dicha variabilidad ha sido relacionado con aspectos fisiol´ogicos como la edad del sujeto, y por tanto con propiedades de maduraci´on del volumen conductor del cerebro. En consecuencia, los estudios actuales se centran en determinar la variabilidad del EEG causada por la edad del sujeto, mediante t´ecnicas de regresi´on por grupos de edades. Los cuales han resultado ´utiles en el caso de Bases Normativas en adultos, al contrario que en el caso de aplicaciones de Bases Normativas en ni˜nos. Este hecho se debe a que en la etapa de ni˜nez la maduraci´on del cerebro var´ıa de forma considerable. A este inconveniente se suma que la variabilidad intrapersonal de EEG no se debe ´unicamente al factor de la edad del sujeto por lo que los par´ametros de regresi´on han de ser modificados mediante la intervenci´on humana. Como alternativa a los m´etodos de regresi´on utilizados actualmente, se propone realizar una estratificaci´on de sujetos de la Base Normativa que tenga en cuenta expl´ıcitamente los cambios significativos en su EEG. La t´ecnica de estratificaci´on requiere de unas herramientas matem´aticas que determinen un criterio de separaci´on. La clusterizaci´on da soporte al uso de estas herramientas matem´aticas en la clasificaci´on de los sujetos de la Base Normativa de EEG, ofreciendo un m´etodo de separaci´on autom´atica. En este cap´ıtulo se detalla el dise˜no de la t´ecnica de clusterizaci´on basada en cambios significativos del EEG para solventar el problema de variabilidad de la se˜nal EEG en Bases Normativas. La secci´on 5.2 introduce el concepto de clusterizaci´on como una t´ecnica de separaci´on de los datos de una distribuci´on, as´ı como los m´etodos estudiados y se justifica la selecci´on del m´etodo jer´arquico como la t´ecnica de clusterizaci´on elegida. Esta secci´on tambi´en recoge estudios relacionados con la medida de la distancia del propio m´etodo de clusterizaci´on, as´ı como el m´etodo para seleccionar el punto de corte del gr´afico dendogrma. 40 5. Clusterizaci´on 5.2 Clusterizaci´on y M´etodos 5.2. Clusterizaci´on y M´etodos La clusterizaci´on permite a partir de una serie de observaciones determinar si existen clases en las que dichas observaciones puedan ser agrupadas. De forma gen´erica la clusterizaci´on es un m´etodo de separaci´on no supervisado en el que no se conoce a priori ni el n´umero de clases ni la pertenencia de las observaciones a ´estos (en el caso de que se conozca el n´umero de clases el problema se simplifica y si adem´as se conocen las etiquetas que caracterizan cada clase el problema es trivial). El m´etodo de clusterizaci´on permite resolver la siguiente situaci´on: dado un conjunto de sujetos (de N elementos) caracterizados por la informaci´on de n variables Xj, (j= 1,2, ..., n), se plantea clasificarlos de manera que los individuos pertenecientes a un grupo (cluster) sean tan similares entre s´ı como sea posible. El problema de clusterizaci´on sobre el desarrollo de una Base Normativa consiste en separar los sujetos de la misma, los cuales est´an formados por N elementos determinados por el n´umero de epochs de EEG de cada uno de ellos. Cada epoch est´a caracterizado por el conjunto de 76 variables de potencias. Por tanto, el problema se presenta mediante una clusterizaci´on de sujetos de acuerdo a 76 variables unidimensionales. Entre las herramientas existentes para desarrollar clusterizaci´on se encuentran: m´etodos basados en particiones, m´etodos jer´arquicos, m´etodos basados en densidades, etc. Entre los m´etodos mencionados se estudiaron aquellos elegidos por su aproximaci´on al presente problema. En primer lugar es estudi´o un m´etodo basado en particiones denominado k-medias, el cual est´a basado en la agrupaci´on de observaciones seg´un la proximidad a unos centros de clusters definidos inicialmente, los cuales var´ıan en el proceso de desarrollo. Esta t´ecnica fue desarrollada en cada variable unidimensional, sin embargo los resultados no fueron favorables, dado que los clusters resultantes no eran consistentes sobre las dimensiones. En un segundo estudio, se desarroll´o el m´etodo de mezclas de gaussianas (en ingl´es Gaussian Mixture Model, GMM). El cual supone que los datos se han generado a partir de una mezcla de k distribuciones multigaussianas. Debido a no poseer variables multigaussianas se realiz´o una prueba de GMM a partir de una variable de potencias (para una banda y canal determinado), la cual es modelada por una gaussiana unidimensional. Como resultado no se obtuvo una soluci´on consistente, es decir, en cada ejecuci´on de la t´ecnica GMM variaba la agrupaci´on de sujetos resultantes. Para solucionar la falta de fiabilidad de los resultados se decidi´o aplicar el m´etodo de remuestreo bootstraping. En la cual se recomputaba la t´ecnica GMM Nveces, y el resultado de la clusterizaci´on vino dada por la grupaci´on con mayor moda de las Nrepeticiones. Esta t´ecnica se repiti´o para el resto de variables de potencias absolutas resultando diferentes clusterizaciones. Para evitar este problema se estudi´o reducir la dimensionalidad de los datos mediante la t´ecnica de s´ıntesis de informaci´on: an´alisis de componentes principales (en ingl´es, Principal Component Analysis PCA). Esta t´ecnica se rechaz´o porque una reducci´on de dimensionalidad produce inevitablemente una p´erdida indeseada de informaci´on. Finalmente, se decidi´o la aplicaci´on del clustering ascendente jer´arquico por su simpli- 41 5. Clusterizaci´on 5.2 Clusterizaci´on y M´etodos cidad y porque no requiere elecci´on del n´umero de clusters a diferencia de los anteriores estudiados. En el clustering ascendente jer´arquico se pretende ir agrupando en cada paso aquellos dos objetos (o conglomerados) m´as cercanos, para de esta forma ir construyendo una estructura conocida como dendograma o ´arbol de clasificaci´on. El dendograma establece una relaci´on ordenada de los grupos previamente definidos y la longitud de sus ramas es una representaci´on de la distancia entre los distintos nodos del mismo (Figura 5.1) El dendograma posibilita la obtenci´on de distintas particiones, simplemente variando el nivel de corte de dicha estructura. El dise˜no del m´etodo de clustering esta unido la elecci´on de una medida de la distancia de similitud entre clusters y a la definici´on de un umbral que permita cortar el dendrograma para determinar los clusters resultantes. Figura 5.1: Ejemplo de un ´arbol de clasificaci´on o dendograma resultante al aplicar un m´etodo de clusterizaci´on jer´arquico. 5.2.1. Distancia de similitud La medida de similitud indica la fuerza de la relaci´on entre dos observaciones. La elecci´on de la medida de similitud proporciona el criterio de agrupaci´on entre sujetos de la Base Normativa. A continuaci´on se detallan los estudios realizados para la selecci´on de la medida de distancia para llevar a cabo el clustering ascendente jer´arquico. El primer estudio se realiz´o en t´erminos de inferencia estad´ıstica, mediante la aplicaci´on de una prueba de hip´otesis. Entre las pruebas de hip´otesis, se barajaron las siguientes pruebas: el Z-test y el contraste para la diferencia de medias de dos poblaciones normales con datos independientes. Para informaci´on m´as detallada consultar el apartado D.2. Sin embargo, ambas m´etricas fueron rechazadas por su car´acter estad´ıstico, cuyos resultados depend´ıan fuertemente de la muestra de poblaci´on utilizada. Un segundo dise˜no contempl´o el uso de la distancia de Bhattacharyya, la cual mide la semejanza de dos distribuciones de probabilidad discretas y es usada generalmente para medir la posibilidad de separaci´on de clases en la clasificaci´on. Este coeficiente es un n´umero comprendido entre cero y uno que expresa la similitud entre dos histogramas. Concretamente es una medida de aproximaci´on de la cantidad de solapamiento entre dos funciones de distribuci´on. Para una descripci´on m´as detallada consultar el anexo D.5. 42 5. Clusterizaci´on 5.3 Resultados Para aplicar esta distancia al problema de clusterizaci´on presente fue necesario definir una distancia global DB, que reuniera la distancia de Bhattacharyya aplicada en cada una de las distribuciones unidimensionales. Dada DBij, la distancia marginal de Bhattacharyya para el conjunto de potencias del canal ien la banda j, la distancia global viene dada por: DB =v u u t 19 X i=1 4 X j=1 DB2 ij (5.1) Se eligi´o esta m´etrica entre otras alternativas como pod´ıan ser la mediana, media, o m´axima de las distancias marginales, por ser robusta ante espurios adem´as de favorecer la contribuci´on de las distancias marginales DBij de mayor valor sobre la distancia global DB. Es decir, las distancias DBij con valores altos contribuir´an en mucha mayor medida que las de valores bajos. 5.3. Resultados A continuaci´on se detallan los resultados obtenidos de aplicar el m´etodo del codo en la selecci´on del punto de corte del dendograma. Adem´as, se describen los resultados obtenidos mediante la clusterizaci´on de los sujetos que conforman la Base Normativa. Se realiz´o un estudio para determinar el punto de corte del dendograma resultante en la aplicaci´on de m´etodo de clusterizaci´on ascendente jer´arquico sobre los sujetos de la Base Normativa, que define las agrupaciones de los sujetos en clusters. El ´unico criterio para determinar el punto de corte sobre el dendograma es la funci´on que representa la distancia m´axima entre clusters frente al n´umero de clusters. ´ Esta es una funci´on mon´otona creciente, sobre la cual se podr´ıa definir una funci´on de coste. Sin embargo, se propone una alternativa que permite calcular el umbral de forma autom´atica. El m´etodo consiste en calcular el codo de dicha funci´on. El codo de una curva se puede definir como el punto con mayor gradiente dentro de una curva mon´otona. Sin embargo, en la pr´actica se trabaja con curvas poco o nada suavizadas, con una gran cantidad de picos y en algunas ocasiones la curva no es mon´otona en todos sus puntos. Determinar el codo de una curva no es un procedimiento trivial, debido a su definici´on ambigua. Partiendo de un m´etodo desarrollado en la literatura [19] se realizaron peque˜nas modificaciones para adaptarlo a funciones crecientes. Este m´etodo tiene en cuenta la tendencia global de la curva y no ´unicamente la diferencia entre valores sucesivos en la curva. Dada la funci´on F(m) el codo se determina siguiendo el siguiente m´etodo: 1. Calcular la diferencia de pendientes consecutivas, teniendo la siguiente funci´on: DiffFun(m)=F(m−1)+F(m+1)−2F(m), donde DiffFun es la funci´on diferencia. 2. Detectar los ncambios locales significantes en DiffFun, seleccionando los nm´aximos locales de la funci´on DiffFun(m). Una vez localizado estos m´aximos, ordenarlos de 43 5. Clusterizaci´on 5.3 Resultados forma descendente en funci´on del valor de DiffFun(m). 3. Calcular el ´angulo que forma la curva F(m) en cada uno de los nm´aximos locales calculados anteriormente. En la figura 5.2 se muestra el c´alculo este ´angulo como la suma de α1yα2siguiendo la expresi´on: Angle(m) = arctan(1/|F(m)−F(m−1)|) + arctan(1/|F(m+ 1) −F(m)|(5.2) 4. El valor del codo viene dado por el primer m´ınimo de la funci´on Angle(m). Figura 5.2: C´alculo del ´angulo que forma la curva F(m) descrita en el punto 3 del m´etodo de c´alculo del codo. La Figura 5.3 muestra el resultado de aplicar el c´alculo del codo sobre la funci´on que describe la distancia entre clusters en funci´on de la distancia entre clusters, resultante en la aplicaci´on del m´etodo de clusterizaci´on ascendente jer´arquico. El codo se sit´ua en el punto [3, 8.107], ´este representa el punto de corte del dendograma en una distancia de BC=8.107 resultando tres clusters sobre la Base Normativa. Definido el punto de corte aBC=8.107 sobre el dendograma la Figura 5.4 representa el dendograma resultante en forma circular. Como soluci´on se han obtenido tres clusters representados por los colores (rojo, azul y verde). En consecuencia mediante la clusterizaci´on jer´arquica ha sido posible realizar una estratificaci´on sobre los sujetos de la Base Normativa desarrollada mediante los cambios significativos en los par´ametros de potencias absolutas de la se˜nal EEG. Como resultado se ha obtenido un total de tres clusters sobre los 61 sujetos que han sido sometidos al m´etodo de clusterizaci´on. 44 5. Clusterizaci´on 5.3 Resultados Figura 5.3: Funci´on distancia entre clusters en funci´on de los clusters resultantes del m´etodo clusterizaci´on ascendente jer´arquica aplicado sobre la Base Normativa. El punto en rojo muestra el codo de dicha funci´on. Figura 5.4: Dendograma circular resultante del m´etodo de clusterizaci´on ascendente jer´arquico aplicado sobre la Base Normativa compuesta por 61 sujetos. El corte de dicho dendograma se ha realizado sobre el punto 8.107 resultante de aplicar el m´etodo del codo. 45