Modelos matem´aticos para el estudio de la activaci´on de

Modelos matemáticos para el estudio de la activación de la corteza cerebral Juan Alberto Escamilla Reyna 1999

Introducci´ on La modelación matemática, en el área de las neurociencias, ha tenido éxito en los distintos niveles de organización del cerebro, baste recordar uno de los éxitos más famosos, el logrado por los neurofisiólogos A. Hodgkin y A. Huxley (H-H), que en el a˜ no 1952 [57], construyeron un modelo matemático para explicar la dinámica del potencial de membrana de una neurona, en función de las corrientes iónicas subyacentes. El modelo se toma todav´ıa como ejemplo para la elaboración de modelos “realistas” de neuronas [21], [99], [100]. El modelo de H-H consta de cuatro ecuaciones diferenciales no lineales –otros modelos constan de más ecuaciones–. El estudio matemático de estos modelos es muy complejo y se ha abordado usando diversas herramientas matemáticas : métodos numéricos, teor´ıa cualitativa de ecuaciones diferenciales, teor´ıa de bifurcaciones, métodos de peque˜ no parámetro, etc. [16], [21], [28], [97], [100]. Otros éxitos famosos, son los logrados en el área de las redes neuronales, desde el trabajo pionero de Mc Culloch-Pitts [80], en el cual se construye una red neuronal, con neuronas artificiales muy simplificadas, pero basada en algunas propiedades de las neuronas reales, hasta los trabajos de Hopfield [59], [60]. Los trabajos de Hopfield han sido muy importantes para el estudio de la memoria asociativa. Actualmente, el área de las redes neuronales es muy fruct´ıfera, se desarrolla fundamentalmente en dos l´ıneas de investigación: la primera, tiene como objetivo construir redes neuronales que resuelvan problemas en diversas áreas, computación, matemáticas, etc., independientemente de que si las propiedades de las neuronas y los parámetros que intervienen en la red tienen alguna relación significativa con las neuronas reales. La segunda, intenta construir redes neuronales “realistas”, que ayuden a explicar fenómenos fisiológicos reales. Para la primera l´ınea de investigación se pueden ver los trabajos [10], [12], [43], para la segunda [24], [40], [41], [75]. En esta área de investigación también se presentan problemas matemáticos complejos que se han abordado con métodos numéricos, ecuaciones diferenciales, teor´ıa de bifurcaciones y técnicas de peque˜ no parámetro. En el área de la electroencefalograf´ıa, se han elaborado modelos considerando al cerebro como un medio conductor de la electricidad, y con una aproximación cuasi-estática, para tratar el problema de localización de I

fuentes, para una distribución de potencial en el cuero cabelludo. Este problema no tiene solución u ńica y usualmente se agregan hipótesis, sugeridas por la anatom´ıa y fisiolog´ıa del cerebro, a las posibles fuentes que originan la distribución de potencial dada sobre el cuero cabelludo. Estos modelos plantean problemas matemáticos que se ubican en el área matemática de los problemas inversos. También se auxilian con los métodos matemáticos de la teor´ıa del potencial. La principal limitante de los modelos de medio conductor cuasi-estáticos es que no pueden informar sobre el origen fisiológico de los patrones temporales r´ıtmicos observados en el electroencefalograma. En cuanto al estudio dinámico de la actividad cortical y su relación con la dinámica de los electroencefalogramas (EEG), se han elaborado diversos modelos. Estos se pueden clasificar, grosso modo, en modelos locales, modelos globales y modelos locales-globales. En los modelos locales se trata de explicar la dinámica de la actividad cortical en función de lo que ocurre a nivel sináptico [33], [34], [75], [76], [77], [95]. En los modelos globales, se considera que la dinámica de la actividad cortical se puede explicar usando propiedades macroscópicas de la propia corteza [87], [88], [89], [90], [92] y en los modelos locales-globales, se mezclan ambos niveles [64], [66], [68], [104], [105]. Entre los modelos globales, se encuentra el elaborado por Nu˜ nez [89], [92]. Es un modelo global lineal de interacción neuronal que describe la propagación de ondas de frecuencias bajas en la corteza cerebral, sus resultados coinciden –cualitativamente– con datos experimentales relacionados con el ritmo alfa y con algunos otros fenómenos electroencelográficos [89], [92]. Sin embargo, su modelo no incluye otros fenómenos fisiológicos, como los correspondientes a los diversos estados del sue˜ no. El modelo relaciona la actividad sináptica –excitatoria e inhibitoria– en cada elemento de corteza x en el tiempo t, con el n´ umero de potenciales de acción que se “generan” en los demás elementos de corteza en tiempos anteriores, teniendo en cuenta las conexiones sinápticas y las velocidades de propagación de los potenciales de acción. El objetivo de esta tesis es elaborar y estudiar modelos que describen la dinámica de la actividad cortical, utilizando variables y parámetros que se puedan medir y as´ı poder correlacionarlos con las mediciones reales de dicha actividad. Los modelos que elaboramos desarrollan las ideas de Nu˜ nez en al II

menos dos direcciones: se introducen nuevos parámetros y se toma en cuenta la no homogeneidad de la corteza cerebral. Además, el enfoque matemático usado –ecuaciones diferenciales parciales ordinarias, técnicas de peque˜ no parámetro, etc.– nos permitió estudiar modelos no lineales de una manera menos complicada. Ver cap´ıtulos 2 y 3 de esta tesis. La tesis fue dividida en tres cap´ıtulos. En el cap´ıtulo 1, se hace un recuento breve sobre la estructura y funcionamiento del cerebro en sus distintos niveles, los diversos métodos y técnicas para obtener información y comentarios sobre los diversos modelos matemáticos que existen a nivel neurona, redes y estructuras completas. En el cap´ıtulo 2, describimos los fundamentos anatómicos y fisiológicos de los modelos matemáticos que construimos. Analizamos dichos modelos y los comparamos con el modelo de Nu˜ nez. En el tercer cap´ıtulo describimos los resultados fundamentales que obtuvimos.

III

Contenido Cap´ıtulo 1 El cerebro : su estructura, funciones y diferentes formas de medir y modelar su actividad 1.1 1.2 1.3 1.4 1.5 1.6 1.7 1.8 1.9

El cerebro : un sistema complejo....................................................... Algunas estructuras importantes del cerebro..................................... La neurona y su fisiolog´ıa.................................................................. Bases f´ısicas del potencial de reposo de la membrana y del potencial de acción............................................................................................ Técnicas para obtener información del cerebro.................................. Algunos modelos de estructuras del cerebro...................................... Redes neuronales................................................................................ Modelos de medio conductor.............................................................. Conclusiones.......................................................................................

1 2 4 8 12 19 22 30 34

Cap´ıtulo 2 Modelaci´ on global de la actividad sin´ aptica de la corteza cerebral 2.1 Introducción......................................................................................... 2.2 Bases fisiológicas para la modelación de la actividad sináptica en la corteza cerebral.................................................................................... 2.3 Descripción y simplificaciones en el modelo global de actividad sináptica en una banda de corteza....................................................... 2.4 Algunos modelos alternativos al modelo de Nu˜ nez.............................. 2.5 Desarrollo y estudio de un modelo de la actividad sináptica de la corteza cerebral.................................................................................... 2.6 Transformación del sistema de ecuaciones integrales en un sistema de ecuaciones diferenciales parciales......................................................... 2.7 La ecuación de activaci´ on para la corteza cerebral.............................. 2.8 Construcción de un modelo compartimentado de activaci´ on de la corteza cerebral....................................................................................

IV

35 36 44 51 56 60 65 68

Cap´ıtulo 3 An´ alisis matem´ atico del modelo de activaci´ on de la corteza cerebral 3.1 Un modelo matemático para una rebanada de corteza cerebral.......... 3.2 Controlabilidad de los estados fisiológicos de la corteza cerebral........ 3.3 Existencia de ondas viajeras en la corteza cerebral............................. 3.4 Ondas estacionarias en la corteza cerebral..........................................

82 89 106 126

Conclusiones .............................................................................................. Bibliograf´ıa ................................................................................................

129 131

V

Cap´ıtulo 1 EL CEREBRO: SU ESTRUCTURA, FUNCIONES Y DIFERENTES FORMAS DE MEDIR Y MODELAR SU ACTIVIDAD 1.1 EL CEREBRO: UN SISTEMA COMPLEJO El cerebro es un sistema extremadamente complejo. Está compuesto de aproximadamente 1011 neuronas densamente interconectadas. Existen alrededor de 1015 conexiones entre las neuronas. Cada neurona se conecta con alrededor de 104 neuronas y recibe conexiones de aproximadamente 104 neuronas. Además las neuronas están rodeadas por un n´ umero igual o mayor de células, llamadas células gliales. Estudiar el diagrama de conexiones del cerebro es una tarea enorme, a pesar de que en muchas áreas del cerebro se ha logrado establecer que las interconexiones entre las neuronas no se hacen al azar, sino que “siguen leyes” bastante definidas. Pero aunque lográramos establecer perfectamente el diagrama de conexiones del cerebro, éste ser´ıa extremadamente dif´ıcil de manejar. Tarea igual de compleja o más es determinar cómo funcionan las distintas estructuras del cerebro y cómo se acoplan entre ellas. El cerebro procesa la información en diferentes escalas espacio-temporales. Se dan interacciones a nivel sináptico, entre grupos de neuronas de diversos tama˜ nos, entre áreas cerebrales, etc. de las cuales surgen fenómenos que no se pueden explicar en 1

una sola escala espacio-temporal y que requieren un tratamiento cuidadoso en varias de estas escalas. 1.2 ALGUNAS ESTRUCTURAS IMPORTANTES DEL CEREBRO La información que se presenta en esta sección, fue obtenida de los libros [17], [71], [106]. Para una información más detallada sobre las distintas estructuras del cerebro y sus funciones recurrir a ellos. Se ha logrado establecer que el cerebro se encuentra dividido en muchas estructuras que interact´ uan entre s´ı. De algunas de ellas se conoce aproximadamente su estructura anatómica y las funciones cerebrales en las que están involucradas. La corteza cerebral, la parte exterior del cerebro y, desde el punto de vista evolutivo, la parte más reciente del cerebro, es un órgano muy complejo de aproximadamente dos cent´ımetros de espesor y una superficie aproximada de 1600 cm2 . En la cavidad craneana, la corteza se encuentra densamente plegada, formando giros y surcos. Se divide en cuatro regiones, los lóbulos frontal, parietal, occipital y temporal. Cada uno de ellos cumple con funciones distintas. Se considera que la corteza se encuentra involucrada en las siguientes actividades : percepción sensorial, actividades motoras voluntarias, almacenamiento de distintos tipos de memoria, actividad conciente, etc. A pesar de su complejidad – gran cantidad de circunvoluciones y aproximadamente 1010 neuronas densamente interconectadas– la corteza tiene un conjunto de caracter´ısticas anatómicas y funcionales que facilitan la elaboración de modelos que nos permiten estudiarla teóricamente. Algunas de estas caracter´ısticas son: 1. Es un órgano estratificado que consta de seis capas horizontales de distintas densidades, además la corteza es bastante uniforme en el plano determinado por cada capa.

2

2. Tiene una organización columnar, es decir las neuronas se agrupan por columnas de diversos tama˜ nos - minicolumnas, columnas y macrocolumnas corticales - que se traslapan y se orientan perpendicularmente a la superficie cortical. 3. La corteza cerebral está mucho más densamente interconectada consigo misma que con otras estructuras del cerebro. La estructura que tiene más conexiones con la corteza es el tálamo y éste contribuye con el 1% de las conexiones que tiene la corteza. En el siguiente cap´ıtulo presentaremos mayor información cualitativa y cuantitativa que justifique los modelos teóricos que estudiaremos. Por ahora terminaremos diciendo que en la corteza cerebral humana se procesa la información en m´ ultiples escalas espacio temporales. Se dan interacciones entre neuronas, entre minicolumnas, entre columnas corticales, entre macrocolumnas, entre áreas corticales e interacciones que involucran todas estas escalas. El tálamo es una estación de relevo y centro integrador de las entradas sensoriales a la corteza. Se cree que también regula el nivel de conciencia y los aspectos emocionales de las sensaciones. El hipotálamo se considera como la interfaz cerebral con los sistemas hormonal y autónomo que controlan la homeóstasis interna del cuerpo. El tálamo y el hipotálamo constituyen conjuntamente el diencéfalo o cerebro intermedio. El sistema l´ımbico ocupa una porción cortical y una porción subcortical. Se encarga de las funciones relacionadas con la emoción, la motivación y ciertos tipos de memoria. Dos subestructuras del sistema l´ımbico son el hipocampo y la am´ıgdala. El hipocampo juega un papel importante en la memoria. La am´ıgdala coordina las acciones de los sistemas autónomo y endocrino. El cerebelo tiene una estructura muy sistematizada y homógenea en todo el órgano. La función más importante del cerebelo es la coordinación de la actividad motora y de la “postura”. Los ganglios basales son cinco n´ ucleos altamente interconectados: caudado, putamen, globo pálido, n´ ucleo subtálamico y sustancia nigra. Los 3

ganglios basales se considera que están implicados en aspectos cognitivos superiores del control motor y en otras muchas funciones. Por u ´ltimo, los hemisferios cerebrales están formados por la corteza y tres estructuras profundas: los ganglios basales, el hipocampo y el n´ ucleo amigdalino. Aunque cada hemisferio tiene funciones especializadas, ambos trabajan en asocación en lo que se refiere a funciones perceptivas, cognitivas y motoras superiores. Existen muchas otras estructuras que no mencionamos en esta sección. Tampoco mencionamos las distintas interrelaciones que existen entre las diversas estructuras. Para mayor información se puede recurrir a los libros mencionados al inicio de esta sección. 1.3 LA NEURONA Y SU FISIOLOGÍA Una exposición detallada de los temas expuestos en esta sección y la siguiente se puede encontrar en [17], [20], [71], [106], [116], [117]. Las neuronas del cerebro presentan una gran variedad de formas y tama˜ nos, (Fig. #1) , pero en todas ellas se pueden identificar cuatro partes (Fig. #2): 1. El cuerpo celular o soma. 2. Las dendritas. 3. El axón. 4. Los terminales axónicas o sinápticas. En el funcionamiento de la neurona, cada una de estas partes tiene un papel espec´ıfico: 1. El cuerpo celular es el centro metabólico de la neurona y el lugar donde se procesa la información. Su diámetro var´ıa de 5 a 100µ.

4

2. Las dendritas son el centro receptor de la información de la neurona. Por lo general, la neurona recibe la información (en su mayor parte de tipo eléctrico) a través de las dendritas, aunque algunas veces la recibe directamente a través del cuerpo celular. 3. El axón es una fibra nerviosa que nace del cuerpo celular y es el encargado de transportar las se˜ nales que genera la neurona, en ocasiones transportando la se˜ nal a distancias considerables. El diámetro de los axones en humano, var´ıa de 0.1µ a 20µ. Frecuentemente los axones están envueltos en una vaina grasa, denominada mielina. La vaina está dividida por espacios regulares sin mielina, llamados nodos de Ranvier. La misión de la vaina de mielina es aumentar la velocidad de propagación de las se˜ nales generadas por la neurona. 4. Los terminales axónicos o sinápticos constituyen los elementos de transmisión de la neurona. A través de ellos una neurona contacta y transmite la información a la zona receptora de otra neurona. La zona de contacto se llama sinápsis y por lo general la zona postsináptica se ubica en las dendritas, aunque también ocurren contactos sinápticos en el cuerpo celular de la neurona y algunas veces en el axón de la neurona. Existen aproximadamente 104 contactos sinápticos por neurona. En condiciones de reposo, toda la neurona (cuerpo celular, soma, dendritas) está polarizada de manera tal que el interior tiene un potencial de -70 mV , considerando que la superficie externa de la membrana es el nivel de referencia, es decir, tiene potencial cero. Este potencial de la neurona es llamado potencial de reposo de la membrana. El potencial de reposo puede variar, dependiendo del tipo de neurona, pero el interior de la membrana de la neurona queda siempre cargado negativamente, con respecto al exterior de la misma. Además este potencial puede cambiar en respuesta a diversos est´ımulos como temperatura, PH, concentraciones extracelulares de iones, corrientes eléctricas, etc. Como ya mencionamos anteriormente la neurona recibe la información (de otras neuronas o información producida por medios artificiales) a través de las dendritas. Esta información se presenta como un cambio en el potencial eléctrico en las dendritas; que llamaremos potenciales postsinápticos. Estos potenciales, tienen una magnitud entre 0.1 − 10 mV , son de naturaleza 5

local, pasivos y disminuyen en intensidad progresivamente y no se detectan más allá de 1 ó 2 mm del sitio de origen. El cuerpo celular los suma espacial y temporalmente y cuando se alcanza un cierto umbral, la neurona “genera” una se˜ nal, llamada potencial de acción. Esta se˜ nal es un impulso eléctrico que se propaga por el axón hasta los terminales axónicos, donde la despolarización produce la liberación de neurotransmisores que se propagan a través del espacio sináptico y afectan a la membrana postsináptica. Esta afectación modifica el estado eléctrico de la membrana y puede ser de tipo excitatorio despolarizándola y as´ı ayudar a la neurona receptora a que “genere” más fácilmente un potencial de acción, o de tipo inhibitorio que tiene el efecto contrario. Estos potenciales postsinápticos serán llamados potenciales excitatorios postsinápticos (PEPS) si excitan a la neurona y potenciales inhibitorios postsinápticos (PIPS) si la inhiben. La se˜ nal que genera una neurona, es decir, el potencial de acción tiene las siguientes propiedades : 1. Se propaga activamente a lo largo del axón. 2. No disminuye su intensidad en función de la distancia. 3. Es de naturaleza todo o nada. 4. Es semejante en todas las neuronas. La amplitud del potencial de acción es de unos 100 mV y dura de 0.5 a 2 milisegundos. Su velocidad de propagación depende del diámetro del axón y de su mielinización. 5. La mayor´ıa de las neuronas generan un potencial de acción cuando el est´ımulo aplicado es despolarizante. La magnitud del est´ımulo debe ser mayor que un cierto valor cr´ıtico, llamado el umbral de disparo de la neurona. Este umbral depende de la neurona y a´ un en la misma neurona puede cambiar. También existen neuronas que generan un potencial de acción cuando el est´ımulo aplicado es hiperpolarizante y de cierta duración. 6. Si el mismo est´ımulo aplicado a la neurona, se mantiene durante un cierto tiempo, el umbral de disparo de la neurona disminuye.

6

7. Una vez que se genera un potencial de acción, la neurona es incapaz, durante cierto tiempo, de generar un nuevo potencial de acción, independientemente de la magnitud del est´ımulo que se le aplique, en este caso se dice que la neurona está en estado refractario. Al tiempo que dura la neurona en estado refractario se le llama per´ıodo refractario absoluto. Existe un intervalo de tiempo, llamado per´ıodo refractario relativo, durante el cual el umbral de disparo de la neurona aumenta y con esto dificultar la generación de nuevos potenciales de acción. El per´ıodo refractario absoluto de una neurona determina la frecuencia de generación de potenciales de acción. 8. Un potencial de acción se puede generar mediante est´ımulos de corriente que var´ıen linealmente con el tiempo, con tal que la pendiente sea mayor que un cierto valor cr´ıtico. Si la pendiente es menor que este valor cr´ıtico, no se genera ning´ un potencial de acción, sin importar qué tan grande sea la despolarización de la membrana. A esta propiedad de la membrana se le llama acomodación. El mecanismo de propagación del potencial de acción se puede describir esquemáticamente de la siguiente manera: el est´ımulo inicial en un punto del axón (una despolarización de la neurona de cierta magnitud) provoca corrientes locales que fluyen de manera pasiva a través de la membrana, causando una variación en el potencial de reposo en zonas cercanas al punto estimulado. Esta despolarización en las zonas cercanas “genera” un mecanismo en la membrana que genera un voltaje muchas veces mayor que el inicial, el cual a su vez produce corrientes locales que causan un cambio en el potencial en las zonas adyacentes y as´ı sucesivamente. Está claramente demostrado que el potencial de acción no se puede propagar por el axón como si éste fuera un cable (en forma pasiva) ya que en este caso el potencial decrecer´ıa en forma exponencial y la rapidez de decrecimiento estar´ıa en función de la resistencia de la membrana axónica y de la resistencia a las corrientes que fluyen longitudinalmente a través del axón

7

1.4 BASES FÍSICAS DEL POTENCIAL DE REPOSO DE LA ´ MEMBRANA Y DEL POTENCIAL DE ACCION Toda neurona contiene en su interior varios tipos de iones, entre los que se encuentra el sodio (N a+ ), el potasio (K + ) y el cloro (Cl− ) en distintas concentraciones, además de algunas otras moléculas como aminoácidos, proteinas, etc. La membrana de una neurona tiene distintas permeabilidades con respecto a cada tipo de ion. La difusión de los iones a través de la membrana celular se realiza a través de zonas especiales restringidas de la membrana, llamados “poros” o “canales” y existe selectividad de estos canales para los distintos tipos de iones. La membrana tiene también mecanismos activos para “bombear” algunos de los tipos de iones (principalmente el potasio y el sodio) del interior al exterior de la membrana y rec´ıprocamente. En el exterior de la neurona también se encuentran distintos tipos de iones, principalmente sodio, potasio y cloro, en distintas concentraciones. Para el caso de la membrana neuronal en reposo del axón gigante del calamar, (ver Tabla 1). Los tipos de iones y las concentraciones de éstos, pueden variar, dependiendo del tipo de neurona que se considere. Tabla 1 L´ıquido Potencial de Especie Citoplasma extracelular equilibrio∗ iónica (mM) (mM) (mV) + K 400 20 −75 N a+ 50 440 +55 Cl− 52 560 −60 A− (iones 385 – – orgánicos) La distribuci´ on de los principales iones a trav´ es de la membrana neuronal en reposo del ax´ on gigante del calamar. ∗

Se denomina as´ı al potencial de membrana en el que no hay flujo neto de iones a través de la membrana

8

Actualmente se considera que el potencial de reposo de la membrana se logra cuando los flujos pasivos y activos de N a+ , K + , Cl− están balanceados (Fig. #3), lo cual se logra cuando las fuerzas de difusión, provocadas por las distintas concentraciones de los iones, en el interior y exterior de la membrana se equilibran con las fuerzas eléctricas, producidas por el hecho de que los iones son part´ıculas cargadas eléctricamente. El potencial de acción se logra cuando la neurona se despolariza hasta un nivel cr´ıtico; en ese momento tiene lugar un cambio s´ ubito, donde las barreras normales que se oponen al flujo de sodio y potasio desaparecen momentáaneamente, estableciéndose un flujo local de iones, de suficiente magnitud como para invertir el potencial de membrana, que alcanza unos 50 mV (positivo en el interior) y después se invierte de nuevo para establecer el potencial de reposo. Pero la primera inversión (hacia el interior positivo) ha producido una se˜ nal potente que se extiende y lleva la región adyacente de la membrana a un nivel cr´ıtico, repitiéndose entonces el proceso (Fig. #3). En la siguiente sección explicaremos más detalladamente esto y describiremos los diversos modelos que existen al respecto. Una ecuación importante en la fisiolog´ıa es la ecuación de Nernst, obtenida en 1888 por W. Nernst. Esta ecuación relaciona el potencial de una membrana permeable a un ion con las concentraciones del ion en el interior y el exterior de la membrana. La ecuación se aplica con la condición de que la corriente total iónica sea igual a cero, es decir, cuando la corriente iónica producida por el potencial es igual a la corriente producida por la diferencia de concentraciones en el interior y el exterior de la neurona. Supongamos que tenemos dos compartimientos A y B, donde se tiene una solución que contiene un ion I, a diferentes concentraciones CA y CB y que los compartimientos están separados por una membrana permeable al ion I (Fig. #4)

Figura 4 A

ion I CA +

ion I −CB V 9

B

Como el ion I se encuentra a diferentes concentraciones en A y B, entonces se empezará a difundir del compartimiento de mayor concentración al de menor concentración pero como el ion lleva una carga, el compartimiento hacia el cual se está difundiendo empieza a aumentar su carga as´ı aumenta su voltaje frenando la difusión del ion. El equilibrio dinámico se logrará cuando el voltaje en el compartimiento que recibe los iones difundidos se equilibra con la fuerza de difusión que surge debido a las distintas concentraciones del ion. Es decir, cuando: µ

¶

RT CB (1.1) V = ln ZF CA donde V es el voltaje entre el compartimiento A y el compartimiento B, R es la constante de los gases ideales, T la temperatura (en grados Kelvin), Z la valencia del ion y F la constante de Faraday (la carga de una mole de part´ıculas univalentes), aproximadamente 9.65 × 104 C/mol. Esta ecuación nos permite calcular el voltaje entre A y B en términos de las concentraciones del ion en los compartimientos A y B, cuando el flujo neto de iones entre los compartimientos es cero. Este voltaje entre A y B es llamado el potencial de equilibrio del ion. En la Tabla 1 están los valores de los potenciales de equilibrio de algunos iones. En el caso que los compartimientos contengan más de un ion, con distintas concentraciones y que la membrana contenga canales independientes para el paso de los iones, se puede calcular el potencial de equilibrio de un solo ion, usando (1.1), con la condición de que la membrana tenga abierto, u ńicamente, el canal de ese ion. Si todos los canales de la membrana están abiertos simultáneamente, entonces el voltaje dependerá de las permeabilidades (facilidad con que un ion atraviesa la membrana (cm/seg)) relativas de los iones y de las concentraciones de los iones. Goldman (1943) calculó el potencial final promedio: ÃP

P

RT P + [I + ]A + Pi− [Ii− ]B ln P i+ +i V = P Zi F Pi [Ii ]B + Pi− [Ci− ]A

!

(1.2)

donde Pi+ , Pi− representan las permeabilidades de las especies iónicas positivas y negativas respectivamente y [Ii+ ]∗ , [Ii− ]∗ son las concentraciones de estos iones en el compartimiento ∗.

10

En el caso de la membrana de una neurona sabemos que en el interior y el exterior de la misma se encuentran varios tipos de iones con distintas concentraciones y permeabilidades. Entre los principales se encuentran el sodio (N a+ ), el potasio (K + ) y el cloro (Cl− ). Aplicando (1.2) para estos tres iones obtenemos: "

KT Pk+ [K + ]e + PN a+ [N a+ ]e + PCla− [Cl− ]i V = ln F Pk+ [K + ]i + PN a+ [N a+ ]i + PCla− [Cl− ]e

#

(1.3)

En la ecuación de Goldman se supone que el campo eléctrico que existe a través de la membrana var´ıa uniformemente, por lo que también se supone que los movimientos de los iones a través de la membrana son independientes unos de otros. La importancia de este resultado es que nos permite hacer predicciones cuantitativas del potencial de membrana cuando se tienen variaciones de las concentraciones de los iones y de las permeabilidades de la membrana a los distintos iones. Usando la ecuación de Nernst y con ciertos experimentos de tipo eléctrico se ha logrado demostrar que la membrana de una neurona en estado de reposo no es permeable u ńicamente al potasio sino que también es permeable al sodio, aunque la permeabilidad del sodio es muy peque˜ na con respecto a la del potasio. Se considera que la permeabilidad con respecto a los iones de sodio es aproximadamente el 1% de la correspondiente a los iones de potasio y por ello debe de haber una “bomba” en la membrana que bombea iones de sodio hacia el exterior de la membrana e iones de potasio al interior de la misma. También se logró establecer lo que ocurre cuando una neurona genera un potencial de acción. El potencial de acción se genera cuando se logra una despolarización de la membrana de aproximadamente 15 mV y en ese momento ocurre que la permeabilidad de la membrana al sodio aumenta, por la apertura de canales de N a+ (que dependen del voltaje), lo que hace que el potencial se eleve. Esto a su vez permite que aumente la permeabilidad de la membrana, por un nuevo aumento de los “canales” de sodio, lo que hace que el potencial se eleve y as´ı sucesivamente, hasta que el potencial casi alcanza el potencial de equilibrio del sodio (se establece en la membrana un 11

mecanismo de retroalimentación positiva entre el potencial de membrana y la permeabilidad de ésta al sodio). Este potencial de equilibrio no se alcanza ya que cuando se inicia el aumento de la permeabilidad de la membrana al sodio, un corto tiempo después empieza a aumentar la permeabilidad de la membrana al potasio y una entrada de cloro a la membrana. Además los canales de sodio empiezan un proceso de inexcitabilidad y también se tiene un flujo pasivo de iones de potasio. En resumen, el potencial de acción se origina cuando la estimulación eléctrica alcanza un cierto umbral. Su dinámica se encuentra determinada por la relación que existe entre las distintas corrientes iónicas que se establecen en la membrana, tanto corrientes pasivas como activas. En la sección 1.6 se mostrarán algunos modelos que describen este proceso de manera cuantitativa. Una manera equivalente de estudiar el potencial de membrana de una neurona consiste en construir un circuito eléctrico que incorpore las tres propiedades básicas que utiliza la membrana para generar se˜ nales eléctricas - los canales iónicos, los gradientes de concentración de los iones, la capacidad de la membrana para almacenar carga eléctrica. Cada canal iónico se puede representar como una resistencia (o como una conductancia que es la inversa de la resistencia). La conductancia para un canal iónico ser´ıa en cierta medida, equivalente a la permeabilidad de la membrana a ese ion. La capacidad de la membrana para almacenar carga se puede representar como una capacitancia y los gradientes de concentración como una bater´ıa (Figs. #5,6) El modelo de circuito equivalente es muy importante desde el punto de vista experimental ya que permite el uso de técnicas electrofisiológicas para el estudio de las se˜ nales neurales, en lugar del uso de isótopos radiactivos que se usan para medir los cambios de permeabilidad y difusión de iones. ´ ´ DEL 1.5 TECNICAS PARA OBTENER INFORMACION CEREBRO En el proceso de obtener información sobre la estructura anatómica y funcional del cerebro en sus diversos niveles –neurona, poblaciones de neu12

ronas y todo el cerebro– se han elaborado, con la ayuda de la tecnolog´ıa y otras ciencias, un conjunto de técnicas e instrumentos que se han ido mejorando con el transcurso del tiempo y que han sido de gran utilidad para el diagnóstico cl´ınico, para los neuroanatomistas y los neurofisiólogos. Para el estudio anatómico se han usado instrumentos como el microscopio óptico, el microscopio electrónico, etc. y algunos métodos de tinción selectiva como el método de Golgi, Nissl, Nauta, etc. El método de Golgi, inventado por Camilo Golgi alrededor de 1875, es un método de tinción selectiva que tiene la propiedad de que al te˜ nir un grupo de neuronas, sólo algunas de ellas, aproximadamente una de cada cien se ti˜ nen por completo (dendritas, cuerpo celular, axón) y las demás quedan intactas. Con este método, Santiago Ramón y Cajal logró establecer que el cerebro está constituido por células separadas, que las interconexiones entre ellas no eran al azar, además de describir la arquitectura de una gran cantidad de estructuras del cerebro. Con el método de Nissl se logra destacar los cuerpos celulares de las neuronas, pero no sus dendritas y axones. Por mucho tiempo los neuroanatomistas usaron estos dos métodos de tinción para estudiar anatómicamente el cerebro. Alrededor de 1952 W. Nauta inventó un método de tinción selectiva, que actualmente lleva su nombre. Seg´ un D. H. Hubel, este método es el primer instrumento eficaz para cartografiar las conexiones entre las distintas estructuras del cerebro. El método de Nauta aprovecha el hecho de que, cuando se destruye una neurona, la fibra nerviosa que procede de ella degenera y, antes de desaparecer por completo, puede ser te˜ nida de manera diferente de sus vecinos normales. Los avances tecnológicos y una mayor comprensión del funcionamiento de las neuronas han permitido crear técnicas con las cuales, no sólo se logra obtener información más precisa sobre las conexiones entre regiones cerebrales, identificar neuronas, sino también sobre las zonas activas del cerebro. Existen técnicas que aprovechan el hecho de que las neuronas absorben y sintetizan varias sustancias – prote´ınas, enzimas, etc. –. Algunas de éstas sustancias son transportadas, a través del axón, desde el soma hasta los terminales axónicos (transporte anterógrado). También existe un transporte de sustancias desde los terminales axónicos al cuerpo celular (transporte retrógrado).

13

La autorradiograf´ıa de Cowan (1972) es una técnica que aprovecha el transporte axónico anterógrado para determinar los terminales axónicos de una neurona. Esta técnica consiste en inyectar, en alguna estructura cerebral, aminoácidos marcados radiactivamente, éstos penetran en el soma de la neurona y se incorporan en las sustancias que son transportadas hasta los terminales axónicos, donde se les puede detectar por un proceso de autorradiograf´ıa. También se cuenta con técnicas que aprovechan el transporte axónico retrógrado. Entre las más conocidas esta la técnica de la peroxidasa de rábano, se inyecta en el sistema nervioso, la enzima peroxidasa de rábano, ésta es absorbida por los terminales axónicos y transportada retrógradamente hacia el cuerpo celular, donde puede ser detectada por diversas técnicas, y as´ı conocer el origen del axón. Para mayor información de otras sustancias que se usan para el estudio de las conexiones entre regiones cerebrales, otros métodos para detectar sustancias (la inmunohistoqu´ımica se basa en la detección de sustancias tisulares por medio de anticuerpos), la hibridación in situ, etc. recurrir a [51] y [107]. En este u ´ltimo se hace una revisión detallada de la organización neuronal y sináptica de la corteza cerebral. Basándose en el hecho de que la glucosa es el combustible de las neuronas, y éstas consumen más glucosa cuando están activas que cuando se encuentran en reposo, Louis Sokoloff (1977) inventó una técnica que le permitió medir el metabolismo de las neuronas. La técnica consiste en inyectar en el torrente sangu´ıneo un análogo qu´ımico de la glucosa (la 2–desoxiglucosa), que es absorbida por las neuronas, mediante el mismo sistema que la glucosa, aunque una vez en el interior de la neurona, ya no puede ser metabolizada. Si la 2–desoxiglucosa es marcada radiactivamente e inyectada en el torrente sangu´ıneo, la velocidad con que se acumula en el interior de las neuronas, nos proporciona un ´ındice de la actividad metabólica de las neuronas. Las técnicas que hemos mencionado hasta ahora requieren en alg´ un momento de la destrucción del sujeto de estudio. Con el desarrollo de la tecnolog´ıa se crearon nuevas técnicas menos invasivas, que permit´ıan obtener información sobre la actividad del cerebro humano en vivo, por ejemplo la to14

mograf´ıa axial computarizada (TAC), la tomograf´ıa por emisión de positrones (TEP) y la formación de imagenes por resonancia magnética (IRM), entre otras. Las herramientas teóricas fundamentales para el análisis de los diferentes tipos de tomograf´ıas están vinculadas a la transformada de Radon y los métodos matemáticos y computacionales de reconstrucción por proyecciones. Existe una vasta bibliograf´ıa de estos temas entre los que vale citar a [53], [54] y [85]. La tomograf´ıa axial computarizada se aprovecha del hecho que los distintos tejidos absorben diferentes cantidades de energ´ıa de rayos X; mientras más denso es el tejido, tanto más absorbe. Por ello un haz de rayos X muy concentrado que atraviese el cuerpo emergerá del mismo a un nivel reducido de intensidad que dependerá de los tejidos y órganos encontrados a su paso. Si se dirigen a través del cuerpo haces planos de rayos X con diferentes ángulos de inclinación, se recogerá información suficiente para reconstruir una imagen de la sección del cuerpo, la cual debidamente procesada permite reconstruir una imagen tridimensional de la intensidad de radiación absorbida y por ende de la estructura interna del objeto estudiado. La tomograf´ıa por emisión de positrones, permite cartografiar las estructuras cerebrales activas ya que permite la medición de diversas actividades del cerebro, tales como el metabolismo de la glucosa, el consumo de ox´ıgeno, el riego sangu´ıneo y la interacción con los fármacos. De todas estas variables, el flujo sangu´ıneo ha demostrado ser el indicador más fiable de la función cerebral en cada momento. Con la TEP se pueden localizar cambios de actividad con precisión milimétrica. La formación de imágenes por resonancia magnética (IRM) es un instrumento bastante com´ un para diagnosticar lesiones en tejidos. Tuvo su origen en una técnica de laboratorio denominada resonancia magnética nuclear (RMN), dise˜ nada con el fin de explorar las caracter´ısticas qu´ımicas de las moléculas. La (RMN) se basa en que muchos átomos se comportan como diminutas agujas imantadas en presencia de un campo magnético. En estas condiciones, la aplicación de impulsos de una onda radioeléctrica a una muestra, perturba los átomos de una manera muy precisa, y éstos en consecuencia emiten se˜ nales de radio detectables que identifican el n´ umero y el estado de los átomos pertenecientes a la muestra. El ajuste del campo magnético y de los impulsos de onda radioeléctrica proporciona información 15

sobre la muestra de estudio. En el caso del cerebro P. C. Lauterbur descubrió que pod´ıan formarse imagenes por RMN mediante la detección de protones. La utilidad de estas part´ıculas, abundantes en el cuerpo humano, reside en su gran sensibilidad a los campos magnéticos que los hace comportarse como min´ usculas agujas imantadas. De esta forma se obtuvieron excelentes imagenes anatómicas con detalle muy superior a las conseguidas mediante TAC. El interés por aplicar la IRM a la obtención de imágenes cerebrales se debe a que esta técnica es capaz de detectar una se˜ nal inaccesible a las exploraciones mediante TEP como por ejemplo, el aumento en la concentración de ox´ıgeno que se produce en una zona donde se ha intensificado la actividad neuronal. Otras ventajas sobre la TAC y otras técnicas son: 1. La se˜ nal procede directamente de cambios inducidos por causas funcionales en el tejido cerebral (el cambio en la concentración de ox´ıgeno venenoso). 2. Para obtener una se˜ nal no hay que inyectar ninguna substancia, radiactiva o de otro tipo. 3. La IRM proporciona información anatómica y funcional sobre cada sujeto, lo cual permite una identificación estructural exacta de las zonas activas. 4. Adecuadamente equipada, puede observar en tiempo real el ritmo de cambio de la se˜ nal de ox´ıgeno inducida en el flujo sangu´ıneo. Sin embargo a pesar de las ventajas que proporciona la combinación de estas técnicas, ninguna de ellas puede determinar la evolución temporal del intercambio de información en diferentes zonas del cerebro, lo cual es fundamental para estudiar cómo se coordinan en red diferentes zonas del cerebro para la realización de una actividad dada. El obstáculo que enfrentan la IRM y TEP para este objetivo es la gran diferencia de velocidades entre la actividad neuronal y el ritmo de variación de los niveles de oxigenación en las neuronas. Las se˜ nales viajan de una parte a otra del cerebro en 0.01 segundos o menos, pero los cambios de flujo sangu´ıneo y oxigenación de la sangre son mucho más lentos y ocurren desde décimas de segundo hasta varios segundos más tarde. La IRM no puede seguir las “conversaciones” entre regiones 16

cerebrales. Los u ńicos métodos que responden con rápidez suficiente son las técnicas de registro eléctrico, es decir, la electroencefalograf´ıa (EEG) que detecta la actividad eléctrica del cerebro a partir de mediciones efectuadas sobre el cuero cabelludo y la magnetoencefalograf´ıa (MEG), la cual mide el comportamiento electromagnético que genera la actividad eléctrica en el interior del cerebro. El defecto de estas técnicas es que su resolución espacial es menor que la de otras técnicas y la localización precisa de las fuentes de actividad cerebral sigue siendo un problema muy complejo. En 1929 el psiquiatra alemán Hans Berg demostró que era posible medir –usando electrodos colocados en el cuero cabelludo– la actividad eléctrica del cerebro humano. Demostró que bajo ciertas condiciones en el sujeto de estudio (ojos cerrados, sin est´ımulos exteriores, relajado) el registro de la actividad eléctrica segu´ıa un patrón temporal r´ıtmico al que llamó ritmo α. Este patrón temporal ten´ıa una frecuencia entre 9 y 10 c/s. Los or´ıgenes fisiológicos del ritmo α, siguen siendo motivo de controversia, algunos fisiólogos consideran que el ritmo α se origina en el tálamo y que el tálamo le impone este ritmo a la corteza, a través de las conexiones que existen entre ambos, mientras que otros consideran que este ritmo es una manifestación de la actividad masiva de la corteza cerebral. Posteriormente se han descubierto otros patrones temporales r´ıtmicos que se producen bajo ciertas condiciones corporales especiales. Por ejemplo, el ritmo δ (0 a 5 Hz) se produce bajo sue˜ no profundo, el ritmo β (15 a 30 Hz) se produce bajo intensa actividad mental y otros más. El EEG (medido en el cuero cabelludo o en otras regiones del cerebro) se ha convertido en una herramienta cl´ınica importante para el estudio de ciertas patolog´ıas del cerebro, como la epilepsia, tumores cerebrales, etc. Estas patolog´ıas aparecen en el EEG como patrones anormales de la actividad eléctrica espontánea del cerebro. Los potenciales evocados (los potenciales que se producen como repuesta a un est´ımulo) se han usado en gran medida para el diagnóstico cl´ınico de las afecciones sensoriales (auditivos, visuales, etc.) y en la detección de problemas cognitivos [14], [92], [105]. Uno de los problemas fundamentales de la electroencefalograf´ıa consiste en correlacionar los patrones espacio-temporales obtenidos mediante el EEG con la actividad eléctrica del cerebro y con conductas observables. Si se quiere 17

elaborar modelos que describan la dinámica de la actividad eléctrica del cerebro y que se correlacionen con los patrones espacio-temporales obtenidos con el EEG; se requiere tomar en cuenta las escalas espacio-temporales en que están promediando las mediciones efectuadas ya que los datos electroencefalográficos son obtenidos promediando la actividad neuronal en distintas escalas espacio-temporales. Existen EEG’s con distintas resoluciones tanto espaciales como temporales. Por ello las variables a considerar en nuestros modelos deben estar en congruencia con las escalas espacio-temporales usadas en el EEG. En el siguiente cap´ıtulo detallaremos los problemas que se presentan en esta área de estudio. Los modelos que más han sido estudiados son los que utilizan los principios de la teor´ıa del campo electromagnético considerando al cerebro como un medio continuo conductor, con una cierta geometr´ıa (que depende de las hipótesis del modelo) y cierta conductividad variable. En el problema directo, se considera dada la ubicación de las fuentes de corrientes impresas (zonas activas), su configuración, la geometr´ıa de la cabeza y las conductividades de los distintos comportamientos (cuero cabelludo, cráneo, corteza, l´ıquido extracraneal, etc.) y se calcula la distribución del potencial generado en estas condiciones, comparándose los resultados teóricos con los obtenidos con técnicas electroencefalográficas. Mediante el aporte de ambos datos se puede encontrar alguna información sobre la respuesta electrica de una zona activa dada. El problema inverso consiste en obtener información sobre la ubicación de las fuentes de corrientes impresas y las conductividades en diferentes partes del cerebro usando los datos discretos que se obtienen con el EEG. Este problema no tiene solución u ńica, es decir, existe una multitud de distribuciones de fuentes de corriente en el cerebro que pueden producir la misma distribución de potencial en el cuero cabelludo. Al no tener solución u ńica el problema inverso, usualmente se agregan restricciones, sugeridas por la anatom´ıa y fisiolog´ıa del cerebro, a las posibles fuentes que originan la distribución de potencial dada sobre el cuero cabelludo. Usualmente estos modelos se restringen a la aproximación electrostática por lo que son u ´tiles para estudiar la distribución espacial de generadores de los EEG pero no as´ı sus propiedades temporales [9], [105]. 18

Para estudiar la dinámica espacio-temporal de la actividad eléctrica se pueden construir ecuaciones plausibles, a partir de datos anatómicos y fisiológicos que involucren variables que sean congruentes con las escalas espacio-temporales usadas en los distintos EEG, para estudiar sus soluciones y su relación con las mediciones obtenidas a través del EEG. Finalmente, otro método usado consiste en elaborar una jerarqu´ıa de escalas espacio-temporales. En cada escala espacio-temporal se deducen las ecuaciones que gobiernan las variables que se consideran idóneas en ese nivel y se estudia la relación que existe entre un nivel y el nivel superior inmediato. En la sección siguiente desarrollaremos con mayor detalle algunos de estos modelos. 1.6 ALGUNOS MODELOS DE ESTRUCTURAS DEL CEREBRO En el estudio del cerebro se han elaborado diversos modelos para los distintos niveles del mismo. Se cuenta con modelos a nivel de neurona, grupos de neuronas, de subestructuras del cerebro y del cerebro completo. Cada uno de estos modelos toma en cuenta ciertos aspectos del cerebro dependiendo de los problemas que se quieren estudiar. Es importante recalcar que los elementos y las relaciones que se toman en cuenta para la elaboración del modelo están intimamente relacionados con el nivel del cerebro que se esté considerando y con los fenómenos que se quieren explicar. Uno de los primeros éxitos importantes en la modelación matemática en el estudio del cerebro fue el logrado por A. Hodking y A. Huxley (H.H) en 1952 [57]. El modelo por ellos propuesto describe, por medio de un sistema de ecuaciones diferenciales no lineales, cómo var´ıa temporalmente el potencial de membrana en función de las corrientes iónicas, que fluyen a través de la membrana en el axón de un calamar. El modelo toma en cuenta tres corrientes iónicas, la de potasio (IK ), la de sodio (IN a ) y una corriente de filtraje (IL ). Considera que las corrientes de sodio y potasio fluyen a través de canales cuya conductancia depende del potencial de membrana y del tiempo.

19

De acuerdo con las leyes de Kirchhoff, se tiene que: Cm

dV = IN a + IK + IL , dt

(1.4)

donde Cm es la capacitancia de la membrana. También se tiene: IN a = gN a (V − VN a ) ,

(1.5)

IK = gK (V − VK ) ,

(1.6)

IL = gL (V − VL ) ,

(1.7)

donde g∗ es la conductancia del ion ∗ y V∗ el potencial de equilibrio de dicho ion. Luego se tiene que Cm

dV = gN a (V − VN a ) + gK (V − VK ) + gL (V − VL ) . dt

(1.8)

En este modelo se considera que gN a y gK son dependientes del potencial V y gL es constante, es decir: Cm

dV = m3 h¯ gN a (V − VN a ) + n4 g¯K (V − VK ) + gL (V − VL ) , dt

(1.9)

donde g¯N a , g¯K son las conductancias máximas del sodio y el potasio y 0 < m, n, h < 1 donde m y n son variables de activación y h una variable de inactivación que obedecen las ecuaciones: dm = αm (1 − m) − βm m , dt

(1.10)

dn = αn (1 − n) − βn n , dt

(1.11)

dh = αh (1 − h) − βh h , dt

(1.12)

con coeficientes αm , βm , αn , βn , αh , βh que son funciones obtenidas empir´ıcamente y que dependen del potencial V . 20

Hodking y Huxley describen un gran n´ umero de experimentos controlados que involucran corrientes, conductancias y potenciales, cuyos resultados “coinciden” con las predicciones del modelo. El modelo (1.9), (1.10), (1.11) y (1.12) se construyó bajo la hipótesis de que el potencial de membrana sea uniforme a lo largo del axón, o se considere u ńicamente un pedazo peque˜ no de membrana. Para estudiar situaciones en las que el potencial de membrana var´ıa temporal y espacialmente a lo largo del axón, este modelo no se puede aplicar. Usando resultados de la teor´ıa del cable se puede construir un sistema de ecuaciones diferenciales parciales no lineales que tomen en cuenta la dependencia espacial y generalicen el sistema (1.9), (1.10), (1.11), (1.12). En este caso el sistema de ecuaciones que se obtiene es 1 ∂ 2V ∂V = C + g¯N a m3 h(V − VN a )2πr R ∂x2 ∂t +{¯ gk n4 (V − Vk ) + g¯e (V − Ve )}2πr ,

(1.13)

∂m = αm (1 − m) − βm m , ∂t

(1.14)

∂h = αh (1 − h) − βh h , ∂t

(1.15)

∂n = αn (1 − n) − βn n , (1.16) ∂t donde r es el radio del axón, R es la resistencia por unidad de longitud en el interior del axón y C la capacitancia de la membrana por unidad de longitud. Para un resumen de los experimentos realizados por Hodking y Huxley y de los resultados obtenidos consultar [21]. También se puede consultar [40]. El modelo es importante ya que describe la evolución temporal del potencial de membrana en función de las corrientes iónicas subyacentes. Sin embargo, no toma en cuenta los eventos moleculares y qu´ımicos subyacentes al proceso, y por lo tanto no es posible predecir cómo influyen en el proceso los cambios en el nivel qu´ımico y molecular. Además no es un modelo universal, es decir, se aplica u ńicamente al potencial de membrana de una neurona del calamar gigante. 21

El modelo de H-H incorpora u ńicamente tres corrientes iónicas, sodio, potasio y una corriente de filtración, y debido a esto, el modelo consiste de un sistema de 4 ecuaciones diferenciales ordinarias. Existen otros modelos donde se incorporan más corrientes iónicas, por ejemplo, Bucholtz, [16] elabora un modelo de la L. P. neurona en el ganglio estomatogástrico del cangrejo. Este modelo consiste de un sistema de 14 ecuaciones diferenciales ordinarias. Otros neurofisiólogos elaboran modelos de neuronas que siguen esencialmente los principios del modelo de H-H, los cuales incorporan más corrientes iónicas o toman en cuenta las variaciones de las concentraciones de ciertos iones en el interior de la membrana. En [21] se presentan otros modelos de neuronas, también se pueden ver otros modelos en [99], [100]. El estudio matemático de estos modelos es bastante complejo debido al n´ umero de ecuaciones diferenciales que aparecen. Por lo general son “intratables” anal´ıticamente y a veces sólo se pueden encontrar soluciones por integración n´ umerica. Se han realizado algunos estudios cualitativos de algunos de estos modelos [58], [111]. Una técnica que se usa frecuentemente es reducir, mediante alg´ un mecanismo, el n´ umero de dimensiones del modelo manteniendo las principales caracter´ısticas del mismo. Por ejemplo, Fitzhugh [28], [29], propone un modelo bidimensional para el modelo de H-H Abbot, Kepler [2] proponen un método de reducción más sistemático que el de Fitzhugh, ya que reducen el modelo de H-H a dos dimensiones, y Golomb, Guckenheimer, Gueron [39], aplican este método para reducir el modelo de Buchholtz a 7 dimensiones. Doya, Selverston [26] elaboran un método de reducción de variables usando redes neuronales artificiales. 1.7 REDES NEURONALES El estudio de la estructura y funcionamiento de una neurona es un paso importante en el estudio del cerebro, pero el cerebro está compuesto de aproximadamente 1011 neuronas densamente interconectadas, organizadas en grupos de diversos tama˜ nos. Las neuronas de un mismo grupo pueden interactuar entre s´ı y también pueden interactuar los grupos entre s´ı y frecuentemente se dan estos tipos de interacción de manera simultánea. La estructura 22

de estos grupos y sus mecanismos de interacción no son muy conocidos. Los estudios experimentales de estos grupos de neuronas son dif´ıciles de realizar. Los modelos de redes neuronales (un conjunto de unidades elementales interconectadas) son una herramienta importante para resolver muchos problemas en el área de la ingenier´ıa de la computación, matemáticas, neurofisiolog´ıa, etc. En cada una de estas áreas se interesan en diferentes aspectos de las redes neuronales, dependiendo de los problemas que se quieran resolver. A los neurofisiólogos les interesa estudiar redes neuronales con la finalidad de estudiar el cerebro humano. Ellos desean construir redes neuronales “realistas” que imiten ciertas capacidades del cerebro como las sensoriales, las de aprendizaje y las relacionadas con la memoria. También desean construir redes que describan los mecanismos subyacentes a los patrones r´ıtmicos que se detectan en el cerebro. Para los ingenieros, computólogos y matemáticos, interesa construir redes neuronales con fines computacionales. A estas u ´ltimas se les conoce como redes neuronales artificiales y frecuentemente tienen poca relación con lo que ocurre en el cerebro. Las redes neuronales se construyen siguiendo ciertas propiedades que tiene el cerebro. Están constituidas por un conjunto de n unidades elementales llamadas neuronas artificiales, unidades de procesamiento o nodos que están conectados entre s´ı. Las propiedades de las redes dependen del n´ umero de neuronas de que consta la red y de la forma en que se interconectan. Cada nodo se modela considerando algunas caracter´ısticas de las neuronas reales. Como mencionamos anteriormente una neurona real recibe la información por las dendritas, la procesa en el cuerpo celular y bajo ciertas condiciones emite potenciales de acción. Si la red tiene n nodos, a cada nodo i se le asignan tres cantidades: la entrada neta Si al nodo i, Ui = Ui (t), el estado interno del nodo en el instante t (potencial de membrana de una neurona real), la se˜ nal de salida yi del nodo i (n´ umero de potenciales de acción que dispara la neurona por unidad de tiempo) que es función del estado interno, es decir yi = fi (Ui ) , 23

donde fi transforma el estado interno Ui en una se˜ nal de salida yi . La entrada Si al nodo i, usualmente se considera de la forma Si =

n X

Wij yi ,

j=1

donde Wij es el tercer tipo de magnitud asociada a un nodo y se llama peso sináptico. Esta cantidad nos indica el grado de conexión entre el nodo i y el nodo j. El peso Wij puede ser constante o dependiente del tiempo y puede ser positivo o negativo seg´ un corresponda a una conexión excitatoria o inhibitoria respectivamente. Se considera que el conocimiento se encuentra representado en los pesos sinápticos y que los procesos de aprendizaje implican cambios en los pesos sinápticos. Los modelos de redes neuronales que se obtienen, dependen de cómo se relacionan las cantidades asignadas a cada nodo. El primer modelo que se construyó fue el modelo de Mc Culloch-Pitts [80] que tiene las siguientes caracter´ısticas. 1. Wij constantes.  n  1 , si P W y ≥ θ (la neurona “dispara”), θ umbral de la neurona i i ij i i 2. yi = fi (Ui ) = i=1 

0,

en el otro caso .

En este caso se dice que el nodo o neurona artificial es binario, es decir está activado o está desactivado dependiendo de si la magnitud Si rebasa un umbral o no. En otros modelos lo que se construye es una ecuación diferencial en Ui que describe la evolución temporal de Ui : dUi = Gi (Ui , Si ), i = 1, 2, ..., n . dt Por ejemplo en el modelo de Hopfield [59], [60]: τi

dUi = −αi Ui + (Si + θi ) , dt 24

Si =

n X

Wij gj (Uj ) con g, una función sigmoide .

j=1

En el modelo de Grossberg [43]: 

τi

n X



dUi = −αUi + (γi − βi Ui )  Wij gi (Uj ) + θj  . dt j=1

Cuando se quiere estudiar redes neuronales que aprendan, es decir que sean capaces de modificar sus pesos sinápticos, se debe agregar un sistema de ecuaciones diferenciales que describan la evolución temporal de los pesos sinápticos. Las redes neuronales se forman con un conjunto de neuronas interconectadas. Dependiendo de cómo se agrupan las neuronas, las redes se clasifican en: 1. Redes neuronales monocapa (1 capa). Una capa o nivel es un conjunto de neuronas cuyas entradas provienen de la misma fuente (que puede ser otra capa de neuronas) y cuyas salidas se dirigen al mismo destino. 2. Redes neuronales multicapa. Con respecto a las conexiones en la red, si ninguna salida de las neuronas es entrada de la misma capa o capas precedentes, se dice que la red es de propagación hacia adelante. Cuando las salidas pueden ser conectadas como entradas de neuronas de niveles previos o del mismo nivel, incluyéndose ellas mismas además de conexiones con neuronas de las capas posteriores, se dice que la red es de propagación hacia atrás. Con respecto al aprendizaje de las redes neuronales, modificación de los pesos sinápticos, las redes neuronales se pueden clasificar como: redes neuronales con aprendizaje supervisado y redes neuronales con aprendizaje no supervisado. En las primeras se precisa de un agente externo (supervisor) que controle el proceso de aprendizaje de la red, y que determina la respuesta de la red a partir de una entrada determinada, en las segundas la propia red controla su propio proceso de aprendizaje. En ambos tipos de redes se requiere conocer las reglas de aprendizaje –criterios que se siguen para modificar los pesos sinápticos–. El aprendizaje 25

supervisado se puede llevar a cabo de tres maneras: aprendizaje por corrección de error, por esfuerzo y estocástico. En el aprendizaje no supervisado se puede llevar a cabo, entre otras maneras, el aprendizaje Hebbiano y el aprendizaje competitivo y cooperativo. Para mayor información sobre este tema recurrir a [11], [12], [31] y [55]. El primer modelo de red neuronal artificial, llamada Perceptron, fue desarrollado por Rosenblatt [98]. El Perceptron es una red neuronal con aprendizaje supervisado, con regla de aprendizaje de corrección de error, de propagación hacia adelante que consiste de dos capas, una de entrada - que tiene n neuronas lineales - y la otra capa, que consta de una sola neurona, es la capa de salida. El Perceptron es capaz de decidir cuándo una entrada presentada a la red pertenece a una de las dos clases que es capaz de reconocer. La u ńica neurona de salida del Perceptron realiza la suma ponderada de las entradas, resta el umbral y pasa el resultado a una función de transferencia de tipo escalón. La regla de decisión es responder +1 si el patrón presentado pertenece a la clase A, o −1 si patrón pertenece a la clase B; A y B son las dos clase que el Perceptron es capaz de reconocer, previo entrenamiento. El Perceptron funciona cuando el conjunto de n−adas, cuyas componentes son las entradas a las n−neuronas de la capa de entrada, pueden ser separadas P por un hiperplano n−dimensional - el correspondiente a N i=1 Wi Xi − Θ = 0; donde Wi es el peso correspondiente a la entrada Xi y Θ el umbral. En general esto no siempre se puede hacer. Posteriormente se desarrolla el Perceptron multicapa [11], [31]. El Adaline y el Madaline son redes neuronales desarrolladas por Widrow y Hoff [118] en 1960. La arquitectura del Adaline y del Madaline son esencialmente la misma que la del Perceptron. Lo que cambia es la regla de aprendizaje ya que Adaline y Madaline usan la llamada regla delta de Widrow - Hoff o regla del m´ınimo error cuadrado (LMS), basada en la b´ usqueda de una expresión del error entre la salida deseada y la salida lineal obtenida antes de aplicarle la función de activación escalón (a diferencia de la salida binaria utilizada en el caso del Perceptron). Sus principales aplicaciones están en el campo del procesamiento de se˜ nales. También existen redes con aprendizaje supervisado del tipo estocástico, 26

la máquina de Boltzman [56], etc., del tipo aprendizaje con refuerzo, Linear Reward-Penalty o LR−P . Con respecto a las redes con aprendizaje no supervisado se tienen: la red de Hopfield, Learning Matrix, etc. que utilizan el aprendizaje Hebbiano (ajuste de los pesos de las conexiones de acuerdo con la correlación de los valores de activación de las dos neuronas conectadas). Para las que utilizan el aprendizaje competitivo - cooperativo (las neuronas compiten o cooperan unas con otras con el fin de llevar a cabo una tarea) ver [5], [7]. Para mayor información sobre la arquitectura, reglas de aprendizaje, aplicación y limitaciones de las redes que hemos mencionado y a otras más recurrir a: [31], [55]. Cuando se intenta construir modelos de redes neuronales para el estudio del cerebro, los nodos de las redes se pueden interpretar como una sola neurona o como un grupo de neuronas. Matemáticamente no existe ninguna distinción entre estos dos casos, pero si se desea construir redes neuronales que describan ciertos fenómenos que ocurren en el cerebro y que los resultados se correlacionen con mediciones experimentales, la distinción s´ı es importante. Se ha detectado actividad r´ıtmica en algunas estructuras del cerebro: oscilaciones en el tálamo, oscilaciones en el sistema bulbo olfativo, corteza, etc. Algunos modelos de redes neuronales se elaboran tomando en cuenta datos anatómicos y fisiológicos con la finalidad de describir estos fenómenos temporales r´ıtmicos. Por ejemplo: Lopes da Silva [76], [77] elabora un modelo de red neuronal que intenta explicar los or´ıgenes del ritmo α en el perro con redes neuronales creadas con neuronas de la corteza y del tálamo. Destexhe, Contreras, Sejnowski, y Steriade [24] elaboran un modelo de red neuronal para neuronas del n´ ucleo reticular talámico y obtienen resultados teóricos de tipo oscilatorio con frecuencias de 7 a 14 Hz y otros en el rango de 0.5 - 1 Hz dependiendo de los parámetros que se consideran.

27

Freeman [33], [34], [35], elabora un modelo que describe el modo en que la información sensorial puede ser procesada en un modo distribuido para una de las formas más simples de corteza. El modelo se elabora con base en estudios anatómicos y fisiológicos sobre las interconexiones neuronales entre el bulbo olfativo y la corteza en gatos y conejos e intenta relacionar los resultados teóricos con datos electroencefalográficos obtenidos en experimentos con estos animales. Menon [81], estudia el fenómeno de las oscilaciones temporales en el cerebro, considerando que éstas se deben a la interacción de neuronas inhibitorias y neuronas excitatorias. Elabora un modelo en el cual considera dos grupos de neuronas, una que tiene u ńicamente conexiones inhibitorias y el otro grupo que tiene u ńicamente conexiones excitatorias. Estudia la dinámica de estos dos grupos de neuronas interconectados elaborando un modelo de dos ecuaciones integro-diferenciales acopladas. Nu˜ nez comenta en [89], [92] que la mayor´ıa de los modelos de redes neuronales no surgen para estudiar las bases fisiológicas de los EEG en el cuero cabelludo sino para estudiar algunas funciones cerebrales, como por ejemplo la memoria de corto plazo. También menciona que si se quiere “explicar” por medio de redes neuronales los patrones r´ıtmicos temporales que se miden con electrodos en el cuero cabelludo (EEG), no basta que se construyan redes (a´ un basados en datos anatómicos y fisiológicos) que produzcan soluciones oscilatorias con la misma frecuencia que se obtienen en los EEG. Por lo tanto, es necesario considerar la estructura espacial de la red para considerar los potenciales correspondientes y compararlos con los que se obtienen en los EEG, ya que la magnitud de un potencial depende de la configuración geométrica de sus fuentes. Nu˜ nez elabora un modelo de interacción neuronal en la corteza a la cual considera como un medio continuo, donde la columna de corteza es la unidad fundamental. El modelo consiste en un par de ecuaciones integrales acopladas que relacionan la actividad sináptica excitatoria e inhibitoria en una columna de corteza en un tiempo t con los potenciales de acción que se “disparan” en las otras columnas en tiempos anteriores, tomando en cuenta las conexiones que existen entre las columnas y las velocidades de propagación de los potenciales de acción, ver [87], [88], [89] y [92].

28

Nu˜ nez elabora el modelo con la finalidad de estudiar las bases fisiológicas de los patrones espacio temporales medidos con electrodos ubicados en el cuero cabelludo. Considera que para el estudio de estos patrones espaciotemporales la columna de corteza es la unidad a considerar, ya que los EEG se obtienen con electrodos que captan la actividad eléctrica promedio de varios cent´ımetros cuadrados de corteza, además de la homogenidad que ofrece la corteza a nivel de cada capa. En el siguiente cap´ıtulo analizaremos con más detalle este modelo.

29

1.8 MODELOS DE MEDIO CONDUCTOR En el funcionamiento del cerebro, la actividad eléctrica es un elemento fundamental; la comunicación entre las neuronas y entre estructuras del cerebro se realiza a través de ella. Otro elemento fundamental en el funcionamiento del cerebro es la actividad neuroqu´ımica. Se conoce muy poco de cómo interact´ uan ambos elementos en el funcionamiento del cerebro. La actividad eléctrica del cerebro puede ser medida en muchas escalas espacio-temporales y con diversos grados de resolución se puede medir la actividad eléctrica de una sola neurona, de grupos de neuronas, y usando electrodos colocados en el cuero cabelludo, la del cerebro en su totalidad. Se conocen diversos patrones espacio-temporales de la actividad eléctrica en distintos niveles, y se desea estudiar las bases fisiológicas subyacentes a ellos. Se han usado diversos modelos: redes neuronales y sistemas cooperativoscompetitivos de interacción neuronal, para describir los posibles mecanismos que originan estos patrones espacio-temporales. Pero estos modelos, deben ir acompa˜ nados de otros modelos que toman en cuenta otras caracter´ısticas del cerebro como son las propiedades electromagnéticas derivadas de la activación de fuentes de corrientes impresas en el cerebro y por ello es necesario estudiar cómo se propagan las corrientes eléctricas a través de los distintos tejidos del cerebro. No basta con determinar que ciertas redes neuronales correspondientes a ciertas estructuras del cerebro producen oscilaciones teóricas que coinciden con las oscilaciones detectadas con los EEG medidos en diversas partes del cerebro. Los potenciales dependen de la distribución de las fuentes y de las propiedades eléctricas del medio en que ellos se propagan. En el estudio del cerebro, como medio conductor, una herramienta importante es la teor´ıa del campo electromagnético que en el nivel microscópico nos describe la relación entre cargas, corrientes eléctricas, campos eléctricos y magnéticos y en el nivel macroscópico debe incluir ciertas relaciones experimentales para la permitividad (para aislantes), la conductividad (para conductores) y la permeabilidad (para materiales magnéticos). Pero nuevamente, si se quieren elaborar modelos cuyos resultados se puedan relacionar 30

con datos experimentales, es necesario determinar en qué escala espacio temporal se está trabajando. En cualquier escala espacio temporal que se use, las mediciones electroencefalográficas se obtienen como promedio de la actividad eléctrica de un cierto volumen de tejido y en cierto per´ıodo de tiempo y por tanto, las conductividades, las permitividades y las permeabilidades magnéticas de los tejidos deben considerarse en la misma escala. Además las variables de nuestro modelo deben estar en la misma escala espacio-temporal que la de las mediciones. El problema que se requiere estudiar es determinar la ubicación y propiedades de las fuentes de corriente generadas por la actividad del cerebro a partir de las mediciones eléctricas y magnéticas efectuadas en el cuero cabelludo o en otras estructuras del cerebro. Este problema no tiene solución u ńica; en el sentido de que la misma medición de potencial o campo magnético puede ser generada por diferentes distribuciones de fuentes de corriente concentradas en la masa cerebral. El problema inverso no tiene solución u ńica, aunque se conozca la distribución de potencial (o campo magnético) en todos los puntos del cuero cabelludo. Pero la situación real es todav´ıa peor ya que las mediciones del potencial se obtienen u ńicamente en un conjunto finito de puntos del cuero cabelludo donde están ubicados los electrodos. El problema es sumamente complejo ya que es necesario considerar la geometr´ıa y propiedades electromagnéticas de las distintas regiones que intervienen –corteza, cráneo, cuero cabelludo, etc.– con el objetivo de relacionar los datos obtenidos en el cuero cabelludo con sus posibles fuentes. A´ un el problema directo, es decir, dada la distribución de fuentes de corriente, la conductividad de los diversos tejidos y la permeabilidad magnética, determinar el campo electromagnético en cualquier punto, ya es bastante complejo. Usualmente lo que se hace, en ambos problemas, es introducir simplificaciones con base en información conocida sobre el cerebro y resolver el problema con esas simplificaciones. Estas simplificaciónes dependen mucho de la escala en que se está trabajando. 31

La primera simplificación importante que se introduce es el considerar que se puede aplicar al cerebro la teor´ıa electromagnética cuasi-estática [92]. Es conocido que en un medio conductor con conductividad σ las cargas se distribuyen muy rápidamente, es decir, la densidad ρ de carga decrece como σ

ρ ∼ e− εo t , donde εo es la constante dieléctrica. Para la conductividad del hueso 1 del orden de 200 y para la piel y la masa cerebral es mucho menor. Por este motivo se tiene que

σ εo

es

¯ ¯ ¯ ∂ρ ¯ σ − εσ t ¯ ¯ ¯ ¯∼ e o ¯ ∂t ¯ εo

y como se conoce experimentalmente que el intervalo de tiempo en que la densidad de corriente J~ en el cerebro var´ıa apreciablemente, es mucho mayor que εσo , entonces en la ecuación de continuidad ∂ρ ∇ · J~ + =0 ∂t se puede despreciar

∂ρ ∂t

y queda ∇ · J~ = 0 .

Además se considera que las corrientes que circulan en el cerebro son de dos tipos: Ohmicas, es decir, generadas por cargas en movimiento a través del fluido extracelular y corrientes primarias o impresas, J~p , que se mueven a través de las membranas celulares y aparecen en las zonas activas. Entonces se tiene ~, J~ = J~p − σ E ~ es el campo eléctrico. donde E Del desacoplamiento entre el campo eléctrico y magnético se tiene que ~ = 0 y, por lo tanto, E ~ = ∇u, donde u es el potencial del campo rot E electrostático.

32

De la ecuación de continuidad simplificada se obtiene que en la masa del cerebro se cumple la ecuación ∇u =

1 ∇ · J~p . σ

Si además se considera un sistema de varias capas donde la región Ω1 ocupada por el cerebro está cubierta por otras regiones disjuntas Ωi (i ≥ 2) que representan el fuido intracraneal, el cráneo, la piel, etc., entonces en cada frontera de separación entre dos regiones Ωi y Ωi+1 , se deben satisfacer condiciones del tipo Ã

σi

∂u ∂ni

!−

Ã

= σi+1

∂u ∂ni

!+

+ f~i ,

u− = u+ , donde σi es la conductividad en cada Ωi , ni es el vector unitario normal exterior a Ωi en la parte com´ un de las fronteras de Ωi y Ωi+1 . El supra´ındice − o + indica que el valor l´ımite a la frontera de la cantidad indicada se toma desde Ωi o desde Ωi+1 respectivamente. Mediante f~i se denota la componente normal de la densidad de corriente superficial, la cual en este modelo puede estar presente u ńicamente en la frontera del cerebro, es decir, en la corteza cerebral. Si se denota por Ωe la región cuya frontera exterior Se coincide con el cuero cabelludo, entonces el problema u ńico consiste en determinar J~p y f~1 a partir de las mediciones de u sobre Se . Este problema está mal planteado y por ello se necesita introducir algoritmos de regularización en su solución. Una de las simplificaciones más importantes en el estudio del problema inverso consiste en representar f~1 y J~p en forma de un n´ umero finito de dipolos puntuales y, en este caso, el problema se reduce a la localización de sus posiciones y momentos, es decir, a encontrar 6N parámetros donde N es el n´ umero de dipolos propuestos para la fuente activa [103].

33

1.9 CONCLUSIONES Como hemos podido ver, existe un gran n´ umero de técnicas y métodos para determinar anomal´ıas de carácter anatómico o funcional en el cerebro, as´ı como para el estudio de los patrones temporales ritmicos de la actividad eléctrica del mismo. Sin embargo, ninguno de estos métodos puede describir la interacción neuronal entre diferentes regiones del cerebro ni puede explicar el origen fisiológico de los patrones temporales r´ıtmicos de la actividad eléctrica. El objetivo de este trabajo es estudiar modelos que son una generalización de los introducidos por Paul Nu˜ nez en [89], [92] que, en conjunto con los modelos de medio conductor, permitan describir la evolución temporal de la interacción neuronal masiva en la corteza cerebral y también dar alguna explicación matemáticamente fundamentada de los or´ıgenes fisiológicos de los patrones temporales r´ıtmicos de la actividad eléctrica. Esperamos además que este tipo de trabajo sirva como fundamento para el desarrollo de software más apropiado para el análisis cl´ınico de los EEG.

34

Cap´ıtulo 2 ´ MODELACION GLOBAL DE LA ACTIVI´ DAD SINAPTICA DE LA CORTEZA CEREBRAL ´ 2.1 INTRODUCCION Para el estudio de las correlaciones existentes entre ciertos estados fisiológicos y las mediciones electroencefalográficas –a nivel de cuero cabelludo o de corteza– se han elaborado diversos modelos de interacción neuronal a nivel de corteza-tálamo. Estos modelos se pueden dividir, grosso modo, en modelos locales, modelos globales y modelos mixtos. En los modelos locales [33], [34], [75], [76], [77], [95], se tratan de explicar las propiedades espacio-temporales de los EEG con base en propiedades a nivel de neurona, redes neuronales discretas o redes neuronales distribuidas en el rango de mil´ımetros; las propiedades temporales de los EEG son determinadas por los retrasos sinápticos como son los tiempos de respuesta y de decrecimiento de los potenciales postsinápticos. En los modelos globales [87], [88], [89], [90], [92], se toman como elemento importante los sistemas de fibras corticorticales (cuyas longitudes se encuentran en el rango de cent´ımetros) y el hecho de que la corteza es una estructura cerrada. En estos modelos se considera que el EEG medido sobre el cuero cabelludo es fundamentalmente el resultado de la interacción de neuronas corticales a través de los potenciales de acción. Se considera que 35

ciertos patrones espacio-temporales de los EEG se pueden explicar tomando en cuenta los retrasos debido a la velocidad finita de propagación de los potenciales de acción que viajan por las fibras corticorticales y por el hecho de que la corteza es una estructura cerrada. Nu˜ nez considera [92] que los modelos locales pueden ser importantes en ciertos estados fisiológicos del cerebro y que en otros casos, como en aquellos que interviene un m´ınimo de actividad cognitiva (sue˜ no, coma, anestesia) se pueden explicar con los modelos globales; pero que dada la extrema complejidad de los interacciones neuronales se requieren modelos que toman en cuenta la interacción neuronal en las diversas escalas espacio-temporales [64], [65], [66], [67], [68]. En este cap´ıtulo analizaremos el modelo global elaborado por Nu˜ nez. Este modelo será el punto de partida para la elaboración de nuestro modelo. El modelo describe la evolución de la actividad sináptica en la corteza y consiste en un par de ecuaciones integrales que relacionan la actividad sináptica en cada columna de corteza en un instante de tiempo t con los potenciales de acción que se generan en tiempos anteriores en las demás columnas de corteza, tomando en cuenta las conexiones que existen entre las columnas de corteza y las velocidades de propagación de los potenciales de acción. También se pueden tomar en cuenta los retrasos sinápticos o incluir circuitos locales, convirtiéndolo en un modelo local-global más completo. Con el modelo global, sin tomar en cuenta los retrasos locales, se estudia bajo qué condiciones se tienen soluciones de tipo ondulatorio y se comparan estas soluciones con los datos obtenidos a través del EEG. Sus resultados teóricos “coinciden” con los obtenidos experimentalmente para el ritmo α pero no incluye otros estados fisiológicos. ´ ´ DE 2.2 BASES FISIOLOGICAS PARA LA MODELACION ´ LA ACTIVIDAD SINAPTICA EN LA CORTEZA CEREBRAL La corteza cerebral es un sistema muy complejo que presenta interacciones en distintos niveles: neurona, grupos de neuronas y entre grupos de neuronas. Aun as´ı la corteza presenta ciertas caracter´ısticas anatómicas y fisiológicas 36

que permiten elaborar modelos de interacción neuronal cuyas variables y parámetros se pueden correlacionar con las medidas electroencefalográficas obtenidas en el cuero cabelludo. El modelo que presentaremos toma como unidad fundamental la columna de corteza y pretende describir las interacciones entre las columnas de corteza, a través de sus interconexiones y la interacción de la corteza con otras estructuras del cerebro (principalmente el tálamo); más precisamente se describirá la evolución espacio-temporal de la actividad sináptica de las columnas de corteza. Los elementos del modelo son: a) h± (x, t) − La densidad de sinapsis activas o de actividad sináptica (+ para las sinapsis de tipo excitatorio, − para las de tipo inhibitorio), en la columna x en el tiempo t. h± (x, t) representa el n´ umero de sinapsis activas por cm3 en la columna x en el tiempo t (1/cm3 ) b) g(x0 , t) − La fracción del total de neuronas de la columna x0 que generan potenciales de acción en el tiempo t (adimensional) c) R± (x0 , x, v) − El n´ umero de conexiones sinápticas por cm3 en la columna x por unidad de área en x0 y por unidad de velocidad, es decir, R± (x0 , x, v)dx0 dv representa el n´ umero de conexiones sinápticas por cm3 en la columna x producidas por las neuronas que se encuentran en las columnas de x0 a x0 + dx0 y por cuyos axones los potenciales de acción “viajan” con velocidad entre v y v + dv (+ para las conexiones de tipo excitatorio, − para las de tipo inhibitorio). d) S± (x) − La densidad de sinapsis en la columna x (+, para las sinapsis excitatorias, − para las inhibitorias) (1/cm3 ). Se tiene que: Z ∞Z 0

S

R± (x0 , x, v)dx0 dv = S± (x) ,

donde S es la superficie de la corteza.

37

e) F± (v) − La fracción del total de neuronas de la corteza por unidad de velocidad (seg/cm), es decir, F± (v)dv representa la fracción del total de neuronas de la corteza que generan su potencial de acción con velocidad entre v y v + dv, (+ para las conexiones sinápticas excitatorias, − para las inhibitorias). f ) h0± (x, t) − La densidad de sinapsis activas o de actividad sináptica ( + para las sinapsis de tipo excitatorio, − para las de tipo inhibitorio), en la columna x en el tiempo t. Actividad sináptica producida por los potenciales de acción que se generan en otras estructuras del cerebro, fundamentalmente el tálamo y que “llegan” a la columna x en el tiempo t. h0± (x, t) representa el efecto producido en la columna x en el tiempo t por las estructuras internas del cerebro que tienen conexiones en la columna x. Con estas definiciones es fácil ver que si S es la superficie de la corteza cerebral, entonces se tienen las ecuaciones:

h+ (x, t) =

ho+ (x, t)

h− (x, t) =

ho− (x, t)

+

+

Z ∞Z 0

Z ∞Z 0

d(x, x0 ) R+ (x , x, v)g(x , t − ) dx0 dv , v S

(2.1)

d(x, x0 ) ) dx0 dv , R− (x , x, v)g(x , t − v S

(2.2)

0

0

0

0

donde d(x, x0 ) es la distancia geodésica entre dos puntos de la corteza S. Este sistema de ecuaciones integrales nos dice que la densidad de actividad sináptica tiene dos componentes. La primera de ellas es debida al “input” que reciben la corteza y las fibras corticorticales proveniente de estructuras interiores durante un estado fisiológico fijo. La componente integral expresa la idea de que los potenciales de acción que se disparan en la columna x0 se trasladan a través de las fibras corticorticales causando actividad sináptica en la columna x en un tiempo posterior que depende de la distancia que separa x de x0 y de la velocidad de propagación de los potenciales de acción. 38

En general esta distancia no tiene que coincidir con la euclideana aunque puede ser considerada como tal para el análisis cualitativo. A continuación enumeraremos las hipótesis de trabajo y las principales consideraciones anatómicas para el trabajo con este modelo. 1. No están considerados los retrasos debidos a la transmisión sináptica –retrasos locales– sino, u ńicamente el debido a la velocidad finita de propagación de los potenciales de acción a lo largo de los axones (fibras corticorticales). 2. En el modelo (2.1), (2.2) se desprecia la existencia de circuitos de retroalimentación entre la corteza y el tálamo ya que el n´ umero de fibras corticorticales es mucho mayor que las fibras talamocorticales. En efecto las primeras representan el 99% del total [72]. No resulta d´ıficil incluir en este modelo global la presencia de estos circuitos de retroalimentación. Un ejemplo importante de estos circuitos es el circuito córtico-estriato-nigro-tálamocortical asociado a la actividad motora. 3. Las neuronas corticales están colocadas en columnas que presentan una gran homogeneidad en su distribución espacial con respecto a la profundidad. Esto quiere decir que es imposible distinguir entre dos cortes perpendiculares realizados en dos secciones diferentes de la corteza cerebral. Por este motivo podemos considerar que la corteza tiene una estructura estratificada uniforme compuesta por seis capas debido a las diferencias en las densidades relativas de diferentes tipos de neuronas localizadas a distintas profundidades corticales y que en una primera aproximación las conexiones en el plano de la corteza son del mismo tipo en todas las direcciones y por ello sólo dependen de la distancia entre las columnas consideradas es decir: R± (x0 , x, v) = R± (d(x, x0 ), v) .

(2.3)

Las neuronas corticales se pueden clasificar de acuerdo con la forma del cuerpo celular, por la longitud de su axón y por la extensión y distribución espacial del árbol dendr´ıtico. Además del tipo de transmisor

39

que libera, el blanco axonal, el tipo de sinapsis, etc. Para un estudio más detallado sobre esta tema recurrir a [82]. En primera aproximación, las neuronas corticales las dividiremos en dos grandes tipos: excitatorias (neuronas piramidales y estrelladas espinosas) e inhibitorias (todas las interneuronas). Las interneuronas inhibitorias se dividen en varios tipos y se consideran que pueden tener varias escalas espacio–temporales. Pero en nuestro caso tomaremos que la inhibición ocurre rápidamente. Las piramidales tienden a ocupar vol´ umenes cil´ındricos con axones y dendritas alineadas perpendicularmente a la superficie cortical y aproximadamente de 2/3 a 3/4 de las neuronas corticales son de este tipo. Las neuronas estrelladas ocupan volumenes más esféricos. El cuerpo celular de las neuronas ocupa alrededor de un 10% de la superficie neuronal. La mayor parte del área superficial de la neurona está compuesta por sus ramificaciones dendr´ıticas. El n´ umero de ramas dendr´ıticas que parten de una neurona y que intersectan a una superficie esférica concéntrica, centrada en el cuerpo celular de la neurona, decrece exponencialmente seg´ un una ley del tipo: N´ umero de ramas ≈ exp(−λ· distancia radial) donde λ−1 (cm) es la escala caracter´ıstica de longitud que determina el decrecimiento en e−1 veces del n´ umero total de conexiones corticales a partir de una neurona. Por ejemplo, para las células estrelladas y las dendritas basales de las células piramidales en las cortezas visual y motora del gato y el ratón se ha comprobado experimentalmente que λ = 350 cm−1 . La relación para el decrecimiento de las ramas generadas a partir de una dendrita es de orden polinomial en lugar de exponencial. 4. Las interconexiones entre neuronas corticales pueden ser divididas en dos grandes grupos: las conexiones intracorticales que son de rango corto y que tienen una longitud promedio de menos de 1 mm en humanos y las conexiones corticorticales que tienen una longitud promedio de hasta varios cent´ımetros en humanos. La materia blanca subcortical está compuesta por las fibras de asociación corticorticales que están ubicadas paralelamente a la superficie 40

cortical (Fig. #7). Mientras que las conexiones intracorticales pueden ser tanto excitatorias como inhibitorias, las fibras de asociación (axones de las células piramidales) que constituyen la mayor parte de la materia blanca son solamente conexiones de tipo excitatorio. Las neuronas piramidales que env´ıan axones a otras regiones interiores del cerebro como el tálamo, terminan u ńicamente en sinapsis excitatorias. En las neuronas piramidales existen otro tipo de neuronas intracorticales que también terminan en sinapsis excitatorias. Sin embargo la cantidad de sinapsis excitatorias en humanos realizadas por las fibras de asociación (corticorticales) es mucho mayor que las que son producidas por axones intracorticales, aproximadamente el 80% del total de las sinapsis excitatorias son por las fibras de asociación. Aunque en ratones es del orden del 50%. Para mayor información sobre las fibras de asociación se puede consultar [72]. 5. La mayor´ıa de las fibras de asociación mielinadas tienen diámetros entre 1 y 2 µ y en ellas las velocidades de los potenciales de acción son proporcionales al diámetro y fluct´ uan entre 6 y 12 m/seg. 6. El área de la superficie total de ambas mitades del cerebro, las cuales están separadas por el cuerpo calloso es del orden de 1600 cm2 . Las 2/3 partes de esta superficie están en las fisuras incluyendo la fisura sagital que separa las dos mitades del cerebro. De estas dos terceras partes la mitad del área está contenida en las peque˜ nas fisuras incluyendo la fisura sagital. Esto nos dice que el área efectiva del medio donde se propagan las ondas es diferente de la superficie cortical. Un problema importante en los modelos globales del EEG consiste en saber cuándo ambas mitades del cerebro deben ser consideradas por separado. La materia blanca está formada por fibras de aferentes y eferentes que van y vienen del interior del cerebro y la médula espinal, junto con fibras de asociación que conectan diferentes áreas corticales

41

y fibras comisurales que conectan las dos mitades del cerebro fundamentalmente a través del cuerpo calloso. En humanos el n´ umero total de fibras comisurales que van de un hemisferio a otro del cerebro a través del cuerpo calloso es del orden de 2 × 108 , la cual es significativamente menor que el n´ umero de fibras de asociación que entran a la corteza, que es del orden de 1010 . Esto hace pensar que las fibras comisurales permiten una conexión muy débil entre los dos hemisferios cerebrales, por lo cual parece ser que el análisis de los fenómenos ondulatorios debidos a la actividad cerebral debe restringirse a un sólo hemisferio que puede identificarse con la cuarta parte de un elipsoide de revolución. El cerebro humano tiene 45cm de circunferencia (sin contar las fisuras) de la superficie: del polo occipital al polo frontal, a lo largo de la l´ınea media. Esta idea resulta consistente con los resultados experimentales obtenidos en pacientes a quienes se les ha medido la frecuencia del ritmo α con los dos hemisferios acoplados y con cuerpo calloso cortado. Si los dos hemisferios acoplados se comportaran como un sistema u ńico se esperar´ıa que al efectuar el corte en el cuerpo calloso hubiera un significativo incremento de la frecuencia del ritmo α, lo cual no ocurre. 7. Los potenciales del EEG medidos en el cuero cabelludo se deben a una actividad sincrónica de neuronas corticales que cubren un área de al menos varios cent´ımetros cuadrados de corteza. Por este motivo la medición del EEG corresponde a un nivel jerárquico macroscópico de modelación donde la homogeneidad e isotrop´ıa son dos consideraciones importantes. En cambio, del hecho que las conexiones corticorticales son inhomógeneas a niveles jerárquicos peque˜ nos se concluye que los modelos de actividad medida con microelectrodos requiere la inclusión de parámetros de conexión preferencial. Como las neuronas piramidales generan grandes momentos dipolares y están orientadas perpendicularmente a la superficie de la corteza, ellas son capaces de producir corrientes sincronizadas de suficiente magnitud para ser medidas sobre el cuero cabelludo. Debido a esto y a la superioridad numérica de las neuronas piramidales junto con el hecho de que 42

son ellas las que producen la mayor parte de sinapsis excitatorias en la corteza a través de sus axones (fibras corticorticales), es que son estas neuronas las principales responsables de la medición del EEG sobre el cuero cabelludo. De aqu´ı concluimos que para una primera modelación de la contribución global de la corteza al EEG se debe considerar la actividad sináptica excitatoria producida por las neuronas piramidales as´ı como la actividad sináptica inhibitoria producida sobre ellas por las neuronas estrelladas. Como mencionamos al inicio de este cap´ıtulo, en esta tesis nos dedicaremos al estudio de modelos teóricos de la actividad global de la corteza cerebral que dan lugar al EEG. Estos modelos tienen en cuenta las consideraciones anatómicas que hemos descrito en los puntos 1-7 anteriores. A diferencia de ellos, los modelos de circuitos locales deben incluir parámetros de conexión preferenciales, en lugar de la homogeneidad que se considera en los modelos de actividad global, lo cual se pone de manifiesto en la selección de las funciones R± (x0 , x, v) para cada v´ıa del circuito. Si los circuitos son peque˜ nos puede considerarse que la velocidad de propagación v es infinita no siendo as´ı en los circuitos en que aparecen interacciones de largo alcance, en los cuales el tiempo que se demora al viajar un potencial de acción de un lugar a otro del circuito puede ser del orden de las constantes de tiempo de las neuronas piramidales. Por otra parte, como ya mencionamos en el punto 1, los modelos locales deben considerar los retrasos temporales debido a la transmisión sináptica, que son más peque˜ nos que los retrasos debidos a la velocidad finita de propagación de los potenciales de acción a través de los axones, y por ello no se consideran en los modelos globales.

43

´ Y SIMPLIFICACIONES EN EL 2.3 DESCRIPCION ´ MODELO GLOBAL DE ACTIVIDAD SINAPTICA EN UNA BANDA DE CORTEZA. Para evitar repeticiones innecesarias, en este cap´ıtulo el sub´ındice + siempre se referirá a las sinapsis de tipo excitatorio y el sub´ındice − a las de tipo inhibitorio Aplicaremos el modelo general (2.1), (2.2) al estudio de la actividad del EEG producida por una banda infinitamente delgada de la corteza que por conveniencia, supondremos infinitamente larga. Para el análisis cualitativo obviaremos la geometr´ıa de la corteza y supondremos que la banda coincide con una recta. En este caso, que llamaremos el caso unidimensional, los elementos del modelo son los mismos, sólo cambian las unidades con que se trabajará. As´ı, h± (x, t) representará el n´ umero de sinapsis activas por cm2 en la columna x en el tiempo t (1/cm2 ). R± (x0 , x, v) representará el n´ umero de sinapsis, por cm2 en la columna x en el tiempo t por unidad de longitud en x0 y por unidad de velocidad (seg/cm4 ). Es decir, R± (x0 , x, v)dx0 dv representará el n´ umero de sinapsis, por cm2 en la columna x en el tiempo t producidas por las neuronas que están ubicadas desde la columna x0 hasta la columna x0 +dx0 y por cuyos axones, los potenciales de acción “viajan” con velocidad entre v y v + dv, h0± (x, t) representará el n´ umero de sinapsis activas en la columna 2 x por cm , esta actividad sináptica es la producida por los potenciales de acción que se generan en estructuras internas del cerebro, principalmente el tálamo y que llegan en el tiempo t a la columna x, g(x0 , t) seguirá siendo una cantidad sin dimensiones. Todas las complicaciones del sistema de conexiones intracorticales e intercorticales están “escondidas” en la representación de las funciones R+ y R− . Introduciremos algunas hipótesis adicionales que nos permitirán obtener una expresión para estas funciones que al menos cualitativamente reune las caracter´ısticas principales de este sistema de conexiones.

44

Hip´ otesis 1 Comenzaremos suponiendo que cada columna está constituida por una cantidad finita N (que no depende de la columna) de sistemas de fibras que se superponen. En primera aproximación se podr´ıa considerar que N = 6 corresponde a la distribución de las neuronas en las 6 capas en que usualmente se divide la corteza.

Hip´ otesis 2 Como las sinapsis inhibitorias se producen sólo a nivel local, consideraremos que la densidad de actividad sináptica inhibitoria en una columna de corteza x en un tiempo t se produce por los potenciales de acción que se generan en las columnas “cercanas” a la columna x en tiempos “anteriores”, es decir, las columnas ubicadas fuera de una vecindad de la columna de corteza x no contribuyen a la actividad sináptica inhibitoria en la columna x. Es pertinente mencionar que Nu˜ nez enfoca esta cuestión de manera diferente: considera que la velocidad de propagación de los potenciales de acción que “viajan” por las fibras que tienen conexión sináptica de tipo inhibitorio es “infinita”. Esta hipótesis la justifica por el hecho de que la distancia que recorren los potenciales de acción por estas fibras es muy peque˜ na. Nosotros consideramos que esta hipótesis no refleja que las sinapsis inhibitorias se producen sólo a nivel local, lo que refleja es que el “efecto” de los potenciales de acción que “viajan” por las fibras que tienen conexión sináptica de tipo inhibitorio es “sentida” instantáneamente en las columnas donde se tienen estas sinapsis inhibitorias. En realidad, para tener en cuenta el hecho de que las sinapsis inhibitorias se producen sólo a nivel local, debemos considerar que los potenciales de acción que “viajan” por las fibras que tienen conexión sináptica de tipo inhibitorio, recorren esas distancias en tiempos más peque˜ nos que los usados por los potenciales de acción que “viajan” por las fibras de largo alcance. Más adelante estudiaremos un modelo que tome en cuenta esta hipótesis. Desde el punto de vista matemático, nuestra hipótesis (el suponer que la actividad sináptica inhibitoria en una columna de corteza es afectada 45

u ńicamente por los potenciales de acción que se generan en columnas cercanas) nos conducirá a estudiar, en el caso l´ımite, un proceso global de propagación de ondas de actividad sináptica inhibitoria modulada por un proceso local de propagación de ondas de actividad sináptica inhibitoria; mientras que en el caso de Nu˜ nez, lo que se estudia es un proceso global de propagación de ondas de actividad sináptica excitatoria acoplada con un proceso estacionario de actividad sináptica inhibitoria. De esta forma definiremos, para cada n = 1, ..., N : F±n (v) - La fracción del total de neuronas de la corteza por unidad de velocidad (seg/cm) cuyos potenciales “viajan” por las fibras del sistema n-ésimo. Es decir, F±n (v)dv representará la fracción del total de neuronas de la corteza cuyos potenciales “viajan” por las fibras del n-ésimo sistema con velocidad entre v y v + dv. Supondremos que esta cantidad no depende de en cuál columna se ubiquen las neuronas. Ahora la función R± (x0 , x, v) se puede representar en la forma R± (x0 , x, v) =

N X

n R± (x0 , x, v) .

(2.4)

n=1 n R± (x0 , x, v) (seg/cm2 ) representa la densidad de conexiones sinápticas en la columna x por unidad de longitud en x0 y por unidad de velocidad ason ciadas con las fibras del n-ésimo sistema. As´ı R± (x0 , x, v)dx0 dv representará el n´ umero de sinapsis por cm2 en la columna x producidas por las neuronas que están ubicadas desde la columna x0 hasta la columna x0 + dx0 y cuyos potenciales de acción “viajan” por el n-ésimo sistema de fibras con velocidad entre v y v + dv. n Supondremos que R± se representa en la forma n n R± (x0 , x, v) = r± (x0 , x)F±n (v) ,

(2.5)

n (x0 , x) (1/cm3 ) será la densidad de conexiones sinápticas en la donde r± columna x, asociadas al n-ésimo sistema de fibras, por unidad de longitud en x0 , as´ı rn (x0 , x)dx0 representará el n´ umero de sinapsis por cm2 en la columna x producidas por neuronas que están ubicadas desde la columna x0 hasta la columna x0 + dx0 y cuyos axones están ubicados en el n-ésimo sistema de fibras.

46

n Ahora r± (x0 , x) se puede escribir como el producto de dos factores n r± (x0 , x) = tn± (x0 , x)S±n (x) ,

(2.6)

donde tn (x0 , x) (1/cm) representa el n´ umero de neuronas por unidad de longitud en x0 , cuyos axones establecen conexión sináptica en la columna x y dichos axones se encuentran en el n-ésimo sistema de fibras, S±n (x) representa la densidad de sinapsis (1/cm2 ) en la columna x, producida por neuronas cuyos axones se encuentran en el n-ésimo sistema de fibras. De (2.6) se tiene n que R± se puede escribir como el producto de tres factores n R± (x0 , x, v) = tn± (x0 , x)S±n (x)F±n (v) .

(2.7)

Ahora, haremos la siguiente suposición Hip´ otesis 3 La densidad de fibras representadas por la función tn± y que conectan neuronas de la columna x0 con la columna x es igual a la densidad de fibras que conectan neuronas de la columna x con x0 . Esta hipótesis es fundamental para la simplificación del modelo (2.1), (2.2) y significa que tn± (x0 , x) = tn± (x, x0 ) .

(2.8)

Ya hemos comentado que la función tn± (x0 , x) decrece exponencialmente con la distancia entre x0 , y x. Para que se cumpla la propiedad (2.8) supondremos además que: Hip´ otesis 4 La función tn± (x0 , x) depende u ńicamente de la distancia entre los puntos x y x en la forma: 0

λn± −λn± |x0 −x| e , (2.9) 2 donde (λn± )−1 (cm) es la escala longitudinal que determina el decrecimiento de conexiones entre dos puntos de la corteza. tn± (x0 , x) =

47

De (2.4) a (2.9) obtenemos R± (x0 , x, v) =

N 1X n 0 S±n (x)λn± F±n (v)e−λ± |x −x| . 2 n=1

(2.10)

Hip´ otesis 5 En primera aproximación supondremos que cada sistema de fibras está n caracterizado por una velocidad v± de propagación de potenciales de acción. En este caso: n ). F±n (v) = δ(v − v±

(2.11)

N 1X 0 n −λn S±n (x)λn± δ(v − v± )e ± |x −x| . 2 n=1

(2.12)

Finalmente obtenemos R± (x0 , x, v) =

Pasemos ahora a imponer las restricciones que exigiremos a la función g(x0 , t). Hip´ otesis 6 Los potenciales de acción que se generan en la columna x0 en el instante t dependen de la densidad de actividad sináptica - potenciales postsinapticos - y de la intensidad de la actividad neuroqu´ımica (concentraciones de neurotransmisores) que se realice en dicha columna en ese instante. En este trabajo sólo consideraremos la dependencia de g con respecto a h+ y h− , es decir g(x0 , t) = G(h+ (x0 , t), h− (x0 , t)) ,

(2.13)

donde la función G es suave y tiene un comportamiento sigmoidal cuya forma exacta es muy d´ıficil de predecir aunque esto no es importante para el análisis cualitativo que llevaremos a cabo posteriormente.

48

Es decir, la función G(h+ , h− ) es creciente con respecto a h+ para cada h− fijo y decreciente con respecto a h− para cada h+ fijo. Además existen curvas de nivel γ(h+ , h− ) = α, γ(h+ , h− ) = β ,

(2.14)

tales que   0,      

si γ(h+ , h− ) < α ,

G(h+ , h− ) =  0 < G < 1 , si α < γ(h+ , h− ) < β ,       1,

(2.15)

si γ(h+ , h− ) > β .

Para cada h− fijo (respectivamente h+ fijo) las funciones G(·, h− ) (resp. G(h+ , ·)) deben tener un comportamiento sigmoidal usual. Más adelante propondremos una expresión expl´ıcita para las curvas de nivel (2.14). Las relaciones sigmoidales son muy frecuentes en los modelos del sistema nervioso. Por ejemplo, también existe una relación de este tipo entre las concentraciones de neurotransmisores y las frecuencias de disparo de las neuronas o también la frecuencia de actividad sináptica. Todas las hipótesis anteriores nos conducen a una simplificación del modelo determinado por las ecuaciones (2.1) y (2.2) h+ (x, t) = ho+ (x, t) + 12

PN

n=1

S+n (x)λn+

R∞ −∞

n

0

e−λ+ |x −x| G(h+ (x0 , t −

|x0 −x| ), h− (x0 , t n v±

−

|x0 −x| ))dx0 vn ±

, (2.16)

h− (x, t) = ho− (x, t) + 12

PN

n=1

S−n (x)λn−

R x+² −λn |x0 −x| 0| |x0 −x| 0 − G(h+ (x0 , t − |x−x ), h (x , t − ))dx0 , n − x−² e v vn −

49

−

(2.17)

donde le parámetro ² > 0 corresponde al carácter local de la actividad inhibitoria. Sustituyendo (2.13) en (2.16) y (2.17) obtenemos h+ (x, t) = ho+ (x, t) +

1 2

R 0 −λn n ∞ n + |x −x| g(x0 , t n=1 S+ (x)λ+ −∞ e

PN

−

|x0 −x| )dx0 n v+

(2.18) ,

h− (x, t) = ho− (x, t) +

1 2

PN

n n n=1 S− (x)λ−

R x+² −λn |x0 −x| 0| − g(x0 , t − |x−x )dx0 . x−² e vn

(2.19)

−

En el sistema (2.18), (2.19), existen un conjunto de simplificaciones no mencionadas expl´ıcitamente con anterioridad. Consideramos que es importante explicar estas simplificaciones, ya que además de provocarnos algunas confusiones, el considerarlas nos permitirá la construcción de un modelo más detallado, el cual también consistirá de un sistema de ecuaciones integrales. Este sistema se podrá transformar en un sistema de ecuaciones diferenciales parciales de tipo hiperbólico en donde podremos analizar con mayor detalle las distintas “fuerzas” que act´ uan en el medio cortical. Es pertinente mencionar que al intentar transformar el sistema (2.18), (2.19), en un sistema de ecuaciones diferenciales hiperbólicas nos surgieron algunas complicaciones que no permitieron dicha transformación. Algunas de las simplificaciones introducidas en el sistema (2.18), (2.19) y no mencionadas por Nu˜ nez son: No distingue por capas la densidad de actividad sináptica ni la generación de los potenciales de acción; este hecho nos obliga a considerar que el efecto producido en la actividad sináptica en una columna x en un instante t, producido por los potenciales de acción que se generan en tiempos anteriores en las demás columnas, se distribuye de manera uniforme sobre cada capa de la columna x (en cada capa el efecto total se dividirá por N , donde N es el n´ umero de capas). A nivel de columna tampoco distingue la contribución de la actividad sináptica de cada capa en la columna x a la generación de 50

los potenciales de acción en esa columna. Otro problema de estos modelos es que tampoco distinguen la actividad sináptica producida por neuronas de diferente clase. En la sección siguiente presentaremos otros modelos que subsanan algunas de estas simplificaciones. 2.4 ALGUNOS MODELOS ALTERNATIVOS AL MODELO ˜ DE NUNEZ En esta sección analizaremos algunos modelos alternativos al modelo de Nu˜ nez. Estos modelos subsanan algunas de las simplificaciones mencionadas en la sección anterior. También analizaremos brevemente la relación de estos modelos con el modelo de Nu˜ nez y las relaciones que tienen entre s´ı. 1. El modelo más general de capas ser´ıa considerar todas las interacciones que se pueden establecer entre todas las capas de cada par de columnas, as´ı obtendr´ıamos el sistema de 2N 2 ecuaciones integrales siguiente para un modelo de N capas. 0j hij ± (x, t) = h± (x, t)

+

R∞R∞ 0

ij 0 ij 0 −∞ R± (x , x, v)g (x , t −

|x−x0 | )dx0 dv vij

(2.20) ,

donde: hij aptica en la j-ésima capa de ± (x, t) − La densidad de actividad sin´ la columna x en el tiempo t producida por potenciales de acción que se generan en las otras columnas en tiempos anteriores y que viajan por el i-ésimo sistema de fibras, junto con la provocada por los potenciales de acción que se disparan en tiempos anteriores en estructuras interiores del cerebro. aptica en la j-ésima capa de h0j ± (x, t) − La densidad de actividad sin´ la columna x producida por potenciales de acción que se generan 51

en estructuras internas del cerebro y que tienen conexión sináptica con la capa j-ésima de la columna x. ij R± (x0 , x, v) − El n´ umero de conexiones sinápticas/cm2 en la j-ésima capa de la columna x por unidad de longitud en la capa i-ésima de la columna x0 y por unidad de velocidad.

g ij (x0 , t) − La fracción del total de neuronas de la capa i-ésima de la columna x0 que generan potenciales de acción en el tiempo t y que tienen contacto sináptico con neuronas de la capa j-ésima de las otras columnas. Con respecto al modelo de Nu˜ nez se tienen las siguientes relaciones R± (x, t) =

N X N X

ij (x, t) , R±

(2.21)

i=1 j=1

gij (x, t) = αij (x)g(x, t) ,

h± (x, t) =

N X N X ij

h± (x, t) ,

(2.22)

(2.23)

i=1 j=1

donde αij (x) es un coeficiente entre 0 y 1 que indica la parte de g(x, t) que se disparan en la capa i-ésima de la columna x y que llegan a la capa j-ésima de las otras columnas. 2. Otra alternativa ser´ıa considerar hj± (x, t) = h0j ± (x, t) +

1 N

R∞R∞ 0

−∞

j R± (x0 , x, v)

donde: 52

PN

i=1

αij g i (x0 , t −

|x0 −x| )dx0 dv i v±

, (2.24)

hj± (x, t) − La densidad total de actividad sináptica en la j-ésima capa de la columna x en el tiempo t. g j (x0 , t) − La fracción del total de neuronas de la capa j-ésima de la columna x0 que generan potenciales de acción en el tiempo t. h0j aptica en la capa j-ésima de ± (x, t) − La densidad de actividad sin´ la columna x producida por potenciales de acción que se generan en tiempos anteriores en estructuras internas del cerebro. j R± umero de conexiones sinápticas/cm2 en la columna (x0 , x, v) − El n´ x que llegan a la capa j-ésima de la columna x, por unidad de longitud en la columna x0 y por unidad de velocidad.

αij − Los pesos que corresponden a la contribución de la capa iésima de la columna x0 a la actividad sináptica en la capa j-ésima de la columna x. Obviamente 0 ≤ αij ≤ 1 y podr´ıa considerarse P αij = αij (x, x0 , t). También se tiene que N j=1 αij ≤ 1. Con respecto al modelo 1 se está eliminando la contribución particular de cada capa de la columna x0 a la actividad sináptica en la capa j-ésima de la columna x, además de que se tiene un problema similar que en el sistema (2.18) (2.19) ya que no se distingue tampoco el efecto de g i (x0 , t) en cada capa de la columna x, asumiendo que se distribuye uniformemente sobre cada capa de x (la contribución a la actividad 0 0 sináptica en la capa j-ésima de la columna x ser´ıa g (xN ,t) ). La ventaja que se tiene en este modelo es que g j (x0 , t) depende u ńicamente de j j 0 0 h+ (x , t) y h− (x , t); es decir de la actividad sináptica en esa capa y no de la actividad sináptica en todas las capas de esa columna. Se puede considerar que se cumple la relación h± (x, t) =

PN

i=1

βi (x)hi± (x, t) ,

donde βi (x) son ciertos pesos positivos tales que

PN

i=1

βi (x) = 1.

3. También se puede considerar el siguiente sistema 53

(2.25)

¯ j± (x, t) = h0 (x, t) h ± +

R∞R∞ 0

j 0 j 0 −∞ R± (x , x, v)g (x , t

−

|x−x0 | )dx0 dv v

(2.26) ,

donde ¯ j± (x, t) − La densidad de actividad sináptica en la columna x en el h tiempo t producida por potenciales de acción que se generan en la capa j-ésima de las otras columnas en tiempos anteriores más la actividad sináptica en la columna x producida por potenciales de acción que se generan en tiempos anteriores en otras estructuras interiores del cerebro. En este caso no se distingue, como en el caso 2, la actividad sináptica en la columna por capas, sino por los potenciales de acción que provienen de cada capa en otras columnas. g¯j (x0 , t) − La fracción del total de neuronas de la columna x0 que generan potenciales de acción en el tiempo t y que viajan por las fibras del j-ésimo sistema. h0± (x, t) − La densidad de actividad sináptica en la columna x en el tiempo t producida por los potenciales de acción que se generan en tiempos anteriores en otras estructuras interiores del cerebro. j (x0 , x, v) − El n´ R± umero de conexiones sinápticas por cm2 en la columna x por unidad de longitud en la capa j-ésima de la columna x0 y por unidad de velocidad.

El problema que se presenta en este modelo es que g¯j (x, t) depende de toda la actividad sináptica que ocurre en la columna x, es decir de ¯ N (x, t)(x, t). ¯ 1 (x, t), . . . , h h ± ± Estos tres modelos y el modelo de Nu˜ nez son las posibilidades más importantes que se pueden encontrar en un modelo de N capas. Estos 4 modelos están relacionados. El modelo 1 es el modelo más detallado de los cuatro, pero es el que incluye más ecuaciones y parámetros que nos dificultar´ıan el 54

estudio del mismo. Estos modelos están relacionados entre s´ı de la siguiente manera: En el modelo 2 se está considerando la densidad de actividad sináptica por capas en cada columna pero sin distinguir la contribución a esta actividad de las distintas capas de las otras columnas, además tampoco se distingue el efecto de una capa de una columna (g j (x, t)) sobre cada capa de las otras columnas, se considera que el efecto (g j (x, t)) se distribuye de manera uniforme sobre las capas de cada columna afectada. La relación entre el modelo 2 y el 1 se da a través de considerar que Ã

g ij

|x − x0 | x0 , t − vij

!

Ã

= αij g i

|x − x0 | x0 , t − vi

!

,

(2.27)

i.e. el efecto de la i-ésima capa de la columna x sobre la j-ésima capa de las columnas afectadas es independiente de la capa que afecta y se distribuye de manera “uniforme” sobre cada capa. También se cumple que N X 1 i 0 ij R± (x , x, v) = (x0 , x, v) R± N j=1

(2.28)

y por u ´ltimo hi± (x, t) =

N X ij

h± (x, t) .

(2.29)

hi± (x, t) ,

(2.30)

j=1

En el modelo 3 se tiene obviamente que

h± (x, t) =

N X i=1

donde h+ (x, t) es la densidad total de actividad sináptica que aparece en el modelo de Nu˜ nez. Nótese, además que g j (x0 , t) = g¯j (x0 , t), pero g j sólo j ¯ i , h− se tiene que g¯j ; depende de h+ y h− con respecto a las variables h + N 1 ¯ . Con más exactitud ¯ ,...,h depende de h ± ± ¯ j (h+ , h− ) = G ¯j g¯ = G j

ÃN X i=1

55

¯ i , h− h +

N X i=1

!

¯i h −

(2.31)

y como el modelo de Nu˜ nez se relaciona con el modelo 2 a través de h± (x, t) =

N X j

h± (x, t) ,

(2.32)

j=1

g j (x, t) =

g(x, t) . N

(2.33)

Entonces también el modelo (2.18) (2.19) se relaciona con el modelo 1. El modelo que eligiremos para estudiar será el modelo 3 ya que es un poco más detallado que el de Nu˜ nez y resulta un modelo con menos ecuaciones y parámetros que los otros dos modelos restantes. En la siguiente sección desarrollaremos el modelo 3. 2.5 DESARROLLO Y ESTUDIO DE UN MODELO DE LA ´ ACTIVIDAD SINAPTICA DE LA CORTEZA CEREBRAL En esta sección desarrollaremos el modelo 3, que como mencionamos en la sección anterior es un poco más preciso que el del sistema (2.18), (2.19), además este sistema se puede transformar en un sistema de ecuaciones diferenciales parciales. El modelo tiene también una conexión muy clara en el sistema (2.18), (2.19). En este modelo distinguiremos la contribución de cada capa de cada columna a la actividad sináptica de una columna de corteza dada en un instante de tiempo y lo que estudiaremos es la evolución espaciotemporal de esa actividad sináptica tomando en cuenta la distribución de conexiones sinápticas, las velocidades de propagación de los potenciales de acción de las neuronas de la capa correspondiente. Para ser más precisos introduciremos los elementos de que consiste el modelo. Sean: ¯ j± (x, t) − La densidad de actividad sináptica (potenciales postsinápticos) h en la columna x en el tiempo t, n´ umero de sinapsis activas/cm2 , producida por potenciales de acción que se generan en las otras columnas y que viajan por las fibras del sistema j-ésimo, más la actividad sináptica 56

en la columna x producida por potenciales de acción que se generan en tiempos anteriores en otras estructuras interiores del cerebro. g¯j (x, t) − La fracción del total de neuronas de la columna x0 que generan potenciales de acción en el tiempo t y que viajan por las fibras del j-ésimo sistema. ¯ 0 (x, t) − La densidad de actividad sináptica en la columna x en el tiempo h ± t producida por los potenciales de acción que se generan en estructuras interiores del cerebro en tiempos anteriores y que tienen conexión sináptica en la columna x. j R± (x0 , x, v) − El n´ umero de conexiones sinápticas/cm2 , del j-ésimo sistema en la columna x por unidad de longitud en x0 y por unidad de velocidad.

¯ j± (x, t), g¯j (x, t) y h± (x, t), g(x, t) es la siguiente La relación entre h N X ¯j

h± (x, t) .

h± (x, t) =

(2.34)

j=1

La misma idea que se usó para obtener la ecuación (2.18) se puede aplicar ¯ j± (x, t). As´ı que para para obtener un sistema de ecuaciones integrales para h j = 1, . . . , N ¯ 0 (x, t) ¯ j+ (x, t) = h h + +

R ∞ −λj |x−x0 | j 0 1 j j e + g¯ (x , t λ S (x) −∞ 2 + +

−

|x0 −x| )dx0 j v+

(2.35) .

¯ j− (x, t) se tiene Para h ¯ j− (x, t) = h ¯ 0 (x, t) + 1 S−j (x)λj− h − 2

Z x+² x−²

j

0

e−λ− |x−x | g¯j (x0 , t −

57

|x − x0 | )dx0 . (2.36) j v−

Ahora, sabemos que los potenciales de acción que se generan en cada capa de cada columna (¯ g j (x, t)) dependen de la actividad sináptica que existe en esa columna (no depende u ńicamente de la actividad sináptica en esa capa ¯ j± (x, t)) y de la actividad neuroqu´ımica en la por la forma en que se definió h misma columna. En este trabajo sólo consideraremos la dependencia de g¯j (x, t) con res¯ j± (x, t), j = 1, . . . , N ) es decir consideraremos pecto a (h N X

¯j ( g¯j (x, t) = G

¯ i (x, t), h +

N X

hi− (x, t)) .

(2.37)

i=1

i=1

La forma exacta de esta dependencia es d´ıficil de predecir, pero cualitativamente consideramos que debe ser de tipo sigmoide en cada variable. As´ı ¯ i (x, t), i 6= k y h ¯ j− (x, t) fijo debemos para cada k = 1, . . . , N y para cada h + ¯ 1 , ..., h ¯ k−1 , h ¯ k , . . . , hN , h ¯1 , . . . , h ¯ N ) es una función creciente de ¯ j (h tener que G + + + + − − tipo sigmoide. ¯1 , . . . , h ¯ N fijo debemos tener que Para cada 1 ≤ k ≤ N y para cada h + + ¯1 , . . . , h ¯ N , ·) es una función decreciente de tipo sigmoide. ¯ j (h G + + Además existen curvas de nivel ¯ 1 , ..., h ¯ j+ , . . . , hN , h− ) = βj , (2.38) ¯ 1 , ..., h ¯ j+ , . . . , hN , h− ) = αj y γj (h γj (h + + + + tales que

¯ N , h− ) = ¯1 , h ¯ j+ . . . h ¯ j (h G + +

  0,             

¯1 , h ¯ N , h− ) < αj , ¯ j+ . . . h si γj (h + +

¯j < 1 , 0
¯1 , h ¯ N , h − ) < βj , ¯ j+ . . . h si αj < γj (h + +

1,

¯ j+ . . . h ¯ N , h− ) > βj . ¯1 , h si γj (h + + (2.39)

Más adelante propondremos una expresión expl´ıcita para las curvas de nivel (2.38). Con estas condiciones el sistema (2.35), (2.36) se transforma en 58

¯ 0 (x, t) ¯ j+ (x, t) = h h + + 12 λj+ S+j (x)

R∞ −∞

j

0 PN ¯ i |x0 −x| PN 0 )dx0 , ), i=1 hi− (x0 , t − |xv−x| i j i=1 h+ (x , t −

¯j ( e−λ+ |x −x| G 0

v+

+

(2.40)

¯ j− (x, t) = h ¯ 0 (x, t) h − + 12 λj− S−j (x)

R x+² −λj |x0 −x| j PN PN ¯ i |x0 −x| 0 ¯ i (x0 , t − |x0 −x| ¯ ( i=1 h − G ), h (x , t − )dx0 ) . j i i=1 − + x−² e v− v −

(2.41)

Como ya hemos mencionado anteriormente consideramos más adecuado para el análisis del modelo, transformar el sistema (2.35), (2.36) en un sistema ¯ j− (x, t). Esto u ¯ j+ (x, t), h ´ltimo lo de ecuaciones diferenciales parciales para h haremos en la siguiente sección.

59

´ DEL SISTEMA DE ECUACIONES 2.6 TRANSFORMACION INTEGRALES EN UN SISTEMA DE ECUACIONES DIFERENCIALES PARCIALES Usando la regla de Leibnitz para derivar bajo el signo de la integral, ¯j ∂2h

¯j ∂h

¯j ∂h

¯j ∂2h

podemos obtener ∂x± , ∂x2± , ∂t± , ∂t2± , no escribiremos estas derivadas ya que además de constar de muchos términos, no son necesarias para los siguientes desarrollos. Escribiremos solamente el sistema de ecuaciones diferenciales ¯ j+ como para h ¯ j− . parciales resultante, tanto para h ¯ j+ se obtiene Para h ¯j ∂2h + ∂t2

j + 2λj+ v+

· j 2 ) (v+

¯j ∂2h + ∂x2

¯j ∂h + ∂t

j 2¯j + (λj+ v+ ) h+ =

j 0 ¯j ) (x) ∂ h 2(S+ + j ∂x S (x)

−

µ

j 00 j j 0 ) (x) (x)(S+ ) (x))2 −S+ 2((S+

+

j 2 (S+ ) (x)

+

h

¶

j j j j g¯ (x, t) + S+ (x) λj+ v+ +λj+ v+

i

∂¯ gj (x, t) ∂t

¯ j+ h

¸

(2.42)

+ F+j (x, t) .

¯ j− se obtiene Para h ¯j ∂2h − ∂t2

j + 2λj− v−

· j 2 ) (v−

¯j ∂2h − ∂x2

−

¯j ∂h − ∂t

j 2¯j ) h− = + (λj− v−

j 0 ¯j ) (x) ∂ h 2(S− − j ∂x S (x)

µ

j 00 j j 0 ) (x) (x)(S− ) (x))2 −S− 2((S−

+

j 2 (S− ) (x)

−

h

j j j j g¯ (x, t) + S− (x) λj+ v− +λj− v−

i

∂¯ gj (x, t) ∂t

¶

¯ j− h

¸

(2.43)

+ F−j (x, t) + F j (x, t, ²) ,

donde ·

F+j (x, t)

=

¯0 ∂2h + ∂t2

−

j 2 (v+ )

+ ·

¯ j ∂h + 2λj+ v+ ∂t

¯0 ∂2h + ∂x2

0

−

+

j 2¯0 (λj+ v+ ) h+

j 0 ¯0 2(S+ ) (x) ∂ h + j ∂x S (x) +

60

µ

+

¸

j 0 j j 00 2((S+ ) (x))2 −S+ (x)(S+ ) (x) j 2 (S+ ) (x)

¶

¸

¯0 , h + (2.44)

·

F−j (x, t) =

¯0 ∂2h − ∂t2

−

j 2 (v− )

j + 2λj− v−

·

¯0 ∂2h − ∂x2

−

¯0 ∂h − ∂t

j 2¯0 + (λj− v− ) h−

j 0 ¯0 2(S− ) (x) ∂ h − j ∂x S (x)

µ

+

¸

j 0 j j 00 2((S− ) (x))2 −S− (x)(S− ) (x) j 2 (S− ) (x)

−

·

j

j j j j −λ− ² j F j (x, t, ²) = − 21 λj− v− S− (x) v− λ− e g¯ (x − ², t − j

j

g + e−λ− ² ∂¯ (x − ², t − ∂t

² j ) v−

j

j −λ− ² j g¯ (x + ², t − e + λj− v−

+

j gj e−λ− ² ∂¯ (x ∂t

+ ², t −

² j ) v−

j

¶

¸

¯0 , h − (2.45)

² j ) v−

j

j −λ− ² ∂¯ g − v− e (x − ², t − ∂x

² j ) v−

² j ) v−

+

j v−

j gj e−λ− ² ∂¯ (x ∂x

¸

+ ², t −

² j ) v−

.

(2.46) Con la finalidad de presentar concisamente y más adecuadamente para nuestros cálculos el sistema (2.42), (2.43), definiremos 4N operadores diferenciales parciales, actuando sobre las funciones de dos variables con derivadas parciales hasta de segundo orden continuas. Denotaremos estos operadores an definidos como sigue con los s´ımbolos Li,j ± , 1 ≤ i ≤ 2, 1 ≤ j ≤ N y estar´ L1,j ± [u] = ·

L2,j ± [u]

=

j 2 ) (v±

∂2u ∂x2

−

∂2u ∂t2

j ∂u + 2λj± v± + (λ± v± )2 u , ∂t

j 0 2(S± ) (x) ∂u j S (x) ∂x ±

µ

+

j 0 j j 00 2((S± ) (x))2 −S± (x)(S± ) (x) j 2 (S± ) (x)

(2.47) ¶ ¸

u ,

(2.48) donde u = u(x, t) es una función con derivadas parciales hasta de segundo orden continuas.

61

Con esta nueva notación, el sistema (2.42), (2.43) se convierte en el sistema 2,j ¯ j ¯j L1,j + [h+ ] = L+ [h+ ] h

j j j j + λj+ v+ S+ (x) λj+ v+ g¯ (x, t) +

i

∂¯ gj (x, t) ∂t

(2.49)

+ F+j (x, t) , 2,j ¯ j ¯j L1,j − [h− ] = L− [h− ]

h

j j j j + λj− v− S− (x) λj− v− g¯ (x, t) +

i

∂¯ gj (x, t) ∂t

(2.50)

+ F−j (x, t) + F (x, t, ²) . Para reducir el sistema a un sistema equivalente más sencillo de estudiar haremos uso de los siguientes resultados los cuales son muy fáciles de probar ¯j Lema 1 Sean Li,j ± los operadores diferenciales definidos en (2.48) y sean L± los operadores diferenciales parciales definidos por ¯ j± [u] = L

∂2u ∂t2

2

j ∂ u j ∂u j 2 ) u − v± + 2λj± v± + (λj± v± . ∂t ∂x2

(2.51)

Entonces ¯ j± [¯ uj± ] = L

1

j S± (x)

h

i

2,j j j L1,j ± [u± ] − L± [u± ] ,

(2.52)

donde u¯j± (x, t) =

uj± (x, t) . S±j (x)

¯ j± como en el Lema 1, entonces Lema 2 Sea L j

j

j

j

¯ j+ [uj+ ] = ∂+2,j [¯ eλ+ v+ t L uj+ ] , ¯ j− [uj− ] = ∂−2,j [¯ eλ− v− t L uj− ] , 62

(2.53) (2.54)

donde j

j

u¯j± (x, t) = eλ± v± t uj± (x, t) y ∂+2,j [u] =

2 ∂ 2u j 2∂ u − (v ) , + ∂t2 ∂x2

(2.55)

∂−2,j [u] =

2 ∂ 2u j 2∂ u − (v ) . − ∂t2 ∂x2

(2.56)

Aplicando el Lema 1 al sistema (2.49), (2.50) ese sistema es equivalente al sistema "

¯ j+ [H+j ] L

"

¯ j− [H−j ] L

#

∂¯ gj 1 j j j g¯ (x, t) + λj+ v+ = j F+j (x, t) + λj+ v+ (x, t) , ∂t S+ (x)

(2.57)

#

∂¯ gj 1 1 j j j g¯ (x, t) + λj− v− = j F−j (x, t)+λj− v− , (x, t) +F j (x, t, ²) j ∂t S− (x) S− (x) (2.58)

donde

H+j (x, t)

¯ j± (x, t) h = j S± (x)

y aplicando el Lema 2 al sistema (2.57), (2.58), este sistema es equivalente a "

¯ +j ] ∂+2,j [H

j λj+ v+ t

=e

"

¯ −j ] ∂−2,j [H

j λj− v− t

=e

##

"

F+j (x, t) ∂¯ gj j j j j j (x, t) + λ v λ v g ¯ (x, t) + + + + + ∂t S+j (x) #

"

,

(2.59) #

F j (x, t, ²) ∂¯ gj F−j (x, t) j j j j j (x, t) + + λ , v λ v g ¯ (x, t) + − − − − ∂t S−j (x) S−j (x) (2.60) 63

donde ¯ ±j (x, t) = eλj± v±j t H±j (x, t) . H El sistema (2.57), (2.58) es un sistema de 4N ecuaciones diferenciales par¯ i PN hi ), ¯ j (PN h ciales no lineales acopladas (recordemos que g¯j (x, t) = G i=1 + , i=1 − ¯ j es que sea una función donde lo u ńico que le hemos pedido hasta ahora a G ¯ i ), además, la ecuación (2.59) desigmoide de cada una de las variables h ± pende de un peque˜ no parámetro ² > 0. Este sistema es sumamente complejo de estudiar, as´ı que para hacer más tratable el sistema, estudiaremos primero el caso N = 1, esto significa que se tiene u ńicamente una sola capa en cada columna. En este caso el sistema a estudiar es "

#

g ¯ + [H+ ] = 1 F+ (x, t) + λ+ v+ λ+ v+ g¯(x, t) + ∂¯ L (x, t) , S+ (x) ∂t "

¯ − [H− ] = L

(2.61)

#

1 ∂¯ g 1 F− (x, t) + λ− v− λ− v− g¯(x, t) + (x, t) + F (x, t, ²) . S− (x) ∂t S− (x) (2.62)

Es decir, tenemos dos ecuaciones diferenciales parciales con la ecuación (2.62) dependiendo de un parámetro peque˜ no ² > 0, la función g¯ depende ¯ ± , de manera sigmoide en cada una de ellas u ńicamente de h ¯ + , h− ) . g¯(x, t) = G(h

(2.63)

En la siguiente sección estudiaremos el sistema (2.61), (2.62) usando técnicas de parámetro peque˜ no, es decir, estudiaremos primeramente el sistema degenerado, − el sistema (2.61), (2.62) con ² = 0 − , para después, con base en este estudio, obtener información sobre el sistema no degenerado. ¯ que dependerá de una sola Además introduciremos una expresión para G, ¯ ¯ − , a esta nueva variable variable, la cual es una combinación lineal de h+ , y h la llamaremos variable de activación, ya que el hecho de que se generen o no potenciales de acción en una columna dependen de esta nueva variable.

64

´ DE ACTIVACION ´ PARA LA CORTEZA 2.7 LA ECUACION CEREBRAL Sabemos que g¯(x, t) depende de una manera muy compleja de la actividad ¯ + , y de la actividad sináptica inhibitoria h ¯ − , i.e. g¯(x, t) sináptica excitatoria h ¯ ¯ ¯ − ) debe ser ¯ ¯ es igual a G(h+ , h− ), cualitativamente consideramos que G(., h ¯ + , para cada h ¯ − fijo y que debe ser una una función sigmoide creciente de h ¯ − para cada h ¯ + fijo. Este planteamiento función sigmoide decreciente de h es sumamente general y dif´ıcil de estudiar, as´ı que estudiaremos un caso ¯+ y h ¯ − . El caso que particular de la dependencia de g¯ con respecto a h estudiaremos está basado en lo siguiente: 1. Si V± representa el valor del potencial postsináptico en una sinapsis, ¯ ± (x, t) representará la densidad de potencial postsináptico entonces V± h en la columna x en el tiempo t. 2. Habrá activación en la columna x (generación de potenciales de acción) cuando V+ H+ (x, t) − V− H− (x, t) ≥ V0 , donde V0 es un cierto umbral promediado de disparo de la columna x, en caso contrario no habrá activación. 3. Si V+ H+ (x, t) − V− H− (x, t) es muy grande, habrá saturación de la columna x en el tiempo t y entonces, no aumentará significativamente la actividad de la columna en ese tiempo. Basados en estos tres puntos consideramos que: ¯ + H+ (x, t) − V− H− (x, t)) , g¯(x, t) = G(V

(2.64)

¯ es una función sigmoide de una sola variable. G ¯ es una función donde G sigmoide de la densidad de potencial postsináptico de la columna. En realidad, como veremos después, se pueden elegir varias expresiones ¯ cuantitativas para G. Llamaremos a la variable u(x, t) = V+ H+ (x, t) − V− H− (x, t) variable de activación de la columna x.

65

La dinámica espacio-temporal de u(x, t) nos permitirá estudiar las zonas de la corteza que están activas en un tiempo determinado. Para cada t, las zonas activas de la corteza estarán caracterizadas por el conjunto E(t) = {x en la corteza|u(x, t) ≥ V0 } . Para obtener la ecuación diferencial parcial que determina la dinámica ¯ + [H+ ] = de la variable de activación u(x, t) multiplicamos la ecuación L ∂¯ g 1 ¯ + [H− ] F (x, t) + λ+ v+ [λ+ v+ g¯(x, t) + ∂t (x, t)] por V+ y le restamos V− L S+ (x) + en ambos lados de la ecuación, as´ı se obtiene "

¯ + [u] = V+ λ+ v+ L

#

∂¯ g F+ (x, t) ¯ + [H− ] . (2.65) λ+ v+ g¯(x, t) + (x, t) +V+ −V− L ∂t S+ (x)

Además tenemos la ecuación para las inhibiciones "

#

g F− (x, t) F (x, t, ²) ¯ − [H− ] = λ− v− λ− v− g¯(x, t) + ∂¯ L (x, t) + + . ∂t S− (x) S− (x)

(2.66)

Como mencionamos en la sección anterior, primeramente estudiaremos el caso degenerado (² = 0) del sistema (2.65), (2.66) "

¯ + [u] = V+ λ+ v+ L

#

∂¯ g F+ (x, t) ¯ + [H− ] . (2.67) λ+ v+ g¯(x, t) + (x, t) +V+ −V− L ∂t S+ (x)

Además tenemos la ecuación para las inhibiciones ¯ − [H− ] = F− (x, t) . L S− (x)

(2.68)

Usando el lema 2, se puede reducir el sistema (2.67), (2.68) al estudio del sistema equivalente "

¯ + [u] = V+ λ+ v+ L

#

∂¯ g ∂u F+ (x, t) ¯ + [H− ] , (2.69) λ+ v+ g¯(u) + (u) + V+ − V− L ∂u ∂t S+ (x)

66

¯ − ) = e−λ− v− t F− (x, t) , ∂−2 (H S− (x)

(2.70)

donde ¯ H(x, t) = eλ− v− t H− (x, t) . ¯ ∂ 2 y calculando L ¯ + [e−λ− v− t H ¯ − (x, t)] y reUtilizando la definición de L, − ordenando los términos, obtenemos el sistema ∂2u ∂t2

h

2

2 ∂ u + 2λ+ v+ ∂u + (λ+ v+ )2 u − v+ = V+ λ+ v+ λ+ v+ g¯(u) + g¯0 (u) ∂u ∂t ∂x2 ∂t

·

−V− e−λ− v− t (1 − + S+V+(x) F+ (x, t) +

2 ¯− v+ ∂2H 2 ) ∂t2 v−

¯

i

− ¯− + 2(λ+ v+ − λ− v− ) ∂ H + (λ+ v+ − λ− v− )2 H ∂t

1 F− (x,t) −λ− v− t 2 S (x) e v− −

¸

, (2.71)

2 ¯ ¯− ∂2H e−λ− v− t 2 ∂ H− = v + F− (x, t) . − ∂t2 ∂x2 S− (x)

(2.72)

Nótese que la ecuación (2.72) es independiente de g¯, as´ı que este sistema se puede interpretar como un proceso global de propagación de ondas de la variable de activación, modulado por un proceso de propagación de ondas de las inhibiciones.

67

´ DE UN MODELO 2.8 CONSTRUCCION ´ DE LA COMPARTIMENTADO DE ACTIVACION CORTEZA CEREBRAL. En el modelo que hemos analizado en las secciones anteriores –con las diversas formas propuestas para la distribución de conexiones– se eliminan los efectos producidos en cada columna por columnas lejanas a ella. La eliminación de los “efectos lejanos”, se debe a la elección de una distribución de conexiones que decrece exponencialmente con las distancias entre columnas con la misma escala de decrecimiento, en todas las direcciones sobre la corteza. Una manera de tomar en cuenta estos “efectos lejanos”, es considerar a la corteza cerebral dividida en n compartimentos disjuntos interconectados, considerando la distribución de conexiones para cada par de compartimentos de decrecimiento exponencial, pero con una escala de decrecimiento que depende del par de compartimentos considerados. As´ı, “los efectos” producidos en una columna de cierto compartimento estarán clasificados dependiendo de los compartimentos con que interact´ ua. En nuestro modelo consideramos que las velocidades de propagación de los potenciales de acción no dependen de los compartimentos que se están considerando. La división de la corteza en n compartimentos, además de ser una mejor aproximación de la corteza real que la propuesta en las secciones anteriores, nos permite correlacionar nuestros resultados teóricos con las mediciones reales de actividad cortical. Además nos permitirá estudiar el diagrama de conexiones funcional que se establece en un estado fisiológico dado. Se podr´ıa estudiar, también la interrelación que existe entre las distintas áreas de la corteza, siguiendo la división propuesta por el anatomista K. Brodmann (Fig. # 9). A continuación describiremos más detalladamente este modelo. Como ya hemos mencionado, en el presente modelo consideramos a la corteza cerebral dividida en n compartimentos ajenos Ωi , i = 1, . . . , n. Por simplicidad y siguiendo las ideas de las secciones anteriores, consideraremos una banda de corteza en la cual, las ondas que se propagan por ella, tienen longitud de onda muy peque˜ na con respecto a la longitud de la misma. Con estas simplificaciones puede pensarse en la corteza como una recta dividida en n intervalos disjuntos, identificando los dos intervalos infinitos como uno solo (considerando que los “efectos ” de los dos intervalos infinitos provocan 68

el mismo efecto sobre los otros intervalos). Lo que deseamos es relacionar la densidad de actividad sináptica – excitatoria e inhibitoria – en una columna de corteza con los potenciales de acción que se generan en las otras columnas en tiempos anteriores, pero como tenemos dividida a la corteza en n compartimentos, aumentará el n´ umero de ecuaciones a considerar. A continuación procederemos a definir las variables y parámetros que intervienen en el modelo, deduciremos las ecuaciones en forma integral y al igual que en las secciones anteriores transformaremos éstas en un sistema de ecuaciones diferenciales parciales. Como en las secciones anteriores, el sub´ındice + se referirá a las sinapsis excitatorias y el − a las sinapsis inhibitorias. Sean : a) hij aptica por ± (x, t) − La densidad de sinapsis activas o de actividad sin´ cm2 en la columna x ∈ Ωj , en el tiempo t, producidos por los potenciales de acción que se generan en el compartimento Ωi en tiempos anteriores hij ± (x, t) representa el efecto de Ωi en la columna x ∈ Ωj en el instante t (1/cm2 ). b) f±j (x, t) − La densidad de sinapsis activas por cm2 en la columna x ∈ Ωj en el instante t, producidas por los potenciales de acción que se generan en tiempos “anteriores” en las estructuras internas del cerebro, principalmente el tálamo, y que tienen conexión sináptica con la columna x ∈ Ωj . Supondremos que f−j (x, t) ≡ 0 para cada j. c) hj± (x, t) =

n X ij

h± (x, t) + f±j (x, t) .

(2.73)

i=1

Representará la densidad de actividad sináptica total en la columna x ∈ Ωj en el instante t, producida por los potenciales de acción generados en todas las regiones de la corteza y las estructuras internas del cerebro, principalmente el tálamo. ij d) R± umero de conexiones sinápticas por cm2 en la columna (x0 , x, v) − El n´ x ∈ Ωj por unidad de longitud en x0 ∈ Ωi y por unidad de velocidad, es ij decir, R± (x0 , x, v)dx0 dv representa el n´ umero de conexiones sinápticas

69

por cm2 en la columna x ∈ Ωj que provienen de las neuronas que están ubicadas desde la columna x0 hasta la columna x0 + dx0 y por cuyos axones los potenciales de acción “viajan” con velocidad entre v y v + dv. e) g(x, t) − La fracción del total de neuronas en la columna x que generan potenciales de acción en el instante t. f ) g j (x, t) − La fracción del total de neuronas en la columna x ∈ Ωj que generan potenciales de acción en el tiempo t. Obviamente g(x, t) = g j (x, t) para x ∈ Ωj .

(2.74)

Lo que queremos conocer es la evolución espacio-temporal de la densidad de actividad sináptica total en cada columna de cada compartimento de la corteza (hij ± (x, t), x ∈ Ωj , j = 1, . . . , n). Para lograr esto, primero obtendremos la evolución espacio-temporal de la densidad de actividad sináptica en cada columna de un compartimento dado de la corteza producida por la interacción con otro compartimento de la misma (hij ± (x, t), x ∈ Ωj , i, j = 1, . . . , n) y le a˜ nadiremos la producida por las estructuras internas del cerebro (f+j (x, t), j = 1, . . . , n). Consideraremos que por cada par de compartimentos Ωi , Ωj se cumplen las siguientes hipótesis : 1. Existe una u ńica velocidad de propagación, v± , de los potenciales de acción que viajan de Ωi a Ωj y que no depende de Ωi ni de Ωj . ij (x0 , x, v)) depende u ńicamente de la dis2. La densidad de conexiones (R± tancia entre las columnas y tienen su decrecimiento exponencial, con olo depende de las regiones cada escala de decrecimiento (λij ± ) que s´ Ωi , Ωj .

R± (x0 , x, v) =

n X

ij R± (x0 , x, v),

i=1

donde

70

x ∈ Ωj ,

(2.75)

ij R± (x0 , x, v)

=

 ij ij 0 −λ |x−x0 |  Sj (x)δ(v − v± ) ,  α± (x )e ±  

0,

si x ∈ Ωj , x0 ∈ Ωi en otro caso (2.76)

0

R± (x , x, v) representará el n´ umero de conexiones sinápticas por cm2 en x, por unidad de longitud en x0 y por unidad de velocidad. Sj (x) es la ij densidad de conexiones sinápticas en x ∈ Ωj (1/cm2 ) y α± (x0 ) (1/cm) es un factor de normalización, sujeto a la condición n Z X i=1 R

ij

0

ij α± (x0 )e−λ± |x−x | dx0 = 1, para cada x ∈ Ωj .

(2.77)

En lo que sigue, identificaremos Ωi (i = 1, . . . , n) con el conjunto de ij ij ij , (x0 ) = α± (x0 ) constante; α± los n´ umeros reales y si suponemos α± entonces de (2.77) obtenemos : 2

ij n X α± i=1

λij ±

= 1.

(2.78)

Debido al carácter local de las inhibiciones, podemos suponer que ij = 0 si i 6= j. α−

As´ı, en lo que sigue consideraremos que : ij = α+

λij λjj + jj = − . , α− 2n 2

(2.79)

3. Los potenciales de acción que se generan en la columna x ∈ Ωj en el instante t, dependen de la actividad sináptica total en x ∈ Ωi en el instante t y de la densidad de la actividad neuroqu´ımica (concentración de neurotransmisores) que se realice en dicha columna en ese instante. Nosotros sólo consideraremos la dependencia de g i con respecto a la actividad total en la columna x en el instante t. Además, de la suposición de que la actividad inhibitoria en la corteza sólo tiene un carácter local, se concluye que hki − (x, t) ≡ 0 si k 6= i. 71

Finalmente como las conexiones aferentes a la corteza son sólo excitatorias, tenemos que f−j (x, t) ≡ 0 y as´ı obtenemos g i (x, t) = g i (

n X

i ii hki + (x, t) + f+ (x, t), h− (x, t)) .

(2.80)

k=1

Siguiendo el mismo razonamiento expuesto en las secciones anteriores, tenemos que :

hij − (x, t)

Z ∞ jj λij |x − x0 | 0 + e−λ+ |x−x | g(x0 , t − )dx0 , = Sj (x) 2n v+ −∞

hjj − (x, t) =

Z x+² jj λjj |x − x0 | 0 − Sj (x) )dx0 , e−λ+ |x−x | g(x0 , t − 2 v− x−²

(2.81)

(2.82)

donde el parámetro ² > 0 corresponde a la localidad de la actividad inhibitoria. De nuevo, debido al carácter local del proceso inhibitorio, consideraremos que hij − (x, t) ≡ 0 si i 6= j. Introduciendo las nuevas funciones incógnitas H+ij (x, t)

nhij (x, t) hjj + (x, t) jj = , H− (x, t) = − Sj (x) Sj (x)

(2.83)

2 y realizando el mismo proceso de aplicar el operador ∂ 2 a cada lado de ∂x (2.81), (2.82) y después pasar al l´ımite cuando ² → 0, se obtiene, debido a (2.74) 2 ij ∂H+ij ∂ 2 H+ij ij ij ij 2 2 2 ∂ H+ v ) v H − v + 2λ + (λ + + + + + + ∂t2 ∂t ∂x2

"

=

λij + v+

#

ij λij + v+ g (x, t)

∂g ij (x, t) , + ∂t

x ∈ Ωj ,

2 jj ∂ 2 H−jj ∂H−jj jj jj 2 2 jj 2 ∂ H− + 2λ− v− + (λ− ) v− H− − v− = 0, ∂t2 ∂t ∂x2

72

x ∈ Ωj ,

(2.84)

(2.85)

donde g ij (x, t) es la fracción del total de neuronas que generan potenciales de acción en x ∈ Ωj debido a las sinapsis que provienen de neuronas en Ωi . En lo que sigue supondremos que g ij (x, t) = µij (x)g j (x, t) ,

(2.86)

donde g j (x, t) está definida en (2.74) y 0 ≤ µij (x) ≤ 1 es tal que n X

µi,j (x) = 1, para toda x ∈ Ωj .

i=1

Las funciones µij (x); i, j = 1, . . . , n determinan el diagrama de conexión funcional en la corteza, correspondiente a un estado fisiológico, caracterizado por una función de activación de neuronas g j (x, t) para cada compartimento Ωj ; j = 1, . . . , n. El hecho de que µij (x) pueda ser diferente de µkj para i 6= k nos permite considerar una distribución no homógenea de activación en x ∈ Ωj correspondiente a los compartimentos Ωk y Ωi que no pod´ıamos considerar en el modelo estudiado en las secciones anteriores. Sea V± el valor absoluto del potencial postsináptico producido por una sinapsis (+ excitatoria, − inhibitoria). Introduciremos las nuevas variables h

uij (x, t) =

V+ (n − 1) H+ij (x, t) +

f j (x,t) Sj (x)

i

− V− H−jj (x, t)

(V+ + V− )

,

(2.87)

i, j = 1, . . . , n , que llamaremos variable de activación parcial en x ∈ Ωj , en el instante t, provocada por la acción de Ωi sobre Ωj . De (2.84), (2.85), (2.86), (2.87) se obtiene el siguiente sistema de ecuaciones para las variables de activación : 2 ∂ 2 uij ∂uij ij ij 2 2 ∂ uij v v ) u − v + 2λ + (λ = + + ij + + + ∂t2 ∂t ∂x2

"

λij ∂g j + v+ V+ (n − 1) ij j µij (x) λ+ v+ g (x, t) + (x, t) (V+ + V− ) ∂t

#

+ Fij (x, t) − Iij (x, t), x ∈ Ωj ; i, j = 1, . . . , n , 73

(2.88)

donde 

³

´

j

2 f (x,t) ∂ V+ (n − 1)  ∂ Sj (x) ij Fij (x, t) = v + 2λ  + + (V+ + V− ) ∂t2

2 2 +(λij + ) v+

j

f (x, t) 2 − v+ Sj (x)

∂2

³

f j (x,t) Sj (x) ∂x2

³

f j (x,t) Sj (x)

´

∂t

´   ,

"

(2.89)

#

2 jj V− ∂ 2 H−ij ∂H−ij ij ij 2 jj 2 ∂ H− Iij (x, t) = + 2λ . v + (λ v ) H − v + + + + − + V+ + V− ∂t2 ∂t ∂x2 (2.90) Nótese que, en (2.90), H−jj (x, t) sólo depende de las condiciones iniciales de inhibición, debido a que satisface la ecuación (2.85). Teniendo en cuenta (2.80) podemos suponer que en (2.88) se tiene

g j (x, t) = gj (uj ) , donde uj (x, t) =

n X

(2.91)

uij (x, t)

i=1

es la llamada variable de activación total en x ∈ Ωj . De Ã ! n X f j (x, t) ij 0≤ ≤ 1, para cada x ∈ Ωj H+ (x, t) + Sj (x) i=1 y 0 ≤ H−jj (x, t) ≤ 1, para cada x ∈ Ωj , se tiene que −1 ≤ uj (x, t) ≤ 1, para cada x ∈ Ωj . Entonces de (2.91) y el comportamiento cualitativo que, a partir del conocimiento fisiológico, debe tener la función de activación gj , podemos suponer que gj (uj ) es una función creciente continuamente diferenciable de la variable uj , definida en todo R y tal que

74

gj =

   0,  

creciente , 1,

si uj ≤ −1 , si − 1 ≤ uj ≤ 1 , si uj ≥ 1 .

(2.92)

Para tener más arbitrariedad en la selección de gj , podemos suponer que es una función sigmoide de la variable uj que está muy cercana al valor 0 en uj = −1 y muy cercana al valor 1 en uj = 1. Haremos ahora la suposición de que las conexiones excitatorias entre compartimentos diferentes son pocas en comparación con las que hay dentro de cada compartimento, es decir, λjj + ¿ 1, si i 6= j; i, j = 1, . . . , n . λij +

(2.93)

Entonces, si ponemos λjj λjj + v+ + = βj , = ², jj λ− v− λij +

i 6= j

(2.94)

y pasamos a las nuevas coordenadas jj λjj + v+ t = τ, λ+ x = y

(2.95)

en (2.88) para todo i, j = 1, . . . , n, obtenemos un sistema con parámetro peque˜ no ² para las variables de activación parcial en coordenadas adimensionales ²2

2 ∂ 2 u˜ij ∂ u˜ij ˜ij 2∂ u + 2² + u ˜ − ² = ij 2 ∂τ ∂τ ∂y 2

"

#

V+ (n − 1) ∂˜ gj µ ˜ij (y) g˜j (y, τ ) + ² (y, τ ) + (V+ + V− ) ∂τ (2.96)

+

h i 1 ˜ij (y, τ ) − I˜ij (y, τ ) ; F 2 (λij + v+ )

75

1 y ∈ Ωj , i, j = 1, . . . , n, λjj +

i 6= j ,

donde u˜ij (y, τ ) = uij (x, t), µ ˜ij (y) = µij (x), g˜j (y, τ ) = g j (x, t), F˜ij (y, τ ) = Fij (x, τ ), I˜ij (y, τ ) = Iij (x, t). Además se tiene ∂ 2 u˜jj ∂ u˜jj ∂ 2 u˜jj + 2 + u ˜ − = jj ∂τ 2 ∂τ ∂y 2 "

#

V+ (n − 1) ∂˜ g µ ˜jj (y) g˜j (y, τ ) + (y, τ ) + (V+ + V− ) ∂τ (2.97) i h 1 ˜jj (y, τ ) − I˜jj (y, τ ) , 1 y ∈ Ωj . F (λ+ v+ )2 λjj +

De (2.91) tenemos que g˜j (y, τ ) = g˜j (˜ uj ) , donde u˜j (y, τ ) =

n X

u˜ij (y, τ )

i=1

es la variable de activación total de la región Ωj en coordenadas adimensionales. 76

En lo que sigue consideraremos que las variables de activación se obtienen en forma de promedios espaciales en cada compartimento Ωj , es decir u˜ij (y, τ ) = u˜ij (τ ), u˜ij (y) = u˜ij , g˜j (y, τ ) = g˜j (τ ), F˜ij (y, τ ) = F˜ij (τ ) , I˜ij (y, τ ) = I˜ij (τ ) . De esta forma, el sistema (2.96), (2.97) se reduce al sistema ²2 "

d2 u˜ij d˜ uij + 2² + u˜ij = dτ 2 dτ #

i h d˜ gj 1 V+ (n − 1) j ˜ij (τ ) − I˜ij (τ ) , g˜ (τ ) + ² F µ ˜ij (τ ) + ij (V+ + V− ) dτ (λ+ v+ )2

i, j = 1, . . . , n,

(2.98)

i 6= j ,

"

#

∂˜ gj d2 u˜jj d˜ ujj V+ (n − 1) j g ˜ (τ ) + + 2 + u ˜ = µ ˜ (τ ) jj ij dτ 2 dτ (V+ + V− ) ∂τ

(2.99)

+

i h 1 ˜jj (τ ) − I˜jj (τ ) . F 2 (λjj + v+ )

Por simplicidad supondremos que la entrada aferente Fij (τ ) a la región Ωj proveniente de estructuras internas del cerebro es peque˜ na, es decir : ˜ F˜ij (τ ) = (λjj + v+ )Gij (τ ), i 6= j , ˜ ij (τ ) está acotada. donde G 77

(2.100)

Por otra parte, de (2.90) se tiene que "

#

2 ˜ jj ˜ −jj (τ ) 1 V− dH 2 d H− (τ ) ˜ij (τ ) = I ² + 2² + H−jj (τ ) , (2.101) 2 2 (V + V ) dτ dτ (λij v ) + − + +

˜ −jj (τ ) = H−jj (τ ) es, seg´ donde H un (2.85), solución de la ecuación 2 ˜ jj 2 d H− βj 2

˜ −jj dH ˜ −jj = 0 , + 2βj +H dτ dτ donde βj viene definido en (2.94). De (2.102) se tiene que

(2.102)

˜ −jj (τ ) = (aτ + b)e−τ /βj , H donde a y b son parámetros que dependen de las condiciones iniciales de inhibición. Finalmente, la ecuación (2.98) se transforma en ²2

d2 u˜ij d˜ uij + 2² + u˜ij = 2 dτ dτ #

"

V+ (n − 1) d˜ gj ˜ ij (τ ) µ ˜ij g˜j (τ ) + ² (τ ) + ²2 G (V+ + V− ) dτ −

(2.103)

V− ˜ jj (τ ) para cada i, j = 1, . . . , n, i 6= j . H (V+ + V− ) −

En la aproximación de orden cero, para valores de τ no cercanos a 0, la solución u˜0ij (τ ) de (2.103) satisface la ecuación u˜0ij =

V+ (n − 1) V− ˜ −jj (τ ), i, j = 1, . . . , n, i 6= j , µ ˜ij g˜j (˜ u0j ) − H (V+ + V− ) (V+ + V− ) (2.104)

donde u˜0j (τ ) =

n X i=1

78

u˜0ij (τ ),

de (2.104) se concluye que

u˜0j =

n X

u˜0ij (τ ) =

i=1

n X V+ (n − 1) (n − 1)V− ˜ jj µ ˜ij + u˜0jj − g˜j (˜ u0j ) H− (τ ) , (V+ + V− ) (V + V ) + − i=1,i6=j (2.105)

de donde se tiene : u˜0j =

V+ (n − 1) (n − 1)V− ˜ jj (1 − µ ˜jj )˜ g (˜ u0j ) − u˜0jj + H (τ ) . (V+ + V− ) (V+ + V− ) −

(2.106)

Si suponemos que gj es una función sigmoide con un solo punto de inflexión, entonces la ecuación (2.106) admite una solución u ńica u0j , para 0 cualquier u˜jj , si se cumple la condición V+ (n − 1) 0 (1 − µ ˜jj )˜ gj (u) < 1, para cada j = 1, . . . , n . (V+ + V− )

(2.107)

En este caso, aplicando al sistema (2.99), (2.103) los resultados de la teor´ıa de perturbaciones singulares para sistemas con parámetro peque˜ no y de la ecuación (2.106), tendremos que para cualquier compacto contenido en {τ > 0} la aproximación de orden cero de la variable de activación total, satisface la ecuación d˜ u0j d2 u˜0j +2 + u˜0j = 2 dτ dτ "

#

d˜ u0 V+ (n − 1) 0 µ ˜jj g˜j (˜ u0j ) + g˜j (˜ u0j ) j + (V+ + V− ) dτ #

"

d˜ gj 0 V+ (n − 1) d2 g˜j 0 u0j ) (1 − µ ˜jj ) (˜ uj ) + 2 (˜ u ) + g˜j (˜ 2 (V+ + V− ) dτ dτ j "

˜ −jj ˜ −jj dH (n − 1)V− d2 H ˜ −jj + 2 +H − (V+ + V− ) dτ 2 dτ 79

#

+

1 2 (λjj + v+ )

h

i

F˜jj (τ ) − I˜jj (τ ) .

˜ −jj y que Teniendo en cuenta la expresión de H 1 V− I˜jj (τ ) = jj 2 (V+ + V− ) (λ+ v+ ) =

Ã

˜ −jj ˜ jj d2 H dH ˜ −jj +2 − +H dτ dτ

!

V− e−τ /βj (βj − 1) [2aβj + (βj − 1)(aτ + b)] . (V+ + V− )βj2

Obtenemos finalmente la ecuación "

#

d2 u˜0j V+ (n − 1) 0 (1 − µ ˜jj )gj (˜ u0j ) + 1− (V+ + V− ) dτ 2

"

#

d˜ u0j V+ (n − 1) 0 2− (2 − µ ˜jj )gj (˜ u0j ) + (V+ + V− ) dτ Ã

d˜ u0j V+ (n − 1) 00 − (1 − µ ˜jj )gj (˜ u0j ) (V+ + V− ) dτ u˜0j 1

˜

−

+

(2.108)

#

"

Fjj (τ ) 2 (λjj + v+ )

!2

V+ (n − 1) g˜j (˜ u0j ) = − (V+ + V− )

nV− (βj − 1)e−τ /βj [2aβj + (βj − 1)(aτ + b)] , (V+ + V− )βj2

donde el primer sumando del miembro derecho está asociado a la entrada aferente a la región Ωj proveniente de estructuras interiores del cerebro y el segundo sumando corresponde al control inhibitorio en Ωj . La ecuación (2.108) describe, en la aproximación de orden cero con respecto a las potencias de ², la evolución temporal de la variable de activación total en cada compartimento Ωj . 80

Nótese que el coeficiente (1 − µ ˜jj ) corresponde a la contribución de los compartimentos Ωj con i 6= j a la activación de Ωj . Nótese además que bajo las suposiciones efectuadas se obtienen ecuaciones desacopladas para la variable de activación total en cada Ωj , aunque como se˜ nalamos en el párrafo anterior el efecto de los restantes compartimentos está encerrado en el coeficiente (1 − µ ˜jj ). Finalmente hemos obtenido el teorema siguiente. Teorema 1 Supongamos que la corteza cerebral está dividida en n compartimentos mutuamente ajenos Ωj , j = 1, . . . , n y que se cumplen las siguientes hipótesis: 1. La función de distribución de densidad de conexiones entre dos localizaciones de la corteza está descrita por (2.75), (2.76), (2.79). 2. La función que describe la fracci´ on del total de neuronas en una localizaci´ on determinada que generan potenciales de acci´ on debido a las sinapsis que provienen de neuronas ubicadas en otros compartimentos viene dada por las expresiones (2.86), (2.91) y (2.92). 3. Las conexiones excitatorias entre compartimentos diferentes son pocos en comparaci´ on con los que hay dentro de cada compartimento, seg´ un se expresa en (2.93). 4. Las entradas sinápticas aferentes a la corteza son peque˜ nas en el sentido de (2.100). 5. Entonces para la aproximaci´ on de orden cero de los promedios espaciales en cada Ωj de la variable de activación total u0j (τ )

=

n X

u0ij (τ )

(ver 2.87),

i=1

se tiene la ecuaci´ on (2.108), si se cumple la condici´ on (2.107).

81

Cap´ıtulo 3 ´ ´ ANALISIS MATEMATICO DEL MODELO ´ DE ACTIVACION DE LA CORTEZA CEREBRAL. ´ 3.1 UN MODELO MATEMATICO PARA UNA REBANADA DE CORTEZA CEREBRAL Un problema importante en fisiolog´ıa es determinar bajo qué est´ımulos externos se puede lograr que una peque˜ na rebanada de corteza cerebral oscile en frecuencias determinadas de antemano. Con el modelo que hemos obtenido es posible estudiar este problema. Como estar´ıamos considerando u ńicamente una peque˜ na rebanada de corteza cerebral, lo que tendr´ıamos que hacer en nuestro modelo, es considerar que no se tiene variación espacial; en este caso, nuestro modelo ser´ıa un sistema de dos ecuaciones diferenciales ordinarias no lineales. Lo que nos proponemos estudiar, es las condiciones que debe cumplir el est´ımulo externo para que nuestro sistema de ecuaciones tenga soluciones periódicas. As´ı que, el sistema de ecuaciones diferenciales ordinarias a estudiar es

82

d2 u dt2

h

+ 2λ+ v+ du + (λ+ v+ )2 u = V+ λ+ v+ λ+ v+ g¯(u) + g¯0 (u) du dt dt ·

−V− e−λ− v− t (1 −

2 ¯− v+ d2 H ) 2 dt2 v−

i

¯

− ¯− + 2(λ+ v+ − λ− v− ) dH + (λ+ v+ − λ− v− )2 H dt

¸

+f+ (t) , (3.1) ¯− d2 H = f− (t) , (3.2) dt2 donde f+ (t), f− (t) son funciones que están relacionadas con los est´ımulos externos a la peque˜ na rebanada de corteza. La ecuación (3.2) es muy fácil de resolver, sus soluciones son de la forma ¯ − (t) = R(t) + at + b , H

(3.3)

donde a, b son constantes y R(t) =

Z t µZ s 0

0

¶

f− (z)dz ds .

(3.4)

Sustituyendo (3.3) en (3.1), obtenemos la ecuación diferencial ordinaria no lineal que necesitamos estudiar d2 u dt2

h

+ 2λ+ v+ du + (λ+ v+ )2 u = V+ λ+ v+ λ+ v+ g¯(u) + g¯0 (u) du dt dt ·

−V− e

−λ− v− t

(1 −

i

2 v+ d2 R 2 ) dt2 v−

i

(3.5)

+2(λ+ v+ − λ− v− )( dR + a) + (λ+ v+ − λ− v− )2 (R(t) + at + b) dt +f+ (t) .

A continuación realizaremos un cambio de variable que nos permitirá reducir la ecuación a una ecuación diferencial equivalente más sencilla de estudiar. 83

Lema 3 Sea τ = λ+ v+ t, entonces la ecuaci´ on (3.5) es equivalente a la ecuaci´ on diferencial d2 u dτ 2

h

+ 2 du + u = V+ g¯(u) + g¯0 (u) du dτ dτ

−V− e

−

λ− v− τ λ+ v+

·

(1 −

2 v+ d2 R 2 ) dτ 2 v−

i

+ 2(λ+ v+ − λ− v− )( (λ+1v+ )2 dR + dτ

) +(λ+ v+ − λ− v− )2 ( (λR(τ 2 + + v+ )

aτ (λ+ v+ )3

+

b ) (λ+ v+ )2

i

a ) (λ+ v+ )2

(3.6)

+f+ (τ ) . Demostraci´ on. Dividimos la ecuación (3.5) por (λ+ v+ )2 y usamos la regla de la cadena para obtener las derivadas necesarias de las funciones que aparecen en la ecuación con respecto a τ . En el siguiente resultado construiremos, a partir de la ecuación (3.6), una familia de ecuaciones diferenciales de primer orden dependiente de un parámetro. Esta familia de ecuaciones diferenciales de primer orden es equivalente a la ecuación diferencial de segundo orden (3.6), en el sentido que se enuncia en la siguiente proposición. Proposici´ on 1 Si u es una solución de la ecuaci´ on diferencial (3.6) con 0 0 condiciones iniciales u(0) = u0 y u (0) = u0 , entonces u es solución de la ecuaci´ on diferencial de primer orden du + u = V+ g¯(u) + e−τ (F+ (τ ) + K) , dτ

84

(3.7)

donde F+ (τ ) =

Ã

Rτ

−

−V− e

0

λ− v− ν λ+ v+

·

(1 −

2 v+ d2 R 2 ) dν 2 v−

+ 2(λ+ v+ − λ− v− )( (λ+1v+ )2 dR + dν +

(λ+ v+ − λ− v− )2 ( (λR(ν) 2 + + v+ )

a ) (λ+ v+ )2

aν (λ+ v+ )3

+

i

b ) (λ+ v+ )2

´

+ f+ (ν) eν dν (3.8)

y K = u00 + u0 − V+ g(u0 ) . Rec´ıprocamente, si u es una solución de (3.7) con u(0) = u0 , entonces, u es una solución de (3.6) con condiciones iniciales u(0) = u0 y u0 (0) = u00 . Demostraci´ on. Si multiplicamos la ecuación (3.6) por eτ , obtenemos eτ

h

d2 u dτ 2

i

h

i

+ 2 du + u = eτ V+ g¯(u) + g¯0 (u) du + dτ dτ

dF+ dτ

Como "

"

#

d2 u d τ du du e ( + u) = eτ +2 +u 2 dτ dτ dτ dτ

#

y V+

d τ du [e g¯(u)] = V+ [¯ g (u) + g¯0 (u) ] dτ dτ

(3.6) se transforma en "

#

d d τ du e ( + u) = [V+ eτ g¯(u) + F+ (τ )] . dτ dτ dτ Entonces 85

.

(3.9)

eτ (

du + u) = V+ eτ g¯(u) + F+ (τ ) + K , dτ

donde K = u00 + u0 + V+ g¯(u0 ) , entonces du + u = V+ g¯(u) + e−τ (F+ (τ ) + K) . dτ Es decir, u es solución de la ecuación (3.7). El rec´ıproco se demuestra de manera similar. Ahora consideraremos que la entrada externa inhibitoria a la corteza cerebral es cero, es decir, tomaremos f− (t) ≡ 0. En este caso la ecuación (3.6) se reduce a "

d2 u du du 0 + 2 + u = V g ¯ (u) + g ¯ (u) + dτ 2 dτ dτ

−V− e

−

λ− v− τ λ+ v+

"

#

Ã

a aτ b 2(λ+ v+ − λ− v− ) + (λ+ v+ − λ− v− )2 + 2 3 (λ+ v+ ) (λ+ v+ ) (λ+ v+ )2

+f+ (τ ) y

F+ (τ ) =

Rτ

Ã −

−V− e

0

Ã 2

+(λ+ v+ − λ− v− )

λ− v− ν λ+ v+

"

2(λ+ v+ − λ− v− )

b aν + 3 (λ+ v+ ) (λ+ v+ )2

!#

a (λ+ v+ )2 !

+ f+ (ν) eν dν

En el siguiente resultado, probaremos que si escogemos el est´ımulo externo de tal manera que e−τ (F+ (τ ) + K) sea una función periódica en τ de per´ıodo 86

!#

T , se pueden proponer condiciones sobre g¯ para que la ecuación diferencial (3.7) tenga una u ńica solución periódica, de per´ıodo T . Observaci´ on Para lograr que e−τ (F+ (τ ) + K) sea una función periódica de per´ıodo T , basta tomar f+ (τ ) = p(τ ) + q(τ ) donde p(τ ) sea una función cotinua periódica de per´ıodo T y limτ →∞ q(τ ) = 0. Además, en nuestro caso se puede demostrar que la función q(τ ) decrece a cero en forma exponencial cuando τ tiende a ∞. En lo que sigue B = B[0, ∞) = {u : [0, ∞) → R| u es continua y acotada en [0, ∞)}, para u ∈ B, ||u||∞ = supt∈[0 ∞] |u(t)| sabemos que B es un espacio de Banach con ||u||∞ . Sea M = {u ∈ B : u(0) = u0 }, M es un subespacio cerrado de B. Teorema 2 Si e−τ (F+ (τ )+K) es una función continua en [0, ∞) y peri´ odica 0 de per´ıodo T , g¯ una función tal que V+ supu∈R |¯ g (u)| < 1, entonces la ecuaci´ on diferencial du + u = V+ g¯(u) + e−τ (F+ (τ ) + K) dτ tiene una u ńica solución acotada en [0 ∞)

(3.10)

Demostraci´ on. Por el método de variación de parámetros, u es una solución de (3.10) si y sólo si −τ

u(τ ) = e u0 +

Z τ 0

(e−τ es [V+ g¯(u) + e−s (F+ (s) + K)])ds .

(3.11)

Vamos a construir un operador T : M → M que tiene un u ńico punto fijo y este punto fijo será una solución de la ecuación (3.11) y por lo tanto de la ecuación (3.10). 87

Definamos T : M → M por −τ

T u(τ ) = e u0 +

Z τ 0

(e−τ es [V+ g¯(u) + e−s (F+ (s) + K)])ds .

(3.12)

1. Vamos a demostrar que T u es un elemento del espacio M |T u(τ )| ≤ |e−τ u0 | + ≤ |u0 | + V+ e−τ

Rτ 0

|e−τ +s |V+ g¯(u(s))|ds +

Rτ 0

(|e−τ (F+ (s) + K)|)ds

R Rτ s −τ τ 0 |F+ (s) + K|ds 0 e ds + e

≤ |u0 | + V+ e−τ (eτ − e0 ) + ||F+ (τ ) + K||∞ e−τ τ , por lo tanto |T u(τ )| ≤ |u0 | + V+ (1 − e−τ ) + ||F+ (τ ) + K||∞ eτ τ ≤ |u0 | + V+ + ||F+ (τ ) + K||∞ e−τ τ . Como limτ →∞ τ e−τ = 0 se tiene que |T u(τ )| ≤ m para cada τ ∈ [0 ∞). As´ı si u ∈ M entonces, T u ∈ M 2. Ahora vamos a probar que T es una contracción. Sean u1 , u2 ∈ M , entonces |T u1 (τ ) − T u2 (τ )| ≤

Rτ 0

|e−τ +s [V+ (¯ g (u1 (s))) − (¯ g (u2 (s)))]|ds

≤ V+ e−τ supu∈R |¯ g 0 (u)|||u1 − u2 ||∞

R∞ s 0 e ds

= supu∈R |¯ g 0 (u)|||u1 − u2 ||∞ V+ [1 − e−τ ] ≤ Sup|¯ g 0 (u)|u∈R ||u1 − u2 ||∞ V+ , por lo tanto ||T u1 − T u2 ||∞ ≤ k1 ||u − 1 − u2 ||∞ , con k1 = V+ supu∈R |¯ g 0 (u)| < 1. El operador T es una contracción y por el teorema del punto fijo de Banach, T tiene un u ńico punto fijo u, el cual es una solución de la ecuación diferencial (3.10) En resumen hemos probado que si e−τ (F+ (τ )+K) es una función periódica, de per´ıodo T , entonces la ecuación (3.10) tiene una u ńica solución acotada 88

en [0, ∞). Por el Teorema 4.11 de la página 116 del libro [39], se tiene que la ecuación diferencial (3.9) tiene una solución periódica de per´ıodo T , por la proposición 1, u es una solución periódica de la ecuación (3.6). Para garantizar que la solución periódica u cumpla con | u(τ ) |≤ max{V+ , V− } , es necesario considerar, en lugar del conjunto M , el conjunto cerrado M 0 , donde M 0 es M 0 = {u ∈ B | u(0) = u0 y exista τ0 ≥ 0 tal que | u(τ ) |≤ max{V+ , V− } para cada τ ≥ τ0 } . Se puede demostrar que el operador T definido en (3.12) tiene la propiedad de que T (M 0 ) ⊂ M 0 y as´ı la prueba del teorema 2 se puede repetir. 3.2 CONTROLABILIDAD DE LOS ESTADOS ´ FISIOLOGICOS DE LA CORTEZA CEREBRAL. En esta sección estudiaremos el problema del control de los estados fisiológicos de la corteza cerebral. Utilizaremos como elemento de control a las inhibiciones. Se considera que el papel de las inhibiciones es controlar el proceso de excitación. Esta hipótesis se refleja en el sistema (3.1), (3.2), ya que la función de densidad de actividad sináptica inhibitoria (H¯− ) aparece en la ecuación para la función de activación (u), pero no al revés. En lo que sigue consideraremos que 1. No existen entradas inhibitorias a la corteza cerebral, esto se traduce en que f− (t) ≡ 0 en la ecuación diferencial (3.1). 2. Los tiempos caracter´ısticos de las inhibiciones son mucho menores que los de las excitaciones. Por esto consideraremos que λ + v+ ¿ 1. λ − v− 89

(3.13)

As´ı, la ecuación diferencial que estudiaremos es d2 u dt2

h

+ (λ+ v+ )2 u = V+ λ+ v+ λ+ v+ g¯(u) + g¯0 (u) du + 2λ+ v+ du dt dt

i

−V− e−λ− v− t [2(λ+ v+ − λ− v− )α + (λ+ v+ − λ− v− )2 (αt + β)]

(3.14)

+f+ (t) . ¯ − debe satisfacer la ecuación (3.2) y por lo tanto, debe ser de la ya que H forma ¯ − (t) = αt + β . H Usando el cambio de variable τ = λ+ v+ t y el lema 3, la ecuación diferencial (3.14) es equivalente a la ecuación d2 u dτ 2

h

+ 2 du + u = V+ g¯(u) + g¯0 (u) du dτ dτ λ

−V− e

v

− λ− v− τ + +

i

h

2(λ+ v+ − λ− v− ) (λ+av+ )2 + (λ+ v+ − λ− v− )2 ( (λ+av+ )3 τ +

+f+ (τ ) . (3.15) As´ı, el problema de control que queremos estudiar es el siguiente: dados u0 , u00 , T, uT y u0T , lo que queremos demostrar es que existen a, b n´ umeros reales y u solución de la ecuación (3.15) tal que u(0) = u0 ,

u(T ) = uT y

u0 (0) = u00 ,

u0 (T ) = u0T .

Es decir queremos ver si el sistema dinámico determinado por la ecuación (3.15) es controlable. Para estudiar este problema, nos auxiliaremos de la proposición 1, la cual nos garantiza que u es una solución de (3.15) que satisface las condiciones iniciales u(0) = u0 , u0 (0) = u00 si, y sólo si, u satisface la ecuación diferencial de primer orden

90

i

b ) (λ+ v+ )2 )

du + u − V+ g¯(u) = (f0 (τ ) + ca,b (τ ) + K)e−τ , dτ

(3.16)

K = u0 (0) + u(0) − V+ g¯(u(0)) ,

(3.17)

donde

f0 (t) es una función relacionada con las entradas excitatorias a la corteza cerebral y ca,b (t) es la función (" ³

ca,b (t) = V−

e

" ³

+

e

λ

λ

v

´ "Ã

1− λ− v− t + +

´

v

1− λ− v− t + +

#Ã

−1

!

#

#

λ − v− t 1 1 −1 − + a λ + v+ λ+ v+ λ+ v+ λ + v+ ! )

λ− v− −1 b . λ+ v+

(3.18)

Las constantes a y b están relacionadas con las condiciones iniciales de la inhibición. El problema de control que estudiaremos es: dado u0 , u00 , T > 0, uT , u0T ¿existirán a, b n´ umeros reales y u solución de la ecuación (3.16) tal que u(0) = u0 ,

u(T ) = uT , ?

0

u (0) =

u00

,

0

u (T ) =

u0T

.

Este problema es equivalente a probar que la función U : R2 −→ R2 definida por Ã !

Ã

!

a u(T, a, b) U = 0 b u (T, a, b)

es sobreyectiva, donde u(t, a, b) es la solución de la ecuación diferencial (3.16) que cumple la condición inicial u(0) = u0 . Nota. Por la importancia que tienen las constantes a, b, es que ahora denotaremos a las soluciones de (3.16) por u(t, a, b). La dependencia de las condiciones iniciales se omitirá, ya que éstas condiciones estarán fijas. 91

Sabemos que U es una función diferenciable en R2 . Para demostrar que U es suprayectiva, probaremos que U es un mapeo abierto y también que U (R2 ) es un conjunto cerrado, y como R2 es conexo, se debe tener que U (R2 ) = R2 . Para simplificar la notación, introduciremos las siguientes funciones v(t) = u(t, a + ², b) − u(t, a, b) , ³

f1 (t) = e

λ

v

´ "Ã

1− λ− v− t + +

Ã ³

f2 (t) = Como

λ+ v+ λ− v−

e

!

(3.19) #

λ − v− t 1 1 −1 − + , (3.20) λ + v+ λ + v+ λ + v+ λ+ v+ λ

´

v

1− λ− v− t + +

!Ã

−1

!

λ − v− −1 . λ + v+

(3.21)

¿ 1 se tiene que f1 (t) > 0 y f2 (t) < 0 para cada t ≥ 0 .

Lema 4 Para cada 0 < τ ≤ T se cumple f1 (τ )f2 (T ) − f2 (τ )f1 (T ) < 0 .

(3.22)

Demostraci´ on. Sea f3 (τ ) =

f1 (τ ) . −f2 (τ )

(3.23)

Probaremos primero que f3 es una función decreciente en (0, T ], demostrando que f30 ≤ 0 en (0, T ]. Multiplicando y dividiendo (3.23) por λ

e

v

−(1− λ− v− )τ + +

y sustituyendo (3.20) y (3.21) en (3.23), se tiene 92





1  f3 (τ ) = λ + v+

τ λ

e

+

v

( λ− v− −1)τ + +

1 λ− v − λ+ v +

−1

−1

.

(3.24)

Derivando (3.24), se tiene ³

f30 (τ )

e

λ− v− −1 λ+ v+

´

³

"

τ

λ

v

´

1− λ− v− τ

=Ã ³ ´ !2 1 − e λ− v− −1 τ −1 e λ+ v+

+ +

Ã

! #

λ− v− − −1 τ . λ+ v+

(3.25)

El primer factor de (3.25) es mayor que cero, as´ı que, para probar que < 0 en (0, T ], necesitamos probar que el segundo factor de (3.25) es menor que cero; para esto, consideremos la función f30

³

f4 (τ ) = 1 − e

λ

´

v

1− λ− v− τ + +

Ã

!

λ − v− − −1 τ. λ + v+

Como f4 (0) = 0, entonces, si probamos que f4 es decreciente en (0, T ], se tendr´ıa que f4 < 0 en (0, T ]. Pero Ã

f40 (τ ) =

! ³

λ− v− −1 e λ+ v+

λ

v

´

1− λ− v− τ + +

Ã

!2 Ã ³

λ − v− − −1 λ + v+

e

λ

v

´

1− λ− v− τ + +

!

− 1 ≤ 0,

por lo tanto, f4 es una función decreciente en (0, T ]. As´ı, f4 (τ ) ≤ 0 para τ ∈ (0, T ]. Esto u ´ltimo implica que f30 ≤ 0 en (0, T ], que es lo que quer´ıamos probar. Ahora, por ser f4 una función decreciente, para cada τ ∈ (0, T ],se tiene f1 (t) f1 (T ) ≥ −f2 (t) −f2 (T ) y como f2 (τ ) < 0 para cada τ , se tiene f1 (τ )f2 (T ) − f2 (τ )f1 (T ) < 0 . En lo que sigue, vamos a necesitar las siguientes derivadas ∂u0 ∂u0 ∂u (T, a, b), ∂u (T, a, b), (T, a, b) y (T, a, b) . ∂b ∂a ∂a ∂b As´ı que vamos a calcular estas derivadas. 93

1.

∂u (T, a, b) . ∂a

Como u(t, a, b), u(t, a + ², b) son soluciones de (3.16) se cumple que du (τ, a dτ

+ ², b) + u(τ, a + ², b) − V+ g¯(u(τ, a + ², b)) = (3.26) −τ

(f0 (τ ) + ca+²,b (τ ) + K)e du (τ, a, b) dτ

,

+ u(τ, a, b) − V+ g¯(u(τ, a, b)) = (3.27) −τ

(f0 (τ ) + ca,b (τ ) + K)e

,

entonces, restando (3.27) a (3.26) y usando (3.19) se tiene dv + v − V+ (¯ g (u(τ, a + ², b)) − g¯(u(τ, a, b))) = (ca+²,b (τ ) − ca,b (τ ))e−τ . dτ Ahora, tomando la parte lineal de g¯(u(τ, a, b)), se tiene dv + v − V+ (¯ g 0 (u(τ, a, b)))v = (ca+²,b (τ ) − ca,b (τ ))e−τ , dτ pero como ca+²,b (τ ) − ca,b (τ ) = V− ²f1 (τ ) y

v(T ) = V− ²

Z T 0

e−t e−

RT t

(1−V+ g¯0 (u(τ,a,b)))dτ

f1 (t) dt ,

(3.28)

entonces, usando (3.19) y dividiendo (3.28) por ², se tiene u(T, a + ², b) − u(T, a, b) v(T ) = ² ² = V− 94

Z T 0

e−t e−

RT t

(1−V+ g¯0 (u(τ,a,b)))dτ

f1 (t) dt .

Tomando el l´ımite cuando ² −→ 0, se obtiene Z T RT ∂u −t − t (1−V+ g¯0 (u(τ,a,b)))dτ (T, a, b) = V− e e f1 (t) dt . (3.29) ∂a 0

De manera similar, se puede calcular 2.

∂u (T, a, b) ∂b

.

∂u (T, a, b). ∂b

Z T RT 0 ∂u (T, a, b) = V− e−t e− t (1−V+ g¯ (u(τ,a,b)))dτ f2 (t) dt. ∂b 0

3.

(3.30)

∂u0 (T, a, b). ∂a

Como u(t, a, b), u(t, a + ², b) son soluciones de (3.16) se cumple

u0 (T, a, b) = V+ g¯(u(T, a, b)) − u(T, a, b) + (f (T ) + c(T ) + K)e−T (3.31) y

u0 (T, a+², b) = V+ g¯(u(T, a+², b))−u(T, a+², b)+(f (T )+c(T )+K)e−T , (3.32) Restando (3.31) a (3.32), se obtiene u0 (T, a + ², b) − u0 (T, a, b) = V+ [¯ g (u(T, a + ², b)) − g¯(u(T, a, b))] − [u(T, a + ², b) − u(T, a, b)] + [ca+²,b (T ) − ca,b (T )] = [V+ g¯0 (u(T, a, b)) − 1] [u(T, a + ², b) − u(T, a, b)] + e−T ²V− f1 (T ) . Entonces u0 (T, a + ², b) − u0 (T, a, b) = ² 95

"

#

u(T, a + ², b) − u(T, a, b) [V+ g¯ (u(T, a, b)) − 1] + e−T V− f1 (T ) , ² 0

Tomando el l´ımite cuando ² −→ 0, obtenemos ∂u0 ∂u (T, a, b) = [V+ g¯0 (u(T, a, b)) − 1] (T, a, b) + V− e−T f1 (T ) . (3.33) ∂a ∂a De la misma manera, se obtiene 4.

∂u0 (T, a, b). ∂b

∂u0 ∂u (T, a, b) = [V+ g¯0 (u(T, a, b)) − 1] (T, a, b) + V− e−T f2 (T ) . (3.34) ∂b ∂b

Proposici´ on 2 Sea U : R2 −→ R2 definida por Ã !

Ã

!

a u(T, a, b) = 0 , U b u (T, a, b) donde u(t, a, b) es la soluci´ on de la ecuaci´ on diferencial (3.16) que cumple la condici´ on inicial u(0) = u0 . Entonces U es un mapeo abierto. Demostraci´ on. Basta probar que localmente U es un mapeo abierto. Para probar que U es localmente un mapeo abierto, demostraremos que  ∂u (T, a, b)  ∂a DU (a, b) = det 

∂u0 (T, a, b) ∂a

para cada a, b elementos de R. 96

∂u (T, a, b) ∂b ∂u0 (T, a, b) ∂b

   6= 0

(3.35)

DU (a, b) =

∂u ∂u0 ∂u0 ∂u (T, a, b) (T, a, b) − (T, a, b) (T, a, b) . ∂a ∂b ∂a ∂b

(3.36)

Sustituyendo (3.33) y (3.34) en (3.36) se obtiene "

#

∂u ∂u DU (a, b) = (T, a, b) (T, a, b) [V+ g¯0 (u(T, a, b)) − 1] + V− e−T f2 (T ) ∂a ∂b "

#

∂u ∂u ∂u − (T, a, b) (T, a, b) [V+ g¯0 (u(T, a, b)) − 1] − V− (T, a, b)f1 (T ) , ∂b ∂a ∂b entonces ( −T

DU (a, b) = V− e

)

∂u ∂u (T, a, b)f2 (T ) − (T, a, b)f1 (T ) ∂a ∂b

.

(3.37)

Sustituyendo (3.29) y (3.30) en (3.37) se obtiene "Z

DU (a, b) = V−2 e−T

T 0

−s −

e e

RT s

0

(1−V+ g¯ (u))dθ

#

{f1 (s)f2 (T ) − f2 (s)f1 (T )} ds . (3.38)

Por el lema 4, se tiene f1 (t)f2 (T ) − f2 (t)f1 (T ) < 0 para 0 < t ≤ T y por lo tanto el integrando en (3.38) es estrictamente menor que 0 en el intervalo [0, T ], y por lo tanto DU (a, b) < 0 para cada a, b ∈ R . As´ı, hemos probado que U es localmente un mapeo abierto y por lo tanto que U es un mapeo abierto en R2 . Ahora demostraremos que U (R2 ) es un conjunto cerrado en R2 . Primeramente probaremos el siguiente lema. 97

Lema 5 Para cada T > 0 sea BT = {(u(T, a, b), u0 (T, a, b)) : (a, b) ∈ A}. Si BT es un conjunto acotado, entonces A es un conjunto acotado. Demostraci´ on. Primero demostraremos que los conjuntos ( −T

α(T, a, b) = e

)

Z T

ca,b (s)ds : (a, b) ∈ A

0

(3.39)

y {β(T, a, b) = ca,b (T ) : (a, b) ∈ A}

(3.40)

son acotados. Como u es solución de la ecuación (3.16), es decir satisface la ecuación du + u − V+ g¯(u) = (f0 (τ ) + ca,b (τ ) + K)e−τ , dτ entonces, por el método de variación de parámetros, u debe satisfacer −τ

−τ

u(τ ) = e u(0)+e V+

Z τ 0

s

−τ

g¯(u(s))e ds+e

Z τ 0

(f (s)+ca,b (s)+K)ds . (3.41)

Evaluando (3.41) en T se tiene −T

u(T, a, b) = u(0)e

+ V+

Z T 0

g¯(u(s))e

s−T

−T

ds + e

Z T 0

(f (s) + ca,b (s) + K)ds . (3.42)

Despejando α(T, a, b) de (3.42) se tiene −T

α(T, a, b) = u(T, a, b)−u(0)e

−V+

Z T 0

s−T

g¯(u(s))e

ds−e

−T

Z T 0

(f0 (s)+K)ds , (3.43)

entonces | α(T, a, b) |≤| u(T, a, b) | + | u(0)e−T | + | V+

RT 0

s−T

g¯(u(s))e

−T

ds | + | e 98

RT 0

(3.44) (f0 (s) + K)ds | .

Ahora, por hipótesis existe N1 (T ) tal que | u(T, a, b) |≤ N1 (T ) para cada (a, b) ∈ A .

(3.45)

Además, como g¯ es acotada por 1 en R y e−τ ≤ 1 para τ ≤ 0, se tiene ¯ ¯ Z T ¯ ¯ ¯ ¯ s−T ¯ ≤ V+ T ¯ V+ g ¯ (u(s))e ds ¯ ¯ 0

(3.46)

y como f0 (τ ) es continua en el intervalo [0, T ], se tiene que, existe M (T ) tal que ¯ ¯Z ¯ T ¯ ¯ ¯ ¯ (f 0 (s) + K)ds ¯ ≤ M (T ) . ¯ 0 ¯

(3.47)

Usando (3.45), (3.46), (3.47) en (3.44) se prueba que | α(T, a, b) |≤ N1 (T ) + u0 + V+ (T ) + M (T ) para cada (a, b) ∈ A

(3.48)

que es lo que quer´ıamos probar. Por (3.41) se tiene u0 (T, a, b) + u(T, a, b) − V+ g¯(u(T, a, b)) = (f (T ) + ca,b (T ) + K)e−T . (3.49) Despejando β(T, a, b) de (3.49) se tiene que β(T, a, b) = u0 (T, a, b) + u(T, a, b) − V+ g¯(u(T, a, b)) − (f (T ) + K)e−T , (3.50) entonces,

| β(T, a, b) |≤| u0 (T, a, b) | + | u(T, a, b) | + (3.51) | V+ g¯(u(T, a, b)) | + | (f (T ) + K)e−T | .

Por la hipótesis existe N2 (T ) tal que | u0 (T, a, b) |≤ N2 (T ), para cada (a, b) ∈ A .

(3.52)

Usando (3.45) y (3.52) se tiene que | β(T, a, b) |≤ N2 (T ) + N1 (T ) + V+ + M1 (T ), para cada (a, b) ∈ A (3.53) 99

y por lo tanto se obtiene que {β(T, a, b) = ca,b (T ) : (a, b) ∈ A} es un conjunto acotado en R. Ahora, demostraremos que a y b deben satisfacer un par de ecuaciones. Z T

α(T, a, b) = ³

−V− e

λ

v

´

1− λ− v− T + +

ca,b (s)ds =

0

³

(

− 1 a + V− (−1)e λ + v+

λ

v

´

1− λ− v− T + +

Ã

!

)

λ − v− +1− −1 T b λ + v+

= A(T )a + B(T )b , donde ³

A(T ) = −V− e ³

(

λ

λ

v

+ +

v

´

− 1− λ− v− T

B(T ) = V− (−1)e

´

1− λ− v− T

+ +

1 , λ + v+ Ã

!

λ − v− −1 T +1− λ + v+

)

y β(T, a, b) = V− f1 (T )a + f2 (T )b . Por lo tanto, tenemos que se cumple A(T )a + B(T )b = α(T, a, b) , (3.54) V− f1 (T )a + V− f2 (T )b = β(T, a, b) . Si probamos que

  

A(T )

B(T )

V− f1 (T ) V− f2 (T ) 100

  

es invertible, se tendr´ıa que   

a b





  =

A(T )

−1

B(T )

Ã

 

V− f1 (T ) V− f2 (T )

α(T, a, b) β(T, a, b)

!

para cada (a, b) ∈ A.

Pero como hemos demostrado que α(T, a, b), β(T, a, b) son conjuntos acotados para cada T > 0, se concluye que A es un conjunto acotado. Probaremos que A(T )f2 (T ) − B(T )f1 (T ) 6= 0 . V−2 ³

A(T )f2 (T ) − B(T )f1 (T ) =e V−2 ³

Ã

λ

´

v

´

v

+ +

Ã

1 f2 (T )− λ + v+

! !

λ− v− +1− − 1 t f1 (T ) . λ+ v+

1− λ− v− T

−e

λ

1− λ− v− T

(3.55)

+ +

Como f2 (T ) < 0, ³

λ

´

v

1− λ− v− T

e

+ +

1 f2 (T ) < 0 . λ+ v+

Sabemos también que f1 (T ) > 0, as´ı si probamos que Ã ³

e

λ

´

v

1− λ− v− T + +

Ã

!

λ − v− −1+ −1 T λ + v+

!

< 0,

(3.56)

demostrar´ıamos la desigualdad A(T )V− f2 (T ) − B(T )V− f1 (T ) < 0. V−2 Para probar (3.56), consideremos la función Ã ³

f5 (t) = e

λ

v

´

1− λ− v− t + +

Ã

!

λ − v− −1+ −1 T λ + v+ 101

!

.

Se tiene que f5 (0) = 0 y su derivada es ! Ã ³

Ã

f50 (t) Como

λ− v− λ+ v+

λ − v− t e = 1− λ + v+

λ

´

v

1− λ− v− t + +

!

−1

.

À 1, se concluye que f50 (t) < 0 para cada t ≥ 0 .

Por lo tanto f5 es un función decreciente en [0, ∞) y como f5 (0) = 0, se tiene f5 (T ) ≤ f5 (0) = 0 . Con esto u ´ltimo probamos (3.56). Por lo tanto hemos demostrado que   

A(T )

B(T )

V− f1 (T )

  

V− f2 (T )

es invertible. Que es lo que quer´ıamos demostrar. En la siguiente proposición, probaremos que U (R2 ) es un conjunto cerrado en R2 . Proposici´ on 3 Si U : R2 −→ R2 definida por  

a

U





 

= 



u(T, a, b)

   ,

u0 (T, a, b)

b

donde u(t, a, b) es una solución de (3.16) con u(0) = u0 . Entonces U (R2 ) es un conjunto cerrado en R2 . 102

Demostraci´ on. Sea (uT , u0T ) ∈ U (R2 ), existe una sucesión {(an , bn ) : n ∈ N } tal que 



 lim  U 

n→∞

an





   = lim  n→∞

bn

u(T, an , bn )





  =

u0 (T, an , bn )

uT

   .

u0T

³ ´ a

Por demostrar que existen b en R2 y u solución de (3.16) tal que 

U 

a





  =

b

u(T, a, b)





  =

u0 (T, a, b)

uT

   .

u0T

Como {u(T, an , bn ), u0 (T, an , bn ) : n ∈ N } es un conjunto acotado en R2 , por el lema 5 se tiene que {(an , bn ) | n ∈ N } es un conjunto acotado en R2 , por el teorema de Bolzano-Weierstrass, existe una subsucesión convergente {(anj , bnj ) | j ∈ N } de {(an , bn ) | n ∈ N }. Sea (a, b) el l´ımite de esta subsucesión. Por la continuidad de las soluciones u(t, anj , bnj ) se tiene que lim u(T, anj , bnj ) = u(T, a, b) ,

j→∞

as´ı u(T, a, b) = uT . Además, como u0 (T, anj , bnj ) = −u(T, anj , bnj ) + V+ g¯(u(T, anj , bnj )) +(f (T ) + canj ,bnj (T ) + K)e−T , lim (−u(T, anj , bnj )) = −u(T, a, b) ,

j→∞

103

lim (V+ g¯(u(T, anj , bnj ))) = V+ g¯(u(T, a, b)) ,

j→∞

y lim ((f (T ) + canj ,bnj (T ) + K)e−T ) = (f (T ) + ca,b (T ) + K)e−T

j→∞

se tiene que lim u0 (T, anj , bnj ) = u0 (T, a, b) .

j→∞

por lo tanto hemos demostrado que U (R2 ) es un conjunto cerrado en R2 . En resumen, hemos probado que el mapeo U : R2 −→ R2 definido por  

a

U





  =

u(T, a, b)

   ,

u0 (T, a, b)

b

donde u(T, a, b) es una solución de (3.16) con condición inicial u(0) = u0 es un mapeo sobreyectivo. Teorema 3 El sistema dinámico d2 u dτ 2

h

+ 2 du + u = V+ g¯(u) + g¯0 (u) du dτ dτ λ

−V− e

v

− λ− v− τ + +

i

h

2(λ+ v+ − λ− v− ) (λ+av+ )2 + (λ+ v+ − λ− v− )2 ( (λ+av+ )3 τ +

+f+ (τ ) (3.57) umeros es controlable en R2 , es decir, dados u0 , u00 y T, uT y u0T , existen a, b n´ reales y u soluci´ on de la ecuaci´ on (3.57) tal que u(0) = u0 ,

u(T ) = uT y u0 (T ) = u0T .

u0 (0) = u00 , 104

i

b ) (λ+ v+ )2 )

Demostraci´ on. Sean u0 , u00 y T, uT y u0T dados. Por la proposición 2, U : R2 → R2 definido por  

a

U





  =

b

u(T, a, b)

   ,

u0 (T, a, b)

donde u(t, a, b) es una solución de du + u − V+ g¯(u) = (f0 (τ ) + ca,b (τ ) + K)e−τ , dτ con K = u0 (0) + u(0) − V+ g¯(u(0)), es un mapeo abierto en R2 . Por la proposición 3, U es un mapeo cerrado en R2 , por lo tanto como R2 es conexo, se tiene que U es suprayectiva, es decir, existen a, b ∈ R y u solución de du + u − V+ g¯(u) = (f0 (τ ) + ca,b (τ ) + K)e−τ , dτ con u(0) = u0 ,

u(T ) = uT ,

u0 (0) = u00 ,

u0 (T ) = u0T

y K = u0 (0) + u(0) − V+ g¯(u(0)) . Pero por la proposición 1, u también es solución de (3.57), y as´ı hemos demostrado el teorema. 105

3.3 EXISTENCIA DE ONDAS VIAJERAS EN LA CORTEZA CEREBRAL En esta sección estudiaremos la existencia de estados fisiológicos que se corresponden con la existencia de ondas viajeras (ondas que se propagan por el medio cortical sin perder su forma) en la ecuación (2.71). En el caso de modelos matemáticos de neuronas se ha demostrado la existencia de ondas viajeras [30], [48], [101]. También se ha estudiado este problema en otro tipo de procesos de reacción difusión [102]. Buscaremos bajo qué condiciones la ecuación diferencial parcial "

#

2 ∂2u ∂u ∂u 2 2 ∂ u 0 + 2λ v + (λ v ) u − v = V λ v λ v g ¯ (u) + g ¯ (u) + + + + + + + + + + ∂t2 ∂t ∂x2 ∂t (3.58) tenga soluciones de la forma: u(x, t) = φ(x − ct) para algunos valores de la constante c. Como en la sección anterior, haremos algunos cambios de variable, que no alteran el resultado buscado, pero que nos permitirán construir una ecuación diferencial parcial equivalente a (3.58) más sencilla de estudiar. Si hacemos el cambio de variable τ = λ+ v+ t y y = λ+ x, la ecuación (3.58) es equivalente a la ecuación

"

#

∂2u ∂u ∂ 2u ∂u +2 + u − 2 = V+ g¯(u) + g¯0 (u) . 2 ∂τ ∂τ ∂y ∂t

(3.59)

Para esta ecuación buscaremos soluciones de la forma u(y, τ ) = φ(y − cτ ) = φ(s) con s = y − cτ y tal que φ(−∞) = s0 = φ(∞) donde s0 es constante. Antes de estudiar este tema, enunciaremos algunos resultados importantes que usaremos más adelante en esta sección. Estos resultados no serán demostrados. Para ver la demostración de los Teoremas 4 y 5, recurrir a [23]. Teorema 4 Sea A una matriz n × n y f : R → Rn una función continua acotada Consideremos la ecuaci´ on diferencial

106

du = Au + f (t) . (3.60) dt La condici´ on necesaria y suficiente para que la ecuaci´ on (3.60) tenga una u ńica solución acotada en R, es que σ(A) no intersecte al eje imaginario. En este caso la solución acotada se expresa en la forma. u(t) =

Z ∞ −∞

G(t − τ )f (τ )dτ ,

(3.61)

donde G(t) es la función de Green principal de A, la cual se define como (

G(t) =

eAt P− , −eAt P+ ,

t>0 t < 0,

(3.62)

donde P− , P+ son, respectivamente los proyectores sobre los subespacios propios con parte real positiva, negativa.

Teorema 5 Consideremos la ecuaci´ on du = Au + f (t, u) . dt Supongamos que se cumplen las condiciones siguientes

(3.63)

1. A es una matriz n × n, tal que σ(A) ∩ ({0} × R) = ∅ .

(3.64)

2. Sea ρ > 0. Existe M > 0 tal que | f (t, u) |≤ M

para t ∈ R y | u |≤ ρ .

(3.65)

3. Existe q > 0 tal que

| f (t, u1 ) − f (t, u2 |≤ q | u1 − u2 | para | u1 |≤ ρ y | u2 |≤ ρ . (3.66)

107

4. Existe N > 0 y r > 0 tal que | G(t) |≤ N e−rt para cada t ∈ R .

(3.67)

5. M y q cumplen con M≤

ρr r y q< . 2N 2N

(3.68)

Entonces, la ecuaci´ on (3.63) tiene una u ńica solución acotada en R, tal que | u(t) |≤ ρ para cada t ∈ R .

(3.69)

Ahora, demostraremos un lema, que usaremos posteriormente. on diferencial Lema 6 Sea f : R2 → R2 continua. Supongamos que la ecuaci´ dx = f (x), x ∈ R2 (3.70) dt tiene un u ńico punto estacionario x0 , el cual es un punto silla. Entonces la condici´ on necesaria y suficiente para que la ecuaci´ on (3.70) tenga una solución homocl´ınica es que tenga una solución no trivial acotada en R. Demostraci´ on. La suficiencia es evidente. Sea φ una solución no trivial acotada de (3.70). Como x0 es el u ńico punto estacionario de la ecuación (3.70), y como además es un punto silla, el teorema de Poincaré-Bendixon, nos asegura que 1. limt→±∞ φ(t) = x0 o 2. ω(φ) = {x0 } = α(φ), donde ω(φ) y α(φ) son, respectivamente los conjuntos ω-l´ımite y α-l´ımite de la órbita generada por la solución φ . 108

En ambos casos, φ ser´ıa la solución homocl´ınica no trivial buscada. En lo que sigue, retomaremos el problema de la existencia de las ondas viajeras para la ecuación (3.59) Ã

!

∂u ∂2u ∂u ∂ 2u + 2 + u − = V+ g¯0 (u) + g¯(u) . 2 2 ∂τ ∂τ ∂y ∂τ Es decir, buscaremos soluciones de la forma u(y, τ ) = φ(y − cτ ) ,

(3.71)

donde φ es una función de una variable que cumple que el lim φ(s) = φ0 ,

s→±∞

(3.72)

donde φ0 es una constante. Ahora veremos que la existencia de soluciones de la ecuación (3.59) que son de la forma (3.71) y que cumplen con la condición (3.72), es equivalente a la existencia de soluciones de una ecuación diferencial ordinaria de segundo orden, que cumplen con la condición (3.72). Si u es una solución de la ecuación (3.59) de la forma (3.71), que cumple con la condición (3.70), entonces usando la regla de la cadena, obtenemos ∂u = −c dφ , ∂τ ds

2 ∂ 2 u = c2 d φ , ∂τ 2 ds2

∂u = dφ , ∂y ds

2 ∂ 2u = d φ . 2 ∂y ds2

Sustituyendo estas derivadas en la ecuación (3.59) y reordenando la ecuación, obtenemos que φ debe ser solución de la ecuación d2 φ dφ (c − 1) 2 + c(V+ g¯0 (φ) − 2) + (φ − V+ g¯(φ)) = 0 . (3.73) ds ds Rec´ıprocamente, si φ es solución de (3.73) que cumple con la condición (3.72) entonces la función 2

u(y, τ ) = φ(y − cτ ) es una solución de la ecuación (3.59). 109

Observaci´ on. En lo que sigue supondremos que 0 < c < 1, ya que esto corresponde a una velocidad real de propagación, cV+ , la cual no puede ser superior a la velocidad de propagación de los potenciales de acción. También tenemos que la ecuación (3.73) es equivalente al sistema de ecuaciones de primer orden. Ã !



d φ =  dt ψ

0

1

1 1−c2

−2c 1−c2

  

Ã !

Ã

!

0 φ V+ , + ψ 1 − c2 c¯ g 0 (φ)ψ − g¯(φ)

(3.74)

donde ψ(s) = dφ ds La existencia de ondas viajeras para la ecuación (3.59) es equivalente a la existencia de soluciones homocl´ınicas para la ecuación (3.74) para alg´ un valor de c ∈ (0, 1). De (3.73) es evidente que para c = 1 no puede haber soluciones homocl´ınicas para esta ecuación. En lo que sigue supondremos que g¯0 (φ) <

1 para cada φ ∈ R . V+

(3.75)

Con esta condición veremos que el sistema (3.74) tiene un u ńico punto estacionario para cualquier valor de la constante c que cumple 0 < c < 1. Además, este punto estacionario es independiente de c y es un punto silla del sistema. Con este resultado, el problema de encontrar soluciones homocl´ınicas para el sistema (3.74) se reduce al de buscar soluciones acotadas en R, diferentes de la trivial de este sistema. ³ ´ φ Como las soluciones estacionarias de (3.74) son de la forma 00 donde φ0 es ra´ız de la ecuación φ − V+ g¯(φ) = 0 . De la condición V+ g¯0 (φ) < V1+ , se concluye que esta ecuación tiene una u ńica solución φ0 , que corresponde a una u ńica solución estacionaria de (3.74). 110

Para determinar si el punto estacionario es un punto silla, estudiaremos la “parte lineal” del sistema en un entorno del punto estacionario. Calculando la derivada del lado derecho de (3.74) en el punto estacionario y haciendo algunos cálculos obtenemos que la “parte lineal” del sistema en un entorno del punto estacionario es Ã !



d φ  = dt ψ

0



1

1 (V+ g¯0 (φ0 ) c2 −1

− 1)

c (2 c2 −1

 

− V+ g¯0 (φ0 ))

Ã !

φ . ψ

(3.76)

Su ecuación caracter´ıstica es : ·

¸

c 1 λ λ+ 2 (V+ g¯0 (φ0 ) − 2) + 2 (V+ g¯0 (φ0 ) − 1) = 0 . c −1 c −1

(3.77)

Las ra´ıces de (3.77) son

λ± =

c (2 c2 −1

− V+ g¯0 (φ0 )) ±

q

( c2c−1 (2 − V+ g¯0 (φ0 ))2 +

1 (V+ g¯0 (φ0 ) c2 −1

2

− 1))

.

(3.78) Ahora, como V+ g¯ (φ0 ) < 1 y 0 < c < 1, obtenemos 0 < λ+ y λ− < 0. 0

³ ´ φ

Por lo tanto 00 es un punto silla de (3.76) para cualquier valor de 0 < c < 1. Seg´ un el lema 6, la existencia de soluciones homocl´ınicas de (3.74) es ³ ´ φ0 equivalente a la existencia de soluciones acotadas diferentes de la trivial 0 . Por el teorema 4, si se denota por 

Ac =  

0

1

1 1−c2

−2c 1−c2

  

(3.79)

y G(t, c) la correspondiente función de Green principal, entonces las soluciones acotadas en R de (3.74) coinciden con las soluciones acotadas de la ecuación integral Ã ! V+ Z ∞ 0 φ G(t − τ, c) dτ . = 1 − c2 −∞ c¯ g 0 (φ(τ )ψ(τ ) − g¯(τ ) ψ

Ã !

111

(3.80)

Ahora estudiaremos bajo qué condiciones la ecuación integral (3.80) tiene soluciones acotadas. Lema 7 La función de Green principal de Ac es     1 − c c2 − 1   t   1 c−1   e    ,  2     −1 1 + c 

si t > 0

    1 + c 1 − c2   t     − 12 e c+1    ,     1 1−c

si t < 0 .

G(t, c) = 

(3.81)

Demostraci´ on. Seg´ un (3.62), requerimos calcular eAc t y los proyectores P+ , P− . Los valores propios de Ac son 1 1 , . (3.82) 1+c c−1 Usando la forma canónica de Jordan de Ac , y realizando algunos cálculos, que no escribiremos aqu´ı, obtenemos 

eAc t

1 t  = e 1+c  2

1 + c 1 − c2 1

1−c





1 t    + e c−1  2

1 − c c2 − 1 −1

   .

(3.83)

1+c

Los proyectores son Ã !

Pc+

x 1 =  2 y

Ã !

Pc−





1 x =  y 2

1 + c 1 − c2 1

112

Ã !  x , 

y

1−c

1 − c c2 − 1 −1



1+c

(3.84)



Ã !  x . 

y

(3.85)

Sustituyendo (3.83), (3.84) y (3.85) en (3.62), se obtiene     1 − c c2 − 1   t   1 c−1    e   ,  2     −1 1 + c 

si t > 0

   1+c   t    1 c+1  −2e      1

si t < 0 .

G(t, c) = 

1 − c2

(3.86)

   ,

1−c

Lema 8

−|t|

| G(t) |≤ (1 + c)e 1+c .

(3.87)

Demostraci´ on. Se concluye del lema 7 y del hecho que la norma del operador   M (x, y) = 

a11 a12



Ã !  x 

a21 a22

y

cumple | M |≤ 2K , donde K=

max

1≤i≤2, 1≤j≤2

| aij | .

Para obtener estimados más simples –sin perder generalidad– del integrando en la ecuación integral, supondremos 1 con α > 0, β > 0 . (3.88) 1 + e−αu+β Obteniendo la primera y segunda derivada de g¯, definida por (3.88) se tiene g¯(u) =

g¯0 (u) = α¯ g (1 − g¯) , 113

(3.89)

g¯00 (u) = α¯ g 0 (1 − 2¯ g) ,

(3.90)

| g¯0 |≤ α ,

(3.91)

| g¯00 |≤ 2α2 .

(3.92)

Lema 9 Sea ρ > 0. Si | φ |2 + | ψ |2 ≤ ρ2 , | φi |2 + | ψi |2 ≤ ρ2 para 1 ≤ i ≤ 2 . Entonces, se tienen las siguientes desigualdades | c¯ g 0 (φ)ψ − g¯(φ) |≤ cαρ + 1 ,

(3.93)

| (c¯ g 0 (φ1 )ψ1 − g¯(φ1 )) − (c¯ g 0 (φ2 )ψ2 − g¯(φ2 )) |≤ 2α(αρc + 1) .

(3.94)

Demostraci´ on. Se sigue inmediatamente de (3.91) y (3.92). Proposici´ on 4 Supongamos que en el sistema (3.74) la función g¯ tiene la forma 1 g¯(φ) = , α > 0, y β > 0 . 1 + e−αφ+β 1 } Supongamos además que V+ α < 14 . Entonces para cada ρ > max{φ0 , 2α existe 0 < c(ρ) < 1 tal que para cada 0 < c < ³c(ρ), la u ńica solución acotada ´ φ de este sistema contenida en la bola cerrada { ψ ∈ R2 :| φ |2 + | ψ |2 ≤ ρ2 },

³

es la solución estacionaria ³

c(ρ) =

1+

1 αρ

+

1 4α2 ρV+

´

´ φ0 as c(ρ) viene dada por 0 . Adem´

+

r³

1+

1 αρ

+

2 114

1 4α2 ρV+

´2

−4

³

1 αρ

−

1 4α2 ρV+

´

. (3.95)

Observaci´ on.

³ ´ φ

La condición V+ α < 14 garantiza que 00 es el u ńico punto estacionario del sistema (3.74), donde φ0 es la u ńica solución de la ecuación φ − V+ g¯(φ) = 0 . Demostraci´ on. Demostraremos que bajo las condiciones de esta proposición, el sistema (3.74) cumple con las condiciones del teorema 4, y que por lo tanto este sistema ńica solución acotada u³ contenida en la bola cerrada ³ ´ debe tener una u ´ φ φ0 2 2 2 2 { ψ ∈ R :| φ | + | ψ | ≤ ρ }, pero como 0 es también una solución acotada del sistema contenida en la misma bola cerrada, se debe tener que ³ ´ φ0 ńica solución del sistema (3.74) en esta bola cerrada. 0 es la u Pasemos a demostrar que el sistema (3.74) cumple con las hipótesis del teorema 5. Los valores propios de la matriz Ac son 1 , 1+c

1 , c−1

por lo tanto Ac cumple con la condición (3.64). Por el lema 9 se tiene V+ 1 − c2

¯Ã ¯Ã !¯ !¯ ¯ ¯ φ ¯ ¯ V+ 0 ¯ ¯ ¯ ¯ (cαρ + 1) para ¯ ¯ ¯≤ ¯ ≤ ρ, 0 2 ¯ c¯ ¯ ¯ g (φ)ψ − g¯(φ) 1−c ψ ¯2

V+ 1 − c2

¯Ã !¯ ¯ ¯ 0 ¯ ¯ ¯ ¯ 0 0 ¯ c¯ g (φ1 ) − g¯(φ1 ) − (c¯ g (φ2 )ψ2 − g¯(φ2 )) ¯ ¯Ã

!

Ã

¯ φ V+ φ2 1 ¯ ≤ 2α(αρc + 1) ¯ − 2 ¯ 1−c ψ1 ψ2

para

¯Ã !¯ ¯Ã !¯ ¯ φ ¯ ¯ φ ¯ 1 ¯ 2 ¯ ¯ ¯ ¯ ¯ ≤ρ y ¯ ¯ ≤ ρ. ¯ ψ1 ¯ ¯ ψ2 ¯ 2

2

115

!¯ ¯ ¯ ¯ ¯

2

2

As´ı, también se cumplen las condiciones (3.66), (3.67), con M = (cαρ + 1)

V+ , 1 − c2

q = 2α(αρc + 1)

V+ . 1 − c2

Por los lemas 7 y 8, se tiene que la función de Green principal, G(t, c), de Ac cumple −|t| | G(t, c) |≤ (1 + c)e 1+c , por lo tanto, también se cumple (3.67) con N =1+c y r =

1 . 1+c

Sólo nos falta probar la condición (3.68), es decir, nos falta probar que bajo las condiciones de la proposición 4 se cumple (cαρ + 1)V+ ρ ≤ , 2 1−c 2(1 + c)2

(3.96)

2α(αρc + 1)V+ 1 < . 2 1−c 2(1 + c)2

(3.97)

Multiplicando (3.96), (3.97) por 1−c2 y simplificando obtenemos que este sistema de desigualdades es equivalente a (cαρ + 1)V+ ≤

ρ(1 − c) , 2(1 + c)

2α(αρc + 1)V+ <

(3.98)

1−c . 2(1 + c)

(3.99)

Realizando algunas operaciones y reordenando obtenemos que el sistema (3.98), (3.99) es equivalente al sistema Ã

!

1 1 1 1 P1 (c) = c2 + 1 + + c+ − ≤ 0, αρ 2V+ α αρ 2V+ α Ã

!

Ã

1 1 1 1 + 2 c+ − 2 P2 (c) = c + 1 + αρ 4α ρV+ αρ 4α ρV+ 2

116

(3.100)

!

< 0.

(3.101)

As´ı que necesitamos probar que P1 (c) ≤ 0 , P2 (c) < 0 .

(3.102) (3.103)

Para probar esto, obtengamos las ra´ıces de los polinomios cuadráticos P1 (c), P2 (c). Las ra´ıces de P1 son ³

λ1± (ρ) =

− 1+

1 αρ

+

1 2V+ α

´

±

r³

1 αρ

1+

+

1 2V+ α

´2

−4

³

1 αρ

−

1 2V+ α

2

´

. (3.104)

Las ra´ıces de P2 son ³

λ2± (ρ) =

− 1+

1 αρ

+

1 2V+ α

´

±

r³

1+

1 αρ

+

1 4ρV+ α2

´2

−4

2

Como por hipótesis V+ α <

1 4

yρ>

1 2α

³

1 αρ

−

´

1

4α2 ρV

+

.

(3.105)

se tiene

1 1 − 2 < 0, αρ 4α ρV+ por lo tanto las ra´ıces de P1 (c) son reales, una negativa, λ1− (ρ), y la otra positiva λ1+ (ρ). Lo mismo ocurre con P2 (c), tiene una ra´ız negativa, λ2− (ρ) y otra ra´ız positiva, λ2+ (ρ). Entonces P1 (c) ≤ 0 para λ1− (ρ) ≤ c ≤ λ1+ (ρ) ,

(3.106)

P2 (c) < 0 para λ1− (ρ) < c < λ2+ (ρ) .

(3.107)

Además λ2+ (ρ) < λ1+ (ρ) , 117

(3.108)

esto u ´ltimo se sigue de la expresión para cada uno de ellos y del hecho que 4α2 ρV+ = 2αρ(2V+ α) > (2V+ α) ya que ρ >

1 . 2α

También es fácil probar que λ1+ (ρ) < 1. As´ı, si tomamos 0 < c(ρ) = < 1 se tiene que, para cualquier 0 < c < c(ρ)

λ1+ (ρ)

P1 (c) ≤ 0 , P2 (c) < 0 , que es lo que deseábamos demostrar. Resumiendo, con las condiciones de la proposición 4, el sistema (3.74) cumple con las condiciones del teorema 5, para cualquier 0 < c < c(ρ). Esto es lo que quer´ıamos demostrar. 1 }, Con la proposición 4, hemos demostrado que si V+ α < 14 , ρ > max{φ0 , 2α entonces el sistema (3.74) no puede tener soluciones homocl´ınicas para 0 < c < c(ρ), por lo tanto, la no existencia de ondas viajeras para la ecuación (3.59), cuando tomamos 0 < c < c(ρ). 1 } Ahora pasaremos a analizar qué sucede cuando tomamos ρ > max{φ0 , 2α y c(ρ) < c < 1. 1 } el espacio β , Para esto, consideremos, para ρ > max{φ0 , 2α ρ (Ã !

βρ =

Ã !

φ φ : R −→ R2 : ψ ψ

)

es continua, k φ

k2∞

+kψ

k2∞ ≤

ρ

2

.

Sabemos que βρ es un de Banach. Consideremos la familia de µ³ espacio ´ ¶ φ operadores integrales F ψ , c definida por ÃÃ !

F

!

Ã !

φ φ ,c = − F1 ψ ψ

ÃÃ !

!

φ ,c , ψ

(3.109)

donde ÃÃ !

F1

! Ã ! V+ Z ∞ φ 0 ,c = G(t − τ, c) dτ . c¯ g 0 (φ(τ ))ψ(τ ) − g¯(φ(τ )) ψ 1 − c2 −∞ (3.110)

118

µ³ ´ ¶ φ ue de βρ a βρ , c debe cumplir con la condición ψ , c act´

Para que F

(3.100). De aqu´ı conclu´ımos que aquellos valores de c para los cuales es posible encontrar soluciones homocl´ınicas del sistema (3.74), son los que cumplen λ2+ (ρ) < c < λ1+ (ρ) . As´ı consideremos F : βρ × (λ2+ (ρ), λ1+ (ρ)) −→ βρ definido por (3.109) . ³

φ

´

Como 00 es un punto estacionario del sistema (3.74) para cada 0 < c < 1, tenemos que ÃÃ !

F

!

φ , c = 0 para cada λ2+ (ρ) < c < λ1+ (ρ) . 0 ³³ ´

´

Entonces, si lograramos probar que DF φ00 , c es invertible, y que se cumplen ³³ las otras condiciones del teorema de la función impl´ıcita ³ en ´ un en´ ´ φ0 φ0 torno de , c , se concluir´ıa que existen U, V vecindades de 0 y c res0 pectivamente, tal que para cada c∗ ∈ V , existe una u ńica L (c∗ ) = tal que ÃÃ !

³ ´ φ ψ

∈U

!

φ , c∗ = 0 , ψ

F ³ ´ φ

pero entonces, ψ ser´ıa una solución acotada en R, no trivial del sistema (3.74) y por lo tanto una solución homocl´ınica del mismo sistema, a la cual le corresponder´ıa una onda viajera de la ecuación (3.59) con velocidad c∗ . As´ı que, en lo que sigue probaremos que el operador integral F con λ1+ (ρ) < c < λ2+ (ρ), cumple con las condiciones del teorema de la función impl´ıcita. ³³ ´ ´ Vamos ahora a calcular la derivada de Frechet de F en ψφ , c con λ1+ (ρ) < c < λ2+ (ρ). Por la forma que tiene F, basta calcular la derivada de Frechet de la función F1 . Además se tendr´ıa que ÃÃ !

DF

!

φ , c = I − DF1 ψ 119

ÃÃ !

!

φ ,c . ψ

(3.111)

³ ´ φ

Para hacer esto, calcularemos las derivadas de Gato de F1 en ψ con res³³ ´ ´ pecto a φ y ψ respectivamente. Denotaremos estas derivadas por δφ F1 ψφ , c , δψ F2

³³ ´ φ ψ

´

, c respectivamente.

Lema 10 Bajo las condiciones de la proposici´ on 4, el³³operador ´ ´ integral F definido en (3.110), es de clase C 1 en un entorno de φ00 , c . Adem´ as la derivada parcial de Frechet de F con respecto a DF

³³ ´ φ ψ

³ ´ φ ψ

´

, la cual denotamos por

, c , viene dada por ÃÃ !

DF con

ÃÃ !

!

φ , c = I − DF1 ψ

!

φ , 0 = δφ F1 ψ

DF1

ÃÃ !

ÃÃ !

!

φ ,c , ψ

!

φ , c + δψ F1 ψ

ÃÃ !

φ ,c ψ

!

(3.112)

y ÃÃ !

!

φ ,c , ψ

δφ F1

ÃÃ !

δψ F1

φ ,c ψ

!

son las derivadas de Gato, con respecto a φ y ψ, respectivamente, las cuales vienen dadas por δφ F1

³³ ´ φ ψ

´

, c (h) =



" ! # Ã t −τ c+1 −V+  e c+1 Z ∞ c+1 e [ c¯ g 00 (φ(τ ))ψ(τ ) − g¯0 (φ(τ ) ]h(τ ) dτ 2 1 + c t 1 "

Ã

!

#



 −τ e c−1 Z t c−1 + e c−1 [ c¯ g 00 (φ(τ ))ψ(τ ) − g¯0 (φ(τ )) ]h(τ ) dτ  , c + 1 −∞ 1 t

ÃÃ !

δψ F1



!

Ã

(3.113)

!

−τ −V+ c  e c+1 Z ∞ c+1 φ c+1 0 , c (k) = e g¯ (φ(τ ))k(τ ) dτ ψ 2 1 + c t 1

+

t c−1

e 1+c

Z t −∞

Ã −τ

e c−1

t

!

 

c−1 0 g¯ (φ(τ ))k(τ ) dτ .  1 120

(3.114)

Demostraci´ on. Observemos que F1

³³ ´ φ ψ

´

,c =

 t " Ã ! # −τ −V+  e c+1 Z ∞ c+1 c+1 0 e c¯ g (φ(τ ))ψ(τ ) − g¯(φ(τ )) dτ 2 1 + c t 1 "

Ã

!



#

 −τ e c−1 Z t c−1 + e c−1 (c¯ g 0 (φ(τ ))ψ(τ ) − g¯(φ(τ )) dτ  . c + 1 −∞ 1 t

Calculemos la derivada de Gato de F1 en el punto

³³ ´ φ ψ

aφ ÃÃ !

δφ F1  



!

F1 φ , c (h) = lim  ²→0 ψ



Ã

−τ c+1 V+ e c+1 Z ∞ c+1  =− lim e  ²→0 2 1+c t 1 t

³³

´

φ+²h ψ

´

, c − F1

Ã

³³ ´ φ ψ



dτ

(¯ g 0 (φ(τ ) + ²h(τ )) − g¯(φ(τ )))ψ(τ ) c ²

g¯(φ(τ )) + ²h(τ ) − g¯(φ(τ )) − ²

δφ F1

´

!#

!" Ã

Ã

Entonces

,c

(¯ g 0 (φ(τ ) + ²h(τ )) − g¯0 (φ(τ )))ψ(τ ) c ²

g¯(φ(τ ) + ²h(τ ) − g¯(φ(τ )) − ² t

φ ψ

!" Ã

Ã

−τ e c+1 Z t c+1 c−1 + e c + 1 −∞ 1

´

, c con respecto

³³ ´

²

(3.115)

´

, c (h) = 121

!#

#)

dτ

.

!

!



Ã ! # " t −τ V+  e c+1 Z ∞ c+1 c+1 00 0 − [c¯ g (φ(τ ))ψ(τ )h(τ ) − g¯ (φ(τ ))h(τ )] dτ e 1 2 1 + c t

+

t 1+c

e 1+c

Z t

"

Ã

De manera similar, se obtiene la derivada de Gato de F1 en respecto a ψ ÃÃ !

δψ F1

 

#

c−1 [c¯ g 00 (φ(τ )ψ(τ )h(τ ) − g¯0 (φ(τ ))h(τ )] dτ  . 1

−τ

e c−1

−∞

!



!

Ã

³³ ´ φ ψ

´

, c con

!

−τ φ V+ c  e c+1 Z ∞ c+1 c+1 0 , c (k) = e g¯ (φ(τ ))k(τ )dτ ψ 2 c + 1 t 1 t



!

Ã

 −τ c−1 0 e c−1 Z t c−1 g¯ (φ(τ ))k(τ )dτ . + e  1 1 + c −∞ t

Además, δφ F1

³³ ´ φ φ

´

, c , δψ F1

³³ ´ φ φ

´

, c son lineales como función de h y k ³ ´ φ

respectivamente y continuas como función de ψ , se tiene que la derivada de ³ ´ φ Frechet de F1 con respecto a ψ existe y se puede calcular con la expresión (3.112). ³³ ´ ´ También es fácil probar que DF ψφ , c es continua en un entorno de

(

³ ´ φ ψ ,c

).

Lema 11 Bajo las condiciones de la proposici´ on 4, DF erador invertible. Demostraci´ on. Como

ÃÃ

DF

!

!

φ0 , c = I − DF1 0 122

ÃÃ

!

!

φ0 ,c , 0

³³ ´ φ0 0

´

, c es un op-

basta probar que

Evaluando (3.112) en ÃÃ

DF1

!

φ0 ,c 0

° ÃÃ ! !° ° ° φ0 ° °DF1 , c °° < 1 . ° ° 0 ³³ ´ ´ φ0 0

!Ã !

h = δφ F1 k

, c obtenemos

ÃÃ

!

!

φ0 , 0 (h) + δψ F1 0

ÃÃ

!

!

φ0 , c (k) , (3.116) 0

entonces ° ÃÃ ! ! Ã !° ° φ0 h °° ° °DF1 ,c ° ° 0 k °

∞

° ° ° ° ÃÃ ! ! ÃÃ ! ! ° ° ° ° φ0 φ0 ° ° ° ≤ °δφ F1 , c (h)° +°δψ F1 , c (k)°° ° ° ° ° 0 0 ∞

. ∞

(3.117)

Obtengamos estimados para ° ÃÃ ! !° ° ° φ0 ° °δφ F1 , c °° ° ° 0

,

∞

Evaluando (3.113) en δφ F1

³³ ´ φ0 0

´

³³ ´ φ0 0

° ÃÃ ! !° ° ° φ0 ° °δψ F1 , c °° ° ° 0

.

∞

´

, c obtenemos

, c (h) =



!

Ã

Ã

!



 −τ −τ e c−1 Z t c−1 c+1 c−1 V+ g¯(φ0 )  e c+1 Z ∞ c+1 h(τ )dτ + h(τ )dτ − e e  c + 1 t 1 1 2 1 + c −∞ 

t

Ã

t

!

Ã

!



 −τ −τ c+1 e c−1 Z t c−1 c−1 V+ c¯ g 0 (φ0 )  e c+1 Z ∞ c+1 k(τ )dτ − k(τ )d(τ ) , + e e c + 1 t  2 1 c + 1 −∞ 1 t

t

entonces ° ° ÃÃ ! ! ° ° φ ° ° , c (h) ° °δφ F1 ° ° 0

∞

√

 t −τ V+ g¯ (φ0 ) √ 2 e c+1 Z ∞ c+1 ≤ c + 2c + 2 e dτ k h k∞  2 c+1 t 0

+ c2 − 2c + 2

t c−1

e 1+c

Z t −∞

123

τ

e c−1 dτ k h k∞

  

.

Calculando las integrales impropias y simplificando ° ° ÃÃ ! ! ° ° φ0 ° °δφ F1 , c (h)°° ° ° 0

≤ ∞

´ √ V+ g¯(φ0 ) ³√ 2 c + 2c + 2 + c2 − 2c + 2 k h k∞ . 2

Para 0 < c < 1 se tiene que √ √

c2 + 2c + 2 es creciente,

c2 + 2c − 2 es decreciente.

As´ı se tiene que ° ° ÃÃ ! ! ° ° φ0 ° °δφ F1 , c (h)°° ° ° 0

√ ´ V+ g¯0 (φ0 ) ³√ 5 + 2 k h k∞ . 2

(3.118)

√ ´ V+ g¯0 (φ0 )c ³√ 5 + 2 k k k∞ . 2

(3.119)

≤

∞

De manera similar se obtiene ° ° ÃÃ ! ! ° ° φ0 ° ° °δψ F1 , c (k)° ° ° 0

≤ ∞

Sustituyendo (3.118), (3.119) en (3.117), obtenemos ° ÃÃ ! ! Ã !° ¶ ° ³√ √ ´ µq φ0 h °° V+ g¯0 (φ0 ) ° 2 2 (1 + c) 5 + 2 k h k∞ + k k k∞ . °DF1 ,c °≤ ° 2 k ° 0

Como 0 < c < 1 y V+ α <

1 4

obtenemos

° ÃÃ ! ! Ã !° ° φ0 h °° ° °DF ,c ° < 1, ° k ° 0

que es lo que quer´ıamos probar. Finalmente, tenemos el siguiente teorema.

124

Teorema 6 Supongamos que en el sistema (3.74), la función g¯ tiene la forma 1 g¯(φ) = , α, β > 0 1 + e−αφ+β 1 y que αV+ < 14 . Entonces para cada ρ > max{φ0 , 2α }, existen c1 (ρ), c2 (ρ) con 0 < c1 (ρ) < c2 (ρ) < 1, tales que, para cada c, que cumpla 0 < c < c1 (ρ), la u ńica solución acotada de³este sistema, contenida en la bola cerrada β[0, ρ], ´ φ0 es la solución estacionaria 0 . Mientras que, para cada c, con c1 (ρ) ≤ c ≤ c2 (ρ), existe una solución acotada, no trivial, del sistema, la cual también es una solución homocl´ınica, y a la que le corresponde una onda viajera de la ecuaci´ on (3.59) con velocidad v. La velocidad real de traslaci´ on de esta onda viajera es cV+ . Adem´ as c1 (ρ), c2 (ρ) son las ra´ıces positivas de los polinomios cuadráticos definidos en (3.100), (3.101).

Demostraci´ on. La existencia de c1 (ρ), fue demostrada en la proposición 4. Por los lemas 10 y 11 el operador integral F definido en (3.101) cumple con ³³las´ condi´ ciones del teorema de la función impl´ıcita en un entorno de φ00 , c , es decir, F

³³ ´ φ ψ

´

, c es de clase C 0 en un entorno de

³³ ´ φ0 0

´

,0 , F

³³ ´ φ0 ³0

´

,c = 0 y

´ φ0 0 ³ y´ c resφ pectivamente tales que, para cada c∗ ∈ V , existe un u ńico elemento ψ ∈ U

DF

³³ ´ φ0 0

´

, c es invertible. Entonces existen entornos U y V de

no trivial tal que

ÃÃ !

F

!

φ , c∗ = 0 . ψ

³ ´ φ

Entonces ψ es una solución no trivial, acotada de (3.80) y por el teorema ³ ´ φ 3, ψ es una solución acotada, no trivial del sistema (3.74), la cual, por el lema 6, es una solución homocl´ınica del sistema (3.74). 3.4 ONDAS ESTACIONARIAS EN LA CORTEZA CEREBRAL.

125

En esta sección estudiaremos las oscilaciones peque˜ nas –con la misma frecuencia angular w para cada x ∈ R– de la ecuación ∂2u ∂t2

2

2 ∂ u + 2λ+ v+ ∂u + (λ+ v+ )2 u − v+ = ∂t ∂x2

h

V+ λ+ v+ λ+ v+ g¯(u) + g¯0 (u) ∂u ∂t

i

(3.120)

alrededor de un estado de reposo de la misma (una solución u(x, t) ≡ u0 ). Este problema corresponde a considerar una banda circular de ancho infinitamente peque˜ no de corteza cerebral. Lo que queremos saber es si en nuestro modelo existen ondas estacionarias que se correspondan con los ritmos básicos del cerebro. Para estudiar este problema consideraremos una solución u(x, t) ≡ u0 de la ecuación (3.120). Tomaremos la parte lineal de esta ecuación alrededor de u0 y en la ecuación lineal buscaremos soluciones de la forma u(x, t) = u0 + L(x)eiwt , con L(x) una función periódica real, de periodo 2π y L(x) peque˜ no.

126

(3.121)

Observaci´ on. En realidad, lo que buscamos son soluciones de la ecuación lineal de la forma u(x, t) = u0 + A(x) cos wt donde A(x) es una función periódica real, de periodo 2π. Pero como se tiene una ecuación lineal, basta encontrar soluciones de la forma (3.121) y después tomar la parte real de estas soluciones. Ahora, si u(x, t) ≡ u0 es una solución de (3.120) la “parte lineal” de esta ecuación alrededor de u0 es 2 ¯ ∂ 2 u¯ ∂ u¯ 0 2 0 2 ∂ u + λ v (2 − V g ¯ (u )) + [(λ v ) (1 − V g ¯ (u ))]¯ u − v =0 + + + o + + + 0 + 2 ∂t ∂t ∂x2 (3.122) donde

u¯ = u − u0

Teorema 7 Sea u0 ∈ (0, V+ ) tal que, g¯0 (u0 ) = V2+ , supongamos además que 1 − λ2+ > 0. Entonces la ecuaci´ on (3.122) tiene soluciones de la forma u¯(x, t) = L(x)eiwt , con L(x) una función peri´ odica real, de periodo 2π. Más concretamente, si µq

w = v+

¶

1−

λ2+

,

(3.123)

u¯(x, t) = A1 cos xeiwt

(3.124)

entonces con A1 constante, es una solución de (3.122).

127

Observaci´ on 2 , v+

Como g¯0 (u0 ) =

esto nos garantiza que

u(x, t) ≡ u0 es solución de la ecuación (3.120). Demostraci´ on. Sea w como en (3.123). Obtengamos las derivadas de (3.124) ∂u ¯ ∂t

= (iw)A1 cos xeiwt ,

∂2u ¯ ∂t2

= −w2 A1 cos xeiwt ,

∂u ¯ ∂x

= −A1 sen xeiwt ,

∂2u ¯ ∂x2

= −A1 cos xeiwt .

Sustituyendo estas derivadas en la ecuación (3.122) obtenemos −w2 A1 cos xeiwt + λ+ v+ (2 − V+ g¯0 (u0 ))iwA1 cos xeiwt 2 +[(λ+ v+ )2 (1 − V+ g¯0 (u0 )]A1 cos xeiwt + v+ A1 cos xeiwt = 2 2 (λ2+ − 1)v+ A1 cos xeiwt + (λ+ v+ )2 A1 cos xeiwt + v+ A1 cos xeiwt = 0 ,

ya que g¯0 (u0 ) =

2 V+

q

y w = v+ 1 − λ2+ .

128

CONCLUSIONES. El principal modelo que hemos estudiado, a pesar de ser un modelo que tiene algunas simplificaciones importantes : modelo unidimensional sin curvatura, considerar una distribución de conexiones que depende u ńicamente de la distancia entre columnas y con escalas de decrecimiento que no dependen de las regiones corticales consideradas, las complicaciones intra e intercorticales están incluidas en las funciones R± , no tomar en cuenta expl´ıcitamente la actividad neuroqu´ımica, etc., nos ha permitido obtener un conjunto de resultados que pueden ser sujetos a ser estudiados experimentalmente y que en principio coinciden cualitativamente con algunos datos experimentales ya conocidos. Se debe destacar la relativa simplicidad del modelo obtenido y la posibilidad inmediata de incluir, aunque sea parcialmente, algunos elementos no considerados en el modelo principal como por ejemplo, escalas de decrecimiento que dependan de las regiones de corteza considerada, que la función de generación de potenciales de acción g¯ dependa de la columna y de la función de activación (ver el modelo de la sección 2.8) además de la posibilidad de construir un modelo matemático bidimensional de la corteza cerebral que además incluya, aunque sea parcialmente, la geometr´ıa real de la corteza cerebral. Obviamente el estudio matemático de este problema resulta mucho más complejo. Consideramos que todo lo anteriormente enunciado, son los elementos más valiosos de la tesis. La inclusión en el modelo de la actividad neuroqu´ımica nos parece que es un problema sumamente complejo pero en primera aproximación, considerando que la liberación de neurotransmisores es proporcional a la actividad sináptica, podr´ıa ser incluida en el mismo. A continuación enlistamos las principales conclusiones de nuestro trabajo. 1. En este trabajo se han continuado desarrollando las ideas de Nu˜ nez, [86], [87], [89], [92], en al menos dos direcciones : se han introducido modificaciones y agregado nuevos elementos al modelo. Además el tratamiento matemático del modelo ha sido diferente y más completo que el que le dió Nu˜ nez. 2. En cuanto a resultados, se ha demostrado, que bajo ciertas condiciones,(teorema 6) existen ondas viajeras de actividad cortical que se propagan por la corteza cerebral con una velocidad de propagación menor que la velocidad de propagación de la actividad excitatoria. 129

También se ha demostrado,(teorema 2) que la actividad cortical de una rebanada de corteza se puede acoplar con la de una entrada aferente periódica sostenida. Esto ha sido demostrado experimentalmente [25], [78]. Además en estos trabajos también se concluye, experimentalmente, el papel controlador de la actividad inhibitoria en la corteza. Lo cual probamos en el teorema 3. 3. Los resultados obtenidos,(teorema 2, 6, 7) corresponden a estados fisiológicos muy especiales, ésto u ´ltimo se refleja en las restricciones impuestas a la función de activación neuronal. Consideramos que es necesario estudiar, posteriormente, la significación fisiológica de estas restricciones. 4. Hemos obtenido un modelo (sección 2.8), que en conjunto con estudios del problema inverso electroencefalográfico, nos permitirá estudiar, posteriormente, el diagrama de conexiones funcionales y de la función de activación neuronal en un estado fisiológico determinado. 5. El modelo compartimentado (2.8) puede ser usado para el estudio de otros sistemas del cerebro, como el sistema de los ganglios basales y el sistema l´ımbico. 6. Posteriormente será necesario realizar simulaciones computacionales de las soluciones del modelo con diversos valores de los parámetros involucrados y compararlos con mediciones experimentales.

130

Bibliograf´ıa [1] Abbott, L. F. [1990] : A network of oscillators, J. Phys. A. Math. Gen. 23 : 3835-3859. [2] Abbot, L. F., F., Kepler, T. B. [1990] : Model neurons from Hodgkin-Huxley to Hopfield In : Garrido (ed) Statistical mechanics of neural networks. Springer, Berlin, Heidelberg, New York. [3] Abeles, M. [1982] : Local Cortical Circuits in Studies in Brain Function 6, Springer Verlag, New York. [4] Amick, C, J, [1991] : A problem in neural networks. Proceedings of the Royal Society of Edinburgh. 118 : 225-236. [5] Amari, S. [1972] : Competitive and cooperative aspects in dynamics of neural excitation and self-organization, IEEE Trans on computers. [6] Amari, S. [1982] : A mathematical theory of self-organizing nerve sytems. In Proceedings of Biomathematics : Current Status and Perspectives. Ricciardi, L. M., Scott, A., ed. [7] Amari, S., Arbib, M. A (eds) [1982] : Competition and Cooperation in Neural Networks, Springer Verlag, New York. [8] Amari, S., Takeuchi, A. [1978] : Mathematical theory on formation of category detecting nerve cells. Biol. Cybernetics 29. 127-136. [9] Amir, A. [1994] : Uniqueness of the generators of brain evoked potential maps, IEEE Trans. on Biomedical-Engineering. Vol. 41,no.1-11. 131

[10] Andersen, P., Andersson, S.A. [1968] : Physiological basis of the alpha rhytm. Appleton-Century-Crofts, New York. [11] Anderson, J. A., Rosenfeld, E. [1988] : Neurocomputing Foundation of Research. Cambridge. The M.I.T. Press. [12] Beale, R.,Jackson, T. [1994] : Neural Computing : An Introduction, Institute of Physics Publishing, Bristol and Philadelphia. [13] Baumgärtel, H. [1985] : Analytic Perturbation Theory for Matrices and Operators. Birkhäuser Verlag Basel. Boston-Stuttgart. [14] Bogacz, J. y colaboradores. [1985] : Los potenciales evocados en el hombre. Ateneo, Buenos Aires. [15] Bower, J. M. [1992] : Modeling the nervous system. Trends in Neuroscience 15 : 411-412. [16] Buchholtz, F., Golowasch, J., Epstein, I. R., Marder, E. [1992] : Mathematical model of an identified stomatogastric neuron. J.Neurophysiol. 67 : 332-340. [17] Cardinali, D. P. [1992] Manual de Neurofisiolog´ıa. Ediciones D´ıaz de Santos, S.A. Madrid, Espa˜ na. [18] Chauvet, G. A. [1993a] : Non-locality in biological systems results from hierarchy. Application to nervous system. Journal Math. Biol. 31. 475-486. [19] Cichocki, A., Unbehauen, R. [1993] : Neural Networks for Optimazation and Signal Processing, John Wiley, New York. [20] Coro, F. [1996] : Fisiolog´ıa Celular, Un enfoque biof´ısico, Textos Cient´ıficos , BUAP, Puebla. [21] Cronin, J. [1987] : Mathematical aspects of Hodgkin-Huxley neural theory, Cambridge University Press. [22] Cuffin, B.N. [1987] : The role of model and computacional experiments in the biomagnetic inverse problem.Phys.Med.Biol. 132

[23] Daleckii, Ju. L., Krein, M. G. [1974] : Stability of Solutions of Differential Equations in Banach Space. American Mathematical Society, Providence, Rhode Island. [24] Destexhe, A., Contreras, D., Sejnowski, T. J., Steriade, M. [1994] : A model of spindle rhytmicity in the isolated thalamic reticular nucleus. Journal of Neurophysiology. Vol. 72, No. 2. 803-1818. [25] Dickson, C. T., Alonso, A [1997] : Muscarinic induction of Synchronous population activity in the entorhinal cortex. Journal of Neuroscience. 17 (17:6729-6744). [26] Doya, K., Selverston, A.I. [1993] : Dimension reduction of biological neuron models by artificial neural networks, Neural Comp. Vol. 6 no. 4 : 697-717. [27] Farkas, Miklos. [1994] : Periodics Motions, Springer Verlag, New York. [28] Fitzhugh, R. [1961] : Impulses and physiological stages in theoretical models of nerve membrane. Biophys, J. 1 : 445-466. [29] Fitzhugh, R. [1969] : Mathematical models of excitation and propagation in nerve. In Biological Engineering H.P. Schwann.Mc Graw Hill Publications, New York. [30] Flores, G. [1991] : Stability Analysis for the Slow Traveling Pulse of the Fitzhugh-Nagumo System, Siam J. Math. Anal. Vol. 22, no. 2, 392-399. [31] Freeman J. A., Skapura, D.M. [1993] : Redes Neuronales. Algoritmos, aplicaciones y técnicas de programación. Addison-Wesley Iberoamericana Wilmington, Delaware. U.S.A. [32] Freeman, W.J.[1975] : Mass Action in the Nervous System, Academic Press, New York. [33] Freeman, W.J. [1979] : Nonlinear dynamics of paleocortex manifested in olfactory EEG, Biological Cybernetics. 35 : 21-37.

133

[34] Freeman, W. J. [1992] : Predictions on neocortical dynamics derived from studies in paleocortex. In : Induced rhythms of the brain, ed. E.Basar, T. H. Bullock. Birkhauser Boston Inc. 183-201. [35] Freeman, W. J. [1994] : Neural mechanisms underlying destabilization of cortex by sensory input. Physica D 75 : 151 - 64. [36] Fukushima, K. [1975] : Cognitron : a self organizing multilayered neural network, Biolog. Cybernetics. Vol. 20 : 121-136. [37] Fukushima, K. [1980] : Neocognitron : a self organizing neural network model for a mechanism of pattern recognition unaffected by shift in position, Biolog. Cybernetics. Vol. 36,193-202. [38] Galusinski, C. [1998] : Existence and Continuity of Uniform Exponential Attractors of the Singularity Perturbed Hodgkin-Huxley System. Journal of Differential Equations. 144 : 99-169. [39] Golomb, D., Guckenheimer, J., Gueron, S. [1993] : Reduction of a Channel Based Model for a Stomatogastric Ganglion LP Neuron. Biol. Cybern. 68 : 129-137. [40] Golomb, D., Rinzel, J. [1993] : Dynamics of globallly coupled inhibitory neurons with heterogeneity., Physical Review. E. Vol. 48,no. 6. [41] Golomb, D., Rinzel, J. [1994] : Clustering in globally coupled inhibitary neurons. Physics D. 72 : 259-282. [42] Golomb, D., Wang, X. J., Rinzel, J. [1994] : Synchronization Properties of Spindle Oscillations in a Thalamic Reticular Nucleus Model. J. Neurophysiology. Vol. 72, no.3, 1109-1126. [43] Grossberg, S. [1976] : Adaptative pattern classification and universal recording. Biolog. Cybernetics. Vol. 23 : 121-134 y 187-202. [44] Guckenheimer, J., Holmes, P. [1983] : Nonlinear oscillations dynamical systems, and bifurcations of vector fields. Springer, Berlin, Heidelberg, New York.

134

[45] Hale, J. K. [1963] : Oscillations in Nonlinear Systems. Dover, New York. [46] Hale, J. K. [1969] : Ordinary Differential Equations. Wiley Interscience. New York, London, Sydney, Toronto. [47] Hale, J. ; Kocak, H. [1991] : Dynamics and Bifurcations, Springer Verlag, New York. [48] Hastings, S. [1975/1976] : On traveling wave solutions of the Hodgkin-Huxley equations, Arch. Rational Mech. Anal. 60 : 229257. [49] Hastings, S. [1976] : On existence of homoclinic and periodic orbits for the Fitzhugh-Nagumo equations,Quart J. Math.27 : 123-134. [50] Hebb, D. O. [1949] : The organization of behavior. Wiley, New York. [51] Heiner, L. ; Záborszky, L. [1989] : Neuroanatomical Tract–Tracing Methods 2. Plenum Press. New York and London. [52] Heller, L. [1990] : Return Current in encephalography Variationals Principles, Biophysical Journal. Vol. 57 : 601-606. [53] Herman G. T. [1980] : Image reconstruction from projections. The fundamental of computerized tomography. Academic Press. [54] Herman G. T., Tuy H. K., Langenberg K. J., Sabatier P.C. [1985] : Basic methods of tomography and inverse problems. Malvern Physics Series. [55] Hilera, J. R., Mart´ınez, V. J. [1995] : Redes neuronales artificiales. Fundamentos, modelos y aplicaciones. Addison-Wesley Iberoamericana. Madrid, Espa˜ na. [56] Hinton, G. Ackley, D. Sejnowski, T. [1984] : Boltzman machines : Constraint satisfaction networks than learn. Carnegie-Mellon University, Departament of Computer Science Technical Report. (CMU-CS-84-119). 135

[57] Hodgkin, A. L. ; Huxley, A.F. [1952] : A quantitative description of membrane current an its application to conductance and excitation in nerve, J. Physiol. Lond.117 : 500-544 [58] Holden, A. V., Hyde, J., and Muhamad, M. A. [1991] : Equilibria, Periodicity, Bursting and Chaos in Neural Activity, in Pasemann, F., and Doebner, H. D. (eds), World Scientific Pub., Singapore : 97-187. [59] Hopfield, J. J. [1982] : Neural networks and physical systems with emergent collective computational abilities, Proc. Natl.Acad.Sci. 79 : 2554-58. [60] Hopfield, J. J [1984] : Neurones with graded response have collective computational propierties like those of two states neurones. Proccedings of National Academy of Science. 81 : 3088-92 [61] Hoppensteadt, F. [1986] : An Introduction to the Mathematics of Neurons, Cambridge University Press. [62] Hoppensteadt, F.C. [1993] : Analysis and Simulation of Chaotic Systems, Springer Verlag, New York. [63] Hoppensteadt, F., Izhikevich, E. M. [1997] : Weakly Connected Neural Networks. Springer Verlag. New York, Inc. [64] Ingber, L. ; Nu˜ nez P. [1990] : Multiple scales of statistical physics of Neocortex. Aplication to the Electroencephalography. Mathematical and Computer Modelling. [65] Ingber, L. [1991] : Statistical mechanics of neurocortical interactions : A scaling paradigm applied to electroencephalography., Physical Review. A 44 : 4017-60. [66] Ingber, L. [1994] : Statistical Mechanics of Multiple Scales of Neocortical Interactions. [67] Ingber, L. [1994] : Statistical Mechanics of Neocortical Interactions. : Path-integral evolution of Short-term memory. Physical Review E 49 : 4652-64. 136

[68] Ingber, L. [1995] : Statistical mechanics of multiple scales of neocortical interactions. In : Neocortical dynamics and human EEG. rhythms, ed. P.L. Nu˜ nez. Oxford University Press. [69] Jahnsen, H. ; Llinas, R. [1984a] : Electrophysiological properties of guinea-pig thalamic neurons : An in vitro studio. J. Physiol. (London). 349 : 205-226. [70] Jahnsen, H. ; Llinas, R. [1984b] : Ionic basis for the electroresponsiveness and oscillatory properties of guinea-pig thalamic neurons in vitro. J. Physiol. (London) 349 : 227-2447. [71] Kandel, E. R., Schwartz, J. H. ; Jessell, T. M. [1997] : Neurociencia y conducta. Prentice Hall. Madrid. [72] Katznelson, R. D. [1981] : Normal modes of the brain : Neuroanatomical basis and a physiological theoretical model. In : Electric fields of the brain : The Neurophysics of EEG. P.L. Nu˜ nez. Oxford University Press. New York. [73] Kepler, T. B., Abbot, L. F., Marder, E. [1992] : Reduction of conductance based neuron models. Biol. Cybern. 66 : 381-387. [74] Kevorkian, J., Cole, J. D. [1996] : Multiple Scale and Singular Perturbation Methods, Springer-Verlag, New York, Inc. [75] Lopes da Silva, F. H. [1991] : Neural mechanisms underlying brain waves : From neural membranes to networks electroencephagr. Clin. Neurophysiol. 79 : 81-93. [76] Lopes da Silva, F. H., van Rotterdam, A., Barts, P., van Heusden, E., Burr, B. [1976] : Models of neurones populations the basic mechanism of rhytmicity. Progress in Brain Research. 45 : 282308. [77] Lopes da Silva, F. A. ; Hoeks, A., Smiths, H. and Zettergerg L. H. [1974] : Model of Brain Rhytmic Activity, Kybernetik. 15 : 27-37. [78] Lukatch, H. S.; Mac Iver, B. [1997] : Physiology, Pharmacology and Topography of cholinergic neocortical oscillations in vitro. Journal of Neuroscience. 2427-2444. 137

[79] Mac. Gregor, R. J. [1987] : Neural and Brain Modeling. Academic Press, Inc. San Diego. Cal. [80] Mc Culloch, W. S. ; Pitts, W. H. [1943] : A logical calculus of ideas immanent in nervous activity, Bull. Math. Biophys. 5 : 115-33. [81] Menon, V. [1991] : Population oscillations in neuronal groups. International Journal of Neural Systems. Vol. 2, no. 3, 237-62. [82] Mischenko, E. F.; Rozov, N. Kh. [1980] : Differential Equations with Small Parameters and Relaxation Oscillations, Plenum, New York. [83] Mountcastle, V. B. [1978] : An organizing principle for cerebral function : the unit module and the distributed system. In : The Mindful Brain : Cortical Organization and the Group Selective theory of Higher Function, ed. G. M. Edelman and V. B. Mountcastle. MIT Press, Cambridge. [84] Mountcastle, V. B. [1998] : Perceptual Neuroscience the Cerebral Cortex. Harvard University Press. Cambridge, Massachusetts, and London, England. [85] Natterer, F. [1986] : The mathematics of Computerized Tomography. Teubner, Stuttgart. [86] Nu˜ nez, P. L. [1981] : A study of origins of the time dependencies of scalp EEG : I. Theoretical basis, IEEE Trans. Bio-Med Eng. 28 : 271-280, 1981b. [87] Nu˜ nez, P. L. [1981] : A study of origins of the time dependencies of scalp EEG: II. Experimental support of theory,IEEE Trans. BioMed Eng. 28 : 281-288, 1981c. [88] Nu˜ nez, P. L. [1988] : Global Contribution to cortical dynamics : theoretical and experimental evidence for standing waves phenomena. In : Dynamics of Sensory and Cognitive Processing by the Brain. (E. Baser, Ed.). Springer Verlag. New York.

138

[89] Nu˜ nez, Paul L; Katznelson, R. D. [1981] : Electric Fields on the Brain, the Neurophysics of EEG, Oxford University Press, New York. [90] Nu˜ nez, P. [1989] : Generation of Human EEG by a combination of long and short range neurocortical interactions, Brain topography. 1 : 199-215. nez, P., Srinivasan, R. [1993] : Implications of recording [91] Nu˜ strategy for estimates of neurocortical dynamics with electroencephalography chaos. 3 : 257-66. [92] Nu˜ nez, P. L. [1995] : Neocortical Dynamics and Human EEG Rythms, Oxford University Press. [93] Plonsey, R. [1969] : Bioelectric Phenomena, Mc Graw Hill, New York. [94] Petsche, H., Sterc, J. [1968] : The significance of the cortex for the traveling phenomenon of brain waves. Electroencephalography and Clinical Neurophysiology. 25 : 1222. [95] Rall, W., Shepherd, G. M. [1968] : Theoretical reconstruction of field potentials and dendrodendritic synaptic interactions in olfactory bulb. J. Neurophysiol. 31 : 884-915. [96] Ramacher, U. [1993] : Hamiltonian Dynamics of Neural Networks. Neural Networks. 6 : 547-557. [97] Rinzel, J. ; Ermentrout, G. B. [1989] : Analysis of neural excitability and oscillations. In Methods of neural Modelling : From Synapses to networks (C. Koch and I. Segev. Eds), Cambridge, MA : MIT Press. [98] Rosenblatt, F. [1957] : The perception a perceiving and recognizing automation, , Cornell Aeronautical Laboratory Report. [99] Rose, R. M. ; Hindmarsh, J. L. [1985] : A model of a thalamic neuron, Proc. R. Soc. Lond.

139

[100] Rush, M. E., Rinzel, J. [1994] : Analysis of bursting in a thalamic neuron model. Biol. Cybern. 71, 281-291. [101] Sánchez, F., Maini, P. K. [1997] : Travelling wave phenomena in non-linear diffusion degenerate Nagumo equations. Journal of mathematical biology, 35 : 713-728. [102] Sánchez, F., Maini, P. K. [1997] : Travelling wave phenomena in some degenerate reaction diffusion equations. Journal of differential equations. Vol. 117, no. 2 : 281-319. [103] Sarvas, J. [1987] : Basic Mathematical and Electromagnetic Concepts of the Biomagnetic Inverse Problem. Phys. Med. Biol. Vol.32, no.1 : 11-22. [104] Silberstein, R. B. [1994] : Neuromodulation of neocortical dynamics. In : Neocortical dynamics and human EEG rhythms. ed. P. L. Nu˜ nez. Oxford University Press. [105] Silberstein, R. B. [1994] : Steady-State Visually Evoked Potentials, Brain Resonances, and Cognitive Processes In : Neocortical dynamics and human EEG rhythms. ed. P. L. Nu˜ nez. Oxford University Press. [106] Smith, C. U. M. [1982] : El cerebro. Sexta ed. Alianza Universidad Madrid. [107] Somogyi, P. ; Tamás, G. ; Lujan, R., Buhl, E. [1998] : Salient features of synaptic organisation in the cerebral cortex. Brain research reviews, 26 : 113-135. [108] Temam, R. [1998] : Infinite-Dimensional Dynamical Systems. 2 ed. Springer Verlag, New York, Berlin, Heidelberg. [109] Tijonov A.N., Goncharsky A.V. [1987] : Problemas mal planteados de las ciencias naturales. (Editorial de la Universidad Estatal de Mosc´ u). [110] Tikhonov, V. S. [1985] : Differential Equations. Springer Verlag. Berlin , Heidelberg. 140

[111] Troy, W. C. [1978] : The bifurcation of periodic solutions in the Hodgkin-Huxley equations. Quart. Appl. Math., Vol. 36, 73-83. [112] Tyson, J. J. ; Keener, J. P. [1988] : Singular Perturbation Theory of Traveling Waves in Excitable Medic (a review). Physica.D 32 : 327-361. [113] van Rotterdam, A., Lópes da Silva, F. H., van den Ende, J., Viergeuer, M. A., Hermans, A. J. [1982] : A model of the spaticaltemporal characteristics of the alpha rhytm. Bull. Math. Biol. 4e : 283-305. [114] Wang, X. J. ; Rinzel, J.[1993] : Spindle rhythmicity in the reticularis thalami nucleus : synchronization among mutually inhibitory neurons. Neuroscience. 53 : 899-904. [115] Wang, X. J. [1994] : Multiple dynamical modes of thalamic relay neurons : Rhythmic bursting and intermittent phase-locking. Neuroscience. Vol. 59, no.1, 21-31. [116] Weiss, T. F. [1996] : Celular biophysics, Vol. 1, Transport. MIT Press, Cambridge, MA. [117] Weiss, T. F. [1996] : Celular biophysics, Vol. 2, Electrical properties. MIT Press, Cambridge, MA. [118] Widrow, B., Hoff, M. [1960] : Adaptative Switching Circuits IREWESCON Convention Record, Part 4, 96-104. [119] Wilson, H. R., Cowan, J. D. [1972] : Excitatory and inhibitory interaction in localized populations of model neurons. Biophys. Journal 12. 1-24. [120] Wilson, H. R., Cowan, J. D. [1973] : A mathematical theory of the functional dynamics of cortical and thalamic nervous tissue. Kybernetik 13. 13-80. [121] Wright, J.J. ; Liley, T. J. [1996] : Dynamics of the Brain at Global and Microscopic scales : Neural Networks and the EEG, Behavioral and the brain sciences. 19 : 285-320. 141

Modelos matem´aticos para el estudio de la activaci´on de

Recommend Documents