CN113797450A - Radiotherapy plan intensity distribution optimization method and device and radiotherapy equipment - Google Patents
Radiotherapy plan intensity distribution optimization method and device and radiotherapy equipment Download PDFInfo
- Publication number
- CN113797450A CN113797450A CN202110586492.7A CN202110586492A CN113797450A CN 113797450 A CN113797450 A CN 113797450A CN 202110586492 A CN202110586492 A CN 202110586492A CN 113797450 A CN113797450 A CN 113797450A
- Authority
- CN
- China
- Prior art keywords
- dose
- intensity distribution
- algorithm
- radiotherapy plan
- neural network
- 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.)
- Pending
Links
- 238000000034 method Methods 0.000 title claims abstract description 73
- 238000005457 optimization Methods 0.000 title claims abstract description 61
- 238000001959 radiotherapy Methods 0.000 title claims abstract description 57
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 53
- 230000008021 deposition Effects 0.000 claims abstract description 47
- 239000011159 matrix material Substances 0.000 claims abstract description 46
- 238000013507 mapping Methods 0.000 claims abstract description 16
- 230000000694 effects Effects 0.000 claims abstract description 14
- 238000013528 artificial neural network Methods 0.000 claims abstract description 13
- 238000013527 convolutional neural network Methods 0.000 claims description 31
- 238000012549 training Methods 0.000 claims description 20
- 238000011176 pooling Methods 0.000 claims description 15
- 238000012360 testing method Methods 0.000 claims description 15
- 230000005855 radiation Effects 0.000 claims description 9
- 210000000056 organ Anatomy 0.000 claims description 8
- 230000004913 activation Effects 0.000 claims description 6
- 230000004927 fusion Effects 0.000 claims description 5
- 238000005070 sampling Methods 0.000 claims description 5
- 230000008569 process Effects 0.000 abstract description 21
- 230000011218 segmentation Effects 0.000 abstract description 5
- 230000006870 function Effects 0.000 description 16
- 238000004364 calculation method Methods 0.000 description 12
- 230000001537 neural effect Effects 0.000 description 8
- 238000010586 diagram Methods 0.000 description 7
- 238000005259 measurement Methods 0.000 description 5
- 238000005516 engineering process Methods 0.000 description 4
- 210000000920 organ at risk Anatomy 0.000 description 4
- 239000002245 particle Substances 0.000 description 4
- 238000012545 processing Methods 0.000 description 4
- 238000012937 correction Methods 0.000 description 3
- 230000004907 flux Effects 0.000 description 3
- 230000003213 activating effect Effects 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 2
- 238000011478 gradient descent method Methods 0.000 description 2
- 230000001788 irregular Effects 0.000 description 2
- 238000004088 simulation Methods 0.000 description 2
- 229910052721 tungsten Inorganic materials 0.000 description 2
- 238000010200 validation analysis Methods 0.000 description 2
- 238000012795 verification Methods 0.000 description 2
- 210000001015 abdomen Anatomy 0.000 description 1
- 230000003321 amplification Effects 0.000 description 1
- 230000001455 anti-clotting effect Effects 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 239000002131 composite material Substances 0.000 description 1
- 238000011109 contamination Methods 0.000 description 1
- 230000003203 everyday effect Effects 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 238000002721 intensity-modulated radiation therapy Methods 0.000 description 1
- 230000003993 interaction Effects 0.000 description 1
- 210000000867 larynx Anatomy 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000003199 nucleic acid amplification method Methods 0.000 description 1
- 238000003825 pressing Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 238000002948 stochastic simulation Methods 0.000 description 1
- 210000003437 trachea Anatomy 0.000 description 1
- 238000012546 transfer Methods 0.000 description 1
- 238000013519 translation Methods 0.000 description 1
- WFKWXMTUELFFGS-UHFFFAOYSA-N tungsten Chemical compound [W] WFKWXMTUELFFGS-UHFFFAOYSA-N 0.000 description 1
- 239000010937 tungsten Substances 0.000 description 1
Images
Classifications
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61N—ELECTROTHERAPY; MAGNETOTHERAPY; RADIATION THERAPY; ULTRASOUND THERAPY
- A61N5/00—Radiation therapy
- A61N5/10—X-ray therapy; Gamma-ray therapy; Particle-irradiation therapy
- A61N5/103—Treatment planning systems
- A61N5/1031—Treatment planning systems using a specific method of dose optimization
-
- A—HUMAN NECESSITIES
- A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
- A61N—ELECTROTHERAPY; MAGNETOTHERAPY; RADIATION THERAPY; ULTRASOUND THERAPY
- A61N5/00—Radiation therapy
- A61N5/10—X-ray therapy; Gamma-ray therapy; Particle-irradiation therapy
- A61N5/1048—Monitoring, verifying, controlling systems and methods
- A61N5/1071—Monitoring, verifying, controlling systems and methods for verifying the dose delivered by the treatment plan
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
- G16H20/00—ICT specially adapted for therapies or health-improving plans, e.g. for handling prescriptions, for steering therapy or for monitoring patient compliance
- G16H20/40—ICT specially adapted for therapies or health-improving plans, e.g. for handling prescriptions, for steering therapy or for monitoring patient compliance relating to mechanical, radiation or invasive therapies, e.g. surgery, laser therapy, dialysis or acupuncture
Landscapes
- Health & Medical Sciences (AREA)
- Engineering & Computer Science (AREA)
- Biomedical Technology (AREA)
- General Health & Medical Sciences (AREA)
- Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
- Public Health (AREA)
- Radiology & Medical Imaging (AREA)
- Pathology (AREA)
- Life Sciences & Earth Sciences (AREA)
- Animal Behavior & Ethology (AREA)
- Veterinary Medicine (AREA)
- Medical Informatics (AREA)
- Primary Health Care (AREA)
- Epidemiology (AREA)
- Urology & Nephrology (AREA)
- Surgery (AREA)
- Radiation-Therapy Devices (AREA)
Abstract
The invention discloses a radiotherapy plan intensity distribution optimization method, a radiotherapy plan intensity distribution optimization device and radiotherapy equipment, wherein the radiotherapy plan intensity distribution optimization method comprises the following steps: mapping the PB dose to the MC dose via a neural network; solving the gradient of the MC dose relative to the PB dose to obtain a deposition matrix based on an MC algorithm; and substituting the deposition matrix based on the MC algorithm into the FMO to realize the optimization of the intensity distribution of the radiotherapy plan. The invention applies the MC algorithm in FMO to ensure that the intensity distribution obtained in the optimization process is real and reliable, the obtained subfield exposure is more reasonable by using the intensity distribution to carry out multi-leaf grating subfield segmentation, the finally obtained planned dose distribution is ensured to be closer to the clinical requirement, and the planning effect is obviously improved.
Description
Technical Field
The invention belongs to the technical field of radiotherapy, and particularly relates to a radiotherapy plan intensity distribution optimization method and device and radiotherapy equipment.
Background
The inverse planning intensity modulated radiation therapy technology is a mainstream technology adopted by radiotherapy departments of medical institutions at present, and the optimization method comprises the following steps: two-step process, one-step process (direct sub-field optimization). The two-step process refers to: optimizing intensity distribution and dividing the grating sub-fields. The intensity distribution optimization is a process of solving the optimal ray intensity distribution which can meet the clinical requirements according to the combination of the patient image and the clinical requirements; raster subfield segmentation is the process of converting the optimal intensity distribution into a series of subfields that can be performed by a physical Multi-layer raster (MLC). In the intensity distribution Optimization (FMO), each continuous beam is divided into a number of beamlets (beams), the intensity values of each beamlet are solved using Optimization methods such as gradient descent, taking into account the dose deposition of each beamlet in the medium, respectively. Finally, for a specific point, the dose contributions of each beam to it are superimposed, resulting in the total dose for that point.
An important step in FMO is to calculate the dose deposition of the radiation intensity in human tissue, and the typical software uses Pencil Beam (PB). The dose value of the PB algorithm is the result of a two-dimensional convolution of the pencil-beam kernel, which describes a beam of energy-dispersive patterns incident in an infinitely homogeneous medium, with the intensity distribution. The method for solving the pencil beam kernel is as follows: MC simulation, direct from measurement data, deconvolution from measured broad beam radiation dose distribution, finite pencil beam method. For the case of inhomogeneous media, 3D density correction of the medium, scattering of the rays before reaching the phantom, and electronic contamination are considered. Although the PB algorithm can better handle the situation of irregular radiation fields, the calculation accuracy on the non-uniform tissues of the human body and the irregular contours of the body surface is still insufficient, and particularly for the cavity tissues such as the trachea and the larynx, the dose difference can reach 100 percent, as shown in figure 2. In order to better simulate the dose deposition of rays in inhomogeneous tissues, many studies propose methods such as primary ray fluence correction, kernel function correction and the like, but these methods are inherently data fitting methods and have far lower accuracy than the MC method.
The MC algorithm, also known as a computer stochastic simulation method, is based on the principle of random number generation, by simulating the process of interaction of a considerable number of particles in a medium, the distribution of flux, energy spectrum, dose in the whole medium is obtained. Since MC is a process that realistically simulates the deposition of energy from a large number of particles (up to a hundred million), the results of the simulation are recognized as being the most accurate.
In the FMO of commercial software, the PB algorithm is mostly used in consideration of the operation speed, and since the PB algorithm itself is not highly accurate, the optimized optimal intensity distribution is pseudo-optimal, and thus the dose distribution calculated by the intensity distribution through the final MC algorithm does not necessarily meet the expected clinical requirement, and the result is contrary to the original purpose of the clinical requirement setting. Therefore, the use of the MC algorithm in planning optimization is an important research direction.
Most researchers have generated MC deposition matrices for all beams in advance using tools such as EGS4, and called values directly in intensity distribution optimization. The efficiency of performing MC dose calculation directly on the CPU is low, and in order to increase the operation speed, the application of MC dose calculation on the GPU is gradually increased, and MC dose calculation software libraries based on the GPU, such as gDPM and GPUMCD, are developed. The MC optimization based on the CPU provides a solution, the application on the GPU can effectively improve the calculation speed, but the speed is still too slow for the intensity distribution optimization, and the requirement of commercial software cannot be met. In addition, the two methods require the deposition matrix to be stored in a computer, occupy a large amount of hard disks and memory space, and each patient needs to be calculated once, so that the time cost is high, and the method cannot meet the requirement of medical institutions with a large number of patients every day.
Disclosure of Invention
In order to solve the technical problems, the invention provides a radiotherapy plan intensity distribution optimization method, a radiotherapy plan intensity distribution optimization device and radiotherapy equipment.
In order to achieve the purpose, the technical scheme of the invention is as follows:
on one hand, the invention discloses a radiotherapy plan intensity distribution optimization method, which comprises the following steps:
s1: mapping the PB dose to the MC dose via a neural network;
s2: solving the gradient of the MC dose relative to the PB dose to obtain a deposition matrix based on an MC algorithm;
s3: and substituting the deposition matrix based on the MC algorithm into the FMO to realize the optimization of the intensity distribution of the radiotherapy plan.
On the basis of the technical scheme, the following improvements can be made:
preferably, S1 specifically includes the following steps:
s1.1: preparing data, namely selecting CT of N1 patients, adding N2 fields to each patient, respectively calculating PB dose and MC dose of each field to form N1N 2 groups of data, wherein each group of data comprises corresponding CT, PB dose and MC dose;
s1.2: dividing data, namely dividing N1N 2 groups of data into a training set test set in proportion;
s1.3: creating a convolutional neural network;
s1.4: training and establishing a good convolutional neural network by adopting training set data;
s1.5: and testing the test set data on the trained convolutional neural network to obtain a predicted dose, and calculating the gamma passing rate of the predicted dose and the MC dose.
Preferably, the convolutional neural network comprises: the device comprises an input layer, a fusion layer, a convolution layer, an activation layer, a pooling layer, an upper sampling layer and an output layer.
As a preferred embodiment of the present invention,
wherein A is a deposition matrix based on the MC algorithm;
DPBand DMCThe dosage based on the PB algorithm and the dosage based on the MC algorithm are respectively;
ajmfor each value of the deposition matrix based on the MC algorithm;
dDPB,k/dxmrepresents the derivative of the kth PB-dose with respect to the mth intensity value;
Ndoserepresenting the number of doses;
from the above equation, a deposition matrix based on the MC algorithm can be obtained by solving the gradient of the convolutional neural network of the PB dose to MC dose mapping.
As a preferred scheme, FMO specifically comprises the steps of:
s3.1: dividing the radiation field into pencil beams;
s3.2: introducing a deposition matrix based on an MC algorithm;
s3.3: setting an initial intensity value as an optimized initial value;
s3.4: solving a derivative of an objective function by using a deposition matrix based on an MC algorithm;
s3.5: solving the optimal step length under the current gradient;
s3.6: updating intensity values, dose and objective function;
s3.7: and repeating the steps S3.4-S3.6 until the iteration is stopped after the stop condition is met.
Preferably, the method further comprises the following steps:
s4: and calculating the dose in the human body and the DVH of each organ according to the radiotherapy plan intensity distribution optimization method, and evaluating the plan effect.
On the other hand, the invention also discloses a radiotherapy plan intensity distribution optimizing device, which comprises:
the neural network mapping module is used for mapping the PB dose to the MC dose through a neural network;
the deposition matrix obtaining module is used for solving the gradient of the MC dose relative to the PB dose to obtain a deposition matrix based on an MC algorithm;
and the intensity distribution optimization module is used for substituting the deposition matrix based on the MC algorithm into the FMO to realize the optimization of the intensity distribution of the radiotherapy plan.
In another aspect, the present invention also discloses a radiotherapy apparatus for performing any one of the radiotherapy plan intensity distribution optimization methods described above, or, comprising: the radiotherapy plan intensity distribution optimizing device.
The radiotherapy plan intensity distribution optimization method, the radiotherapy plan intensity distribution optimization device and the radiotherapy equipment have the following beneficial effects:
the invention applies the MC algorithm in FMO to ensure that the intensity distribution obtained in the optimization process is real and reliable, the obtained subfield exposure is more reasonable by using the intensity distribution to carry out multi-leaf grating subfield segmentation, the finally obtained planned dose distribution is ensured to be closer to the clinical requirement, and the planning effect is obviously improved.
Drawings
In order to more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings needed to be used in the embodiments will be briefly described below, it should be understood that the following drawings only illustrate some embodiments of the present invention and therefore should not be considered as limiting the scope, and for those skilled in the art, other related drawings can be obtained according to the drawings without inventive efforts.
Fig. 1 is a flowchart of a radiotherapy plan intensity distribution optimization method according to an embodiment of the present invention.
FIG. 2 is a graph comparing MC dose to PB dose;
FIG. 2(a) is a MC dose map;
FIG. 2(b) is a PB dose map.
FIG. 3 is a diagram of a layer of data in a training set according to an embodiment of the present invention;
FIG. 3(a) is a CT image;
FIG. 3(b) is a PB dose data plot;
FIG. 3(c) is a chart of MC dose data.
Fig. 4 is a flowchart of an intensity profile optimization method according to an embodiment of the present invention.
Fig. 5 is a flowchart of a method for mapping a PB dose to an MC dose according to an embodiment of the present invention.
Fig. 6 is a schematic structural diagram of a convolutional neural network according to an embodiment of the present invention.
FIG. 7 is a graph comparing neural dose gamma passage rate with PB dose and gamma passage rate predicted by the convolutional neural network provided in the embodiment of the present invention;
FIG. 7(a) is a MC dose map;
FIG. 7(b) is an MC dose distribution plot;
FIG. 7(c) is a GroudTruth diagram;
FIG. 7(d) is a PB dose map;
FIG. 7(e) is a PB dose profile;
FIG. 7(f) is a graph of PB dose gamma passage;
FIG. 7(g) is a neural dose map predicted by a convolutional neural network;
FIG. 7(h) is a neural dose distribution plot predicted by a convolutional neural network;
FIG. 7(i) is a graph of neural dose gamma passage predicted by the convolutional neural network.
Fig. 8 is a diagram of a convolution process of the 3 × 3 convolution kernel and the 4 × 4 matrix (outer 0 is complemented) according to an embodiment of the present invention.
Fig. 9 is a diagram of an inverse pooling upsampling process provided by an embodiment of the present invention.
Fig. 10 is a diagram of an interpolation upsampling process according to an embodiment of the present invention.
Fig. 11 is a diagram of a deconvolution upsampling process provided by an embodiment of the present invention.
Fig. 12 is a comparison graph of optimized DVH based on convolutional neural network and optimized DVH based on MC according to the embodiment of the present invention.
Detailed Description
Preferred embodiments of the present invention will be described in detail below with reference to the accompanying drawings.
The technical solutions in the embodiments of the present invention will be clearly and completely described below with reference to the drawings in the embodiments of the present invention, and it is obvious that the described embodiments are only a part of the embodiments of the present invention, and not all of the embodiments. All other embodiments, which can be derived by a person skilled in the art from the embodiments given herein without making any creative effort, shall fall within the protection scope of the present invention.
In order to achieve the object of the present invention, in some embodiments of a radiotherapy planning intensity distribution optimization method, a radiotherapy planning intensity distribution optimization device and a radiotherapy apparatus, as shown in fig. 1, the radiotherapy planning intensity distribution optimization method includes the following steps:
s1: mapping the PB dose to the MC dose via a neural network;
s2: solving the gradient of the MC dose relative to the PB dose to obtain a deposition matrix based on an MC algorithm;
s3: and substituting the deposition matrix based on the MC algorithm into the FMO to realize the optimization of the intensity distribution of the radiotherapy plan.
Notably, in S1, the PB and MC doses may be calculated, but are not limited to, by the planning system TPS.
The calculation method of the PB dose can calculate the dose deposition of the main rays and the scattered rays in the tissues according to a formula or a test result, and the calculation method of the MC dose is used for simulating the dose deposition of billions of particles in the tissues.
In order to further optimize the implementation effect of the present invention, in other embodiments, the remaining features are the same, except that S1 specifically includes the following steps:
s1.1: preparing data, namely selecting CT of N1 patients, adding N2 fields to each patient, respectively calculating PB dose and MC dose of each field to form N1N 2 groups of data, wherein each group of data comprises corresponding CT, PB dose and MC dose;
s1.2: dividing data, namely dividing N1N 2 groups of data into a training set test set in proportion;
s1.3: creating a convolutional neural network;
s1.4: training and establishing a good convolutional neural network by adopting training set data;
s1.5: and testing the test set data on the trained convolutional neural network to obtain a predicted dose, and calculating the gamma passing rate of the predicted dose and the MC dose.
Further, the convolutional neural network includes: input layer, fusion layer, convolution layer, activation layer, pooling layer, upsampling layer, output layer, as shown in fig. 6.
In a specific embodiment, the method for mapping the PB dose to the MC dose can be specifically divided into the following steps, as shown in fig. 5;
s1.1: preparing data:
selecting 50 patients, wherein the CT of the patients covers the head and neck, the chest and abdomen, the pelvic cavity and other parts; each patient is added with 72 fields (with the interval of 0.5 degrees), and the isocenter positions of the fields are distributed in different CT layers as far as possible, so that the fields can cover more CT ranges;
respectively calculating the PB dose and the MC dose of each radiation field, wherein the dose network is kept uniform and can adopt a 3mm grid;
overall, there were 3600 sets of CT, PB, MC dose data, as shown in fig. 3;
since the sampling pitch, size of the CT data does not coincide with the dose grid, the CT data needs to be interpolated to a size and spacing that coincides with the dose grid.
S1.2: data division: and (3) pressing 3600 groups of data as 7: 2: the scale of 1 is divided into a training set, a validation set, and a test set.
S1.3: creating a convolutional neural network: the convolutional neural network consists of an input layer, a fusion layer, a convolutional layer, an activation layer, a pooling layer, an upper sampling layer and an output layer.
And performing random cropping (crop) on the CT and PB dose images to improve the generalization capability of a neural network, wherein the size of the processed image is 128 × 128, the CT and PB dose images form a group of input data (128 × 128 × 2) through a fusion layer, then forming a 128 × 128 × 64 tensor through 64 convolution kernel processing, continuing convolution processing to form a 128 × 128 × 64 tensor, forming a 64 × 64 × 64 tensor through maximum pooling processing, and so on until the tensor of 128 × 128 × 64 is obtained after upsampling and convolution processing, and obtaining a 128 × 128 output image after convolution again.
The true output value of the convolutional neural network is the MC dose, and the loss function is defined as the Mean Square Error (MSE) between the MC dose and the predicted dose, which is the Mean of the sum of the squared distances between the target variable and the predicted value:
wherein: n is 128 × 128.
S1.4: training a convolutional neural network;
the network training process is to continuously iterate and solve the convolution kernel, and is a feature extraction process;
the specific steps of training are as follows:
i. randomly selecting data (batch size is 512), and performing online data amplification on an input image, wherein the method comprises translation and 2D rotation, and aims to increase the data volume and improve the generalization capability of a convolutional neural network;
using random gradient descent method to back-propagate to obtain all the weight and offset gradients of the network, using a learning rate of 10-5Updating the weight and the offset, obtaining a new predicted value on the basis, comparing the predicted value with the MC dose to obtain the mean square error, and repeating the process for 200 times to finish an epoch;
repeating steps i), ii) until the mean square error is less than a threshold or a maximum epoch is reached (e.g.: 5000) Stopping iteration, and finishing the training of the convolutional neural network;
the above calculation process is performed on nvidirtx 2080 Ti.
In this embodiment, the data is divided into a training set, a validation set, and a test set. In the training process, generally, during training, a verification set is run once after several epochs are finished to see the effect. Therefore, the problems of the model or the parameters can be found in time, the training is stopped in time, and the model is readjusted or adjusted. In addition, the generalization capability of the model can be verified, and if the effect on the verification set is much worse than on the training set, whether the model is over-fitted or not should be considered.
S1.5: testing the effect;
testing on a trained convolutional neural network, inputting CT and PB dose data of a test set to obtain predicted dose, storing all the predicted dose data in a local area, and calculating gamma passing rates (3mm, 3 percent and 10 percent) of the predicted dose and the MC dose by using dose verification software (such as DVS and LinaTech), wherein the passing rate is more than 95 percent, namely the pass rate is qualified.
gamma passage is a quantitative measure of the Difference between "measure" and "calculate" and is a combination of Distance To Agent (DTA) and Dose Difference (Dose Difference, DD):
r(rm,rc)＝|rc-rm|；
δ(rm,rc)＝Dc(rc)-Dm(rm)；
wherein:
r(rm,rc) Representing for the measurement point rmAnd calculating the point rcThe distance between them;
δ(rm,rc) Represents the dose difference between the two;
ΔdMrepresentative distance threshold (which may be but is not limited to taking 3 mm);
ΔDMa threshold representing a dose (which may be but is not limited to taking 3% of the maximum dose);
when gamma (r)m) When the value is less than or equal to 1, the point is passed, otherwise, the point is not passed.
If only the distance factor (delta (r) is consideredm,rc) 0), then the basis for the passing of the measurement point is: within a distance threshold (e.g., 3mm), there is a calculation point where the dose value equals the measurement point;
if only the dosage factor (r) is consideredm,rc) 0), then the basis for the passing of the measurement point is: at the same location, the difference between the calculated dose and the measured dose is no greater than a dose threshold (e.g., 3% of maximum).
As shown in fig. 7, which shows the gamma-passing rate of PB dose and neural dose (based on neural network dose, i.e. predicted dose of convolutional neural network) in a certain layer, the MC dose is true, it can be seen that the neural dose is closer to the MC dose, and the statistical result shows that the gamma-passing rate of neural dose reaches 99%, while the gamma-passing rate of PB dose is only 62%.
In order to further optimize the implementation effect of the present invention, in other embodiments, the rest of the feature technologies are the same, except that a deposition matrix based on the MC algorithm is obtained, and the principle is as follows:
wherein A is a deposition matrix based on the MC algorithm;
DPBand DMCThe dosage based on the PB algorithm and the dosage based on the MC algorithm are respectively;
ajmfor each value of the deposition matrix based on the MC algorithm;
dDPB,k/dxmrepresents the derivative of the kth PB-dose with respect to the mth intensity value;
Ndoserepresenting the number of doses;
from the above equation, a deposition matrix based on the MC algorithm can be obtained by solving the gradient of the convolutional neural network of the PB dose to MC dose mapping.
Notably, the above formula transforms the problem "calculate deposit matrix using MC algorithm" into "solve derivative of MC dose to PB dose", where: derivative dD of PB dose to intensity value xPB,k/dxmThe method is characterized in that the deposition matrix obtained by the original algorithm is a matrix calculated by using a method of superposing main rays and scattered rays in TPS software.
The essence of solving the derivative of the MC dose to the PB dose is the gradient solution of a convolutional neural network mapping the PB dose to the MC dose, and the solution method can refer to a back propagation algorithm.
1) The gradient of each type of network is solved, and the gradient solving of the convolutional layer (with active layer), the pooling layer and the upsampling layer is described below.
i. Gradient of convolutional layer with active layer: as shown in FIG. 8, let input be X, weight and bias for convolutional layers be W and b, respectively, output be Y, and activation function be
Then there are:
conv2(W, X) represents a discrete two-dimensional convolution of W and X, which can also be written as:
thus is provided with
Wherein:
k reflects the convolution kernel size (2 × k + 1) ═ kernel size), and when a 3 × 3 convolution kernel is used, k is 1;
activating a function
If it is a ReLU function, its derivative is:
activating a function
If it is a Sigmoid function, its derivative is:
the derivative of the convolutional layer (with activation function) can be obtained according to the above formula.
Gradient of pooling layer: taking Max Pooling (Max Pooling), a 2 × 2 convolution kernel as an example, the output of the layer is equal to the maximum of the corresponding input, i.e., y (i, j) ═ Max { x (p, q), p ═ i-1, i, q ═ j-1, j }, then the gradient can be expressed as:
gradient of the upsampled layer: there are three methods for upsampling (Up sampling): interpolation, deconvolution, and anticlotting;
a) as shown in fig. 10, the interpolation method: enlarging an input image to a certain size, calculating a pixel value of each point by using an interpolation method, taking bilinear interpolation as an example, and expressing the pixel value of a certain point of the upsampled image as the following formula, wherein the position of the upsampled image I, J corresponds to I, J of the original image, and outputting the gradient of the input image as interpolation weight w;
b) as shown in fig. 11, deconvolution: convolution can be expressed as the product of two matrixes, deconvolution is the multiplication of the transpose of a convolution kernel matrix and an original matrix, and is a conversion method from small size to large size, and is also convolution in nature, and the gradient calculation method can refer to the convolution method.
c) As shown in fig. 9, anti-pooling: taking the maximum pooling as an example, recording the coordinates of max-pooling in the corresponding kernel, amplifying an element according to the kernel in the anti-pooling process, filling the element according to the previous coordinates, and filling 0 in other positions.
Its gradient calculation is similar to the pooling layer method: the gradient at the maximum passes directly, and the other positions have a gradient of 0.
2) Calculating the gradient reversely in turn: considering a network comprising an input layer A, a hidden layer B, an output layer C, the information transfer direction is A → B → C, C for each point CijRelative to each point A in AstThe gradient of (c) can be passed through B, written as:
considering the whole convolutional neural network, starting from the last output layer, calculating the gradients relative to the previous layer, the previous two layers and the previous three layers in sequence, and continuously recursing forwards until the gradient relative to the input layer is obtained through calculation.
In order to further optimize the implementation effect of the present invention, in other embodiments, the remaining features are the same, except that the FMO specifically includes the following steps:
s3.1: dividing the radiation field into pencil beams;
s3.2: introducing a deposition matrix based on an MC algorithm;
s3.3: setting an initial intensity value as an optimized initial value;
s3.4: solving a derivative of an objective function by using a deposition matrix based on an MC algorithm;
s3.5: solving the optimal step length under the current gradient;
s3.6: updating intensity values, dose and objective function;
s3.7: and repeating the steps S3.4-S3.6 until the iteration is stopped after the stop condition is met.
In a specific embodiment, the intensity profile optimization is the first step of a two-step process of planning optimization, and the step of performing this optimization using a gradient descent method is shown in FIG. 4, and includes:
s3.1: field segmentation into pencil beams (beamlets): dividing the set range of the field tungsten gate into N multiplied by M rectangular networks, wherein N represents the width direction of the grating, M represents the movement direction of the grating, the grid in the N direction is divided according to the width of the grating under normal conditions, so that the next grating sub-field division is conveniently executed, and the M direction is divided according to equal intervals;
s3.2: introducing a deposition matrix A based on an MC algorithm;
s3.3: setting an initial strength value: as an initial value for optimization;
s3.4: solving the derivative of the objective function: the optimization objective function defined here is:
wherein:
NOARthe number of organs at risk;
Nithe number of voxel points of the ith crisis organ;
wjis the weight of the jth voxel point;
δjthe calculated dose is a penalty factor, and is 0 when the calculated dose meets the preset constraint, otherwise, the calculated dose is 1;
pjis the desired dose value;
in the same way, the method for preparing the composite material,
is a quadratic error of the target area, where NPTVThe number of voxel points of the target area;
ajmfor the dose contribution of the mth pencil beam to the jth voxel, the resulting matrix is the deposition matrix a based on the MC algorithm.
Thus, the derivative of the objective function is:
s3.5: solving the optimal step length: solving the optimal step length under the current gradient;
s3.6: updating intensity values, dose and objective function;
s3.7: judging whether a stop condition is met: when the argument (intensity value x) or the objective function Jobj(x) Is smaller than the threshold epsilon, the iteration is stopped, otherwise steps S3.4-S3.4 are repeated.
In order to further optimize the implementation effect of the present invention, in other embodiments, the rest of feature technologies are the same, except that the radiotherapy plan intensity distribution optimization method further comprises the following steps:
s4: and calculating the dose in the human body and the DVH of each organ according to the radiotherapy plan intensity distribution optimization method, and evaluating the plan effect.
In a specific embodiment, the gDPM is used directly to compute the deposition matrix (total particle number about 2 × 10) for all pencil beams (beams)9) And substituting the optimized intensity distribution into an optimization process to obtain the optimal intensity distribution, and further obtaining the dose distribution as a benchmark (baseline) for intensity distribution optimization.
The method disclosed by the invention is used for obtaining a deposition matrix based on the MC algorithm, adding the deposition matrix into an optimization flow to obtain optimized intensity distribution, calculating the dose in a human body by using the MC algorithm, and counting the DVH of organs.
The DVH of the two samples are compared, and the planning effect of the invention is verified.
The results of both are shown in fig. 12, in fig. 12: neural dose represents the result based on the neural network (namely the result of the invention), and MC dose represents the result based on random calculation, and the results are almost completely consistent, thereby verifying the feasibility of the invention.
On the other hand, the embodiment of the invention also discloses a radiotherapy plan intensity distribution optimizing device, which comprises:
the neural network mapping module is used for mapping the PB dose to the MC dose through a neural network;
the deposition matrix obtaining module is used for solving the gradient of the MC dose relative to the PB dose to obtain a deposition matrix based on an MC algorithm;
and the intensity distribution optimization module is used for substituting the deposition matrix based on the MC algorithm into the FMO to realize the optimization of the intensity distribution of the radiotherapy plan.
On the other hand, the embodiment of the present invention further discloses a radiotherapy apparatus, configured to perform the radiotherapy plan intensity distribution optimization method disclosed in any of the above embodiments, or, the radiotherapy apparatus includes: the radiotherapy planning intensity distribution optimizing device disclosed in the above embodiment.
The radiotherapy plan intensity distribution optimization method, the radiotherapy plan intensity distribution optimization device and the radiotherapy equipment have the following beneficial effects:
the invention applies the MC algorithm in FMO to ensure that the intensity distribution obtained in the optimization process is real and reliable, the obtained subfield exposure is more reasonable by using the intensity distribution to carry out multi-leaf grating subfield segmentation, the finally obtained planned dose distribution is ensured to be closer to the clinical requirement, and the planning effect is obviously improved.
The above embodiments are merely illustrative of the technical concept and features of the present invention, and the purpose thereof is to enable those skilled in the art to understand the content of the present invention and implement the present invention, and not to limit the scope of the present invention, and all equivalent changes or modifications made according to the spirit of the present invention should be covered in the scope of the present invention.
Claims (8)
1. A radiotherapy plan intensity distribution optimization method is characterized by comprising the following steps:
s1: mapping the PB dose to the MC dose via a neural network;
s2: solving the gradient of the MC dose relative to the PB dose to obtain a deposition matrix based on an MC algorithm;
s3: and substituting the deposition matrix based on the MC algorithm into the FMO to realize the optimization of the intensity distribution of the radiotherapy plan.
2. The radiotherapy plan intensity distribution optimization method of claim 1, wherein S1 specifically comprises the following steps:
s1.1: preparing data, namely selecting CT of N1 patients, adding N2 fields to each patient, respectively calculating PB dose and MC dose of each field to form N1N 2 groups of data, wherein each group of data comprises corresponding CT, PB dose and MC dose;
s1.2: data division, namely dividing N1N 2 groups of data into a training set and a test set in proportion;
s1.3: creating a convolutional neural network;
s1.4: training and establishing a good convolutional neural network by adopting training set data;
s1.5: and testing the test set data on the trained convolutional neural network to obtain a predicted dose, and calculating the gamma passing rate of the predicted dose and the MC dose.
3. The radiotherapy plan intensity distribution optimization method of claim 2, wherein the convolutional neural network comprises: the device comprises an input layer, a fusion layer, a convolution layer, an activation layer, a pooling layer, an upper sampling layer and an output layer.
4. The radiotherapy plan intensity distribution optimization method of claim 3,
wherein A is a deposition matrix based on the MC algorithm;
DPBand DMCThe dosage based on the PB algorithm and the dosage based on the MC algorithm are respectively;
ajmfor each value of the deposition matrix based on the MC algorithm;
dDPB,k/dxmrepresents the derivative of the kth PB-dose with respect to the mth intensity value;
Ndoserepresenting the number of doses;
from the above equation, a deposition matrix based on the MC algorithm can be obtained by solving the gradient of the convolutional neural network of the PB dose to MC dose mapping.
5. The radiotherapy plan intensity distribution optimization method of claim 4,
the FMO specifically comprises the following steps:
s3.1: dividing the radiation field into pencil beams;
s3.2: introducing a deposition matrix based on an MC algorithm;
s3.3: setting an initial intensity value as an optimized initial value;
s3.4: solving a derivative of an objective function by using a deposition matrix based on an MC algorithm;
s3.5: solving the optimal step length under the current gradient;
s3.6: updating intensity values, dose and objective function;
s3.7: and repeating the steps S3.4-S3.6 until the iteration is stopped after the stop condition is met.
6. The radiotherapy plan intensity distribution optimization method of any one of claims 1 to 5, further comprising the steps of:
s4: and calculating the dose in the human body and the DVH of each organ according to the radiotherapy plan intensity distribution optimization method, and evaluating the plan effect.
7. Radiotherapy plan intensity distribution optimizing device characterized by, includes:
the neural network mapping module is used for mapping the PB dose to the MC dose through a neural network;
the deposition matrix obtaining module is used for solving the gradient of the MC dose relative to the PB dose to obtain a deposition matrix based on an MC algorithm;
and the intensity distribution optimization module is used for substituting the deposition matrix based on the MC algorithm into the FMO to realize the optimization of the intensity distribution of the radiotherapy plan.
8. Radiotherapy apparatus for performing a radiotherapy plan intensity distribution optimization method according to any one of claims 1 to 6, or comprising: the radiotherapy plan intensity distribution optimizing device of claim 7.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110586492.7A CN113797450A (en) | 2021-05-27 | 2021-05-27 | Radiotherapy plan intensity distribution optimization method and device and radiotherapy equipment |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202110586492.7A CN113797450A (en) | 2021-05-27 | 2021-05-27 | Radiotherapy plan intensity distribution optimization method and device and radiotherapy equipment |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN113797450A true CN113797450A (en) | 2021-12-17 |
Family
ID=78942413
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202110586492.7A Pending CN113797450A (en) | 2021-05-27 | 2021-05-27 | Radiotherapy plan intensity distribution optimization method and device and radiotherapy equipment |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN113797450A (en) |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101954148A (en) * | 2010-09-15 | 2011-01-26 | 四川大学 | Method for accelerating dosage calculation in radiotherapy based on GPU (Graphics Processing Unit) |
| CN110327554A (en) * | 2019-07-08 | 2019-10-15 | 南方医科大学 | Intensity modulated radiation therapy plan optimization method and application based on predicted dose distribution guidance |
| US20200118306A1 (en) * | 2018-10-12 | 2020-04-16 | Korea Advanced Institute Of Science And Technology | Method for processing unmatched low-dose x-ray computed tomography image using neural network and apparatus therefor |
| CN112233200A (en) * | 2020-11-05 | 2021-01-15 | 福建自贸试验区厦门片区Manteia数据科技有限公司 | Dose determination method and device |
-
2021
- 2021-05-27 CN CN202110586492.7A patent/CN113797450A/en active Pending
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101954148A (en) * | 2010-09-15 | 2011-01-26 | 四川大学 | Method for accelerating dosage calculation in radiotherapy based on GPU (Graphics Processing Unit) |
| US20200118306A1 (en) * | 2018-10-12 | 2020-04-16 | Korea Advanced Institute Of Science And Technology | Method for processing unmatched low-dose x-ray computed tomography image using neural network and apparatus therefor |
| CN110327554A (en) * | 2019-07-08 | 2019-10-15 | 南方医科大学 | Intensity modulated radiation therapy plan optimization method and application based on predicted dose distribution guidance |
| CN112233200A (en) * | 2020-11-05 | 2021-01-15 | 福建自贸试验区厦门片区Manteia数据科技有限公司 | Dose determination method and device |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| US11413474B2 (en) | System and method for modelling of dose calculation in radiotherapy treatment planning | |
| Mao et al. | RapidBrachyDL: rapid radiation dose calculations in brachytherapy via deep learning | |
| Knöös et al. | Comparison of dose calculation algorithms for treatment planning in external photon beam therapy for clinical situations | |
| US8401148B2 (en) | Non-voxel-based broad-beam (NVBB) algorithm for intensity modulated radiation therapy dose calculation and plan optimization | |
| US7450687B2 (en) | Method for verification of intensity modulated radiation therapy | |
| Yepes et al. | A GPU implementation of a track-repeating algorithm for proton radiotherapy dose calculations | |
| US20060259282A1 (en) | Deterministic computation of radiation transport for radiotherapy dose calculations and scatter correction for image reconstruction | |
| US20080091388A1 (en) | Method for calculation radiation doses from acquired image data | |
| US20070201614A1 (en) | Method and system for optimizing dose delivery of radiation | |
| JP2013526883A (en) | A method for calculating the load on which ionizing radiation is deposited. | |
| CN110570923B (en) | Mixed Monte Carlo radiotherapy reverse optimization method, equipment and storage medium | |
| US9495513B2 (en) | GPU-based fast dose calculator for cancer therapy | |
| Da Silva et al. | Sub-second pencil beam dose calculation on GPU for adaptive proton therapy | |
| CN107551411B (en) | Proton heavy ion intensity modulated radiotherapy robust optimization method aiming at range uncertainty | |
| Pawlicki et al. | Monte Carlo simulation for MLC-based intensity-modulated radiotherapy | |
| Penfold | Image reconstruction and Monte Carlo simulations in the development of proton computed tomography for applications in proton radiation therapy | |
| US9486644B2 (en) | Method and system for dose determination of radiation therapy | |
| Maleike et al. | Simulation and visualization of dose uncertainties due to interfractional organ motion | |
| CN113178242B (en) | Automatic plan optimization system based on coupled generation countermeasure network | |
| Vazquez et al. | A deep learning-based approach for statistical robustness evaluation in proton therapy treatment planning: a feasibility study | |
| Vasseur et al. | Dose calculations using artificial neural networks: a feasibility study for photon beams | |
| CN113797450A (en) | Radiotherapy plan intensity distribution optimization method and device and radiotherapy equipment | |
| Yeo et al. | Dose reconstruction for intensity-modulated radiation therapy using a non-iterative method and portal dose image | |
| CN115270588A (en) | Fluence spectral distribution calculation modeling and fluence spectral distribution calculation method and device | |
| Tseng | Dose Calculation Strategies Toward Efficient and Accurate Treatment Planning for Modern Radiotherapy |
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 |