Issue 
A&A
Volume 683, March 2024



Article Number  A103  
Number of page(s)  18  
Section  Cosmology (including clusters of galaxies)  
DOI  https://doi.org/10.1051/00046361/202347986  
Published online  12 March 2024 
KiDS1000 cosmology: Combined second and thirdorder shear statistics
^{1}
University of Bonn, ArgelanderInstitut für Astronomie, Auf dem Hügel 71, 53121 Bonn, Germany
email: pburger@astro.unibonn.de
^{2}
Universität Innsbruck, Institut für Astro und Teilchenphysik, Technikerstr. 25/8, 6020 Innsbruck, Austria
^{3}
Department of Astronomy and Astrophysics, University of California, Santa Cruz, 1156 High Street, Santa Cruz, CA 95064, USA
^{4}
E.A Milne Centre, University of Hull, Cottingham Road, Hull HU6 7RX, UK
^{5}
INAF – Osservatorio Astronomico di Trieste, via Tiepolo 11, 34131 Trieste, Italy
^{6}
INFN – Sezione di Trieste, 34100 Trieste, Italy
^{7}
IFPU – Institute for Fundamental Physics of the Universe, via Beirut 2, 34151 Trieste, Italy
^{8}
UniversitätsSternwarte, Fakultät für Physik, LudwigMaximiliansUniversität München, Scheinerstr.1, 81679 München, Germany
^{9}
MaxPlanckInstitut für Astrophysik, KarlSchwarzschildStraße 1, 85741 Garching, Germany
^{10}
School of Mathematics, Statistics and Physics, Newcastle University, Newcastle upon Tyne NE1 7RU, UK
^{11}
Ruhr University Bochum, Faculty of Physics and Astronomy, Astronomical Institute (AIRUB), German Centre for Cosmological Lensing, 44780 Bochum, Germany
^{12}
Leiden Observatory, Leiden University, PO Box 9513 2300 RA Leiden, The Netherlands
^{13}
AixMarseille Université, CNRS, CNES, LAM, Marseille, France
Received:
15
September
2023
Accepted:
6
December
2023
Aims. In this work, we perform the first cosmological parameter analysis of the fourth release of Kilo Degree Survey (KiDS1000) data with second and thirdorder shear statistics. This paper builds on a series of studies aimed at describing the roadmap to thirdorder shear statistics.
Methods. We derived and tested a combined model of the secondorder shear statistic, namely, the COSEBIs and the thirdorder aperture mass statistics 〈ℳ_{ap}^{3}〉 in a tomographic setup. We validated our pipeline with Nbody mock simulations of the KiDS1000 data release. To model the second and thirdorder statistics, we used the latest version of HMCODE2020 for the power spectrum and BIHALOFIT for the bispectrum. Furthermore, we used an analytic description to model intrinsic alignments and hydrodynamical simulations to model the effect of baryonic feedback processes. Lastly, we decreased the dimension of the data vector significantly by considering only equal smoothing radii for the 〈ℳ_{ap}^{3}〉 part of the data vector. This makes it possible to carry out a data analysis of the KiDS1000 data release using a combined analysis of COSEBIs and thirdorder shear statistics.
Results. We first validated the accuracy of our modelling by analysing a noisefree mock data vector, assuming the KiDS1000 error budget, finding a shift in the maximum of the posterior distribution of the matter density parameter, ΔΩ_{m} < 0.02 σ_{Ωm}, and of the structure growth parameter, ΔS_{8} < 0.05 σ_{S8}. Lastly, we performed the first KiDS1000 cosmological analysis using a combined analysis of second and thirdorder shear statistics, where we constrained Ω_{m} = 0.248_{−0.055}^{+0.062} and S_{8} = σ_{8}√(Ω_{m}/0.3 )= 0.772 ± 0.022. The geometric average on the errors of Ω_{m} and S_{8} of the combined statistics decreases, compared to the secondorder statistic, by a factor of 2.2.
Key words: gravitation / gravitational lensing: weak / methods: analytical / methods: numerical / cosmological parameters / largescale structure of Universe
© The Authors 2024
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.
This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.
1. Introduction
Gravitational lensing describes the deflection of light by massive objects. It is sensitive to baryonic and dark matter and, therefore, ideal for probing the total matter distribution in the Universe. Since the distribution of matter is highly sensitive to cosmological parameters, it is excellent to test and probe the standard model of cosmology, called the Λ cold dark matter model (ΛCDM). Although the ΛCDM model can describe observations of the early Universe, such as the cosmic microwave background (CMB; e.g., Planck Collaboration VI 2020), or the Local Universe, such as the observed largescale structure (LSS) of matter and galaxies (Sánchez et al. 2017), with remarkable accuracy, it is being put under stress due to tension observed between early and local probes. A ∼2σ tension is in the structure growth parameter , where σ_{8} is the normalisation of the power spectrum and Ω_{m} is the total matter density parameter (Hildebrandt et al. 2017; Planck Collaboration VI 2020; Joudaki et al. 2020; Heymans et al. 2021; DES Collaboration 2022; Di Valentino et al. 2021; Dalal et al. 2023), suggests that the Local Universe is less clustered than what is expected from earlytime measurements when extrapolated under the ΛCDM model. If these tensions are not due to systematics, extensions to the ΛCDM model are necessary, which will be tested with the next generation of cosmic shear surveys such as Euclid (Laureijs et al. 2011) or the Vera Rubin Observatory Legacy Survey of Space and Time (LSST, Ivezic et al. 2008).
Commonly, twopoint statistics for weak lensing and galaxy positions are used to infer cosmological parameters since they can be modelled accurately and systematic inaccuracies are well understood (Schneider et al. 1998; Troxel et al. 2018; Hildebrandt et al. 2017, 2020; Hikage et al. 2019; Asgari et al. 2019, 2020). Twopoint statistics are excellent for capturing the entire information content of a Gaussian random field. However, nonlinear gravitational instabilities created a significant amount of nonGaussian features during the evolution of the Universe, such that the local matter distribution departed strongly from a Gaussian field. Therefore, higher order statistics are needed to extract all the available information in the local LSS of matter and galaxies. Furthermore, such higher order statistics usually depend differently on cosmological parameters and systematic effects such as intrinsic alignment (IA), meaning that a joint investigation of second and higher order statistics tightens cosmological parameter constraints (see, e.g., Kilbinger & Schneider 2005; Bergé et al. 2010; Pires et al. 2012; Fu et al. 2014; Pyne & Joachimi 2021). Recently used examples of higher order statistics are the peak count statistics (Martinet et al. 2018; HarnoisDéraps et al. 2021), persistent homology (Heydenreich et al. 2021), density split statistics (Gruen et al. 2018; Burger et al. 2023), and the integrated threepoint correlation function used in (Halder et al. 2021; Halder & Barreira 2022), along with a second and thirdorder convergence moment analysis (Gatti et al. 2022).
This work considers second and thirdorder shear statistics, where the former probes the variance and the latter the skewness of the LSS at various scales. Our chosen secondorder statistic is the E_{n}modes of the COSEBIs (Schneider et al. 2010), and the thirdorder statistic is described in Schneider et al. (2005) and recently measured in the Dark Energy Survey Year (DES) 3 Results (Secco et al. 2022). Furthermore, this work belongs to a series of papers that aim for cosmological parameter analyses using thirdorder shear statistics. In the first paper, Heydenreich et al. (2023, hereafter H23), we validated the analytical fitting formulae for a nontomographic analysis and the conversion from threepoint correlation function (3PCF) of cosmic shear to . We found that , even though they are combined from the shear 3PCF at different scales, contain a similar amount of information on Ω_{m} and S_{8} as the 3PCF itself. The fact that can be measured from 3PCF is very convenient since this allows unbiased estimates for any survey geometry. A fast computation method of the aperture mass statistics by measuring the shear 3PCF is tackled in Porth et al. (2023). Lastly, adding thirdorder statistics increases the dimension of the data vector significantly, and an analytical expression for the covariance is preferred, which is derived and validated for in a nontomographic setup in Linke et al. (2023). However, as this analysis considers combining secondorder statistics with in a tomographic setup, and a joint covariance matrix has not been derived yet (Wielders et al., in prep.), we must still rely on numerical simulations to determine the covariance matrix.
This article presents the first cosmological parameter analysis using the fourth data release of KiDS (KiDS1000), combining second and thirdorder shear statistics. We show the cosmological results, preceded by validation of several extensions of the analytical model to allow for a tomographic analysis and to include astrophysical effects such as IA (e.g., Joachimi et al. 2015) or baryonic feedback processes (Chisari et al. 2015). We aim to find the smallest set of scales that retain most of the cosmological information.
The paper is structured as follows: In Sect. 2, we review the basics of the second and thirdorder shear statistics, extend the modelling to a tomographic analysis and describe our method to model the IA analytically. In Sect. 4, we describe the KiDS1000 data and in Sect. 5, we introduce the simulation data used to validate our model and to correct it for baryonic feedback processes. In Sect. 6, we briefly review our method to perform cosmological parameter interference. In Sects. 7 and 8, we validate our analysis pipeline against several systematics. The final cosmological results are presented in Sect. 9, and we present our conclusions in Sect. 10.
2. Theoretical background
This section offers a review of the basics of weak gravitational lensing formalism and aperture statistics. For more detailed reviews, we refer to Bartelmann & Schneider (2001, Hoekstra & Jain (2008), Munshi et al. (2008), Bartelmann (2010), and Kilbinger (2015). In this work, we assume a spatially flat universe, such that the comoving angulardiameter distance is expressed as f_{K}[χ(z)] = χ(z), where χ(z) is the comoving distance at redshift z. Given the matter density, ρ(x, z), at comoving position, x, and redshift, z, the density contrast is , where is the average matter density at redshift, z. The dimensionless surface mass density or convergence, κ, for sources at redshift, z, is determined by the lineofsight integration as
where ϑ is the angular position on the sky, H_{0} is the Hubble constant, and a is the scale factor. The second argument of δ simultaneously describes the radial direction and the cosmological epoch, related through the lightcone condition c dt=a(z) dχ.
2.1. Limber projections of power and bispectrum
Given the Fourier transform of the matter density contrast, the matter power spectrum, P_{δδ}(k, z), and bispectrum, B_{δδδ}(k_{1}, k_{2}, k_{3}, z), are:
where δ_{D} is the Diracdelta distribution. The statistical isotropy of the Universe implies that the power and bispectrum only depend on the moduli of the kvectors.
The projected power and bispectrum can then be computed using the Limber approximation (Limber 1954; Kaiser & Jaffe 1997; Bernardeau et al. 1997; Schneider et al. 1998; LoVerde & Afshordi 2008),
where ℓ_{3} = ℓ_{1} + ℓ_{2}, and g^{(i)}(χ) denotes the lensing efficiency and is defined as
with n^{(i)}(χ) being the redshift probability distribution of the ith tomographic z_{ph}bin. Since the Limber approximation breaks down for small values of ℓ (Kilbinger et al. 2017), we consider only scales below 4° to model . In principle, we could have also used the modified Limber approximation and shift ℓ_{1, 2, 3} → ℓ_{1, 2, 3} + 1/2, but since we found a difference at maximum for the largest filter radii of 1.5%, we neglected it here. To model the nonlinear matter power spectrum P_{δδ}, we use the revised HMCODE2020 model of Mead et al. (2021) and for the matter bispectrum B_{δδδ} the BIHALOFIT of Takahashi et al. (2020).
2.2. Nonlinear alignment model
The impact of galaxy IA is a known contaminating signal to the cosmic shear measurements that must be accounted for in all weak lensing studies. To model the effects of intrinsic alignment, we use the nonlinear alignment (NLA) model (Bridle & King 2007), which is a oneparameter model described as
where D_{+} is the linear growth factor at redshift z, Mpc^{3}, as calibrated in Brown et al. (2002), and A_{IA} captures the coupling strength between the matter density and the tidal field.
By considering the tidal alignment field, δ_{I}, as a biased tracer of the matter density contrast field δ, neglecting all higher order bias terms, we get:
With this, we find that:
where P_{δδ}(k, z) is the nonlinear matter power spectrum. The projected power spectra then follow to
and the total projected power spectrum becomes
The term describes the actual lensing signal, which is given in Eq. (4) with the weighting kernel in Eq. (6). The II contribution describes how two galaxies spatially close together tend to be aligned. The term GI(z_{i} > z_{j}) describes the fact that high matter density regions align the lower redshift galaxies but also affect the shear of the background galaxies. While the II term is dominant if galaxies of the same tomographic bin are considered, GI and IG start to dominate if galaxies of separated tomographic bins are considered. With z_{i} < z_{j}, the term GI is expected to vanish for two tomographic bins with no significant redshift overlap.
Following the ansatz for modelling IA in the power spectrum, we get
where B_{δδδ} is the nonlinear matter bispectrum, which we calculated with BIHALOFIT. The projected bispectra are
The total projected bispectrum is:
where the tuple is . The interpretation of all these terms is analogous to the ones from the power spectrum, where is the actual pure lensing signal given by Eq. (5) with the weighting kernel in Eq. (6).
2.3. Aperture mass statistics
One of the major problems of weak lensing mass reconstruction techniques is the masssheet degeneracy (hereafter MSD; Falco et al. 1985; Schneider & Seitz 1995), which corresponds to adding a uniform surface mass density without affecting lensing observables such as the shear. However, it is possible to define quantities invariant under the MSD, one example being the aperture mass statistics (Schneider 1996; Bartelmann & Schneider 2001). Another advantage of aperture mass statistics is that they separate the signal into socalled E and Bmodes (Schneider et al. 2002), where, to leading order, the weak gravitational lensing effect cannot create Bmodes. Lastly, as shown in H23, aperture statistics are an excellent strategy to compress thousands of bins of 3PCF into a few hundred bins.
The aperture mass ℳ_{ap} at position ϑ with filter radius θ_{ap} is defined through the convergence κ, as follows:
where U_{θap}(ϑ′) is a compensated filter such that ∫dϑ′ ϑ′ U_{θap}(ϑ′) = 0. The tangential shear component γ_{t} of the complex shear in Cartesian coordinates γ = γ_{1} + iγ_{2}, is defined as , where ϕ is the polar angle of ϑ′. Given the tangential shear γ_{t}, the aperture mass ℳ_{ap} can also be calculated as
where Q_{θap} is related to U_{θap} via
We define , denote by the Fourier transform of u and use the filter function introduced in Crittenden et al. (2002),
2.4. Modelling aperture mass moments
The expectation value of the aperture mass ⟨ℳ_{ap}⟩(θ_{ap}), which approximates the ensemble average over all positions ϑ, vanishes by construction. However, the secondorder (variance) of the aperture mass is nonzero and can be calculated as
Equivalently, the thirdorder moment of the aperture statistics , can be computed from the convergence bispectrum via (Jarvis et al. 2004; Schneider et al. 2005)
where ℓ_{3} = ℓ_{1} + ℓ_{2}. Later, we differentiate between equal filter radii, for which θ_{ap, 1} = θ_{ap, 2} = θ_{ap, 3} and nonequal filter radii, where the values of θ_{ap, i} are all allowed to vary.
2.5. Modelling COSEBIs
Schneider et al. (2010) introduced the complete orthogonal sets of E/Bintegrals (COSEBIs), which are defined via the twopoint shear statistics ξ_{±}(ϑ) = ⟨γ_{t}γ_{t}⟩(ϑ)±⟨γ_{×}γ_{×}⟩(ϑ) on a finite angular range
where T_{±n}(ϑ) are filter functions with support in [ϑ_{min}, ϑ_{max}] (Schneider et al. 2010). If
where J_{0, 4} are the zeroth and fourth order Bessel function, then the COSEBIs offer the advantage of cleanly separating all welldefined E and Bmodes within the range [ϑ_{min}, ϑ_{max}]. This is not given, for instance, for secondorder aperture mass statistics, as this would require information of ξ_{±} over the full space. Analytically, the COSEBIs can be calculated from the Emode power spectrum defined in Eq. (4) and a Bmode angular power spectra as:
Since Bmodes cannot be created by gravitational lensing directly, we neglected the modelling of B_{n} for this work. Therefore, we focus, henceforth, only on the E_{n} modes from the COSEBIs.
3. Measuring shear statistics
As discussed in H23 can be estimated from data in three ways, which we review below. The first uses the convergence field, κ, the second uses the shear field, and the third uses correlation functions. The latter is used to measure the real data and compute the covariance matrix used for the real data analysis.
3.1. Measuring shear statistics from convergence maps
If the convergence field, κ, is available (e.g., from simulations), the easiest way to measure ℳ_{ap} is to use Eq. (21). If the convergence field does not contain masks, this estimator is also unbiased as long as a border of the size of the filter function is removed from the aperture mass field before calculating the spatial average. However, no boundaries need to be removed for unmasked fullsky convergence fields to get an unbiased estimator. To estimate the aperture mass from (fullsky) convergence maps, the maps are smoothed with the HEALPY^{1} (Zonca et al. 2019) function SMOOTHING. This function needs a beam window function created by the function BEAM2BL, which is determined from the corresponding U_{θap} filter.
In the absence of Bmodes the E_{n} modes can also be calculated by using convolutions as (Schneider et al. 2010):
where
with T_{+n} introduced in Eq. (29).
3.2. Measuring shear statistics from galaxy catalogues
In a real survey, the convergence is not observable and can be inferred only from the measured shear. However, estimating ℳ_{ap} from these reconstructed convergence maps is not accurate as the reconstructed convergence is necessarily smoothed and potentially also contains other systematic effects caused by masks or boundaries of the survey (Seitz & Schneider 1997).
However, as motivated by Eq. (22), the aperture mass can also be estimated from an unmasked shear field (e.g., for simulations). Given a galaxy catalogue, where galaxies are only found at specific positions such that the number of galaxies within an aperture varies on the sky, Eq. (22) needs to get modified to
where the sum over the filter functions serves as a normalisation, ε_{t} are the observed galaxy ellipticities converted into their tangential component, and ϑ_{i} are their respective positions. We note that we sampled the galaxies on a grid using a cloudincell method^{2}, so that convolutions can determine the aperture masses. Since we used this approach solely for finite fields, we applied the same cutoff of to all aperture mass maps, where is the largest filter radius. For both approaches, where the aperture mass is determined, the second and thirdorder aperture statistics follow by multiplying the respective aperture mass maps pixelwise and then taking the average of all pixel values.
3.3. Measuring shear statistics from correlation functions
The first two methods to estimate unbiased aperture statistics are not applicable to observed data with masks. Another disadvantage is the removal of the edges, which can be significant for a complex survey footprint, leading to decreased statistical certainty. The best method to estimate aperture shear statistics or the E_{n} is by measuring them from two and threepoint shear correlation functions (Schneider et al. 2002, 2005; Jarvis et al. 2004). The advantage of correlation functions is that they can be estimated for any survey geometry. The measurement of the twopoint shear correlation functions ξ_{±}(θ_{ap}) is unbiased and easily performed by TREECORR (Jarvis et al. 2004) and the conversion to COSEBIs are given in Eqs. (27) and (28).
For the measurement of the aperture statistics, , we refer to our companion paper, namely, Porth et al. (2023). It describes an efficient estimation procedure of the natural components of the shear 3PCF (Schneider & Lombardi 2003), which is then transformed into aperture statistics. The corresponding equations for a tomographic setup can be found in their Sect. 5.2.
4. Observed data
This analysis explores the fourth data release of KiDS (Kuijken et al. 2015, 2019; de Jong et al. 2015, 2017). The weak lensing data observed with the highquality VSTOmegaCAM camera is public^{3}. It is collectively called ‘KiDS1000’ as it covers ∼1000 deg^{2} of images, which is then reduced to an effective area of 777.4 deg^{2} after masking. The significant advantage, compared to previous weak lensing surveys and data releases, is its overlap with its partner survey VIKING (VISTA Kilodegree Infrared Galaxy survey, Edge et al. 2013), which observes galaxy images at infrared wavelength. Therefore, galaxies were observed in nine optical and nearinfrared bands, u, g, r, i, Z, Y, J, H, K_{s}, allowing for better control over redshift uncertainties (Hildebrandt et al. 2021, hereafter H21).
The KiDS1000 cosmic shear catalogue is divided into five tomographic z_{ph}bins, whose redshifts are calibrated using the selforganising map (SOM) method^{4} described in Wright et al. (2020). The redshift distributions of all five tomographic bins are shown in Fig. 1, and were initially presented in H21. The residual systematic uncertainties on the redshift distributions are listed in Table 1 and are included in this work as nuisance parameters which we marginalise over. We note that they are correlated, which we account for by using their correlation matrix for the marginalisation. The galaxy shear ellipticities and their corresponding weights, w, are estimated by the lensfit tool (Miller et al. 2013; Fenech Conti et al. 2017; Kannawadi et al. 2019) and are described in more detail in Giblin et al. (2021). These shearrelated systematic effects shift the S_{8} parameter by (at most) 0.1σ when measured by cosmic shear twopoint functions. The resulting systematics are stated in Table 1, where we marginalise over the shear multiplicative mbias correction in the resulting posteriors.
Fig. 1. Redshift distribution of the five tomographic z_{ph}bins of the KiDS1000 data. 
Overview of the observational KiDS1000 data.
5. Simulated data
We use several simulated data sets created to resemble the observed KiDS1000 data, to validate our inference pipeline, study the impact of key systematic uncertainties, and forecast the expected KiDS1000 analysis.
In particular, we used the fullsky gravitational lensing simulations described in Takahashi et al. (2017, hereafter T17) to generate data vectors and numerical covariance matrices. The cosmoSLICS+IA simulations, described in HarnoisDéraps et al. (2022), were used to test the modelling of IA; while the Magneticum lensing simulations, first introduced in Hirschmann et al. (2014), were used to infuse Baryon feedback on the model whose strength we regulate with a free parameter in the posterior estimation.
5.1. Takahashi simulations
Since the T17 simulations are used in this series of previous works, we only reiterate the essential details here. The T17 simulations follow the nonlinear evolution of 2048^{3} particles evolved in a large series of nested cosmological volumes with side length starting at L = 450 Mpc h at low redshift, and increasing at higher redshift, resulting in 108 different fullsky realisations. These were produced by the Gadget3 Nbody code (Springel 2005) and are publicly available^{5}. The cosmological parameters of the matter and vacuum energy density are fixed to Ω_{m} = 1 − Ω_{Λ} = 0.279, the baryon density parameter to Ω_{b} = 0.046, the dimensionless Hubble constant to h = 0.7, the normalisation of the power spectrum to σ_{8} = 0.82, and the spectral index to n_{s} = 0.97. The shear information of the T17 simulations is given in terms of 108 γ and κ fullsky realisations, where each realisation is divided into 38 ascending redshift slices. To reproduce the KiDS1000 data for a given tomographic z_{ph}bin shown in Fig. 1, we built the weighted average of the first 30 γ and κ redshift slices. The weight for each redshift slice is measured by integrating the n(z) over the corresponding width of the redshift slice.
We decided on two approaches (convergence maps vs. galaxy ellipticity) to measure the data vectors and covariance matrices from the T17 simulations. While the first approach (Sect. 5.1.1) has the advantage that a large number of mock data is available to measure a reliable covariance matrix even for large data vectors, the second approach (Sect. 5.1.2) has the advantage that the exact galaxy positions are used, which implies that the holes (masks) match the data.
5.1.1. Convergence mock data
The first approach is to convolve the full sky convergence maps with the T_{±, n}(θ_{ap}) or U(θ_{ap}) filters, then multiply with corresponding other smoothed maps and then take the spatial average. For the covariance matrix, we divided the smoothed convergence maps into 48 subpatches, where each patch has an area of 860 deg^{2}. Since the κ maps are first convolved with the filter functions and the subpatches have common edges, the individual patches are not fully independent from each other, which decreases the covariance matrix. However, we found in our testing that selecting only 18 subpatches, that have no common edges, gives almost identical results. Shape noise is added to the convergence maps by drawing random numbers from a Gaussian distribution with a vanishing mean and a standard deviation, as follows:
with the pixel area as A_{pix} = 0.74 arcmin^{2} (NSIDE = 4096), n_{eff} as the effective galaxy number density that included the lensfit weights, and σ_{ϵ} as the shape noise contribution given in Table 1. With 48 subpatches and 108 realisations, the covariance matrix is measured from 5184 mock data. As the actual KiDS1000 area is roughly 777.4 deg^{2}, we rescaled the covariance by 860 deg^{2}/777.4 deg^{2} ≈ 1.1. The reference data vector is measured from one T17 realisation with a resolution of A_{pix} = 0.18 arcmin^{2} (NSIDE = 8192) without shape noise.
5.1.2. Galaxy shear catalogue mock data
For real surveys with a complex topology, the convolved maps give biased results (Seitz & Schneider 1997). Therefore, we decided on our second approach to measure our statistics from second and thirdorder shear correlation functions. For this, we created galaxy shear catalogues from projected γ and κ fields by extracting the shear information only at the true positions of the observed galaxies in KiDS1000. Since the correlation functions need to be measured to very small scales, we used realisations with a pixel resolution of A_{pix} = 0.18 arcmin^{2} (NSIDE = 8192). We checked that this resolution is sufficient by comparing the data vectors obtained from catalogues constructed from the same initial conditions on the reference resolution (NSIDE = 8192) and a higher resolution (NSIDE = 16 384), finding a maximum deviation of 1% for an aperture filter scale of 4′. To increase the number of mock data, we shifted the galaxy positions 18 times by 20° along the lines of constant declination and extracted the shear information at the new positions. Afterwards, the shifted galaxy positions and the corresponding lensing information are shifted back to the original footprint. The back shifting is done only for simplicity and has no physical reason. With this procedure, 1944 almost independent mock data were created, from which the covariance and the reference data vector were measured. To add shape noise, we combined the twocomponent reduced shear g = γ/(1 − κ) of each object with a shape noise contribution, ϵ^{s}, to create observed ellipticities ϵ^{obs} (Seitz & Schneider 1997), as follows:
The quantities here are all complex numbers and the asterisk ‘’ indicates complex conjugation. To mock the ϵ^{s}, we use the observed ellipticities of the KiDS1000 data and randomly rotate them to erase the underlying correlated shear signal. This procedure offers the advantage that the resulting distribution of ϵ^{s} matches the distribution of . The estimated mean and dispersion σ_{ϵ} are given in Table 1. Furthermore, when computing our statistics via correlation functions, we need to consider the corresponding weight w for each shape measurement, which ensures that we use the correct effective number density. Since the lensing weights and the intrinsic ellipticities of source galaxies are correlated, adding the rotated ellipticities to the shear signal preserves this correlation.
5.1.3. Data vector measurement and modelling
We follow the analysis of Asgari et al. (2021, hereafter A21) as it is the fiducial cosmic shear analysis of the KiDS1000^{6}. We only use the first five E_{n}moments determined from twopoint correlation functions, measured from to 300′ in 400 radial bins. The resulting E_{n} are shown in Fig. 2, where we checked that the B_{n} are consistent with zero. The E_{n} data vector for γ^{7} and g is the mean of all 1944 T17 mock data; and for the κ data vector, we used one realisation with resolution A_{pix} = 0.18 arcmin^{2} (NSIDE = 8192). The analytical model was computed using OSMOSIS (Zuntz et al. 2015a).
Fig. 2. Measured and modelled E_{n} vector for the first five moments in the T17 mock data. The green and blue dots are the mean of all 1944 KiDS1000 mock data that are also used to compute the covariance matrix with a resolution A_{pix} = 0.18 arcmin^{2} and shape noise. The red dots represent the data vector measured from one fullsky T17 realisation with a resolution A_{pix} = 0.18 arcmin^{2} and no shape noise. The grey band indicates the expected uncertainty from the KiDS1000 survey. 
The data vectors are measured from the same mock data that are also used to compute the E_{n}. The analytical model pipeline is described in H23. In Fig. 3, we show for aperture filter radii θ_{ap} ∈ {4′,6′,8′,10′,14′,18′,22′,26′,32′,36′}, where we used only equal scales. We show them here for some specific z_{ph}bin combinations. In Fig. B.2, we show also for nonequal aperture filter radii θ_{ap} ∈ {4′,8′,16′,32′}.
Fig. 3. Same as Fig. 2, but here the measured data and model are the vector for equalscale aperture filter radii θ_{ap} ∈ {4′,6′,8′,10′,14′,18′,22′,26′,32′,36′} in the T17 mock data for some selected z_{ph}bin combinations. The full data set is shown in Fig. B.1. 
5.2. CosmoSLICS+IA
The cosmoSLICS+IA simulations are cosmic shear mock galaxy catalogues infused with the nonlinear alignment model of Bridle & King (2007), which is ideally suited for testing and validating our analytical IA model. They are based on Nbody simulations of identical box size and particle density as the SLICS (HarnoisDéraps et al. 2018), which were already used in previous works of this series; we therefore only summarise the essential details.
The cosmology corresponds to the fiducial cosmoSLICS model presented in HarnoisDéraps et al. (2019), with Ω_{m} = 0.2905, Ω_{Λ} = 0.7095, Ω_{b} = 0.0473, h = 0.6898, σ_{8} = 0.836, w_{0} = −1.0 and n_{s} = 0.969. We use the full set of 50 simulated galaxy catalogues, each covering 100 deg^{2} and reproducing the KiDS1000 n(z) (see Fig. 1) and galaxy number densities n_{eff} specified in H21. As we use these mock data only to infuse IA effects into the T17 data vector, shape noise is not included. Following the methods described in HarnoisDéraps et al. (2022), the intrinsic ellipticity components of these galaxies are computed as
where f_{IA} is defined in Eq. (7), and ϕ is the gravitational potential. The partial derivatives of the gravitational potential describe the Cartesian components of the projected tidal field tensors. The terms were then combined with the noisefree cosmic shear signal using Eq. (36), resulting in an IAcontaminated weak lensing sample that is consistent with the NLA model. As these mock data cover a square patch without any masking, we used the methods described in Sect. 3.2 to estimate the aperture statistics. Furthermore, they were used only to validate the IA modelling, which we discuss in Appendix A.
5.3. Magneticum simulations
The feedback processes due to baryonic matter significantly affect the distribution of the LSS, such that the clustering of the matter is reduced on intracluster scales by up to tens of per cent (van Daalen et al. 2011). However, in a quantitative sense, this suppression is not well understood (Chisari et al. 2015). We use the Magneticum simulations (Hirschmann et al. 2014) to investigate the impact of baryonic feedback processes. Magneticum was run using GADGET3 code which is a more efficient version of Gadget 2 (Springel 2005) that includes modern smoothed particle hydrodynamics (Beck et al. 2016). The dark matter particle mass is 6.9 × 10^{8} h^{−1} M_{⊙} and gas particle mass 1.4 × 10^{8} h^{−1} M_{⊙}. The underlying matter fields were constructed from the Magneticum Pathfinder simulations^{8} in the Run2 with a comoving volume of side 352 h^{−1} Mpc and Run2b with a comoving volume of side 640 h^{−1} Mpc, and are described in Hirschmann et al. (2014) and Ragagnin et al. (2017), respectively. The Magneticum simulations account for radiative cooling, star formation, supernovae and active galactic nuclei.
The Magneticum shear maps used in this work were first presented in Castro et al. (2018). Although the underlying creation of these shear maps is similar to the production of the cosmoSLICS, meaning that the mock data follow the same redshift distribution and number density as the KiDS1000 data, there are three main differences. First, only ten instead of 50 pseudoindependent light cones are available. Second, the galaxies are placed at the exact galaxy positions as in the observed data rather than randomly. Third, as these mock data are flat sky simulations and cover only an area of 100 deg^{2}, the KiDS1000 footprint is divided into 18 regions. This results in 18 catalogues for each of the ten pseudoindependent light cones. Due to the nontrivial geometry of the masks, we used correlation functions to estimate the shear statistics.
To incorporate baryonic feedback processes into our modelling pipeline using the Magneticum simulations, we calculated the data vector, d^{DM}, using dark matter only and the data vector, d^{DM + BA}, where dark matter and baryons are included. With them we modified the model vector, m, to
where we introduce the parameter A_{ba}. A value of zero for A_{ba} means no baryonic feedback, and a value of one is exactly the strength of baryonic feedback as in the Magneticum. Furthermore, we interpolate between zero and one and extrapolate to five, which we chose arbitrarily. We note that for the joint analysis of second and thirdorder statistics, we need at least extrapolation to A_{ba} = 2 (at 95% confidence), which approximately agrees with the baryonic strength of other hydrodynamical simulations (Martinet et al. 2021). The thirdorder analysis alone is bound by the prior A_{ba} = 5. Although extrapolating to these large values is probably incorrect, we introduced A_{ba} only to give the model some flexibility to account for baryonic feedback effects without having a direct physical meaning such as the ejected gas.
6. Cosmological parameter inference methodology
In the following sections, we determine multiple posterior distributions using Markov chain Monte Carlo (MCMC) samplings for different model ingredients.
We varied the two cosmological parameters Ω_{m} and S_{8} while fixing the remaining parameters to the T17 cosmology. Additionally, we varied the intrinsic alignment amplitude, A_{IA}, and the parameter A_{ba} that accounts for the strengths of the response to baryonic feedback^{9} and was described in Sect. 5.3. Lastly, we always marginalise over the shifts of the redshift distributions δ⟨z⟩ given in H21 and the mbias (Giblin et al. 2021) of all five tomographic bins shown in Table 1. For these nuisance parameters, we assume Gaussian priors and account for the fact that the uncertainties of shifts of the redshift distributions are correlated using its correlation matrix. We give an overview of the priors in Table 2.
All varied parameters and their flat prior knowledge.
If the covariance matrix, , is measured from simulations, it is a random variable. We followed the method from Percival et al. (2022), which leads to credible intervals that can also be interpreted as confidence intervals with approximately the same coverage probability. The posterior distribution of a model vector, m, that depends on n_{p} parameters, p, if the covariance matrix is measured from n_{r} mock data, is given by:
where d is the measured data vector and
The powerlaw index m is
with n_{d} being the number of data points and
To give a quantified value for the comparison of different modelling choices, we use the Figure of Merit (FoM), which we calculate as
where C is the parameter covariance matrix of the S_{8} − Ω_{m} plane resulting from the MCMC process.
Depending on the used filter combinations, the modelling of takes around 10–20 min, which stands as an obstacle to directly running an MCMC using the model. To circumvent this issue, we use the emulation tool contained in CosmoPower (Spurio Mancini et al. 2022). We trained the emulator on 3000 model points for the analysis based on mock data in the parameter space {Ω_{m}, S_{8}, A_{IA}} and 5000 model points for the real data analysis in the parameter space {Ω_{m}, S_{8}, A_{IA}, δ⟨z⟩}, which we distributed in a Latin hypercube for the given prior range, and fixed all other parameters to the ones used in the T17 simulations. The mbias and the A_{ba} do not need to be emulated as they follow a simple correction formalism that is applied during the MCMC sampling. The δ⟨z⟩ of each z_{ph}bin for the training and testing are distributed uncorrelated between [ − 0.06, 0.06]. The accuracy of the emulator is tested by comparing the emulator prediction with the model at 500 independent points in the same parameter space. The fractional error of each vector element is smaller than 2% for E_{n} and smaller than 5% for , which is well within the accuracy of the HMCODE2020 or the BIHALOFIT model itself. We used a MetropolisHastings sampler for MCMC^{10}, where we used 1000 walkers running each 20 000 steps and cut the first 2000 steps away to ensure that the posteriors are not biased by the burnin phase.
7. Reduction of data vector combinations and filter radii
Generally, as long as the covariance is converged, the more information is used, the higher the constraining power. However, for thirdorder statistics with four different filter radii and five tomographic bins, the data vector contains 700 elements. To obtain a reliable covariance matrix, roughly ∼10 times the dimension of the data vector is needed. This implies that modelling and measuring the model/data vector and covariance matrix are timeconsuming and require thousands of simulations. Therefore, it is inevitable for this and for future analyses to compress the data vector without substantial information loss. We note that this section is only concerned with reducing the dimension of the data vector and is not meant as a proper forecast, which we do in the following sections. The spatial resolution of the mock data used to compute the covariance induces an error smaller than < 1% in the E_{n} covariance matrix. However, as we only reduce the elements of , while using all combinations of the E_{n}, this is unproblematic for the purpose of this section. We note again that the data vector was measured from one realisation with a resolution of A_{pix} = 0.18 arcmin^{2}, and that we use the κ based covariance for this section to ensure the covariance matrix is converged even for the data vector with the largest dimension. As a first check, we show in Fig. 4 multiple element choices that discard a significant part of the data vector. Using only the five autotomographic bins is not a good choice for reduction, as seen in the fact that the FoM is reduced by 32%. However, using only equalscale aperture radii is a better reduction, reducing the data vector to 28% while reducing the FoM by only ≈8%. Using only equalscale aperture radii also has the advantage of being modelled faster. Therefore, we continue using only the equalscale aperture radii for the rest of this work, which reduces the dimension of the data vector from 775 to 215.
Fig. 4. Posterior distribution for data vectors and covariance measured from the T17 convergence maps catalogues. Here only specific parts of the data vector are used. Here ‘only equal z_{ph}bins’ means that only the auto tomographic bins are used. Lastly, all nonequal filter radii are discarded for ‘only equal filter radii’. 
Next, we investigated the choice of aperture filter radii. The aperture filter radii under consideration are θ_{ap} ∈ {4′,6′,8′,10′,14′,18′,22′,26′,32′,36′}. The resulting posteriors are displayed in Fig. 5, which reveals that using only the large filter radii above > 22′ have the worst constraining power. However, using only the four smallest filter radii is also not the best choice, as they are largely correlated and miss the largescale information. The best choice of filter function is to use a filter from each range such as θ_{ap} ∈ {4′,8′,14′,32′}.
Fig. 5. Posterior distribution for different filter radii combinations if the data vector and covariances are measured from the T17 convergence maps. 
Next, we considered the full data vector, meaning the , and the E_{n} for the next compression strategy. The idea is to decrease the number of elements by considering only those with the highest constraining power on S_{8}. We start with the element with the highest S/N, then consecutively add those vector elements that maximise the S_{8} Fisher information content. This Fisher information is calculated as (Tegmark et al. 1997)
where the partial derivatives are computed with a fivepoint stencil beam (Fornberg 1988),
where and . For this analysis, we used only the equalscale filter radii θ_{ap} ∈ {4′,6′,8′,10′,14′,18′,22′,26′,32′,36′}. The first ∼200 elements are sufficient to get converged posteriors. It is also interesting to see which elements help increase the constraining power. Unsurprisingly, the first elements are all COSEBIs, but among the first 100, approximately half are elements. Furthermore, it is interesting to see that crosstomographic bin elements are more likely to be selected by our method, which is expected because they have a higher S/N than autotomographic bin elements. Nevertheless, we also observe that using equalscale filter radii θ_{ap} ∈ {4′,8′,14′,32′} results in a better FoM then using the best 215 elements. This is likely because we optimised only the S_{8} parameter here while fixing all others. Therefore, we also checked for the equal scale filter case if a principal component analysis (PCA) applied to the covariance matrix performs better in compressing the data. However, it needs more elements to get the same constraining power as our FoM maximiser. This is probably because a PCA considers only the covariance matrix and ignores the derivatives, meaning that a PCA is not necessarily sensitive to cosmology.
Finally, we give in Table 3 an overview of the FoM for the Ω_{m}–S_{8} plane and the size of the data vector for some more element choices. As the covariance matrix is measured from 5184 mock data, we can assume that for all data vector dimensions under consideration, the covariance matrix is converged. Nevertheless, for the next sections, we use the covariance matrix measured from 1944 galaxy shear catalogues, limiting us to a maximum dimension of our data vector ∼200 elements. Given the investigations in this section, we restrict the further sections to the case where we use all E_{n} elements and all elements with the equalscale filter radii θ_{ap} ∈ {4′,8′,14′,32′}, as this resulted in the best FoM for the restricted dimension of the data vector.
Overview of the data vector choices and their relative constraining power.
8. Validation of data vector estimator
In a real survey analysis, two further difficulties arise. The first issue is that the lensing information is not given in terms of convergence maps but by point estimates (galaxies). The finite area where galaxies are measured implies that no lensing information is available outside that area. It is, therefore, necessary to measure the statistic from correlation functions described in Sect. 3.3. Second, these point estimates are the reduced shear g = γ/(1 − κ), which increases the signal and therefore needs to be accounted for. To correct these effects, we measured data vectors without shape noise and with the largest available resolution A_{pix} = 0.05 arcmin^{2}. Since these data vectors result only from one realisation without shape noise, we can expect that the ratio is a good approximation for the reduced shear effect. Although the deviations are small (as seen in Figs. A.3 and A.4), we decided to scale the model vectors by the ratio of g and γ data vector.
The resulting posteriors are shown in Fig. 6. Our first observation is that similar to the finding of H23, combining secondorder with thirdorder shear statistics significantly improves the constraining power on S_{8} and Ω_{m}. Compared to the E_{n}only, the S_{8}Ω_{m} FoM increases by 93% and for the only by 233%. For these improvements, we have not considered that the posteriors of the individual statistics are bound by the priors on Ω_{m}, meaning that the improvements are lower bounds. Compared to the analysis based on κ mock data (see Fig. 5), the posteriors are broader because the covariance based on κ maps also uses information outside the patches since the boundaries were not removed. This decreases the variance between the patches. A further difference is that the κ analysis is not subject to masks and uses a slightly larger effective area. Although we rescaled the covariance to correct for the different effective areas, we found in Linke et al. (2023) that this rescaling is not necessarily accurate. Lastly, we notice that our modelling within the KiDS1000 uncertainty is accurate, which we quantified by measuring the shift of the maximum of the posterior (MP) from the true values in the matter density parameter ΔΩ_{m} < 0.02 σ_{Ωm} and in the clustering amplitude ΔS_{8} < 0.05 σ_{S8} with respect to the averaged noisy mock data vector and the KiDS1000 uncertainty. We define the MP as the maximum of the onedimensional marginal distributions.
Fig. 6. Posterior distribution for data vectors and covariance measured from the T17 galaxy shear catalogues. The covariance matrix and reference data vector were measured from 1944 noisy g mock data. 
9. Cosmological results
Finally, we are ready to present the first cosmological constraints from the KiDS1000 data using second and thirdorder shear statistics, displayed in Figs. 7 and 8. The data vectors are described in more detail in A21 for the E_{n}, and the in Porth et al. (2023). Given the fact that the covariance matrix is measured only from 1944 realisations using the reduced shear g, we decide to build our data and model vector from all E_{n} modes and with all equal filter radii θ_{ap} ∈ {4′,8′,14′,32′}. To control whether the model accurately fits the data, we minimised the from the real data and the from each of the 1944 T17 mock data used to compute the covariance matrix. To estimate the probability of measuring a (pvalue), we counted the number of that are greater than divided by 1944. The resulting pvalue for the E_{n}only, the only, and the combination are given in Table 4. The resulting pvalues for all combinations are better than 0.1, indicating that our covariance is well matched to the observed data and our model is accurate enough to describe the data. Our maximum posterior χ^{2} value increases to 81 if we swap our E_{n} covariance matrix to the analytical expression in Joachimi et al. (2021) calculated at the MP parameter values in A21. This is unsurprising as the numerical covariance is computed at the T17 cosmology, giving a signal larger than that computed at the MP of A21. Furthermore, our covariance matrix is measured from reduced shear mock data that slightly increases the covariance matrix. We also notice that the model and the data seem inconsistent for some z_{ph}bin combinations. However, adjacent COSEBI modes are highly correlated, so visually inspecting the model’s goodness of fit to the data is misleading. We refer to A21 for further details on this discrepancy.
Fig. 7. Measured and modelled E_{n} for the first five components. The blue points show the measurements from the KiDS1000 data (A21). The blue error bars indicate the KiDS1000 uncertainty. The different dashed lines show analytical descriptions at the MP if the combination of or only E_{n} is used. 
Fig. 8. Measured and modelled with filter radii θ_{ap} ∈ {4′,8′,14′,32′}. The blue error bars indicate the KiDS1000 uncertainty. The different dashed lines show analytical descriptions at the MP if the combination of or only is used. 
MP values with their marginal 68% credible intervals.
We show the resulting posteriors in Fig. 9, where we marginalised over the shift in the redshift distribution and multiplicative shear correction, both stated in Table 1. We improve the constraints on S_{8} by 23% and on Ω_{m} by 47% compared to the E_{n}only case. This shows how powerful a combined analysis of the second and thirdorder shear statistics is. The constraints on A_{IA} and A_{ba} are basically untouched, showing that is not helpful for constraining these nuisance parameters.
Fig. 9. Posterior distribution for the real KiDS1000 data vector while the covariance is measured from the 1944 T17 galaxy shear catalogues. Here the E_{n} are compared to using all available combinations for the aperture filter radii θ_{ap} ∈ {4′,8′,14′,32′} and redshift bin combinations. We note that the prior (see Table 2) range on Ω_{m} is enlarged compared to the previous validation plots. The FoM values of the E_{n}only or only case should not be compared to Fig. 6, as priors of Ω_{m} bind their contours. 
Compared to the maximum of the onedimensional marginal distributions constraints of E_{n} measurements given in A21 ( and ), we have slightly larger constraints when using only the E_{n}. This is because we use a numerical covariance, which is larger than the analytical one due to the underlying cosmology and the fact that we model the reduced shear effect. Furthermore, we note that A21 varied Ω_{cdm}h^{2}, Ω_{b}h^{2} and h, while we fixed h and Ω_{b} and varied only Ω_{m}. We also find that if we use the same pipeline as A21 but allow larger Ω_{cdm}h^{2}, the posteriors increase towards larger Ω_{m} and therefore get more consistent with our results. Furthermore, we use a different sampler compared to A21 and different baryon feedback process modelling. Nevertheless, our results from the E_{n} analysis are consistent with A21 within 0.03 σ in Ω_{m} and within 0.14 σ in S_{8}. Similar to A21, we also perform an internal consistency check, removing one z_{ph}bin at a time and finding consistent results. We discuss this in more detail in Appendix A, and we find that (as expected) the fifth z_{ph}bin is most important to constrain Ω_{m} and S_{8}. We further discuss in Appendix A modelling checks regarding the infusion of IA, baryon feedback and the reduced shear correction. The baryon feedback and reduced shear corrections are always subdominant compared to a shift in S_{8}. The IA, however, is important, especially if lower z_{ph}bins are included and must be accurately modelled.
Finally, we notice that all statistics are consistent with A_{IA} = 0, although both statistics alone seem to favour positive A_{IA}. Interestingly, the joint analysis is shifted to lower A_{IA} and is slightly more constraining, which we do not observe in the validation in Sect. 8 where all posteriors accurately peak at the input value and all constraining power comes from the E_{n}. This might indicate that for real data with nonzero A_{IA}, thirdorder shear statistics can contribute (at least a bit) to constraining IA, in line with predictions by, for instance, Pyne & Joachimi (2021) and Troxel & Ishak (2012). Furthermore, the baryonic parameter can be confined to A_{ba} < 1.4 at 68% confidence with our statistic. Here, we should note that does not help in constraining A_{ba}, which is probably due to the fact that changes in S_{8} absorb all A_{ba} effects.
10. Conclusions
This work validates the combined modelling of second and thirdorder shear statistics, showing that its accuracy is wellsuited for a KiDS1000 analysis. Our secondorder shear statistic of choice is the E_{n}modes of the COSEBIs (Schneider et al. 2010) and the thirdorder statistic is (Schneider et al. 2005). In particular, we incorporated intrinsic alignment modelling based on the nonlinear alignment model of Bridle & King (2007) and validated its accuracy against simulations infused with IA effects. This test is also interesting for other simulationbased analyses for which IA cannot be modelled analytically. We incorporated the impact of the baryonic feedback process by measuring a response function using the Magneticum simulations. Since the amplitude of this response function has no physical meaning, it is considered a nuisance parameter, which does not bias our cosmological parameter predictions.
We investigate which parts of the data vector can be neglected without losing too much cosmological information. This is important because the data vector for thirdorder shear statistics, due to the possibility it offers of combining three different filters with three different redshift bins, inflating the data vector to several hundred elements easily. Therefore, both a numerical and an analytical covariance are difficult to compute. We find that crosstomographic redshift bins contain a large amount of cosmological information. Using only equalscale filter radii but all available tomographic bin combinations was the best data compression strategy, which comes with the advantage that equalscale filters are faster to compute analytically as they require a lower integral accuracy. Next, we investigated the chosen filter radii. For this we measured and modelled the for θ_{ap} ∈ {4′,6′,8′,10′,14′,18′,22′,26′,32′,36′}. We find that filter radii above 30′ do not contribute much constraining power and that having many small filter radii is unnecessary. The best option is to use filter radii θ_{ap} ∈ {4′,8′,14′,32′}. We also tested whether selecting elements that maximise the Fisher information matrix on S_{8} is more optimal but find no difference compared to selecting all equalscale filter radii. However, this method can be used to speed up the modelling and measurement of future analyses by discarding irrelevant elements.
Next, we validated the assumption that the correlation function estimator gives accurate results. This is important since only correlation function estimators have the potential to give unbiased results for real data. We used realistic mock data created from the T17 simulations. In particular, we created several galaxy catalogues where the positions of the galaxies are exactly at the KiDS1000 galaxy positions. We had to rely on correlation functions to measure the second and thirdorder statistics, which give unbiased results also if the data has a complex topology. Our first finding is that our modelling and measurement result in unbiased cosmological parameters given the KiDS1000 uncertainty. Second, we find that using the reduced shear, or the shear itself, does not change the results and can, therefore, be ignored for the KiDS1000 data analysis.
We conclude this paper with an analysis of the real KiDS1000 data. Overall, our E_{n} constraints are less informative than the original KiDS1000 analysis. This is mostly because we used a numerical covariance matrix from T17 simulations. However, the chosen sampler, the priors on the cosmological parameters, and the modelling strategy of the baryonic feedback processes impact our constraints, too. We find an S_{8} = 0.772 ± 0.022 and an , which are improved compared to the E_{n}only case by 23% and by 47%, respectively. With a pvalue of 0.25, we also find a good agreement of model and data given the KiDS1000 uncertainty. This demonstrates that combining second and thirdorder statistics is powerful in constraining cosmological parameters. The gain in constraining power in Ω_{m} is also interesting for combined weak lensing and galaxy clustering analysis because the constraining power in Ω_{m} for clustering analysis comes with the issue of further nuisance parameters such as galaxy bias. However, since second and thirdorder shear statistics constrain Ω_{m} quite well, combining it with clustering statistics might enable us to learn more about these nuisance parameters.
We leave the optimisation of the IA and baryonic feedback modelling for future analysis. An especially interesting improvement would be a more physically motivated description of the baryon feedback processes, which are identical for power and bispectrum. As all baryon feedback models rely on hydro simulations, we have to use the same simulations for power and bispectrum. Furthermore, although we found that the reduced shear and limber approximation is sufficient for a KiDS1000 analysis, we probably have to model these effects for future Stage IV surveys. Lastly, we ignore the effect of source clustering for this work. Although Gatti et al. (2024) found it relevant for thirdorder weak lensing statistics based on convergence mass maps, we expect it to be less critical for our analysis, which uses shear catalogues and no mass map reconstruction. For future Stage IV surveys, this needs to be investigated.
When placing the properties (shear) of a galaxy in a pixel grid, instead of shifting it to the pixel centre, we assume that the galaxy itself has the size of a pixel and distribute its properties to all neighbouring pixel centres weighted by their relative distance to the galaxy. This improves the accuracy on scales smaller than the pixel size.
The KiDS data products are public and available through http://kids.strw.leidenuniv.nl/DR4
T17 simulations: http://cosmo.phys.hirosakiu.ac.jp/takahasi/allsky_raytracing/
For updated KiDS1000 cosmic shear analysis, we refer the reader to DES & KiDS Collaboration (2023) and Li et al. (2023).
To compute the individual γ data vectors with shape noise we have used Eq. (36) and replaced
The code we made use of can be find here: https://github.com/justinalsing/affine
Acknowledgments
This paper went through the KiDS review process, and we want to thank the KiDS internal reviewer for the fruitful comments to improve this paper. We thank Mike Jarvis for maintaining treecorr and his constructive and fruitful comments. Furthermore, we thank Alessio Mancini for developing the CosmoPower emulator, which made the analysis pipeline of this work significantly faster. Some of the results in this paper have been derived using the HEALPY and HEALIX package (currently http://healpix.sourceforge.net; Górski et al. 2005; Zonca et al. 2019). The figures in this work were created with MATPLOTLIB (Hunter 2007) and CHAINCONSUMER (Hinton 2016). We further make use of COSMOSIS (Zuntz et al. 2015a), NUMPY (Harris et al. 2020) and SCIPY (Jones et al. 2001) software packages. L.P. acknowledges support from the DLR grant 50QE2002. SH is supported by the US Department of Energy, Office of Science, Office of High Energy Physics under Award Number DESC0019301. We acknowledge use of the lux supercomputer at UC Santa Cruz, funded by NSF MRI grant AST 1828315. L.L. is supported by the Austrian Science Fund (FWF) [ESP 357N]. T.C. is supported by the INFN INDARK PD51 grant and by the FARE MIUR grant ‘ClustersXEuclid’ R165SBKTMA. K.D. acknowledges support by the COMPLEX project from the European Research Council (ERC) under the European Union’s Horizon 2020 research and innovation program grant agreement ERC2019AdG 882679 as well as support by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) under Germany’s Excellence Strategy  EXC2094  390783311. J.H.D. acknowledges support from an STFC Ernest Rutherford Fellowship (project reference ST/S004858/1). H.H. is supported by a DFG Heisenberg grant (Hi 1495/51), the DFG Collaborative Research Center SFB1491, as well as an ERC Consolidator Grant (No. 770935). K.K. acknowledges support from the Royal Society and Imperial College. N.M. acknowledges support from the Centre National d’Etudes Spatiales (CNES) fellowship. Author contributions: All authors contributed to the development and writing of this paper. The authorship list is given in three groups: the lead authors (P.A.B., L.P., S.H., L.L., N.W., P.S.), followed by an alphabetical group. This alphabetical group includes those who have either made a significant contribution to the data products or to the scientific analysis. The KiDS results in this paper are based on observations made with ESO Telescopes at the La Silla Paranal Observatory under programme IDs 177.A 3016, 177.A3017, 177.A3018 and 179.A2004, and on data products produced by the KiDS consortium. The KiDS production team acknowledges support from: Deutsche Forschungsgemeinschaft, ERC, NOVA and NWOM grants; Target; the University of Padova, and the University Federico II (Naples). Data processing for VIKING has been contributed by the VISTA Data Flow System at CASU, Cambridge and WFAU, Edinburgh. The calculations for the hydrodynamical simulations were carried out at the Leibniz Supercomputer Center (LRZ) under the project pr83li (Magneticum). J.H.D. acknowledges support from an STFC Ernest Rutherford Fellowship (project reference ST/S004858/1).
References
 Asgari, M., Heymans, C., Hildebrandt, H., et al. 2019, A&A, 624, A134 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Asgari, M., Tröster, T., Heymans, C., et al. 2020, A&A, 634, A127 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Asgari, M., Lin, C.A., Joachimi, B., et al. 2021, A&A, 645, A104 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Bartelmann, M. 2010, Class. Quant. Grav., 27, 233001 [CrossRef] [Google Scholar]
 Bartelmann, M., & Schneider, P. 2001, Phys. Rep., 340, 291 [Google Scholar]
 Beck, A. M., Murante, G., Arth, A., et al. 2016, MNRAS, 455, 2110 [Google Scholar]
 Bergé, J., Amara, A., & Réfrégier, A. 2010, ApJ, 712, 992 [CrossRef] [Google Scholar]
 Bernardeau, F., van Waerbeke, L., & Mellier, Y. 1997, A&A, 322, 1 [NASA ADS] [Google Scholar]
 Bridle, S., & King, L. 2007, New J. Phys., 9, 444 [Google Scholar]
 Brown, M. L., Taylor, A. N., Hambly, N. C., & Dye, S. 2002, MNRAS, 333, 501 [NASA ADS] [CrossRef] [Google Scholar]
 Burger, P. A., Friedrich, O., HarnoisDéraps, J., et al. 2023, A&A, 669, A69 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Castro, T., Quartin, M., Giocoli, C., Borgani, S., & Dolag, K. 2018, MNRAS, 478, 1305 [Google Scholar]
 Chisari, N., Codis, S., Laigle, C., et al. 2015, MNRAS, 454, 2736 [NASA ADS] [CrossRef] [Google Scholar]
 Crittenden, R. G., Natarajan, P., Pen, U.L., & Theuns, T. 2002, ApJ, 568, 20 [NASA ADS] [CrossRef] [Google Scholar]
 Dalal, R., Li, X., Nicola, A., et al. 2023, arXiv eprints [arXiv:2304.00701] [Google Scholar]
 de Jong, J. T. A., Verdoes Kleijn, G. A., Boxhoorn, D. R., et al. 2015, A&A, 582, A62 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 de Jong, J. T. A., Verdoes Kleijn, G. A., Erben, T., et al. 2017, A&A, 604, A134 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 DES Collaboration (Abbott, T. M. C., et al.) 2022, Phys. Rev. D, 105, 023520 [CrossRef] [Google Scholar]
 DES & KiDS Collaboration (Abbott, T. M. C., et al.) 2023, Open J. Astrophys., 6, 36 [NASA ADS] [Google Scholar]
 Di Valentino, E., Anchordoqui, L. A., Akarsu, Ö., et al. 2021, Astropart. Phys., 131, 102604 [NASA ADS] [CrossRef] [Google Scholar]
 Edge, A., Sutherland, W., Kuijken, K., et al. 2013, The Messenger, 154, 32 [NASA ADS] [Google Scholar]
 Falco, E. E., Gorenstein, M. V., & Shapiro, I. I. 1985, ApJ, 289, L1 [Google Scholar]
 Fenech Conti, I., Herbonnet, R., Hoekstra, H., et al. 2017, MNRAS, 467, 1627 [NASA ADS] [Google Scholar]
 Fornberg, B. 1988, Math. Comp., 184, 699 [CrossRef] [Google Scholar]
 Fu, L., Kilbinger, M., Erben, T., et al. 2014, MNRAS, 441, 2725 [Google Scholar]
 Gatti, M., Jain, B., Chang, C., et al. 2022, Phys. Rev. D, 106, 083509P [NASA ADS] [CrossRef] [Google Scholar]
 Gatti, M., Jeffrey, N., Whiteway, L., et al. 2024, MNRAS, 527, L115 [Google Scholar]
 Giblin, B., Heymans, C., Asgari, M., et al. 2021, A&A, 645, A105 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Górski, K. M., Hivon, E., Banday, A. J., et al. 2005, ApJ, 622, 759 [Google Scholar]
 Gruen, D., Friedrich, O., Krause, E., et al. 2018, Phys. Rev. D, 98, 023507 [Google Scholar]
 Halder, A., & Barreira, A. 2022, MNRAS, 515, 4639 [NASA ADS] [CrossRef] [Google Scholar]
 Halder, A., Friedrich, O., Seitz, S., & Varga, T. N. 2021, MNRAS, 506, 2780 [CrossRef] [Google Scholar]
 HarnoisDéraps, J., Amon, A., Choi, A., et al. 2018, MNRAS, 481, 1337 [Google Scholar]
 HarnoisDéraps, J., Giblin, B., & Joachimi, B. 2019, A&A, 631, A160 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 HarnoisDéraps, J., Martinet, N., Castro, T., et al. 2021, MNRAS, 506, 1623 [CrossRef] [Google Scholar]
 HarnoisDéraps, J., Martinet, N., & Reischke, R. 2022, MNRAS, 509, 3868 [Google Scholar]
 Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357 [Google Scholar]
 Heydenreich, S., Brück, B., & HarnoisDéraps, J. 2021, A&A, 648, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Heydenreich, S., Linke, L., Burger, P., & Schneider, P. 2023, A&A, 672, A44 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Heymans, C., Tröster, T., Asgari, M., et al. 2021, A&A, 646, A140 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hikage, C., Oguri, M., Hamana, T., et al. 2019, PASJ, 71, 43 [Google Scholar]
 Hildebrandt, H., Viola, M., Heymans, C., et al. 2017, MNRAS, 465, 1454 [Google Scholar]
 Hildebrandt, H., Köhlinger, F., van den Busch, J. L., et al. 2020, A&A, 633, A69 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Hildebrandt, H., van den Busch, J. L., Wright, A. H., et al. 2021, A&A, 647, A124 [EDP Sciences] [Google Scholar]
 Hinton, S. R. 2016, J. Open Source Softw., 1, 00045 [NASA ADS] [CrossRef] [Google Scholar]
 Hirschmann, M., Dolag, K., Saro, A., et al. 2014, MNRAS, 442, 2304 [Google Scholar]
 Hoekstra, H., & Jain, B. 2008, Ann. Rev. Nucl. Part. Sci., 58, 99 [NASA ADS] [CrossRef] [Google Scholar]
 Hunter, J. D. 2007, Comput. Sci. Eng., 9, 90 [Google Scholar]
 Ivezic, Z., Axelrod, T., Brandt, W. N., et al. 2008, Serb. Astron. J., 176, 1 [Google Scholar]
 Jarvis, M., Bernstein, G., & Jain, B. 2004, MNRAS, 352, 338 [Google Scholar]
 Joachimi, B., Cacciato, M., Kitching, T. D., et al. 2015, Space. Sci. Rev., 193, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Joachimi, B., Lin, C. A., Asgari, M., et al. 2021, A&A, 646, A129 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Jones, E., Oliphant, T., Peterson, P., et al. 2001, SciPy: Open Source Scientific Tools for Python, http://www.scipy.org/ [Google Scholar]
 Joudaki, S., Hildebrandt, H., Traykova, D., et al. 2020, A&A, 638, L1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kaiser, N., & Jaffe, A. 1997, ApJ, 484, 545 [NASA ADS] [CrossRef] [Google Scholar]
 Kannawadi, A., Hoekstra, H., Miller, L., et al. 2019, A&A, 624, A92 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kilbinger, M. 2015, Rep. Progr. Phys., 78, 086901 [Google Scholar]
 Kilbinger, M., & Schneider, P. 2005, A&A, 442, 69 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kilbinger, M., Heymans, C., Asgari, M., et al. 2017, MNRAS, 472, 2126 [Google Scholar]
 Kuijken, K., Heymans, C., Hildebrandt, H., et al. 2015, MNRAS, 454, 3500 [Google Scholar]
 Kuijken, K., Heymans, C., Dvornik, A., et al. 2019, A&A, 625, A2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Laureijs, R., Amiaux, J., Arduini, S., et al. 2011, arXiv eprints [arXiv:1110.3193] [Google Scholar]
 Li, S. S., Hoekstra, H., Kuijken, K., et al. 2023, A&A, 679, A133 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Limber, D. N. 1954, ApJ, 119, 655 [NASA ADS] [CrossRef] [Google Scholar]
 Linke, L., Heydenreich, S., Burger, P. A., & Schneider, P. 2023, A&A, 672, A185 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 LoVerde, M., & Afshordi, N. 2008, Phys. Rev. D, 78, 123506 [NASA ADS] [CrossRef] [Google Scholar]
 Martinet, N., Schneider, P., Hildebrandt, H., et al. 2018, MNRAS, 474, 712 [Google Scholar]
 Martinet, N., Castro, T., HarnoisDéraps, J., et al. 2021, A&A, 648, A115 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Mead, A. J., Brieden, S., Tröster, T., & Heymans, C. 2021, MNRAS, 502, 1401 [Google Scholar]
 Miller, L., Heymans, C., Kitching, T. D., et al. 2013, MNRAS, 429, 2858 [Google Scholar]
 Munshi, D., Valageas, P., van Waerbeke, L., & Heavens, A. 2008, Phys. Rep., 462, 67 [NASA ADS] [CrossRef] [Google Scholar]
 Percival, W. J., Friedrich, O., Sellentin, E., & Heavens, A. 2022, MNRAS, 510, 3207 [NASA ADS] [CrossRef] [Google Scholar]
 Pires, S., Leonard, A., & Starck, J.L. 2012, MNRAS, 423, 983 [NASA ADS] [CrossRef] [Google Scholar]
 Planck Collaboration VI. 2020, A&A, 641, A6 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Porth, L., Heydenreich, S., Burger, P., Linke, L., & Schneider, P. 2023, arXiv eprints [arXiv:2309.08601] [Google Scholar]
 Pyne, S., & Joachimi, B. 2021, MNRAS, 503, 2300 [NASA ADS] [CrossRef] [Google Scholar]
 Ragagnin, A., Dolag, K., Biffi, V., et al. 2017, Astron. Comput., 20, 52 [Google Scholar]
 Sánchez, A. G., Scoccimarro, R., Crocce, M., et al. 2017, MNRAS, 464, 1640 [Google Scholar]
 Schneider, P. 1996, MNRAS, 283, 837 [Google Scholar]
 Schneider, P., & Lombardi, M. 2003, A&A, 397, 809 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Schneider, P., & Seitz, C. 1995, A&A, 294, 411 [NASA ADS] [Google Scholar]
 Schneider, P., van Waerbeke, L., Jain, B., & Kruse, G. 1998, MNRAS, 296, 873 [NASA ADS] [CrossRef] [Google Scholar]
 Schneider, P., van Waerbeke, L., & Mellier, Y. 2002, A&A, 389, 729 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Schneider, P., Kilbinger, M., & Lombardi, M. 2005, A&A, 431, 9 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Schneider, P., Eifler, T., & Krause, E. 2010, A&A, 520, A116 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Secco, L. F., Jarvis, M., Jain, B., et al. 2022, Phys. Rev. D, 105, 103537 [NASA ADS] [CrossRef] [Google Scholar]
 Seitz, C., & Schneider, P. 1997, A&A, 318, 687 [NASA ADS] [Google Scholar]
 Springel, V. 2005, MNRAS, 364, 1105 [Google Scholar]
 Spurio Mancini, A., Piras, D., Alsing, J., Joachimi, B., & Hobson, M. P. 2022, MNRAS, 511, 1771 [NASA ADS] [CrossRef] [Google Scholar]
 Takahashi, R., Hamana, T., Shirasaki, M., et al. 2017, ApJ, 850, 24 [Google Scholar]
 Takahashi, R., Nishimichi, T., Namikawa, T., et al. 2020, ApJ, 895, 113 [NASA ADS] [CrossRef] [Google Scholar]
 Tegmark, M., Taylor, A. N., & Heavens, A. F. 1997, ApJ, 480, 22 [NASA ADS] [CrossRef] [Google Scholar]
 Troxel, M. A., & Ishak, M. 2012, MNRAS, 419, 1804 [CrossRef] [Google Scholar]
 Troxel, M. A., MacCrann, N., Zuntz, J., et al. 2018, Phys. Rev. D, 98, 043528 [Google Scholar]
 van Daalen, M. P., Schaye, J., Booth, C. M., & Dalla Vecchia, C. 2011, MNRAS, 415, 3649 [Google Scholar]
 van den Busch, J. L., Wright, A. H., Hildebrandt, H., et al. 2022, A&A, 664, A170 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Wright, A. H., Hildebrandt, H., van den Busch, J. L., et al. 2020, A&A, 640, L14 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Zonca, A., Singer, L., Lenz, D., et al. 2019, J. Open Source Softw., 4, 1298 [Google Scholar]
 Zuntz, J., Paterno, M., Jennings, E., et al. 2015, Astron. Comput., 12, 45 [NASA ADS] [CrossRef] [Google Scholar]
Appendix A: Consistency checks and modelling choices
As a first consistency check, we again inferred the joint analysis posteriors but removed one z_{ph}bin at a time. The results are shown in Fig. A.1 and reveal (most importantly) that all z_{ph}bins are internally consistent with each other. However, we further observe that z_{ph}bin one and z_{ph}bin two are not important for inferring cosmological parameters, as expected, given their S/N as shown in Porth et al. (2023). For the A_{IA} parameter, in turn, z_{ph}bin one is very important, as the low redshift is most sensitive to IA. Next, we observe that the third z_{ph}bin results are basically the same as for all five z_{ph}bins, meaning that the third bin can eventually be discarded in future KiDS analyses if the dimension needs to be decreased further. Lastly, if either the fourth or fifth z_{ph}bin is removed, the constraining power on Ω_{m} and S_{8} drastically decreases. For the fifth z_{ph}bin, this also affects the constraining power on A_{IA}.
Fig. A.2. Posterior distribution for data vectors and covariance measured from the 1944 T17 noisy g mock data. Here the T17 data vector is infused with IA measured from the cosmoSLICS+IA. The modelling is described in Sect. 2.2. 
Now we investigate the accuracy of the modelling strategy of IA described in Sect. 2.2. To validate our IA modelling, we measured the data vector from the cosmoSLICS+IA for A_{IA} = { − 1, 0, 1} and then took the ratio of the data vector at A_{IA} ∈ { − 1, 1} and divided it by the data vector at A_{IA} = 0. This ratio is multiplied with the averaged T17 γ mock data vector to infuse IA. For this analysis, we used only equalscale filter radii θ_{ap} ∈ {4′,8′,14′,32′}. As shown in Fig. , our model nicely recovers the input IA amplitude without shifting the other parameter posteriors. Interestingly, the larger A_{IA} the larger the FoM, which is due to the increased degeneracy breaking seen in the Ω_{m}A_{IA} panel.
As a further modelling check, we investigate in Fig. A.3 and Fig. A.4 the impact on our model vector of several ingredients, namely the reduced shear correction, the infusion of baryonic feedback processes and intrinsic alignments (IA). As a reference, we always use a model without reduced shear, A_{IA} = 0, A_{ba} = 0, S_{8} = 0.78, and Ω_{m} = 0.25. We then only change one of these effects while fixing the others. We also show the impact on the model vector when we change the cosmology to S_{8} = 0.8. The change of the model vectors is divided by the square root of the covariance diagonal, indicating the relevant importance given the KiDS1000 uncertainty. The reduced shear and baryon feedback correction are always subdominant compared to the impact of S_{8} or IA. For the higher z_{ph}bins, the S_{8} change dominates, and the IA dominates if z_{ph}bin one or two are included. This result is unsurprising given that the lower z_{ph}bins are more affected by IA as the overall shear signal is much lower. However, we observe especially for the that the largest IA S/N emerge if lower z_{ph}bins are combined with large z_{ph}bins, which is due to the fact that the large z_{ph}bins drastically deplete the noise. The takeaway from this investigation is that for the current KiDS1000 analysis, our model correction for the reduced shear and baryon feedback is sufficient. The biggest impact, especially if lower z_{ph}bins are included, comes from IA, which is the biggest concern of current weak lensing analyses.
Fig. A.3. Illustration of several effects that change the E_{n} model vector. The fiducial model is without reduced shear, A_{IA} = 0, A_{ba} = 0, S_{8} = 0.78 and Ω_{m} = 0.25. For the different points, one of the effects is changed. We scaled the model differences by the KiDS1000 uncertainty, indicating the relevance of the change in the context of this work. 
Appendix B: Additional material
This section offers a collection of figures to support the analysis in the main text. In Fig. B.1, we show the full data vector for equal filter radii, from which we showed a subpart in Fig. 3. In Fig. B.2, we show the data vector if nonequal filter radii are used and measured from the convergence T17 maps.
Fig. B.1. Same as Fig. 3, but here the measured data and model are the vector for equalscale aperture filter radii θ_{ap} ∈ {4′,6′,8′,10′,14′,18′,22′,26′,32′,36′} in the T17 mock data for z_{ph}bin combinations. 
Fig. B.2. Measured data and model vector for all combinations of aperture filter radii θ_{ap} ∈ {4′,8′,16′,32′}. The data vector results from one fullsky convergence realisation without shape noise. The grey band shows the expected uncertainty of KiDS1000 estimated with the convergence maps. The is scaled by the third root of the product of the corresponding filter radii. 
All Tables
All Figures
Fig. 1. Redshift distribution of the five tomographic z_{ph}bins of the KiDS1000 data. 

In the text 
Fig. 2. Measured and modelled E_{n} vector for the first five moments in the T17 mock data. The green and blue dots are the mean of all 1944 KiDS1000 mock data that are also used to compute the covariance matrix with a resolution A_{pix} = 0.18 arcmin^{2} and shape noise. The red dots represent the data vector measured from one fullsky T17 realisation with a resolution A_{pix} = 0.18 arcmin^{2} and no shape noise. The grey band indicates the expected uncertainty from the KiDS1000 survey. 

In the text 
Fig. 3. Same as Fig. 2, but here the measured data and model are the vector for equalscale aperture filter radii θ_{ap} ∈ {4′,6′,8′,10′,14′,18′,22′,26′,32′,36′} in the T17 mock data for some selected z_{ph}bin combinations. The full data set is shown in Fig. B.1. 

In the text 
Fig. 4. Posterior distribution for data vectors and covariance measured from the T17 convergence maps catalogues. Here only specific parts of the data vector are used. Here ‘only equal z_{ph}bins’ means that only the auto tomographic bins are used. Lastly, all nonequal filter radii are discarded for ‘only equal filter radii’. 

In the text 
Fig. 5. Posterior distribution for different filter radii combinations if the data vector and covariances are measured from the T17 convergence maps. 

In the text 
Fig. 6. Posterior distribution for data vectors and covariance measured from the T17 galaxy shear catalogues. The covariance matrix and reference data vector were measured from 1944 noisy g mock data. 

In the text 
Fig. 7. Measured and modelled E_{n} for the first five components. The blue points show the measurements from the KiDS1000 data (A21). The blue error bars indicate the KiDS1000 uncertainty. The different dashed lines show analytical descriptions at the MP if the combination of or only E_{n} is used. 

In the text 
Fig. 8. Measured and modelled with filter radii θ_{ap} ∈ {4′,8′,14′,32′}. The blue error bars indicate the KiDS1000 uncertainty. The different dashed lines show analytical descriptions at the MP if the combination of or only is used. 

In the text 
Fig. 9. Posterior distribution for the real KiDS1000 data vector while the covariance is measured from the 1944 T17 galaxy shear catalogues. Here the E_{n} are compared to using all available combinations for the aperture filter radii θ_{ap} ∈ {4′,8′,14′,32′} and redshift bin combinations. We note that the prior (see Table 2) range on Ω_{m} is enlarged compared to the previous validation plots. The FoM values of the E_{n}only or only case should not be compared to Fig. 6, as priors of Ω_{m} bind their contours. 

In the text 
Fig. A.1. Same as black contours in Fig. 9 but while removing one z_{ph}bin at a time. 

In the text 
Fig. A.2. Posterior distribution for data vectors and covariance measured from the 1944 T17 noisy g mock data. Here the T17 data vector is infused with IA measured from the cosmoSLICS+IA. The modelling is described in Sect. 2.2. 

In the text 
Fig. A.3. Illustration of several effects that change the E_{n} model vector. The fiducial model is without reduced shear, A_{IA} = 0, A_{ba} = 0, S_{8} = 0.78 and Ω_{m} = 0.25. For the different points, one of the effects is changed. We scaled the model differences by the KiDS1000 uncertainty, indicating the relevance of the change in the context of this work. 

In the text 
Fig. A.4. Same as Fig. A.3 but for the model vector. 

In the text 
Fig. B.1. Same as Fig. 3, but here the measured data and model are the vector for equalscale aperture filter radii θ_{ap} ∈ {4′,6′,8′,10′,14′,18′,22′,26′,32′,36′} in the T17 mock data for z_{ph}bin combinations. 

In the text 
Fig. B.2. Measured data and model vector for all combinations of aperture filter radii θ_{ap} ∈ {4′,8′,16′,32′}. The data vector results from one fullsky convergence realisation without shape noise. The grey band shows the expected uncertainty of KiDS1000 estimated with the convergence maps. The is scaled by the third root of the product of the corresponding filter radii. 

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.