A method and system for simulating the dose of heavy ion fixed head radiotherapy

Abstract

The application belongs to the technical field of radiotherapy dose simulation measurement, and provides a heavy ion fixed head radiotherapy dose simulation measurement method and system. The method first performs limited point and / or slice dose measurement in the solid water phantom according to the layered and key area encryption sampling principle to obtain sparse dose data and sampling positions; then inputs the sampling positions into the coding-decoding structure and the mask perception neural network with attention mechanism in the decoding stage, and outputs the complete three-dimensional dose distribution; finally, the simulation results containing the dose distribution and the uncertainty distribution are output. The application improves the dose reconstruction accuracy and the clinical quality control efficiency, and is suitable for heavy ion radiotherapy dose verification.

Description

TECHNICAL FIELD

[0001] The present application relates to the technical field of radiotherapy dose simulation measurement, in particular to a heavy ion fixed head radiotherapy dose simulation measurement method and system. BACKGROUND

[0002] Heavy ion radiotherapy has important advantages in tumor treatment due to its high biological effect ratio and clear Bragg peak characteristics. In order to ensure clinical safety and treatment accuracy, the beam dose distribution must be simulated and verified before treatment. The existing method generally uses solid water phantom as tissue equivalent, and obtains the dose distribution in the phantom through the detector. However, under the condition of fixed head, the beam direction is limited, and the small field dose changes steeply, so it is difficult to completely restore the three-dimensional dose distribution by relying on limited sparse measurement points, especially in the Bragg peak edge and high gradient area, and the existing technology has the problem of insufficient accuracy. SUMMARY

[0003] To this end, the present application provides a heavy ion fixed head radiotherapy dose simulation measurement method, device, electronic equipment, medium and computer program product, to at least partially solve the above technical problems.

[0004] The present application provides a heavy ion fixed head radiotherapy dose simulation measurement method, comprising the following method steps:

[0005] Performing limited point and / or slice dose measurement in the solid water phantom, using the layered and key area encryption sampling principle to sample sparse dose data and sampling positions, wherein the key area includes the depth range of the Bragg peak peak value, the rising edge and falling edge of the Bragg peak, the field edge area and the small field lateral boundary;

[0006] Inputting the sparse dose data and the sampling position mask into the neural network with mask perception, and outputting the complete three-dimensional dose distribution through the neural network, wherein the neural network adopts an encoding-decoding structure, and the decoding stage is provided with an attention mechanism;

[0007] Outputting the simulation results, the simulation results including dose distribution and uncertainty distribution, and the uncertainty distribution being constructed based on uncertainty index, wherein the uncertainty evaluation index includes sampling data density and network reconstruction error.

[0008] Another aspect of the present application also provides a heavy ion fixed head radiotherapy dose simulation measurement device, comprising:

[0009] The measurement data acquisition module is used for performing limited point and / or slice dose measurement in the solid water phantom, and sparse dose data and sampling positions are obtained by using the layered and key area encryption sampling principle, wherein the key area includes a depth range where a Bragg peak peak value is located, a rising edge and a falling edge of the Bragg peak, a field edge area and a small field lateral boundary.

[0010] The three-dimensional dose distribution acquisition module is used for inputting the sparse dose data and the sampling position mask into a mask-aware neural network, and outputting a complete three-dimensional dose distribution through the neural network, wherein the neural network adopts an encoding-decoding structure, and an attention mechanism is arranged in the decoding stage.

[0011] The simulation result output module is used for outputting the simulation result including a dose distribution and an uncertainty distribution, and the uncertainty distribution is constructed based on an uncertainty index, wherein the uncertainty evaluation index includes a sampling data density and a network reconstruction error.

[0012] In another aspect of the present application, a computer readable storage medium is provided, which stores computer program instructions, and the computer program instructions can be executed by a processor to implement the heavy ion fixed head radiotherapy dose simulation measurement method as described above.

[0013] In another aspect of the present application, a computer program product is provided, which includes a computer program, and the computer program is executed by a processor to implement the heavy ion fixed head radiotherapy dose simulation measurement method as described above.

[0014] The heavy ion fixed head radiotherapy dose simulation measurement method and system can accurately obtain sparse dose data in key areas by layered and key area encryption sampling, reduce the measurement workload while ensuring the integrity of the core area information, realize high-precision reconstruction of three-dimensional dose distribution by using an encoding-decoding neural network with an attention mechanism combined with joint loss function optimization, reduce the Bragg peak position error and improve the gamma index pass rate, output simulation results containing uncertainty distribution, intuitively indicate the reliability of dose estimation, and assist clinical quality control decision-making. The whole process is adapted to the clinical automation process, improves the radiotherapy dose verification precision and efficiency, and reduces the treatment risk. BRIEF DESCRIPTION OF DRAWINGS

[0015] In order to more clearly illustrate the technical solutions of the embodiments of the present application, the following will briefly introduce the drawings needed to be used in the embodiments. It should be understood that the following drawings only show some embodiments of the present application, and therefore should not be regarded as a limitation on the scope. For those skilled in the art, other related drawings can also be obtained without creative labor on the basis of these drawings.

[0016] Other features, objects, and advantages of the application will become more apparent from the following detailed description of non-limiting embodiments thereof, when read in connection with the following accompanying drawings:

[0017] Figure 1 A heavy ion fixed head radiotherapy dose simulation measurement method provided by an embodiment of the application.

[0018] Figure 2 A position mask construction method provided by an embodiment of the application.

[0019] Figure 3 A neural network architecture provided by an embodiment of the application.

[0020] Figure 4 A heavy ion fixed head radiotherapy dose simulation measurement system provided by an embodiment of the application.

[0021] Figure 5 A device provided by an embodiment of the application. DETAILED DESCRIPTION

[0022] To make the objectives, technical solutions, and advantages of the embodiments of the application clearer, the following will be combined with the accompanying drawings for the embodiments of the application to clearly and completely describe the technical solutions in the embodiments of the application. Obviously, the described embodiments are some but not all of the embodiments of the application. Based on the embodiments in the application, all other embodiments obtained by a person of ordinary skill in the art without creative work fall within the protection scope of the application.

[0023] To solve the above problems, the application provides a heavy ion fixed head radiotherapy dose simulation measurement method, and the technical solutions of the application will be described in detail in combination with various embodiments.

[0024] As shown in Figure 1 The application discloses a heavy ion fixed head radiotherapy dose simulation measurement method 100, comprising the following method steps:

[0025] S101, performing dose measurement of limited points and / or slice layers in a solid water phantom, and obtaining sparse dose data and sampling positions by adopting a layered and key area encryption sampling principle, wherein the key area includes a depth range where a Bragg peak peak value is located, a Bragg peak rising edge and a falling edge, an irradiation field edge angle area, and a small field lateral boundary;

[0026] S102, inputting the sparse dose data and the sampling position mask into a neural network with mask perception, and outputting a complete three-dimensional dose distribution through the neural network, wherein the neural network adopts an encoding-decoding structure, and an attention mechanism is arranged in a decoding stage.

[0027] S103, outputting the simulation result, the simulation result including a dose distribution and an uncertainty distribution, the uncertainty distribution being constructed based on an uncertainty index, wherein the uncertainty index includes a sampling data density and a network reconstruction error.

[0028] In one embodiment, for step S101, a solid water phantom is selected, which meets the quality control standards of radiotherapy, such as a material with a density close to human soft tissue and a ray attenuation coefficient equivalent to human tissue. The size needs to cover the depth and lateral range that may be involved in clinical treatment of the heavy ion beam (for example, the length is not less than the maximum Bragg peak depth of the beam + 5 cm, and the lateral cross section is not less than the maximum irradiation field diameter of the beam + 3 cm), to avoid dose measurement deviation caused by beam edge exceeding the range of the water phantom.

[0029] The solid water phantom is fixed on the positioning guide rail of the radiotherapy treatment bed, and the position of the water phantom is calibrated through a laser positioning system, to ensure that the central axis of the water phantom coincides with the central axis of the heavy ion beam, and the "depth direction" (beam incidence direction) of the water phantom is consistent with the lifting direction of the treatment bed, to ensure that the depth coordinate of the subsequent measurement corresponds to the "tumor depth" concept in clinical treatment.

[0030] Point dose detectors, such as semiconductor detectors, thermoluminescent dosimeters (TLDs), or ionization chambers, are set to collect "finite point" dose data. Surface dose measurement devices, such as dose films or two-dimensional array detectors, are set to collect "slice" dose data.

[0031] Optionally, a layered and key area encryption sampling principle is adopted (based on the dose distribution physical properties of the heavy ion beam, key areas are measured with smaller sampling intervals and more sampling points), the measurement layers are divided according to the depth direction (Z-axis) of the beam, and the sampling points are divided according to the lateral direction (X-Y plane), which specifically includes, for example:

[0032] Depth direction layering: starting from the beam incidence surface (the surface of the solid water phantom, Z=0), measurement depths are set at regular intervals and in the Bragg peak area. For example, in the plateau region (a region with gentle dose change) before the Bragg peak, one measurement depth is set every 5 cm; in the depth range of the Bragg peak (pre-calculated according to the energy of the heavy ion, for example, the Bragg peak depth of carbon ions with an energy of 150 MeV / u is about 10 cm), one measurement depth is set every 1 mm; in the tail region (a region with rapid dose decline) after the Bragg peak, one measurement depth is set every 2 mm, to ensure complete coverage of the rising edge, peak point, and falling edge of the Bragg peak.

[0033] Lateral plane sampling point arrangement: Within the XY plane at each measurement depth, sampling points are set according to the gradient region at the center and edge. For a circular irradiation field, with the beam center as the origin, sampling points are set in the radial direction at "0cm (center), 0.5cm, 1cm, 2cm, irradiation field boundary (e.g., 3cm), and 1cm outside the boundary". For a rectangular irradiation field, sampling points are set in the major and minor axes at the above intervals, and additional sampling points are added in the "corner region" of the irradiation field (where the dose gradient is steeper) (e.g., one point every 0.3cm) to avoid loss of dose information at the corners of the small field.

[0034] Optionally, select 3-5 key depth planes as cutting layers. For example, the following areas are preferentially covered:

[0035] Typical depth of the plateau region after beam incident (e.g., Z=2cm, dose stability region).

[0036] The depth of the Bragg peak (e.g., Z=10cm, the region of highest dose).

[0037] The middle depth of the Bragg peak's falling edge (e.g., Z=10.5cm, the region of maximum dose gradient).

[0038] The depths of the upper and lower margins of the tumor target area (simulating the key layers where the tumor is located in clinical treatment).

[0039] For example, measurements are performed using corresponding equipment for each slice, taking a two-dimensional array detector as an example.

[0040] Align the detection surface of the two-dimensional array detector with the slice depth plane, and ensure that the center of the detector coincides with the beam center, so that the effective detection area of ​​the detector completely covers the irradiation field and edge areas (e.g., the transverse dimension of the detector is 2 cm larger than the irradiation field). After the beam irradiation is started, the detector collects the dose value of each detection unit in real time, directly outputs the two-dimensional dose data of the slice, and records the three-dimensional coordinates of each unit of the detector (X and Y coordinates are the position of the unit in the array, and Z coordinate is the slice depth) as the sampling position information of the slice.

[0041] Optionally, after the measurement is completed, the measurement data is standardized and organized. For data with limited points, for example, it is grouped by depth-lateral coordinates, and the dose value (unit: Gy) corresponding to the three-dimensional coordinates (X, Z, Y) of each sampling point is recorded, and outliers (such as data where the deviation of a single measurement value from the average value exceeds 5%) are removed.

[0042] For slice data, for example, the two-dimensional dose distribution is converted into discrete point data. For example, the XY plane of the slice is divided into grids at fixed intervals (e.g., 0.2 cm), the dose value of each grid node is extracted, and the discrete point dose set of the slice is formed. This set is then merged with the finite point data to form sparse dose data (which is called "sparse" because the total number of sampling points is much less than the total number of voxels in the three-dimensional dose field).

[0043] For sampling location recording, for example, a three-dimensional rectangular coordinate system is established with the center of the incident surface of the solid water model's beam as the origin (X=0, Y=0, Z=0), as exemplarily:

[0044] X-axis: Horizontal direction (consistent with the left-right direction of the treatment bed);

[0045] Y-axis: Vertical direction (consistent with the up-and-down direction of the treatment bed);

[0046] Z-axis: Depth direction (consistent with the incident direction of the beam, the Z value increases with increasing depth).

[0047] Record the (X,Y,Z) coordinates corresponding to each sparse dose data to form a list of sampling locations, and label the measurement type of each sampling point (e.g., point detector - plateau area, film slice - Bragg peak) to facilitate the subsequent neural network to distinguish sampling data from different sources and optimize reconstruction accuracy.

[0048] In one embodiment, for step S102, before inputting the model, the sparse dose data is preprocessed, for example, by removing outlier data that exceeds a threshold and normalizing the dose values. For example, the sparse dose data is mapped to the [0,1] interval using a linear normalization method with reference to the dose range of the high-resolution dose field generated by Monte Carlo calculation.

[0049] Simultaneously, the sparse dose data is converted from a discrete list format of coordinates-dose to a three-dimensional mesh sparse matrix that matches the input dimension of the neural network. For example, a three-dimensional mesh is constructed based on the physical dimensions of the solid water phantom (e.g., the mesh voxel size is set to 0.1mm × 0.1mm × 0.1mm, consistent with the accuracy requirements for clinical dose validation). Standardized dose values ​​are only entered at voxel locations with sampled data; voxel locations without sampled data are temporarily left empty, forming a preliminary sparse three-dimensional matrix, preparing for subsequent combination with a sampling location mask.

[0050] In this embodiment, the sampling location mask explicitly indicates the true sampling location of sparse dose data to the network, helping the network to prioritize areas with measurement data during the decoding stage and reasonably infer areas without measurement data based on information from these areas. For example, such as... Figure 2 As shown, it can be constructed through the following steps:

[0051] S201, 3D mesh coordinate system alignment, ensuring that the 3D mesh of the mask is completely consistent with the 3D mesh of the sparse dose data. Specifically, with the center of the solid water phantom beam incident surface as the origin (X=0, Y=0, Z=0), the X-axis corresponds to the horizontal direction of the water phantom, the Y-axis corresponds to the vertical direction, and the Z-axis corresponds to the beam depth direction. The size and number of mesh voxels are exactly the same as the 3D mesh of the sparse dose data, ensuring the correspondence between dose data and mask positions.

[0052] S202, Mask value definition and padding: The mask value is defined using a binary marking method. Specifically, for voxels with sparse dose data in the 3D mesh (i.e., coordinate points of finite point measurements and discretized mesh nodes of slice measurements), the mask value is set to "1", indicating that there is real measurement data at this location, and the network can directly use this data for feature learning; for voxels without sparse dose data, the mask value is set to "0", indicating that the dose inference at this location needs to be performed by the network based on the surrounding areas with measurement data.

[0053] For continuous planar data from slice measurements, the two-dimensional dose distribution of the slice needs to be discretized into voxel data in a three-dimensional grid (e.g., if the slice depth is Z=10cm, the mask value corresponding to all voxels with dose data (X,Y,10cm) in the plane is set to "1") to ensure that the continuous sampling area of ​​the slice is completely marked in the mask and to avoid the network misjudging the area as a region without data.

[0054] S203, Mask edge smoothing: For example, the boundary between the "1" value region (with data) and the "0" value region (without data) in the mask is smoothed using Gaussian blur (e.g., setting the Gaussian kernel size to 3×3×3). This avoids hard transitions at the boundary, causing the network to focus on boundary noise rather than the true dose distribution gradient during training. The smoothed mask can still clearly distinguish between data regions and non-data regions, while better reflecting the continuous characteristics of gradient changes in the heavy ion dose distribution.

[0055] In one embodiment, the input layer design enables the collaborative input of sparse dose data and sampling location masks, ensuring that the network can simultaneously read dose value information and sampling location information.

[0056] Specifically, a dual-input channel structure is adopted. For example, the normalized sparse dose data three-dimensional matrix is ​​used as channel 1 (dose channel), and the constructed sampling position mask three-dimensional matrix is ​​used as channel 2 (position channel). The two types of data are completely aligned in spatial dimension (the number of voxels and the coordinate position are consistent), and together they constitute the input tensor of the network.

[0057] The dimensions of the input tensor can be represented as (C×D×H×W), where (C=2) (corresponding to two channels), D (depth direction), H (lateral direction of Y-axis), and W (lateral direction of X-axis) are the number of voxels in the 3D mesh, respectively.

[0058] The dual-channel design enables the network to extract both dose magnitude and data reliability features during the encoding stage. For example, for regions with high dose values ​​in channel 1 (which may be Bragg peak regions), the network can confirm the reliability of the data in that region by combining the information that the mask value of that region is "1" in channel 2, and prioritize the construction of high-resolution features based on that data. For regions with null values ​​in channel 1 and mask values ​​of "0" in channel 2, the network explicitly needs to make inferences based on the surrounding reliable data.

[0059] In one embodiment, such as Figure 3 As shown, the neural network adopts an encoder-decoder structure, including an input layer, an encoder layer, a skip connection layer, a decoder layer, and an output layer. The decoder layer has an attention mechanism to enhance the reconstruction accuracy in high-dose gradient regions.

[0060] The encoding stage involves progressively compressing the spatial dimension and extracting more abstract multi-scale dose features (such as global beam trends and local high-dose region features) from the dual-channel input of sparse dose data and sampling location masks. This provides feature support for the reconstruction in the decoding stage. For example, it consists of four convolutional layers.

[0061] Specifically, the functions and parameters of each level are designed as follows:

[0062] The input layer and the first convolutional layer (shallow feature extraction) receive, for example, a preprocessed dual-channel input tensor, and use a (3×3×3) three-dimensional convolutional kernel (adapting to the spatial characteristics of the three-dimensional dose field). The number of convolutional kernels is set to 64, the stride is 1, and the padding is set to Same (to ensure that the size of the feature map after convolution is consistent with the input). The ReLU activation function is used to perform shallow feature extraction on the input data.

[0063] The 2nd to 4th convolutional layers (deep feature compression and abstraction) use a combination of convolution and pooling to gradually compress spatial dimensions and improve feature abstraction.

[0064] Specifically, for convolution operations, a (3×3×3) three-dimensional convolution kernel is used, with the number of kernels increasing sequentially (128 in the second layer, 256 in the third layer, and 512 in the fourth layer). The stride is 1 and the padding is Same. By increasing the number of kernels, the feature representation capability is improved, and more complex dose distribution patterns are captured (such as the correlation between the depth position of the Bragg peak and the lateral diffusion range).

[0065] For pooling operations, each convolutional layer is followed by a (2×2×2) three-dimensional max pooling (stride 2), which compresses the spatial dimension to half of its original size (e.g., (D×H×W) becomes (D / 2×H / 2×W / 2)). At the same time, it retains the largest feature value (i.e. the most significant dose feature, such as high dose peak and obvious gradient change) within each pooling window, reducing the amount of computation while strengthening the core features.

[0066] The final encoding stage outputs a 512-channel deep feature map (dimension (512×D / 8×H / 8×W / 8), which is compressed to 1 / 8 of the input due to 3 pooling operations).

[0067] To address the loss of detail caused by spatial dimension compression during the encoding stage, the network incorporates symmetrical skip connections between the encoding and decoding layers. This allows feature maps from each layer in the encoding stage to be directly passed to the corresponding layers in the decoding stage, supplementing the detail information. The specific implementation is as follows:

[0068] Skip connections employ a channel concatenation method. For example, the feature maps output from layers 1-3 of the encoding stage (the original convolutional feature maps without pooling) are concatenated with the feature maps from the corresponding layers (layers 3-1) of the decoding stage after transposed convolution in the channel dimension (instead of simple addition), ensuring that the detailed features of the encoding stage are not covered by the restoration operation in the decoding stage.

[0069] Among them, the first layer of encoding (shallow, high resolution) corresponds to the third layer of decoding (the spatial dimension has been partially restored during the restoration process, and fine-grained details are supplemented); the third layer of encoding (deep, low resolution) corresponds to the first layer of decoding (the spatial dimension is close to the final output, and global feature associations are supplemented), forming a complementary relationship between encoding details and decoding restoration.

[0070] The decoding stage involves taking the deep abstract features output from the encoding stage, gradually restoring the spatial dimension through transposed convolutions, and combining this with the detailed features supplemented by skip connections to reconstruct a complete 3D dose distribution with the same size as the input grid. This stage consists of four transposed convolutional layers, each with the following functions:

[0071] The first to third transposed convolutional layers (spatial dimension restoration and feature fusion) exemplarily employ a (2×2×2) three-dimensional transposed convolutional kernel (stride 2). Each transposed convolution expands the spatial dimension to twice its original size (e.g., (D / 8×H / 8×W / 8) is gradually restored to (D / 4×H / 4×W / 4), (D / 2×H / 2×W / 2), and (D×H×W)). The number of convolutional kernels decreases sequentially (256 in the first layer, 128 in the second layer, and 64 in the third layer), symmetrical to the number of convolutional kernels in the encoding stage, ensuring that the feature dimension gradually adapts to the final output.

[0072] After each transposed convolutional layer, the feature maps passed through the corresponding layers in the encoding stage are concatenated by channels. For example, the 256-channel feature map of the decoding layer 1 and the 256-channel feature map of the encoding layer 3 are concatenated to form 512 channels. Then, feature fusion is performed through (1×1×1) convolutional kernels (the number of which is the same as the output channels of the transposed convolutional layer) to eliminate feature redundancy after channel concatenation and retain effective details.

[0073] After the spatial dimension is restored to be consistent with the input ((D×H×W)), the 4th transposed convolutional layer ((3×3×3) convolutional kernel, 1 output channel, stride 1, padding is Same) converts the 64-channel feature map into a single-channel initial three-dimensional dose distribution matrix. This matrix has initially contained the complete dose distribution outline, but the detail accuracy of the high dose gradient region still needs to be further enhanced through the attention mechanism.

[0074] In this embodiment, the attention mechanism is embedded between the third and fourth transposed convolutional layers in the decoding stage: at this point, the spatial dimension has been fully restored ((D×H×W)), the feature map has been fused with the detailed information from the encoding stage, and the final dose distribution has not yet been output. Applying attention weights here can directly enhance the accuracy of key regions in the dose distribution that is about to be generated, thus avoiding the inability to correct it later.

[0075] Specifically, the 64-channel feature map output from the third transposed convolutional layer is received and decoded. First, a global average pooling operation is performed (averaging all voxels in each channel) to obtain a (64×1×1×1) channel feature vector. This vector can characterize the contribution of each channel to the key features of the dose distribution (e.g., a channel is better at capturing high-dose regions, and its average pooling value is higher).

[0076] The 64-dimensional channel feature vector obtained by global average pooling is input into a two-layer fully connected network (32 neurons in the first layer and 64 neurons in the second layer), and the 64 channel weight values ​​(range 0-1) are output through the Sigmoid activation function:

[0077] The channel weights are multiplied by the original 64-channel feature map in a channel-wise manner (the feature map of each channel is multiplied by the corresponding channel weight) to obtain the weighted feature map.

[0078] For the weighted feature map, the spatial gradient value of each voxel is calculated. Specifically, the dose change rate of each voxel in the depth direction (Z-axis) and the lateral direction (X-axis, Y-axis) is calculated using the (3×3×3) Sobel operator. The square root of the sum of the squares of the change rates in the three directions is taken to obtain the spatial gradient value. The higher the value, the more the voxel belongs to a high dose gradient region (such as the edge of the Bragg peak). The spatial gradient value is normalized to the range of 0-1 using the Sigmoid activation function to obtain the spatial weight value of each voxel, thereby realizing the localization of key regions at the spatial level.

[0079] The spatial weights are element-wise multiplied with the weighted feature map to obtain the final attention-weighted feature map. It can be understood that the voxels in high-dose gradient regions (such as the rising and falling edges of the Bragg peak) have high spatial weights, so their feature values ​​are amplified. The network will prioritize learning the details of this region (such as subtle changes in dose value and clear outlines of edges) during subsequent transposed convolutions. On the other hand, the voxels in regions with flat doses have low spatial weights, so their feature values ​​are weakened, thus avoiding the network from over-focusing on non-critical regions and wasting resources.

[0080] During network training, the weight parameters of the attention mechanism (weights of fully connected layers, coefficients of Sobel operators) are optimized together with the parameters of the encoder-decoder layers through backpropagation using a joint loss function.

[0081] In some embodiments, during the training of the neural network, a high-resolution dose distribution generated by Monte Carlo calculation is used as the ground truth, sparse sampled data is used as input, and a joint loss function of voxel error, Bragg peak position error and energy conservation constraint is used for optimization.

[0082] The Monte-Carlo method has been selected as the gold standard source of truth for high-resolution dose distribution because it can accurately simulate the interaction between heavy ions and matter (such as ionization and nuclear reactions). The generation process is consistent with the parameters of clinical radiotherapy scenarios, including beam parameters, tissue equivalent parameters, and spatial resolution.

[0083] The energy deposition process of heavy ion beams in a solid water model is simulated using Monte Carlo programs (such as FLUKA and GEANT4), and a three-dimensional dose distribution matrix (i.e., high-resolution true value) is output.

[0084] The generated true values ​​are physically validated: the shape of the Bragg peak is checked (whether it conforms to the typical characteristics of plateau-rising edge-peak-falling edge) and energy conservation (e.g., the deviation between the total dose in the water model and the total incident energy of the beam is ≤0.5%). Abnormal true values ​​that do not conform to physical laws are eliminated to ensure the reliability of each true value sample.

[0085] Sparse sampling data is extracted from the Monte Carlo ground truth and used as the training input for the network. The construction logic is completely consistent with the actual clinical measurement process. The sampling is carried out according to finite point sampling and slice sampling. The dose value is extracted from the ground truth for each sampling point to form finite point data of coordinate-dose. The process of merging finite point and slice data will not be described in detail here.

[0086] For each sparse sampled data, a sampling location mask is generated synchronously (consistent with the mask format in the inference stage): in the same 3D grid as the ground truth, voxels with sparse sampled data are marked as 1, and voxels without data are marked as 0, forming a sample triplet of sparse input-mask-ground truth, ensuring that the format of the training input and the inference input are completely matched.

[0087] The dose data (sparse input and true value) of all samples are linearly normalized to [0,1] (based on the maximum and minimum dose values ​​of each sample's true value) to eliminate the difference in absolute dose values ​​caused by different beam energies;

[0088] Outlier filtering is performed on sparse input data (data that exceeds ±10% of the true dose range due to sampling error is removed) to ensure the reliability of input data and avoid noise interference during the training process.

[0089] In one embodiment, the loss function is optimized using a joint loss function of voxel error, Bragg peak position error, and energy conservation constraint.

[0090] Specifically, regarding voxel error loss It is used to measure the dose deviation between the three-dimensional dose distribution output by the network and the Monte Carlo true value at each voxel, ensuring that the overall distribution trend is consistent; the mean square error (MSE) is used to calculate the dose distribution output by the network and the true value at all voxels;

[0091] Bragg peak position error loss Based on the requirements for Bragg peak position error, the depth position accuracy of Bragg peak is optimized to avoid insufficient dose to the clinical treatment target area or damage to normal tissue due to position deviation.

[0092] For example, the calculation method includes:

[0093] Bragg peak location: for network output With truth value The maximum dose value of each XY section in the depth direction (Z-axis) is calculated, and the average depth value corresponding to the maximum dose of all sections is taken to obtain the result. (Predicted Bragg peak depth) and (True Bragg peak depth);

[0094] Position error calculation: Calculate the absolute position error to obtain... ;

[0095] Loss value definition: A piecewise function is used to enhance the penalty for exceeding tolerances, and the formula is as follows:

[0096]

[0097] When the error is ≤1mm, the loss is 0; when the error exceeds the limit, the network is forced to correct the position deviation by applying a penalty based on the square of the error.

[0098] Energy conservation constraints prevent the energy of the heavy ion beam from being completely deposited within the solid water phantom. The total dose output by the network must be consistent with the true total dose (energy conservation) to avoid deviations in the total dose due to overfitting of local data by the network, which could affect the accuracy of clinical dosing.

[0099] For example, the calculation method includes:

[0100] To each and Summing all voxel doses yields the total dose. and ;

[0101] Calculate relative deviation ;

[0102] The loss value is also defined using a piecewise function, and the formula is as follows:

[0103]

[0104] When the relative deviation is ≤2% (physically reasonable range), the loss is 0; when it exceeds the limit, the deviation value is penalized to ensure that the network output conforms to the law of energy conservation.

[0105] Based on clinical priorities (e.g., Bragg peak position accuracy > overall dose accuracy > energy conservation accuracy), different weights are assigned to the three loss terms. For example, the formula is as follows:

[0106]

[0107] As an example, the position of the Bragg peak directly determines whether the tumor target area is accurately irradiated, so it is given the highest weight (0.5); voxel error affects the reliability of the overall dose distribution, so it is given the second highest weight (0.3); energy conservation is a physical constraint, so it is given a weight (0.2). The sum of the weights of the three is 1, which ensures that the loss value is quantifiable and stable.

[0108] Based on the prepared dataset and joint loss function, the network parameters are iteratively adjusted through backpropagation, enabling the network to gradually learn the mapping relationship from sparse input to high-resolution ground truth. For example, the specific process includes:

[0109] Network parameter initialization and hyperparameter settings, specifically including:

[0110] He normal initialization is used for the convolution kernels and transposed convolution kernels of the encoder-decoder structure (to adapt to the ReLU activation function and avoid gradient vanishing); Xavier uniform initialization is used for the parameters of the fully connected layer of the attention mechanism to ensure that the initial parameter distribution is reasonable and to accelerate training convergence.

[0111] Choose the Adam optimizer (fast convergence speed, low sensitivity to learning rate), set the initial learning rate to 0.001, and dynamically decrease it with each training round (decreasing to 0.9 times the original value every 10 rounds) to avoid overfitting in the later stages;

[0112] Set the batch size to 8-16 (to balance computational efficiency and memory usage, and ensure that each batch contains multiple sample types).

[0113] The initial training rounds are set to 50 rounds. The validation set loss is used as the monitoring metric. If the validation set loss does not decrease for 5 consecutive rounds (overfitting occurs), training is stopped early.

[0114] The iterative training process specifically includes,

[0115] Forward propagation: The prediction results are generated by randomly sampling a batch of samples (including sparse input and mask) from the training set and inputting them into a 3D super-resolution neural network; the network extracts features through an encoder-decoder structure and strengthens high-gradient regions through an attention mechanism, and finally outputs the predicted 3D dose distribution. .

[0116] Loss Calculation and Backpropagation: Adjusting Network Parameters

[0117] based on Monte Carlo true value of this batch of samples Calculate joint loss ;

[0118] Through backpropagation, the parameters of each layer (convolutional kernel, attention weights, etc.) are calculated progressively from the output layer to the input layer. The gradient;

[0119] The optimizer updates the network parameters according to the gradient direction and the set learning rate, reducing the joint loss and making the network prediction results closer to the true value.

[0120] Validation and adjustment to avoid overfitting, specifically including:

[0121] After each training round, the network performance is evaluated using a validation set: The validation set is calculated... Bragg peak position error, gamma index (2% / 2mm) pass rate;

[0122] If the validation set If the loss continues to rise (overfitting), an early stopping strategy is adopted to stop training and load the network parameters (optimal model) with the minimum loss on the validation set.

[0123] In one embodiment, for step S103, the three-dimensional dose distribution obtained after neural network optimization undergoes format conversion, physical quantity restoration, and clinical annotation to form an output file that meets radiotherapy quality control standards.

[0124] Specifically, during network training and optimization, the dose data undergoes standardization processing within the [0,1] interval (referencing the dose range of the Monte Carlo true value). In the output stage, the standardized dose value is first reverse-mapped to the real physical dose (unit: Gy). The restoration method corresponds to the normalization method. For example, D_real = D_norm × (D_max - D_min) + D_min), where D_norm is the standardized dose value output by the network, and D_max and D_min are the maximum and minimum real dose values ​​of the corresponding Monte Carlo high-resolution dose field, respectively.

[0125] The restored three-dimensional dose distribution matrix (voxel size 0.1mm×0.1mm×0.1mm) is converted into the DICOM-RTDose format (radiotherapy dose-specific format) commonly used in the radiotherapy field. Further details are omitted here.

[0126] In the output three-dimensional dose distribution, dose parameters for key clinically relevant areas are automatically labeled, facilitating rapid assessment of efficacy and risk. For example,

[0127] Bragg peak region: Marks the center depth, peak dose, and peak width of the Bragg peak;

[0128] Tumor target simulation area: If the input data includes target location information (such as the pre-set target equivalent layer in the solid water model during measurement), then the average dose, maximum dose and minimum dose within the target area are marked.

[0129] Normal tissue simulation area: The maximum dose in the non-target area along the injection beam path, providing a reference for clinical assessment of the risk of damage to normal tissue.

[0130] In this embodiment, the uncertainty evaluation index is used to indicate the reliability of dose estimation at different voxel locations, and its core is based on the calculation of sampling data density and network reconstruction error.

[0131] Specifically, the magnitude of uncertainty is positively correlated with two factors: one is the sampling data density around the voxel (the sparser the sampling points, the higher the uncertainty); the other is the prediction error when the network reconstructs the voxel (the greater the deviation from the Monte-Carlo true value, the higher the uncertainty). Before calculation, two basic matrices need to be constructed:

[0132] Sampling density matrix: Count the number of sampling points within the 3×3×3 voxel range around each voxel (the number of voxels with a mask value of "1"). The more sampling points, the higher the density value (range 0 - 27);

[0133] Reconstruction error matrix: Calculate the absolute error between the restored true dose and the Monte-Carlo true value for each voxel. The greater the error, the lower the network reconstruction accuracy.

[0134] Normalize the sampling density matrix and the reconstruction error matrix to [0,1] respectively, and then fuse them according to the weights to obtain the uncertainty value U of each voxel. Exemplarily, the calculation formula can be simplified as:

[0135]

[0136] where is the normalized sampling density value (the higher the sampling density, the smaller it is, and the lower the contribution to uncertainty);

[0137] is the normalized reconstruction error value (the greater the error, the higher the contribution to uncertainty);

[0138] The weight setting (sampling density 0.6, reconstruction error 0.4) is based on clinical priorities: the reliability of sampling data has a greater impact on uncertainty, so a higher weight is given. It can be understood that the specific numerical value of the weight can be set according to actual applications, and the present invention does not limit it.

[0139] Optionally, after calculation, the uncertainty value is divided into three risk levels,

[0140] Low uncertainty (U≤0.3): The dose estimate is reliable and can be used for treatment plan verification;

[0141] Medium uncertainty (0.3 < U≤0.6): There is a certain deviation in the dose estimate, and other verification methods (such as additional point dose measurement) need to be combined to confirm; <OOO0321>

[0142] High uncertainty (U>0.6): The dose estimate is unreliable, and it is necessary to trace back the measurement or optimize the network reconstruction process to check for problems such as insufficient sampling or model deviation.

[0143] Figure 4A heavy ion fixation radiotherapy dose simulation measurement system 400 is shown. Embodiments of this system are... Figure 1 Corresponding to the method embodiments shown, this system can be specifically applied to various electronic devices, including:

[0144] The measurement data acquisition module 401 is used to perform dose measurements at limited points and / or in sections within a solid water model. It uses a layered and key area densified sampling principle to obtain sparse dose data and sampling locations. The key areas include the depth range of the Bragg peak, the rising and falling edges of the Bragg peak, the corner areas of the irradiation field, and the lateral boundaries of the small field.

[0145] The three-dimensional dose distribution acquisition module 402 is used to input the sparse dose data and the sampling position mask into a mask-aware neural network, and output a complete three-dimensional dose distribution through the neural network. The neural network adopts an encoder-decoder structure, and the decoding stage is equipped with an attention mechanism.

[0146] The simulation result output module 403 is used to output the simulation results, which include dose distribution and uncertainty distribution. The uncertainty distribution is constructed based on uncertainty indicators, including sampling data density and network reconstruction error.

[0147] Based on the same inventive concept, this application also provides an electronic device. The method corresponding to the electronic device can be the method in the foregoing embodiments, and its problem-solving principle is similar to that method. The electronic device provided in this application includes: at least one processor; and a memory communicatively connected to the at least one processor; wherein the memory stores instructions executable by the at least one processor, and the instructions are executed by the at least one processor to enable the at least one processor to execute the methods and / or technical solutions of the foregoing embodiments of this application.

[0148] The electronic device can be a user device, or a device formed by integrating user devices and network devices through a network, or it can be an application running on the aforementioned devices. The user device includes, but is not limited to, various terminal devices such as computers, mobile phones, tablets, smartwatches, and smart bands. The network device includes, but is not limited to, network hosts, single network servers, multiple network server sets, or cloud computing-based computer sets, and can be used to implement some processing functions when setting an alarm clock. Here, the cloud consists of a large number of hosts or network servers based on cloud computing. Cloud computing is a type of distributed computing, consisting of a virtual computer composed of a group of loosely coupled computer sets.

[0149] Figure 5The diagram illustrates the structure of an apparatus suitable for implementing the methods and / or technical solutions in the embodiments of this application. The apparatus 500 includes a central processing unit (CPU) 501, which can perform various appropriate actions and processes based on a program stored in a read-only memory (ROM) 502 or a program loaded from a storage portion 508 into a random access memory (RAM) 503. The RAM 503 also stores various programs and data required for system operation. The CPU 501, ROM 502, and RAM 503 are interconnected via a bus 504. An input / output (I / O) interface 505 is also connected to the bus 504.

[0150] The following components are connected to I / O interface 505: input section 506 including keyboard, mouse, touch screen, microphone, infrared sensor, etc.; output section 507 including cathode ray tube (CRT), liquid crystal display (LCD), LED display, OLED display, etc., and speakers, etc.; storage section 508 including one or more computer-readable media such as hard disk, optical disk, magnetic disk, semiconductor memory, etc.; and communication section 509 including network interface card such as LAN (local area network) card, modem, etc. Communication section 509 performs communication processing via a network such as the Internet.

[0151] In particular, the methods and / or embodiments in this application can be implemented as computer software programs. For example, the embodiments disclosed in this application include a computer program product comprising a computer program carried on a computer-readable medium, the computer program containing program code for performing the methods shown in the flowchart. When the computer program is executed by the central processing unit (CPU) 501, it performs the functions defined in the methods of this application.

[0152] Another embodiment of this application provides a computer-readable storage medium having computer program instructions stored thereon, which can be executed by a processor to implement the methods and / or technical solutions of any one or more embodiments of this application described above.

[0153] Specifically, this embodiment may employ any combination of one or more computer-readable media. A computer-readable medium may be a computer-readable signal medium or a computer-readable storage medium. A computer-readable storage medium may be, for example, but not limited to, an electrical, magnetic, optical, electromagnetic, infrared, or semiconductor system, apparatus, or device, or any combination thereof. More specific examples of computer-readable storage media (a non-exhaustive list) include: an electrical connection having one or more wires, a portable computer disk, a hard disk, a random access memory, a read-only memory, an erasable programmable read-only memory, an optical fiber, a portable compact disk read-only memory, an optical storage device, a magnetic storage device, or any suitable combination thereof. In this document, a computer-readable storage medium may be any tangible medium containing or storing a program that can be used by or in connection with an instruction execution system, apparatus, or device.

[0154] Computer-readable signal media may include data signals propagated in baseband or as part of a carrier wave, carrying computer-readable program code. Such propagated data signals may take various forms, including but not limited to electromagnetic signals, optical signals, or any suitable combination thereof. Computer-readable signal media may also be any computer-readable medium other than computer-readable storage media, capable of sending, propagating, or transmitting programs for use by or in connection with an instruction execution system, apparatus, or device.

[0155] Program code contained on a computer-readable medium may be transmitted using any suitable medium, including but not limited to wireless, wire, optical fiber, RF, etc., or any suitable combination thereof.

[0156] Computer program code for performing the operations of this application can be written in one or more programming languages ​​or a combination thereof, including object-oriented programming languages ​​such as Java, Smalltalk, and C++, as well as conventional procedural programming languages. The program code can be executed entirely on the user's computer, partially on the user's computer, as a standalone software package, partially on the user's computer and partially on a remote computer, or entirely on a remote computer or server. In cases involving remote computers, the remote computer can be connected to the user's computer via any type of network, including a local area network (LAN) or a wide area network (WAN), or it can be connected to an external computer (e.g., via the Internet using an Internet service provider).

[0157] The flowcharts or block diagrams in the accompanying drawings illustrate the architecture, functionality, and operation of possible implementations of devices, methods, and computer program products according to various embodiments of this application. In this regard, each block in a flowchart or block diagram may represent a module, segment, or portion of code containing one or more executable instructions for implementing a specified logical function. It should also be noted that in some alternative implementations, the functions indicated in the blocks may occur in a different order than those indicated in the drawings. For example, two consecutively indicated blocks may actually be executed substantially in parallel, and they may sometimes be executed in reverse order, depending on the functions involved. It should also be noted that each block in the block diagrams and / or flowcharts, and combinations of blocks in the block diagrams and / or flowcharts, can be implemented using a dedicated hardware-specific system that performs the specified function or operation, or using a combination of dedicated hardware and computer instructions.

[0158] Furthermore, it is clear that the word "comprising" does not exclude other units or steps, and the singular does not exclude the plural. Multiple units or devices recited in a device claim may also be implemented by a single unit or device through software or hardware. Terms such as "first," "second," etc., are used to indicate names and do not indicate any specific order.

Claims

1. A method for simulating and measuring the dose of a heavy ion fixation radiotherapy device, characterized in that, Includes the following steps: Dose measurements were performed at limited points and / or in sections within a solid water phantom. Sparse dose data were obtained by using a layered and key area densified sampling principle, as well as the sampling locations. Key areas included the depth range of the Bragg peak, the rising and falling edges of the Bragg peak, the corner areas of the irradiation field, and the lateral boundaries of the small field. The sparse dose data and sampling position mask are input into a mask-aware neural network, which outputs a complete three-dimensional dose distribution. The neural network adopts an encoder-decoder structure, and an attention mechanism is provided in the decoding stage. During the training of the neural network, the high-resolution dose distribution generated by Monte Carlo calculation is used as the ground truth, and the sparse sampling data is used as the input. The network is optimized using a joint loss function of voxel error, Bragg peak position error and energy conservation constraint. The simulation results are output, including dose distribution and uncertainty distribution. The uncertainty distribution is constructed based on uncertainty indices, which include sampling data density and network reconstruction error. The sampling location mask construction includes: Align the 3D grid coordinate system to ensure that the 3D grid of the mask is completely consistent with the 3D grid of the sparse dose data; For voxels with sparse dose data in the 3D mesh, set the mask value to 1; for voxels without sparse dose data, set the mask value to 0. The mask edges are smoothed using Gaussian blur. The alignment of the three-dimensional grid coordinate system ensures that the three-dimensional grid of the mask is completely consistent with the three-dimensional grid of the sparse dose data, including: With the center of the solid water model beam incident surface as the origin (X=0, Y=0, Z=0), the X-axis corresponds to the horizontal direction of the water model, the Y-axis corresponds to the vertical direction, and the Z-axis corresponds to the beam depth direction. The size and number of grid voxels are exactly the same as the three-dimensional grid of the sparse dose data to ensure the correspondence between dose data and mask position.

2. The method for simulating and measuring the dose of a heavy ion fixation radiotherapy head according to claim 1, characterized in that, The neural network has an input layer based on a dual-input channel structure, with each channel corresponding to a three-dimensional matrix of sparse dose data and a three-dimensional matrix of sampling position mask, respectively. A symmetrical skip connection is set between its encoding and decoding layers, and the unpooled original feature maps output sequentially from multiple layers in the encoding stage are concatenated with the feature maps of the corresponding layers in the decoding stage after transposed convolution in the channel dimension.

3. A heavy ion fixation radiotherapy dose simulation and measurement system, characterized in that, include: The measurement data acquisition module is used to perform dose measurements at limited points and / or in sections within a solid water model. It uses a layered and key area densified sampling principle to obtain sparse dose data and sampling locations. The key areas include the depth range of the Bragg peak, the rising and falling edges of the Bragg peak, the corner areas of the irradiation field, and the lateral boundaries of the small field. The three-dimensional dose distribution acquisition module is used to input the sparse dose data and the sampling position mask into the mask-aware neural network, and output the complete three-dimensional dose distribution through the neural network. The neural network adopts an encoder-decoder structure, and the decoding stage is equipped with an attention mechanism. During the training process of the neural network, the high-resolution dose distribution generated by Monte Carlo calculation is used as the ground truth, the sparse sampling data is used as the input, and the joint loss function of voxel error, Bragg peak position error and energy conservation constraint is used for optimization. The simulation results output module is used to output simulation results including dose distribution and uncertainty distribution, wherein the uncertainty distribution is constructed based on uncertainty indices, and wherein the uncertainty indices include sampling data density and network reconstruction error. The three-dimensional dose distribution acquisition module is specifically used for aligning the three-dimensional grid coordinate system to ensure that the three-dimensional grid of the mask is completely consistent with the three-dimensional grid of the sparse dose data. For voxels with sparse dose data in the 3D mesh, set the mask value to 1; for voxels without sparse dose data, set the mask value to 0. The mask edges are smoothed using Gaussian blur. The three-dimensional dose distribution acquisition module is further configured to take the center of the solid water model beam incident surface as the origin (X=0, Y=0, Z=0), with the X-axis corresponding to the horizontal direction of the water model, the Y-axis corresponding to the vertical direction of the water model, and the Z-axis corresponding to the beam depth direction. The size and number of the grid voxels are exactly the same as the three-dimensional grid of the sparse dose data to ensure the correspondence between dose data and mask position.

4. The heavy ion fixation head radiotherapy dose simulation measurement system according to claim 3, characterized in that, It also includes, The neural network has an input layer based on a dual-input channel structure, with each channel corresponding to a three-dimensional matrix of sparse dose data and a three-dimensional matrix of sampling position mask, respectively. A symmetrical skip connection is set between its encoding and decoding layers, and the unpooled original feature maps output sequentially from multiple layers in the encoding stage are concatenated with the feature maps of the corresponding layers in the decoding stage after transposed convolution in the channel dimension.

5. An electronic device, characterized in that, include: At least one processor; and a memory communicatively connected to the at least one processor; wherein, The memory stores instructions that can be executed by the at least one processor to enable the at least one processor to perform the method of any one of claims 1-2.

6. A computer-readable medium having computer program instructions stored thereon, characterized in that, The computer program instructions can be executed by a processor to implement the method as described in any one of claims 1-2.

7. A computer program product, comprising a computer program, characterized in that, When the computer program is executed by a processor, it implements the method as described in any one of claims 1-2.

Images