XMM-Newton Science Analysis System


imweightadd (tools-1.68.1) [xmmsas_20230412_1735-21.0.0]

Null-hypothesis probability distribution for a weighted sum of Poissonian variates

The approach followed in the present section is essentially that of Fay and Feuer [1]. See also Stewart []

Let $X$ be a random deviate which follows a Poissonian probability distribution about the expectation value $a = \langle X \rangle$. Although the Poisson distribution itself is only defined for integer $X$, one can find the following continuous `envelope function' to the Poisson values:

  $$e(X;a) = \frac{a^X \exp(-a)}{\Gamma(X+1)},$$ (3)

where $\Gamma$ is the gamma function

$$\Gamma(a) = \int_0^\infty dt \ e^{-1} \ t^{a-1}.$$

Note that both the expectation value $\langle X^2 \rangle$ and variance $\langle X \rangle-\langle X \rangle^2$ of the function $e$ are identical to those of the corresponding discrete Poisson distribution, namely both equal to $a$. Now, given $N$ random deviates $X_i$, each of which follows a distinct Poisson distribution about its average $a_i$, let us form the weighted sum

$$Y = \sum_{i=1}^N w_i X_i.$$

The expectation value $\mu$ of $Y$ is

$$\mu = \langle Y \rangle = \sum_{i=1}^N w_i \langle X_i \rangle = \sum_{i=1}^N w_i a_i;$$

the variance $v$ can in similar fashion be shown to equal

$$v = \sum_{i=1}^N w_i^2 a_i.$$

The probability function $p(Y)$ only has values where all the $X_i$ are integer and generally speaking may be expected to be a messy-looking and intractable function. However, recall that for purposes of source detection we are not interested in $p$ but in the integral $P$ of $p$ from a particular sample $y$ of $Y$ to infinity. $P$ is stepwise continuous, the steps becoming smaller and denser as $y$ increases. The envelope function with the same expectation value and variance as $Y$ is given by

  $$e(Y;\mu,v) = \frac{(k\mu)^{(kY)} exp(-k\mu)}{\Gamma(kY+1)},$$ (4)

where

$$k = \mu / v.$$

In plain English, what comparison of equations 3 and 4 suggests is that a weighted sum of Poisson variates behaves approximately like a single Poisson distribution with the same average and variance. Therefore since, for a single Poisson variate $X$, the probability $P$ for $X$ to be greater than some sampled value $x$ is given by

$$P(X\ge x;a) = 1 - Q(x,a)$$

(where $Q$ is the incomplete gamma function we met with in section 3.2) we postulate that the equivalent expression for $Y$ is approximately given as follows:

  $$P(Y\ge y;\mu,v) \sim 1 - Q(ky,k\mu).$$ (5)

Equation 5 is used both in the present program and in boxdetect to estimate the null-hypothesis probability distribution of the weighted sum of Poissonian images.