Radiotherapy robot multi-dimensional radiotherapy path planning method based on RRT algorithm

Abstract

The invention provides a radiotherapy robot multi-dimensional radiotherapy path planning method based on an RRT algorithm. The method comprises the steps of firstly modeling a radiotherapy path planning problem of a dynamic tumor, then fitting a quasi-periodic motion model of the tumor by utilizing an ant colony algorithm, constructing a target function of radiotherapy path point selection by fusing the model, and finally providing a radiotherapy path point search strategy based on an RRT algorithm. And a dose distribution evaluation index is formulated to evaluate the optimization effect of the radiotherapy path points. Compared with the prior art, the method has the advantages that the tumor position, the radiotherapy path point and the radiation dose are unified into a multi-dimensional variable, factors influencing the radiotherapy effect are optimized at the same time, and the radiotherapy precision and the applicability to dynamic tumors in the human body breathing movement process can be improved.

Description

Technical Field

[0001] The present invention relates to the technical field of medical robots, and in particular to a multi-dimensional radiotherapy path planning method for a radiotherapy robot based on an RRT algorithm, which can solve the beam direction optimization problem of a precision radiotherapy robot facing a dynamic tumor. Background Art

[0002] The beam angle orientation optimization problem (BAO) is a combinatorial optimization problem. The traditional approach is to iteratively search for a set of beam combinations through a heuristic algorithm. However, due to the large number of spatial combinations of beam angles, the amount of calculation is too large, especially when considering multiple constraints at the same time, it is easy to fall into the local optimum.

[0003] With the advancement of science and technology, the traditional method of manually selecting beam directions based on physician experience or fixed templates has been basically replaced by automatic beam angle optimization algorithms. Traditional heuristic algorithms work well in generating simple coplanar spot radiotherapy plans, but with the development of radiotherapy equipment, non-coplanar radiotherapy and VMAT continuous radiotherapy technology have gradually dominated radiotherapy operations. When faced with complex beam solution spaces, heuristic algorithms have too much computation and are prone to fall into local optimality. Although beam optimization methods based on deep learning and reinforcement learning have good adaptability, they require a large amount of training data and computing resources, and do not take into account the impact of respiratory motion on radiotherapy accuracy.

[0004] Although many effective solutions have been proposed for BAO, most of them are for static tumors. The few solutions for dynamic tumors do not consider the radiation source, tumor, and dose as a whole, resulting in insufficient adaptability of radiotherapy solutions to dynamic tumors. Therefore, studying a radiotherapy path planning method that can take tumor respiratory motion into account will help significantly improve radiotherapy accuracy and reduce the radiation damage to healthy human tissues caused by radiotherapy point errors. Summary of the invention

[0005] The purpose of the present invention is to provide a multi-dimensional radiotherapy path planning method for a radiotherapy robot based on an RRT algorithm, so as to solve the beam direction optimization problem of a precise radiotherapy robot when facing a dynamic tumor.

[0006] The technical solution of the present invention is: a multi-dimensional radiotherapy path planning method for a radiotherapy robot based on an RRT algorithm, comprising the following steps:

[0007] Step 1: The tumor position, radiation source position, and radiation dose of the corresponding radiation point are uniformly optimized as a multidimensional variable. According to the spatial accessibility of the radiotherapy arm, the radiotherapy sphere is selected as the upper hemisphere centered on the tumor, and the radiotherapy path points are optimized;

[0008] Step 2: According to the fluctuation of the tumor during human respiratory movement, a quasi-periodic respiratory motion model of the tumor is fitted based on the ant colony algorithm;

[0009] Step 3: Construct an objective function for selecting radiotherapy path points based on the dose constraints of risk organs, radiotherapy target areas, and healthy tissues;

[0010] Step 4: Optimize and select the radiotherapy path point sequence based on the RRT algorithm;

[0011] In the above technical solution, in step 1, according to the physiological displacement of the tumor caused by human respiratory movement, the tumor position, radiation source position and radiation dose are considered uniformly and constructed as a multidimensional variable for optimization.

[0012] In the above technical solution, in order to avoid collision between the robot arm and other equipment, the radiotherapy sphere is selected as the upper hemisphere, and the direction angle is determined by the medical method. represents the position of the path point on the radiotherapy sphere.

[0013] In the above technical solution, in step 2, the effect of tumor respiratory motion on dose delivery is considered, the periodic characteristics of tumor respiratory motion are fitted based on the ant colony algorithm (ACO), and the amplitude, phase and frequency of tumor motion are analyzed.

[0014] Since the influence of human respiratory movement will cause dynamic displacement of tumors in related parts, in order to improve the accuracy of radiotherapy, it is necessary to integrate the respiratory movement of the tumor into the optimization problem of radiotherapy path points, and plan the radiotherapy path points and the corresponding tumor positions at the same time. In the radiotherapy simulation platform, the model data of the target area and healthy tissue are composed of voxel units. The voxel units are transformed in the same position as the voxel units near the tumor according to the quasi-periodic motion model of the tumor to simulate dynamic radiotherapy and calculate the dynamic dose distribution. Therefore, a quasi-periodic model of tumor motion is required. Considering the influence of baseline drift, the quasi-periodic respiratory motion model of the tumor can be expressed as:

[0015]

[0016] Where τ(t) is the coordinate value of the tumor principal element in the Z-axis direction, the posture parameter C0 represents the baseline position, C1 represents the baseline drift speed, A is the oscillation amplitude, and t represents time. is the phase shift, f i Represents the oscillation frequency, all posture parameters have the same oscillation frequency, N represents the number of sin functions of different frequencies and phases contained in the tumor quasi-periodic respiratory motion model, and its value represents the steepness and flatness of the model shape. The parameters in the above formula can be obtained by fitting preoperative data.

[0017] In the above technical solution, in step 3, a multidimensional radiotherapy path point sequence X selecting objective function integrating respiratory motion tracking is constructed according to the dose constraints of the tumor target area, risk organs, and healthy tissues, while taking into account the goal of optimizing radiotherapy time.

[0018] The target area tissue optimization function is defined as:

[0019]

[0020] The risk area optimization function is defined as:

[0021]

[0022] The general organization optimization function is defined as follows:

[0023]

[0024] where d pre is the prescription dose of the target area, d max is the tolerance dose for the risk organ, is the upper dose penalty factor for risk organs. PTV 、N OAR 、N nor are the number of voxel units in the target area, risk organ, and normal tissue, respectively. L is the number of path points, d ij is the dose impact matrix, which indicates the radiation dose received by the i-th voxel unit in the model from the j-th beam. This matrix can be used to optimize the radiotherapy plan so that the side effects on the surrounding normal tissues are minimized while treating the tumor tissue. j It is the dose weight of the path point, which is used for the subsequent control point dose allocation.

[0025] (x) + is a step function:

[0026]

[0027] Where x is the actual irradiation dose. For the target area, when the actual irradiation dose reaches or exceeds the prescribed dose, the step function value is 0, and the increase in the irradiation dose objective function value will not be affected. When the actual irradiation received by the target area is lower than the prescribed dose, the step function value is 1. At this time, the objective function value is the variance of the actual dose and the prescribed dose, so the actual irradiation dose will be adjusted in the direction of increasing the target area irradiation dose. For risk organs, when the actual irradiation dose exceeds the maximum dose constraint, the step function is 1. At this time, the greater the deviation between the actual dose and the maximum dose constraint, the greater the objective function value, so that the actual dose is adjusted in the direction of reducing the dose of the risk organ. When the actual irradiation dose is lower than the maximum dose constraint, the step function is 0 and no longer affects the overall objective function value.

[0028] The time objective function is defined as:

[0029] f time (X) = min(d(ω n-1 ,r i ) / ω) (6)

[0030] where r i is the i-th path point sampled, ω n-1 is the path point at the previous moment, d(ω n-1 ,r i ) represents the angle difference between the current sampling path point and the path point at the previous moment, and ω is the angular velocity of the radiotherapy arm.

[0031] The overall objective function is:

[0032] F(X)=∑ r∈R δ r f r (X) (7)

[0033] Where X represents the multidimensional radiotherapy path point sequence, R represents the real number set, and f r (X) is the rth objective function, δ r is the weight of the rth objective function.

[0034] After the BAO problem is determined, different irradiation doses need to be assigned to the selected radiotherapy points to achieve targeted irradiation. This problem is called flux map optimization (FMO). The problem solved by FMO is to optimize the dose weight ω of the beam at each radiotherapy point so that the deviation between the actual cumulative radiation dose and the expected dose is minimized while satisfying the dose constraint of the organ at risk. The present invention solves the flux map optimization problem by iteratively optimizing the total objective function composed of the weighted objective functions of the various types of regional objective functions mentioned above through the optimizer of the treatment planning platform.

[0035] In the above technical solution, the RRT algorithm is used in step 4 for radiotherapy path planning, and the optimal radiotherapy path point is selected from the sampling points based on the fused radiotherapy path point objective function and added to the radiotherapy sequence. The specific steps are as follows:

[0036] (1): Randomly generate a starting radiotherapy path point and update it to the current radiotherapy path point

[0037] (2): With the current radiotherapy path point as the center, the direction angle corresponding to the current radiotherapy path point Randomly sample b path points within a certain threshold range and remove the sampling points that do not meet the constraints.

[0038] (3): Calculate the overall objective function value of the remaining sampling points combined with the previous path points, select the sampling point with the minimum objective function value and add it to the radiotherapy path sequence, and update the current path point at the same time.

[0039] (4): Determine whether the number of radiotherapy path points is greater than the set threshold, and whether the reduction of the objective function caused by the added radiotherapy path points is less than the given threshold. If both are satisfied, stop the optimization; otherwise, continue searching for additional path points.

[0040] (5): Calculate the dose weights {ω1,ω2,...,ω M}, and assign radiation dose to each path point according to the dose weight, as shown in the following formula:

[0041]

[0042] where γ i is the radiation dose at the ith path point, D pre is the prescribed dose, ω i is the dose weight of the ith path point calculated after flux optimization, and M is the number of path points.

[0043] (6): Generate radiotherapy path point sequence X.

[0044] The advantages of the present invention are:

[0045] The present invention provides a radiotherapy path planning method for precise radiotherapy robots for dynamic tumors, which treats the radiation source position, tumor position, and radiation dose as a multidimensional variable for unified optimization, and formulates a dynamic radiotherapy path point selection objective function that integrates tumor respiratory motion; in order to solve the NP-hard combinatorial optimization problem, the present invention proposes a radiotherapy path point search strategy based on the RRT algorithm, and finally constructs a complete optimization framework for radiotherapy path points for dynamic tumors. Compared with the existing methods, the present invention unifies the tumor position, radiotherapy path points, and radiation dose into a multidimensional variable, and simultaneously optimizes the factors that affect the radiotherapy effect, which can improve the accuracy of radiotherapy and the applicability to dynamic tumors during human respiratory motion. BRIEF DESCRIPTION OF THE DRAWINGS

[0046] The present invention will be further described below in conjunction with the accompanying drawings and embodiments:

[0047] Figure 1 This is a flow chart of a method according to Embodiment 1 of the present invention.

[0048] Figure 2 This is a schematic diagram of screening radiotherapy path points according to the method proposed in the present invention.

[0049] Figure 3 A set of beam directions is selected based on the method proposed in the present invention. DETAILED DESCRIPTION

[0050] In order to achieve the purpose of the present invention, Figure 1 As shown, in one embodiment of the present invention, a multi-dimensional radiotherapy path planning method of a radiotherapy robot based on an RRT algorithm comprises the following steps:

[0051] In step 1, considering that the tumor will undergo physiological displacement due to human respiratory movement, the tumor position, radiation source position and radiation dose are considered uniformly and constructed as a multidimensional variable for optimization.

[0052] In order to avoid collision between the robot arm and other equipment, the radiotherapy sphere is selected as the upper hemisphere, and the direction angle is determined by the medical method. represents the position of the path point on the radiotherapy sphere.

[0053] In step 2, in order to consider the impact of tumor respiratory motion on dose delivery, the periodic characteristics of tumor respiratory motion were fitted based on the ant colony algorithm (ACO), and the amplitude, phase and frequency of tumor motion were analyzed.

[0054] Since the influence of human respiratory movement will cause dynamic displacement of tumors in related parts, in order to improve the accuracy of radiotherapy, it is necessary to integrate the respiratory movement of the tumor into the optimization problem of radiotherapy path points, and plan the radiotherapy path points and the corresponding tumor positions at the same time. In the radiotherapy simulation platform, the model data of the target area and healthy tissue are composed of voxel units. The voxel units are transformed in the same position as the voxel units near the tumor according to the quasi-periodic motion model of the tumor to simulate dynamic radiotherapy and calculate the dynamic dose distribution. Therefore, a quasi-periodic model of tumor motion is required. Considering the influence of baseline drift, the quasi-periodic respiratory motion model of the tumor can be expressed as:

[0055]

[0056] Where τ(t) is the coordinate value of the tumor principal element in the Z-axis direction, and the pose parameter C0 represents the baseline position.

[0057] C1 represents the baseline drift speed, A is the oscillation amplitude, t represents the time, is the phase shift, f i Represents the oscillation frequency, all posture parameters have the same oscillation frequency, N represents the number of sin functions of different frequencies and phases contained in the tumor quasi-periodic respiratory motion model, and its value represents the steepness and flatness of the model shape. The parameters in the above formula can be obtained by fitting preoperative data.

[0058] In step 3, based on the dose constraints of the tumor target, risk organs, and healthy tissues, and considering the goal of optimizing the radiotherapy time, a multidimensional radiotherapy path point sequence X that integrates respiratory motion tracking is constructed to select the objective function.

[0059] The target area tissue optimization function is defined as:

[0060]

[0061] The risk area optimization function is defined as:

[0062]

[0063] The general organization optimization function is defined as follows:

[0064]

[0065] where d pre is the prescription dose of the target area, d max is the tolerance dose of the risk organ, and μ is the upper dose penalty factor of the risk organ. PTV 、N OAR 、N nor are the number of voxel units in the target area, risk organ, and normal tissue, respectively. L is the number of path points, d ijis the dose impact matrix, which indicates the radiation dose received by the i-th voxel unit in the model from the j-th beam. This matrix can be used to optimize the radiotherapy plan so that the side effects on the surrounding normal tissues are minimized while treating the tumor tissue. j It is the dose weight of the path point, which is used for the subsequent control point dose allocation.

[0066] (x) + is a step function:

[0067]

[0068] Where x is the actual irradiation dose. For the target area, when the actual irradiation dose reaches or exceeds the prescribed dose, the step function value is 0, and the increase in the irradiation dose objective function value will not be affected. When the actual irradiation received by the target area is lower than the prescribed dose, the step function value is 1. At this time, the objective function value is the variance of the actual dose and the prescribed dose, so the actual irradiation dose will be adjusted in the direction of increasing the target area irradiation dose. For risk organs, when the actual irradiation dose exceeds the maximum dose constraint, the step function is 1. At this time, the greater the deviation between the actual dose and the maximum dose constraint, the greater the objective function value, so that the actual dose is adjusted in the direction of reducing the dose of the risk organ. When the actual irradiation dose is lower than the maximum dose constraint, the step function is 0 and no longer affects the overall objective function value.

[0069] The time objective function is defined as:

[0070] f time (X) = min(d(ω n-1 ,r i ) / ω) (6)

[0071] where r i is the i-th path point sampled, ω n-1 is the path point at the previous moment, d(ω n-1 ,r i ) represents the angle difference between the current sampling path point and the path point at the previous moment, and ω is the angular velocity of the radiotherapy arm.

[0072] The overall objective function is:

[0073] F(X)=∑ r∈R δ r f r (X) (7)

[0074] Where X represents the multidimensional radiotherapy path point sequence, R represents the real number set, and f r (X) is the rth objective function, δ r is the weight of the rth objective function.

[0075] After the beam angle orientation (BAO) is determined, different irradiation doses need to be assigned to the selected radiotherapy points to achieve targeted irradiation. This problem is called the flux map optimization (FMO) problem. The problem solved by FMO is to optimize the dose weight ω of the beam at each radiotherapy point so that the deviation between the actual radiation cumulative dose and the expected dose is minimized while satisfying the dose constraint of the endangered organ. The present invention solves the flux map optimization problem by iteratively optimizing the total objective function composed of the weighted objective functions of various types of regional objective functions through the optimizer of the radiotherapy simulation platform.

[0076] In step 4, the RRT algorithm is used for radiotherapy path planning. Based on the fused radiotherapy path point objective function, the optimal radiotherapy path point is selected from the sampling points and added to the radiotherapy sequence. To calculate the combined effect between radiotherapy path points. Starting from the empty beam set, by calculating the objective function value after combining the sampling points with the path point sequence in the beam set, radiotherapy path points are iteratively added one by one to the beam set until the set upper limit of the number of radiotherapy points is reached or the added radiotherapy points do not significantly reduce the objective function value. The specific steps are as follows:

[0077] (1): Randomly generate a starting radiotherapy path point and update it as the current radiotherapy path point.

[0078] (2): With the current radiotherapy path point as the center, the direction angle corresponding to the current radiotherapy path point Randomly sample b path points within a certain threshold range, where b is a positive integer, and remove sampling points that do not meet the constraints.

[0079] (3): Calculate the overall objective function value of the remaining sampling points combined with the previous path points, select the sampling point with the minimum objective function value and add it to the radiotherapy path point sequence, and update the current path point at the same time.

[0080] (4): Determine whether the number of radiotherapy path points is greater than the set threshold, and whether the reduction of the objective function caused by the added radiotherapy path points is less than the given threshold. If both are satisfied, stop the optimization; otherwise, continue searching for additional path points. (5):

[0082] Calculate the dose weights {ω1,ω2,...,ω M}, and assign radiation dose to each path point according to the dose weight, as shown in the following formula:

[0083]

[0084] where γ i is the radiation dose at the ith path point, D preis the prescribed dose, ω i is the dose weight of the ith path point calculated after flux optimization, and M is the number of path points.

[0085] (6): Generate radiotherapy path point sequence X.

[0086] Of course, the above embodiments are only for illustrating the technical concept and features of the present invention, and their purpose is to enable people familiar with the technology to understand the content of the present invention and implement it accordingly, and they cannot be used to limit the protection scope of the present invention. Any modifications made according to the spirit of the main technical solution of the present invention should be included in the protection scope of the present invention.

Claims

1. A multi-dimensional radiotherapy path planning method for a radiotherapy robot based on an RRT algorithm, characterized in that: The steps include: Step 1: The tumor position, radiation source position, and radiation dose of the corresponding radiation point are uniformly optimized as a multidimensional variable. According to the spatial accessibility of the radiotherapy arm, the radiotherapy sphere is selected as the upper hemisphere centered on the tumor, and the radiotherapy path points are optimized; Step 2: According to the fluctuation of the tumor during the human respiratory movement, a quasi-periodic respiratory motion model of the tumor is fitted based on the ant colony algorithm; Step 3: Construct an objective function for selecting radiotherapy path points based on the dose constraints of risk organs, radiotherapy target areas, and healthy tissues; Step 4: Optimize and select the radiotherapy path point sequence based on the RRT algorithm.

2. The multi-dimensional radiotherapy path planning method of a radiotherapy robot based on the RRT algorithm according to claim 1, characterized in that: In step 1, based on the physiological displacement of the tumor caused by human respiratory movement, the tumor position, radiation source position and radiation dose are considered uniformly and constructed as a multidimensional variable for optimization.

3. The multi-dimensional radiotherapy path planning method of the radiotherapy robot based on the RRT algorithm according to claim 1 is characterized in that: In step 1, the radiotherapy sphere is selected as the upper hemisphere, and the direction angle is determined by the medical method. represents the position of the path point on the radiotherapy sphere; θ∈[0°,360°).

4. The multi-dimensional radiotherapy path planning method of a radiotherapy robot based on the RRT algorithm according to claim 1, characterized in that: In step 2, the respiratory motion of the tumor is integrated into the radiotherapy path point optimization problem, and the radiotherapy path points and the corresponding tumor positions are planned at the same time; according to the influence of baseline drift, the quasi-periodic respiratory motion model of the tumor is expressed as: Where τ(t) is the coordinate value of the tumor principal element in the Z-axis direction, the posture parameter C0 represents the baseline position, C1 represents the baseline drift speed, A is the oscillation amplitude, and t represents time. is the phase shift, f i Represents the oscillation frequency, all posture parameters have the same oscillation frequency, N represents the number of sin functions of different frequencies and phases contained in the quasi-periodic respiratory motion model of the tumor, and its value represents the steepness and flatness of the model shape.

5. The multi-dimensional radiotherapy path planning method of the radiotherapy robot based on the RRT algorithm according to claim 1, characterized in that: In step 2, the fitted tumor quasi-periodic respiratory motion model is used as the input of the radiotherapy simulation platform, and the voxel units near the tumor are subjected to the same position transformation according to the quasi-periodic respiratory motion model of the tumor to realize the simulation of dynamic radiotherapy and calculate the dynamic dose distribution.

6. The multi-dimensional radiotherapy path planning method of a radiotherapy robot based on the RRT algorithm according to claim 1, characterized in that: In step 3, the objective function for selecting the radiotherapy path point is constructed based on the dose constraints of the tumor target area, risk organs, and healthy tissues and the fusion of the time optimization objectives. F(X)=∑ r∈R δ r f r (X) (2) Where X represents the multidimensional radiotherapy path point sequence, R represents the real number set, and f r (X) is the rth objective function, δ r is the weight of the rth objective function.

7. The multi-dimensional radiotherapy path planning method of the radiotherapy robot based on the RRT algorithm according to claim 1, characterized in that: In step 4, a radiotherapy path point search strategy is formulated based on the RRT algorithm, and radiotherapy path points are selected from the sampling points according to the objective function.

8. The multi-dimensional radiotherapy path planning method of the radiotherapy robot based on the RRT algorithm according to claim 7, characterized in that: The step 4 comprises: starting from an empty beam set, calculating the objective function value after the sampling point is combined with the path point sequence in the beam set, iteratively adding radiotherapy path points to the beam set one by one, until the set upper limit of the number of radiotherapy points is reached, and the reduction of the objective function by the added radiotherapy points is less than a given threshold, then stopping the optimization.

9. The multi-dimensional radiotherapy path planning method of the radiotherapy robot based on the RRT algorithm according to claim 7, characterized in that: The step 4 comprises: (1): Randomly generate a starting radiotherapy path point and update it as the current radiotherapy path point; (2): With the current radiotherapy path point as the center, randomly sample at least one path point in the interval within a certain threshold range of the direction angle corresponding to the current radiotherapy path point, and remove the sampling points that do not meet the constraints; (3): Calculate the overall objective function value of the remaining sampling points combined with the previous path points, select the sampling point with the minimum objective function value and add it to the radiotherapy path point sequence, and update the current path point at the same time; (4): Determine whether the number of radiotherapy path points is greater than the set threshold, and whether the reduction of the objective function caused by the added radiotherapy path points is less than the given threshold. If so, stop the optimization; otherwise, continue searching for additional path points; (5): Calculate the dose weights {ω1,ω2,...,ω M }, and assign radiation dose to each path point according to the dose weight, as shown in the following formula: where γ i is the radiation dose at the ith path point, D pre is the prescribed dose, ω i is the dose weight calculated by flux optimization for the i-th path point, and M is the number of path points; (6): Generate a sequence of radiotherapy path points.

Images