Prostate MRI-TRUS deformable image registration method based on biomechanical generation

Abstract

The invention belongs to the technical field of medical image registration, and particularly relates to a prostate MRI-TRUS deformable image registration method based on biomechanical generation. According to the method, biomechanics obtained through finite element analysis is learned into a low-dimensional deformation manifold in a priori mode through a generative model, and the learned manifold is optimized through shape similarity so as to find an optimal deformation field. According to the method, the domain-invariant prostate deformation rule is directly learned from biomechanical prior, so that the influence of modal difference of the image and corresponding structure deficiency on registration can be effectively overcome; in the reasoning stage, a deformation field with implicit biomechanical constraints can be directly output, and the limitation of explicit regularization is effectively overcome; and the precision of MRI-TRUS deformable image registration and the biomechanical rationality of a deformation field can be remarkably improved, and excellent performance is shown in the aspect of keeping the internal structure of the prostate. The method does not depend on image texture features to establish a spatial corresponding relation, and has high generalization and generalizability.

Description

Technical Field

[0001] The invention belongs to the technical field of medical image registration, and in particular relates to a prostate MRI-TRUS image registration method. Background Art

[0002] Prostate cancer is one of the most prevalent tumor diseases in middle-aged and elderly men, and its diagnosis usually relies on 12-core biopsy guided by transrectal ultrasound (TRUS). However, the resolution of TRUS images is low, making it difficult to clearly display the lesion area in the prostate, which can easily lead to inaccurate positioning of the treatment target and increase the false negative rate. Magnetic resonance imaging (MRI) has higher clarity and contrast, which can improve the detection rate and positioning accuracy of prostate cancer, but it cannot be used for real-time imaging due to its speed and cost. Therefore, the fusion of preoperative MRI and TRUS has become the preferred method for prostate targeted intervention. However, due to changes in patient position, the compression of the prostate gland by the ultrasound probe in the rectum, and physiological activities such as breathing, the prostate often deforms and shifts during surgery, causing the target area on the TRUS image to deviate from the preoperative imaging position, which poses a huge challenge to the accurate fusion of MRI and TRUS.

[0003] The key to the fusion of preoperative MRI and intraoperative TRUS images lies in their accurate registration. Compared with other medical image registration tasks, MRI-TRUS deformable image registration of the prostate faces the following challenges: First, MRI and TRUS differ significantly due to different imaging principles, and ultrasound images have low resolution and relatively limited field of view, which leads to the loss of some anatomical structures visible in MRI, making it difficult for the registration method to obtain the corresponding contextual position information between organs and the internal features of the prostate. In addition, the compression of the prostate by the TRUS probe will cause large deformation. Soft tissue deformable registration is an ill-posed inverse problem with many possible spatial transformations.

[0004] With the rapid development of artificial intelligence, many MRI-TRUS deformable image registration methods based on deep learning have emerged. In order to overcome the difficulties caused by domain differences and structural missing, they either use weak supervision labels such as anatomical annotations and segmentation to drive registration, or simplify multimodal problems into single-modal problems through image translation. In addition, these methods usually use explicit regularization terms (such as L2 norm) to restrict the solution in order to improve the rationality of the output deformation field. However, they still have some limitations: First, most existing methods rely on inferring the deformation field through image features, which is easily affected by noise and artifacts in ultrasound. Second, explicit regularization introduces parameter setting problems. If the regularization weight is too small, the deformation field may appear physically unreasonable folding; if the weight is too large, it will limit the model's ability to capture complex deformations. Summary of the invention

[0005] The purpose of the present invention is to provide a prostate MRI-TRUS deformable image registration method based on biomechanical generation to overcome the influence of noise and artifacts in ultrasound and provide flexible and effective implicit biomechanical constraints for the output deformation field, thereby improving the registration accuracy and generalization of the model.

[0006] The present invention provides a prostate MRI-TRUS deformable image registration method based on biomechanical generation, which uses a generative model to learn the domain-invariant prostate biomechanical deformation law, and realizes registration through optimization in latent space; including: embedding the deformation law obtained in the finite element biomechanical prior into the low-dimensional manifold in the latent space through a VAE-based generative model, so that the generative model can generate biomechanically reasonable deformation without explicit regularization terms; then optimizing on the learned manifold through shape matching to find the deformation field that best suits the image to be registered, thereby effectively overcoming the influence of noise and artifacts in ultrasound and improving the registration accuracy and generalization of the model. The specific steps are:

[0007] Step 1: Preprocess the image. This includes first normalizing the image intensity and performing histogram equalization; then resampling all voxels to a uniform size; then cropping the image to the same size centered on the region of interest; and finally performing preliminary rigid alignment.

[0008] Step 2, obtain biomechanical priors. This includes first using TetGen [5] The prostate segmentation results of MRI and TRUS images are gridded and then non-rigid coherent point drift (CPD) is used [1]The algorithm aligns the surface points of the MRI grid to the surface points of the TRUS. Specifically, each point m in the MR prostate surface point set is regarded as the centroid of a Gaussian Mixture Model (GMM), while the surface point of the US prostate represents the observation x of the GMM. All components of the GMM are set to follow the mean y m (the position of the MR point) and the variance is σ 2 Gaussian distribution with equal membership probability Then the probability density function of GMM is expressed as:

[0009]

[0010] Where M represents the number of MRI surface points, ω is the weight of the uniform distribution weight added to represent noise, 0≤1≤ω.

[0011] Then, the iterative expectation maximization algorithm is used to maximize the likelihood E of the GMM observations, as shown in the following formula:

[0012]

[0013] Where n represents the points in the TRUS point set, and its number is N; in solving the final θ and σ 2 After that, the deformed position of the moving point set can be calculated.

[0014] Next, the neo-Hookean hyperelastic model is used [6] To model the prostate, the strain energy function of the neo-Hookean material is defined as follows:

[0015]

[0016]

[0017] Where λ is the first Lamé parameter and μ is the shear modulus or second Lamé parameter, both of which are material properties and can be calculated from Young's modulus E and Poisson's ratio v; I1 = trace(FF T ) and J = det(F) are invariants, where F is the deformation gradient tensor, through Calculated, where X and x represent the material point positions before and after deformation, respectively. Then, using FEBio 4.5.0 [2] Finite element analysis was performed and the surface point displacement field obtained from MRI to TRUS was used as the boundary strip. In the calculation, the load on the prostate was considered to be applied smoothly and the boundary conditions were applied in five steps to obtain an incremental iterative solution, i.e., the amplitude of each step increased by one fifth of the total boundary displacement in order to promote convergence and avoid unreasonable deformation.

[0018] Finally, the results of the finite element analysis are linearly interpolated to obtain a voxel-level dense deformation field with the same size and pixel spacing as the image pair to be registered. The complete biomechanical prior acquisition process can be found in Figure 1 (a) shown.

[0019] Step 3: Use the generative model to learn the biomechanical prior as a deformation manifold in the latent space. A generative model based on a variational autoencoder (VAE) is used to implicitly embed the deformation field in step 2 into the deformation manifold in the latent space of the generative model, thereby compactly representing the deformation during prostate surgery in a low-dimensional space. Specifically, the generative model consists of an encoder and a multi-scale decoder, where the encoder extracts features from the input deformation field φ and embeds it into a low-dimensional latent space, which is then represented by a Gaussian distribution constructed by mean μ and variance σ through reparameterization; the decoder reconstructs the latent vector z sampled from the distribution into a deformation field. The skip connection and upsampling connection can make full use of the multi-scale information in the process. In addition, the first-order spatial derivative of the deformation field is added to the reconstruction term of the decoder. and To reduce the impact of rigid translation, the training loss function of the generated model is as follows:

[0020]

[0021] Among them, p(z)~N(0,I) represents the prior Gaussian distribution, q θ represents the encoder with parameter θ, D KL Denotes KL divergence. The first two terms of the function are used to penalize the reconstruction error between the input and output in all directions, where α is a hyperparameter that controls the weights between different reconstruction terms; the latter term forces the learned distribution to be as similar as possible to the prior distribution, and β is a hyperparameter that weighs the quality of reconstruction and the degree of regularity of the latent space.

[0022] After training, the VAE decoder The z sampled from the latent space can be reconstructed into the corresponding deformation field, that is, In this way, the biomechanical deformation law of the prostate can be implicitly embedded into the output of the model, making It can generate biomechanically plausible deformations without explicit regularization terms. The training process of the generative model is as follows Figure 1 (b) as shown.

[0023] Step 4: Implement deformable registration through optimization in latent space. First, by Initialize z by random sampling from , and convert it into the corresponding initial deformation field by the trained VAE decoder, and then use the Spatial Transformer Networks (STN) [4] Apply it to the MRI prostate segmentation image to obtain the deformed MRI prostate segmentation image. Then calculate the Dice loss of MR and US prostate segmentation, and drive the optimization process of z in the latent space through back propagation. The loss function of the image in this process is:

[0024]

[0025] Among them, s f and m Denote the prostate segmentation of TRUS and MRI, respectively, and ° denotes the spatial transformation completed by the spatial transformer network STN. In addition, in order to ensure that z is within the range of the learned manifold when searching in the latent space, we add an L2 regularization term (i.e., the L2 norm (Euclidean norm) square of the parameter) to encourage z to be in the high probability area under the corresponding Gaussian prior, and the final loss function is:

[0026]

[0027] Among them, ρ controls the probability that z falls within the Gaussian prior range. By minimizing the above formula, we can find the latent vector z that can be mapped to the deformation field that best suits the input image along the manifold in the latent space. Until two consecutive calculations of L reg When it is less than the set threshold ∈, the optimization process stops, that is, Where i represents the optimization iteration step. The image registration process based on latent space optimization is as follows: Figure 1 (c) as shown.

[0028] In the above formula, the parameters involved are specifically set as follows: ω=0, α=10, β=1, μ=0.001, ∈=0.0001.

[0029] Compared with the prior art, the biomechanically generated prostate MRI-TRUS deformable image registration method of the present invention has the following advantages:

[0030] (1) The prostate deformation pattern is learned directly from the biomechanical prior through the VAE latent space, without relying on image features to establish spatial correspondence, thus reducing the impact of modality differences and missing corresponding structures;

[0031] (2) Implicit biomechanical constraints are imposed on the output deformation field, overcoming the limitations of traditional explicit regularization terms and improving the registration accuracy and rationality of deformation fields;

[0032] (3) By optimizing the latent space to generate the deformation field and achieve registration, the impact of ultrasound artifacts and noise is reduced and the generalization of the model is improved. BRIEF DESCRIPTION OF THE DRAWINGS

[0033] Figure 1 The overall framework diagram of the present invention. Among them, (a) is the biomechanical prior acquisition process, (b) is the training process of the generative model, and (c) is the image registration process based on latent space optimization.

[0034] Figure 2 Prostate MRI-TRUS image registration results. The first column is the input original MRI and TRUS, and the next five columns are the registration results of the method of the present invention and different comparison methods. Columns 1 to 5 are marked with white outlines of the original MRI prostate label, the original TRUS prostate label, and the MRI prostate label after deformation by different methods. The second row shows the deformation field of the corresponding method represented by a grid.

[0035] Figure 3 Visualization of the internal prostate tumor before and after registration. The original TRUS tumor label is marked with a black outline, and the MRI tumor label before and after deformation is marked with a white outline. After registration, the misalignment of the internal tumor is significantly improved, and it has better alignment in both margins and spatial positions. DETAILED DESCRIPTION

[0036] The present invention is further described below by way of embodiments in conjunction with the accompanying drawings.

[0037] Example 1: The method of the present invention is used to perform deformable prostate MRI-TRUS image registration. The specific implementation process is as follows:

[0038] Step 1: A total of 575 training data were collected from the prostate MRI-TRUS biopsy dataset. First, the image was intensity normalized and histogram equalized, and then all voxels were resampled to a uniform size. The image was cropped to the same size of 128×128×128 centered on the region of interest, and finally a preliminary rigid alignment was performed.

[0039] Step 2, first use TetGen to grid the prostate segmentation results of MRI and TRUS images, and then use the non-rigid CPD algorithm to align the surface points of the MR grid to the surface points of the US. Among them, ω = 0, which represents the weight of the noise distribution, β = 100, which represents the width of the smoothing Gaussian filter to enhance consistency, and γ = 3.3, which controls the maximum likelihood fitting and deformation field regularization. Then, the neo-Hookean hyperelastic model with Poisson's ratio v = 0.49 and Young's modulus E = 5kPa is used to model the prostate, and finite element analysis is performed using FEBio 4.5.0. Finally, the finite element results are linearly interpolated to obtain a voxel-level dense deformation field with the same size and pixel spacing as the image pair to be registered, thereby obtaining a biomechanical prior.

[0040] Step 3, use the VAE-based generative model to learn the biomechanical prior as a deformation manifold in the latent space. The encoder consists of a 4-layer convolutional neural network, whose output is flattened into a one-dimensional vector by the Flatten layer, and then reparameterized into a Gaussian distribution constructed by mean μ and variance σ. The dimension of the latent space where the distribution is located is set to 64. After that, the latent vector z is restored to its original shape by the Unflatten layer and then input into the decoder to reconstruct the corresponding deformation field. The decoder consists of 4 layers of deconvolution and contains two jump connections and upsampling to make full use of multi-scale information. We set the relevant parameters to α = 10, β = 1, and use the Adam optimizer with a learning rate of 0.0001 for iterative optimization. The weights are continuously updated, and the loss value gradually decreases. Finally, the VAE-based prostate deformation field generation model is trained.

[0041] Step 4, first pass The trained VAE decoder converts it into the corresponding initial deformation field, and then uses the spatial transformer network STN to act on the MRI prostate segmentation image to obtain the deformed MRI prostate segmentation image. Then the Dice similarity loss of MR and US prostate segmentation is calculated, and the optimization process of z in the latent space is driven by back propagation until the L calculated twice in a row is reg When the value is less than the set threshold ∈=0.0001, the optimization process stops and the final deformation field and the MRI image deformed by the deformation field are output. The method of the present invention is compared with other existing methods. The results are as follows: Figure 2 The test set has a total of 122 cases. The dice coefficient values, 95% Hausdorf distance (HDE), average distance (ASD) and Jacobian determinant (JD) of different methods are shown in Table 1 (mean ± variance).

[0042] Table 1

[0043]

[0044] The experimental results of this embodiment show that the method of the present invention has higher registration accuracy, better generalization, and can maintain the biomechanical rationality of internal deformation of organs; the method we proposed does not require parameter adjustment and can adapt to different deformation conditions, proving that it has higher flexibility than the explicit regularization method and is therefore more practical.

[0045] References:

[0046] [1]Myronenko A,Song X.Point Set Registration:Coherent Point Drift[J].IEEE transactions on pattern analysis and machine intelligence,2010,32(12):2262-2275.

[0047] [2]Maas SA, Ellis BJ, Ateshian GA, et al. FEBio:finite elements for biomechanics[J].2012.

[0048] [3]Kingma D P.Auto-encoding variational bayes[J].arXiv preprintarXiv:1312.6114,2013.

[0049] [4]Jaderberg M,Simonyan K,Zisserman A,et al.Spatial transformernetworks[J].Advances in neural information processing systems,2015,28.

[0050] [5]Si H,TetGen A.Aquality tetrahedral mesh generator and a 3ddelaunay triangulator[J].Cited on,2009,61.

[0052] [6]Bonet J,Wood R D.Nonlinear continuum mechanics for finite elementanalysis[M].Cambridge university press,1997。

Claims

1. A biomechanically generated prostate MRI-TRUS deformable image registration method, characterized in that: Generative models are used to learn domain-invariant prostate biomechanical deformation laws, and registration is achieved through optimization in latent space; including: embedding the deformation laws obtained in finite element biomechanical priors into low-dimensional manifolds in latent space through VAE-based generative models, so that the generative model can generate biomechanically reasonable deformations without explicit regularization terms; then, through shape matching, optimization is performed on the learned manifold to find the deformation field that best suits the image to be registered. The specific steps are as follows: Step 1: preprocess the image, including normalizing the image intensity and performing histogram equalization; resampling all voxels to a uniform size; cropping the image to the same size centered on the region of interest; and finally performing preliminary rigid alignment. Step 2, obtain biomechanical priors; including, using TetGen to grid the prostate segmentation results of MRI and TRUS images, and then using the non-rigid coherent point drift (CPD) algorithm to align the surface points of the MR grid to the surface points of the US; specifically, each point m in the MR prostate surface point set is regarded as the centroid of the Gaussian mixture model, and the surface point of the US prostate represents the observation value x of the GMM; all components of the GMM are set to follow the mean y m That is, the position and variance of the MR point are σ 2 Gaussian distribution with equal membership probability Then the probability density function of GMM is expressed as: Where M represents the number of MRI surface points, ω is the weight of the uniform distribution weight added to represent noise, 0≤1≤ω; Then, the iterative expectation maximization algorithm is used to maximize the likelihood E of the GMM observations, as shown below: Where n represents the points in the TRUS point set, and its number is N; in solving the final θ and σ 2 After that, the deformed position of the moving point set can be calculated; Next, the neo-Hookean hyperelastic model is used to model the prostate. The strain energy function of the neo-Hookean material is defined as follows: Where λ is the first Lamé parameter and μ is the shear modulus or second Lamé parameter, both of which are material properties and can be calculated from Young's modulus E and Poisson's ratio v; I1 = trace(FF T ) and J = det(F) are invariants, where F is the deformation gradient tensor, through The calculation results show that, where X and x represent the material point positions before and after deformation, respectively; FEBio 4.5.0 was used to perform finite element analysis, and the surface point displacement field obtained from MRI to TRUS was used as the boundary strip; in the calculation, the load on the prostate was considered to be applied smoothly, and the boundary conditions were applied in five steps to obtain an incremental iterative solution, that is, the amplitude of each step was increased by one-fifth of the total boundary displacement in order to promote convergence and avoid unreasonable deformation; Finally, the results of the finite element analysis are linearly interpolated to obtain a voxel-level dense deformation field with the same size and pixel spacing as the image pair to be registered; Step 3, using a generative model to learn the biomechanical prior as a deformation manifold in the latent space; using a generative model based on a variational autoencoder (VAE), the deformation field in step 2 is implicitly embedded into the deformation manifold in the latent space of the generative model, thereby compactly representing the deformation during prostate surgery in a low-dimensional space; specifically, the generative model consists of an encoder and a multi-scale decoder; the encoder extracts features from the input deformation field φ and embeds it into a low-dimensional latent space, which is then represented by a Gaussian distribution constructed by a mean μ and a variance σ through reparameterization; the decoder reconstructs the latent vector z sampled from the distribution into a deformation field The skip connection and upsampling connection can make full use of the multi-scale information in the process; the first-order spatial derivative of the deformation field is added to the reconstruction term of the decoder. and To reduce the impact of rigid translation, the training loss function of the generated model is as follows: Among them, p(z)~N(0,I) represents the prior Gaussian distribution, q θ represents the encoder with parameter θ, D KL represents the KL divergence; the first two terms of the function are used to penalize the reconstruction error between the input and output in all directions, where α is a hyperparameter that controls the weight between different reconstruction terms; the latter term forces the learned distribution to be as similar as possible to the prior distribution, and β is a hyperparameter that weighs the reconstruction quality and the degree of regularity of the latent space; After training, the VAE decoder Reconstruct the z sampled from the latent space into the corresponding deformation field, that is, Therefore, the biomechanical deformation law of the prostate is implicitly embedded in the output of the model, making Ability to generate biomechanically plausible deformations without explicit regularization terms; Step 4: Deformable registration is achieved through optimization in latent space; first, by z is initialized by random sampling in , and converted into the corresponding initial deformation field by the trained VAE decoder, and then the spatial transformer network (STN) is used to act on the prostate segmentation image of MRI to obtain the deformed MRI prostate segmentation image; then the Dice similarity loss of MR and US prostate segmentation is calculated, and the optimization process of z in the latent space is driven by back propagation; the loss function of the image in this process is: Among them, s f and m Represents TRUS and MRI prostate segmentation, ° represents the spatial transformer network (STN) ] Completed spatial transformation; in addition, in order to ensure that z is within the range of the learned manifold when searching in the latent space, an L2 regularization term is added to encourage z to be in a high probability area under the corresponding Gaussian prior, and the final loss function is: Among them, ρ controls the possibility that z falls within the Gaussian prior range; by minimizing the above formula, we can find the potential vector z that can be mapped to the deformation field that best suits the input image along the manifold in the latent space; until two consecutive calculations of L reg When it is less than the set threshold ∈, the optimization process stops, that is, where i represents the optimization iteration step.

2. The prostate MRI-TRUS deformable image registration method according to claim 1, characterized in that: In the above formula, the parameters involved are specifically set as follows: ω=0, α=10, β=1, μ=0.001, ∈=0.0001.

3. The prostate MRI-TRUS deformable image registration method according to claim 1, characterized in that: In the generative model, the encoder consists of a 4-layer convolutional neural network, whose output is flattened into a one-dimensional vector by the Flatten layer and then reparameterized into a Gaussian distribution constructed by a mean μ and a variance σ; then, the latent vector z is restored to its original shape by the Unflatten layer and then input into the decoder to reconstruct the corresponding deformation field; the decoder consists of 4 layers of deconvolution and contains two stages of skip connections and upsampling to make full use of multi-scale information.

Images