https://orcid.org/0000-0002-6274-1424
Abstract
Foreground removal techniques for CMB analyses make specific assumptions about the properties of foregrounds in temperature and in polarization. By investigating the statistics of foreground components more understanding about the degree to which these assumptions are valid can be obtained. In this work we investigate E- and B-mode maps of the two strongest polarized foregrounds, synchrotron and thermal dust emission, with regards to their similarity with Gaussian processes, their spectral variations, and cross-correlations. We perform tests in patches of ∼3.7° size collectively covering the full sky and find most of them to conform to their Gaussian expectation according to the statistics in use. Correlations exhibit distinct differences in E- and B-mode signals, which point towards necessities in foreground removal methods. We discuss potential consequences and possible further directions.
1 INTRODUCTION
The polarization of the Cosmic Microwave Background (CMB) offers unique insights into the dynamics of light and matter around redshift z ∼ 1100, the epoch of re-ionization around z ∼ 6, and possibly primordial gravitational waves. (For a recent review of the corresponding physics and experimental achievements, see Staggs, Dunkley & Page 2018.) This latter contribution is conveniently characterized by the tensor-to-scalar ratio r, whose detection is yet outstanding. While the most recent analyses result in upper bounds of r < 0.07 (BICEP2/Keck Collaboration et al. 2015; Planck Collaboration VI 2018a), future surveys forecast sensitivities down to r ∼ 10−3 or even 10−4 (e.g. Abazajian et al. 2016), based on our current knowledge of those sources interfering with a clean measurement. Galactic foregrounds and systematics, especially in polarization, seem to offer the greatest challenges, and the corresponding methods for their removal (Planck Collaboration XLVI 2016c; Planck Collaboration IV 2018c) are perpetually tailored to our continuously improving understanding of their influences.
The Galactic foreground components contributing most to the polarized microwave and millimetre sky are synchrotron radiation, at low, and thermal dust emission, at high frequencies; others are presumably only polarized at the per cent level. Many properties of the polarized synchrotron and thermal dust skies are already under thorough investigation (Krachmalnicoff et al. 2018; Planck Collaboration XIX 2015a; Planck Collaboration XXII 2015b; Planck Collaboration XI 2018d). For instance, out of immediate relevance to the extraction of cosmological parameters from the CMB, those properties related to their power spectra have been of great interest (Planck Collaboration XXX 2016b), of which the surprisingly constant ratio of approximately 2 of the dust’s E- to B-mode power (Planck Collaboration XXXVIII 2016a) at intermediate to high multipoles has sparked considerable attention (Caldwell, Hirata & Kamionkowski 2017; Kandel, Lazarian & Pogosyan 2018; Kritsuk, Flauger & Ustyugov 2018). While also of astrophysical interest, statistical investigations of Galactic foregrounds can also assist their very removal to obtain a sufficiently cleaned CMB map.
Searches for and detection of residual foreground emission in CMB products have been practiced already since the first release of the WMAP data (e.g. Naselsky, Doroshkevich & Verkhodanov 2003, 2004; Dineen & Coles 2004) and more recently on Planck data (e.g. von Hausegger et al. 2016), highlighting potential downsides of assumptions made in foreground removal algorithms. Besides the general hope for reducing residual contamination in the final product, meticulous searches for statistical peculiarities in the CMB (Planck Collaboration XVI 2016e), e.g. primordial non-Gaussianity (Bartolo, Matarrese & Riotto 2010), must be safe from distortion of these very distributions from contaminants, such as systematic noise and/or foregrounds. Only recently more focus has been placed on understanding the relationship of foreground’s statistics with those which are used to search for non-Gaussianity in the final CMB products, see e.g. Hill (2018); Jung, Racine & van Tent 2018; Coulton & Spergel 2019. In this work we would like to take a step back and investigate the statistical behaviour of foregrounds in a general manner, independent of studies of specific CMB non-Gaussianities, in order to, among others, finally make connections to general issues in foreground removal methods. Also this has been of recent interest (Ben-David, von Hausegger & Jackson 2015; Rana et al. 2018) regarding temperature maps of Galactic synchrotron radiation, which calls for similar studies of maps of polarization.
At present, foreground separation techniques place requirements on certain properties of the foreground components – some explicitly, such as parametric foreground fitting algorithms (Eriksen et al. 2008), and some implicitly, such as so-called blind foreground removal techniques (Bennett et al. 2003; Cardoso et al. 2008; Delabrouille et al. 2009). By the example of the latter, we shall, in Section 2 of this paper, provide motivation to analyse statistics of these foreground components regarding their Gaussianity and their relations to each other. The evaluation of all such algorithms is most effectively assisted by a set of simulated foregrounds, which in turn depend on assumptions about their statistics: The simulation of certain foregrounds’ properties is often underpinned by assuming Gaussian behaviour (Tegmark & Efstathiou 1996; Tegmark et al. 2000), and small-scale structure is commonly added into existing maps by simply generating Gaussian fluctuations according to extrapolations of the available power spectra (Remazeilles et al. 2015; Hervías-Caimapo, Bonaldi & Brown 2016). As accurate simulations are instrumental in making realistic predictions for forthcoming surveys, exploiting possible Gaussian patterns in foregrounds’ statistics will prove itself useful. In other words, understanding foregrounds’ statistics will enable to justify (and expand on) assumptions in current foreground simulation procedures, such as those of Thorne et al. (2017) or Delabrouille et al. (2013). Some recent effort has already gone into the determination at which scales the assumption of Gaussianity in foregrounds holds (Ben-David et al. 2015; Rana et al. 2018). These findings are to be extended to studies of polarization, which for the first time is attempted in this work, by the investigation of Galactic foreground E- and B-mode maps.
In short, the aim of this paper is to elucidate assumptions inherent to current foreground removal techniques and to, once convinced that foreground residuals will be present in CMB maps at some level, generally characterize their statistics.
After a general motivation in Section 2, our procedure for classifying Gaussianity is briefly reviewed in Section 3. In Section 4 we then present the results from E- and B-mode maps at 23 and 353 GHz, corresponding to the polarization signal from synchrotron and from thermal dust emission, respectively. We conclude in Section 5.
2 MOTIVATION
To illustrate the necessity of studying the statistics of foregrounds, we provide the following motivation. A commonly used framework to understand so-called blind foreground subtraction algorithms is that of Internal Linear Combination, or ILC, which has been introduced to CMB science by the WMAP collaboration (Bennett et al. 2003). Sky maps of different frequency are linearly combined with the aim to arrive at a map of the CMB anisotropy. More sophisticated algorithms are in use today, see e.g. Planck Collaboration IX (2016d), albeit they produce CMB products consistent with an ILC solution. In this section we shall briefly review essential elements and assumptions of the ILC method, as a representative of such blind foreground removal techniques.
In this rough exploration we have not considered noise terms. However, their inclusion would only lead to the addition of the noise covariance matrix, and the corresponding cross-covariance matrices in equation (2), and not change the arguments about the foreground properties made above. In particular, equations (8)–(9) look equivalent for non-correlated noise terms. Yet there is an essential difference between uncorrelated foreground and uncorrelated noise: By the continuing improvement in detector technology, one may obtain lower noise levels and thereby reduce the influence of this term, the contribution from uncorrelated foregrounds, however, remains the same. In the case of Gaussian, uncorrelated foregrounds (or noise) equation (8) completely describes their properties, and in the light of the previous discussion we find Gaussianity in a region Ω to arise from either the foreground ‘template’ F(p) itself, or from the direction-dependent term Δν(p), while the respective other remains nearly constant. In particular, any residual of the sort described in equation (5) will itself be Gaussian. The later distinction of such contamination from the also Gaussian CMB will be challenging.
It is these considerations which motivate a spatially resolved investigation of foregrounds with regards to their statistics – in specific, we shall investigate their similarity (or dissimilarity) to Gaussian variables in small patches distributed over the sky, in addition to studying correlations between foreground maps on the same scales. In a previous study (Ben-David et al. 2015), we have performed such tests on a full-sky temperature map at 408 MHz, representative for Galactic synchrotron radiation; we here extend the analysis to polarized foregrounds. Even though formulated for temperature fields, above considerations hold also for the removal of polarization foregrounds.3 In the remainder of this paper we will investigate the E- and B-modes of full-sky polarized foreground maps.
3 HOW TO CLASSIFY GAUSSIANITY
There is a multitude of ways in which distributions can be investigated with regards to their shape. In particular Gaussian distributions (or deviations thereof) have been in the prime focus of such methods ( D’Agostino 1986). As motivated above, we here are interested in the (dis)similarity of foregrounds with such, that follow Gaussian distributions. While no estimator can detect any sort of non-Gaussianity with equally high efficiency, and we do not have any preconception with regards to a specific type of distribution, we restrict ourselves to the third and fourth normalized statistical moments, the skewness and kurtosis, as done in our earlier work.
Example demonstrating the method. Top panel: The skewness map of the K-band’s E-mode. Middle panel: Histogram of both skewness, γ1, and kurtosis, γ2, from the patches of 1000 simulated sky maps based on the power spectrum of the K-band’s E-mode. The contours delineate 1, 2, and 3σ. Bottom panel: The skewness–kurtosis histogram for the K-band’s E-mode map. The contours from the histogram in the middle panel are used to discern Gaussian from less Gaussian patches.
Spatially localized structures, i.e. correlations in the pixel domain, will generally lead to correlations in the Fourier phases, and even a Gaussian, statistically homogeneous and isotropic signal, i.e. where the phases are drawn at random from a uniform distribution, can still have ‘preferred clustering’, as, for example, is expected for the CMB.6 As we wish to perform our tests in the pixel domain, we must establish a fair measure by which to estimate the degree of (or departure from) Gaussianity, given the computed values of skewness and kurtosis. We do this by generating simulated sky maps based on the power spectrum of the investigated sky map, i.e. they all have the same power spectrum as the input map, while their phases are selected from a uniform distribution as noted above. A set of 1000 such realizations provides us with |$12000\times n_s^2$| skewness and kurtosis values, which we cast into a 2-dimensional histogram (cf. equation 11) to approximate the skewness–kurtosis probability distribution from which Gaussian samples would be drawn. Contours on this distribution will serve as guidelines when deciding on the degree of Gaussianity of the input sky map’s patches, as shown in the middle and bottom panels of Fig. 1. More detail on the computations and background information can be found in Ben-David et al. (2015).
4 RESULTS
As described above, we are interested in extending our previous formalism of Gaussianity checks to the E- and B- polarization maps, as well as performing correlation studies on them. The maps under investigation here are the E and B maps of the WMAP K-band map (Bennett et al. 2013) and the Planck 353 GHz map (Planck Collaboration III 2018b) primarily, which in polarization well describe synchrotron and thermal dust polarization, respectively.7 For the correlation studies we also utilize the WMAP Ka-band polarization maps and those of the Planck 217 GHz map. All maps are smoothed with a 1° beam, prior to all calculations. The results we obtain in this manner can point out departures from Gaussian statistics of the foregrounds for the scales investigated and, in addition, locate these very departures on the sky. Also the correlations are performed locally via mosaic correlations.
In this context it is important to emphasize that finding any one region on the E- or B-mode sky to be particularly Gaussian (or non-Gaussian) does not equate to stating that physical processes in the same direction contribute to this finding. This is due to the non-local definition of the E- and B-modes: Each point on the E/B sky receives contributions from all other points on the Q/U sky, except its own direction (and its antipole; for a nice explanation, see Rotti & Huffenberger 2019). One might argue, that therefore the many contributions are expected to render the signal in any particular direction Gaussian by the central limit theorem. None the less, the introductory arguments and searches for non-Gaussianity in the CMB’s E- or B-modes motivate a characterization of the foregrounds along similar lines.
As in our previous work, Ben-David et al. (2015), we present results from investigating maps at an NSide = 512 in patches of NSide = 16, which corresponds to scales of approximately 3.7°. However, this methodology can be applied to patches of different size (and shape) for the analysis of different scales. Despite smoothing, the presence of noise in the polarization maps of both WMAP and Planck influences the signal at high Galactic latitudes. While we demonstrate in appendix B that this does not severely affect the results in skewness and kurtosis (subSection 4.1; also Figs A3 and A4), we shall see the clear consequence of this in subSection 4.2.
4.1 Skewness and kurtosis
We calculate skewness and kurtosis values for the E and the B maps in NSide = 16 pixels as shown by the example in the middle panel of Fig. 1. As described in the previous section, we quantify the degree of non-Gaussianity of the respective patches on the sky as their deviations from an ensemble of Gaussian simulations. As demonstrated in the bottom panel of Fig. 1 for the E map, we hereto overlay the resulting 2-dimensional skewness–kurtosis histograms of the E or B maps with contours from the corresponding simulations. Using the simulations’ probability distribution functions we can assign a probability to each pair of skewness and kurtosis values, and thereby to each patch. By doing this we compute maps of departures from the Gaussian expectation in units of standard deviations for both E- and B-modes. (In Figs A3 and A4 we show the corresponding maps for Q- and U-modes for completeness.)
4.1.1 Synchrotron polarization
We show maps of the standard deviations computed for the E- and B-modes of WMAP’s K-band map in Fig. 2. Most patches deliver values below 2σ and therefore, according to our estimator, qualify as consistent with arising from a Gaussian sample. Apart from hints towards peculiarities along the Galactic plane the distribution of values further does not seem to prefer certain regions on the sky. This finding resembles the one for the 408 MHz Haslam map, shown in Ben-David et al. (2015).
Deviations of the skewness and kurtosis values from their Gaussian expectation in units of standard deviations, σ, quantified by the Gaussian realizations as described in the text. The values were computed from the WMAP K-band E- (top panel) and B-mode (bottom panel) maps in patches of NSide = 16. The standard deviations are calculated using the probability distribution functions obtained from the simulations.
4.1.2 Dust polarization
Equivalent to the synchrotron results in the previous section we compute those of Planck’s 353 GHz polarization maps for the same patch sizes. These appear to be very similar, see Fig. 3, apart from a slight surplus of high-skewness patches, which lead to a correspondingly higher abundance of >2σ patches. This tendency is also observed in dust temperature maps and can be traced back to the presence of pronounced filaments and point sources in the dust maps. While point sources are not expected to be intrinsically polarized, filaments can be visible on polarization maps. In particular, the E-mode will inherit power from polarized filaments, due to its sensitivity to elongated structures (see e.g. Liu et al. 2018), and indeed we find the 353 GHz E-mode map to give |${\approx }26{{\ \rm per\ cent}}$| more patches with values above 2σ than the B-mode map.
Same as Fig. 2, just for the polarization of the Planck 353 GHz map.
4.2 Correlations
In Section 2 we discussed the possibility of spatially varying spectral indices, or, more generally, spatial variation in the scaling coefficient of a ‘template’ per foreground component. Another point of discussion was formed by regions Ω on the sky within which foreground or noise signal are uncorrelated from band to band. Both lead to terms in the |$\mathcal {S}$| map’s variance which prevent a minimization which is independent of sky location, and thereby constitute potential contamination in the final CMB product. We shall here investigate the same patches with regards to these properties which previously were tested for Gaussianity.
As before, we at first focus on low- and high-frequency foregrounds separately. We each consider correlations between the signals measured in two bands, which are dominated by the same foreground mechanism. To make comparison to the results of the previous subsection, also here we compute Pearson correlations in patches of NSide = 16 – a method known as mosaic correlation (as presented by Verkhodanov, Khabibullina & Majorova 2009; and used in von Hausegger & Liu 2015, for the case of correlations among foreground temperature maps). Any correlation short from perfect, i.e. any change in a signal from one frequency to one nearby can only be induced by spatial variations of the spectral index (or of an equivalent scaling of the foreground), or by instrumental noise; for the frequencies considered, no other polarized foreground components nor CMB should play a significant role besides synchrotron or thermal dust.
4.2.1 Synchrotron polarization
To detect changes in the synchrotron sky in nearby frequencies we compare the WMAP K-band (23 GHz) with the Ka-band (33 GHz) polarization maps. Neither of free–free emission, spinning dust emission, nor molecular line emission, is expected to be polarized at a level comparable to synchrotron emission over most of the sky, such that we can assume both maps to carry signal from either synchrotron emission or noise. We present the resulting mosaic correlation maps, |$K_{\mathrm{K}\times \mathrm{Ka}}^E$| and |$K_{\mathrm{K}\times \mathrm{Ka}}^B$|, in the left-hand panels of Fig. 4. Distinct is the area of most pronounced correlations along the Galactic plane, slightly extending also along the North Polar Spur, especially for the E-mode maps. While positive correlations overweigh, at intermediate to high Galactic latitudes the correlations weaken and expose many patches with correlations close to zero, see also Fig. A1. In order to understand the origin of this reduced correlation we investigated WMAP K- and Ka-band single-year maps, which showed that the dominating contributor at high Galactic latitudes is most likely noise. However, recall that it is irrelevant to the point presented in Section 2, whether this loss in correlation arises from changes in the foregrounds’ morphology or from uncorrelated instrumental noise in the two bands – both present challenges for ILC-like methods.
Mosaic correlations between E-mode (upper panels) and B-mode maps (lower panels). Left-hand panels: Correlations between WMAP K- and Ka-band maps. Middle panels: Correlations between Planck 217 and 353 GHz maps. Right-hand panels: Correlations between WMAP K-band and Planck 353 GHz maps.
4.2.2 Dust polarization
Applying the same method to the 217 and 353 GHz Planck maps to obtain results for the polarized thermal dust sky, leads to the correlations presented in the middle panels of Fig. 4. Similar to the case at low frequencies, high correlation can be observed close to the Galactic plane and lower to vanishing correlation towards high Galactic latitudes. The larger abundance of high correlation patches, see also Fig. A2, is largely due to higher signal-to-noise ratios in the maps (a conclusion also supported by analysing half-mission maps), but could also indicate generally lower amount of spectral index variation in thermal dust emission.
4.2.3 Correlations between synchrotron and dust
Before concluding we still need to consider correlations between synchrotron and dust emission. In the case of more than a single foreground, the variance |$\langle \mathcal {S}^2\rangle _\Omega$| will also contain cross-terms between the different foreground components, such that the weights wν need to be determined according to high or low cross-correlation. We show the corresponding mosaic correlation maps in the right-hand panels of Fig. 4 for both E- and B-modes of the WMAP K-band and the Planck 353 GHz map. The tendency in both is higher correlations along the Galactic plane. However, away from the Galactic plane no strong correlation between the polarized emission of synchrotron and thermal dust emission is established on the scales given by the patch size investigated here. We further point out that while the high Galactic latitudes are most likely, as before, affected by instrumental noise, noticeable differences between the correlations at low Galactic latitudes from E- and B-modes can be seen – those latitudes where instrumental noise is subdominant. We return to this point in the discussion. However, higher sensitivity observations will be needed to characterize the foregrounds’ polarized emission better both at low and at high frequencies.
5 DISCUSSION AND CONCLUSION
While also of astrophysical interest statistical investigations of Galactic foregrounds also have great relevance in regards to sufficient foreground removal, as we showed in this paper. In particular distinct statistical properties of E- and B-mode signals require an adequate, perhaps separate treatment in foreground separation algorithms.
In this work we investigated statistics of the two strongest CMB foregrounds in polarization, synchrotron radiation and thermal dust emission, in order to draw conclusions on the feasibility of obtaining a clean CMB map. In particular variations of foreground properties across the sky were of motivated interest wherefore we employed two methods which both work in pre-defined patches on the sky – the skewness-kurtosis method for classification of foreground maps as Gaussian processes (Ben-David et al. 2015), and the mosaic correlations (von Hausegger & Liu 2015). The E- and B-mode maps under investigation were smoothed to 1° and the selected patches were of extent ∼3.7°, which can be roughly translated into the multipole range ℓ ∈ [50, 180]. We summarize our findings as follows:
On scales between approximately 1° and 3.7°, the E- and B-mode maps of both synchrotron and thermal dust polarization exhibit distributions consistent with those expected from Gaussian realizations over most of the sky, with preferred regions of departure only along the Galactic plane, cf. Figs 2 and 3.
On the same scales we find spatial variation of the frequency spectra of both synchrotron and thermal dust polarization to be negligible along the Galactic plane, cf. Fig. 4.
At intermediate and higher Galactic latitudes increased instrumental noise prevents us from drawing conclusions about the spectral properties. However, the observed decrease in correlations also in the case of noise will impede ILC-like foreground removal algorithms, as elaborated in Section 2.
Mosaic correlations between synchrotron and dust polarization maps were distinctly different for E- and B-mode maps along the Galactic plane.
It should not be failed to mention that correlations between synchrotron and dust polarization have been subject to a many of recent studies that assess the potential level of foreground contribution at those frequencies where the CMB signal is strongest (Choi & Page 2015; Krachmalnicoff et al. 2016, 2018; Planck Collaboration XI 2018d). Their findings underline the existence of these correlations among the components also at high Galactic latitudes, and correlations were quantified for different scales by computing the corresponding cross-power spectra. However, they do not consider spatial variation of the correlation other than imposing Galactic masks of different extent. For the scales chosen in this work, we pointed out these spatial variations here, and, in addition to the differences between those from E- and B-modes, showed their relation to principles in current foreground separation algorithms.
In light of the increasingly precise measurements at microwave frequencies performed for the observation of the CMB, an equally precise understanding of the behaviour of those sources interfering with a clean measurement is required. Given that foreground separation algorithms and the subsequent inclusion of the CMB products in a combined framework for determining the cosmological parameters focus on the less contaminated regions of the sky, away from the Galactic plane, these results will need to be taken into account in further analyses of the results. In other words, the supposedly cleanest regions of the sky might contain the most stubborn foregrounds.
ACKNOWLEDGEMENTS
We thank Pavel Naselsky for helpful and exciting discussions. This work was funded in part by the Danish National Research Foundation (DNRF) and by Villum Fonden through the Deep Space project. HL is supported by the Youth Innovation Promotion Association, CAS. Some of the results in this paper have been derived using the healpix (Gorski et al. 2005) package implemented into healpy. We thank both the referee and the editor for their constructive comments.
Footnotes
It should be clear that in the case where equation (3; or its extension to n foregrounds) is fulfilled, the minimization of the variance becomes redundant. In the case of a single foreground only two observations at different frequencies are required to solve for the weights wν:
Indeed, with the onset of more precise measurements of the radio sky, it will become increasingly clear that the spectral index of synchrotron emission varies across the sky (Taylor (Krachmalnicoff et al. 2018; Taylor 2018), and also for thermal dust emission such variation has been observed and is known as de-correlation (Planck Collaboration L 2017; see however, Sheehy & Slosar 2018; Planck Collaboration XI 2018d).
In addition to the assumptions about constant spectral indices or the like, the ILC approach in polarization requires the polarization angles of each component to be constant across frequencies. It is not clear whether the ILC in polarization is best performed on the Stokes parameters Q and U, or their non-local transformations E and B. Depending on the particular method, both are in use. Also combinations of both in decompositions like those proposed in Liu, Creswell & Naselsky (2018) and Rotti & Huffenberger (2019) might be of interest. Also for this reason we show corresponding results for Q and U in the appendix, Figs A3 and A4.
In fact, the observed synchrotron sky arises from an interplay between many influencing factors, such as supernova rate and their spatial distribution, the cosmic ray electrons’ energy and their spatial distribution, the alignment and strength of Galactic magnetic fields, absorption effects of the interstellar medium, etc. In polarization one must further consider depolarization effects along the line of sight, such as the frequency-dependent Faraday rotation, which depends on most of the previously mentioned quantities in a non-trivial manner.
Further note that, any signal observed through an instrument is automatically convolved with the instrument’s beam function, which produces named clustering to below the beam’s size.
REFERENCES
APPENDIX A: ADDITIONAL FIGURES
Histograms of the mosaic correlations between WMAP K- and Ka-band E-mode (upper panel) and B-mode maps (lower panel), as shown in Fig. 4. The inset shows in blue those regions from which the blue distribution is plotted, corresponding to the WMAP polarization mask. The grey histogram contains the values from the full sky.
Same as Fig. A1, but for mosaic correlations between the Planck 353 GHz and the 217 GHz maps. The blue mask here is the Planck polarization mask. In addition we select patches via a more aggressive mask, the Planck Gal20 mask, shown in orange.
Deviations of the skewness and kurtosis values from their Gaussian expectation in units of standard deviations, σ, quantified by the Gaussian realizations as described in the text. The values were computed from the WMAP K-band Q- (top panel) and U-mode (bottom panel) maps in patches of NSide = 16. The standard deviations are calculated using the probability distribution functions obtained from the simulations.
Same as Fig. A3, just for the polarization of the Planck 353 GHz map.
