Noise Estimation for Single-Slice Sinogram of Low-Dose X-Ray Computed Tomography Using Homogenous Patch

We present a new method to estimate noise for a single-slice sinogram of low-dose CT based on the homogenous patches centered at a special pixel, called center point, which has the smallest variance among all sinogram pixels. The homogenous patch, composed by homogenous points, is formed by the points similar to the center point using similarity sorting, similarity decreasing searching, and variance analysis in a very large neighborhood VLN to avoid manual selection of parameter for similarity measures.Homogenous pixels in the VLN allow us find the largest number of samples, who have the highest similarities to the center point, for noise estimation, and the noise level can be estimated according to unbiased estimation. Experimental results show that for the simulated noisy sinograms, the method proposed in this paper can obtain satisfied noise estimation results, especially for sinograms with relatively serious noises.


Introduction
With continued technology advancement and wider applications, use of computed Tomography CT is increasing.However radiation exposure and associated risk of cancer for patients receiving CT examination have been an increasing concern in recent years.Thus, minimizing X-ray exposure to patients has been one of the major efforts in the CT field 1-3 .
A simple and cost-effective means to achieve low-dose CT applications is to lower X-ray tube current mA as low as achievable.However, dose reduction generally leads to a nonlocal neighborhood with size 21 × 21 and very recently is used in low-dose CT imaging 25-32 .In the very large neighborhood, the homogenous patches are determined by finding similar points to the center point through similarity sorting, decreasing similarity searching and variance analysis.The noise parameters will be estimated on this homogenous patch using unbiased estimate.Since very large neighborhood provides more reliable estimation, proposed method can get satisfied noise estimation results.
The remainder of this paper is arranged as the noise models will be discussed in Section 2; the motivation and method will be given in Sections 3 and 4, respectively; and then the experimental results will be given in Section 5; finally, it is the conclusion, future works, and acknowledgment.

Noise Models
For clinical X-ray system, its detected that X-ray intensity follows a compound Poisson distribution 9-11 .However, for its complex expression, it has no analytic formula for the likelihood function.Various approximations have been proposed for it.Two common models include Poisson distribution and Gaussian distribution.
In this section, we will introduce signal-independent Gaussian noise SIGN , Poisson noise, and signal-dependent Gaussian noise, as well as their relations and the reason for addressing SIGN.

Signal-Independent Gaussian Noise (SIGN)
SIGN is a common noise for imaging system.Poisson noise and signal-dependent Gaussian noise can be converted to SIGN using scale transforms which will be discussed in Section 2. 3.
Let the original projection data be {x i }, i 1, . . ., m, where i is the index of the ith bin.The signal has been corrupted by additive noise {n i }, i 1, . . ., m, and one noisy observation where n i is an observation for the random variable N i as normal N 0, σ 2 N and independent to the Gaussian random variable X i where the uppercase letters denote the random variables and the lower-case letters denote the observations for respective variables.

Poisson Model and Signal-Dependent Gaussian Model
The photon noise is due to the limited number of photons collected by the detector 42 .For a given attenuating path in the imaged subject, N 0 i, α and N i, α denote the incident and the penetrated photon numbers, respectively.Here, i denote the index of detector channel or bin and α is the index of projection angle.In the presence of noises, the sinogram should be considered as a random process and the attenuating path is given by where N 0 i, α is a constant and N i, α is Poisson distribution with mean N.

Mathematical Problems in Engineering
Thus we have Both its mean value and variance are N. Gaussian distributions of ployenergetic systems were assumed based on limited theorem for high-flux levels and followed many repeated experiments in 6 ; we have where μ i is the mean and σ 2 i is the variance of the projection data at detector channel or bin i, γ is a scaling parameter and f i is a parameter adaptive to different detector bins.
The most common conclusion for the relation between Poisson distribution and Gaussian distribution is that the photon count will obey Gaussian distribution for the case with large incident intensity and Poisson distribution with feeble intensity 6 .In addition, in 42 , the authors deduce the equivalency between Poisson model and Gaussian model.Therefore, both theories indicate that these two noises have similar statistical properties and can be unified into a whole framework.

Scale Transformations
As a statistical method for treating nonnormally distributed or signal dependent noise, scale transformations are widely-used to stabilize the variance 6, 31 .In 4 , authors indicate that the variance of their signal-dependent Gaussian model represented in 2.4 can be stabilized by where c is the expected constant standard deviation for transformed data and c 1 is an arbitrary constant or where c 2 also is an arbitrary constant.Since both Poisson noise and signal-dependent Gaussian noise can be converted to SIGN, the noise estimation for low-dose CT can start from SIGN to focus on our new method itself.

Backgrounds and Motivation
Following the SIGN model discussed in Section 2.1, we will introduce backgrounds and motivation for proposed method.
Noise estimation for low-dose CT projection data is an estimation problem.Using the terms in Section 2.1, the goal for SIGN estimation is to get unbiased estimated value σ N of σ N from the noisy observation y i .Here σ N is the parameter for population distribution of the noise.
However, only one noisy observation y i is provided.Thus we have to share statistical information each other.One reasonable way for sharing statistical information is to assume independent identical distributions iid for nearby similar pixels among noisy projection data sinogram .In this paper, the similar point is the point with high block similarity to the center point.This section will discuss the motivation for finding samples of unbiased estimator.

Unbiased Estimation
Suppose we have a statistical model parameterized by θ, and a statistic θ which serves as an estimator of θ based on any samples x.That is, we assume that population data follows Gaussian distribution with a fixed and unknown constant θ, and then we construct some estimator θ that maps samples to values that we hope are close to θ.Then the bias of this estimator is defined to be where E • denotes expected value over the distribution, that is, averaging overall possible samples x.An estimator is said to be unbiased if its bias is equal to zero for parameter θ.Suppose x {x 1 , . . ., x n } are observations for a Gaussian random variable with expectation μ and variance σ, the unbiased estimators for these samples are

3.2
In theory, the sample estimators shown in 3.2 are asymptotically unbiased and efficient for the sample size is moderate or large.Thus unbiased estimation requires providing enough samples to approximate the population distribution.

Measure Similarity
In order to measure similarity between two different points x i,j and x s,t of noisy sinogram where i / s or j / t, we should measure two 7 × 7 blocks centered at x i,j and x s,t , respectively.The similarity is measured as Generally, we can predefine a threshold T to find the similar points of x i,j in a nonlocal neighborhood with 21 × 21 centered at x i,j x s,t is a similar point of x i,j if S i, j , s, t ≤ T, x s,t is not a similar point of x i,j if S i, j , s, t > T,

Motivation
Since only samples with moderate or large size can approximate the population distribution well, we must find enough samples to estimate the noise correctly.Moreover, iid assumption requires to share statistical information among similar points of noisy sinogram.Both requirements coincide with the start point of existing methods.That is, noise estimation should be performed at smoother versions of the images to suppress the influence of the complex structures of images.
However, estimation noise using smoother versions of the images has high computation burden and overestimate noise levels 24 .In this paper, we propose a method to estimate noise in a VLN by similarity sorting and variance analysis.
Just as discussed in this section, the key objective is to find enough similar points for noise estimation.Unlike existing global methods, which estimate noise levels using whole images, the proposed method try to estimate noise levels in a more local way to reduce computation burden and increase flexibility for estimation.
In order to accomplish the above objectives enough samples with a more local structure , VLN is used to find similar points samples .Moreover, the VLN should be put in a suitable position to ensure the enough samples can be found.That is, VLN should be put in a homogenous region.Since the center of VLN can determine the position of VLN, how to put VLN in a suitable position can be converted to how to locate the center point in a homogenous region.
The noise level of low-dose CT is low.In this situation, only variance is enough for describing the local homogeneity roughly.That is, large variance relates to a square near singularities while small variance relates to a homogenous square.Therefore, by comparing variances of squares with fixed size for all pixels in sinogram, the center point is defined as the pixel with the smallest variance among all pixels.
After determining the center point, we propose a new method for noise estimation based on similarity sorting, similarity decreasing searching, and variance analysis.Its main start point is from how to avoid threshold setting in finding similar points to the center point since the threshold selection depends on the noise level which makes it a "chicken and egg" problem.
Similarity sorting provides a similarity decreasing sequence SDS .Thus the samples must be formed by the up at the front points in the SDS since points with smaller similarities to the center point maybe the outliers for the estimation.
In order to find the largest number samples with the highest similarities in the VLN, similarity decreasing searching combining variance analysis is used.That is, the samples are formed by adding 50 points each time according to the order of the SDS.Thus it ensures that each addition adds the points with the highest similarities in residual SDS.
Moreover, in order to find the largest number samples used in estimation, we must find when the outliers are added.It can be achieved by variance analysis for each addition.When the variance for an addition becomes large suddenly, it means some outliers are added.Therefore, the real samples are the last before this addition.The detailed algorithm will be introduced in Section 4 and two-number examples can be found in Section 5.3.

The Method
In this section, we will give the framework for noise estimation.Just as shown in Figure 1, it includes three steps.
1 Find the center point: It is the key step to locate the VLN, which should be put in a homogenous region.Thus the center point is defined as the point with the smallest local variance in all sinogram pixels.
2 Determine the homogenous patch: The homogenous patch is composed by similar points to the center point samples , and these similar points are searched using similarity sorting, similarity decreasing searching and variance analysis.
3 Estimate noise: The parameters are estimated using unbiased estimator by the samples on the homogenous patch.
By these steps, the largest number of samples with the highest similarities can be obtained.These samples form a very large homogenous patch for noise estimation.Thus even using simple unbiased estimate, satisfied results can be obtained.The details for these three steps will be given in the remainder of this section.

Find Center Pixel
The center pixel should locate at the center of a homogenous region.However, in noise, finding a large regular homogenous region correctly and putting a pixel of sinogram in the center of this region is not a trivial task.
Motivated by 27 , irregular pixel patch composed by similar pixels to the center point is a reasonable choice.Thus the choice for the center point becomes an easy task which only needs to ensure it is the center of a small homogenous square.
Since the variance can describe how far the numbers lie from the mean, small variances will relate to homogenous squares while large variances will relate to squares near singularities.Thus the center point is chosen as the pixel with the smallest variance in a 7 × 7 where X is the sinogram and is the mean.By this way, we can find a point centered at a homogenous square and then extend it to an irregular homogenous patch.

Determine Similar Points (Homogenous Patch)
After finding the center point x i,j , the irregular homogenous patch P composed by the similar points to x i,j should be determined.
Each similar point is found by computing the similarity between itself and x i,j using 3.3 .However, unlike existing methods which use a predefined threshold value T shown in Section 4.1, our method proposes a different scheme by similarities sorting and variance analysis in a VLN.
The main advantage for this scheme is that it avoids threshold setting in finding similar points.Since the threshold should be set according to different noise level, noise estimation and threshold setting become a "chickens and eggs" question.
The most straight motivation for the proposed method is that the small variance indicates only the homogenous samples, while the large variance indicates not only the homogenous samples but also outliers.Thus the variances can be considered as a sign for outliers.It reminds us that if we sort the points according to the their similarities to x i,j in a VLN and then compute the variances for increasing adding points until the big variance is reached, the most reliable noise estimation can be obtained by the as large as possible number for the samples.
Thus we can add samples according to the fact that samples with bigger similarities will be added earlier, and compute the variance after each addition.When the added variance become big suddenly, it means that there are some outliers that are mixed in estimation.Thus the real noise level is estimated using the last before this addition.
In summary, the steps for finding similar points are as follows.
1 Compute the similarity between the center point X i,j and each point in the 41 × 41 square using 3.3 .
2 Sort the similarities to form a similarity descending sequence SDS .
3 Add 50 samples each time to form collections of samples and then estimate the standard deviations SDs of the collection of samples to form an SDs sequence.
4 Compute the variance difference using the i 1th and ith, i 1, . . ., 31, elements of the SDs sequence to form a difference variances sequence and from the 21 addition, find the variance becomes big suddenly is bigger than 5 , and denote its index k.
5 The estimated variance of the noise is k − 1th element of the SDs sequence obtained on the step 3.

Experimental Results
In this section, we will compare our method to some well-known noise estimation methods.Just as discussed in Section 1, most of existing methods for noise estimation are not suitable for sinogram for their complex prior, segmentation for images or practical setting parameters.Only methods which estimate noise by filtering images to suppress singularities previously can be used for noise estimation of sinogram 16-19, 24 .We will introduce two of these methods firstly in Section 5.2 and then compare them with the proposed method in Section 5.3.

Data
Five test images are used in this paper: a thorax phantom acquired from a GE HiSpeed multislice CT scanner see Figure 2  The SIG noises are added on the projection data in Figure 3 directly.Since the interesting SIG noises are with small variances SIG noises with large variances are not accepted in clinic situations , the variances are form 25 to 225 the standard deviations are from 5 to 15 .

Comparison Methods
Noise estimation is not a trivial task in image processing because of the complex structures for images.Some researchers suggest that the noise levels should be estimated by filtering images previously to suppress the image structures 17, 18, 24 .The methods compared with the proposed method are as follows.

Wavelets
This method suppresses image structures using wavelets 18 .Since coefficients in the HH subbands of wavelets only preserving high frequency energy for images, the image structures, which mainly reside in low frequency coefficients, can be considered to be suppressed in the HH subbands.Thus based on robust median estimator, the level of the noise can be estimated as where X is the image, W is the wavelet operator, HH indicates the coefficients in the HH subbands, | • | is the absolute value, and median is the median estimator.

Fast Estimation (FE)
This method is proposed in 17 .It only has two steps for noise estimation: the first step is filtering the image using Laplacian operator and the noise level is estimated using where w and h are the width and height of the image, respectively, ⊗ is the convolution operator, L is the Laplacian operator:

Experimental Results
In this section, three methods, proposed method PM , wavelets, and Fast estimation FE , are compared.In order to compare these three methods on a fair stage, the parameters used in these methods are fixed for all noise levels and all test images.According to this standard, the method proposed recently in 24 is not included because one of its parameter must be adjusted and selected with the assumption that the real noise levels are known.The wavelet used for wavelets is Symlets with support 4. In summary, the parameters for PM are the squares used for finding the center point which are 7 × 7 squares see Section 4.1 ; the size of the square computing similarity between the considering pixel and the center point is 7 × 7 see Section 3.2 ; the size for very large neighborhood is 41 × 41 see Section 4.2 step 1 ; adding 50 samples each time to form the SDs sequence see Section 4.2 step 3 ; the threshold of absolute difference between adjacent standard deviations for finding the variances become ssuddenly 5; in order to avoid the influence of the outliers, the comparison for the variance is from the 21th addition for samples see Section 4.2 step 4 .
Firstly, in order to show the role of variance analysis for the SDs sequence, we will use two Figures , Figures 2 c and 2 d , for analyzing the SDs sequence see Table 1 , where "Real SD" is the real standard deviation SD for the added noise, "No" is the times adding samples, "ESD" is the estimated standard deviation for each sample adding, "Diff" is the difference between the ith and i−1th estimated variances, the bold digits indicate the difference between two adjoint variances beyond the 5 bigger diff , and the italic digits are the final estimated SDs.
According to the steps in Subsection 4.2, if the index for bigger diff is i, the final estimated SD is the i − 1th SD.For example, in the second column, which is Figure 2 c added noise with real SD 5, we can find that the bigger diff is 19.68 and its index is 26 see the last no blank row of the first and second columns .Thus the index for final estimated SD is 25 represented by an italic digit, and the SD is 4.90 see the first and second columns with no. 25 .
The final estimated SDs can be get through four steps see Section 4.2 .The following two tables Tables 2 and 3 give the final estimated SDs for sinograms in Figure 3 added with noise whose SDs are from 5 to 15.
Table 2 gives the final estimated SDs for phantoms see Figures 3 a -3 c , "Real SD" is the real standard deviation SD for the added noise, and "PM," "W," and "FE" are the abs. of the proposed method, wavelets 18 , and fast estimation 17 , respectively.The bold digits represent the best estimate using three estimators while the digits with * are the worst estimate.In order to ignore small differences between two estimators, two digits will not be signed if the difference between two biases difference of estimated SDs and the real SD is smaller than 0.05.For example, in b , "PM" and "W" have very similar estimated results, they are not signed.
From Table 2, we can see that "PM" and "W" have very similar estimated results for b , while "FE" has the worst estimated results in all noise levels; estimated results using "PM" are the best in most of noise levels, especially in relatively serious noises of c and the worst estimator is "FE" also, which overestimates noise levels for both subfigures.
Table 3 gives the final estimated SDs for real projection data see Figures 3 d -3 e .The all signs on Table 3 are the same as they are on Table 2.
From Table 3, we also can find that "PM" and "W" have very similar estimated results, while "FE" is the worst estimator in all noise levels for its overestimated SDs.For e , "PM" has better performance in relatively serious noises, while "W" has better performance in relatively low noises.
It should be indicated that the estimated results using "PM" on a are not satisfied especially in low noises.It reminds us, that there are also future works which can be done for improving the estimated results for the proposed method which will be discussed in Section 7.

Conclusion
In this paper, we propose a new method to estimate noise for sinograms of low-dose CT.The proposed method can obtain estimated results both for and real projection data, especially in relatively serious noises, which demonstrate its potential for noise estimation of sinograms of low-dose CT.
Based on the similarity sorting and variance analysis in a very large neighborhood whose scale is 41 × 41, we can find enough similar samples to obtain reliable estimated results which make this proposed local method have very similar estimated results to the best existing global methods.
In addition, avoiding convolution for suppressing image structures and relatively homogenous local structure makes the proposed method also be easily generalized to the more complex noises, such as, Possion noise and Gaussian compound noise.Thus the proposed method is also a promising method for real sinograms of low-dose CT.

Future Works
Although this paper proposes a new powerful method for simulated sinogram noise estimation, it may be improved as follows.
1 How to find a center point in a large homogenous patch to ensure that there are enough points to obtain reliable estimation.We try to use multiresolution method to solve this problem.

Figure 1 :
Figure 1: The flow chart for the proposed method.

Figure 2 :
Figure 2: Test images.a Thorax phantom data acquired at 120 kVp, 150 mA.b A phantom data produced by Matlab.c Modified Shepp-Logan head phantom whose size is 256 × 256.d A CT image of an abdomen of a 58-year-old man with 120 kVp, 150 mAs.e A CT image of a chest of a 62-year-old woman with 120 kVp, 150 mAs.

Figure 3 :
Figure 3: The relative projection data for the test images.

Table 1 :
Variance analysis to estimate the noise level.

Table 2 :
Estimated noise levels for the phantoms.

Table 3 :
Estimated noise levels for the real projection data.
nos. 2010CB732501, 2011CB302801/2011CB302802 , China Postdoctoral Science Foundation of China no.20070410843 , and Open Foundation of Visual Computing and Virtual Reality Key Laboratory Of Sichuan Province no.J2010N03 .