Resting-state brain activity data assimilation method based on twin brain model

Abstract

The application discloses a resting-state brain activity data assimilation method based on a twin brain model; the method steps are as follows: 1) using the pulse neuron equation set depicting the basic working activity process of the brain as a process model, coupling with the real biological structure of the brain, and designing an observation model for observing the output results of the whole brain process model at different time points; 2) initializing all neuron parameters of the twin brain model with different distributions in different regions; 3) taking the thalamic region as the only observable region and assimilation region, and representing the assimilation process through particle set evolution and correction. The application obtains specific neuron parameters by assimilating the implicit distributed hyperparameters, and only uses the thalamic region activity signal as the observation signal, so that the assimilation efficiency is greatly improved under the condition of ensuring the precision and condition of data assimilation. The application standardizes the resting-state data assimilation process of the twin brain model, and is convenient and easy to operate and realize.

Description

Technical Field

[0001] This invention belongs to the field of neuroscience technology, specifically, it relates to a method for assimilating resting-state brain activity data based on a twin brain model. Background Technology

[0002] The human brain is one of the most complex information processing systems in nature. Understanding its dynamic plasticity and complex network structure is crucial. The structure and working mechanism of the human brain are important frontier topics in current scientific research. Studies have found that the realization of the resting-state function of the human brain depends on nonlinear interactions and complex feedback activities at different spatiotemporal scales, including local neuronal population activity, the coordinated work of neural circuits, and the dynamic connectivity patterns of the whole brain. Whole-brain network models, through detailed modeling of neurons, synapses, and neural circuits, can integrate experimental data from multiple modes into a common framework, thus providing an experimental platform for testing neural activity mechanisms and studying brain function.

[0003] Neuron-based whole-brain models heavily rely on large-scale biological data for both modeling and validation. A fundamental problem in the modeling process is how to effectively estimate network parameters using recorded data from biological experiments, enabling the model to generate activity patterns corresponding to physiological data.

[0004] The goal of whole-brain model parameter estimation is to fit physiological signals at timescales such as fMRI, BOLD, or electrophysiological (EEG, MEG, sEEG). Blood oxygenation response signals collected in fMRI experiments are typically used directly to calculate functional connectivity (FC) for that task, also known as empirical functional connectivity (eFC). Electrical signals are preprocessed and subjected to bandpass filtering and Hilbert transform to extract the amplitude envelope, which is then used to calculate functional connectivity under that state. Functional connectivity is usually assumed to be fixed, i.e., unchanging over time, reflecting the specific characteristics of interregional cooperation under that task state. For example... Figure 1 Traditional resting-state fitting methods determine a final set of parameters that best fit the FC by traversing the parameter space of the free variables and the output of the simulation model. 1，2 .

[0005] Besides FC, other empirical features of biological experimental data can also be selected for model fitting. For example, Jirsa 3Parameters are estimated by fitting the spectral distribution of stereotactic (SEEG) electrodes. When the research problem requires certain specific physiological characteristics, the model can estimate parameters by fitting these characteristics, resulting in a framework model adapted to the problem. There are many methods for fitting experimental data or other data features, such as Bayesian methods and backpropagation algorithms based on loss functions, when a generative model is known. Summary of the Invention

[0006] To address the shortcomings of the existing technologies, the present invention aims to propose a resting-state brain activity data assimilation method based on a well-established twin brain model. This method, based on an initially constructed whole-brain network model of spiking neurons, seeks to explore key parameters of the neuron model to reproduce the real resting-state activity of the human brain. The present invention uses a reasonable and effective data assimilation method to fit the resting-state physiological activity of the real human brain and restores the high biological relevance of activity data in functional magnetic resonance imaging (fMRI) mode.

[0007] The technical solution of the present invention is described in detail below.

[0008] This invention provides a method for assimilating resting-state brain activity data based on a twin brain model, comprising the following steps: S1: Using the established twin brain model as the process model, the spiking neuron equation, which describes the basic working activities of the brain, is used as the basic dynamic process. Coupled with the real biological structure of the brain, the state evolution of all neurons at different times is obtained. S2: Construct an observation model for the output results of the whole-brain process model at different times based on the simulation and observation results of the process model at different times in the study area; S3: Initialize all neuron parameters of the twin brain model using different distributions in different regions, with different hyperparameters distributed across different regions; S4: Using the Thalamus region of the thalamus as both the observed and assimilation regions, the assimilation process is represented through particle ensemble evolution and correction. The process model incorporates hyperparameter walk error simulation to obtain the state evolution results of different brain regions at different times. Based on ensemble Kalman filtering (EnKF), the assimilation amount of this region is calculated at different times. Real-world thalamic region activity signals are used to assimilate the state variables of the thalamus region at each time point, outputting a hyperparameter evolution sequence showing the distribution of neuronal parameters in the assimilated thalamus region; where: In step S2, the method for constructing observation models at different times based on the simulation and observation results of process models in different regions at different times is as follows: ; In the formula, the whole brain is divided into common... Each brain region, This represents the observation of region j at time t, where the neurons corresponding to brain region j are... ,in, , Represents neurons The issuance rate at time t, This indicates the number of neurons in brain region j. OM represents the observation noise in the observation model of region j; In step S3, in the twin brain model, Maximum synaptic conductance of neuron # As free variables, the point-directed parameters of each brain region in the whole-brain model are set to follow a distribution using prior information. The specific hyperparameters for each region are: Specifically .

[0009] In this invention, in step S1, the twin brain model is a twin brain model based on leaky integration-ignition LIF neurons with instantaneous jumps and exponentially decaying synapses.

[0010] In this invention, step S1 involves constructing process models at different times to obtain the state evolution of neurons, specifically as follows: ; In the formula, t represents the discrete process time, n represents the number of nodes in the twin brain model, i.e., the number of neurons, T represents the maximum time, and PM represents the process model. , These represent the neuronal activity of the whole brain at times t+1 and t, respectively. , These represent the states of the 1st and nth neurons in the study region at time t+1, respectively. , These represent the states of the 1st and nth neurons at time t, respectively. , Let represent the conductance parameters of the 1st and nth neurons at time t, respectively. , denoted as the noise variables of the 1st and nth neurons at time t, respectively, represent the errors in simulating state evolution.

[0011] In this invention, step S4 is specifically as follows: (1) Consider the dependence of the conductance variable on hyperparameters of the distribution By assigning random walk behavior to brain region j at time t, the process model re-simulates the evolution of each state variable of the study region at different times, specifically: ; ; In the formula, Let J be the hyperparameter of brain region j at time t-1. The random white noise signal in brain region j at time t-1; The observation model then transforms the simulation results of each state variable of each grid at different times in the study area into the observation space to obtain the observation output; (2) Under resting assimilation, the thalamus region is used as the sole research area, and the set of brain region labels under this region is assumed to be: Let the firing rate state of the neurons in this region be denoted as . The state matrix corresponding to this region is ,Right now ; in Let be the neuronal firing rate state vector of the thalamus region; similarly, the observation matrix is... And the hyperparameter matrix of this region is The area has a total of For each particle, the assimilation quantity is as follows: ; in, noise vector The covariance matrix, Neuron state and observation signals The Kalman gain matrix, Hyperparameters for the distribution of conductivity variables and observation signals The Kalman gain matrix; (3) By using real experimental signals from the thalamus region as real observations Assimilation at the current moment Study of thalamic regions The Neuron firing rate state vector of each simulated particle and hyperparameter vector ,in Numbering the simulated particles using the EnKF algorithm; the specific process is as follows: ; in, Indicates the current time thalamic region The Observations from simulated particles Indicates the current time thalamic region The Noise simulated by a single particle; The system's state and unknown hyperparameters are then assimilated, with the particle mean value used as the current time point. A reasonable estimate after assimilation, i.e. ; For the current assimilated particle states and hyperparameter status This serves as the initial state for the next moment, continuing into the next cycle, thereby obtaining the system state and hyperparameter sequence at each time point.

[0012] Compared with existing technologies, the advantages of this invention are as follows: This invention obtains specific neuronal parameters by assimilating implicitly distributed hyperparameters, and utilizes only thalamic region activity signals as observation signals, achieving a significant improvement in assimilation efficiency while ensuring the accuracy of data assimilation. This invention standardizes the resting-state data assimilation process of twin brain models, making it convenient and easy to implement, resulting in twin brain models that are sufficiently similar to real resting-state biological brains. This has significant theoretical and practical implications for brain computation in the field of neuroscience. Attached Figure Description

[0013] Figure 1 Comparison of whole-brain network models in restoring resting-state activity.

[0014] Figure 2 Resting-state structural assimilation results of digital twin brain. Detailed Implementation

[0015] This invention provides a method for assimilating resting-state brain activity data based on a twin brain model, comprising the following steps: 1) Using the spiking neuron equations that characterize the basic working processes of the brain as a process model, coupled with the real biological structure of the brain, an observation model is designed to observe the output of the whole brain process model at different times. 2) Initialize all neuron parameters of the twin brain model using different distributions in different regions, with different hyperparameters distributed across different regions; 3) Using the thalamus region as the only observable and assimilation region, the process model incorporates hyperparameter walk error to simulate the state evolution of different brain regions at different times. Based on ensemble Kalman filtering (EnKF), the assimilation amount of this region is calculated at different times. The state variables of the thalamus region at each time time are assimilated using real thalamus region activity signals. The hyperparameter evolution sequence of the distribution of neuronal parameters in the assimilated thalamus region is output.

[0016] The following will, in conjunction with embodiments, provide a complete and clear explanation of the technical solutions and basic principles of the embodiments of the present invention, based on a twin brain model established using leaky integration-ignition (LIF) neurons with instantaneous jumps and exponentially decaying synapses.

[0017] Example 1

[0018] The spiking neural network model is the foundation of the method in this invention, consisting of two parts: basic computational units and network structure. The neuron model studied here is a leaky integration-ignition (LIF) neuron with instantaneous jumps and exponentially decaying synapses. Each neuron model receives postsynaptic current as input and describes the generation of action potential time points as output. The dynamic representation of the subthreshold membrane potential Vi of an excitatory (inhibitory) neuron is as follows:

[0019] Where u∈AMPA, NMDA, GABAa, GABAb represent different synaptic neurotransmitters. , These represent synaptic conduction leaks with equilibrium potentials VL and Vu, respectively. This represents the timescale of the synaptic delay effect, with AMPA channels exhibiting a faster decay constant and NMDA channels decaying much more slowly. The dimensionless parameter... Synaptic weights, gating variables This represents the fraction of open channels in neurons. The neuronal capacitance of a neuron; when the membrane potential reaches the firing threshold. At this time, the neuron emits a spike, and the membrane potential is at... Reset, and during the refractory period Keep it at that value.

[0020] The aforementioned neuron equations, under custom structural coupling, form a process model, and the output of this model... It is fed into a classic Balloon model to output the observed signal. This corresponds to the BOLD signal in biological observation data.

[0021] Many previous studies have suggested that the resting state of the human brain can be considered to be driven solely by somatosensory input as external energy. 4 The method of this invention considers a major relay point for sensory information—the thalamus. 5 This serves as a region for signal assimilation. The actual observed signals in whole-brain network models consist only of the resting-state BOLD signal in the thalamus. This invention assimilates the input current in the thalamus. Other state signals enable the model to exhibit biologically similar resting-state whole-brain activity patterns.

[0022] In an embodiment of the present invention, the resting-state BOLD signal recorded experimentally in the thalamus is used as assimilated data, and it is assumed that... Obedience to hyperparameters Controlled distribution The external current parameters of the region are estimated using the method described above. This allows for the fitting of the actual observed BOLD signal. After assimilation, the process model and observation model are re-simulated using only the current parameters assimilated in the model. The correlation between the observed BOLD obtained from this simulation and the actual BOLD is used to evaluate the results of the twin brain resting-state assimilation.

[0023] Figure 2 The results of deploying the present invention in the model are presented. The BOLD of each voxel can be matched well with the experimentally measured BOLD signal and converge to the vicinity of the true value in a short time. From the perspective of similarity, the average correlation of the whole brain can reach 0.79.

[0024] The present invention, based on the hierarchical ensemble Kalman filter algorithm and combined with multilayer Bayesian estimation theory, assimilates the external stimulus currents of the thalamus (resting-state driving source) in a neural network. Here, we assume that the external current parameters of neurons in the same region (e.g., the thalamus) of the network follow the same probability distribution; therefore, only the hyperparameters of the distribution need to be estimated, and the parameters are updated accordingly when the hyperparameters are updated. In this embodiment, data assimilation is used to fit the external currents of the thalamus, thereby driving the operation of the whole brain and reconstructing the resting-state brain. For the blood oxygenation level-dependent BOLD signal obtained from resting-state human brain fMRI scans, the external input currents of neurons in the thalamic region of the whole-brain network model are assimilated, and then the brain operation is driven through the connection between the thalamus and other brain regions, thus realizing the resting-state reproduction and simulation of the whole-brain model.

[0025] References

[0026] 1.HIPP JF, HAWELLEK DJ, CORBETTA M, et al. Large-scale corticalcorrelation structure of spontaneous oscillatory activity[J]. Natureneuroscience, 2012, 15(6): 884-890. 2.DECO G, CABRAL J, WOOLRICH M W, et al. Single or multiplefrequencygenerators in on-going brain activity: A mechanistic whole-brainmodel of empiricalmeg data[J]. Neuroimage, 2017, 152: 538-550. 3.JIRSA V K, PROIX T, PERDIKIS D, et al. The virtual epilepticpatient: individualizedwhole-brain models of epilepsy spread[J]. Neuroimage,2017, 145:377-388. 4.BERNTSON G G, GIANAROS P J, TSAKIRIS M. Interoception and theautonomicnervous system: Bottom-up meets top-down[J]. The interoceptivemind:From homeostasis to awareness, 2019, 9: 3-23. 5.PHILLIPS J W, SCHULMANN A, HARA E, et al. A repeated moleculararchitectureacross thalamic pathways[J]. Nature neuroscience, 2019, 22(11):1925-1935。

Claims

1. A method for assimilating resting-state brain activity data based on a twin brain model, characterized in that, Includes the following steps: S1: Using the established twin brain model as the process model, the spiking neuron equation, which describes the basic working activities of the brain, is used as the basic dynamic process. Coupled with the real biological structure of the brain, the state evolution of all neurons at different times is obtained. S2: Construct an observation model for the output results of the whole-brain process model at different times based on the simulation and observation results of the process model at different times in the study area; S3: Initialize all neuron parameters of the twin brain model using different distributions in different regions, with different hyperparameters distributed across different regions; S4: Using the Thalamus region of the thalamus as both the observed and assimilation regions, the assimilation process is represented through particle ensemble evolution and correction. The process model incorporates hyperparameter walk error simulation to obtain the state evolution results of different brain regions at different times. Based on ensemble Kalman filtering (EnKF), the assimilation amount of this region is calculated at different times. Real-world thalamic region activity signals are used to assimilate the state variables of the thalamus region at each time point, outputting a hyperparameter evolution sequence showing the distribution of neuronal parameters in the assimilated thalamus region; where: In step S2, the method for constructing observation models at different times based on the simulation and observation results of process models in different regions at different times is as follows: ； In the formula, the whole brain is divided into common... Each brain region, This represents the observation of region j at time t, where the neurons corresponding to brain region j are... ,in, , Represents neurons The issuance rate at time t, This indicates the number of neurons in brain region j. OM represents the observation noise in the observation model of region j; In step S3, in the twin brain model, Maximum synaptic conductance of neuron # As free variables, the point-directed parameters of each brain region in the whole-brain model are set to follow a distribution using prior information. The specific hyperparameters for each region are: Specifically 。 2. The method for assimilation of resting-state brain activity data based on a twin brain model as described in claim 1, characterized in that, In step S1, the twin brain model is a twin brain model based on leaky integration-ignition LIF neurons with instantaneous jumps and exponentially decaying synapses.

3. The method for assimilation of resting-state brain activity data based on a twin brain model as described in claim 1, characterized in that, In step S1, process models at different time points are constructed to obtain the state evolution of neurons, specifically as follows: ； In the formula, t represents the discrete process time, n represents the number of nodes in the twin brain model, i.e., the number of neurons, T represents the maximum time, and PM represents the process model. , These represent the neuronal activity of the whole brain at times t+1 and t, respectively. , These represent the states of the 1st and nth neurons in the study region at time t+1, respectively. , These represent the states of the 1st and nth neurons at time t, respectively. , Let represent the conductance parameters of the 1st and nth neurons at time t, respectively. , denoted as the noise variables of the 1st and nth neurons at time t, respectively, represent the errors in simulating state evolution.

4. The method for assimilation of resting-state brain activity data based on a twin brain model as described in claim 1, characterized in that, The specific steps of step S4 are as follows: (1) Consider the dependence of the conductance variable on hyperparameters of the distribution By assigning random walk behavior to brain region j at time t, the process model re-simulates the evolution of each state variable of the study region at different times, specifically: ； ； In the formula, Let J be the hyperparameter of brain region j at time t-1. The random white noise signal in brain region j at time t-1; The observation model then transforms the simulation results of each state variable of each grid at different times in the study area into the observation space to obtain the observation output; (2) Under resting assimilation, the thalamic region is used as the sole research area, and the set of brain region labels under this region is assumed to be: Let the firing rate state of the neurons in this region be denoted as . The state matrix corresponding to this region is ,Right now ； in Let be the neuronal firing rate state vector of the thalamus region; similarly, the observation matrix is... And the hyperparameter matrix of this region is The area has a total of For each particle, the assimilation quantity is as follows: ； in, noise vector The covariance matrix, Neuron state and observation signals The Kalman gain matrix, Hyperparameters for the distribution of conductivity variables and observed signals The Kalman gain matrix; (3) By using real experimental signals from the thalamus region as real observations Assimilation at the current moment Study of thalamic regions The Neuron firing rate state vector of a simulated particle and hyperparameter vector ,in Numbering the simulated particles using the EnKF algorithm; the specific process is as follows: ； in, Indicates the current time thalamic region The Observations from simulated particles Indicates the current time thalamic region The Noise simulated by a single particle; The system's state and unknown hyperparameters are then assimilated, with the particle mean value used as the current time point. A reasonable estimate after assimilation, i.e. ； For the current assimilated particle states and hyperparameter status This serves as the initial state for the next moment, continuing into the next cycle, thereby obtaining the system state and hyperparameter sequence at each time point.

Images