Download Estudio de dinámica de una red neuronal con

Estudio de dinámica de una red neuronal con modelo de
Hindmarsh-Rose
Carlos Garcı́a Argos (carlos.garciaargos@estudiante.uam.es)
22 de noviembre de 2008
Resumen
Blah
Índice
1. Introducción
1
2. Modelo de Hindmarsh-Rose
2.1. Implementación del modelo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1
2
3. Resultados del estudio mediante simulaciones
3.1. Simulaciones con una neurona aislada . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.2. Simulaciones con una red de 3 neuronas acopladas . . . . . . . . . . . . . . . . . . . . . . . . .
4
4
4
4. Conclusiones
6
1.
Introducción
Esta práctica está basada en el trabajo previo discutido en [1], y consiste en la implementación de un modelo de
red neuronal simple, formada por 3 neuronas interconectadas por conductancias que simulan una sinapsis eléctrica.
En este trabajo se presenta la implementación llevada a cabo para realizar las simulaciones en software de
cálculo numérico. Concretamente, el software aquı́ incluido se ejecuta tanto en GNU Octave [2] como en MATLAB
[3].
Las simulaciones se realizan en base al modelo de Hindmarsh-Rose [4], el cual se ha resuelto mediante métodos
numéricos para EDOs de Runge-Kutta.
2.
Modelo de Hindmarsh-Rose
Las ecuaciones que gobiernan la dinámica del modelo que se ha simulado son:
dxi (t)
dt
dyi (t)
dt
dzi (t)
dt
= yi (t) + 3x2i (t) − x3i (t) − zi (t) + ei −
X
gij (xi (t) − xj (t))
j6=i
=
1 − 5x2i (t) − yi (t)
= µ (−zi (t) + S (xi (t) + 1,6))
Se trata de unas ecuaciones diferenciales acopladas, que forman un sistema de ecuaciones no lineales que
puede ser resuelto mediante técnicas numéricas. En este caso, se ha optado por un método de Runge-Kutta de
quinto orden con 6 evaluaciones de la función (RK56). Existen dos variantes de dicho método, que únicamente se
diferencian por los pesos de los valores de la función evaluados.
La decisión de optar por un método de Runge-Kutta de orden elevado, frente a un método más simple como
el método de Euler o un ODE23 se debe a la inestabilidad del primero y a la falta de precisión para el segundo,
1
que requiere valores de paso demasiado pequeños para ser utilizado de forma práctica en un lote de muchas
simulaciones.
En la siguiente sección se incluye el código fuente de los programas implementados que pueden utilizarse para
la verificación de los resultados que se presentan en la sección 3.
Los valores de los parámetros µ y S se fijan a:
µ = 0,0021
S=4
En todos los casos. Por otro lado, el parámetro ei gobierna el régimen de la neurona:
Régimen regular: ei = 3
Régimen caótico: ei = 3,281
2.1.
Implementación del modelo
2.1.1.
hr red.m
Este archivo contiene las ecuaciones de Hindmarsh-Rose y los parámetros de cada una de las neuronas. Es
posible adaptarlo para simular una red más grande o más pequeña, indicando el número de neuronas en la variable
N y modificando los array de los parámetros µ, S y ei (que pueden configurarse diferentes para cada neurona),
ası́ como el array de conductancias entre las neuronas.
Para llevar a cabo las diferentes simulaciones de la práctica, es necesario modificar los valores de los vectores
E y G. Además, para simular una neurona aislada puede modificarse el valor de N a 1 en este archivo y en el
siguiente.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Neurociencia Computacional 1
% Practica 2: Modelo de red - Ecuaciones del modelo HR
% Carlos Garcia Argos (carlos.garciaargos@estudiante.uam.es)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function dydt = hr_red(t,x)
S = 4;
mu = [0.0021, 0.0021, 0.0021];
E = [3.28100,3.28100,3.28100];
% N = numero de neuronas
N = 3;
g1 = 0.03;
g2 = 0.03;
G = [0
g2
g1
g2
0
g1
g1;
g1;
0];
for i=0:N-1
% Suma de los estimulos que vienen de las otras neuronas de la red
suma = 0;
for j=1:N
suma = suma + G(i+1,j)*(x(i+1) - x(j));
end;
% Las tres ecuaciones del modelo
dydt(3*i+1,1) = x(3*i+2) + 3*x(3*i+1)ˆ2 - x(3*i+1)ˆ3 - x(3*i+3)
+ E(i+1) + suma;
dydt(3*i+2,1) = 1 - 5* x(3*i+1)ˆ2 - x(3*i+2);
dydt(3*i+3,1) = mu(i+1)*(-x(3*i+3) + S*(x(3*i+1) + 1.6));
end;
2
2.1.2.
prac2hr red.m
Este programa simplemente hace la llamada al método de Runge-Kutta con los vectores de tiempo y de condiciones iniciales correspondientes. Se ha tomado un vector de condiciones iniciales nulo, a pesar de que serı́a
posible obtener mejores resultados con uno mejor adaptado, por carecer de información sobre las posibles condiciones iniciales en este escenario.
En este archivo también se incluye el parámetro N que indica el número de neuronas a simular, solamente para
establecer el vector de condiciones iniciales.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Neurociencia Computacional 1
% Practica 2: Modelo de red
% Carlos Garcia Argos (carlos.garciaargos@estudiante.uam.es)
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Parametros de la simulacion: paso y tiempo final
dt = 0.05;
fin = 10000;
N = 3;
% Condiciones iniciales
yinit = zeros(1,3*N);
% Vector de tiempo
tspan = 0:dt:fin;
% Llamada al metodo de resolucion del SEDO
[t_plot,y_plot] = rk56a(’hr_red’, 0, fin, yinit’, dt);
2.1.3.
%
%
%
%
%
rk56a.m
Implementacion del Metodo de Runge-Kutta de quinto orden
Metodo A de orden 5 con 6 evaluaciones
Uso: [t,y,dy]=rk56a(funcion,inicio,fin,alfa,paso)
(C) - Carlos Garcia Argos, 20/VIII/2001
function [ts,ys,fs]=rk56a(f,inicio,fin,alfa,h)
if nargin<5
disp(’Faltan argumentos de entrada’);
return;
else
ts=inicio:h:fin;
ys(:,1)=alfa;
n=(fin-inicio)/h;
for i=1:n
k1=feval(f,ts(i),ys(:,i));
k2=feval(f,ts(i)+h,ys(:,i)+h*k1);
k3=feval(f,ts(i)+h,ys(:,i)+(k1+k2)*h/2);
k4=feval(f,ts(i)+h/4,ys(:,i)+h*(14*k1+5*k2-3*k3)/64);
k5=feval(f,ts(i)+h/3,ys(:,i)+h/96*(-12*k1-12*k2+8*k3+64*k4));
k6=feval(f,ts(i)+h*(3/4),ys(:,i)+h/64*(-9*k2+5*k3+16*k4+36*k5));
fs(:,i)=k1;
ys(:,i+1)=ys(:,i)+(7*k1+7*k3+32*k4+12*k5+32*k6)*h/90;
end;
3
end;
3.
Resultados del estudio mediante simulaciones
3.1.
Simulaciones con una neurona aislada
3.1.1.
Definición de las simulaciones
Las simulaciones para una neurona aislada que se han llevado a cabo han consistido en dos casuı́sticas:
1. Una sola neurona en régimen regular.
2. Una sola neurona en régimen caótico.
La longitud de las simulaciones se ha establecido de manera que sea posible observar varios ciclos de funcionamiento de las neuronas y que el tiempo de simulación esté aceptablemente acotado alrededor de una hora
en un procesador actual. Se ha optado por una simulación hasta 20000 unidades de tiempo (que en el modelo de
Hindmarsh-Rose es adimensional).
A continuación se presentan las gráficas de ambos casos, tanto el tiempo completo de simulación como una
muestra en un intervalo más reducido, de 5000 unidades de tiempo.
3.1.2.
Gráficas de resultados
La Figura 1 se muestra el resultado obtenido para el caso del régimen regular. Se observa que las ráfagas siguen
una distribución uniforme, con una distribución claramente periódica de dichas ráfagas.
(a) Tiempo total de la simulación
(b) Zoom
Figura 1: Resultado de la simulación para una neurona regular
Por el contrario, como se observa en la Figura 2 se ilustra el caso de una neurona en régimen caótico. En
este caso, las ráfagas no muestran ningún comportamiento periódico, sino que cada ráfaga dispara en un instante
diferente e independiente de las anteriores, además de tener duraciones diferentes.
3.2.
Simulaciones con una red de 3 neuronas acopladas
3.2.1.
Definición de las simulaciones
Para el caso de tener 3 neuronas acopladas mediante sinapsis eléctricas, se han simulado los mismos casos que
en [1], que se ilustran en la Tabla 1.
Para el caso de tener las neuronas desacopladas, ...
4
(a) Tiempo total de la simulación
(b) Zoom
Figura 2: Resultado de la simulación de una neurona caótica
Simulación
A
B
C
D
E
F
g1
0
0,01
0,03
0
0,01
0
g2
0
0,01
0,03
0,02
0,01
0,02
Régimen
Caótico
Caótico
Caótico
Caótico
Regular
Regular
Cuadro 1: Parámetros para las simulaciones con la red de 3 neuronas
5
3.2.2.
4.
Gráficas de resultados
Conclusiones
En esta práctica se ha hecho blah blah blah.
Referencias
[1] F.B. Rodriguez, P. Varona, R. Huerta, M.I. Rabinovich, Henry D.I. Abarbanel. Richer Network Dynamics of
Intrinsically Non-regular Neurons Measured through Mutual Information. International Work-Conference on
Artificial and Natural Neural Networks, 2001.
[2] GNU/Octave. http://www.octave.org.
[3] MATLAB. http://www.mathworks.es/products/matlab.
[4] J.L. Hindmarsh, R.M. Rose. A Model of Neuronal Bursting Using Three Coupled First Order Differential
Equations. Proceedings of the Royal society of London, 1984.
6