Photon counting CT image noise reduction method based on multi-channel 3D U-net

Abstract

The invention discloses a photon counting CT image noise reduction method based on multi-channel 3D U-net, relates to the technical field of medical image processing, and aims to solve the problems of insufficient spatial information utilization of a traditional 2D CNN, complex calculation of a 3DCNN and poor noise reduction effect of a traditional method. According to the method, high-energy, low-energy and all-energy images obtained by photon counting CT are used as three-dimensional multi-channel input, and a matrix containing space and energy spectrum dimensions is formed after preprocessing. Multi-scale three-dimensional features are extracted through a 3D U-net encoder, deep semantic information is aggregated through a bottleneck layer, a decoder carries out deconvolution up-sampling and is fused with an encoder feature map through jump connection, and meanwhile, an attention module is embedded in the jump connection to generate an attention mask so as to strengthen key area features; a mixed loss function training network including mean square error loss and structural similarity loss is adopted, and pixel-level precision and structural retention are balanced.

Description

Technical Field

[0001] The present invention relates to the field of medical image processing technology, and in particular to a photon counting CT image denoising method based on a multi-channel 3D U-net. Background Art

[0002] Traditional 2D convolutional neural networks (CNNs) have been widely used in CT image denoising. However, because 2D CNNs can only process single slices and cannot fully utilize three-dimensional spatial information, they have certain limitations in terms of detail preservation and structural consistency in the denoised images. 3D convolutional neural networks (3D CNNs), a recent breakthrough technology, have demonstrated strong spatial information learning capabilities in medical image processing, but they still face challenges such as high network training complexity and high computational resource requirements.

[0003] Medical CT images are usually stored in the form of three-dimensional volume data. 2D CNN can only process a single slice and cannot capture inter-layer information, resulting in poor consistency of the denoised image in the Z-axis direction.

[0004] Traditional noise reduction methods have shortcomings. For example, BM3D and TV filtering can easily introduce blurring effects while removing noise, resulting in edge loss.

[0005] Existing methods do not have the ability to combine multi-energy spectral image analysis because traditional CT cannot distinguish multiple energies at once and has poor noise reduction capabilities.

[0006] Based on the above problems, the present invention proposes a photon counting CT image denoising method based on multi-channel 3D U-net. Summary of the Invention

[0007] The purpose of the present invention is to provide a photon counting CT image denoising method based on multi-channel 3D U-net to solve the problems raised in the prior art.

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

[0009] The photon counting CT image denoising method based on multi-channel 3D U-net includes the following steps:

[0010] S1. Use high-energy, low-energy, and full-energy images acquired by photon counting CT as three-dimensional multi-channel input data to form an image matrix containing spatial dimensions and energy spectrum channels, providing multi-spectral raw data for subsequent three-dimensional noise reduction processing;

[0011] S2, the input data is processed layer by layer through the 3D U-net encoder, and multi-scale three-dimensional spatial features are extracted using 3D convolutional layers, batch normalization layers and ReLU activation functions, and the spatial dimension is gradually reduced in combination with the maximum pooling operation;

[0012] S3, the bottleneck layer at the end of the encoder aggregates deep features to extract abstract features representing the core structure of the CT image, providing a high-level semantic information foundation for the decoder to restore the image;

[0013] S4, upsampling the bottleneck layer abstract features through the decoder's deconvolution operation, gradually restoring the image's spatial dimensions, and completing the initial conversion from high-level features to image dimensions;

[0014] S5. Concatenate the high-resolution feature maps of the corresponding layers in the encoder with the upsampled feature maps of the decoder through skip connections. Before concatenation, align the spatial dimensions through padding operations to fuse high- and low-level features.

[0015] S6. Embed an attention module in the skip connection path, perform linear transformation, nonlinear activation and normalization on the input feature map, generate an attention mask and multiply it element-wise with the original feature map;

[0016] S7. A hybrid loss function consisting of structural similarity loss and mean square error loss is used to train the network. The parameters are adjusted through back propagation, and the high-energy, low-energy, and full-energy denoised images corresponding to the input channels are finally output.

[0017] 9. The photon counting CT image denoising method based on multi-channel 3D U-net according to claim 1, wherein said S1 further comprises the following contents:

[0018] S1.1: Use photon counting CT equipment to obtain high-energy H, low-energy L, and full-energy T images of the same scan object, denoted as I H , I L , I T , all three are three-dimensional volume data, with a unified spatial size of H×W×D; H is height, W is width, and D is depth; the original image is then preprocessed, including removing detector background noise, correcting energy response differences, and normalizing pixel values ​​to the interval [0,1];

[0019] S1.2: The pre-processed high-energy, low-energy, and full-energy images are superimposed according to the channel dimension to form an input matrix X with dimensions H×W×D×3, where the first three dimensions correspond to the spatial dimensions and the fourth dimension is the energy spectrum channel in the order of high energy, low energy, and full energy, i.e., X=[I H , I L , I T ], as the input data of the 3D U-net network, enables the network to simultaneously learn the noise distribution characteristics of different energy spectrum channels, and provide a multi-dimensional data foundation for subsequent three-dimensional feature learning.

[0020] Said S2 further comprises the following contents:

[0021] S2.1: The input matrix X is subjected to layer-by-layer feature extraction through a 3D U-net encoder. Each layer first performs a 3×3×3 convolution operation (Convolution 3×3) with a stride of 1 and padding of 1 to maintain the spatial dimension. Batch normalization (BN) is then used to reduce internal covariate shift. The ReLU activation function is then used to introduce a nonlinear transformation (Convolution+ReLU+BN) module, forming the basic processing unit of the encoder. The details are as follows:

[0022] For 3D convolution operations, a convolution kernel of size 3×3×3 is used, with a stride of 1 and a padding of 1 to ensure that the spatial dimension remains unchanged. The calculation formula is as follows:

[0023]

[0024] Where W (l) is the convolution kernel of the lth layer, b (l) is the bias term, l represents the level of the network, F (l) Represents the output feature map after the l-th layer convolution operation in the 3D U-net encoder. i, j, and k correspond to the position index of the convolution kernel in the three spatial dimensions of depth D, height H, and width W, respectively, and the value range is 0, 1, and 2; X (l-1) Represents the output of the (l-1)th layer; (x, y, z) is the feature map coordinate, and three-dimensional local features are extracted through convolution operations across space and channels;

[0025] For batch normalization, the convolution output is normalized to reduce the internal covariate offset. The formula is as follows:

[0026]

[0027] in, μ is the feature map after batch normalization of convolution output; β is the batch mean, is the batch variance, γ (l) and β (l) is a learnable parameter, ε=10 -5 A small constant to prevent division by zero;

[0028] The ReLU activation function introduces a nonlinear function through a linear rectification function, and performs a nonlinear transformation on the feature map after batch normalization. The calculation formula is as follows:

[0029]

[0030] Among them, A (l) is the feature map after activation;

[0031] S2.2: After each convolution activation, the feature map is downsampled through a 2×2×2 max pooling layer with a stride of 2. The calculation formula is:

[0032]

[0033] P (l) Represents the feature map after pooling of the first layer; 2x, 2y, 2z represent the starting coordinates of the pooling window, the stride is 2, the windows do not overlap, and the output size is halved;

[0034] Then the feature map size is halved to extract multi-scale features: shallow features retain edge and texture details, and deep features capture the global structure of tissues and organs to form encoder features {P1, P2, P3, P bot}, where P1, P2, and P3 are feature maps generated by the encoder through multi-layer maximum pooling operations, corresponding to different downsampling stages. The subscripts 1, 2, and 3 in the numbers represent the pooling levels. bot is the deepest feature map before the bottleneck layer.

[0035] The S3 further includes the following:

[0036] S3.1: After the encoder undergoes multi-layer 3D convolution and maximum pooling processing, the deepest feature map is formed at the end, denoted as P pre-bot , dimension is H pre-bot ×W pre-bot ×D pre-bot ×C pre-bot ;in, k is the total number of encoder layers, C pre-bot is the number of channel layers; the bottleneck layer performs continuous 3D convolution processing, including 2-3 convolution layers, and the configuration of each layer is consistent with the encoder to ensure the spatial dimension H pre-bot ×W pre-bot ×D pre-bot The number of channels will gradually increase to aggregate deep semantic information.

[0037] S3.2: The bottleneck layer performs convolution on the encoder’s last pooled output and outputs an abstract feature map F. core , whose dimension is H bot ×W bot ×D bot ×C core ; Among them, C core >C pre-bot , C core Indicates the number of channels in the bottleneck layer; the feature map contains the core semantic information of the CT image, which is different from the local details such as edges and textures in the shallow layer of the encoder.

[0038] Said S4 further comprises the following contents:

[0039] S4.1: The decoder receives the abstract feature map F output by the bottleneck layer core , and restore its spatial dimension to 2H through deconvolution operation bot ×2W bot ×2D bot , the deconvolution calculation formula is:

[0040] U (1) =TransposeConv(F core ,ker nelsize=2×2×2, stride=2)

[0041] Among them, U (1) represents the number of output channels after the first layer of deconvolution of the decoder; TransposeConv(·) represents the three-dimensional deconvolution operation;

[0042] The deconvolution operation expands the spatial dimension by inserting zero values ​​between pixels of the input feature map and performing convolution, converting the low-resolution abstract features of the bottleneck layer into high-resolution features layer by layer, providing size-matched input for subsequent fusion with the encoder feature map.

[0043] Said S5 further comprises the following contents:

[0044] S5.1: Feature map U generated by the decoder after deconvolution upsampling (m) The high-resolution feature map P at the level corresponding to the encoder n , there are differences in spatial dimensions and number of channels; the decoder aligns the spatial dimensions by padding after the deconvolution operation. Specifically, the deconvolution output feature map size is H'×W'×D', and the corresponding encoder level feature map size is H"×W"×D", then symmetrical padding is performed in the height, width, and depth directions to make H'=H", W'=W", D'=D"; then the number of channels is adjusted by 1×1×1 convolution so that U (m) The number of channels C' and P n The number of channels C” is consistent; the formula is as follows:

[0045] U align =Conv3D(U (m) ,ker nelsize=1×1×1, stride=1, outchannels=C″);

[0046] Among them, U align is the decoder feature map after channel number alignment, Conv3D represents the three-dimensional convolution operation;

[0047] S5.2: After dimension alignment, the decoder feature map U is concatenated through the concatenation layer. alignThe feature map P corresponding to the encoder level n Perform element-by-element addition to fuse the low-level spatial details extracted by the encoder with the high-level semantic features generated by the decoder deconvolution to form a fused feature map F containing multi-scale information. merge ; The formula is:

[0048] F merge =P n +U align .

[0049] Said S6 further comprises the following contents:

[0050] S6.1: In the skip connection path, the attention module receives the output feature map U from the previous layer of the decoder (m-1) and the feature map P from the corresponding level of the encoder n , for the decoder feature map U (m-1) Perform linear transformation to reduce feature dimensions. The formula is:

[0051]

[0052] Among them, Conv3D represents a three-dimensional convolution operation, and the output channels are reduced to C / 2, keeping the spatial dimension unchanged;

[0053] S6.2: Encoder feature map P n Perform a 1×1×1 convolution to make the number of channels consistent with the features after dimensionality reduction of the decoder. The formula is as follows:

[0054]

[0055] Then add the two elements one by one, and then perform nonlinear activation through the ReLU function. The formula is as follows:

[0056] δ=ReLU(φ X (P n )+φ G (U (m-1) ));

[0057] S6.3: Perform a second 1×1×1 convolution on the activated feature δ to restore the number of channels to C, and then normalize it through the Sigmoid function to generate the attention mask S. The calculation formula is as follows:

[0058] S=σ(Conv3D(δ, kernelsize=1×1×1, stride=1, outchannels=C));

[0059] Among them, o' represents the Sigmoid function, the dimension of S is H×W×D×C, and the element value range is [0,1]. Then the attention mask S is combined with the original encoder feature map Pn Multiply element by element and output the weighted feature map Y.

[0060] The S7 further includes the following:

[0061] S7.1: The multi-channel 3D Unet network enhances the U-Net model by introducing an attention module and a hybrid loss function. The network starts with an input CT image of 256×256×10×3 and undergoes continuous convolution and maximum pooling operations to reduce the spatial dimension while extracting progressively more complex features; the hybrid loss function is composed of the mean squared error loss L MSE and structural similarity loss L SSIM The linear combination is as follows:

[0062] L hybrid =αL SSIM +βL MSE ;

[0063] Among them, L MSE Measures the voxel-level pixel difference between the predicted image and the true noise-free image. The calculation formula is as follows:

[0064]

[0065] Among them, Y is the three-channel denoised image output by the network, is the corresponding noise-free reference image, B is the batch size, H, W, D are the image space dimensions, and C = 3 corresponds to the high energy, low energy, and full energy channels;

[0066] L SSIM Calculated based on the structural similarity index SSIM, the formula is as follows:

[0067]

[0068] SSIM compares the brightness, contrast, and structural information of images, with a value range of [-1, 1], where 1 indicates complete similarity. α and β are weighting parameters used to balance the importance of pixel-level accuracy and structure preservation.

[0069] S7.2: The Adam optimizer is used to iteratively optimize the network parameters. The specific process is as follows: first, the preprocessed batch multi-channel CT images and their noise-free reference images are input into the network; then forward propagation is performed, and the data is sequentially passed through the encoder's feature extraction and downsampling, and the decoder's feature fusion and upsampling to generate the denoised image Y, and the hybrid loss L is calculated based on the formula in sub-step S7.1 hybrid; Then backpropagation is performed, and the gradient of the loss function with respect to the parameters of each layer is calculated through automatic differentiation; finally, the Adam optimizer updates the parameters according to the gradient; by repeating this process until the loss converges, after processing by the decoder, the output layer generates high-energy, low-energy, and full-energy denoised images corresponding to the input channel, completing noise suppression and structure preservation, and the output image dimension is consistent with the input.

[0070] Compared with the prior art, the present invention has the following beneficial effects:

[0071] 1. 3D convolution improves spatial consistency: Through the 3D U-net encoder and decoder structure, a 3×3×3 convolution kernel and maximum pooling operations are used to extract multi-scale 3D spatial features layer by layer. This fully utilizes the inter-layer information of CT images and avoids the poor inter-layer consistency caused by 2D CNN processing only a single slice. This effectively reduces the blurring effect of the denoised image and improves the integrity of the 3D spatial structure.

[0072] 2. Multi-channel input enhances noise pattern recognition: Photon counting CT high-energy, low-energy, and full-energy images (IH, IL, and IT) are stacked channel-wise into an H×W×D×3 input matrix, enabling the network to simultaneously learn the noise distribution characteristics of different energy spectrum channels. Compared to single-channel input, multi-channel data provides richer energy information interaction, helping the network more accurately identify noise patterns in each channel, specifically suppressing noise while retaining the effective signal of the corresponding energy spectrum, thereby improving noise reduction effectiveness.

[0073] 3. Skip connections fuse high- and low-level features: The encoder's multi-scale feature maps (P1, P2, P3, Pbot) generated through max pooling are fused element-by-element with the decoder's upsampled feature maps after deconvolution. Padding is used to align spatial dimensions before concatenation, and 1×1×1 convolutions are used to adjust the number of channels. This effectively combines low-level spatial details with high-level semantic features (tissue and organ structure), avoiding the dimensionality mismatch that can occur with direct concatenation in traditional U-nets and improving feature fusion efficiency.

[0074] 4. Hybrid Loss Function Co-Optimization: A linear combination of LMSE and LSSIM (Lhybrid = αLSSIM + βLMSE) is used to balance the importance of pixel-level accuracy and structure preservation. Parameters α and β (e.g., α = 0.006 and β = 0.008) are optimized through trial and error to ensure that the network effectively removes noise while preserving diagnostic details, thereby enhancing the clinical value of CT images. BRIEF DESCRIPTION OF THE DRAWINGS

[0075] Figure 1It is a multi-channel 3D U-net network structure of the photon counting CT image denoising method based on multi-channel 3D U-net of the present invention;

[0076] Figure 2 This is a comparison chart of different methods of the photon counting CT image denoising method based on multi-channel 3D U-net in the mouse denoising task;

[0077] Figure 3 This is an enlarged comparison of different methods of the photon counting CT image denoising method based on multi-channel 3D U-net in the mouse denoising task;

[0078] Figure 4 This is a comparison chart of different methods of the photon counting CT image denoising method based on multi-channel 3D U-net in the walnut denoising task;

[0079] Figure 5 This is an enlarged comparison of different methods of the photon counting CT image denoising method based on multi-channel 3D U-net in the walnut denoising task. DETAILED DESCRIPTION

[0080] The following will clearly and completely describe the technical solutions in the embodiments of the present invention in conjunction with the accompanying drawings. 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 making creative efforts are within the scope of protection of the present invention.

[0081] Example: Figure 1-Figure 5 As shown, the present invention provides a technical solution.

[0082] The photon counting CT image denoising method based on multi-channel 3D U-net includes the following steps:

[0083] S1. Use high-energy, low-energy, and full-energy images acquired by photon counting CT as three-dimensional multi-channel input data to form an image matrix containing spatial dimensions and energy spectrum channels, providing multi-spectral raw data for subsequent three-dimensional noise reduction processing;

[0084] S2, the input data is processed layer by layer through the 3D U-net encoder, and multi-scale three-dimensional spatial features are extracted using 3D convolutional layers, batch normalization layers and ReLU activation functions, and the spatial dimension is gradually reduced in combination with the maximum pooling operation;

[0085] S3, the bottleneck layer at the end of the encoder aggregates deep features to extract abstract features representing the core structure of the CT image, providing a high-level semantic information foundation for the decoder to restore the image;

[0086] S4, upsampling the bottleneck layer abstract features through the decoder's deconvolution operation, gradually restoring the image's spatial dimensions, and completing the initial conversion from high-level features to image dimensions;

[0087] S5. Concatenate the high-resolution feature maps of the corresponding layers in the encoder with the upsampled feature maps of the decoder through skip connections. Before concatenation, align the spatial dimensions through padding operations to fuse high- and low-level features.

[0088] S6. Embed an attention module in the skip connection path, perform linear transformation, nonlinear activation and normalization on the input feature map, generate an attention mask and multiply it element-wise with the original feature map;

[0089] S7. A hybrid loss function consisting of structural similarity loss and mean square error loss is used to train the network. The parameters are adjusted through back propagation, and the high-energy, low-energy, and full-energy denoised images corresponding to the input channels are finally output.

[0090] 10. The photon counting CT image denoising method based on multi-channel 3D U-net according to claim 1, wherein said S1 further comprises the following contents:

[0091] S1.1: Use photon counting CT equipment to obtain high-energy H, low-energy L, and full-energy T images of the same scan object, denoted as I H , I L , I T , all three are three-dimensional volume data, with a unified spatial size of H×W×D; H is height, W is width, and D is depth; the original image is then preprocessed, including removing detector background noise, correcting energy response differences, and normalizing pixel values ​​to the interval [0,1];

[0092] S1.2: The pre-processed high-energy, low-energy, and full-energy images are superimposed according to the channel dimension to form an input matrix X with dimensions H×W×D×3. This part corresponds to the Input layer in the network structure diagram, where the first three dimensions correspond to the spatial dimensions and the fourth dimension is the energy spectrum channel, in the order of high energy, low energy, and full energy, i.e., X=[I H , I L , I T ], as the input data of the 3D U-net network, enables the network to simultaneously learn the noise distribution characteristics of different energy spectrum channels, and provide a multi-dimensional data foundation for subsequent three-dimensional feature learning.

[0093] Said S2 further comprises the following contents:

[0094] S2.1: The input matrix X is subjected to layer-by-layer feature extraction through a 3D U-net encoder. Each layer first performs a 3×3×3 convolution operation (Convolution 3×3) with a stride of 1 and padding of 1 to maintain the spatial dimension. Batch normalization (BN) is then used to reduce internal covariate shift. The ReLU activation function is then used to introduce a nonlinear transformation (Convolution+ReLU+BN) module, forming the basic processing unit of the encoder. The details are as follows:

[0095] For 3D convolution operations, a convolution kernel of size 3×3×3 is used, with a stride of 1 and a padding of 1 to ensure that the spatial dimension remains unchanged. The calculation formula is as follows:

[0096]

[0097] Where W (l) is the convolution kernel of the lth layer, b (l) is the bias term, l represents the level of the network, F (l) Represents the output feature map after the l-th layer convolution operation in the 3D U-net encoder. i, j, and k correspond to the position index of the convolution kernel in the three spatial dimensions of depth D, height H, and width W, respectively, and the value range is 0, 1, and 2; X (l-1) Represents the output of the (l-1)th layer; (x, y, z) is the feature map coordinate, and three-dimensional local features are extracted through convolution operations across space and channels;

[0098] For batch normalization, the convolution output is normalized to reduce the internal covariate offset. The formula is as follows:

[0099]

[0100] in, μ is the feature map after batch normalization of convolution output; β is the batch mean, is the batch variance, γ (l) and β (l) is a learnable parameter, ε=10 -5 A small constant to prevent division by zero;

[0101] The ReLU activation function introduces a nonlinear function through a linear rectification function, and performs a nonlinear transformation on the feature map after batch normalization. The calculation formula is as follows:

[0102]

[0103] Among them, A (l) is the feature map after activation;

[0104] S2.2: After each convolution activation, the feature map is downsampled by a 2×2×2 max pooling layer with a stride of 2. The calculation formula is:

[0105]

[0106] P (l) Represents the feature map after pooling of the first layer; 2x, 2y, 2z represent the starting coordinates of the pooling window, the stride is 2, the windows do not overlap, and the output size is halved;

[0107] Then the feature map size is halved to extract multi-scale features: shallow features retain edge and texture details, and deep features capture the global structure of tissues and organs to form encoder features {P1, P2, P3, P bot}, where P1, P2, and P3 are feature maps generated by the encoder through multi-layer maximum pooling operations, corresponding to different downsampling stages. The subscripts 1, 2, and 3 in the numbers represent the pooling levels. bot is the deepest feature map before the bottleneck layer.

[0108] The S3 further includes the following:

[0109] S3.1: After the encoder undergoes multi-layer 3D convolution and maximum pooling processing, the deepest feature map is formed at the end, denoted as P pre-bot , dimension is H pre-bot ×W pre-bot ×D pre-bot ×C pre-bot ;in, k is the total number of encoder layers, C pre-bot is the number of channel layers; the bottleneck layer performs continuous 3D convolution processing, including 2-3 convolution layers, and the configuration of each layer is consistent with the encoder to ensure the spatial dimension H pre-bot ×W pre-bot ×D pre-bot The number of channels will gradually increase to aggregate deep semantic information.

[0110] S3.2: The bottleneck layer performs convolution on the encoder’s last pooled output and outputs an abstract feature map F. core , whose dimension is H bot ×W bot ×D bot ×C core ; Among them, C core >C pre-bot , C core Indicates the number of channels in the bottleneck layer; the feature map contains the core semantic information of the CT image, which is different from the local details such as edges and textures in the shallow layer of the encoder.

[0111] Said S4 further comprises the following contents:

[0112] S4.1: The decoder receives the abstract feature map F output by the bottleneck layer core , through the deconvolution operation Up Sampling2x2, its spatial dimension is restored to 2H bot ×2W bot ×2D bot , the deconvolution calculation formula for restoring the image size layer by layer is:

[0113] U (1) =TransposeConv(F core ,ker nelsize=2×2×2, stride=2);

[0114] Among them, U (1) represents the number of output channels after the first layer of deconvolution of the decoder; TransposeConv(·) represents the three-dimensional deconvolution operation;

[0115] The deconvolution operation expands the spatial dimension by inserting zero values ​​between pixels of the input feature map and performing convolution, converting the low-resolution abstract features of the bottleneck layer into high-resolution features layer by layer, providing size-matched input for subsequent fusion with the encoder feature map.

[0116] Said S5 further comprises the following contents:

[0117] S5.1: Feature map U generated by the decoder after deconvolution upsampling (m) The high-resolution feature map P at the level corresponding to the encoder n , fused through Skip connection; first align the spatial dimensions through symmetrical padding operation, and then use 1×1×1 convolution Convolution 1×1 to adjust the number of channels, corresponding to the Convolution 1×1 module in the figure, so that the two dimensions match. Specifically, the deconvolution output feature map size is H'×W'×D', and the corresponding encoder level feature map size is H"×W"×D", then symmetrical padding is performed in the height, width, and depth directions to make H'=H", W'=W", D'=D"; then adjust the number of channels through 1×1×1 convolution so that U (m) The number of channels C' and P n The number of channels C” is consistent; the formula is as follows:

[0118] U align =Conv3D(U (m) ,kernelsize=1×1×1,stride=1,outchannels=C");

[0119] Among them, Ualign is the decoder feature map after channel number alignment, Conv3D represents the three-dimensional convolution operation;

[0120] S5.2: After dimension alignment, the decoder feature map U is concatenated through the concatenation layer. align The feature map P corresponding to the encoder level n Perform element-by-element addition to fuse the low-level spatial details extracted by the encoder with the high-level semantic features generated by the decoder deconvolution to form a fused feature map F containing multi-scale information. merge ; The formula is:

[0121] F merge =P n +U align .

[0122] Said S6 further comprises the following contents:

[0123] S6.1: Embed the attention module Attention Block in the skip connection path. The attention module receives the output feature map U from the previous layer of the decoder (m-1) and the feature map P from the corresponding level of the encoder n , for the decoder feature map U (m-1) Perform linear transformation to reduce feature dimensions. The formula is:

[0124]

[0125] Among them, Conv3D represents a three-dimensional convolution operation, and the output channels are reduced to C / 2, keeping the spatial dimension unchanged;

[0126] S6.2: Encoder feature map P n Perform a 1×1×1 convolution to make the number of channels consistent with the features after dimensionality reduction of the decoder. The formula is as follows:

[0127]

[0128] Then add the two elements one by one, and then perform nonlinear activation through the ReLU function. The formula is as follows:

[0129] δ=ReLU(φ X (P n )+φ G (U (m-1) ));

[0130] S6.3: Perform a second 1×1×1 convolution on the activated feature δ to restore the number of channels to C, and then normalize it through the Sigmoid function to generate the attention mask S. The calculation formula is as follows:

[0131] S=σ(Conv3D(δ, kernelsize=1×1×1, stride=1, outchannels=C));

[0132] Among them, o' represents the Sigmoid function, the dimension of S is H×W×D×C, and the element value range is [0,1]. Then the attention mask S is combined with the original encoder feature map P n Multiply element by element and output the weighted feature map Y.

[0133] The S7 further includes the following:

[0134] S7.1: The multi-channel 3D Unet network enhances the U-Net model by introducing an attention module and a hybrid loss function. The network starts with an input CT image of 256×256×10×3 and undergoes continuous convolution and maximum pooling operations to reduce the spatial dimension while extracting progressively more complex features; the hybrid loss function is composed of the mean squared error loss L MSE and structural similarity loss L SSIM The linear combination is as follows:

[0135] L hybriad =αL SSM +βL MSE ;

[0136] Among them, L MSE Measures the voxel-level pixel difference between the predicted image and the true noise-free image. The calculation formula is as follows:

[0137]

[0138] Among them, Y is the three-channel denoised image output by the network, is the corresponding noise-free reference image, B is the batch size, H, W, D are the image space dimensions, and C = 3 corresponds to the high energy, low energy, and full energy channels;

[0139] L SSIM Calculated based on the structural similarity index SSIM, the formula is as follows:

[0140]

[0141] SSIM compares the brightness, contrast, and structural information of images, with a value range of [-1, 1], where 1 indicates complete similarity. α and β are weighting parameters used to balance the importance of pixel-level accuracy and structure preservation.

[0142] S7.2: The Adam optimizer is used to iteratively optimize the network parameters. The specific process is as follows: first, the preprocessed batch multi-channel CT images and their noise-free reference images are input into the network; then forward propagation is performed, and the data is sequentially passed through the encoder's feature extraction and downsampling, and the decoder's feature fusion and upsampling to generate the denoised image Y, and the hybrid loss L is calculated based on the formula in sub-step S7.1 hybrid ; Then backpropagation is performed, and the gradient of the loss function with respect to the parameters of each layer is calculated through automatic differentiation; finally, the Adam optimizer updates the parameters according to the gradient; by repeating this process until the loss converges, after processing by the decoder, the output layer generates high-energy, low-energy, and full-energy denoised images corresponding to the input channel, completing noise suppression and structure preservation, and the output image dimension is consistent with the input.

[0143] The following is a detailed example of a multi-channel 3D U-net-based photon counting CT image denoising method, using CT scan data of a mouse brain (size 256×256×32) and a walnut (size 256×256×20), including specific data and results:

[0144] First, high-energy, low-energy, and full-energy three-channel images of both images were acquired using a photon-counting CT device. The initial pixel range for the mouse data was 0-4096, and for the walnut data was 0-4096. In the preprocessing stage, detector background noise was removed (e.g., the mouse data was subtracted from the background mean of 50), energy response differences were corrected (a gain correction factor of 1.2 for the low-energy channel), and pixel values ​​were normalized to [0, 1]. Subsequently, the pixels were superimposed channel by channel to form input matrices with dimensions of 256×256×32×3 for the mouse and 256×256×20×3 for the walnut.

[0145] The 3D U-net encoder consists of four downsampling stages. Each layer undergoes 3×3×3 convolutions, batch normalization, and ReLU activations before being halved in size using 2×2×2 max pooling. The first stage increases the number of channels from 3 to 32, resulting in an output of 128×128×16×32. The second stage has 64 channels, resulting in an output of 64×64×8×64. The third stage has 128 channels, resulting in an output of 32×32×4×128. The fourth stage has 256 channels, which, after pooling, produces the deepest feature layer of 16×16×2×256. The bottleneck layer uses two convolutions to increase the number of channels to 512, extracting core abstract features such as the overall structure of the mouse brain and the distribution of walnut cavities.

[0146] The decoder restores the dimensions layer by layer through deconvolution. For example, the first deconvolution layer converts 16×16×2×512 to 32×32×4×256. To address the dimensionality difference between the decoder and encoder features, the spatial dimensions are symmetrically padded, followed by 1×1 convolution to match the number of channels (e.g., 256 to 128). Skip connections are used to fuse high- and low-level features. The attention module reduces the dimensionality of the decoder features (128 channels) and the encoder features (128 channels) to 64 channels. After adding activations, an attention mask is generated to highlight key areas such as mouse brain blood vessels and walnut seed coats while suppressing noise.

[0147] The training adopts the mixed loss function (α=0.006,β=0.008), Adam optimizer (learning rate 1×10 -4 ), with a batch size of 16 and training for 150 epochs. The initial mean squared error loss on the mouse training set was 0.012, which converged to 0.003; the structural similarity loss decreased from 0.35 to 0.08 (SSIM increased from 0.65 to 0.92).

[0148] The final noise reduction results: The peak signal-to-noise ratio (PSNR) of the mouse brain image increased from 28.5dB to 34.48dB, the structural similarity (SSIM) increased from 0.72 to 0.7806, the background noise standard deviation decreased from 200 to 80, the cerebral vascular details were clear, close to the noise-free reference image (NDCT); the PSNR of the walnut image increased from 32.2dB to 37.64dB, the SSIM increased from 0.85 to 0.8885, the blurriness of the seed coat membrane edge was reduced by 40%, and the background noise level was the lowest. Figure 2-5 Comparison shows that this method is significantly better than traditional methods such as BM3D and TV in preserving details in the densely vascular areas of mice and the walnut seed coat area, effectively improving the quality of low-dose CT images.

[0149] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above and that the invention can be embodied in other specific forms without departing from the spirit or essential characteristics of the invention. Therefore, the embodiments should be considered in all respects as illustrative and non-restrictive, and the scope of the invention is defined by the appended claims, not the foregoing description, and all variations within the meaning and range of equivalents of the claims are intended to be included therein. Any reference sign in a claim should not be construed as limiting the claim to which it relates.

Claims

1. A photon counting CT image denoising method based on a multi-channel 3D U-net, characterized by: The following steps are involved: S1. Use high-energy, low-energy, and full-energy images acquired by photon counting CT as three-dimensional multi-channel input data to form an image matrix containing spatial dimensions and energy spectrum channels, providing multi-spectral raw data for subsequent three-dimensional noise reduction processing; S2, the input data is processed layer by layer through the 3D U-net encoder, and multi-scale three-dimensional spatial features are extracted using 3D convolutional layers, batch normalization layers and ReLU activation functions, and the spatial dimension is gradually reduced in combination with the maximum pooling operation; S3, the bottleneck layer at the end of the encoder aggregates deep features to extract abstract features representing the core structure of the CT image, providing a high-level semantic information foundation for the decoder to restore the image; S4, upsampling the bottleneck layer abstract features through the decoder's deconvolution operation, gradually restoring the image's spatial dimensions, and completing the initial conversion from high-level features to image dimensions; S5. Concatenate the high-resolution feature maps of the corresponding layers in the encoder with the upsampled feature maps of the decoder through skip connections. Before concatenation, align the spatial dimensions through padding operations to fuse high- and low-level features. S6. Embed an attention module in the skip connection path, perform linear transformation, nonlinear activation and normalization on the input feature map, generate an attention mask and multiply it element-wise with the original feature map; S7. A hybrid loss function consisting of structural similarity loss and mean square error loss is used to train the network. The parameters are adjusted through back propagation, and the high-energy, low-energy, and full-energy denoised images corresponding to the input channels are finally output.

2. The photon counting CT image denoising method based on multi-channel 3D U-net according to claim 1 is characterized in that: Said S1 further comprises the following contents: S1.1: Use photon counting CT equipment to obtain high-energy H, low-energy L, and full-energy T images of the same scan object, denoted as I H , I L , I T , all three are three-dimensional volume data, with a unified spatial size of H×W×D; H is height, W is width, and D is depth; the original image is then preprocessed, including removing detector background noise, correcting energy response differences, and normalizing pixel values ​​to the interval [0,1]; S1.2: The pre-processed high-energy, low-energy, and full-energy images are superimposed according to the channel dimension to form an input matrix X with dimensions H×W×D×3, where the first three dimensions correspond to the spatial dimensions and the fourth dimension is the energy spectrum channel in the order of high energy, low energy, and full energy, i.e., X=[I H , I L , I T ], as the input data of the 3D U-net network, enables the network to simultaneously learn the noise distribution characteristics of different energy spectrum channels, and provide a multi-dimensional data foundation for subsequent three-dimensional feature learning.

3. The photon counting CT image denoising method based on multi-channel 3D U-net according to claim 1, characterized in that: Said S2 further comprises the following contents: S2.1: The input matrix X is subjected to layer-by-layer feature extraction using a 3D U-net encoder. Each layer first performs a 3×3×3 convolution operation with a stride of 1 and padding of 1 to maintain the spatial dimension. Batch normalization (BN) is then used to reduce internal covariate shifts. Finally, a ReLU activation function is used to introduce a nonlinear transformation module, forming the basic processing unit of the encoder. The details are as follows: For 3D convolution operations, a convolution kernel of size 3×3×3 is used, with a stride of 1 and a padding of 1 to ensure that the spatial dimension remains unchanged. The calculation formula is as follows: Where W (l) is the convolution kernel of the lth layer, b (l) is the bias term, l represents the level of the network, F (l) Represents the output feature map after the l-th layer convolution operation in the 3D U-net encoder. i, j, and k correspond to the position index of the convolution kernel in the three spatial dimensions of depth D, height H, and width W, respectively, and the value range is 0, 1, and 2; X (l-1) Represents the output of the (l-1)th layer; (x, y, z) is the feature map coordinate, and three-dimensional local features are extracted through convolution operations across space and channels; For batch normalization, the convolution output is normalized to reduce the internal covariate offset. The formula is as follows: in, μ is the feature map after batch normalization of convolution output; β is the batch mean, is the batch variance, γ (l) and β (l) is a learnable parameter, ε=10 -5 A small constant to prevent division by zero; The ReLU activation function introduces a nonlinear function through a linear rectification function, and performs a nonlinear transformation on the feature map after batch normalization. The calculation formula is as follows: Among them, A (l) is the feature map after activation; S2.2: After each convolution activation, the feature map is downsampled through a 2×2×2 max pooling layer with a stride of 2. The calculation formula is: P (l) Represents the feature map after pooling of the first layer; 2x, 2y, 2z represent the starting coordinates of the pooling window, the stride is 2, the windows do not overlap, and the output size is halved; Then the feature map size is halved to extract multi-scale features: shallow features retain edge and texture details, and deep features capture the global structure of tissues and organs to form encoder features {P1, P2, P3, P bot }, where P1, P2, and P3 are feature maps generated by the encoder through multi-layer maximum pooling operations, corresponding to different downsampling stages. The subscripts 1, 2, and 3 in the numbers represent the pooling levels. bot is the deepest feature map before the bottleneck layer.

4. The photon counting CT image denoising method based on multi-channel 3D U-net according to claim 1, characterized in that: The S3 further includes the following: S3.1: After the encoder undergoes multi-layer 3D convolution and maximum pooling processing, the deepest feature map is formed at the end, denoted as P pre-bot , dimension is H pre-bot ×W pre-bot ×D pre-bot ×C pre-bot ;in, k is the total number of encoder layers, C pre-bot is the number of channel layers; the bottleneck layer performs continuous 3D convolution processing, including 2-3 convolution layers, and the configuration of each layer is consistent with the encoder to ensure the spatial dimension H pre-bot ×W pre-bot ×D pre-bot The number of channels will gradually increase to aggregate deep semantic information. S3.2: The bottleneck layer performs convolution on the encoder’s last pooled output and outputs an abstract feature map F. core , whose dimension is H bot ×W bot ×D bot ×C core ; Among them, C core >C pre-bot , C core Indicates the number of channels in the bottleneck layer; the feature map contains the core semantic information of the CT image, which is different from the local details such as edges and textures in the shallow layer of the encoder.

5. The photon counting CT image denoising method based on multi-channel 3D U-net according to claim 1, characterized in that: Said S4 further comprises the following contents: S4.1: The decoder receives the abstract feature map F output by the bottleneck layer core , and restore its spatial dimension to 2H through deconvolution operation bot ×2W bot ×2D bot , the deconvolution calculation formula is: U (1) =TransposeConv(F core ,ker nelsize=2×2×2,stride=2); Among them, U (1) represents the number of output channels after the first layer of deconvolution of the decoder; TransposeConv(·) represents the three-dimensional deconvolution operation; The deconvolution operation expands the spatial dimension by inserting zero values ​​between pixels of the input feature map and performing convolution, converting the low-resolution abstract features of the bottleneck layer into high-resolution features layer by layer, providing size-matched input for subsequent fusion with the encoder feature map.

6. The photon counting CT image denoising method based on multi-channel 3D U-net according to claim 1, characterized in that: Said S5 further comprises the following contents: S5.1: Feature map U generated by the decoder after deconvolution upsampling (m) The high-resolution feature map P at the level corresponding to the encoder n , there are differences in spatial dimensions and number of channels; the decoder aligns the spatial dimensions by padding after the deconvolution operation. Specifically, the deconvolution output feature map size is H'×W'×D', and the corresponding encoder level feature map size is H"×W"×D", then symmetrical padding is performed in the height, width, and depth directions to make H'=H", W'=W", D'=D"; then the number of channels is adjusted by 1×1×1 convolution so that U (m) The number of channels C' and P n The number of channels C” is consistent; the formula is as follows: U align =Conv3D(U (m) ,ker nelsize=1×1×1,stride=1,outchannels=C″); Among them, U align is the decoder feature map after channel number alignment, Conv3D represents the three-dimensional convolution operation; S5.2: After dimension alignment, the decoder feature map U is concatenated through the concatenation layer. align The feature map P corresponding to the encoder level n Perform element-by-element addition to fuse the low-level spatial details extracted by the encoder with the high-level semantic features generated by the decoder deconvolution to form a fused feature map F containing multi-scale information. merge ; The formula is: F merge =P n +U align 。 7. The photon counting CT image denoising method based on multi-channel 3D U-net according to claim 1, characterized in that: Said S6 further comprises the following contents: S6.1: In the skip connection path, the attention module receives the output feature map U from the previous layer of the decoder (m-1) and the feature map P from the corresponding level of the encoder n , for the decoder feature map U (m-1) Perform linear transformation to reduce feature dimensions. The formula is: Among them, Conv3D represents a three-dimensional convolution operation, and the output channels are reduced to C / 2, keeping the spatial dimension unchanged; S6.2: Encoder feature map P n Perform a 1×1×1 convolution to make the number of channels consistent with the features after dimensionality reduction of the decoder. The formula is as follows: Then add the two elements one by one, and then perform nonlinear activation through the ReLU function. The formula is as follows: δ=ReLU(φX(P n )+φ G (The (m-1) ); S6.3: Perform a second 1×1×1 convolution on the activated feature δ to restore the number of channels to C, and then normalize it through the Sigmoid function to generate the attention mask S. The calculation formula is as follows: S=σ(Conv3D(δ,ker nelsize=1×1×1, strid=1, outchannels=C)); Among them, o' represents the Sigmoid function, the dimension of S is H×W×D×C, and the element value range is [0,1]. Then the attention mask S is combined with the original encoder feature map P n Multiply element by element and output the weighted feature map Y.

8. The photon counting CT image denoising method based on multi-channel 3D U-net according to claim 1, characterized in that: The S7 further includes the following: S7.1: The multi-channel 3D Unet network enhances the U-Net model by introducing an attention module and a hybrid loss function. The network starts with an input CT image of 256×256×10×3 and undergoes continuous convolution and maximum pooling operations to reduce the spatial dimension while extracting progressively more complex features; the hybrid loss function is composed of the mean squared error loss L MSE and structural similarity loss L SSIM The linear combination is as follows: L hybrid =αL SSIM +βL MSE ; Among them, L MSE Measures the voxel-level pixel difference between the predicted image and the true noise-free image. The calculation formula is as follows: Among them, Y is the three-channel denoised image output by the network, is the corresponding noise-free reference image, B is the batch size, H, W, D are the image space dimensions, and C = 3 corresponds to the high energy, low energy, and full energy channels; L SSIM Calculated based on the structural similarity index SSIM, the formula is as follows: SSIM compares the brightness, contrast, and structural information of images, with a value range of [-1, 1], where 1 indicates complete similarity. α and β are weighting parameters used to balance the importance of pixel-level accuracy and structure preservation. S7.2: The Adam optimizer is used to iteratively optimize the network parameters. The specific process is as follows: first, the preprocessed batch multi-channel CT images and their noise-free reference images are input into the network; then forward propagation is performed, and the data is sequentially passed through the encoder's feature extraction and downsampling, and the decoder's feature fusion and upsampling to generate the denoised image Y, and the hybrid loss L is calculated based on the formula in sub-step S7.1 hybrid ; Then backpropagation is performed, and the gradient of the loss function with respect to the parameters of each layer is calculated through automatic differentiation; finally, the Adam optimizer updates the parameters according to the gradient; by repeating this process until the loss converges, after processing by the decoder, the high-energy, low-energy, and full-energy denoised images corresponding to the input channel are generated through the output layer, completing noise suppression and structure preservation, and the output image dimension is consistent with the input.

Images