Automatic Noise Estimation from the Multiresolution Support

We will say that a multiresolution support (Starck et al. 1995) of an image describes in a logical or Boolean way if an image I contains information at a given scale j and at a given position (x,y). If M^(I)(j,x,y) = 1 (or = true), then I contains information at scale j and at the position (x,y). M depends on several parameters:

• The input image;
• The algorithm used for the multiresolution decomposition;
• The noise;
• All constraints that we want the support to satisfy in addition.

Such a support results from the data and the treatment (noise estimation, etc.), and from our knowledge of the objects contained in the data (size of objects, linearity, etc.). In the most general case, a priori information is not available to us.

The multiresolution support of an image is computed in several steps:

• 1.  
  The wavelet transform of the image is computed.
• 2.  
  Booleanization of each scale leads to the multiresolution support.
• 3.  
  A priori knowledge can be introduced by modifying the support.

The last step depends on the knowledge we have of our images. For instance, if we know there is no interesting object smaller or larger than a given size in our image, we can suppress, in the support, anything that is due to that kind of object. This can often be done conveniently by the use of mathematical morphology. In the most general setting, we naturally have no information to add to the multiresolution support.

The wavelet transform of an image by an algorithm such as the à trous method (Shensa 1992) produces, at each scale j, a set {w_j}. This has the same number of pixels as the image, and thus this wavelet transform is a redundant one. The original image c₀ can be expressed as the sum of all the wavelet scales and the smoothed array c_p:

A pixel at position (x,y) can also be expressed as the sum all the wavelet coefficients at this position, plus the smoothed array:

The multiresolution support will be obtained by detecting at each scale the significant coefficients. The multiresolution support is defined by

Given stationary Gaussian noise, in order to define if w_j is significant it suffices to compare w_j(x,y) to kσ_j. Often k is chosen as 3. If w_j(x,y) is small, it is not significant and could be due to noise. If w_j(x,y) is large, it is significant. That is,

So we need to estimate, in the case of Gaussian noise models, the noise standard deviation at each scale.

The appropriate value of σ_j in the succession of wavelet planes is assessed from the standard deviation of the noise σ_I in the original image and from study of the noise in the wavelet space. This study consists of simulating an image containing Gaussian noise with a standard deviation equal to 1 and taking the wavelet transform of this image. Then we compute the standard deviation σ^e_j at each scale. We obtain a curve σ^e_j as a function of j, which gives the behavior of the noise in the wavelet space. (Note that if we had used an orthogonal wavelet transform, this curve would be linear.) As a result of the properties of the wavelet transform, we have σ_j = σ_Iσ^e_j. The standard deviation of the noise at a scale j of the image is equal to the standard deviation of the noise of the image multiplied by the standard deviation of the noise of the scale j of the wavelet transform.

The algorithm to create the multiresolution support is as follows:

• 1.  
  We compute the wavelet transform of the image.
• 2.  
  We estimate the noise standard deviation at each scale. We deduce the statistically significant level at each scale.
• 3.  
  Booleanization of each scale leads to the multiresolution support.
• 4.  
  Modification using a priori knowledge (if desired).

Step 4 is optional. A typical use of a priori knowledge is the suppression of isolated pixels in the multiresolution support in the case in which the image is obtained with a point‐spread function (PSF) of more than 1 pixel. Then we can be sure that isolated pixels are residual noise that has been detected as significant coefficients. If we use a pyramidal wavelet transform method (Starck et al. 1996) or a nonlinear multiresolution transform, the same method can be used.

In order to visualize the support, we can create an image S defined by

Figure 1 shows such a multiresolution support visualization of an image of galaxy NGC 2997.

Zoom In Zoom Out Reset image size

We simulate an image containing around 100 objects (stars and galaxies) and a background. Then we derive a set of images by adding a Gaussian noise with different standard deviation. The original image had a mean value of 1431 ADU (analog‐to‐digital units), minimum and maximum values of 1236 and 18,365, and a σ (standard deviation) of 546. The added Gaussian values, using the original image's mean value, had standard deviations ranging from 2 (an almost imperceptible change) to 5000 (so that noise totally dominates).

Results presented in Table 1 show the robustness of the standard deviation estimation from the multiresolution support. The error percentages are defined as

Note how, with imperceptible change (σ small, close to 2), most methods perform badly relative to this admittedly tiny value of σ. Broadly speaking, for sufficiently large noise, the various methods do well. We see, however, that the multiresolution support does remarkably well, for all values of the noise.

Figure 2 shows a simulated image mostly of galaxies, as used by Freudling & Caulet (1993) and Caulet & Freudling (1993). We added Gaussian noise, which leads to varying signal‐to‐noise ratios (SNR) (defined as 10×log ₁₀[σ(input image)/σ(added noise)]). Figure 3 shows an image with SNR = −0.002. This image was produced by adding noise of standard deviation arbitrarily chosen to be equal to the overall standard deviation of the original image. For varying SNRs, the following three methods were assessed. First, a median‐smoothed version of the image (which provides a crude estimate of the image background) is subtracted, and the noise is estimated from this difference image. This noise estimate is made more robust by 3 σ clipping. Second, the multiresolution support approach described in this article was used, in which the specification of the multiresolution support was iteratively refined. Third, the noise statistics in a running block or square window are used, and the averages of these values are returned.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 2 gives results obtained. The first column gives the noise as a standard deviation. The second column expresses this as SNR. The third, fourth, and fifth columns give automatic noise estimates using the median/3 σ clipping method, the iterative multiresolution support method, and the block method. The excellent estimation capability of the iterative multiresolution support method can be seen: the fourth column values are very close to the first column values.

The goal of this section is to estimate to what extent an error in noise evaluation is significant for the detection of objects in astronomical images. The methodology used is to take an image with objects of known positions, to add noise, and then to compare the true object list with that found by an object detection program.

Starting from a simulated image containing 116 objects, their positions being known, we add a Gaussian noise with a standard deviation σ_N. Using an astronomical package for source detection, in which the noise standard deviation is a free parameter, it is possible to measure the impact of an error on the standard deviation in the results. Introducing an error between −20% and +20%, we compare the results of the detection using two criteria:

• Percentage of real detected objects:

  
• Percentage of false detections:

  

  Since extended sources are sometimes detected as several smaller objects, we can have also multiple detections for one object.

A positive error is equivalent to an increase in σ_N, which implies detection with a higher threshold. Fewer objects are detected, but they are more likely to be true detections. Conversely, a negative error means decreasing the threshold, thus detecting more objects but making more false detections.

We define a parameter, called ratio, which links together ρ(o) and ρ(e):

Figure 4 shows values of this ratio for various values of standard deviations (σ) of added noise. From Figure 4, it seems that an error between −10% and +10% is acceptable. In this case, most methods are satisfactory. This result is important because it means that for a large survey, in which the amount of data is very great, we do not need a high‐precision estimate of the standard deviation of the noise. In such a situation, the computation time may well be a more important issue than the quality. Obviously, in this case, the multiresolution support, which requires somewhat greater computation time, should not be used.

Zoom In Zoom Out Reset image size

So far we have targeted small‐valued nuisance data. Now we turn our attention to very large‐valued nuisance or artifact data.

The ISOCAM infrared camera is one of the four instruments on board ISO, which was launched successfully on 1995 November 17. It operates in the 2.5–17 μm range and was developed by the French Service d'Astrophysique of CEA Saclay. For one observation, the ISOCAM instrument delivers a set of frames, so several measurements are available for a sky position. This is important because cosmic‐ray hits seriously affect the detector. Figure 5 (top) shows an example of what a detector pixel sees. Data averaging without taking into account these glitches (cosmic‐ray hits) would lead to a large error. Figure 5 (middle) shows the result using a deglitching algorithm (Starck et al. 1997). A parameter of this method is the standard deviation of the noise. As can be seen in Figure 5 (middle), many weak glitches have not been eliminated because of overestimation of noise variance. The effect of the error is a large number of false source detections in the final calibrated image. In Figure 5 (bottom), the same deglitching algorithm was used, but we used as input the noise variance found by the multiresolution method. The result is better, and all false detections have disappeared in the calibrated image.

Zoom In Zoom Out Reset image size

If the noise in the data I is Poisson, the transform

acts as if the data arose from a Gaussian white noise model (Anscombe 1948), with σ = l, under the assumption that the mean value of I is large. The arrival of photons, and their expression by electron counts, on CCD detectors may be modeled by a Poisson distribution. In addition, there is additive Gaussian readout noise. The Anscombe transformation has been extended to take this combined noise into account. The generalization of the variance stabilizing Anscombe formula is derived as (Murtagh, Starck, & Bijaoui 1995)

where α is the gain, σ is the standard deviation, and g is the mean of the readout noise.

If the readout noise standard deviation is not well known, we can use the variance stabilization transform to estimate it. This is done by using the fact that we know that the standard deviation in t(I) should be equal to 1. The algorithm is as follows:

• 1.  
  Set n to 0, S_min to 0, and S_max to σ(I).
• 2.  
  Set r_n to (S_min + S_max)/2.
• 3.  
  Compute the transform of I with a readout noise standard deviation equal to r_n.
• 4.  
  Estimate the standard deviation of the noise σ_S in t(I) by the method described in § 3.
• 5.  
  If σ_S<1, then S_min = r_n; otherwise S_max = r_n.
• 6.  
  If S_max - S_min>, then n = n + l and go to step 2.
• 7.  
  r_n = (S_min + S_max)/2 is a good estimation of the readout noise standard deviation.

Starting from the same image as in § 4.1, we simulate a set of images by adding Poisson noise and Gaussian noise. Using the algorithm described above, the Gaussian component was sought. The gain was taken as equal to one in all simulations. Results of these simulations are presented in Table 3. As we can see, the standard deviation of the readout noise can be found with good accuracy (less than 1%) when the readout noise is greater than 45. For very low readout noise, the error can become relatively high.