Remember me
ORIGINAL ARTICLE
Year : 2023 | Volume
: 38
| Issue : 1 | Page : 8-15
99m-Tc MDP bone scan image enhancement using pipeline application of dynamic stochastic resonance algorithm and block-matching 3D filter
Anil Kumar Pandey1, Gagandeep Kaur1, Jagrati Chaudhary1, Angel Hemrom1, Jasim Jaleel1, Param Dev Sharma2, Chetan Patel1, Rakesh Kumar1
1 Department of Nuclear Medicine, All India Institute of Medical Sciences, New Delhi, India
2 Department of Computer Science, SGTB Khalsa College, University of Delhi, New Delhi, India
Correspondence Address:
Dr. Anil Kumar Pandey
Department of Nuclear Medicine, All India Institute of Medical Sciences, New Delhi - 110 029
India
Source of Support: None, Conflict of Interest: None
CheckDOI: 10.4103/ijnm.ijnm_78_22
Abstract
Introduction: In this pilot study, we have proposed and evaluated pipelined application of the dynamic stochastic resonance (DSR) algorithm and block-matching 3D (BM3D) filter for the enhancement of nuclear medicine images. The enhanced images out of the pipeline were compared with the corresponding enhanced images obtained using individual applications of DSR and BM3D algorithm. Materials and Methods: Twenty 99m-Tc MDP bone scan images acquired on SymbiaT6 SPECT/CT gamma camera system fitted with low-energy high-resolution collimators were exported in DICOM format to a personal computer and converted into PNG format. These PNG images were processed using the proposed algorithm in MATLAB. Two nuclear medicine physicians visually compared each input and its corresponding three enhanced images to select the best-enhanced image. The image quality metrics (Brightness, Global Contrast Factor (GCF), Contrast per pixel (CPP), and Blur) were used to assess the image quality objectively. The Wilcoxon signed test was applied to find a statistically significant difference in Brightness, GCF, CPP, and Blur of enhanced and its input images at a level of significance. Results: Images enhanced using the pipelined application of SR and BM3D were selected as the best images by both nuclear medicine physicians. Based on Brightness, Global Contrast Factor (GCF), CPP, and Blur, the image quality of our proposed pipeline was significantly better than enhanced images obtained using individual applications of DSR and BM3D algorithm. The proposed method was found to be very successful in enhancing details in the low count region of input images. The enhanced images were bright, smooth, and had better target-to-background ratio compared to input images. Conclusion: The pipelined application of DSR and BM3D algorithm produced enhancement in nuclear medicine images having following characteristics: bright, smooth, better target-to-background ratio, and improved visibility of details in the low count regions of the input image, as compared to individual enhancements by application of DSR or BM3D algorithm.
Keywords: 99m-Tc MDP bone scan, block-matching 3D algorithm, dynamic stochastic resonance, nuclear medicine image
Introduction 
Nuclear medicine is the gold standard for the functional assessment of the human body.[1] However, the nuclear medicine images are noisy and have low contrast, which leads to difficulty in unequivocal differentiation and delineation surrounding structures. Despite technological advancements in both hardware and software, the observed images in routine NM practice have remained noisy and low contrast [Figure 1]. Specifically, confident interpretation of planar nuclear medicine images is still a challenging task.[2]
The collection of the limited number of gamma photons per unit area from the field of view of the imaging system in the reasonable imaging period is the major reason for the poor image quality. Health risk associated with the radiation exposure limits the amount of radioisotope that can be administered, and so we have to deal with the images obtained from acceptable dose limits only.[3]
There are reports available in the literature, in which dark images have been enhanced using the principle of dynamic stochastic resonance (DSR). In this article,[4] the inherent randomness associated with pixel counts was considered as noise to be added to the corresponding pixel using the mathematical formulation described in Equation 1, and this addition improved the dark images significantly. The same algorithm can be applied on the nuclear medicine images to improve the counts per pixel. In the spatial domain, the image itself has inherent noise and which can be added to the image to improve dark areas. The random noise associated with a pixel count is equal to the square root of the pixel counts, and hence, the counts added to a pixel are proportional to the pixel counts. Visually, the output images processed using DSR resemble the image acquired with an extended period, although it might not be exactly similar based on quantification. However, the increase in pixel counts is also associated with the increase in random noise within the pixel which necessitates the subsequent denoising.[5]
Block-matching with 3D filtering (BM3D) is one of the most efficient images denoizing algorithms available currently. In recent years, it has gained much interest within the image processing community due to its performance advantage, both in terms of reducing noise from images and maintaining visual details. BM3D algorithm cleans the image without disturbing the edges or sharp boundaries.[6] In this study, we propose and evaluate a pipelined application of DSR and BM3D denoizing algorithm on a 99m-Tc MDP bone scan for image enhancement.
Materials and Methods Dynamic stochastic resonance algorithm
The idea “how addition of noise can make signal detectable” can be illustrated with the help of schematic diagram. In [Figure 2], a one-dimensional subthreshold signal (A) after the addition of noise crosses the threshold at different points of time (B) and makes the signal detected (C). There are several ways, in which the right amount of noise can be added to the image to make the details visible both in the dark and bright regions of the image.
Figure 2: Schematic diagram illustrating how stochastic resonance can help in detection of sub threshold signalDSR adds the optimum amount of noise iteratively, to maximize signal-to-noise ratio. The process increases the mean and variance of the image by increasing the brightness and contrast, respectively.
Based on analogy to Benzi et al.'s[7] double well model for the global climate, mathematical formulation of the theory of DSR for contrast enhancement of dark image has been proposed.[4] In the context of image processing, the image pixel and its intensity value are treated as a discrete particles in Brownian motion More Details and its energy, respectively. Benzi's et al.[7] double-well theory suggests two states of image contrast (low and high) and at optimum intrinsic noise density (or an optimum number of oscillations), the image pixel (particle) makes a transition into the other well (from low-to-high contrast, and vice versa). DSR refers to the transition of low-intensity pixel value into high intensity, and vice versa, after the addition of the optimum amount of noise to the image. The optimum amount of noise is reached by stochastic approximation using a discrete iterative equation (Eq. 1).
Where 
denotes the sequence of input signal and noise, with the initial condition being. x(0) = 0 Here, a = 0.5, b = 0.185, △t is the sampling time, taken as 0.15. The Input is the nuclear medicine image. The nuclear medicine image can be viewed as containing signal (image formation) as well as noise (due to the limited number of counts per pixel in the image).
The following simulated example illustrates how DSR improves the low count image (dark image):
Consider a one-dimensional low count image of size 1×11, denoted by X [Table 1] where Xi represents the value of pixel or pixel counts at index i; for example in [Table 1], X1 = 2, and X5 = 11. The pixel value of X lies in the range, (2–31) a grayscale image has a minimum intensity of 0 and maximum intensity 255. The dark image X was processed using the equation (1), and the output after eight iterations is denoted by XSR, the value of XSR was scaled between 0 and 255 and the resultant image have been denoted by XSR_Normalized.
Table 1: Simulated example showing transformation of pixel intensities of the original (input) image by dynamic stochastic resonance, and local contrastWhen the discrete iterative equation Eq.(1) is visualized as a discretized initial value problem (IVP) with original pixel intensity X at a point as the initial value and rate of its change during iterations of DSR the corresponding temporal differential equation, application of classic Runge–Kutta method estimates the solution of the IVP. Vectorized solution for all points in X was used to compute local contrast in the enhanced image. [Table 1] shows the computed values of the local contrast of the enhanced image from the local contrast of the input (original) image. The table also presents the normalized computed values (scaled between 0 and 255).
An inspection of the data provided in [Table 1] leads to the observations:
Dark image [pixel intensity range (2–31) has been transformed into a bright image (145–255) (that is increased counts per pixel) and also there is increase in global contrast of the image (from 29 to 110; contrast = max– min) as result of the application of DSR algorithmThere is also increase in local contrast of the image as a result of the application of DSR algorithm; however, there exists no relationship between the local contrast in the enhanced image and local contrast in the input image. In other words, the increase in contrast at each location in the images is random. The local contrast of the input image [denoted by lc_X, [Table 1]] and local contrast of the enhanced image [denoted by lc_XSR_Normalized, [Table 1]] was calculated as the difference between the pixel value of the current location minus the pixel value at the location immediately right to itThe intensity value of the output image is concentrated toward the higher side of gray scale that is in the range (145–255).However, in order to perceive all the details in the image, the intensity values should be distributed in the entire range of the gray scale that is from dark to bright, not concentrated in either of the dark or gray or bright intensity range. Thus, it is important to optimize the number of iterations of the noise to be added.
99m-Tc MDP bone scan images were processing using the DSR algorithm, we kept the number of iterations equal to 25 for processing all 20 images. After processing the image with DSR algorithm, the image can be assumed to be contaminated with Additive White Gaussian Noise (AWGN); the pixel counts might have increased to be above 30 and in this case the pixel counts can be assumed to follow the Gaussian distribution. As AWGN is random in nature and it corrupts almost all areas of images, it is challenging to remove AWGN from images. It becomes increasingly difficult to preserve the small details of an image as the error level increases.
Block-matching 3D filtering
Readers interested in mathematical details may refer to the article written by the author of the BM3D algorithm.[6] Following the description of BM3D algorithm for image denoizing is in brief without mathematical details:
Block-matching and 3D filtering,[6] exploits nonlocal image modeling[8] with a method termed grouping and collaborative filtering. Proposed in 2007 by Dabov et al. BM3D is currently the state-of-the-art method for image denoizing and outperforms all other algorithms when it comes to removing AWGN at a reasonable computational cost.
There are two principal properties of natural scene images: the existence of mutually similar patches within the close neighborhood and the local correlation of pixel values within a single patch; the BM3D filter utilizes these properties for image denoizing. Similar 2D patches are grouped into 3D data arrays that allow to exploit of both intra-patch correlations and inter-patch correlations. As a result, the group enjoys correlation in all three dimensions, allowing a highly sparse representation of the image in the transform domain. The term sparse refers to linear transformation domain representations which have few high magnitude coefficients and many low magnitude coefficients. Thus, sparseness refers to the energy compaction property (i.e., most of the image details are representation by few large coefficients while noise is spread across a range of small coefficients). This sparse representation allows us to separate the noise from the true signal by applying shrinkage on the coefficients in the transform domain. This approach of exploiting similarity and estimating the original signal is collaborative filtering which has three steps:
For each reference patch, find similar patches from the input image by classifying them according to some similarity criteria and transform them into a 3D data array by grouping the matched 2D blocksApply shrinkage of the coefficients in the transformed 3D spectrum to attenuate the noiseApply inverse 3D transform to the shrunken coefficients and return the obtained 2D estimates of the grouped blocks to their original positions.The collaborative filtering process gives a 3D estimation of the jointly filtered 2D blocks. As the grouped 2D blocks are similar, the transformation can achieve a very high sparse representation of the original signal. This process reveals the finest details shared by grouped blocks by preserving the unique features of each individual block.
The BM3D algorithm can be distinctively decomposed into two almost identical steps. The two steps in BM3D work as follows:
Step 1: In this step, an intermediate (i.e., basic estimate) denoized image is estimated using hard thresholding during the collaborative filtering process by taking the original image (noisy image) as the input:
Grouping and collaborative filtering: the input image is processed block by block. These blocks are called reference blocks. For each reference block, similar blocks to the currently processed one are found using a similarity measure (block matching). A 3D group (array) is then built by stacking the matched blocks and a collaborative filter is applied to the grouped blocks. In this step, hard thresholding is applied during shrinkage of the coefficients in the transform domainAggregation: After collaborative filtering, we get an estimate of each block and a variable number of estimates for each pixel due to the overlapping of blocks. The output of the first step is obtained by weighted averaging of all the achieved block-wise estimates that have overlapped.Step 2: The step produces the final denoized estimate based on both the input noisy image and the basic estimate obtained from step 1. Here, instead of hard thresholding, Wiener filtering is used as the shrinkage method. It applies Wiener filtering to the original noisy input image by using the basic estimate obtained from step 1 as an oracle.
The authors of this article[6] have already provided the data and MATLAB codes for validation of their algorithm on their official website. The MATLAB program provided by them has been used in this study. The important input parameter that which user needs to supply is the value of estimated noise in the image (sigma). The value of sigma controls the amount of smoothness which was kept equal to 25 irrespective of the characteristics of the input image.
Image acquisition and processing
Twenty 99m-Tc MDP bone scan images were acquired using the following protocol: 555–740 MBq of 99m-Tc MDP was administered three to 4 h before acquisition. Anterior and posterior images of the whole body were acquired at table scan speed of 20 cm/min on Siemens SymbiaT6 SPECT/CT scanner. The scanner was equipped with low-energy high-resolution collimator and the image acquisition matrix size was . These images were exported in DICOM format and transferred to a personal computer. An in-house MATLAB program read these DICOM images and converted them into PNG format. The pixel depth of the image was 8-bit in PNG format. The conversion was done so that after processing they can be faithfully displayed on the personal computer display monitor.
These images were processed individually using our proposed method, that is, applying the DSR algorithm and BM3D algorithm in pipelined fashion.
Qualitative assessment of image quality
Two nuclear medicine physicians visually compared each output image (obtained after the application of the DSR algorithm, application of BM3D denoizing algorithm, and pipelined combination of DSR and BM3D algorithm) with the corresponding input image and selected one image from among the four images. Their selection was biased with the following characteristics: the image is bright enough so that even small structures were visible unequivocally; had good contrast; had sufficient anatomical details visible and; had good spatial resolution.
Quantitative assessment of image quality
The quality of images was assessed using Brightness, Global Contrast Factor (GCF),[9]Contrast per pixel (CPP),[10] and Blur.[11]Brightness was estimated as mean intensity of the image which measures the perceived brightness of the image. GCF is the average local contrast of smaller image fractions. Higher the GCF, higher is the detail in the image. CPP is defined as the average absolute difference of luminance value with the adjacent pixels. Higher the CPP value, more the CPP. Generally, if the difference with the neighboring pixels is high, distinguishing object details is easy. Blur is the nonreference perceptual blur metric with values ranging from 0 to 1. Zero value of blur indicates no smoothing while 1 indicates heavy smoothing. Higher the value of all these factors (Brightness, GCF, CPP, and Blur), the more favorable will be the image quality except in the case of blur. Blur value should be optimally high, and should not be close to 1 as this will result in undesirable over smoothening of the image.
Statistical analysis
Shapiro–Wilk normality test was applied on the brightness, GCF, CPP, and Blur of all the four categories of data. We found P value less than the significance level of 0.05 which indicates that normality is violated. Thus, we applied Kruskal–Wallis rank sum test and pairwise comparison using the Wilcoxon rank-sum test with continuity correction.
R open-source statistical software[12] was used for the statistical analysis (for Wilcoxon signed–rank test) and boxplots. The EBImage package[13] developed to be used with R software was used for arranging the images. All experiments were performed on a personal computer having Windows 7 Home Basic (copyright© 2023 Microsoft Corporation) 64-bit operating system, 2 GB RAM, and Intel (R) core (TM) i3-2120 CPU @ 3.30 GHz processor.
Results In all 20 cases, both the nuclear medicine physicians selected the image obtained from the pipelined application of DSR and the BM3D algorithm. These images were smooth, and had sufficient brightness enabling effortless perception of details of low count areas of the input image. The benefit of our proposed method is more clearly observable on low count areas of the image where structures and boundaries of the image became clearly visible [[Figure 3] (bottom right), [Figure 4] (bottom right), [Figure 5] (bottom left)] in the enhanced image compared with its input image.
Figure 3: Pipelined combination of stochastic resonance and BM3D algorithm Top Row: (left) original input image, (right) BM3D denoised input image. Bottom row: (left) SR preprocessed input, (right) BM3D denoised SR preprocessed input. BM3D: Block-matching 3D
Figure 4: Pipelined combination of stochastic resonance and BM3D algorithm Top Row: (left) original input image, (right) BM3D denoised input image. Bottom row: (left) SR preprocessed input, (right) BM3D denoised SR preprocessed input. BM3D: Block-matching 3D
Figure 5: Pipelined combination of stochastic resonance and BM3D algorithm Top Row: (left) original input image, (right) BM3D denoised input image. Bottom row: (left) SR preprocessed input, (right) BM3D denoised SR preprocessed input. BM3D: Block-matching 3DAlthough images processed with DSR were noisy, these images were brighter than other images and details were also easily perceptible [[Figure 3] (bottom left), [Figure 4] (bottom left), [Figure 5] (bottom left)].
When BM3D is directly applied on the input image, the output image is comparatively less bright than the input image, target-to-background ratio is better, background nearly disappears because of excessive smoothing; however, the effect of smoothening on edges/boundaries is minimal. The ribs and vertebrae are much better visualized/clearly delineated in the output image as compared to the input image [[Figure 3] (top right), [Figure 4] (top right), [Figure 5] (top right)]; however, there is no improvement in the dark region of the image, that is, details in the dark region are identical in both input and output images.
In both images, processed using DSR and images processed using the proposed pipeline (Inp-SR-BM3D), increase in the size of bright lesions, and merging of two or more very closely spaced bright lesions were noticed [[Figure 5] (bottom left), [Figure 5] (bottom right)].
Results of objective image quality assessment
Images processed with our proposed pipeline were found brightest among all the four categories [[Figure 6] (bottom right); median value of D is highest [Table 2]].
Figure 6: Box-plot for-Top Left: Blur, Top Right: Global Contrast Factor, Bottom left: Contrast per pixel, and Bottom right: Brightness; A: Input Image, B: Input-SR, C: Input-BM3D, D: Input-SR-BM3D. BM3D: Block-matching 3D
Table 2: Pairwise comparison using Wilcoxon rank sum test with continuity correctionWhile inspecting the boxplot [[Figure 6] Top Right], the reader can see that the GCF of B was the highest, however, the difference between B (highest) and D (second highest) was statistically insignificant [Table 2].
The CPP of the D was found superior as compared to all other four categories [[Figure 6] bottom left]. Although the median value of B and D were quite close, statistically, there was no statistically significant difference between the two [Table 2].
The blur value of the D was found the highest among all four categories [[Figure 6] top left]. This indicates the efficient smoothening or denoizing of the image by our proposed pipeline [Table 2].
From the above analysis, it is evident that the images enhanced by our proposed pipeline were brightest, smooth and having much better contrast.
Discussion In this study, we have proposed and evaluated an image enhancement method that uses the pipelined application of DSR and BM3D algorithms for nuclear medicine images. Twenty 99m-Tc MDP bone scan images were enhanced using DSR, BM3D algorithm, and pipelined application of DSR and BM3D algorithm. Two NM physicians compared enhanced images with the corresponding input image and selected the best images out of four images (one input and three enhanced images). In all 20 cases, the NMPs selected enhanced images obtained from the result of the application of our proposed method, that is, pipelined application of DSR and BM3D algorithm. All the above visual observations are supported by the quantitative analysis using image quality metrics: Brightness, Blur, GCF, and CPP.
Increase in size of the bright lesions was noticed in the enhanced image obtained using the proposed method. This might be due to the addition of more than the sufficient number of counts in the bright area of the image, leading to oversaturation of pixel counts at pixels containing bright lesions and its surroundings.
Because of this effect, some of the closely spaced bright lesions were found to be merged and visualized to be as single bright lesions. It is to be noted that a fixed number of iterations (equal to 25) was used independently of the counting statistics of the input image to be enhanced. In future, we would like to take up the issue of optimization of the number of iterations to avoid the increase in the size of the bright lesion and oversaturation of pixels surrounding it. The input parameters of the BM3D algorithm were also not tailored based on the amount of noise present in the input image (the value of sigma was equal to 25 for all images). Thus, in the future study, we would like to develop a scheme to fine-tune both numbers of iterations of DSR and sigma of BM3D in order to get an enhanced image without oversaturation of pixels.
Several image enhancements have been proposed and evaluated for nuclear medicine images, and some of these have been implemented by vendors and come as a part of image processing software package with SPECT/CT or PET/CT scanner.[14],[15],[16],[17],[18] NM Physicians/NM Technologists use some of these in nuclear medicine practice as per their comfort levels. We have not compared the result of our study with any of such routinely used software tool for image enhancement. In such pipelined application, there are limited data in the literature related to the enhancement of nuclear medicine images.[19]
The highlight of this study is that we have used personal computers for all implementation independent of costly commercial vendor hardware and software available with the SPECT, SPECT/CT, and PET/CT systems.
Conclusion The pipelined application of the DSR and BM3D algorithm to 99m-Tc MDP bone scan images produces a visually brighter, smoother, and better target-to-background ratio image compared to the corresponding input image.
Financial support and sponsorship
Nil.
Conflicts of interest
There are no conflicts of interest.
References 
Comments (0)