Tumor radiotherapy dose prediction method based on diffusion model

Abstract

The invention discloses a tumor radiotherapy dose prediction method based on a diffusion model. The implementation scheme is as follows: 1) data collection; 2) data preprocessing; 3) constructing a dose prediction model; 4) training a prediction model; and 5) dose prediction. According to the dose prediction model constructed by the invention, a three-branch encoder is adopted to efficiently extract key features of input information, and accurate fusion of multi-source input features is realized by means of a feature fusion module MFF; the designed full-scale jump connection based on deformable convolution can make full use of the feature distribution of input conditions on different scales, thereby significantly improving the accuracy and robustness of dose prediction. The method can quickly generate a dose distribution result rich in high-frequency information, and has important scientific value and application prospect in the field of precise tumor radiotherapy.

Description

Technical Field

[0001] The present invention belongs to the technical field of radiotherapy dose prediction, and particularly relates to a tumor radiotherapy dose prediction method based on a diffusion model. Background Art

[0002] Radiotherapy is a widely used tumor treatment method in clinical practice at present, which can significantly improve the survival rate and quality of life of patients. The design of radiotherapy plans is the core of the entire radiotherapy process, and its goal is to ensure that the planning target volume (PTV) receives sufficient dose irradiation while minimizing the radiation damage to organs at risk (OAR) to the greatest extent. However, the traditional radiotherapy plan design highly depends on the professional level and clinical experience of physicists, and requires repeated adjustment of the relevant parameters of the target area and organs at risk, which is not only time-consuming and laborious, but also there is a certain degree of uncertainty in the plan quality. Therefore, researching automatic radiotherapy plans and developing efficient and accurate tumor dose prediction algorithms are of great significance for improving radiotherapy efficiency and optimizing treatment plans.

[0003] In recent years, with the rapid development of deep learning technology, many deep learning-based dose prediction methods have been successively proposed and have made remarkable progress. However, the existing methods still have certain limitations, such as fuzzy prediction results and loss of high-frequency information, and this high-frequency information is crucial for accurately depicting the ray path and dose attenuation characteristics. In addition, the existing methods fail to fully integrate multi-source information such as CT images, PTV contour maps, OAR contour maps, and non-intensity-modulated dose images, and still lack attention and utilization of different-scale features. Summary of the Invention

[0004] The present invention fully considers the problems existing in the prior art, and its purpose is to propose a tumor radiotherapy dose prediction method based on a diffusion model. By optimizing the model structure and parameter configuration, high-precision prediction of tumor radiotherapy dose is achieved.

[0005] I. Technical Principle

[0006] Existing dose prediction algorithms generally have problems such as fuzzy prediction results and loss of high-frequency information, which affect the accuracy and quality of dose distribution. To solve these problems and generate high-quality dose distribution results, the present invention constructs a tumor radiotherapy dose prediction method based on a diffusion model. To fully extract and utilize the complementary features of multi-source inputs, the present invention designs a three-branch encoder and a feature fusion module MFF. First, the key features of multi-source inputs are extracted through the three-branch encoder, and then the deep fusion of features is achieved through the MFF module, thereby enhancing the model's attention to key regions and improving the accuracy of dose prediction. To effectively utilize the feature information of input data at different scales, the present invention designs a full-scale skip connection based on deformable convolution, enabling the full integration of high- and low-frequency information at different scales, thereby significantly improving the accuracy and precision of dose prediction.

[0007] II. According to the above principle, the present invention is implemented through the following solutions:

[0008] A tumor radiotherapy dose prediction method based on a diffusion model, comprising the following steps:

[0009] (1) Collect tumor radiotherapy plan data to form a database.

[0010] (2) Data preprocessing, specifically including the following steps:

[0011] (2-a) Obtain the patient's CT image, the contour image of the target PTV, the contour image of the organ at risk OAR, the unmodulated dose image, and the clinical dose image from the database described in step (1).

[0012] (2-b) Calculate the distance image based on the target PTV contour image and the organ at risk OAR contour image obtained in step (2-a) The calculation formula is:

[0013]

[0014] where (i,j,k) is any voxel point in the CT image, represents the set of voxel points on the surface of the c-th PTV or OAR contour, the subscript c represents the total number of PTVs and OARs in the CT image, s c represents any voxel point within, ‖(i,j,k),s c ‖ represents the distance from the voxel point (i,j,k) to the point s c and min represents taking the minimum value.

[0015] (2-c) Sample all the images obtained in step (2-a) and the distance image obtained in (2-b) to the same size and perform normalization processing.

[0016] (3) Construct a dose prediction model, which specifically includes the following steps:

[0017] (3-a) Noise addition process: At each time step t, add Gaussian noise ε to the clinical dose image obtained in step (2) t to obtain the corresponding intermediate feature map y t , where t ∈ [0, T], and T is a positive integer representing the total number of time steps.

[0018] (3-b) Denoising process: Predict the noise added at each time step through a noise predictor and subtract it from y t to finally generate the predicted dose distribution map; The construction of the noise predictor specifically includes the following steps: Ⅰ. Construct a three-branch encoder: Each branch includes five modules. The first module consists of a convolutional layer, and the second to fifth modules each consist of two residual blocks and a downsampling layer; The three branches are respectively used to extract the features of x

[0019] , y CT , y t , x other , where x CT is the CT image, and x other includes the PTV contour and its distance image, the OAR contour and its distance image, and the non-intensity-modulated dose image; x CT , y t , x other are respectively processed by the first two modules of the three branches and then matrix-added to obtain the feature map Similarly, x CT , y t , x other are respectively processed by the first three modules of the three branches and then matrix-added to obtain the feature map are respectively processed by the first four modules of the three branches and then matrix-added to obtain the feature map x CT , y t , x other respectively obtain the feature maps f1, f3, and f2 after being processed by the three branches;

[0020] Ⅱ. Construct a feature fusion module MFF: The feature fusion module consists of an interactive attention block, a gated attention block, and a Mamba block; After f1 and f2 are processed by the interactive attention block, the result f IA and f3 are input into the gated attention block together, and the obtained result f GA is then processed by the Mamba block to obtain the fusion result f out of f1, f2, and f3;

[0021] Ⅲ. Construct the decoder: The decoder consists of five modules. The first to fourth modules are each composed of two residual blocks and an upsampling layer, and the fifth module is a convolutional layer; the f obtained from step (3-b)-Ⅱ out After being processed by the first module of the decoder, a feature map is obtained

[0022] Ⅳ. Full-scale skip connection: Resample the obtained in step (3-b)-Ⅰ to the same scale as the obtained in step (3-b)-Ⅲ, process them separately using deformable convolutional blocks, and then concatenate them. The concatenated result is successively processed by a deformable convolutional block and the second decoding block to obtain a feature map Similarly, resample the to the same scale as the and process them separately using deformable convolutional blocks and then concatenate them. The concatenated result is successively processed by a deformable convolutional block and the third decoding block to obtain a feature map Resample the to the same scale as the and process them separately using deformable convolutional blocks and then concatenate them. The concatenated result is successively processed by a deformable convolutional block and the second decoding block to obtain a feature map Then, after being processed by the fifth decoder, a prediction result is obtained Complete the construction of the noise predictor.

[0023] (4) Train the prediction model:

[0024] Train the prediction model. Use the AdamW optimizer to update the parameters of the model until the loss no longer decreases, and obtain the trained prediction model. The loss function used is:

[0025]

[0026] where M is the total number of voxel points.

[0027] (5) Dose prediction:

[0028] Use the dose prediction model trained in step (4) to perform dose prediction on the data to be tested, and obtain the dose prediction result.

[0029] The residual block in step (3-b) includes two convolutional blocks and a residual connection. The convolutional block consists of a normalization layer, an activation layer, and a convolutional layer.

[0030] The deformable convolutional block in step (3-b) consists of a normalization layer, a deformable convolutional layer, a SiLU activation function layer, and a residual connection.

[0031] Compared with the prior art, the present invention has the following advantages:

[0032] (1) The present invention constructs a tumor radiotherapy dose prediction method based on a diffusion model, effectively solving the problem of overly smooth prediction results in existing dose prediction models. By generating dose distribution images rich in high-frequency information, the accuracy of dose prediction is improved, and the efficiency of automatic plan formulation is significantly enhanced.

[0033] (2) The designed feature fusion module MFF of the present invention can fully fuse the features of multi-source inputs such as CT images, PTV contour maps, and OAR contour maps, realizing the effective utilization of complementary features, enabling the model to focus on important regions, and thus significantly improving the dose prediction accuracy in the PTV and OAR regions.

[0034] (3) The designed full-scale skip connection based on deformable convolution of the present invention can make full use of the different-scale features extracted by the encoder, effectively integrate its high- and low-frequency information, and thus significantly improve the accuracy of dose prediction. BRIEF DESCRIPTION OF THE DRAWINGS

[0035] Figure 1 Flowchart of a tumor radiotherapy dose prediction method based on a diffusion model according to an embodiment of the present invention;

[0036] Figure 2 Structure diagram of the dose prediction model according to an embodiment of the present invention;

[0037] Figure 3 Structure diagram of the noise predictor according to an embodiment of the present invention;

[0038] Figure 4 Structure diagram of the feature fusion module MFF according to an embodiment of the present invention;

[0039] Figure 5 Example diagram of the full-scale skip connection according to an embodiment of the present invention;

[0040] Figure 6 Comparison diagram of the dose prediction results according to an embodiment of the present invention and the prediction results of other methods. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0041] The following describes the specific embodiments of the present invention:

[0042] Example 1

[0043] Figure 1 As shown is a flowchart of a tumor radiotherapy dose prediction method based on a diffusion model according to an embodiment of the present invention, and the specific steps are as follows:

[0044] Step 1, collect tumor radiotherapy plan data to form a database.

[0045] Step 2, data preprocessing, specifically including the following steps:

[0046] (2-a) Obtain the patient's CT image, the contour image of the target PTV, the contour image of the organ at risk OAR, the non-intensity-modulated dose image, and the clinical dose image from the database described in Step 1.

[0047] (2-b) Calculate the distance image based on the target PTV contour image and the organ at risk OAR contour image obtained in step (2-a). The calculation formula is:

[0048]

[0049] where (i,j,k) is any voxel point in the CT image. represents the set of voxel points on the surface of the c-th PTV or OAR contour, and the subscript c represents the total number of PTVs and OARs in the CT image. s c represents any voxel point in, ‖(i,j,k),s c ‖ represents the distance from the voxel point (i,j,k) to the point s c and min represents taking the minimum value.

[0050] (2-c) Sample all the images obtained in step (2-a) and the distance image obtained in (2-b) to the same size and perform normalization processing.

[0051] Step 3, construct a dose prediction model, specifically including the following steps:

[0052] Figure 2 The following shows the structure diagram of the dose prediction model according to the embodiment of the present invention. This model mainly includes a noise addition process and a denoising process.

[0053] Noise addition process: At each time step t, add Gaussian noise ε t to the clinical dose image obtained in Step 2 t to obtain the corresponding intermediate feature map y

[0054] where t ∈ [0,T], and T represents the total number of time steps, taking T = 1000. Denoising process: Predict the noise added at each time step through a noise predictor t and subtract it from y to finally generate the predicted dose distribution map.

[0055] Figure 3The following is a structural diagram of the noise predictor according to an embodiment of the present invention, which is composed of an encoder, a feature fusion module MFF, a decoder, and a full-scale skip connection. The construction of the noise predictor specifically includes the following steps:

[0056] Ⅰ. Construct a three-branch encoder: Each branch includes five modules. The first module consists of a convolutional layer, and the second to fifth modules each consist of two residual blocks and a downsampling layer; the three branches are respectively used to extract features of x with a size of 1×128×128 CT , y with a size of 1×128×128 t , and x with a size of 13×128×128 other , where x CT is a CT image, and x other includes the PTV contour and its distance image, the OAR contour and its distance image, and the non-intensity-modulated dose image; x CT , y t , and x other are respectively processed by the first two modules of the three branches and then matrix-added to obtain a feature map with a size of 32×64×64 Similarly, x CT , y t , and x other are respectively processed by the first three modules of the three branches and then matrix-added to obtain a feature map with a size of 64×32×32 are respectively processed by the first four modules of the three branches and then matrix-added to obtain a feature map with a size of 128×16×16 x CT , y t , and x other respectively obtain feature maps f1, f3, and f2 with a size of 256×8×8 after being processed by the three branches;

[0057] Ⅱ. Figure 4 The following is a structural diagram of the feature fusion module MFF according to an embodiment of the present invention. Its construction steps are as follows: The feature fusion module consists of an interactive attention block, a gated attention block, and a Mamba block; after f1 and f2 are processed by the interactive attention block, the result f IA and f3 are input into the gated attention block together, and the obtained result f GA is further processed by the Mamba block to obtain the fusion result f out of f1, f2, and f3, with a size of 256×8×8;

[0058] Ⅲ. Construct a decoder: The decoder includes five modules. The first to fourth modules each consist of two residual blocks and an upsampling layer, and the fifth template is a convolutional layer; f obtained from step 3-Ⅱ outThe feature map is obtained after being processed by the first module of the decoder with a size of 128×16×16;

[0059] Ⅳ. Full-scale skip connection: Resample the obtained in step 3-Ⅰ to the same scale as the obtained in step 3-Ⅲ, process them separately using deformable convolutional blocks, then splice them. The spliced result is successively processed by a deformable convolutional block and the second decoding block to obtain a feature map with a size of 64×32×32 Similarly, resample the to the same scale as the and process them separately using deformable convolutional blocks, then splice them. The spliced result is successively processed by a deformable convolutional block and the third decoding block to obtain a feature map with a size of 32×64×64 Resample the to the same scale as the and process them separately using deformable convolutional blocks, then splice them. The spliced result is successively processed by a deformable convolutional block and the second decoding block to obtain a feature map with a size of 16×128×128 Then, it is processed by the fifth decoder to obtain a prediction result with a size of 1×128×128 The construction of the noise predictor is completed;

[0060] Figure 5 The following shows an example diagram of the full-scale skip connection in the embodiment of the present invention, which demonstrates the process of obtaining and to obtain

[0061] Step 4, training the prediction model

[0062] Train the prediction model, and use the AdamW optimizer to update the parameters of the model until the loss no longer decreases, obtaining the trained prediction model. The loss function used is:

[0063]

[0064] where M is the total number of voxel points.

[0065] Step 5, dose prediction

[0066] Use the dose prediction model trained in step 4 to perform dose prediction on the data to be tested, obtaining the dose prediction result.

[0067] Example 2

[0068] ​The method in Example 1 was used to conduct a dose prediction experiment on the collected liver tumor database. This experiment used the PyTorch 2.1.0 framework and was trained and tested in a hardware environment equipped with four NVIDIA GeForce 2080Ti GPUs.

[0069] In this implementation case, the Mean Absolute Error (MAE), D x , D mean , D max , V x , and the Homogeneity Index (HI), six indicators, were used to conduct experiments on five dose prediction networks, namely DCNN, DeepLabV3+, C3D, MD-Dose, DoseDiff, and the present invention, in the liver tumor database.

[0070] MAE is defined as the average of the sum of the absolute values of the differences between all predicted values and true values, and the calculation formula is:

[0071]

[0072] where m represents the total number of voxels, y c represents the true dose value of voxel i, and y p represents the predicted dose value of pixel i.

[0073] D x represents the minimum absorbed dose covering x% of the PTV or OAR volume. D mean represents the average dose absorbed by the PTV or OAR. D max represents the highest dose absorbed by any voxel in the PTV or OAR. V x represents the percentage of the volume within the PTV or OAR where the absorbed dose is greater than or equal to a specific dose threshold x in the total volume of this region.

[0074] HI is used to measure the dose distribution uniformity of the target region, and the calculation formula is as follows:

[0075]

[0076] The results obtained from the experiment are shown in Tables 1 and 2. Table 1 shows the comparison of the MAE index of different methods, and Table 2 shows the comparison of the indexes D x , D mean , D max , V x , and HI of different methods. It can be found that, compared with other methods, most of the evaluation indexes of the present invention on the PTV and OARs reach the optimal level.

[0077] Figure 6 The figure shows a comparison chart of the dose prediction results of the embodiments of the present invention and the dose prediction results of other methods. The trained prediction model is used to predict the test set, and the test example results are as follows Figure 6 shown, where Figure 6 (g) is the clinical dose distribution of the test example and also the label of the test example, Figure 6 (a)-(f) are the dose prediction results of DCNN, DeepLabV3+, C3D, MD-Dose, DoseDiff, and the present invention respectively; the first and third rows are the dose distribution results obtained from the test data; the second and fourth rows represent the difference between the predicted dose distribution result and the clinical dose distribution. The legend on the right shows the dose change range in Gy, and the brighter the area, the higher the dose value. The area with the highest brightness in the figure is the area where the tumor is located. It can be seen that the prediction error of the method proposed in the present invention is the smallest, and the beam shape in the prediction result is clear, and the high-frequency details are predicted more accurately.

[0078] The above-described embodiments are only the preferred embodiments of the present invention, and do not limit the scope of implementation of the present invention. Therefore, any changes made according to the structure and principle of the present invention should be covered within the protection scope of the present invention.

[0079] Table 1

[0080]

[0081] Table 2

[0082]

Claims

1. A tumor radiotherapy dose prediction method based on a diffusion model, characterized in that It includes the following steps: (1) Collect tumor radiotherapy plan data to form a database; (2) Data preprocessing, specifically including the following steps: (2-a) Obtain the patient's CT image, the contour image of the target volume PTV, the contour image of the organs at risk OAR, the non-intensity-modulated dose image, and the clinical dose image from the database described in step (1); (2-b) Calculate the distance image based on the target PTV contour image and the organ at risk OAR contour image obtained in step (2-a). The calculation formula is as follows: where (i, j, k) is any voxel point in the CT image, represents the set of surface voxel points of the c-th PTV or OAR contour, and the subscript c represents the total number of PTVs and OARs in the CT image, s c represents any voxel point in, ‖(i, j, k), s c ‖ represents the distance from the voxel point (i, j, k) to the point s c and min represents finding the minimum value; (2-c) Sample all the images obtained in step (2-a) and the distance image obtained in (2-b) to the same size and perform normalization processing; (3) Construct a dose prediction model, specifically including the following steps: (3-a) Noise addition process. At each time step t, Gaussian noise ε is added to the clinical dose image obtained in step (2) t to obtain the corresponding intermediate feature map y t , where t ∈ [0, T] and T is a positive integer representing the total number of time steps; (3-b) Denoising process, the noise added at each time step is predicted by the noise predictor And from y t Minus Finally, the predicted dose distribution map is generated; the construction of the noise predictor specifically includes the following steps: Ⅰ. Construct a three-branch encoder: Each branch consists of five modules. The first module is composed of a convolutional layer, and the second to fifth modules are each composed of two residual blocks and a downsampling layer; the three branches are respectively used to extract the features of x CT , y t , x other , where x CT is a CT image, and x other includes the PTV contour and its distance image, the OAR contour and its distance image, and the non-intensity-modulated dose image; x CT , y t , x other are respectively processed by the first two modules of the three branches and then added matrix-wise to obtain a feature map Similarly, x CT , y t , x other are respectively processed by the first three modules of the three branches and then added matrix-wise to obtain a feature map are respectively processed by the first four modules of the three branches and then added matrix-wise to obtain a feature map x CT , y t , x other are respectively processed by the three branches to obtain feature maps f1, f3, and f2; II. Construct the Feature Fusion Module MFF: The feature fusion module consists of an interactive attention block, a gated attention block, and a Mamba block; after f1 and f2 are processed by the interactive attention block, the result f IA is input into the gated attention block together with f3, and the resulting f GA is then processed by the Mamba block to obtain the fusion result f of f1, f2, and f3 out ; Ⅲ. Construct the decoder: The decoder contains five modules. The first to the fourth modules are each composed of two residual blocks and an upsampling layer, and the fifth module is a convolutional layer; the f obtained from step (3-b)-Ⅱ out After being processed by the first module of the decoder, a feature map is obtained Ⅳ. Full-scale skip connection: Resample the result obtained in step (3-b)-Ⅰ to the same scale as that obtained in step (3-b)-Ⅲ, process them separately using deformable convolutional blocks, and then splice them. The spliced result is successively processed by a deformable convolutional block and a second decoding block to obtain a feature map Similarly, resample to the same scale as that obtained in step (3-b)-Ⅲ, process them separately using deformable convolutional blocks, and then splice them. The spliced result is successively processed by a deformable convolutional block and a third decoding block to obtain a feature map Resample to the same scale as process them separately using deformable convolutional blocks, and then splice them. The spliced result is successively processed by a deformable convolutional block and a third decoding block to obtain a feature map Resample to the same scale as process them separately using deformable convolutional blocks, and then splice them. The spliced result is successively processed by a deformable convolutional block and a second decoding block to obtain a feature map Then, it is processed by a fifth decoder to obtain a prediction result Complete the construction of the noise predictor; (4) Train the prediction model: Train the prediction model. Use the AdamW optimizer to update the parameters of the model until the loss no longer decreases, and obtain the trained prediction model. The loss function used is: where M is the total number of voxel points; (5) Dose prediction: Use the trained dose prediction model in step (4) to perform dose prediction on the data to be tested, and obtain the dose prediction result.

2. The method for predicting tumor radiotherapy dose based on a diffusion model according to claim 1, wherein The residual block in step (3-b) includes two convolutional blocks and a residual connection. The convolutional block consists of a normalization layer, an activation layer, and a convolutional layer.

3. The tumor radiotherapy dose prediction method based on a diffusion model according to claim 1, characterized in that, The deformable convolutional block in step (3-b) is composed of a normalization layer, a deformable convolutional layer, a SiLU activation function layer, and a residual connection.

Images