Weighted Poisson Uncertainties

I recently ran across a rather awkward mathematical problem. I'm trying to make a histogram of photon arrival phases, complete with an uncertainty on the number of photons in each bin. Normally this is done by just taking the square root of the number of photons, which is at least approximately right based on Poisson statistics. But in my problem — data from the Fermi space telescope, which is not very good at localizing low-energy gamma rays — the photons are weighted: for each photon I have a probability that it really came from the source. So the values in the histogram should be the total probability. But what should the uncertainty be? The short version is: the square root of the sum of the squares of the weights.

ETA: This is in the literature, without justification as far as I can tell. See below.
I've worked out the details as best I can in an ipython notebook. It's a pain to include ipython notebooks directly into blogs at this point, so I've posted it on the ipython notebook viewer. I'll just add a few comments here.

The math works out, at least based on what I think is a reasonable mathematical model for the problem. I also ran statistical tests based on the same model, and they confirm that it works well enough for reasonably large numbers of photons, provided that the weight isn't dominated by a small number of heavily-weighted photons. And the formula is convenient and easy to remember.

Edited to add: I found a paper, Espinoza et al. 2012, that uses this formula, and it cites Pletsch et al. 2012, which simply states that this is the right formula. There's also Guillemot et al. 2012, which uses a Monte Carlo approach - essentially using all the bins as a population from which to draw the photon weights. But since Lucas Guillemot is the second author on the Espinoza et al. paper, I'm guessing they've figured out that the simple formula is the way to go.

No comments: