Monopole and dipole estimation for multi-frequency sky maps by linear regression
I. K. Wehus, U. Fuskeland, H. K. Eriksen, A. J. Banday, C. Dickinson, T Ghosh, K. M. Gorski, C. R. Lawrence, J. P. Leahy, D. Maino, P. Reich, W. Reich
AAstronomy & Astrophysics manuscript no. PIP˙104˙Wehus˙arxiv c (cid:13)
ESO 2018October 16, 2018
Monopole and dipole estimation for multi-frequency sky maps bylinear regression
I. K. Wehus , U. Fuskeland , H. K. Eriksen , A. J. Banday , , C. Dickinson , T. Ghosh , K. M. G´orski , , C. R. Lawrence , J. P. Leahy ,D. Maino , , P. Reich , and W. Reich
1. Jet Propulsion Laboratory, California Institute of Technology, 4800 Oak Grove Drive, Pasadena, California, U.S.A., ∗
2. Institute of Theoretical Astrophysics, University of Oslo, Blindern, Oslo, Norway3. Universit´e de Toulouse, UPS-OMP, IRAP, F-31028 Toulouse cedex 4, France4. CNRS, IRAP, 9 Av. colonel Roche, BP 44346, F-31028 Toulouse cedex 4, France5. Jodrell Bank Centre for Astrophysics, Alan Turing Building, School of Physics and Astronomy,The University of Manchester, Oxford Road, Manchester, M13 9PL, U.K.6. Institut d’Astrophysique Spatiale, CNRS (UMR8617) Universit´e Paris-Sud 11, Bˆatiment 121, Orsay, France7. Warsaw University Observatory, Aleje Ujazdowskie 4, 00-478 Warszawa, Poland8. Dipartimento di Fisica, Universit`a degli Studi di Milano, Via Celoria, 16, Milano, Italy9. INAF / IASF Milano, Via E. Bassini 15, Milano, Italy10. Max-Planck-Institut f¨ur Radioastronomie, Auf dem H¨ugel 69, 53121 Bonn, GermanyPreprint online version: October 16, 2018
ABSTRACT
We describe a simple but e ffi cient method for deriving a consistent set of monopole and dipole corrections for multi-frequency skymap data sets, allowing robust parametric component separation with the same data set. The computational core of this method islinear regression between pairs of frequency maps, often called “T-T plots”. Individual contributions from monopole and dipole termsare determined by performing the regression locally in patches on the sky, while the degeneracy between di ff erent frequencies is liftedwhenever the dominant foreground component exhibits a significant spatial spectral index variation. Based on this method, we presenttwo di ff erent, but each internally consistent, sets of monopole and dipole coe ffi cients for the 9-year WMAP, Planck µ m,Haslam 408 MHz and Reich & Reich 1420 MHz maps. The two sets have been derived with di ff erent analysis assumptions and dataselection, and provides an estimate of residual systematic uncertainties. In general, our values are in good agreement with previouslypublished results. Among the most notable results are a relative dipole between the WMAP and Planck experiments of 10–15 µ K(depending on frequency), an estimate of the 408 MHz map monopole of 8 . ± . . ± .
03 K pointing towards Galactic coordinates ( l , b ) = (308 ◦ , − ◦ ) ± ◦ . These values represent the sum of any instrumentaland data processing o ff sets, as well as any Galactic or extra-Galactic component that is spectrally uniform over the full sky. Key words. methods: statistical – cosmology: observations – Galaxy: general – radio continuum: general
1. Introduction
The cosmic microwave background (CMB) fluctuations consistof small variations with a root-mean-square (RMS) of 70 µ K im-printed on top of a mean temperature of 2.73 K and a Doppler-induced dipole of ∼ ∼ − relativeto the total signal. To minimize systematic uncertainties, mod-ern CMB anisotropy experiments are therefore forced to employsome form of di ff erential measuring technique, eliminating thelarge 2.73 K o ff set already at the instrument level. Both COBE -DMR (Smoot et al. 1992) and WMAP (Bennett et al. 2013) em-ployed coupled di ff erencing assemblies that only recorded tem-perature di ff erences between two positions on the sky, while for Planck the instrumental o ff set is large and unknown, and can-not be used to constrain the monopole (Planck Coll. I 2013).However, while necessary for systematics suppression, this alsoimplies that these experiments are intrinsically unable to mea-sure the true absolute zero-point (or monopole) of their finalmaps. In addition, the dipole is also associated with a large rela- ∗ e-mail: [email protected] tive uncertainty because of the large numerical value of the CMBDoppler dipole; a small relative error in the determination of theCMB dipole direction can induce a dipole error of many mi-crokelvins in a CMB map. Typically, this will be strongly corre-lated among frequencies within a single experiment, though, andso informative priors can be imposed among frequency chan-nels within a given experiment. For foreground-dominated fre-quency channels, instrumental systematics, such as gain fluctua-tions, may induce dipole errors.Removing the monopole and dipole from CMB data setsdoes not constitute a major limitation in terms of CMB-basedcosmological parameter analysis, since losing a handful har-monic modes out of many thousands only negligibly reducesthe total amount of available information. However, it doeshave a significant indirect impact because of the presence ofnon-cosmological foreground contamination from Galactic andextra-Galactic sources. In order to obtain a clean image of thecosmological CMB fluctuations, this foreground contaminationhas to be removed from the raw sky maps through some form ofcomponent separation prior to power spectrum and parameter es-timation. A wide range of such methods have been already been a r X i v : . [ a s t r o - ph . C O ] N ov . K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression proposed (e.g., Planck Coll. XII 2013, and references therein),ranging from the simplest of template fitting and internal linearcombination approaches through blind- or semi-blind image pro-cessing techniques to full-blown parametric Bayesian methodsemploying physical models. Except for the very simplest meth-ods, all of these exploit the fact that all foreground frequencyspectra are qualitatively di ff erent from the CMB spectrum. Forinstance, while the CMB spectrum is that of a perfect black-body, thermal dust emission can be well approximated by thatof a one- or two-component greybody. However, such relationsclearly only hold if there are no arbitrary o ff sets between the dif-ferent frequency maps. In other words, spurious monopole anddipole errors can bias any estimation algorithm that exploits fre-quency dependencies, and this can in turn lead to leakage be-tween various components, and eventually contamination in theCMB estimate.A number of methods have already been proposed in the lit-erature for estimating monopoles, while fewer have addressedresidual dipoles. One example of the former is the co-secantmethod adopted by the WMAP team (Bennett et al. 2003). Inthis case, the Galactic signal is approximated as plane-parallel inGalactic coordinates, with an amplitude falling roughly propor-tionally with the co-secant of the latitude. The major weaknessof this method is that the Galaxy is neither plane-parallel nor fol-lows a co-secant, and the method does also not account for resid-ual dipole terms. A second approach was proposed by Eriksenet al. (2008), who include the monopole and dipole terms asadditional free parameters within a global Bayesian parametricframework. The major weakness with this method is a large de-generacy between the monopole and dipole terms relative to theunknown zero-point of each foreground; it is possible to add aconstant to each foreground amplitude, and then subtract a corre-sponding frequency-scaled o ff set from each monopole, leavingthe net sum unchanged.A third widely used method for setting the zero-level of ra-dio maps is that of linear regression, or through so-called “T-Tplots”. This method has a long and prominent history in radio as-tronomy (see, e.g., Turtle et al. 1962; Davies et al. 1996; Reich &Reich 1988; Reich et al. 2004; Wehus et al. 2013, and referencestherein), as it provides a highly robust estimate of the spectralindex of a single signal component given observations at twodi ff erent frequencies. When plotting the measured pixel valuesat one frequency as a function of the measured pixel values at theother frequency, the spectral index is given (up to a constant fac-tor) by the slope of the resulting T-T plot, which is easily foundby linear regression. The main virtue of this spectral index es-timate is that it is completely insensitive to any constant o ff setin either of the two frequency maps, since these only a ff ect theregression intercept, not the slope. Intuitively, o ff sets only shiftthe scatter plot horizontally or vertically, but they do not deformor rotate it.In the following, we exploit the same idea to estimate bothspurious monopoles and dipoles by noting that when the truevalues are correctly determined, the regression intercept has tobe zero: If the foreground signal at one frequency is exactly zero,it also has to be exactly zero at the other frequency . This impliesthat a single T–T plot constrains both the spectral index of thecomponent and the relative o ff sets of the two maps, m and m ,to m = am + b , where a and b are the slope and intercept of the Line emission processes, such as that arising from carbon-monoxide (CO), are clearly exceptions to this, and regions with sig-nificant line emission must be masked before applying the method pre-sented in this paper.
30 GHz amplitude ( µ K) GH z a m p lit ud e ( µ K ) β = -2 β = -3 β = -4
30 GHz offset GH z o ff s e t Fig. 1: Schematic illustration of a single T–T plot evaluatedfor an ideal low-frequency
Planck simulation (30 and 44 GHz).Each dot represents the observed values at two frequencies for asingle pixel. The slope of the distribution is given by the spectralindex of the signal component in the field, as indicated by thedashed lines. A constant o ff set in either frequency map simplytranslates the entire T–T plot either horizontally or vertically. Formaps without spurious o ff sets, the best-fit straight line shouldpass through the origin; the fundamental idea of the algorithmpresented in this paper is to ensure that this is the case for anysu ffi ciently small patch of the sky.T–T plot. The main goal of the present paper is to develop thissimple idea into a complete and robust method for determiningboth monopoles and dipoles from a set of multi-frequency skymaps.
2. Method
Before describing the main method, we note that robust linear re-gression is in general di ffi cult in the low signal-to-noise regime,as both the slope and intercept are associated with large uncer-tainties. We therefore adopt a two-step process in which we firstuse the main T–T plot method for high signal-to-noise frequencychannels, followed by a direct template fit method for low signal-to-noise channels. The second stage, however, is both concep-tually and implementationally straightforward and well estab-lished in the literature; the di ffi cult task is to set the o ff sets cor-rectly for the high signal-to-noise components, and this is ourprimary concern in the following. We start by considering a basic data model consisting of a singlesignal component on the form d ν ( p ) = F ( ν ) s ( p ) + n ν ( p ) + m ν , (1)where d ν ( p ) denotes the observed sky map value at frequency ν and pixel p , s represents the true sky signal (at some reference
2. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression frequency, ν ), F ( ν ) denotes the (normalized) frequency spec-trum of the sky signal (often called mixing matrix in the com-ponent separation literature), n ν ( p ) is instrumental noise, and m ν is the unknown spurious o ff set we want to estimate and remove.For now, m ν is taken to be a pure constant (monopole), and F ( ν )is assumed constant over the entire observed field.For many radio astronomy applications, the signal spectrumcan be approximated by a power-law, F ( ν ) = ( ν/ν ) β , at leastover some limited frequency range. The e ff ective spectral index, β , can be computed between any two frequencies, ν and ν , us-ing the so-called “T–T plot” technique. Let us first consider thenoiseless case without o ff sets. For this case, the spectral index isgiven by d = (cid:32) ν ν (cid:33) β d ⇔ β = log (cid:16) d d (cid:17) log (cid:16) ν ν (cid:17) . (2)Defining the slope as a ≡ ( ν /ν ) β , and including noise and spu-rious o ff sets, this expression reads d − ( n + m ) = a [ d − ( n + m )] (3) d = ad + ( m − a m ) + ( n − a n ) (4) ≡ ad + b + n (5)Thus, the observed signal at one frequency is related linearly tosignal at the other frequency, with a slope given uniquely by thespectral index, and an intercept given by b ≡ ( m − a m ).This relation is often conveniently visualized in terms ofscatter (or “T–T”) plots, as illustrated in Fig. 1. Each dot indi-cates the observed data values at the two frequencies for onepixel, while the dashed lines indicate three models with dif-ferent values of β . A spurious o ff set in either frequency mapcorresponds directly to a vertical or horizontal translation ofthe entire scatter plot, respectively, but does not change theslope. Therefore, the spectral index is fully insensitive to spuri-ous constant o ff sets when estimated by this T–T plot technique.However, as shown below spurious dipoles do bias the spectralindex, because they introduce a gradient across the field, result-ing in a net additional tilt. The appropriate correction for this isdiscussed below.As long as F ( ν ) (cid:44) ν , it is clear that if the data vanishin one frequency, it also have to disappear in the other. Thus, fordata free of any spurious o ff sets the best-fit line through the T–Tplot must pass through the origin. Correspondingly, a non-zeroregression intercept can only be due to the o ff set term in Eq. 5,implying that m and m must be related by m − a m = b . (6)In other words, the true pair of o ff sets have to lie somewherealong the best-fit regression line in the T–T plot. From a singlescatter plot it is impossible to determine the precise location, butthis simple linear relation nevertheless forms the core unit of ouralgorithm.Before proceeding with the algorithm, we make a few com-ments concerning the implementation of the fitting procedurefor a and b . First, since both d and d are noisy quantities, thestandard method based on the normal equations does not strictlyapply, as that estimator is known to su ff er from so-called atten-uation bias ; noise in the descriptor variable, d , biases a low(e.g., Draper & Smith 1998). Second, in real data sets strongoutliers occur quite frequently; a typical example in the CMBsetting is unmasked point sources. Since we will only apply thismethod to high signal-to-noise data sets in the following, thesecond problem is definitely more pressing for our purposes. We accordingly adopt the non-parametric and highly robust Theil-Sen estimator in the following (Theil 1950): We estimate theslope as ˆ a = median[( d ( p ) − d ( p )) / ( d ( p ) − d ( p ))] eval-uated over all pixel pairs ( p , p ), and the intercept as ˆ b = median[ d ( p ) − ˆ ad ( p )]. If one wants to apply the same method tolow signal-to-noise observations, a more appropriate estimator isthe Deming estimator (Deming 1948), which properly accountsfor uncertainties in both directions. We have implemented bothof these estimators in our codes, and find fully consistent resultsfor the cases considered in this paper. Above we assumed that all pixels in the region of interest havethe same spectral index, and the entire map may therefore beanalyzed within a single T–T plot. In reality, the true spectral in-dex varies across the sky to some extent. For instance, the spec-tral index of synchrotron emission typically ranges between, say, β s = − . − .
1, while the thermal dust emissivity ranges be-tween, say, β d = . (cid:48) pixels (corresponding to a HEALPixresolution parameter of N side = ◦ N side = ×
16 high-resolution pixels.Such partitioning into smaller regions is not only useful inorder to ensure nearly constant spectral indices, but it allows infact an absolute determination of the individual o ff sets, m and m . Let the full data set be partitioned into N regions, each with a(nearly) constant spectral index β i , and perform an independentlinear regression procedure for each region, as outlined above.In this case, one will obtain an independent linear constraint on m and m , given by Eq. 6, from each region . These may becombined into the following (over-determined) linear system, − a − a ... − a N (cid:34) m m (cid:35) = b b ... b N (7)Writing this linear system in a matrix form, Am = b , ithas a unique solution given by the normal equations, ˆ m = ( A t A ) − A t b . Thus, any spatial variation in the spectral indexbreaks the degeneracy between the o ff sets at the two frequen-cies, and at least formally allows absolute determination of both. Partitioning the sky into sub-regions further allows us to esti-mate additional degrees of freedom. Specifically, suppose thatthe total o ff set parameter space may be spanned by some setof basis vectors, T k ( p ), each with an unknown amplitude, z k .The archetypal example is the space spanned by a monopole and
3. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression three dipoles, such that the total e ff ective o ff set for a given fre-quency map is˜ m = Tz (8) = φ sin θ sin θ sin φ cos θ φ sin θ sin θ sin φ cos θ ... ... ... ... φ n sin θ n sin θ n sin φ n cos θ n z z z z , (9)where subscripts indicate pixel number. To account for these newdegrees of freedom within a single region i , Eq. 6 generalizes to (cid:88) k = T k ( i ) m k , − a (cid:88) k = T k ( i ) m k , = b , (10)where we approximate the net impact of the additional templatesas constant over the region, i.e., T k ( i ) = / N pix (cid:80) p ∈ i T k ( p ). Thefull joint all-regions linear system, corresponding to Eq. 7, is cor-respondingly generalized into Ax = b , where A now is an N × − a i T k ( i ) for region i and template component k in the first four columns, and T k ( i ) in the four last columns. The x vector has 8 elements containing the template amplitudes forthe first map in the first four entries, and the template amplitudesfor the second map in the four last entries. The right-hand side, b , is identical to that in Eq. 7. Again, this system is solved by thenormal equations, ˆ x = ( A t A ) − A t b .While spurious monopoles do not change the net slope of aT–T plot, spurious dipoles do. And the larger area the consideredregion covers, the larger the e ff ect is. We account for this e ff ectby iteration. That is, we 1) solve Eq. 10 as described above, 2)subtract the derived monopole and dipole estimates from the rawinput maps, and 3) iterate until all o ff set updates are smaller than,say, 1% of the total value. Typically three or four iterations areneeded to reach convergence.In Fig. 2 four di ff erent simulations illustrate the various casesdiscussed so far. Each line corresponds to the best-fit linear fitto the T–T plot of a synchrotron-only sky evaluated at 30 and44 GHz, corresponding to the two lowest Planck frequencies.Each sky map has an angular resolution of 1 ◦ FWHM, and ispixelized on an N side =
256 HEALPix grid. The spectral indexis chosen to be constant within each N side = β = − ± .
2, resulting in 768independent T–T plots.In the top left panel, we show the ideal case with neitherspurious o ff sets nor instrumental noise. In this case, we see thatall lines truly converge on the origin, as they should. In thetop right panel we have added a − µ K o ff set to the 30 GHzchannel, and a + µ K o ff set to the 44 GHz channel. The en-tire plot is translated accordingly, now focusing on the point( d , d ) = ( − , µ K. Intuitively, the main goal of monopolecorrection is to re-center the focal point on the origin.Next, the bottom left panel shows the e ff ect of adding a spu-rious dipole to each frequency band. For a single region, this isalmost equivalent to a simple translation, just like a monopole;however, each scatter plot is translated di ff erently, depending onits position on the sky. When considering all T–T plots simulta-neously, the overall distribution therefore appears smeared out,and possibly o ff set from the correct position, depending on therelative orientation between the two dipoles and the dominantGalactic signal. The intuitive goal of dipole correction is to makethe focus point of this plot as sharp as possible.Finally, the bottom right panel illustrates the e ff ect of instru-mental noise. This simply smears out each individual line, mak-ing it harder to assess where the lines converge to a single point. -40-2002040-100 -50 0 50 100-40-2002040 -50 0 50 100 Ideal MonopoleDipole Noise
30 GHz amplitude ( µ K) GH z a m p lit ud e ( µ K ) Fig. 2: T–T summary plot of four di ff erent simulations. Eachcoloured line corresponds to the best-fit to a single HEALPix N side = N side = ff sets nor instrumental noise; alllines converge perfectly on the origin. In the top right panel spu-rious o ff sets of −
30 and + µ K are added to the two frequen-cies, resulting in a simple corresponding translation of the entireplot. The bottom left panel shows the e ff ect of a spurious dipole;the lines no longer converge on a single point, as each regionis e ff ectively translated by a di ff erent o ff set. Finally, the bottomright panel illustrates the e ff ect of instrumental noise.When the instrumental noise becomes comparable to the signal,robustly estimating the slopes and intercepts becomes very di ffi -cult, and we choose for now not to be aggressive in this respect;for low signal-to-noise cases, we find that simpler template fit-ting methods yield more robust results. The algorithm presented in the previous section constitutes thecentral engine in our method, and is already at this stage a self-contained and complete method for ideal data sets. However, realdata are seldom ideal, and several adjustments and extensionsare usually required before the method becomes practical. In thissection we present a list of these issues, as well as their solutions.
First, most multi-frequency data sets typically have di ff erent an-gular resolutions at di ff erent frequencies. In order to estimatethe spectral indices (i.e., slopes) reliably across frequencies, itis therefore necessary to smooth all bands to a common angularresolution. To do so, we decompose each sky map into spheri-cal harmonics, d ( p ) = (cid:80) (cid:96) m a (cid:96) m Y (cid:96) m ( p ). According to the spher-ical convolution theorem, a convolution in pixel domain trans-lates into a multiplication in harmonic domain by the convolu-tion theorem, such that if b (cid:96) is the Legendre transform of the
4. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression intrinsic instrumental beam, and b (cid:96) is the Legendre transformof the desired common beam, the smoothed map is given byˆ d ( p ) = (cid:80) (cid:96) m ( b (cid:96) / b (cid:96) ) a (cid:96) m Y (cid:96) m ( p ).This degradation does not change the monopole or dipolesof the original map, and the derived low-resolution o ff set cor-rections can therefore be applied also to the full-resolution dataset. However, it is important to note that information is lost inthis process, in the sense that the T–T plots exhibits a smallerdynamic range after smoothing, e ff ectively making it harder topinpoint the optimal solution. One should therefore not smoothmore than necessary to bring the frequency maps to a commonresolution. Second, many recent data sets have more than two frequencychannels, whereas the T–T plot method intrinsically only in-volves pairs of maps. To deal with multiple maps, we order themaps according to frequency, and derive Eq. 7 for each pair ofneighbouring frequencies. Considering the simplest case withonly a monopole degree-of-freedom for each of k frequencybands, this results in the following joint system, − a . . . − a . . .
00 0 − a . . . ... ... ... ... ... . . . m m ... m k = b b ... b k (11)Generalization to dipole estimation and multiple regions isstraightforward, although somewhat notationally involved. Notealso, of course, that nothing prevents including frequency pairsbeyond neighboring in the fit; the only limitation is that anygiven frequency pair should be well described by a single domi-nant signal component, which typically sets an e ff ective limit onthe allowed frequency range. For frequencies between 30 and 143 GHz, the
Planck andWMAP sky maps are dominated by CMB fluctuations ratherthan di ff use Galactic foregrounds. Over the cleanest regions ofthe sky, these fluctuations can therefore in principle serve asthe signal for evaluating the scatter plot slopes and intercepts.However, this is non-trivial for at least two reasons. First, theCMB variance is typically of the order of 70 µ K on degree angu-lar scales, while both the desired o ff set precision and the instru-mental noise are typically just a few microkelvins. Second, sincethe CMB frequency spectrum follows a blackbody, it has by def-inition a constant spectral index (equal to 0) at all frequenciesand positions on the sky. O ff set estimation on CMB fluctuationstherefore give at most relative results, and even those are associ-ated with relatively large uncertainties.The most straightforward solution, and the one adopted inthe following, is to subtract an estimate of the CMB sky (e.g.,Planck Coll. XII 2013; Bennett et al. 2013) from all sky mapsprior to o ff set estimation. The accuracy of this estimate is notcritical, as the nature of CMB fluctuations are very close toGaussian, and they therefore mostly add random noise to the T–Tplots. As long as the CMB uncertainties are significantly smallerthan the absolute foreground amplitude of the relevant channel,which usually is the case with current CMB experiments, the fitis stable. However, in order to ensure that the resulting o ff set correction estimates are directly applicable to the original skymaps, one should ensure that whatever CMB estimate is sub-tracted is orthogonal to all basis vectors, T k , which in practiceimplies making sure that it does not have any monopole or dipolecomponents. However, an uncertainty in this estimate will in theend translate directly into an uncertainty with the very specificCMB frequency spectrum, and will therefore not confuse com-ponent separation algorithms, as long as these also allow for amonopole and dipole component in the CMB fit. Related to the previous issue, the algorithm intrinsically assumesthe presence of only one dominant foreground component per re-gion. For data sets such as the low-frequency and synchrotrondominated 408 MHz Haslam and 1420 MHz Reich & Reichmaps, this is for most parts of the sky not an issue. Neither is itfor the high-frequency
Planck channels above 143 GHz, whichare dominated by thermal dust emission. However, for the low-foreground
Planck and WMAP CMB channels between, say, 30and 100 GHz somewhat greater care is warranted. For these fre-quencies, the overall signal budget is made up by a combina-tion of synchrotron, free-free, anomalous microwave emission(AME; spinning dust), CO and thermal dust emission (PlanckColl. XII 2013). For some frequencies, this complication can besolved by masking, by removing spatially localized componentssuch as CO and free-free. For other frequencies, say, between 20and 40 GHz, where both synchrotron and AME are significantand not spatially localized, more sophisticated component sep-aration methods should be used, simultaneously accounting formultiple components.
Near the foreground minimum around 70 GHz, low-frequencyforegrounds (AME, free-free and synchrotron) and thermal dustcontribute (by definition) equally, and the scatter plot techniqueis therefore intrinsically unreliable. In addition, the foregroundsignal-to-noise level is low, further destabilizing the T–T plottechnique. On the other hand, the absolute foreground levels arealso correspondingly low, and even simple methods are typicallyable to derive quite accurate monopole and dipole estimates. Wetherefore adopt the following template fitting technique for fre-quencies between 33 (WMAP Ka-band) and 100 GHz: we firstderive absolute o ff set corrections for all high-foreground fre-quencies, and apply these to the respective (CMB-subtracted)channel maps. We then adopt the nearest neighbors on either sideof the foreground minimum (i.e., the Planck
30 and 143 GHzmaps) as low-frequency foreground and thermal dust templates,respectively. We then fit these together with the monopole andthree dipoles to the low-foreground channels by solving the gen-eralized normal equations, z = ( T T N − T ) − T T N − d , (12)where T is the N pix × N is the (assumed diagonal) pixel noise covariance ma-trix. A very conservative mask excluding both point sources andresidual Galactic contamination is applied in the fit, as shown inFig. 4.We typically find that the uncertainties from this fit are on theorder of a few microkelvins. However, even such small uncer-tainties can in principle be detected by detailed χ based analy-ses, for instance as implemented in parametric foreground meth-ods like Commander (Eriksen et al. 2008). For such applications,
5. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression -20020253035-65-60-5590950510159101112202224 1 2 3 4
Iteration number -17-16-15-14 O ff s e t a m p lit ud e ( µ K ) Monopole @ 30 GHzX dipole @ 30 GHzY dipole @ 30 GHzZ dipole @ 30 GHzMonopole @ 44 GHzX dipole @ 44 GHzY dipole @ 44 GHzZ dipole @ 44 GHz
Fig. 3: Recovered o ff set coe ffi cients as a function of analysis iter-ation for a two-band simulation. The true input values are shownas dashed horizontal lines. The values change between iterationsbecause spurious dipoles bias the slope of the T–T plots; itera-tive dipole correction remove this bias. Overall, the method re-produces the true input values to (cid:46) µ K. However, the adoptedbootstrap uncertainties tend to underestimate the uncertaintiesin the sub-dominant channels (44 GHz in this case), and are onlyintended to give an indication of the true uncertainties. End-to-end simulations are required for fully reliable uncertainty esti-mation.we recommend that one fits (or at least verifies) the o ff set am-plitudes for these channels within the joint analyses itself, ac-counting simultaneously for foreground amplitudes and overallo ff sets. Fig. 4: Joint mask used for the analysis of the WMAP, Planck and 100 µ m maps. To improve the rigidity and physicality of the fit, it can be ad-vantageous to impose a positivity prior on the foreground am-plitudes, by requiring that the post-correction frequency map isnon-negative everywhere. We implement this as an optional fea-ture in our codes as follows. We locate the coldest set of pixelson the sky, separated by at least 10 ◦ on the sky. Each pixel valuedefines an independent inequality constraint on the o ff set coe ffi -cients given by d ( p ) − (cid:88) i = T i ( p ) z i ≥ . (13)For cases with significant noise contribution, the inequalityshould be relaxed by adding N σ ( p ) to the right-hand side, where σ ( p ) is the instrumental noise RMS in pixel p , and N is a thresh-old level in units of σ ; we adopt 4 σ as our positivity threshold.Note that optimization of, say, Eq. 11 under these constraints isa substantially computationally more complicated problem thansolving the linear normal equations, and must be performed us-ing non-linear methods. Finally, we make a short note on estimation of statistical uncer-tainties, but emphasize that this topic is by far the most compli-cated part of the entire procedure, and the method outlined hereis only intended to give a rough estimate of the uncertainties. Thefundamental problem is that for most current experiments, themonopole and dipole coe ffi cient uncertainties are vastly dom-inated by systematic e ff ects (foreground modelling, optical im-perfections etc.), rather than instrumental noise. Noise-based un-certainties are therefore virtually meaningless for describing trueuncertainties. For this reason, we adopt a Monte Carlo basedbootstrap method for now, aiming to capture some of theseintrinsic systematic uncertainties in a non-parametric manner.From a data set consisting of m disjoint sky regions, we selectrandomly a sub-sample of m regions (i.e., one region may be in-cluded several times), and perform the full analysis on this sub-sample in the same manner as for the original data set. This pro-cess is repeated typically 100 times, and the resulting varianceamong those 100 resamples is taken as the bootstrap uncertainty.Further, any tunable parameter, such as whether to perform theanalysis on N side =
4, 8 or 16 regions, are also drawn randomlywithin their allowed ranges between each resample.The main advantage of this bootstrap approach is that it doesto some extent account for foreground modelling uncertainties
6. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression −
500 500WMAP K - Planck 70 [ µ K] −
500 500WMAP K - Planck 70 [ µ K] −
50 50Planck 100 - Planck 70 [ µ K] −
50 50Planck 100 - Planck 70 [ µ K] −
200 200Planck 30 - Planck 70 [ µ K] −
200 200Planck 30 - Planck 70 [ µ K] −
100 100Planck 143 - Planck 70 [ µ K] −
100 100Planck 143 - Planck 70 [ µ K] −
100 100WMAP Ka - Planck 70 [ µ K] −
100 100WMAP Ka - Planck 70 [ µ K] −
500 500Planck 217 - Planck 70 [ µ K] −
500 500Planck 217 - Planck 70 [ µ K] −
50 50WMAP Q - Planck 70 [ µ K] −
50 50WMAP Q - Planck 70 [ µ K] − µ K] − µ K] −
50 50Planck 44 - Planck 70 [ µ K] −
50 50Planck 44 - Planck 70 [ µ K] − µ K] − µ K] −
50 50WMAP V - Planck 70 [ µ K] −
50 50WMAP V - Planck 70 [ µ K] − µ K] − µ K] −
50 50WMAP W - Planck 70 [ µ K] −
50 50WMAP W - Planck 70 [ µ K] − µ m - Planck 70 [ µ K] × − µ m - Planck 70 [ µ K] × −
250 250Planck 70 [ µ K] −
250 250Planck 70 [ µ K] Fig. 5: Di ff erence maps between the CMB and high-frequency channels and the 70 GHz map, before (first and third columns) andafter (second and fourth columns) o ff set corrections. The bottom row shows the pre- and post-correction 70 GHz map itself.and algorithmic choices. However, it is still determined by theactually realized sky, and can therefore only account for thosevariables that vary within our actual data set; not those that arefixed for a given realization, but in principle stochastic. Themethod will therefore necessarily underestimate the true uncer- tainties, and we caution against interpreting these as Gaussian68% confidence levels. For fully reliable uncertainty estimation,proper end-to-end simulations (including di ff erent foregroundrealizations in each simulation) are likely to be the only truly
7. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression satisfactory solution, which is outside the scope of the currentpaper.
Before applying our method to real observations, we analyzea simple well-controlled simulation, to check that the methodproduces sensible results in this case. Specifically, we onceagain consider a synchrotron-plus-noise simulation, with the sig-nal component evaluated at the two lowest
Planck frequencies,30 and 44 GHz, but this time adopting the spatial spectral in-dex distribution derived from the 408 and 1420 MHz maps inSect. 4, and summarized in Fig. 9. In addition, we add spuriousmonopoles and dipoles to both maps ranging between −
60 and + µ K.The results from this simulation are summarized in Fig. 3for each o ff set parameter, as a function of analysis iteration. Thedashed lines show the true value, and the uncertainties indicatethe bootstrap errors described above. The general behaviour seenin this plot is typical for all cases we have analyzed, and there-fore serves as a useful tool for building intuition about the per-formance of the method. First, and most importantly, we see thatthe method overall reproduces the true input values in terms ofabsolute values to a precision of at worst 1–2 µ K.Second, we see that the largest changes are observed be-tween the first and the second iteration. This is due to the al-ready mentioned fact that a spurious dipole introduces a bias inthe e ff ective slope (or spectral index) of the T–T plots, and thisin turn leads to a leakage of dipole power into the monopole.However, even the first-order dipole correction leads to a vastlyimproved monopole estimate, and subsequently nearly stable re-sults; the 30 GHz monopole jumps directly from + µ K to thetrue value of − µ K in the second iteration. Additional itera-tions only change the results with small amounts.Third, as already stressed in Sect. 2.4, we see that the boot-strap uncertainties do not adequately describe the true uncertain-ties in the fit for all coe ffi cients. While they do a reasonable jobfor the dominant 30 GHz channel, they underestimate the uncer-tainties at 44 GHz by up to a factor of four or five. However, itis again important to note that the absolute uncertainties for thesame coe ffi cients are small. In general, the bootstrap errors tendto produce a reasonable estimate for the dominant channels, butunderestimate the uncertainties in the sub-dominant channels. Inthis paper, we will never quote uncertainties smaller than 1 µ Kfor any component, even if the formal bootstrap uncertainty fora few cases is as low as 0 . µ K.
3. Data and analysis setup
We now turn our attention to a set of 17 publicly available full-sky maps of the radio, millimeter and sub-millimeter sky, withthe goal of establishing a consistent set of o ff set coe ffi cientsthat can be used for multi-experiment CMB component sepa-ration analysis. We include in the following 1) the nine Planck µ m map by Schlegel et al. (1998). For Planck and WMAP we use the non-foreground cleaned co-addedfrequency maps; for the 408 MHz map, we use the version avail- able on LAMBDA that has no filtering applied; and the 100 µ mmap accounts for the DIRBE and IRAS bandpasses. All mapsare downgraded to a common resolution of 1 ◦ , and repixelizedat a HEALPix resolution of N side = Planck webase this on the provided high-resolution noise variance maps.For the two low-frequency and the 100 µ m channels we enforcea strict positivity prior, and simply demand that there should beno negative pixels at all.The 408 and 1420 MHz maps are analyzed together, and sep-arately from the other 15 maps, which are all analyzed jointly.For the low-frequency analysis, we adopt the WMAP KQ85mask (Bennett et al. 2013), and for the main analysis we addi-tionally apply the Planck
GAL60 Galactic and PCCS + SZ pointsource mask (Planck Coll. XII 2013). This joint mask is firstsmoothed by a 10 (cid:48)
Gaussian beam, and thresholded at a value of0.5, to remove the very smallest source holes from the
Planck mask. Then it is smoothed by a 1 ◦ Gaussian beam, and thresh-olded at 0.99, to exclude residual foregrounds leaking out of themask after smoothing each frequency map to the common reso-lution of 1 ◦ . The resulting mask is shown in Fig. 4, and permitsa total of 38 % of the sky.We use the main T–T scatter plot technique for the 408–1420 MHz combination, as well as for the combination ofWMAP K-band and Planck
30 GHz and for all frequenciesabove 143 GHz. For the frequencies between WMAP Ka-bandand
Planck
100 GHz we use the template fit technique describedin Sect. 2.4, and adopt the
Planck
30 and 143 GHz channels asforeground tracers for low-frequency foregrounds and thermaldust, respectively. For CMB suppression, we adopt the
Planck
Commander solution from which the monopole and dipole hasbeen removed outside the
Commander analysis mask (PlanckColl. XII 2013). For robustness, we also performed an anal-ysis using the low-foreground 9-year WMAP ILC CMB map(Bennett et al. 2013), obtaining nearly identical results.To explore overall stability with respect to analysis choices,we additionally analyze a subset of the above. As reported byPlanck Coll. XII (2013), the spectral index of thermal dust be-low 353 GHz is found to be lower than the expected value of1.3–1.7 over extended regions of the sky. This may be explainedeither in terms of systematic uncertainties in the maps, or a breakin the spectral index around 353 GHz, or simply a general break-down of the simple one-component greybody model. When in-cluding the high-frequency channels, as done above, this featureis regularized by high signal-to-noise measurements, resultingin a far more stable fit. However, this stability comes at a cost;if there happens to exist a second thermal dust component, thescatter plot technique is not well defined. In this case, it wouldbe better for the low-frequency channels to exclude the high-frequency channels, and thereby reduce the overall sensitivityto this second component. The second data set considered heretherefore comprises the 14 frequency bands up to (and includ-ing) the
Planck
353 GHz channel. However, as noted by PlanckColl. XII (2013), this does result in a large uncertainty for the353 GHz monopole. Therefore, we additionally fix this num-ber to 308 µ K, a value determined by Planck Coll. VII (2013)through cross-calibration with HI observations.We make three additional changes for this particular anal-ysis. First, we fix the WMAP dipoles to zero, recognizing thatthe WMAP scanning strategy should be well suited to measure http: // lambda.gsfc.nasa.gov8. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression Table 1: Monopole and dipole estimates for two di ff erent analysis configurations. The top section summarizes the results for themain analysis of all 17 frequency maps considered in this paper, and the bottom section shows the results from the reduced 14frequency data set employing external priors and di ff erent masks. Uncertainties are defined to be the maximum of the Monte Carlo-based bootstrap error described in Section 2.4 and 1 µ K. Conversion between Galactic Cartesian and polar coordinates is givenby ( l , b ) = (90 ◦ − acos[ Z / A ] , atan2[ Y / A , X / A ]), where A = √ X + Y + Z is the dipole amplitude. Temperatures are given inthermodynamic units. Sky map Method Partner(s) Monopole X dipole Y dipole Z dipole Unit17 band combination; tuned mask; no external priorsHaslam 408 MHz . . . . T–T fit 1420 MHz 8 . ± . . ± . . ± . − . ± . . ± .
02 0 . ± . − . ± . − . ± .
03 KWMAP K-band . . . . . T–T fit 30 GHz − ± ± ± ± µ KWMAP Ka-band . . . . Template fit 30 +
143 GHz 3 ± ± ± ± µ KWMAP Q-band . . . . . Template fit 30 +
143 GHz 1 ± ± ± ± µ KWMAP V-band . . . . . Template fit 30 +
143 GHz 1 ± ± − ± ± µ KWMAP W-band . . . . Template fit 30 +
143 GHz 2 ± ± − ± ± µ K Planck
30 GHz . . . . . T–T fit K-band 10 ± − ± ± ± µ K Planck
44 GHz . . . . . Template fit 30 +
143 GHz 3 ± ± ± − ± µ K Planck
70 GHz . . . . . Template fit 30 +
143 GHz 14 ± ± ± − ± µ K Planck
100 GHz . . . . Template fit 30 +
143 GHz 15 ± ± ± − ± µ K Planck
143 GHz . . . . T–T fit 217 GHz 34 ± ± ± − ± µ K Planck
217 GHz . . . . T–T fit 143 +
353 GHz 84 ± ± ± − ± µ K Planck
353 GHz . . . . T–T fit 217 +
545 GHz 315 ± ±
16 22 ± − ± µ K Planck
545 GHz . . . . T–T fit 353 +
857 GHz 0 . ± .
01 0 . ± . − . ± .
02 0 . ± .
01 MJy / sr Planck
857 GHz . . . . T–T fit 545 GHz + µ m 0 . ± .
03 0 . ± . − . ± .
04 0 . ± .
03 MJy / srSFD 100 µ m . . . . . . T–T fit 857GHz 0 . ± .
02 0 . ± . − . ± .
01 0 . ± .
01 MJy / sr14 band combination; WMAP KQ85 mask + | b | > ◦ ; fixed 353 GHz monopole and WMAP dipolesHaslam 408 MHz . . . . T–T fit 1420 MHz 7 . . . − . .
28 0 . − . − .
12 KWMAP K-band . . . . . T–T fit 30 GHz 27 0 a a a µ KWMAP Ka-band . . . . T–T fit 30 GHz + Q-band 16 0 a a a µ KWMAP Q-band . . . . . T–T fit Ka-band +
44 GHz 10 0 a a a µ KWMAP V-band . . . . . T–T fit 30 +
143 GHz 6 0 a a a µ KWMAP W-band . . . . T–T fit 30 +
143 GHz 7 0 a a a µ K Planck
30 GHz . . . . . T–T fit K-band 29 − − µ K Planck
44 GHz . . . . . T–T fit Q-band + V-band 12 1 7 − µ K Planck
70 GHz . . . . . T–T fit V-band + W-band 20 1 7 − µ K Planck
100 GHz . . . . T–T fit W-band +
143 GHz 17 3 8 − µ K Planck
143 GHz . . . . T–T fit 100 GHz +
217 GHz 33 2 13 − µ K Planck
217 GHz . . . . T–T fit 143 +
353 GHz 73 2 13 − µ K Planck
353 GHz . . . . T–T fit 217 GHz 308 b − µ K a Dipoles fixed to o ffi cial WMAP values. b
353 GHz monopole fixed by HI cross-correlation (Planck Coll. VII 2013). this particular mode. Second, to probe sensitivity to sky cover-age, we impose a less conservative sky cut consisting only of theunion of the WMAP KQ85 mask and a straight | b | > ◦ mask,in total allowing 54% of the sky. In this case, we use the T–Tplot technique for all frequencies.Analogous to fixing the WMAP dipole to zero in the con-sistency run, it is of course simple to impose additional externalconstraints whenever available, and these will always improvethe rigidity of the overall fit. For example, the dominant sourceof dipole uncertainty for the CMB-dominated channels is theCMB dipole itself, as indeed is demonstrated in the following. Inthese cases, one may therefore impose a sharp prior on the dipoledirection by including only a single dipole template in Eq. 10,with a direction equal to the CMB dipole. A second example arethe high-frequency channels above 353 GHz, which are stronglysignal dominated and may be assumed to have lower relativedipole errors than the CMB channels. Imposing a zero prior onthese components may be justifiable. However, in this paper,we fit explicitly for all components without such priors, both todemonstrate the method and to derive minimal-assumption andconservative estimates for all channels.
4. Results
A complete summary of our main results is presented in Table 1,listing the monopole and dipole coe ffi cients for all sky maps con-sidered in the analysis. The top section shows the results fromthe main 17 frequency analysis, and the bottom section the re-sults from the reduced 14 frequency analysis. The uncertaintiesare taken to be the maximum of the bootstrap errors discussedin Sect. 2.4 and 1 µ K (see Section 2.5); no uncertainties are re-ported for the 14 frequency analysis, as the presence of externalpriors make these very hard to assess.When interpreting these results, it is important to remem-ber that the algorithm has been explicitly designed to leavefeatures that are spatially varying (i.e., Galactic structures) un-changed, while any isotropic signal is fitted out. Thus, the fit-ted monopole consists of the sum of any instrumental and dataprocessing o ff sets and any Galactic or extra-Galactic componentthat is spectrally uniform over the full sky. Three typical exam-ples are the CMB monopole of 2.73 K, the mean value of theCosmic Infrared Background (CIB), and the mean generated byextra-Galactic point sources. Note, though, that the zero-level ofthe Galactic foregrounds are not necessarily fitted out, becausethese components have a spatially varying spectral index, and amonopole at one frequency therefore does not correspond to amonopole at other frequencies. Only components that are well
9. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression
30 100 300 1000 3000
Frequency (GHz) A n t e nn a t e m p e r a t u r e ( µ K ) Without offset correctionsWith offset corrections
Fig. 6: Mean antenna temperature as a function of frequency,measured outside the union mask before (red points) and af-ter (black points) o ff set corrections. The dashed lines shows thebest-fit sum of a low-frequency power-law (with free spectralindex) and a one-component greybody (with fixed emissivity of β d = . T K =
18 K; Planck Collaboration XI2013; Planck Collaboration XII 2013; Planck Collaboration Int.XVII 2014; Planck Collaboration Int. XXII 2014) to each of thetwo cases.approximated as monopoles at all frequencies are removed byour fit.In Fig. 5 we show the di ff erence between each of the CMBand high-frequency channels with the 70 GHz channel beforeand after applying the monopole and dipole corrections. The70 GHz channel is chosen as a reference because it has the low-est foreground contamination, and each di ff erence map shouldtherefore ideally be dominated by red colors. Several interestingpoints can be seen by eye in this plot: First, we see significantrelative dipoles in many of the pre-correction maps. Just a fewexamples are WMAP Q-band, Planck
LFI 44 GHz, and
Planck
HFI 100 GHz. After applying our dipole corrections, no suchclear features are seen any more. Second, we see that severalof the maps are dominated by blue colors, suggesting a rela-tive monopole o ff set. Again, after applying our monopole cor-rections, such features are no longer visible.In Fig. 6 we plot the mean temperature of the same maps out-side the union mask adopted in this analysis, adopting antennatemperature units, both before and after o ff set corrections. Thedashed lines show the best-fit sum of a low-frequency power-lawand a high-frequency one-component greybody component toeach of the two data sets. While the pre-correction mean exhibitsquite random behaviour between channels, the post-correctionmean follows quite well the expected physical behaviour. Thebest-fit low-frequency power-law indices are β s = − .
36 and − .
87, respectively.As described in Sect. 2.4, we adopt a straight template-fitprocedure for the low-foreground frequencies between 33 and100 GHz. A potential worry is therefore that a spatially vary-ing foreground spectral index can leak low-multipole power intothe monopole and dipole coe ffi cients. To get an intuitive feelingfor the magnitude of this e ff ect, we show in Fig. 7 the resid-ual maps obtained by subtracting the best-fit templates (CMB,monopole, dipole, low-frequency foregrounds /
30 GHz and ther-mal dust /
143 GHz) from each frequency map. While correlated large-scale features are indeed seen in these figures, indicatingthe presence of spatial variations, it is important to note that thecolor scale ranges between −
10 and 10 µ K, and the magnitudesof these features are therefore small.Next, in Fig. 8 we compare the pre- and post-correction408 and 1420 MHz maps. The single most striking feature inthis plot is a clear dipole extending from south to north in the1420 MHz map. Converting the Cartesian dipole coe ffi cients forthe 1420 MHz map in Table 1 into spherical coordinates, wefind that the best-fit dipole is 0 . ± . l , b ) = (308 ◦ , − ◦ ) ± ◦ . Interestingly,this direction is consistent with the Equatorial south pole, ( l , b ) = (303 ◦ , − ◦ ), possibly suggesting that the observed dipole mightbe interpreted e ff ectively in terms of an declination dependento ff set.After o ff set corrections, one can still see hints of an east-west type dipole in both the 408 and 1420 MHz maps. It is notdirectly obvious whether this feature is physical or not, as thedipole X and Y dipole coe ffi cients for these maps are quite large(as well as correlated), and if for instance the Y-dipole in the408 MHz is shifted by only one standard deviation, from 0.7 to2.1K, the visual impression becomes far more symmetric. Onthe other hand, it is worth noting that the low-frequency WMAPpolarization maps shows a similar asymmetry (see, e.g., Fig. 4 inBennett et al. 2013), due to the presence of strong synchrotronradiation near the Fan region.In Fig. 9 we compare the spectral index between the 408 and1420 MHz maps before and after o ff set corrections, as estimatedfrom the T–T plot distributions. Note that since monopoles donot change the index at all, these di ff erences are all due to thedipole correction. Thus, the typical dipole-induced bias on themean spectral index as computed over N side = ± . σ when including16% more sky area, with the biggest di ff erence is seen for the408 MHz monopole at 1 . σ .On the high frequency side, we see that the 353 GHzmonopole of 315 ± µ K derived within the 17 band solution is inexcellent agreement with HI cross-correlation result of 308 µ K(Planck Coll. VII 2013). The 353 GHz dipoles shows larger vari-ations at the 2 σ level. Imposing an HI prior also on this compo-nent should prove useful for the 14-band analysis.For the CMB-dominated frequencies, we see that the dipolesare overall in good agreement in terms of absolute numbers, de-spite the fact that the WMAP dipoles are forced to zero. Thisrobustness gives added credibility to the derived Planck dipoles.On the other hand, it is also clear that the reported uncertaintiesare too small, as already seen in simulations.Finally, we note relatively large di ff erences in the monopolevalues for frequencies between 23 (K-band) and 94 GHz (W-band), both for WMAP and Planck . This may be linked to therelatively large K-band dipole of 16 ± µ K in the 17-band analy-sis; when forcing this value to be zero, the monopole at the samefrequency increases from −
14 to 27 µ K, and this must in turn beaccommodated by increased monopoles at higher frequencies.
5. Summary and conclusions
One of the most di ffi cult tasks in physical CMB component sep-aration is the determination of absolute o ff sets, i.e., spurious
10. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression −
10 10WMAP Ka-band [ µ K] −
10 10WMAP Q-band [ µ K] −
10 10Planck 44 GHz [ µ K] −
10 10WMAP V [ µ K] −
10 10Planck 70 GHz [ µ K] −
10 10WMAP W-band [ µ K] −
10 10Planck 100 GHz [ µ K] Fig. 7: Map residuals after subtracting the best-fit template set (CMB, monopole, dipole, synchrotron /
30 GHz and thermal dust / ff sets maybias the estimation of spectral parameters significantly, and thiscan in turn lead to errors in the actual CMB map. Ideally, the op-timal approach would be to fit the o ff sets and foreground modeljointly, as for instance implemented by Gibbs sampling (Eriksenet al. 2008); however, if the data model contains a large numberof spectral indices, say, one per square degree pixel, the system isoften not su ffi ciently rigid to uniquely determine the optimal so-lution. The o ff sets are nearly degenerate with the e ff ective zero-level of each foreground component, and the only feature thatbreaks this degeneracy is spatial variation in the spectral param-eters.In this paper, we have presented an alternative method forestimating spurious monopoles and dipoles in multi-frequencydata sets. This method builds on a well-established methodol-ogy from the radio astronomy literature called T–T plots. Themain advantages of this method over the Bayesian approach arethat 1) it makes minimal assumptions about the nature of thesignal components, 2) it is computationally cheap, and 3) it istrivial to tune the number of regions to the point that the degener-acy between the spectral indices and the o ff sets are broken. Thelatter can of course also be implemented within the Bayesianframework, at which point we expect the two methods to per-form similarly. The main disadvantage of the current method isa relatively large systematic uncertainty when no signal com-ponent dominates, which for CMB purposes typically happensnear the foreground minimum at 70 GHz. In the present paper,we have adopted a straight template-fit approach for these fre- quencies, but note that a true multi-component fit is certainlypreferable. The most likely application for the current method istherefore to set the o ff sets at the foreground-dominated frequen-cies, which will then serve as an “anchor” for the Bayesian fit,e ff ectively breaking the degeneracies between the o ff sets and thespectral parameters.The main products presented in the paper are two di ff erent,but each internally consistent, sets of monopole and dipole coef-ficients. Overall, the two sets agree well with each other, exceptfor a single common monopole extending from 23 to 94 GHz.This is largely due to the significant systematic uncertainty asso-ciated with the WMAP K-band and Planck
30 GHz o ff sets. Werecommend that methods employing our o ff set values for sub-sequent analyses consider both sets for systematic uncertaintyassessment.An early version of the method presented in this paper wasalready adopted by the 2013 Planck release to determine thezero-levels for the physical component separation results, butonly applied to the
Planck frequencies between 30 and 353 GHz(Planck Coll. XII 2013). In this paper we have extended the totaldata set to also include the WMAP frequencies, the
Planck µ m SFD map. Overall, theresults presented here are in good agreement with the earlier re-sults, with a largest relative monopole di ff erence of 2 µ K for allchannels between 30 and 217 GHz. The only significant outlier isthe 353 GHz channel, for which we derive a value of 315 ± µ K,whereas Planck Coll. XII (2013) obtained a value of 414 µ K. Athird and fully independent estimate of this number was pro-
11. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression
Input Haslam 408 MHz - 8.9K[K]
Corrected Haslam 408 MHz [K]
Input Reich & Reich 1420 MHz - 3.2K[K]
Corrected Reich & Reich 1420 MHz [K]
Fig. 8: Comparison of the 408 (top) and 1420 (bottom) MHz low-frequency maps before (left) and after (right) o ff set corrections.vided by Planck Coll. VII (2013) based on cross-correlationwith HI observations, who reported a value of 308 ± µ K, in ex-cellent agreement with our result. Similarly, at 857 GHz PlanckColl. VII (2013) reports a value of 0 . ± . / sr fromHI measurements, which is to be compared with our fully in-ternal estimate of 0 . ± .
03 mK; for the 545 GHz the corre-sponding numbers are 0 . ± .
017 to be compared with ourvalue of 0 . ± .
01. Other channels are also in good agreement.However, for the SFD 100 µ m map, we note that the X-dipolehas both a large value and a large uncertainty, suggesting a lessconstrained fit overall. This is not unexpected, as this particularchannel is only coupled to the Planck
857 GHz map through along frequency extrapolation. The o ff set values for this channelare clearly less robust than for the HFI channels, and its roleis primarily to stabilize the 857 GHz results, rather than deriveindependent and robust o ff sets for the SFD map itself.On the low-frequency side, the 408 MHz and 1420 MHzmonopoles have already been subject of considerable discussionin the literature. Haslam et al. (1982) estimated that the zero-level uncertainty of the 408 MHz survey was ± ± . . − .
13 K for the 1420MHz survey, where the uniform background is 2 . − . ff set for the 408 MHz survey, when assuming that the zero-level at 1420 MHz is correct. Later, Reich et al. (2004) usedan improved source-count correction, which results in a uniformbackground at 408 MHz of 5 . ± . − . ff sets of − . ± . − . ± .
14 K at1420 MHz. All together these studies show that the zero-level ofthe 408 MHz survey is too low and requires corrections between + . + . + .
13 K. All these studies assume that the remain-ing extended emission in the survey maps is of Galactic origin.Spectral index maps (Reich & Reich 1988; Lawson et al. 1987,e.g.,) do not contradict this assumption. However, the minimain the survey maps at 408 MHz and 1420 MHz, with the zero-level, isotropic source-count and CMB corrections by the Reichet al. (2004) values, are about 9 . . ◦ angular res-olution. This means that there is room for a Galactic contribu-tion to a monopole and higher order components, but also fora larger isotropic extra-galactic component than calculated fromsource counts, as discussed in Sun & Reich (2010). Without anycorrection the survey minima are 12 . . ◦ reso-lution, respectively, which constrain the monopole componentswhen derived directly from the survey data. Recently, Bennettet al. (2013) used the same co-secant method as they applied tothe WMAP CMB frequencies to derive a background value of7 . . ± . . ± .
025 K at 1420 MHz. They concludethat their results agree with the extra-galactic component from
12. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression -4 -2
Before dipole correction -4 -2
After dipole correction -0.4 0.4
Difference
Fig. 9: Estimates of the spectral index derived from the combina-tion of the 408 and 1420 MHz maps using T–T plots, both before(top) and after (middle) applying the o ff set corrections. The bot-tom panel shows the di ff erence. These di ff erences show directlythe magnitude of the bias in T–T plot based spectral indices dueto spurious dipoles.ARCADE 2 (Fixsen et al. 2011). Subramanian & Cowsik (2013)also modelled the Galactic disk and halo emission and showedthat the simple co-secant method used to fit Galactic emission isthe reason for the ARCADE 2 excess and that there is no need foran isotropic extra-galactic component beyond that from sourcecounts. The currently deepest source count data by Vernstrom etal. (2014), which also include di ff use extra-galactic emission at1.75 GHz, give a strong indication that the excessive temperaturefound by ARCADE 2 is not of extra-galactic origin. Our current monopole values of 8 . ± . . ± .
02 K at1420 MHz are in the possible range allowed by the survey data.While the method presented in this paper is very general,and can deal with spurious monopole and dipoles of any origin,it is worth noting that by far the dominant source of dipole un-certainty in current CMB maps comes from estimating the solarCMB dipole. For instance, even small calibration uncertaintiescan lead to a significant uncertainties in the recovered dipole.However, this uncertainty is perfectly correlated between chan-nels within a given experiment, and it is therefore possible toimpose the prior that the corrections should be identical acrossfrequencies.Finally, we conclude with two important caveats. First andforemost, it is important to realize that while the method pre-sented here is extremely e ffi cient at establishing relative o ff setsbetween channels (which by far is the most important problemfor most component separation algorithms), it requires both highsignal-to-noise observations and significant spatial spectral vari-ations across the sky to determine absolute o ff sets. With the datasets studied in this paper, it appears that these criteria hold forboth the low-frequency 408 and 1420 MHz maps and the high-frequency channels above 100 GHz, but not for the intermedi-ate CMB channels between 23 and 94 GHz. Again, as shown inSection 3 the systematic uncertainty on the WMAP K-band o ff -sets is large. More conservatively, we caution against adoptingany of the derived o ff sets for frequencies between 23 and 44 GHzwithout further cross-checks because of the presence of multiplesignificant foreground components (i.e., synchrotron, free-freeand AME). Instead a full parametric fit is should be used forthese channels. Second, we emphasize that accurate uncertaintyestimation is di ffi cult, because the dominant sources of uncer-tainties are generally systematic in nature. In the present paper,we have adopted an internal bootstrap approach, which is able tocapture some, but not all, of these uncertainties. For more real-istic systematic uncertainty assessment, we recommend propa-gating both o ff set sets provided in this paper through subsequentanalyses. Acknowledgements.
We thank Greg Dobler for useful discussions. The de-velopment of Planck has been supported by: ESA; CNES and CNRS / INSU-IN2P3-INP (France); ASI, CNR, and INAF (Italy); NASA and DoE(USA); STFC and UKSA (UK); CSIC, MICINN, JA and RES (Spain);Tekes, AoF and CSC (Finland); DLR and MPG (Germany); CSA(Canada); DTU Space (Denmark); SER / SSO (Switzerland); RCN (Norway);SFI (Ireland); FCT / MCTES (Portugal); PRACE (EU). A description ofthe Planck Collaboration and a list of its members, including thetechnical or scientific activities in which they have been involved,can be found at . Part of the research was carried outat the Jet Propulsion Laboratory, California Institute of Technology, under acontract with NASA. This project was supported by the ERC Starting GrantStG2010-257080. Some of the results in this paper have been derived using the
HEALPix package.
References
Bennett, C. L., Hill, R. S., Hinshaw, G., et al. 2003, ApJS, 148, 97Bennett, C. L., Larson, D., Weiland, J. L., et al. 2013, ApJS, 208, 20Davies, R. D., Watson, R. A., & Gutierrez, C. M. 1996, MNRAS, 278, 925Deming, W. E. 1943. “Statistical adjustment of data.” Wiley, NY, ISBN 0-486-64685-8Draper, N. R., & Smith, H. 1998 “Applied Regression Analysis”, 3rd edition,Wiley, NY, ISBN 0-471-17082-8Eriksen, H. K., Jewell, J. B., Dickinson, C., et al. 2008, ApJ, 676, 10Fixsen, D. J., Kogut, A., Levin, S., et al. 2011, ApJ, 734, 5Fornengo, N., Lineros, R.A., Regis, M., Taoso, M., 2014, arXiv:1402.2218v2G´orski, K. M., Hivon, E., Banday, A. J.,Wandelt, B. D., Hansen, F. K., Reinecke,M., Bartelman, M. 2005, ApJ, 622, 759Haslam, C. G. T., Salter, C. J., Sto ff el, H., & Wilson, W. E. 1982, A&AS, 47, 1
13. K. Wehus et al.: Monopole and dipole estimation for multi-frequency sky maps by linear regression