Personalized real-time regulation and control method and device for nerve and blood vessel coupling model based on transcranial direct current stimulation response
By constructing a personalized transcranial direct current stimulation neurovascular coupling model, and combining multimodal data and neurotransmitter dynamics, the tDCS parameters were dynamically adjusted, which solved the problems of individual differences and inaccurate real-time feedback in post-stroke aphasia patients, and achieved precise control of neurovascular coupling and optimization of treatment effects.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- BEIJING UNIV OF TECH
- Filing Date
- 2026-02-10
- Publication Date
- 2026-05-01
Smart Images

Figure CN121944375A_ABST
Abstract
Description
Technical Field
[0001] This application relates to the field of neuromodulation and medical technology for post-stroke aphasia, specifically to a personalized real-time modulation method and device based on a neurovascular coupling model of transcranial direct current stimulation response. Background Technology
[0002] Post-stroke aphasia (PSA) is an acquired language disorder caused by damage to the dominant language area of the brain hemisphere due to stroke. It manifests as varying degrees of impairment in six areas of language function: spontaneous speech, auditory comprehension, repetition, naming, and writing. PSA is one of the most common functional impairments after stroke, with an incidence rate of up to 32% after the first stroke. Although some spontaneous recovery occurs in patients with acute PSA, most will still experience some residual speech dysfunction.
[0003] Transcranial direct current stimulation (tDCS), a non-invasive brain stimulation technique, influences brain activity by applying low-intensity direct current directly to the scalp, altering neuronal excitability. Clinically, tDCS is used to modulate various neuropsychiatric disorders, such as depression, anxiety, and post-stroke aphasia. Existing research indicates that the neural effects of tDCS extend beyond changes in neural activity; it may also indirectly affect cerebral blood flow by altering neurovascular coupling mechanisms. Therefore, optimizing the impact of tDCS stimulation on neurovascular coupling can help improve model performance and provide new directions for early research on post-stroke aphasia. However, the effectiveness of tDCS varies significantly among individuals, and traditional tDCS stimulation protocols are difficult to optimize for individual neurovascular coupling characteristics. Therefore, achieving individualized tDCS stimulation modulation remains a significant technical challenge in the field of neuromodulation.
[0004] Neurovascular coupling (NVC) refers to the relationship between neuronal activity and local cerebral blood flow. Under normal circumstances, when neuronal activity increases in a brain region, local blood flow in that region also increases to meet the increased metabolic demands; this phenomenon is crucial for the normal functioning of the brain. Dysregulation of neurovascular coupling is often closely associated with various neurological diseases, such as stroke and Alzheimer's disease.
[0005] Research on tDCS began in the 1990s and has gradually developed into a commonly used neuromodulation tool, achieving certain application results in areas such as neuromodulation, cognitive enhancement, and pain management. Its main principle involves placing electrodes on the scalp surface and applying a low-intensity direct current (typically 1-2 mA) to the cerebral cortex, altering the membrane potential of neurons and thus regulating neural activity. However, the effectiveness of tDCS varies considerably among individuals. Factors such as brain anatomy, neural network activity, and cerebral blood flow all influence the results. These differences mean that a universal tDCS stimulation protocol cannot meet the needs of different patients, resulting in a high degree of uncertainty in the modulatory effect.
[0006] To address this issue, recent studies have attempted to enhance neuromodulation effects by personalizing and real-time adjusting tDCS stimulation protocols. For example, some studies have employed neuroimaging techniques (such as functional magnetic resonance imaging (fMRI) and functional near-infrared spectroscopy (fNIRS)) to monitor brain activity and cerebral blood flow responses, thereby adjusting the tDCS current intensity and stimulation protocol based on individual neurovascular coupling characteristics. Neuroimaging techniques have demonstrated that brain stimulation can controllably alter brain networks, creating a closed-loop neuromodulation mechanism. Through this closed-loop real-time feedback mechanism, researchers can dynamically adjust tDCS parameters based on real-time neural activity and blood flow data to achieve optimal modulation results.
[0007] Furthermore, some studies have established computational models to simulate the neurovascular coupling process and the neural effects of tDCS, thereby predicting individual responses to tDCS and further optimizing stimulation protocols. Therefore, this invention promotes NVC improvement through personalized tDCS regulation. tDCS stimulation can increase the excitability of neurons in the damaged area, activate glutamate release, thereby triggering vasodilation responses through astrocytes and improving local CBF; tDCS promotes endothelial nitric oxide (NO) synthesis, enhances vasodilation capacity, and improves the reactivity of microvessels damaged after stroke; tDCS regulates Ca²⁺ signaling activity within astrocytes, optimizing their control over local blood flow and achieving dynamic synchronization between neural activity and blood flow; long-term, periodic stimulation promotes the plasticity of the brain microvascular network, increases capillary density, and expands the physiological reserve of coupling capacity. Summary of the Invention
[0008] Therefore, this application provides a personalized real-time control method and device based on a neurovascular coupling model of transcranial direct current stimulation response, in order to solve the problems of individual differences, inaccurate real-time feedback and complex equipment in the prior art.
[0009] To achieve the above objectives, this application provides the following technical solution:
[0010] In a first aspect, a personalized real-time modulation method based on a neurovascular coupling model of transcranial direct current stimulation response is characterized by comprising:
[0011] Step S1: Acquire multimodal data composed of magnetic resonance imaging (T1w) data, diffusion tensor imaging (DTI) data, and resting-state functional magnetic resonance imaging (rs-fMRI) data from healthy individuals and patients with post-stroke aphasia, respectively, and preprocess the multimodal data, then extract the structural and functional connectivity matrices respectively.
[0012] Step S2: Based on the preprocessed multimodal data, perform finite element simulation modeling to simulate the electric field distribution of tDCS in the brains of stroke aphasia patients.
[0013] Step S3: Construct a large-scale dynamic mean-field model (DMF) that integrates neurotransmitter dynamics and synaptic plasticity mechanisms, and obtain neural activity;
[0014] Step S4: Establish a large-scale neurovascular coupling model of the whole brain that links neural activity and hemodynamic response in brain regions;
[0015] Step S5: Combine finite element simulation and neurovascular coupling model to construct an individualized tDCS real-time control model.
[0016] Preferably, in step S1, the multimodal data, composed of magnetic resonance imaging (T1w) data, diffusion tensor imaging (DTI) data, and resting-state functional magnetic resonance imaging (rs-fMRI) data of patients with post-stroke aphasia, is obtained from the OpenNeuro database. The structural connectivity matrix and functional connectivity matrix are extracted using the AAL-90 brain atlas in Freesurfer.
[0017] Preferably, in step S1, when extracting the structural connectivity matrix and functional connectivity matrix using the AAL-90 brain atlas in Freesurfer, the cerebral cortex is divided into 90 regions.
[0018] Preferably, in step S2, the multimodal data of healthy individuals and post-stroke aphasia patients are preprocessed, including visualization, artifact removal, temporal registration, head movement correction, structural-functional registration, and spatial normalization. Finite element simulation modeling is then used to simulate the electric field strength and current distribution in the brains of post-stroke aphasia patients using tDCS.
[0019] Preferably, in step S3, a large-scale dynamic mean field (DMF) model integrating neurotransmitter dynamics and synaptic plasticity mechanisms is constructed to describe the whole-brain neural activity of patients with post-stroke aphasia. Its nodal dynamic equations are expressed as follows:
[0020]
[0021]
[0022]
[0023] Where i = 1, 2, ..., N, and N = 90, representing the i-th cortical or subcortical brain region delineated based on the AAL-90 brain atlas. In the dynamic mean-field model, NMDA-mediated synaptic gating variables... The dynamics are determined by the time constant. Control, synaptic activation rate parameter set to . This represents the average synaptic gating variable mediated by NMDA in the i-th brain region. Essentially, it is a synaptic-level state variable, with a value range of 0 ≤ 0. ≤1. Formula (1) includes neurodynamics and synaptic parameters. This represents the synaptic time constant, with a value of 100 ms; This represents the synaptic activation rate constant, with a value of 0.641 / 1000ms. -1 σ represents the noise amplitude, with a value of 0.001 nA. This represents Gaussian white noise with zero mean and unit variance. and The dynamics of NMDA are determined by σ=0.001, which is the key setting for resting-state oscillations.
[0024] It is an input-output function, in the form of a nonlinear sigmoid, representing the average firing rate of a local neuron population. Formula (2) includes the parameters of the input-output function, where a, b, and d represent the parameters of the input-output function, with values of 270 nC. −1 , 108Hz, 0.154s, where nC −1 =(10 -9 C) −1 =10 9 C −1, , 1 Hz = 1 s⁻¹. The input current is mapped to the average discharge rate (Hz), which is approximately linear in the low input region and saturates in the high input region, equivalent to the population behavior of LIF neurons.
[0025] This represents the total synaptic input current received by the i-th brain region, in A = Cs -1 10 9 C −1×(C s⁻¹) = 10 9 s⁻¹. Formula (3) includes synaptic coupling and network parameters. This represents the intensity of local excitatory self-feedback, with a value of 0.9; The value represents the unit synaptic coupling strength, which is 0.2609 nA; G represents the global coupling strength, with a scanning parameter range of 0-3 and an optimal coupling strength of 0.69. This represents the normalized structural connectivity matrix extracted based on DTI. This represents the external input current, with a value of 0.3nA, used to maintain a low discharge steady state. =0.9 is to ensure the stability of the local loop; It is obtained by the spiking → mean-field mapping; G is the most important adjustable parameter, used to fit the full-field (FC). This is to regulate the system to a low discharge (2–3 Hz) resting state.
[0026] To further understand the modulating effects of transcranial direct current stimulation (tDCS) and changes in the neurochemical environment under stroke pathological conditions on local neurodynamics, a neurotransmitter regulation mechanism was introduced based on the nodal dynamics model to simulate the excitatory neurotransmitter glutamate and the inhibitory neurotransmitter... The effect of changes in γ-aminobutyric acid (GABA) concentration on synaptic input and excitation-inhibition balance.
[0027] The modulation mechanism of neurotransmitters on synaptic gain is crucial. Physiologically, glutamate primarily mediates excitatory synaptic transmission, while GABA primarily mediates inhibitory synaptic transmission. Therefore, in the model, the effect of changes in neurotransmitter concentration is reflected by adjusting the synaptic gain parameters of different types of synapses. Increased glutamate concentration mainly enhances excitatory synaptic gain (J). EE J EI Increased GABA concentration primarily enhances inhibitory synaptic gain (J). IE J II ).
[0028] The excitatory-inhibitory neuronal population dynamics model refers to further dividing a local neuronal population within each brain region into excitatory (E) and inhibitory (I) subpopulations. The kinetic equations are based on chemical synaptic input and membrane potential changes; these equations are commonly used to describe the dynamic behavior of neuronal populations, for example:
[0029]
[0030]
[0031] in: These represent the average activity levels of excitatory and inhibitory neuronal populations, respectively, and are usually expressed as average firing rates; The time constant of the excitability (NMDA) gate variable is 100 ms, indicating that the population responds slowly to changes in input. This represents the time constant of the inhibition (GABA-A) gated variable, with a value of 10 ms, indicating that the population responds quickly to changes in input. These are functions The derivative with respect to time t represents the rate of change of the average activity level (usually expressed as the firing rate of neurons) of the excitatory and inhibitory neuronal populations over time. , : represent the nonlinear activation functions of excitatory and inhibitory neuron populations, respectively, used to describe how input current is converted into population firing rate, typically exhibiting S-shaped characteristics. Here... It is the formula (2) Replace with The result, substituted into the following:
[0032]
[0033] Where a = 270nC −1 b=108Hz, d=0.154s.
[0034] The total input currents received by the excitatory and inhibitory neuronal populations are expressed as follows:
[0035]
[0036]
[0037]
[0038] Among them, J EE J EI J IE J II : Represents the effective coupling strength of different types of synapses, reflecting the neurotransmitter regulation effect, with values of 1.5, 1.0, 1.2, and 0.5 respectively; , : Represents the external input current, used to simulate the external input current of tDCS modulation, with values of 0.3nA and 0.3nA respectively.
[0039] 6. The personalized real-time control method for a neurovascular coupling model based on transcranial direct current stimulation response according to claim 5, wherein the kinetic equation is expressed as in step S4, is characterized in that the Balloon–Windkessel hemodynamic model proposed by Friston et al. is used to analyze local neural synaptic activity. This is converted into an observable BOLD signal. The model describes the biophysical coupling relationship between changes in vasodilation, cerebral blood flow, blood volume, and deoxyhemoglobin content induced by neural activity, as follows:
[0040]
[0041]
[0042]
[0043]
[0044] The meaning of each variable in formulas (10)-(13) is as follows: Indicates a vasodilatory signal. Representation function The derivative with respect to time t; Indicates the first The neural synaptic activity input of each brain region (from the DMF model) is generally taken as follows: = It will trigger vasodilatory signals. The increase in this signal is influenced by self-regulatory feedback, which indirectly determines subsequent BOLD-related variables such as vasodilation signals, blood flow, and blood volume. Indicates cerebral blood flow. Representation function The derivative with respect to time t, at rest =1; Indicates cerebral blood volume. Representation function The derivative with respect to time t, at rest =1; This indicates the deoxyhemoglobin content. Representation function The derivative with respect to time t, at rest =1.
[0045] The values for each parameter in formulas (10)-(13) are assigned as follows: Indicates the first The signal attenuation rate for each region is set to 0.65s. -1 ; This represents the strength of the self-adjusting feedback, with a value of 0.41s. -1 ; This represents the hemodynamic mean transit time, with a value of 0.98 s; This represents the Grubb index (blood volume-to-blood flow index), with a value of 0.32. This represents the resting oxygen extraction fraction, with a value of 0.34.
[0046] Based on the above hemodynamic variables, the first The BOLD signal in each brain region is represented as follows:
[0047]
[0048] in, Indicates resting blood volume fraction. =0.02; It is the resting oxygen extraction fraction, obtained from formula (11). =0.34; Indicates the deoxyhemoglobin content; Indicates cerebral blood volume; , and This represents a set of weighting coefficients related to blood oxygen dynamics. =2.38, =2, =0.88. This BOLD model is a local model, which only maps the synaptic activity of the corresponding brain region to the BOLD signal of that brain region, and does not directly introduce cross-brain region blood flow coupling.
[0049] Preferably, in step S5, to maximize the intervention effect of transcranial direct current stimulation (tDCS) on patients with post-stroke aphasia, the model parameters and stimulation parameters are jointly optimized based on the aforementioned neurodynamic model and neurovascular coupling model, and a personalized stimulation plan is formulated accordingly. The goal of optimization is to enhance the coupling strength between neural activity and blood flow response in the target brain region, thereby improving local cerebral blood flow and metabolic response, and thus promoting functional recovery.
[0050]
[0051]
[0052]
[0053] in, Indicates the first The total synaptic input current of the excitatory neuronal population received by a brain region at time t reflects the level of neural activity or the strength of synaptic drive in that brain region. Indicates the first The blood oxygen metabolism level of a brain region at time t is used to characterize the intensity of the hemodynamic-metabolic response of that brain region and can be correlated with blood oxygen consumption rate (CMRO2) or hemodynamic state. = ; This indicates that within the selected time window or stimulation cycle, the first... The maximum value of the excitatory synaptic input current in each brain region is used to normalize neural activity. This indicates that within the corresponding time window, the first... The maximum value of the blood oxygen metabolism response in each brain region was used to normalize the blood flow / metabolic response. (Known) and Parameters used: w=0.9 =0.2609nA, =0.3nA, a=270nC−1, b=108Hz, d=0.154s, =100ms =0.641 / 1000ms⁻¹, σ=0.001nA, therefore The above parameter values can be used to determine this. Formula (15) quantitatively characterizes the coupling strength between neural activity and blood flow response in the target brain region by calculating the normalized product of excitatory neural activity and blood flow metabolic response, and is used to characterize the neurovascular coupling (NVC) state. Coupling degree The larger the value, the stronger the neural activity and hemodynamic response. These goals are achieved by minimizing a specific objective function, which is typically the difference in BOLD signal, the degree of coupling between neural activity and hemodynamics, and the degree of recovery of motor function.
[0054] Indicates the first BOLD fMRI signals of brain regions before transcranial direct current stimulation (tDCS) intervention; Indicates the first BOLD fMRI signal of a brain region after transcranial direct current stimulation (tDCS) intervention; |·| represents the absolute value of the change in BOLD signal before and after stimulation, used to measure the degree of functional change caused by stimulation. The parameters used for outputting BOLD are known: =0.02, =0.34, =2.38, =2, 0.88, =0.98, =0.32, =0.65 =0.41, therefore as well as All of these can be determined through the above parameter values. Formula (16) is used to quantitatively describe the degree of change in the oxygen level dependent signal (BOLD) of the target brain region before and after transcranial direct current stimulation, so as to reflect the regulatory effect of stimulation on brain function.
[0055] This represents the comprehensive optimization objective function, used to guide the optimization of personalized transcranial direct current stimulation parameters; These are weighting coefficients, with values of 0.5 and 0.5 respectively, and satisfying the following conditions: + =1, used to balance the relative importance of neurovascular coupling indicators and BOLD functional response indicators in the optimization process. It can be adjusted according to the actual situation to balance the relative importance between different targets. Formula (17) optimizes the neurovascular coupling strength and BOLD functional response changes together, which can improve the functional observability and applicability of the stimulation scheme while ensuring the rationality of the neural regulation mechanism.
[0056] During the stimulation parameter optimization phase, the system's core objective is to enhance the neurovascular coupling strength between neural activity and hemodynamic response in the target brain region, while simultaneously inhibiting abnormal activation in non-target brain regions. After each cycle of transcranial direct current stimulation, the system evaluates the objective function based on real-time or offline acquired neuroimaging and hemodynamic data, and dynamically updates stimulation parameters such as current intensity and electrode position. These stimulation parameter updates follow the following relationship:
[0057]
[0058]
[0059] in, The stimulation parameters for the k-th stimulation cycle are given, and the parameter vector includes at least the current intensity, electrode position, and stimulation duration. To control the adjustment range of the updated stimulus parameters each time; is the parameter adjustment coefficient, with a value range of [0.01, 0.3], and here it is set to 0.1. It is used to control the magnitude of stimulus parameter update and avoid over-adjustment leading to stimulus instability; ∇J is the change (or approximate gradient) of the objective function in adjacent stimulus cycles. The value range is 1-2mA. =1mA; The value range is 10-20 min. =20min; Corresponding to the positions of the anode and cathode of the electrode plate, ={F7, F8}. Through the above method, the stimulation parameters are adaptively adjusted, thereby gradually optimizing the neurovascular coupling state and promoting functional recovery.
[0060] Secondly, a personalized real-time control device based on a neurovascular coupling model of transcranial direct current stimulation response includes:
[0061] Data acquisition module: used to extract multimodal data composed of magnetic resonance imaging (T1w) data, diffusion tensor imaging (DTI) data, and resting-state functional magnetic resonance imaging (rs-fMRI) data from healthy individuals and patients with post-stroke aphasia;
[0062] Preprocessing and extracting the connectivity matrix module: This module is used to preprocess multimodal data and extract the structural and functional connectivity matrices respectively.
[0063] Finite element simulation module: used to simulate the electric field distribution of tDCS in the brains of stroke patients with aphasia;
[0064] Large-scale dynamic mean-field model building module: used to construct large-scale dynamic mean-field models that integrate neurotransmitter dynamics and synaptic plasticity mechanisms;
[0065] Neural activity output module: used for neural activity obtained from large-scale dynamic mean-field models;
[0066] Whole-brain large-scale neurovascular coupling model establishment module: used to link neural activity and hemodynamic responses in brain regions;
[0067] tDCS Real-Time Control Module: Used to combine finite element simulation and neurovascular coupling model to construct individualized tDCS real-time control model.
[0068] Thirdly, a personalized real-time control device includes a computer program / instructions, characterized in that the computer program / instructions, when executed by a processor, implement the steps of the method according to any one of claims 1 to 6.
[0069] Compared with the prior art, this application has at least the following beneficial effects:
[0070] This application provides a personalized real-time control method and device based on a neurovascular coupling model of transcranial direct current stimulation response. It acquires multimodal data composed of magnetic resonance imaging (T1w), diffusion tensor imaging (DTI), and resting-state functional magnetic resonance imaging (rs-fMRI) data from healthy individuals and post-stroke aphasia patients. The multimodal data is preprocessed, and structural and functional connectivity matrices are extracted separately. Based on the preprocessed multimodal data, finite element simulation modeling is performed to simulate the electric field distribution of tDCS in the brains of post-stroke aphasia patients. A large-scale dynamic mean-field model (DMF) integrating neurotransmitter dynamics and synaptic plasticity mechanisms is constructed to obtain neural activity. A whole-brain large-scale neurovascular coupling model linking neural activity and hemodynamic responses in brain regions is established. Finally, a personalized real-time control model of tDCS is constructed by combining finite element simulation and the neurovascular coupling model.
[0071] This technical solution precisely assesses an individual's neurovascular coupling characteristics and, combined with monitoring data, dynamically adjusts tDCS stimulation parameters using a closed-loop real-time feedback mechanism. This application overcomes the shortcomings of existing technologies, such as individual variability, inaccurate real-time feedback, and complex equipment. It allows for personalized real-time control of neural activity and hemodynamic parameters in post-stroke aphasia patients using a neurovascular coupling model, thereby obtaining real-time tDCS stimulation effects. Through multiple rounds of adjustments and evaluations, the system continuously optimizes the stimulation protocol and improves the stimulation effect until it reaches its optimal state. This helps physicians understand the patient's optimal stimulation location, intensity, and other parameters. After the control is completed, the subsequent control protocol is further adjusted based on the patient's recovery progress and symptom improvement to ensure long-term efficacy. Attached Figure Description
[0072] To more intuitively illustrate the prior art and this application, exemplary drawings are provided below. It should be understood that the specific shapes and structures shown in the drawings should not generally be regarded as limiting conditions for implementing this application; for example, based on the technical concept disclosed in this application and the exemplary drawings, those skilled in the art are able to easily make conventional adjustments or further optimizations to the addition / reduction / classification, specific shapes, positional relationships, connection methods, size ratios, etc. of certain units (components).
[0073] Figure 1 A system structure diagram of personalized real-time control based on a neurovascular coupling model of transcranial direct current stimulation response provided in Embodiment 1 of this application;
[0074] Figure 2 A flowchart illustrating the personalized real-time control of a neurovascular coupling model based on transcranial direct current stimulation response provided in Embodiment 1 of this application;
[0075] Figure 3 This refers to all 45 ROIs in each cerebral hemisphere of the AAL map provided in Embodiment 1 of this application;
[0076] Figure 4 This is a flowchart of the tDCS stimulation process provided in Embodiment 1 of this application;
[0077] Figure 5 The result graph of principal component analysis (PCA) provided in Embodiment 1 of this application;
[0078] Figure 6 This is a diagram of the BOLD signal results provided in Embodiment 1 of this application;
[0079] Figure 7 A graph showing the results of excitatory and inhibitory synaptic discharge rates in various brain regions provided in Embodiment 1 of this application;
[0080] Figure 8 The functional connectivity values between the damaged brain region and other brain regions in the whole brain provided in Embodiment 1 of this application;
[0081] Figure 9 The distance between the functional connection FCs and the correlation between the empirical FC and the simulated FC are provided for Embodiment 1 of this application. Detailed Implementation
[0082] The present application will be further described in detail below with reference to the accompanying drawings and specific embodiments.
[0083] In the description of this application: unless otherwise stated, "a plurality of" means two or more. The terms "first," "second," "third," etc., in this application are intended to distinguish the objects referred to and do not have any special meaning in terms of technical connotation (e.g., they should not be construed as an emphasis on importance or order). Expressions such as "including," "comprising," and "having" also mean "not limited to" (certain units, components, materials, steps, etc.).
[0084] The terms used in this application, such as "upper," "lower," "left," "right," and "middle," are generally used to indicate the general relative positional relationship for the purpose of intuitive understanding by referring to the accompanying drawings, and are not absolute limitations on the positional relationship in the actual product.
[0085] Example 1:
[0086] Please see Figure 1 and Figure 2 This embodiment provides a personalized real-time modulation method for a neurovascular coupling model based on transcranial direct current stimulation response, including:
[0087] S1: Acquire multimodal data composed of magnetic resonance imaging (T1w) data, diffusion tensor imaging (DTI) data, and resting-state functional magnetic resonance imaging (rs-fMRI) data from healthy individuals and patients with post-stroke aphasia, respectively. Preprocess the multimodal data and then extract the structural and functional connectivity matrices.
[0088] Specifically, this embodiment obtains magnetic resonance imaging (T1w), diffusion tensor imaging (DTI), and resting-state functional magnetic resonance imaging (rs-fMRI) data from the OpenNeuro database for healthy individuals and patients with post-stroke aphasia.
[0089] The OpenNeuro database is dedicated to analyzing and sharing neuroimaging data. Based on the Brain Imaging Data Structure (BIDS) specification, it provides researchers with a unified environment for uploading, processing, and sharing neuroimaging datasets. This post-stroke aphasia dataset includes anonymized images, behavioral measurements, and demographic details from a cohort of individuals with chronic stroke and speech and language impairments. The dataset also includes multiple modalities, including anatomical (T1w, T2w, FLAIR), functional (fMRI), and resting-state data. Furthermore, the dataset is fully anonymized, allowing for public sharing, which is crucial for education and research.
[0090] Multimodal data from healthy individuals and post-stroke aphasia patients underwent preprocessing including visualization, artifact removal, temporal registration, head movement correction, structural-functional registration, and spatial normalization. For diffusion tensor imaging (DTI) data from healthy individuals and post-stroke aphasia patients, DSI_Studio software was used for quality control and resampling. The AAL-90 brain atlas in Freesurfer was used for data extraction, dividing the cerebral cortex into 90 regions to obtain a 90×90 structural connectivity matrix (SC).
[0091] In the preprocessing of magnetic resonance imaging (T1w) data and resting-state functional magnetic resonance imaging (rs-fMRI) data of healthy individuals and patients with post-stroke aphasia, time registration, head movement correction, structural-functional registration and spatial normalization were performed using SPM12 in MATLAB. Similarly, the AAL-90 brain atlas in Freesurfer was used for extraction, dividing the cerebral cortex into 90 regions to obtain a 90×90 functional connectivity matrix (FC).
[0092] Pearson correlation analysis was used to reconstruct the functional connectome FC(i,j) to accurately characterize the strength and patterns of connections between neurons, providing strong support for in-depth exploration of the intrinsic structure of brain functional networks. The calculation formula is as follows:
[0093]
[0094] in, and These represent the time-series signal values of brain regions i and j at time point t, respectively; and Then, these represent the average values of the time series of brain regions i and j, respectively.
[0095] S2: Based on the preprocessed multimodal data, finite element simulation modeling is performed to simulate the electric field distribution of tDCS in the brains of stroke aphasia patients;
[0096] Specifically, this step is based on finite element simulation modeling to simulate the electric field strength and current distribution of tDCS in each node of the brain of a stroke-induced aphasia patient.
[0097] An individualized head model was established by preprocessing MRI and DTI image data, including image segmentation, denoising, and standardization, to generate a clear head geometry model for subsequent finite element analysis.
[0098] Finite element mesh construction: Based on the individualized head geometry model, a finite element mesh is generated. The head is divided into multiple small elements (such as tetrahedral or hexahedral elements), and the size and shape of each element are adjusted according to the complexity of the anatomical structure. The mesh density is ensured to be sufficiently high in key regions (such as electrode placement areas and cerebral cortex regions) to improve computational accuracy.
[0099] Setting Material Parameters: The head structure is automatically segmented using the mri2mesh tool in the SimNIBS open-source software. Based on the conductivity and dielectric constant of different tissues (such as scalp, skull, cerebrospinal fluid, gray matter, and white matter), corresponding material parameters are assigned to each element. These parameters are typically obtained from literature or through experimental measurements. Ensuring the accuracy of the material parameters is crucial, as they directly affect the calculation of current and electric field distributions.
[0100] Apply boundary conditions and excitation sources: In the finite element model, set the position and shape of the electrodes as the excitation source of the current. Specify the current intensity and polarity of the electrodes according to the experimental design for tDCS (transcranial direct current stimulation). Set boundary conditions, such as zero-potential boundary conditions on the scalp surface, to simulate the actual electrical stimulation environment.
[0101] The potential distribution under transcranial direct current stimulation was solved using the finite element method. The potential distribution satisfies the Poisson equation of the quasi-static conductivity model.
[0102]
[0103] in, To improve electrical conductivity, The potential distribution is shown. Tissue conductivity was determined using the literature and SimNIBS default parameter settings. The conductivity of scalp, skull, cerebrospinal fluid, gray matter, and white matter was set to 0.465, 0.010, 1.654, 0.275, and 0.126 S / m, respectively.
[0104] Model validation and optimization: Validate the calculation results, checking whether the current density and electric field distribution conform to physical laws and experimental expectations. Adjust the finite element mesh density, material parameters, or boundary conditions as needed to optimize the model's accuracy and reliability.
[0105] By following the steps above, a finite element model can be constructed to calculate the current density and electric field intensity distribution of each node, providing a foundation for building an individualized tDCS real-time control model.
[0106] Set electrode parameters and placement. Based on the 10-10 system or personalized brain region localization (such as Broca's area), define the position (F7, F8), shape (usually rectangular or circular), size (such as 5 cm × 4 cm), and current intensity (such as ±1 mA) of the anode and cathode electrodes.
[0107] After obtaining the potential distribution, the electric field strength and current density can be further calculated:
[0108]
[0109]
[0110] When the amplitude of the stimulation current changes, since the Poisson equation responds linearly to the injected current, the electric field can be expressed as:
[0111]
[0112] in This represents the electric field distribution when the injected current is 1 mA, and is used to quickly predict the electric field changes under any current intensity.
[0113] S3: Construct a large-scale dynamic mean-field model (DMF) integrating neurotransmitter dynamics and synaptic plasticity mechanisms to obtain neural activity; divide brain regions based on brain atlases and construct neural network topologies. Incorporate neurotransmitter regulation mechanisms into the node dynamics model to simulate the effects of changes in neurotransmitter concentrations such as glutamate and GABA on synaptic input.
[0114] Specifically, this step is based on neurotransmitter dynamics and synaptic plasticity mechanisms. Brain region segmentation and network structure setting are performed using brain atlases (such as AAL-90) to divide brain regions and construct a neural network topology. Connection weights can be derived from DTI, SC data, or the default structure connection matrix. Through node dynamics modeling, each brain region is modeled as a mean-field unit containing a population of excitatory (E) and inhibitory (I) neurons, thus forming a mean-field model of the neuronal population.
[0115] A large-scale dynamic mean-field (DMF) model integrating neurotransmitter dynamics and synaptic plasticity mechanisms was constructed to describe whole-brain neural activity in patients with post-stroke aphasia. Its nodal dynamic equations are expressed as follows:
[0116]
[0117]
[0118]
[0119] Where i = 1, 2, ..., N, and N = 90, representing the i-th cortical or subcortical brain region delineated based on the AAL-90 brain atlas. In the dynamic mean-field model, NMDA-mediated synaptic gating variables... The dynamics are determined by the time constant. Control, synaptic activation rate parameter set to . This represents the average synaptic gating variable mediated by NMDA in the i-th brain region. Essentially, it is a synaptic-level state variable, with a value range of 0 ≤ 0. ≤1. Formula (6) includes neurodynamics and synaptic parameters. This represents the synaptic time constant, with a value of 100 ms; This represents the synaptic activation rate constant, with a value of 0.641 / 1000ms. -1 σ represents the noise amplitude, with a value of 0.001 nA. This represents Gaussian white noise with zero mean and unit variance. and The dynamics of NMDA are determined by σ=0.001, which is the key setting for resting-state oscillations.
[0120] It is an input-output function, in the form of a nonlinear sigmoid, representing the average firing rate of a local neuron population. Formula (7) includes the parameters of the input-output function, where a, b, and d represent the parameters of the input-output function, with values of 270 nC. −1 , 108Hz, 0.154s, where nC−1 =(10 -9 C) −1 =10 9 C −1, , 1 Hz = 1 s⁻¹. The input current is mapped to the average discharge rate (Hz), which is approximately linear in the low input region and saturates in the high input region, equivalent to the population behavior of LIF neurons.
[0121] This represents the total synaptic input current received by the i-th brain region, in A = Cs -1 10 9 C −1 ×(C s⁻¹) = 10 9 s⁻¹. Formula (8) includes synaptic coupling and network parameters. This represents the intensity of local excitatory self-feedback, with a value of 0.9; The value represents the unit synaptic coupling strength, which is 0.2609 nA; G represents the global coupling strength, with a scanning parameter range of 0-3 and an optimal coupling strength of 0.69. This represents the normalized structural connectivity matrix extracted based on DTI. This represents the external input current, with a value of 0.3nA, used to maintain a low discharge steady state. =0.9 is to ensure the stability of the local loop; It is obtained by the spiking → mean-field mapping; G is the most important adjustable parameter, used to fit the full-field (FC). This is to adjust the system to a low-discharge (2–3 Hz) resting state. An initial value [w] is set. G σ]=[0.9 0.3 0.690.001] T It is used to simulate FC and SC.
[0122] Neurotransmitter regulation mechanisms are incorporated into node dynamics models to simulate the effects of changes in neurotransmitter concentrations, such as glutamate (an excitatory neurotransmitter) and GABA (an inhibitory neurotransmitter), on synaptic input. Excitatory-inhibitory neuronal population dynamics models further divide local neuronal populations within each brain region into excitatory (E) and inhibitory (I) subpopulations. The kinetic equations are based on chemical synaptic input and membrane potential changes; these equations are commonly used to describe the dynamic behavior of neuronal populations, for example:
[0123]
[0124]
[0125] in: These represent the average activity levels of excitatory and inhibitory neuronal populations, respectively, and are usually expressed as average firing rates; The time constant of the excitability (NMDA) gate variable is 100 ms, indicating that the population responds slowly to changes in input. This represents the time constant of the inhibition (GABA-A) gated variable, with a value of 10 ms, indicating that the population responds quickly to changes in input. These are functions The derivative with respect to time t represents the rate of change of the average activity level (usually expressed as the firing rate of neurons) of the excitatory and inhibitory neuronal populations over time. , : represent the nonlinear activation functions of excitatory and inhibitory neuron populations, respectively, used to describe how input current is converted into population firing rate, typically exhibiting S-shaped characteristics. Here... It is the formula (7) Replace with The result, substituted into the following:
[0126]
[0127] Where a = 270nC −1 b=108Hz, d=0.154s.
[0128] The total input currents received by the excitatory and inhibitory neuronal populations are expressed as follows:
[0129]
[0130]
[0131]
[0132] Among them, J EE J EI J IE J II : Represents the effective coupling strength of different types of synapses, reflecting the neurotransmitter regulation effect, with values of 1.5, 1.0, 1.2, and 0.5 respectively; , : Represents the external input current, used to simulate the external input current of tDCS modulation, with values of 0.3nA and 0.3nA respectively.
[0133] These equations describe how neuronal populations change over time based on received input currents and their current activity levels. These equations allow for the simulation and study of how neuronal populations in different brain regions interact and influence each other's activity. By implementing the modulation of synaptic plasticity using the tDCS electric field, the electric field strength is mapped to changes in synaptic weights (e.g., LTP / LTD).
[0134] The Pearson correlation coefficient is a statistic that measures the degree of linear correlation between two variables, with values ranging from -1 to 1. Here we use the Pearson correlation coefficient to quantify the simulated FC (full-force coefficient). ) and empirical FC ( The fit between ) ).
[0135]
[0136] Where n is a vector and The number of elements in the middle. and They are vectors and The i-th element in. and They are vectors and The average value.
[0137] The Expectation-Maximization (EM) algorithm is an iterative algorithm used to maximize the likelihood function in the presence of latent variables. In neural network simulations, the EM algorithm is used to optimize model parameters to maximize the goodness of fit between simulated functional connectivity (FC) and empirical functional connectivity (FC). In the EM algorithm, we typically maximize the goodness of fit by iterating through the following steps: :
[0138] E-step (Expectation step): Estimates the expected value of the latent variables based on the current parameters.
[0139] M-step (Maximization step): Maximizes the expected log-likelihood function, thereby updating the parameter estimates.
[0140] In each iteration, we calculate To evaluate the goodness of fit between simulated FC and empirical FC, and based on The model parameters are adjusted based on the changes in [the parameters]. This process is repeated until [the changes are completed]. Convergence or reaching the preset number of iterations, thus obtaining w, G and The optimal values for these four parameters.
[0141] S4: Establish a large-scale neurovascular coupling model of the entire brain that links neural activity in brain regions with hemodynamic responses. The core objective is to integrate the neural activity (excitatory activity) of brain regions obtained in S3 with the data from S3. Using this as a driving signal, the Balloon-Windkessel model describes the changes in hemodynamics, further generating a BOLD signal associated with fMRI. The Balloon-Windkessel model is a commonly used mathematical model in neurovascular coupling (NVC) research, which can link neural activity in brain regions with hemodynamic responses (such as local blood flow and oxygen consumption), thereby affecting the BOLD signal.
[0142] A complete physiological mechanism refers to the close coordination between neuronal activity and local cerebral blood flow (CBF), ensuring that active brain regions receive an adequate supply of oxygen and glucose. In patients with post-stroke aphasia, the NVC mechanism is often impaired due to neuronal damage, vascular dysfunction, or astrocyte abnormalities, leading to a mismatch between neural activity and blood flow supply, thereby limiting the regulatory effect. The core idea of this model is to drive blood flow and oxygen exchange processes through changes in metabolic demand induced by neural activity. The Balloon-Windkessel model integrates neural networks, metabolic processes, and hemodynamics, and is a key tool for studying the interaction between cerebral blood flow and neural activity.
[0143] Specifically, the Balloon–Windkessel hemodynamic model proposed by Friston et al. was used to convert local neural synaptic activity into observable BOLD signals. This model describes the biophysical coupling between changes in vasodilation, cerebral blood flow, blood volume, and deoxyhemoglobin levels induced by neural activity. For the i-th brain region, its synaptic activity level is denoted as […]. Directly take it as the NMDA-mediated average synaptic gating variable Vasodilatory signals are recorded as However, it still represents the total synaptic input current in formula (3), which is simply through... = It indirectly determines subsequent BOLD-related variables such as vasodilation signals, blood flow, and blood volume.
[0144]
[0145]
[0146]
[0147]
[0148] The meaning of each variable in formulas (16)-(19) is as follows: Indicates a vasodilatory signal. Representation function The derivative with respect to time t; Indicates the first The neural synaptic activity input of each brain region (from the DMF model) is generally taken as follows: = It will trigger vasodilatory signals. The increase in this signal is influenced by self-regulatory feedback, which indirectly determines subsequent BOLD-related variables such as vasodilation signals, blood flow, and blood volume. Indicates cerebral blood flow. Representation function The derivative with respect to time t, at rest =1; Indicates cerebral blood volume. Representation function The derivative with respect to time t, at rest =1; This indicates the deoxyhemoglobin content. Representation function The derivative with respect to time t, at rest =1.
[0149] The values for each parameter in formulas (16)-(19) are assigned as follows: Indicates the first The signal attenuation rate for each region is set to 0.65s. -1 ; This represents the strength of the self-adjusting feedback, with a value of 0.41s. -1 ; This represents the hemodynamic mean transit time, with a value of 0.98 s; This represents the Grubb index (blood volume-to-blood flow index), with a value of 0.32. This represents the resting oxygen extraction fraction, with a value of 0.34.
[0150] Based on the above hemodynamic variables, the first The BOLD signal in each brain region is represented as follows:
[0151]
[0152] in, Indicates resting blood volume fraction. =0.02; It is the resting oxygen extraction fraction, obtained from formula (19). =0.34; Indicates the deoxyhemoglobin content; Indicates cerebral blood volume; , and This represents a set of weighting coefficients related to blood oxygen dynamics. =2.38, =2, =0.88. This BOLD model is a local model, which only maps the synaptic activity of the corresponding brain region to the BOLD signal of that brain region, and does not directly introduce cross-brain region blood flow coupling.
[0153] In summary, this model describes the biophysical coupling relationship between changes in vasodilation, cerebral blood flow, blood volume, and deoxyhemoglobin content induced by neural activity.
[0154] S5: Combining finite element simulation and neurovascular coupling model, construct an individualized tDCS real-time control model.
[0155] Specifically, the electric field distribution (finite element simulation), neurodynamic response (DMF model), and blood flow changes (balloon model) obtained in the previous steps are integrated to establish a system model capable of simulating and predicting individual brain responses and used for feedback regulation. By monitoring changes and implementing effective regulation, during tDCS stimulation, the system adjusts stimulation parameters by monitoring real-time changes in neural activity and cerebral blood flow, ensuring that each cycle of tDCS stimulation adapts to the individual's brain response, thereby optimizing the regulatory effect and avoiding the problem of mismatch between stimulation parameters and neurovascular coupling.
[0156] An electric field-neurodynamic coupling interface was established, and the individualized electric field vector distribution obtained from finite element simulation was utilized. The normal component (scalar value) of the electric field vector at each cortical grid point, normE, is used as an external input to the excitatory units in the DMF model. An excitatory unit is a neural population composed of excitatory neurons (usually glutamatergic cells) that generates positive feedback responses to input signals, promoting the propagation of neural activity.
[0157] The formula is , For the first The electrical modulation input received by the excitatory units in each brain region; This represents the normal component (average value) of the electric field in this brain region. This is a coupling strength parameter used to adjust the magnitude of the effect of the electric field on the neuron population.
[0158] By using a personalized tDCS real-time control model, the personalized stimulation current intensity is obtained, thereby obtaining the optimal stimulation parameters in the personalized tDCS real-time control scheme for patients with post-stroke aphasia (PSA).
[0159] To maximize the intervention effect of transcranial direct current stimulation (tDCS) on patients with post-stroke aphasia, the model parameters and stimulation parameters were jointly optimized based on the aforementioned neurodynamic model and neurovascular coupling model, and personalized stimulation protocols were developed accordingly. The optimization aimed to enhance the coupling strength between neural activity and blood flow response in the target brain region, thereby improving local cerebral blood flow and metabolic response and promoting functional recovery.
[0160]
[0161]
[0162]
[0163] in, Indicates the first The total synaptic input current of the excitatory neuronal population received by a brain region at time t reflects the level of neural activity or the strength of synaptic drive in that brain region. Indicates the first The blood oxygen metabolism level of a brain region at time t is used to characterize the intensity of the hemodynamic-metabolic response of that brain region and can be correlated with blood oxygen consumption rate (CMRO2) or hemodynamic state. = ; This indicates that within the selected time window or stimulation cycle, the first... The maximum value of the excitatory synaptic input current in each brain region is used to normalize neural activity. This indicates that within the corresponding time window, the first... The maximum value of the blood oxygen metabolism response in each brain region was used to normalize the blood flow / metabolic response. (Known) and Parameters used: w=0.9 =0.2609nA, =0.3nA, a=270nC−1, b=108Hz, d=0.154s, =100ms =0.641 / 1000ms⁻¹, σ=0.001nA, therefore The above parameter values can be used to determine this. Formula (21) quantitatively characterizes the coupling strength between neural activity and blood flow response in the target brain region by calculating the normalized product of excitatory neural activity and blood flow metabolic response, and is used to characterize the neurovascular coupling (NVC) state. Coupling degree The larger the value, the stronger the neural activity and hemodynamic response. These goals are achieved by minimizing a specific objective function, which is typically the difference in BOLD signal, the degree of coupling between neural activity and hemodynamics, and the degree of recovery of motor function.
[0164] Indicates the first BOLD fMRI signals of brain regions before transcranial direct current stimulation (tDCS) intervention; Indicates the first BOLD fMRI signal of a brain region after transcranial direct current stimulation (tDCS) intervention; |·| represents the absolute value of the change in BOLD signal before and after stimulation, used to measure the degree of functional change caused by stimulation. The parameters used for outputting BOLD are known: =0.02, =0.34, =2.38, =2, 0.88, =0.98, =0.32, =0.65 =0.41, therefore as well as All of these can be determined through the above parameter values. Formula (22) is used to quantitatively describe the degree of change in the oxygen-dependent signal (BOLD) of the target brain region before and after transcranial direct current stimulation, so as to reflect the regulatory effect of stimulation on brain function.
[0165] This represents the comprehensive optimization objective function, used to guide the optimization of personalized transcranial direct current stimulation parameters; These are weighting coefficients, with values of 0.5 and 0.5 respectively, and satisfying the following conditions: + =1, used to balance the relative importance of neurovascular coupling indicators and BOLD functional response indicators in the optimization process. It can be adjusted according to the actual situation to balance the relative importance between different objectives. Formula (23) optimizes the neurovascular coupling strength and BOLD functional response changes together, which can improve the functional observability and applicability of the stimulation scheme while ensuring the rationality of the neural regulation mechanism.
[0166] During the stimulation parameter optimization phase, the system's core objective is to enhance the neurovascular coupling strength between neural activity and hemodynamic response in the target brain region, while simultaneously inhibiting abnormal activation in non-target brain regions. After each cycle of transcranial direct current stimulation, the system evaluates the objective function based on real-time or offline acquired neuroimaging and hemodynamic data, and dynamically updates stimulation parameters such as current intensity and electrode position. These stimulation parameter updates follow the following relationship:
[0167]
[0168]
[0169] in, The stimulation parameters for the k-th stimulation cycle are given, and the parameter vector includes at least the current intensity, electrode position, and stimulation duration. To control the adjustment range of the updated stimulus parameters each time; is the parameter adjustment coefficient, with a value range of [0.01, 0.3], and here it is set to 0.1. It is used to control the magnitude of stimulus parameter update and avoid over-adjustment leading to stimulus instability; ∇J is the change (or approximate gradient) of the objective function in adjacent stimulus cycles. The value range is 1-2mA. =1mA; The value range is 10-20 min. =20min; Corresponding to the positions of the anode and cathode of the electrode plate, ={F7, F8}. Through the above method, the stimulation parameters are adaptively adjusted, thereby gradually optimizing the neurovascular coupling state and promoting functional recovery.
[0170] Calculate the FCD (Frechet Class Distance) to determine the size N of the dataset, i.e., the number of samples in the dataset. For each sample i (from 1 to N), calculate the empirical cumulative distribution function (FCD). ) and simulated cumulative distribution function ( ).
[0171] For each pair (i,j), calculate the square of the difference between the empirical cumulative distribution function and the simulated cumulative distribution function, i.e. Sum the squares of the differences between all i and j: Take the square root of the result from the previous step and divide it by N to obtain the FCD value:
[0172]
[0173] The smaller the FCD value, the smaller the difference between the empirical cumulative distribution function and the simulated cumulative distribution function, and the better the model fits.
[0174] In summary, this scheme realizes a closed-loop neural regulation logic of "monitoring-regulation-optimization", ensuring that stimulation parameters are dynamically adjusted according to individual neurophysiological feedback, thereby improving the efficiency of brain function recovery and the regulatory effect.
[0175] Example 2
[0176] This embodiment provides a personalized real-time control device based on a neurovascular coupling model of transcranial direct current stimulation response, including:
[0177] Data acquisition module: used to extract multimodal data composed of magnetic resonance imaging (T1w) data, diffusion tensor imaging (DTI) data, and resting-state functional magnetic resonance imaging (rs-fMRI) data from healthy individuals and patients with post-stroke aphasia;
[0178] Preprocessing and extracting the connectivity matrix module: This module is used to preprocess multimodal data and extract the structural and functional connectivity matrices respectively.
[0179] Finite element simulation module: used to simulate the electric field distribution of tDCS in the brains of stroke patients with aphasia;
[0180] Large-scale dynamic mean-field model building module: used to construct large-scale dynamic mean-field models that integrate neurotransmitter dynamics and synaptic plasticity mechanisms;
[0181] Neural activity output module: used for neural activity obtained from large-scale dynamic mean-field models;
[0182] Whole-brain large-scale neurovascular coupling model establishment module: used to link neural activity and hemodynamic responses in brain regions;
[0183] tDCS Real-Time Control Module: Used to combine finite element simulation and neurovascular coupling model to construct individualized tDCS real-time control model.
[0184] Example 3
[0185] This embodiment provides a computer program product, including a computer program / instruction, characterized in that, when the computer program / instruction is executed by a processor, it implements a personalized real-time control method based on a neurovascular coupling model of transcranial direct current stimulation response.
[0186] The technical features of the above embodiments can be combined in any way (as long as there is no contradiction in the combination of these technical features). For the sake of brevity, not all possible combinations of the technical features in the above embodiments are described; these embodiments not explicitly written should also be considered to be within the scope of this specification.
Claims
1. A personalized real-time control method based on a neurovascular coupling model of transcranial direct current stimulation response, characterized in that, include: Step S1: Acquire multimodal data composed of magnetic resonance imaging (T1w) data, diffusion tensor imaging (DTI) data, and resting-state functional magnetic resonance imaging (rs-fMRI) data from healthy individuals and patients with post-stroke aphasia, respectively, and preprocess the multimodal data, then extract the structural and functional connectivity matrices respectively. Step S2: Based on the preprocessed multimodal data, perform finite element simulation modeling to simulate the electric field distribution of tDCS in the brains of stroke aphasia patients. Step S3: Construct a large-scale dynamic mean-field model (DMF) that integrates neurotransmitter dynamics and synaptic plasticity mechanisms, and obtain neural activity; Step S4: Establish a large-scale neurovascular coupling model of the whole brain that links neural activity and hemodynamic response in brain regions; Step S5: Combine finite element simulation and neurovascular coupling model to construct an individualized tDCS real-time control model.
2. The personalized real-time control method for a neurovascular coupling model based on transcranial direct current stimulation response according to claim 1, characterized in that, In step S1, the multimodal data composed of magnetic resonance imaging (T1w) data, diffusion tensor imaging (DTI) data, and resting-state functional magnetic resonance imaging (rs-fMRI) data of patients with post-stroke aphasia was obtained from the OpenNeuro database; the structural connectivity matrix and functional connectivity matrix were extracted using the AAL-90 brain atlas in Freesurfer.
3. The personalized real-time control method for a neurovascular coupling model based on transcranial direct current stimulation response according to claim 2, characterized in that, In step S1, when extracting the structural connectivity matrix and functional connectivity matrix using the AAL-90 brain atlas in Freesurfer, the cerebral cortex is divided into 90 regions.
4. The personalized real-time control method for a neurovascular coupling model based on transcranial direct current stimulation response according to claim 1, characterized in that, In step S2, the multimodal data of healthy individuals and patients with post-stroke aphasia are preprocessed, including visualization, artifact removal, temporal registration, head movement correction, structural-functional registration, and spatial normalization. Finite element simulation modeling was used to simulate the electric field strength and current distribution in the brain of a stroke-induced aphasia patient using tDCS.
5. The personalized real-time control method for a neurovascular coupling model based on transcranial direct current stimulation response according to claim 1, characterized in that, In step S3, a large-scale dynamic mean field (DMF) model integrating neurotransmitter dynamics and synaptic plasticity mechanisms is constructed to describe the whole-brain neural activity of patients with post-stroke aphasia. Its nodal dynamic equations are expressed as follows: Where i = 1, 2, ..., N, and N = 90, representing the i-th cortical or subcortical brain region classified based on the AAL-90 brain atlas; in the dynamic mean-field model, NMDA-mediated synaptic gating variables The dynamics are determined by the time constant. Control, synaptic activation rate parameter set to This represents the average synaptic gating variable mediated by NMDA in the i-th brain region. Essentially, it is a synaptic-level state variable, with a value range of 0 ≤ 0. ≤1; Formula (1) includes neurodynamics and synaptic parameters, This represents the synaptic time constant, with a value of 100 ms; This represents the synaptic activation rate constant, with a value of 0.641 / 1000ms. -1 σ represents the noise amplitude, with a value of 0.001 nA. This represents Gaussian white noise with zero mean and unit variance. and The dynamics of NMDA are determined by σ=0.001, which is the key setting for resting-state oscillations; It is an input-output function, in the form of a nonlinear sigmoid, representing the average firing rate of a local neuron population; Formula (2) includes the input-output function parameters, where a, b, and d represent the parameters of the input-output function, with values of 270 nC respectively. −1 , 108Hz, 0.154s, where nC −1 =(10 -9 C) −1 =10 9 C −1, 1 Hz = 1 s⁻¹ The input current is mapped to the average discharge rate (Hz), which is approximately linear in the low input region and saturates in the high input region, equivalent to the population behavior of LIF neurons; This represents the total synaptic input current received by the i-th brain region, in A = Cs -1 10 9 C −1 ×(C s⁻¹) =10 9 s⁻¹; Formula (3) includes synaptic coupling and network parameters, This represents the intensity of local excitatory self-feedback, with a value of 0.9; The value represents the unit synaptic coupling strength, which is 0.2609 nA; G represents the global coupling strength, with a scanning parameter range of 0-3 and an optimal coupling strength of 0.
69. This represents the normalized structural connectivity matrix extracted based on DTI. This represents the external input current, with a value of 0.3nA. =0.9; It is obtained by the spiking → mean-field mapping; G is an adjustable parameter used to fit the full-field (FC). This is to regulate the system to a low discharge state, i.e., a resting state of 2–3 Hz; To further understand the modulating effects of transcranial direct current stimulation (tDCS) and changes in the neurochemical environment under stroke pathological conditions on local neurodynamics, a neurotransmitter regulation mechanism was introduced based on the nodal dynamics model to simulate the excitatory neurotransmitter glutamate and the inhibitory neurotransmitter... The effect of changes in γ-aminobutyric acid (GABA) concentration on synaptic input and excitation-inhibition balance; The modulation mechanism of neurotransmitters on synaptic gain is crucial. Physiologically, glutamate mainly mediates excitatory synaptic transmission, while GABA mainly mediates inhibitory synaptic transmission. Therefore, in the model, the effect of neurotransmitter concentration changes is reflected by adjusting the synaptic gain parameters of different types of synapses. Increased glutamate concentration mainly enhances excitatory synaptic gain (J). EE J EI Increased GABA concentration primarily enhances inhibitory synaptic gain (J). IE J II ); The excitatory-inhibitory neuronal population dynamics model refers to further dividing the local neuronal population within each brain region into two subpopulations: excitatory (E) and inhibitory (I). The dynamic equations are based on chemical synaptic input and membrane potential changes, and these equations are commonly used to describe the dynamic behavior of neuronal populations. in: These represent the average activity levels of excitatory and inhibitory neuronal populations, respectively, and are usually expressed as average firing rates; The time constant of the excitability (NMDA) gate variable is 100 ms, indicating that the population responds slowly to changes in input. This represents the time constant of the inhibition (GABA-A) gated variable, with a value of 10 ms, indicating that the population responds quickly to changes in input. These are functions The derivative with respect to time t represents the average activity level of the excitatory and inhibitory neuronal populations, and is expressed as the rate of change of the neuron's firing rate over time. , : These represent the nonlinear activation functions of excitatory and inhibitory neuronal populations, respectively, used to describe how input current is converted into population firing rate, and typically exhibit S-shaped characteristics; here It is the formula (2) Replace with The result, substituted into the following: Where a = 270nC −1 b=108Hz, d=0.154s; The total input currents received by the excitatory and inhibitory neuronal populations are expressed as follows: Among them, J EE J EI J IE J II : Represents the effective coupling strength of different types of synapses, reflecting the neurotransmitter regulation effect, with values of 1.5, 1.0, 1.2, and 0.5 respectively; , : Represents the external input current, used to simulate the external input current of tDCS modulation, with values of 0.3nA and 0.3nA respectively.
6. The personalized real-time control method for a neurovascular coupling model based on transcranial direct current stimulation response according to claim 5, wherein the kinetic equation is expressed as in step S4, characterized in that, Using the Balloon–Windkessel hemodynamic model, local neural synaptic activity was analyzed. The signal is converted into an observable BOLD signal. The model describes the biophysical coupling relationship between changes in vasodilation, cerebral blood flow, blood volume, and deoxyhemoglobin content induced by neural activity, as follows: The meaning of each variable in formulas (10)-(13) is as follows: Indicates a vasodilatory signal. Representation function The derivative with respect to time t; Indicates the first The neural synaptic activity input of each brain region is taken = It will trigger vasodilatory signals. The increase in this signal is influenced by self-regulatory feedback, which indirectly determines subsequent BOLD-related variables such as vasodilation signals, blood flow, and blood volume. Indicates cerebral blood flow. Representation function The derivative with respect to time t, at rest =1; Indicates cerebral blood volume. Representation function The derivative with respect to time t, at rest =1; This indicates the deoxyhemoglobin content. Representation function The derivative with respect to time t, at rest =1; The values for each parameter in formulas (10)-(13) are assigned as follows: Indicates the first The signal attenuation rate for each region is set to 0.65s. -1 ; This represents the strength of the self-adjusting feedback, with a value of 0.41s. -1 ; This represents the hemodynamic mean transit time, with a value of 0.98 s; The Grubb index represents the relationship between blood volume and blood flow, with a value of 0.
32. This represents the resting oxygen extraction fraction, with a value of 0.
34. Based on the above hemodynamic variables, the first The BOLD signal in each brain region is represented as follows: in, Indicates resting blood volume fraction. =0.02; It is the resting oxygen extraction fraction, obtained from formula (11). =0.34; Indicates the deoxyhemoglobin content; Indicates cerebral blood volume; , and This represents a set of weighting coefficients related to blood oxygen dynamics. =2.38, =2, =0.88; This BOLD model is a local model, which only maps the synaptic activity of the corresponding brain region to the BOLD signal of that brain region, and does not directly introduce cross-brain region blood flow coupling.
7. The personalized real-time control method for a neurovascular coupling model based on transcranial direct current stimulation response according to claim 6, characterized in that, In step S5, to maximize the intervention effect of transcranial direct current stimulation (tDCS) on patients with post-stroke aphasia, the model parameters and stimulation parameters are jointly optimized based on the aforementioned neurodynamic model and neurovascular coupling model, and a personalized stimulation plan is formulated accordingly. The goal of the optimization is to enhance the coupling strength between neural activity and blood flow response in the target brain region, thereby improving local cerebral blood flow and metabolic response and promoting functional recovery. in, Indicates the first The total synaptic input current of the excitatory neuronal population received by a brain region at time t reflects the level of neural activity or the strength of synaptic drive in that brain region. Indicates the first The blood oxygen metabolism level of a brain region at time t is used to characterize the intensity of the hemodynamic-metabolic response in that brain region and is related to the blood oxygen consumption rate or hemodynamic state. = ; This indicates that within the selected time window or stimulation cycle, the first... The maximum value of the excitatory synaptic input current in each brain region is used to normalize neural activity. This indicates that within the corresponding time window, the first... The maximum value of the blood oxygen metabolism response in each brain region was used to normalize the blood flow / metabolic response; known and Parameters used: w=0.9 =0.2609nA, =0.3nA, a=270nC−1, b=108Hz, d=0.154s, =100ms =0.641 / 1000ms⁻¹, σ=0.001nA, therefore The above parameter values can be used to determine the coupling strength between neural activity and blood flow response in the target brain region by calculating the normalized product of excitatory neural activity and blood flow metabolic response, which is used to characterize the neurovascular coupling (NVC) state; coupling degree The larger the value, the stronger the neural activity and hemodynamic response; these goals are achieved by minimizing a specific objective function, which is the difference in BOLD signal, the coupling degree of neural activity and hemodynamics, and the degree of recovery of motor function. Indicates the first BOLD fMRI signals of individual brain regions before transcranial direct current stimulation (tDCS) intervention; Indicates the first BOLD fMRI signals of several brain regions after transcranial direct current stimulation (tDCS) intervention; |·| represents the absolute value of the change in BOLD signal before and after stimulation, used to measure the degree of functional change caused by stimulation; the parameters used to output BOLD are known: =0.02, =0.34, =2.38, =2, 0.88, =0.98, =0.32, =0.65 =0.41, therefore as well as All of these can be determined through the above parameter values; Formula (16) is used to quantitatively describe the degree of change in the oxygenation level dependent signal (BOLD) of the target brain region before and after transcranial direct current stimulation; This represents the comprehensive optimization objective function, used to guide the optimization of personalized transcranial direct current stimulation parameters; These are weighting coefficients, with values of 0.5 and 0.5 respectively, and satisfying the following conditions: + =1, used to balance the relative importance of neurovascular coupling index and BOLD functional response index in the optimization process; Formula (17) optimizes the neurovascular coupling strength and BOLD functional response changes together; During the stimulation parameter optimization phase, the system's core objective is to enhance the neurovascular coupling strength between neural activity and hemodynamic response in the target brain region, while simultaneously inhibiting abnormal activation in non-target brain regions. After each cycle of transcranial direct current stimulation, the system evaluates the objective function based on real-time or offline acquired neuroimaging and hemodynamic data, and dynamically updates the stimulation parameters, including current intensity and electrode position. These updates follow the following relationship: in, The stimulation parameters for the k-th stimulation cycle are given, and the parameter vector includes at least the current intensity, electrode position, and stimulation duration. To control the adjustment range of the updated stimulus parameters each time; is the parameter adjustment coefficient, with a value range of [0.01, 0.3], and here it is set to 0.
1. It is used to control the magnitude of stimulus parameter updates and avoid over-adjustment leading to stimulus instability; ∇J is the change (or approximate gradient) of the objective function in adjacent stimulus cycles. The value range is 1-2mA. =1mA; The value range is 10-20 min. =20min; Corresponding to the positions of the anode and cathode of the electrode plate, ={F7, F8}.
8. A personalized real-time control device based on a neurovascular coupling model of transcranial direct current stimulation response, characterized in that, include: Data acquisition module: used to extract multimodal data composed of magnetic resonance imaging (T1w) data, diffusion tensor imaging (DTI) data, and resting-state functional magnetic resonance imaging (rs-fMRI) data from healthy individuals and patients with post-stroke aphasia; Preprocessing and extracting the connectivity matrix module: This module is used to preprocess multimodal data and extract the structural and functional connectivity matrices respectively. Finite element simulation module: used to simulate the electric field distribution of tDCS in the brains of stroke patients with aphasia; Large-scale dynamic mean-field model building module: used to construct large-scale dynamic mean-field models that integrate neurotransmitter dynamics and synaptic plasticity mechanisms; Neural activity output module: used for neural activity obtained from large-scale dynamic mean-field models; Whole-brain large-scale neurovascular coupling model establishment module: used to link neural activity and hemodynamic responses in brain regions; tDCS Real-Time Control Module: Used to combine finite element simulation and neurovascular coupling model to construct individualized tDCS real-time control model.
9. A computer program product, comprising a computer program / instructions, characterized in that, When the computer program / instructions are executed by the processor, they implement the steps of the method according to any one of claims 1 to 7.