A Multi-Temporal Medical Image Registration Method Based on Multiple Inputs and Spatial Attention

Abstract

This invention discloses a multi-temporal medical image registration method based on multi-branch input and spatial attention, belonging to the field of medical image registration technology. This method employs multi-branch input to ensure that while learning the features of the image pair (floating and fixed images), it fully utilizes the spatial semantic information of each image pair, overcoming the limitations of deep learning-based methods in understanding spatial semantic information during multi-temporal medical image registration. The method introduces an LK module to increase the effective receptive field, enabling it to handle not only long-distance spatial information but also capture and fuse spatial information at a finer scale. The method also creates a VSAW module, assigning different weights to pixels in different regions of the medical image, which can more accurately capture structural features in the image and improve the registration accuracy of multi-temporal malignant tumor CT images with large-scale deformation.

Description

Technical Field

[0001] This invention relates to the field of medical image registration technology, specifically to a multi-temporal medical image registration method based on multiple inputs and spatial attention. Background Technology

[0002] The incidence of malignant tumors is increasing year by year and showing a trend towards affecting younger people. For patients with intermediate and advanced cancer, the current common treatment regimen is external beam radiation therapy combined with brachytherapy and concurrent chemoradiotherapy. Brachytherapy uses a very small radiation source placed directly on the surface of the tumor or inserted into the tumor to deliver a high dose of radiation. Multiple radiation sessions are performed during brachytherapy, with a CT scan taken before each session. To ensure that the tumor receives a high dose of radiation while surrounding normal organs receive only a low dose, doctors need to assess changes in the tumor and the radiation dose already received by organs at risk. Due to differences in patient positioning and organ movement at different examination times, significant shape differences (large-scale deformation) occur in the CT images of the same tissue or organ from different examination times, leading to errors in the superposition of radiation doses. Therefore, achieving registration between multiple temporal images during brachytherapy is of practical significance for doctors.

[0003] Existing medical image registration methods are mainly divided into traditional methods and deep learning-based methods. Traditional medical image registration methods typically optimize an objective function for each image pair, continuously updating the transformation parameters required by the objective function based on a defined similarity metric until the similarity metric reaches its optimum. For images with large datasets and significant organ deformation, this is computationally time-consuming and requires substantial computational resources. Currently, deep learning-based methods are mostly unsupervised registration methods based on U-shaped structures and weakly supervised registration methods incorporating anatomical structures. These methods quickly calculate the deformation field by directly evaluating the function, resulting in higher computational efficiency compared to traditional medical image registration methods. However, deep learning-based methods typically use a pair of images stitched together as input. This input method makes it difficult to fully represent the independent anatomical information within each image pair, limiting the ability of deep learning-based methods to understand the spatial semantic information of multi-temporal medical image registration. Furthermore, unsupervised registration methods based on U-shaped structures struggle to capture global deformation information when dealing with images with large-scale deformation, leading to uneven spatial information distribution and making it difficult to maintain structural consistency between images. Summary of the Invention

[0004] In view of the shortcomings of the prior art, the purpose of this invention is to provide a multi-temporal medical image registration method based on multiple inputs and spatial attention, which aims to improve the registration accuracy of multi-temporal malignant tumor CT images with large-scale deformation while reducing registration time and resource costs.

[0005] The technical solution of this invention is:

[0006] A multi-temporal medical image registration method based on multiple inputs and spatial attention, comprising the following steps:

[0007] Step 1: Read the fixed image OX, floating image OY, label mask X corresponding to the fixed image and label mask Y corresponding to the floating image from the training set, and use Lanczos interpolation to adjust the resolution of these four images to H×W, with T channels, represented as OX(T×H×W), OY(T×H×W), X(T×H×W), Y(T×H×W);

[0008] Step 2: Multiply the read label mask X with the fixed image OX, and multiply the label mask Y with the floating image OY to obtain the region of interest images M(T×H×W) and N(T×H×W);

[0009] Step 3: First, stitch together images M and N to obtain image MN. Then, extract features from these three branches through the encoding stage to obtain feature maps E1 to E5.

[0010] Step 4: Perform feature extraction and activation operations on E5 to obtain feature map D5;

[0011] Step 5: Combining the feature maps E1 to E4 obtained in Step 3, D5 is decoded to recover the features and obtain two velocity fields F'. xy (2×H×W) and F' yx (2×H×W);

[0012] Step 6: Create the VSAW module and process the velocity field F' obtained in Step 5. xy and F' yx The velocity field F was obtained by applying the VSAW module respectively. xy (2×H×W) and F yx (2×H×W);

[0013] Step 7: Based on the velocity field F xy With F yx Synthesizing the deformation field from M to N and the deformation field from N to M pass To distort X into a deformed image And the image that distorts OX into a deformed form. pass Distort Y into a deformed image And the image that distorts OY into a deformed form.

[0014] Furthermore, according to the multi-temporal medical image registration method based on multiple inputs and spatial attention, the feature extraction of the images of the three branches M, N, and MN after the encoding stage includes the following steps:

[0015] Step 3.1: Perform feature extraction operations on M, N and MN respectively to obtain feature map S1 (C×H×W) of M, feature map S2 (C×W×H) of N and feature map S3 (P×W×H) of MN. Then concatenate feature maps S1, S2 and S3 to obtain P1. Perform feature extraction operations on P1 to obtain feature map E1.

[0016] Step 3.2: Perform feature extraction on E1 to obtain feature map P2, and then pass P2 through the existing LK module to obtain feature map E2;

[0017] Step 3.3: Following the method in Step 3.2, extract features from E2 to obtain feature map E3;

[0018] Step 3.4: Following the method in Step 3.2, extract features from E3 to obtain feature map E4;

[0019] Step 3.5: Perform feature extraction and activation operations on E4 to obtain feature map E5.

[0020] Furthermore, according to the multi-temporal medical image registration method based on multiple inputs and spatial attention, step 5 involves performing feature reconstruction on D5 through a decoding stage to obtain two velocity fields F. xy and F yx The method includes: performing feature restoration on D5 to obtain feature map D4; performing feature restoration on D4 to obtain feature map D3; restoring features from D3 to obtain feature map D2; restoring features from D2 to obtain feature map D1; and performing two convolution operations on D1 to obtain the velocity field F'. xy (2×H×W) and F' yx (2×H×W).

[0021] Furthermore, according to the multi-temporal medical image registration method based on multiple inputs and spatial attention, in step 6, the velocity field F' is... xy Applying the VSAW module to obtain F xy (2×H×W) includes the following steps:

[0022] Step 6.1.1: For the velocity field F' xy We calculate the weight of one channel dx to obtain the weight map W6, and multiply W6 with dx to obtain the new feature representation W7 (1×H×W);

[0023] Step 6.1.2: For the velocity field F' xyThe other channel dy is used to calculate the weight to obtain the weight map H6, and H6 is multiplied by dy to obtain the new feature representation H7;

[0024] Step 6.1.3: Combine W7 and H7 to obtain a new velocity field F. xy (2×H×W).

[0025] Furthermore, according to the multi-temporal medical image registration method based on multiple inputs and spatial attention, the specific process of step 6.1.1 is as follows:

[0026] (1) Perform two parallel convolution operations on dx with kernel size of 3×3 and stride of 1, change the number of feature map channels in the dx direction, and obtain feature maps W1 and W2 respectively;

[0027] (2) Perform two parallel convolution operations on W1 and W2 with a kernel size of 3×3 and a stride of 1 to perform feature restoration operations to obtain W3 and W4;

[0028] (3) In order to enhance the representation ability of the features, the restored feature maps W3 and W4 are added to dx respectively to obtain feature maps W'3 and W'4, and then feature maps W'3 and W'4 are added to obtain feature map W5;

[0029] (4) Pass W5 through the Sigmoid activation function to obtain a weight map W6 with dx ranging from 0 to 1;

[0030] (5) Multiply dx by W6 to form a new channel representation W7.

[0031] Furthermore, based on the aforementioned multi-temporal medical image registration method based on multiple inputs and spatial attention, the obtained velocity field F is... xy With F yx The deformation field from M to N is synthesized using the following formula. and the deformation field from N to M

[0032]

[0033]

[0034] in, By F xy The scaling and squaring layers are used to obtain a representation of transforming M into the average shape of M and N; By F yx The scaling and squaring layers are used to obtain a representation of transforming N into the average shape of M and N.

[0035] Compared with the prior art, the present invention has the following beneficial effects:

[0036] 1. The method of the present invention adopts multiple inputs, including three branches: fixed image X, floating image Y, and image XY after stitching the fixed image and the floating image. This ensures that while learning the features of the image pair of floating image and fixed image, the spatial semantic information of each image pair is fully utilized, thus solving the limitation of deep learning-based methods in understanding the spatial semantic information of multi-temporal medical image registration.

[0037] 2. The method of this invention introduces an LK (Large kernel) module. The LK module uses four parallel paths to increase the effective receptive field, which can not only process long-distance spatial information, but also capture and fuse spatial information at a finer scale.

[0038] 3. The method of this invention creates a VSAW (Velocity-based Space Attention Weight) module, which assigns different weights to pixels in different regions of medical images, enabling more accurate capture of structural features in the images and improving the registration accuracy of multi-temporal malignant tumor CT images with large-scale deformation. Attached Figure Description

[0039] Figure 1 This is a flowchart of the multi-temporal medical image registration method based on multiple inputs and spatial attention in this embodiment;

[0040] Figure 2 This is a schematic diagram of the overall structure of the multi-temporal medical image registration method based on multiple inputs and spatial attention in this embodiment;

[0041] Figure 3 This is a diagram of the velocity field generation network structure in this embodiment;

[0042] Figure 4 This is a schematic diagram of the structure of the LK module introduced in this embodiment;

[0043] Figure 5 This is a schematic diagram of the VSAW module used in this embodiment;

[0044] Figure 6(A) is a floating image and its corresponding label mask; (B) is a display of the registration result of (A) using the existing medical image registration method Affine; (C) is a display of the registration result of (A) using the existing medical image registration method Ac-RgNet; (D) is a display of the registration result of (A) using the existing medical image registration method DIRNet; (E) is a display of the registration result of (A) using the existing medical image registration method VoxelMorph; (F) is a display of the registration result of (A) using a multi-temporal medical image registration method based on multiple inputs and spatial attention; (G) is a fixed image and its corresponding label mask. Detailed Implementation

[0045] To facilitate understanding of this application, a more comprehensive description of this application will be provided below with reference to the accompanying drawings.

[0046] like Figure 1 and Figure 2 As shown, the multi-temporal medical image registration method based on multiple inputs and spatial attention in this embodiment includes the following steps:

[0047] Step 1: Read the fixed image OX, floating image OY, label mask X corresponding to the fixed image and label mask Y corresponding to the floating image from the training set, and use Lanczos interpolation to adjust the resolution of these four images to H×W, with T channels, represented as OX(T×H×W), OY(T×H×W), X(T×H×W), Y(T×H×W);

[0048] In this embodiment, the fixed image OX, the floating image OY, the label mask X corresponding to the fixed image, and the label mask Y corresponding to the floating image are all read using Lanczos interpolation to adjust the resolution of the images to 288×368, with 1 channel, represented as 1×288×368. The first dimension is 1, representing the number of channels for OX, OY, X, and Y; the second dimension is 288, representing the height of these four images; and the third dimension is 368, representing the width of these four feature maps.

[0049] Step 2: Multiply the read label mask X with the fixed image OX, and multiply the label mask Y with the floating image OY to obtain the region of interest images M(T×H×W) and N(T×H×W);

[0050] In this embodiment, the image obtained by multiplying the tag mask X with the fixed image OX is M (1×288×368), and the image obtained by multiplying the tag mask Y with the floating image OY is N (1×288×368).

[0051] Step 3: First, stitch images M and N together to obtain image MN. Then, perform feature extraction on these three branches of images through the encoding stage to obtain feature maps E1 to E5.

[0052] The feature extraction process for the images of branches M, N, and MN after the encoding stage includes the following steps:

[0053] Step 3.1: As Figure 3 As shown, feature extraction operations are performed on M, N and MN respectively to obtain feature map S1 (C×H×W) of M, feature map S2 (C×W×H) of N and feature map S3 (P×W×H) of MN. Then, feature maps S1, S2 and S3 are concatenated to obtain P1. The concatenated feature map is then subjected to feature extraction operations to obtain feature map E1.

[0054] In this example, M undergoes a convolution operation with a kernel size of 3×3 and a stride of 1, followed by an activation operation, to obtain feature maps S1 (4×288×368), S2 (4×288×368), and S3 (8×288×368). Feature maps S1, S2, and S3 are then concatenated. The concatenated feature map P1 is then subjected to a convolution operation with a kernel size of 3×3 and a stride of 1, followed by an activation operation, to obtain feature map E1 (16×288×368).

[0055] Step 3.2: Perform feature extraction on E1 to obtain feature map P2, and then pass P2 through the existing LK module to obtain feature map E2;

[0056] The LK module mentioned above is as follows: Figure 4 The diagram shows four parallel layers: a k×k LK convolutional layer, a 3×3 convolutional layer, a 1×1 convolutional layer, and an identity shortcut. The LK module in this example includes a 7×7 LK convolutional layer, a 3×3 convolutional layer, a 1×1 convolutional layer, and an identity shortcut.

[0057] In this embodiment, E1 undergoes a convolution operation with a kernel size of 3×3 and a stride of 2, followed by an activation operation, to obtain feature map P2 (32×144×184). P2 is then processed by the LK module for feature extraction to obtain feature map E2 (32×144×184).

[0058] Step 3.3: Following the method in Step 3.2, extract features from E2 to obtain feature map E3;

[0059] Following the same method as in step 3.2, E2 is first subjected to convolution and activation operations to obtain feature map P3 (64×72×92). Then, P3 is subjected to LK module for feature extraction to obtain feature map E3 (64×72×92).

[0060] Step 3.4: Following the method in step 3.2, extract features from E3 to obtain feature map E4.

[0061] First, perform convolution and activation operations on E3 to obtain feature map P4 (128×36×46). Then, perform feature extraction on P4 through the LK module to obtain feature map E4 (128×36×46).

[0062] Step 3.5: Perform feature extraction and activation operations on E4 to obtain feature map E5.

[0063] E4 undergoes a convolution operation with a kernel size of 3×3 and a stride of 2, followed by an activation operation, resulting in feature map E5 (256×18×23).

[0064] Step 4: Perform feature extraction and activation operations on E5 to obtain the feature map.

[0065] E5 undergoes a convolution operation with a kernel size of 3×3 and a stride of 1, followed by an activation operation, resulting in feature map D5 (128×18×23).

[0066] Step 5: Combining the feature maps E1 to E4 obtained in Step 3, D5 is decoded to recover the features and obtain two velocity fields F'. xy (2×H×W) and F' yx (2×H×W);

[0067] Step 5.1: Perform feature reconstruction on D5 to obtain feature map D4;

[0068] In this embodiment, the detailed process of feature restoration of D5 to obtain D4 is as follows: First, D5 is upsampled to change the feature map size, resulting in feature map D'4 (128×36×46); then D'4 is concatenated with E4 to obtain feature map Q4 (256×36×46); then Q4 is subjected to two sets of convolution operations with kernel size of 3×3 and stride of 1, and activation operations to change the number of feature map channels, resulting in feature map D4 (64×36×46).

[0069] Step 5.2: Perform feature restoration on D4 to obtain feature map D3.

[0070] The feature map D4 is upsampled to change its size, resulting in feature map D'3 (64×72×92). D'3 is then concatenated with E3 to obtain feature map Q3 (128×72×92). Q3 is then subjected to two sets of convolution operations with kernels of size 3×3 and stride of 1, and activation operations to change the number of channels in the feature map, resulting in feature map D3 (32×72×92).

[0071] Step 5.3: Restore the features of D3 to obtain feature map D2.

[0072] The feature map D3 is upsampled to change its size, resulting in feature map D'2 (32×144×184). D'2 is then concatenated with E2 to obtain feature map Q2 (64×144×184). Q2 is then subjected to two sets of convolution operations with kernels of size 3×3 and stride of 1, and activation operations to change the number of channels in the feature map, resulting in feature map D2 (32×144×184).

[0073] Step 5.4: Restore the features of D2 to obtain feature map D1.

[0074] The feature map D2 is upsampled to change its size, resulting in feature map D'1 (32×288×368). D'1 is then concatenated with E1 to obtain feature map Q1 (48×288×368). Q1 is then subjected to two sets of convolution operations with kernels of size 3×3 and stride of 1, and activation operations to change the number of channels in the feature map, resulting in feature map D1 (32×288×368).

[0075] Step 5.5: Perform two convolution operations on D1 to obtain the velocity field F'. xy (2×H×W) and F' yx (2×H×W);

[0076] D1 is subjected to two convolution operations with a kernel size of 5×5 and a stride of 1 to generate the velocity field F'. xy (2×H×W) and F' yx (2×H×W), using the softsign activation function to limit the range of the velocity field to [-1, 1], and multiplying by a constant c to normalize the velocity field in the range [-c, c].

[0077] In this embodiment, D1 is subjected to two convolution operations with a kernel size of 5×5 and a stride of 1 to generate a velocity field F'. xy (2×288×368) and F' yx (2×288×368), the range of the velocity field is limited to [-1, 1] using the softsign activation function, and the velocity field is normalized in the range [-1, 1] by multiplying by a constant 1.

[0078] Step 6: Create the VSAW module and process the velocity field F' obtained in Step 5. xy and F' yx The velocity field F was obtained by applying the VSAW module respectively. xy (2×H×W) and F yx (2×H×W);

[0079] In this embodiment, the velocity field F' xy (2×288×368) and F' yx The first dimension in (2×288×368) is 2, representing F' xy F' yx The number of channels is given, and two channels are represented as dx and dy. dx and dy are processed by convolutional layers to obtain their respective weight maps, where the weight mask values ​​range from [0,1]. Then, dx and dy are multiplied by their respective weight mask values ​​to generate a new velocity field F. xy (2×H×W) and F yx (2×H×W).

[0080] Step 6.1: Set F' xy F is obtained through the VSAW module. xy (2×H×W).

[0081] Step 6.1.1: Calculate the weights of channel dx to obtain the weight map W6, and multiply it with dx to obtain a new feature representation W7 (1×H×W);

[0082] In this embodiment, as Figure 5 As shown, the specific process of obtaining W7 (1×288×368) from dx through the VSAW module is as follows:

[0083] (1) F' xy The first channel is represented as dx (1×288×368). Two convolution operations with a kernel size of 3×3 and a stride of 1 are performed on dx. The number of feature map channels in the dx direction is changed to obtain feature maps W1 (16×288×368) and W2 (16×288×368) respectively.

[0084] (2) Perform two parallel convolution operations on W1 and W2 with a kernel size of 3×3 and a stride of 1 to perform feature restoration operations to obtain W3 (1×288×368) and W4 (1×288×368);

[0085] (3) In order to enhance the representation ability of features, the restored feature maps W3 and W4 are added to dx respectively to obtain feature maps W'3 (1×288×368) and W'4 (1×288×368), and then feature maps W'3 and W'4 are added to obtain feature map W5 (1×288×368), which improves the connection and dependency between features;

[0086] (4) Pass W5 through the Sigmoid activation function to obtain a weighted map W6 (1×288×368) with dx ranging from 0 to 1;

[0087] (5) Multiply dx by W6 to form a new channel representation W7 (1×288×368).

[0088] Step 6.1.2: Calculate the weights for channel dy to obtain the weight map H6, and multiply it by dy to obtain the new feature representation H7; following the same method as in step 5.1.1, calculate the weights for F'. xy The second channel is represented as dy(1×288×368). First, calculate dy...

[0089] Perform two parallel convolution operations to obtain feature maps H1 (16×288×368) and H2 (16×288×368); then perform convolution operations on H1 and H2...

[0090] H2 is subjected to feature reduction operations to obtain H3 (1×288×368) and H4 (1×288×368); the reduced feature maps H3 and H4 are added to dy to obtain feature maps H'3 (1×288×368) and H'4 (1×288×368), respectively; feature maps H'3 and H'4 are then added to obtain feature map H5 (1×288×368), which improves the connection and dependency between features; H5 is processed by the Sigmoid activation function to obtain the weight map H6 (1×288×368) with dy ranging from 0 to 1; dy is multiplied by H6 to form a new channel representation W7 (1×288×368).

[0091] Step 6.1.3: Combine W7 and H7 to obtain a new velocity field F. xy (2×288×368).

[0092] Step 6.2: Set F' yx F is obtained through the VSAW module. yx .

[0093] The new velocity field F is obtained using the same method as in step 5.1. yx (2×288×368).

[0094] Step 7: Analyze the obtained velocity field F xy With F yxThe deformation field from M to N is synthesized using the following formula. and the deformation field from N to M

[0095]

[0096]

[0097] in, By F xy The scaling and squaring layers are used to obtain a representation of transforming M into the average shape of M and N; By F yx The scaling and squaring layers are used to obtain a representation of transforming N into the average shape of M and N.

[0098] The existing Spatial Transformation Network (STN) is used to apply deformation operations to the X, OX, Y, and OY axes of the image. To distort X into a deformed image To distort OX into a deformed image pass Distort Y into a deformed image OY is distorted into a deformed image.

[0099] When verifying the method of the present invention, the calculation is first performed. and Y and The Dice coefficient between X and and OY and The similarity between the model and OX, as well as the smoothness of the deformation field, are used to update the parameters of the method of the present invention.

[0100] In this embodiment, and Y and The Dice coefficient between X and X is expressed as and OY and The similarity between and OX is represented as The smoothness of the deformation field is expressed as Therefore, the loss function of the method of the present invention is expressed as:

[0101]

[0102] Where α represents the regularization parameter for the smoothness of the deformation field, and in this embodiment, the effect is best when α = 0.3. Φ is used to uniformly represent... and

[0103] for In its specific application, this invention uses normalized cross-correlation (NCC) to represent the local cross-correlation between the fixed image OX and the floating image OY. The higher the NCC, the better the alignment effect between OX and OY. The calculation formula is as follows:

[0104]

[0105]

[0106] Among them, OX i This represents the value of the i-th pixel in the fixed image OX. Let σ(OX) represent the mean of the fixed image OX, and let σ(OX) represent the standard deviation of the fixed image. j This represents the j-th pixel value of the floating image OY. Let σ(OY) represent the mean of the floating image OY, and let σ(OY) represent the standard deviation of the floating image.

[0107] The smoothing loss function is specifically expressed as follows:

[0108]

[0109] The goal of this loss function is to minimize the smoothness of the deformation field in both the horizontal and vertical directions, achieved by squaring the smoothness in each direction and then taking the mean. i dy represents the difference between two adjacent columns of pixels in the horizontal direction of the deformation field. i This represents the difference in the deformation field between two adjacent rows of pixels in the vertical direction. x n represents the total number of pixels used to calculate the smoothness of the deformation field in the horizontal direction. y This represents the total number of pixels used when calculating the smoothness of the deformation field in the vertical direction.

[0110] Finally, by calculating the above loss function, it is determined whether the loss value has reached the optimal value. The gradient value is backpropagated to update the parameters of the method of this invention, and the updated training weights are saved in the specified location.

[0111] Compared with existing medical image registration methods such as Affine, Ac-RgNet, DIRNet, and VoxelMorph, the method of this invention performs better on medical images with large-scale deformation, and shows superior results in Dice similarity coefficient (DSC), 95% Hausdorff distance (HD95), and surface overlap (SO) evaluation metrics.

[0112] Table 1. DSC, HD95, and SO values ​​for different methods

[0113]

[0114] As shown in Table 1, the method of the present invention has significant improvements over other methods in terms of DSC and HD95 indices. Specifically, the DSC index reaches 0.8286, the HD95 index reaches 3.7700, and the SO index reaches 0.7611.

[0115] Next, according to Figure 6 Conduct qualitative analysis. Figure 6 In the diagram, (A) shows the floating image and its corresponding label mask; (B) shows the registration result of (A) using the existing medical image registration method Affine; (C) shows the registration result of (A) using the existing medical image registration method Ac-RgNet; (D) shows the registration result of (A) using the existing medical image registration method DIRNet; (E) shows the registration result of (A) using the existing medical image registration method VoxelMorph; (F) shows the registration result of (A) using a multi-branch input and spatial attention-based multi-temporal medical image registration method; and (G) shows the fixed image and its corresponding label mask. It can be seen that the registration result using the method of this invention shows a higher similarity to image (G).

[0116] In summary, the multi-temporal medical image registration method based on multi-branch input and spatial attention in this invention achieves more accurate results in the registration of multi-temporal medical images. This is because other methods do not have sufficient deformation capability to register images with large organ deformation, while this invention, through multi-branch input and LK module, can make full use of the spatial semantic information of each image pair and enable the network to accept a sufficiently large receptive field. In addition, to ensure the existence of the inverse deformation field, this invention adds a VSAW module and applies spatial attention weights to the velocity field, which can obtain more accurate registration results.

[0117] It should be understood that, inspired by the technical concept of this invention, those skilled in the art can make various improvements or modifications based on the above content without departing from the scope of this invention, and these modifications still fall within the protection scope of this invention.

Claims

1. A multi-temporal medical image registration method based on multiple inputs and spatial attention, characterized in that, The method includes the following steps: Step 1: Read the fixed image OX, floating image OY, label mask X corresponding to the fixed image and label mask Y corresponding to the floating image from the training set, and use Lanczos interpolation to adjust the resolution of these four images to H×W, with T channels, represented as OX(T×H×W), OY(T×H×W), x(T×H×W), Y(T×H×W); Step 2: Multiply the read label mask x with the fixed image OX, and multiply the label mask Y with the floating image OY to obtain the region of interest images M(T×H×W) and N(T×H×W); Step 3: First, stitch together images M and N to obtain image MN. Then, extract features from these three branches through the encoding stage to obtain feature maps E1 to E5. Step 4: Perform feature extraction and activation operations on E5 to obtain feature map D5; Step 5: Combining the feature maps E1 to E4 obtained in Step 3, D5 is decoded to restore the features and obtain two velocity fields F′. xy (2×H×W) and F′ yx (2×H×W); Step 6: Create the VSAW module and process the velocity field F′ obtained in Step 5. xy and F′ yx The velocity field F was obtained by applying the VSAW module respectively. xy (2×H×W) and F yx (2×H×W); Step 7: Based on the velocity field F xy With F yx Synthesizing the deformation field from M to N and the deformation field from N to M pass To distort X into a deformed image And the image that distorts OX into a deformed form. pass Distort Y into a deformed image And the image that distorts OY into a deformed form.

2. The multi-temporal medical image registration method based on multiple inputs and spatial attention as described in claim 1, characterized in that, The feature extraction process for the images of branches M, N, and MN after the encoding stage includes the following steps: Step 3.1: Perform feature extraction operations on M, N and MN respectively to obtain feature map S1 (C×H×W) of M, feature map S2 (C×W×H) of N and feature map S3 (P×W×H) of MN. Then concatenate feature maps S1, S2 and S3 to obtain P1. Perform feature extraction operations on P1 to obtain feature map E1. Step 3.2: Perform feature extraction on E1 to obtain feature map P2, and then pass P2 through the existing LK module to obtain feature map E2; Step 3.3: Following the method in Step 3.2, extract features from E2 to obtain feature map E3; Step 3.4: Following the method in Step 3.2, extract features from E3 to obtain feature map E4; Step 3.5: Perform feature extraction and activation operations on E4 to obtain feature map E5.

3. The multi-temporal medical image registration method based on multiple inputs and spatial attention as described in claim 1, characterized in that, In step 5, D5 is decoded to obtain two velocity fields F. xy and F yx The method includes: performing feature restoration on D5 to obtain feature map D4; performing feature restoration on D4 to obtain feature map D3; restoring features from D3 to obtain feature map D2; restoring features from D2 to obtain feature map D1; and performing two convolution operations on D1 to obtain the velocity field F′. xy (2×H×W) and F′ yx (2×H×W).

4. The multi-temporal medical image registration method based on multiple inputs and spatial attention as described in claim 1, characterized in that, In step 6, the velocity field F′ xy Applying the VSAW module to obtain F xy (2×H×W) includes the following steps: Step 6.1.1: For the velocity field F′ xy We calculate the weight of one channel dx to obtain the weight map W6, and multiply W6 with dx to obtain the new feature representation W7 (1×H×W); Step 6.1.2: For the velocity field F′ xy The other channel dy is used to calculate the weight to obtain the weight map H6, and H6 is multiplied by dy to obtain the new feature representation H7; Step 6.1.3: Combine W7 and H7 to obtain a new velocity field F. xy (2×H×W).

5. The multi-temporal medical image registration method based on multiple inputs and spatial attention as described in claim 4, characterized in that, The specific process of step 6.1.1 is as follows: (1) Perform two parallel convolution operations on dx with kernel size of 3×3 and stride of 1, change the number of feature map channels in the dx direction, and obtain feature maps W1 and W2 respectively; (2) Perform two parallel convolution operations on W1 and W2 with a kernel size of 3×3 and a stride of l to perform feature restoration operations to obtain W3 and W4; (3) In order to enhance the representation ability of the features, the restored feature maps W3 and W4 are added to dx respectively to obtain feature maps W′3 and W′4, and then feature maps W′3 and W′4 are added to obtain feature map W5; (4) Pass W5 through the Sigmoid activation function to obtain a weight map W6 with dx ranging from 0 to 1; (5) Multiply dx by W6 to form a new channel representation W7.

6. The multi-temporal medical image registration method based on multiple inputs and spatial attention as described in claim 1, characterized in that, The obtained velocity field F xy With F yx The deformation field from M to N is synthesized using the following formula. and the deformation field from N to M in, By F xy The scaling and squaring layers are used to obtain a representation of transforming M into the average shape of M and N; By F yx The scaling and squaring layers are used to obtain a representation of transforming N into the average shape of M and N.

Images