The AstroStat Slog

As mentioned before, background subtraction plays a big role in astrophysical analyses. For a variety of reasons, it is not a good idea to subtract out background counts from source counts, especially in the low-counts Poisson regime. What Bayesians recommend instead is to set up a model for the intensity of the source and the background and to infer these intensities given the data.

For instance, suppose as before, that C counts are observed in a region of the image that overlaps a putative source, and B counts in an adjacent, non-overlapping region that is mostly devoid of the source and which is r times larger in area and exposure than the source region. Further suppose that a fraction f of the source falls in the so-called source region (typically, f~0.9) and a fraction g falls in the background region (we strive to make g~0). Then the observed counts can be written as Poisson realizations of intensities,

C = Poisson(φ_S) ≡ Poisson(f ­ θ_S + θ_B) , and
B = Poisson(φ_B) ≡ Poisson(g ­ θ_S + r ­ θ_B) ,

where the subscripts denote the model intensities in the source (S) or background (B) regions.

The joint probability distribution of the model intensities,

p(φ_S φ_B | C B) dφ_S dφ_B

can be rewritten in terms of the interesting parameters by transforming the variables,

≡ p(θ_S θ_B | C B) J(φ_S φ_B ; θ_S θ_B) d θ_S d θ_B

where J(φ_S φ_B ; θ_S θ_B) is the Jacobian of the coordinate transformation, and thus

= p(θ_S θ_B | C B) (r ­ f – g) d θ_S d θ_B .

The posterior probability distribution of the source intensity, θ_S, can be derived by marginalizing this over the background intensity parameter, θ_B. A number of people have done this calculation in the case f=1,g=0 (e.g., Loredo 1992, SCMA II, p275; see also van Dyk et al. 2001, ApJ 584, 224). The general case is slightly more convoluted, but is still a straightforward calculation (Kashyap et al. 2008, AAS-HEAD 9, 03.02); but more on that another time.