CN108478936B - The method and apparatus of proton therapeutic dosage and range is determined by proton-induced thermoacoustic signal - Google Patents

Abstract

The invention discloses a method and equipment for determining the dose and range of proton therapy through proton-induced thermoacoustic signals. Induce the propagation path of the thermoacoustic signal in the body to perform ray tracing, and calculate the minimum travel time of the first wave; use the perturbation method to determine the partial derivative matrix of the travel time with respect to the position of the thermoacoustic source, and establish the inversion equation; obtain the thermoacoustic signal by iterative fitting Source position and velocity model; solve the double coupling problem of thermoacoustic source position disturbance and velocity model disturbance by parameter separation method; further correct the position of thermoacoustic source and the velocity model of proton-induced thermoacoustic signal propagation in the human body through the disturbance quantity ; Determine the dose and range of proton therapy through the energy distribution relationship of the proton beam. The experimental result obtained by the invention not only saves running time but also ensures positioning accuracy, and can be widely used in proton therapy.

Description

Method and method for determining dose and range of proton therapy by proton-induced thermoacoustic signal equipment

technical field

The present invention relates to the field of medical physics, in particular to a method and equipment for determining the dose and range of proton therapy through proton-induced thermoacoustic signals.

Background technique

For a long time, tumors have seriously threatened people's health. At present, radiation therapy is one of the commonly used methods to control and treat tumors, but in the process of clinical radiation therapy, as the radiation enters the human body, the dose and the depth of the radiation decay exponentially, which makes the normal tissues or organs before and after the tumor Suffer varying degrees of impact, resulting in radiation reactions and damage. Proton therapy for tumors can effectively improve this situation, because the physical characteristics of the Bragg (Bragg) peak are formed during proton beam irradiation, so that a large amount of energy is deposited at the position of the Bragg peak and quickly deposited. When the Bragg peak is located in the treatment target area, the tumor surrounding Normal tissues and organs receive very little dose, which greatly reduces treatment damage. In addition, the depth of the Bragg peak is energy-dependent. The outside world can properly control the peak width according to the tumor volume by adjusting the energy of the proton beam incident on the human body, so that the high-energy areas can be concentrated in tumor areas of different depths and sizes, thereby achieving purpose of treatment.

One of the difficulties in proton therapy is the positioning of the Bragg peak. Due to its physical characteristics, once the Bragg peak is not positioned correctly, a large amount of energy will be deposited in the non-lesional area, which will cause greater damage to normal tissues and organs. The analysis of proton dose and positioning by the KC Jones team in the United States is aimed at the propagation of proton beams in water. In their experiments, it was found that the proton beams would induce ultrasonic signals called α and γ waves during the propagation process, and α waves are Bragg peaks. The source of the γ wave is the Bragg peak of the proton, and the arrival time of the α and γ peaks at the sensor reflects the distance from the proton beam propagation axis and the center of the Bragg peak, respectively. Therefore, through the α, γ peaks The location of Bragg Peak can be known by the time difference. However, water is a homogeneous medium, and the propagation speed of sound waves in water is known, while the human body medium is inhomogeneous, and the speed model of sound waves propagating in the human body is unknown, so the method of the KC Jones team cannot be used in actual proton therapy used in .

Contents of the invention

The technical problem to be solved by the present invention is to provide a method and equipment for determining the dose and range of proton therapy through proton-induced thermoacoustic signals in view of the defects in the prior art.

The technical solution adopted by the present invention to solve its technical problems is:

The present invention provides a method for determining the dose and range of proton therapy through proton-induced thermoacoustic signals, comprising the following steps:

Step 1. Set up ultrasonic sensors at intervals near the irradiation path of the proton beam, receive proton-induced thermoacoustic signals through the ultrasonic sensors, and analyze them to extract travel time and waveform information;

Step 2. Use the linear travel time interpolation algorithm of wave front expansion to perform ray tracing on the propagation path of the proton-induced thermoacoustic signal in the body, and calculate the minimum travel time of the first arrival wave;

Step 3, using the perturbation method to determine the partial derivative matrix of the travel time with respect to the position of the thermal sound source, and establish the inversion equation;

Step 4. Solve the nonlinear least squares problem composed of the inversion equation, and obtain the position and velocity model of the thermal sound source through iterative fitting;

Step 5. Solve the problem of double coupling of the disturbance quantity of the thermoacoustic source position and the disturbance quantity of the velocity model through the parameter separation method;

Step 6, further modifying the position of the thermoacoustic source and the velocity model of the proton-induced thermoacoustic signal propagating in the human body through the disturbance amount;

Step 7. Determine the dose and range of proton therapy according to the energy distribution relationship of the proton beam.

Further, the thermal acoustic source position in the method of the present invention includes the Bragg peak position.

Further, in step 4 of the present invention, the Levenberg-Marquardt algorithm is used to solve the non-linear least squares problem constituted by the inversion equation.

Further, in step 5 of the present invention, the method for solving the dual coupling problem of the disturbance quantity of the thermoacoustic source position and the disturbance quantity of the velocity model through the parameter separation method is specifically as follows:

The inversion equation is: δt=Aδx+Bδv;

Among them, A is the partial derivative matrix of the travel time residual to the position of the thermal sound source, B is the partial derivative matrix of the travel time residual to the velocity, x is the position of the thermal sound source, v is the velocity model, and δv is the disturbance of the velocity model, δx is the disturbance amount at the location of the thermal sound source, and δt is the travel time residual;

First perform singular value decomposition A=UΛV ^T on A to obtain matrix U;

Then the matrix U is partitioned U=[U ₂ U ₁ U ₀ ], assuming that the matrix A has r singular values, p thermal sound source positions, and the number of sensors is d, then U ₂ is the first r columns of U, and U ₁ is the pr column, U ₀ is the dp column. Its purpose is to separate the parameters. After the parameters are separated by U, we get:

Calculate the disturbance amount δv of the velocity model according to the above formula, and then calculate the correction amount δ(δx) of the disturbance amount at the position of the thermal sound source by A ⁺ Aδ(δx)=-A ⁺ Bδv.

Further, in step 6 of the present invention, the method of further correcting the position of the thermoacoustic source and the velocity model of the proton-induced thermoacoustic signal propagating in the human body through the disturbance amount is specifically:

The disturbance amount δx of the thermal sound source position is corrected by the obtained correction amount δ(δx) of the disturbance amount of the thermal sound source position, that is, δx=δx+δ(δx), and is corrected by the corrected disturbance amount δx of the sound source position The initial sound source position is x ₀ , that is, x=x ₀ +δx, and the initial velocity model v ₀ is corrected by the obtained disturbance amount δv of the velocity model, that is, v=v ₀ +δv.

The present invention provides a device for determining the dose and range of proton therapy through proton-induced thermoacoustic signals, including:

A plurality of ultrasonic sensors are arranged at intervals near the irradiation path of the proton beam, and receive proton-induced thermoacoustic signals through the ultrasonic sensors;

The data processor, connected with the ultrasonic sensor, is used to analyze and extract the travel time and waveform information of the proton-induced thermoacoustic signal; use the linear travel time interpolation algorithm of wavefront expansion to perform ray tracing on the propagation path of the proton-induced thermoacoustic signal in the body, And calculate the minimum travel time of the first arrival wave; use the perturbation method to determine the partial derivative matrix of the travel time to the position of the thermal sound source, and establish the inversion equation; solve the nonlinear least squares problem composed of the inversion equation, and obtain the thermal sound source through iterative fitting Sound source position and velocity model; use parameter separation method to solve the double coupling problem of thermoacoustic source position disturbance and velocity model disturbance; further correct the position of thermoacoustic source and the velocity of proton-induced thermoacoustic signal propagation in the human body through disturbance quantity Model; determine the dose and range of proton therapy through the energy distribution relationship of the proton beam.

The beneficial effects produced by the present invention are: the method and equipment for determining the dosage and range of proton therapy by means of proton-induced thermoacoustic signals of the present invention use the ultrasonic signals excited when proton beams enter the human body (i.e., proton-induced thermoacoustic signals) to determine the proton therapy dose and range. The proton energy distribution range and therapeutic dose in the human body during treatment, the proton-induced thermoacoustic signal is received by an acoustic sensor installed on the human body, and the proton energy distribution is obtained by analyzing the travel time information, waveform information, and spectrum information of the thermoacoustic signal Parameters related to the scope and treatment dose; fully considering the uncertainty of proton energy distribution positioning under the unknown situation of the velocity model, and solving the problem of double coupling of position and velocity by using the method of parameter separation, the results obtained by the present invention will be More scientific and more precise.

Description of drawings

The present invention will be further described below in conjunction with accompanying drawing and embodiment, in the accompanying drawing:

Fig. 1 is a schematic structural diagram of an embodiment of the present invention.

Detailed ways

In order to make the object, technical solution and advantages of the present invention clearer, the present invention will be further described in detail below in conjunction with the accompanying drawings and embodiments. It should be understood that the specific embodiments described here are only used to explain the present invention, not to limit the present invention.

In order to solve the problem of positioning the Bragg peak in the process of proton therapy for tumors, the propagation speed of the unknown thermoacoustic signal in the human body. The invention discloses a method and equipment for determining the dosage and range of proton therapy by proton-induced thermoacoustic signals. Proton energy distribution range and therapeutic dose in the body. In this invention, it is necessary to use sensors placed at intervals along the proton irradiation direction close to the body surface of the human body to receive the thermoacoustic signal induced by the proton, and obtain the arrival time of the thermoacoustic signal. Then, the computer program simulates the propagation path of the thermoacoustic signal in the body to obtain the simulated travel time. Then use the perturbation method to determine the inversion equation and fit the initial position of the Bragg peak through the Levenberg-Marquardt algorithm. In order to correct the position information and the velocity model, the present invention adopts the velocity model and the disturbance amount of the position information separately obtained by parameter separation, and repeats the above steps several times until the travel time residual meets certain conditions, and outputs the final Bragg peak position and velocity model . Finally, the dose distribution at the Bragg peak is obtained through the energy distribution relationship. The experimental results obtained by the invention not only save running time but also ensure positioning accuracy, so it can be widely used in proton therapy.

As shown in Figure 1, the present invention is mainly based on proton-induced thermoacoustic signaling, considering that proton beams will induce alpha and gamma waves during propagation, and propose a method for determining the dose and range of proton therapy by proton-induced thermoacoustic signals Method and equipment, it uses the ultrasonic signal excited when the proton beam enters the human body to determine the proton energy distribution range and therapeutic dose in the human body during proton therapy, the proton-induced thermoacoustic signal is received by the acoustic sensor installed on the human body surface, and the thermal The travel time information, waveform information, and spectrum information of the acoustic signal are analyzed to obtain the relevant parameters of the proton energy distribution range and therapeutic dose. (1) Proton-excited ultrasound signals are used as the information source; (2) The required therapeutic monitoring parameters are obtained by analyzing the ultrasonic signals; (3) The obtained therapeutic monitoring parameters include but not limited to the position of the Bragg peak, etc. The invention fully considers the uncertainty of proton energy distribution positioning under the condition that the velocity model is unknown, and solves the problem of double coupling of position and velocity by using a method of parameter separation. The results obtained by the invention will be more scientific and more accurate.

In the embodiment of the present invention, x is defined as the position of the sound source (including but not limited to Bragg peak), v is defined as the velocity model, x ₀ is the initial sound source position, v ₀ is the initial velocity model, t _obs is the sensor received Travel time, t _cal is the minimum travel time calculated by using the ray tracing method, δx is the disturbance of the sound source position, δv is the disturbance of the velocity model, A is the partial derivative matrix of the travel time residual with respect to the sound source position, and B is the travel time The partial derivative matrix of the residual to the velocity, δt is the travel time residual, δt=t _obs -t _cal , and δ(δx) is the correction amount of the sound source position disturbance.

Step 1. Place sensors at intervals near the proton beam irradiation path on the patient so that they are close to the surface of the human body to reduce the interference of noise in the air to the signal. The arranged sensor is used to receive the thermoacoustic signal, record the signal waveform and extract the travel time t _obs .

Step 2: Set the initial velocity distribution v ₀ and the initial sound source position x ₀ through the prior knowledge and the patient's previous examination report, and grid the initial velocity distribution model, and the velocity in each grid is consistent. The linear traveltime interpolation (Linear Traveltime Interpolation, LTI) algorithm of wavefront expansion was used to trace the propagation path of the proton-induced thermoacoustic signal in the body, and the minimum traveltime t _cal of the first wave was calculated. The method of wavefront expansion can reduce the running time of the program while maintaining the calculation accuracy.

Step 3, use the perturbation method to determine the partial derivative matrix A of the travel time with respect to the position, and establish an inversion equation δt=Aδx that has nothing to do with the velocity.

Step 4, use the Levenberg-Marquardt algorithm to solve the nonlinear least squares function formed in step 3, and obtain the disturbance δx of the sound source position through multiple iterations of fitting.

Step 5, solve the problem of double coupling of the disturbance of the sound source position and the disturbance of the velocity model through the parameter separation method. The inversion equation is δt=Aδx+Bδv. Using the parameter separation method can reduce the calculation amount of the program and improve the stability of the program. In this part, the singular value decomposition A=UΛV ^T is performed on A first, and then the matrix U is partitioned U=[U ₂ U ₁ U ₀ ]. Suppose the matrix A has r singular values, p sound source positions, and the number of sensors is d. Then U ₂ is the first r column of U, U ₁ is the pr column, and U ₀ is the dp column. Its purpose is to separate parameters. After parameter separation by U, we get:

Calculate the disturbance amount δv of the velocity model according to the above formula, and then calculate the correction amount δ(δx) of the disturbance amount at the sound source position through A ⁺ Aδ(δx)=-A ⁺ Bδv.

Step 6, the correction amount δ(δx) of the disturbance amount of the sound source position obtained in step 5 is corrected Step 4 obtains the disturbance amount δx of the sound source position, that is, δx=δx+δ(δx), and the corrected sound source position is used The disturbance amount δx of , modifies the initial sound source position x ₀ , that is, x=x ₀ +δx. And modify the initial velocity model v ₀ by the disturbance δv of the velocity model obtained in step 5, that is, v=v ₀ +δv.

Step 7, repeat steps 2-6 until the objective function converges or reaches the preset number of cycles, at this time the obtained sound source position x is the position of our final desired thermoacoustic signal source (or Bragg peak). Through the energy distribution relationship of the proton beam, the proton energy distribution range and therapeutic dose in the human body during proton therapy are determined.

It should be understood that those skilled in the art can make improvements or changes based on the above description, and all these improvements and changes should belong to the protection scope of the appended claims of the present invention.

Claims (5)

1. A method for determining proton therapy dose and scope by proton-induced thermoacoustic signal, is characterized in that, comprises the following steps: Step 1. Set up ultrasonic sensors at intervals near the irradiation path of the proton beam, receive proton-induced thermoacoustic signals through the ultrasonic sensors, and analyze them to extract travel time and waveform information; Step 2. Use the linear travel time interpolation algorithm of wave front expansion to perform ray tracing on the propagation path of the proton-induced thermoacoustic signal in the body, and calculate the minimum travel time of the first arrival wave; Step 3, using the perturbation method to determine the partial derivative matrix of the travel time with respect to the position of the thermal sound source, and establish the inversion equation; Step 4. Solve the nonlinear least squares problem composed of the inversion equation, and obtain the position and velocity model of the thermal sound source through iterative fitting; Step 5. Solve the problem of double coupling of the disturbance quantity of the thermoacoustic source position and the disturbance quantity of the velocity model through the parameter separation method; In step 5, the method of solving the dual coupling problem of the disturbance quantity of the thermoacoustic source position and the disturbance quantity of the velocity model through the parameter separation method is as follows: The inversion equation is: δt=Aδx+Bδv; Among them, A is the partial derivative matrix of the travel time residual to the position of the thermal sound source, B is the partial derivative matrix of the travel time residual to the velocity, x is the position of the thermal sound source, v is the velocity model, and δv is the disturbance of the velocity model, δx is the disturbance amount at the location of the thermal sound source, and δt is the travel time residual; First perform singular value decomposition A=UΛV ^T on A to obtain matrix U; Then the matrix U is partitioned U=[U ₂ U ₁ U ₀ ], assuming that the matrix A has r singular values, p thermal sound source positions, and the number of sensors is d, then U ₂ is the first r columns of U, and U ₁ is the pr column, and U ₀ is the dp column. The purpose is to separate the parameters. After the parameters are separated by U, we get: Calculate the disturbance quantity δv of the velocity model according to the above formula, and then calculate the correction quantity δ(δx) of the disturbance quantity at the position of the thermoacoustic source through A ⁺ Aδ(δx)=-A ⁺ Bδv; Step 6, further modifying the position of the thermoacoustic source and the velocity model of the proton-induced thermoacoustic signal propagating in the human body through the disturbance amount; Step 7. Determine the dose and range of proton therapy according to the energy distribution relationship of the proton beam. 2. The method for determining the dose and range of proton therapy by proton-induced thermoacoustic signal according to claim 1, characterized in that, the position of the thermoacoustic source in the method comprises the Bragg peak position. 3. The method for determining the dose and range of proton therapy through proton-induced thermoacoustic signals according to claim 1, characterized in that in step 4, the Levenberg-Marquardt algorithm is used to solve the nonlinear least squares problem composed of inversion equations. 4. The method for determining the dose and range of proton therapy by proton-induced thermoacoustic signal according to claim 1, characterized in that, in step 6, the position of the thermoacoustic source and the position of the proton-induced thermoacoustic signal are further corrected by the disturbance amount in the human body. The method of velocity model of propagation is specifically as follows: The disturbance amount δx of the thermal sound source position is corrected by the obtained correction amount δ(δx) of the disturbance amount of the thermal sound source position, that is, δx=δx+δ(δx), and is corrected by the corrected disturbance amount δx of the sound source position The initial sound source position is x ₀ , that is, x=x ₀ +δx, and the initial velocity model v ₀ is corrected by the obtained disturbance amount δv of the velocity model, that is, v=v ₀ +δv. 5. A device for determining the dose and range of proton therapy through proton-induced thermoacoustic signals, characterized in that it includes: A plurality of ultrasonic sensors are arranged at intervals near the irradiation path of the proton beam, and receive proton-induced thermoacoustic signals through the ultrasonic sensors; The data processor, connected with the ultrasonic sensor, is used to analyze and extract the travel time and waveform information of the proton-induced thermoacoustic signal; use the linear travel time interpolation algorithm of wavefront expansion to perform ray tracing on the propagation path of the proton-induced thermoacoustic signal in the body, And calculate the minimum travel time of the first arrival wave; use the perturbation method to determine the partial derivative matrix of the travel time to the position of the thermal sound source, and establish the inversion equation; solve the nonlinear least squares problem composed of the inversion equation, and obtain the thermal sound source through iterative fitting Sound source position and velocity model; use parameter separation method to solve the double coupling problem of thermoacoustic source position disturbance and velocity model disturbance; further correct the position of thermoacoustic source and the velocity of proton-induced thermoacoustic signal propagation in the human body through disturbance quantity Model; determine the dose and range of proton therapy through the energy distribution relationship of the proton beam; The method to solve the double coupling problem of the disturbance quantity of the thermoacoustic source position and the disturbance quantity of the velocity model through the parameter separation method is as follows: The inversion equation is: δt=Aδx+Bδv; Among them, A is the partial derivative matrix of the travel time residual to the position of the thermal sound source, B is the partial derivative matrix of the travel time residual to the velocity, x is the position of the thermal sound source, v is the velocity model, and δv is the disturbance of the velocity model, δx is the disturbance amount at the location of the thermal sound source, and δt is the travel time residual; First perform singular value decomposition A=UΛV ^T on A to obtain matrix U; Then the matrix U is partitioned U=[U ₂ U ₁ U ₀ ], assuming that the matrix A has r singular values, p thermal sound source positions, and the number of sensors is d, then U ₂ is the first r columns of U, and U ₁ is the pr column, and U ₀ is the dp column. The purpose is to separate the parameters. After the parameters are separated by U, we get: Calculate the disturbance amount δv of the velocity model according to the above formula, and then calculate the correction amount δ(δx) of the disturbance amount at the position of the thermal sound source by A ⁺ Aδ(δx)=-A ⁺ Bδv.