CN103606165A - Medical image processing method based on network phase synchronization - Google Patents

Abstract

The invention discloses a medical image processing method based on network phase synchronization. The method mainly solves the problem that an existing medical image processing method is not high in diagnosis recognition rate and poor in local detail processing. The method comprises the implementation steps that (1) an original medical image is read in, the characteristics of pixel points are extracted, and watershed cutting is carried out; (2) the similarity between areas after cutting is calculated; (3) the areas are used as nodes, the similarity is used as a network weight, and a network is built; (4) the initial phase of the network is generated randomly, an iteration is carried out on the phase until the phase is stable, and the stable node phase is extracted; (5) original classification is carried out on the nodes through the stable phase, the slope of the nodes is calculated, and disaggregated classification is carried out on the nodes according to the slope; (6) the misclassified nodes are redistributed, and category serial numbers are signed to the pixel points; (7) different grey values are assigned to the inhomogeneous pixel points, and the split image is obtained. Compared with an existing method, the diagnosis recognition rate of the medical image is improved.

Description

Phase locked medical image processing method Network Based

Technical field

The invention belongs to image processing field, relate to the processing of medical image, can be used for monitoring of diseases and distribute, study pathogenesis and disease auxiliary diagnosis.

Background technology

Medical image comprises the shooting of CT, positron emission chromatography imaging technique PET, single photon radial faults SPECT, Magnetic resonance imaging MR, and the image that obtains of other medical imaging device.Complicacy and diversity due to medical image, and the property difference of the image-forming principle of medical image and tissue itself, and the formation of image is subject to the impact such as noise, field offset effect, local bulk effect and histokinesis etc., medical image is compared with normal image, has the features such as fuzzy, unevenness.In addition, human dissection institutional framework and complex-shaped, the difference between Different Individual is large, also makes the problem of medical image processing more complicated.

In recent years, along with medical imaging is in clinical medical successful Application, the emerging technology of some medical image processing, as fuzzy mathematics, mathematical morphology, digital topology, artificial intelligence etc., makes remarkable progress, and new Medical Imaging Technology emerges in an endless stream.Wherein, the watershed algorithm based on mathematical morphology that the people such as Digabel proposes because of have computing velocity soon, cut apart the characteristics such as accurate and accurate positioning image edge, and be subject to extensive concern.The ultimate principle of watershed algorithm is that image is regarded as to the topological landforms in geodesy, and gray level image is considered as landform appearance, and grey scale pixel value is to the sea level elevation that should put.Each local minizing point, and its range of influence is called as reception basin, the border of reception basin is known as watershed divide.Watershed algorithm is by finding watershed divide to Image Segmentation Using.Yet the defect of this method is to have over-segmentation phenomenon.

Summary of the invention

The object of the invention is to for the problem that medical image is fuzzy, inhomogeneous, individual difference is large, propose a kind of phase locked medical image processing method Network Based, the over-segmentation problem existing to solve existing method, the diagnosis and distinguish rate of raising medical image.

Technical thought of the present invention is: the Image Segmentation Using with dividing ridge method after to feature extraction, the zonule after cutting apart is corresponded to the node in network, and set up network model; By oscillator formula, node is carried out to phase place iteration, complete the merging of overdivided region, improve the effect that image is cut apart.Implementation step comprises as follows:

(1) extract the pixel feature of primitive medicine image, and the pixel feature of extracting is carried out to watershed segmentation, obtain N region, N >=1000;

(2) calculate the similarity between every two regions, and each region and the similarity of himself, N * (N+1)/2 similarity value obtained;

(3) use N region as each node in network, use N * (N+1)/2 similarity value as the weights in network;

(4) the random N number that produces in [0,2 π], as the initial phase of each node in network, adopts oscillator formula to carry out phase place iteration to the node of the N in network, when all node N meet formula θ _i(k)-θ _i(k-3) during < 0.01, finishing iteration, θ _i(k) represent the phase place of the k time iteration of node i, iteration total degree is designated as m;

(5) [0,2 π] is divided into n ₁individual sub-range, statistics phase place is fallen the node number in each sub-range, with n ₁individual sub-range sequence number, as horizontal ordinate, is usingd and is dropped on sub-range interior nodes number as ordinate, draws PHASE DISTRIBUTION curve, statistics on curve between adjacent trough all nodes, be classified as a class, N node is divided into L altogether ₁class, L ₁numerical value identical with trough number, the classification sequence number of N node is stored in the vectorial F that length is N, if the arbitrary node i in N node belongs to j class, i in F value f _i=j, j=1 ..., L ₁;

(6) the m generation of use node and m-3, for the slope of all nodes of phase calculation, obtain N slope value, to L ₁each class in class is carried out following operation: [0.01,0.01] is divided into n ₂individual sub-range, statistics slope value is fallen the node number in each sub-range, with n ₂individual sub-range sequence number, as horizontal ordinate, is usingd and is dropped on sub-range node number as ordinate, draws L ₁bar slope distribution curve, adds up on each curve all nodes between adjacent trough, is classified as a class, and all nodes are divided into L altogether ₂class, L ₂with L ₁the trough sum of bar slope distribution curve is identical, upgrades the classification sequence number of each node in vectorial F;

(7) detect L ₂the nodes of each class in individual class, if the nodes of certain class is less than threshold value T=N/200, redistributes such all nodes, and the similarity of the arbitrary node i in such and class j is designated as to v _ij, the L of calculating i and all classes ₂individual similarity, to L ₂individual similarity sorts, and maximizing is given class corresponding to maximal value again by node i, upgrades the classification sequence number of node i in vectorial F;

(8) by the classification sequence number of vectorial F storage, i.e. classification sequence number under the region of pixel, mark is to each pixel;

(9) compose different gray-scale values to inhomogeneous pixel, the image after being cut apart, and output.

The present invention compared with prior art has the following advantages:

1, the present invention adopts the region after oscillator formula is cut apart watershed algorithm to carry out phase place iteration, and operand is little, and computing velocity is fast;

2, the present invention merges the region after cutting apart by extracting phase place and slope, has made up the over-segmentation defect of watershed algorithm, has improved the diagnosis and distinguish rate of medical image;

3, the present invention fully takes into account the complicacy of medical image, adopts the method for sequencing of similarity to reclassify sorting out wrong region, can accurately locate focus region.

Simulation result shows, mammogram figure is adopted to medical image processing method of the present invention, has improved the discrimination of image, contributes to the Accurate Diagnosis of medical image.

Accompanying drawing explanation

Fig. 1 is process flow diagram of the present invention;

Fig. 2 is the original mammogram figure that emulation of the present invention adopts;

Fig. 3 is that the present invention carries out the mammogram figure after watershed segmentation to Fig. 2;

Fig. 4 carries out the mammogram figure after the merging of region by the present invention and existing two kinds of methods to Fig. 3.

Embodiment

With reference to Fig. 1, the present invention is based on the synchronous medical image processing method of network phase, comprise the steps:

Step 1, extracts the feature of primitive medicine image pixel point and pixel feature is carried out to watershed segmentation.

(1a) input primitive medicine image, this example is selected original mammogram figure as shown in Figure 2, and size is 1024 * 1024;

(1b) extract the feature based on gray level co-occurrence matrixes, for any pixel i, extract contrast, consistance, energy [0,45,90,135] 12 dimensional features on four direction, then extract 10 dimensional features based on non-sampling wavelet transformation, above-mentioned 12 dimensional features and this two feature of 10 dimensional features are merged to composition 22 dimensional vectors, as the feature of i pixel;

(1c), to all pixel repeating steps (1b) in image, obtain the feature of all pixels of primitive medicine image;

(1d) the pixel feature of extracting is carried out to watershed segmentation, obtain N region, N >=1000, the mammogram figure after cutting apart, as shown in Figure 3.

Step 2, the similarity between zoning.

(2a) calculate the feature of All Ranges, obtain the feature x in N region ₁... x _n, wherein, the feature x in i region _iequal the mean value of all pixel features in i region;

(2b) use the feature x in N region ₁... x _n, calculate the similarity between every two regions, and each region and the similarity of himself, the computing formula of similarity is:

$w_{\mathrm{ij}}=e^{-\frac{{(x_i-x_j)}^2}{{2\mathrm{\&sigma;}}^2}},i=1,...N,j=1,...N,$

Wherein, w _ijrepresent the similarity between i region and j region, x _ithe feature in i region, x _jbe the feature in j region, σ is scale parameter, and span is (0,1).

Step 3, as each node in network, as the weights in network, sets up network by N * (N+1)/2 similarity value with N region.

Step 4, adopts oscillator phase synchronization formula to carry out phase place iteration.

(4a) arithmetic mean of all weights N * (N+1)/2 in computational grid

generator matrix E, the value e of the capable j row of i in this matrix E _ijfor:

$\left\{\begin{array}{c}{e}_{\mathrm{ij}}=1,w_{\mathrm{ij}}\&GreaterEqual;\stackrel{\&OverBar;}{w}\\ e_{\mathrm{ij}}=0,w_{\mathrm{ij}}<\stackrel{\&OverBar;}{w}\end{array}\right.,i=1,...N,j=1,...N;$

(4b) in [0,2 π], produce the random number of N, as the initial phase of N node, be designated as θ (0), θ (0)={ θ ₁(0) ... θ _i(0) ..., θ _n(0) }, wherein, θ _i(0) be the initial phase of i node;

(4c) according to the initial phase θ of i node _i(0), obtaining oscillator phase synchronization formula is:

${\mathrm{\&theta;}}_i(k+1)={\mathrm{\&theta;}}_i\left(k\right)+\underset{j=1}{\overset{N}{\mathrm{\&Sigma;}}}\frac{\mathrm{\&Gamma;}}{{\mathrm{\&Gamma;}}_{\max }}\sin ({\mathrm{\&theta;}}_j\left(k\right)-{\mathrm{\&theta;}}_i\left(k\right)),i=1,...N,$

Wherein, the nodes that N is network, k is iterations, θ _i(k) be the phase place of the k time iteration of node i, Γ is stiffness of coupling, works as e _ij=1 o'clock, Γ=3, worked as e _ij=0 o'clock, Γ=-1,

(4d), from k=0, by the above-mentioned oscillator phase synchronization formula of initial phase θ (0) substitution of N node, N node carried out to phase place iteration;

(4e) when iterations k>3, it is stable whether decision node reaches, if all N node all meets formula θ _i(k)-θ _i(k-3) < 0.01, and node phase place reaches stable, finishing iteration, and iteration total degree is designated as m;

(4f) when iterations k≤3, or be not that all nodes all meet formula θ _i(k)-θ _i(k-3) < 0.01, the k of N node for phase theta ₁(k) ... θ _n(k) in the above-mentioned oscillator phase synchronization iterative formula of substitution, calculate the k+1 of N node for phase place, and iterations k is added to 1, return to step (4e).

Step 5, classifies to node according to phase place.

[0,2 π] is divided into n ₁individual sub-range, statistics phase place is fallen the node number in each sub-range, with n ₁individual sub-range sequence number, as horizontal ordinate, is usingd and is dropped on sub-range interior nodes number as ordinate, draws PHASE DISTRIBUTION curve, and statistics all nodes between adjacent trough on curve, are classified as a class, and N node is divided into L altogether ₁class, L ₁numerical value identical with trough number, the classification sequence number of N node is stored in the vectorial F that length is N, if the arbitrary node i in N node belongs to j class, i in F value f _i=j, j=1 ... L ₁.

Step 6, classifies to node according to slope.

(6a) slope S for all node N of phase calculation with m generation of node and m-3:

S={s ₁,...s _i,...s _N}，

Wherein, the slope s of i node _ifor:

$s_i=\frac{{\mathrm{\&theta;}}_i\left(m\right)-{\mathrm{\&theta;}}_i(m-3)}{3},i=1,...N;$

(6b) to L ₁each class in class is carried out following operation: [0.01,0.01] is divided into n ₂individual sub-range, statistics slope value is fallen the node number in each sub-range, with n ₂individual sub-range sequence number, as horizontal ordinate, is usingd and is dropped on sub-range node number as ordinate, draws L ₁bar slope distribution curve, adds up on each curve all nodes between adjacent trough, is classified as a class, and all nodes are divided into L altogether ₂class, L ₂with L ₁the trough sum of bar slope distribution curve is identical, upgrades the classification sequence number of each node in vectorial F.

Step 7, redistribute may misclassification node.

(7a) detect L ₂the nodes of any class in individual class, if such nodes is less than threshold value T=N/200, calculates the L of i node and all classes in such ₂individual similarity v _ij:

$v_{\mathrm{ij}}=\underset{n=1,f_n=j}{\overset{N}{\mathrm{\&Sigma;}}}{w}_{\mathrm{in}},j=1,...L_2;$

(7b) to L ₂individual similarity sorts, and maximizing is given class corresponding to maximal value again by i node, and upgrades the node classification sequence number in vectorial F;

(7c) to all L ₂individual class repeating step (7a)-(7b), is less than all nodes node in the class of threshold value T and reclassifies.

Step 8, by the classification sequence number of vectorial F storage, i.e. classification sequence number under the region of pixel, mark is to each pixel.

Step 9, composes different gray-scale values to inhomogeneous pixel, the image after being cut apart, and output.

(9a) calculate the gray-scale value h of i class _i:

$h_i=255\&times;\frac{L_2-i+1}{L_2},i=1,...L_2,$

By h _ibe assigned to all pixels of i class;

(9b) make i=i+1, repeating step (9a), until i=L ₂, L now ₂all assignment is complete for all pixels of individual class, and the gray level image that all pixels after assignment are formed is as the final image of cutting apart;

(9c) export the final image of cutting apart.

Effect of the present invention can be by further illustrating the emulation experiment of mammogram figure below:

1, simulated conditions

Emulation of the present invention is at windowsXP, SPI, and CPUInterCore2Duo, basic frequency 2.33Ghz, software platform is Matlab2007b operation.The original mammogram that emulation is selected derives from common data sets MIAS, obtains altogether 330 original mammograms, and wherein Fig. 2 is a width canceration mammogram wherein.

2, emulation content and result

Emulation 1, carries out watershed segmentation with the present invention to the mammogram of canceration shown in Fig. 2, the mammogram after obtaining Fig. 3 and cutting apart.

Emulation 2, by the present invention and existing Kmeas method and existing Fcm method, the image after Fig. 3 is cut apart carries out region classification, and result is as Fig. 4.Wherein Fig. 4 a adopts network phase synchronous method of the present invention Fig. 3 to be carried out to the result of region classification; Fig. 4 b adopts existing Kmeas method Fig. 3 to be carried out to the result of region classification; Fig. 4 c adopts existing Fcm method Fig. 3 to be carried out to the result of region classification.

From Fig. 4 a, can find out, the present invention can detect the lesion region of mammary gland exactly, and compares with Fig. 4 c with Fig. 4 b, and the detected lesion region of the present invention is without assorted point, and region consistance is higher, and aspect the local detail processing to tumor region, advantage is comparatively obvious.

Above-mentioned simulation result shows, the present invention can improve compared to traditional method the order of accuarcy that image is cut apart, and obtains better segmentation result, is a kind of effective medical image processing method.

Above example is implemented take technical solution of the present invention under prerequisite, provided detailed embodiment and concrete operating process, but protection scope of the present invention is not limited to above-described embodiment.

Claims (6)

1. a phase locked medical image processing method Network Based, comprises the steps:

(1) extract the pixel feature of primitive medicine image, and the pixel feature of extracting is carried out to watershed segmentation, obtain N region, N >=1000;

(2) calculate the similarity between every two regions, and each region and the similarity of himself, N * (N+1)/2 similarity value obtained;

(3) use N region as each node in network, use N * (N+1)/2 similarity value as the weights in network;

(4) the random N number that produces in [0,2 π], as the initial phase of each node in network, adopts oscillator formula to carry out phase place iteration to the node of the N in network, when all node N meet formula θ _i(k)-θ _i(k-3) during < 0.01, finishing iteration, θ _i(k) represent the phase place of the k time iteration of node i, iteration total degree is designated as m;

(5) [0,2 π] is divided into n ₁individual sub-range, statistics phase place is fallen the node number in each sub-range, with n ₁individual sub-range sequence number, as horizontal ordinate, is usingd and is dropped on sub-range interior nodes number as ordinate, draws PHASE DISTRIBUTION curve, statistics on curve between adjacent trough all nodes, be classified as a class, N node is divided into L altogether ₁class, L ₁numerical value identical with trough number, the classification sequence number of N node is stored in the vectorial F that length is N, if the arbitrary node i in N node belongs to j class, i in F value f _i=j, j=1 ..., L ₁;

(6) the m generation of use node and m-3, for the slope of all nodes of phase calculation, obtain N slope value, to L ₁each class in class is carried out following operation: [0.01,0.01] is divided into n ₂individual sub-range, statistics slope value is fallen the node number in each sub-range, with n ₂individual sub-range sequence number, as horizontal ordinate, is usingd and is dropped on sub-range node number as ordinate, draws L ₁bar slope distribution curve, adds up on each curve all nodes between adjacent trough, is classified as a class, and all nodes are divided into L altogether ₂class, L ₂with L ₁the trough sum of bar slope distribution curve is identical, upgrades the classification sequence number of each node in vectorial F;

(7) detect L ₂the nodes of each class in individual class, if the nodes of certain class is less than threshold value T=N/200, redistributes such all nodes, and the similarity of the arbitrary node i in such and class j is designated as to v _ij, the L of calculating i and all classes ₂individual similarity, to L ₂individual similarity sorts, and maximizing is given class corresponding to maximal value again by node i, upgrades the classification sequence number of node i in vectorial F;

(8) by the classification sequence number of vectorial F storage, i.e. classification sequence number under the region of pixel, mark is to each pixel;

(9) compose different gray-scale values to inhomogeneous pixel, the image after being cut apart, and output.

2. according to the Network Based phase locked medical image processing method described in claims 1, the pixel feature of the described extraction primitive medicine image of step (1) wherein, according to following steps, carry out:

(2a) extract the feature based on gray level co-occurrence matrixes, for any pixel i, extract contrast, consistance, 12 dimensional features of energy on [0,45,90,135] four direction;

(2b), for any pixel i, extract 10 dimensional features based on non-sampling wavelet transformation;

(2c) above-mentioned 12 dimensional features and this two feature of 10 dimensional features are merged to composition 22 dimensional vectors, as the feature of i pixel;

(2d) to all pixel repeating steps (2a) in image-(2c), obtain the feature of all pixels of primitive medicine image.

3. according to the Network Based phase locked medical image processing method described in claims 1, the similarity between the described zoning of step (2) wherein, carry out as follows:

(3a) calculate the feature of All Ranges, obtain the feature x in N region ₁... x _n, wherein, the feature x in i region _iequal the mean value of all pixel features in i region;

(3b) use the feature x in N region ₁... x _n, calculate the similarity between every two regions, and each region and the similarity of himself, the computing formula of similarity is:

Wherein, w _ijrepresent the similarity between i region and j region, x _ithe feature in i region, x _jbe the feature in j region, σ is scale parameter, and span is (0,1).

4. according to the Network Based phase locked medical image processing method described in claims 1, wherein the described employing oscillator phase synchronization formula of step (4) carries out phase place iteration to the node of the N in network, carries out as follows:

(4a) arithmetic mean of all weights N * (N+1)/2 in computational grid

generator matrix E, the value e of the capable j row of i in this matrix E _ijfor:

(4b) in [0,2 π], produce the random number of N, as the initial phase of N node, be designated as θ (0), θ (0)={ θ ₁(0) ... θ _i(0) ..., θ _n(0) }, wherein, θ _i(0) be the initial phase of i node;

(4c) according to the initial phase θ of i node _i(0), obtaining oscillator phase synchronization iterative formula is:

Wherein, the nodes that N is network, k is iterations, mistake! Do not find Reference source.θ _i(k) be the phase place of the k time iteration of i node, Γ is stiffness of coupling, works as e _ij=1 o'clock, Γ=3, worked as e _ij=0 o'clock, Γ=-1,

(4d), from k=0, by the above-mentioned oscillator phase synchronization iterative formula of initial phase θ (0) substitution of N node, N node carried out to phase place iteration.

5. according to the Network Based phase locked medical image processing method described in claims 1, the slope of the described all nodes of calculating of step (6) wherein, carries out according to following formula:

Wherein, s _iit is the slope of i node.

6. according to the Network Based phase locked medical image processing method described in claims 1, the wherein described computing node i of step (7) and the L of all classes ₂individual similarity, according to following formula, carry out:

Wherein, v _ijit is the similarity of node i and class j.

Images