CN106780518A - A kind of MR image three-dimensional interactive segmentation methods of the movable contour model cut based on random walk and figure - Google Patents

Abstract

The present invention discloses a kind of MR image three-dimensional interactive segmentation methods of movable contour model cut based on random walk and figure, and it includes step：S1, selected seed point partial 3 d MR view data of the interception comprising hypophysoma, and then obtain the initial boundary curved surface of segmentation；S2, based on initial boundary curved surface, sets up hybrid activity skeleton pattern, obtains model energy function；S3, discretization is carried out to energy function；S4, each tissue points for the partial 3 d MR images that will be intercepted carry out the structure of figure as the node of figure using 6 neighborhoods of each tissue points；Initial value is assigned respectively to the tissue points inside and outside initial boundary curved surface；And according to the energy function after discrete between node, between node and source point, the side connected and between node and meeting point assigns corresponding weights respectively；S5, carries out figure and cuts calculating based on the figure for building, and obtains segmentation result, and the boundary surface extracted in segmentation result is segmentation contour；S6, the initial boundary curved surface in S1 is substituted with current segmentation contour, and repeat step S2 to S5, until segmentation result is restrained, obtains final segmentation contour.Using the three-dimensional segmentation of the achievable hypophysoma MR images of the present invention, the accuracy of hypophysoma image segmentation is improved.

Description

3D Interaction of MR Images Based on Active Contour Model Based on Random Walk and Graph Cut segmentation method

technical field

The invention relates to the technical field of image segmentation, in particular to a three-dimensional interactive segmentation method of an MR image based on an active contour model of random walk and graph cut for area segmentation of a magnetic resonance imaging human brain image.

Background technique

Magnetic resonance imaging (MRI) has excellent soft tissue resolution and can image various parts of the human body from multiple angles and planes. The application of various imaging methods and pulse sequence technology is helpful for the localization and qualitative diagnosis of pituitary tumors . Usually, tumor segmentation in clinical medical images is done by doctors manually drawing the boundaries of pituitary tumors sequentially on two-dimensional images. According to the Response Evaluation in Solid Tumors (RECIST) standard, the radiologist will identify the high-brightness/low-brightness tissue of the brain on the two-dimensional image, and measure the diameter of the pituitary tumor lesion in one or two directions. However, due to the variable shape of the pituitary tumor , a single diameter does not represent the volume of pituitary tumors. In addition, due to the huge amount of data in the brain scan sequence images, manual segmentation is time-consuming and laborious. Moreover, the results of manual segmentation and image information analysis are affected by the professional knowledge and operational proficiency of doctors, which often vary from person to person, are highly subjective, and are prone to unpredictable results. Multimodal pituitary tumor magnetic resonance (Magnetic resonance, MR) image processing and analysis, comprehensive use of modern computer technologies such as image processing, image analysis, pattern recognition, machine learning, etc., can greatly improve the current status of clinical medical image processing and analysis of pituitary tumors, Provide more efficient, accurate and objective analysis results.

Pituitary tumors belong to one type of brain tumors. In clinical practice, it is found that the medical imaging of multimodal pituitary tumors has complex characteristics. The characteristics of its MR images are similar to those of other brain tumors. The main manifestations are: the structure of the pituitary tumor is complex; its space occupation, shape and size are uncertain; There is aliasing, and the boundaries are blurred; there are large differences in the organizational structure of different individuals; and due to the influence of various factors in the principle of MRI technology and the image acquisition process, there are noises and artifacts in the image. Therefore, it is very difficult to achieve the correct segmentation of pituitary tumors. The segmentation tasks of brain tumors include normal brain tissue (white matter, gray matter and cerebrospinal fluid), diseased tissue (edema, tumor and its internal necrosis and cystic degeneration, etc.). Since the 1970s, various types of brain tumor segmentation methods have emerged, which can be mainly divided into region-based segmentation methods, edge-based segmentation methods, and segmentation methods that combine regions and boundaries, as well as neural networks, fuzzy aggregation, etc. classes, Markov random field models, graph theory-based and graph-based other segmentation methods.

Due to the complexity of brain medical images, various segmentation methods can only solve specific problems, and there is no general segmentation method. Current research trends present the following characteristics:

(1) Combination of multiple segmentation algorithms that combine their respective advantages. For example, the image segmentation method combined with fuzzy clustering and deformation model, the segmentation method based on region growing and neural network; and the segmentation algorithm based on graph theory developed in recent years, such as joint segmentation algorithm of graph cut and active shape model, graph Cutting and deformation model joint segmentation algorithm, etc.

(2) Semi-automatic segmentation technology is still the key direction of research. Fully automatic segmentation algorithm can greatly improve the segmentation efficiency of brain tumors, but in clinical practice, it is necessary to allow doctors to perform special segmentation operations according to specific cases. Based on the automatic segmentation results, the interactive semi-automatic segmentation method enables doctors to evaluate the segmentation results and manually intervene in the segmentation, making up for the shortcomings of the automatic segmentation algorithm in performing special segmentation tasks, and achieving accurate segmentation and quantitative analysis. For example, brain tumor and lung tumor segmentation algorithms using random walk; also developed some interactive segmentation software, such as 3D Slicer and so on.

The automatic and semi-automatic interactive brain tumor segmentation algorithms that have emerged in recent years can be roughly divided into two categories based on machine learning: generative algorithms and discriminative algorithms. The generative algorithm learns the gray-scale spatial distribution of healthy brain tissue types from a series of prior images (images with known tissue classifications), and uses the learned distribution function to "predict" healthy tissue and find boundaries to segment tumors. The discriminative algorithm directly learns the grayscale difference between tumor and normal brain tissue, and uses the learned discriminant function to classify unlabeled pixels/voxels. The learning step can be performed either offline using some manually labeled patient training images, or online by the user marking seed regions on the current patient image. Generative algorithms have advantages for normal brain tissue segmentation, but have difficulty in segmenting tumors with arbitrary shapes. The advantage of the discriminant algorithm is that it combines the gray distribution of normal brain tissue and tumor, which is beneficial to the segmentation of tumors; but it only uses local image features, and the requirements for image gray levels of different tissues in multimodal image data are relatively consistent. This is difficult to achieve, especially for MR images. Combining the two methods, fusing background anatomical knowledge of brain structures and grayscale characteristics of lesions, may improve the accuracy of segmentation.

In recent years, a multimodal image segmentation method for brain tumors based on deformation models has also emerged. This type of method obtains prior knowledge such as constraint information, target position, size, and shape from the image, which not only utilizes the underlying information of the image itself but also combines the upper information of the human anatomical structure, and is suitable for the three-dimensional data of multimodal images of brain tumors. segmentation. In the deformation model, the geometric deformation model based on the level set theory has been widely valued and developed because it can solve the problem of topology change when the curve evolves, and other characteristics are superior to the parametric deformation model. A single level set function can only divide the image into two parts, the target area and the background area, that is, the image two-phase segmentation. In order to meet the needs of multiple target area segmentation in the image, a multiphase level set (Multiphase Level Set) method was produced. After nearly two decades of development, the idea of multiphase level sets has been gradually applied to the Mumford-Shah model and the Chan-Vese model. The LBF model proposed by Li et al. introduces the Gaussian kernel function, which enhances the ability of the level set method to capture edges. Only one level set function is used on some multiphase images, which solves the problem that the multiphase level set algorithm needs to use multiple levels. The computational complexity problem brought by the set function; and developed into a normalized LBF model by Peng et al. But at present, the solution of level set method still has the problems of complex method and large amount of calculation.

The research team represented by J.Egger proposed the use of graph theory-based graph cut and graph search methods, and designed a template-based semi-automatic pituitary tumor segmentation method. However, judging from their published research papers, the main cases used are regular and approximately elliptical pituitary structures, and the segmentation difficulty is relatively small. Adding a shape model to the energy function of the graph cut algorithm and the deformation model is also a common segmentation method for segmenting brain tumors with relatively fixed shape and little change. However, most pituitary tumors have irregular shapes and may infiltrate other brain tissues. Therefore, templates or shape models play a limited role in the segmentation of irregularly shaped pituitary tumors.

Contents of the invention

The technical problem to be solved by the present invention is to propose a three-dimensional interactive MR image based on a random walk (Random walk, RW) and graph cut (Graph cut) active contour model (Random walks and graph cuts based active contourmodel, RGCACM) The segmentation method realizes the three-dimensional segmentation of pituitary tumor MR images and improves the accuracy of pituitary tumor image segmentation.

The technical solution adopted by the present invention is specifically: a method for three-dimensional interactive segmentation of MR images based on an active contour model of random walk and graph cut, comprising the following steps:

S1, to obtain the initial boundary surface, including steps:

S11, selecting the seed point required by the random walk algorithm, using the seed point as the three-dimensional center point, and intercepting the local three-dimensional MR image data including the pituitary tumor from the three-dimensional brain MR data;

S12, using a three-dimensional random walk algorithm to process the intercepted local three-dimensional MR image data to obtain a three-dimensional map of the probability that each voxel point reaches the seed point;

S13, using the maximum expected value algorithm to select a probability threshold to segment the suspected target area, using the voxels in the suspected target area as the foreground point of the active contour model algorithm, and the boundary surface of the suspected target area is the initial boundary surface;

S2. Based on the initial boundary surface, a hybrid active contour model is established: the standard geodesic model is used for the boundary energy item in the energy function of the standard active contour model, and the gray maximum posterior probability described by the local Gaussian model is used for the regional energy item, and then obtained The energy function of the hybrid active contour model;

S3, model discretization, including:

S31, define a binary variable for each voxel point, and assign values 1 and 0 to represent foreground points and background points respectively, so as to discretize the regional energy term function;

S32, using the cut measure to discretize the boundary energy term function, that is, the geodesic model;

S33, obtaining a discretized energy function;

S4, Constructing a graph: For the intercepted local 3D MR image data, each voxel point is used as a node of the graph, and the 6 neighborhoods of each voxel point are used to construct the graph; for the initial boundary surface and the initial boundary surface The voxel points of are assigned initial values; and according to the energy function after discretization, assign corresponding weights to the edges connecting nodes, the edges connecting nodes to source points, and the edges connecting nodes to sink points;

S5, performing graph cut calculation based on the graph constructed in step S4: using the maximum flow and minimum cut algorithms to obtain the segmentation result, and extracting the boundary surface in the segmentation result as the segmentation contour;

S6, using the current segmentation contour as the initial boundary surface, iterating step S2 to step S5 repeatedly until the segmentation result converges, and outputting the final segmentation result.

Further, the present invention also includes step S7, using three-dimensional median filtering to process the evolution surface obtained in step S5, to obtain a filtered segmentation result; step S6 uses the current filtered segmentation contour as the initial boundary surface.

Preferably, in step S11, the step of selecting the seed point is to select a slice at the central position in the sagittal two-dimensional display image, and select a point at the approximate central position of the pituitary tumor in the slice as the seed point. The approximate central position is near the central position. This selection step is manual selection, so judging whether the selected seed point is close to the central position depends on the selector, and is based on the position as close as possible to the central position. The general error range is 30% of the diameter of the pituitary tumor. % within.

Preferably, in step S11, the size of the local three-dimensional MR image data including the pituitary tumor is 71×71×41 voxels.

In step S13, it is determined according to the gray value distribution of the probability map, and the value range is 0-1.

Preferably, in step S2, the given image domain I(x): Ω→R is a two-dimensional image; the closed curve C divides the image into two independent areas, and L(C) represents the length of C; let the area inside the curve C be Ω ₁ , and the area outside the curve C be Ω ₂ , for each pixel p in the image domain Ω, consider its circular neighborhood with radius ρ, defined as O _x ＝{y:|xy|≤ρ}, x is the central pixel of the circular neighborhood, and y is the circle Any pixel in the shape neighborhood except x;

Then the pixel gray probability density function in the area O _x ∩Ω _i is expressed as:

where u _i (x) and σ _i (x) are the local gray mean and variance of pixel x in area O _x ∩Ω _i respectively, I(y) is the gray value of pixel y,

The boundary energy item E _Edge adopts the standard geodesic model, and the regional energy item E _Reigon adopts the energy function E in the form of the gray maximum posterior probability described by the local Gaussian model as follows:

where β is an arbitrary constant, and the function w is a weight function.

Further, in step S2 of the present invention, the weight function w adopts a truncated Gaussian kernel form, expressed as:

In the formula, α is a constant that makes ∫w(x, y)=1. That is, the present invention assigns a large weight value to the gray value I(y) of the pixel y close to the center point x in the neighborhood O _x .

The discretization of the model in step S3 of the present invention is to integrate with the energy function of the graph cut algorithm. Preferably, in step S31, a binary variable x _p is defined for each voxel point p to discretize the area energy term of the energy function:

Object represents the foreground point, which is the segmentation target, and Background represents the background point;

The discretization of the regional energy term function is expressed as:

p, q refer to different voxels in the image, x _p , x _q refer to binary values (0 or 1) of different voxels, and N(p) refers to the neighborhood of voxel p.

u _s (p), u _t (p) and σ _s (p), σ _t (p) are the discrete forms of u _i (x) and σ _i (x) respectively, expressed as follows:

In step S32, the boundary energy term function is discretized and expressed as:

ω _pq is a non-negative edge weight, expressed in a three-dimensional form:

ΔΦ _pq = Δφ _pq · Δψ _pq = π ² /4, corresponding to the angular interval of the unit circle in the 6-neighborhood, δ is the size of the unit cell, and the possible value δ=1, |e _pq | is the unit length of the boundary e _pq , the value |e _pq |=1 in the present invention, D is a constant Riemann measure. Preferably, ω _pq =π/4 in the present invention.

Based on the discretized area items and boundary items, the discretized hybrid activity model energy function can be obtained.

Preferably, in order to avoid the balance problem of the area and boundary items in the standard energy function of the graph cut algorithm, and the area item will affect the graph cut algorithm tends to segment small areas. The energy function after discretization adopts the form of multiplication, that is, the regional energy item is used as the weight item of the boundary energy item, and the energy function after discretization is expressed as:

E＝E _Edge (p,q)·E' _Region (p) (10)

If the neighborhood voxel p, q are located inside or outside the segmentation target at the same time, there are |E _s (p)+E _t (q)|＞0 and |E _s (q)+E _t (p)|＞0 hold ; If p, q are located on both sides of the target boundary, then |E _s (p)+E _t (q)|≈0 or |E _s (q)+E _t (p)|≈0 holds true. So only when the neighborhood voxels p, q are located on both sides of the target boundary, that is, when p, q belong to the target and background regions respectively, the formula (10) can get the minimum value.

Preferably, in step S4, define the edge n-link weights connected to each other between voxel points as E _Edge (p,q)·E' _Region (p); if |E _s (p)|＞|E _t (p)|and x _p =1, then the weight of the edge t-link connecting the voxel point and the sink point T (referring to the segmentation target ₎ is set to infinity; if |E _s (p)|<|E _t (p)|and x _p =0, then the weight of the edge t-link connecting the voxel point to the source point S (referring to the background) is set to infinity. Finally, the maximum flow/minimum cut algorithm can be used to solve the segmentation result, that is, the energy function of the graph is solved to obtain the label (0 or 1) of each node in the graph, which is the prior art.

Beneficial effect

The method of the present invention is based on the improvement and expansion of the Graph cuts based active contour model (GCACM) algorithm, and utilizes the active contour model to solve the stepping defect of the segmentation edge of the graph cut algorithm, and also utilizes the maximum flow/ The minimum cut algorithm solves the problem of computational complexity of the level set method of the active contour model and the disadvantage of being easily trapped in local points. In order to solve the problem of boundary seed point selection in the graph cut algorithm, the random walk algorithm is used to roughly determine the target area, and the relatively fixed anatomical position of the pituitary tumor (behind the nose, in the middle of the eyes, below the optic chiasm, etc.) is used as a random walk algorithm Manual interaction provides seed points. Combining position prior knowledge with graph cut and active contour models to achieve 3D segmentation of pituitary tumors can improve the effectiveness of quantitative analysis such as pituitary tumor localization, pathological association with surrounding tissues, and volume measurement.

That is, the present invention provides a feasible and effective semi-automatic pituitary tumor MR image segmentation algorithm for the first time, combining two typical image segmentation algorithms based on graph theory and active contour model, which can realize regular and irregular pituitary tumors Segmentation and volume estimation, especially for invasive pituitary tumors, also have good accuracy. The interactive operation is simple and easy, and the remaining segmentation tasks can be automatically completed only by selecting two points with the user's mouse. In the algorithm, the local Gaussian distribution is used to describe the gray level of the image area, which can effectively solve the problem of uneven gray level in the pituitary tumor lesion area of MR images.

Description of drawings

Fig. 1 shows the schematic flow sheet of the method of the present invention;

Fig. 2 shows the segmentation process and result schematic diagram of the present invention, wherein a-user selects a foreground point (yellow point) and a background point (green point) to be superimposed on the original sagittal two-dimensional image, and the red box represents the intercepted Region of interest; b- the initial contour superimposed on the probability map; c- the segmentation contour superimposed on the original image; d- the 3D rendering model of the segmentation result;

FIG. 3 is a schematic diagram of the circular neighborhood of the pixel in step S2, and the dotted circle represents the neighborhood of the pixel x.

detailed description

It will be further described below in conjunction with the accompanying drawings and specific embodiments.

The present invention is based on the MR image three-dimensional interactive segmentation method of the active contour model of random walk and graph cut, comprises the following steps:

S1, to obtain the initial boundary surface, including steps:

S11, selecting the seed point required by the random walk algorithm, using the seed point as the three-dimensional center point, and intercepting the local three-dimensional MR image data including the pituitary tumor from the three-dimensional brain MR data;

S12, using a three-dimensional random walk algorithm to process the intercepted local three-dimensional MR image data to obtain a three-dimensional map of the probability that each voxel point reaches the seed point;

S13, using the maximum expected value algorithm to select a probability threshold to segment the suspected target area, using the voxels in the suspected target area as the foreground point of the active contour model algorithm, and the boundary surface of the suspected target area is the initial boundary surface;

S2. Based on the initial boundary surface, a hybrid active contour model is established: the standard geodesic model is used for the boundary energy item in the energy function of the standard active contour model, and the gray maximum posterior probability described by the local Gaussian model is used for the regional energy item, and then obtained The energy function of the hybrid active contour model;

S3, model discretization, including:

S31, define a binary variable for each voxel point, and assign values 1 and 0 to represent foreground points and background points respectively, so as to discretize the regional energy term function;

S32, using the cut measure to discretize the boundary energy term function, that is, the geodesic model;

S33, obtaining a discretized energy function;

S4, Constructing a graph: For the intercepted local 3D MR image data, each voxel point is used as a node of the graph, and the 6 neighborhoods of each voxel point are used to construct the graph; for the initial boundary surface and the initial boundary surface The voxel points of are assigned initial values; and according to the energy function after discretization, assign corresponding weights to the edges connecting nodes, the edges connecting nodes to source points, and the edges connecting nodes to sink points;

S5, performing graph cut calculation based on the graph constructed in step S4: using the maximum flow and minimum cut algorithms to obtain the segmentation result, and extracting the boundary surface in the segmentation result as the segmentation contour;

S6, using the current segmentation contour as the initial boundary surface, iterating step S2 to step S5 repeatedly until the segmentation result converges, and outputting the final segmentation result.

Further, the present invention can also include step S7, which is implemented between S5 and S6, that is, the evolution surface obtained in step S5 is processed by three-dimensional median filtering to obtain the filtered segmentation result; step S6 uses the current filtered segmentation contour as the initial boundary surface.

Example 1

Referring to Fig. 1, the steps of the method in this embodiment are briefly described as follows: image segmentation seed point selection and initial boundary surface acquisition, iterative segmentation of the improved GCACM algorithm, and post-processing of the three-dimensional median filter of the segmentation results. That is, the following initialization steps, segmentation steps and post-processing steps are described in detail as follows.

1. Initialization steps. It is mainly for the user to interactively manually select the seed points required by the Random walk algorithm, and to obtain the initial boundary surface required by the segmentation algorithm.

In order to reduce the computational complexity of the data and make use of the user's medical background knowledge, the data cube containing the pituitary tumor was intercepted from the original three-dimensional brain MR data. The specific method is to select a slice at the central position in the sagittal two-dimensional display image, and then select a point in the approximate central position of the pituitary tumor with the mouse. Take this point as the center to intercept the data cube (in this experiment, the size of 71×71×41 is selected), and this point is also used as the foreground point of the Random walk algorithm, and then randomly select a point in the background area of the slice image as the random walk algorithm. background point.

The intercepted data cube is processed by the three-dimensional random walk algorithm, and the probability three-dimensional map of each voxel reaching the foreground seed point is obtained. The EM algorithm is used to select the probability threshold, and the area of the suspected target is initially segmented. The voxels in the area will be used as the foreground points of GCACM, and the boundary surface of the area is the initial boundary surface.

The process of obtaining the initial boundary surface is shown in Figure 2. Taking a sagittal two-dimensional image approximately at the center of the pituitary tumor as an example, first a- the user selects a foreground point (yellow point) and a background point (green point) to superimpose On the original sagittal two-dimensional image, the red box represents the intercepted region of interest; in b is the initial contour superimposed on the probability map, in c is the segmentation contour superimposed on the original image, and in d is the segmentation The resulting 3D rendered model.

2. Segmentation steps, including the design of the hybrid active contour model, model discretization, graph construction and solution.

(1) Establishment of hybrid active contour model

Referring to Figure 3, a given image domain I(x):Ω→R is a two-dimensional image. A closed curve C divides the image into two independent regions, and L(C) represents the length of C. Let the area inside the curve C be Ω ₁ , and the area outside the curve C be Ω ₂ . For each pixel p in the image domain Ω, consider its circular neighborhood with radius ρ, defined as O _x ={y:xy≤ρ}.

The pixel gray probability density function in the area O _x ∩Ω _i can be expressed as follows:

In the formula, u _i (x) and σ _i (x) are the local gray mean and variance of the pixel x in the area O _x ∩Ω _i , respectively.

The energy function combines the boundary energy item and the regional energy item, where the boundary item E _Edge adopts a standard geodesic model, and the regional item E _Reigon adopts a gray-scale local Gaussian model and is expressed as a maximum a posteriori probability (Maximum a posteriori probability, MAP) form. The energy function can be written as follows:

where β is an arbitrary constant. The function w is a weight function, which assigns a large weight value to the gray value I(y) of the pixel y close to the center point x in the neighborhood O _x . In this algorithm, w in the form of a truncated Gaussian kernel is used, expressed as follows:

where α is a constant such that ∫w(x,y)=1 holds.

(2) Model discretization

In order to integrate with the energy function of the graph cut algorithm, a binary variable x _p is defined for each voxel point p in the graph, and the values 1 and 0 are assigned to represent the foreground point (segmentation target) and background point (background), respectively, and the energy function of the area item Discretization, as shown in formula (4).

Using the cut metric proposed by Boykov et al. to discretize the geodesic model, the boundary term of the energy function is expressed as:

where ω _pq is a non-negative edge weight, which can be expressed in three-dimensional form as follows:

In the above formula, ΔΦ _pq = Δφ _pq · Δψ _pq = π ² /4 corresponds to the angular interval of the unit circle in the 6 neighborhoods, δ is the size of the unit cell, δ=1, |e _pq | is the unit length of the boundary e _pq , |e _pq |=1, D is a constant Riemann measure, and ω _pq =π/4 is taken in the present invention.

The discrete form of the area term of the energy function can be expressed as follows:

In the formula, p and q refer to different voxels in the image, x _p and x _q refer to the binary values (0 or 1) of different voxels, N(p) refers to the neighborhood of voxel p; u _s (p ), u _t (p) and σ _s (p), σ _t (p) are the discrete forms of u _i (x) and σ _i (x) respectively, as follows:

In order to avoid the balance problem of the area and boundary items in the standard energy function of the graph cut algorithm, and the area item will affect the graph cut algorithm tends to segment small areas. This algorithm uses the energy function in the form of multiplication, that is, the area item is used as the weight item of the boundary item, as shown in formulas (10) and (11):

E＝E _Edge (p,q)·E' _Region (p) (10)

If the neighboring pixels p and q are located inside or outside the segmentation target (ie pituitary tumor area) at the same time, there are |E _s (p)+E _t (q)|＞0 and |E _s (q)+E _t (p )|＞0 is established; if p, q are located on both sides of the target boundary, then |E _s (p)+E _t (q)|≈0 or |E _s (q)+E _t (p)|≈0 is established . So formula (10) can get the minimum value only when the neighboring pixels p, q are located on both sides of the target boundary. The minimum value of formula (10) means that the energy of the image is the minimum value, that is, the energy function achieves the optimal result, and when the energy reaches the optimal value, the segmentation is completed. At this time, the voxel points belonging to the target and background areas are assigned different tag value.

(3) Construction and solution of graph

The present invention uses each voxel in the MR data volume as a graph node, and uses six neighborhoods to construct a three-dimensional graph. The n-link weight of the connected edges between voxels in the graph is set to E _Edge (p,q)·E' _Regi o _n (p); if |E _s (p)|＞|E _t (p)| And x _p =1, then the t-link weight of the edge connecting the voxel to the sink T (corresponding target) is set to infinity; if |E _s (p)|<|E _t (p)| and x _p = 0, the t-link weight of the edge connecting the voxel to the source point S (corresponding to the background) is set to infinity. Finally, the maximum flow/minimum cut algorithm is used to solve the segmentation result, which is the prior art here and will not be described in detail.

3. Post-processing

The evolution surface obtained in the segmentation step is post-processed with a three-dimensional median filter, and then iteratively segmented according to the segmentation result until the result converges to obtain the final segmentation result

The iterative segmentation process of the present invention is briefly described as follows:

(1) Use the Random walk algorithm to obtain the initial boundary surface S _ur ; the binary variable x _p corresponding to the voxel p inside the initial assignment S _ur is 1, and the x _p corresponding to the voxel p outside S _ur is 0;

(2) Calculate the mean and variance u _s (p), u _t (p) and σ _s (p), σ _t (p) of the local gray scale using formulas (8) and (9);

(3) construct graph with the energy function represented by formula (10);

(4) Use the maximum flow/minimum cut algorithm to solve the segmentation result, and then update the value of the binary variable x _p corresponding to the voxel p according to the segmentation result. meta variable value;

(5) Carry out smooth post-processing to the segmentation result with a three-dimensional median filter;

(6) Replace the initial contour with the current segmentation contour, that is, the initial boundary surface in S2, and repeat steps (2)-(5) until the segmentation result converges, and the obtained segmentation result does not change, which means convergence.

4. Experimental results

Using the present invention to conduct a segmentation experiment on the MR T1W data of 23 pituitary tumor patients, the schematic diagrams of the results are shown in c and d in FIG. 2 . Table 1 shows the minimum value, maximum value, mean value and standard deviation value of Dice coefficient value (Dice similarity coefficient, DSC) of 23 case segmentation results. The DSC value is defined as the relative overlapping volume ratio of the experimental result V _E and the gold standard V _G , as shown in formula (12):

In addition, in order to show the superiority of the segmentation effect of the method of the present invention, this algorithm and two other published segmentation algorithms (GrowCut algorithm based on region growth and general GCACM algorithm loaded in the software 3DSlice) were compared for the DSC value of the segmentation result . Refer to Table 1 for the comparison results. It can be seen that the average DSC value of the segmentation result of this algorithm is 88.36%, which is better than 81.61% of the GrowCut algorithm and 70.88% of the general GCACM algorithm. It shows that this algorithm has high validity and accuracy in the application of pituitary tumor segmentation.

Table 1. DSC values of segmentation results of this algorithm and two published segmentation algorithms

Results Proposed method Grow Cut GCACM min 74.86% 67.62% 49.45% max 93.35% 89.07% 78.55% mean±SD 88.36%±5.61% 81.61%±7.37% 70.88%±8.08%

The above is only a preferred embodiment of the present invention, it should be pointed out that for those of ordinary skill in the art, without departing from the technical principle of the present invention, some improvements and modifications can also be made. It should also be regarded as the protection scope of the present invention.

Claims (10)

1. a kind of MR image three-dimensional interactive segmentation method based on the active contour model of random walk and graph cut, it is characterized in that, comprises the following steps: S1, to obtain the initial boundary surface, including steps: S11, selecting the seed point required by the random walk algorithm, using the seed point as the three-dimensional center point, and intercepting the local three-dimensional MR image data including the pituitary tumor from the three-dimensional brain MR data; S12, using a three-dimensional random walk algorithm to process the intercepted local three-dimensional MR image data to obtain a three-dimensional map of the probability that each voxel point reaches the seed point; S13, using the maximum expected value algorithm to select a probability threshold to segment the suspected target area, using the voxels in the suspected target area as the foreground point of the active contour model algorithm, and the boundary surface of the suspected target area is the initial boundary surface; S2. Based on the initial boundary surface, a hybrid active contour model is established: the standard geodesic model is used for the boundary energy item in the energy function of the standard active contour model, and the gray maximum posterior probability described by the local Gaussian model is used for the regional energy item, and then obtained The energy function of the hybrid active contour model; S3, model discretization, including: S31, define a binary variable for each voxel point, and assign values 1 and 0 to represent foreground points and background points respectively, so as to discretize the regional energy term function; S32, using the cut measure to discretize the boundary energy term function, that is, the geodesic model; S33, obtaining a discretized energy function; S4, Constructing a graph: For the intercepted local 3D MR image data, each voxel point is used as a node of the graph, and the 6 neighborhoods of each voxel point are used to construct the graph; for the initial boundary surface and the initial boundary surface The voxel points of are assigned initial values; and according to the energy function after discretization, assign corresponding weights to the edges connecting nodes, the edges connecting nodes to source points, and the edges connecting nodes to sink points; S5, performing graph cut calculation based on the graph constructed in step S4: using the maximum flow and minimum cut algorithms to obtain the segmentation result, and extracting the boundary surface in the segmentation result as the segmentation contour; S6, using the current segmentation contour as the initial boundary surface, repeating iteration steps S2 to S5 until the segmentation results converge, and outputting the final segmentation results. 2. The method according to claim 1, further comprising step S7, utilizing the three-dimensional median filter to process the evolution surface obtained in step S5, to obtain the filtered segmentation result; step S6 to obtain the segmented contour after the current filtering as the initial boundary surface. 3. The method according to claim 1 or 2, characterized in that, in step S11, the selection step of the seed point is, in the two-dimensional display image of the sagittal position, select a slice at the central position, within the pituitary tumor in the slice Select a point at the approximate center position as the seed point. 4. The method according to claim 1 or 2, wherein in step S11, the size of the data cube containing the pituitary tumor is 71×71×41 voxels. 5. The method according to claim 1 or 2, characterized in that, in step S2, the given image domain I(x): Ω→R is a two-dimensional image; the closed curve C divides the image into two independent areas, and L(C) represents the length of C; let the area inside the curve C be Ω ₁ , and the area outside the curve C be Ω ₂ , for each pixel p in the image domain Ω, consider its circular neighborhood with radius ρ, defined as O _x ＝{y:|xy|≤ρ}, x is the central pixel of the circular neighborhood, and y is the circle Any pixel in the shape neighborhood except x; Then the pixel gray probability density function in the area O _x ∩Ω _i is expressed as: ${\mathrm{<span\; class="notranslate">\; p</span>}}_{{\mathrm{<span\; class="notranslate">\; \&Omega;</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{\mathrm{<span\; class="notranslate">\; 1</span>}}{\sqrt{\mathrm{<span\; class="notranslate">\; 2</span>}\mathrm{<span\; class="notranslate">\; \&pi;</span>}}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; )</span>}\mathrm{<span\; class="notranslate">\; exp</span>}<span\; class="notranslate">\; (</span><span\; class="notranslate">\; -</span>\frac{{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; u</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}{\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; )</span>$ where u _i (x) and σ _i (x) are the local gray mean and variance of pixel x in area O _x ∩Ω _i respectively, I(y) is the gray value of pixel y, The boundary energy item E _Edge adopts the standard geodesic model, and the regional energy item E _Reigon adopts the energy function E in the form of the gray maximum posterior probability described by the local Gaussian model as follows: $\begin{array}{c}\mathrm{<span\; class="notranslate">\; E.</span>}<span\; class="notranslate">\; =</span>{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; E.</span>}\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; g</span>}\mathrm{<span\; class="notranslate">\; e</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; Re</span>}\mathrm{<span\; class="notranslate">\; g</span>}\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; o</span>}\mathrm{<span\; class="notranslate">\; no</span>}}<span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; 0</span>}}{\overset{\mathrm{<span\; class="notranslate">\; L</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; c</span>}<span\; class="notranslate">\; )</span>}{<span\; class="notranslate">\; \&Integral;</span>}}\mathrm{<span\; class="notranslate">\; g</span>}<span\; class="notranslate">\; (</span><span\; class="notranslate">\; \&dtri;</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; c</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; +</span>\underset{\mathrm{<span\; class="notranslate">\; i</span>}<span\; class="notranslate">\; =</span>\mathrm{<span\; class="notranslate">\; 1</span>}}{\overset{\mathrm{<span\; class="notranslate">\; 2</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}}\underset{{\mathrm{<span\; class="notranslate">\; \&Omega;</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}}{<span\; class="notranslate">\; \&Integral;</span><span\; class="notranslate">\; \&Integral;</span>}<span\; class="notranslate">\; \&lsqb;</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; log</span>}{\mathrm{<span\; class="notranslate">\; p</span>}}_{{\mathrm{<span\; class="notranslate">\; \&Omega;</span>}}_{\mathrm{<span\; class="notranslate">\; i</span>}}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; \&rsqb;</span>\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; x</span>}\\ <span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; 0</span>}}{\overset{\mathrm{<span\; class="notranslate">\; L</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; c</span>}<span\; class="notranslate">\; )</span>}{<span\; class="notranslate">\; \&Integral;</span>}}\frac{\mathrm{<span\; class="notranslate">\; 1</span>}}{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; \&beta;</span>}<span\; class="notranslate">\; |</span><span\; class="notranslate">\; \&dtri;</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; c</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span><span\; class="notranslate">\; |</span>}\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; the\; s</span>}<span\; class="notranslate">\; -</span>\underset{{\mathrm{<span\; class="notranslate">\; \&Omega;</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}{<span\; class="notranslate">\; \&Integral;</span><span\; class="notranslate">\; \&Integral;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; log</span>}{\mathrm{<span\; class="notranslate">\; p</span>}}_{{\mathrm{<span\; class="notranslate">\; \&Omega;</span>}}_{\mathrm{<span\; class="notranslate">\; 1</span>}}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; x</span>}\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; -</span>\underset{{\mathrm{<span\; class="notranslate">\; \&Omega;</span>}}_{\mathrm{<span\; class="notranslate">\; 2</span>}}}{<span\; class="notranslate">\; \&Integral;</span><span\; class="notranslate">\; \&Integral;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; log</span>}{\mathrm{<span\; class="notranslate">\; p</span>}}_{{\mathrm{<span\; class="notranslate">\; \&Omega;</span>}}_{\mathrm{<span\; class="notranslate">\; 2</span>}}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; x</span>}\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; the\; y</span>}\end{array}<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 2</span>}<span\; class="notranslate">\; )</span>$ where β is an arbitrary constant, and the function w is a weight function. 6. The method according to claim 5, characterized in that, in step S2, the weight function w adopts a truncated Gaussian kernel form, expressed as: $\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\left\{\begin{array}{ccc}\frac{\mathrm{<span\; class="notranslate">\; 1</span>}}{\mathrm{<span\; class="notranslate">\; \&alpha;</span>}}\mathrm{<span\; class="notranslate">\; exp</span>}<span\; class="notranslate">\; (</span><span\; class="notranslate">\; -</span>\frac{{<span\; class="notranslate">\; |</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; |</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}{\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; ,</span>& \mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; f</span>}& <span\; class="notranslate">\; |</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; |</span><span\; class="notranslate">\; \&le;</span>\mathrm{<span\; class="notranslate">\; \&rho;</span>}\\ \mathrm{<span\; class="notranslate">\; 0</span>}<span\; class="notranslate">\; ,</span>& \mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; f</span>}& <span\; class="notranslate">\; |</span>\mathrm{<span\; class="notranslate">\; x</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; the\; y</span>}<span\; class="notranslate">\; |</span><span\; class="notranslate">\; ></span>\mathrm{<span\; class="notranslate">\; \&rho;</span>}\end{array}\right.<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 3</span>}<span\; class="notranslate">\; )</span>$ In the formula, α is a constant that makes ∫w(x, y)=1. 7. The method according to claim 5, characterized in that, in step S31, a binary variable x _{p is defined for each voxel point p} : ${\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}}<span\; class="notranslate">\; =</span>\left\{\begin{array}{cc}\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; ,</span>& \mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; \&Element;</span>\mathrm{<span\; class="notranslate">\; o</span>}\mathrm{<span\; class="notranslate">\; b</span>}\mathrm{<span\; class="notranslate">\; j</span>}\mathrm{<span\; class="notranslate">\; e</span>}\mathrm{<span\; class="notranslate">\; c</span>}\mathrm{<span\; class="notranslate">\; t</span>}\\ \mathrm{<span\; class="notranslate">\; 0</span>}<span\; class="notranslate">\; ,</span>& \mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; \&Element;</span>\mathrm{<span\; class="notranslate">\; B</span>}\mathrm{<span\; class="notranslate">\; a</span>}\mathrm{<span\; class="notranslate">\; c</span>}\mathrm{<span\; class="notranslate">\; k</span>}\mathrm{<span\; class="notranslate">\; g</span>}\mathrm{<span\; class="notranslate">\; r</span>}\mathrm{<span\; class="notranslate">\; o</span>}\mathrm{<span\; class="notranslate">\; u</span>}\mathrm{<span\; class="notranslate">\; no</span>}\mathrm{<span\; class="notranslate">\; d</span>}\end{array}\right.<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 4</span>}<span\; class="notranslate">\; )</span>$ Object represents the foreground point, which is the segmentation target, and Background represents the background point; The regional energy term of the energy function is discretized, and the discretized regional energy term function is expressed as: $\begin{array}{c}{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; Re</span>}\mathrm{<span\; class="notranslate">\; g</span>}\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; o</span>}\mathrm{<span\; class="notranslate">\; no</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span>\\ <span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; p</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\underset{\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; \&Element;</span>\mathrm{<span\; class="notranslate">\; N</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span>}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; log\&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; +</span>\frac{{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; u</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}{\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; )</span>\\ <span\; class="notranslate">\; +</span>\underset{\mathrm{<span\; class="notranslate">\; p</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\underset{\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; \&Element;</span>\mathrm{<span\; class="notranslate">\; N</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span>}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; log\&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; +</span>\frac{{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; u</span>}}_{\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}{\mathrm{<span\; class="notranslate">\; 2</span>}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; t</span>}}{<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}\end{array}<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 5</span>}<span\; class="notranslate">\; )</span>$ P and q refer to different voxels in the image, x _p and x _q refer to binary values (0 or 1) of different voxels, and N(p) refers to the neighborhood of voxel p; u _s (p), u _t (p) and σ _s (p), σ _t (p) are the discrete forms of u _i (x) and σ _i (x) respectively as follows: $\left\{\begin{array}{c}{\mathrm{<span\; class="notranslate">\; u</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{\underset{\mathrm{<span\; class="notranslate">\; q</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; )</span>}{\underset{\mathrm{<span\; class="notranslate">\; q</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; )</span>}\\ {\mathrm{<span\; class="notranslate">\; u</span>}}_{\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{\underset{\mathrm{<span\; class="notranslate">\; q</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}}{\underset{\mathrm{<span\; class="notranslate">\; q</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}}\end{array}\right.<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 6</span>}<span\; class="notranslate">\; )</span>$ $\left\{\begin{array}{c}{\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{\underset{\mathrm{<span\; class="notranslate">\; q</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; u</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; )</span>}{\underset{\mathrm{<span\; class="notranslate">\; q</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; )</span>}\\ {\mathrm{<span\; class="notranslate">\; \&sigma;</span>}}_{\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\frac{\underset{\mathrm{<span\; class="notranslate">\; q</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span>{<span\; class="notranslate">\; (</span>{\mathrm{<span\; class="notranslate">\; u</span>}}_{\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}}{\underset{\mathrm{<span\; class="notranslate">\; q</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; w</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}}\end{array}\right.<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 7</span>}<span\; class="notranslate">\; )</span>$ In step S32, the boundary energy term function is discretized and expressed as: ${\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; E.</span>}\mathrm{<span\; class="notranslate">\; d</span>}\mathrm{<span\; class="notranslate">\; g</span>}\mathrm{<span\; class="notranslate">\; e</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; ,</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; p</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\underset{\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; \&Element;</span>\mathrm{<span\; class="notranslate">\; N</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span>}{<span\; class="notranslate">\; \&Sigma;</span>}\frac{{\mathrm{<span\; class="notranslate">\; \&omega;</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; (</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}}<span\; class="notranslate">\; )</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; -</span>{\mathrm{<span\; class="notranslate">\; x</span>}}_{\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; )</span>}{\mathrm{<span\; class="notranslate">\; 1</span>}<span\; class="notranslate">\; +</span>\mathrm{<span\; class="notranslate">\; \&beta;</span>}<span\; class="notranslate">\; |</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; -</span>\mathrm{<span\; class="notranslate">\; I</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; |</span>}<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 8</span>}<span\; class="notranslate">\; )</span>$ ω _pq is a non-negative edge weight, expressed in a three-dimensional form: ${\mathrm{<span\; class="notranslate">\; \&omega;</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; =</span>\frac{{\mathrm{<span\; class="notranslate">\; \&delta;</span>}}^{\mathrm{<span\; class="notranslate">\; 3</span>}}<span\; class="notranslate">\; \&CenterDot;</span>{<span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; e</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; |</span>}^{\mathrm{<span\; class="notranslate">\; 3</span>}}<span\; class="notranslate">\; \&CenterDot;</span>{\mathrm{<span\; class="notranslate">\; \&Delta;\&Phi;</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; \&CenterDot;</span>\mathrm{<span\; class="notranslate">\; det</span>}\mathrm{<span\; class="notranslate">\; D.</span>}}{\mathrm{<span\; class="notranslate">\; \&pi;</span>}<span\; class="notranslate">\; \&CenterDot;</span>{<span\; class="notranslate">\; (</span>{{\mathrm{<span\; class="notranslate">\; e</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}\mathrm{<span\; class="notranslate">\; q</span>}}}^{\mathrm{<span\; class="notranslate">\; T</span>}}<span\; class="notranslate">\; \&CenterDot;</span>\mathrm{<span\; class="notranslate">\; D.</span>}<span\; class="notranslate">\; \&CenterDot;</span>{\mathrm{<span\; class="notranslate">\; e</span>}}_{\mathrm{<span\; class="notranslate">\; p</span>}\mathrm{<span\; class="notranslate">\; q</span>}}<span\; class="notranslate">\; )</span>}^{\mathrm{<span\; class="notranslate">\; 2</span>}}}<span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 9</span>}<span\; class="notranslate">\; )</span>$ ΔΦ _pq = Δφ _pq · Δψ _pq = π ² /4, corresponding to the angular interval of the unit circle in the 6-neighborhood, δ is the size of the unit cell, |e _pq | is the unit length of the boundary e _pq , D is the constant Riemannian measure . 8. The method according to claim 7, characterized in that ω _pq =π/4. 9. The method according to claim 7, characterized in that, the energy function after the discretization adopts a multiplication form, that is, the regional energy item is used as the weight item of the boundary energy item, and the energy function after the discretization is expressed as: E＝E _Edge (p,q)·E' _Region (p) (10) ${\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; Re</span>}\mathrm{<span\; class="notranslate">\; g</span>}\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; o</span>}\mathrm{<span\; class="notranslate">\; no</span>}}^{<span\; class="notranslate">\; \&prime;</span>}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; =</span>\underset{\mathrm{<span\; class="notranslate">\; p</span>}}{<span\; class="notranslate">\; \&Sigma;</span>}\mathrm{<span\; class="notranslate">\; m</span>}\mathrm{<span\; class="notranslate">\; i</span>}\mathrm{<span\; class="notranslate">\; no</span>}<span\; class="notranslate">\; (</span><span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; |</span><span\; class="notranslate">\; ,</span><span\; class="notranslate">\; |</span>{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; the\; s</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; q</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; +</span>{\mathrm{<span\; class="notranslate">\; E.</span>}}_{\mathrm{<span\; class="notranslate">\; t</span>}}<span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; p</span>}<span\; class="notranslate">\; )</span><span\; class="notranslate">\; |</span><span\; class="notranslate">\; )</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; -</span><span\; class="notranslate">\; (</span>\mathrm{<span\; class="notranslate">\; 11</span>}<span\; class="notranslate">\; )</span>$ If the neighboring pixels p and q are located inside or outside the segmentation target at the same time, |E _s (p)+E _t (q)|＞0 and |E _s (q)+E _t (p)|＞0 are established; If p, q are located on both sides of the target boundary, then |E _s (p)+E _t (q)|≈0 or |E _s (q)+E _t (p)|≈0. 10. method according to claim 8, is characterized in that, define the edge n-link weight that interconnects between voxel point to be E _Edge (p, q)·E' _Region (p); If | E _s (p)|＞|E _t (p)| and x _p =1, then the weight of the edge t-link connecting the voxel point and the sink point T (referring to the segmentation target) is set to infinity; if |E _s (p)|<|E _t (p)| and x _p =0, then the weight of the edge t-link connecting the voxel point to the source point S (referring to the background) is set to infinity.