The 2-point correlation function covariance with fewer mocks

We present FitCov an approach for accurate estimation of the covariance of 2-point correlation functions that requires fewer mocks than the standard mock-based covariance. This can be achieved by dividing a set of mocks into jackknife regions and fitting the correction term first introduced in Mohammad & Percival (2022), such that the mean of the jackknife covariances corresponds to the one from the mocks. This extends the model beyond the shot-noise limited regime, allowing it to be used for denser samples of galaxies. We test the performance of our fitted jackknife approach, both in terms of accuracy and precision, using lognormal mocks with varying densities and approximate EZmocks mimicking the DESI LRG and ELG samples in the redshift range of 𝑧 = [ 0 . 8 , 1 . 1 ] . We find that the Mohammad-Percival correction produces a bias in the 2-point correlation function covariance matrix that grows with number density and that our fitted jackknife approach does not. We also study the effect of the covariance on the uncertainty of cosmological parameters by performing a full-shape analysis. We demonstrate that our fitted jackknife approach based on 25 mocks can recover unbiased and as precise cosmological parameters as the ones obtained from a covariance matrix based on 1000 or 1500 mocks, while the Mohammad-Percival correction produces uncertainties that are twice as large. The number of mocks required to obtain an accurate estimation of the covariance for the 2-point correlation function is therefore reduced by a factor of 40-60. The FitCov code that accompanies this paper is available at this GitHub repository.


INTRODUCTION
A new generation of cosmological surveys such as Dark Energy Spectroscopic Instrument (DESI Collaboration et al. 2016, DESI ★ E-mail: strusov@lpnhe.in2p3.frCollaboration et al. 2022) have started taking data and even more will in the coming years with e.g. the start of operations of Euclid (Laureĳs et al. 2011) and the Vera Rubin Observatory (Ivezić et al. 2019).Therefore, it is becoming vital to develop methods for deriving covariance matrices in order to estimate the uncertainties on the cosmological parameters of interest.
Existing methods of evaluating the covariance matrix that quantifies the errors on the galaxy 2-point correlation function of galaxy redshift surveys can be separated into three different categories: mock-based, analytic and internal, each best suited to different scenarios.Mock-based covariance matrices are built from a large suite of numerical simulations, 'mock' catalogues, that mimic the properties of the cosmological surveys with high fidelity.These mocks need to be i) accurate in the sense that they have to reproduce the two-and higher-point statistics with limited biases and ii) numerous in order to avoid sample variance, which introduces noise in the covariance matrices that could bias the inferred parameter uncertainties (e.g.Dawson et al. 2013, Percival et al. 2014).
Analytic approaches provide expectation values of the large-scale structure statistics directly and are much less computationally expensive.However, that requires a description of the non-Gaussian terms that enter the four-point correlation function, which is needed to compute the covariance of the two-point correlation function.Accurate modelling of the non-linear gravitational evolution, galaxy bias, redshift-space distortions and shot noise is thus a challenge to compute analytic covariance matrices.The modelling usually relies on Perturbation Theory (PT) which limits the domain of accuracy to the quasi-linear regime when the density perturbations remain small compared to unity.Moreover, one also needs to account for survey geometry and window function effects.Recent progress in this direction has been made to develop codes for the power spectrum (CovaPT, Wadekar et al. 2020).Additionally, we can mention semianalytic approaches, which use the data to calibrate themselves, for example, RascalC code (Philcox et al. 2019, O'Connell et al. 2016).
Finally, data-based or internal methods, such as jackknife and bootstrap, are often used especially when large sets of mocks are not available.They consist in resampling the survey data by slicing the original data into sub-samples and weighting these sub-samples following specific prescriptions.In the standard jackknife approach, for a given jackknife realisation , the sub-samples have unit weight except the sub-sample indexed  that is weighted 0, hence this approach is also called 'delete-one' jackknife resampling.Internal resampling methods do not rely on any assumption about the underlying gravity model and are thus less sensitive to unknown physics.However, they can lack precision and suffer from biases, as discussed in Norberg et al. (2009), Friedrich et al. (2016) and Favole et al. (2021).One fundamental deficiency of all internal covariance estimators is the large-scale bias coming from a lack of a proper estimation of the super-sample covariance (Lacasa, Fabien & Kunz, Martin 2017), which is due to the lack of modes larger than the survey size.Recently, a correction to the standard jackknife resampling method was proposed in Mohammad & Percival (2022) which consists in introducing a different weighting scheme for the cross-pairs than for the auto-pairs, where the auto-pairs are made up of objects that lie in the same sub-sample and cross-pairs of two objects that reside in two distinct sub-samples.Indeed, the choice of assigning weights to pairs of objects is arbitrary and Mohammad & Percival (2022) tested different prescriptions.They found that by adjusting the weighting of the pairs that compose the estimates of the two-point correlation function, they were able to provide more accurate estimates of the variance than the standard jackknife.However, it remains an internal estimator, with the associated characteristic fundamental problems such as super-sample covariance.
In this work, we follow a similar methodology but propose to go beyond that work by i) considering some cross-pairs that were neglected in both the standard jackknife and the jackknife method with Mohammad-Percival correction, ii) fitting the appropriate weighting scheme to a mock-based covariance built from a smaller number of mocks than for traditional mock-based approach.The paper's outline is as follows: in Section 2 we review the formalism associated with the standard jackknife resampling method and the correction proposed in Mohammad & Percival (2022).We introduce there the formalism of our proposed hybrid approach, which performance on mocks is presented in Section 3 and compared with the original correction for jackknife and with mock-based method for estimating the covariance matrix.We conclude and discuss further prospects in Section 4.

COVARIANCE ESTIMATORS
In the present paper, we work in configuration space.We use the Landy-Szalay estimator with double-counting assumed, (Landy & Szalay 1993), which can be written as: where  is the redshift space separation of a pair of galaxies,  is the cosine of the angle between the separation vector and the line of sight,  (, ) is the 2-point correlation function in redshift space,  (, ) are the binned auto-pair counts of the data catalogue, (, ) are the binned pair counts computed from a matching random catalogue, and  (, ) are the binned cross-pair counts between the random and the data catalogue.All pair counts are assumed to be suitably normalised in Eq. 1.
The 2-point correlation function can be decomposed into Legendre multipoles defined as: where ℓ is the order of the multipole, and  ℓ () are the Legendre polynomials.

Covariance from data or data-like mocks
Cosmological simulations can be divided into two categories: i) precise and expensive computationally N-body simulations, which are known to treat properly non-linear gravitational evolution; ii) less accurate approximate mock methods, such as BAM (Balaguera-Antolínez et al. 2018), COLA (Tassev et al. 2013), EZmock (Chuang et al. 2014), (Zhao et al. 2021), FastPM (Feng et al. 2016), GLAM (Klypin & Prada 2018), lognormal, PATCHY (Kitaura et al. 2016) etc.They can provide a good covariance for scales > 10 ℎ −1 Mpc, but small-scale clustering is not properly resolved.Assuming a survey with  m mocks, the covariance matrix of the 2-point correlation function is defined as: where [  ]  is the  th bin of correlation function of the  th mock, and ⟨  ⟩ is the mean over the  m mocks of the  th bin of the correlation function.
However, for some subsets of modern surveys, like the DESI Bright Galaxy Survey (BGS), the number of galaxies and their number density sometimes becomes so large, that even these approximate methods become expensive computationally, posing a problem.

Jackknife covariance
Jackknife is a data resampling approach that involves creating multiple sub-samples of the same dataset by systematically excluding regions of the data.When applied to the cosmological surveys, the footprint is divided into regions of similar area and it is these that are systematically excluded to make the multiple sub-samples.
This approach has the advantage of making no assumptions regarding non-linear evolution and non-standard physics, and at the same time is extremely cheap from the computational perspective, as it does not require expensive production of thousands of mocks.Assuming we have cut our dataset into  jk pieces, the covariance matrix is: where [  ] is the  th bin of the correlation function of the  th jackknife region, and ⟨  ⟩ is its mean over all the  jk jackknife regions.The coefficient on the right-hand side is larger than the corresponding factor in Eq. 3 as it compensates for the reduction in the covariance due to the overlap between the subsamples.
In practice, we consider the galaxy 2-point correlation function and the ,   and  pair counts mentioned in the Landy-Szalay estimator defined in Eq. (1).

Standard approach
We will assume the number of sub-samples is  jk and work in terms of pair counts rather than correlation functions.For simplicity, we will denote as   the auto-counts that are contributed by pairs of galaxies that both reside in the  th area of the survey (the areas that are systematically excluded to form the jackknife sub-samples) and   the cross-counts between galaxies in this  th area and those in the jackknife sub-sample that is made by excluding this area.The counts in the jackknife sub-sample   are related to the overall number of counts in the full survey  tot and the above quantities by (5) where in defining each of these pair counts we count each unique pair only once.The total number of auto-and cross-pairs can be related to their means over the jackknife samples by and, as we account for double counting with the cross-pairs only while looking at the full sample, we need to divide the obtained estimate by 2 to be consistent with the auto-pairs: where =1   .Following (Mohammad & Percival 2022), we choose to define an estimator of the normalised auto-pairs  a, in a specific realisation, such that  a =  by and the estimator of the normalised cross-pairs  , such that  c =  by where it was taken into account that the cross-pairs contribute to the total estimate twice, while the auto-pairs only once.We can then further compute for each jackknife realization the deviation from the mean value of the auto paircounts and cross paircounts We can now express how the covariance of each type of pair count can be represented in terms of the estimators above, if we assume the following definition for the covariance, where   are just some pair counts of type : By replacing ( 1 ,  2 ) by ( , ) or (, ) or (, ) in Eq. 12 and using Eqs. 10 and 11, one obtains: 2 jk ( jk − 1) This gives all the components needed to compute the covariance of , using its definition in Eq. 5: Note how the terms scale differently with the number of the jackknife regions.Mohammad & Percival (2022) argue that this inconsistent scaling is the source of the bias that arises with the standard jackknife approach.In the next sections, we will see how adjusting this scaling can enable one to recover an unbiased covariance estimator and demonstrate the need for going beyond the Mohammad-Percival correction to get unbiased covariance estimators in all regimes of galaxy number density.

Mohammad-Percival correction
Mohammad & Percival (2022) proposed to weight the cross-pairs  in order to fix the mismatch in the scaling, as seen in Eq. 16.With this weight  multiplying all the  pair counts, the expression for   becomes The definition of  , is then generalised to: which also changes slightly the mean of this quantity as   () = .
Following the steps from equations ( 9), ( 11) and ( 14), the modified expression for the covariance of the  paircounts weighted by  is We see that for  = 1 we recover the ordinary jackknife, as it will remove the cross-pairs in the same way as it removes the autopairs.Alternatively, by choosing  =  jk / 2 + √ 2( jk − 1) we can achieve equal scaling for the first two terms.Therefore, under the assumption of cov(, ) = 0 we indeed have all the terms scaling with  jk in same manner, which can be seen by rewriting the expression for cov( (),  ()) as In order to illustrate the effect of introducing the  weighting of Mohammad & Percival (2022), we create 1000 Poisson random catalogues in a box with a size of 1 Gpc/h, divide them into 125 cubic regions and then compute the covariance matrices of the real-space correlation function.We do this for both the standard jackknife and jackknife with the Mohammad-Percival correction.The results are presented in Fig. 1.We show the ratio of the mean of the diagonal elements,  2 ≡   , of the covariance matrix between jackknifebased  jk and mock-based  (estimated directly using Eq.3), where the blue curve uses the standard jackknife and the orange one includes the Mohammad-Percival correction.The standard jackknife is overestimating the covariance with respect to that from the mocks, while introducing the  weighting of Mohammad & Percival (2022) for the cross-pairs removes this bias.

Hybrid approach
The real galaxy density has physical correlations and so galaxy distributions are not Poisson distributions.Therefore, the assumption of cov(, ) = 0 is not valid.With the  weighting of the cross-pairs that was introduced in Section 2.2.2, Eq. 15 becomes Comparison of the accuracy in the estimate of the diagonal elements of the covariance matrix for the real-space correlation functions as a function of scale obtained from 1000 cubic box independent mock catalogues.The ratio is the mean of the diagonal elements obtained using different jackknife approaches to those obtained directly from the ensemble of mocks.The noticeable scale-dependent bias that is visible for the standard jackknife estimate is absent when the Mohammad-Percival correction is employed.
We can see that adopting any general fixed value of  unfortunately leaves the scaling of cov(, ) different from those of cov( , ) and cov(, ), so, in order to try to recover the benefits of the Mohammad-Percival approach, we are treating  as a free parameter.We propose therefore to augment the jackknife method with  weighting where the value of  is tuned by fitting the covariance estimate from a limited number of mocks.A scheme that represents the approach is shown in Fig. 2. First, let us assume we have a set of   mocks  = { 1 ...   }.Then, /  denotes the set of mocks with the  th mock removed.Then, we refer to the mock covariance from such a set /  as    [/  ].We also introduce the -dependent jackknife covariance obtained from a mock   with a chosen  weighting as    [  ] (), from correlation functions constructed with counts following eq.( 17).
Having that in our possession, we are able to estimate the uncertainty on the diagonal elements of the covariance Ξ   (diag()).First, we resample the given set of mocks and produce   covariances    [/  ].Then we compute the covariance matrix of the diagonals Ξ   (diag()), where we limit ourselves to the diagonal elements as there are not enough degrees of freedom to build a covariance of matrices (Wishart 1928): In general  m should be greater than the number of elements in the fitted part of the covariance.However, in the case of a small  m , one can restrict this to just the diagonal elements of Ξ   , to ensure that covariance matrix stays non-singular.The next step consists of finding which specific  is needed to obtain a realisation of the covariance matrix to describe    [𝑆].First, we can write the  dependent estimator of the covariance    () based on the mean of  m  dependent jackknife covariances: Then, the  2 of the   () describing the   [] can be written as: Following that, we minimise  2  by varying , such that we obtain  2  ( min ) = min(  2  ()).To justify using the Gaussian likelihood in this procedure, we first notice that we are using only the diagonals of the covariance matrix.That allows us, with sufficiently large  m , to approximate the distribution of the separate bins of the diagonals   with a Gaussian.Therefore, our proposed estimator of the  dependent covariance matrix  (fit)   can be defined as: While only the diagonal of   that is found.In the original Mohammad-Percival approach, the contribution of the cross-pairs to the covariance is adjusted to match that of the auto-pairs.Our hybrid approach allows us to adjust the cross-pair contribution on the  weighted covariance so that the covariance matches the one obtained from the limited set of mocks.We will show in the next section that by doing so, we can greatly reduce the bias that can appear for dense samples when using the fixed  weighting of Mohammad & Percival (2022).However, the hybrid approach does require more than a single mock to create a covariance estimate, but in the next section we will also show that the number of mocks needed is significantly reduced compared to a purely mock-based approach.

TESTS ON MOCKS
We test the performance of the fitted jackknife method with respect to other covariance matrix estimation methods on different sets of mocks that include RSD and some geometrical effects that we will describe in subsequent sections.For each specific set of mocks we also generate a set of matching random synthetic catalogues.
In section 3.1 we present the methodology of the tests that we perform on our mocks.In section 3.2, a set of tests is performed on lognormal mocks produced by the MockFactory code1 with three number densities to explore shot noise-dominated and sample variancedominated regimes, but also to mimic the DESI LRG and ELG samples.In section 3.3 approximate EZmocks mimicking the DESI LRG and ELG samples are used to provide a mock-based covariance matrix which has the level of statistical precision of expected from the DESI Year-5 data.The corresponding number densities can be seen in Fig. 3 for LRG EZmocks in red, ELG EZmocks in purple and the different lognormal mocks at n = (2, 5, 15) × 10 −4 [Mpc/h] −3 in blue, orange and green respectively.We use 1500 lognormal mocks for each space density, and 1000 ELG and LRG EZ mocks respectively.

Methodology
Both the random and data samples are divided into  jk = 196 jackknife regions (the results, shown in Sec.3.2, are not sensitive to  jk ) and FKP weights   for each point  in the dataset are assigned as follows: where  0 = 10 4 ℎ 3 Mpc −3 is the power spectrum estimate at the given redshift.The FKP weights (Feldman et al. 1994) minimize the variance of the power spectrum estimate for samples that have a number density that varies with redshift.Then, the correlation functions are computed using pycorr2 for both the samples and the jackknife realisations, which allows us to obtain    ,    () and , defined in eq. ( 3), eq. ( 23) and eq.( 25).
In order to test the robustness and precision of different covariance estimators using our set of lognormal mocks we have used the procedure described in Fig. 4.There we have created a set of 30 fitted and jackknife covariances and inferred cosmological parameters from 50 randomly selected mocks.In total we have 1500 pairs of covariances and mocks, which give us a set of 1500 fits for both the jackknife and fitted covariance approaches.As the mocks have the same cosmological parameters, and the covariances are considered estimators of the same underlying "true" covariance matrix, we then compared the spread of parameter values and their uncertainties to the one obtained from fitting separately each of the 1500 mocks to the conventional mock-based covariance matrix.The same is then repeated for the approximate mocks, with the difference that this time we have only 1000 mocks, bringing us to the sets of 20 fitted and jackknife covariances.

Lognormal mocks
In order to quickly test our approach with different parameters, such as number density, we produce a set of lognormal mocks which are often used as a simple approximation to the non-linear density field that evolves from Gaussian initial conditions.The lognormal distributed density contrast (ì ) is related to a Gaussian field The two-point correlation function  () is related to the correlation function of the Gaussian field   () as: So, a fiducial power spectrum () can be transformed into the correlation function  (), which is then converted to the correlation function of the Gaussian field using eq.( 28).We Fourier transform it to the power spectrum   () and eventually generate the Fourier space Gaussian field  () as: where   ,   are Gaussian random variables with unit variance and zero mean, and  is the volume of the simulation.After simulating the Fourier Gaussian field  () on the grid, we then use Fast Fourier Transform (FFT) to transform it and obtain the regular configuration space Gaussian field  ().This is then transformed into the overdensity field using eq.( 27).The expectation value for the number of galaxies in a particular cell is computed given a fixed mean number density n, and galaxies are then drawn using the Poisson distribution and placed randomly the cell.Velocities are then assigned using the linearised continuity equation: where () is a scale factor, and which is solved using Zeldovich approximation (Zel'dovich 1970).
Eventually, the RSD effect is modelled at a chosen redshift using the velocity information by affecting the coordinates of the galaxy   as: where   rsd are the redshift-distorted coordinates,  is the linear growth rate of structure, ì  is the velocity of the galaxy, and   is the line of sight.

Dependence on number density
We create 3 sets of lognormal mocks, each set containing 1500 realisations, for number densities n = 2 × 10 −4 , 5 × 10 −4 and 15 × 10 −4 ℎ 3 Mpc −3 at  = 1.Each of the realisations is made from a cubic box with a volume of (2Gpc/ℎ) 3 with grid of size 384 3 and fiducial cosmology with ℎ = 0.674,  8 = 0.816 and Ω (0) m = 0.31.The CLASS code (Blas et al. 2011) is used to generate the initial power spectrum.Redshift space distortions are then added, and each box is cut to have a footprint that covers 15% of the full sky.Each mock is then analysed in the redshift range from 0.8 to 1.2, and the corresponding randoms are generated, which are about 4 times denser than the data mocks.The procedure to obtain the fitted jackknife covariance is summarised in Fig. 2 and explained in the previous section.Here, we use  m = 50 mocks.We measure correlation functions from the mocks in bins of 5ℎ −1 Mpc.Fig. 5 presents the  parameter value distribution, obtained from the fits of the covariances.
Fig. 6 shows a measure of the relative bias Δ 2 ( ℓ )/( 2 Mock ) between a jackknife-based covariance matrix and the mock-based covariance as a function of pair separation .For simplicity we only consider the diagonal elements of each covariance matrix estimate.This relative bias is defined as where ( ℓ ) is the variance on a given multipole  obtained from the jackknife method,  Mock ( ℓ ) is the variance on the same multipole obtained from the 1500 lognormal mocks and ( 2 Mock ) is the uncertainty on the mock-based error bar, determined by applying the classical jackknife delete-one mock estimator to the set of mocks from which the covariance is estimated.
The left panel of Fig. 6 shows this relative bias of the jackknife method with the Mohammad-Percival correction while the right panel shows the result for our fitted jackknife method.In both cases, the monopole,  0 , is displayed in the top panel, the quadrupole,  2 , in the middle and the hexadecapole,  4 , in the bottom.The coloured lines show different number densities and the solid lines are the baseline configuration of 196 jackknife regions while the dashed lines show the test of using 100 jackknife regions instead.As expected, the underestimation slightly worsens with the increase in the number of jackknife regions, as predicted by eq. ( 15).
However, as the number density n increases, the underestimation of the jackknife method with the Mohammad-Percival correction becomes more and more significant, especially for n = 15 × 10 −4 ℎ 3 Mpc −3 .This underestimation is not visible on the jackknife covariance matrix estimates produced from the random catalogues as shown in Fig. 1.As explained in the previous section, the clustering of the data leads to higher covariance due to additional covariance coming from cross-correlations between  and  pair counts.
Additionally, there is no strong dependence on the number density for the fitted jackknife method which makes it more robust whatever the density regime of the galaxy sample of interest.It should be noted that for low-density regimes optimal  seems to be closer to the default value of Mohammad-Percival approach, and its fitted estimation in our method introduces additional uncertainty, which makes our method more imprecise as () decreases.

Effect on the cosmological parameters
To test the performance of different covariance estimation techniques we infer   8 ,  ∥ and  ⊥ by fitting the theoretical predictions for the multipoles to the ones from the mocks using covariances from estimators reviewed earlier.The fit is performed using a 5-parameter model, which is based on Lagrangian Perturbation Theory and includes the linear growth rate   8 , Alcock-Paczynski parameters (Alcock & Paczynski 1979)  ∥ and  ⊥ , first-and second-order biases  1 ,  2 and the effective Fingers Of God parameter (FOG)  FOG .The theoretical power spectrum  FOG is obtained using the MomentumExpansion module of thevelocileptors package (for more details, see Chen et al. 2020Chen et al. , 2021)).The Fingers-Of-God effect is modelled following Taruya et al. (2010), as where ( ) is the power spectrum without the FOG effect obtained with velocileptors,  FOG is the one-dimensional velocity dispersion and n is the LOS direction unit vector.The power spectrum  FOG ( ) is then transformed into the 2-point correlation function  th (, ) using a Fast-Fourier-Transform and from that we compute the theoretical correlation function multipoles   ℎ ℓ (, ).Once we have the correlation function multipoles  ℓ (  ), and covariance matrix , we can obtain the likelihood  (  1 , ...,   ): where  Table 1.For each of the estimation methods we tabulate the standard deviation  of (   8 −   8 )/  (   8 ), over independent fits, .For the mock covariance method  ≈ 1 (as expected when all the fits are performed consistently with the same covariance), for the fitted covariance method it is also quite close to unity, but for the jackknife method  > 1.4, which shows a much higher degree of deviation from the truth.(Hartlap, J. et al. 2007) on the inverse of the covariance matrix such that the original uncorrected covariance matrix denoted as  (orig)   and the corrected inverse covariance matrix Σ −1   are related by: where  is the number of discrete samples, and  is the number of entries in the data vector (number of bins used).We use a likelihood maximisation method to find the  2 minima using iminuit (Dembinski & et al. 2020).Errors are estimated from the region of Δ 2 = 1 of the marginalized  2 distribution, and they are allowed to be asymmetric.The choice of a frequentist method of analysis is motivated by its low computational cost.In Fig. 7, the first row shows the distributions of reduced  2 for different choices of n, and the other rows show the marginalised 2D-distributions of parameters and their uncertainties for respectively   8 ,  ∥ and  ⊥ .The distributions of reduced  2 show the goodness of the individual fits for the three methods.The contours in the bottom panel show how, for all the parameters, the spread from the Mohammad-Percival jackknife in green is in general much wider than the one from the mock covariance in red both in terms of uncertainty and parameter values.While in case of the fitted jackknife covariance, the blue contours are very similar to the mock covariance ones.Presumably, this improvement comes from using 50 realisations rather than one.In Fig. 8 we also show in the same form the performance from the standard jackknife in comparison with the Mohammad-Percival corrected jackknife and mock-based covariance.As expected, the standard jackknife produces slightly larger contours, which are noticeably shifted with respect to the mock covariance, especially for   8 .To additionally test the validity of our inference approaches, we will define the quantity where  is an inferred parameter from a specific fit, η is the mean from all the fits, and () is the error estimation from a specific fit.The distribution of quantity  is called a pull distribution.If  follows a Gaussian distribution, the distribution of  will form a normal distribution with x = 0 and () = 1.
For the mock covariance, we fit the 1500 available samples, while for the Mohammad-Percival jackknife and for the fitted jackknife 50 random mocks are fitted using 30 realisations of the covariance, under the assumption that all of the covariance estimators are probing the same underlying likelihood.
Pull distributions for   8 ,  ∥ and  ⊥ are presented in Fig. 9 for each number density of the lognormal mocks.The fitted jackknife and mock covariance pull distributions have Gaussian-shape with  = 1 normal distributions as expected, while the pull distributions obtained when using Mohammad-Percival jackknife are slightly wider, which is due the covariance being less precise.We can see it quantitatively in the Table 1, where the standard deviations of the distributions from Fig. 9 are presented.That is due to various shifts of the distributions obtained from fitting to different jackknife covariances.This is not the case for the fitted approach, however.
Overall, for all the number densities, the performance of the fitted jackknife method using 50 mocks is much better than that of the standard jackknife with the Mohammad-Percival correction, and, most importantly, it gives similar performance as the covariance matrix created from 1500 mocks.
We also test the performance of the approach when varying the number of mocks used for producing the fitted covariance.We test using 10, 25 and 50 mocks and report the results on the cosmological fits in Fig. 10, following the same methodology as explained before for 50 mocks.The precision on the marginalised 2D contours of the cosmological parameters of interest starts to drop noticeably when 10 mocks are used, while it remains stable between 25 and 50 mocks.

Approximate mocks
Approximate mocks based on the extended Zeldovich approximation described in Zhao et al. (2021) are used to mimic the DESI LRG and ELG samples.These mocks are expected also to reproduce the clustering in the quasi-linear regime, although they are less accurate than N-body simulations.They provide a better representation of the real survey and better reproduce the non-Gaussian effects, which are not present in the lognormal mocks.The EZmocks used here are built using a 4-parameter model that is calibrated to match the clustering of N-body simulations, the 25 AbacusSummit simulations designed to meet the DESI requirements (Maksimova et al. 2021).The 4 model parameters are: (1)  c -critical density required to overcome the background expansion; (2)  exp -responsible for the exponential cut-off of the halo bias relation; (3)  -argument in the power law probability distribution function () =  ×   of having  galaxies in the limited volume; (4)  is the standard deviation for the distribution modelling peculiar velocities.
In this work, we use a set of 1000 EZmocks generated from Nbody simulations with 6 Gpc/h box size.The fiducial cosmology employed is Planck 2018 (Aghanim et al. 2020), and the boxes are generated at  = 0.8 for the LRGs and  = 1.1 for the ELGs.We use the redshift range of  = [0.8,1.1] and the mocks are cut to a footprint that reproduces that planned for the 5-year DESI data in order to match the expected final precision of the mock-based   32) representing the bias of the specific covariance estimation approach plotted for three multipoles of LRG and ELG EZmocks (left and right panels respectively).Solid lines are with Mohammad-Percival correction and dashed lines for the fitted jackknife.covariance matrix.The comparison of the difference with the mock covariance for the single realisation of the jackknife covariance and the fitted covariance is presented in Fig. 11.
On Fig. 12 the relative bias of the diagonals of jackknife-based vs mock-based covariances as defined by eq.32 are shown for the LRG sample on the left and for the ELG sample on the right, in a similar way to Fig. 6.First, The same trend is seen for the Mohammad-Percival jackknife as we found with the lognormal mocks: the bias of the jackknife method with the Mohammad-Percival correction tends to increase with number density, so from LRG to ELG, and the fitted jackknife is still able to mitigate it.However, we can also notice that   ), where  is a separate fit for each of the methods.We can see, that for the mock covariance, it is close to 1 (as it is supposed to be when all of the fits share the same covariance.),for fitted covariance it is closing on it, and for jackknife usually takes values > 1.4, which shows a much higher degree of deviation from what we assumed to be the truth.
the differences are less pronounced in the case of the EZmocks which is due to a bigger volume being probed by the same number density.
In Appendix A, we test the impact of the size of the footprint on the diagonal elements of the covariance matrices by considering the North Galactic Cap, South Galactic Cap and full footprint separately.As in the previous section, we also infer the values of the cosmological parameters   8 ,  ∥ and  ⊥ , using the same methodology as for the lognormal mocks.The results of the fits are shown in Fig. 13 where the first row shows the  2 /dof distribution and the other rows show the marginalised 2D contours for best-fit values and uncertainties on the cosmological parameters.We confirm the findings with the lognormal mocks that the fitted jackknife method provides results which are in much better agreement with the mock-based method while the jackknife method with the Mohammad-Percival correction over-estimates clearly the uncertainties on all the cosmological parameters.The effect is also stronger as the number density of the galaxy sample increases.Moreover, as we have fewer mocks than for the tests with the lognormal mocks, we can notice that the fitted covariance based on 50 mocks actually produces smaller contours overall than the mock covariance which uses 1000 EZmocks.
In Fig. 14, we show the pull distribution as defined by eq.36 for the cosmological parameters and the standard deviations of the   8 distribution, which is taken as an example, are presented on Table 2.The results are also similar to the ones obtained with the lognormal mocks: both the fitted jackknife and mock covariances produce a Gaussian shape with  = 1, while the standard deviation of the pull distribution obtained using the Mohammad-Percival correction for the jackknife method is larger (=1.5, 1.8 for LRG and ELG respectively).This quantitative test thus demonstrates that the fitted jackknife method performs better in estimating an unbiased and accurate covariance matrix for the two-point correlation function.
Overall, throughout all of the tests for varying number densities, different types of mocks and number of fitted mocks, the fitted jackknife approach shows a considerable improvement over the correction for standard jackknife proposed by Mohammad & Percival (2022).The fitted jackknife approach can achieve an unbiased estimate of the covariance matrix with similar precision to a mock-based covariance but with the major advantage of requiring a much smaller number of mocks.

CONCLUSIONS
Obtaining an accurate covariance matrix is a key ingredient for any cosmological analysis, and raises a significant challenge due to the limitations in computing power for mock-based methods or in the assumptions used in the analytical approaches.Additionally, as was shown in a series of reviews comparing different approximate methods, they still have problems reproducing exactly the results of more computationally intensive codes, especially in the non-linear regime (Lippich et al. 2018, Blot et al. 2019, Colavincenzo et al. 2018).Some works also focused on decreasing the number of simulations needed to obtain a precise covariance matrix (Chartier et al. 2021), for example combining the results from N-body and approximate simulations.
In this work we have attempted to tackle this challenge with the use of internal resampling methods.In Section 2, we review the basics of the jackknife formalism for two-point correlation function covariance estimation and perform a test on a toy model which confirms the improvement brought by a correction to the standard jackknife approach proposed by Mohammad & Percival (2022).Instead of using an analytically fixed correction to some terms that enter the jackknife covariance matrix, we propose to fit the correction to a mock-based covariance obtained from a small number of mocks.Moreover, we also noticed an unconstrained term in the different pairs that comprise the jackknife estimate of the covariance matrix, which we propose to account for by the same fitted jackknife procedure.In Section 3, we have tested this fitted jackknife covariance method and compared its performance with respect to the jackknife method with Mohammad-Percival correction and to a mock-based approach using lognormal mocks and approximate EZ mocks.We showed that the underestimation of the covariance obtained when using the Mohammad-Percival correction increases with galaxy number density while the fitted jackknife covariance remains unbiased.Performing the cosmological inference showed that the fitted jackknife covariance based on 50 mocks performs with the same accuracy as the covariance created from 1000-1500 mocks, both in terms of precision (unbiased constraints) and accuracy (similar uncertainties).There is also a significant decrease in computational power needed and we also stress that the method is simple to implement on top of the standard jackknife covariance computation.We provide a Python package that contains the implementation of the fitted jackknife method: https://github.com/theonefromnowhere/FitCovFuture work may include further tests of such a fitted jackknife covariance estimation technique when applied to scales smaller than ∼ 20ℎ −1 Mpc.We plan to investigate the small scales in another work that aims at fitting the clustering of DESI Early Data with this method and mock-based covariances in order to estimate the galaxy-halo connection for different galaxy samples.A similar technique could also be developed in Fourier space, however, it would require a proper treatment of the window function effects when splitting the footprint into subsamples, together with a significant computational effort.We leave for future work the application of such techniques to other statistics, such as 3-point correlation function.Such a fitted jackknife covariance method can also be beneficial for multi-tracer analysis where it could accommodate all the degrees of freedom needed without requiring too many additional mocks.We plan to continue this work and apply the multi-tracer technique on the upcoming DESI Bright Galaxy Survey (Zarrouk et al. 2021, Hahn et al. 2022) whose high-density sampling make it a challenging test of the performance of the fitted jackknife covariance method.

Figure 2 .
Figure 2. Schematic describing the procedure to obtain the fitted covariance    fit as defined in Eq. (25) and discussed in Section 2.3.

Figure 3 .
Figure 3. Number density dependence on redshift for different datasets used.The lognormal mock samples were chosen to have a constant density selection function, to simplify the matters, while LRG and ELG mock samples follow the expected values from the corresponding DESI survey subsets.

Figure 4 .
Figure 4. Schematic view of the procedure to test different covariance matrix estimators, as described in section 3.1.

Figure 7 .
Figure 7.The summary of the results from the cosmological fits from the lognormal mocks with varying density (one for each column and with density in (Mpc/h) 3 indicated at the top) for the three covariance matrix estimation methods: jackknife covariance with Mohammad-Percival correction in green, fitted jackknife covariance in blue and mock covariance in red.The top panels shows the histograms of the reduced  2 , while the three bottom ones show the marginalised 2D-distributions of parameters and their uncertainties for   8 ,  ∥ and  ⊥ , obtained from the set of fits described in the Section 3.1

Figure 8 .
Figure 8.The summary of the cosmological fits from the lognormal mocks with a varying density.Similar to Fig.7but for a slightly different set of covariance matrix methods: jackknife covariance with the Mohammad-Percival correction in green, mock covariance in red (the same contours as on Fig.7) and standard jackknife covariance in blue.

Figure 9 .
Figure 9. Pull distributions for different covariance estimation techniques with results from fits on various lognormal mocks, shown for 3 different number densities indicated at the top in (Mpc/h) 3 .Line colors follow those in Fig. 7

Figure 10 .
Figure 10.The summary of the cosmological fits when using different numbers of mocks to obtain the fitted jackknife covariance: the default number of 50 mocks in red, 25 mocks in blue and 10 mocks in green.The figure is organised like Fig. 7.

Figure 11 .
Figure 11.Comparison of the deviation of jackknife and fit covariances from the mock covariance multiplied by a square of separation for multipoles ℓ = 0, 2, 4 for the EZ LRG mocks.

Figure 12 .
Figure12.The quantity defined in Eq. (32) representing the bias of the specific covariance estimation approach plotted for three multipoles of LRG and ELG EZmocks (left and right panels respectively).Solid lines are with Mohammad-Percival correction and dashed lines for the fitted jackknife.

Figure 13 .Figure 14 .Survey
Figure 13.The summary of the cosmological fits for the EZ mocks for LRGs and ELGs (left and right column respectively), similar to Fig. 7 layout.