Issue 
A&A
Volume 648, April 2021



Article Number  A115  
Number of page(s)  12  
Section  Cosmology (including clusters of galaxies)  
DOI  https://doi.org/10.1051/00046361/202040155  
Published online  23 April 2021 
Impact of baryons in cosmic shear analyses with tomographic aperture mass statistics
^{1}
AixMarseille Univ., CNRS, CNES, LAM, Marseille, France
email: nicolas.martinet@lam.fr
^{2}
Dipartimento di Fisica, Sezione di Astronomia, Università di Trieste, Via Tiepolo 11, 34143 Trieste, Italy
^{3}
INAF – Osservatorio Astronomico di Trieste, Via Tiepolo 11, 34131 Trieste, Italy
^{4}
IFPU – Institute for Fundamental Physics of the Universe, Via Beirut 2, 34151 Trieste, Italy
^{5}
INFN – Sezione di Trieste, 34100 Trieste, Italy
^{6}
Astrophysics Research Institute, Liverpool John Moores University, 146 Brownlow Hill, Liverpool L3 5RF, UK
^{7}
Institute for Astronomy, University of Edinburgh, Royal Observatory, Blackford Hill, Edinburgh EH9 3HJ, UK
^{8}
INAF – Osservatorio di Astrofisica e Scienza dello Spazio di Bologna, Via Gobetti 93/3, 40129 Bologna, Italy
^{9}
INFN – Sezione di Bologna, Viale Berti Pichat 6/2, 40127 Bologna, Italy
^{10}
University Observatory Munich, Scheinerstr. 1, 81679 Munchen, Germany
^{11}
MaxPlanck Institut fur Astrophysik, KarlSchwarzschild Str. 1, 85741 Garching, Germany
Received:
17
December
2020
Accepted:
15
February
2021
NonGaussian cosmic shear statistics based on weaklensing aperture mass (M_{ap}) maps can outperform the classical shear twopoint correlation function (γ2PCF) in terms of cosmological constraining power. However, reaching the full potential of these new estimators requires accurate modeling of the physics of baryons as the extra nonGaussian information mostly resides at small scales. We present one such modeling based on the Magneticum hydrodynamical simulation for the KiDS450 and DESY1 surveys and a Euclidlike survey. We compute the bias due to baryons on the lensing PDF and the distribution of peaks and voids in M_{ap} maps and propagate it to the cosmological forecasts on the structure growth parameter S_{8}, the matter density parameter Ω_{m}, and the dark energy equation of state w_{0} using the SLICS and cosmoSLICS sets of darkmatteronly simulations. We report a negative bias of a few percent on S_{8} and Ω_{m} and also measure a positive bias of the same level on w_{0} when including a tomographic decomposition. These biases reach ∼5% when combining M_{ap} statistics with the γ2PCF as these estimators show similar dependency on the AGN feedback. We verify that these biases constitute a less than 1σ shift on the probed cosmological parameters for current cosmic shear surveys. However, baryons need to be accounted for at the percentage level for future Stage IV surveys and we propose to include the uncertainty on the AGN feedback amplitude by marginalizing over this parameter using multiple simulations such as those presented in this paper. Finally, we explore the possibility of mitigating the impact of baryons by filtering the M_{ap} map but find that this process would require suppressing the smallscale information to a point where the constraints would no longer be competitive.
Key words: cosmology: observations / dark energy / largescale structure of Universe / dark matter / gravitational lensing: weak / surveys
© N. Martinet et al. 2021
Open Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
1. Introduction
Weak lensing cosmic shear is one of the most powerful cosmological probes of the latetime Universe. So far, most analyses have focused on studying the correlation between the shape distortions of pairs of galaxies as a function of their separation: the socalled shear twopoint correlation function (γ2PCF; e.g., Troxel et al. 2018; Hikage et al. 2019; Asgari et al. 2021). However, it is becoming clear that other estimators, and in particular those based on weaklensing mass maps, outperform the standard γ2PCF (e.g., Dietrich & Hartlap 2010; Ajani et al. 2020; Zürcher et al. 2021; Coulton et al. 2020; Martinet et al. 2021). Indeed, mass maps are highly sensitive to the nonGaussian part of the matter distribution that arises from the nonlinear growth of structures, which contains information that is overlooked by twopoint estimators. Consequently, the combination of both probes yields even tighter constraints, as seen in recent applications to observational data (e.g., Martinet et al. 2018; HarnoisDéraps et al. 2020). With future Stage IV surveys, this combination is also expected to improve not only our measurement of the growth of structure parameter S_{8} by a factor of two, but also that of the dark energy equation of state w_{0} (Martinet et al. 2021) and of the sum of neutrino masses Σm_{ν} (Li et al. 2019), by factors of three and two, respectively.
Nevertheless, these nonGaussian estimators are difficult to predict theoretically because of limits in our understanding of the nonlinear growth of structures (see e.g., Fan et al. 2010; Lin & Kilbinger 2015; Shan et al. 2018; Giocoli et al. 2018a; Barthelemy et al. 2021, for some attempts) and are instead modeled with Nbody simulations. This can significantly increase the computational cost of such analyses, but resorting to Nbody simulations is also necessary to accurately model the γ2PCF at scales affected by nonlinearities (e.g., Euclid Collaboration 2020). Moreover, many public simulation suites can be exploited for this purpose, including the Scinet LIghtCone Simulations (HarnoisDéraps et al. 2015, SLICS), the cosmoSLICS (HarnoisDéraps et al. 2019), or the MassiveNuS (Liu et al. 2018).
As mass map estimators focus on small scales of typically around a few arcminutes, they can be severely affected by baryonic feedback, which can bias the inferred cosmological constraints. To avoid this, an effective approach is to quantify the impact of baryons on the estimator with hydrodynamical simulations, and to apply a correction factor to the model extracted from dark matter (DM) only simulations. For mass maps, this approach was pioneered in Yang et al. (2013) and Osato et al. (2015). However, feedback from active galactic nuclei (AGN) was not included in the first analysis and only with a low amplitude of the feedback in the second. As a result, these two studies found a mild impact of the baryonic physics on the distribution of peaks in mass maps and therefore underestimated the bias on cosmological parameters (e.g., Weiss et al. 2019; Coulton et al. 2020).
Stateoftheart hydrodynamical simulations with box lengths of hundreds of megaparsecs and including realistic AGN feedback later enabled refinement of the measure of the bias on mass map estimators due to baryons. Fong et al. (2019) measured a ∼10% reduction in the number of high signaltonoiseratio (S/N) M_{ap} peaks in the BAHAMAS simulations (McCarthy et al. 2017), with a particular look at possible degeneracies between the effect of baryons and that of massive neutrinos. Osato et al. (2021) measured biases of the same order of magnitude on the number of peaks, the number of minima, and the lensing probability distribution function (PDF) in the TNG 300 simulation (Pillepich et al. 2018), also highlighting a less pronounced bias at higher redshift.
An interesting alternative to hydrodynamical simulations is the baryonification method described in Schneider & Teyssier (2015) where particle positions are shifted in DMonly simulations to mimic the impact of baryons. This method now accurately reproduces AGN feedback and star formation compared to hydrodynamical simulations (Aricò et al. 2020) and could offer an efficient way of decreasing the computational resources needed to include baryonic effects in cosmological models. This technique is however not tested in the present article as it requires the particle positions, which are generally not stored for a posteriori applications. Weiss et al. (2019) applied this baryonification method to model the impact of baryons on peak statistics and also noted that the latter could be mitigated by smoothing the small scales by applying a Gaussian filter to the mass map. This smoothing is nevertheless likely to also reduce the statistical power of the mass map estimators, a hypothesis that can only be verified by propagating this bias to the cosmological parameter inference.
Recently, Coulton et al. (2020) performed this propagation and measured the impact of baryons directly on the forecasts of the matter density Ω_{m}, the amplitude of fluctuations A_{S}, and the sum of the neutrino masses for peaks and minima using the BAHAMAS simulations. In a nontomographic approach, these latter authors found larger biases for peaks than for minima in their LSSTlike mock data, concluding that the latter is potentially more robust against baryons.
Building on these previous analyses, and exploiting the cosmological analysis pipeline introduced in Martinet et al. (2021), we measure for the first time the effect of baryons on the dark matter and dark energy cosmological parameters in a tomographic StageIV lensing survey setup. We focus on the particular case of aperturemass (M_{ap}; Schneider 1996) maps, which are particularly well suited for cosmological analyses. We model the effect of baryons with the Magneticum^{1} hydrodynamical simulation suite (e.g., Castro et al. 2020), which includes all key ingredients about the physics of baryons, such as AGN feedback (Springel et al. 2005; Fabjan et al. 2010; Hirschmann et al. 2014) and star formation (Springel & Hernquist 2003). Both effects redistribute the matter in and around DM haloes, but the exact amplitude of the feedback varies between simulations (see Chisari et al. 2019, for a recent review on feedback in cosmological hydrodynamical simulations). Extending the simulation suites of Martinet et al. (2021) based on the SLICS and the cosmoSLICS, we construct Euclidlike mocks from the Magneticum hydrodynamical simulation. We measure the impact of baryons on the γ2PCF, on the lensing peaks, minima, and on the lensing PDF, in the form of a multiplicative baryon bias correction factor. We then propagate this correction into the cosmological inference pipeline described in Martinet et al. (2021), and investigate the impact on the parameter forecasts. Finally, we explore various mitigation schemes to decrease the baryondominated smallscale contribution to the M_{ap} computation.
We introduce the Magneticum simulation and compare it to other stateoftheart hydrodynamical simulations in Sect. 2. We present the methodology that we employ to measure the bias due to baryons in Sect. 3. We then measure their impact on different data vectors (DV) in Sect. 4.1 and propagate the effect to the cosmological forecasts in Sect. 4.2. We test different mitigation setups in Sect. 5 and conclude in Sect. 6. Finally, we adapt our mocks to the KiDS and DES surveys in Appendix A and measure the impact of baryons on the cosmological constraints by Martinet et al. (2018), HarnoisDéraps et al. (2020) in these two surveys.
2. Modeling baryonic effects
2.1. The Magneticum hydrodynamical simulation
The Magneticum suite (Biffi et al. 2013; Saro et al. 2014; Steinborn et al. 2015, 2016; Dolag et al. 2015, 2016; Teklu et al. 2015; Remus et al. 2017a; Castro et al. 2020) is a compilation of Nbody and hydrodynamical simulations describing the cosmic evolution of the Universe. In total, the suite follows up to 2 × 10^{11} particles divided in DM, gas, stars, and black holes. The simulations were performed with the TreePM+SPH code PGadget3 – a higher performance and more efficient version of the publicly available Gadget2 code (Springel 2005) developed concomitantly to its successor Gadget4 (Springel et al. 2020). Notably, the smoothedparticlehydrodynamics (SPH) solver implements the improved model of Beck et al. (2016). The particle dynamics is coupled to different astrophysical effects such as radiative cooling, heating by a uniform evolving UV background, star formation (Springel & Hernquist 2003), stellar evolution, and chemical enrichment processes (Tornatore et al. 2007). Cooling is implemented following the metallicitydependent formulation presented in Wiersma et al. (2009) using cooling tables produced by the publicly available CLOUDY photoionization code (Ferland et al. 1998). Lastly, AGN feedback and black hole growth are modeled as described in Hirschmann et al. (2014).
The subgrid physics model of our simulation reproduces a long list of observations, from galaxy properties (Teklu et al. 2015; Remus et al. 2017b) and AGN population (Hirschmann et al. 2014; Steinborn et al. 2016), to the intergalactic and intercluster medium (Dolag et al. 2016; Bocquet et al. 2016; Gupta et al. 2017). Of particular importance for this work is the robustness of the AGN feedback of our model, which has a clear footprint on the matter powerspectrum, suppressing its amplitude with respect to the DMonly case on k ∼ 10 Mpc h^{−1} by ∼15%. The specific suppression range and amplitude depend strongly on the AGN model. In Fig. 1 we compare the ratio of the matter powerspectrum computed from Box 2 (see Table 1 for details on the simulations) Hydro and DMonly with other simulations presented in Chisari et al. (2019): BAHAMAS (McCarthy et al. 2017), EAGLE (Schaye et al. 2015), Horizon (Dubois et al. 2014), Illustris (Vogelsberger et al. 2014), IllustrisTNG (Pillepich et al. 2018), and OWLS (Schaye et al. 2010). As can be observed, Magneticum provides an AGN feedback suppression consistent with other simulations, in particular with BAHAMAS and to a lesser extent with TNG300, which were used in other recent mass map analyses. We note that the AGN feedback also controls the gas fraction of halos (and the conversion efficiency into stars) such that the simulations presented in Fig. 1 also differ in these quantities although the matter power spectrum is a good indicator with which to assess the current theoretical uncertainties in modeling baryonic physics.
Fig. 1.
Ratio of power spectra at z = 0 between various hydrodynamical simulations and their corresponding DMonly run. The Magneticum used in this analysis is representative of an average behavior with a loss of power of about 15% at k ∼ 10 Mpc h^{−1} due to AGN feedback. 
Subset of the Magneticum simulation suite used in this work.
2.2. Past lightcone reconstruction
The Magneticum past lightcone reconstructions are performed using the Simulation LIght conE buildeR code – SLICER^{2} – closely following the pipeline presented in Castro et al. (2018) and forked from a MapSim branch (Giocoli et al. 2015, 2018b; Hilbert et al. 2020). Briefly, lightcones are built in postprocessing, assigning particles to predetermined 2D mass maps according to the triangularshaped cloud mass assignment scheme. The geometry of the past lightcones is a squarebased pyramid (in angular coordinates) where the observer is located at the z = 0 vertex, and which extends to z_{max}. The opening angle is chosen to be 10° and the angular resolution of the lightcone mass planes is 3.6 arcsec. SLICER allows mass maps thicknesses that are a rational fraction of the box size, which we have chosen to be half of the box size. Every two mass maps, particles are randomly shifted and reflected with respect to one of the box axes (accounting for periodic boundary conditions) in order to avoid the repetition of the same cosmic structure along the line of sight.
Subsequently, mass maps are converted into surface density maps Σ(x, y) as
where n indicates the number of particles, m_{j} the interpolated contribution of the jth particle to the pixel at position x, y, and L_{p} is the physical size of the pixel in units of Mpc h^{−1}. Given a source plane, convergence maps κ(x, y) are created by weighting these maps by the corresponding critical surface mass density,
and integrating over the past lightcone. Here, c indicates the speed of light, G the Newton’s constant, and D_{l}, D_{s}, and D_{ls} are the angular diameter distances between observer–lens, observer–source, and lens–source, respectively. The shear components (γ_{1}, γ_{2}) are obtained from the standard inversion technique of Kaiser & Squires (1993).
As Box 2b Hydro has not been run down to z = 0, our past lightcones are built from grafting Box 2 in the range z = [0, 0.248] with Box 2b for z = [0.248, 3.44]. We used the DMonly runs to validate our approach, and verified that for the source redshift distribution of interest in this paper, grafting only marginally affects the lensing statistics for z_{s} ≳ 1.0 with respect to the same statistics inferred from a contiguous lightcone^{3}. In total, we produced 20 pseudoindependent lightcones – 10 Hydro and 10 DMonly – with properties summarized in Table 1.
3. Methodology
In this section, we first give a brief description of our cosmological forecasting pipeline in Sect. 3.1, and then detail the calculation of baryonic bias in Sect. 3.2.
3.1. From shear to cosmology
Cosmological forecasts are computed with the same pipeline as in Martinet et al. (2021). While we recapitulate the salient points of this analysis here, we refer the reader to this publication for more details.
– wCDM Simulations: the cosmology dependence of M_{ap} statistics is emulated with radial basis functions based on measurements from StageIV mock data covering 25 cosmologies in S_{8} − Ω_{m} − w_{0} − h, organised in a latin hypercube. These were built from the cosmoSLICS Nbody simulations (HarnoisDéraps et al. 2019) and contain 50 mocks per cosmology, covering ten lightcones and five shape noise realizations to further increase our precision on the model. The covariance is measured from a separate suite of N_{s} = 928 fully independent ΛCDM mocks extracted from the SLICS Nbody simulations (HarnoisDéraps et al. 2018) with the same survey properties. These StageIV mocks match the expected 30 arcmin^{−2} galaxy density of Euclid with a realistic redshift distribution in the range 0 < z < 3 (Laureijs et al. 2011), and an opening angle of 100 deg^{2} each. While galaxy positions and intrinsic ellipticities are chosen randomly in the covariance mocks, they are fixed across the different cosmologies in cosmoSLICS mocks so as to reduce the impact of shape noise on the model.
– Measurements: From each mock, we compute a 1024 × 1024 pixel M_{ap} map (Schneider 1996) as a convolution between the tangential ellipticity ϵ_{t} around a pixel θ_{0}, and the Schirmer et al. (2007) compensated Q filter adapted to the detection of matter halos:
with
This filter is a simpler form of the NFW (Navarro et al. 1997) profile with an additional exponential attenuation in the center and at the edge of the aperture. x_{c} controls the tilt between the core and the edge of the profile and we introduce an optional inner radius parameter θ_{in} to govern the exponential cutoff in the inner part. The sum in Eq. (3) is carried out over galaxies at position θ_{i} within an aperture of radius θ_{ap} = 10′ by default (corresponding to an effective smoothing scale of θ_{ap}x_{c} = 1.5′) and centered on θ_{0}, and is normalized by the galaxy density n_{gal} within the aperture. We estimate the noise due to intrinsic ellipticities in the mocks to define several M_{ap} map statistics based on the S/N distribution of pixels of these maps: peaks and voids (pixels with values greater or smaller than their 8 neighbors), as well as the full distribution of pixels (1D M_{ap}, often referred to as the lensing PDF). Following Martinet et al. (2021), these distributions are organised in eight bins between −2.5 < S/N < 5.5 for peaks, eight bins between −5 < S/N < 3 for voids, and in nine bins between −4 < S/N < 5 for the 1D M_{ap}. For each mock, we also compute the γ2PCF ξ_{±}(ϑ) using the Athena software (Kilbinger et al. 2014), with the estimators defined as (Schneider et al. 2002): , where the sum is over galaxy pairs ij separated by an angle ϑ and potentially in different tomographic bins. The results are binned in eight angular bins logarithmically spaced between 0.1′ and 60.5′ for ξ_{+} and between 0.5′ and 300′ for ξ_{−}. The quantities ϵ_{t, ×} are the tangential and “cross” components of the ellipticities.
– Tomography: We separate the mock data into five redshift bins: 0 < z_{1} < 0.47, 0.47 < z_{2} < 0.72, 0.72 < z_{3} < 0.96, 0.96 < z_{4} < 1.33, and 1.33 < z_{5} < 3. We reconstruct M_{ap} maps from the galaxy samples in individual redshift slices (autoM_{ap}) and from the combination of multiple redshift bins from two and up to five (the crossM_{ap} terms introduced in Martinet et al. 2021). The latter showed that including the crossmaps yields a significant improvement in the forecast precision. For the γ2PCF, we include both the auto and the crosstomographic correlations.
– Likelihood: Finally, we compute the likelihood of the data given the cosmology p(xπ) on a fourdimensional grid, taking the “fiducial” cosmoSLICS simulations as our observation DV x. We use a Studentt likelihood (Sellentin & Heavens 2016),
which generalizes the multivariate Gaussian, and we account for the noise through the covariance matrix Σ computed at a fixed cosmology π_{0}. In the above equation, x_{m}(π) refers to the DV modeled from the cosmoSLICS at cosmology π and both x and x_{m}(π) correspond to the mean DV over the different noise realizations. The posterior likelihood is obtained from Bayes’ theorem with the cosmoSLICS parameter space as a fixed prior (p(π) = 1 within that space and 0 outside). We also improve the accuracy of the emulator used in Martinet et al. (2021), initially limited by the number of points in the interpolation grid, by performing a second interpolation in the parameter space restricted to the hypervolume where the likelihood is nonzero. We measure the 1σ uncertainty on each cosmological parameter by finding the range of values enclosing 68% of the posterior likelihood previously marginalized over all other parameters. As the posterior distribution on h extends up to the prior limit imposed by the simulation space, we only discuss the forecasts for the other three probed parameters^{4}: S_{8}, w_{0}, and Ω_{m}. All predictions are computed for a 100 deg^{2} area to ensure that the few percent uncertainties in the model (primarily due to the accuracy of the Nbody simulations) are lower than the precision of the forecasts.
3.2. Measuring biases
In this article, we use the pipeline described above to study the impact of the baryonic bias on a StageIV cosmic shear data analysis. However, we note that we could instead propagate biases due to shear measurement uncertainty, mean photometric redshift inaccuracy, galaxy intrinsic alignments, or sourcelens coupling (see e.g., Kacprzak et al. 2016; HarnoisDéraps et al. 2020), but we leave this to future work. The strategy is to bias the observation DV and compare the positions of maximum likelihood to the nobias case.
We use the DMonly and equivalent hydrodynamical runs of the Magneticum simulations to measure the cosmological bias due to baryonic effects. We create galaxy mocks that reproduce the same StageIV survey properties as the cosmoSLICS. In particular, galaxy redshifts, positions, and intrinsic ellipticities are identical to those of the model. We generate 50 realizations of the hydro and DM lightcones: ten different lines of sight to lower the sample variance, each populated with five different realizations of intrinsic ellipticities to converge on an average shape noise contribution. We compute the M_{ap} map in every mock, extract the DV, and measure the ratio between the average over the 50 DVs measured in the hydrodynamical and in the DMonly mocks.
This multiplicative correction factor is computed for each S/N bin and serves to infuse our Magneticum baryon model in our (DMonly) observation data. Infusing baryons into the model avoids being affected by small residual differences between the Magneticum DMonly and the cosmoSLICS DMonly simulations, such as the finite box effect (e.g., HarnoisDéraps & van Waerbeke 2015; Euclid Collaboration 2020), the nonlinear physics modeled by the different Poisson solvers, or the chosen distances between lensing planes (e.g., Takahashi et al. 2017; Zorrilla Matilla et al. 2020). In other words, it ensures that the differences in the likelihood maxima are only due to baryons. We neglect here the possible dependence of baryons on cosmology, a hypothesis well supported by the recent findings of van Daalen et al. (2020).
4. Impact of baryons
We examine in this section the impact of baryons on the different DVs (Sect. 4.1) and on the cosmology inferred by a StageIV survey (Sect. 4.2).
4.1. Impact on computed statistics
Figure 2 shows the impact of baryons on the γ2PCF, presented as the fractional difference between ξ_{±} measured in the hydrodynamical and in the equivalent DMonly runs. This figure shows the baryonic bias in the absence of tomography for visual purposes; however similar curves are observed in the tomographic case. We measure a decrease in the amplitude of ξ_{±} of up to 10−15% at scales below a few arcminutes, which is fully consistent with the suppression in the matter power spectrum seen in Fig. 1.
Fig. 2.
Relative change in the γ2PCF due to baryons. The orange and purple curves represent the ξ_{+} and ξ_{−} estimators, respectively. The error bars correspond to the diagonal elements of the SLICS covariance matrix rescaled to 15 000 deg^{2}. 
In Fig. 3 we show the equivalent measurement on the M_{ap} statistics: voids, peaks, and 1D M_{ap}. The impact of baryons on these DVs is more complex than in the case of the γ2PCF given the different physical origins of each part of the M_{ap} distributions. The most plausible scenario is described in Osato et al. (2021): Overdense regions are diluted due to AGN feedback, leading to highS/N structures with smaller amplitudes. This expelled material can be deposited in lowdensity regions, which likely explains the reduction in pixels with highly negative S/N. We note that the cause of the latter effect is not yet fully understood and could involve other baryonic processes such as interaction of internal and accretion shocks (e.g., Zhang et al. 2020). The increase of the distribution at S/N close to zero accounts for the density smoothing due to this redistribution of matter which leads to a higher number of small S/N structures. This reasoning is deducted from the lensing PDF behavior but holds for peaks and voids. We do not probe the impact of radiative cooling seen at very high S/N in Osato et al. (2021) as we focus here on a lower, more conservative S/N range. This would also necessitate a finer pixel scale than our fiducial 0.59′. Indeed these effects appear to be significant only at scales lower than 0.5′ in the M_{ap} map according to the location of the ξ_{+} upturn in Fig. 2 and consistent with the propagation of the radiative cooling scale of ∼15 Mpc h^{−1} into the lensing PDF using Eq. (11) of Castro et al. (2018).
Fig. 3.
Relative change in the M_{ap} estimators due to baryons. The blue, green, and red curves represent void counts, peak counts, and the lensing PDF, respectively. The error bars correspond to the diagonal elements of the SLICS covariance matrix rescaled to 15 000 deg^{2}. 
Quantitatively, we find a decrease of between 5% and 13% in the extremal values of our DVs depending on the considered DV and S/N bin. In particular, for peaks of S/N = 4 we measure a depletion of ∼5% due to baryons, in perfect agreement with measurements from Fong et al. (2019) and Coulton et al. (2020) on the BAHAMAS simulations, which present a similar amplitude of the baryonic feedback (see Fig. 1). This is also of the same order of magnitude as the results from the TNG simulations (Osato et al. 2021) and the baryonification method (Weiss et al. 2019).
When applying tomography, we also note a decrease in the impact of baryons in each slice compared to the combined case. Using thin source slices, Osato et al. (2021) also found that the baryon bias is lower at higher redshift. However, the lower effect in our case is more likely due to an increase in the noise due to lower galaxy densities as noted by HarnoisDéraps et al. (2020) in their tomographic peak count analysis of DES data. This is supported by the fact that we find similar baryonic effects in all five redshift slices, which we designed to have identical galaxy densities.
In Figs. 2 and 3 the error bars are computed for the expected 15 000 deg^{2} of the Euclid survey by arearescaling the SLICS covariance matrix as Cov→Cov (100/15 000). These figures show that the impact of baryons on all tested estimators will be significant with respect to the statistical precision of Stage IV cosmic shear surveys. For the present 100 deg^{2} analysis, the changes are below the statistical noise, which is why it is necessary to fix the noise in the cosmoSLICS model and the Magneticum mocks.
4.2. Propagation to cosmological constraints
We propagate the impact of baryons on the cosmological parameter forecasts following the method described in Sect. 3. We focus on the peak statistics as it is the most widely used mass map estimator in the literature, but we also present results for voids and the lensing PDF. All values reported in this section are summarized in Table 2.
Forecasts on the biases due to baryons in 100 deg^{2}Euclidlike mocks.
Figures 4–6 show the 1 and 2σ contours of the marginalized 2D and 1D likelihoods for S_{8}, Ω_{m}, and w_{0} for various configurations of the peak statistics: without tomography, with tomography, and combined with the γ2PCF with tomography. In all these figures, the blue contours and curves correspond to the forecasts for the DMonly observations, while the data yielding the green contours are infused with the baryon bias. As already noted in Martinet et al. (2021), we accurately recover the input parameter values for the DMonly case, with only a ∼1σ bias on Ω_{m} in some cases, which is due to the sampling of the initial parameter space (see the reference above for more details).
Fig. 4.
Forecast of cosmological parameters from peak counts in a 100 deg^{2} survey at Euclid depth without tomography. Marginalized 2D (1 and 2σ contours) and 1D (full likelihood) constraints are displayed in blue for the DMonly case and in green when including baryon physics. Dashed lines correspond to the 1σ constraints in the 1D marginalized likelihood. The blue and green crosses indicate the best estimate in the 2D constraints, while the red crosses and lines indicate the input cosmological parameter values. Gray crosses correspond to parameters of the 25 cosmoSLICS simulations that are used to estimate the cosmology dependence of the number of peaks. 
Fig. 5.
Same as Fig. 4 but for a tomographic analysis with five redshift slices and including auto and crossM_{ap} terms between redshift slices. 
Fig. 6.
Same as Fig. 4 but for the combination of peak counts and γ2PCF for a tomographic analysis with five redshift slices and including auto and crossterms between redshift slices. 
We first consider the nontomographic peak results in Fig. 4. We see a shift towards smaller values of S_{8} (%)) and Ω_{m} (ΔΩ_{m} = −0.018 (−6.1%)) when including baryons. As expected, including baryonic feedback mimics the effect of having lower matter density/structure growth, the two being highly degenerate for cosmic shear. These results are fully consistent with the Ω_{m} shift recently found by Coulton et al. (2020) for peaks measured in the BAHAMAS simulations, but are larger than those found in the earlier work of Osato et al. (2015), who reported a −1.5% shift in Ω_{m} from simulations with weaker AGN feedback. Finally, we see negligible changes in w_{0}, in part because the sensitivity to that parameter is quite low in the nontomographic case.
When including tomography (see Fig. 5), the constraining power increases significantly. We measure a very similar effect to that seen in the nontomographic case for S_{8} and Ω_{m} but with slightly smaller shifts (ΔS_{8} = −0.024 (−2.9%) and ΔΩ_{m} = −0.005 (−1.7%)). This behavior is also found in the DESY1 Magneticum mocks (see HarnoisDéraps et al. 2020, and Appendix A): with tomography, the noise in each mass map increases and tends to wash out the effect of baryons. We also find a small positive shift in the maximum of the 1D likelihood for w_{0} (Δw_{0} = 0.035 (3.5%)). However, we note that this result is only supported by a small distortion of the likelihood, which otherwise agrees fairly well with the DMonly case. Although this effect could be physically motivated, it is likely due to degeneracies in the parameter space notably between S_{8} and w_{0}, and it will be interesting to see if it persists in future analyses.
If we now combine peaks with the γ2PCF in the tomographic setup (Fig. 6), the results are similar to those from peaks alone, but accentuated in amplitude. We find ΔS_{8} = −0.037 (−4.4%), ΔΩ_{m} = − 0.019 (−6.5%), and Δw_{0} = 0.047 (4.7%). Both peaks and γ2PCF are affected in a similar manner by baryons and thus their combination shows a larger effect. Although one could have hoped that the impact of baryons would diminish when adding the information of the γ2PCF that partly comes from larger scales, this study demonstrates on contrary that baryonic physics do not vanish in this combination and need to be accounted for in future surveys.
Table 2 shows the forecast biases for all the different tested DVs and their combination with the γ2PCF. Overall, baryons impact the different M_{ap} estimators in a similar manner, as expected from the similarities in how they affect each DV in Fig. 3. We report smaller biases for voids than for peaks, as already noted in Coulton et al. (2020) in the nontomographical case, and confirm this trend when including tomography, with possibly a very low bias on w_{0}. When combined with the γ2PCF we also note a difference of sign in the w_{0} bias which could highlight a different sensitivity of voids to this parameter, but is more likely due to a residual small peak in the likelihood from the interpolation of the DV in this case. The lensing PDF presents a comparable bias to peaks with slightly lower shifts as well. Although voids are the least affected by baryons, this analysis shows that all estimators present biases of a few percent on at least one of the probed cosmological parameters. Considering the degeneracies between cosmological parameters, this highlights the necessity for accounting for baryons when modeling the dependence of nonGaussian M_{ap} estimators on cosmology. Because these results depend on the particular implementation of baryonic feedback processes in the Magneticum simulation, we run an additional test where we only infuse half of the DV baryonic bias to compute the cosmological forecasts. Although we cannot accurately model the response of the M_{ap} statistics to AGN feedback with only one simulation, this case likely corresponds to a much lower feedback amplitude. We find a bias on cosmological parameters which approaches the percent value in the tomographic case. This strong dependency of the cosmological parameters on the amplitude of the infused baryonic effect suggests that the impact of baryons could be mitigated by integrating a modeling of AGN feedback in the likelihood, a possibility which is not studied here.
Although these biases of a few percent are worrisome for future Stage IV cosmic shear surveys, we note that they remain fairly small compared to the statistical precision of current surveys. Our 100 deg^{2} mocks at Euclid depth include about ten million galaxies, a similar number to Stage III surveys. In this case, the biases are always below 1σ for every parameter and configuration, except for the S_{8} parameter when combining M_{ap} estimators and γ2PCF where the bias can reach up to 2σ. In Appendix A, we tailor our mocks to the KiDS450 and DESY1 surveys to verify the impact of baryons on the peak statistics analyses conducted in Martinet et al. (2018) and HarnoisDéraps et al. (2020), respectively. Our findings validate the choice of neglecting baryons in current Stage III peak count analyses.
5. Mitigating baryons with small scale cut
In the case of γ2PCF, the impact of baryons can be mitigated by discarding small scales which are the most affected, as seen from Fig. 2. Although this process decreases the statistical precision as it removes part of the signal, it provides a gain in accuracy without needing to run computationally expensive hydrodynamical simulations. This tradeoff was notably chosen by the DES collaboration in the first year data release (Troxel et al. 2018). Recently, Taylor et al. (2018, 2021) proposed an efficient way to cut the smallscale contribution to cosmic shear twopoint estimators based on the nulling scheme presented in Bernardeau et al. (2014). Alternatively, baryons can be modeled with halobased codes (HMCode, Mead et al. 2021) or from libraries of power spectra (van Daalen et al. 2020), allowing the inclusion of smaller angular scales and increasing the statistical power as in Hildebrandt et al. (2017) and Huang et al. (2021).
As we do not yet have access to a library of hydrodynamical simulations to estimate the baryonic bias, we follow the work of Weiss et al. (2019) and explore different possibilities to mitigate the effect of baryons on M_{ap} statistics. As suggested by their analysis, we vary the size of the aperture θ_{ap} but we also investigate possible variations to the Schirmer et al. (2007)Q filter shape entering Eq. (3). In particular, we vary the aperture filter size θ_{ap}, the tilt parameter x_{c}, and the inner filter radius parameter θ_{in}.
The top part of Fig. 7 shows the impact of varying θ_{ap}, θ_{in}, and x_{c} on the shape of the Q filter. The fiducial values that we use in the rest of the paper are θ_{ap} = 10′, θ_{in} = 0.4′, and x_{c} = 0.15. In this figure, we only present results when varying one parameter at a time for θ_{ap} = 20′, θ_{in} = 2.4′, and x_{c} = 0.45 to better highlight the effect of each parameter. However, we investigate multiple values in the ranges 3.3′ ≤ θ_{ap} ≤ 106.7′, 0.4′≤θ_{in} ≤ 3′, and 0.05 ≤ x_{c} ≤ 1 and explore variations of all three parameters together. Increasing any of these three parameter values reduces the relative importance of the small scales in different manners. Increasing θ_{ap} adds galaxies further away from the aperture center; increasing θ_{in} removes galaxies close to the aperture center without distorting the general weighting; finally, increasing x_{c} upweights distant galaxies without modifying the number of galaxies captured by the aperture. We recompute the baryonic bias from M_{ap} constructed with the new Q filters and observe, in the bottom part of Fig. 7, that none of these methods is completely efficient at removing the impact of baryons. This result is somewhat contradictory to our measurement that the effect of baryons in our simulations concentrates on the central few arcminutes of the M_{ap} profile around massive halos. This is because with mass map statistics we cannot discard close galaxies as is done for γ2PCF: we only reduce their weight when centered on them which does not prevent them from entering other apertures on a θ_{ap} scale. We note a mild improvement when increasing θ_{in} and x_{c}. Increasing θ_{ap} to 20.0′ still decreases the impact of baryons by ∼50% at the extreme S/N values, and only the largest θ_{ap} > 100′ can bring it to zero, but the cost in precision is high, as we show next.
Fig. 7.
Top: profiles of the Q filters defined via Eq. (4) and used to mitigate the impact of baryons by reducing the weight of small scales. Bottom: relative change in the M_{ap} peak counts due to baryons for the various mitigation setups. 
We apply these filter modifications to the measurements from the cosmology and covariance mocks, and carry out a full cosmological forecast for each of these. In Table 3, we present the variations of the bias on the inferred cosmology for these different configurations relative to the fiducial Q filter. We also report the change in the forecast precision due to the reduction in smallscale information. We show these results for peaks in a configuration without tomography but we find similar behavior for other M_{ap} estimators and including tomography: the bias due to baryons is reduced (e.g., ), but at the cost of a reduction in the statistical precision (). Quantitatively, we find with θ_{ap} = 20′ that the bias is reduced by a factor of almost two and three on S_{8} and Ω_{m} respectively, but with a loss of respectively 12% and 17% on the forecast precision, and a loss of 15% on w_{0}. The two other variations are less efficient but also retain more of the cosmological information. The cut at θ_{in} = 2.4′ decreases the bias by ∼30% at the cost of a 5%, ∼10%, and ∼20% wider statistical precision on S_{8}, w_{0}, and Ω_{m}, respectively. We find similar results for all the Q filter configurations using higher parameter values than the fiducial: the gain in accuracy is always balanced by a significant loss in precision. However, when using smaller values of θ_{ap}, θ_{in}, or x_{c}, the impact of baryons is increased because of the larger contribution of the smallscale baryonic features but the constraining power is also degraded as we chose the fiducial Q filter to maximize the forecast precision in Martinet et al. (2021).
Changes in the measured bias due to baryons and in the associated forecast precision for various mitigation schemes in the case of M_{ap} peaks without tomography.
Overall, none of the smallscale cuts we applied to mitigate baryonic effects are able to decrease the bias below a few percent accuracy whilst preserving strong statistical precision. Using the effective scale of θ_{ap}x_{c} = 16′ recommended in the analysis of Weiss et al. (2019) for Euclidlike mocks, we confirm that the bias becomes consistent with zero. However, such large smoothing scale results in a decrease in the constraining power by a factor of more than 2, 1.5, and 3 on S_{8}, w_{0}, and Ω_{m}, respectively, motivating a full forwardmodel approach of the baryonic bias.
6. Conclusion
In this paper we investigate the impact of baryons on various M_{ap} statistics: peaks, voids, and the lensing PDF. Baryonic physics is modeled with the stateoftheart Magneticum hydrodynamical simulations, and its impact on the data vector is propagated into full cosmological forecasts on S_{8}, w_{0}, and Ω_{m}, for a StageIV lensing survey. The likelihood sampling exploits the cosmological pipeline of Martinet et al. (2021), which is based on the SLICS and cosmoSLICS DMonly simulations. Our results are summarized below:
– Baryons bias the measured M_{ap} estimators by about 5−10% on most S/N bins, notably decreasing the number counts at extreme S/N values, while increasing the number of intermediate S/N features. This is a direct consequence of strong baryonic feedback, which dilutes the density profile of massive halos and decreases their S/N in the M_{ap} map.
– In our StageIV survey setup without tomography, the baryonic feedback propagates into a negative bias of about −3% on S_{8} for every estimator. The bias on Ω_{m} depends on the estimator and ranges from zero for voids to −6% for peaks. When including a tomographic decomposition with five redshift slices including crosstomographic bins, these biases are slightly lowered, likely because of the increased shape noise in each tomographic slice, but remain of the order of a few percent. We observe positive bias in our tomographic setup of the order of 3−5% on w_{0}, although it is not clear at the moment whether this has a physical origin or is caused by parameter degeneracies that are not fully captured by our likelihood.
– Biases on all parameters are increased to ∼5% when combining any M_{ap} statistics with γ2PCF in the tomographic analysis. This combined analysis is maximally affected as the contributions of baryons are in the same direction and of similar amplitude for the individual probes.
– After investigating a range of scale cuts on the M_{ap} statistics, we find that it is difficult to efficiently lower the impact of baryons without significantly degrading the statistical power. In line with Weiss et al. (2019), we find that only an overly large aperture size could lower the bias to subpercent level. A large portion of the signal is lost with such filtering, leading to less competitive cosmological constraints.
We built Magneticum mocks to measure the impact of baryons on peak counts analyses of current Stage III surveys, namely in KiDS450 (Martinet et al. 2018) and DESY1 (HarnoisDéraps et al. 2020), and show that it remains below the statistical uncertainties associated to these surveys. However, this will not be the case for future Stage IV surveys for which baryons need to be accounted for in order to reach percentagelevel precision whilst remaining unbiased.
In this article, we present a correction scheme based on the DV: we measure a corrective factor from hydrodynamical simulations and apply it to mock observations in order to estimate how biased would become cosmological constraints in cases where this step was omitted. However, we note that these results fully depend on the amplitude of the baryonic feedback modeled by the Magneticum simulations, which is still uncertain as seen from the scatter between the different stateoftheart simulations. A more accurate correction would consist in modeling the impact of baryons on M_{ap} statistics from a set of hydrodynamical simulations with various feedback amplitude values, and to marginalize over the extra free parameters when computing the cosmological constraints. Coulton et al. (2020) recently showed the feasibility of this approach using the BAHAMAS simulations run with three different amplitudes of the AGN feedback. The baryonification method described in Schneider & Teyssier (2015) is particularly suited to such an analysis and is therefore a promising tool to design future sets of Nbody simulations that explore both the cosmology dependence and the response to baryons of nonGaussian statistics.
We note that a similar grafting strategy has been used for the highprecision “Clone” data used by the CFHTLenS team (HarnoisDéraps et al. 2012).
Acknowledgments
We thank our KiDS and Euclid Collaborators for useful discussions. NM acknowledges support from a fellowship of the Centre National d’Etudes Spatiales (CNES). TC is supported by the INFN INDARK PD51 grant and by the PRINMIUR 2015W7KAWC grant. JHD acknowledges support from an STFC Ernest Rutherford Fellowship (project reference ST/S004858/1). CG acknowledges support from PRIN MIUR 2017 WSCC32 “Zooming into dark matter and protogalaxies with massive lensing clusters”. KD acknowledges support by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under Germany’s Excellence Strategy – EXC2094 – 390783311. Magneticum has been run using the “LeibnizRechenzentrum” with CPU time assigned to the Project ‘pr86re’ and ‘pr83li’. The SLICS numerical simulations can be found at http://slics.roe.ac.uk/, while the cosmoSLICS can be made available upon request. Computations for these Nbody simulations were enabled by Compute Ontario (www.computeontario.ca), Westgrid (www.westgrid.ca) and Compute Canada (www.computecanada.ca). This work has been carried out thanks to the support of the OCEVU Labex (ANR11LABX0060) and of the Excellence Initiative of AixMarseille University – A*MIDEX, part of the French “Investissements d’Avenir” program. All authors contributed to the development of this paper. NM (lead) conducted the analysis. TC and JHD (both coleads) equally participated by producing the Magneticum and the SLICS and cosmoSLICS mocks respectively. EJ contributed to the analysis pipeline while CG and KD prepared the initial configuration of the Magneticum mocks.
References
 Ajani, V., Peel, A., Pettorino, V., et al. 2020, Phys. Rev. D, 102, 103531 [CrossRef] [Google Scholar]
 Aricò, G., Angulo, R. E., HernándezMonteagudo, C., et al. 2020, MNRAS, 495, 4800 [Google Scholar]
 Asgari, M., Lin, C.A., Joachimi, B., et al. 2021, A&A, 645, A104 [CrossRef] [EDP Sciences] [Google Scholar]
 Barthelemy, A., Codis, S., & Bernardeau, F. 2021, MNRAS, 503, 5204 [Google Scholar]
 Beck, A. M., Murante, G., Arth, A., et al. 2016, MNRAS, 455, 2110 [Google Scholar]
 Bernardeau, F., Nishimichi, T., & Taruya, A. 2014, MNRAS, 445, 1526 [Google Scholar]
 Biffi, V., Dolag, K., & Böhringer, H. 2013, MNRAS, 428, 1395 [Google Scholar]
 Bocquet, S., Saro, A., Dolag, K., & Mohr, J. J. 2016, MNRAS, 456, 2361 [Google Scholar]
 Castro, T., Quartin, M., Giocoli, C., Borgani, S., & Dolag, K. 2018, MNRAS, 478, 1305 [Google Scholar]
 Castro, T., Borgani, S., Dolag, K., et al. 2020, MNRAS, 500, 2316 [Google Scholar]
 Chisari, N. E., Mead, A. J., Joudaki, S., et al. 2019, Open J. Astrophys., 2, 4 [Google Scholar]
 Coulton, W. R., Liu, J., McCarthy, I. G., & Osato, K. 2020, MNRAS, 495, 2531 [CrossRef] [Google Scholar]
 Dietrich, J. P., & Hartlap, J. 2010, MNRAS, 402, 1049 [NASA ADS] [CrossRef] [Google Scholar]
 Dolag, K., Gaensler, B., Beck, A., & Beck, M. 2015, MNRAS, 451, 4277 [Google Scholar]
 Dolag, K., Komatsu, E., & Sunyaev, R. 2016, MNRAS, 463, 1797 [Google Scholar]
 Dubois, Y., Pichon, C., Welker, C., et al. 2014, MNRAS, 444, 1453 [NASA ADS] [CrossRef] [Google Scholar]
 Euclid Collaboration (Knabenhans, M., et al.) 2020, MNRAS, submitted [arXiv:2010.11288] [Google Scholar]
 Fabjan, D., Borgani, S., Tornatore, L., et al. 2010, MNRAS, 401, 1670 [NASA ADS] [CrossRef] [Google Scholar]
 Fan, Z., Shan, H., & Liu, J. 2010, ApJ, 719, 1408 [NASA ADS] [CrossRef] [Google Scholar]
 Ferland, G. J., Korista, K. T., Verner, D. A., et al. 1998, PASP, 110, 761 [Google Scholar]
 Fong, M., Choi, M., Catlett, V., et al. 2019, MNRAS, 488, 3340 [Google Scholar]
 Giocoli, C., Metcalf, R. B., Baldi, M., et al. 2015, MNRAS, 452, 2757 [Google Scholar]
 Giocoli, C., Moscardini, L., Baldi, M., Meneghetti, M., & Metcalf, R. B. 2018a, MNRAS, 478, 5436 [NASA ADS] [CrossRef] [Google Scholar]
 Giocoli, C., Baldi, M., & Moscardini, L. 2018b, MNRAS, 481, 2813 [Google Scholar]
 Gupta, N., Saro, A., Mohr, J., Dolag, K., & Liu, J. 2017, MNRAS, 469, 3069 [Google Scholar]
 HarnoisDéraps, J., & van Waerbeke, L. 2015, MNRAS, 450, 2857 [Google Scholar]
 HarnoisDéraps, J., Vafaei, S., & Van Waerbeke, L. 2012, MNRAS, 426, 1262 [Google Scholar]
 HarnoisDéraps, J., van Waerbeke, L., Viola, M., & Heymans, C. 2015, MNRAS, 450, 1212 [NASA ADS] [CrossRef] [Google Scholar]
 HarnoisDéraps, J., Amon, A., Choi, A., et al. 2018, MNRAS, 481, 1337 [NASA ADS] [CrossRef] [Google Scholar]
 HarnoisDéraps, J., Giblin, B., & Joachimi, B. 2019, A&A, 631, A160 [CrossRef] [EDP Sciences] [Google Scholar]
 HarnoisDéraps, J., Martinet, N., & Castro, T. 2020, ArXiv eprints [arXiv:2012.02777] [Google Scholar]
 Hikage, C., Oguri, M., Hamana, T., et al. 2019, PASJ, 71, 43 [NASA ADS] [CrossRef] [Google Scholar]
 Hilbert, S., Barreira, A., Fabbian, G., et al. 2020, MNRAS, 493, 305 [Google Scholar]
 Hildebrandt, H., Viola, M., Heymans, C., et al. 2017, MNRAS, 465, 1454 [Google Scholar]
 Hirschmann, M., Dolag, K., Saro, A., et al. 2014, MNRAS, 442, 2304 [Google Scholar]
 Huang, H.J., Eifler, T., Mandelbaum, R., et al. 2021, MNRAS, 502, 6010 [Google Scholar]
 Joudaki, S., Hildebrandt, H., Traykova, D., et al. 2020, A&A, 638, L1 [CrossRef] [EDP Sciences] [Google Scholar]
 Kacprzak, T., Kirk, D., Friedrich, O., et al. 2016, MNRAS, 463, 3653 [NASA ADS] [CrossRef] [Google Scholar]
 Kaiser, N., & Squires, G. 1993, ApJ, 404, 441 [NASA ADS] [CrossRef] [Google Scholar]
 Kilbinger, M., Bonnett, C., & Coupon, J. 2014, Astrophysics Source Code Library [record ascl:1402.026] [Google Scholar]
 Laureijs, R., Amiaux, J., Arduini, S., et al. 2011, ArXiv eprints [arXiv:1110.3193] [Google Scholar]
 Li, Z., Liu, J., Matilla, J. M. Z., & Coulton, W. R. 2019, Phys. Rev. D, 99 [Google Scholar]
 Lin, C.A., & Kilbinger, M. 2015, A&A, 576, A24 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Liu, J., Bird, S., Zorrilla Matilla, J. M., et al. 2018, JCAP, 2018, 049 [Google Scholar]
 Martinet, N., Schneider, P., Hildebrandt, H., et al. 2018, MNRAS, 474, 712 [NASA ADS] [CrossRef] [Google Scholar]
 Martinet, N., HarnoisDéraps, J., Jullo, E., & Schneider, P. 2021, A&A, 646, A62 [EDP Sciences] [Google Scholar]
 McCarthy, I. G., Schaye, J., Bird, S., & Le Brun, A. M. C. 2017, MNRAS, 465, 2936 [NASA ADS] [CrossRef] [Google Scholar]
 Mead, A. J., Brieden, S., Tröster, T., & Heymans, C. 2021, MNRAS, 502, 1401 [Google Scholar]
 Navarro, J. F., Frenk, C. S., & White, S. D. M. 1997, ApJ, 490, 493 [NASA ADS] [CrossRef] [Google Scholar]
 Osato, K., Shirasaki, M., & Yoshida, N. 2015, ApJ, 806, 186 [Google Scholar]
 Osato, K., Liu, J., & Haiman, Z. 2021, MNRAS, 502, 5593 [Google Scholar]
 Pillepich, A., Springel, V., Nelson, D., et al. 2018, MNRAS, 473, 4077 [Google Scholar]
 Planck Collaboration XIII. 2016, A&A, 594, A13 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Remus, R.S., Dolag, K., & Hoffmann, T. L. 2017a, Galaxies, 5, 49 [Google Scholar]
 Remus, R.S., Dolag, K., Naab, T., et al. 2017b, MNRAS, 464, 3742 [Google Scholar]
 Saro, A., Liu, J., Mohr, J. J., et al. 2014, MNRAS, 440, 2610 [Google Scholar]
 Schaye, J., Dalla Vecchia, C., Booth, C. M., et al. 2010, MNRAS, 402, 1536 [Google Scholar]
 Schaye, J., Crain, R. A., Bower, R. G., et al. 2015, MNRAS, 446, 521 [Google Scholar]
 Schirmer, M., Erben, T., Hetterscheidt, M., & Schneider, P. 2007, A&A, 462, 875 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Schneider, P. 1996, MNRAS, 283, 837 [Google Scholar]
 Schneider, A., & Teyssier, R. 2015, JCAP, 2015, 049 [Google Scholar]
 Schneider, P., van Waerbeke, L., Kilbinger, M., & Mellier, Y. 2002, A&A, 396, 1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Sellentin, E., & Heavens, A. F. 2016, MNRAS, 456, L132 [Google Scholar]
 Shan, H., Liu, X., Hildebrandt, H., et al. 2018, MNRAS, 474, 1116 [Google Scholar]
 Springel, V. 2005, MNRAS, 364, 1105 [Google Scholar]
 Springel, V., & Hernquist, L. 2003, MNRAS, 339, 289 [Google Scholar]
 Springel, V., Di Matteo, T., & Hernquist, L. 2005, MNRAS, 361, 776 [Google Scholar]
 Springel, V., Pakmor, R., Zier, O., & Reinecke, M. 2020, MNRAS, submitted [arXiv:2010.03567] [Google Scholar]
 Steinborn, L. K., Dolag, K., Hirschmann, M., Prieto, M. A., & Remus, R.S. 2015, MNRAS, 448, 1504 [Google Scholar]
 Steinborn, L. K., Dolag, K., Comerford, J. M., et al. 2016, MNRAS, 458, 1013 [Google Scholar]
 Takahashi, R., Hamana, T., Shirasaki, M., et al. 2017, ApJ, 850, 24 [Google Scholar]
 Taylor, P. L., Bernardeau, F., & Kitching, T. D. 2018, Phys. Rev. D, 98, 083514 [Google Scholar]
 Taylor, P. L., Bernardeau, F., & Huff, E. 2021, Phys. Rev. D, 103, 043531 [Google Scholar]
 Teklu, A. F., Remus, R.S., Dolag, K., et al. 2015, ApJ, 812, 29 [Google Scholar]
 Tornatore, L., Borgani, S., Dolag, K., & Matteucci, F. 2007, MNRAS, 382, 1050 [NASA ADS] [CrossRef] [Google Scholar]
 Troxel, M. A., MacCrann, N., Zuntz, J., et al. 2018, Phys. Rev. D, 98, 043528 [Google Scholar]
 van Daalen, M. P., McCarthy, I. G., & Schaye, J. 2020, MNRAS, 491, 2424 [Google Scholar]
 Vogelsberger, M., Genel, S., Springel, V., et al. 2014, MNRAS, 444, 1518 [NASA ADS] [CrossRef] [Google Scholar]
 Weiss, A. J., Schneider, A., Sgier, R., et al. 2019, JCAP, 2019, 011 [Google Scholar]
 Wiersma, R. P. C., Schaye, J., & Smith, B. D. 2009, MNRAS, 393, 99 [Google Scholar]
 Yang, X., Kratochvil, J. M., Huffenberger, K., Haiman, Z., & May, M. 2013, Phys. Rev. D, 87, 023511 [Google Scholar]
 Zhang, C., Churazov, E., Dolag, K., Forman, W. R., & Zhuravleva, I. 2020, MNRAS, 494, 4539 [Google Scholar]
 Zorrilla Matilla, J. M., Waterval, S., & Haiman, Z. 2020, AJ, 159, 284 [Google Scholar]
 Zürcher, D., Fluri, J., Sgier, R., Kacprzak, T., & Refregier, A. 2021, JCAP, 2021, 028 [Google Scholar]
Appendix A: Correcting for baryons in KiDS and DES
In this section we build Magneticum mocks for current Stage III surveys to investigate the effect of baryons on peak counts in recent analyses. We first revisit the results of Martinet et al. (2018) on the KiDS450 data, and then present a comparison with those obtained for the DESY1 analysis of HarnoisDéraps et al. (2020).
A.1. Estimation of the bias for KiDS450 (Martinet et al. 2018)
We follow the same procedure as in Martinet et al. (2018) to create KiDS mocks, but now we additionally include baryonic physics and examine the impact on the inferred cosmology. In short, we use the KiDS450 redshift distribution calibrated with the direct spectroscopic method (DIR; Hildebrandt et al. 2017) for a single tomographic slice including all galaxies with photometric redshifts 0.1 < z < 0.9 (no tomography is applied). The galaxy density of the mocks is approximately 7.5 arcmin^{−2} and we tile the full 450 deg^{2} by repeating our 100 deg^{2} simulation mocks. In this process we set the positions and the amplitude of the intrinsic ellipticities of the mock galaxies to that of the observed data. We use independent simulated shear values from the ten lines of sight reconstructed in the Magneticum and apply the same five different random rotations of the intrinsic ellipticities as in the KiDS mocks used to compute the model. We build M_{ap} maps with the same Schirmer et al. (2007) filter with an aperture size θ_{ap} = 12.5′. The effect of baryons is measured as the ratio between the mean distribution of peaks in the M_{ap} maps built from the 50 hydrodynamical and DMonly mocks, in 12 bins of S/N ranging from 0 to 4. It is then infused as a multiplicative factor applied to the observed data. In this approach we aim to remove the effect of baryons from the data rather than modifying the model (both methods are equivalent). The correction is then propagated to the cosmological constraints using the same pipeline as in Martinet et al. (2018), which we recall is built from the 157 Dietrich & Hartlap (2010)Nbody simulations that pave the σ_{8} − Ω_{m} plane, plus 35 at the fiducial cosmology for the covariance matrix estimation. When including various lightcones and shape noise realizations, the model is constructed from a total of 3925 mocks, with an additional 175 pseudoindependent mocks for the covariance matrix.
The effect of baryons on the KiDS450 peak counts analysis is shown in Fig. A.1. The blue contours correspond to the 1 and 2σ constraints presented in Fig. 7 of Martinet et al. (2018) and which do not include systematic errors. The black lines show the constraints including the effect of baryons. We see a positive shift in both σ_{8} and Ω_{m} resulting in a change of the structure growth parameter ΔS_{8} = 0.021 (2.8%). We note that the bias is positive in this case: given a fixed observed DV which already includes baryons, including them in the model results in a higher S_{8} cosmology, i.e. the number of large S/N peaks in the model is reduced because of the baryons. This is in contrast with the simulationbased approach in the rest of the article where we estimate the bias from infusing baryons to the DMonly observation. The effect in KiDS is lower than the 3.4% measured in the Euclidlike mocks, likely because of the larger noise with the lower galaxy density of KiDS. It is nevertheless larger than the 2.3% correction used in Martinet et al. (2018) and derived from the simulations of Osato et al. (2015), as expected from the weaker AGN feedback implemented in the latter study. When compared to the statistical precision of the KiDS450 analysis, the S_{8} shift due to baryons remains small with a value of 0.27σ. With this updated baryon model, we revise the S_{8} constraints from Martinet et al. (2018), including the other sources of systematic error described therein: multiplicative shear bias, mean photometric redshift, intrinsic alignment, and shearposition coupling. As the impact of baryons is now fully modeled, and because the constraining power is identical with or without baryons, we no longer need to inflate the uncertainty on S_{8}. We find , a result slightly closer to the Planck estimate (Planck Collaboration XIII 2016, the red contours in Fig. A.1) than the previously reported.
Fig. A.1.
Effect of baryons on the KiDS450 M_{ap} peaks cosmological constraints of Martinet et al. (2018). Blue and black contours correspond to the 1 and 2σ constraints without and with accounting for baryons respectively. Green and red contours represent the best KiDS450 tomographic γ2PCF constraints (Hildebrandt et al. 2017) and Planck cosmic microwave background constraints (Planck Collaboration XIII 2016) that were available at the time. 
A.2. Estimation of the bias for DESY1 (HarnoisDéraps et al. 2020)
The Magneticum simulations have recently been used in HarnoisDéraps et al. (2020) to investigate the impact of baryons on the peak statistics analysis of the DESY1 data. These mocks use the DIRcalibrated redshift distributions used in the cosmic shear reanalysis by Joudaki et al. (2020), carried out in four photometric redshift bins between 0.2 and 1.3. We review the measurement of the baryons bias here and compare the results to those obtained in the previous sections.
As described in HarnoisDéraps et al. (2020), we tile the full DESY1 survey (1321 deg^{2}) with our 100 deg^{2}Magneticum mocks and fix the galaxy positions and intrinsic ellipticity amplitudes to that of the data. We create 100 mock surveys for the DMonly and the hydrodynamical Magneticum simulations to lower the sample variance and shape noise. In this analysis, the model and the covariance matrix are evaluated from the same cosmoSLICS and SLICS Nbody simulations used in the main part of this article respectively, improving in accuracy from the Dietrich & Hartlap (2010)Nbody simulations used in the KiDS450 peak count analysis. The M_{ap} are computed with an aperture size of θ_{ap} = 12.5′ in the four autotomographic bins and in the crossbins. The DV is the concatenation of the peak distributions in 12 bins between 0 < S/N < 4 for all tomographic configurations. A total of 2600 mock surveys are used for the wCDM peak count model, and 1240 for the covariance matrix.
The impact of baryons on the tomographic constraints from M_{ap} peak statistics in DESY1 are shown in Fig. A.2. The blue contours correspond to the 1 and 2σ constraints on the S_{8} and Ω_{m} parameters after marginalisation over the photometric redshift and the shear calibration uncertainties, while the black contours further include the correction due to baryonic physics. As expected, we find again a positive shift towards larger values of S_{8} and Ω_{m}. Quantitatively, the shift is of ΔS_{8} = 0.013 (1.8%), lower than the 2.9% found with the Euclidlike mocks, again due to the lower galaxy density in the DESY1 data, which is ∼5 arcmin^{−2}. In terms of statistical precision, this bias corresponds to a 0.32σ shift and can be safely ignored in the DESY1 analysis. HarnoisDéraps et al. (2020) explored other sources of systematics in their analysis, and show in their Fig. 17 that baryons and intrinsic alignments are the two most important effects, while sourcelens clustering and uncertainty in the nonlinear growth of structure could be safely ignored.
Fig. A.2.
Effect of baryons on the DESY1 M_{ap} peaks cosmological constraints of HarnoisDéraps et al. (2020). Blue and black contours correspond to the 1 and 2σ constraints without and with accounting for baryons respectively. 
All Tables
Changes in the measured bias due to baryons and in the associated forecast precision for various mitigation schemes in the case of M_{ap} peaks without tomography.
All Figures
Fig. 1.
Ratio of power spectra at z = 0 between various hydrodynamical simulations and their corresponding DMonly run. The Magneticum used in this analysis is representative of an average behavior with a loss of power of about 15% at k ∼ 10 Mpc h^{−1} due to AGN feedback. 

In the text 
Fig. 2.
Relative change in the γ2PCF due to baryons. The orange and purple curves represent the ξ_{+} and ξ_{−} estimators, respectively. The error bars correspond to the diagonal elements of the SLICS covariance matrix rescaled to 15 000 deg^{2}. 

In the text 
Fig. 3.
Relative change in the M_{ap} estimators due to baryons. The blue, green, and red curves represent void counts, peak counts, and the lensing PDF, respectively. The error bars correspond to the diagonal elements of the SLICS covariance matrix rescaled to 15 000 deg^{2}. 

In the text 
Fig. 4.
Forecast of cosmological parameters from peak counts in a 100 deg^{2} survey at Euclid depth without tomography. Marginalized 2D (1 and 2σ contours) and 1D (full likelihood) constraints are displayed in blue for the DMonly case and in green when including baryon physics. Dashed lines correspond to the 1σ constraints in the 1D marginalized likelihood. The blue and green crosses indicate the best estimate in the 2D constraints, while the red crosses and lines indicate the input cosmological parameter values. Gray crosses correspond to parameters of the 25 cosmoSLICS simulations that are used to estimate the cosmology dependence of the number of peaks. 

In the text 
Fig. 5.
Same as Fig. 4 but for a tomographic analysis with five redshift slices and including auto and crossM_{ap} terms between redshift slices. 

In the text 
Fig. 6.
Same as Fig. 4 but for the combination of peak counts and γ2PCF for a tomographic analysis with five redshift slices and including auto and crossterms between redshift slices. 

In the text 
Fig. 7.
Top: profiles of the Q filters defined via Eq. (4) and used to mitigate the impact of baryons by reducing the weight of small scales. Bottom: relative change in the M_{ap} peak counts due to baryons for the various mitigation setups. 

In the text 
Fig. A.1.
Effect of baryons on the KiDS450 M_{ap} peaks cosmological constraints of Martinet et al. (2018). Blue and black contours correspond to the 1 and 2σ constraints without and with accounting for baryons respectively. Green and red contours represent the best KiDS450 tomographic γ2PCF constraints (Hildebrandt et al. 2017) and Planck cosmic microwave background constraints (Planck Collaboration XIII 2016) that were available at the time. 

In the text 
Fig. A.2.
Effect of baryons on the DESY1 M_{ap} peaks cosmological constraints of HarnoisDéraps et al. (2020). Blue and black contours correspond to the 1 and 2σ constraints without and with accounting for baryons respectively. 

In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.