A neurosurgical operation path intelligent planning method and system

Abstract

The application discloses a kind of neurosurgery operation path intelligent planning method and system, it is related to neurosurgical navigation technical field.The method comprises: obtaining preoperative multi-modal image, segmenting brain function area and constructing individualized brain function network, calculating the comprehensive network topological weight of each function area;Based on the weight, function area sensitivity and anisotropic distance along fiber bundle running, construct continuous function potential field;The function potential field is coupled with physical obstacle field, and the optimal operation path is planned using the path search algorithm guided by potential field;Intraoperative displacement data is updated to update the potential field and carry out local re-planning;The potential field and path are visualized and displayed through augmented reality equipment.The application converts the topological properties of brain function network into continuous potential field to drive path planning, which helps to actively bypass high network value area while avoiding physical obstacles, and improves the protection of functional connectivity integrity in operation path planning.

Description

Technical Field

[0001] This invention relates to the field of neurosurgical navigation technology, and more specifically, to a method and system for intelligent planning of neurosurgical surgical pathways. Background Technology

[0002] Neurosurgical pathway planning needs to avoid physical obstacles such as important blood vessels and ventricles while protecting functional brain areas and nerve fiber bundles as much as possible to reduce the risk of postoperative functional impairment. Existing pathway planning methods are mainly based on three-dimensional reconstruction and geometric analysis of preoperative images, searching for feasible paths from the entry point to the lesion by setting safe distance constraints. These methods treat functional brain areas as discrete, rigid obstacles, relying primarily on geometric obstacle avoidance strategies. They fail to fully consider the complex connections between functional areas and their importance in the whole-brain functional network, and have limited ability to protect pivotal areas connecting multiple functional areas.

[0003] Chinese patent CN114948199A discloses a surgical assistance system and a surgical path planning method. This system reconstructs a 3D model using an image processing module and employs a spatial trajectory planning algorithm to find the optimal surgical path within the surgical constraint area. This method focuses on path search under static geometric constraints, primarily protecting brain functional areas by setting a safety distance, and does not incorporate the network topology attributes of the functional areas into the cost function of path planning.

[0004] Chinese patent CN121242727A discloses a method and system for precise planning of neurosurgical brain puncture paths. This method calculates functional interaction risks by extracting the volume of brain functional areas and the direction of nerve fiber bundles, and constructs a risk gradient map to score and screen candidate paths. While this method introduces more dimensions of assessment in terms of functional protection, the path generation method mainly relies on the scoring and screening of a pre-set set of candidate paths. The initiative of path planning and its adaptability to dynamic changes such as intraoperative brain drift need further improvement. Therefore, a method and system for intelligent planning of neurosurgical surgical pathways is proposed to address the above problems. Summary of the Invention

[0005] To overcome the aforementioned deficiencies of the prior art, embodiments of the present invention provide an intelligent planning method and system for neurosurgical surgical pathways, aiming to address the problems in existing neurosurgical surgical pathway planning methods, such as the reliance on geometric obstacle avoidance for the protection of brain functional areas, the failure to integrate the topological attributes of the brain functional network into the pathway planning process, resulting in insufficient protection of network hub areas, and the need to improve pathway adaptability during intraoperative brain drift.

[0006] To achieve the above objectives, the present invention provides an intelligent planning method for neurosurgical surgical pathways, comprising the following steps: S1. Acquire the patient's preoperative multimodal image data, segment brain functional areas based on the preoperative multimodal image data and obtain a binary three-dimensional spatial mask for each brain functional area, construct an individualized brain functional network, and calculate the comprehensive network topology weight of each brain functional area based on the brain functional network. The comprehensive network topology weight is used to quantify the hub status of the functional area in the whole brain functional network. S2. Determine the surgical area based on the target lesion point and the entry point into the skull. Based on the comprehensive network topology weights, the sensitivity of each brain functional area, and the anisotropic distance from any point in space to the boundary of each functional area, construct a continuous functional potential field covering the surgical area. The anisotropic distance is the geodesic distance from the point to the boundary of the functional area along the direction of the white matter fiber bundle. The continuous functional potential field reflects the range of influence of the functional area propagating along the fiber bundle and its network importance. S3. Based on the continuous functional potential field and the physical obstacle field, a coupled field is constructed. In the coupled field, a path search algorithm is used to plan the optimal surgical path from the entry point to the target lesion point. The coupled field integrates the constraints of functional protection and physical obstacle avoidance.

[0007] Furthermore, the construction of the individualized brain functional network includes: Based on diffusion tensor imaging data, fiber connections between brain functional areas are extracted, and a structural connectivity matrix is ​​constructed. Based on functional magnetic resonance imaging data, time-series correlations of brain functional regions were extracted, and a functional connectivity matrix was constructed. By fusing the structural connectivity matrix and the functional connectivity matrix, a weighted adjacency matrix is ​​obtained, which characterizes the overall structural and functional connectivity strength between brain functional areas.

[0008] Furthermore, the comprehensive network topology weights are obtained by weighted fusion of multiple graph theory indicators for each brain functional region node. These multiple graph theory indicators include weighted degree centrality, betweenness centrality, eigenvector centrality, and intra-module connectivity. These indicators reflect the topological importance of nodes in the whole-brain network from different dimensions.

[0009] Furthermore, the anisotropic distance is obtained by solving the equation, and the propagation speed in the equation is positively correlated with the anisotropic component value in the diffusion tensor imaging data, which makes the functional influence propagate faster and farther along the fiber bundle direction. The potential energy value at any point in the continuous functional potential field is the superposition of the anisotropic decaying potential energy generated by all brain functional areas at that point. The anisotropic decaying potential energy is calculated by using the anisotropic distance as the independent variable and the comprehensive network topology weight as the proportional coefficient through a decay function. The decay function is a monotonically decreasing function of distance, which makes the network hub region generate higher potential energy.

[0010] Furthermore, the path search algorithm is a potential field-guided fast exploration random tree algorithm, and its sampling probability is determined according to the potential energy value of the continuous functional potential field. The lower the potential energy value, the higher the sampling probability of the region, thus guiding the search to concentrate on low-risk regions. Its node expansion direction is a weighted sum of the random sampling direction and the negative gradient direction of potential energy, so that path growth takes into account both random exploration and the direction of potential energy decrease.

[0011] Furthermore, the optimal surgical path is obtained by minimizing a preset cost function, which includes a path length integral term, a physical obstacle risk integral term along the path, a functional potential energy integral term along the path, and an angle penalty term between the path direction and the potential energy gradient direction. The angle penalty term encourages the path to move along the equipotential line direction to avoid crossing the potential energy abrupt change region.

[0012] Furthermore, the method also includes an intraoperative path update step: T1, to obtain real-time intraoperative imaging data or key point displacement data collected by intraoperative sensors, in order to obtain current brain tissue deformation information; T2, based on the intraoperative real-time image data, the displacement field of the brain tissue is calculated by non-rigid image registration, or based on the key point displacement data, the displacement field of the brain tissue is solved by a biomechanical model, the displacement field describing the spatial deformation of the brain tissue from preoperative to intraoperative. T3, based on the displacement field, the preoperatively constructed continuous functional potential field is mapped to the current spatial position to obtain the intraoperative real-time functional potential field, so that the functional protection information is aligned with the current anatomical state; T4, starting from the current position of the surgical instrument, performs local replanning of the remaining path based on the real-time functional potential field during the operation and the current position of the surgical instrument, in order to compensate for the path deviation caused by brain drift.

[0013] Furthermore, the local replanning employs a dynamic fast exploratory random tree algorithm, which only prunes and re-expands the paths in the regions affected by the displacement field, reducing computational overhead and maintaining the continuity of the original paths.

[0014] Furthermore, the method also includes a visualization guidance step: the continuous functional potential field is superimposed on the patient's actual anatomical structure in the form of a semi-transparent color cloud map using an augmented reality device, and the optimal surgical path is displayed in real time, providing doctors with intuitive functional risk distribution and path guidance.

[0015] To achieve the above objectives, the present invention also provides an intelligent neurosurgical path planning system for executing the intelligent neurosurgical path planning method described in any of the above claims, the system comprising: The preoperative modeling module is used to perform operations such as acquiring preoperative multimodal image data, segmenting brain functional areas and obtaining binary three-dimensional spatial masks for each brain functional area, constructing individualized brain functional networks, calculating the topological weights of the integrated network, and constructing a continuous functional potential field. The path planning module is used to perform the operation of searching for the optimal surgical path in the coupled field formed by the continuous functional potential field and the physical obstacle field; The intraoperative update module is used to perform operations such as acquiring intraoperative displacement data, updating the continuous functional potential field according to the displacement field, and triggering local path replanning. The navigation guidance module is used to perform operations that visualize the continuous functional potential field and the optimal surgical path through an augmented reality device.

[0016] The technical effects and advantages of this invention are as follows: This invention constructs an individualized brain functional network and calculates the comprehensive network topological weights for each functional region. These weights integrate multiple graph theory metrics, such as degree centrality and betweenness centrality, which helps quantify the pivotal role of each functional region in the whole-brain functional network. Introducing these weights into the functional potential field construction enables network hub regions to generate higher functional potential energy. During path planning, this guides the surgical path to actively bypass high-network-value regions, thus improving the protection of the overall integrity of functional connectivity.

[0017] This invention employs anisotropic distance measurement based on diffusion tensor imaging to determine the functional influence range. By solving the equation, the geodesic distance along the fiber bundle is obtained, allowing the influence of the functional region to propagate along the white matter fiber bundle direction. Combining anisotropic distance and integrated network topology weights to construct a continuous functional potential field enables a more comprehensive reflection of the potential risk of functional region damage during path planning, achieving a transition from discrete obstacle avoidance to continuous potential field guidance.

[0018] This invention couples a continuous functional potential field with a physical obstacle field and employs a potential field-guided fast exploration random tree algorithm for path search. This algorithm concentrates sampling points in low-potential-energy regions through potential energy bias sampling, and integrates random directions with the negative gradient direction of the potential energy during expansion. This helps generate surgical paths that balance physical safety and functional protection while ensuring search efficiency.

[0019] This invention calculates the brain tissue displacement field using real-time intraoperative imaging or key point displacement data. Based on this displacement field, it dynamically distorts the preoperatively constructed functional potential field to obtain the real-time intraoperative functional potential field. Local path replanning within the displacement-affected area helps compensate for spatial errors caused by brain drift, ensuring the planned path remains consistent with the patient's current functional anatomy.

[0020] This invention uses augmented reality to overlay the real-time functional potential field during surgery onto the patient's actual anatomical structure as a color cloud map, and displays the optimal surgical path and instrument position in real time. When instruments deviate from the path or approach high-potential areas, a warning signal is issued, helping to provide doctors with intuitive risk perception and operational guidance, and improving the human-computer interaction process. Attached Figure Description

[0021] Figure 1 This is a flowchart illustrating the overall implementation of the method of the present invention; Figure 2 This is a flowchart illustrating the construction process of the continuous functional potential field of the present invention. Figure 3 This is the potential field-guided path planning and optimization branch diagram of the present invention; Figure 4 This is a flowchart of the intraoperative potential field update and local replanning process of the present invention. Detailed Implementation

[0022] 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.

[0023] As attached Figures 1 to 4 The intelligent planning method and system for neurosurgical surgical pathways shown herein are implemented as follows: Example 1: Preoperative Data Acquisition and Functional Network Construction This embodiment uses brain tumor puncture surgery as an example to illustrate the preoperative data processing flow.

[0024] Step S1: Acquire the patient's preoperative multimodal imaging data using a 3.0T MRI scanner: T1-weighted structural image: MPRAGE sequence, TR=2300ms, TE=2.98ms, voxel size 1mm×1mm×1mm.

[0025] Diffusion tensor imaging: EPI sequence, TR=8500ms, TE=89ms, b=1000s / mm 2 There are 30 gradient directions and the voxel size is 2mm×2mm×2mm.

[0026] Resting-state functional magnetic resonance imaging: EPI sequence, TR=2000ms, TE=30ms, voxel size 3mm×3mm×3mm.

[0027] Use FSLFLIRT to linearly register all images to the T1 image space to ensure spatial coordinate system consistency.

[0028] Brain functional region segmentation and 3D spatial mask acquisition: Based on T1 images, a U-Net convolutional neural network is used to segment brain functional regions. The input to this network is a 1×256×256×256 T1 image, with the image grayscale values ​​pre-normalized to the [0, 1] interval. The output is a probability map for each functional region, with the same size as the input image.

[0029] The network is pre-trained on the HCP dataset and adopts an encoder-decoder structure, which includes 4 downsampling and 4 upsampling operations, and skip connections to fuse high and low layer features.

[0030] The output probability map is thresholded; in this embodiment, the threshold is set to 0.5, resulting in a binarized three-dimensional spatial mask for each functional area. .

[0031] Perform connected component analysis on the probabilistic graph and remove components with a volume smaller than 100 mm. 3 The threshold for isolated areas is based on clinical experience and can be appropriately adjusted according to the anatomical volume of different functional areas. This embodiment adopts a unified threshold for simplification to ensure the anatomical rationality of the mask.

[0032] The segmented functional areas include: Broca's area, Wernicke's area, primary motor cortex, primary sensory cortex, supplementary motor area, thalamus, and basal ganglia, totaling seven functional areas.

[0033] Simultaneously, a lesion area mask is obtained based on clinical annotations, and the target points of the lesions are calculated. , Obtain the centroid coordinates of the lesion region. Entry point. The three-dimensional coordinates are pre-defined by neurosurgeons based on their experience with surgical approaches.

[0034] Construction of the structural connection matrix: Fiber tracing was performed using MRtrix3 software with the FACT algorithm (Fiber Assignment Continuous Tracing Algorithm). The algorithm takes diffusion tensor imaging data and FA maps as input and outputs a sequence of coordinates of whole-brain fiber tracts in .trk format.

[0035] The tracking parameters were set as follows: angle threshold 45°, FA threshold 0.2, step size 1mm, and seed points were evenly distributed throughout the whole brain with an interval of 2mm.

[0036] After fiber tracking is completed, artifact fibers with a length of less than 10mm are removed. Based on DTI resolution, fibers with a length of less than 10mm are mostly artifacts. If short fibers are located in important functional pathways such as the arcuate bundle, they need to be retained. In this embodiment, the process is simplified, and only fibers with FA values ​​greater than 0.2 are retained for subsequent analysis.

[0037] Count the number of fiber connections between each pair of brain functional areas. Average fiber length Average FA value Construct the structural connection matrix. Its elements Defined as: ; in This represents the maximum number of fibers connected across all functional areas. This represents the maximum FA value in the corresponding connection. For functional area pairs without fiber connections... .

[0038] Functional connection matrix construction: Extract the mask for each functional area The time series of all voxels were calculated, and the average time series for each functional region was calculated. FSLMCFLIT was used for head motion correction, and time points with head motion displacement greater than 2 mm were removed.

[0039] The average time series is preprocessed: linear drift is removed, and bandpass filtering is applied to retain the signal in the 0.01-0.1 Hz frequency band. The Pearson correlation coefficient is calculated for each pair of functional area average time series. The result obtained by Fisher Z transform This transformation makes the correlation coefficient satisfy the normal distribution assumption, which facilitates subsequent statistical analysis.

[0040] Functional connection matrix elements ,when Set to 0 to remove weak correlations.

[0041] Weighted brain functional network fusion: Fusion Structure Connection Matrix Functional connection matrix The weighted adjacency matrix is ​​obtained. Its elements Defined as: ; in This is the fusion coefficient, ranging from 0.3 to 0.7. It is dynamically adjusted based on the signal-to-noise ratio statistical analysis of structural and functional connections. This indicates a greater reliance on functional connectivity. This indicates a greater reliance on structural connections; in this embodiment, we take... Using brain functional areas as nodes, To construct an individualized weighted brain function network based on the edge weights.

[0042] Comprehensive network topology weight calculation: For each brain functional region node Calculate the following graph theory indices: Weighted Degree Centrality : It reflects the total strength of the connection between a node and other nodes.

[0043] Betweenness centrality The Brandes algorithm is used to calculate the number of nodes passed through. The ratio of the number of shortest paths to the total number of shortest paths.

[0044] eigenvector centrality Solving the characteristic equation ,in Adjacency matrix The largest eigenvalue is found by the power iteration method. It reflects the overall connectivity strength of the network.

[0045] Intra-module connectivity The Louvain algorithm is used to divide the network into functional modules. The input of the algorithm is a weighted adjacency matrix. The output is the module label for each functional area. For nodes The average connection strength with other nodes within its module.

[0046] Will , , , Normalized to the interval [0, 1] respectively, and then using the minimum-maximum normalization method, we obtain... , , , Comprehensive network topology weights The weighted sum of the above indicators: ; in ~ The weighting coefficients have values ​​ranging from 0 to 1 and satisfy the following conditions: Based on clinical expert scores and principal component analysis, the sample used in this embodiment was determined to be... , , , . Quantified functional areas Its pivotal role in the whole-brain functional network The larger the value, the more critical the functional area is in network information transmission; damage to this area may lead to more widespread functional impairments.

[0047] Step S2: Construct a continuous functional potential field Based on the clinical surgical approach, with the lesion target point and entry point Determine the surgical area by extending 20mm outwards from the connecting line as the center. This serves as the spatial domain for subsequent potential field calculations.

[0048] Shortest anisotropic distance calculation: for any point inside Calculate its value to each functional area Shortest anisotropic distance at the boundary This distance is defined as the distance from point [point name missing] along the direction of travel of the white matter fiber bundle. To the functional area The shortest path length at the boundary is obtained by solving the equation: ; in From point Arrival time to the functional area boundary For point The local propagation speed at that location. Anisotropy fraction in diffusion tensor imaging Related, defined as: ; The base velocity is set to 1 mm / unit time in this embodiment; The amplification factor is set to 2 in this embodiment. To avoid the speed being zero, let... Minimum value mm / unit time, boundary condition is within the functional area. .

[0049] The equation is solved using a fast-moving algorithm, which calculates the arrival time field by constructing a narrow band (with a width of 5 grid points, which can achieve a balance between computational accuracy and speed) and moving outward from the boundary. This algorithm is applicable to anisotropic velocity fields.

[0050] When narrowband propulsion covers the entire The algorithm terminates when the region is reached. The input to the algorithm is the velocity field. (Defined on grid points, the velocity at any point is obtained through trilinear interpolation) and functional region boundary conditions, the output is the arrival time field. . That is what we are looking for. .

[0051] Calculation of anisotropic decay potential energy: Each functional area At point Anisotropic decay potential energy generated at the location Defined as: ; in The comprehensive network topology weights obtained in step S1, Functional area sensitivity, functional area sensitivity Calculated based on volume and fiber coverage ratio: , Volume of functional area (unit: mm) 3 ), Fiber coverage ratio, i.e., the proportion of voxels within a functional region that are traversed by fiber bundles. Attenuation function. Using Gaussian form: ; The attenuation scale is measured in the range of 5-15 mm. mm indicates a limited scope of influence. mm indicates a wide range of influence, which is determined in this embodiment based on the statistical distribution range of the functional regions in the white matter fiber bundles: For language functional areas (Broca area, Wernicke area) mm, for the motor functional area (primary motor cortex) mm, for other functional areas take mm. For points within the function area, set directly. This is to ensure that the planned route does not cross the core of the functional area.

[0052] Total functional potential field: point Total functional potential energy at the point The sum of the potential energy generated by all functional zones at this point: ; exist The internal sampling is performed in a grid pattern with 1mm intervals to generate a three-dimensional array for storage. This forms a continuous functional potential field. The potential energy value at the boundary is extrapolated outward by 5mm (using nearest neighbor extrapolation), and then Gaussian smoothed. mm can effectively smooth boundary artifacts without changing the internal potential field distribution, thus ensuring continuity at the boundary.

[0053] Step S3: Path planning based on coupled field Physical obstacle field construction: Based on preoperative MRI images, key structures such as cerebral blood vessels, ventricles, and dura mater are segmented, and a binary obstacle mask is generated. For any point inside The Euclidean distance to the nearest obstacle is calculated through distance transformation. Physical obstacle risk Defined as: ; In this embodiment, the following is taken mm is used to avoid the denominator being zero.

[0054] Coupled field construction: functional potential field and physical obstacle risk field Normalized to the interval [0, 20], and then obtained by minimum-maximum normalization. and Coupled field Defined as the weighted sum of the two: ; in This is a weighting coefficient, ranging from 0.3 to 0.8, used to balance functional protection and physical obstacle avoidance. It is adjusted according to the lesion location: a larger value is used when the lesion is closer to the functional area; in this embodiment, it is [value missing]. .

[0055] Potential-guided fast exploratory random tree algorithm (PG-RRT*) With entry point Initialize a random tree for the root node Each node in the tree stores its coordinates, parent node index, and the cost from the root node to that node. Initially, the tree contains only the root node. .

[0056] Algorithm parameter settings: Maximum number of iterations in this embodiment Step length mm.

[0057] 1. Sampling steps: Generate random points .

[0058] Rejection sampling method is used: Candidate points are generated uniformly within the area. With probability Accept this point, where This is a temperature parameter, with a value range of 0.1-0.5. The sampling is highly concentrated in the low potential energy region. When sampling is nearly uniform, in this embodiment, a value of 0.2 (relative to...) is used. The scale [0, 20], 0.2 makes the sampling moderately concentrated. If accepted, then Otherwise, the sampling points are generated repeatedly. This strategy concentrates sampling points in low-potential regions.

[0059] 2. Nearest Neighbor Search: In a tree Search distance The nearest node is denoted as .

[0060] 3. Direction Fusion: Calculate random direction vectors The central difference method is used to calculate gradient at (One-sided difference is used at the boundary), negative gradient direction .

[0061] Integration direction ,in These are adaptive coefficients that vary with the number of iterations. The value increased linearly from 0.2 to 0.8 (based on experience, with an initial focus on exploration and a later focus on utilization): .

[0062] 4. Extended steps: along Direction by step Move to get a new node .

[0063] examine and Whether the line passes through an obstacle: Sample evenly at 0.5mm intervals along the line. 0.5mm is smaller than the minimum obstacle size (usually >2mm) to ensure detection accuracy.

[0064] If any sampling point is located within the obstacle mask, it is considered a collision; otherwise, it is considered a non-collision. Add to tree And calculate from the root node to The cost .

[0065] 5. Reconnect: In neighborhood radius Search for potential parent nodes within the node and select the one that makes them... The smallest node is selected as the parent node, and the tree structure is updated.

[0066] 6. Iteration check: Repeat steps 1-5 until the maximum number of iterations is reached. or With the target point The distance is less than .

[0067] Cost function definition: For any path (and ), its cost Defined as the integral along the path: ; in: ds is the differential of the path arc length. The trapezoidal rule is used for numerical integration, and the integration result is the path length. , unit mm.

[0068] For the path at point Tangential direction and potential energy gradient direction at the point The included angle, .when hour, .

[0069] ~ These are weighting coefficients, ranging from 0.1 to 5.0. The relative magnitude of each weight determines the optimization priority, and normalization can be used to make the sum equal to 1. In this embodiment, the weights are obtained through grid search on the validation set. , , , .

[0070] This cost function encourages the path to avoid regions of high physical risk and high functional potential, while also moving along equipotential lines (perpendicular to the gradient) to avoid crossing regions of abrupt potential changes.

[0071] After the algorithm converges, from the target point Backtracking to the root node This yields an initial path composed of polyline segments.

[0072] To facilitate smooth movement of surgical instruments along the path and avoid tissue damage caused by sharp turns, a cubic B-spline curve was used to smooth the path. The number of control points was set to 1 / 3 of the number of path nodes (an empirical value to balance smoothness and path fidelity), and a maximum deviation distance was set. The constraint is 0.5 mm (0.5 mm is less than the positioning accuracy of surgical instruments, which is usually 1 mm).

[0073] The smoothed path is resampled at 0.5mm intervals to verify whether path points enter high-risk areas (in this embodiment, normalized paths are assumed to be...). Or normalized High risk, and (Normalized to the [0, 20] interval). If high-risk points exist, fine-tune the control points and smooth again until the safety conditions are met.

[0074] Intraoperative pathway update: Step T1: Acquire intraoperative data. There are two ways to acquire intraoperative data: Method 1: Intraoperative MRI scan to obtain three-dimensional images of the current brain tissue, with the same scanning parameters as before the operation.

[0075] Method 2: Use the NDIPolaris optical positioning system to collect key point displacement data, including cortical vascular markers, lesion boundary points, etc., with a sampling frequency of 20Hz.

[0076] When both data sources are available, the intraoperative image registration result should be used first.

[0077] Step T2: Calculate the displacement field Different methods are used depending on the data source: If real-time intraoperative images are available: Non-rigid image registration is performed using ANTs software, employing a B-spline free deformation model, with mutual information as the similarity measure, 200 iterations, and a grid spacing of 10 mm. The algorithm inputs preoperative and intraoperative images, and the output is the whole-brain displacement field. .

[0078] If only the displacements of key points are obtained: using these points as boundary conditions, solve the linear elastic biomechanical model. The governing equations are: ; in , The value is the Lamé constant (the subscript m represents the material parameter), which is taken as in this embodiment. , , For vector Laplacian operators.

[0079] The brain tissue segmentation results were converted into tetrahedral meshes, generated using TetGen with a cell size of 2 mm, and solved using FEBio software. The inputs were boundary displacements and material parameters, and the output was the whole-brain displacement field. .

[0080] Calculate displacement field Then, check the mapping. Does the Jacobian determinant exist? If the region is folded, it indicates that the mapping has folded, and the displacement field needs to be Gaussian smoothed. mm) until the entire area The maximum number of iterations is 10.

[0081] Step T3: Update the potential field.

[0082] According to the displacement field The preoperative functional potential field is mapped to the current spatial location to obtain the real-time functional potential field during surgery. .

[0083] For the current spatial point The corresponding position before surgery is ,but ,like If not at a grid point, trilinear interpolation is used to obtain the result. .

[0084] Physical obstacle field Update in the same way, after the update and The data is transmitted in real time to the path planning module for subsequent replanning.

[0085] Step T4: Local replanning. Based on the current position of the surgical instruments. The root node is the target point of the lesion. For the target point, at the displacement amplitude Local replanning is performed within the specified area. In this embodiment, a displacement threshold of 1 mm is used.

[0086] Employing a dynamic fast-exploration random tree algorithm: preserving regions in the original tree located at low displacements ( For nodes and their subtrees, only prune path branches within the displacement-affected region, and within the affected region... Rerun the PG-RRT* algorithm for the root node to generate a new path fragment, and smoothly concatenate it with the preserved path (local smoothing is performed at the concatenation point to ensure continuity) to obtain the updated complete path.

[0087] After replanning, samples are taken along the new path at 0.5mm intervals (the same interval used for the smoothed path verification) to re-verify whether all points meet the requirements. and If the safety conditions are not met, return to step T4 to continue optimization until they are met.

[0088] Visual guidance: The following information is displayed in real time using the Microsoft HoloLens 2 augmented reality device: Functional potential field Displayed as a semi-transparent color cloud map, using a red-yellow-green color mapping and linear interpolation to generate continuous color transitions, with red corresponding to the normalized potential energy value. In the region where green corresponds to the normalized potential energy value The area.

[0089] Optimal surgical path: shown by blue lines, with distance information marked at key points along the path.

[0090] Current surgical instrument position: displayed by a green crosshair. In this embodiment, a deviation warning threshold is set: when the instrument deviates from the planned path by more than 1mm, the cursor turns red, and the beep frequency increases linearly with the deviation distance. ,in Distance (mm) The beep frequency (Hz) indicates that when the instrument is less than 2mm away from the high potential energy zone, the cursor will flash and the beep will continue, while an avoidance direction arrow will be displayed, pointing in the direction where the potential energy decreases the fastest.

[0091] System Implementation Example: A neurosurgical surgical path planning system, comprising: Preoperative modeling module: Used to perform steps S1-S2, taking DICOM format preoperative multimodal image data as input and outputting a continuous functional potential field. and physical obstacle risk field It is stored in a three-dimensional array format.

[0092] Path planning module: Used to execute step S3, with the following input: , Entry point and lesion target point The output is the coordinate sequence of the optimal surgical path.

[0093] Intraoperative update module: Used to perform operations T1-T4, connected to intraoperative imaging equipment or optical positioning system, to acquire intraoperative data in real time, dynamically update potential field and path, and output updated path coordinate sequence.

[0094] Navigation and guidance module: Used to perform visual guidance operations, connects to augmented reality devices, receives path data and instrument position data, generates a visual interface and overlays it in the doctor's field of vision.

[0095] The modules communicate with each other via a local area network, and the data formats used are NIfTI medical image format and custom path point sequence format.

[0096] The above description is only a specific embodiment of the present invention, but the protection scope of the present invention is not limited thereto. Any changes or substitutions that can be easily conceived by those skilled in the art within the technical scope disclosed in the present invention should be included within the protection scope of the present invention.

Claims

1. A method for intelligent planning of neurosurgical surgical pathways, characterized in that, Includes the following steps: S1: Acquire the patient's preoperative multimodal image data, segment brain functional areas based on the preoperative multimodal image data and obtain a binary three-dimensional spatial mask for each brain functional area, construct an individualized brain functional network, and calculate the comprehensive network topology weight of each brain functional area based on the brain functional network. S2: Determine the surgical area based on the target lesion point and the entry point into the skull. Based on the integrated network topology weights, the sensitivity of each brain functional area, and the anisotropic distance from the spatial point to the boundary of each functional area, construct a continuous functional potential field covering the surgical area. The anisotropic distance is the geodesic distance from the point to the boundary of the functional area along the direction of the white matter fiber bundle. S3: Construct a coupled field based on the continuous functional potential field and the physical obstacle field, and use a path search algorithm in the coupled field to plan the optimal surgical path from the entry point to the target lesion point.

2. The intelligent planning method for neurosurgical surgical pathways according to claim 1, characterized in that, The construction of personalized brain functional networks includes: Based on diffusion tensor imaging data, fiber connections between brain functional areas are extracted, and a structural connectivity matrix is ​​constructed. Based on functional magnetic resonance imaging data, time-series correlations of brain functional regions were extracted, and a functional connectivity matrix was constructed. By fusing the structural connection matrix and the functional connection matrix, a weighted adjacency matrix is ​​obtained.

3. The intelligent planning method for neurosurgical surgical pathways according to claim 2, characterized in that, The comprehensive network topology weights are obtained by weighted fusion of multiple graph theory indicators for each brain functional region node, including weighted degree centrality, betweenness centrality, eigenvector centrality, and intra-module connectivity.

4. The intelligent planning method for neurosurgical surgical pathways according to claim 1, characterized in that, The anisotropic distance is obtained by solving the equation, and the propagation speed in the equation is positively correlated with the anisotropic component value in the diffusion tensor imaging data. The potential energy value at any point in the continuous functional potential field is the superposition of the anisotropic decaying potential energy generated by all brain functional areas at that point. The anisotropic decaying potential energy is calculated by using the anisotropic distance as the independent variable and the comprehensive network topology weight as the proportional coefficient through a decay function.

5. The intelligent planning method for neurosurgical surgical pathways according to claim 1, characterized in that, The path search algorithm is a fast exploratory random tree algorithm guided by the potential field. Its sampling probability is determined according to the potential energy value of the continuous functional potential field. The lower the potential energy value, the higher the sampling probability of the region. Its node expansion direction is a weighted sum of the random sampling direction and the negative gradient direction of the potential energy.

6. The intelligent planning method for neurosurgical surgical pathways according to claim 5, characterized in that, The optimal surgical path is obtained by minimizing a preset cost function, which includes a path length integral term, a physical obstacle risk integral term along the path, a functional potential energy integral term along the path, and an angle penalty term between the path direction and the potential energy gradient direction.

7. The intelligent planning method for neurosurgical surgical pathways according to claim 1, characterized in that, It also includes the intraoperative pathway update step: T1: Acquire real-time intraoperative imaging data or key point displacement data collected by intraoperative sensors; T2: Calculate the displacement field of brain tissue by image registration based on the intraoperative real-time image data, or solve the displacement field of brain tissue by biomechanical model based on the key point displacement data. T3: Based on the displacement field, the continuous functional potential field constructed before surgery is mapped to the current spatial position to obtain the real-time functional potential field during surgery; T4: Based on the real-time functional potential field during the operation and the current position of the surgical instruments, perform local replanning of the remaining path.

8. The intelligent planning method for neurosurgical surgical pathways according to claim 7, characterized in that, The local replanning uses a dynamic fast exploratory random tree algorithm, which only prunes and re-expands the paths in the region affected by the displacement field.

9. The intelligent planning method for neurosurgical surgical pathways according to claim 1, characterized in that, It also includes a visualization guidance step: using an augmented reality device to overlay the continuous functional potential field as a semi-transparent color cloud map onto the patient's actual anatomical structure, and display the optimal surgical path in real time.

10. A neurosurgical surgical path planning system, characterized in that, The system is used to execute the intelligent planning method for neurosurgical surgical pathways according to any one of claims 1 to 9, the system comprising: The preoperative modeling module is used to perform operations such as acquiring preoperative multimodal image data, segmenting brain functional areas and obtaining binary three-dimensional spatial masks for each brain functional area, constructing individualized brain functional networks, calculating the topological weights of the integrated network, and constructing a continuous functional potential field. The path planning module is used to perform the operation of searching for the optimal surgical path in the coupled field formed by the continuous functional potential field and the physical obstacle field; The intraoperative update module is used to perform operations such as acquiring intraoperative displacement data, updating the continuous functional potential field according to the displacement field, and triggering local path replanning. The navigation guidance module is used to perform operations that visualize the continuous functional potential field and the optimal surgical path through an augmented reality device.

Images