CN117649416A - Robust chest CT image registration method - Google Patents
Robust chest CT image registration method Download PDFInfo
- Publication number
- CN117649416A CN117649416A CN202410124244.4A CN202410124244A CN117649416A CN 117649416 A CN117649416 A CN 117649416A CN 202410124244 A CN202410124244 A CN 202410124244A CN 117649416 A CN117649416 A CN 117649416A
- Authority
- CN
- China
- Prior art keywords
- image
- mask
- lung
- bone
- spine
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 53
- 238000006073 displacement reaction Methods 0.000 claims abstract description 95
- 238000005070 sampling Methods 0.000 claims abstract description 3
- 210000004072 lung Anatomy 0.000 claims description 95
- 210000000988 bone and bone Anatomy 0.000 claims description 87
- 230000006870 function Effects 0.000 claims description 63
- 230000008569 process Effects 0.000 claims description 15
- 238000006467 substitution reaction Methods 0.000 claims description 7
- 230000000694 effects Effects 0.000 claims description 5
- 230000009466 transformation Effects 0.000 claims description 4
- 238000013136 deep learning model Methods 0.000 claims description 3
- 230000000877 morphologic effect Effects 0.000 claims description 3
- 230000011218 segmentation Effects 0.000 claims description 3
- 238000012804 iterative process Methods 0.000 claims description 2
- 230000035945 sensitivity Effects 0.000 abstract description 2
- 238000005457 optimization Methods 0.000 description 12
- 230000033001 locomotion Effects 0.000 description 5
- 238000013135 deep learning Methods 0.000 description 4
- 210000000056 organ Anatomy 0.000 description 4
- 210000003484 anatomy Anatomy 0.000 description 3
- 210000004204 blood vessel Anatomy 0.000 description 2
- 230000007423 decrease Effects 0.000 description 2
- 238000013461 design Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 238000003384 imaging method Methods 0.000 description 2
- 230000004807 localization Effects 0.000 description 2
- 238000012986 modification Methods 0.000 description 2
- 230000001902 propagating effect Effects 0.000 description 2
- 206010028980 Neoplasm Diseases 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000008859 change Effects 0.000 description 1
- 238000013170 computed tomography imaging Methods 0.000 description 1
- 238000007408 cone-beam computed tomography Methods 0.000 description 1
- 238000002786 image-guided radiation therapy Methods 0.000 description 1
- 238000002675 image-guided surgery Methods 0.000 description 1
- 238000002513 implantation Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 239000003550 marker Substances 0.000 description 1
- 230000007246 mechanism Effects 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
- 238000012887 quadratic function Methods 0.000 description 1
- 210000003437 trachea Anatomy 0.000 description 1
Abstract
The invention discloses a robust chest CT image registration method, which is characterized in that different dissimilarity measures and regularization are constructed for each region in a divide-and-conquer way, gaussian pyramids are embedded into different resolutions in a sampling mode, then a control lattice is used for controlling a displacement field, a smooth Welsch function is used for constructing an objective function, an alternative function is used for accelerating iteration efficiency, and parameters are dynamically adjusted for each iteration to reduce sensitivity to noise and increase robustness. The invention can effectively solve the problem of low registration precision of the high-noise CT image.
Description
Technical Field
The invention relates to the field of medical image processing, in particular to a robust chest CT image registration method.
Background
CT image registration is the basis for a variety of clinical applications such as image guided surgery, surgical planning, tumor monitoring, motion tracking, and the like. Chest CT images are often subject to artifacts and noise due to limitations in imaging equipment, imaging dose, and scanning conditions. For example, cone beam CT imaging produces high noise due to low dose, and implantation of metallic objects in a patient can cause significant star or streak artifacts. In many cases, registration of high noise CT images is required, as in image-guided radiation therapy, where pre-processed low noise CT images are registered to high noise real-time CBCT images acquired at a fast, low dose to clearly show the target anatomy. These noises present great difficulties for chest CT image registration.
Registration methods can generally be divided into two categories: optimization-based and deep learning-based methods. Optimization-based methods include non-parametric methods that directly optimize the displacement vector for each voxel, and parametric methods that use a small set of parameters to control the entire displacement field. With the progress of deep learning, recent registration studies are focusing more and more on a deep learning-based method. Although the deep learning based registration method improves efficiency, there is room for improvement in terms of accuracy, generalization, and robustness over the optimization based method.
Current chest CT image registration methods are numerous, some of which design specific dissimilarity measures to guide registration using more representative features, however they are global and do not fully account for the various regions of different features within the chest CT image. Some methods designNorms or +.>The displacement field of norms regularizes to handle sliding interfaces of some organ surfaces, however, optimization of the objective function is difficult due to non-smoothness. Furthermore, these methods typically do not directly take into account the effects of noise, resulting in robustness of the algorithm being neglected.
Recently, researchers have proposed using Welsch functions to achieve robust registration in the field of rigid and non-rigid registration studies of three-dimensional meshes. However, three-dimensional grid registration is very different from image registration (in particular three-dimensional chest CT image registration), the distance metric and regularization mode of the three-dimensional grid registration are very different, the form of the optimized target variable is also different, and the three-dimensional grid registration cannot be used for CT image registration unless a great deal of labor force is used to creatively modify the algorithm.
Therefore, how to comprehensively consider different intensity values and different regions of motion features and to eliminate noise influence as much as possible is a problem to be solved in chest CT image registration.
Disclosure of Invention
The invention provides a robust chest CT image registration method, which can consider different intensity values and different areas of motion characteristics, and uses a dynamic noise tolerance mechanism to eliminate noise influence as much as possible, so as to enhance robustness.
In order to achieve the above purpose, the invention adopts the following technical scheme: a robust chest CT image registration method specifically comprises the following steps:
s1: respectively dividing a lung mask, a bone mask and a spine mask from the floating image and the fixed image; constructing a Gaussian pyramid, and embedding floating/fixed images of the lung, the bone and other areas except the lung and the bone, fixed masks of the spine and fixed masks of other areas except the spine into the Gaussian pyramid according to the resolution corresponding to the layer number;
s2: controlling a dense displacement field of the whole floating image by using a control lattice displacement field, using different dissimilarity measures and regularization for different areas, and constructing an objective function by using a Welsch function and a substitution function thereof;
s3: the pyramid layer number is 1 initially, all floating images, fixed masks and Welsch parameters are updated according to the Gaussian pyramid layer number, an AMSGrad optimizer is used for iterating and minimizing an objective function, the Welsch parameters are dynamically updated in the iteration process, and an optimal control lattice displacement field is solved;
s4: judging whether the layer is a Gaussian pyramid final layer, if not, adding one layer number, up-sampling a control lattice displacement field to the size of the next layer, and returning to the step S3; if yes, go to step S5;
s5: and twisting the floating image according to the dense displacement field of the final layer to obtain a registered floating image.
Further, the steps of separating the lung mask, bone mask, and spine mask from the floating image and the fixed image, respectively, include the following: and respectively obtaining the lung mask from the floating image and the fixed image through a deep learning model, obtaining the bone mask through threshold segmentation, obtaining the expanded lung mask through morphological expansion operation on the lung mask, performing AND operation on the expanded lung mask and the bone mask, extracting the largest connected sub-block to obtain a partially missing spine mask, and performing OR operation on the partially missing spine mask and the middle part of the bone mask to obtain the complete spine mask.
Further, in S1, constructing a gaussian pyramid, and embedding the floating/fixed image of the lung, the bone, the rest of the area except the lung and the bone, the fixed mask of the spine, and the fixed mask of the rest of the area except the spine into the gaussian pyramid according to the resolution corresponding to the number of layers, the steps include: firstly, constructing the Gaussian pyramid, namely, registering and setting different layers, setting corresponding resolution for each layer, wherein the higher the layer number is, the higher the resolution is, so that the registration is carried out from the lower layer to the higher layer from the thick layer to the thin layer; then subtracting the lung mask and the bone mask from the whole image mask of the floating image and the fixed image to obtain a rest region mask except for the lung and the bone, and extracting the floating/fixed image of the lung, the bone and the rest region except for the lung and the bone through the lung mask, the bone mask and the rest region mask except for the lung and the bone; the spine mask of the fixed image is the spine fixed mask, and the spine fixed mask is subtracted from the whole image mask of the fixed image to obtain the fixed masks of the rest areas except the spine; finally, floating/fixed images of the lung, bone, remaining regions other than the lung, bone, the fixed mask of the spine, the fixed mask of the remaining regions other than the spine are embedded into each layer of the gaussian pyramid at a resolution sampled to correspond to.
Further, in S2, the dense displacement field step of controlling the floating image using the control lattice displacement field includes the following: embedding uniform control lattices into the floating image according to preset intervals of height, width and depth, so that the displacement vector of each voxel of the floating image is controlled by a control point on the control lattice, namely, the registered target is converted into search of the displacement vector of the optimal sparse control point, and the process of controlling the displacement vector of the nearby voxels by the control point is realized by using free deformation based on B-spline.
Further, in S2, using different dissimilarity measures and regularization steps for different regions includes the following: first, using floating/fixed images of the lung, bone, and the rest of the region except the lung, bone embedded in the gaussian pyramid in S1, using different dissimilarity measures for the lung, bone, and the rest of the region except the lung, bone; specifically, for bone regions, measuring dissimilarity using normalized cross-correlation NCC is defined as NCC dissimilarity; for the lungs and the rest of the area except the lungs, bones, using normalized gradient field NGF to measure dissimilarity is defined as NGF dissimilarity, NCC dissimilarity is expressed as:
，
wherein,representing a domain of a floating image or a fixed image, wherein the floating image has the same domain as the fixed image,/->Is the coordinates of the voxel in the floating or fixed image, < >>Representing floating image +.>A fixed image is represented and a fixed image is represented,representing the average intensity value of the fixed image, +.>Representing the average intensity value of the floating image, +.>Representing the passing of dense displacement field->Distorted floating image +.>At->Intensity values at that point, warping refers to moving each voxel in a dense displacement field and then interpolating the entire image into a new uniform image, with NGF dissimilarity expressed as:
，
wherein,is used here only to illustrate the effect of the symbols on the internal parameters,/->Is the number of voxels in the floating image, +.>Is a smaller value for avoiding errors, < >>Representing the gradient, then using different control lattice displacement field positive for the spine, the rest of the region except the spine, using the fixation mask of the spine embedded with the Gaussian pyramid, the fixation mask of the rest of the region except the spine in step S1Then the method comprises the steps of (1) dissolving; in particular, for the spinal region, the square of the control lattice displacement field is used +.>Norms to regularize the control lattice displacement field; for the rest of the region except the spine, use is made of +.>Norms are used to regularize the control lattice displacement field.
Further, in S2, for NGF dissimilarity and control of lattice displacement field gradientsNorm regularization using Welsch function +.>To be modified, wherein->To control the parameters of the Welsch function sparsity, the modified parameters are called WelschNGF dissimilarity and Welsch regularization, and the WelschNGF dissimilarity is expressed as:
，
welsch regularization is expressed as:
，
wherein,representing a control lattice domain>Representing the control lattice displacement field +.>Representing the position of a control point on the control lattice, +.>Representing gradient->And->The WelschNGF dissimilarity and the parameters controlling the sparsity of the Welsch function in the Welsch regularization, respectively. Then, the WelschNGF dissimilarity and the Welsch regularization are modified by using a substitution function Sur, and the modified functions are respectively called SurWelschNGF dissimilarity and SurWelsch regularization, and the SurWelschNGF dissimilarity is expressed as:
，
wherein,representing the optimal displacement field for all voxels currently in the iterative process,
the survivinsch regularization is expressed as:
，
wherein,and finally, obtaining an objective function of the registration problem:
，
wherein,weight for bone region dissimilarity measure, +.>Weight for lung region dissimilarity measure, +.>For the dissimilarity measure weight of the rest of the region except the bone and the lung,/for>Regularized weights for spinal region, +.>Regularizing weights for the remaining regions except the spine,/->Fix the image for the bone region,>floating image for bone region->Fixing the image for the lung region,>floating image for lung area->Fixing the image for the remaining areas except for bones, lungs,>floating the image for the remaining area except bone, lung,>mask for fixing the spinal region,/->Fixing mask for the rest of the area except the spine area,/->Representing multiplication by element.
Further, in S3, the Gaussian pyramid layer number is initially 1, according to the Gaussian pyramid layer numberThe steps of updating the floating image, the fixed mask and the Welsch parameters include the following: the Gaussian pyramid layer is initially a first layer, floating/fixed images of the lung, the bone and other areas except the lung and the bone of the current Gaussian pyramid layer and fixed masks of the spine and other areas except the spine are input into the objective function, and parameters in the objective function are input into the objective functionSet to->Where i is the current number of layers, s is the lattice spacing, and l is the average voxel spacing.
Further, in S3, using the AMSGrad optimizer to iterate the minimization objective function, dynamically updating the Welsch parameter during the iteration process, and solving the optimal control lattice displacement field includes the following steps: for the objective function, using an AMSGrad optimizer to perform gradient descent, searching for a control lattice displacement field which minimizes an objective function value, and after each iteration of gradient descent, calculating the average NGF dissimilarity between each region of the floating image after the dense displacement field transformation obtained by the previous iteration and the region of the fixed imageAnd set it to +.f for the term corresponding to that region in the objective function at the next iteration>Values.
Further, in S4, the step of upsampling the control lattice displacement field to the next layer size includes the following: and upsampling the control lattice displacement field to the resolution corresponding to the next layer of the Gaussian pyramid to serve as an initial control lattice displacement field optimized by the objective function of the next layer.
The beneficial effects of the invention are as follows: compared with the prior art, the method provided by the invention uses a divide-and-conquer method to process different metrics and displacement field regularization of different anatomical structures, utilizes the Welsch function to strengthen smoothness and sparsity, and introduces a dynamic parameter updating scheme to adaptively adjust the tolerance to noise and enhance the robustness. In a word, the invention can effectively solve the problem of low registration accuracy of the high-noise CT image.
Drawings
Fig. 1 is a general flow chart of a robust chest CT image registration method provided by an embodiment of the present invention.
Fig. 2 is a specific flowchart of a robust chest CT image registration method according to an embodiment of the present invention.
Fig. 3 is a determination parameter of a robust chest CT image registration method according to an embodiment of the present inventionSchematic diagram of initial value setting scheme.
Fig. 4 is a schematic diagram of a registration result of a robust chest CT image registration method according to an embodiment of the present invention.
Detailed Description
In order to make the objects, technical solutions and advantages of the present invention more apparent, embodiments of the present invention will be described in detail below with reference to the accompanying drawings and examples, it being apparent that the described examples are only some, but not all, examples of the present invention.
An embodiment of the present invention provides a robust chest CT image registration method, as shown in fig. 1, which mainly includes:
step S1, respectively dividing a lung mask, a bone mask and a spine mask from a floating image and a fixed image; a gaussian pyramid is constructed, and floating/fixed images and spine fixing masks of the lung, the bone and the rest of the area except the lung and the bone and the rest of the area except the spine are embedded into the gaussian pyramid according to the resolution corresponding to the layer number.
In the embodiment of the invention, the floating image and the fixed image are three-dimensional chest CT sequence images. Since chest CT images exhibit large scale distortions, direct registration at the original resolution may lead to undesirable local minima. Thus, we integrate different anatomies onto one unified control grid, enabling coarse to fine registration.
The specific process is as follows:
1) Firstly, respectively obtaining a lung mask from a floating image and a fixed image through a deep learning model, and obtaining the bone mask through threshold segmentation.
2) Further, morphological expansion operation is carried out on the lung mask to obtain an expanded lung mask, AND operation is carried out on the expanded lung mask and the bone mask, the largest connected sub-block is extracted, a partially missing spinal mask is obtained, and then OR operation is carried out on the partially missing spinal mask and the middle part of the bone mask, so that the spinal mask is obtained.
3) Further, the lung mask and the bone mask are subtracted from the entire image mask of the floating image and the fixed image to obtain the remaining region mask excluding the lung and the bone, and the floating/fixed image of the lung, the bone, and the remaining region excluding the lung and the bone is extracted from these masks. The spine mask is subtracted from the entire image mask of the fixed image to obtain the remaining area fixed mask except the spine.
4) Further, a gaussian pyramid is constructed, as shown in fig. 2, 5 levels are set for registration, corresponding resolutions are set for each level, and the higher the number of levels is, the higher the resolution is, so that registration is performed from a lower level to a higher level from thick to thin. The resolution of the last layer is the original resolution.
5) Further, floating/fixed images of the lung, bone, and the remaining regions other than the lung, bone, and fixed masks of the spine, and the remaining regions other than the spine, are embedded in each layer of the gaussian pyramid at the resolution sampled to correspond. When registration proceeds to a layer, dissimilarity measures and regularization are performed using the image and mask embedded in the layer. With the additional guidance provided by the mask, higher accuracy can be achieved at lower levels, preventing errors from propagating and accumulating between layers.
And S2, controlling a dense displacement field of the whole floating image by using a control lattice displacement field, using different dissimilarity measures and regularization on different areas, and constructing an objective function by using a Welsch function and a substitution function thereof.
The specific process is as follows:
1) The registration efficiency can be improved by controlling the lattice, i.e. controlling the entire dense displacement field with a set of sparse and evenly distributed control points. Therefore, a uniform control lattice is embedded in the image according to preset axial intervals, so that the displacement vector of each voxel of the image is controlled by the control points on the control lattice, namely, the registered target is converted into the search of the displacement vector of the optimal sparse control point. The control point spacing of the embodiment of the invention is set as. The process of the control point controlling the displacement vector of nearby voxels is implemented using B-spline based free deformation.
2) Further, using floating/fixed images of the lungs, bones, and other areas except the lungs, bones embedded in the gaussian pyramid, different dissimilarity measures are used for the lungs, bones, and other areas except the lungs, bones. Specifically, since bone regions are mainly rigidly displaced as a whole, measuring dissimilarity using normalized cross-correlation NCC for bone regions is defined as NCC dissimilarity; for the lungs and the rest of the areas except the lungs, bones, the normalized gradient field NGF is used to measure what is defined as NGF dissimilarity, since the air and organ internal intensity values in the lung area are mostly uniform and not suitable for guiding localization, while the blood vessel rich structures inside the lungs and the boundaries of the organs can be used for guiding localization. NCC dissimilarity is expressed as:
，
wherein,representing a floating image or a fixed image domain (floating image has the same domain as fixed image),is that the voxels float in the image or solidCoordinates of the fixed image>Representing floating image +.>Representing a fixed image +.>Representing the average intensity value of the fixed image, +.>Representing the average intensity value of the floating image, +.>Representing the passing of dense displacement field->Distorted floating image +.>At->Intensity values at. "warping" refers to moving each voxel in a dense displacement field and then interpolating the entire image into a new uniform image. NGF dissimilarity is expressed as:
，
wherein, ，
f, g are used herein only to illustrate the effect of the symbol on the internal parameters,is the number of voxels in the floating image. />Is used for avoiding errorA smaller value of the difference.
3) Further, using a fixed mask embedded in the spine of the gaussian pyramid, the rest of the area except the spine, different displacement fields are regularized to the spine, the rest of the area except the spine. Specifically, since the spinal region is hardly displaced, the square of the control lattice displacement field is used for the spinal regionRegularizing the control lattice displacement field by a norm to punish the displacement; for the rest of the area except the spine, the control of the gradient of the lattice displacement field is used because of the displacement, because the displacement of the most area is smoother>Norms regularize the control lattice displacement field to penalize differences in the displacements of adjacent control points.
4) Further, a smooth Welsch function is usedModification of NGF dissimilarity term and control of lattice displacement field gradient>Normals, wherein->To control the parameters of the Welsch function sparsity, the modified parameters are called WelschNGF dissimilarity and Welsch regularization, and the WelschNGF dissimilarity is expressed as:
，
welsch regularization is expressed as:
，
wherein,representing a control lattice domain>Representing the control lattice displacement field +.>Representing the coordinates of a control point on the control lattice, +.>Representing the gradient.
5) Further, since the Welsch function has been demonstrated to exist as a convex quadratic upper bound function at any point y，/>Time->Is greater than->Otherwise->Equal to->. Thus, pair->Can be converted into p +.>And due to the minimization ofIs a convex quadratic function, and has lower optimization complexity. In the embodiment of the invention, the optimal dense displacement field searched by the current optimizer is used for each iteration>And controlling the lattice shift field->As->The WelschNGF dissimilarity and the Welsch regularization term are modified as substitution functions, and the modified term is respectively called SurWelschNGF dissimilarity and SurWelsch regularization. "transformation" refers to nesting original functions using substitution functions and removing constant terms that do not affect optimization. SurWelschNGF dissimilarity is expressed as:
，
the survivinsch regularization is expressed as:
，
6) Further, the dissimilarity measure terms of bones, lungs, other regions except the lungs and bones are multiplied by the respective weights and summed, so that the objective function of the registration problem is obtained:
，
wherein each ofAnd->Weights in the objective function are normalized for dissimilarity and regularization of regions.
And S3, updating the image and the Welsch parameters according to the Gaussian pyramid layer number (1 initially), using an AMSGrad optimizer to iterate and minimize an objective function, dynamically updating the Welsch parameters in the iteration process, and solving the optimal control lattice displacement field.
The specific process is as follows:
1) The floating/fixed images of the lungs, bones, the rest of the area except the lungs, bones, and the fixed mask of the spine, the rest of the area except the spine of the current gaussian pyramid layer (initially the first layer) are input into the objective function.
2) Parameters (parameters)And->Plays a key role in the robustness of the method: observations show that the WelschNGF term and the Welsch regularized weight term are ++according to their own losses>And->Dynamic updating of canonical Gaussian functions, +.>And->Is the standard deviation of the gaussian function. This means that the larger the loss value of the corresponding point, the smaller the weight it takes in the objective function loss. In the final phase of the optimization, < > the>And->Should be small enough that the gaussian weights effectively attenuate the effects of large errors caused by noise. If it is to->And->Fixing at a smaller value at all times reduces the optimisationSince at the beginning many corresponding points are lost so much that if the gaussian weights are small, then there is not enough penalty to be ignored in the optimization process. Therefore, we have done the optimization process on +.>And->And dynamically updating.
3) Regularization parameters for WelschInitially, since the resolution is low, the number of control points is small, and the slip needs to be strictly controlled to prevent errors from propagating and accumulating between layers. Since the greatest sliding occurs at the lung boundary, it is assumed that the lung boundary can move up to 100 mm along each axis (y-axis and z-axis) parallel to the sliding plane (the plane made up of the y-axis and z-axis) while the spine remains stationary, as shown in FIG. 3, when the gradient of displacement of the boundary-most lung voxel (i.e., the middle voxel in FIG. 3) in the x-axis is +.>The gradient of the displacement in the y-axis is +.>Gradient in z-axis is +.>Thus total displacement gradient +.>The norm penalty is +.>Where i is the current number of layers, s is the lattice spacing, and l is the average voxel spacing. We will->The initial value is set to +.>Sufficient to severely penalize the sliding. When the Gaussian pyramid enters the next layer, the Gaussian pyramid is addedUpdated to 1/4 of the own. Thus, when initializing the objective function per layer, it is +_based on the current layer number>Parameter in the objective function +.>Set to->。
4) Further, the objective function is minimized by using an AMSGrad optimizer, and the basic idea is to gradient the objective functionAnd then back-propagating to adjust the control lattice displacement field. For each iteration we will next iterate +.>Setting the average NGF dissimilarity between the floating image and the fixed image after the dense displacement field transformation obtained by the previous iteration>I.e. "new welsch parameter" in fig. 3. It should be noted that the above procedure requires separate treatment of the lung region using the WelschNGF dissimilarity measure and the rest of the region except the lung, bone. According to the well known->Principle, when normalized gradient field NGF spatial distance is less than +.>When the weight is large, the optimization is realizedContains enough terms that when normalized gradient field NGF spatial distance is greater than +.>When the weight is small, the noisy corresponding point can be eliminated. As the optimization proceeds, the average dissimilarity gradually decreases, +.>And also decreases.
5) Further, when the iteration reaches the limit of the number of times, or the difference value of each variable of two adjacent iterations is smaller than the threshold value, the iteration of the current layer is ended. In the present embodiment, the iteration number limit of 5 layers is set to 800, 800, 800, 100, 50, and the threshold is set to 800。
And S4, judging whether the layer is a Gaussian pyramid final layer. If the control lattice displacement field is not the final layer, the optimal control lattice displacement field is up-sampled to the control lattice displacement field size required by the next layer, and the control lattice displacement field is used as the initial control lattice displacement field of the next layer. If the layer is the last layer, the optimal dense displacement field is the dense displacement field obtained by registration.
Step S5: and twisting the floating image according to the dense displacement field obtained by registration to obtain a registered floating image, and finishing registration.
Fig. 4 is an example of registration of an embodiment of the present invention, visualizing the displacement of 300 marker points of two chest CTs with different displacements resulting from the registration.
In the embodiment of the invention, different dissimilarity measures are designed for different regions by utilizing the abundant trachea and blood vessels in the lung region, the rigid movement of the bone region and the obvious intensity change between organs. Furthermore, embodiments of the present invention customize the proper displacement field regularization for each region, taking into account the minimal displacement of the spine and the smooth and sliding motion that occurs in other regions.
In addition, the embodiment of the invention also utilizes a Welsch function to treat NGF incompatibilitiesSimilarity measureSparsity is induced in regularization of the displacement field of norms. The embodiments of the present invention then process the Welsch term using an optimization-Minimization (optimization-Minimization) algorithm. The algorithm builds a proxy function of the objective function based on the current variable values and minimizes the proxy function to update the variables, thereby reducing the complexity of the optimization process. In addition, the embodiment of the invention introduces a dynamic parameter updating scheme, so that sparsity is dynamically updated, the tolerance to noise is adjusted in a self-adaptive mode, the robustness is enhanced, and the high-noise chest CT image registration is realized.
The technical scheme of the embodiment can be seen that: the method of the invention constructs different dissimilarity measure and regularization for each region by dividing and controlling, samples into different resolutions, embeds Gaussian pyramid, then uses control lattice to control displacement field, constructs objective function by smooth Welsch function, and uses alternative function to accelerate iteration efficiency, and dynamically adjusts parameters each iteration to reduce sensitivity to noise and increase robustness. In summary, the invention can effectively solve the problem of low registration accuracy of the high-noise CT image.
The foregoing description of the preferred embodiments of the invention is not intended to limit the invention to the precise form disclosed, and any such modifications, equivalents, and alternatives falling within the spirit and scope of the invention are intended to be included within the scope of the invention.
Claims (9)
1. A robust chest CT image registration method, comprising the steps of:
s1: respectively dividing a lung mask, a bone mask and a spine mask from the floating image and the fixed image; constructing a Gaussian pyramid, and embedding floating/fixed images of the lung, the bone and other areas except the lung and the bone, fixed masks of the spine and fixed masks of other areas except the spine into the Gaussian pyramid according to the resolution corresponding to the layer number;
s2: controlling a dense displacement field of the whole floating image by using a control lattice displacement field, using different dissimilarity measures and regularization for different areas, and constructing an objective function by using a Welsch function and a substitution function thereof;
s3: the pyramid layer number is 1 initially, all floating images, fixed masks and Welsch parameters are updated according to the Gaussian pyramid layer number, an AMSGrad optimizer is used for iterating and minimizing an objective function, the Welsch parameters are dynamically updated in the iteration process, and an optimal control lattice displacement field is solved;
s4: judging whether the layer is a Gaussian pyramid final layer, if not, adding one layer number, up-sampling a control lattice displacement field to the size of the next layer, and returning to the step S3; if yes, go to step S5;
s5: and twisting the floating image according to the dense displacement field of the final layer to obtain a registered floating image.
2. The method of registration of robust chest CT images of claim 1, wherein in S1, the step of separately segmenting the lung mask, bone mask, spine mask from the floating image and the fixed image comprises: and respectively obtaining the lung mask from the floating image and the fixed image through a deep learning model, obtaining the bone mask through threshold segmentation, obtaining the expanded lung mask through morphological expansion operation on the lung mask, performing AND operation on the expanded lung mask and the bone mask, extracting the largest connected sub-block to obtain a partially missing spine mask, and performing OR operation on the partially missing spine mask and the middle part of the bone mask to obtain the complete spine mask.
3. The robust chest CT image registration method according to claim 2, wherein in S1, the step of constructing a gaussian pyramid, embedding floating/fixed images of the lung, bone, remaining regions other than the lung, bone, and fixed masks of the spine, fixed masks of remaining regions other than the spine, into the gaussian pyramid with a resolution corresponding to the number of layers comprises: firstly, constructing the Gaussian pyramid, namely setting different layers for registration, setting corresponding resolution for each layer, and enabling the registration to be carried out from low layers to high layers from thick to thin as the layer number is higher and the resolution is higher; then subtracting the lung mask and the bone mask from the whole image mask of the floating image and the fixed image to obtain a rest region mask except for the lung and the bone, and extracting the floating/fixed image of the lung, the bone and the rest region except for the lung and the bone through the lung mask, the bone mask and the rest region mask except for the lung and the bone; the spine mask of the fixed image is the spine fixed mask, and the spine fixed mask is subtracted from the whole image mask of the fixed image to obtain the fixed masks of the rest areas except the spine; finally, floating/fixed images of the lung, bone, remaining regions other than the lung, bone, the fixed mask of the spine, the fixed mask of the remaining regions other than the spine are embedded into each layer of the gaussian pyramid at a resolution sampled to correspond to.
4. A robust chest CT image registration method according to claim 3, wherein in S2 the dense displacement field step of controlling the floating image using the control lattice displacement field comprises the following: embedding uniform control lattices into the floating image according to preset intervals of height, width and depth, so that the displacement vector of each voxel of the floating image is controlled by a control point on the control lattice, namely, the registered target is converted into search of the displacement vector of the optimal sparse control point, and the process of controlling the displacement vector of the nearby voxels by the control point is realized by using free deformation based on B-spline.
5. The robust chest CT image registration method according to claim 4, wherein in S2, using different dissimilarity measures and regularization steps for different regions comprises: first, using floating/fixed images of the lung, bone, and the rest of the region except the lung, bone embedded in the gaussian pyramid in S1, using different dissimilarity measures for the lung, bone, and the rest of the region except the lung, bone; specifically, for bone regions, we use normalized cross-correlation NCC to measure dissimilarity, defined as NCC dissimilarity; for the lungs and the rest of the area except the lungs, bones, normalized gradient field NGF was used to measure dissimilarity, defined as NGF dissimilarity, NCC dissimilarity expressed as:
，
wherein,representing a floating image or a domain of a fixed image, wherein the floating image has the same domain as the fixed image,is the coordinates of the voxel in the floating or fixed image, < >>Representing floating image +.>Representing a fixed image +.>Representing the average intensity value of the fixed image, +.>Representing the average intensity value of the floating image, +.>Representing the passing of dense displacement field->Distorted floating image +.>At->Intensity values at that point, warping refers to moving each voxel in a dense displacement field and then interpolating the entire image into a new uniform image, with NGF dissimilarity expressed as:
，
wherein,is used here only to illustrate the effect of the symbols on the internal parameters,/->Is the number of voxels in the floating image, +.>Is a smaller value for avoiding errors, < >>Representing the gradient, and then regularizing the spine and the rest of the regions except the spine by using different control lattice displacement fields by using the fixation mask of the spine embedded with the Gaussian pyramid and the fixation mask of the rest of the regions except the spine in the step S1; in particular, for the spinal region, the square of the control lattice displacement field is used +.>Norms to regularize the control lattice displacement field; for the rest of the region except the spine, use is made of +.>Norms are used to regularize the control lattice displacement field.
6. The robust chest CT image registration method as recited in claim 5, whichCharacterized in that in S2, for NGF dissimilarity and control of lattice displacement field gradientNorm regularization using Welsch function +.>To be modified, wherein->To control the parameters of the Welsch function sparsity, the modified parameters are called WelschNGF dissimilarity and Welsch regularization, and the WelschNGF dissimilarity is expressed as:
，
welsch regularization is expressed as:
，
wherein the method comprises the steps ofRepresenting a control lattice domain>Representing the control lattice displacement field +.>Representing the position of a control point on the control lattice, +.>Representing gradient->And->Parameters for controlling the sparsity of the Welsch function in the WelschNGF dissimilarity and the Welsch regularization are respectively; the WelschNGF dissimilarity and Welsch regularization are then modified using a substitution function Sur, and the modified dissimilarity and WelschNGF regularization are referred to as SurWelschNGF dissimilarity and SurWelsch regularization, respectively, and the SurWelschNGF dissimilarity is expressed as:
，
wherein,representing the optimal displacement field for all voxels currently in the iterative process,
the survivinsch regularization is expressed as:
，
wherein,and finally, obtaining an objective function of the registration problem:
，
wherein,weight for bone region dissimilarity measure, +.>Weight for lung region dissimilarity measure, +.>For the dissimilarity measure weight of the rest of the region except the bone and the lung,/for>Regularized weights for spinal region, +.>Regularizing weights for the remaining regions except the spine,/->Fix the image for the bone region,>floating image for bone region->An image is fixed for the lung region,floating image for lung area->Fixing the image for the remaining areas except for bones, lungs,>floating the image for the remaining area except bone, lung,>mask for fixing the spinal region,/->Fixing mask for the rest of the area except the spine area,/->Representing multiplication by element.
7. The robust chest CT image registration method of claim 6, wherein SIn the step 3, the number of layers of the Gaussian pyramid is initially 1, and the steps of updating the floating image, the fixed mask and the Welsch parameters according to the number of layers of the Gaussian pyramid comprise the following contents: the Gaussian pyramid layer is initially a first layer, floating/fixed images of the lung, the bone and other areas except the lung and the bone of the current Gaussian pyramid layer and fixed masks of the spine and other areas except the spine are input into the objective function, and parameters in the objective function are input into the objective functionSet to->Where i is the current number of layers, s is the lattice spacing, and l is the average voxel spacing.
8. The robust chest CT image registration method of claim 7, wherein in S3, the step of iteratively minimizing the objective function using an AMSGrad optimizer, dynamically updating Welsch parameters during the iteration, and solving for the optimal control lattice displacement field comprises: for the objective function, performing gradient descent by using an AMSGrad optimizer, searching for a control lattice displacement field which minimizes an objective function value, and after each iteration of gradient descent, calculating the average NGF dissimilarity between each region of the floating image after the dense displacement field transformation obtained by the previous iteration and the fixed image corresponding to each regionAnd set it to +_ of the term corresponding to each region in the objective function at the next iteration>Values.
9. The method of robust chest CT image registration of claim 8, wherein in S4, the step of upsampling the control lattice displacement field to the next layer size comprises: and upsampling the control lattice displacement field to the resolution corresponding to the next layer of the Gaussian pyramid to serve as an initial control lattice displacement field optimized by the objective function of the next layer.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202410124244.4A CN117649416B (en) | 2024-01-30 | 2024-01-30 | Robust chest CT image registration method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202410124244.4A CN117649416B (en) | 2024-01-30 | 2024-01-30 | Robust chest CT image registration method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN117649416A true CN117649416A (en) | 2024-03-05 |
| CN117649416B CN117649416B (en) | 2024-04-12 |
Family
ID=90050010
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202410124244.4A Active CN117649416B (en) | 2024-01-30 | 2024-01-30 | Robust chest CT image registration method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN117649416B (en) |
Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20090046951A1 (en) * | 2007-07-16 | 2009-02-19 | Nikos Paragios | System and Method for Dense Image Registration Using Markov Random Fields and Efficient Linear Programming |
| US20140192046A1 (en) * | 2013-01-07 | 2014-07-10 | Ecole Centrale Paris | Method and device for elastic registration between a two-dimensional digital image and a slice of a three-dimensional volume with overlapping content |
| CN105303547A (en) * | 2014-07-11 | 2016-02-03 | 东北大学 | Multiphase CT image registration method based on grid matching Demons algorithm |
| US20160335777A1 (en) * | 2015-05-13 | 2016-11-17 | Anja Borsdorf | Method for 2D/3D Registration, Computational Apparatus, and Computer Program |
| CN106204561A (en) * | 2016-07-04 | 2016-12-07 | 西安电子科技大学 | Prostate multi-modality images non-rigid registration method based on mixed model |
| CN107230223A (en) * | 2017-06-09 | 2017-10-03 | 中国科学院苏州生物医学工程技术研究所 | Liver's three-dimensional multimode state method for registering images based on discontinuous fluid |
| CN111798500A (en) * | 2020-07-20 | 2020-10-20 | 陕西科技大学 | Differential homoembryo non-rigid registration algorithm based on hierarchical neighborhood spectrum features |
| CN113822796A (en) * | 2021-09-18 | 2021-12-21 | 长春理工大学 | Multi-modal brain image registration method based on improved image pyramid |
| CN116883467A (en) * | 2023-07-17 | 2023-10-13 | 大连理工大学 | Non-rigid registration method for medical image |
-
2024
- 2024-01-30 CN CN202410124244.4A patent/CN117649416B/en active Active
Patent Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20090046951A1 (en) * | 2007-07-16 | 2009-02-19 | Nikos Paragios | System and Method for Dense Image Registration Using Markov Random Fields and Efficient Linear Programming |
| US20140192046A1 (en) * | 2013-01-07 | 2014-07-10 | Ecole Centrale Paris | Method and device for elastic registration between a two-dimensional digital image and a slice of a three-dimensional volume with overlapping content |
| CN105303547A (en) * | 2014-07-11 | 2016-02-03 | 东北大学 | Multiphase CT image registration method based on grid matching Demons algorithm |
| US20160335777A1 (en) * | 2015-05-13 | 2016-11-17 | Anja Borsdorf | Method for 2D/3D Registration, Computational Apparatus, and Computer Program |
| CN106204561A (en) * | 2016-07-04 | 2016-12-07 | 西安电子科技大学 | Prostate multi-modality images non-rigid registration method based on mixed model |
| CN107230223A (en) * | 2017-06-09 | 2017-10-03 | 中国科学院苏州生物医学工程技术研究所 | Liver's three-dimensional multimode state method for registering images based on discontinuous fluid |
| CN111798500A (en) * | 2020-07-20 | 2020-10-20 | 陕西科技大学 | Differential homoembryo non-rigid registration algorithm based on hierarchical neighborhood spectrum features |
| CN113822796A (en) * | 2021-09-18 | 2021-12-21 | 长春理工大学 | Multi-modal brain image registration method based on improved image pyramid |
| CN116883467A (en) * | 2023-07-17 | 2023-10-13 | 大连理工大学 | Non-rigid registration method for medical image |
Non-Patent Citations (4)
| Title |
|---|
| ROMAN SCHAFFERT, ET AL.: "Robust Multi-View 2-D/3-D Registration Using Point-To-Plane Correspondence Model", 《IEEE TRANSACTIONS ON MEDICAL IMAGING》, vol. 39, no. 1, 13 June 2019 (2019-06-13) * |
| ZONGTAO HU, ET AL.: "Non-rigid Registration of White Matter Tractography Using Coherent Point Drift Algorithm", 《MULTIMODAL BRAIN IMAGE ANALYSIS AND MATHEMATICAL FOUNDATIONS OF COMPUTATIONAL ANATOMY》, 10 October 2019 (2019-10-10) * |
| 杭利华;蒋佩钊;邓冲;: "基于仿射变换与B样条自由形变的医学图像配准", 兰州交通大学学报, no. 03, 15 June 2013 (2013-06-15) * |
| 胡宗涛等: "针对颅脑放疗规划的海马体自动勾画平台及其验证", 《北京生物医学工程》, vol. 39, no. 4, 31 August 2020 (2020-08-31) * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN117649416B (en) | 2024-04-12 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN104835112B (en) | A kind of liver multiphase phase CT image interfusion methods | |
| Klein et al. | Evaluation of optimization methods for nonrigid medical image registration using mutual information and B-splines | |
| GB2559013B (en) | A method and system for image processing | |
| Tagare | Deformable 2-D template matching using orthogonal curves | |
| CN104424629B (en) | A kind of x-ray chest radiograph lung segmentation method and apparatus | |
| US8582854B2 (en) | Method and system for automatic coronary artery detection | |
| JP4885138B2 (en) | Method and system for motion correction in a sequence of images | |
| CN109064476B (en) | CT chest radiography lung tissue image segmentation method based on level set | |
| US20090190809A1 (en) | Method and Apparatus for Efficient Automated Re-Contouring of Four-Dimensional Medical Imagery Using Surface Displacement Fields | |
| CN107451995B (en) | Method, device, equipment and storage medium for extracting cardiovascular from CTA image | |
| CN112602099A (en) | Deep learning based registration | |
| US9269139B2 (en) | Rib suppression in radiographic images | |
| CN112614169B (en) | 2D/3D spine CT (computed tomography) level registration method based on deep learning network | |
| CN102243759A (en) | Three-dimensional lung vessel image segmentation method based on geometric deformation model | |
| US9042619B2 (en) | Method and system for automatic native and bypass coronary ostia detection in cardiac computed tomography volumes | |
| Zeng et al. | Liver segmentation in magnetic resonance imaging via mean shape fitting with fully convolutional neural networks | |
| Zhou et al. | Interactive medical image segmentation using snake and multiscale curve editing | |
| CN106846330A (en) | Human liver's feature modeling and vascular pattern space normalizing method | |
| CN111815663A (en) | Hepatic vessel segmentation system based on Hessian matrix and gray scale method | |
| JP2006506153A (en) | Image registration method and apparatus | |
| CN110288581A (en) | A kind of dividing method based on holding shape convexity Level Set Models | |
| CN117649416B (en) | Robust chest CT image registration method | |
| Wei et al. | Morphology-preserving smoothing on polygonized isosurfaces of inhomogeneous binary volumes | |
| CN113327274B (en) | Lung CT image registration method and system integrating segmentation function | |
| CN111260546B (en) | Image processing method and device and electronic equipment |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant |