Method for field optimization based on deep neural networks and evolutionary algorithms

Abstract

The application discloses a kind of based on deep neural network and evolution algorithm's irradiation field optimization method, comprising: step one, data is preprocessed;Step two, model training, obtain dose prediction map and determine constraint equation accordingly;Step three, determine objective function and complete population initialization;Step four, optimize irradiation field set by evolution algorithm, after the cross variation and environment selection operation of population maximum iteration number times, the irradiation field scheme corresponding to the three variables of the number of sub-field, irradiation intensity and the leaf shape of multileaf collimator that the Pareto optimal solution of P are obtained.The application can realize the optimization of irradiation field by neural network and evolution algorithm, so as to improve the degree of automation and efficiency of intensity modulation radiotherapy.

Description

Technical Field

[0001] This invention relates to the field of illumination field optimization, specifically to an illumination field optimization method based on deep neural networks and evolutionary algorithms. Background Technology

[0002] Intensity-modulated radiation (IMRT) systems are mature and widely used core software for controlling electron linear accelerators and multi-leaf collimators. By changing the positions of multiple pairs of independently moving parallel blades, the multi-leaf collimator can form an irregular radiation field, effectively irradiating the target area while avoiding damage to normal areas. To determine the optimal blade positions, physicists need to establish a suitable optimization model based on the data requirements and solve it using optimization methods. Irradiation field optimization requires simultaneously optimizing the positions of thousands of multi-leaf collimator blades to meet the dose requirements of dozens of regions, while considering the dynamic errors caused by the actual conditions of the target. Therefore, irradiation field optimization is a complex optimization problem with challenges such as large-scale decision variables, numerous targets, discrete coding, and noise, while also demanding high efficiency from the optimization algorithm. However, currently commercially available optimization systems are mainly based on mathematical programming optimization methods, which can only obtain a local optimum each time. This forces physicists to repeatedly adjust the weights between targets and continuously change the parameters within the targets based on experience to obtain diverse solutions and ultimately determine the optimal irradiation scheme. In short, how to improve the optimization effect and reduce the difficulty of irradiation has always been a key focus and challenge in this field.

[0003] Radiation field optimization aims to search for the optimal leaf positions and irradiation intensities of a multi-leaf collimator by meshing the rectangular irradiation area and pre-calculating the dose deposition on each voxel for each mesh. This partially masks the mesh to form an irregular target region, effectively treating the target region while avoiding damage to normal areas. Fluorescence map optimization aims to search for the optimal irradiation intensities for each mesh based on the data to be optimized; a zero irradiation intensity indicates that the mesh is masked. This step can typically be implemented using analytical methods or metaheuristic algorithms, but is limited by the physical constraints of the electron linear accelerator and the multi-leaf collimator. Although existing planning methods can achieve relatively fast convergence speeds with the support of mathematical programming, they can only obtain a single solution in a single optimization process, resulting in low automation and inefficiency in existing planning systems. Summary of the Invention

[0004] To overcome the shortcomings of existing methods, this invention proposes an irradiation field optimization method based on deep neural networks and evolutionary algorithms, aiming to obtain the optimal irradiation field set scheme with little or no human intervention, thereby improving the effect and efficiency of intensity-modulated radiation.

[0005] To achieve the above-mentioned objectives, the present invention adopts the following technical solution:

[0006] The invention discloses a irradiation field optimization method based on deep neural network and evolutionary algorithm, characterized in that the method is performed according to the following steps:

[0007] Step 1: preprocess data;

[0008] Obtain one CT image and H-1 organ contour images and perform normalization processing on them respectively to obtain a preprocessed CT image and an organ contour image set; wherein H is the total number of images; obtain an optimal dose map GT;

[0009] Step 2: build a dose prediction network based on 3D-Unet network, which is used to predict the dose of each target region, and establish a constraint equation according to the predicted dose;

[0010] Step 2.1: the dose prediction network consists of an upsampling path, a downsampling path, and a convolution layer with a dimension of b×b×b, and skip connections are used to introduce underlying features in the downsampling path into the upsampling path. The dimension of an input image is defined as Y×Y×Y, where Y is the size of the input image; b represents the size of the convolution kernel;

[0011] Step 2.1.1: the downsampling path is composed of T identical downsampling modules, and each downsampling module is sequentially composed of a repeated downsampling layer, a convolution layer, an instance normalization layer and a rectified linear unit. After an input feature map is processed by one downsampling module, the resolution is reduced by half and the number of channels is doubled;

[0012] The CT image and the organ contour image set are input into the dose prediction network, and after encoding processing by the downsampling path, 2 T ·H (Y·2 -T )×(Y·2 -T )×(Y·2 -T ) feature maps are obtained;

[0013] Step 2.1.2: the upsampling path is composed of T identical upsampling modules, and each upsampling module is sequentially composed of a repeated upsampling layer, a skip cascade, a convolution layer, an instance normalization layer and a rectified linear unit;

[0014] Step 2.1.3: 2 T ·H (Y·2 -T )×(Y·2 -T )×(Y·2 -T ) feature maps are input into the upsampling path for decoding processing, and then H×Y×Y×Y feature maps are obtained;

[0015] Step 2.1.4: The H×Y×Y×Y feature map is processed by a convolutional layer of dimension b×b×b to generate a predicted dose map. The dose values ​​of all voxels are extracted from the predicted dose map to obtain the dose value prediction set d = {d1, d2…d...}. n …d N}; where d n The dose value predicted for the nth voxel; N is the number of voxels used in the target region of the CT image;

[0016] Step 2.2: Construct the loss L(N,GT) between the predicted dose map d and the optimal dose map GT using equation (1), and use it as the supervision loss for training the dose prediction network. Update the weight parameters of the network through the backpropagation algorithm, and stop network training when the maximum number of iterations is reached, thereby obtaining the optimal dose prediction network.

[0017]

[0018] In equation (2), GT n This represents the optimal dose value for the nth voxel in the optimal dose map GT;

[0019] Step 2.3: Input the CT image to be predicted into the optimal dose prediction network for processing, and obtain the dose value d of the voxel in the optimal dose prediction map. pre ;

[0020] Step 2.4: Construct constraint equations using equation (2);

[0021]

[0022] In equation (2), ε is the optimization threshold; and we have:

[0023]

[0024] y g =a j h j，g (4)

[0025] In equations (3) and (4), B represents the number of irradiations, and D represents the number of irradiations. g，n Let D be the dose precipitation coefficient at the g-th sampling point of the n-th voxel, and D g，n ∈D, where D is the dose-precipitation matrix, C j Let y be the set of sampling points for the j-th irradiation. g Let a be the irradiation intensity at the g-th sampling point in the irradiation field scheme to be optimized. j h is the irradiation intensity of the j-th irradiation. j，g h represents the leaf shading at the g-th sampling point during the i-th irradiation. j，g=0 indicates that the g-th sampling point is blocked by the blades of the multi-leaf collimator during the j-th irradiation;

[0026] Step 3: Determine the objective function equation and initialize N. pop A scheme for assembling illumination fields;

[0027] Step 3.1: Construct the overall objective function F(x):

[0028] Step 3.1.1: Construct the minimum dose equation MinDose(x, V, d) for each voxel using equation (5). min ):

[0029]

[0030] In equation (5), x represents the irradiation field scheme using a three-layer coding method, where the three-layer coding refers to coding the irradiation intensity, the number of subfields, and the leaf shape respectively; d min The minimum dose for setting dose constraints;

[0031] Step 3.1.2: Construct the maximum dose constraint equation MaxDose(x, V, d) for each voxel using equation (6). max ):

[0032]

[0033] In equation (6), d max The maximum dose for the set dose constraint;

[0034] Step 3.1.3: Construct the dose-volume constraint equation DVC(x, V, d) to restrict all voxels using equation (7). dvc ,r):

[0035]

[0036] In equation (7), d dvc For the maximum dose under the set dose-volume constraint, d′ r This represents the predicted dose value of the r-th voxel after sorting d in descending order;

[0037] Step 3.1.4: Construct the first objective function f1(x) using equation (8):

[0038]

[0039] Step 3.1.5: Construct the second objective function f2(x) using equation (9):

[0040]

[0041] In equation (9), PTV is the planned volume of the target cells in the CT image;

[0042] Step 3.1.6: Construct the third objective function f3(x) using equation (10):

[0043]

[0044] In equation (10), OAR is the normal region in the CT image;

[0045] Step 3.1.7: Construct the third objective function f4(x) using equation (11):

[0046]

[0047] Step 3.1.8: Construct the overall objective function F(x) using equation (12):

[0048] Minimize F(x)=(f1(x),…,f4(x)) (12)

[0049] Step 3.2: Calculate the fitness of the subfield quantity variable on the objective function in each dimension: FN = {FN1, FN2, ..., FN} s …FN s The set of fitness values ​​of each dimension of the blade shape variable on the objective function is FH = {FH1, FH2, ..., FH}. s …FH s}, the set of fitness values ​​FH for each row of the leaf shape variable in the s-th dimension. s ={FH s，1 FH s，2 …FH s，k …FH s，U}, the set of fitness values ​​FH for each sampling point in the k-th row of the s-th dimension of the two-dimensional matrix of blade shape variables. s，k ={FH s，k，1 FH s，k，2 …FH s，k，w …FH s，U，W}, where FN s FHs represents the fitness of the s-th dimension subfield quantity variable on the objective function, and FH represents the fitness matrix of the s-th dimension leaf shape variable on the objective function. sk FH represents the set of fitness scores of the sampling points in the k-th row of the s-th dimension of a two-dimensional matrix representing blade shape variables on the objective function. skw The fitness of the sampling point in the k-th row and w-th column of the s-th dimension of the two-dimensional matrix representing the leaf shape variable on the objective function is given by S, where S is the total number of illumination rays from different angles as the upper limit of the dimension, and the leaf shape is a matrix of size U×W, where U is the number of rows of the matrix and W is the number of columns of the matrix.

[0050] Step 3.3: Define the maximum number of iterations as G. max The current iteration number is G, and G is initialized to 1;

[0051] Define the population of generation G Let N represent the i-th irradiation field scheme of the G-th generation, where each irradiation field scheme is composed of the variables of subfield quantity, leaf shape, and irradiation intensity; pop The number of schemes in the population;

[0052] Step 3.3.1: Define and randomly initialize the population P of generation G. G The irradiation intensity variable, where the P of the Gth generation population G The i-th illumination field scheme Irradiation intensity variable in, P represents the population of generation G. G The i-th illumination field scheme The s-th real number encoding in the middle, and a max This represents the maximum irradiation intensity.

[0053] Step 3.3.2: Define and randomly initialize the population P of generation G. G The variable of the number of subfields, where the population P in the Gth generation. G The i-th illumination field scheme Binary encoding of the subfield quantity variable Head. P represents the population of generation G. G The i-th illumination field scheme The s-th binary code of the subfield quantity variable, and 0 represents the population P of generation G. G The s-th subfield of the i-th illumination field scheme does not exist, and 1 indicates that the population P of the G-th generation does not exist. G The s-th subfield of the i-th illumination field scheme exists;

[0054] Step 3.3.3: Define and randomly generate the population P in the Gth generation. G Leaf shape variables, where the Gth generation population P G The i-th illumination field scheme Leaf shape variable in, For the i-th illumination field scheme The s-th dimension of the two-dimensional matrix of the blade shape variable. Represents the i-th illumination field scheme The k-th row of the s-th dimension of the two-dimensional matrix of the leaf shape variable. Represents the i-th illumination field scheme Is the sampling point located in the k-th row and w-th column of the s-th dimension of the two-dimensional matrix of the blade shape variable obscured by the blades of the multi-leaf collimator? A value of 0 indicates the i-th illumination field scheme. The sampling point in the kth row and wth column of the s-th dimension of the two-dimensional matrix of the blade shape variable is blocked by the blade of the multi-leaf collimator.

[0055] Step 4: Obtain an optimal set of illumination fields using an evolutionary algorithm;

[0056] Step 4.1: Define the population P of generation G. G The number N of offspring irradiation field schemes off and initialize N off =1;

[0057] Step 4.2: For the Gth generation population P G The three variables undergo crossover and mutation operations, resulting in N variables. pop The Gth generation population P G offspring population O G ;

[0058] Step 4.3: Transfer the Gth generation population P G and its offspring population O G After merging, a new population P′ of generation G is obtained. G ;

[0059] Step 4.4: Based on the overall objective function f(x), calculate the new population P in the Gth generation. G Perform non-dominated sorting to obtain M non-dominated front surfaces, and denote the total number of illumination field schemes located on the first M-1 front surfaces as z;

[0060] Step 4.5: Calculate the crowding distance of each illumination field scheme in the M non-dominated front surfaces, which is used to sort the illumination field schemes of the Mth front surface in descending order, and take the top N of the descending order. pop The -z illumination field schemes are merged with the first M-1 illumination field schemes of the leading edge to obtain the population P of the G+1 generation. G+1 ;

[0061] Step 4.6: Assign G+1 to G, and repeat step 4.1 sequentially until G > G. max Until then, we obtain the Gth... max Generation population

[0062] Step 4.7, put the Gth maxGeneration population The illumination field schemes located on the first Pareto front are selected as the optimal illumination field set and output.

[0063] The illumination field optimization method based on deep neural networks and evolutionary algorithms described in this invention is characterized in that step 4.2 is performed as follows:

[0064] Step 4.2.1: Transfer the Gth generation population P G N pop Each illumination field scheme is randomly shuffled and placed into a mating pool;

[0065] Step 4.2.2: Randomly select two schemes from the mating pool. and As two parent generations;

[0066] Step 4.2.3, for and The radiation intensity variable in the data is cross-referenced and mutated to produce the Nth variant. off The radiation intensity variables of each offspring: where the crossover method is simulated binary crossover, and the mutation method is polynomial mutation;

[0067] Step 4.2.4, for and The number of subfields in the variable is cross-crossed and mutated to produce the Nth mutated subfield. off The variable is the number of offspring fields; the crossover method is simulated binary crossover; the mutation method is to randomly select the Nth offspring field after crossover based on the fitness FN obtained in step 3.2. off The number of offspring fields in each generation is determined by two binary codes that are 1. The binary code with the larger fitness value is flipped, and the Nth crossover is randomly selected. off In the variable representing the number of offspring fields of each offspring, two binary codes with a value of 0 are flipped, and the binary code with the smaller fitness value is flipped.

[0068] Step 4.2.5, for and The leaf shape variables in the model are cross-crossed and mutated to produce the Nth mutated version. off The leaf shape variable of each offspring; where the crossing method is to cross the Nth generation. off The binary encoding of each child is set to the same as that of the parent. Similarly, based on the fitness FH obtained in step 3.2, from Two binary codes with a value of 1 are randomly selected from the leaf shape mask, and the binary code with the larger fitness value is flipped. Two binary codes with a value of 0 are randomly selected from the leaf shape mask, and the binary code with the smaller fitness value is flipped; where, Indicates the parent generation after inversion.

[0069] The mutation method is the same as that in 4.2.4;

[0070] Step 4.2.6, N off +1 is assigned to N eff Repeat step 4.2.2 sequentially until N is reached. off >N pop until.

[0071] The present invention provides an electronic device, including a memory and a processor, wherein the memory is used to store a program that supports the processor in executing the illumination field optimization method, and the processor is configured to execute the program stored in the memory.

[0072] The present invention discloses a computer-readable storage medium on which a computer program is stored, wherein the computer program is executed by a processor to perform the steps of the illumination field optimization method.

[0073] Compared with existing technologies, the beneficial effects of this invention are reflected in:

[0074] 1. Unlike other methods, this invention optimizes the number of light irradiations, replacing the traditional method of having physicists define the number of irradiations themselves. This improves the automation of the entire process. At the same time, the multi-objective evolutionary algorithm provides more solutions instead of a single solution to choose from, which greatly saves the time for solution development.

[0075] 2. The dose prediction model built in this invention enhances feature extraction, can make full use of image context information for prediction, and the upsampled restored image can locate the target object, which meets the needs of medical image segmentation and dose prediction.

[0076] 3. This invention employs a unique crossover mutation operator to generate offspring illumination fields, which is more suitable for illumination field optimization compared to traditional methods.

[0077] 4. This invention fully considers the relationship between the shape of the blades of the multi-leaf collimator and the position of the irradiated light, making the optimization of irradiation intensity and the optimization of blade shape more in line with reality. Attached Figure Description

[0078] Figure 1 This is a schematic diagram of the network architecture used.

[0079] Figure 2 This is a flowchart of the present invention;

[0080] Figure 3 is a schematic diagram of the encoding mode of the present invention;

[0081] Figure 4 is a schematic diagram of the two-step method. DETAILED DESCRIPTION

[0082] In this embodiment, for an irradiation field optimization method based on a deep neural network and an evolutionary algorithm, dose prediction is implemented through a deep neural network first, and then segmental segmentation is performed on the predicted dose prediction map through an evolutionary algorithm to obtain a plurality of optimal irradiation field set schemes for physicists to select, which effectively improves the speed and accuracy of the generated optimal irradiation field scheme set. For the specific process, refer to Figure 2 , which includes:

[0083] Step 1: preprocess data;

[0084] acquire one CT image and H-1 organ contour images respectively, and perform normalization processing on them separately to obtain a preprocessed CT image and an organ contour image set; wherein H is the total number of images; acquire an optimal dose map GT;

[0085] Step 2: build a dose prediction network based on a 3D-Unet network, which is used to predict the dose of each target area, and establish a constraint equation according to the predicted dose;

[0086] Step 2.1: the dose prediction network is composed of an upsampling path, a downsampling path, and a convolutional layer with a dimension of b×b×b, and skip connections are used to introduce underlying features in the downsampling path into the upsampling path. The dimension of the input image is defined as Y×Y×Y, wherein Y is the size of the input image; b represents the size of the convolution kernel; the weights of the convolutional layer are initialized according to a normal distribution according to the MSRA method, the batch size of the Adam optimizer is set to α, the weight decay is set to β, the initial learning rate is set to γ, and the maximum number of iterations is defined as E max , and the current number of iterations is E.

[0087] Step 2.1.1: the downsampling path is composed of T identical downsampling modules, and each downsampling module is sequentially composed of a repeated downsampling layer, a convolutional layer, an instance normalization layer and a rectified linear unit. After the input feature map is processed by each downsampling module, the resolution is reduced by half and the number of channels is doubled; for example Figure 1 as shown, if the number of channels of the first module of the downsampling path is 32, the number of channels of the second module is 64;

[0088] the CT image and the organ contour image set are input into the dose prediction network, and after the encoding processing of the downsampling path, 2 T ·H sheets (Y·2 -T )×(Y·2-T ) × (Y·2 -T ) feature maps;

[0089] Step 2.1.2, the upsampling path consists of T identical upsampling modules, and each upsampling module is sequentially composed of a repeated upsampling layer, a skip concatenation, a convolution layer, an instance normalization layer and a rectified linear unit;

[0090] Step 2.1.3, 2 T ·H (Y·2 -T ) × (Y·2 -T ) × (Y·2 -T ) feature maps are input into the upsampling path for decoding processing, and then an H×Y×Y×Y feature map is obtained;

[0091] Step 2.1.4, the H×Y×Y×Y feature map is processed by a convolution layer with dimension b×b×b, and a predicted dose map is generated, and dose values of all voxels are extracted from the predicted dose map, thereby obtaining a dose value prediction set d={d₁, d₂ … d n …d N}; wherein, d n is a dose value predicted for the nth voxel; N is the number of voxels used for a target region in the CT image;

[0092] Step 2.2, construct the loss L(N, GT) between the predicted dose map d and the optimal dose map GT using formula (1), use it as a supervision loss for training the dose prediction network, update the weight parameters of the network through a back propagation algorithm, and stop network training when the maximum number of iterations is reached, thereby obtaining an optimal dose prediction network;

[0093]

[0094] In formula (2), GT n is the optimal dose value of the nth voxel in the optimal dose map GT;

[0095] Step 2.3, input a CT image to be predicted into the optimal dose prediction network for processing, and obtain the dose value d of voxels of the optimal dose prediction map pre ;

[0096] Step 2.4, construct a constraint equation using formula (2);

[0097]

[0098] In formula (2), ε is an optimization threshold; and:

[0099]

[0100] yg =a j h j，g (4)

[0101] In equations (3) and (4), B represents the number of irradiations, and D represents the number of irradiations. g，n Let D be the dose precipitation coefficient at the g-th sampling point of the n-th voxel, and D g，n ∈D, where D is the dose-precipitation matrix, C j Let y be the set of sampling points for the j-th irradiation. g Let a be the irradiation intensity at the g-th sampling point in the irradiation field scheme to be optimized. j h is the irradiation intensity of the j-th irradiation. j，g h represents the leaf shading at the g-th sampling point during the i-th irradiation. j，g =0 indicates that the j-th irradiation at the g-th sampling point is blocked by the blades of the multi-leaf collimator; for example Figure 4 As shown, this is a two-step process for dose prediction and generating an irradiation field by dividing the subfield using an evolutionary algorithm.

[0102] Step 3: Determine the objective function equation and initialize N. pop A scheme for assembling illumination fields;

[0103] Step 3.1: Construct the overall objective function F(x):

[0104] Step 3.1.1: Construct the minimum dose equation MinDose(x, V, d) for each voxel using equation (5). min ):

[0105]

[0106] In equation (5), x represents the irradiation field scheme using three-layer coding, where three-layer coding refers to coding the irradiation intensity, the number of subfields, and the leaf shape separately; d min The minimum dose for setting dose constraints;

[0107] Step 3.1.2: Construct the maximum dose constraint equation MaxDose(x, V, d) for each voxel using equation (6). max ):

[0108]

[0109] In equation (6), d max The maximum dose for the set dose constraint;

[0110] Step 3.1.3: Construct the dose-volume constraint equation DVC(x, V, ddvc, r) to restrict all voxels using equation (7):

[0111]

[0112] In equation (7), d dvc For the maximum dose under the set dose-volume constraint, d′ r This represents the predicted dose value of the r-th voxel after sorting d in descending order;

[0113] Step 3.1.4: Construct the first objective function f1(x) using equation (8):

[0114]

[0115] Step 3.1.5: Construct the second objective function f2(x) using equation (9):

[0116]

[0117] In equation (9), PTV is the planned volume of the target cells in the CT image;

[0118] Step 3.1.6: Construct the third objective function f3(x) using equation (10):

[0119]

[0120] In equation (10), OAR is the normal region in the CT image;

[0121] Step 3.1.7: Construct the third objective function f4(x) using equation (11):

[0122]

[0123] Step 3.1.8: Construct the overall objective function F(x) using equation (12):

[0124] Minimize F(x)=(f1(x),…,f4(x)) (12)

[0125] Step 3.2: Calculate the fitness of the subfield quantity variable on the objective function in each dimension: FN = {FN1, FN2, ..., FN} s …FN s The set of fitness values ​​of each dimension of the blade shape variable on the objective function is FH = {FH1, FH2, ..., FH}. s …FH s}, the set of fitness values ​​FH for each row of the leaf shape variable in the s-th dimension. s ={FH s，1 FH s，2 …FH s，k …FH s，U}, the set of fitness values ​​FH for each sampling point in the k-th row of the s-th dimension of the two-dimensional matrix of blade shape variables. s，k ={FH s，k，1 FH s，k，2 …FH s，k，w …FH s，U，W}, where FN s FH represents the fitness of the s-th dimension subfield quantity variable on the objective function. s FH represents the matrix representing the fitness of the s-th dimension blade shape variable on the objective function. sk FH represents the set of fitness scores of the sampling points in the k-th row of the s-th dimension of a two-dimensional matrix representing blade shape variables on the objective function. skw The fitness of the sampling point in the k-th row and w-th column of the s-th dimension of the two-dimensional matrix representing the leaf shape variable is on the objective function; S is the total number of illumination rays from different angles as the upper limit of the dimension; the leaf shape is a matrix of size U×W, where U is the number of rows of the matrix and W is the number of columns of the matrix; the fitness value reflects the importance of the dimension variable on the objective function, and the smaller the fitness value, the more important the dimension variable is.

[0126] Step 3.3: Define the maximum number of iterations as G. max The current iteration number is G, and G is initialized to 1;

[0127] Define the population of generation G Let N represent the i-th irradiation field scheme of the G-th generation, where each irradiation field scheme is composed of the variables of subfield quantity, leaf shape, and irradiation intensity; pop The number of schemes in the population;

[0128] Step 3.3.1: Define and randomly initialize the irradiation intensity variable of the population, where the irradiation intensity variable of the i-th irradiation field scheme in the G-th generation is... in, Let represent the s-th real number encoding in the i-th illumination field scheme, and a max This is the maximum irradiation intensity; such as Figure 3 The irradiance map of part b in the example matrix is ​​shown. This indicates that the irradiation intensity on the first subfield in the i-th irradiation field scheme of generation G is 0.1;

[0129] Step 3.3.2: Define and randomly initialize the population P of generation G. G The variable of the number of subfields, where the population P in the Gth generation. G The i-th illumination field scheme Binary encoding of the subfield quantity variable and P represents the population of generation G. G The i-th illumination field scheme The s-th binary code of the subfield quantity variable, and 0 represents the population P of generation G. G The s-th subfield of the i-th illumination field scheme does not exist, and 1 indicates that the population P of the G-th generation does not exist. G The s-th subfield exists in the i-th illumination field scheme; such as Figure 3 The number of subfields in part a is shown in the example matrix. This indicates that the first subfield does not exist in the i-th illumination field scheme of generation G;

[0130] Step 3.3.3: Define and randomly generate the population P in the Gth generation. G Leaf shape variables, where the Gth generation population P G The i-th illumination field scheme Leaf shape variable in, For the i-th illumination field scheme The s-th dimension of the two-dimensional matrix of the blade shape variable. Represents the i-th illumination field scheme The k-th row of the s-th dimension of the two-dimensional matrix of the leaf shape variable. Represents the i-th illumination field scheme Is the sampling point located in the k-th row and w-th column of the s-th dimension of the two-dimensional matrix of the blade shape variable obscured by the blades of the multi-leaf collimator? A value of 0 indicates the i-th illumination field scheme. The sampling point located in the k-th row and w-th column of the s-th dimension of the two-dimensional matrix of the blade shape variable is obscured by the blades of the multi-leaf collimator; such as Figure 3 The blade shape diagram in part c is shown in the example matrix. This indicates that the location of the first sampling point in the k-th row of the s-dimensional two-dimensional matrix of the i-th illumination field scheme of the G-th generation is blocked by the blades of the multi-leaf collimator.

[0131] Step 4: Obtain an optimal set of illumination fields using an evolutionary algorithm;

[0132] Step 4.1: Define the maximum number of iterations as G. max The current iteration number is G, and G is initialized to 1;

[0133] Step 4.2: Define the population P of generation G. G The number N of offspring irradiation field schemes off and initialize N off =1;

[0134] Step 4.3: For the Gth generation population P G The three variables undergo crossover and mutation operations, resulting in N variables. pop The Gth generation population P G offspring population O G ;

[0135] Step 4.3.1: Transfer the Gth generation population P G N pop Each illumination field scheme is randomly shuffled and placed into a mating pool;

[0136] Step 4.3.2: Randomly select two schemes from the mating pool. and As two parent generations;

[0137] Step 4.3.3, for and The radiation intensity variable in the data is cross-referenced and mutated to produce the Nth variant. off The radiation intensity variable for each offspring:

[0138] The crossover method is simulated binary crossover, and the mutation method is polynomial mutation;

[0139] Step 4.3.4: Follow the process in step 4.2.2. and The number of subfields in the variable is cross-crossed and mutated to produce the Nth mutated subfield. off The variable representing the number of offspring in each generation:

[0140] The crossover method is a simulated binary crossover, and the mutation method is to mutate according to the fitness FN obtained in step 3.2, and to randomly select the Nth person with the same probability. off For each offspring's child field number variable, if two dimensions are equal to 1, the dimension with the larger fitness value is flipped. Then, with equal probability, two dimensions with zero values ​​are randomly selected, and the dimension with the smaller fitness value is flipped. For example, if the offspring... The offspring number variable is (1, 1, 0, 0, 0). If the dimension selected as 1 is flipped, then dimensions 1 and 2 will be selected. Also, assuming the fitness of dimension 1 is greater than that of dimension 2, the mutated offspring individuals... The number of subfields is (0, 1, 0, 0, 0);

[0141] Step 4.3.5, for and The leaf shape variables in the model are cross-crossed and mutated to produce the Nth mutated version. off Leaf shape variables for each offspring:

[0142] The crossover method is to cross over the offspring The binary encoding is set to be the same as Similarly, based on the fitness FH obtained in step 3.4, the same probability is used to execute the process from... From the leaf shape mask, randomly select two dimensions with a fitness value of 1, and choose the dimension with the larger fitness value for the Nth dimension. off The leaf shape variable of each offspring is flipped, and values ​​of 1 are set to 0. This is repeated with equal probability. In the mask, randomly select two dimensions that are 0, and choose the dimension with the smaller fitness value to pair with the Nth dimension. off The leaf shape variable of each offspring is inverted, setting values ​​of 0 to 1. Among these, Indicates the parent generation after inversion.

[0143] The mutation method is the same as that in 4.3.4; for example, if the parent... and The leaf shape variables are (0, 1, 0, 1, 0) and (1, 0, 1, 0, 0), respectively. When dimension 1 is flipped, dimensions 2 and 4 will be selected. Also, assuming the fitness value of dimension 2 is less than the fitness value of dimension 4, then the offspring... The leaf shape variable is (0, 1, 0, 0, 0); when the offspring selects dimension 0 for flipping, dimensions 1 and 3 will be selected. Also, assuming the fitness value of dimension 1 is less than the fitness value of dimension 3, then the offspring... The leaf shape variable is (1, 1, 0, 1, 0);

[0144] Step 4.3.6, N off +1 is assigned to N off Repeat step 4.3.2 sequentially until N. off >N pop until;

[0145] Step 4.4: Transfer the Gth generation population P G and its offspring population O G After merging, a new population P′ of generation G is obtained. G ;

[0146] Step 4.5: Based on the overall objective function f(x), calculate the new population P′ in the Gth generation. G Perform non-dominated sorting to obtain M non-dominated front surfaces, and denote the total number of illumination field schemes located on the first M-1 front surfaces as Z;

[0147] Step 4.6: Calculate the crowding distance of each illumination field scheme in the M non-dominated front surfaces, which is used to sort the illumination field schemes of the Mth front surface in descending order, and take the top N of the descending order. popThe Z illumination field schemes are merged with the first M-1 illumination field schemes of the leading edge to obtain the population P of the G+1 generation. G+1 ;

[0148] Step 4.7: Assign G+1 to G, and repeat step 4.2 sequentially until G > G. max Until then, we obtain the Gth... max Generation population

[0149] Step 4.8, put the Gth max Generation population The illumination field schemes located on the first Pareto front are selected as the optimal illumination field set and output.

[0150] In this embodiment, an electronic device includes a memory and a processor. The memory stores a program that supports the processor in executing the above-described method, and the processor is configured to execute the program stored in the memory.

[0151] In this embodiment, a computer-readable storage medium stores a computer program, which is executed by a processor to perform the steps of the above method.

Claims

1. A method for optimizing illumination field based on deep neural networks and evolutionary algorithms, characterized in that, It is done in the following steps: Step 1: Preprocess the data; Acquire a CT image and The organ contour images were processed and normalized to obtain preprocessed CT images and organ contour image sets; among them... Total number of images; obtain the optimal dose map ; Step 2: Build a dose prediction network based on the 3D-Unet network to predict the dose in each target region, and establish constraint equations based on the predicted dose; Step 2.1: The dose prediction network consists of an upsampling path, a downsampling path, and a dimension of... It consists of convolutional layers, and uses jumpers to introduce low-level features from the downsampling path into the upsampling path, defining the dimension of the input image as... Where Y is the size of the input image; Indicates the size of the convolution kernel; Step 2.1.1, the downsampling path is... The system consists of several identical downsampling modules, each of which is composed of a repeated downsampling layer, a convolutional layer, an instance normalization layer, and a rectified linear unit. After each downsampling module processes the input feature map, the resolution is reduced by half and the number of channels is doubled. The CT images and organ contour image set are input into the dose prediction network, and after encoding processing by the downsampling path, the following is obtained: open( ×( ×( Feature map; Step 2.1.2, the upsampling path is... It consists of identical upsampling modules, and each upsampling module is composed of repeating upsampling layers, skipping cascades, convolutional layers, instance normalization layers and rectified linear units in sequence; Step 2.1.3 open( ×( ×( After the feature map is input into the upsampling path for decoding, the result is obtained. Feature map; Step 2.1.4 The feature map is processed by a dimension of The convolutional layer is processed to generate a predicted dose map, and the dose values ​​of all voxels are extracted from the predicted dose map to obtain a dose value prediction set. ;in, For the first Dosage values ​​predicted by individual factors; The number of voxels used in the target region of the CT image; Step 2.2: Construct the predicted dose map using equation (1). With optimal dose diagram Losses between The weight parameters of the network are updated through backpropagation and used as a supervisory loss for training the dose prediction network. The network training is stopped when the maximum number of iterations is reached, thereby obtaining the optimal dose prediction network. = (1) In equation (1), Optimal Dosage Chart Middle The optimal dose of individual pigments; Step 2.3: Input the CT image to be predicted into the optimal dose prediction network for processing, and obtain the dose values ​​of the voxels in the optimal dose prediction map. ; Step 2.4: Construct constraint equations using equation (2); (2) In equation (2), To optimize the threshold; and we have: = (3) (4) In equations (3) and (4), For the number of irradiations, For the first The first individual element The dose deposition coefficient at each sampling point, and , This is a dose-precipitation matrix. For the first The set of sampling points for the second irradiation. The first irradiation field scheme to be optimized Irradiation intensity at each sampling point For the first The intensity of the irradiation. For the first During the second irradiation Leaf occlusion at each sampling point =0 indicates the first The second irradiation One sampling point was blocked by the blades of the multi-leaf collimator; Step 3: Determine the objective function equation and initialize it. A scheme for assembling illumination fields; Step 3.1: Construct the overall objective function : Step 3.1.1: Construct the minimum dose equation for each voxel using equation (5). : (5) In equation (5), The three-layer coding scheme refers to coding the irradiation intensity, the number of subfields, and the leaf shape respectively. The minimum dose for setting dose constraints; Step 3.1.2: Construct the maximum dose constraint equation for each voxel using equation (6). : (6) In equation (6), The maximum dose for the set dose constraint; Step 3.1.3: Construct dose-volume constraint equations to restrict all voxels using equation (7). : (7) In equation (7), The maximum dose for the set dose-volume constraint. Indicates to The dose value predicted for the r-th voxel after sorting in descending order; Step 3.1.4: Construct the first objective function using equation (8). : (8) Step 3.1.5: Construct the second objective function using equation (9). : (9) In equation (9), The planned volume of the target cells in the CT image; Step 3.1.6: Construct the third objective function using equation (10). : (10) In equation (10), This refers to a normal area in a CT image. Step 3.1.7: Construct the third objective function using equation (11). : (11) Step 3.1.8: Construct the overall objective function using equation (12). : (12) Step 3.2: Calculate the fitness of the subfield quantity variable on the objective function in each dimension. { , The set of fitness values ​​of each dimension of the leaf shape variable on the objective function. { , The set of fitness values ​​for each row of the leaf shape variable in the s-th dimension. ={ The s-th dimension of the two-dimensional matrix of blade shape variables The set of fitness values ​​for each sampling point in the row ={ ,in, Indicates the first The fitness of the vidin field quantity variable on the objective function. Indicates the first The fitness matrix of the blade shape variable on the objective function. The s-th dimension of the two-dimensional matrix representing the blade shape variable is... The set of fitness scores of the sampled points in a row on the objective function. The fitness of the sampling point in the k-th row and w-th column of the s-th dimension of the two-dimensional matrix representing the blade shape variable on the objective function; The upper limit of the dimension is the total number of rays from different angles. The shape of the leaf is a matrix of size U×W, where U is the number of rows of the matrix and W is the number of columns of the matrix. Step 3.3: Define the maximum number of iterations as... The current iteration number is and initialize ; Define the population of generation G ={ }, Indicates the first The generation Each irradiation field scheme consists of a subfield number variable, a leaf shape variable, and an irradiation intensity variable. The number of schemes in the population; Step 3.3.1: Define and randomly initialize the population of generation G. The irradiation intensity variable, where the Gth generation population Middle One illumination field scheme Irradiation intensity variable ={ ， ,in, Represents the Gth generation population Middle One illumination field scheme Middle Encoding of real numbers, and , This represents the maximum irradiation intensity. Step 3.3.2: Define and randomly initialize the population of generation G. The variable of the number of offspring, where the population in generation G. Middle One illumination field scheme Binary encoding of the subfield quantity variable ,and ={ ， , Represents the Gth generation population Middle One illumination field scheme The first of the subfield quantity variables A binary code, and 0 indicates the Gth generation population Middle The first illumination field scheme The field does not exist; 1 indicates the population of generation G. Middle The first illumination field scheme A wild individual exists; Step 3.3.3: Define and randomly generate the population for generation G. Leaf shape variables, of which, the Gth generation population Middle One illumination field scheme Leaf shape variable ={ ， } ,in, For the first One illumination field scheme The s-th dimension of the two-dimensional matrix of the blade shape variable. ={ , Indicates the first One illumination field scheme The k-th row of the s-th dimension of the two-dimensional matrix of the blade shape variable. ={ , Indicates the first One illumination field scheme Is the sampling point located in the k-th row and w-th column of the s-th dimension of the two-dimensional matrix of the blade shape variable obscured by the blades of the multi-leaf collimator? A value of 0 indicates the first point. One illumination field scheme The sampling point in the kth row and wth column of the s-th dimension of the two-dimensional matrix of the blade shape variable is blocked by the blade of the multi-leaf collimator. Step 4: Obtain an optimal set of illumination fields using an evolutionary algorithm; Step 4.1, Define the first Generation population Number of offspring irradiation field schemes and initialize =1; Step 4.2, for the first Generation population The three variables undergo crossover and mutation operations, resulting in a number of... The Generation population offspring population ; Step 4.3, place the first Generation population and its offspring population After merging, we get the first New populations of the generation ; Step 4.4: Based on the overall objective function , for the New populations of the generation Perform a non-dominated sort to obtain A non-dominated frontier surface, and will be located in the front The total number of illumination field schemes for each front surface is denoted as ; Step 4.5, Calculation The crowding distance of each illumination field scheme in each non-dominated frontal surface is used to determine the first... The illumination field schemes of the leading edge are sorted in descending order, and the schemes in descending order are selected. The illumination field scheme and the previous The illumination field schemes of the first leading edge are merged to obtain the second... Generation population ; Step 4.6, +1 is assigned to Repeat step 4.1 sequentially until... > Until then, thus obtaining the first Generation population ; Step 4.7, place the first Generation population The illumination field schemes located on the first Pareto front are selected as the optimal illumination field set and output.

2. The illumination field optimization method based on deep neural networks and evolutionary algorithms according to claim 1, characterized in that, Step 4.2 is performed as follows: Step 4.2.1, place the first Generation population of Each illumination field scheme is randomly shuffled and placed into a mating pool; Step 4.2.2: Randomly select two schemes from the mating pool. and As two parent generations; Step 4.2.3, for and The irradiation intensity variable in the data is cross-referenced and mutated to produce the mutated result. The radiation intensity variables of each offspring: where the crossover method is simulated binary crossover, and the mutation method is polynomial mutation; Step 4.2.4, for and The number of subfields in the variable is cross-crossed and mutated to produce the mutated number of subfields. The variable is the number of offspring fields; the crossover method is simulated binary crossover; the mutation method is based on the fitness obtained in step 3.

2. Randomly select the first crossover. The number of offspring fields in each generation is determined by two binary codes that are 1. The binary code with the larger fitness value is flipped, and the resulting crossover is randomly selected. In the variable representing the number of offspring fields of each offspring, two binary codes with a value of 0 are flipped, and the binary code with the smaller fitness value is flipped. Step 4.2.5, for and The leaf shape variable in the model is cross-crossed and mutated to produce the mutated first... The leaf shape variable of each offspring; where the crossing method is to cross the first generation. The binary encoding of each child is set to the same as that of the parent. Similarly, based on the fitness obtained in step 3.2 ,from Two binary codes with a value of 1 are randomly selected from the leaf shape mask, and the binary code with the larger fitness value is flipped. Two binary codes with a value of 0 are randomly selected from the leaf shape mask, and the binary code with the smaller fitness value is flipped; where, Indicates the parent generation after inversion. ; The mutation method is the same as that in 4.2.4; Step 4.2.6, will +1 is assigned to Repeat step 4.2.2 sequentially until... > until.

3. An electronic device, comprising a memory and a processor, characterized in that, The memory is used to store a program that supports the processor in executing the irradiation field optimization method of claim 1 or 2, and the processor is configured to execute the program stored in the memory.

4. A computer-readable storage medium storing a computer program thereon, characterized in that, The computer program is executed by the processor to perform the steps of the illumination field optimization method according to claim 1 or 2.

Images