Histogram and fit to Poisson distribution
kronion
Thu, 05/24/2012 - 01:05 am
my task is as follows: For calibration, I take microscopic images of a fluorescent solution. Now we want to know the intensity I. The intensity values should follow a Poisson distribution around the expectation value (especially in the case of low intensity, above that it is OK approximated with Gaussian). Since it is not well approximated by a Gaussian distribution in the case of low intensities, I want to find the expectation value using the Poisson distribution.
My question is now: How do I fit a Poisson distribution to the histogram?
I made a cumulative histogram from the image, normalized to 1 (Histogram/CUM/P/B=4 myimage,myimage_Hist) and tried to fit the CDF (http://en.wikipedia.org/wiki/Poisson_distribution):
<br />
Function Poisson_CUM(w,k) : FitFunc<br />
Wave w<br />
Variable k<br />
<br />
return GammaInc(floor(k+1),w[0])/factorial(floor(k))<br />
End
Function Poisson_CUM(w,k) : FitFunc<br />
Wave w<br />
Variable k<br />
<br />
return GammaInc(floor(k+1),w[0])/factorial(floor(k))<br />
End
But this does not work, so I am wondering whether there is another way of doing it?
What I have done is to generate approximate artificial Poisson data using Knuth's method (cited in your reference) for a specified mean value. Then I created a histogram of that data. Finally that was followed by a curve fit to the Poisson distribution:
Wave w
Variable x
//CurveFitDialog/ These comments were created by the Curve Fitting dialog. Altering them will
//CurveFitDialog/ make the function less convenient to work with in the Curve Fitting dialog.
//CurveFitDialog/ Equation:
//CurveFitDialog/ f(x) = exp(-mu)*(mu^x) / factorial(x)
//CurveFitDialog/ End of Equation
//CurveFitDialog/ Independent Variables 1
//CurveFitDialog/ x
//CurveFitDialog/ Coefficients 1
//CurveFitDialog/ w[0] = mu
return exp(-w[0])*(w[0]^x) / factorial(x)
End
You have to be careful that the fitting takes place only on the integer values of the data and of the Histogram. Do not use the /CUM flag on the Histogram. I attach a figure showing an analytic Poisson PDF with mean value 6, the histogram of 200 data variates simulated using that mean, and the curve fit of the Histogram data to the above function. The fit shown gave mu = 6.1284 . I used a limited number of "measurement" samples to demonstrate that the fit is not exactly along the analytic result, and that the histogram can deviate a lot from the ideal case.
If I have mis-interpreted your question, please provide more details.
May 24, 2012 at 08:45 am - Permalink
I have a feeling that Kronion is misinterpreting something, or has not described his problem well. Light intensity values should be Poisson-distributed because measuring light intensity is effectively a process of counting photons, and counting processes have Poisson-distributed uncertainty. But that would be the uncertainty for the intensity value measured at any one pixel in the image, not the histogram for the image as a whole.
So perhaps Kronion could try again to describe exactly what the problem is.
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
May 24, 2012 at 09:15 am - Permalink
thanks a lot for the replies!
To clarify the problem: I take a confocal image of a fluorophore solution. Thus, every pixel should have the same intensity, more precisely they should be Poisson-distributed around the same expectation value EXP(I). To determine EXP(I), I want to treat the 1024x1024 pixels as 1mio independent measurements of EXP(I). I look at the histogram of the image and want to find EXP(I) by fitting a Poisson PDF of CDF. Attached are images of the histogram and the cumulative histogram.
However, the best estimator should actually be the mean value (I indicated it with the black line) - but it does not look like that to me --> I might have a wrong understanding of the maths.
My fitting problem is probably related to making the histogram in a good way (choose bins so that I have integers, as you pointed out). I am used to Matlab and there you can specify an array containing the bins - have to find the equivalent now in Igor :) or at least a way to see the scale (so far I can only change the scale, but can only guess how the automatically generated histogram-wave is scaled.
May 24, 2012 at 11:01 am - Permalink
Igor has both StatsPoissonCDF and StatsPoissonPDF that can be wrapped in a user-defined fitting function. As stated before, you will need to be careful that the histogram bin values fall on integer values, or StatsPoissonXXX will give bad values.
I wonder if you need some scaling factor, or an X offset. It seem like a Poisson curve with peak at such large numbers is more nearly Gaussian.
Your picture looks very much like a Gamma distribution, which could be fit using a user-define fitting function that wraps our StatsGammaPDF or StatsGammaCDF functions. Since they are continuous distributions, fitting will be more straightforward.
John Weeks
WaveMetrics, Inc.
support@wavemetrics.com
May 25, 2012 at 09:06 am - Permalink
John has beat me to it in replying, but I will add my two cents. I interpret your last posting to mean that you are essentially trying to determine the statistics of the number of fluorophores per pixel. You clearly will have to establish a scaling factor between intensity and fluorophore density. From the shape of your pdf histogram it clearly is closest to a Poisson distribution with a rather low mean, probably around 4. As John noted it is also qualitatively not a great fit to the Poisson distribution. Since your intensity measurements give smoothly distributed values, it will be an interesting exercise for you to adapt them to any discrete density function, like the Poisson. My guess is that the measurement scaling factor should be treated as another fitting parameter.
I think the Igor binning parameters are unambiguous, as they are set by your choice. If you use a discrete pdf that has integer arguments, I think it is most convenient to set your scaled histogram bin size to one.
May 25, 2012 at 09:29 am - Permalink