# Evidence of cross-correlation between the CMB lensing and the -ray sky

###### Abstract

We report the measurement of the angular power spectrum of cross-correlation between the unresolved component of the Fermi-LAT -ray sky-maps and the CMB lensing potential map reconstructed by the Planck satellite. The matter distribution in the Universe determines the bending of light coming from the last scattering surface. At the same time, the matter density drives the growth history of astrophysical objects, including their capability at generating non-thermal phenomena, which in turn give rise to -ray emissions. The Planck lensing map provides information on the integrated distribution of matter, while the integrated history of -ray emitters is imprinted in the Fermi-LAT sky maps. We report here the first evidence of their correlation. We find that the multipole dependence of the cross-correlation measurement is in agreement with current models of the -ray luminosity function for AGN and star forming galaxies, with a statistical evidence of 3.0. Moreover, its amplitude can in general be matched only assuming that these extra-galactic emitters are also the bulk contribution of the measured isotopic -ray background (IGRB) intensity. This leaves little room for a big contribution from galactic sources to the IGRB measured by Fermi-LAT, pointing toward a direct evidence of the extragalactic origin of the IGRB.

98.80.-k, 98.80.Es, 98.70.Rz, 98.62.Sb

## 1. Introduction

The weak gravitational lensing by large scale structures imprints the integrated dark matter distribution onto the cosmic microwave background (CMB) anisotropies. It results in a remapping of the CMB observables, which depends on the line-of-sight integral of the gravitational potential, with a broad kernel peaking at a redshift , and which is referred to as the lensing potential (Blanchard & Schneider, 1987) (see also (Lewis & Challinor, 2006) for a review). This process perturbates the statistical properties of the CMB observables, which are primarily very close to Gaussian fields. This non-Gaussian signature can be exploited to extract the lensing potential from the CMB maps (Okamoto & Hu, 2003). Using such a technique, the Planck Collaboration obtained nearly all-sky maps of the lensing potential reconstructed from an undusted CMB temperature map (Ade et al., 2013) and from both foregorund-cleaned CMB temperature and polarisation maps (Ade et al., 2015). They provide us with an estimate of the matter distribution, mainly sensitive to halos located at .

On the other hand, the accretion of baryonic matter in halos also creates active astrophysical objects. They can host violent phenomena, such as, e.g., supernova explosions and relativistic outflows, which are able to accelerate particles to high-energies. Particles with GeV-TeV energy interacting with the ambient medium emit -ray radiation, mostly by means of production and decay of neutral pions, inverse Compton scattering, and non-thermal bremsstrahlung. In addition, the same dark matter which forms the halos could produce -rays, through its self-annihilation or decay. In the past few decades, the all-sky diffuse -ray emission has been measured, but its origin and composition remain key open questions in high-energy astrophysics. The featureless energy spectrum of the isotropic -ray background (IGRB) (Abdo et al., 2010; Ackermann et al., 2014) and its flat angular power spectrum (APS) (Ackermann et al., 2012) make the IGRB identification a complex task. The cross-correlation of the IGRB with large scale structure tracers is a very valuable technique for understanding its composition (Xia et al., 2011; Camera et al., 2013; Fornengo & Regis, 2013; Ando, Benoit-Lèvy & Komatsu, 2014; Shirasaki, Horiuchi & Yoshida, 2014; Ando, 2014).

In this work, we first show that the lensing potential map estimated by Planck and the -ray sky observed by Fermi-LAT do correlate, by reporting a measurement of their cross-correlation APS. This stems from their common origin associated to extragalactic structures, and we discuss the extragalactic -ray background (EGB) properties which can explain the measurement.

The adopted cosmological model throughout this paper is the six-parameter CDM Planck best fitting model reported in (Ade et al., 2013).

## 2. Data and Analysis

We use the -ray measurements obtained by the Fermi-LAT in its first 68 months of operations, from early August 2008 to late April 2014. We have processed the data with the Fermi Science Tools version v9r32p5, using the Pass7-reprocessed instrument response functions for the clean event class (P7REP CLEAN V15) for both front and back conversion types of events, which have been taken together. We have selected photon counts from 700 MeV to 300 GeV, subdivided into 70 energy bins (uniform in log scale) and mapped with a pixel size of 0.125 (suitable for subsequent HEALPIX projection with ). The Fermi-LAT exposure maps have been derived on the same energy grid and resolution, and we adopted a step size , in order to have sufficiently refined exposures. From the count and exposure map cubes, we have finally derived the full-sky flux maps. For the cross-correlation analysis, we have grouped the energy sections in 6 bins (with boundaries at: 0.7, 0.99, 2.0, 5.1, 10.2, 48.7, 300 GeV).

The maps are contaminated by the galactic foreground: since we are interested in the extragalactic signal only, the maps have been cleaned by subtracting the Fermi-LAT galactic model gll_iem_v05, which can be obtained from the Fermi-LAT website^{1}^{1}1http://fermi.gsfc.nasa.gov/ssc/data/access/lat/ BackgroundModels.html.

We account for Fermi-LAT PSF attenuation (which can be relevant at the angular scales of interest) by correcting the measured APS through a beam window function built as described in (Ackermann et al., 2012).

As part of its public data releases of 2015 (2013), the Planck Collaboration provided a CMB lensing convergence (potential) map that was employed in our analysis. We use the convergence harmonic coefficients instead of the potential to reduce the steepness of the APS as a function of multipoles. The methodology followed to derive a convergence map from the 2013 potential map is described in (Ade et al., 2013).

We then mask regions contaminated by galactic foreground and extragalactic sources. For the lensing maps, we adopt the publicly released masks, that both preverse about 70% of the sky. We stress that these masks largely account for the galactic dust emission and the carbon-monoxide lines (that may correlate to -ray foreground). In order to mitigate against multipole mixing, we further use an apodisation over 5 degrees for the Planck 2013 analysis. For the -ray maps we prepare two masks combining the Planck 2013 lensing mask, a cut for galactic latitudes and excluding a angular radius around each source in the 2-year Fermi-LAT catalog (2FGL) (Nolan et al., 2012) and the 4-year Fermi-LAT catalog (3FGL) (Acero et al., 2015), respectively. The 2FGL and 3FGL masks are apodized over and , respectively (the first choice being meant to provide us with a more conservative test), and the resulting effective sky fraction available is about 24% (2FGL) and 23% (3FGL). We explored different apodizations and sets of galactic masks (including larger galactic cuts and an additonal mask for the region of the so-called “Fermi Bubbles” (Su et al., 2010)), finding consistent results.

The cross-correlation APS between the Planck lensing map and the Fermi-LAT -ray map is estimated using a pseudo- approach (Hivon et al., 2001). To this aim, we make use of the publicly available tool PolSpice (Chon et al., 2004; Szapudi, Prunet & Colombi, 2001). Although the PolSpice algorithm properly deconvolves the signal APS from mask effects, it is known not to be a minimum variance algorithm (Efstathiou, 2003). Thus the associated covariance matrix is likely to be an overestimation of the actual uncertainty, and the significances reported throughout the paper can in turn be considered as conservative.

We band-pass filter the cross-correlation APS in the multipole range in order to reduce possible contamination from systematic effects. This multipole range was defined in (Ade et al., 2013) as a confidence interval retaining of the lensing information, with multipoles requiring large mean-field bias corrections. Similarly, multipoles above few hundreds are hardly accessible with the Fermi-LAT sensitivity and angular resolution, and low multipoles correspond to the scales where the foreground cleaning has the largest impact (Ackermann et al., 2012). Since the expected signal is predicted to scale as (see next Section), the analysis is performed in terms of . In order to further mitigate mask mode-mixing, we average the APS in six linear bins of width .

First, we measure the cross-correlation APS between the CMB lensing and a single -ray map derived from the integrated counts at GeV. A hint of a signal in the low- range (with a peak at ) is present, while the larger- bins are compatible with no deviation from a null signal. We estimate the global significance of this low- peak evaluating the ratio of the measured APS over its error, , in three multipole-bins: , , . The errors are the diagonal elements of the covariance matrix obtained from PolSpice. We found the off-diagonal terms of the binned covariance matrix to be negligible. Considering the Planck 2015 map with 3FGL mask (which is our reference analysis), the significances in the three multipole bins amount to: , and , respectively. The significance of the first bin for the four analyses arising from the combination of the Planck maps (2013 and 2015 releases) and the -ray point-source masks (2FGL and 3FGL) is reported in the first line of Table 1.

In order to better exploit all the available information encoded in the maps, we can combine the cross-correlation from the different -ray energy bins introduced above. Since the EGB spectrum roughly scales with (see inset in Fig. 1), low energy bins have larger statistics. On the other hand, the Fermi-LAT point spread function significantly improves at high energy (Ackermann et al., 2012). We therefore expect an information gain by splitting the signal in different energy bins. A minimum variance combination of the 6 single -bin measurements in a given multipole bin can be defined as:

(1) |

where is the sub-matrix for the covariance in the bin , and . Note that, after having checked the stability of our results against the inclusion of the correlation among different multipole bins , we choose not to include them for simplicity. We normalize by means of the factor (with and ) to make it approximately flat in energy.

The computation of the full covariance matrix including correlation among different -bins is not straightforward. Whereas the correlation terms among different multipole bins within an -bin are provided by PolSpice, we estimate the off-diagonal correlation between the and bins (with ) using a two-step process. We first derive a semi-analytic Gaussian approximation (averaged in the multipole bin ):

(2) |

where is the cross-correlation APS, estimated using a benchmark theoretical prediction discussed in the next Section. (Note that this term is in any case subdominant in Eq. (2)). and are the auto-correlation APS that we estimate from the corresponding maps using PolSpice and is the cross-correlation APS between two energy bins and . As a sanity test, we checked that the noise-subtracted estimate (where is the power spectrum of the shot noise and is the beam function) agrees well with the autocorrelation APS reported by the Fermi-LAT Collaboration (Ackermann et al., 2012). Similarly, our is consistent with theoretical expectations, once corrected for the noise APS provided in the Planck public data release (Ade et al., 2013). The factor corrects for the effective available fraction of the sky, but Eq. (2) might actually underestimate the impact of masks. To have a more conservative error estimate we derive a scaling coefficient from , where is obtained from PolSpice and from Eq. (2), and then define the off-diagonal terms of the covariance matrix as . The reliability of this scaling is further supported by the fact we are using the same mask for all the -ray maps.

The combined APS of Eq. (1) is shown in Fig. 1 for the four cases considered. Error bars are given by . The different analyses are in excellent agreement with each other. As for the analysis with gamma-rays integrated above 1 GeV, we estimate the significance of the cross-correlation signal in the multipole-bins , , . The significances now amount to: , and , respectively. A comparison with the results of the previous analysis shows that by adding spectral information increases the significance of the signal in the low- sector, while in the larger- bins the cross-correlation are still compatible with zero. The results obtained so far therefore show an evidence of correlation for multipoles below .

Energy | Multipole | Statistical | Significance | |||||

test | P15-3FGL | P15-2FGL | P13-3FGL | P13-2FGL | ||||

Single -bin | GeV | Single -bin | ||||||

6 -bins | GeV | Single -bin | ||||||

6 -bins | GeV | 6 -bins, | Model fitting |

As a cross-check for the stability of -ray data, we repeat the analysis considering the data from the first 150 weeks and subsequent 150 weeks separately. The obtained APS are compatible and, once combined together, very closely resemble the APS of the full period presented above.

The subtraction of the galactic foreground in the -ray maps has a significant systematic uncertainty related to the modeling of the galactic diffuse emission, which can affect anisotropies on large scales (Ackermann et al., 2012). The foreground residuals in the lensing map are instead thought to be more under control since they do not show up in the autocorrelation studies (Ade et al., 2015). Assuming the lensing map to be free from galactic contaminations, the presence of a gamma-ray galactic foreground in the maps would not provide a cross-correlation signal. Rather it would only act as a noise term. To test this, we performed the same analysis discussed above but employing -ray maps where the foreground was not subtracted. We found the same central values for the cross-correlation APS points, but with larger errors (and so lower statistical significance), consistent with the fact that the galactic foregrounds contribute to the error budget but not to the signal. This suggests that possible contaminations of the APS from a galactic foreground bispectrum are small.

In the next Section we will show that the derived APS can be explained in terms of gamma-rays emission from astrophysical sources emitting mostly at intermediate redshifts.

## 3. Interpretation

We now move on to discuss the agreement between theoretical models and the measurements reported in Fig. 1.

In the Limber approximation (Limber, 1953), the theoretical two-point cross-correlation APS can be computed as:

(3) |

where denotes the radial comoving distance, and are the window functions for lensing and rays, and is the three-dimensional (3D) power-spectrum (PS) of the cross-correlation. For the latter we follow the halo model approach (see, e.g., (Cooray & Sheth, 2002) for a review), where can be split in the one-halo and two-halo components as (see (Fornengo & Regis, 2013) for their expressions).

The CMB lensing window function is given by (Bartelmann, 2010):

(4) |

where is the Hubble constant, is the matter-density parameter, and is the comoving distance to the last-scattering surface.

The window function for a -ray emitter is (see, e.g., (Camera et al., 2013)):

(5) |

where is the -ray luminosity per unit energy range, is -ray luminosity-function (GLF), and is the optical depth for absorption (Stecker et al., 2007).

We consider four different extragalactic -ray populations: star forming galaxies (SFG), misaligned AGN (mAGN), and two subclasses of blazars, BL Lacertae (BL Lac) and flat spectrum radio quasars (FSRQs). The GLFs of the last three source classes are taken from the best-fit models of, respectively, (Di Mauro et al., 2014), (Ajello et al., 2014), and (Ajello et al., 2012). In the case of SFG we consider the infrared luminosity function from (Gruppioni et al., 2013) (adding up spiral, starburst, and SF-AGN populations of their Table 8), and linking and infrared luminosities by means of the relation derived in (Ackermann et al., 2012). The energy spectrum is assumed to be a power-law with spectral indexes (SFG), (mAGN), (BL Lac), and (FSRQ). The model fairly reproduces Fermi-LAT measurements for both the EGB (see the upper-right inset of Fig. 1) and the -ray autocorrelation APS. For the latter, we found a flat APS (given by the 1-halo term and dominated by BL Lac contribution) with for GeV.

The cross-correlation power spectrum at the intermediate scales considered here is mostly set by the linear part of the clustering, , which is similar in the various cases (i.e., it is related to the linear total matter PS ) except for the specific bias term, with negligible contribution from . In other words, we approximately have , where the “effective” bias of a -ray population is: with being the bias between the -ray source and matter, as a function of luminosity and redshift. To estimate the latter, we use the halo bias (Sheth & Tormen, 1999) (setting and the relation (setting the mass of the halo hosting astrophysical objects with a certain luminosity ), as described in (Camera et al., 2014). Comparing the bias of blazars obtained in our analysis with the bias derived in (Allevato et al., 2014), we found the latter to be somewhat larger than our estimates. We obtain a mean mass hosting the object of (BL LAC) and (FSRQ) contrary to of (Allevato et al., 2014). However, our measurement probes the unresolved (individually fainter) component which reside in less massive halos than the brightest blazar sub-sample considered in (Allevato et al., 2014). Thus the two results are not in contradiction.

The cross-correlation APS predicted in the models of the four -ray emitters described above and their collective contribution are shown in Fig. 1.

With the theoretical model at hand we can fit its overall amplitude by minimizing the , which is computed by means of the full covariance matrix introduced above. The statistical significance of the model is derived computing the between null signal and best-fit model. We obtain with 3.0 significance which shows a statistically significant preference for a signal with the correct features expected from the extragalactic gamma-ray emission.

The window functions of the considered -ray populations are all peaked at . To explore in a more general way the kind of -ray model preferred by the data, we compute in Fig. 1 the signals from two Gaussian window functions , one peaked at low redshift (model G0.1 with ), and one peaked at high redshift (model G2 with and ), both normalized to match the Fermi-LAT EGB measurement above 1 GeV (and bias modelled as for mAGNs). We found () and (). For peaked at the relative contribution of small (more distant) objects with respect to larger objects increases, while no power is detected at small scales (above ). This slightly reduces the statistical significance (although with the current data accuracy we cannot exclude this possibility). On the contrary, peaked at low would provide the right bump at low , increasing the statistical significance. However, the large value of the overall amplitude translates into , which is typically way too large for a low- population (see e.g., (Cooray & Sheth, 2002)). Note also that, since the window function of the CMB lensing peaks at moderately high redshift, as mentioned in the Introduction, its overlapping is more effective with high- -populations rather than low- emitters. Therefore, in the latter case, the required becomes slightly larger.

Above arguments seem to suggest that, in order to reproduce the observed cross-correlation, the bulk of -ray contribution to the EGB have to reside at intermediate redshift.

Fig. 2 shows the measured cross-correlation APS for different energy bins and averaged in the multipole bin . The spectrum is consistent with the benchmark model and similar to the Fermi-LAT EGB spectrum (having spectral index close to ), although possibly slightly softer.

## 4. Discussion and Conclusions

We reported the first indication of a cross-correlation between the unresolved -ray sky and CMB lensing. The analysis also points towards a direct evidence that the IGRB is of extragalactic origin. The analysis has been based on the -ray data of the first 68 month of operation of the Fermi-LAT and on the 2013 public release by the Planck Collaboration of the CMB lensing potential map. Current models of AGN and SFG can fit well the amplitude, angular dependence and energy spectrum of the observed APS. The size of the signal appears to be robust against variations of the analysis assumptions. Data exhibit a preference for a signal with the correct features expected from the extragalactic gamma-ray emission with a significance.

The forthcoming Fermi-LAT Pass-8 reprocessed events will allow for a more refined assessment of the signal. Moreover, the technique presented in this work can be also applied for cross-correlating the -ray sky with probes of the large scale structure of the Universe at different redshifts (such as galaxy catalogues and weak-lensing surveys). Such a tomographic analysis of the EGB will provide invaluable information about its composition.

Contaminations from foreground, either real (e.g., a dust or point sources bi-spectrum) or spurious, cannot, at present, be totally excluded, but they are significantly disfavoured. A model of -ray populations built to explain the EGB, and not tuned to the measurement presented here, matches well the data both in features and normalization. More generically, a population of extragalactic -ray emitters following matter clustering at large scales with GLF peaked at intermediate and with agrees well with the data, once the associated EGB is normalized to fit the Fermi-LAT measurement of the IGRB. On the contrary, if, for example, the contribution to the IGRB is reduced to 50%, the required bias would become , which is likely unrealistically large. This implies that the presented results can be considered as a first direct proof that the majority of the IGRB is emitted by extragalactic structures.

## References

- Abdo et al. (2010) Abdo, A. A., et al. [Fermi-LAT Collab.] 2010, Phys. Rev. Lett., 104, 101101.
- Acero et al. (2015) Acero, F., et al. [Fermi-LAT Collab] 2015, arXiv:1501.02003.
- Ackermann et al. (2012) Ackermann, M., et al. [Fermi-LAT Collab] 2012, Phys. Rev. D, 85, 083007.
- Ackermann et al. (2012) Ackermann, M., et al. [Fermi-LAT Collab] 2012, ApJS, 203, 4.
- Ackermann et al. (2012) Ackermann, M. et al. [Fermi-LAT Collab] 2012, ApJ, 755, 164.
- Ackermann et al. (2014) Ackermann, M. et al. [Fermi-LAT Collab] 2014, arXiv:1410.3696.
- Ade et al. (2013) Ade, P. A. R., et al. [Planck Collab] 2014, A&A, 571, AA16.
- Ade et al. (2013) Ade, P. A. R., et al. [Planck Collab] 2014, A&A, 571, AA17.
- Ade et al. (2015) Ade, P. A. R., et al. [Planck Collab] 2015, arXiv:1502.01591.
- Ajello et al. (2012) Ajello, M., et al. 2012, ApJ, 751, 108.
- Ajello et al. (2014) Ajello, M. et al. 2014, ApJ, 780, 73.
- Allevato et al. (2014) Allevato, V., Finoguenov, A., & Cappelluti, N. 2014, ApJ, 797, 96.
- Ando (2014) Ando, S. arXiv:1407.8502 [astro-ph.CO].
- Ando, Benoit-Lèvy & Komatsu (2014) Ando, S., Benoit-Lèvy, A., & Komatsu, E. 2014, Phys. Rev. D, 90, 023514.
- Bartelmann (2010) Bartelmann, M. 2010, Class. Quant. Grav., 27, 233011
- Blanchard & Schneider (1987) Blanchard, A., & Schneider, J., 1987, A&A, 184, 1.
- Camera et al. (2013) Camera, S., et al. 2013, ApJ, 771, L5.
- Camera et al. (2014) Camera, S., et al. arXiv:1411.4651 [astro-ph.CO].
- Chon et al. (2004) Chon,G., et al. 2004, MNRAS, 350, 914.
- Cooray & Sheth (2002) Cooray, A., & Sheth, R. 2002, Phys. Rep., 372, 1.
- Di Mauro et al. (2014) Di Mauro, M. et al. 2014, ApJ, 780, 161.
- Efstathiou (2003) Efstathiou G., 2004, MNRAS, 349, 603.
- Fornengo & Regis (2013) Fornengo, N., & Regis, M. 2014, Front. Physics, 2, 6.
- Gruppioni et al. (2013) Gruppioni, C., et al. 2013, MNRAS, 432, 23.
- Hivon et al. (2001) Hivon, E., et al. 2002, ApJ, 567, 2.
- Lewis & Challinor (2006) Lewis, A., & Challinor, A. 2006, Phys. Rept., 429, 1.
- Limber (1953) Limber, D. N. 1953,ApJ, 117, 134L.
- Nolan et al. (2012) Nolan,P. L., et al. 2012, ApJS, 199, 31.
- Okamoto & Hu (2003) Okamoto, T., & W. Hu, W. 2003, Phys. Rev. D, 67, 083002.
- Sheth & Tormen (1999) Sheth, R. K., & Tormen, G. 1999, MNRAS, 308, 119.
- Shirasaki, Horiuchi & Yoshida (2014) Shirasaki, M., Horiuchi, S., & Yoshida, N. 2014, Phys. Rev. D, 90, 063502.
- Stecker et al. (2007) Stecker, F. W., Malkan, M. A., & Scully, S. T. 2007, ApJ, 658, 1392.
- Su et al. (2010) Su, M., Slatyer, T. R., & Finkbeiner, D. P. 2010, ApJ, 724, 1044.
- Szapudi, Prunet & Colombi (2001) Szapudi, I., Prunet, S., & Colombi, S. 2001, ApJ, 561, L11.
- Xia et al. (2011) Xia, J.-Q., et al. 2011, MNRAS, 416, 2247.