ON COMPUTING UPPER LIMITS TO SOURCE INTENSITIES

Our study is carried out in the context of background-contaminated detection of point sources in photon counting detectors, as in X-ray astronomical data. We set up the problem for the case of uncomplicated source detection (i.e., ignoring source confusion, intrinsic background variations, and instrumental effects such as vignetting, detector efficiency, PSF structure, bad pixels, etc). However, the methodology we develop is sufficiently general to apply in complex situations.

There is an important, subtle, and often overlooked distinction between an upper limit and the upper bound of a confidence interval, and the primary goal of this paper is to illuminate this difference. The confidence interval is the result of inference on the source intensity, while the upper limit is a measure of the power of the detection process. We can precisely state this difference in the context of an example. Suppose that we have a typical case of a source detection problem, where counts are collected in a region containing a putative or possible source and are compared with counts from a source-free region that defines the background. If the source counts exceed the threshold for detection, the source is considered to be detected. This detection threshold is usually determined by limiting the probability of a false detection. If the threshold were lower, there would be more false detections. Given this setup, we might ask how bright must a source be in order to ensure detection. Although statistically there are no guarantees, the upper limit is the minimum brightness that ensures a certain probability of detection. Critically, this value can be computed before the source counts are observed. It is based on two probabilities, (1) the probability of a false detection which determines the detection threshold and (2) the minimum probability that a bright source is detected. Although the upper limit is primarily of interest when the observed counts are less than the detection threshold, it does not depend on the observed counts. This is in sharp contrast to a confidence interval for the source intensity, that is typically of the form "source intensity estimate plus or minus an error bar," where the estimate and error bars depend directly on the observed source counts in the putative source region. Of course, the functional form of the confidence interval may be more complicated than in this example, especially in low count settings, but any reasonable confidence interval depends on the source counts, unlike upper limits. The fact that upper limits do not depend on the source counts while confidence intervals do should not be viewed as an advantage of one quantity or the other. Rather it reflects their differing goals. Upper limits quantify the power of the detection procedure and confidence intervals describe likely values of the source intensity. These distinctions are highlighted with illustrative examples in Section 3.1.1.

To formalize discussion, suppose that a known source has an intrinsic intensity in a given passband of λ_S and that the background intensity under the source is λ_B. Further suppose that the source is observed for a duration τ_S and that n_S counts are collected, and similarly, a separate measurement of the background could be made over a duration τ_B and n_B counts are collected. If the background counts are collected in an area r times the source area,⁷ we can relate the observed counts to the expected intensities,

$$\begin{eqnarray} \nonumber n_B| (\lambda _B,r,\tau _B) &\sim & {\rm Poisson}(r\tau _B\lambda _B) \\ n_S| (\lambda _S,\lambda _B,\tau _S) &\sim & {\rm Poisson}\Big (\tau _S(\lambda _S+\lambda _B)\Big) \,, \end{eqnarray} \tag{ 1 }$$

in the background and source regions respectively, where n_B and n_S are independent. For simplicity, we begin by assuming that λ_B is known.

Confidence intervals give a set of values for the source intensity that are consistent with the observed data. They are typically part of the inference problem for the source intensity. The basic strategy is to compute an interval of parameter values so that on repeated observations a certain proportion of the intervals contain the true value of the parameter. It is in this "repeated-observation" sense that a classical confidence interval has a given probability of containing the true parameter. Bayesian credible intervals have a more direct probabilistic interpretation. They are computed by deriving the posterior probability distribution of the source intensity parameter given the observed counts and finding an interval with nominal probability of containing the true rate (see, e.g., Loredo 1992; van Dyk et al. 2001; Kashyap et al. 2008). In summary, confidence intervals are frequentist in nature meaning that they are interpreted in terms of repeated observations of the source. Credible intervals, on the other hand, are Bayesian in nature meaning that they represent an interval of a certain posterior (or other Bayesian) probability.

2.2.1. Confidence Intervals

From a frequentist point of view randomness stems only from data collection—it is the data, not the parameters that are random. Often we use a 95% interval, but intervals may be at any level and we more generally refer to an L% confidence interval. Thus, the proper interpretation of a given interval is

• L% of experiments (i.e., observations) with intervals computed in this way will result in intervals that contain the true value of the source intensity.

In frequentist terms, this means that in any given experiment one cannot know whether the true source intensity is contained in the interval but if the experiment is repeated a large number of times, about L% of the resulting intervals will contain the true value. Strictly speaking, the more colloquial understanding that there is an L% chance that the "source intensity is contained in the reported confidence interval," is incorrect.

Put another way, a confidence interval for the source intensity gives values of λ_S that are plausible given the observed counts. Suppose that λ_B = 3 and that for each value of λ_S, we construct an interval ${\cal I}(\lambda _S)$ of possible values of the source counts that has at least an L% chance: ${\rm Pr}(n_S\in {\cal I}(\lambda _S)|\lambda _S,\lambda _B,\tau _S) \ge L\%$. Once the source count is observed and assuming λ_B is known, a confidence interval can be constructed as the set of values of λ_S for which the observed count is contained in ${\cal I}(\lambda _S)$,

$$\begin{equation} \lbrace \lambda _S: n_S\in {\cal I}(\lambda _S)\rbrace \,. \end{equation} \tag{ 2 }$$

In repeated observations, at least L% of intervals computed in this way cover the true value of λ_S.

The frequency coverage of confidence intervals is illustrated in Figure 1, where the confidence intervals of Garwood (1936) for a Poisson mean are plotted as boxes of width equal to the interval for various cases of observed counts (see van Dyk's discussion in Mandelkern 2002). For given values of λ_S and λ_B, the probability of the possible values of the observed counts can be computed using Equation (1); for each of these possible values, there will be a different confidence interval. Thus, the confidence intervals themselves have their own probabilities, which are represented as the heights of the boxes in Figure 1. Because these are 95% confidence intervals, the cumulative heights of the boxes that contain the true value of λ_S in their horizontal range must be at least 0.95. It is common practice to only report confidence intervals for detected sources so that only intervals corresponding to n_S above some threshold are reported. Unfortunately, this upsets the probability that the interval contains λ_S upon repeated observations. Standard confidence intervals are designed to contain the true value of the parameter (say) 95% of time, i.e., in 95% of data sets. If some of the confidence intervals are taken away (i.e., are not reported because, e.g., the counts are too small), there is no reason to expect that 95% of those remaining will contain the true value of the parameter. This is because instead of summing over all the values of n_S to get a probability that exceeds 95%, we are summing over only those values that are greater than the detection threshold. This results in a form of Eddington bias and is discussed in detail in Section 3.4.

Zoom In Zoom Out Reset image size

2.2.2. Credible Intervals

In a Bayesian setting, probability is used to quantify uncertainty in knowledge and in this regard parameters are typically viewed as random quantities. This distinction leads to a more intuitive interpretation of the credible interval. A credible interval at the L% level, for example, is any interval that contains the true value of the parameter L% of the time according to its posterior distribution. (See Park et al. 2008 for a discussion on interval selection.) Thus, from a Bayesian perspective, it is proper to say that there is an L% chance that the source intensity is contained in the reported credible interval. The corresponding credible intervals look similar to the confidence intervals in Figure 1, at least in high count scenarios.⁸

So far we have considered a very simple problem with only one unknown parameter, λ_S. The situation is more complicated if there are unknown nuisance parameters, such as λ_B. In this case, frequency based intervals typically are constructed using asymptotic arguments and/or by conditioning on ancillary statistics that yield a conditional sampling distribution that does not depend on the nuisance parameter. Identifying ancillary statistics can be a subtle task and the resulting intervals may not be unique. Bayesian intervals, can be constructed using a simple and clear principle known as marginalization. For example, if λ_B is unknown, the marginal posterior distribution of λ_S is simply

$$\begin{eqnarray} &&p(\lambda _S| n_S, n_B,\tau _S, \tau _B, r) = \int p(\lambda _S, \lambda _B| n_S, n_B,\tau _S,\tau _B, r) d\lambda _B\, \nonumber\\ &&\quad = \int p(\lambda _S| \lambda _B, n_S, \tau _S) p(\lambda _B| n_B, \tau _B, r) d\lambda _B\,. \end{eqnarray} \tag{ 3 }$$

Credible intervals for λ_S are computed just as before, but using the marginal posterior distribution. ⁹

We emphasize that neither confidence nor credible intervals directly quantify the detection sensitivity of an experiment. To do this, we consider the detection problem in detail, which from a statistical point of view is a test of the hypothesis that there is no source emission in the given energy band,¹⁰ i.e., a test of λ_S = 0. Formally, we test the null hypothesis that λ_S = 0 against the alternative hypothesis that λ_S>0. The test is conducted using a test statistic that we denote by ${{ {\cal S}}}$. An obvious choice for ${ {\cal S}}$ is the counts in the source region, n_S; larger values of n_S are indicative of a detection of a source, since they become increasingly less likely to have been obtained as a random fluctuation from the background. Other choices for ${ {\cal S}}$ are the signal-to-noise ratio (as in the case of sliding-cell local-detect algorithms; see Section 4) or the value of the correlation of a counts image with a basis function (as in the case of wavelet-based algorithms) or a suitably calibrated likelihood-ratio test statistic (as in the case of γ-ray detectors like Fermi; see Mattox et al. 1996). The count in the source region is an example of a test statistic that is stochastically increasing¹¹ in λ_S. For any fixed ${ {\cal S}}^\star$, λ_B, τ_S, τ_B, and r, the probability that the test statistic ${ {\cal S}}$ is less than the threshold ${ {\cal S}}^\star$ decreases as λ_S increases, i.e., ${\rm Pr}({{ {\cal S}}} \le { {\cal S}}^\star | \lambda _S,\lambda _B,\tau _S,r)$ decreases as λ_S increases. We assume that ${ {\cal S}}$ is stochastically increasing in λ_S throughout.¹²

Because larger values of ${ {\cal S}}$ indicate a source, we need to determine how large ${ {\cal S}}$ must be before we can declare a source detection. This is done by limiting the probability of a false detection, also known as a Type I error. Thus, the detection threshold ${ {\cal S}}^\star$ is the smallest value such that

$$\begin{equation} {\rm Pr}({ {\cal S}}> { {\cal S}}^\star | \lambda _S=0,\lambda _B,\tau _S,\tau _B,r)\le \alpha, \end{equation} \tag{ 4 }$$

where α is the maximum allowed probability of a false detection¹³ and we declare a detection if the observed value of ${ {\cal S}}$ is strictly greater than ${ {\cal S}}^\star$:

• If ${ {\cal S}}\le { {\cal S}}^\star$ we conclude there is insufficient evidence to declare a source detection.
• If ${ {\cal S}}> { {\cal S}}^\star$ we conclude there is sufficient evidence to declare a source detection.

We call ${ {\cal S}}^\star$ the α-level detection threshold and sometimes write ${ {\cal S}}^\star (\alpha)$ to emphasize its dependance on α (see Figure 2). Note that α is a bound on the probability of a Type I error; the actual probability of a Type I error is given by the probability on the left-hand side of Equation (4). Due to the discrete nature of the Poisson distribution, the bound is generally not achieved and the actual probability of a Type-I error is less than α.

Zoom In Zoom Out Reset image size

Although its role in the definition of the detection threshold indicates that it is viewed as the more important concern, a false detection, also known as a "false positive," is not the only type of error. A false negative, or Type II error, occurs when a real source goes undetected (see Figure 3). The probability of a false negative is quantified through the power of the test to detect a source as a function of its intensity,

$$\begin{equation} \beta (\lambda _S) = {\rm Pr}({ {\cal S}}> { {\cal S}}^\star | \lambda _S,\lambda _B, \tau _S, \tau _B, r) \,. \end{equation} \tag{ 5 }$$

Equation (5) gives the probability of a detection. For any λ_S>0, this is the power of the test or one minus the probability of a false negative.¹⁴ For λ_S = 0, Equation (5) gives the probability of a false detection (cf. Equation (4)) and consequently, β(0) ⩽ α. This reflects the trade-off in any detection algorithm: the compromise between minimizing the number of false detections against maximizing the number of true detections. That is, if the detection threshold is set low enough to detect weaker sources, the algorithm will also produce a larger number of false positives that are actually background fluctuations. Conversely, the more stringent the criterion for detection, the smaller the probability of detecting a real source (this is illustrated by the location of the threshold ${ {\cal S}}^\star$ that defines both α and β in Figure 3). Note that although our notation emphasizes the dependence of the power on λ_S, it also depends on λ_B, τ_S, τ_B, and r.

Zoom In Zoom Out Reset image size

The power calculation is shown for the simple Poisson case in Figure 4, where β(λ_S) is plotted for different instances of λ_B and for different levels of the detection threshold ${ {\cal S}}^\star$. As expected, stronger sources are invariably detected. For a given source intensity, an increase in the background or a larger detection threshold (i.e., lower α) both cause the detection probability to decrease. In a typical observation, the background and the detection threshold are already known, and thus it is possible to state precisely the intensity λ_S at which the source will be detected at a certain probability. We may set a certain minimum probability, β_min of detecting a "bright" source by setting the exposure time long enough so that any source with intensity greater than a certain pre-specified cutoff has probability β_min or more of being detected. Conversely, we can determine how bright a source must be in order to have probability β_min or more of being detected with a given exposure time. This allows us to define an upper limit on the source intensity by setting a minimum probability of detecting the source. This latter calculation is the topic of Section 3.1 and the basis of our definition of an upper limit.

Zoom In Zoom Out Reset image size

Power calculations are generally used to determine the minimum exposure time required to ensure a minimum probability of source detection (see Appendix B). In Section 3, we use them to construct upper limits.

In astronomy, upper limits are inextricably bound to source detection: by an upper limit, an astronomer means

• The maximum intensity that a source can have without having at least a probability of β_min of being detected under an α-level detection threshold.

or conversely,

• The smallest intensity that a source can have with at least a probability of β_min of being detected under an α-level detection threshold.

Unlike a confidence interval, the upper limit depends directly on the detection process and in particular on the maximum probability of a false detection and the minimum power of the test, that is on α and β_min, respectively. In this way, an upper limit incorporates both the probabilities of a Type I and a Type II error. Formally, we define the upper limit, ${{ {\cal U}}}(\alpha, \beta _{\rm min})$ to be the smallest λ_S such that

$$\begin{equation} {\rm Pr}({{ {\cal S}}} > { {\cal S}}^\star (\alpha) | \lambda _S,\lambda _B, \tau _S, \tau _B, r)\ge \beta _{\rm min}. \end{equation} \tag{ 6 }$$

Commonly used values for β_min throughout statistics are 0.8 and 0.9. If β_min ≈ 1, ${ {\cal U}}(\alpha, \beta _{\rm min})$ represents the intensity of a source that is unlikely to go undetected, and we can conclude that an undetected source is unlikely to have intensity greater than ${ {\cal U}}(\alpha, \beta _{\rm min})$.

The simplest example occurs in the hypothetical situation when λ_B is known to be zero and there is no background observation. In this case, we set ${ {\cal S}}=n_S$ and note that Pr(n_S>0|λ_S = 0, λ_B = 0, τ_S) = 0 so the detection threshold is zero counts and we declare a detection if there is even a single count. (Recall, we declare a detection only if ${ {\cal S}}$ is strictly greater than ${ {\cal S}}^\star$.) The upper limit in this case is the smallest value of λ_S with probability of detection greater than β_min. Figure 5 plots Pr(n_S>0|λ_S, λ_B = 0, τ_S) as a function of τ_Sλ_S, thus giving ${ {\cal U}}(\alpha, \beta _{\rm min})$ for any given τ_S and every value of β_min. Note that the upper limit decreases in inverse proportion to τ_S.

Zoom In Zoom Out Reset image size

When λ_B is greater than zero but well determined and can be considered to be known, the detection threshold using ${ {\cal S}}=n_S$ is given in Equation (4). With this threshold in hand we can determine the maximum intensity a source can have with significant probability of not producing a large enough fluctuation above the background for detection. This is the upper limit.

In particular, ${ {\cal U}}(\alpha, \beta _{\rm min})$ is the largest value of λ_S such that ${\rm Pr}(n_S\le { {\cal S}}^\star (\alpha) | \lambda _S,\lambda _B, \tau _S) > 1-\beta _{\rm min}$. This is illustrated for three different values of β_min (panels) and three different values of α (line types) in Figure 6. Note that the upper limit increases as β_min increases and as α decreases.

Zoom In Zoom Out Reset image size

3.1.1. Illustrative Examples

So far our definition of an upper limit assumes that there are no unknown nuisance parameters, and in particular that λ_B is known. Unfortunately, the probabilities in Equations (4) and (6) cannot be computed if λ_B is unknown. In this section, we describe several strategies that can be used in the more realistic situation when λ_B is not known precisely.

The most conservative procedure ensures that the detection probability of the upper limit is greater than β_min for any possible value of λ_B. Generally speaking, the larger λ_B is, the larger λ_S must be in order to be detected with a given probability, and thus the larger the upper limit. Thus, a useful upper limit requires a finite range, Λ_B, to be specified for λ_B. Given this range, a conservative upper limit can be defined as the smallest λ_S that satisfies

$$\begin{equation} \inf _{\lambda _B\in \Lambda _B} {\rm Pr}({{ {\cal S}}} > { {\cal S}}^\star (\alpha) | \lambda _S,\lambda _B, \tau _S, \tau _B, r) \ \ge \ \beta _{\rm min}. \end{equation} \tag{ 7 }$$

(We use the term infimum ($\inf$) rather than minimum to allow for the case when the minimum may be on the boundary of, but outside, the range of interest. It is the largest number that is smaller than all the numbers in the range. For instance, the minimum of the range {x>0} is undefined, but the infimum is 0.) Unfortunately, unless the range of values Λ_B is relatively precise, this upper limit will often be too large to be useful.

In practice, there is a better solution. The background count provides information on the likely values of λ_B that should be used when computing the upper limit. In particular, the distribution of λ_B given n_B can be computed using standard Bayesian procedures¹⁵ and used to evaluate the expected detection probability,

$${\selectfont\fontsize{9.4}{10}{\begin{equation} \beta (\lambda _S) = \int {\rm Pr}({{ {\cal S}}} > { {\cal S}}^\star (\alpha) | \lambda _S,\lambda _B, \tau _S,\tau _B, r) p(\lambda _B| n_B, \tau _B,r)d\lambda _B, \end{equation}}} \tag{ 8 }$$

where ${ {\cal S}}^\star (\alpha)$ is the smallest value such that

$$\begin{eqnarray} \int {\rm Pr}({ {\cal S}}> { {\cal S}}^\star (\alpha) | \lambda _S &=& 0,\lambda _B,\tau _S,\tau _B, r) p(\lambda _B| n_B, \tau _B, r)d\lambda _B\nonumber\\ &\le& \ \alpha. \end{eqnarray} \tag{ 9 }$$

The upper limit is then computed as the smallest λ_S that satisfies β(λ_S) ⩾ β_min. Unlike the upper limit described in Section 3.1, these calculations require data, in particular, n_B. For this reason, we call the smallest λ_S that satisfies Equation (8) the background count conditional upper limit or bcc upper limit. An intermediate approach that is more practical than using Equation (7) but more conservative than using Equation (8) is to simply compute a high percentile of p(λ_B|n_B), perhaps its 95th percentile. The procedure for known λ_B can then be used with this percentile treated as the known value of λ_B. This is a conservative strategy in that it assumes a nearly worst case scenario for the level of background contamination.

As an illustration, suppose the uncertainty in λ_B given the observed background counts can be summarized in the posterior distribution plotted in the left panel of Figure 7. This is a gamma posterior distribution of the sort that typically arises when data are sampled from a Poisson distribution. Using this distribution, we can compute ${ {\cal S}}^\star$ for any given value of α as the smallest value that satisfies Equation (9); the results are given for three values of α in the legend of Figure 7. We then use Equation (8) to compute β(λ_S) as plotted in the right panel of Figure 7. The upper limit can be computed for any β_min using these curves just as in Figure 6.

Zoom In Zoom Out Reset image size

Although the form of a confidence interval makes it tempting to use its upper bound in place of an upper limit, this is misleading and blurs the distinction between the power of the detection procedure and the confidence with which the flux is measured. As an illustration, we have computed the upper limit for each of the nine panels in Figure 1 (using the true values of λ_B reported in each panel). The results are plotted as solid vertical lines in Figure 8. Note that unlike the upper bound of the confidence interval the upper limit does not depend on n_S and takes on values that only depend on λ_B. Using the upper bound of the confidence interval sometimes overestimates and sometimes underestimates the upper limit. In all but one of the nine cases with the highest λ_S/λ_B (i.e., λ_S = 5, λ_B = 1) the upper limit for λ_S is larger than λ_S. Of the nine cases, this is the one that is mostly likely to result is a source detection and is the only one with a probability of detection greater than β_min = 0.8.

Zoom In Zoom Out Reset image size

Alternatively, we can compute the value of β_min required for the upper bound of Garwood's confidence interval to be interpreted as an upper limit. Figure 9 does this for the three values of λ_B used in the three columns of Figures 1 and 8. Consider how the upper bound of Garwood's confidence interval increases with n_S in Figure 8. Each of these upper bounds can be interpreted as an upper limit, but with an increasing minimum probability of a source detection, β_min. The three panels of Figure 8 plot how the required β_min increases with n_S for three values of λ_B. Note that a source with intensity equal to the upper bound of Garwood's confidence interval can have a detection probability as low as 20% or essentially as high as 100%. Thus, the upper bound does not calibrate the maximum intensity that a source can have with appreciable probability of going undetected in any meaningful way.

Zoom In Zoom Out Reset image size

As mentioned in Section 2.2.1, it is common practice to only report a confidence interval for detected sources. Selectively deciding when to report a confidence interval in this way can dramatically bias the coverage probability of the reported confidence interval.¹⁶ We call this bias a statistical selection bias. Note that this is similar to the Eddington bias (Eddington 1913) that occurs when intensities are measured for sources close to the detection threshold. For sources whose intrinsic intensity is exactly equal to the detection threshold, the average of the intensities of the detections will be overestimated because downward statistical fluctuations result in non-detections and thus no intensity measurements. In extreme cases, this selection bias can lead to a nominal 95% confidence interval having a true coverage rate of well below 25%, meaning that only a small percentage of intervals computed in this way actually contain λ_S. As an illustration, Figure 10 plots the actual coverage of the nominal 95% intervals of Garwood (1936) for a Poisson mean when the confidence intervals are only reported if a source is detected with α = 0.05. These intervals are derived under the assumption that they will be reported regardless of the observed value of n_S. Although alternative intervals could in principle be derived to have proper coverage when only reported for detected sources, judging from Figure 10 such intervals would have to be wider than the intervals plotted in Figure 1. It is critical that if standard confidence intervals are reported they be reported regardless of the observed value of n_S and regardless of whether a source is detected.¹⁷

Zoom In Zoom Out Reset image size

As discussed in Section 1, the detection threshold is sometimes used as an upper limit. Under certain circumstances, this can be justified under our definition of an upper limit. Suppose that some invertible function $f({ {\cal S}})$ can be used as an estimate of λ_S and that for any λ_B, the sampling distribution of $f({ {\cal S}})$ is continuous with median equal to λ_S. That is,

$$\begin{equation} {\rm Pr}(f({ {\cal S}}) > \lambda _S | \lambda _S, \lambda _B, \tau _S, \tau _B, r) = 0.5. \end{equation} \tag{ 10 }$$

Because Equation (10) holds for any λ_S, it holds for $\lambda _S= f({ {\cal S}}^\star)$. That is,

$$\begin{equation} {\rm Pr}(f({ {\cal S}}) > f({ {\cal S}}^\star) | \lambda _S=f({ {\cal S}}^\star), \lambda _B, \tau _S, \tau _B, r) = 0.5. \end{equation} \tag{ 11 }$$

Inverting f and integrating over p(λ_B|n_B, τ_B, r), we have

$$\begin{eqnarray} \int {\rm Pr} ({ {\cal S}}> { {\cal S}}^\star | \lambda _S &=& f({ {\cal S}}^\star), \lambda _B, \tau _S, \tau _B, r) p(\lambda _B|n_B,\tau _B,r)d\lambda _B\nonumber\\ &=& 0.5. \end{eqnarray} \tag{ 12 }$$

Comparing Equation (12) with Equation (8), we see that $f({ {\cal S}}^\star) = { {\cal U}}(\alpha, \beta _{\rm min} = 0.5)$. Thus, if f is an identity function the detection threshold is an upper limit. Although the assumption that the sampling distribution of $f({ {\cal S}})$ has median λ_S for every λ_B, is unrealistic in the Poisson case, it is quite reasonable with Gaussian statistics. Even if this assumption holds, $f({ {\cal S}}^\star)$ is a weak upper limit in that half the time a source with this intensity would go undetected and there is a significant chance that sources with intrinsic intensity larger than $f({ {\cal S}}^\star)$ would remain undetected.

It should be emphasized that even when the detection threshold is used as an upper limit, it is not an "upper limit on the counts," but an upper limit on the intrinsic intensity of the source. The counts are an observed, not an unknown quantity. There is no need to compute upper bounds, error bars, or confidence intervals on known quantities. It is for the unknown source intensity that these measures of uncertainty are useful.

Our analysis of upper limits and confidence interval assumes that the observables are photon counts that we model using the Poisson distribution. However, the machinery we have developed is applicable to any process that uses a significance-based detection threshold. Here, we briefly set out a general recipe to use in more complicated cases. For complex detection algorithms, some of the steps may require Monte Carlo methods.

• 1.  
  Define a probability model for the observable source and background data set given the intrinsic source and background strengths, λ_S and λ_B, respectively. For the simple Poisson case, this is defined in Equation (1). In many applications, these could be approximated using Gaussian distributions. It is typically required that a background data set be observed but in some cases λ_B may be known a priori. The source could be a spectral line, or an extra model component in a spectrum, or possibly even more complex quantities that are not directly related to the intensity of a source.
• 2.  
  Define a test statistic ${ {\cal S}}$ for measuring the strength of the possible source signal. In the simple Poisson case, we set ${ {\cal S}}=n_S$. The "source" could be a spectral line or any extra model component in a spectrum, or a more complex quantity that is not related to the intensity of a source.
• 3.  
  Set the maximum probability of a false detection, α, and compute the corresponding α-level detection threshold, ${ {\cal S}}^\star$. Although ${ {\cal S}}^\star$ depends on λ_B, we can compute the expected ${ {\cal S}}^\star$ by marginalizing over λ_B using p(λ_B|n_B) when λ_B is not known exactly. Likewise, if λ_S is defined as a function of several parameters, the same marginalization procedure can be used to marginalize over any nuisance parameters. In this case, we typically marginalize over p(η|n_B) or perhaps p(η|n_S, n_B), where η is the set of nuisance parameters. In this regard, we are setting α to be a quantile of the posterior predictive distribution of ${ {\cal S}}$, under the constraint that λ_S = 0; see Gelman et al. (1996) and Protassov et al. (2002).
• 4.  
  Compute the probability of detection, β(λ_S), for the adopted detection threshold ${ {\cal S}}^\star$.
• 5.  
  Define the minimum probability of detection at the upper limit, β_min. Traditionally, β_min = 0.5 has been used in conjunction with α = 0.003 in astronomical analysis (see Section 3.5).
• 6.  
  Compute the smallest value of λ_S such that β(λ_S) ⩾ β_min. This is the upper limit.