Help identifying image noise

[Danny Uncanny7]

I am a bit out of my field here, so please excuse me if I use the wrong terminology or don't know obvious theory. This is part of a problem I am working on to implement optimization algorithms on automated star detection for sky surveys. Specifically for dim, clustered stars in noisy fields. That is somewhat beside the point. What I need help with specifically is with the objective function.

So I have a grayscale sky image of arbitrary resolution. It has some noise in it for various reasons. I come up (or rather the optimizer does) with a model of the star locations, brightnesses, and point spread function. If this model is accurate and subtracted from the image, I should end up with a flat background field with noise. I am assuming for the moment that the noise is mostly Gaussian, additive and unaffected by mean pixel value. I am also ignoring pixel saturation. I might have to revisit those assumptions later.

Anyways, so the point is, I want to quantify the accuracy of the model based on the residual flat field it leaves behind. Basically I want to get a residual where the pixel probability function is the same for every pixel independent of position. I thought about taking the Gaussian assumption and just getting the mean square error of the histogram against a fitted Gaussian. But the histogram doesn't account for any spatial randomness. I guess what I mean is that I am trying to find some measure of spatial entropy... if that is the right term.

So long story short, does anyone have any suggestions or ideas of what kind of measurements I can use to tell if my image pixels are spatially uncorrelated. I want to tell if what I am looking at is just noise and no signal. Computationally faster calculations are always better too.

[eta oin shrdlu]

This is a little out of my field too, but here are my comments anyway. The most obvious thing to check first is the noise image's autocorrelation function. Point (m,n) in the autocorrelation image estimates the degree of correlation between points (x,y) and points (x+m,y+n) in your residual image; if all of the pixels are zero-mean uncorrelated random variables then all points except for (0,0) should be small. Note that this only tells you about correlations between uniform translations of the image; correlations involving reflections, rotations, and other more complicated transformations are not reliably detected this way, though you can use similar ideas to test these hypotheses. Also note that you can compute the autocorrelation quickly using FFTs.

I think you're also wise to question the assumptions of additive Gaussian noise; this may easily be wrong for a number of reasons: e.g., noise levels may depend on pixel value, and saturation may chop off one tail of the noise. But it's probably a reasonable place to start.

[Danny Uncanny7]

Well I imagine I don't need to autocorrelate the whole image, just neighbouring pixels. Maybe I could do a moving 5x5 autocorrelation for each pixel. But what specific value could I use as a measure for the whole image? Should I look for any average correlation value above a threshold?

[mfb]

The autocorrelations are a global property for the image, you get a single value for each possible shift of pixels. If you know the size of your stars (let's say they do not occupy more than a nxn grid of pixels), maybe you can ignore the autocorrelation for larger shifts. But that is something you can do later. The autocorrelation for small shifts tells you something about the idea that each pixel should be uncorrelated from its neighbours.

[eta oin shrdlu]

Yeah, like mfb says, the autocorrelation is a function of the entire image. Here is a simple example. I start with a 100x100 random Gaussian noise image, each pixel independent, zero mean, and identically distributed, and compute its autocovariance (since the means are zero, this is the same as the autocorrelation):

autocovariance for pixel-independent Gaussian noise

autocovariance uncorrelated.png (8.45 KiB) Viewed 1099 times

I've shifted the image so that the center pixel represents zero offset; that's the bright one-pixel spot in the center. (I've also adjusted the brightness so that you can see the fluctuations in the rest of the image.) Note that the bright spot is a single pixel wide, indicating no significant correlations between neighboring pixels.

If I blur the image to model a nontrivial point-spread function (and so introducing spatial correlations between neighbors), the resulting autocovariance looks like this:

autocovariance for Gaussian-blurred noise

autocovariance correlated.png (3.75 KiB) Viewed 1099 times

Now the central peak is several pixels wide (you may have to zoom in to see this), giving a measure of the size of the neighboring-pixel spatial correlations (in this case, introduced by the circular Gaussian blur I applied, with standard deviation of 1 pixel).

It is also useful to look at the Fourier transform of these autocovariances; these are what's called the spectral densities. For the original (unblurred) image, the spectral density is flat across the whole image:

spectral density for pixel-independent Gaussian noise

spectral-density uncorrelated.png (9.6 KiB) Viewed 1099 times

For the blurred image, the spectral density shows power only at low spatial frequencies, indicating correlation between neighbors:

spectral density for Gaussian-blurred noise

spectral-density correlated.png (2.92 KiB) Viewed 1099 times

[Danny Uncanny7]

Okay, so what measurement would you recommend distilling from the autocorrelation image? A mean value cutoff threshold maybe?

[eta oin shrdlu]

The most important number would probably be the width of the central peak. The height of the central peak is related to the noise variance; its width is related to the spatial correlation length, which is what you're interested in. Other large peaks or obvious structure in the autocorrelation image or the spectral density may hint at other properties of the correlations.

[Danny Uncanny7]

Hmmm, that seems kind of vague, especially considering it might be a single smudge in a noisy field, not necessarily whole image blurring, the correlation over the whole image might be very small. I'll work on it though.

One of the other approaches I was thinking of was some sort of entropy measure. Shannon entropy applied to an image typically seems to be taken as sum of p*log(p) over the histogram. This would however give the same value for a clear image and the same image with the pixel locations randomly scattered. I was wondering if there is a way to use a spatial equivalent.

The other brute force approach I came up with was to take the standard deviation of the standard deviations of a series of image sub samples. These would be 9 pixel region taken about each pixel. I think that it would work, but it doesn't feel mathematically sound. I generally don't like algorithms with integer values (other than 2) in them.

[starslayer]

I want to tell if what I am looking at is just noise and no signal.

This question might be obvious and dumb, but are you doing dark subtraction, bias frames, and flat fielding for these images? The biggest sources of noise in astronomical images are dark current, readout noise, and the chip not being evenly illuminated/evenly sensitive. Poisson noise is below all of those, and is what's left after you subtract out the first three things.

Also, what are you planning to do with these images? Just counting stars? Or some kind of photometry with them? Your S/N ratio requirements are going to be different for those tasks, and when reducing/analyzing your data, you should keep those requirements in mind.

[Danny Uncanny7]

I am just picking out the stars and giving that information back to my astrophysics friend. I am not dealing with any calibration or hardware, so I assume that this data has been properly conditioned to remove all of the noise that is easily removable. In fact most of the data I am dealing with right now is made-up data for validation purposes. It's easier to tell if your algorithm is working if you know what result it should find.

The neat thing about researching this kind of problem is that it has so many applications. Identifying image noise is import not only in astrophysics, but face recognition, autofocus design, QR code reading, and probably tons of other fields. And noise in general probably has some importance to every single field of research.

[D.B.]

Can I just ask, do you have any more info on what the "optimizer" you refer to in the original post doing? Can it give you probabilistic data of star location, numbers, etc, or is it essentially a black box that spits out a list of 'stars found here, here and here'?

See, as you've explained it, you're trying to determine whether or not the remaining pixels in an image (following subtraction) is best described by all noise or by the presence of another star. You've got a probabilistic model for the noise (gaussian at the moment, possibly getting improved later). If you had access to a probabilistic model of what a star looks like then this would be really really useful information - it's going to be much easier to find a faint signal in noise if you have a good idea what that signal should be. You could also think about testing some standard methods of model selection to see how they do. It's going to be easier to defend using some custom method of finding remaining stars that you come up with yourself later if you've already considered stuff like BIC, Laplace's method, etc, and can show why yours is better.

(EDIT: Of course, perhaps what your friends want is an entirely separate method of spotting missed stars that doesn't use their model, so that systematic errors caused by inaccuracies in said model can be eliminated...)

Also, what are you planning to do with these images? Just counting stars? Or some kind of photometry with them? Your S/N ratio requirements are going to be different for those tasks, and when reducing/analyzing your data, you should keep those requirements in mind.

I think this would be a valuable thing to find out.

[Danny Uncanny7]

Can I just ask, do you have any more info on what the "optimizer" you refer to in the original post doing? Can it give you probabilistic data of star location, numbers, etc, or is it essentially a black box that spits out a list of 'stars found here, here and here'?

See, as you've explained it, you're trying to determine whether or not the remaining pixels in an image (following subtraction) is best described by all noise or by the presence of another star. You've got a probabilistic model for the noise (gaussian at the moment, possibly getting improved later). If you had access to a probabilistic model of what a star looks like then this would be really really useful information - it's going to be much easier to find a faint signal in noise if you have a good idea what that signal should be. You could also think about testing some standard methods of model selection to see how they do. It's going to be easier to defend using some custom method of finding remaining stars that you come up with yourself later if you've already considered stuff like BIC, Laplace's method, etc, and can show why yours is better.

(EDIT: Of course, perhaps what your friends want is an entirely separate method of spotting missed stars that doesn't use their model, so that systematic errors caused by inaccuracies in said model can be eliminated...)

Also, what are you planning to do with these images? Just counting stars? Or some kind of photometry with them? Your S/N ratio requirements are going to be different for those tasks, and when reducing/analyzing your data, you should keep those requirements in mind.

I think this would be a valuable thing to find out.

I'm not sure what you mean by a probabilistic star model. Probabilistic in the general sense that the PSF and each location, and brightness have a probability function? In that case no. I'm not entirely sure how I would treat that kind of model, or evaluate it. I treat the star model as deterministic. And the noise is also deterministic in the sense that there is a measured image. I'm trying to see if it comes from a probabilistic model, but the specific pixels have fixed values. I'm not sure I understand what the benefit of this kind of model would be. In my understanding at least, the stars location and brightness are deterministic values which are measured with some probabilistic reading. So I do have a model of what a star looks like in that sense. One star would be a linear variation of the point spread function. But suppose a cluster of stars is misrepresented by the model, or a star is off location a few pixels. The residual smudge or error would certainly not resemble a missed star. I want to find some objective function that just knows that is bad, and the optimizer takes that as feedback.

The optimizers are just generic numerical algorithms for analysing complex highly dimensional problems and minimizing some resulting value. Usually using a combination of gradient information, point comparison, and space exploration. So in this case, for k stars, I have k, intensities, 2*k subpixel coordinates, and a point spread function of several dimensions. k itself is also a variable that would have to optimized in another outer iteration. I suppose I could apply BIC to the algorithm here. It's the first time I've read about it. I had a few other ideas to prevent overcounting stars (the obvious bad solution would be to have as many stars as pixels with a 1 pixel point spread function). The problem could easily have over a thousand dimensions of optimization space. The optimizer moves through that space sampling the objective function and "learning" until it hits some criteria that tells it that it has found a good enough solution. However, I am guessing that are a lot of short cuts available for this problem. The results are mostly independent between stars. The problem is mainly with clusters of stars, which conventional star finding algorithms do very poorly with.

The potential applications would be for cataloguing data from sky surveys like the Sloan Digital Sky Survey. I think astronomers use that kind of data for things like models of how galaxies form, but I'm not sure what else.

Edit: without letting the problem sprawl out of hand, I guess what I really want to figure out is this: for an arbitrary image (or any data set really), what is the probability that each pixel was independently selected from the same probability distribution? I think that the answer will in some way come down to the autocorrelation, or some form of it.

[Yakk]

Edit: without letting the problem sprawl out of hand, I guess what I really want to figure out is this: for an arbitrary image (or any data set really), what is the probability that each pixel was independently selected from the same probability distribution? I think that the answer will in some way come down to the autocorrelation, or some form of it.

A1) Reading this as "what is the probability that a given probability distribution would produce this arbitrary image?"

Nearly zero, assuming your probability distribution isn't boring.

Any fixed image is very, very unlikely to be produced by a particular probability distribution. If your per-pixel distribution was uniform, then the odds of a pure white image are equal to the odds that your exact image would be produced.

A2) Lets say you actually want an answer to the question asked -- then it depends on your priors. Suppose you have a belief that a certain kind of systematic error could occur -- local smears, creases, circles or what have you. You build a model where a family of those have a certain chance of occurring, and there is a certain chance of a uniform uncorrelated image. Then from this, we could work out what the probability that your image was produced by each model in particular.

However, this depends hugely on what your priors are. And, in fact, you need to build a complete set of models of alternative explanations you want to examine (or a statistical model that effectively contains them) if you want an exact answer.

Your priors, however, are highly important. Usually these priors are implicit in the statistical analysis you end up doing -- maybe you sort your images implicitly into "more ordered" and "less ordered" categories, use that to implicitly inform your priors...

A3) You actually don't want an answer to the question. Rather, you want a tool that you can use to say "well, this image is pretty darn random" or "this image seems to have some kind of pattern", and it outputs some value that gets bigger if there seems to be some kind of pattern that the tool recognizes.

[Danny Uncanny7]

Well I ended up thinking long and hard about A2 and getting something like A3.
Here's the little snippet of Python code that I wrote:

Code: Select all

def objective(image):
    dim=shape(image)
    measure=0
    for i in xrange(1,dim[0]):
        for ii in xrange(1,dim[1]):
            measure+=(image[i,ii]*image[i-1,ii-1]+image[i-1,ii]*image[i,ii-1])**2
        #end
    #end
    return measure
#end


Image in this case should have a zero mean. Measure is an unnormalized value that represents the variance in the covariance of neighboring pixels. Neighbours in this case are taken as those diagonal from each other mostly to simplify the calculation. If the whole thing doesn't run too sluggishly, I might add in the vertical and horizontal neighbours which I think will give a cleaner result, but doubles the number of calculations.

My basic thought is that any model that affects pixels spatially will have a similar effect on neighbouring pixels. I can't think of any model applicable to observations that will similarly affect distant pixels, but not neighbours. So the correlations of neighbouring pixels give an indication of the presence of anything other than noise. The distribution of those correlations gives a pretty good indication of whether they are all linearly independent or not. And rather than summing the squares, I could even save an array of those values to give location specific feedback to the optimizer.