Sergio Valadez-Godínez, Universidad Politécnica de Pénjamo

Resumen

Este estudio aborda la posibilidad de simular el cerebro humano en una computadora, centrándose en la transmisión de información mediante potenciales de acción a nivel celular. Para esto, se analizó el tiempo necesario para simular un milímetro cúbico de la corteza cerebral y la corteza completa utilizando modelos neuronales de generación de pulsos eléctricos. Los modelos neuronales se describen mediante ecuaciones diferenciales ordinarias, y la simulación se realizó con el método numérico más simple. Se observó una relación de ley de potencias entre el tiempo de ejecución de la CPU y el tamaño del paso temporal. Los resultados indican que incluso con pasos temporales eficientes, simular la corteza cerebral completa en una computadora sería prohibitivo, llevando a tardar incluso varios siglos. La simulación a niveles más detallados se vuelve más desafiante, y se concluye que, por el momento, simular el cerebro humano en su totalidad es más restrictivo que la simulación de la corteza cerebral.

Palabras Clave

Simulación; Cerebro; Corteza cerebral; Modelo neuronal; Hodgkin-Huxley; Izhikevich; Leaky Integrate-and-Fire; Método Numérico; Paso temporal;  Ley de potencias;

1. Introducción

Para poder contestar a la pregunta que da título a este trabajo, se necesita comprender el nivel de organización y de estudio del cerebro. Además, a partir de lo anterior, se requiere delimitar los aspectos funcionales mínimos que se simularán y posteriormente realizar un cálculo del tiempo que se requerirá para realizar dicha simulación .

Los niveles de organización cerebral abarcan genes, moléculas, sinapsis, neuronas, redes, capas, mapas y sistemas (Churchland & Sejnowski, 1992). Cada nivel de organización puede ser estudiado mediante el uso de tres enfoques: análisis, organización y procesamiento (Churchland & Sejnowski, 1988). El enfoque de análisis comprende el nivel de teoría computacional, el nivel algorítmico y el nivel de implementación física (Marr & Poggio, 1977). Para facilitar nuestro trabajo, este se centrará en el enfoque de análisis y en el nivel organizativo de la neurona individual, que es una unidad indivisible y funcional básica dentro de una red neuronal (Zador, 2000). Además, nos centraremos únicamente en el nivel macroscópico de la electricidad y la bioquímica dentro de la neurona individual, omitiendo el nivel microscópico de los genes y moléculas. Una neurona individual puede ser abstraída en varios niveles, como lo son los modelos detallados de compartimentos, modelos de compartimentos reducidos, modelos de un solo compartimento, modelos en cascada y modelos de caja negra (Herz et al., 2006). En este sentido, utilizaremos los modelos de un solo compartimento ya que poseen las características mínimas que dotan a una neurona de funcionalidad básica (Abbott et al., 2016).

Una de las funcionalidades elementales dentro de una red es la transmisión de información de una neurona a otra a través de potenciales de acción, también conocidos como pulsos eléctricos (Bohte, 2004). Los potenciales de acción se generan cuando las señales eléctricas entrantes de otras neuronas provocan que el potencial de membrana intracelular supere un umbral de voltaje específico, determinado por las propiedades bioeléctricas y bioquímicas de la neurona (Debanne et al., 2011). El potencial de acción resultante representa un valor binario que será enviado a otra neurona. El envío de un potencial de acción puede ser representado por 1, y el no envío por 0 (Maass, 1997). En la literatura existen varios modelos de un solo compartimento para simular potenciales de acción a través de la simulación de la inyección de una corriente eléctrica. Uno de los más reconocidos es el modelo de Hodgkin y Huxley (HH) (Hodgkin y Huxley, 1952). Este modelo es el más realista y computacionalmente costoso de simular en una computadora (Izhikevich, 2004). Hay otros modelos con comportamientos menos realistas y menos costosos de simular, como el modelo Leaky Integrate-and-Fire (LIF) (Lapicque, 1907; Stein, 1965; Knight, 1972), (en español, Integración y disparo con fuga), y modelos intermedios, como el modelo de Izhikevich (IZH) (Izhikevich, 2003). Estos modelos se expresan como ecuaciones diferenciales ordinarias, como veremos más adelante.

Este estudio está limitado, entonces, a la simulación de una de las funciones básicas que realiza el cerebro a nivel celular, el envío de información a través de potenciales de acción o pulsos eléctricos. Debido a la enorme cantidad de neuronas que posee el cerebro, limitaremos nuestro análisis al cálculo del tiempo que requiere una computadora para simular un milímetro cúbico de la corteza cerebral y a la corteza en su totalidad. Un milímetro cúbico de la corteza consta de 10(4) neuronas (Gerstner et al., (2014) y la corteza en su totalidad de 10(11) neuronas (Izhikevich, 2007) aproximadamente. De esta manera, intentaremos responder a la pregunta de ¿Es posible simular el cerebro humano en una computadora?

Es imposible saber directamente el tiempo que se podría tardar una computadora en simular el cerebro humano, ya que tendríamos que esperar a que terminara la simulación, la cual podría tardar días, meses, años o siglos. Por ese motivo, utilizaremos una ecuación matemática que nos permita predecir el tiempo de la simulación.

2. Metodología

Se realizaron simulaciones de pulsos eléctricos con tres modelos (LIF, IZH, HH), y se analizó el tiempo de ejecución en relación con el tamaño del paso temporal, ajustándose a una función de potencia para predecir el tiempo de simulación en una computadora personal. Los tres modelos fueron simulados numéricamente usando el método de Euler (Euler, 1768), que es el método de solución numérica más simple (Butcher, 2016).

El modelo LIF se define mediante la siguiente ecuación diferencial (Lapicque, 1907; Stein, 1965; Knight, 1972):

donde V_m es el voltaje transmembrana de la neurona, C_m es la capacitancia de la membrana celular, R_m es la resistencia de la membrana celular al flujo de corriente e I_m (t) es la corriente eléctrica de estimulación. Cuando el voltaje transmembrana excede un valor de umbral fijo V_umbral, entonces se reinicia a un nuevo valor V_reposo, es decir,

A continuación, se representa el modelo HH (Hodgkin & Huxley, 1952). Este modelo les permitió a Hodgkin y Huxley ganar el premio Nobel en fisiología y medicina en 1963 por describir la dinámica de generación y envío del pulso eléctrico entre neuronas. El modelo HH consta de varias ecuaciones diferenciales ordinarias acopladas:

en donde V_m es el voltaje transmembrana, C_m=1 µF/cm² es la capacitancia de la membrana celular e I_m (t)) es la corriente de estimulación. Los términos en el lado derecho de la Ec. (3) son las corrientes iónicas de potasio (K^+), sodio (〖Na〗^+)y cloro (〖Cl〗^-), respectivamente. E_K=-72.1 mV, E_Na=52.4 mV y E_L=-49.2 mV son los potenciales de reversa o de reposo de las corrientes iónicas. g ̅_K=36 mS/cm², g ̅_Na=120 mS/cm² y g_L=0.3 mS/cm² son las conductancias máximas de las corrientes iónicas del potasio, sodio y fuga, respectivamente. n^4 y m^3 h son las fracciones de canales iónicos de potasio y sodio abiertos en la membrana celular, respectivamente. La Ec. (4) describe la cinética de las partículas, donde x es la partícula n, m o h, según corresponda. α_x es la constante de velocidad que describe la cinética de la partícula x desde el estado cerrado hasta el estado abierto. Mientras tanto, β_x es la constante de velocidad de la cinética de una partícula x desde el estado abierto hasta el estado cerrado. (1 – x) es la probabilidad de que una partícula x esté en posición cerrada. Las constantes de velocidad α_x y β_x fueron ajustadas a datos experimentales por Hodgkin y Huxley como:

El modelo de Izhikevich (IZH) es una reducción dimensional del modelo de Hodgkin-Huxley (HH). Esto se logra mediante el análisis del plano de fase de la dinámica del modelo HH (Izhikevich, 2007). El resultado final es un modelo que combina la plausibilidad biológica del modelo HH y la eficiencia computacional del modelo LIF. El modelo IZH se define como (Izhikevich, 2003):

en donde V_m es el potencial transmembrana, I_m (t) es la corriente de estimulación, n es una variable de recuperación que representa la corriente de activación del potasio y la corriente de inactivación de los iones de sodio al mismo tiempo, a es la escala de tiempo y b es la sensibilidad de n a las fluctuaciones subumbrales del potencial transmembrana. Cuando el voltaje transmembrana alcanza un valor máximo de 30 mV, se reinician tanto el voltaje transmembrana V_m como la variable de recuperación n, es decir:

Cada uno de los modelos fue configurado con los parámetros que se mencionan en (Valadez-Godínez et al., 2017a). El LIF fue estimulado por una corriente de entrada constante de I_m (t)=1.3 nA para producir un tren de pulsos de aproximadamente 70 pulsos por segundo (vea la Fig. 1a). Los demás parámetros fueron: C_m=10nF, R_m=1 MΩ , V_umbral=1 mV y V_reposo=0 mV. El modelo IZH fue estimulado por una corriente constante de I_m (t)=12.5 μA/〖cm〗^2 para producir aproximadamente 70 potenciales de acción por segundo (vea la Fig. 1b). Los valores de las constantes fueron a=0.02, b=0.2, c=-65 y d=2. El modelo HH fue estimulado por una corriente constante de I_m (t)=11 μA/〖cm〗^2 para generar aproximadamente 70 pulsos por segundo (vea la Fig. 1c).

Figura 1. Potenciales de acción generados por los modelos a) LIF, b) IZH y c) HH. Se aplicó una corriente de estimulación constante para generar una frecuencia de pulsos de aproximadamente 70 picos por segundo.

Como se puede observar, los modelos para generar pulsos eléctricos están representados mediante ecuaciones diferenciales ordinarias. Es más, no se conoce una solución analítica exacta para los modelos IZH y HH (Izhikevich, 2007, Gerstner et al., 2014), mientras que el modelo LIF si posee dicho tipo de solución (Lapicque, 1907).  Por ese motivo, las ecuaciones de los modelos deben resolverse utilizando algún método numérico a través de una secuencia discreta de instantes o pasos de tiempo. Este es el punto clave para predecir la simulación, ya que el método numérico y el tamaño del paso de tiempo impactan en el tiempo y en la precisión de la simulación (Butcher, 2016). Existen muchas familias de métodos numéricos que se pueden utilizar para la simulación. No obstante, los más comúnmente utilizados son los métodos explícitos de Runge-Kutta debido a su facilidad de implementación (Butcher, 2016). Dentro de la familia de Runge-Kutta, cada método tiene un orden, que es el número de evaluaciones de derivadas utilizadas para estimar la pendiente y prever la solución de la neurona en cada paso temporal. Por lo tanto, los métodos de orden superior aumentan el número de evaluaciones de derivadas, produciendo resultados más precisos, pero también aumentando el tiempo de simulación (Butcher, 2016). Hoy en día no se conoce cuál es el método numérico y paso temporal óptimos que permitan simular los pulsos neuronales en el menor tiempo y mayor precisión posibles, como lo hemos intentado hacer en (Valadez-Godínez et al., 2017a,b, 2020).

El tiempo de simulación tiende a variar entre diferentes computadoras debido a sus propias características, como el tipo de procesador, la velocidad del procesador, el bus del procesador, la memoria, etc. Para calcular una ecuación matemática que permita realizar una predicción del tiempo de la simulación en función del paso temporal e independiente de cualquier computadora, se utilizó el método numérico de más bajo orden, el método de Euler (Euler, 1768), que pertenece a la familia de los métodos explícitos de Runge-Kutta. El método de Euler se define como:

donde Δt es el tamaño del paso temporal. Convirtiendo una ecuación diferencial en una función arbitraria, es decir, dy/dt = f(y,t), la solución para y en la siguiente actualización de tiempo y^(i+1) puede extrapolarse a partir de un valor previo y^i.

Los experimentos se realizaron en una computadora personal con el sistema operativo Linux Ubuntu 14.04 LTS de 64 bits, utilizando un procesador Intel Core i7-2600 (3.4 GHz con ocho núcleos). La memoria de la computadora era de 8 GB. Los modelos neuronales fueron implementados en el lenguaje de programación Python 2.6.5 y se calculó el tiempo utilizando la función time.clock() proporcionada por la biblioteca time (https://docs.python.org/2/library/time.html).

Se realizaron simulaciones de pulsos eléctricos en una ventana temporal de 1,000 ms en la computadora personal. Para cada modelo, se registró el tiempo de ejecución de la CPU para diferentes tamaños de paso temporal Δt (Valadez-Godínez et al., 2017a). La Fig. 2 muestra la relación entre el tiempo de ejecución de la CPU y el paso temporal Δt para cada modelo.

Figura 2. Tiempo de ejecución de la CPU en segundos para los modelos LIF, IZH y HH. Adaptado de  (Valadez-Godínez et al., 2017a).

Por inspección visual, se observó que el tiempo de simulación en función del tamaño del paso de tiempo sigue una distribución de ley de potencias (Valadez-Godínez et al., 2017a). Utilizando el método de mínimos cuadrados para el análisis de regresión lineal (Anderson, 2019), se ajustó una función de potencia a los datos observados. Como resultado, esta función altamente confiable (Valadez-Godínez et al., 2017a) nos permite predecir el tiempo de ejecución de un modelo neuronal según el tamaño del paso de tiempo en una computadora personal. Las ecuaciones de ley de potencias que permiten obtener el tiempo que tardaría una simulación en segundos en función del paso temporal para los modelos LIF, IZH y HH se muestran en las Ecs. (10), (11) y (12), respectivamente (Valadez-Godínez et al., 2017a).

3. Resultados

Una vez que se obtuvieron las funciones de ley de potencias, se procedió a utilizarlas para predecir el tiempo de ejecución de los modelos en diferentes pasos temporales. Se calculó la estimación del tiempo de simulación para un milímetro cúbico de la corteza cerebral (10^4 neuronas) y la corteza cerebral completa (10^11 neuronas). Para predecir el tiempo para 10^4 y 10^11 neuronas, simplemente se multiplicó el tiempo predicho para una sola neurona por estos números, respectivamente. Así, el tiempo estimado corresponde a un tiempo de simulación secuencial en nuestra computadora personal. Recuerde que estos valores corresponden al tiempo de simulación necesario para simular un segundo de tiempo biológico (1,000 ms). Hemos supuesto que la computadora tiene suficiente capacidad de memoria para almacenar 10^4 y 10^11 neuronas. Almacenar 10^11 neuronas en una computadora personal es prácticamente imposible hoy en día.

Independientemente del método numérico utilizado para resolver las ecuaciones diferenciales de los modelos, si el tamaño del paso de tiempo de la ecuación diferencial es más pequeño que el tiempo de simulación, la simulación tiende a ser más extensa. Por lo tanto, la tarea de definir el tamaño del paso de tiempo es crucial para reducir el tiempo de simulación. Por ejemplo, si el tamaño del paso de tiempo para el modelo HH es muy pequeño (∆t=0.000000001), entonces nuestra computadora tardará 86 días en simular un segundo biológico de pulsos de una sola neurona usando el modelo HH. Con este mismo tamaño de paso, se simularía un milímetro cúbico de la corteza cerebral en 2,356 años y la corteza cerebral completa en 23,561,643,836 años de manera secuencial.

Las ecuaciones de ley de potencias para predecir el tiempo de simulación tienden a infinito cuando el valor del tamaño del paso de tiempo tiende a cero. Cuando los valores del paso de tiempo son muy cercanos a cero, el tiempo de simulación predicho se extiende a decenas de días para simular un segundo biológico de una sola neurona, cientos de años para simular un milímetro cúbico de la corteza cerebral y decenas de milenios para simular toda la corteza cerebral en una computadora personal de manera secuencial (Valadez-Godínez et al., 2017a). Previamente realizamos un estudio para encontrar pasos temporales eficientes que conservaran un tiempo de simulación y una precisión de los pulsos adecuada y encontramos que el paso temporal para los tres modelos debe de ser de ∆t=0.01 (Valadez-Godínez et al., 2020). Con ese tamaño de paso temporal sustituido en las Ecs. (10), (11) y (12), se encontró que se necesitan esperar 219, 452 y 2096 años para simular la corteza cerebral con cada modelo LIF, IZH y HH, respectivamente. En otras palabras, la simulación de la corteza cerebral completa en una computadora personal resulta prohibitiva. Por lo tanto, simular todo el cerebro humano es más caro computacionalmente.

4. Discusión

Simular la simple generación de pulsos eléctricos a nivel cerebral resulta restrictivo. Sin embargo, sería más prohibitivo si se consideraran simulaciones funcionales al nivel genético, molecular, sináptico, o al nivel de mapas o sistemas. Sería más desafiante computacionalmente si las neuronas se simularan utilizando modelos más detallados donde se considerase la distribución geométrica de su morfológica. Realizar esto último es realmente un reto para la ciencia (Markram et al., 2015). El Human Brain Project (Proyecto Cerebro Humano en español) anunció que, por el momento, no es posible modelar el cerebro humano (https://www.humanbrainproject.eu/en/). Afirmaron que la simulación del cerebro humano está en una etapa muy temprana.

5. Conclusiones

Este estudio aborda la pregunta de si es posible simular el cerebro humano en una computadora, centrándose en la simulación de la transmisión de información a través de potenciales de acción o pulsos eléctricos en el nivel celular. Se limita al análisis del tiempo requerido para simular un milímetro cúbico de la corteza cerebral y la corteza completa, utilizando modelos de neuronas como LIF, HH e IZH.

Los modelos neuronales se representan mediante ecuaciones diferenciales ordinarias, y la simulación de estas se realiza mediante métodos numéricos, en este caso, el método de Euler. Se estableció una relación entre el tiempo de ejecución de la CPU y el tamaño del paso temporal para cada modelo, observando que sigue una distribución de ley de potencias. A partir de esta distribución, se crearon ecuaciones de ley de potencias para predecir el tiempo de simulación en función del paso temporal, y se calculó el tiempo estimado para simular un milímetro cúbico y la corteza cerebral completa. Los resultados mostraron que, incluso con pasos temporales eficientes, la simulación de la corteza cerebral completa en una computadora personal resultó prohibitivo, llevando a tardar varios siglos.

Simular la simple generación de pulsos eléctricos a nivel cerebral ya es costoso computacionalmente, y sería más desafiante si se consideraran simulaciones a niveles más detallados.

Agradecimientos

S. Valadez-Godínez expresa su agradecimiento a la Universidad Politécnica de Pénjamo y al CONAHCYT por su apoyo financiero en la realización de este trabajo de divulgación.

Referencias

Abbott, L. F., DePasquale, B., & Memmesheimer, R. M. (2016). Building functional networks of spiking model neurons. Nat Neurosci, 19(3):350–355. https://doi.org/10.1038/nn.4241

Anderson, D. R., Sweeney, D. J., Williams, T. A., Camm, J. D. & Cochran, J. J. (2019). Statistics for Business & Economics. CENGAGE, fourteenth edition.

Bohte, S. M. (2004). The evidence for neural information processing with precise spike-times: A survey. Natural Computing, 3(2):195–206. https://doi.org/10.1023/b:naco.0000027755.02868.60

Butcher, J. C. (2016). Numerical methods for ordinary differential equations. John Wiley & Sons, Ltd, third edition. https://doi.org/10.1002/9781119121534

Churchland, P. & Sejnowski, T. J. (1988). Perspectives on cognitive neuroscience. Science, 242(4879):741–745. https://doi.org/10.1126/science.3055294

Churchland, P. S. & Sejnowski, T. J. (1992). The computational brain. A Bradford Book.

Debanne, D., Campanac, E., Bialowas, A., Carlier, E., & Alcaraz, G. (2011). Axon physiology. Physiol Rev, 91(2):555–602. https://doi.org/10.1152/physrev.00048.2009

Euler, L. (1768). Institutionum calculi integralis. Volumen primum. Petropoli: Impenfis Academiae Imperialis Scientiarum.

Gerstner, W., Kistler, W. M., Naud, R., & Paninski, L. (2014). Neuronal dynamics: From single neurons to networks and models of cognition. Cambridge University Press, Cambridge.

Herz, A. V. M., Gollisch, T., Machens, C. K., & Jaeger, D. (2006). Modeling single-neuron dynamics and computations: A balance of detail and abstraction. Science, 314(5796):80–85. https://doi.org/10.1126/science.1127240

Hodgkin, A. L. & Huxley, A. F. (1952). A quantitative description of membrane current and its application to conduction and excitation in nerve. J Physiol, 117(4):500–544. https://doi.org/10.1007/bf02459568

Izhikevich, E. M. (2003). Simple model of spiking neurons. IEEE Trans Neural Netw, 14(6):1569–1572. https://doi.org/10.1109/tnn.2003.820440

Izhikevich, E. M. (2004). Which model to use for cortical spiking neurons? IEEE Trans Neural Netw, 15(5):1063–1070. https://doi.org/10.1109/tnn.2004.832719

Izhikevich, E. M. (2007). Dynamical systems in neuroscience: The geometry of excitability and bursting. The MIT Press, Cambridge, Massachusetts.

Knight, B. W. (1972). Dynamics of encoding in a population of neurons. J Gen Physiol, 59(6):734–766. https://doi.org/10.1085/jgp.59.6.734

Lapicque, L. (1907). Recherches quantitatives sur l’excitation électrique des nerfs traitée comme une polarization. J Physiol Pathol Gen, 9:620–635.

Maass, W. (1997). Networks of spiking neurons: The third generation of neural network models. Neural Netw, 10(9):1659–1671. https://doi.org/10.1016/s0893-6080(97)00011-7

Markram, H., Muller, E., Ramaswamy, S., Reimann, M. W., Abdellah, M., Sanchez, C. A., Ailamaki, A., AlonsoNanclares, L., Antille, N., Arsever, S., Kahou, G. A., Berger, T. K., Bilgili, A., Buncic, N., Chalimourda, A., Chindemi, G., Courcol, J. D., Delalondre, F., Delattre, V., Druckmann, S., Dumusc, R., Dynes, J., Eilemann, S., Gal, E., Gevaert, M. E., Ghobril, J. P., Gidon, A., Graham, J. W., Gupta, A., Haenel, V., Hay, E., Heinis, T., Hernando, J. B., Hines, M., Kanari, L., Keller, D., Kenyon, J., Khazen, G., Kim, Y., King, J. G., Kisvarday, Z., Kumbhar, P., Lasserre, S., Le Bé, J. V., Magalhães, B. R., Merchán-Pérez, A., Meystre, J., Morrice, B. R., Muller, J., Muñoz Céspedes, A., Muralidhar, S., Muthurasa, K., Nachbaur, D., Newton, T. H., Nolte, M., Ovcharenko, A., Palacios, J., Pastor, L., Perin, R., Ranjan, R., Riachi, I., Rodríguez, J. R., Riquelme, J. L., Rössert, C., Sfyrakis, K., Shi, Y., Shillcock, J. C., Silberberg, G., Silva, R., Tauheed, F., Telefont, M., Toledo-Rodriguez, M., Tränkler, T., Van Geit, W., Díaz, J. V., Walker, R., Wang, Y., Zaninetta, S. M., DeFelipe, J., Hill, S. L., Segev, I., & Schürmann, F. (2015). Reconstruction and simulation of neocortical microcircuitry. Cell, 163(2):456–492. https://doi.org/10.1016/j.cell.2015.09.029

Marr, D. C. & Poggio, T. (1977). From understanding computation to understanding neural circuitry. Neurosci Res Program Bull, 15(3):470–488.

Stein, R. B. (1965). A theoretical analysis of neuronal variability. Biophys J, 5(2):173–194. https://doi.org/10.1016/s0006-3495(65)86709-1

Valadez-Godínez, S., Sossa, H. & Santiago-Montero, R. (2017a). The step size impact on the computational cost of spiking neuron simulation, in: Proceedings of the 2017 SAI Computing Conference (SAI), London, UK, 2017a, pp. 722–728. https://doi.org/10.1109/SAI.2017.8252176

Valadez-Godínez, S., Sossa, H., & Santiago-Montero, R. (2017b). How the accuracy and computational cost of spiking neuron simulation are affected by the time span and firing rate. Computación y Sistemas, 21(4):841–861. https://doi.org/10.13053/CyS-21-4-2787

Valadez-Godínez, S., Sossa, H. & Santiago-Montero, R. (2020). On the accuracy and computational cost of spiking neuron implementation, Neural Netw, 122:196–217. https://doi.org/10.1016/j.neunet.2019.09.026 Zador, A. M. (2000). The basic unit of computation. Nat Neurosci, 3(Supp):1167–1167. https://doi.org/10.1038/81432

---