The Galactic Faraday depth sky revisited

TL;DR

This paper presents an improved method for reconstructing the Galactic Faraday depth sky using advanced imaging algorithms and additional emission data, revealing new insights into Galactic magnetic structures.

Contribution

It introduces a novel inference algorithm incorporating free-free emission data, enhancing the accuracy of the Galactic Faraday depth map reconstruction.

Findings

Reconstructed an amplitude field resembling the free-free emission map.

Identified regions with deviations indicating different Galactic physics.

Found evidence of magnetic field alignment in the Orion arm.

Abstract

The Galactic Faraday depth sky is a tracer for both the Galactic magnetic field and the thermal electron distribution. It has been previously reconstructed from polarimetric measurements of extra-galactic point sources. Here, we improve on these works by using an updated inference algorithm as well as by taking into account the free-free emission measure map from the Planck survey. In the future, the data situation will improve drastically with the next generation Faraday rotation measurements from SKA and its pathfinders. Anticipating this, the aim of this paper is to update the map reconstruction method with the latest development in imaging based on information field theory. We demonstrate the validity of the new algorithm by applying it to the Oppermann et al. (2012) data compilation and compare our results to the previous map.\\ Despite using exactly the previous data set, a number of novel findings are made: A non-parametric reconstruction of an overall amplitude field resembles the free-free emission measure map of the Galaxy. Folding this free-free map into the analysis allows for more detailed predictions. The joint inference enables us to identify regions with deviations from the assumed correlations between the free-free and Faraday data, thereby pointing us to Galactic structures with distinguishably different physics. We e.g. find evidence for an alignment of the magnetic field within the line of sights along both directions of the Orion arm.

Figures

Equations

$${\Delta }_{\lambda }=\varphi \left(n_{\mathrm{th}},B_{\mathrm{LOS}}\right){\lambda }^2,$$

$${\varphi }_{\mathrm{gal}}=\frac{e^3}{2\pi m_e^2{c}^4}{\int }_{\mathrm{LOS}}dl{n}_{\mathrm{th}}{B}_{\mathrm{gal},\mathrm{LOS}},$$

$$\varphi ={\varphi }_{\mathrm{gal}}+{\varphi }_{\mathrm{etc}}.$$

$$d_{\varphi }={\mathcal{R}}_{\varphi }\varphi +n_{\varphi },$$

$${\left[\mathcal{R}\varphi \right]}_i={\int }_{S^2}d\hat{e}\varphi \left(\hat{e}\right)\delta (\hat{e}-{\hat{e}}_i)$$

$$\mathcal{P}\left(n_{\varphi }\right)=\mathcal{G}\left(n_{\varphi },N_{\varphi }\right)=\frac{1}{\sqrt{|2\pi N_{\varphi }|}}{e}^{-\frac{1}{2}{n}_{\varphi }{N}_{\varphi }^{-1}{n}_{\varphi }^{†}},$$

$$\mathcal{P}\left(d_{\varphi }|\varphi \right)=\mathcal{G}\left(d_{\varphi }-{\mathcal{R}}_{\varphi }\varphi ,N_{\varphi }\right)$$

$$\varphi \equiv \chi e^{\rho }\frac{\mathrm{rad}}{{\mathrm{m}}^2},$$

$$\mathcal{P}\left(\chi ,\rho |S_{\chi },S_{\rho }\right)=\mathcal{G}\left(\chi ,S_{\chi }\right)\mathcal{G}\left(\rho ,S_{\rho }\right).$$

$$d_{\varphi }={\mathcal{R}}_{\varphi }({\varphi }_{\mathrm{gal}}+{\varphi }_{\mathrm{etc}})+n_{\varphi }\equiv R_{\varphi }{\varphi }_{\mathrm{gal}}+{\tilde{n}}_{\varphi }$$

$${\tilde{N}}_{\varphi }=\mathrm{diag}\left({\eta }_{\varphi }{\sigma }_{\varphi }^2\right),$$

$$\mathcal{P}\left({\eta }_{\varphi }|{\alpha }_{\varphi },{\beta }_{\varphi }\right)=\frac{{\beta }_{\varphi }^{{\alpha }_{\varphi }}}{\Gamma ({\alpha }_{\varphi })}{\eta }_{\varphi }^{-{\alpha }_{\varphi }-1}\exp \left(-\frac{{\beta }_{\varphi }}{{\eta }_{\varphi }}\right)$$

$${\mathrm{EM}}_{\mathrm{ff}}\left(n_{\mathrm{th}}\right)={\int }_{\mathrm{LOS}}dl{n}_{\mathrm{th}}^2\equiv {\int }_{\mathrm{LOS}}dl{\mathrm{em}}_{\mathrm{ff}},$$

$$B_{\mathrm{LOS}}{n}_{\mathrm{th}}$$

$$\propto n_{\mathrm{th}}^{1+p_B}\cos (\theta )\propto {\mathrm{em}}_{\mathrm{ff}}^{\frac{1+p_B}{2}}\cos (\theta ).$$

$$\varphi$$

$$\approx L^2⟨|B|n_{\mathrm{th}}{⟩}_L⟨\cos (\theta ){⟩}_L,$$

$$d_{\mathrm{ff}}={\mathcal{R}}_{\mathrm{ff}}{\mathrm{EM}}_{\mathrm{ff}}+n_{\mathrm{ff},\mathrm{obs}}+n_{\mathrm{ff},\mathrm{sys}}\equiv {\mathcal{R}}_{\mathrm{ff}}{\mathrm{EM}}_{\mathrm{ff}}+{\tilde{n}}_{\mathrm{ff}}.$$

$$\mathcal{P}\left(d_{\mathrm{ff}}\right)=\mathcal{G}\left(d_{\mathrm{ff}}-{\mathcal{R}}_{\mathrm{ff}}{\mathrm{EM}}_{\mathrm{ff}}|{\tilde{N}}_{\mathrm{ff}}\right),$$

$${\mathrm{EM}}_{\mathrm{ff}}\equiv e^{ϵ}\cdot \frac{\mathrm{pc}}{\mathrm{c}{\mathrm{m}}^6}$$

$$\varphi$$

$${\varphi }_{\mathrm{Orion}}=e^{\delta }\chi$$

$$\mathcal{P}\left(s|d\right)=\frac{\mathcal{P}\left(d|s\right)\mathcal{P}\left(s\right)}{\mathcal{P}\left(d\right)}$$

$$\mathcal{H}\left(s|d\right)=-\ln \mathcal{P}\left(s|d\right)$$

$$\mathcal{H}\left(\mathbf{s}|\mathbf{d}\right)\equiv -\log \left(\mathcal{P}\left(\mathbf{s}|\mathbf{d}\right)\right)$$

$$\sum _{d\in \mathbf{d}}\mathcal{H}\left(d|\mathbf{s}\right)+\sum _{f\in \mathbf{s}}\mathcal{H}\left(f\right)\text{stackon}[1.5pt]=\text{stretchto}\text{scalerel}\ast [width("=")]\wedge 0.5ex$$

$$\sum _{d\in \mathbf{d}}\frac{1}{2}\left(d-{\mathcal{R}}_d\mathbf{s}\right){\tilde{N}}_d^{-1}{\left(d-{\mathcal{R}}_d\mathbf{s}\right)}^{†}+$$

$$+\frac{{\beta }_d}{{\eta }_d}+({\alpha }_d+1)\ln ({\eta }_d)+\sum _{f\in \mathbf{s}}{f}^{†}{S}_f^{-1}f.$$

$$C_{\ell }^{(f)}={\mathbb{P}}_{\ell }{⟨\hat{f}{\hat{f}}^{\ast }⟩}_{{\mathcal{P}}_{(x)}}={\mathbb{P}}_{\ell }\frac{1}{2\ell +1}\sum _m|{\hat{f}}_{\ell m}{|}^2.$$

$$f=\mathbb{Y}{A}_f{\xi }_f.$$

$$S_f=\mathbb{Y}{A}_f{A}_f^{†}{\mathbb{Y}}^{†}.$$

$$f^{†}{S}_f^{-1}f={\xi }_f^{†}{\xi }_f,$$

Full text

11institutetext: Max Planck Institute for Astrophysics, Karl-Schwarzschildstr.1, 85741 Garching, Germany 22institutetext: Ludwig-Maximilians-Universität München, Geschwister-Scholl-Platz 1, 80539 Munich, Germany

The Galactic Faraday depth sky revisited

Sebastian Hutschenreuter

  

Torsten A. Enßlin

The Galactic Faraday depth sky is a tracer for both the Galactic magnetic field and the thermal electron distribution. It has been previously reconstructed from polarimetric measurements of extra-galactic point sources. Here, we improve on these works by using an updated inference algorithm as well as by taking into account the free-free emission measure map from the Planck survey. In the future, the data situation will improve drastically with the next generation Faraday rotation measurements from SKA and its pathfinders. Anticipating this, the aim of this paper is to update the map reconstruction method with the latest development in imaging based on information field theory. We demonstrate the validity of the new algorithm by applying it to the Oppermann et al. (2012) data compilation and compare our results to the previous map.

Despite using exactly the previous data set, a number of novel findings are made: A non-parametric reconstruction of an overall amplitude field resembles the free-free emission measure map of the Galaxy. Folding this free-free map into the analysis allows for more detailed predictions. The joint inference enables us to identify regions with deviations from the assumed correlations between the free-free and Faraday data, thereby pointing us to Galactic structures with distinguishably different physics. We e.g. find evidence for an alignment of the magnetic field within the line of sights along both directions of the Orion arm.

1 Introduction

The Faraday rotation effect is one of the primary sources of information on astrophysical and cosmological magnetic fields. This includes the fields of planets, stars, other galaxies and galaxy clusters, as well as more curious objects such as radio jets, lobes, and relics. In particular, the study of the Faraday rotation induced by the Milky Way magnetic field is of twofold importance. First of all it is an interesting research topic on its own, due its connection to the formation and structure of our home galaxy, but furthermore it also constitutes a non-negligible foreground component. All polarized light stemming from cosmological sources passes through the Galaxy and interacts with the Galactic magnetic field, which affects its polarization direction. Accurate and highly resolved Galactic foreground templates are therefore a necessary condition for any precision measurement of extra-galactic magnetic fields via Faraday rotation.

Multiple efforts have been made to to map the Galactic Faraday sky, quite a few of them already to remarkable accuracy, e.g. in Frick et al. (2001); Johnston-Hollitt et al. (2004); Dineen & Coles (2005); Xu et al. (2006); Oppermann et al. (2012) and also in Oppermann et al. (2015). Especially the approach in Oppermann et al. (2012) is important to us, as we will rely on the same theoretical framework for our inference algorithms.

Our goal in this work is to further sharpen our knowledge on the Faraday sky by exploiting correlations to bremsstrahlung measured from electron proton interaction in the interstellar medium, commonly known as the Galactic free-free emission. This is well motivated by observation as well as physical considerations, as we will outline later on. In a more general context, we hope to demonstrate an interesting test case of multi-messenger astronomy, and show that the synthesis of long existing data sets still can yield undiscovered information, under the condition that our a priori physical knowledge is employed and the data is treated in a consistent way. Our main result, the revised Faraday sky including the free-free data is shown in Fig. 1(a) and the revised map without free-free data in Fig. 1(b). The results of our inference are publicly available 111https://wwwmpa.mpa-garching.mpg.de/~ensslin/research/data/faraday_revisited.html.

The paper is structured as follows. In Sec. 2 we develop our updated reconstruction method. This involves an amplitude field aiming at representing the overall Galactic structure and a sign field representing magnetic reversals and smaller scaled structures. It turns out that the reconstructed amplitude field resembles the free-free emission measure (EM) of the Galaxy, which is physically plausible. In order to take advantage of this, we extend the method to fold in the measurements of the free-free sky in Sec. 3. We conclude in Sec. 4. Mathematical details of the method can be found in the Appendix.

2 Improving the Inference of the Faraday sky

2.1 The physics

The rotation angle $\Delta_{\lambda}$ of the polarization plane of linearly polarized light due to Faraday rotation can be written in most astrophysical settings as (see e.g. Oppermann et al. 2012, 2015)

[TABLE]

where $\lambda$ is the wavelength at which the light is observed. The Faraday depth $\phi$ depends on the thermal electron density $n_{\mathrm{th}}$ and the line of sight (LOS) component of the magnetic field $B_{\mathrm{LOS}}$. More specifically for a polarized extra-galactic source shining through the interstellar medium, the Galactic part of the Faraday depth $\phi_{\mathrm{gal}}$ can be written as (Oppermann et al. 2012)

[TABLE]

where $B_{\mathrm{gal}}$ is the galactic magnetic field, $e$ is the electron charge, $m_{e}$ is the electron mass and $c$ the speed of light. Unfortunately, contributions to the observed rotation angle not only come from the Galactic magnetic field but also contain extra-galactic, potentially source internal contributions. Similar to the approach taken by Oppermann et al. (2012), it will be helpful to separate the Faraday depth signal into

[TABLE]

The terms $\phi_{\mathrm{gal}}$ and $\phi_{\mathrm{etc}}$ describe the aforementioned Galactic and other e.g. extra-galactic, source intrinsic, or ionospheric contributions, respectively. Our goal in this work will be the estimation of $\phi_{\mathrm{gal}}$ as a correlated field on the sky. This implies the necessity for an at least stochastic estimation of the other contributions. The connection of $\phi$ to the observed data $d_{\phi}$ is the topic of the next section.

2.2 The data

A measurement of the Galactic Faraday depth can be modeled via

[TABLE]

where we ignore for the moment the potential non-galactic contributions mentioned in Eq. (3). The operator $\mathcal{R}_{\phi}$ is called the response and describes the (noiseless) measurement process translating the signal into data space. In this case, the evaluation of the response at the location of source $i$ is sufficiently described via

[TABLE]

with $\widehat{e}$ being the radial unit vector on the 2-sphere. The quantity $n_{\phi}$ in Eq. (4) is the measurement noise, which in the easiest case can be assumed to be independent and Gaussian distributed, according to

[TABLE]

with zero mean and known noise covariance $N_{\phi}$. The noise will play an important role for our estimation of non-galactic contributions to $\phi$, which will be discussed in more detail in Sec. 2.3.2. For the moment we assume $N_{\phi}$ to be simply given by the observational measurement uncertainties, therefore being a diagonal matrix in data space.

As we are conducting a Bayesian inference, our main task is the evaluation of the posterior probability distribution. An integral part of this distribution is the likelihood, which in our case is best described by the following Gaussian:

[TABLE]

The likelihood will be further determined by the the specific modeling of $\phi$, which we will present in Sec. 4. There we will also specify the different prior terms.

The data set used in this work is exactly the same compilation of Faraday rotation measures as used in Oppermann et al. (2012) and is mostly publicly accessible 222https://wwwmpa.mpa-garching.mpg.de/ift/faraday/. This compilation contains rotation measures, corresponding standard deviations and source positions for 41632 sources 333(Bonafede et al. 2010; Broten et al. 1988; Brown et al. 2003; Taylor et al. 2003; Brown et al. 2007; Haverkorn et al. 2006; McClure-Griffiths et al. 2005; Clarke et al. 2001; Clarke 2004; Clegg et al. 1992; Feain et al. 2009, 2011; Gaensler et al. 2001; Hammond 2011; Heald et al. 2009; Braun et al. 2007; Hennessy et al. 1989; Johnston-Hollitt 2011, 2003; Johnston-Hollitt & Ekers 2004; Kato et al. 1987; Kim et al. 1991; Klein et al. 2003; Gregorini et al. 1998; Vigotti et al. 1999; Lawler & Dennison 1982; Dennison 1979; Mao et al. 2010, 2012; Gaensler et al. 2005; Mao et al. 2008; Minter & Spangler 1996; Oren & Wolfe 1995; O’Sullivan 2011; Roy et al. 2005; Rudnick & Jones 1983; Schnitzeler 2011; Simard-Normandin et al. 1981; Tabara & Inoue 1980; Taylor et al. 2009; Condon et al. 1998; Van Eck et al. 2011; Wrobel 1993). Similar to Oppermann et al. (2015) we multiply the uncertainties of the Taylor et al. (sub-)data set (Taylor et al. 2009) with 1.22 according to Stil et al. (2011). No additional polarimetric information such as e.g. Faraday spectra was used in this work. The data and the corresponding standard deviations are shown in Figs. 3(a) and 3(b).

Next, we specify the modeling of $\phi$ in Eq. (4). We will furthermore explain the treatment of the non-galactic components of $d_{\phi}$.

2.3 The model

2.3.1 The sky model

The results of the previous reconstruction (Oppermann et al. 2015) indicates that the Faraday sky has two characteristics that necessarily need to be modeled for an optimal inference. Firstly, the sign of the sky field changes rapidly, which reflects the rather abrupt directional changes of the LOS magnetic field component. Secondly, the absolute value of the signal varies over two orders of magnitude, most notably the Galactic disk will have a much stronger signal than the Galactic poles. For these reasons, we chose to parameterize $\phi:S^{2}\rightarrow\mathbb{R}$ as a point-wise product between two fields living on the sky, or explicitly

[TABLE]

where $\chi:S^{2}\rightarrow\mathbb{R}$ will be called the sign field, as it is able to capture the information on the sign of the Faraday depth, and $\rho:S^{2}\rightarrow\mathbb{R}$ will be called the amplitude field, as $e^{\rho}$ is strictly positive after the exponentiation and is able to model the large amplitude variations of the Faraday sky over orders of magnitude. Both fields are unitless. This is a more generic approach as in Oppermann et al. (2012), where instead of the amplitude field a profile function was used to capture the latitudinal dependence of the overall Galactic Faraday dispersion profile. Note that the sign field is in no way constrained to only contain information on the sign of the signal, but also can capture morphological features, leading to a degeneracy between the two fields. This could be broken by either imposing constraints on the correlation structure for (one of) the fields or by introducing new data that informs the algorithm on one of the two fields. However, it turns out that the symmetry is sufficiently broken by the Gaussian process priors we impose on both fields separately in combination with the specific functional form of $\phi$ as given by Eq. (8). Sign variations can only be captured by the sign field, ensuring that it is structured. The overall amplitudes of these fluctuations, which change as a function of position on the sky, are, however, more easily represented by the amplitude field, which therefore preferentially absorbs those. Both fields are assumed to have independent and Gaussian isotropic statistics. Their joint prior probability function is

[TABLE]

However, in contrast to the distribution in Eq. (6), the covariances $S_{\chi}$ and $S_{\rho}$ are unknown. If this were not the case, the above model in a Bayesian setting would directly lead to a specific application of the well known (non-linear) Wiener filter (Wiener 1966; Enßlin et al. 2009). This is not the case unfortunately, as in our case we have little a priori knowledge on the correlation structure. Similar inference problems have been solved in the past, e.g. in the previous reconstructions of the Faraday sky by (Oppermann et al. 2012, 2015), but also e.g. in (Selig et al. 2015; Pumpe et al. 2018; Leike & Enßlin 2019).

All of these inferences relied on a framework laid out by information field theory (IFT) (Enßlin et al. 2009; Enßlin 2018). IFT connects Bayesian statistics with methods from statistical and quantum field theory, joining them into a inference scheme that connects noisy, incomplete data with the underlying continuous field(s). Considering the problem mentioned above, the requirement of simultaneously inferring the map and the correlation structure of a field leads to the critical filter formalism, which was first formulated in IFT by Enßlin & Weig (2010) and Enßlin & Frommert (2011). We will follow the most modern formulation, as outlined in Knollmüller et al. (2017). We will shortly summarize the modeling of the respective covariance matrices in the Appendix.

2.3.2 The noise model

As mentioned before, we also need to find an estimate for the non-galactic contributions to the Faraday rotation. In this part of the inference, we follow the approach in Oppermann et al. (2012) very closely. We will start by inserting Eq. (3) into the measurement equation (4) and define an effective noise term $\tilde{n}_{\phi}$:

[TABLE]

We have no reason to drop the assumption of Gaussianity for the new noise term $\tilde{n}_{\phi}$, but have to adapt the covariance to include the increase in uncertainty. As we have no a priori guesses on the specific systematics of the different sources, we will infer this increase in uncertainty as well. This is implemented by introducing the following model

[TABLE]

for the new noise covariance $\widetilde{N}$ where $\eta_{\phi}$ are the parameters that need to be inferred and $\sigma^{2}_{\phi}$ is the reported measurement uncertainty. We will assume the $\eta_{\phi}$ to follow an inverse gamma distribution:

[TABLE]

The hyper-parameters $\alpha_{\phi}$ and $\beta_{\phi}$ need to be specified for the inference. In contrast to Oppermann et al. (2012) we do not set $\alpha_{\phi}=1$ and $\beta_{\phi}$ such that the prior expectation value of $\ln(\eta_{\phi})$ is zero, but demand a steeper slope with $\alpha_{\phi}=2$ and the prior expectation value of $\eta_{\phi}$ to be unity. The reason for this is that we are allowing for more degrees of freedom in our new sky model as compared to Oppermann et al. (2012), as we here have a two dimenensional amplitude field, where as Oppermann et al. (2012) only had a one dimensional profile function that was empirically determined without any prior. The field prior we have to impose on this amplitude field works against the data, increasing the tendency to classify data structures as noise. The larger value of $\alpha_{\phi}$ just counter balances this tendency to a degree that our resulting initially revised Faraday map (Fig. 1(b)) looks similar to that of Oppermann et al. (2015) in Fig. 1(c).

2.4 The inference

The hierarchical Bayesian model described in the previous sections can be symbolically depicted as a hierarchical tree, of which the lower branches are shown in Fig. 4. The higher branches of the tree contain the hyper-parameters and are discussed in the appendix together with the more explicit likelihood and prior functions.

The evaluation of the posterior is a non-trivial task in our case, as the non-linearities in the models lead to highly non-Gaussian distributions for which no analytical description is known to exist. Furthermore, the fact that we aim for an simultaneous inference of the fields determining the Faraday map as well as their correlation structures leads to a strong interdependence between all degrees of freedoms. Together with the high dimensionality of the problem with of the order of $10^{6}$ unknown parameters for the desired resolution of $\approx 10$ arcmin this leads to considerable numerical complexity.

IFT solves such problems by providing efficient and well tested algorithms for posterior evaluation. In contrast to Oppermann et al. (2012), we employ here the Metric Gaussian Variational Bayes method (MGVI), which is better suited for simultaneous map and correlation structure inference compared to the maximum a posteriori algorithm used in Oppermann et al. (2012), see Knollmüller et al. (2017) for details. We furthermore implement the aforementioned hierarchical tree as a forward model, which aims among other things for a more efficient numerical coupling between different parts of the model and allows the employment of sampling methods (Knollmüller et al. 2017).

On a technical side, the inference was conducted using the newest version (version 5) of the NIFTy python package (Selig et al. 2013; Steininger et al. 2017; The NIFTy5 team et al. 2019) for numerical information field theory. NIFTy is specifically designed for IFT applications and provides libraries for implementing the above mentioned models as well as for inferring their parameters and the uncertainties of those.

2.5 The results

2.5.1 Comparison to older results

The initially revised Faraday depth map relying on rotation measure (RM) data only and inferred according to the model depicted in Fig. 4 is shown in Fig. 1(b), the previous one of Oppermann et al. (2015) in Fig. 1(c), with corresponding standard deviations in Figs. 5(b) and 5(a), respectively. The difference of both maps are shown in Fig. 2(b). A visual comparison of the two maps reveals satisfying agreement. The difference map, however, shows that there are indeed significant small scaled deviations of the maps in the area of the Galactic disk. These may be a result of the different noise estimation parameters in the updated inference and/or the different treatment of the amplitude field (see the discussion in Sec. 2.3.2). One should note that the revised reconstruction was conducted with a doubled resolution, which may also be responsible for some of the very small scaled differences.

The uncertainty maps nicely demonstrate the impact of the sky modeling on the inference. Fig. 5(a) shows that the uncertainty of the previous reconstruction, which exhibits imprints of the longitudinal galactic profile used in the amplitude modeling. Our updated uncertainty map in Fig. 5(b) shows a decrease in uncertainty in certain regions, especially in the galactic disk. This simply reflects the fact that our updated inference scheme can use the correlation information encoded in the data much more efficiently.

A further quantity which needs to be inspected are the power spectra depicted in Fig. 6, as they have been inferred with different algorithms in the previous and initially revised approaches. We note that we can confirm the slope of the spectrum. The revised map has significantly less power on small scales. Again, the parameter change in the noise estimation may be the reason for this discrepancy, as we know that it can have a considerable impact on the smaller scales, as we change the ability of the algorithm to increase uncertainties to accommodate outliers in the data.

With respect to the error estimation, we show in Fig. 7(a) the comparison of the reported observational noise standard deviation $\sigma$ and our estimate for the noise $\widetilde{\sigma}$ according to Eq. (10). The uncertainties of some sources have been multiplied by factors of up to $1000$, indicating strong internal or extra-galactic contributions to these particular RM’s. Nonetheless, we confirm or even decrease the uncertainties for most measurements, as shown in Fig. 7(b), where we show a normalized density plot in the $\sigma-\widetilde{\sigma}$ plane.

With this we conclude the comparison between the previous and the initially revised results and turn our face to the analysis of the underlying components of our sky model defined in Eq. (8).

2.5.2 The sign and amplitude fields

The component fields $\chi$ and $\rho$ are shown in Figs. 8(b) and 9(b), respectively. The sign field captured most of the small scale structure, but has more or less completely lost all information on the Galactic disc profile, which was absorbed by the amplitude field, as intended by our modeling. $\chi$ shows some similarity to the corresponding signal field inferred by Oppermann et al. (2015), which is shown in Fig. 8(a). The exponentiated amplitude field appears to be relatively smooth with few distinguished features apart from the large scaled disk profile. These features, however, show remarkable similarity to galactic free-free emission maps, which was precisely measured by missions investigating the Cosmic Microwave Background (CMB) such as Planck (Planck Collaboration 2016a), as it is a important foreground component in the microwave sky. We show the Planck free-free EM map in Fig. 9(a). The color scales were chosen to highlight the structures in the maps of which we think that they originate from the same galactic structures. This apparent resemblances in both maps have motivated us in including the free-free map in our inference. A investigation of the physical plausibility of this correlation is given in the next section.

3 Including the Free-Free emission

3.1 The physics

In addition to the visual similarities between the free-free map and our inferred amplitude field, we want to further motivate this apparent congruence (and deviations from it) by investigating the physics of the free-free emission and the Faraday rotation a bit more closely. The $\mathrm{EM}$ traced by the bremsstrahlung emitted by thermal electrons on protons in the warm interstellar medium can be quantified via (see e.g. Smoot (1998))

[TABLE]

where we defined $\mathrm{em}_{\mathrm{ff}}\equiv n^{2}_{\mathrm{th}}$ as a source term for $\mathrm{EM}_{\mathrm{ff}}$.

If we turn to the Faraday depth and want to study Eq. (2) with the intention of an comparison to $\mathrm{EM}_{\mathrm{ff}}$, we must quantify the dependence of $B_{\mathrm{LOS}}$ on the thermal electron density. In general, we can assume the magnetic flux to be frozen into the ISM (de Gouveia Dal Pino 2006). The exact nature of the $B-n_{\mathrm{th}}$ relationship, however, is complicated as it strongly depends on the morphology and dynamics of the plasma under consideration. In the simplest case of an isotropically collapsing structure, flux freezing leads to an $|B|\propto n_{\mathrm{th}}^{\frac{2}{3}}$ dependence (Mestel 1966). For a realistic description the Galactic interstellar medium (ISM), however, one needs to consider non-linear magneto-hydrodynamics (MHD). Under certain simplifying assumptions as shown in Passot & Vázquez-Semadeni (2003) one may decompose the ISM flow into several wave types corresponding to different modes, leading to different density dependencies of the magnetic field strength. Most notably in the low density regime, the dominating mode approaches a constant in $n_{\mathrm{th}}$, implying no correlation between $|B|$ and $n_{\mathrm{th}}$ in this regime. Furthermore, as waves corresponding to different modes may very well travel trough the same regions within a short time span, the $B-n_{\mathrm{th}}$ correlation may be subject to large variability in all regimes of measured electron density (Vázquez-Semadeni 2015). We therefore assume the dependence of the magnetic field strength on $n_{\mathrm{th}}$ to follow a power law with an unknown and spatially dependent coefficient $|B|\propto n_{\mathrm{th}}^{p_{B}}$. We can then rewrite the proportionality of the absolute value of integrand of the Faraday depth in Eq. (2) on $n_{\mathrm{th}}$

[TABLE]

Here we introduced the angle $\theta$ indicating the orientation of the magnetic field with respect to the LOS. Its impact will be discussed later on. If $p_{B}\approx 1$ for all locations along a line of sight, the absolute value of the Faraday integrand does strongly depend on the free-free emissivity and the additional electron density dependence has completely canceled. In the unrealistic case that $p_{B}=-1$, the absolute value of the Faraday integrand would be completely uncorrelated to the free-free source term. Furthermore as discussed before, in some cases the electron density will not be a good tracer for the magnetic field strength and the LOS projection might partly mask the relation in the observables, which again limits the above relationship. Nonetheless, we think that the above considerations motivate the inclusion of the free-free map as a proxy for the Faraday amplitude field in our inference, as for most realistic cases of the power law index $p_{B}$ the integral over both sides of Eq. (3.1) will lead to a correlation between the free-free and the Faraday sky. On the other hand this analysis also indicates that the effects of ISM dynamics on the magnetic field strength have to be considered for a reliable inference.

A second effect that needs to be taken into account are sign reversals of the magnetic field along the line of sight, which could strongly affect the value of $\phi$ in contrast to $\mathrm{EM}_{\mathrm{ff}}$. To consider these, we approximate the line of sight integral via

[TABLE]

where $L$ is the length of the line of sight and $\langle\ldots\rangle_{L}\equiv\frac{1}{L}\int_{\mathrm{LOS}}dl\ldots$ indicates an averaging process along the line of sight. The last approximation implies statistical independence between the magnetic field strength and $\theta$. This shows that to good approximation, the geometry of the magnetic field can lead to a multiplicative term, which is not constrained by the free-free data alone and has to be considered independently. In our setup, it can only be further determined by discrepancies between the Faraday and free-free data sets.

Similar to the first part of this paper, we will first discuss properties of the EM data and then show the modeling of the respective sky maps in Sec. 3.3.

3.2 The data

In this work, we use the free-free EM map of the Planck 2015 release (Planck Collaboration 2016a). This map is a side product of a foreground subtraction process of the Planck team in their pursuit of the CMB sky reconstruction. The free-free sky is a dominant component of this foreground in the low frequency regime, along side with synchrotron and spinning dust emission (Planck Collaboration 2016b). Most components in the microwave sky are hardly distinguishable, as they result from similar structures in the sky. Planck therefore employs a sophisticated component separation algorithm called $\mathtt{COMMANDER}$ (Eriksen et al. 2008), mostly exploiting the different energy spectra of the components. This algorithm uses the Gibbs sampling technique to approximate the posterior distribution for all components, also making use of the WMAP data (Bennett et al. 2013) and a 408 MHz survey map (Remazeilles et al. 2015) for the separation. An updated version of $\mathtt{COMMANDER}$ employed in the 2018 release of Planck takes the correlation structure of the components into account. However, Planck has not released an updated free-free map in 2018.

As a side product of the estimation of the posterior via sampling, Planck provides us not only with the posterior mean of the free-free component, but also with uncertainties, which are a valuable source of information and lets us treat the free-free data as a slightly noisy measurement of the true free-free $\mathrm{EM}$ with known uncertainties. Unfortunately, upon inspecting the Planck data, we found some pixels to have zero variances. As this is rather surprising given the large degeneracy between the components, we expect this to be a numerical artifact of the sampling procedure or of the marginalization process that results in the uncertainties for specific components. Our algorithm relies on accurate uncertainty measures for all pixels, as in case of a severe underestimation of the errors in some specific area, this region is over-weighted by the algorithm in the inference of the correlation structure, thereby immediately affecting other parts of the sky and subsequently all other maps of the inference. This has to be accounted for if one is not completely sure that the noise values are free of systematics. We will therefore perform a similar noise estimation procedure as for the Faraday data.

With this in mind, we write the measurement equation for the free-free data as

[TABLE]

In case of Planck, the response $\mathcal{R}_{\mathrm{ff}}$ is just a unit matrix, as the data is already provided as a full sky map on the resolution we desire. Here, we decompose the noise $n_{\mathrm{ff}}$ into two components namely a observed one obtained by COMMANDER and an unobserved one, which contains systematics. As we do not consider the correlation structure of the noise in for any data set, the procedure of noise estimation is the same from here on as in 2.3.2, apart from different hyper priors. In this case, we choose $\alpha_{\mathrm{ff}}=1$ and $\beta_{\mathrm{ff}}$ such that the expectation value over $\ln(\eta_{\mathrm{ff}})$ is zero. This is Jeffreys prior for the inverse gamma distribution (Oppermann et al. 2012). The lower branches of the free free sky model are shown in Fig. 12. The full extended model for the hierarchical Bayesian model is shown in Fig. 13. Using the considerations above, the likelihood for the free-free sky is again Gaussian in analogy to Eq. (7)

[TABLE]

with the noise covariance $\widetilde{N}_{\mathrm{ff}}$ defined in analogy to Eq. (10). Again, we still need to explain the detailed model for $\mathrm{EM}_{\mathrm{ff}}$ and its priors. This is presented in the next section.

3.3 The free free model

3.3.1 The sky model

We now turn our head to the modeling of the free-free sky and its connection to the amplitude field of the Faraday sky under the constraint of the physical limitations considered in Sec. 3.1. We model $\mathrm{EM}_{\mathrm{ff}}$ with a log-normal approach, again to enforce positivity and to capture the expected large variability of $\mathrm{EM}_{\mathrm{ff}}$ on the sky:

[TABLE]

If the true free-free sky would be a perfect estimator for the amplitude field of the Faraday sky, we could simply replace $\rho$ in Eq. (8) with $\epsilon$ to simultaneously infer both quantities from both data sets. However, as discussed in Section 3.1 a necessary condition for a combination of the information of free-free EM’s and Faraday rotation measures is to account for magnetic field reversals, different scalings with the thermal electron density and potential systematics in the data sets. All issues can be modeled by introducing further terms to the model.

We start with our approach to model the differences in the electron density scaling and the geometry factor introduced by the sign reversals. As we can only work with the LOS integrated quantities, the most natural way of inferring the three dimensional effective scaling factor $\frac{1+p_{B}}{2}$ in Eq. (3.1) is not possible, as this would necessarily imply setting up the full 3D problem, for which we are lacking the necessary radial information. We therefore choose to introduce two scaling fields $\gamma:S^{2}\rightarrow\mathbb{R}$ and $\psi:S^{2}\rightarrow\mathbb{R}$ into the modeling according to Eq. (19) which essentially sets up a linear fit for the logarithm of the Faraday amplitude field in terms of the reconstructed logarithmic free-free sky:

[TABLE]

This model somewhat mimics the approximation derived in Eq. (3.1), with $\gamma$ approximating the role of the power law coefficients and vice versa $e^{\psi}$ the geometry factor. Unfortunately, we loose the interpretation of the $\gamma$ field in terms of $\frac{1+p_{B}}{2}$ in Eq.(3.1) due to the aforementioned lack of consistent 3D modeling. Furthermore the large scale modes of the $\psi$ and $\gamma$ fields are to a degree degenerated, making a quantitative interpretation of the $\psi$ field very hard, although the smaller scaled morphology may be informative. Altogether, these terms should be able to capture sign reversal effects as well as strong deviations from the proportionality of the magnetic field strength on the thermal electron density presented in Eq. (3.1).

We further have to consider the case where the Faraday sky contains amplitude structure not present in the $\mathrm{EM}_{\mathrm{ff}}$ sky. This is done by expanding the model according to

[TABLE]

Here, the second amplitude field $\delta:S^{2}\rightarrow\mathbb{R}$ is not connected to the free-free map, giving it the possibility to capture structures that are in the Faraday data but not in the free-free data. All newly introduced fields are assumed to follow Gaussian statistics, so their prior is similar to Eq. (9), again with unknown correlation structure. One should note at this point that the obvious a priori degeneracies between the newly introduced helper fields are usually broken during the inference, if sufficiently supported by data. This is possible due to the joint correlation structure inference, which constrains the respective fields just to the degree it is demanded by the inference problem.

Although there is motivation for each term, the above modeling might be viewed as ad hoc, as one may certainly argue for the addition or omission of certain terms or a completely different modeling altogether. There are certain limitations for potential models, as they all have to be able to capture the physical characteristics (say e.g. the sign reversal of the Faraday map), and should be kept as simple as possible. Unfortunately these requirements already allow for a range of different parametrizations. We have varied the model during our analysis by e.g. adding or omitting further terms in Eq. (20). We found little variation in the overall morphology of the Faraday sky, apart form the disk, where e.g omission of the $\gamma$ and $\psi$ field lead to an extremely pronounced features, with RM values reaching regimes a magnitude higher than previous reconstructions. As we think that the above presented physical considerations are too important to be neglected, we deem the resulting maps in these cases as unrealistic. Given that all newly introduced helper fields have captured distinct structures as we will show in Sec. 3.4, we view the model in Eq. (20) as one of the simplest cases which can capture all relevant effects, but nonetheless we deem the systematic errors higher than the statistical ones in this inference.

To gain some insights on the robustness of the obtained free-free map, we will further perform an inference just for the free free sky with the model in Eq. (18), unconstrained by Faraday data. The comparison of the result to the outcome of the joint Faraday and free-free inference will be presented in the results section. We can implement this extended model numerically with similar algorithms as described in Sec. 2.4.

3.4 The Results

3.4.1 The revised Faraday map

In the following we show the results of the model depicted in Fig. 13. We show the mean, the differences to the previous and initially revised reconstructions and the power spectrum in Figs. 1(a), 2(a), 2(c) and 6, respectively. The inclusion of the free-free sky has lead to a notable difference in our newly revised estimate of the Faraday sky. Driven by the strong disc in the free-free data, the Faraday sky now has a similarly strongly pronounced disc feature as well. Other features mostly farther away from the Galactic plane remained relatively stable, such as the northern arc of the Gum nebula, which is also strongly visible in the free-free data. In Fig. 5(c), we show the updated uncertainties. These have again narrowed considerably compared to our initial free-free data free reconstruction as shown in Fig. 5(b). The innermost part of the disc, however, is still quite uncertain. Furthermore one should be aware that any model uncertainties are not considered in this plot. The new power spectrum of the Faraday map is shown as the red line in Fig. 6. It is very similar to that of the previous reconstruction, with an notable offset towards smaller scales, implying a slightly steeper power law.

3.4.2 The components

We now show some of the components of the Faraday model in Eq. (20), apart from the $\epsilon$ field in Fig. 10(a), which will be more closely discussed in the following section dedicated to the inference results of the free-free sky.

The revised Faraday amplitude is shown in Fig. 9(d). The comparison to the initially revised amplitude field in Fig. 9(b) and the free-free data in Fig. 9(a) demonstrated the influence of both the Faraday and the free-free data on the new field. The field is enhanced by roughly a factor of 4 compared to the old amplitude field. This may be a result of the symmetry breaking process described in Sec. 2.3, as now the amplitude field is partially constrained by the free-free sky. The comparison of the sign fields will demonstrate that this factor was mostly absorbed there.

The second amplitude field $\delta$ shown in Fig. 10(b), its influence on the total logarithmic Faraday amplitude is demonstrated by Figs. 10(c) and 10(d). The field is dominated by two diffuse structures of enhanced $\delta$, centered in the Galactic plane and separated by roughly [math] longitude. Remembering the discussion in Sec. 3.1, the $\delta$ field is supposed predominately capture amplitude variations between the free-free and the Faraday sky. Larger values of $\delta$ can be caused by systematic errors in either the free-free or the Faraday data and/or a high alignment of the magnetic field with the LOS in this direction. The approximated longitudinal difference of the two $\delta$ enhancements of about [math] as well as their positions seem to correspond to the angular positions of the Orion arm on the sky, which are depicted in the plot as red crosses. In the first quadrant of the Milky Way, at about [math] longitude, this minor arm is known to extend over several kiloparsec up to the Sagittarius arm, with which it might merge at last (Xu et al. 2009). In the other direction towards the third quadrant of the Galaxy, the structure of the arm is much more complex. Even a bifurcation or a crossing with the Perseus arm seems possible (Vázquez et al. 2008). We have indicated possible continuations of the arm in the plot at [math] and [math] longitude, respectively. The $\delta$ field seems to reflect not only the position of the arm, but also its morphology, as the structure in the third quadrant is more diffuse, as opposed to the relatively pronounced structure at [math] longitude. All in all, this would argue for a high magnetic field alignment as a viable possibility for the interpretation of the second amplitude field. As the Sun is positioned directly in the Orion arm, we expect more coherent magnetic fields along the LOS directed pointed directly to the dominant direction of the arm. Remembering the discussion in Sec. 3.1, this then may lead to a disproportionate increase in the observed Faraday measures as the free-free map would indicate at these positions. Furthermore the Faraday sky in Fig. 1(a) indeed shows opposite magnetic signs for those regions, as one would expect for a field traversing the solar location. We show the field

[TABLE]

in Fig. 11(a), which corresponds to our estimate of the enhancement of the Orion Faraday depth, should the above considerations be true and if the $\chi$ field really only captured the sign. The (in terms of amplitude contributions) complementary Faraday component is shown in Fig. 11(b). The considerable differences to Fig. 1(a) in the sky regions corresponding to the Orion arm reveal again the influence of the $\delta$ field there. Of course we can not rule out any of the other explanations for the observed $\delta$ enhancements, the plot in Fig. 11(a) should therefore be taken with care.

The comparison of the amplitude field in Fig. 10(c) without the $\delta$ contributions with the $\epsilon$ field alone (see Fig. 10(a)) reveal the strong rescaling and reweighting of the free-free contribution due to the $\gamma$ and $\psi$ fields.

At last we discuss the sign field $\chi$ in Fig. 8(c). Little has changed in this component morphologically compared to Fig. 8(b). The absolute value of the field has decreased by approximately a factor of 4. This factor was mostly absorbed into the new amplitude field, as discussed before. Apart from that, the field has again fully served its intended purpose, namely capturing the sign of the Faraday sky.

We now discuss the results of the inference of the free-free EM sky.

3.4.3 The revised free-free map

The revised free-free EM sky constrained by both the Planck and the Faraday data is shown in Fig. 9(c), the logarithm of it is shown in Fig. 10(a). The difference to the map inferred using only Planck data is shown in Fig. 14. There is little difference between the maps, indicating that the systematic differences of the data sets were mostly absorbed by the second amplitude field $\delta$ and the rescaling fields $\gamma$ and $\epsilon$. Notable differences appear only towards the Galactic poles, where Planck has measured little to no free-free component, while considerable Faraday rotation exists. The visual comparison of Fig. 9(c) to the data shown in Fig. 9(a) demonstrates good accordance. The main difference is the lack of small structures far away from the disc, which where mostly absorbed into the noise. Closer inspection reveals some small scaled deviations in the disc, which we want to investigate by showing an excerpt of the disc ranging from [math] to [math] longitude in Fig. 15. This shows approximately the same region as the first plot of Fig. 7 in Planck Collaboration (2016b), where the Planck team compares their free-free temperature map with independent measurements via the radio recombination line survey (RRL). While our reconstruction reproduced the overall morphology, some peaks in the picture experience deviations around $20\,\%$. In general our results seem to underestimate the amplitudes compared to Planck. It is noteworthy that the comparison between Planck and RRL resulted in a similar statement, as the Planck features were generally much more pronounced. We will refrain from drawing physical conclusions from these discrepancies as these can only be resolved by a detailed discussion of the Planck data analysis, which is beyond the scope of this work. Note again that the Faraday data has little influence on this result. The changes are mainly driven by the noise estimation and the full consideration of the correlation structure of the free-free EM sky. We think that this deems the inclusion of the correlation structure in the inference an important step for any component separation algorithm, as also shown in Knollmüller & Enßlin (2017).

4 Summary and conclusion

We infer the Galactic Faraday sky map using the same data as Oppermann et al. (2012) using an updated algorithm. This involves the separation of the Faraday sky into the point wise product of a sign and amplitude field, which both are supposed to mainly capture the sign and the Galactic profile of the Faraday sky, respectively. In the first part of the paper, we can largely confirm the old results in map and power spectrum. Smaller differences can be explained by the algorithmic advancements. The biggest improvement here is the updated uncertainty map, where we were able to significantly reduce the error bars over the whole sky, especially in some locations in the galactic disc.

In the further proceedings we notice that the amplitude field has notable similarities to the Galactic free-free EM map. In the second part of the paper we therefore incorporate free free data as a potential proxy for the Faraday amplitude field into our inference. This leads to a new Faraday map which now is notably different to previous results. The Galactic disc is much more pronounced and uncertainties are largely reduced. Since a pronounced Galactic is very sensible, as well as required by the data, which now is more satisfied as in any of the other approaches, we regard this map as the most reliable and therefore to be the main result of this work. Moreover, the component fields introduced in our model to resolve discrepancies between the two data sets reveal indicators for more LOS aligned magnetic field structures in the direction of the Orion arm. We also produce a denoised free-free map, which agrees with the Planck data to large extent.

Acknowledgements.

We acknowledge fruitful discussions with J. Knollmüller and others from the information field theory group at the MPI for Astrophysics, Garching.

Appendix A Generic IFT concepts

These appendices are intended to provide the details on the inference scheme used in this work. We will start with introducing some core concepts of IFT. After that, the details of the correlation structure inference are explained. In this part, we will also show the equations for the full problem posterior of our model above. We close with some details on the used posterior evaluation scheme. Extended and more mathematically rigorous discussions on these topics can be found in Enßlin et al. (2009); Knollmüller et al. (2017); Enßlin (2018); Knollmüller & Enßlin (2019).

A generic expression for the posterior probability distribution for a inference problem can be written using Bayes theorem:

[TABLE]

In our case, the signal $s$ contains all fields that describe the aforementioned properties of the desired sky maps. The left hand sight of Eq. (22) is commonly called the posterior. The right hand side consists of the prior $\mathcal{P}\left(s\right)$, the likelihood $\mathcal{P}\left(d|s\right)$ and the evidence $\mathcal{P}\left(d\right)$. The latter is important for model comparison tests, but negligible in an inference algorithm that aims to find the optimal field realizations of $s$ for a given model, as it does not depend on $s$, but only serves as a normalization constant.

In the further proceeding, we will consider the information Hamiltonian or negative log-probabilities of the above equation,

[TABLE]

This has mostly convenience reasons, as the resulting sums of log-probabilities are usually numerically better behaved as the product of probabilities. Furthermore there is a nice analogy to statistical mechanics and field theory (hence the use of ‘Hamiltonian’), which results in a rich variety of techniques IFT can make use of.

Already using the expressions for the likelihoods from Eqs. (7) and (17), the respective Gaussian priors similar to Eq. (9) and the terms stemming from the noise estimation models in Eq. (11), we can write down the posterior Hamiltonian assuming independent measurements of $\mathbf{d}=\left(d_{\phi},d_{\mathrm{ff}}\right)^{T}$:

[TABLE]

The vector $\mathbf{s}\equiv\left(...\right)^{T}$ contains all Gaussian random fields introduced in the previous sections, which are assumed to be a priori independent from each other. The symbol $\mathrel{\stackon[1.5pt]{=}{\stretchto{\scalerel*[width("=")]{\wedge}{\rule{4.30554pt}{12.91663pt}}}{0.5ex}}}$ indicates the omission of unnecessary additive constants. This equation already contains most of the model. A crucial part is missing, however, namely we have not specified the prior covariances $S_{f}$. The modeling and inference of those is explained in the next section.

Appendix B Correlation structure model

Again we will not follow the details in Oppermann et al. (2012), but use an updated algorithm for the correlation structure inference. A short discussion on the differences of the two approaches and potential implications on the spectrum follows at the end of this section. We have verified our results via a comparison to the previous map only using the data used in Oppermann et al. (2012) in Sec. 2.5.

The assumption of statistical homogeneity and isotropy together with Gaussianity for a field $f$ implies via the Wiener-Khintchin theorem that the respective covariance is diagonal in harmonic space. We will make use of this theorem, as it allows us to describe the correlation structure via a field with at most $N_{\mathrm{pix}}$ dimensions, as opposed to the $N^{2}_{\mathrm{pix}}$ dimensions in the non-harmonic space. The correlation structure of a Gaussian field $f$ is then fully described via its power spectrum defined as

[TABLE]

where $\mathbb{P}_{\ell}$ is a projection operator onto spectral bands $\mathbf{\ell}$ of the power spectrum. This operator simply gives us the possibility to reduce the degrees of freedom of the spectrum, it may very well be set to a unit matrix in case one e.g. expects rather sharp features such as characteristic frequencies.

For the inference, we follow the reparametrization trick laid out in Knollmüller et al. (2017). We decompose the sky fields into

[TABLE]

$\mathbb{Y}$ is the spherical harmonic transform. The amplitude fields $A_{f}$ contain the information on the correlation structure on the respective field, while the $\xi_{f}$ are a priori white Gaussian fields with zero mean and unit covariance which contain the information of the specific field configuration of the field $f$ given the correlation structure in $A_{f}$. Both fields live in harmonic space. One can show that the amplitude fields are related to the covariance via

[TABLE]

This formulation has immediate consequences for the Hamiltonian, as the prior terms are now independent of the correlation structure as

[TABLE]

using the above equations. For a more general discussion of the advantages of this approach we refer the reader to Knollmüller et al. (2017). Now we are left with the inference of the amplitude fields. Those are modeled via a log-normal field $\tau_{f}$ to ensure positivity:

[TABLE]

The inference then follows the approach of Arras et al. (2019), in which $\tau_{f}$ is again transformed according to the reparametrisation trick above in to a white field $\xi_{f}$ and hyper parameters $\sigma_{f}$ and $l_{f}$ which steer the expected variability in the spectrum.

In Oppermann et al. (2012) the power spectrum was modeled via a log-normal field together with the product of an inverse gamma and a regularizing smoothness likelihood. The above demonstrated flattening of the hierarchical structure was not performed.

We show a comparison in the section 2.5 to demonstrate that the revised algorithm reproduces the inference of the correlation structure in Oppermann et al. (2012), and also comment on the minor differences occurring between the two spectra. The full Hamiltonian for the simultaneous inference of Faraday sky, the free-free sky, the correlation structure of their respective component fields and noise corrections terms is then:

[TABLE]

where the last term contains the hyper priors for the $\tau_{f}$ inference and is described in detail in Arras et al. (2019).

Appendix C Evaluating the posterior

The Hamiltonian in Eq. (B) fully describes our information on the problem. Unfortunately, the non-linearities in this expression imply that this posterior distribution is non-Gaussian, making the evaluation a non-trivial task as a general analytic solution is not known to exist. Therefore, we need to choose an appropriate numerical evaluation scheme. For this end, variational Bayes approaches have proven to combine satisfying accurateness in their uncertainty estimates with speed and stability during the minimization (Knollmüller et al. 2017). In a mathematical language, this implies that we will minimize the Kullback-Leibler Divergence (KL)

[TABLE]

where $\widetilde{\mathcal{P}}$ is a Gaussian probability distribution and $\widetilde{\mathcal{H}}$ is the respective Hamiltonian. $\left\langle...\right\rangle_{\widetilde{\mathcal{P}}}$ indicates an averaging procedure of the quantity in the brackets with respect to $\widetilde{\mathcal{P}}$. Minimizing the above expression can be interpreted as a fit of $\widetilde{\mathcal{P}}$ to $\mathcal{P}$. This approximated posterior is then easily evaluated via sampling due to its Gaussian nature.

Both the Hamiltonian and the KL functional can be conveniently implemented into NIFTy (Selig et al. 2013; Steininger et al. 2017; The NIFTy5 team et al. 2019), which in its newest version then takes care of calculating the correct gradients via auto-differentiation and provides a library of suitable minimization schemes for optimizing Eq. (C). Convergence is assessed via predefined convergence criteria and the stability is tested via starting the inference from different (both deliberately and randomly chosen) initial conditions. For further details on the here used MGVI technique see Knollmüller & Enßlin (2019).