May 18, 2019 - Statistics of Poisson fields - analytics
Here we look at statistics of Poisson fields. Essentially there is not anything interesting - the result is that the skewness and higher order moments diverge for a Poisson field without any pixelization. This is not surprising; in a realistic experiment not every [math]\displaystyle{ \ell }[/math] is measured, which truncates the sums and leads to a finite result. This is also analogous to the statement that total variance due to perfect white noise is infinite, which is not surprising.
Definitions
Looked up the terms.
It turns out there is Flux density, labeled as [math]\displaystyle{ S }[/math], units Jy = energy/area/time/frequency. Useful, if you observe a point source - how many Watts do you observe in your 1 cm2 detector in a given frequency range?
Then there is Spectral radiance, labeled as [math]\displaystyle{ B }[/math], units Jy/sr. On top of above you add the fact that your detector is observing only a finite solid angle. This is the quality that is related to temperature as in Mar 31, 2019 - Brightness temperature
Full sky
We will work in the limit where each source is much smaller than the pixel. Let's label the size of each pixel [math]\displaystyle{ A }[/math].
Starting from the single source distribution
For starters, let's imagine we have one source randomly on the sky, with its flux density [math]\displaystyle{ S_1 }[/math] independently chosen from a random distribution [math]\displaystyle{ P(S) }[/math], normalized such that [math]\displaystyle{ \int dS P(S) = 1 . }[/math]
The spectral radiance (can be converted to temperature through the conversion factor) is then nonzero only in a single pixel, where it is
[math]\displaystyle{ B(x) = \frac{S_1}{A} . }[/math]
Notice that the units match. The full sky average of [math]\displaystyle{ B^k }[/math] is
[math]\displaystyle{ \langle B^k \rangle_\mathrm{sky} = \frac{1}{4\pi} \int d^2x\, B^k(x) = \frac{A}{4\pi} \times \frac{S_1^k}{A^k} }[/math]
Expectation value of this, now averaged also over the source flux, is
[math]\displaystyle{ \langle B^k \rangle = \frac{A}{4\pi} \times \frac{1}{A^k} \times \langle S_1^k \rangle = \frac{A^{1-k}}{4\pi} \int dS P(S) S^k }[/math]
Two things are notable:
- The mean spectral radiance does not depend on pixelization
- For [math]\displaystyle{ k \gt 1 }[/math], the expectation value diverges as [math]\displaystyle{ A \rightarrow 0 }[/math]. The same result is obtained in the continuum limit with delta function sources - the spectral radiance increases above all bounds. This makes sense, as we have a constant flux source being imagined into a smaller and smaller angular area.
In reality, over the whole sky we expect to have total of [math]\displaystyle{ \langle N_\mathrm{sources} \rangle }[/math] sources. If the pixels are so small that we do not expect two sources in the same pixel, signal from various pixels is uncorrelated and this leads to
[math]\displaystyle{ \langle B^k \rangle = \frac{\langle N_\mathrm{sources}\rangle A^{1-k}}{4\pi} \int dS P(S) S^k }[/math]
Defining now the number of sources (over the whole sky) having the flux in [math]\displaystyle{ [S, S+dS] }[/math],
[math]\displaystyle{ \frac{dN}{dS} = \langle N_\mathrm{sources} \rangle P(S) , }[/math]
we can finally write
[math]\displaystyle{ \langle B^k \rangle = \frac{A^{1-k}}{4\pi} \int dS \frac{dN}{dS} S^k . }[/math]
Distribution in a pixel
We can also start by focusing on the pixel as the basic unit, instead of the sources. To do that, let's define another probability measure [math]\displaystyle{ F(S, A) }[/math], that now say what is the probability that a pixel of size [math]\displaystyle{ A }[/math] contains sources with total flux density in [math]\displaystyle{ [S, S+dS] }[/math]. It can be related to [math]\displaystyle{ P(S) }[/math] by
[math]\displaystyle{ F(S, A) = \tilde P(1\, \mathrm{source}) P(S) + \frac{1}{2}\tilde P(2\, \mathrm{sources}) \int d\Delta P(S-\Delta) P(\Delta) + \frac{1}{6}\tilde P(3\, \mathrm{sources}) \int d\Delta d\Delta' P(S-\Delta) P(\Delta-\Delta') P(\Delta') + \dots , }[/math]
where [math]\displaystyle{ \tilde P }[/math] denotes the probability that pixel in question contains a given number of sources (so some Poisson process). Notice that [math]\displaystyle{ F(S, A) }[/math] now depends on the size of the pixel, unlike [math]\displaystyle{ P(S) }[/math].
We can repeat the calculation from the previous section and get
[math]\displaystyle{ \langle B^k \rangle = \int dS F(S, A) \frac{S^k}{A^k} . }[/math]
In the limit of small pixels, only the first term is important and
[math]\displaystyle{ \tilde P(1\, \mathrm{source}) = \frac{\langle N_\mathrm{sources} \rangle A}{4\pi} , }[/math]
which reproduces the previous result.
Starting from the map
In simulations, we have [math]\displaystyle{ N_\mathrm{pix} }[/math] of pixels with values [math]\displaystyle{ B_i }[/math]. From these we can estimate
[math]\displaystyle{ F(S,A) = \frac{1}{N_\mathrm{pix}} \sum_\mathrm{pix} \delta(S - B_i A) . }[/math]
Plugging in above, this leads to
[math]\displaystyle{ \langle B^k \rangle = \int dS \frac{1}{N_\mathrm{pix}} \sum_\mathrm{pix} \delta (S-B_i A) \frac{S^k}{A^k} = \frac{1}{N_\mathrm{pix}} \sum_i B_i^k . }[/math]
Notice that there is still the dependence on the size of the pixel - if we had 2000 pixels and a single source in half of them (corresponding to unity surface brightness), we would have
[math]\displaystyle{ \langle B^k \rangle = \frac{1000 \times 1^k}{2000} = \frac{1}{2} }[/math]
Making the pixels half the size, there would be now 4000 pixels with a source in 1000 of them, with two units of surface brightness. This would correspond to
[math]\displaystyle{ \langle B^k \rangle = \frac{1000 \times 2^k}{4000} = \frac{2^k}{4} }[/math]