CN103593842B - Based on intersect climb the mountain memetic quantum evolution calculate medical image registration method - Google Patents

Abstract

The invention discloses the medical image registration method of a kind of memetic quantum evolution calculating of climbing the mountain based on intersection, the problem mainly solving prior art registration weak effect.Implementation step is: read in image, initial phase related parameter；Initial population is produced by chaotic method；Calculate each individual corresponding registration parameter in population；According to registration parameter, floating image is carried out image conversion, the image after being converted；Calculate the similarity between image and the reference picture after conversion, and find this generation optimum registration parameter, according to optimum registration parameter, individualities all in population are carried out Quantum rotating gate renewal；Optimum registration parameter carries out intersection climb the mountain the study of memetic local；Judge whether optimum individual meets and forget condition, perform to forget operator to optimum individual；Control cycling condition, if loop ends, then export registration result.The present invention can obtain more preferable registration result, can be used for the registration of medical image.

Description

Based on intersect climb the mountain memetic quantum evolution calculate medical image registration method

Technical field

The invention belongs to technical field of image processing, relate to medical figure registration, a kind of new climbing based on intersection The medical image registration method that mountain memetic quantum evolution calculates, can be used for follow-up medical image and understands and process.

Background technology

Medical figure registration is a kind of method finding spatial transformation parameter by optimized algorithm, different to two width by this parameter Medical image does spatial alternation, this two width medical image can be made to reach one_to_one corresponding on locus, or have hands Coupling is reached on the point of art meaning.Participate in two width images of registration can be different imaging mode, different patient's or The medical image that same patient shoots at different time.

Medical image registration method generally comprises feature space, search volume, similarity measure and searching algorithm.Wherein:

Feature space is the characteristics of image that method for registering images is used, including gray scale, flex point, boundary line etc.；Search volume It is mode and the scope that in method for registering images, image is carried out deformation, including rigid transformation, non-rigid transformation etc.；Similarity Estimate the quantized value being tolerance registration image with the direct similarity of reference picture, including correlation coefficient, mutual information, normalizing Change mutual information etc.；Searching algorithm is to find the optimized algorithm that registration parameter is used in method for registering images, declines including steepest Method, genetic algorithm etc..According to the difference of the searching algorithm that medical image registration method uses, can be classified as based on tradition The method for registering optimized and method for registering based on intelligent optimization.Different optimized algorithms has for the result of whole method for registering Vital impact.Existing method for registering mainly uses tradition optimized algorithm, and this kind of method utilizes the mathematics of object function Characteristic constitution optimization method, although speed, but registration degree is not high enough, is easily absorbed in local extremum, especially for The non-rigid registration that parameters optimization is more.The most also occur in that medical figure registration based on genetic algorithm, but genetic algorithm is compiled Code length is longer, and still can converge to local extremum, impact registration effect.

Summary of the invention

Present invention aims to the deficiency of existing method for registering, propose a kind of new climb the mountain memetic amount based on intersection The medical image registration method of sub-evolutionary computation, converges to local extremum to reduce, and improves registration effect.

For achieving the above object, technical solution of the present invention comprises the steps:

(1) two images subject to registration are read in, by a wherein changeless image as reference picture I_R, another width is carried out The image of deformation is as floating image I_F, and by spacing L=35 pixel, to floating image uniform sampling, obtain by m Elementary composition uniform grid；

(2) producing one group of individuality by chaotic method, as population, the most each individuality is the vector of a 1 × m；

(3) individuality each in population is converted into candidate solution, then superposes with uniform grid, obtain the control corresponding with each individuality Grid processed, as registration parameter；

(4) according to registration parameter, floating image is done cubic B-spline function to convert, and calculate the floating image after conversion and ginseng Examine the normalized mutual information value between image, find out the registration ginseng of normalized mutual information value maximum in current iteration and correspondence thereof Several and individual, and this individuality is stored as this generation optimum individual p_best；

(5) according to this generation optimum individual p_best, utilize Quantum rotating gate operator that individualities all in this generation population are carried out quantum rotation Revolving door updates, and makes image registration parameter move towards the direction that similarity is bigger；

(6) any one in addition to optimum individual is selected individual, and with optimum individual together as intersecting the memetic that climbs the mountain The locally initial point of learning operator, carries out local and learns, optimize the optimum individual p of current iteration further_best；

(7) optimum individual p is judged_bestWhether change, if p after continuous F iteration_bestDo not change, then optimum Individuality meets forgets condition, and employing is forgotten operator and forgotten by optimum individual, it is to avoid registration parameter is absorbed in local extremum, otherwise, Performing step (8), F is the positive integer more than 2 set；

(8) update iterations g, as iterations reaches maximum, then terminate registration, the Controling network that optimum individual is corresponding Lattice are registration parameter, by this parameter, floating image are done cubic B-spline function conversion, obtain registrating image, otherwise, return Return step (3) and carry out next iteration.

The present invention has the advantage that compared with prior art

1. present invention adds and forget operator, searching algorithm can be prevented effectively from and be absorbed in local extremum, improve registration image and reference The similarity of image, improves registration effect；

2. the present invention adds the local study stage after quantum evolution calculates, calculated to quantum evolution in an iterative process Optimum individual carries out a local optimum again, can global search be combined with Local Search well, makes searching algorithm more Effectively, registration effect is improved.

Accompanying drawing explanation

Fig. 1 is the general flow chart of the present invention；

Fig. 2 is the reference picture of emulation experiment of the present invention；

Fig. 3 is the floating image of emulation experiment of the present invention.

Detailed description of the invention

With reference to Fig. 1, the present invention to implement step as follows:

Step 1, reads in image, initial phase related parameter.

Reading in two images subject to registration, piece image immobilizes, and another width is image to be converted；To wherein immobilize Image as reference picture I_R, another width is carried out the image that converts as floating image I_F, and by spacing L=35 picture The number of element, carries out uniform sampling to floating image, obtains the uniform grid elementary composition by m；

Set maximum iteration time max_g=100, set local study number of times X=10 time, and current iteration number of times g is put Zero.

Step 2, produces initial population by chaotic method.

(2a) a random value θ is generated₁, 0 ＜ θ₁＜ 1, θ₁≠ 0.5, pass through θ₁Calculate whole chaotic individual P=[θ₁,θ₂,...,θ_m], wherein:

θ_i+1=4 θ_i·(1-θ_i),

Wherein, i=1,2 ... m-1, m are the element numbers of uniform grid；

(2b) repeat step (2a) and produce one group of individuality P=[p₁,p₂,...,p_N], composition number is the population of N.

Step 3, calculates each individual corresponding control grid in population.

(3a) selected population P=[p₁,p₂,...,p_NBody p one by one in]_n, n=1,2 ... N；

(3b) by individuality p_nIt is mapped to registration parameter solution space, obtains solution S with m component_n=[S_n1,S_n2,...,S_nm], Each weight expression is as follows:

S_ni=a_i·θ_ni+b_i·(1-θ_ni), i=1,2 ... m,

Wherein, S_niIt is S_nI-th component, a_iAnd b_iIt is value lower limit and the upper limit of i-th component respectively；

(3c) by above-mentioned S_nIt is integrated into and the matrix of uniform grid formed objects, and superposes with uniform grid and obtain individual p_nRight The control grid Φ answered_n；

(3d) repeat step (3a)-(3c), obtain population P=[p₁,p₂,...,p_N] corresponding control grid Φ=[Φ₁,Φ₂,...,Φ_N]。

Step 4, does cubic B-spline function conversion, the image after being converted according to controlling grid, and counts floating image Calculate the similarity of the image after conversion and floating image.

(4a) control grid Φ=[Φ is taken₁,Φ₂,...,Φ_NA grid Φ in]_n, n=1,2 ... N；

(4b) take one pixel of floating image, be designated as (x, y), according to control grid Φ_nCalculate the shift value of this pixel D (x, y), it may be assumed that

Wherein, WithRepresent respectively and x and y is taken It is whole,It is to control grid Φ_nIn this pixel (x, y) value of around 4 × 4 Controling network lattice points, B_q(e) and B_lF () is Cubic B-spline basic function B_dOne in (t), d=0,1,2,3,0≤t ＜ 1, cubic B-spline basic function B_dT () represents such as Under:

$\left\{\begin{array}{c}{B}_0\left(t\right)=\frac{{(1-t)}^3}{6}\\ B_1\left(t\right)=\frac{{3t}^3-{6t}^3+4}{6}\\ B_2\left(t\right)=\frac{{-3t}^3+{3t}^3+3t+1}{6}\\ B_3\left(t\right)=\frac{t^3}{6}\end{array}\right.;$

(4c) repeat step (4b) and calculate the shift value of each pixel of floating image, obtain Φ_nThe corresponding image after conversion I_n；

(4d) reference picture I is calculated_RWith the image I after conversion_nBetween normalized mutual information value NMI_nFor:

${\mathrm{NMI}}_n=\frac{H\left(I_R\right)+H\left(I_n\right)}{H(I_R,I_n)},$

Wherein, H (I_R) it is reference picture I_REntropy, H (I_n) it is the image I after converting_nEntropy, H (I_R,I_n) it is with reference to figure As I_RWith the floating image I after conversion_nCombination entropy；

(4e) repeat step (4a)-(4d) and obtain population P=[p₁,p₂,...,p_N] corresponding normalized mutual information value NMI=[NMI₁,NMI₂,...,NMI_N]；

(4f) NMI=[NMI is found out₁,NMI₂,...,NMI_NMaximum in] and the individuality of correspondence thereof, and this individuality is stored Optimum individual p for current iteration_best, by p_bestThe corresponding grid that controls saves as the optimum registration parameter Φ of current iteration_best, By p_bestCorresponding normalized mutual information value is stored as maximum normalized mutual information NMI_best。

Step 5, according to this generation optimum individual, individualities all to this generation carry out Quantum rotating gate renewal.

(5a) according to size delta θ of the iterations g structure anglec of rotation₀:

Δθ₀=0.015 π × (g/max_g)+0.15 π × exp (-mod (g, 100)/10),

Wherein, max_g is the maximum iteration time set, and mod (g, 100) is to take the remainder that g is divided by 100, Exp (-mod (g, 100)/10) is to take e-mod (g, 100)/10 powers；

(5b) g is taken out for p individual in population；

(5c) the i-th component of individual p is taken outIf this component has crossed the i-th component of optimum individual, thenOtherwise,The i-th component of the individual p in g+1 generationMeet:

${\mathrm{\θ}}_i^{g+1}=\arccos (\cos \left(\mathrm{\Δ}{\mathrm{\θ}}_i^g\right)\·\cos \left({\mathrm{\θ}}_i^g\right)+(-\sin \left(\mathrm{\Δ}{\mathrm{\θ}}_i^g\right))\·\sin \left({\mathrm{\θ}}_i^g\right));$

(5d) step (5c) is repeated until the important update all of individual p；

(5e) step (5b)-(5d) is repeated until all individual update alls in population.

Step 6, selects this generation optimum individual and the arbitrary individuality in addition to optimum individual, carries out intersection and climbs the mountain memetic locally Study.

(6a) any one individual p in addition to optimum individual in population is selected_r；

(6b) described individual p is generated with BLX crossover operator_rWith optimum individual p_bestOne group of offspring individual Off=[off₁,off₂,...,off_R], R is offspring individual sum, wherein kth offspring individual off_kFor:

off_k=(1-ξ_k)·p_r+ξ_k·p_best,

Wherein, k=1,2 ... R, ξ_k=U [-α, 1+ α], U is for being uniformly distributed, and α is adjustable parameter, 0 ＜ α ＜ 1；

(6c) filial generation off=[off is calculated according to step (3) with the method described in step (4)₁,off₂,...,off_R] corresponding normalization is mutual Value of information NMI_off=[NMI_off1,NMI_off2,...,NMI_offR], find out NMI_offIn maximum NMI_offmax, and NMI_offmaxCorresponding optimum offspring individual off_best；

(6d) by optimum offspring individual off_bestCorresponding normalized mutual information value NMI_offmaxWith described individual p_rCorresponding normalizing Change association relationship NMI_rCompare, if NMI_offmax＞ NMI_r, then described individual p is made_rIt is optimum offspring individual off_best, i.e. p_r=off_best, and make described individual p_rCorresponding normalized mutual information value NMI_rIt is optimum offspring individual off_bestCorresponding returns One changes association relationship NMI_offmax, i.e. NMI_r=NMI_offmax, perform step (6f), otherwise, perform step (6e)；

(6e) by optimum offspring individual off_bestCorresponding normalized mutual information value NMI_offmaxWith optimum individual p_bestCorresponding returns One changes association relationship NMI_bestCompare, if NMI_offbest＞ NMI_best, then optimum individual p is made_bestIt it is optimum offspring individual off_best, i.e. p_best=off_best, and make optimum individual p_bestCorresponding normalized mutual information value NMI_bestIt it is optimum offspring individual off_bestCorresponding normalized mutual information value NMI_offmax, i.e. NMI_best=NMI_offmax, otherwise, perform step (6f)；

(6f) judge whether local study number of times reaches X=10 time, if reaching, then stop local study, otherwise, return step (6b)。

Step 7, it may be judged whether meet and forget condition, performs to forget operator.

(7a) setting positive integer F more than 2, this example selects F=10；

(7b) optimum individual p is judged_bestWhether change, if p after continuous F iteration_bestDo not change, then optimum Individual p_bestMeet and forget condition, perform step (7c), otherwise, perform step (8)；

(7c) optimum individual p is deleted_best, and by p_bestCorresponding maximum normalized mutual information NMI_bestZero setting.

Step 8, updates iterations, it may be judged whether reach termination condition.

Update iterations g=g+1, if iterations reaches maximum max_g=100, then terminate registration, by optimum Body p_bestCorresponding control grid Φ_bestAs registration parameter, by this parameter, floating image is done cubic B-spline function and converts, Image after being registrated, maximum normalized mutual information value NMI of output_best, otherwise, return step (3) and carry out next iteration.

The effect of the present invention can be further illustrated by following experiment emulation:

1, experiment condition and method

Hardware platform is: Intel (R) Core (TM) i5-2450M@2.50GHz, 3.91GBRAM.；

Software platform is: MATLAB R2012b；

Experimental technique: respectively by method based on steepest descent method, method based on genetic algorithm and the Realization of Simulation of the present invention Medical figure registration.

2, emulation content

As shown in Figures 2 and 3, Fig. 2 is the brain MRI image of normal person to image used by emulation experiment, and Fig. 3 is patient Brain MRI image.

In order to verify effectiveness of the invention, respectively by method based on steepest descent method, method based on genetic algorithm and Experimental image shown in Fig. 2 and Fig. 3 is registrated by the present invention.Every kind of method is independently emulated 100 times, obtains normalization The meansigma methods of association relationship and two results of the maximum of normalized mutual information value.

3, simulation result

Two results that emulation obtains are as shown in table 1, and wherein method based on steepest descent method is method 1, calculate based on heredity The method of method is method 2.

Table 1 three all method for registering Comparative result table

	Method 1	Method 2	The present invention
				Average normalized mutual information	1.1422	1.1386	1.1490
Maximum normalized mutual information	1.1422	1.1429	1.1544

From table 1, the inventive method compared to the method for method based on steepest descent method and genetic algorithm at average normalizing Changing mutual information and maximum two aspects of normalized mutual information are the most outstanding, improve registration effect, algorithm performance is more excellent, is A kind of effective medical image registration method.

Claims (3)

1. based on intersect climb the mountain memetic quantum evolution calculate a medical image registration method, comprise the steps:

(1) two images subject to registration are read in, by a wherein changeless image as reference picture I_R, another width is carried out The image of deformation is as floating image I_F, and by spacing L=35 pixel, to floating image uniform sampling, obtain by m Elementary composition uniform grid；

(2) producing one group of individuality by chaotic method, as population, the most each individuality is the vector of a 1 × m；

(3) individuality each in population is converted into candidate solution, then superposes with uniform grid, obtain the control corresponding with each individuality Grid processed, as registration parameter；

Described individuality each in population is converted into candidate solution, is by individuality p=[θ₁,θ₂,...,θ_m] it is mapped to registration parameter solution sky Between, obtain the solution S=[S with m component₁,S₂,...,S_m], each weight expression is as follows:

S_i=a_i·θ_i+b_i·(1-θ_i), i=1,2 ... m,

Wherein, S_iIt is the i-th component of S, a_iAnd b_iIt is value lower limit and the upper limit of i-th component respectively；

(4) according to registration parameter, floating image is done cubic B-spline function to convert, and calculate the floating image after conversion and ginseng Examine the normalized mutual information value between image, find out normalized mutual information value maximum in current iteration and the registration of correspondence thereof Parameter and individuality, and this individuality is stored as this generation optimum individual p_best；

(5) according to this generation optimum individual p_best, utilize Quantum rotating gate operator that individualities all in this generation population are carried out quantum rotation Revolving door updates, and makes image registration parameter move towards the direction that similarity is bigger:

(5a) according to size delta θ of the iterations g structure anglec of rotation₀:

Δθ₀=0.015 π × (g/max_g)+0.15 π × exp (-mod (g, 100)/10),

Wherein, max_g is the maximum iteration time set, and mod (g, 100) is to take the remainder that g is divided by 100, Exp (-mod (g, 100)/10) is to take e-mod (g, 100)/10 powers；

(5b) g i-th component for p individual in population is taken outJudge whether this component crosses the i-th of optimum individual Individual component, if crossing, then the anglec of rotationOtherwise, the anglec of rotationThen g+1 generation The i-th component of individual pMeet:

${\mathrm{\θ}}_i^{g+1}=\arccos \left(\cos \left({\mathrm{\Δ\θ}}_i^g\right)\·\cos \left({\mathrm{\θ}}_i^g\right)+\left(-\sin \left({\mathrm{\Δ\θ}}_i^g\right)\right)\·\sin \left({\mathrm{\θ}}_i^g\right)\right);$

(5c) step (5b) is repeated until the important update all of individual p；

(6) any one in addition to optimum individual is selected individual, and with optimum individual together as intersecting the memetic that climbs the mountain The locally initial point of learning operator, carries out local and learns, optimize the optimum individual p of current iteration further_best；

(7) optimum individual p is judged_bestWhether change, if p after continuous F iteration_bestDo not change, then optimum Individuality meets forgets condition, and employing is forgotten operator and forgotten by optimum individual, it is to avoid registration parameter is absorbed in local extremum, no Then, performing step (8), F is the positive integer more than 2 set；

(8) update iterations g, as iterations reaches maximum, then terminate registration, the Controling network that optimum individual is corresponding Lattice are registration parameter, by this parameter, floating image do cubic B-spline function conversion, obtain registrating image, otherwise, Return step (3) and carry out next iteration.

2. according to the medical figure registration side of the memetic quantum evolution calculating of climbing the mountain based on intersection described in claims 1 Method, wherein described in step (2) with chaotic method produce one group individual, carry out as follows:

(2a) a random value θ is first generated₁, 0 < θ₁< 1, θ₁≠ 0.5, then pass through θ₁Calculate whole chaotic individual P=[θ₁,θ₂,...,θ_m], wherein:

θ_i+1=4 θ_i·(1-θ_i),

Wherein, i=1,2 ... m-1, m are the element numbers of uniform grid；

(2b) repeat step (2a) and produce one group of individuality P=[p₁,p₂,...,p_N], composition number is the population of N.

3. according to the medical figure registration side of the memetic quantum evolution calculating of climbing the mountain based on intersection described in claims 1 Method, wherein any one individuality selected in addition to optimum individual described in step (6), and with optimum individual together as friendship Pitch the initial point of memetic local learning operator of climbing the mountain, carry out local and learn, be to carry out as follows:

(6a) any one individual p in addition to optimum individual in population is selected_r, set local study number of times X=10 time；

(6b) described individual p is generated with BLX crossover operator_rWith optimum individual p_bestOne group of offspring individual Off=[off₁,off₂,...,off_R], R is offspring individual sum, wherein kth offspring individual off_kFor:

off_k=(1-ξ_k)·p_r+ξ_k·p_best,

Wherein, k=1,2 ... R, ξ_k=U [-α, 1+ α], U is for being uniformly distributed, and α is adjustable parameter, and 0 < α < 1；

(6c) filial generation off=[off is calculated according to step (3) with the method described in step (4)₁,off₂,...,off_R] corresponding normalization Association relationshipFind out NMI_offIn maximum NMI_offmax, and NMI_offmaxCorresponding optimum offspring individual off_best；

(6d) by optimum offspring individual off_bestCorresponding normalized mutual information value NMI_offmaxWith described individual p_rCorresponding returns One changes association relationship NMI_rCompare, if NMI_offmax>NMI_r, then described individual p is made_rIt is optimum offspring individual off_best, I.e. p_r=off_best, and make described individual p_rCorresponding normalized mutual information value NMI_rIt is optimum offspring individual off_bestCorresponding Normalized mutual information value NMI_offmax, i.e. NMI_r=NMI_offmax, perform step (6f), otherwise, perform step (6e)；

(6e) by optimum offspring individual off_bestCorresponding normalized mutual information value NMI_offmaxWith optimum individual p_bestCorresponding Normalized mutual information value NMI_bestCompare, if NMI_offbest>NMI_best, then optimum individual p is made_bestIt it is optimum offspring individual off_best, i.e. p_best=off_best, and make optimum individual p_bestCorresponding normalized mutual information value NMI_bestIt it is optimum offspring individual off_bestCorresponding normalized mutual information value NMI_offmax, i.e. NMI_best=NMI_offmax, otherwise, perform step (6f)；

(6f) judge whether local study number of times reaches X=10 time, if reaching, then stop local study, otherwise, return step Suddenly (6b).