PureB by Messenger Method
Michael Ray, Colin Bischoff, 2019-09-xx
This posting describes an alternate pure-B estimator applied to the 95 GHz DC 04.00 maps. It uses the "messenger method" (Elsner & Wandelt 2012) to Wiener filter the maps. We treat the E-mode part of the signal covariance as if it was noise covariance (Bunn & Wandelt 2017), so filter equation becomes:
mapwf = (SB/(SB+SE+N)) * map,
where SE,B are the E-mode and B-mode parts of the signal covariance.
For these simulations with unfiltered skies, the signal covariance is diagonal in the spherical harmonic basis. We treat the noise covariance as diagonal in the pixel basis (not strictly true because of the white + 1/ℓ noise spectrum in the data challenge maps). The messenger method iteratively solves fo the Wiener-filtered map by introducing a messenger field, t, with covariance matrix, T, that is proportional to the identity matrix, so it can be written easily in any basis. For each step of the algorithm, we transform back and forth between pixel space and spherical harmonic space, with the messenger field carrying information between the two bases.
As a first demonstration, we applied this technique to the 95 GHz DC 04.00 maps (Gaussian foregrounds, circular 3% fsky mask). The signal covariance matrices (in spherical harmonic space) are set to the standard lensed-LCDM E and B-mode spectra (purely diagonal). We did not include any foreground contribution to the signal covariance. The noise covariance is diagonal in pixel space, with values equal to 0.01 μK2 times the inverse of the mask (checked against a set of noise maps). For unobserved pixels, the noise covariance is set to a very large but finite value.
The messenger method includes parameter, λ, which scales the covariance of the messenger field. Convergence proceeds more quickly for large values of λ, but the final solution is obtained for λ = 1. This means that we want to implement a "cooling schedule" for efficient filtering. In our implementation, λ starts at a value of 1300, is decreased to 100 on the second iteration, and then is reduced by factor η = 0.825 on each subsequent iteration. When λ reaches 1, we iterate an additional 5 times. We have experimented fairly extensively with other cooling schedules and found that this method gives good convergence in a reasonable amount of time. (We are able to filter a single set of Q and U maps in 34 minutes on a single processor.)
Figures 1 and 2 show a map before and after filtering. The input maps is dominated by LCDM E modes. The output map contains B modes from foregrounds and lensing. The output map color scale has been zoomed in so that we can see that the pure-B filtering operation has extended the map into the unobserved region.
After applying the filter, we calculate the BB spectrum with anafast and bin with Δℓ=35, to match the pure-S2hat analysis. Noise bias is measured from 100 noise-only simulations. Bandpower window functions (Figures 4 and 5) and suppression factors are calculated by running ℓ-by-ℓ simulations (but only 6 realizations per ℓ value so far). Figure 3 shows noise-debiased and suppression-factor-corrected power spectra for 100 LCDM+foregrounds+noise simulations. The expectation value line is calculated by applying the bandpower window functions to a LCDM+foregrounds model that matches the sims.
There is some noticeable bias in the last bin especially. Something is clearly off with the calculations for this bin (see bandpower window functions) and needs to be investigated.
As one can see, the expectation value for our estimator runs directly through the middle of our simulations and is very close to our mean. It is worth noting that there is definitely something out of the ordinary happening in our last bin. In this bin, we routinely get results that do not agree with the rest of the data. Because of this, we generally ignore results from the last bin given that we are not too worried with what is happening at that ell range anyway. Shown below are plots of BB to BB and EE to BB band power window functions for both the messenger method estimator and the S2hat estimator. Band power window functions for the messenger estimator were calculated with an input signal C_ell = 1 on the sky. This means that the actual input power for deriving band power window functions is 1 * B_ell squared, where B_ell squared contains the effect of the beam. We also ran 6 simulations at each ell with this input power and then took the mean of the output across these six simulations to get the final band power window functions. Summing the window functions across ell values gives us the total suppression factor.
Figures 4 and 5 show bandpower window functions for input B modes → output BB
Again, it is visually apparent that there is something odd going on in the last bin for the messenger method estimator.
Below are plots the variance in each bin across the 100 realizations which were filtered. The exact same simulations were used as input to the two methods compared (signal plus noise simulations at 95 GHz using the 04.00 experiment configuration).
As shown above, the messenger method beats the S2hat estimator in terms of variance per bin. We believe the main reason for this is that the S2hat estimator uses a simple weighting of 1/N to filter the maps. Since the signal plus noise simulations have a sizeable amount of BB lensing signal contained within them, the fact that we are doing a weiner filter (and thus taking into account the signal covariance) means that we would expect to get a large improvement over an estimator which only knows about noise covariance. For delensed simulations, we would expect to see less of an improvement from the messenger method over the S2hat filtering.