An individualized brain functional area boundary positioning method and system

Abstract

The application belongs to the technical field of image recognition, and discloses a kind of individualized brain function area boundary positioning method and system;Obtain structural MRI data, functional MRI data and diffusion tensor imaging data under the same voxel space, generate functional activation map under resting state condition, and determine coarse functional area mask in target anatomical region according to functional activation map;Coarse functional area mask is expanded in space and edge band is extracted, so as to demarcate a candidate boundary band voxel set around coarse functional area mask;The boundary feature vector corresponding to each voxel is constructed;Boundary feature vector is input into pre-trained boundary discriminant network, and the boundary probability value of each voxel is obtained, and the boundary probability graph of individualized brain function area is generated;Based on boundary probability graph, construct boundary energy function, combine energy minimization segmentation mechanism, divide candidate boundary band voxel set, three-dimensional reconstruction is carried out in individual structure space, and the boundary surface of individualized brain function area is obtained.

Description

Technical Field

[0001] This invention relates to the field of image recognition technology, and more specifically, to a method and system for individualized brain functional area boundary localization. Background Technology

[0002] With the widespread adoption of medical imaging technologies such as functional magnetic resonance imaging (fMRI), diffusion tensor imaging (DTI), and high-resolution structural MRI, brain functional area localization based on image recognition and segmentation is playing an increasingly important role in preoperative neurosurgical planning, epilepsy focus resection, brain tumor resection design, brain stimulation target planning (such as DBS and TMS), and cognitive neuroscience research. In clinical and research practice, it is often necessary to accurately obtain the spatial boundaries of functional areas such as motor, language, and visual cortices at the individual level in order to maximize lesion resection or optimize stimulation targets while ensuring functional preservation.

[0003] Most current mainstream brain functional area localization procedures are based on standard templates and population brain functional atlases. While these methods have some merit in population statistical analysis, they have significant limitations in accurately identifying individualized functional boundaries.

[0004] While existing registration algorithms (such as nonlinear deformation and feature-point-based registration) can align cortical gyri and sulci morphology well at the structural level, significant positional drift and morphological differences exist between brain functional areas across individuals. Even within the same anatomical gyrus, the true activation center and boundaries of functional areas often differ from the population template. Template maps can only provide an "average" functional distribution, lacking a detailed depiction of the transition zones at the boundaries of functional areas. They cannot accurately reflect the true functional boundaries of individual subjects through simple spatial transformations, leading to problems such as "functional area positional shift, blurred boundaries, and even cross-area misjudgment" at the individual level.

[0005] Therefore, an individualized method and system for locating the boundaries of brain functional areas are designed. Summary of the Invention

[0006] To overcome the aforementioned deficiencies of the prior art and to achieve the above objectives, the present invention provides the following technical solution: a method for individualized brain functional area boundary localization, comprising:

[0007] Step S1: Collect multimodal brain imaging data of the target individual, and register and preprocess them to obtain structural MRI data, functional MRI data and diffusion tensor imaging data in the same voxel space;

[0008] Step S2: Based on functional MRI data, generate a functional activation map under resting conditions, and determine a coarse functional area mask within the target anatomical region according to the functional activation map;

[0009] Step S3: Spatial expansion and edge band extraction are performed on the coarse functional region mask to delineate a set of candidate boundary band voxels around the coarse functional region mask;

[0010] Step S4: For each voxel in the candidate boundary zone voxel set, extract multi-scale image features from structural MRI data, functional MRI data, and diffusion tensor imaging data to construct the boundary feature vector corresponding to each voxel.

[0011] Step S5: Input the boundary feature vector into the pre-trained boundary discrimination network to obtain the boundary probability value of each voxel and generate a boundary probability map of individualized brain functional areas.

[0012] Step S6: Construct a boundary energy function based on the boundary probability map, and combine it with the energy minimization segmentation mechanism to divide the candidate boundary zone voxel set to obtain the inner region voxel set, the outer region voxel set, and the boundary voxel set; and based on this, perform three-dimensional reconstruction in the individual structural space to obtain the individualized brain functional area boundary surface, and output this as the final individualized brain functional area boundary localization result.

[0013] Preferably, the method for acquiring multimodal brain imaging data of the target individual and registering and preprocessing it includes:

[0014] For a single target individual, a magnetic resonance imaging (MRI) device is used to scan and acquire multimodal brain imaging data;

[0015] Skull removal, noise suppression, and intensity non-uniformity correction were performed on the structural MRI data to obtain preprocessed structural MRI data.

[0016] The functional MRI data were subjected to time-layer correction, head motion correction and time drift removal, and then registered to the voxel space of the preprocessed structural MRI data.

[0017] Motion eddy current correction and tensor fitting were performed on the diffusion tensor imaging data, and the data were registered to the voxel space of the preprocessed structural MRI data.

[0018] Intensity normalization and spatial resampling were performed on the registered structural MRI data, functional MRI data, and diffusion tensor imaging data to ensure they had the same voxel resolution and coordinate system.

[0019] Preferably, the method for determining the coarse functional area mask within the target anatomical region based on the functional activation map includes:

[0020] In the functional activation graph, an adaptive threshold segmentation is applied to the statistical parameter graph corresponding to the target functional paradigm to obtain the initial activation region.

[0021] Based on preprocessed structural MRI data, the cortical surface of the target anatomical region is extracted through cortical surface reconstruction, and the initial activated region is projected onto the cortical surface;

[0022] Morphological closing operations and connected region filtering are performed on the initial activated region after projection to remove discrete patches and obtain the functionally activated main region with a coherent topological structure.

[0023] The functional activation subject area is transformed from the cortical surface space to the voxel space to form a binary mask, which serves as a coarse functional area mask.

[0024] Preferably, the method for defining a candidate boundary band voxel set around the coarse functional region mask includes:

[0025] Perform expansion and erosion with a preset radius on the coarse functional area mask to obtain an expanded functional area mask and a contracted functional area mask, respectively;

[0026] Calculate the difference between the expanded functional region mask and the contracted functional region mask, and use the voxels in the difference set as the candidate boundary zone voxel set;

[0027] In the candidate boundary zone voxel set, based on the preprocessed structural MRI data, voxels that are less than a preset safe distance from the skull boundary are removed, and these are recorded as the final candidate boundary zone voxel set.

[0028] Preferably, the method for constructing the boundary feature vector corresponding to each voxel includes:

[0029] In the preprocessed structural MRI data, the gray-level gradient magnitude, cortical thickness estimate, and local curvature features of the neighborhood of each candidate boundary zone voxel in the candidate boundary zone voxel set are extracted.

[0030] From functional MRI data registered to voxel space, calculate the activation intensity and spatial gradient magnitude of each candidate boundary zone voxel, as well as the average correlation coefficient with the time series of neighboring voxels.

[0031] From the diffusion tensor imaging data registered to voxel space, calculate the anisotropy fraction, principal fiber orientation, and fiber-tracking-based structural connectivity strength characteristics of each candidate boundary zone voxel.

[0032] The extracted multi-scale image features are then concatenated in a preset order to form a boundary feature vector.

[0033] Preferably, the method in step S5 includes:

[0034] A boundary discrimination network is pre-built based on a deep learning network. The input is defined as a boundary feature vector, and the output is a ternary probability distribution, including the probability value of the inner region, the probability value of the boundary, and the probability value of the outer region.

[0035] The boundary discrimination network is pre-trained using training samples;

[0036] For the target individual, the boundary feature vector corresponding to each voxel in the candidate boundary band voxel set is input into the boundary discrimination network, and the corresponding ternary probability distribution is output.

[0037] Extract the boundary probability values ​​and rearrange them according to their voxel spatial positions to form a three-dimensional boundary probability map.

[0038] Preferably, the method for pre-training the boundary discrimination network using training samples includes:

[0039] The voxels containing the true boundaries of brain functional areas are marked as boundaries, the voxels extending a certain distance toward the functional areas are marked as inner regions, and the voxels extending a certain distance toward non-functional areas are marked as outer regions, thus forming training samples with three types of labels.

[0040] During the training phase, the network is trained using the boundary feature vector as input and the supervision signal constructed from the corresponding three-class labels as the target; at the same time, the boundary class is given higher weight in the loss function.

[0041] After the training process converges, the final boundary discrimination network model parameters are saved.

[0042] Preferably, the method for obtaining the inner region voxel set, the outer region voxel set, and the boundary voxel set includes:

[0043] A voxel graph is constructed using the set of voxels in the candidate boundary zone as graph nodes and neighborhood relationships as adjacent edges.

[0044] Based on the boundary probability map, the boundary probability value of each voxel is converted into a data item cost belonging to the boundary category; based on the preprocessed structural MRI data, the gray level difference and spatial distance between adjacent voxels are calculated, which together constitute the smoothing term cost; in this way, a boundary energy function defined on the three categories of labels: inner region, boundary, and outer region is constructed.

[0045] A graph cut-α extension is used to solve the global energy minimization problem of the boundary energy function, and an optimal class label is assigned to each node in the graph.

[0046] The set of voxels with the optimal label of the inner region is defined as the inner region voxel set, the set of voxels with the optimal label of the outer region is defined as the outer region voxel set, and the set of voxels with the optimal label of the boundary is defined as the boundary voxel set.

[0047] Preferably, the method for obtaining the individualized brain functional area boundary surface includes:

[0048] In voxel space, the set of boundary voxels is taken as the set of spatial target points;

[0049] The moving cube algorithm is used to extract the initial boundary surface from the set of spatial target points;

[0050] The initial boundary surface is smoothed using Laplace to obtain the individualized brain functional area boundary surface.

[0051] The individualized brain functional area boundary surfaces are mapped to the surgical navigation coordinate system or the standard medical imaging coordinate system for three-dimensional visualization and data export.

[0052] A personalized brain functional area boundary localization system includes:

[0053] The multimodal brain image preprocessing module is used to acquire multimodal brain image data of the target individual and perform registration and preprocessing.

[0054] The coarse functional area generation module is used to calculate functional activation maps based on functional MRI data and generate coarse functional area masks within the target anatomical region according to the functional activation maps.

[0055] The candidate boundary band construction module is used to spatially expand and extract edge bands from the coarse functional area mask, and construct a set of candidate boundary band voxels.

[0056] The boundary feature extraction module is used to extract multi-scale image features for each voxel in the candidate boundary band voxel set and construct a boundary feature vector.

[0057] The boundary probability estimation module is used to input the boundary feature vector into the pre-trained boundary discrimination network and generate individualized boundary probability maps of brain functional areas based on the output.

[0058] The boundary surface reconstruction module is used to construct the boundary energy function based on the boundary probability map, and combined with the energy minimization segmentation mechanism, obtain the inner region voxel set, the outer region voxel set and the boundary voxel set, and reconstruct the boundary surface of the individualized brain functional area.

[0059] The results output module is used to output the individualized brain functional area boundary surface as the individualized brain functional area boundary localization result.

[0060] The technical effects and advantages of the individualized brain functional area boundary localization method of this invention are as follows:

[0061] By integrating multimodal data fusion and intelligent boundary discrimination, this study systematically and for the first time unifies structural features such as grayscale gradient, cortical thickness, and curvature with features such as functional activation, temporal correlation, and white matter connectivity. A specially designed boundary discrimination network is used to directly learn the feature patterns of boundary voxels through supervised learning, achieving true individualized functional boundary recognition. Based on the energy minimization segmentation mechanism of the boundary probability map, sub-voxel-level localization accuracy can be achieved, far exceeding that of traditional threshold segmentation methods. Functional activation information is projected onto the individual cortical surface to ensure a high degree of conformity between functional boundaries and anatomical structures.

[0062] By combining features from three dimensions—structure, function, and connectivity—sufficient evidence is provided for boundary discrimination; imaging artifacts and non-brain tissue interference are effectively eliminated through skull boundary distance screening and morphological processing; and energy optimization based on graph cut ensures the topological coherence of the boundary in three-dimensional space.

[0063] The provided precise three-dimensional boundary surface can be directly integrated into the neuronavigation system, providing key guidance for surgeries such as tumor resection and epilepsy focus resection, and providing individualized target localization for brain stimulation therapies such as DBS and TMS, thereby improving treatment efficacy and reducing the risk of side effects. Attached Figure Description

[0064] Figure 1 This is a schematic diagram of the individualized brain functional area boundary localization method of the present invention;

[0065] Figure 2 This is a schematic diagram of the individualized brain functional area boundary localization system in this invention. Detailed Implementation

[0066] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0067] Please see Figure 1 and Figure 2 In this embodiment of the invention, a method for locating the boundaries of individualized brain functional areas includes:

[0068] Step S1: Collect multimodal brain imaging data of the target individual, and register and preprocess them to obtain structural MRI data, functional MRI data and diffusion tensor imaging data in the same voxel space;

[0069] Step S2: Based on functional MRI data, generate a functional activation map under resting state conditions, and determine a coarse functional area mask within the target anatomical region according to the functional activation map; wherein, generating a functional activation map based on functional MRI data under resting state conditions is a mature existing technology that can be implemented using existing conventional techniques such as seed point correlation analysis or independent component analysis.

[0070] Step S3: Spatial expansion and edge band extraction are performed on the coarse functional region mask to delineate a set of candidate boundary band voxels around the coarse functional region mask;

[0071] Step S4: For each voxel in the candidate boundary zone voxel set, extract multi-scale image features from structural MRI data, functional MRI data, and diffusion tensor imaging data to construct the boundary feature vector corresponding to each voxel.

[0072] Step S5: Input the boundary feature vector into the pre-trained boundary discrimination network to obtain the boundary probability value of each voxel and generate a boundary probability map of individualized brain functional areas.

[0073] Step S6: Construct a boundary energy function based on the boundary probability map, and combine it with the energy minimization segmentation mechanism to divide the candidate boundary zone voxel set to obtain the inner region voxel set, the outer region voxel set, and the boundary voxel set; and based on this, perform three-dimensional reconstruction in the individual structural space to obtain the individualized brain functional area boundary surface, and output this as the final individualized brain functional area boundary localization result.

[0074] Methods for acquiring multimodal brain imaging data of target individuals and performing registration and preprocessing include:

[0075] For a single target individual, a magnetic resonance imaging (MRI) device is used to scan and acquire multimodal brain imaging data. Specifically, the scanning bed is placed in the center of the magnet so that the target individual's head is in the area with the most uniform magnetic field. The operator needs to select and execute the preset scanning sequence on the control panel to acquire structural MRI data, functional MRI data, and diffusion tensor imaging data.

[0076] Among them, structural MRI data are three-dimensional volume data, with a voxel size of 1 mm³, used to clearly present the anatomical structures of brain tissue such as gray matter, white matter, and cerebrospinal fluid.

[0077] Functional MRI data is a series of 3D brain snapshots acquired continuously over time. It is a four-dimensional dataset, which adds a time dimension to the three-dimensional space. Its voxel size is usually 2-4 mm³. It is used to record the time series of blood oxygenation level dependent signals of the whole brain in a task state or resting state.

[0078] Diffusion tensor imaging data is a three-dimensional vector field map with directional information, used to fit the diffusion tensor of water molecules within each voxel, which is typically 2-2.5 mm³ in size.

[0079] Structural MRI data underwent cranial removal, noise suppression, and intensity non-uniformity correction to obtain preprocessed structural MRI data. Specifically, cranial removal used the FSL tool's BET function to generate an intracranial surface mask defining the boundary between the brain and skull. Combined with prior knowledge from a standard brain template, non-brain tissue voxels were removed from the original image, resulting in an intracranial brain image containing only gray matter, white matter, and cerebrospinal fluid. Noise suppression smoothed random noise in the image by calculating a weighted average of similar neighboring pixels while preserving tissue boundaries. Intensity non-uniformity correction modeled the intensity non-uniformity in the image as a bias field and iteratively optimized to remove it from the original image, ensuring grayscale consistency across the entire brain. The final preprocessed structural MRI data was obtained.

[0080] Functional MRI data undergoes temporal correction, head movement correction, and temporal drift removal, and is then registered to the voxel space of the preprocessed structural MRI data. In detail, because functional MRI data is scanned layer by layer, images from different layers at the same time point are actually acquired at different times within milliseconds. Temporal correction uses interpolation algorithms to align the time series of all layers to the same reference time point (usually an intermediate layer), eliminating the confusion caused by this acquisition time difference in functional activation analysis. Head movement correction involves registering the three-dimensional volumetric image at each time point with the first time point or a reference image, correcting for minute head movements generated by the individual during scanning, ensuring spatial alignment of images at all time points. Functional MRI data often contains slow, non-physiological baseline drift. Temporal drift removal removes low-frequency drift signals as covariates from the original time series to highlight task-relevant neural activity signals.

[0081] The processed functional MRI data is registered with the preprocessed structural MRI data to maximize their mutual information, and the spatial transformation matrix is ​​obtained. This matrix is ​​then applied to all time-point data to accurately map them to voxel space.

[0082] Motion eddy current correction and tensor fitting are performed on diffusion tensor imaging data, registering them to the voxel space of the preprocessed structural MRI data. In detail, motion eddy current correction involves registering each diffusion tensor imaging data point to a reference image while estimating and correcting image geometric distortion and signal distortion caused by subject head movement and magnetic field eddy current effects. Tensor fitting involves solving the diffusion tensor within each voxel using the least squares method on the corrected data, and calculating derived scalar maps such as the anisotropy fraction map and the principal diffusion direction map.

[0083] Subsequently, it was registered with the preprocessed structural MRI data to obtain the spatial transformation matrix, and this matrix was applied to all derived scalar maps to accurately align them to voxel space.

[0084] Intensity normalization (calculating the mean and standard deviation of the intensity of all voxels in each modality image, then subtracting the mean from the intensity value of each voxel and dividing by the standard deviation) and spatial resampling (resampling all modal data to the same voxel size and 3D mesh space as the structural MRI data) are performed on the registered structural MRI data, functional MRI data, and diffusion tensor imaging data in the same voxel space, forming a multimodal data cube with unified space, intensity, and coordinate system, giving it the same voxel resolution and coordinate system.

[0085] In summary, voxel space refers to the image coordinate system organized as a discrete three-dimensional mesh after reconstruction. Each voxel has its unique (x, y, z) coordinates, used to determine its position in the three-dimensional mesh. Each voxel also stores data information from structural MRI, functional MRI, and diffusion tensor imaging. A voxel is a small three-dimensional cube and the basic unit of voxel space.

[0086] Methods for determining coarse functional area masks within a target anatomical region based on functional activation maps include:

[0087] In the functional activation map, an adaptive threshold segmentation is applied to the statistical parameter map (an image generated after statistical calculation, in which each voxel contains a statistical value representing the signal change in the target anatomical region during task execution compared to when the brain is not performing a specific task) corresponding to the target functional paradigm (the task performed by the target individual during the scan). This yields the initial activation region. Specifically, based on the time series of functional MRI data, the time signal of each voxel in the whole brain is linearly fitted to the reference waveform of the target functional paradigm to calculate a statistical value Qt representing the activation intensity. The Qt values ​​of all voxels are counted to generate a statistical parameter map. Then, an intensity threshold is determined, and all voxels in the statistical parameter map that are above the intensity threshold are marked as 1 (indicating activation), and those below are marked as 0 (indicating non-activation). This yields the initial activation region (in the form of a binary map).

[0088] Based on preprocessed structural MRI data, cortical surface reconstruction (aiming to automatically identify and extract the inner and outer surfaces of the brain's gray matter cortex from structural MRI data and model them as a smooth, continuous two-dimensional curved surface mesh, such as FreeSurfer's recon-all workflow) is used to extract the cortical surface of the target anatomical region and project the initial activation region onto the cortical surface. Specifically, the FreeSurfer recon-all workflow is run to automatically process the structural MRI data of the target individual, extracting a precise two-dimensional cortical surface mesh (in the form of a curved surface) of the target anatomical region. A spatial correspondence between voxel space and cortical surface vertices is pre-established. Then, the binary information (1 or 0) in the initial activation region is precisely assigned from voxel space to the corresponding cortical surface vertices. After this processing, the information carried by functional activation is transformed into information spread across the cortical surface. This allows for intuitive observation of functional activities on the anatomical structure and makes subsequent processing more consistent with the brain's topology.

[0089] Morphological closing operations are performed on the projected initial activation regions (to smooth the boundaries of the initial activation regions and connect neighboring small regions that originally belonged to the same functional area but were broken due to noise, without changing the area and shape of the original region) and connected region filtering (this step is because after morphological closing operations, there may still be multiple mutually separated activation regions in the image. Therefore, the purpose of this process is to retain only the largest one or two connected regions, and to treat the remaining small, discrete patches as noise and remove them). Discrete patches are removed to obtain a functionally activated main region with a coherent topological structure. Specifically, morphological closing operations are performed on the initial activation regions on the grid on the cortical surface. Then, all interconnected clumps are found, sorted by clump size, and the largest clump (or the top N clumps, depending on the specific situation) is retained. The remaining small clumps are deleted, and finally a complete functionally activated main region is obtained.

[0090] The functional activation region is transformed from the cortical surface space to the voxel space to form a binary mask, which serves as the coarse functional region mask. Specifically, the functional activation region existing on the cortical surface is mapped back to the voxel space. In the voxel space, all voxels belonging to this region are marked as 1, and the rest are marked as 0, generating a binary mask file, which is the coarse functional region mask.

[0091] Methods for defining a candidate boundary band set of voxels around a coarse functional region mask include:

[0092] Dilation and erosion with a preset radius are performed on the coarse functional region mask to obtain the expanded functional region mask and the contracted functional region mask, respectively. Since the exact location of the candidate boundary cannot be known in advance, the dilation and erosion operations construct a candidate region that completely covers the real functional boundary.

[0093] Specifically, a preset radius (e.g., 2-3 voxels) is used, which determines the size of the neighborhood (the set of other voxels surrounding a given voxel). An expansion operation is then performed on the binarized coarse functional region mask; that is, all voxels with a value of 1 in the mask are extended to their neighborhoods to obtain an expanded functional region mask. This process ensures that the candidate boundary bands in subsequent analysis completely cover the true functional region boundaries.

[0094] Simultaneously, an erosion operation is performed on it using the same radius (removing voxels with a value of 1 on the boundary of the expanded functional region mask, causing it to shrink inward), to obtain the shrunken functional region mask;

[0095] Calculate the difference between the expanded functional region mask and the contracted functional region mask (the purpose is to obtain a shell that only contains the boundary region, difference = expanded functional region mask - contracted functional region mask), and use the voxels in the difference set as the candidate boundary zone voxel set; the voxels contained in this set constitute a three-dimensional annular region surrounding the coarse functional region.

[0096] In the candidate boundary zone voxel set, based on the preprocessed structural MRI data, voxels whose distance to the skull boundary is less than a preset safety distance are removed, and this is recorded as the final candidate boundary zone voxel set. The purpose is to eliminate imaging artifacts or non-brain tissue interference that may be introduced due to proximity to extracranial structures such as the skull. Specifically, based on the obtained intracranial surface mask, the shortest Euclidean distance from each voxel in the candidate boundary zone set to the intracranial surface mask is calculated. A preset safety distance (e.g., 3 mm) is set, and each voxel in the candidate boundary zone set is traversed to obtain the final candidate boundary zone voxel set.

[0097] The aforementioned design process primarily addresses the problem that template maps in the background techniques cannot accurately reflect the location and morphological differences of individual brain functional areas. They only provide an average functional distribution and lack detailed depiction of transitional zones at the boundaries of functional areas, leading to blurred boundaries. The core of this design is based on the target individual's own MRI data. All operations are performed within the individual's space, and the results reflect only the unique anatomical and functional characteristics of the individual, fundamentally eliminating the systematic bias introduced by group templates. Furthermore, it transforms a raw, noisy statistical map into a core functional area with accurate anatomical location, laying a solid and reliable foundation for subsequent individualized, high-precision boundary localization.

[0098] The goal of constructing a boundary feature vector for each voxel is to extract a set of multi-dimensional features from each voxel in the candidate boundary zone voxel set, comprehensively reflecting its anatomical structure, functional activity, and fiber connectivity characteristics, and combine them into a high-dimensional boundary feature vector to provide an information foundation for subsequent intelligent boundary discrimination. The methods include:

[0099] In the preprocessed structural MRI data, the gray-level gradient magnitude, cortical thickness estimate, and local curvature features of the neighborhood of each candidate boundary zone voxel in the candidate boundary zone voxel set are extracted. The gray-level gradient magnitude measures the intensity change of the image at the voxel's location. It is extracted by pre-defining a neighborhood centered on the current candidate boundary zone voxel, using gradient operators to calculate its gradient in various directions of the coordinate system, and obtaining its comprehensive gradient magnitude. The cortical thickness estimate quantifies the gray matter thickness of the cortex where the voxel is located. It is extracted by calculating the shortest distance between the white matter surface (the interface between gray matter and white matter) and the pia mater surface (the medial interface of cerebrospinal fluid), and this distance is mapped back to voxel space to obtain the voxel thickness estimate. The local curvature features describe the geometric morphology of the cortical surface (such as gyral bulges or sulcal depressions). This feature is extracted because functional boundaries are often related to specific anatomical structural transitions. It is extracted by calculating the Gaussian curvature or average curvature of each vertex on the cortical surface and then mapping these curvature values ​​back to voxel space. The above features are used to represent the microscopic anatomical structure of brain tissue at the location of voxels.

[0100] From functional MRI data registered to voxel space, the activation intensity and spatial gradient magnitude of activation intensity of each candidate boundary zone voxel, as well as the average correlation coefficient with the time series of neighboring voxels, are calculated. Activation intensity reflects the strength of neural activity of the voxel during task execution, i.e., the Qt value. The spatial gradient magnitude of activation intensity measures the spatial rate of change of activation intensity. It is calculated by calculating the spatial gradient of the voxel in its three-dimensional neighborhood on the functional activation map and obtaining the magnitude of the gradient. The average correlation coefficient characterizes the synchronicity of the voxel with surrounding brain regions in functional connectivity. It is calculated by extracting the time series of the candidate voxel and its neighborhood from the functional MRI data, calculating the Pearson correlation coefficient between the voxel and all voxels in the neighborhood, and recording the average of the coefficients as the average correlation coefficient.

[0101] From diffusion tensor imaging data registered to voxel space, the anisotropy fraction, principal fiber orientation, and fiber-tracing-based structural connectivity strength features of each candidate boundary zone voxel are calculated. The anisotropy fraction reflects the fiber integrity and density of white matter, and the value corresponding to the voxel is directly read from the anisotropy fraction map. The principal fiber orientation describes the fiber direction of the most abundant white matter within the voxel. It is calculated by performing eigenvalue decomposition on the diffusion tensor within each voxel to obtain its principal feature vector, which is a three-dimensional unit vector whose direction is the principal fiber orientation of the voxel. The fiber-tracing-based structural connectivity strength features quantify the tightness of white matter fiber connections between the voxel and a specific region of interest (such as the core region within a coarse functional area mask). It is calculated by using the coarse functional area mask as a seed point, performing deterministic fiber tracing throughout the brain, and calculating how many reconstructed fiber bundles pass through the candidate voxel. The number of fiber bundles passing through the voxel is used as the structural connectivity strength feature.

[0102] The extracted multi-scale image features are then arranged in a preset order to form a boundary feature vector.

[0103] To address the problem of background technologies relying solely on functional activation intensity and resulting in vague and singular boundary determination, this invention systematically integrates features from three dimensions: anatomical structure (grayscale gradient, cortical thickness, curvature), functional activity (activation intensity, gradient, correlation), and white matter connectivity (anisotropy fraction, fiber orientation, connectivity strength). This provides each candidate voxel with a three-dimensional "identity profile," enabling boundary determination to move beyond threshold cutting based on a single indicator and instead rely on comprehensive intelligent decision-making based on multiple pieces of evidence. This allows for clear identification of the true functional boundary transition zone, ensuring that boundary determination is entirely based on the unique brain structure and functional connectivity patterns of each individual, accurately capturing individual differences.

[0104] The method in step S5 includes:

[0105] A boundary discrimination network is pre-built based on a deep learning network. The input is defined as a boundary feature vector, and the output is a ternary probability distribution, including the probability value of the inner region, the probability value of the boundary, and the probability value of the outer region.

[0106] The boundary discrimination network is pre-trained using training samples;

[0107] For the target individual, the boundary feature vector corresponding to each voxel in the candidate boundary band voxel set is input into the boundary discrimination network, and the corresponding ternary probability distribution is output.

[0108] Extract the boundary probability values ​​and rearrange them according to their voxel spatial positions to form a three-dimensional boundary probability map.

[0109] Methods for pre-training boundary discrimination networks using training samples include:

[0110] Based on fine manual annotation, the voxels containing the true boundaries of brain functional areas are marked as boundaries, voxels extending a certain distance towards the functional areas (the brain regions that respond most strongly, are most stable, and have the most concentrated range for specific tasks) are marked as inner regions, and voxels extending a certain distance towards non-functional areas are marked as outer regions, forming training samples with three types of labels (e.g., 0 represents inner region, 1 represents boundary, and 2 represents outer region). During training, this label is often converted into one-hot encoding, such as inner region = [1,0,0], boundary = [0,1,0], and outer region = [0,0,1].

[0111] During the training phase, the network is trained using the boundary feature vector as input and the supervision signal constructed from the corresponding three-class labels as the target; at the same time, the boundary class is given higher weight in the loss function.

[0112] After the training process converges, the final boundary discrimination network model parameters are saved. The resulting network has the ability to map the input boundary feature vectors to a ternary probability distribution.

[0113] For example, a five-layer fully connected neural network is constructed as a boundary discrimination network, and the specific network design is as follows;

[0114] The number of neurons in the input layer is set to be the same as the dimension of the boundary feature vector, which is used to receive all the multimodal features of a single voxel; there are three hidden layers with the number of neurons being 512, 256 and 128 respectively. Each layer is followed by a ReLU activation function to introduce a non-linear transformation, and a Dropout mechanism (dropout rate set to 0.3) is applied to prevent overfitting.

[0115] The output layer consists of three neurons, corresponding to the three categories of inner region, boundary and outer region. The softmax activation function is used to normalize the output value of the neurons into a ternary probability distribution [P_in, P_bou, P_out], where P_in, P_bou and P_out represent the probability value of inner region, boundary probability value and outer region probability value respectively, and the sum of the three is 1.

[0116] Construction of the training sample set;

[0117] Collect multimodal brain imaging data from DR (e.g., DR=150) training individuals.

[0118] The surface of the true boundary of the brain functional area is manually marked by at least two experienced neuroradiologists or neurosurgeons. All voxels within d millimeters (e.g., d=1mm) of the true boundary surface are marked as the boundary; voxels extending L1 millimeters (e.g., L1=3mm) from the boundary into the functional area (i.e., the functional core area) are marked as the inner region; voxels extending L2 millimeters (e.g., L2=3mm) from the boundary into the non-functional area are marked as the outer region. The category labels are converted into one-hot codes: inner region=[1,0,0], boundary=[0,1,0], outer region=[0,0,1].

[0119] The Adam optimizer (with an initial learning rate of 1e-4) and a weighted cross-entropy loss function were used, with weights set to inner region = 1.0, boundary = 3.0, and outer region = 1.0 to improve the sensitivity to boundary recognition. The batch size was set to 64, and the training epochs were set to 200. After each epoch, the performance was evaluated on an independent validation set, and training was terminated early when the loss function no longer decreased significantly.

[0120] The boundary feature vector corresponding to each voxel in the candidate boundary band voxel set of the target individual is sequentially input into the pre-trained deep learning MLP, and a ternary probability distribution is output for each voxel. From the ternary probability distribution of each voxel, the second probability value, namely the boundary probability value (P_bou), is extracted separately. Based on the original coordinates of each voxel in the three-dimensional voxel space, all the extracted P_bou values ​​are rearranged into a three-dimensional matrix with the same size as the original brain image. This three-dimensional matrix is ​​the final output boundary probability map. The value of each voxel in the map (between 0 and 1) represents the probability that the position belongs to the boundary of the functional area.

[0121] Methods for obtaining the inner voxel set, outer voxel set, and boundary voxel set include:

[0122] A voxel graph is constructed using the set of voxels in the candidate boundary zone as graph nodes and neighborhood relationships as adjacent edges.

[0123] Based on the boundary probability map, the boundary probability value of each voxel is converted into a data term cost belonging to the boundary category; based on the preprocessed structural MRI data, the gray level difference and spatial distance between adjacent voxels are calculated, which together constitute the smoothing term cost; in this way, a boundary energy function defined on the three labels of inner region, boundary, and outer region is constructed.

[0124] A graph cut-α extension is used to solve the global energy minimization problem of the boundary energy function, and an optimal class label is assigned to each node in the graph.

[0125] The set of voxels with the optimal label of the inner region is defined as the inner region voxel set, the set of voxels with the optimal label of the outer region is defined as the outer region voxel set, and the set of voxels with the optimal label of the boundary is defined as the boundary voxel set.

[0126] For example: Define each voxel in the candidate boundary band voxel set as a graph node, and establish connection edges between each node and all neighboring nodes in the three-dimensional twenty-six neighborhood;

[0127] The boundary energy function E(L) is composed of a weighted sum of the data term cost E_data and the smoothing term cost E_smooth (the weights are customized according to the actual generation requirements). Its goal is to find the label configuration L that minimizes E(L).

[0128] For any voxel v, the data term cost E_data(v=boundary)=1-P_bou(v) is used to measure the penalty for assigning a boundary label to a voxel, and the smoothing term cost E_smooth(v_i,v_j)=λ×exp(-(||I(v_i)-I(v_j)||²÷2σ²))÷dist(v_i,v_j) is used to penalize the case where adjacent voxels are assigned different labels, thereby ensuring the regional smoothness of the segmentation result.

[0129] Where P_bou(v) represents the boundary probability value of voxel v, v_i and v_j represent a pair of adjacent voxels, λ is a preset weighting coefficient (e.g., λ=3.0), σ is a normalization parameter, usually taken as the standard deviation of the voxel gray intensity within the entire candidate boundary band, ||I(v_i)-I(v_j)||² represents the square of the gray difference between adjacent voxel pairs, and I(v_i) and I(v_j) represent the gray intensity of a pair of adjacent voxels.

[0130] Graph cut-α expansion involves multiple iterations, each targeting a specific boundary label. By solving a binary graph cut problem, it determines which nodes should be expanded to that label, thus systematically reducing global energy. After several iterations (usually equal to the number of categories, i.e., 3), the algorithm converges, assigning an optimal category label (inner region, boundary, or outer region) to each node in the graph (i.e., each voxel in the candidate boundary band).

[0131] Traverse all graph nodes and, based on their optimal assigned labels, group them into three separate sets:

[0132] Inner region voxel set: All voxels labeled as inner region;

[0133] Outer region voxel set: All voxels labeled as outer region;

[0134] Boundary voxel set: All voxels labeled as boundary;

[0135] This set of boundary voxels represents the final, precisely defined, and individualized functional area boundaries, which will serve as the direct input for the next step of 3D surface reconstruction.

[0136] The above design process combines the functional boundary possibilities (data terms) provided by the boundary probability map with the anatomical continuity constraints (smoothing terms) provided by structural MRI. Through global energy optimization, it realizes the transformation from fuzzy probability to precise and definite boundary, completely solves the boundary ambiguity problem caused by single threshold segmentation, effectively eliminates discrete and isolated misclassification points, and generates a topologically coherent, spatially consistent and smooth three-dimensional boundary voxel set, laying a solid foundation for subsequent reconstruction of high-quality surfaces.

[0137] Methods for obtaining individualized brain functional area boundary surfaces include:

[0138] In voxel space, the set of boundary voxels is taken as the set of spatial target points;

[0139] The moving cube algorithm is used to extract the initial boundary surface from the set of spatial target points;

[0140] The initial boundary surface is smoothed using Laplace to obtain the individualized brain functional area boundary surface.

[0141] The individualized brain functional area boundary surfaces are mapped to the surgical navigation coordinate system or the standard medical imaging coordinate system for three-dimensional visualization and data export.

[0142] The boundary voxel set is a series of discrete voxel coordinates in voxel space that are marked as "boundaries". This discrete voxel set is converted into a three-dimensional binary volume data, that is, all voxels belonging to the boundary voxel set are assigned a value of 1, and the rest of the voxels are assigned a value of 0.

[0143] The moving cube algorithm is a standard algorithm for extracting isosurfaces from a 3D scalar field. It iterates through each cube (composed of 8 adjacent voxels) of the binary volume data, and determines whether the isosurface (the isosurface threshold is usually set to 0.5) passes through the cube based on the values ​​(0 or 1) of the cube's 8 vertices. It also calculates the intersection points of the isosurface and the cube's edges. Specifically, the prepared binary volume data is input into the algorithm, which automatically traverses the entire volume. Within each cube traversed by the isosurface, one or more triangular patches are generated through linear interpolation. All these triangular patches are connected to form a connected, closed triangular mesh surface, which serves as the initial boundary surface.

[0144] Laplace smoothing is a commonly used mesh smoothing algorithm. Its core idea is to move each vertex on the surface towards the geometric center of its adjacent vertices, gradually smoothing out local undulations through iteration. This is because although the initial surface is close to the final boundary in shape, its surface often has stepped jagged edges due to voxel discretization, making it not smooth enough. Therefore, the specific implementation process is as follows: for each vertex v_i on the initial boundary surface, its new position v_i' is determined by the positions of all its adjacent vertices P(v_i).

[0145] The calculation formula is: v_j is a vertex among all adjacent vertices P(v_i), where u is a relaxation factor (usually set to 0.5) used to control the step size of each movement; P is the number of adjacent vertices; this smoothing process is iterated multiple times (e.g., 10-20 times) until the surface shape is stable and meets the smoothness requirements.

[0146] Example 2, please refer to Figure 2 As shown, for parts not described in detail in this embodiment, please refer to the description in Embodiment 1. A personalized brain functional area boundary localization system is provided, including:

[0147] The multimodal brain image preprocessing module is used to acquire multimodal brain image data of the target individual and perform registration and preprocessing.

[0148] The coarse functional area generation module is used to calculate functional activation maps based on functional MRI data and generate coarse functional area masks within the target anatomical region according to the functional activation maps.

[0149] The candidate boundary band construction module is used to spatially expand and extract edge bands from the coarse functional area mask, and construct a set of candidate boundary band voxels.

[0150] The boundary feature extraction module is used to extract multi-scale image features for each voxel in the candidate boundary band voxel set and construct a boundary feature vector.

[0151] The boundary probability estimation module is used to input the boundary feature vector into the pre-trained boundary discrimination network and generate individualized boundary probability maps of brain functional areas based on the output.

[0152] The boundary surface reconstruction module is used to construct the boundary energy function based on the boundary probability map, and combined with the energy minimization segmentation mechanism, obtain the inner region voxel set, the outer region voxel set and the boundary voxel set, and reconstruct the boundary surface of the individualized brain functional area.

[0153] The results output module is used to output the individualized brain functional area boundary surface as the individualized brain functional area boundary localization result.

[0154] The above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art can still modify the technical solutions described in the foregoing embodiments or make equivalent substitutions for some of the technical features. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

[0155] It should be noted that all formulas in this manual are calculated by removing dimensions and taking their numerical values. The formulas are derived from software simulations based on a large amount of collected data to obtain the most recent real-world results. The preset parameters and thresholds in the formulas are set by those skilled in the art according to the actual situation.

[0156] Although embodiments of the invention have been shown and described, those skilled in the art will understand that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the claims and their equivalents.

Claims

1. A method for locating the boundaries of individualized brain functional areas, characterized in that, include: Step S1: Collect multimodal brain imaging data of the target individual, and register and preprocess them to obtain structural MRI data, functional MRI data and diffusion tensor imaging data in the same voxel space; Step S2: Based on functional MRI data, generate a functional activation map under resting conditions, and determine a coarse functional area mask within the target anatomical region according to the functional activation map; Step S3: Spatial expansion and edge band extraction are performed on the coarse functional region mask to delineate a set of candidate boundary band voxels around the coarse functional region mask; Step S4: For each voxel in the candidate boundary zone voxel set, extract multi-scale image features from structural MRI data, functional MRI data, and diffusion tensor imaging data to construct the boundary feature vector corresponding to each voxel. Step S5: Input the boundary feature vector into the pre-trained boundary discrimination network to obtain the boundary probability value of each voxel and generate a boundary probability map of individualized brain functional areas. Step S6: Construct a boundary energy function based on the boundary probability map, and combine it with the energy minimization segmentation mechanism to divide the candidate boundary zone voxel set to obtain the inner region voxel set, the outer region voxel set, and the boundary voxel set; and based on this, perform three-dimensional reconstruction in the individual structural space to obtain the individualized brain functional area boundary surface, and output this as the final individualized brain functional area boundary localization result.

2. The personalized brain functional area boundary localization method according to claim 1, characterized in that, The method for collecting multimodal brain imaging data of the target individual and registering and preprocessing it includes: For a single target individual, a magnetic resonance imaging (MRI) device is used to scan and acquire multimodal brain imaging data; the multimodal brain imaging data includes structural MRI data, functional MRI data, and diffusion tensor imaging data; Skull removal, noise suppression, and intensity non-uniformity correction were performed on the structural MRI data to obtain preprocessed structural MRI data. The functional MRI data were subjected to time-layer correction, head motion correction and time drift removal, and then registered to the voxel space of the preprocessed structural MRI data. Motion eddy current correction and tensor fitting were performed on the diffusion tensor imaging data, and the data were registered to the voxel space of the preprocessed structural MRI data. Intensity normalization and spatial resampling were performed on the registered structural MRI data, functional MRI data, and diffusion tensor imaging data to ensure they had the same voxel resolution and coordinate system.

3. The individualized brain functional area boundary localization method according to claim 2, characterized in that, The method for determining a coarse functional area mask within a target anatomical region based on a functional activation map includes: In the functional activation graph, an adaptive threshold segmentation is applied to the statistical parameter graph corresponding to the target functional paradigm to obtain the initial activation region. Based on preprocessed structural MRI data, the cortical surface of the target anatomical region is extracted through cortical surface reconstruction, and the initial activated region is projected onto the cortical surface; Morphological closing operations and connected region filtering are performed on the initial activated region after projection to remove discrete patches and obtain the functionally activated main region with a coherent topological structure. The functional activation subject area is transformed from the cortical surface space to the voxel space to form a binary mask, which serves as a coarse functional area mask.

4. The personalized brain functional area boundary localization method according to claim 3, characterized in that, The method for defining a set of candidate boundary band voxels around a coarse functional region mask includes: Perform expansion and erosion with a preset radius on the coarse functional area mask to obtain an expanded functional area mask and a contracted functional area mask, respectively; Calculate the difference between the expanded functional region mask and the contracted functional region mask, and use the voxels in the difference set as the candidate boundary zone voxel set; In the candidate boundary zone voxel set, based on the preprocessed structural MRI data, voxels that are less than a preset safe distance from the skull boundary are removed, and these are recorded as the final candidate boundary zone voxel set.

5. The personalized brain functional area boundary localization method according to claim 4, characterized in that, The method for constructing the boundary feature vector corresponding to each voxel includes: In the preprocessed structural MRI data, the gray-level gradient magnitude, cortical thickness estimate, and local curvature features of the neighborhood of each candidate boundary zone voxel in the candidate boundary zone voxel set are extracted. From functional MRI data registered to voxel space, calculate the activation intensity and spatial gradient magnitude of each candidate boundary zone voxel, as well as the average correlation coefficient with the time series of neighboring voxels. From the diffusion tensor imaging data registered to voxel space, calculate the anisotropy fraction, principal fiber orientation, and fiber-tracking-based structural connectivity strength characteristics of each candidate boundary zone voxel. The extracted multi-scale image features are then concatenated in a preset order to form a boundary feature vector.

6. The personalized brain functional area boundary localization method according to claim 5, characterized in that, The method in step S5 includes: A boundary discrimination network is pre-built based on a deep learning network. The input is defined as a boundary feature vector, and the output is a ternary probability distribution, including the probability value of the inner region, the probability value of the boundary, and the probability value of the outer region. The boundary discrimination network is pre-trained using training samples; For the target individual, the boundary feature vector corresponding to each voxel in the candidate boundary band voxel set is input into the boundary discrimination network, and the corresponding ternary probability distribution is output. Extract the boundary probability values ​​and rearrange them according to their voxel spatial positions to form a three-dimensional boundary probability map.

7. The individualized brain functional area boundary localization method according to claim 6, characterized in that, The method for pre-training the boundary discrimination network using training samples includes: The voxels containing the true boundaries of brain functional areas are marked as boundaries, the voxels extending a certain distance toward the functional areas are marked as inner regions, and the voxels extending a certain distance toward non-functional areas are marked as outer regions, thus forming training samples with three types of labels. During the training phase, the network is trained using the boundary feature vector as input and the supervision signal constructed from the corresponding three-class labels as the target; at the same time, the boundary class is given higher weight in the loss function. After the training process converges, the final boundary discrimination network model parameters are saved.

8. The personalized brain functional area boundary localization method according to claim 7, characterized in that, The method for obtaining the inner region voxel set, the outer region voxel set, and the boundary voxel set includes: A voxel graph is constructed using the set of voxels in the candidate boundary zone as graph nodes and neighborhood relationships as adjacent edges. Based on the boundary probability map, the boundary probability value of each voxel is converted into a data item cost belonging to the boundary category; based on the preprocessed structural MRI data, the gray level difference and spatial distance between adjacent voxels are calculated, which together constitute the smoothing term cost; in this way, a boundary energy function defined on the three categories of labels: inner region, boundary, and outer region is constructed. A graph cut-α extension is used to solve the global energy minimization problem of the boundary energy function, and an optimal class label is assigned to each node in the graph. The set of voxels with the optimal label of the inner region is defined as the inner region voxel set, the set of voxels with the optimal label of the outer region is defined as the outer region voxel set, and the set of voxels with the optimal label of the boundary is defined as the boundary voxel set.

9. The personalized brain functional area boundary localization method according to claim 8, characterized in that, The method for obtaining the individualized brain functional area boundary surface includes: In voxel space, the set of boundary voxels is taken as the set of spatial target points; The moving cube algorithm is used to extract the initial boundary surface from the set of spatial target points; The initial boundary surface is smoothed using Laplace to obtain the individualized brain functional area boundary surface. The individualized brain functional area boundary surfaces are mapped to the surgical navigation coordinate system or the standard medical imaging coordinate system for three-dimensional visualization and data export.

10. A personalized brain functional area boundary localization system, used to implement the personalized brain functional area boundary localization method according to any one of claims 1 to 9, characterized in that, include: The multimodal brain image preprocessing module is used to acquire multimodal brain image data of the target individual and perform registration and preprocessing. The coarse functional area generation module is used to calculate functional activation maps based on functional MRI data and generate coarse functional area masks within the target anatomical region according to the functional activation maps. The candidate boundary band construction module is used to spatially expand and extract edge bands from the coarse functional area mask, and construct a set of candidate boundary band voxels. The boundary feature extraction module is used to extract multi-scale image features for each voxel in the candidate boundary band voxel set and construct a boundary feature vector. The boundary probability estimation module is used to input the boundary feature vector into the pre-trained boundary discrimination network and generate individualized boundary probability maps of brain functional areas based on the output. The boundary surface reconstruction module is used to construct the boundary energy function based on the boundary probability map, and combined with the energy minimization segmentation mechanism, obtain the inner region voxel set, the outer region voxel set and the boundary voxel set, and reconstruct the boundary surface of the individualized brain functional area. The results output module is used to output the individualized brain functional area boundary surface as the individualized brain functional area boundary localization result.

Images