A method and system for removing reflection artifacts in endoscopic photoacoustic tomography images

Abstract

The present invention provides a method and system for removing reflection artifacts in endoscopic photoacoustic tomography images, which relate to the field of medical imaging technology. The method comprises the following steps: obtaining a photoacoustic data group of a plane to be imaged; inputting the photoacoustic data group into a reflection artifact removal model to obtain a light absorption energy distribution diagram after removing the reflection artifacts at the plane to be imaged; using the output of a deep gradient descent network as the input of a deep learning segmentation network to obtain a RAR‑Net neural network, training the RAR‑Net neural network to obtain a reflection artifact removal model, and being able to remove reflection artifacts in endoscopic photoacoustic tomography images and improve the imaging accuracy of endoscopic photoacoustic tomography images.

Description

Technical Field

[0001] The present invention relates to the field of medical imaging technology, and in particular to a method and system for removing reflection artifacts in endoscopic photoacoustic tomography images. Background Art

[0002] Endoscopic photoacoustic tomography (EPAT) is a new type of non-ionizing interventional imaging method that combines the advantages of ultrasound imaging and optical imaging, and can provide good optical contrast and ultrasound resolution. Its imaging principle is: a special catheter is inserted into the lumen to be tested and pushed to the distal end. During the process of slowly withdrawing the catheter, the probe at the top of the catheter emits a non-focused short-pulse laser to irradiate the surrounding tissue. The tissue expands due to heat, and the light energy is quickly converted into heat energy and mechanical stress, thereby generating broadband ultrasound (i.e., photoacoustic waves) and propagating to the tissue surface. The ultrasound detector collects the photoacoustic signal and sends it to the computer. By solving the inverse problem, the initial sound pressure distribution map (structural imaging) and optical characteristic parameter distribution map (quantitative imaging) on ​​the cross section of the cavity are reconstructed.

[0003] The quality of photoacoustic images is closely related to the beamforming process. Beamforming technology is based on the measurement of ultrasonic propagation time. However, this measurement does not take into account that ultrasonic waves may be reflected multiple times when propagating in the medium and each reflected wave may be detected by the ultrasonic transducer. As a result, the reverberation signal is mapped to the wrong position during image reconstruction, resulting in artifacts in the image, namely reflection artifacts, and also reducing the imaging depth. Reflection artifacts are the most common artifacts in photoacoustic images. If reflection artifacts appear in the imaging area of ​​interest, they may be misjudged as useful targets, reducing the accuracy of image recognition and clinical diagnosis. In addition, the acoustic properties of biological tissues are spatially inhomogeneous. Tissue density and the propagation speed of ultrasonic waves in tissues vary with spatial position. When ultrasonic waves propagate in tissues of different components, there are also absorption, scattering and diffusion phenomena, which lead to the attenuation of photoacoustic signals. The phase deviation caused by the change in sound speed and the signal broadening related to acoustic attenuation will reduce the spatial resolution of the image. Summary of the invention

[0004] The object of the present invention is to provide a method and system for removing reflection artifacts in endoscopic photoacoustic tomography images, which can remove reflection artifacts in endoscopic photoacoustic tomography images and improve the imaging accuracy of endoscopic photoacoustic tomography images.

[0005] To achieve the above object, the present invention provides the following solutions:

[0006] A method for removing reflection artifacts in endoscopic photoacoustic tomography images, comprising:

[0007] Acquiring a photoacoustic data set of a plane to be imaged;

[0008] The photoacoustic data set is input into a reflection artifact removal model to obtain a light absorption energy distribution diagram after the reflection artifacts are removed from the imaging plane; the artifact removal model is obtained by training a RAR-Net neural network based on a historical photoacoustic data set; the RAR-Net neural network is obtained by using the output of a deep gradient descent network as the input of a deep learning segmentation network.

[0009] Optionally, before acquiring the photoacoustic signal of the plane to be imaged, the method further includes:

[0010] Acquire historical photoacoustic data groups of multiple imaging planes and historical light absorption energy distribution diagrams of multiple imaging planes after removing reflection artifacts as training data;

[0011] Dividing the training data into a plurality of training sets;

[0012] Using multiple training sets, the RAR-Net neural network was trained with historical photoacoustic data sets as input and historical light absorption energy distribution maps after removing reflection artifacts as the expected output to obtain a reflection artifact removal model.

[0013] Optionally, the method uses multiple training sets, takes the historical photoacoustic data set as input, takes the historical light absorption energy distribution map after removing the reflection artifacts as the expected output, trains the RAR-Net neural network to obtain the reflection artifact removal model, specifically including:

[0014] Determine that the RAR-Net neural network before training is the k×epoch×Ith order RAR-Net neural network;

[0015] Let the iteration round k = 0;

[0016] Set the number of training set traversal epochs to 0;

[0017] Determine any training set as the current training set;

[0018] Let the index of the current training set I = 1;

[0019] Inputting the current training set into a k×epoch×Ith order RAR-Net neural network, adjusting the parameters of the k×epoch×Ith order RAR-Net neural network, and obtaining a k×epoch×(I+1)th order RAR-Net neural network;

[0020] Update the current training set, increase the value of I by 1, and return to the step of "inputting the current training set into the k×epoch×Ith order RAR-Net neural network, adjusting the parameters of the k×epoch×Ith order RAR-Net neural network, and obtaining the k×epoch×(I+1)th order RAR-Net neural network" until the training set is traversed and the k×(epoch+1)th order RAR-Net neural network is obtained;

[0021] Let the value of epoch increase by 1 and return to the step "set the index of the current training set I = 1" until epoch is greater than the threshold of the number of training set traversal times, and the (k+1)×0×1th order RAR-Net neural network is obtained;

[0022] Increase the value of k by 1 and return to the step “set the number of training set traversal epochs = 1” until k is greater than the iteration round threshold, and obtain a k×epoch×I-order RAR-Net neural network as a reflection artifact removal model.

[0023] Optionally, inputting the current training set into a k×epoch×Ith order RAR-Net neural network, adjusting parameters of the k×epoch×Ith order RAR-Net neural network, and obtaining a k×epoch×(I+1)th order RAR-Net neural network specifically includes:

[0024] Inputting the historical photoacoustic data set in the current training set into the k×epoch×I-th order RAR-Net neural network to obtain a historical prediction distribution diagram of light absorption energy at the k×epoch×I-th iteration;

[0025] According to the historical light absorption energy distribution map and the historical predicted light absorption energy distribution map at the k×epoch×Ith iteration, the loss of the deep gradient descent network is determined as the first loss and the loss of the deep learning segmentation network is determined as the second loss;

[0026] Determine the linear activation response of the output layer of the deep gradient descent network as a first linear activation response;

[0027] Determine the linear activation response of the output layer of the deep learning segmentation network as the second linear activation response;

[0028] Determine an error vector of an output layer of a deep gradient descent network as a first output error vector according to the first linear activation response;

[0029] Determine, according to the second linear activation response, an error vector of an output layer of the deep learning segmentation network as a second output error vector;

[0030] According to the first loss and the first output error vector, using a first recursive formula to determine the error vector of each layer of the deep gradient descent network as a first error vector group;

[0031] According to the second loss and the second output error vector, using a second recursive formula to determine the error vector of each layer of the deep learning segmentation network as a second error vector group;

[0032] Determine the gradient of each layer of the k×epoch×I-order RAR-Net neural network according to the first error vector group and the second error vector group;

[0033] According to the gradient of each layer of the k×epoch×I-th order RAR-Net neural network, the Adam algorithm is used to update multiple gradients to obtain the k×epoch×(I+1)-th order RAR-Net neural network.

[0034] Optionally, the first loss is:

[0035] The second loss is:

[0036] The first output error vector is:

[0037] The second output error vector is:

[0038] The first recursive formula is: G,c =(W G,c+1 ) T δ G,c+1 ⊙σ'(Z G,c );

[0039] The second recursive formula is: U,c =(W U,c+1 ) T δ U,c+1 ⊙σ'(Z U,c );

[0040] The gradient of each layer of the RAR-Net neural network is:

[0041] Among them, L G (W G ,b G ) and L U (W U ,b U ) are the first loss and the second loss respectively; W G and b Gare the weight parameter set and bias parameter set of the deep gradient descent network respectively; W U and b U are the weight parameter set and bias parameter set of the deep learning segmentation network respectively; i is the sample number in the training set, i = (1, 2, ..., M); M is the sample capacity of the current training set; x true,i is the historical light absorption energy distribution diagram; x k+1,i and x k,i are the predicted distribution diagrams of the light absorption energy history output at the kth and k+1th iterations respectively; θ k is the parameter of the kth network structure unit; are the same structure with different learning parameters θ k RAR-Net i is the historical photoacoustic data set, x i is the predicted distribution map of light absorption energy history output by the deep learning segmentation network; H(x k,i ) is the same as x k,i The operator concerned; δ G,9 and δ U,37 are the first output error vector and the second output error vector respectively; σ'(·) represents the first-order derivative of the ReLU function, Z G,9 and Z U,37 are the first linear activation response and the second linear activation response respectively; δ G,c , δ U,c , δ G,c+1 and δ U,c+1 are the error vectors of the cth layer and the c+1th layer of the RAR-Net neural network, respectively, where the subscript G represents a deep gradient descent network, and the subscript U represents a deep learning segmentation network; W G,c+1 and W U,c+1 are the weight parameter matrices of the c+1th layer of the RAR-Net neural network; Z G,c and Z U,c are the linear activation responses of the cth layer of the RAR-Net neural network; L(W G ,b G ,W U ,b U ) is the total loss of the RAR-Net neural network; L(W G ,b G ,W U ,b U )=αL G (W G ,b G )+L U (W U ,b U ), α is the scale factor: determined by the convergence speed and reconstruction quality; Wc 、b c , Z c and δ c are the weight parameter matrix, bias parameter vector, linear activation response and error vector of the cth layer of the RAR-Net neural network; in the deep gradient descent network, that is, when c=1,2,3,…,9, W c =W G 、b c =b G , Z c =Z G ,δ c =δ G,c ; In the deep learning segmentation network, that is, when c=10,11,...,37, W c =W U 、b c =b U , Z c =Z U ,δ c =δ U,c ;X c-1 is the output of the c-1th layer network in the forward propagation.

[0042] A system for removing reflection artifacts in endoscopic photoacoustic tomography images, comprising:

[0043] A photoacoustic data set acquisition module, used to acquire a photoacoustic data set of a plane to be imaged;

[0044] A reflection artifact removal imaging module is used to input the photoacoustic data set into a reflection artifact removal model to obtain a light absorption energy distribution diagram after removing reflection artifacts at the plane to be imaged; the artifact removal model is obtained by training a RAR-Net neural network based on a historical photoacoustic data set; the RAR-Net neural network is obtained by using the output of a deep gradient descent network as an input to a deep learning segmentation network.

[0045] Optionally, the system further includes:

[0046] A training data acquisition module, used to acquire historical photoacoustic data groups of multiple imaging planes and historical light absorption energy distribution diagrams of multiple imaging planes after removing reflection artifacts as training data;

[0047] A training set division module, used for dividing the training data into multiple training sets;

[0048] The reflection artifact removal model determination module is used to use multiple training sets, take the historical photoacoustic data group as input, and take the historical light absorption energy distribution map after removing the reflection artifact as the expected output, to train the RAR-Net neural network and obtain the reflection artifact removal model.

[0049] Optionally, the reflection artifact removal model determination module specifically includes:

[0050] An initial neural network determination unit, used to determine that the RAR-Net neural network before training is a k×epoch×I-th order RAR-Net neural network;

[0051] An iteration round assignment unit, used to set the iteration round k=0;

[0052] The training set traversal times assignment unit is used to set the training set traversal times epoch = 0;

[0053] A current training set determining unit, used to determine any training set as the current training set;

[0054] An index assignment unit, used to set the index I of the current training set to 1;

[0055] A neural network updating unit, used for inputting the current training set into a k×epoch×I-th order RAR-Net neural network, adjusting the parameters of the k×epoch×I-th order RAR-Net neural network, and obtaining a k×epoch×(I+1)-th order RAR-Net neural network;

[0056] The current training set updating unit is used to update the current training set, increase the value of I by 1, and call the neural network updating unit until the training set is traversed to obtain the k×(epoch+1)×1th order RAR-Net neural network;

[0057] A training set traversal times updating unit is used to increase the value of epoch by 1 and call the index assignment unit until epoch is greater than the training set traversal times threshold, thereby obtaining a (k+1)×0×1th order RAR-Net neural network;

[0058] The iteration round updating unit is used to increase the value of k by 1 and return to the step of calling the training set traversal times assignment unit until k is greater than the iteration round threshold, and obtain a k×epoch×I-order RAR-Net neural network as a reflection artifact removal model.

[0059] Optionally, the neural network updating unit specifically includes:

[0060] A light absorption energy history prediction distribution map determining subunit is used to input the historical photoacoustic data group in the current training set into the k×epoch×I-th order RAR-Net neural network to obtain the light absorption energy history prediction distribution map at the k×epoch×I-th iteration;

[0061] a loss determination subunit, for determining the loss of the deep gradient descent network as a first loss and determining the loss of the deep learning segmentation network as a second loss according to the historical light absorption energy distribution map and the historical predicted light absorption energy distribution map at the k×epoch×Ith iteration;

[0062] A first linear activation response determination subunit, used to determine that the linear activation response of the output layer of the deep gradient descent network is a first linear activation response;

[0063] A second linear activation response determining subunit, used to determine that the linear activation response of the output layer of the deep learning segmentation network is a second linear activation response;

[0064] A first output error vector determining subunit, configured to determine an error vector of an output layer of a deep gradient descent network as a first output error vector according to the first linear activation response;

[0065] A second output error vector determining subunit, configured to determine, according to the second linear activation response, an error vector of an output layer of the deep learning segmentation network as a second output error vector;

[0066] A first error vector group determining subunit, configured to determine, according to the first loss and the first output error vector, the error vector of each layer of the deep gradient descent network as a first error vector group using a first recursive formula;

[0067] A second error vector group determining subunit, configured to determine, according to the second loss and the second output error vector, the error vector of each layer of the deep learning segmentation network as a second error vector group using a second recursive formula;

[0068] A gradient determination subunit, configured to determine the gradient of each layer of the k×epoch×I-order RAR-Net neural network according to the first error vector group and the second error vector group;

[0069] The gradient updating subunit is used to update multiple gradients according to the gradient of each layer of the k×epoch×I-th order RAR-Net neural network using the Adam algorithm to obtain the k×epoch×(I+1)-th order RAR-Net neural network.

[0070] Optionally, the first loss is:

[0071] The second loss is:

[0072] The first output error vector is:

[0073] The second output error vector is:

[0074] The first recursive formula is: G,c =(W G,c+1 ) T δ G,c+1 ⊙σ'(Z G,c );

[0075] The second recursive formula is: U,c =(W U,c+1 ) T δ U,c+1 ⊙σ'(Z U,c );

[0076] The gradient of each layer of the RAR-Net neural network is:

[0077] Among them, L G (W G ,b G ) and L U (W U ,b U ) are the first loss and the second loss respectively; W G and b G are the weight parameter set and bias parameter set of the deep gradient descent network respectively; W U and b U are the weight parameter set and bias parameter set of the deep learning segmentation network respectively; i is the sample number in the training set, i = (1, 2, ..., M); M is the sample capacity of the current training set; x true,i is the historical light absorption energy distribution diagram; x k,i and x k+1,i are the predicted distribution diagrams of light absorption energy history obtained at the kth and k+1th iterations respectively; θ k is the parameter of the kth network structure unit; are the same structure with different learning parameters θ k RAR-Net i is the historical photoacoustic data set, x i is the predicted distribution map of light absorption energy history output by the deep learning segmentation network; H(x k,i ) is the same as x k,i The operator concerned; δ G,9 and δ U,37 are the first output error vector and the second output error vector respectively; σ'(·) represents the first-order derivative of the ReLU function, Z G,9 and Z U,37 are the first linear activation response and the second linear activation response respectively; δ G,c , δ U,c, δ G,c+1 and δ U,c+1 are the error vectors of the cth layer and the c+1th layer of the RAR-Net neural network, respectively, where the subscript G represents a deep gradient descent network, and the subscript U represents a deep learning segmentation network; W G,c+1 and W U,c+1 are the weight parameter matrices of the c+1th layer of the RAR-Net neural network; Z G,c and Z U,c are the linear activation responses of the cth layer of the RAR-Net neural network; L(W G ,b G ,W U ,b U ) is the total loss of the RAR-Net neural network: L(W G ,b G ,W U ,b U )=αL G (W G ,b G )+L U (W U ,b U ), α is the scale factor: determined by the convergence speed and reconstruction quality; W c 、b c , Z c and δ c are the weight parameter matrix, bias parameter vector, linear activation response and error vector of the cth layer of the RAR-Net neural network; in the deep gradient descent network, that is, when c=1,2,3,…,9, W c =W G 、b c =b G , Z c =Z G , δ c =δ G,c ; In the deep learning segmentation network, that is, when c=10,11,...,37, W c =W U 、b c =b U , Z c =Z U , δ c =δ U,c ;X c-1 is the output of the c-1th layer network in the forward propagation.

[0078] According to the specific embodiments provided by the present invention, the present invention discloses the following technical effects:

[0079] The present invention provides a method and system for removing reflection artifacts in endoscopic photoacoustic tomography images, which obtain a photoacoustic data group of a plane to be imaged; input the photoacoustic data group into a reflection artifact removal model to obtain a light absorption energy distribution diagram after the reflection artifacts are removed at the plane to be imaged; the present invention uses the output of a deep gradient descent network as the input of a deep learning segmentation network, constructs a RAR-Net neural network, trains the RAR-Net neural network to obtain a reflection artifact removal model, and can remove reflection artifacts in endoscopic photoacoustic tomography images, thereby improving the imaging accuracy of endoscopic photoacoustic tomography images. BRIEF DESCRIPTION OF THE DRAWINGS

[0080] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the figures required for use in the embodiments will be briefly introduced below. Obviously, the figures described below are only some embodiments of the present invention. For ordinary technicians in this field, other figures can be obtained based on these figures without paying creative labor.

[0081] Figure 1 Flow chart of a method for removing reflection artifacts in endoscopic photoacoustic tomography images according to an embodiment of the present invention;

[0082] Figure 2 Schematic diagram of endoscopic photoacoustic tomography in an embodiment of the present invention;

[0083] Figure 3 is a schematic diagram of a cross section of a cavity in an embodiment of the present invention;

[0084] Figure 4 This is a diagram of the RAR-Net structure built by the method in the embodiment of the present invention;

[0085] Figure 5 A flow chart of designing a forward propagation RAR-Net in an embodiment of the present invention;

[0086] Figure 6 Flowchart of training and testing RAR-Net in an embodiment of the present invention. DETAILED DESCRIPTION

[0087] The following will be combined with the figures in the embodiments of the present invention to clearly and completely describe the technical solutions in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without creative work are within the scope of protection of the present invention.

[0088] The object of the present invention is to provide a method and system for removing reflection artifacts in endoscopic photoacoustic tomography images, which can remove reflection artifacts in endoscopic photoacoustic tomography images and improve the imaging accuracy of endoscopic photoacoustic tomography images.

[0089] In order to make the above-mentioned objects, features and advantages of the present invention more obvious and easy to understand, the present invention is further described in detail below with reference to the drawings and specific implementation methods.

[0090] like Figure 1 The present invention provides a method for removing reflection artifacts in endoscopic photoacoustic tomography images, comprising:

[0091] Acquiring a photoacoustic data set of a plane to be imaged;

[0092] The photoacoustic data set is input into the reflection artifact removal model to obtain the light absorption energy distribution map after the reflection artifact is removed at the imaging plane; the artifact removal model is obtained by training the RAR-Net neural network based on the historical photoacoustic data set; the RAR-Net neural network is obtained by using the output of the deep gradient descent network as the input of the deep learning segmentation network.

[0093] The method for removing reflection artifacts in endoscopic photoacoustic tomography images provided by the present invention further comprises, before obtaining the photoacoustic signal of the plane to be imaged:

[0094] Acquire historical photoacoustic data groups of multiple imaging planes and historical light absorption energy distribution diagrams of multiple imaging planes after removing reflection artifacts as training data;

[0095] Divide the training data into multiple training sets;

[0096] Using multiple training sets, the RAR-Net neural network was trained with historical photoacoustic data sets as input and historical light absorption energy distribution maps after removing reflection artifacts as the expected output to obtain a reflection artifact removal model.

[0097] Specifically, using multiple training sets, taking the historical photoacoustic data set as input and the historical light absorption energy distribution map after removing the reflection artifacts as the expected output, the RAR-Net neural network is trained to obtain a reflection artifact removal model, which specifically includes:

[0098] Determine that the RAR-Net neural network before training is the k×epoch×Ith order RAR-Net neural network;

[0099] Let the iteration round k = 0;

[0100] Set the number of training set traversal epochs to 0;

[0101] Determine any training set as the current training set;

[0102] Let the index of the current training set I = 1;

[0103] Input the current training set into the k×epoch×Ith order RAR-Net neural network, adjust the parameters of the k×epoch×Ith order RAR-Net neural network, and obtain the k×epoch×(I+1)th order RAR-Net neural network;

[0104] Update the current training set, increase the value of I by 1, and return to the step "input the current training set into the k×epoch×Ith order RAR-Net neural network, adjust the parameters of the k×epoch×Ith order RAR-Net neural network, and obtain the k×epoch×(I+1)th order RAR-Net neural network" until the training set is traversed and the k×(epoch+1)th order RAR-Net neural network is obtained;

[0105] Let the value of epoch increase by 1 and return to the step "set the index of the current training set I = 1" until epoch is greater than the threshold of the number of training set traversal times, and the (k+1)×0×1th order RAR-Net neural network is obtained;

[0106] Increase the value of k by 1 and return to the step “set the number of training set traversal epochs = 1” until k is greater than the iteration round threshold, and obtain a k×epoch×I-order RAR-Net neural network as a reflection artifact removal model.

[0107] Furthermore, the current training set is input into the k×epoch×Ith order RAR-Net neural network, and the parameters of the k×epoch×Ith order RAR-Net neural network are adjusted to obtain the k×epoch×(I+1)th order RAR-Net neural network, which specifically includes:

[0108] Input the historical photoacoustic data set in the current training set into the k×epoch×I-th order RAR-Net neural network to obtain the historical prediction distribution diagram of light absorption energy at the k×epoch×I-th iteration;

[0109] According to the historical light absorption energy distribution map and the historical predicted light absorption energy distribution map at the k×epoch×Ith iteration, the loss of the deep gradient descent network is determined as the first loss and the loss of the deep learning segmentation network is determined as the second loss;

[0110] Determine the linear activation response of the output layer of the deep gradient descent network as a first linear activation response;

[0111] Determine the linear activation response of the output layer of the deep learning segmentation network as the second linear activation response;

[0112] Determine an error vector of an output layer of the deep gradient descent network as a first output error vector according to the first linear activation response;

[0113] Determine, according to the second linear activation response, an error vector of an output layer of the deep learning segmentation network as a second output error vector;

[0114] According to the first loss and the first output error vector, the error vector of each layer of the deep gradient descent network is determined as a first error vector group using a first recursive formula;

[0115] According to the second loss and the second output error vector, the error vector of each layer of the deep learning segmentation network is determined as a second error vector group using a second recursive formula;

[0116] Determine the gradient of each layer of the k×epoch×I-order RAR-Net neural network according to the first error vector group and the second error vector group;

[0117] According to the gradient of each layer of the k×epoch×I-th order RAR-Net neural network, multiple gradients are updated using the Adam algorithm to obtain the k×epoch×(I+1)-th order RAR-Net neural network.

[0118] in,

[0119] The first loss is:

[0120] The second loss is:

[0121] The first output error vector is:

[0122] The second output error vector is:

[0123] The first recursive formula is: G,c =(W G,c+1 ) T δ G,c+1 ⊙σ'(Z G,c );

[0124] The second recursive formula is: U,c =(W U,c+1 ) T δ U,c+1 ⊙σ'(Z U,c );

[0125] The gradient of each layer of the RAR-Net neural network is:

[0126]

[0127] Among them, L G(W G ,b G ) and L U (W U ,b U ) are the first loss and the second loss respectively; W G and b G are the weight parameter set and bias parameter set of the deep gradient descent network respectively; W U and b U are the weight parameter set and bias parameter set of the deep learning segmentation network respectively; i is the sample number in the training set, i = (1, 2, ..., M); M is the sample capacity of the current training set; x true,i is the historical light absorption energy distribution diagram; x k,i and x k+1,i are the predicted distribution diagrams of light absorption energy history obtained at the kth and k+1th iterations respectively; θ k is the parameter of the kth network structure unit; are the same structure with different learning parameters θ k RAR-Net i is the historical photoacoustic data set, x i is the predicted distribution map of light absorption energy history output by the deep learning segmentation network; H(x k,i ) is the same as x k,i The operator concerned; δ G,9 and δ U,37 are the first output error vector and the second output error vector respectively; σ'(·) represents the first-order derivative of the ReLU function, Z G,9 and Z U,37 are the first linear activation response and the second linear activation response respectively; δ G,c , δ U,c , δ G,c+1 and δ U,c+1 are the error vectors of the cth layer and the c+1th layer of the RAR-Net neural network, respectively, where the subscript G represents a deep gradient descent network, and the subscript U represents a deep learning segmentation network; W G,c+1 and W U,c+1 are the weight parameter matrices of the c+1th layer of the RAR-Net neural network; Z G,c and Z U,c are the linear activation responses of the cth layer of the RAR-Net neural network; L(W G ,b G ,W U ,b U ) is the total loss of the RAR-Net neural network, L(W G ,b G ,W U ,b U )=αLG (W G ,b G )+L U (W U ,b U ), α is the scale factor: determined by the convergence speed and reconstruction quality; W c 、b c , Z c and δ c are the weight parameter matrix, bias parameter vector, linear activation response and error vector of the cth layer of the RAR-Net neural network; in the deep gradient descent network, that is, when c=1,2,3,…,9, W c =W G 、b c =b G , Z c =Z G , δ c =δ G,c ; In the deep gradient descent network, that is, when c=10,11,...,37, W c =W U 、b c =b U , Z c =Z U , δ c =δ U,c ;X c-1 is the output of the c-1th layer network in the forward propagation.

[0128] The object of the present invention is to provide a convolutional neural network (CNN) for removing reflection artifacts from endoscopic photoacoustic tomography images, namely RAR-Net (Reflection Artifacts Removal network). The photoacoustic signal with reflection clutter collected by an ultrasonic detector in a biological cavity is input into the trained RAR-Net, and finally a light absorption energy distribution map with suppressed reflection artifacts is output. First, a data set suitable for removing acoustic reflection artifacts from endoscopic photoacoustic tomography images is constructed by a computer simulation method, and the photoacoustic signals at various locations in the imaging plane obtained by forward simulation are used as the input of samples in the data set, and the simulated light absorption energy distribution map is used as the expected output of the sample. Then, RAR-Net is built, and each layer of the network is pre-trained to evaluate the network performance, and then the hyperparameters of each layer of the network are optimized; finally, the measured value of the photoacoustic signal with acoustic reflection clutter is input into the trained RAR-Net, and a high-quality light absorption energy distribution image is output through feature extraction and feature fusion.

[0129] 1. Constructing a deep learning dataset for removing reflection artifacts from endoscopic photoacoustic tomography images

[0130] (1) Establishing a computer simulation model of blood vessel cross section

[0131] like Figure 2 As shown, a computer simulation model of a biological cavity cross section is established, wherein the imaging catheter is located at the center of the cross section, and the surrounding parts are the lumen, the diseased tissue and the tube wall from the inside to the outside.

[0132] (2) Forward simulation of the imaging process

[0133] The Monte Carlo method is used to simulate the transmission process of light in tissues to obtain the light radiation flux at each location in the imaging plane. Combined with the preset tissue light absorption coefficient, the theoretical value of light absorption energy density is obtained:

[0134] A(r)=μ a (r)Φ(r) (1)

[0135] Then the theoretical value of initial sound pressure is obtained:

[0136] p0(r)=Γ(r)A(r) (2)

[0137] In the formula, μ a (r) is the tissue light absorption coefficient at position r, Φ(r), A(r), p0(r) and Γ(r) are the light radiation flux, light absorption energy density, initial acoustic pressure and thermal expansion coefficient at position r in the imaging plane, respectively.

[0138] The following wave equation is used to describe the propagation of photoacoustic waves in lossy media:

[0139]

[0140] Where p(r,t) is the sound pressure at position r and time t, u(r,t) is the vibration velocity of the medium at position r and time t, c(r) is the sound velocity at position r, ρ(r,t) is the acoustic density of the inhomogeneous medium, ρ0 is the density of the homogeneous medium, τ(r) and η(r) are the acoustic absorption coefficient and dispersion coefficient, respectively:

[0141]

[0142] Where a is the power law exponent (for biological tissues, its value range is 1≤a≤1.5), μ=(10 -7 / 2π)cm - 1 rad -1 s is the absorption coefficient.

[0143] The initial condition for formula (3) to be valid is:

[0144]

[0145] The k-space pseudo-spectral method is used to discretize equation (3) to obtain:

[0146]

[0147] in

[0148] Where j is an imaginary unit, ξ = (x, y) is a point in the two-dimensional imaging plane, and k ξ yes The spatial wave number components in the ξ direction, F{…} and F -1 {…} are the two-dimensional Fourier transform and inverse transform, Δt is the time step, ρ ξ (r,t) and u ξ (r, t) are the components of the acoustic density and vibration velocity of the medium in the directions of the coordinate axes of the XY rectangular coordinate system, abs and disp represent the absorption term and dispersion term, respectively. The discrete representation is:

[0149]

[0150]

[0151] The k-space method is used to calculate the time step Δt:

[0152]

[0153] Where c0 is the speed of sound in a homogeneous medium, Δx is the grid distance, and CFL = 0.3 is the trade-off between simulation accuracy and computational speed. The initial acoustic density is

[0154]

[0155] (3) Construct and divide the data set

[0156] The photoacoustic signal distribution map and light absorption energy distribution map obtained by simulation are used as the input and expected output of the samples in the data set, respectively, and a data set suitable for removing acoustic reflection artifacts in EPAT images is obtained.

[0157] All samples in the data set are randomly shuffled and divided into training set, validation set and test set in a ratio of 6:2:2. The validation set is used to select the best network model. The training set is expanded by using data enhancement methods such as rotation, scaling, mirror flipping and translation, and the samples are evenly divided into N groups, each of which is a small batch training set containing M samples.

[0158] 2. Design of RAR-Net for forward propagation

[0159] like Figure 3 As shown, the RAR-Net constructed by the method of the present invention is a combination of a deep gradient descent (DGD) network and a U-net. The DGD network realizes the transition from the signal domain to the image domain, and the U-net (deep learning segmentation network) optimizes low-quality images and outputs high-quality images with suppressed reflection artifacts.

[0160] The construction process of RAR-Net is as follows Figure 4 As shown, the specific steps 2- are as follows:

[0161] Step 2-1 Determine the network architecture and parameters

[0162] The convolution operation in the DGD network uses a filter kernel size of 3×3×T, where T is the number of feature maps in the current layer input. The same convolution with a step size of 1 is used. The initial number of filter kernels (i.e., the number of feature channels) of the two input images of the network unit is 16, and the activation function uses the linear rectification function ReLU. The pixel values ​​of the initial input and output images of the network unit are summed through the far jump connection to obtain the output image of the current layer. For example Figure 3 In the example, the initial input of the first network structure unit is x0, and the output image is the feature image with a quantity of 1. The output image of the current layer is obtained by simply summing the pixel values ​​of the two.

[0163] The convolution operation in U-net uses a filter kernel size of 3×3×T. The last convolution operation uses a 1×1 convolution, and its filter kernel size is 1×1×T. The convolution method uses the same convolution with a step size of 1, where the initial number of filter kernels, that is, the number of feature channels, is u =64. The activation function uses the linear rectification function ReLU, the filter kernel size used in the pooling operation is 2×2×T, and the pooling method uses the maximum pooling with a step size of 2. The filter kernel size used in the upsampling operation is 2×2×T, and the step size is 2.

[0164] Step 2-2 Initialize network parameters

[0165] Maximum number of network layers k max =2, layer number k=0, feature layer number c=1.

[0166] Step 2-3 calculates the gradient of the k-th layer network structure unit:

[0167]

[0168] Where H is the forward imaging operator (i.e., the photoacoustic signal generated by the tissue is obtained from the light absorption energy density), H *is the adjoint operator of H, ▽d(·) is the data fitting information, y i is the input photoacoustic signal matrix of the i-th sample in the training set, x k,i is the light absorption energy distribution diagram output by the i-th sample in the data set after the k-th iteration, where i is the sample number in the training set. k,i and ▽d(y i ,H(x k,i )) as the input of the network.

[0169] In step 2-4, set the number of convolutions conv=0 and the number of feature channels channels=16.

[0170] Step 2-5: Input the M samples in the first group of small batch training sets in the training set into RAR-Net as feature map c , which is also the feature layer generated by the input photoacoustic signal matrix after the c-1th convolution.

[0171] Step 2-6 Map c The convolution operation is performed in parallel through multiple feature channels to generate multiple sets of linear activation responses, and the ReLU function values ​​of the linear activation responses are calculated to obtain the feature layer of the convolution operation.

[0172] Step 2-7: Let c←c+1, conv←conv+1, channels←channels×2.

[0173] In step 2-8, execute step 2-6, then go to step 2-9.

[0174] Step 2-9 Let c←c+1, conv←conv+1, map c =map c1 +map c2 , where map c1 and map c2 They are x k,i and ▽d(y i ,H(x k,i ))Feature map after conv convolution.

[0175] Step 2-10 halves the number of feature channels, i.e., channels←channels / 2.

[0176] Step 2-11 executes step 2-6, and then executes step 2-12.

[0177] In step 2-12, let c←c+1, conv←conv+1, channels=1.

[0178] Step 2-13 executes step 2-6, and then executes step 2-14.

[0179] Step 2-14 Let c←c+1, conv←conv+1, map c ←x k,i +map c , then calculate the ReLU function value of the linear activation response to obtain the feature layer of the convolution operation.

[0180] Step 2-15: Let k←k+1. If k<2, execute step 2-16; otherwise, return to step 2-3.

[0181] Step 2-16 outputs a predicted image of M samples, where c=9.

[0182] Step 2-17 Initialize parameters: c = 10, pool = 0, channels u =64, upsample=0. Where pool is the number of pooling operations, and upsample is the number of upsampling operations.

[0183] Step 2-18 Input the output of the DGD network from the contraction path into U-net as the feature map mapu c , namely mapu c It is the feature layer generated by the input image after the c-1th convolution, pooling and upsampling.

[0184] Step 2-19: Let U-net have the number of convolution operations conv u =0.

[0185] Step 2-20 for mapu c The convolution operation is performed in parallel through multiple feature channels to generate multiple sets of linear activation responses. The ReLU function value of the linear activation response is then calculated to obtain the feature layer of the convolution operation.

[0186] Step 2-21 Let c←c+1, conv u ←conv u +1. If conv u ≥2, then execute steps 2-22; if conv u <2, go to step 2-20.

[0187] Step 2-22: channels u ←channels u ×2, for mapu c Perform a pooling operation and mark the feature layer before the pool operation as map pool .

[0188] Step 2-23: Let c←c+1, pool←pool+1. If pool≥4, let pool=3 and execute step 2-24; otherwise, go to step 2-19.

[0189] Step 2-24: Let conv u =0.

[0190] Step 2-25 for mapu c The convolution operation is performed in parallel through multiple feature channels to generate multiple sets of linear activation responses. The ReLU function value of the linear activation response is then calculated to obtain the feature layer of the convolution operation.

[0191] Step 2-26 Let c←c+1, conv u ←conv u +1. If conv u ≥2, let channels u ←channels u / 2, execute step 2-27; otherwise, go to step 2-25.

[0192] Step 2-27 for mapu c Perform an upsampling operation to obtain the upsampled feature layer.

[0193] Step 2-28 Map pool Perform the splicing operation, map pool Correspondingly stacked on the feature layer obtained by the upsampling operation.

[0194] Step 2-29 Let c←c+1, pool←pool–1, upsample←upsample+1. If upsample≥4, let conv u =0, execute step 2-30; otherwise, go to step 2-24.

[0195] Step 2-30 for mapu c The convolution operation is performed in parallel through multiple feature channels to generate multiple sets of linear activation responses. The ReLU function value of the linear activation response is then calculated to obtain the feature layer of the convolution operation.

[0196] Step 2-31 Let c←c+1, conv u ←conv u +1. If conv u ≥2, let channels u =1, go to step 2-32; otherwise, go to step 2-30.

[0197] Step 2-32 for mapu cPerform a 1×1 convolution operation to generate a linear activation response. Then use the ReLU function to calculate the function value corresponding to the linear activation response to obtain the optimized image. At this time, c=37.

[0198] 3. Optimize network parameters layer by layer

[0199] like Figure 5-6 , the specific steps of RAR-Net training process are as follows:

[0200] In step 3-1, let k=0 and the cumulative number of corrections of the network parameter gradient be t=0.

[0201] In step 3-2, let epoch=0 for traversing all small batch training sets in training set A, and let the index I of the small batch training set in A=1.

[0202] Step 3-3 respectively convert the first group of small batch training set A in training set A I The M samples in are input into the built RAR-Net for forward propagation.

[0203] Step 3-4 Set the network's loss function

[0204] The loss function of the DGD network is:

[0205]

[0206] The loss function of U-net is:

[0207]

[0208] The loss function of RAR-Net is:

[0209] L(W G ,b G ,W U ,b U )=αL G (W G ,b G )+L U (W U ,b U ) (15)

[0210] in:

[0211] x k+1,i =G θk (▽d(y i ,H(x k,i )),x k,i ) (16)

[0212] ▽d(y i ,H(x k,i ))=H* [H(x k,i )-y i ] (17)

[0213] Where, L G (W G ,b G ) is the loss function of the DGD network, L U (W U ,b U ) is the loss function of U-net, L(W G ,b G ,W U ,b U ) is the loss function of the entire RAR-Net, W G and b G are the weight parameter set and bias parameter set of the DGD network, W U and b U They are the weight parameter set and bias parameter set of U-net, x k+1,i and x k,i are the light absorption energy distribution diagrams outputted at the k+1th and kth iterations, respectively, x true,i is the simulated light absorption energy distribution diagram, θ k is the parameter of the kth network structure unit, are the same structure with different learning parameters θ k RAR-Net, y i is the photoacoustic signal matrix measured by the detector, x i is the light absorption energy distribution map output by U-net, H is the forward imaging operator (i.e., the photoacoustic signal generated by the tissue is obtained from the light absorption energy density), and H * is the adjoint operator of H, and ▽d(·) is the data fitting information.

[0214] Since RAR-Net aims to learn more accurate domain transformation rather than relying on image post-processing, the scale factor α is set (α>1), and the specific value of α is determined by the convergence speed and reconstruction quality.

[0215] Step 3-5 According to the recursive relation Z c =W c X c-1 +b c , respectively calculate the gradients of the parameters of the DGD network output layer (i.e., the 9th layer) and the U-net output layer (i.e., the 37th layer).

[0216] The gradients of the parameters of the DGD network output layer (i.e., the 9th layer) are:

[0217]

[0218] In the formula, Z G,9 is the linear activation response generated by the 9th layer in the forward propagation of the DGD network, W G,9 is the weight parameter matrix of the 9th layer of the DGD network, b G,9 is the bias parameter vector of the 9th layer of the DGD network, X8 is the output of the 8th layer of the DGD network in the forward propagation, “⊙” represents the Hadamard product of the matrix, σ'(·) is the first-order derivative of the ReLU function, and the superscript “T” means to find the transpose of the matrix.

[0219] Calculate the gradient of each parameter of the U-net output layer (i.e. layer 37):

[0220]

[0221] Among them, Z U,37 is the linear activation response generated by the 37th layer in the forward propagation of U-net, W U,37 is the weight parameter matrix of the 37th layer of U-net, b U,37 is the bias parameter vector of the 37th layer of U-net, M is the sample size of the current training set, and X 36 It is the output of the 36th layer network in the U-net forward propagation.

[0222] Steps 3-6 calculate the error vectors of each layer of the DGD network and U-net respectively.

[0223] The error vector of the DGD network output layer (i.e., layer 9) is:

[0224]

[0225] In the formula, δ G,9 is the error vector of the output layer (i.e., layer 9) of the DGD network. G,9 The error vector from the 8th layer to the 2nd layer is obtained recursively:

[0226] δ G,c =(W G,c+1 ) T δ G,c+1 ⊙σ'(Z G,c ) (twenty one)

[0227] Among them, δ G,c and δ G,c+1 They are the error vectors of the cth and c+1th layers of the DGD network, c=8,7,...,3,2.

[0228] Calculate the error vector of the U-net output layer (i.e. layer 37):

[0229]

[0230] In the formula, δ U,37 is the error vector of the U-net output layer (i.e., layer 37). According to δ U,37 The error vectors from the 37th layer to the 10th layer are obtained recursively:

[0231] δ U,c =(W U,c+1 ) T δ U,c+1 ⊙σ'(Z U,c ) (twenty three)

[0232] Among them, δ U,c and δ U,c+1 They are the error vectors of the c-th layer and the c+1-th layer of U-net, c=36,35,...,11,10.

[0233] Steps 3-7 calculate the gradient of the parameters of each layer of RAR-Net:

[0234]

[0235] Among them, W c 、b c , Z c and δ c are the weight parameter matrix, bias parameter vector, linear activation response and error vector of the cth layer of RAR-Net. When c=2,3,……,9, W c =W G 、b c =b G , Z c =Z G ,δ c =δ G,c When c=10,11,......,37,W c =W U 、b c =b U , Z c =Z U ,δ c =δ U,c , X c-1 is the output of the c-1th layer network in the forward propagation.

[0236] Steps 3-8 use the Adam algorithm to optimize the network, and the updated network gradient is:

[0237]

[0238]

[0239] Here, μ is the learning rate, which is set to 0.003 in the present invention.

[0240] Step 3-9 Let \(I\leftarrow I + 1\). If \(I\lt N + 1\), then go to Step 3-2; otherwise, let \(epoch\leftarrow epoch+1\), and execute Step 3-10.

[0241] Step 3-10 Determine whether to stop traversing: If \(epoch\lt epoch\) max , \(epoch\) max is the maximum value of \(epoch\), then let \(epoch\leftarrow epoch + 1\), and go to Step 3-2; otherwise, execute Step 3-11.

[0242] Step 3-11 Let \(k\leftarrow k + 1\). If \(k\lt2\), then return to Step 3-2; otherwise, execute Step 3-12.

[0243] Step 3-12 Determine the final RAR-Net model.

[0244] IV. Reconstruct the light absorption energy distribution map using the optimized RAR-Net

[0245] Input the measured photoacoustic signal data containing reflection clutter actually collected by the ultrasonic detector into the trained RAR-Net, and output a high-quality light absorption energy distribution image.

[0246] In addition, the present invention also provides a system for removing reflection artifacts in endoscopic photoacoustic tomography images, including:

[0247] A photoacoustic data group acquisition module, configured to acquire a photoacoustic data group of a plane to be imaged;

[0248] A reflection artifact removal imaging module, configured to input the photoacoustic data group into a reflection artifact removal model to obtain a light absorption energy distribution map with reflection artifacts removed at the plane to be imaged; the artifact removal model is obtained by training the RAR-Net neural network according to historical photoacoustic data groups; the RAR-Net neural network is obtained by taking the output of the deep gradient descent network as the input of the deep learning segmentation network.

[0249] A training data acquisition module, configured to acquire historical photoacoustic data groups of multiple imaging planes and historical light absorption energy distribution maps with reflection artifacts removed from multiple imaging planes as training data;

[0250] A training set division module, configured to divide the training data into multiple training sets;

[0251] A reflection artifact removal model determination module, configured to use multiple training sets, with the historical photoacoustic data group as the input and the historical light absorption energy distribution map with reflection artifacts removed as the expected output, to train the RAR-Net neural network to obtain a reflection artifact removal model.

[0252] Among them, the reflection artifact removal model determination module specifically includes:

[0253] An initial neural network determination unit, used to determine that the RAR-Net neural network before training is a k×epoch×I-th order RAR-Net neural network;

[0254] An iteration round assignment unit, used to set the iteration round k=0;

[0255] The training set traversal times assignment unit is used to set the training set traversal times epoch = 0;

[0256] A current training set determining unit, used to determine any training set as the current training set;

[0257] An index assignment unit, used to set the index I of the current training set to 1;

[0258] A neural network updating unit is used to input the current training set into the k×epoch×Ith order RAR-Net neural network, adjust the parameters of the k×epoch×Ith order RAR-Net neural network, and obtain the k×epoch×(I+1)th order RAR-Net neural network;

[0259] The current training set update unit is used to update the current training set, increase the value of I by 1, and call the neural network update unit until the training set is traversed to obtain the k×(epoch+1)×1th order RAR-Net neural network;

[0260] The training set traversal times updating unit is used to increase the value of epoch by 1 and call the index assignment unit until epoch is greater than the training set traversal times threshold, thus obtaining the (k+1)×0×1th order RAR-Net neural network;

[0261] The iteration round updating unit is used to increase the value of k by 1 and return to the step to call the training set traversal times assignment unit until k is greater than the iteration round threshold, and obtain a k×epoch×I-order RAR-Net neural network as a reflection artifact removal model.

[0262] Specifically, the neural network update unit includes:

[0263] A light absorption energy history prediction distribution map determining subunit is used to input the historical photoacoustic data group in the current training set into the k×epoch×I-th order RAR-Net neural network to obtain the light absorption energy history prediction distribution map at the k×epoch×I-th iteration;

[0264] a loss determination subunit, for determining the loss of the deep gradient descent network as a first loss and determining the loss of the deep learning segmentation network as a second loss according to the historical light absorption energy distribution map and the historical predicted light absorption energy distribution map at the k×epoch×Ith iteration;

[0265] A first linear activation response determination subunit, used to determine that the linear activation response of the output layer of the deep gradient descent network is a first linear activation response;

[0266] A second linear activation response determining subunit, used to determine that the linear activation response of the output layer of the deep learning segmentation network is a second linear activation response;

[0267] A first output error vector determining subunit, configured to determine an error vector of an output layer of a deep gradient descent network as a first output error vector according to a first linear activation response;

[0268] A second output error vector determining subunit, configured to determine an error vector of an output layer of the deep learning segmentation network as a second output error vector according to a second linear activation response;

[0269] A first error vector group determining subunit, configured to determine, according to the first loss and the first output error vector, the error vector of each layer of the deep gradient descent network as a first error vector group using a first recursive formula;

[0270] A second error vector group determining subunit is used to determine the error vector of each layer of the deep learning segmentation network as a second error vector group using a second recursive formula according to the second loss and the second output error vector;

[0271] A gradient determination subunit, for determining the gradient of each layer of a k×epoch×I-order RAR-Net neural network according to the first error vector group and the second error vector group;

[0272] The gradient updating subunit is used to update multiple gradients according to the gradient of each layer of the k×epoch×I-th order RAR-Net neural network using the Adam algorithm to obtain the k×epoch×(I+1)-th order RAR-Net neural network.

[0273] in,

[0274] The first loss is:

[0275] The second loss is:

[0276] The first output error vector is:

[0277] The second output error vector is:

[0278] The first recursive formula is: G,c =(W G,c+1 ) T δ G,c+1 ⊙σ'(Z G,c );

[0279] The second recursive formula is: U,c =(W U,c+1 ) T δ U,c+1 ⊙σ'(Z U,c );

[0280] The gradient of each layer of the RAR-Net neural network is:

[0281] Among them, L G (W G ,b G ) and L U (W U ,b U ) are the first loss and the second loss respectively; W G and b G are the weight parameter set and bias parameter set of the deep gradient descent network respectively; W U and b U are the weight parameter set and bias parameter set of the deep learning segmentation network respectively; i is the sample number in the training set, i = (1, 2, ..., M); M is the sample capacity of the current training set; x true,i is the historical light absorption energy distribution diagram; x k,i and x k+1,i are the predicted distribution diagrams of light absorption energy history obtained at the kth and k+1th iterations respectively; θ k is the parameter of the kth network structure unit; are the same structure with different learning parameters θ k RAR-Net i is the historical photoacoustic data set; x i is the predicted distribution map of light absorption energy history output by the deep learning segmentation network; H(x k,i ) is the same as x k,i The operator concerned; δ G,c ,δ U,c ,δ G,c+1 and δ U,c+1 are the error vectors of the cth layer and the c+1th layer of the RAR-Net neural network, respectively. The subscript G represents the deep gradient descent network, and the subscript U represents the deep learning segmentation network. σ'(·) represents the first-order derivative of the ReLU function, and Z G,9 and ZU,37 are the first linear activation response and the second linear activation response respectively; W G,c+1 and W U,c+1 are the weight parameter matrices of the cth layer of the RAR-Net neural network; Z G,c and Z U,c are the linear activation responses of the cth layer of the RAR-Net neural network; L(W G ,b G ,W U ,b U ) is the total loss of the RAR-Net neural network, L(W G ,b G ,W U ,b U )=αL G (W G ,b G )+L U (W U ,b U ), α is the scale factor: determined by the convergence speed and reconstruction quality; W c 、b c , Z c and δ c are the weight parameter matrix, bias parameter vector, linear activation response and error vector of the cth layer of the RAR-Net neural network; in the deep gradient descent network, that is, when c=1,2,3,…,9, W c =W G 、b c =b G , Z c =Z G ,δ c =δ G,c ; In the deep learning segmentation network, that is, when c=10,11,...,37, W c =W U 、b c =b U , Z c =Z U ,δ c =δ U,c ;X c-1 is the output of the c-1th layer network in the forward propagation.

[0282] The various embodiments in this specification are described in a progressive manner, and each embodiment focuses on the differences from other embodiments. The same or similar parts between the various embodiments can be referenced to each other.

[0283] This article uses specific examples to illustrate the principles and implementation methods of the present invention. The above examples are only used to help understand the method and core ideas of the present invention. At the same time, for those skilled in the art, according to the ideas of the present invention, there will be changes in the specific implementation methods and application scope. In summary, the content of this specification should not be understood as limiting the present invention.

Claims

1. A method for removing reflection artifacts in endoscopic photoacoustic tomography images, characterized in that: The method comprises: Acquiring a photoacoustic data set of a plane to be imaged; The photoacoustic data set is input into a reflection artifact removal model to obtain a light absorption energy distribution diagram after the reflection artifact is removed at the imaging plane; the artifact removal model is obtained by training a RAR-Net neural network based on the historical photoacoustic data set; the RAR-Net neural network is obtained by using the output of a deep gradient descent network as the input of a deep learning segmentation network; The training process of the RAR-Net neural network includes: inputting the current training set into the k×epoch×Ith order RAR-Net neural network, adjusting the parameters of the k×epoch×Ith order RAR-Net neural network, and obtaining the k×epoch×(I+1)th order RAR-Net neural network; k is the iteration round; epoch is the number of training set traversals; I is the index of the current training set; Input the current training set into the k×epoch×Ith order RAR-Net neural network, adjust the parameters of the k×epoch×Ith order RAR-Net neural network, and obtain the k×epoch×(I+1)th order RAR-Net neural network, specifically including: Inputting the historical photoacoustic data set in the current training set into the k×epoch×I-th order RAR-Net neural network to obtain a historical prediction distribution diagram of light absorption energy at the k×epoch×I-th iteration; According to the historical light absorption energy distribution map and the historical predicted light absorption energy distribution map at the k×epoch×Ith iteration, the loss of the deep gradient descent network is determined as the first loss and the loss of the deep learning segmentation network is determined as the second loss; Determine the linear activation response of the output layer of the deep gradient descent network as a first linear activation response; Determine the linear activation response of the output layer of the deep learning segmentation network as the second linear activation response; Determine an error vector of an output layer of a deep gradient descent network as a first output error vector according to the first linear activation response; Determine, according to the second linear activation response, an error vector of an output layer of the deep learning segmentation network as a second output error vector; According to the first loss and the first output error vector, using a first recursive formula to determine the error vector of each layer of the deep gradient descent network as a first error vector group; According to the second loss and the second output error vector, using a second recursive formula to determine the error vector of each layer of the deep learning segmentation network as a second error vector group; Determine the gradient of each layer of the k×epoch×I-order RAR-Net neural network according to the first error vector group and the second error vector group; According to the gradient of each layer of the k×epoch×I-th order RAR-Net neural network, the Adam algorithm is used to update the multiple gradients to obtain the k×epoch×(I+1)-th order RAR-Net neural network; in, The first loss is: The second loss is: The first output error vector is: The second output error vector is: The first recursive formula is: G,c =(W G,c+1 ) T δ G,c+1 ⊙σ'(Z G,c ); The second recursive formula is: U,c =(W U,c+1 ) T δ U,c+1 ⊙σ'(Z U,c ); The gradient of each layer of the RAR-Net neural network is: Among them, L G (W G ,b G ) and L U (W U ,b U ) are the first loss and the second loss respectively; W G and b G are the weight parameter set and bias parameter set of the deep gradient descent network respectively; W U and b U are the weight parameter set and bias parameter set of the deep learning segmentation network respectively; i is the sample number in the training set, i = (1, 2, ..., M); M is the sample capacity of the current training set; x true,i is the historical light absorption energy distribution diagram; x k,i and x k+1,i are the predicted distribution diagrams of light absorption energy history obtained at the kth and k+1th iterations respectively; θ k is the parameter of the kth network structure unit; are the same structure with different learning parameters θ k RAR-Net i is the historical photoacoustic data set; x i is the predicted distribution map of light absorption energy history output by the deep learning segmentation network; H(x k,i ) is the same as x k,i The operator concerned; δ G,c , δ U,c , δ G,c+1 and δ U,c+1 are the error vectors of the cth layer and the c+1th layer of the RAR-Net neural network, respectively, where the subscript G represents a deep gradient descent network and the subscript U represents a deep learning segmentation network; δ G,9 and δ U,37 are the first output error vector and the second output error vector respectively; σ'(·) represents the first-order derivative of the ReLU function, Z G,9 and Z U,37 are the first linear activation response and the second linear activation response respectively; W G,c+1 and W U,c+1 are the weight parameter matrices of the c+1th layer of the RAR-Net neural network; Z G,c and Z U,c are the linear activation responses of the cth layer of the RAR-Net neural network; L(W G ,b G ,W U ,b U ) is the total loss of the RAR-Net neural network, L(W G ,b G ,W U ,b U )=αL G (W G ,b G )+L U (W U ,b U ), α is the scale factor: determined by the convergence speed and reconstruction quality; W c 、b c , Z c and δ c are the weight parameter matrix, bias parameter vector, linear activation response and error vector of the cth layer of the RAR-Net neural network; in the deep gradient descent network, that is, when c=1,2,3,…,9, W c =W G 、b c =b G , Z c =Z G , δ c =δ G,c ; In the deep learning segmentation network, that is, when c=10,11,...,37, W c =W U 、b c =b U , Z c =Z U , δ c =δ U,c ;X c-1 is the output of the c-1th layer network in the forward propagation.

2. The method for removing reflection artifacts in endoscopic photoacoustic tomography images according to claim 1, characterized in that: Before acquiring the photoacoustic signal of the plane to be imaged, the method further includes: Acquire historical photoacoustic data groups of multiple imaging planes and historical light absorption energy distribution diagrams of multiple imaging planes after removing reflection artifacts as training data; Dividing the training data into a plurality of training sets; Using multiple training sets, the RAR-Net neural network was trained with historical photoacoustic data sets as input and historical light absorption energy distribution maps after removing reflection artifacts as the expected output to obtain a reflection artifact removal model.

3. The method for removing reflection artifacts in endoscopic photoacoustic tomography images according to claim 2, characterized in that: The method uses multiple training sets, takes the historical photoacoustic data set as input, takes the historical light absorption energy distribution map after removing the reflection artifacts as the expected output, trains the RAR-Net neural network, and obtains the reflection artifact removal model, specifically including: Determine that the RAR-Net neural network before training is the k×epoch×Ith order RAR-Net neural network; Let the iteration round k = 0; Set the number of training set traversal epochs to 0; Determine any training set as the current training set; Let the index of the current training set I = 1; Inputting the current training set into a k×epoch×Ith order RAR-Net neural network, adjusting the parameters of the k×epoch×Ith order RAR-Net neural network, and obtaining a k×epoch×(I+1)th order RAR-Net neural network; Update the current training set, increase the value of I by 1, and return to the step "input the current training set into the k×epoch×Ith order RAR-Net neural network, adjust the parameters of the k×epoch×Ith order RAR-Net neural network, and obtain the k×epoch×(I+1)th order RAR-Net neural network" until the training set is traversed to obtain the k×(epoch+1)th order RAR-Net neural network; Let the value of epoch increase by 1 and return to step "set the index of the current training set I = 1" until epoch is greater than the threshold of the number of training set traversal times, and the (k+1)×0×1th order RAR-Net neural network is obtained; Increase the value of k by 1 and return to step "set the number of training set traversal epochs = 1" until k is greater than the iteration round threshold, and obtain a k×epoch×I-order RAR-Net neural network as a reflection artifact removal model.

4. A system for removing reflection artifacts in endoscopic photoacoustic tomography images, characterized in that: The system comprises: A photoacoustic data set acquisition module, used to acquire a photoacoustic data set of a plane to be imaged; A reflection artifact removal imaging module is used to input the photoacoustic data set into a reflection artifact removal model to obtain a light absorption energy distribution diagram after removing reflection artifacts at the imaging plane; the artifact removal model is obtained by training a RAR-Net neural network based on the historical photoacoustic data set; the RAR-Net neural network uses the output of a deep gradient descent network as an input to a deep learning segmentation network; The neural network updating unit in the reflection artifact removal model determination module specifically includes: The light absorption energy history prediction distribution map determining subunit is used to input the historical photoacoustic data group in the current training set into the k×epoch×I-th order RAR-Net neural network to obtain the light absorption energy history prediction distribution map at the k×epoch×I-th iteration; k is the iteration round; epoch is the number of training set traversals; I is the index of the current training set; a loss determination subunit, for determining the loss of the deep gradient descent network as a first loss and determining the loss of the deep learning segmentation network as a second loss according to the historical light absorption energy distribution map and the historical predicted light absorption energy distribution map at the k×epoch×Ith iteration; A first linear activation response determination subunit, used to determine that the linear activation response of the output layer of the deep gradient descent network is a first linear activation response; A second linear activation response determining subunit, used to determine that the linear activation response of the output layer of the deep learning segmentation network is a second linear activation response; A first output error vector determining subunit, configured to determine an error vector of an output layer of a deep gradient descent network as a first output error vector according to the first linear activation response; A second output error vector determining subunit, configured to determine, according to the second linear activation response, an error vector of an output layer of the deep learning segmentation network as a second output error vector; A first error vector group determining subunit, configured to determine, according to the first loss and the first output error vector, the error vector of each layer of the deep gradient descent network as a first error vector group using a first recursive formula; A second error vector group determining subunit, configured to determine, according to the second loss and the second output error vector, the error vector of each layer of the deep learning segmentation network as a second error vector group using a second recursive formula; A gradient determination subunit, configured to determine the gradient of each layer of the k×epoch×I-order RAR-Net neural network according to the first error vector group and the second error vector group; A gradient updating subunit, used to update a plurality of the gradients according to the gradient of each layer of the k×epoch×I-th order RAR-Net neural network using the Adam algorithm to obtain the k×epoch×(I+1)-th order RAR-Net neural network; The first loss is: The second loss is: The first output error vector is: The second output error vector is: The first recursive formula is: G,c =(W G,c+1 ) T δ G,c+1 ⊙σ'(Z G,c ); The second recursive formula is: U,c =(W U,c+1 ) T δ U,c+1 ⊙σ'(Z U,c ); The gradient of each layer of the RAR-Net neural network is: Among them, L G (W G ,b G ) and L U (W U ,b U ) are the first loss and the second loss respectively; W G and b G are the weight parameter set and bias parameter set of the deep gradient descent network respectively; W U and b U are the weight parameter set and bias parameter set of the deep learning segmentation network respectively; i is the sample number in the training set, i = (1, 2, ..., M); M is the sample capacity of the current training set; x true,i is the historical light absorption energy distribution diagram; x k,i and x k+1,i are the predicted distribution diagrams of light absorption energy history obtained at the kth and k+1th iterations respectively; θ k is the parameter of the kth network structure unit; are the same structure with different learning parameters θ k RAR-Net i is the historical photoacoustic data set; x i is the predicted distribution map of light absorption energy history output by the deep learning segmentation network; H(x k,i ) is the same as x k,i The operator concerned; δ G,c , δ U,c , δ G,c+1 and δ U,c+1 are the error vectors of the cth layer and the c+1th layer of the RAR-Net neural network, respectively, where the subscript G represents a deep gradient descent network and the subscript U represents a deep learning segmentation network; δ G,9 and δ U,37 are the first output error vector and the second output error vector respectively; σ'(·) represents the first-order derivative of the ReLU function, Z G,9 and Z U,37 are the first linear activation response and the second linear activation response respectively; W G,c+1 and W U,c+1 are the weight parameter matrices of the cth layer of the RAR-Net neural network; Z G,c and Z U,c are the linear activation responses of the cth layer of the RAR-Net neural network; L(W G ,b G ,W U ,b U ) is the total loss of the RAR-Net neural network, L(W G ,b G ,W U ,b U )=αL G (W G ,b G )+L U (W U ,b U ), α is the scale factor, which is determined by the convergence speed and reconstruction quality; W c 、b c , Z c and δ c are the weight parameter matrix, bias parameter vector, linear activation response and error vector of the cth layer of the RAR-Net neural network; in the deep gradient descent network, that is, when c=1,2,3,…,9, W c =W G 、b c =b G , Z c =Z G , δ c =δ G,c ; In the deep learning segmentation network, that is, when c=10,11,...,37, W c =W U 、b c =b U , Z c =Z U , δ c =δ U,c ;X c-1 is the output of the c-1th layer network in the forward propagation.

5. The system for removing reflection artifacts in endoscopic photoacoustic tomography images according to claim 4, characterized in that: The system further comprises: A training data acquisition module, used to acquire historical photoacoustic data groups of multiple imaging planes and historical light absorption energy distribution diagrams of multiple imaging planes after removing reflection artifacts as training data; A training set division module, used for dividing the training data into multiple training sets; The reflection artifact removal model determination module is used to use multiple training sets, take the historical photoacoustic data group as input, and take the historical light absorption energy distribution map after removing the reflection artifact as the expected output, to train the RAR-Net neural network and obtain the reflection artifact removal model.

6. The system for removing reflection artifacts in endoscopic photoacoustic tomography images according to claim 5, characterized in that: The reflection artifact removal model determination module specifically includes: An initial neural network determination unit, used to determine that the RAR-Net neural network before training is a k×epoch×I-th order RAR-Net neural network; An iteration round assignment unit, used to set the iteration round k=0; The training set traversal times assignment unit is used to set the training set traversal times epoch = 0; A current training set determination unit, used to determine any training set as a current training set; An index assignment unit, used to set the index I of the current training set to 1; A neural network updating unit, used for inputting the current training set into a k×epoch×I-th order RAR-Net neural network, adjusting the parameters of the k×epoch×I-th order RAR-Net neural network, and obtaining a k×epoch×(I+1)-th order RAR-Net neural network; The current training set updating unit is used to update the current training set, increase the value of I by 1, and call the neural network updating unit until the training set is traversed to obtain the k×(epoch+1)×1th order RAR-Net neural network; A training set traversal times updating unit is used to increase the value of epoch by 1 and call the index assignment unit until epoch is greater than the training set traversal times threshold, thereby obtaining a (k+1)×0×1th order RAR-Net neural network; The iteration round updating unit is used to increase the value of k by 1 and return to the step of calling the training set traversal times assignment unit until k is greater than the iteration round threshold, and obtain a k×epoch×I-order RAR-Net neural network as a reflection artifact removal model.

Images