# COMPARING AND COMBINING CMB DATASETS

## 1. INTRODUCTION

Undetected systematic errors are one of the main obstacles to using cosmic microwave background (CMB) measurements to constrain cosmological models. One of the best ways to address this problem is to compare experiments whose data sets overlap both in sky coverage and angular scale, to see whether they are consistent.<sup>3</sup><sup>3</sup>3 The best way to address the problem is clearly to design CMB experiments to be more immune to systematic errors in the first place. The next best thing to do is carefully examine the raw time-ordered data from an experiment for specific forms of systematic errors that may be expected, and for general signs of systematic errors for data removal or correction. The cross-check between experiments that are discussed in this paper are by no means a substitute for this, but rather a way of catching additional systematic errors that have slipped through the cracks and not been detected by the team that reduced the raw data.  If they are consistent, a useful second step is to combine them into a single sky map retaining all their cosmological information.

Both of these steps are simple for maps with identical resolution and beam shape. For e.g. the 53 and 90 GHz COBE DMR maps (Bennett et al. 1996), comparing (1) was done by subtracting the maps and checking whether the difference was consistent with pure noise, whereas combining (2) was done by simply averaging the two maps, weighting pixels by their inverse variance. Unfortunately, both steps are usually more complicated. MAP, Planck and most current experiments have different angular resolution in different channels. Many current experiments probe the sky in an even more complicated way, with e.g. double beams, triple beams, interferometric beams or complicated elongated software-modulated beams. Correlated noise further complicates the problem.

Despite these difficulties, precision comparisons between different experiments are crucial. Some of the best evidence so far for detection of CMB fluctuations comes from the success of such comparisons in the past — between FIRS and DMR (Ganga et al. 1993), Tenerife and DMR (Lineweaver et al. 1995), MSAM and Saskatoon (Knox et al. 1998; hereafter K98), two years of Python data (Ruhl et al. 1995), three years of Saskatoon data (Tegmark et al. 1996a) and two flights of MSAM (Inman et al. 1997). The fact that many of these non-COBE data sets were contaminated by systematic errors made the success of these cross-checks even more encouraging.

Data sets are currently growing rapidly in number, size and quality, often overlapping. It is therefore quite timely to develop methods that generalize both steps (1) and (2) to arbitrary experiments. This is the purpose of the present Letter. In the larger context of CMB data analysis, this is important between the steps of mapmaking (Wright et al. 1996; Wright 1996; Tegmark 1997a) and power spectrum estimation (Tegmark 1997b; Bond et al. 1998) in the pipeline.

## 2. Notation

Let us first establish some notation that will be used throughout this paper. Consider a pixelized CMB sky map at some resolution consisting of $`m`$ numbers $`x_1,\mathrm{},x_m`$, where $`x_i`$ is the temperature in the $`i^{th}`$ pixel. Suppose two experiments $`i=1,2`$ have measured $`n_i`$ numbers $`y_1,\mathrm{},y_{n_i}`$, each probing some linear combination of the sky temperatures $`x_i`$. Grouping these numbers into vectors $`𝐱`$, $`𝐲_1`$ and $`𝐲_2`$ of length $`m`$, $`n_1`$ and $`n_2`$, we can generally write<sup>4</sup><sup>4</sup>4 This assumes both that the experimental data have perfect linearity, and that there are no non-zero experimental offsets. Ideally, experimentalists should model and remove both nonlinearities and offsets as part of their data reduction, thereby making equation (1) applicable to their final data product. If offsets of an unknown amplitude remain nonetheless, multiplication of a data set $`𝐲=\mathrm{𝐀𝐱}+𝐧+offsets`$ by a matrix $`𝐏`$ that projects out these offsets will produce a new data set $`𝐲_1=\mathrm{𝐏𝐲}`$ that satisfies equation (1), merely with the slightly more complicated noise vector $`𝐧_1=\mathrm{𝐏𝐧}`$ and with $`𝐀_1=\mathrm{𝐏𝐀}`$. A detailed example of this procedure can be found in de Oliveira et al. (1998). 

$$𝐲_1=𝐀_1𝐱+𝐧_1,𝐲_2=𝐀_2𝐱+𝐧_2$$
(1)   
for some known matrices $`𝐀_i`$ incorporating the beam shapes and some random noise vectors $`𝐧_i`$ with zero mean $`(𝐧_i=\mathrm{𝟎})`$. We will refer to $`𝐱`$ as the “true sky”. $`𝐲_1`$ and $`𝐲_2`$ can be either time-ordered data or some linear combination thereof, for instance pixelized maps. It is sometimes useful to define the larger matrices and vectors

$$𝐀\left(\genfrac{}{}{0pt}{}{𝐀_1}{𝐀_2}\right),𝐲\left(\genfrac{}{}{0pt}{}{𝐲_1}{𝐲_2}\right),𝐧\left(\genfrac{}{}{0pt}{}{𝐧_1}{𝐧_2}\right),$$
(2)   
which allows us to rewrite equation (1) as

$$𝐲=\mathrm{𝐀𝐱}+𝐧.$$
(3)   
Let us write the noise covariance matrix as

$$𝐍\mathrm{𝐧𝐧}^t=\left(\begin{array}{cc}𝐍_1& 𝐍_{12}\\ 𝐍_{12}^t& 𝐍_2\end{array}\right),$$
(4)   
where $`𝐍_1𝐧_1𝐧_1^t`$, $`𝐍_2𝐧_2𝐧_2^t`$ and $`𝐍_{12}𝐧_1𝐧_2^t`$.

We derive useful consistency tests in §2 both for the special case of identical observations $`(𝐀_1=𝐀_2)`$ and for the general case $`𝐀_1𝐀_2`$, then show how to combine data sets without destroying information in §3.

## 3. Comparing data sets: are they consistent?

### 3.1. The “null-buster” test

Let us first consider the simple case where the two data sets measure the same thing, i.e., $`𝐀_1=𝐀_2`$. This often applies for two different channels of the same experiment at the same frequency. We can then form a difference map $`𝐳𝐱_1𝐱_2`$, which in the absence of systematic errors should consist of pure noise.

Consider the null hypothesis $`H_0`$ that such a data set $`𝐳`$ consists of pure noise, i.e., $`𝐳=\mathrm{𝟎}`$, $`\mathrm{𝐳𝐳}^t=𝐍`$ for some noise covariance matrix $`𝐍`$. Suppose we have reason to suspect that the alternative hypothesis $`H_1`$ is true, where $`𝐳=\mathrm{𝟎}`$, $`\mathrm{𝐳𝐳}^t=𝐍+\mathrm{SS}`$ for some signal covariance matrix $`\mathrm{SS}`$, and want to try to rule out $`H_0`$ by using a test statistic $`q`$ that is a quadratic function of the data:

$$q𝐳^t\mathrm{𝐄𝐳}^t=\text{tr}\{\mathrm{𝐄𝐳𝐳}^t\}.$$
(5)   
Depending on whether $`H_0`$ or $`H_1`$ is true, the mean of $`q`$ will be $`q_0=\text{tr}\mathrm{𝐄𝐍}`$ or $`q_1=\text{tr}\mathrm{𝐄𝐍}+\text{tr}𝐄\mathrm{SS}`$, respectively. If $`𝐳`$ has a multivariate Gaussian probability distribution<sup>5</sup><sup>5</sup>5 We will only make the assumption of Gaussianity for the CMB and the detector noise, not for the systematic errors. In fact, systematics, such as foreground signals, data spikes and atmospherics seldom have a Gaussian probability distribution, but they vanish under the null hypothesis that we are trying to rule out. , then $`q`$ will have a variance $`(\mathrm{\Delta }q)^2=2\text{tr}\mathrm{𝐄𝐍𝐄𝐍}`$ if $`H_0`$ is true. Therefore the quantity $`\nu (qq_0)/\mathrm{\Delta }q`$ gives the number of standard deviations (“sigmas”) by which the observed $`q`$-value exceeds the mean expected under the null hypothesis. If we observe $`\nu 1`$, we can thus conclude that $`H_0`$ is ruled out at high significance. Which choice of $`𝐄`$ has the greatest statistical power to reject $`H_0`$ if $`H_1`$ is true, i.e., which $`𝐄`$ maximizes the expectation value

$$\nu \frac{q_1q_0}{\mathrm{\Delta }q}=\frac{\text{tr}𝐄\mathrm{SS}}{\left[2\text{tr}\mathrm{𝐄𝐍𝐄𝐍}\right]^{1/2}}\mathrm{?}$$
(6)   
Since rescaling $`𝐄`$ by a constant leaves $`\nu `$ invariant, let us for simplicity normalize $`𝐄`$ so that the denominator equals unity. We thus want to maximize $`\text{tr}𝐄\mathrm{SS}`$ subject to the constraint that $`\text{tr}\mathrm{𝐄𝐍𝐄𝐍}=1/2`$. Using the method of Lagrange multipliers with $`L=\text{tr}𝐄\mathrm{SS}\lambda \text{tr}\mathrm{𝐄𝐍𝐄𝐍}/2`$ and differentiating $`L`$ with respect to the components of $`𝐄`$, this gives the solution $`𝐄𝐍^1\mathrm{SS}𝐍^1`$. This leaves us with our optimal “null-buster” statistic<sup>6</sup><sup>6</sup>6 K98 discuss a test using the likelihood ratio, which is the best solution to a slightly different problem. Translated into our notation, it corresponds to the choice $`𝐄=𝐍^1[𝐍+\mathrm{SS}]^1=𝐍^1\mathrm{SS}[𝐍+\mathrm{SS}]^1`$. Note that whereas equation (7) is independent of the normalization of $`\mathrm{SS}`$ (the “shape” of the signal matters, but not its amplitude), the likelihood ratio test requires an assumed amplitude. 

$$\nu \frac{𝐳^t𝐍^1\mathrm{SS}𝐍^1𝐳\text{tr}𝐍^1\mathrm{SS}}{\left[2\text{tr}\left\{𝐍^1\mathrm{SS}𝐍^1\mathrm{SS}\right\}\right]^{1/2}},$$
(7)   
which will rule out the null hypothesis $`H_0`$ with the largest average significance if $`H_1`$ is true. For the special case $`\mathrm{SS}=𝐍`$, we see that this reduces to a standard $`\chi ^2`$-test with $`\nu =(\chi ^2n)/\sqrt{2n}`$, $`\chi ^2𝐳^t𝐍^1𝐳`$. Whenever we have reason to suspect systematic errors of a certain form (producing a signal $`\mathrm{SS}`$), the null-buster test will thus be more sensitive to these systematic errors than the $`\chi ^2`$-test, which is a general-purpose tool. This issue is elaborated in K98, which also provides a useful general discussion of consistency tests.

### 3.2. An example: calibration errors

The null-buster is useful for comparing two CMB maps $`𝐲_1`$ and $`𝐲_2`$ that have the same shape, beam size and pixelization. Let $`\mathrm{SS}`$ denote the expected contribution to $`𝐲_i𝐲_{i}^{}{}_{}{}^{t}`$ from CMB fluctuations, i.e., $`\mathrm{SS}=𝐀_1\mathrm{𝐂𝐀}_1^t=𝐀_2\mathrm{𝐂𝐀}_2^t`$, where $`𝐂`$ is the map covariance matrix

$$𝐂_{ij}=\underset{ł=2}{\overset{\mathrm{}}{}}\left(\frac{2ł+1}{4\pi }\right)P_ł(\widehat{𝐫}_i\widehat{𝐫}_j)C_ł;$$
(8)   
the unit vector $`\widehat{𝐫}_i`$ gives the direction towards the $`i^{th}`$ pixel and $`C_ł`$ is the expected or observed CMB power spectrum (the normalization of $`C_ł`$ is irrelevant here — only the shape matters). Now suppose one or both of the data sets $`𝐲_1`$ and $`𝐲_2`$ have a linear calibration error, i.e., are off by some constant multiplicative factors. Consider a difference map of the form $`𝐳𝐲_1r𝐲_2`$ for some factor $`r`$, and plot $`\nu `$ as a function of $`r`$ using equation (7) with the null hypothesis being that $`𝐳`$ is pure noise, $`i.e.`$, that $`𝐍𝐍_1+r[𝐍_{12}+𝐍_{12}^t]+r^2𝐍_2`$. An example is shown in Figure 1. If $`\nu 1`$ for $`r=1`$, then we have a significant detection of signal not common to the two maps. If $`\nu (1)1`$ but $`\nu (r)\mathrm{}<1`$ for some other $`r`$-value, this would show that there is a relative error of $`r`$ in the normalization between the two maps. If the maps are at different frequencies, this could also indicate that they are dominated by foreground contamination with a frequency dependence different from the CMB.

### 3.3. If the beam shape differs

Above we found the null-buster to be a useful consistency test when applied to $`𝐳𝐲_1r𝐲_2`$, since $`𝐀_1=𝐀_2`$ implied that $`𝐲_1𝐲_2`$ should consist of mere noise and be independent of the (a priori unknown) signal. If $`𝐀_1𝐀_2`$, this is no longer true, since the two data sets are not measuring the same thing. We therefore perform our null-buster test with the difference map redefined to be

$$𝐳𝐀_3𝐲_1r𝐀_4𝐲_2$$
(9)   
for some matrices $`𝐀_3`$ and $`𝐀_4`$, and choose these matrices so that the new data vectors $`𝐀_3𝐲_1`$ and $`𝐀_4𝐲_2`$ measure at least approximately the same sky signal. Substituting equation (1) into equation (9) shows that this corresponds to the requirement

$$𝐀_3𝐀_1𝐀_4𝐀_2,$$
(10)   
which would make $`𝐳`$ approximately signal-independent for $`r=1`$. In addition to equation (10), we clearly want $`𝐀_3𝐀_1`$ and $`𝐀_4𝐀_2`$ to have as large rank as possible, to avoid destroying more information than necessary (otherwise say $`𝐀_3=𝐀_4=\mathrm{𝟎}`$ would do the trick).

If one data set can be written as a linear combination of the other (say $`𝐀_2=\mathrm{𝐌𝐀}_1`$ for some matrix $`𝐌`$), then the best choice is clearly $`𝐀_3=𝐌`$, $`𝐀_4=𝐈`$. This is the case for identically sampled maps where $`𝐲_1`$ has higher resolution than $`𝐲_2`$, as well as for cases where a map (say DMR or QMAP) is compared with more complicated weighted averages (say Saskatoon) of the same sky region at the same or lower resolution. This approach was adopted by, e.g., Lineweaver et al. (1995), for comparing the (smoothed) Tenerife data to DMR.

We will now tackle the problem for the general case. Our solution (there may be others) involves performing a signal-to-noise eigenmode analysis (Bond 1994; Bunn & Sugiyama 1995; Tegmark et al. 1996b) three times, in an unusual way. We will first present this procedure with no derivation, then show that it solves our problem. We start by solving the generalized eigenvalue problems 

$`\left[𝐀_1𝐂_1𝐀_1^t\right]𝐁_1`$ $`=`$ $`𝐍_1𝐁_1𝚲_1,`$ (11)   
$`\left[𝐀_2𝐂_2𝐀_2^t\right]𝐁_2`$ $`=`$ $`𝐍_2𝐁_2𝚲_2,`$ (12)   
where the eigenvectors are the columns of the matrices $`𝐁_1`$ and $`𝐁_2`$, normalized so that $`𝐁_i^t𝐍_i𝐁_i=𝐈`$, and the corresponding eigenvectors are elements of the diagonal matrices $`𝚲_1`$ and $`𝚲_2`$, sorted in decreasing order. We then reduce the width of $`𝐁_1`$ and $`𝐁_2`$ by throwing away all eigenvectors with eigenvalues below some cutoff $`\lambda _{min}`$, and define a new smaller data set

$$\stackrel{~}{𝐲}=\left(\genfrac{}{}{0pt}{}{𝐁_1^t𝐲_1}{𝐁_2^t𝐲_2}\right).$$
(13)   
This will have the covariance matrix $`\stackrel{~}{𝐲}\stackrel{~}{𝐲}^t=\stackrel{~}{\mathrm{SS}}+\stackrel{~}{𝐍}`$, where

$`\stackrel{~}{\mathrm{SS}}`$ $`=`$ $`\left(\begin{array}{cc}𝚲_1& 𝐁_1^t𝐀_1\mathrm{𝐂𝐀}_2^t𝐁_2\\ 𝐁_2^t𝐀_2\mathrm{𝐂𝐀}_1^t𝐁_1& 𝚲_2\end{array}\right),`$ (16)   
$`\stackrel{~}{𝐍}`$ $`=`$ $`\left(\begin{array}{cc}𝐈& 𝐁_1^t𝐍_{12}𝐁_2\\ 𝐁_2^t𝐍_{12}^t𝐁_1& 𝐈\end{array}\right).`$ (19)   
We then solve the generalized eigenvalue problem

$$\stackrel{~}{\mathrm{SS}}𝐁=\stackrel{~}{𝐍}𝐁𝚲,$$
(20)   
where the eigenvectors are the columns of the matrix $`𝐁`$, normalized so that $`𝐁^t\stackrel{~}{𝐍}𝐁=𝐈`$, and the corresponding eigenvectors are elements of the diagonal matrix $`𝚲`$, sorted in increasing order. Finally, we reduce the width of $`𝐁`$ by throwing away all eigenvectors with eigenvalues above some cutoff $`\lambda _{max}`$, leaving us with a matrix of the form

$$𝐁=\left(\genfrac{}{}{0pt}{}{𝐁_3}{𝐁_4}\right).$$
(21)   
By choosing a tiny cutoff (say $`\lambda _{max}=10^2`$), we ensure that the transformed data vector $`𝐁^t\stackrel{~}{𝐲}`$ is completely dominated by detector noise, with a for all practical purposes negligible CMB signal. This means that

$$𝐁^t\stackrel{~}{𝐲}=[𝐁_3^t𝐁_1^t𝐀_1+𝐁_4^t𝐁_2^t𝐀_2]𝐱\mathrm{𝟎},$$
(22)   
so we can solve our original problem by defining

$`𝐀_3`$ $``$ $`𝐁_3^t𝐁_1^t,`$ (23)   
$`𝐀_4`$ $``$ $`𝐁_4^t𝐁_2^t.`$ (24)   
Why was the first eigenmode step necessary? We went through the extra trouble of solving equations (11) and (12) and applying the cutoff $`\lambda _{min}`$ because otherwise, a number $`(𝐁^t\stackrel{~}{𝐲})_i`$ in our final data vector could be noise-dominated for two different reasons:

1. Because the signal contribution from $`𝐲_1`$ approximately cancels that from $`𝐲_2`$.
2. Because it is a noise-dominated mode from $`𝐲_1`$ or $`𝐲_2`$ (or some combination thereof).

It is clearly only the first case that interests us when comparing data sets. In the latter case, applying the null-buster only tests for systematic errors internally, within each data set, and this is best done before comparing it with other data sets. We therefore throw away all noise-dominated modes from the individual data sets, choosing say $`\lambda _{min}=1`$, before proceeding to the final eigenvalue problem (20). A lower threshold may be appropriate as well — as long as we choose $`\lambda _{min}\lambda _{max}`$, we know that the lack of signal in $`𝐁^t\stackrel{~}{𝐲}`$ will be due mainly to subtracting the data sets.

## 4. Combining data sets

### 4.1. Combining maps of different beam shape

Suppose that we have performed all the tests described above and conclude that the data sets $`𝐲_1`$ and $`𝐲_2`$ are consistent. We then wish to simplify future calculations by combining the two data sets into a single map $`\stackrel{~}{𝐱}`$, inverting the (usually over-determined) system of linear equations (3). A physically different but mathematically identical problem was solved in Tegmark (1997a), showing that the minimum-variance choice

$$\stackrel{~}{𝐱}=[𝐀^t𝐍^1𝐀]^1𝐀^t𝐍^1𝐲$$
(25)   
retains all the cosmological information of the original data sets. Substituting this into equation (3) shows that the combined map is unbiased ($`\stackrel{~}{𝐱}=𝐱`$) and that the pixel noise $`\text{}\stackrel{~}{𝐱}𝐱`$ has the covariance matrix

$$𝚺\text{}\text{}^t=[𝐀^t𝐍^1𝐀]^1.$$
(26)   
A common special case is that where the two data sets have uncorrelated noise $`(𝐍_{12}=\mathrm{𝟎})`$, simplifying the solution to

$`\stackrel{~}{𝐱}`$ $`=`$ $`𝚺\left[𝐀_1^t𝐍_1^1𝐲_1+𝐀_2^t𝐍_2^1𝐲_2\right],`$ (27)   
$`𝚺`$ $`=`$ $`[𝐀_1^t𝐍_1^1𝐀_1+𝐀_2^t𝐍_2^1𝐀_2]^1.`$ (28)   
An even simpler case occurs if the first data set is already a sky map (say the QMAP map), so that $`𝐀_1=𝐈`$, and we wish to combine it with a more complicated data set covering the same sky region (say the Saskatoon observations). For this case, equation (27) can be rewritten as

$$\stackrel{~}{𝐱}=𝐲_1+𝚺𝐀_2^t𝐍_2^1(𝐲_2𝐀_2𝐲_1),$$
(29)   
which has a simple interpretation. The vector $`𝐀_2𝐲_1`$ is just map 1 convolved with the observing strategy of experiment 2, so the factor $`(𝐲_2𝐀_2𝐲_1)`$ contains only noise. The map $`\stackrel{~}{𝐱}`$ is thus obtained by correcting $`𝐲_1`$ with a pure noise term that partially cancels some of its noisiest modes.

This important case applies also to combining two maps with different angular resolution. For instance, if $`𝐲_1`$ and $`𝐲_2`$ have narrow Gaussian beams of resolution $`\theta _1`$ and $`\theta _2`$, with $`\theta _2>\theta _1`$, then we define $`𝐱`$ to be the true sky map smoothed by $`\theta `$ and have $`(𝐀_2)_{ij}=\mathrm{exp}[\theta _{ij}^2/2\theta ^2]/2\pi \theta ^2`$ where $`\theta _{ij}=\mathrm{cos}^1(\widehat{𝐫}_i\widehat{𝐫}_j)`$ is the angular separation between pixels $`i`$ and $`j`$ and $`\theta =(\theta _2^2\theta _1^2)^{1/2}`$ is the extra smoothing in map 2. Despite occasional claims to the contrary, this shows that two maps at different resolution can be combined without destroying any information, without first degrading the higher resolution down to the lower one by smoothing. Instead, equation (25) will use the lower resolution map $`𝐲_2`$ to improve the accuracy of the large-scale fluctuations in $`𝐲_1`$ (as was done by Schlegel et al. 1998 when combining the DIRBE and IRAS maps), retaining all the information from both maps.

Note that even the case of two identical maps ($`𝐀_1=𝐀_2=𝐈`$) can be non-trivial. Equation (27) shows that the optimal combination is $`\stackrel{~}{𝐱}=𝚺\left[𝐍_1^1𝐲_1+𝐍_2^1𝐲_2\right]`$. This was used in the combined QMAP analysis (de Oliveira-Costa et al. 1998), and reduces to separate averaging for each pixel only if the two maps have vanishing or identical noise correlations.

### 4.2. Combining maps at different frequencies to remove foregrounds

As described in Tegmark (1998) and further elaborated in White (1998), foregrounds can be treated as simply an additional source of noise that is correlated between channels $`(𝐍_{12}\mathrm{𝟎})`$. This means that if we have $`d`$ data sets measured at different frequencies $`\nu _\alpha `$, $`\alpha =1,2,\mathrm{},d`$, each defined by their own matrix $`𝐀_\alpha `$, the best<sup>7</sup><sup>7</sup>7 Specifically, under the assumption that foregrounds and detector noise have multivariate Gaussian probability distributions, one one can show (Tegmark 1997) that this foreground removal method retains all the cosmological information present in the multifrequency set of input maps, i.e., constitutes an information-theoretically lossless form of data compression reducing all data down to a single CMB map. One does not need to assume that the CMB itself is Gaussian. Systematics such as foreground signals, data spikes, and atmospherics seldom have a Gaussian probability distribution – in this more general case, the removal method is no longer strictly lossless, but retains the feature that it minimizes the total rms of foregrounds and noise assuming only their finite second moment.  way to combine them is still given by equation (25) — we simply need to include more physics in the noise covariance matrix $`𝐍`$. Let us be more explicit about this. The noise covariance matrix $`𝐍`$ will be of size $`n\times n`$, where $`n=_{\alpha =1}^dn_\alpha `$, i.e., the total number of numbers in all data sets combined. We will therefore write the elements of $`𝐍`$ as $`𝐍_{\alpha i\beta j}`$, where $`\alpha `$ and $`\beta `$ determine the data sets and $`i`$ and $`j`$ the numbers therein. If there are $`f`$ foreground components, this noise matrix becomes a sum $`𝐍=_{k=0}^f𝐍^{(k)}`$, where $`𝐍^{(0)}`$ is the contribution from instrumental noise and the other terms are the contributions from foregrounds (synchrotron emission, bremsstrahlung, dust, point sources, etc.). Each of these foreground matrices will be of the form

$$𝐍_{\alpha i\beta j}^{(k)}=𝐅_{\alpha \beta }^{(k)}[𝐀_\alpha 𝐂^{(k)}𝐀_\beta ^t]_{ij},$$
(30)   
where the $`d\times d`$ matrix $`𝐅^{(k)}`$ gives the covariance of the $`i^{th}`$ foreground between frequencies, and the $`m\times m`$ matrix $`𝐂^{(k)}`$ gives its spatial covariance between pixels. For example, if the data sets are identical maps ($`𝐀_\alpha =𝐈`$) at $`d=2`$ different frequencies, we obtain the $`(2m)\times (2m)`$ block matrix

$$𝐍_{\alpha i\beta j}^{(k)}=𝐅^{(k)}𝐂^{(k)}=\left(\begin{array}{cc}F^{(k)}_{11}𝐂^{(k)}& F^{(k)}_{12}𝐂^{(k)}\\ F^{(k)}_{21}𝐂^{(k)}& F^{(k)}_{22}𝐂^{(k)}\end{array}\right).$$
(31)   
Explicit models specifying the foreground dependence on frequency and position (the matrices $`𝐅`$ and $`𝐂`$) can be found in Tegmark (1998). Although it has been common to characterize $`𝐂`$ by a power spectrum as in equation (8), the above formalism clearly works even if we break the isotropy by introducing an additional dependence on galactic latitude, say.

This shows how to best combine data sets at different frequencies. When combining different multifrequency experiments, it is desirable to first combine corresponding maps at the same frequency and apply the null-buster test for systematic errors. The different frequencies can then be merged afterwards, as a second step.

## 5. CONCLUSIONS

We have derived a set of useful tools for comparing and combining CMB data sets, all based on simple matrix operations, and drawn the following conclusions:

1. The “null-buster” test is better at detecting systematic errors than a simple $`\chi ^2`$-test.
2. Such a consistency test can be performed even between two experiments with quite different beam shape and observing strategy.
3. When combining two maps of different angular resolution, one need not smooth the higher resolution down to the lower one.
4. When combining two identical maps, one should generally not do the averaging separately for each pixel.
5. Our map merging method also handles the case of e.g. Planck, where the beams are elliptical rather than round and the different detectors have different beam orientations.
6. The foreground removal method of Tegmark (1998) is a special case of the combination technique we derived, and can be carried out even for foreground models involving frequency dependence, scale dependence, latitude dependence and spectral index variations in combination.

As the available CMB data sets continue to increase in quantity and quality, it will be useful to perform such cross-checks against systematic errors and incrementally combine all consistent data sets into a single state-of-the art map containing our entire knowledge of the CMB sky.<sup>8</sup><sup>8</sup>8 The signal and noise covariance matrices $`\mathrm{SS}`$ and $`𝐍`$ of such a map could of course be rather complicated. However, the same complication is encountered if one maintains separate experiment-specific maps and only attempts to combine the power spectrum measurements. This is because the sample variance of the band power measurements from different experiments are generally not independent, and computation of the correlation requires knowledge of these matrices. A disadvantage of combining only at the power spectrum level is that this generally sacrifices useful information on relative phases, whereas the merged map retains all the cosmological information from both data sets. For instance, near degeneracies in the noise covariance matrices of two maps can often be broken by combining them. 

Both theories and new observations can then be tested against this combined map as it gradually grows in size, quality and resolution.

The author wishes to thank Angelica de Oliveira-Costa and Martin White for helpful comments. Support for this work was provided by NASA though grant NAG5-6034 and Hubble Fellowship HF-01084.01-96A from STScI, operated by AURA, Inc. under NASA contract NAS5-26555.

This paper is available at h t t p://www.sns.ias.edu/$`\stackrel{~}{}`$max/comparing.html